Deconvolução de processo sísmico não-estacionário

Deconvolution of non-stationary seismic process

Research Article

Deconvolução de processo sísmico não-estacionário

L. W. B. Leite & M. P. C. da Rocha

Universidade Federal do Pará

Curso de Pós-Graduação em Geofísica

R. Augusto Corrêa, S/N

Campus Universitário do Guamá

Fone: (0XX91) 211-1692 - Fax: (0XX91) 211-1693

e-mail: [e-mail]

O presente trabalho trata da aplicação do filtro Kalman-Bucy (FKB), organizado como uma deconvolução (FKBD), para extração da função refletividade a partir de dados sísmicos. Isto significa que o processo é descrito como estocástico não-estacionário, e corresponde a uma generalização da teoria de Wiener-Kolmogorov. A descrição matemática do FKB conserva a relação com a do filtro Wiener-Hopf (FWH) que trata da contra-parte com um processo estocástico estacionário. A estratégia de ataque ao problema é estruturada em partes: (a) Critério de otimização; (b) Conhecimento a priori; (c) Algoritmo; e (d) Qualidade. O conhecimento a priori inclui o modelo convolucional, e estabelece estatísticas para as suas componentes do modelo (pulso-fonte efetivo, função refletividade, ruídos geológico e local). Para demostrar a versatilidade, a aplicabilidade e limitações do método, elaboramos experimentos sistemáticos de deconvolução sob várias situações de nível de ruídos aditivos e de pulso-fonte efetivo. Demonstramos, em primeiro lugar, a necessidade de filtros equalizadores e, em segundo lugar, que o fator de coerência espectral é uma boa medida numérica da qualidade do processo. Justificamos também o presente estudo para a aplicação em dados reais, como exemplificado.

Palavras-chave: Filtro de Kalman; Deconvolução no domínio do tempo; Processos não-estacionários.

Deconvolution of non-stationary seismic process - The present paper treats the application of the Kalman-Bucy filter (KBF), organized as a deconvolution (KBDF), for the extraction of the reflectivity function from seismic data. This means that the process is described as non-stationary, and corresponds to a generalization of the Wiener-Kolmogorov theory. The mathematical description of the KBF preserves its relationship to the Wiener-Hopf filter (WHF) that deals with the counterpart stationary stochastic process. The strategy to solve the problem is structured in parts: (a) The optimization criterion; (b) The a priori knowledge; (c) The algorithm; and (d) The quality. The a priori knowledge includes the convolutional model, and establishes statistics to its components (effective source wavelet, reflectivity function, and geological and local noises). To demonstrate the versatility, applicability and limitations of the method, we performed systematic deconvolution experiments under several situations of additive noise levels and effective source wavelet. First, we demonstrate the necessity of equalizer filters, and second that the spectral coherence factor is a good measure of the quality of the process. We also justify the present study for its application in real data, as exemplified.

Key words: Kalman filter; Deconvolution in time domain; Non-stationary processes.

INTRODUÇÃO

A finalidade central deste trabalho é o estudo e a aplicação do filtro adaptativo, tão publicamente celebrado sob o nome Kalman-Bucy (FKB), para extrair informação da função refletividade de seções de reflexão sísmica através da deconvolução do pulso-fonte efetivo. O método faz uso da propriedade de que dados reais representam realizações de um processo estocástico não-estacionário, e dados de campos dinâmicos (sísmicos e eletromagnéticos) e de campos potenciais (gravimétricos e magnéticos) são exemplos de dados geofísicos com estas características. Uma dificuldade sempre encontrada é de como descrever a componente ruído presente com a informação desejada, e a simplificação mais consistente é o conceito de série branca (passa-tudo, ou banda-passante). O interesse numa etapa posterior é dirigido à recuperação de amplitudes verdadeiras, objetivando medir propriedades petrofísicas do meio.

A solução fundamental do problema Wiener-Hopf (WH) generalizado, não-estacionário, é devido a Kalman & Bucy (1961), e referências básicas de aplicação na sísmica são Bayless & Brigham (1970), Crump (1974) e Mendel (1983), entre outros. O problema e a sua solução são organizados em duas etapas: a primeira consta da geração do sinal, e a segunda da sua avaliação. A nossa descrição trata de forma independente os modelos convolucional e o de estado, para depois fazer a correspondência entre os conceitos e os termos.

A necessidade de técnicas de deconvolução justificadas pelos objetivos traçados no processamento de seções sísmicas de reflexão, e que permitem a interpretação dos dados registrados. A representação detalhada do sinal sísmico requer um modelo relativamente complicado, e o processamento usa um conjunto de técnicas baseadas em propriedades geométricas do levantamento e estocásticas da informação.

Existem dois princípios básicos para tratar de dados observados, e em particular dados representados no sismograma: um é o determinístico, e o outro é o estocástico. O método determinístico consiste na utilização de teorias físicas de propagação de ondas envolvendo soluções de equações integral e diferenciais satisfazendo condições de contorno e iniciais. O método estocástico consiste na utilização da teoria estatística de séries temporais para serem usadas nas expressões das leis da dinâmicas como fato estatístico. Ambos os métodos são interligados através do principio da decomposição como descrito pelo teorema de Wold (Wold, 1953; Robinson, 1967).

COMPONENTES DO MODELO CONVOLUCIONAL EM 1D

O teorema de Wold estabelece que um processo aleatório não-branco pode ser decomposto sucessivamente em um determinístico e em um não-determinístico. Deste fundamento, segue que um sismograma (sinal) é modelado pela convolução de um pulso-fonte (componente determinística) com a resposta do meio ao impulso (componente aleatória), sendo a ele ainda acrescido ruídos. A descrição é baseada em Silva & Robinson (1979), onde o meio geológico é matematicamente representado pela função refletividade, s(t) , e o sinal-fonte original, u(t), escolhido aqui para as simulações é o denominado de Berlage (Aldridge, 1990). O modelo convolucional é uma condição a priori muito forte para representar o modelo do traço sísmico, e numa forma compacta as equações são:

Um modelo geológico conveniente para representar a função refletividade, s(t), é denominado de Poisson-Gaussiano, e ele é descrito pelas seguintes propriedades:

(1) A subsuperfície é relativamente horizontal, perfeitamente elástica e formada por camadas homogêneas e isotrópicas.

(2) A diferença da impedância acústica nas interfaces entre as camadas é tal que os coeficientes de reflexão correspondentes são não-correlacionáveis, média zero, e as amplitudes têm distribuição Gaussiana.

(3) Estes coeficientes de reflexão são suficientemente pequenos, de modo que reflexões múltiplas e efeitos de transmissão são ignorados.

Obviamente estas considerações não são totalmente válidas; no entanto, elas são relativamente fortes e não são incomuns na geofísica. (Mendel, 1983).

A função s(t) consiste de uma distribuição de impulsos simples, causal, descrita por:

onde a_i et _i referem-se a amplitude e localização de cada impulso, respectivamente. A partir das propriedades estabelecidas, a_i é uma variável aleatória com distribuição Gaussiana com média zero, m = 0, e variância = ; e t_i é descrito como um processo de Poisson de parâmetro l (Papoulis, 1965). A fundamentação desta distribuição no processamento está relacionada ao ganho dinâmico, e as simulações aqui realizadas têm como pressuposto, como já mencionado acima, que a distribuição dos coeficientes de reflexão é puramente aleatória (série branca), e sendo demonstrado matematicamente que, nesta condição especial e restrita, os efeitos de múltiplas são cancelados (Silva & Robinson, 1979). No entanto, os cálculos para gerar séries brancas resultam em aproximações, e consequentemente a reverberação gerada é considerada como ruído aditivo (Ver Fig.1).

Existem várias técnicas para estimar s(t), com o objetivo de resolver esta ambigüidade que não é removível, e o problema fundamental na aplicação de técnicas de deconvolução em dados reais é a incerteza no conhecimento do pulso-fonte efetivo, u(t). No WH convencional, a função autocorrelação pode ser submetida a uma janela de truncamento, e são várias as considerações: (a) Que o pulso-fonte efetivo é estacionário e de fase-mínima; (b) Que o espectro da função refletividade é branco; e (c) Que o espectro dos dados não é limitado a uma banda, ou sofreu distorções. A fundamentação da deconvolução ótima está no teorema que demonstra ser o operador WH, h(t), fortemente de fase-mínima (Robinson & Wold, 1962) e, portanto, especializada (Ver Fig. 2).

O pulso-fonte efetivo pode ser descrito como resultante de várias contribuições convolucionais independentes ao longo da sua trajetória, e que podem ser classificadas em tempo-invariante e tempo-variante, fase-mínima e fase-não-mínima, sendo elas: a assinatura da fonte (pulso-fonte original), fantasmas, anelasticidade, divergência esférica, sensor, instrumentação. A estes ainda podemos adicionar os efeitos de dispersão na propagação em camadas finas. (Burridge et al., 1988). Intencionalmente, na presente análise com ilustrações, colocamos estas componentes em conjunto, e estes efeitos são compostos numa forma global tempo-variante (Ver Fig. 3).

A componente ruído, r(t), pode ainda ser decomposta em componentes na forma:

onde u(t) representa a componente geológica relacionada à distribuição de s(t), e h(t) é a componente relacionada ao sistema de medidas e a efeitos locais e externos. Uma das características da deconvolução aqui abordada é a amplificação de ruído fora das bandas de freqüência do pulso-fonte original (e podemos também adicionar como fora da banda-passante do sensor) que rege a geração do sinal (e o registro), e isto representa uma condição indesejável. Como conseqüência, a banda de interesse da deconvolução deve ser limitada por um filtro banda-passante equalizador que pode ser incluído no processo como p(t) , definindo o pulso-fonte desejado por:

Com isto, o modelo convolucional fica submetido a uma nova descrição na forma:

O PROBLEMA WH GENERALIZADO

O caso paralelo para comparação é o filtro WH, onde o desejado é um operador ótimo, tempo-invariante, h(t), que sobre o sinal medido, g(t), minimiza um erro médio-quadrático entre a saída real, , e a saída desejada, d(t). O critério é, então, a minimização da variância dos resíduos dada por:

que resulta nas equações normais entre os desvios e as observações, e que tem a forma:

A equação da operação básica da saída é dada através da integral da convolução:

É demostrado que a resposta do filtro ao impulso é reduzida a uma equação integral denominada WH, dada por:

f_dg(t) e f_gg(t) são respectivamente as funções correlação cruzada e autocorrelação estocásticas teóricas, sendo g(t) e d(t), portanto, sinais aleatórios estacionários. A condição de ótimo estabelece uma condição a priori em que a distribuição que rege os dados é a Gaussiana. Para especificar o FWH é necessário resolver a equação integral acima, e a desvantagem considerada é a de pressupor a estacionariedade do processo, o que não satisfaz rigorosamente as características do problema em questão, que tem por descrição a não-estacionariedade. O conceito de estacionariedade utilizado é o restrito, sendo, portanto, mais uma outra condição a priori, e isto significa que ela é definida sobre a espectância e sobre a covariância como constantes no caso WH. (Leite & Alves: "Efeito de janelas de suavização e de tipo de fase nos filtros ótimos", manuscrito submetido a RBGf).

As generalidades do problema significa: (1) extensão à não-estacionariedade, h(t,t); (2) limitação dos dados a uma janela móvel, T £ t £ t₀ e t₀ £ s £ T; e (3) forma matricial para multicanal. Estas condições não são satisfeitas pela integral da convolução, sendo a saída real, , reescrita na forma de uma média-móvel de acordo com a teoria denominada de Wiener-Kolmogorov que estabelece a relação simples:

De forma semelhante ao caso WH, esta integral é relacionada às correlações f_dg(t,s) e f_gg(t,s) através da expressão:

Se a componente desejada é a refletividade, isto estabelece a correspondência vista adiante, em que s(t) = d(t). A equação integral acima é do primeiro tipo, e ela representa as dificuldades inerentes à solução geral dos problemas de deconvolução e de inversão geofísica (Polyanin & Manzhirov, 1998; Wing, 1991). Por outro lado, ela é útil na representação de processos aleatórios multidimencionais e não-estacionários, onde inclui observações finitas e estimativas tempo-variante. A solução festejada por Kalman & Bucy (1961) foi a conversão desta equação integral para equações diferenciais ordinárias lineares e não-lineares adaptáveis ao cálculo numérico, e a seguir fazemos uma apresentação resumida.

O ALGORITMO KB

A solução do problema não-estacionário é obtida através do formalismo da teoria de variáveis de estado, que na forma matricial, tempo-variante, as equações de geração são dadas por:

x(t) é a função-vetor variável de estado,F(t) , G(t) e H(t) são funções-matrizes com elementos variáveis em t; w(t) é a função-vetor de geração do estado x(t); y(t) é a função-vetor da saída selecionada através da forma de H(t); e v(t) é a função-vetor ruído aditivo no processo (Ogata, 1990).

Características gerais a priori para os processos w(t) e v(t) envolvidos são definidas na segunda parte do desenvolvimento da solução KB, através do valor médio, da autocorrelação e da correlação cruzada que, para e fixos, são dadas por:

Na nomenclatura de f(t,t), t e t são as variáveis para o cálculo posterior de integração com deslocamento, e de forma a ficar coerente com a definição do lado direito. d(t) é o delta de Dirac, que multiplicado por Q(t) e R(t) define matrizes diagonais para as autocorrelações teóricas dos canais.

O longo algebrismo para mostrar as relações entre as equações integrais WH generalizadas, as equações de estado e as premissas estocásticas, e que resultam na transformação para equações diferenciais lineares e não-lineares apropriadas ao cálculo numérico, foram analisadas anteriormente, por exemplo em Rocha (1998). A solução na forma continua está aqui resumida na Tab. 1, e a Tab. 2 contém a forma discretizada.

A aplicação do FKB a um sismograma com ruído, g_p(k), é um procedimento adaptativo, consistindo de uma seqüência de operações ponto-a-ponto, organizadas em etapas resumidas na Tab. 3. A iteração começa depois das condições iniciais, e retorna ao estágio (2) para cálculo da próxima amostra (k+1). Os superescritos se referem à extrapolação ao redor dos pontos de amostragem (Gelb et al., 1974; Lewis, 1986).

O FKB DE DECONVOLUÇÃO (FKBD)

O FKB como descrito na seção anterior está na forma geral, e deve ser especificado para cada problema. No presente caso, é necessário identificar as componentes e variáveis com o modelo de deconvolução sísmica, e para isto organizamos a descrição em 9 etapas numeradas e dadas a seguir.

(1) O principio básico é reconhecer a relação entre a média-móvel e as equações dinâmicas do FKB. Para isto, um traço sísmico discretizado, tempo-variante, é expresso na forma canônica:

A equação acima é estendida para uma estrutura multicanal, com o vetor dos traços g_j(k), e uma matriz com o pulso tempo-variante, então:

Desta relação fazemos a correspondência:

(2) A matriz do pulso-fonte efetivo, tempo-variante, H_ji(k) é construída na forma:

(3) O vetor de estado é definido como sendo a função refletividade:

O valor inicial das iterações pode ser definido por um vetor nulo.

(4) Para completar as equações do sistema sísmico dinâmico é necessário estabelecer o processo recursivo de geração do vetor de estado. A fórmula de Crump (1974) propõe que:

onde w(k) é um processo estocástico branco já definido anteriormente, e êle representa a geração do estado s(k). Esta equação extrapola os coeficientes s(k) através dos L coeficientes anteriores, mas esta forma não tem um significado físico, sendo necessário definir os coeficientes b_i(k) através de estratégias especiais. Uma proposta é a de inclusão de um modelo direto no processo, como o do tipo Goupillaud, para ser usado na extrapolação de valores b_i(k), o que representa o elo com a teoria de propagação de ondas elásticas. Nas simulações aqui realizadas, os valores de b(k) foram escolhidos a partir de uma seqüência de experimentos, e o valor desta função foi mantido constante e igual a 1/40.

(5) A equação de estado é escrita como:

onde:

As expressões acima constituem as equações de estado do sistema dinâmico para o modelo convolucional.

(6) A matriz de covariância do erro, P(k), é definida como uma matriz identidade para iniciar o processo.

(7) A matriz diagonal representa a variância do ruído associado com a saída do processo, e consta como uma informação a priori. De forma prática, podemos relacionar diretamente com a análise de ruído no sismograma.

(8) A matriz diagonal representa a variância da componente aleatória associada com o processo de geração dos coeficientes de reflexão, e é considerada como uma condição a priori.

(9) Para completar as etapas acima é necessário informar o pulso-fonte efetivo. Neste caso, ele é admitido como sendo de fase-mínima, e estimado a partir da autocorrelação do sismograma, e da aplicação da transformada de Hilbert (TH). O fluxograma deste cálculo numérico é descrito em mais outras 6 etapas descritas imediatamente abaixo.

(a) A autocorrelação f_gg(t) é calculada em janelas especificadas.

(b) A transformada de Fourier (TF) é aplicada a

(c) O logaritmo natural de é calculado para formar ln A(f).

(d) A fase-mínima é calculada pela TH:

(e) O pulso-fonte efetivo é composto por:

(f) A TF inversa gera o pulso-efetivo:

Como a descrição do modelo convolucional é independente do de variáveis de estado, reunimos as equações na Tab. 4, e na Tab. 5 escrevemos a equivalência entre os símbolos matemáticos.

LIMITAÇÕES DA FUNÇÃO REFLETIVIDADE

Seções de reflexão sísmica de afastamento-nulo contabilizam deformações espectrais devido ao processamento de correção de sobre-tempo normal, que é uma variável do modelo de macro-velocidades (estiramento do pulso), e adicionalmente do empilhamento (média sobre os traços), tudo com base na teoria do raio. Outras componentes canônicas, tempo-variante, dependentes do parâmetro horizontal do raio, podem ser adicionadas à descrição do pulso-fonte efetivo para representar estes efeitos ao longo da trajetória. Esta etapa do processamento é muito importante, e a refletividade na seção afastamento-nulo deve ser considerada como uma média dentro de um intervalo de ângulos de incidência.

A aplicação do operador de média-móvel é canonicamente expressa por:

que resulta na representação por resíduos na forma:

O operador FK, h(t,t), inclui na sua formulação uma estimativa independente do pulso-fonte efetivo, u_p(n). Isto implica numa deconvolução com resíduos u_r(t) e r_r(t) devido à diferença entre os sismogramas de dados reais daqueles teóricos abstratos.

A qualidade da imagem de saída do filtro é julgada por um interprete com base na resolução da imagem, e se ela é plausível. Um procedimento mais sofisticado, e que demanda mais esforços, é o que busca a extração de parâmetros físicos das rochas, o que é mais complicado se a distribuição for espacial 3D.

APLICAÇÃO EM DADOS SINTÉTICOS

A Fig. 1 contém a função refletividade simples em duas formas: uma sem ruído, e a outra com ruído geológico aleatório. A Fig. 2 contém o pulso-fonte original usado para gerar os sismogramas através da convolução linear com as funções refletividade e, neste exercício, ele tem a propriedade, mas não restritiva, de ser de fase-mínima. As outras duas partes da Fig. 2 contém os pulsos-efetivos para deconvolução, como extraídos dos sismogramas pela transformada de Hilbert e pela deconvolução homomórfica.

A Fig. 3 contém sismogramas escolhidos para exemplificar a performance do FKBD na presença de três níveis diferentes de ruído como informados nas figuras. O valor da relação sinal/ruído (S/R) foi calculada pela razão simples das variâncias do tipo:

onde e são as médias aritméticas simples.

Como o processo requer o conhecimento do pulso-fonte efetivo, montamos três casos distintos de pulsos deconvolutivos para comparação. O primeiro caso considera o sinal-fonte original, com os resultados mostrados na Fig. 4(a, b, c) para os três níveis de ruído. O segundo caso considera o sinal-fonte estimado pela TH (Fig. 2), com os resultados mostrados na Fig. 5(a, b, c) repetido para os três níveis de ruídos. O terceiro caso considera o sinal-fonte estimado pela deconvolução homomórfica (Gomes, 1998), com os resultados mostrados na Fig. 6(a, b, c) também para os três níveis de ruído. Observamos em todos estes exemplos que o FKBD realiza a compressão do pulso sísmico muito bem, e que ele apresenta uma perda de resolução com o aumento de ruído no sismograma, o que é coerente.

A Fig. 7 é um exemplo selecionado para mostrar o efeito do FKBD nos espectros de amplitudes dos sinais de entrada, de saída e a sua função transferência, sendo todos eles referentes ao sismograma com a S/R=99,86. Observamos que o FKBD age muito bem caracterizado como um filtro de ganho linear ascendente em direção das altas freqüências. Estas observações sugerem a aplicação de um filtro equalizador, e para isto selecionamos o passa-baixa de Ormsby antes do FKBD, com a finalidade de atenuar as altas freqüências. Esta aplicação resulta num aumento de qualidade da resolução da saída do FKBD como pode ser visto na Fig. 8.

Além da avaliação visual, optamos em utilizar a medida de coerência espectral entre dois canais de sinal, para analisar a qualidade da saída do FKBD com relação a entrada, uma vez que este conceito tem por base a correlação cruzada. Um resultado típico é mostrado nas Figs. 9 e 10, com o sismograma ruidoso sendo um dos canais, e a refletividade recuperada como sendo o outro canal. As medidas de coerência espectral complexa são calculados pela seguinte expressão:

A_m e B_m são as transformadas de Fourier dos canais, e A^* e B^* são os seus respectivos complexos conjugados. Desta função são calculados a amplitude e a fase, dados respectivamente por:

Esta medida estabelece o limite máximo da potência comum, quando os dois sinais estão perfeitamente coerentes em fase numa determinada freqüência. A amplitude e a fase variam teoricamente de forma aleatória. (Leite, 1998).

Para analisarmos a aplicabilidade da função coerência descrita acima, geramos 180 sismogramas variando a razão S/R, processamos com o FKBD, e em seguida calculamos e para cada par de entrada-saída. Se observa nas Figs. 9 e 10 que o fator de coerência é bastante oscilatório na faixa de alto ruído (baixo S/R), enquanto que na faixa de baixo ruído (alto S/R, entre 75 e 100), a amplitude e a fase tem uma resposta quase plana, indicando que a razão S/R admissível nos dados é de no mínimo S/R=75. Este nível mínimo serve para medir a conservação de forma entre o sinal de saída com relação ao de entrada, o que é importante para caracterizar o processamento.

APLICAÇÃO EM DADOS REAIS

O dado real utilizado neste trabalho foi o publicado por Buhl et al. (1974), e os resultados obtidos com o FKBD descrito neste trabalho estão nas Fig. 11, 12 e 13. A seção é montada com o efeito visual da repetição, e o sismograma é descrito como tipicamente não-estacionário, e na Fig. 13 foi aplicado um ganho dinâmico. O pulso-fonte efetivo foi calculado através da TH. Neste exemplo foi aplicado um filtro passa-baixa do tipo Ormsby para e qualização. Podemos notar claramente o efeito compressivo do FKBD, possibilitando identificar melhor as reflexões individuais ("interfaces") superiores, e as dificuldades introduzidas com o ganho automático.

CONCLUSÕES

A comparação com estudos paralelos mostram que a implementação do FKBD pode ser mais simples do que a do FWH convencional, e a aplicação do FKB dispensa medidas de estacionariedade. A versatilidade do FKB está diretamente relacionada à capacidade de generalização do problema WH quanto a não-estacionariedade na janela dos dados, condição natural de dados geofísicos. O problema reescrito na forma de média-móvel matricial permite a transformação para equações diferenciais ordinárias lineares e não-lineares, através do sistema de variáveis de estado, com premissas sobre os processos estocásticos envolvidos.

Os experimentos realizados e os resultados obtidos, com exemplos aqui apresentados, nos permite concluir que o FKBD resolve bem a compressão do pulso-fonte efetivo em modelos sintéticos e reais.

O FKBD tem, interessantemente, uma resposta linear ascendente com a freqüência, o que não é desastroso em vez de ser, por exemplo, exponencial. Devido a este fato, descoberto numericamente, demonstramos a necessidade de complementar o processamento com um filtro equalizador, e que o passa-baixa de Ormsby é uma opção.

Poderíamos considerar como uma desvantagem do FKBD a incapacidade simultânea de extrair o pulso-fonte efetivo diretamente dos dados e, portanto, de realizar a deconvolução denominada de dupla por alguns autores, e de cega por outros. Isto significa que é necessário usar técnicas independentes e complementares para extrair esta informação bem como de leis para descrever as componentes tempo-variantes (anelasticidade, divergência esférica e espalhamento).

A deconvolução baseada na teoria WH convencional se limita ao operador teoricamente estacionário e de fase-mínima, negligenciando a componente de fase-não-mínima do pulso-fonte efetivo. Consideramos como vantagem o operador do FKBD ser não-estacionário, e a matriz do pulso-fonte efetivo poder incluir qualquer forma de pulso com uma distribuição qualquer de fase, o que torna o método potencialmente mais versátil, porém com mais engenharia na sua descrição teórica e na aplicação prática.

AGRADECIMENTOS

Os autores agradecem ao Prof. Dr. Enders A. Robinson pelas sugestões e comentários para o desenvolvimento deste artigo, e estendem os agradecimentos aos "referees" ocultos pela sugestões que lhes foram úteis. O primeiro autor agradece a CAPES pela bolsa de pós-doutorado no exterior (Universidade de Karlsruhe, Alemanha, Instituto de Geofísica), e ao WIT-Consorcium, na pessoa do Prof. Dr. Peter H. W. Hubral, pelas discussões cientificas e pelo apoio neste ano sabático acadêmico, o que contribuiu para a conclusão deste trabalho.

Manuscript submitted February 17, 1999

Revised version accepted July 18, 2000

DECONVOLUTION OF NON-STATIONARY SEISMIC PROCESS

The present paper treats the application of the Kalman-Bucy filter (KBF), organized as a deconvolution (KBDF), for the extraction of the reflectivity function from seismic data. The attack strategy to the problem is structured in parts: (a) The optimization criterion; (b) The a priori knowledge; (c) The algorithm; and (d) The quality. The a priori knowledge includes the convolutional model, and established statistics to its components. First, we demonstrate the necessity of equalizer filters, and second that the spectral coherence factor is a good measure of the quality of the process. We also justify the present study for its application in real data.

Wold theorem establishes that non-white process can be decomposed in a deterministic and in a random component. Based on this principle, the deterministic part of a seismogram (signal), g(t) is given by the convolution of an effective-pulse, u(t), with the medium impulse response, s(t); and the random component represents the noise, r(t). The convolutional model is a very strong a priori condition to represent the seismic trace, and in a canonical form:

The effective source-pulse can be constructed as the result from several independent convolutional contributions along its trajectory, and that can be classified as time-variant, time-invariant, minimum-phase, and non-minimum-phase. These contributions are the source-time, ghosts, anelasticity, spherical divergence, dispersion, sensor, instrumentation, etc. These effects are summarized in a global time-variant effective-pulse, u(t,t). One strong characteristic of the deconvolution here studied is the amplification outside the frequency band of the original source-pulse, and of the sensor. As a result, a band-pass equalizer, p(t), must limit the frequency band in the form: and The convolutional model is then subject to a new description as:

The medium, s(t), is conveniently represented by the so-called Poisson-Gaussian model, with the following properties: (1) The subsurface is relatively horizontal perfectly elastic, and formed by homogeneous and isotropic layers; (2) The distribution of impedance contrast is non-correlated; (3) The reflection coefficients are sufficiently small, such that multiple reflection and transmissions effects are ignored. Obviously, these conditions are not totally valid, but they are relatively strong, and they are common in geophysics. The function s(t) consists of a causal distribution of impulses, and it is described by: , with a priori conditions on its statistics.

The generalized multichannel Wiener-Kolmogorov problem is organized in the form of matrix equations, and it corresponds to finding a time varying operator h(t,t), under a condition set by the first-order integral equation:

for the input g(t), resulting in the actual output :

We keep the description of the convolutional model independent of the state model.

The solution is to transform the previous integral equation to linear and non-linear differential equations through the method of state variables. In generalized matrix form:

Where x(t) is the state variable vector function, F(t), G(t) and H(t)are matrix functions of time variable elements; w(t) is a vector function that generates the state x(t); y(t) is a vector function of the selected output by the structure of H(t); and v(t) is a vector function of additive noise in the process. The a priori conditions for the processes w(t) and v(t) are defined by the expected value of the autocorrelation, and of the crosscorrelation, for _t and _t fixed, and are given by:

The long length solution is summarized in Tab.1.

The KBF as described previously is in a general form, and it must be specified for each problem. For the present case, we identify its components and variables with de description of the seismic deconvolution model, and we summarize this correspondence with the equations organized in Tab.4.

We used the spectral coherence measure between 2 channels to systematically analyze the relationship between the output and the input of the KBDF. These measurements are calculated by the following expression:

A_m and B_m are the Fourier transforms of the channels (A^* and B^* are their respective complex conjugate). From this functions the amplitude and the phase are calculated respectively by:

This measure establishes the maximum limit of common power between the two signals, if they are perfectly coherent in phase in a certain frequency.

Zero-offsets sections account for spectral deformations due to normal moveout corrections, which is a variable in the macro-velocity model (pulse stretching), and additionally in the stacking procedures (average over the traces), and all are based on ray theory. Other canonical components, time-variant, dependent on the horizontal ray parameter, can be added to the description of the effective-pulse to represent these effects along the trajectory. This step of the processing is very important, and the reflectivity in the zero-offset section must be considered as an average in the interval of incidence angles.

The KBDF has an interesting response linear ascendant with frequency. We show the necessity for complementing the processing with an equalizer filter, and the low-pass Ormsby is an option. The performed experiments and the obtained results, allow us to conclude that the KBDF resolves well the compression of the effective-pulse in synthetic and in real data.

NOTE ABOUT THE AUTHOR

Marcus Pinto da Costa Rocha

Graduado em Engenharia Civil pela Universidade Federal do Pará (1992; Mestre em Geofísica (1998) pelo CPGf ( Curso de Pós Graduação em Geofísica da UFPa); Atualmente Cursando o doutorado em Geofísica no CPGf. Professor Assistente da UFPa, lotado no Departamento de Matemática. Membro da SBGf e SBMAC e SBPC.

Lourenildo Williane Barbosa Leite

Bacharelado em Geologia pela Universitity, Federal de Pernambuco (1996); Mestrado e Doutorado pela Saint Louis University, Missouri, USA (1983); Professor do departamento de Geofísica da Universidade Federal do Pará, a partir de 1979; Coordenador do Curso de Pós Graduação em Geofísica, CG/UFPa; Orientador de tese de Doutorado e dissertação de Mestrado. Leciona disciplina voltada à sísmica, sismologia, tratamento e inversão de dados geofísicos; Membro da SBMAC, SBGf, SSA e SBPC.

• ALDRIDGE, D. F. -1990- The Berlage Wavelet. Geophysics, 55(1):1508-1511.
• BAYLESS, J. W. & BRIGHAM, E. O. ž1970-. Application of the Kalman filter to continuous signal restoration. Geophysics, 35(1): 2-23.
• BURRIDGE, R. PAPANICOLAUO, G. S. & WHITE, B. S. -1998- One-dimensional wave propagation in a highly discontinuous medium. Wave Motion, v. 10, p. 19-44.
• BUHL, P, STOPHA, P. L. & BRYAN, G. M. ž1974- . The application of homomorphic deconvolution to shallow-water marine seismology ž Part II: Real Data. Geophysics, 39(4): 417-420.
• CRUMP, N. ž1974- A Kalman filter approach to the deconvolution of seismic signals. Geophysical Research, 39(1): 1-13.
• GELB, A., KASPER, J. F., NASH, R. A., PRICE, C. F. & SUTHERLAND, A. A. -1974- Applied Optimal Estimation. The MIT Press. Massachusetts, USA.
• GOMES, M. V. N. ž1998- Aplicaçăo de deconvoluçăo homomórfica a dados geofísicos. Dissertaçăo de Mestrado, UFPA. Curso de Pós-Graduaçăo em Geofísica, (4): 83-92.
• KALMAN, R. E. & BUCY, R. E. -1961- New results in linear filtering and prediction theory. Trans. ASME, Series D, Journal of Basic Engineering, 83:95-107.
• LEITE, L. W. B. ž1998- Introduçăo a Análise Espectral em Geofísica. FADESP. Belém, Pará.
• LEWIS, F. L. -1986- Optimal Estimation. John Wiley & Sons. New York, USA.
• MENDEL, J. M., -1983- Optimal Seismic Deconvolution. Academic Pres. New York, USA, (4): 35-130.
• OGATA, K. -1990- Modern Control Engineering. Prentice-Hall Englewood Cliffs. New Jersey, USA.
• PAPOULIS, A. ž1965- Probability, Random Variables, and Stochastic Processes. MacGraw-Hill. New York, USA, (6): 160-180.
• POLYANIN, A. D. & MANZHIROV, A. V.-1998- Handbook of Integral Equations. CRC Press. New York, USA. (10): 491-525.
• ROBINSON, E. A. ž1967 žPredictive Deconvolution. Geophysics, (4): 418-484.
• ROBINSON, E. A. & WOLD, H. ž1962- Structural properties of Stationary Stochastic Processes with Applications. Em M. Rosenblatt (editor), Proceedings of the Symposium on Time Series Analysis at Brown University. John Wiley and Sons, Inc. New York, USA. (11): 170-196.
• ROCHA, M. P. C. -1998- Aplicaçăo do método de Kalman a dados geofísicos. Dissertaçăo de Mestrado, UFPA. Curso de Pós-Graduaçăo em Geofísica, (3): 18-47.
• SILVA, M. T. & ROBINSON, E. A. ž1979- Deconvolution of Geophysical Time Series in the Exploration for Oil and Natural Gas. Elsevier Scientific. Amsterdam, Netherlands, (4): 113-168.
• WING, G. M. -1991- A Primer on Integral Equations of the First Kind: The problem of Deconvolution and Unfolding. SIAM. Philadelphia, USA. : 1-133.
• WOLD, H. -1953- A Study in the Analysis of Stationary Time Series (2^nd edition). Almqwist & Wiksell (Eds). Stockholm, Sweden.

Datas de Publicação

• Publicação nesta coleção
  10 Out 2001
• Data do Fascículo
  Mar 2000

Histórico

• Recebido
  17 Fev 1999
• Aceito
  18 Jul 2000