Inversão 1D de SEV: como funciona o ajuste de camadas

O que a inversão faz

Inverter uma SEV é procurar as resistividades e espessuras de um meio em camadas cuja curva teórica reproduz a curva de campo dentro do erro de medida. É um ajuste iterativo: a cada rodada o programa mede o desajuste entre a curva calculada e a medida e mexe nos parâmetros para reduzi-lo. O resultado é um modelo possível, não uma fotografia do subsolo. Ele depende do número de camadas escolhido, do modelo de partida e da qualidade dos dados, e convive com a equivalência.

Se a SEV ainda é novidade, comece pelo guia completo de SEV e pela introdução à eletrorresistividade. Aqui o foco é o que acontece dentro do programa de inversão.

Problema direto e problema inverso

O problema direto vai do modelo para a curva: dadas as camadas, calcula a resistividade aparente em cada AB/2. O problema inverso vai da curva para o modelo, e é resolvido chamando o problema direto muitas vezes.

• Direto (modelagem direta): tem uma única resposta. Um modelo produz exatamente uma curva.
• Inverso (inversão): não tem resposta única. Muitos modelos produzem curvas que diferem menos que o erro das leituras, e o dado não tem como escolher entre eles.

Por isso a inversão de SEV sempre combina três coisas: um modelo direto rápido e exato, um método de otimização e alguma informação de fora da curva (número de camadas, modelo inicial, regularização, poços) para escolher entre as soluções possíveis.

O modelo de camadas horizontais

A SEV 1D supõe que o subsolo é uma pilha de camadas planas, horizontais, homogêneas e isotrópicas, que se estendem lateralmente sob toda a linha. A última é um semiespaço, de espessura infinita.

Um modelo de L camadas tem 2L − 1 parâmetros: L resistividades (ρ₁ … ρ_L) e L − 1 espessuras (h₁ … h_L−1). As profundidades das interfaces são somas das espessuras. Essa hipótese é o que torna o problema 1D: a resistividade só varia com a profundidade. Onde há variação lateral forte (contatos, diques, encostas, paleocanais), a curva de campo deixa de ser a de um meio em camadas e a inversão 1D devolve um modelo distorcido, às vezes com RMS alto, às vezes não. Nesses casos, o caminho é o caminhamento elétrico ou a tomografia 2D.

Como a curva sintética é calculada

A curva sintética sai de duas etapas: a função transformada de resistividade, calculada de baixo para cima a partir das camadas, e uma integral de Hankel que a converte em resistividade aparente. Desde Ghosh (1971), essa integral é calculada por um filtro linear, uma soma ponderada de poucos valores da transformada.

A transformada de resistividade

A transformada T(λ) resume o efeito de todas as camadas para cada “frequência espacial” λ. Ela começa no semiespaço, com T = ρ_L, e sobe camada por camada pela recorrência de Pekeris:

Para λ grande (pequenas distâncias), tanh → 1 e T tende a ρ₁: a curva começa na resistividade da primeira camada. Para λ pequeno (grandes distâncias), T tende a ρ_L: o fim da curva tende à resistividade do semiespaço.

O filtro linear

Para o arranjo Schlumberger, a resistividade aparente é uma integral de T(λ) com a função de Bessel J₁. Ghosh (1971) mostrou que, com uma mudança de variável para escala logarítmica, essa integral vira uma convolução, que pode ser calculada por um filtro digital: um conjunto fixo de abscissas Xⱼ e pesos Wⱼ.

Cada ponto da curva custa só algumas dezenas de avaliações da recorrência, o que torna a modelagem instantânea e permite que a inversão a repita centenas de vezes. Desde então foram publicados muitos filtros. Guptasarma (1982), por exemplo, otimizou filtros curtos para ganhar precisão com poucos coeficientes. Um filtro curto demais, porém, pode dar resultados ruins em contrastes fortes: a precisão do filtro precisa ser testada contra soluções exatas.

Exemplo: a curva de um modelo H

A tabela abaixo foi calculada com o mesmo modelo direto do Rhoa para três camadas: 100 Ω·m com 2 m, 20,0 Ω·m com 8 m e um substrato de 500 Ω·m.

Repare em duas coisas. O mínimo da curva (cerca de 30,1 Ω·m) fica bem acima dos 20,0 Ω·m da camada condutora: a resistividade aparente nunca chega à verdadeira de uma camada intermediária. E, em AB/2 = 100 m, a curva ainda está em 168 Ω·m, longe dos 500 Ω·m do substrato. É por isso que a inversão é necessária: ler as resistividades direto da curva subestima os contrastes.

Ajuste por mínimos quadrados em escala log

A inversão escolhe os parâmetros que minimizam a soma dos quadrados das diferenças entre os logaritmos da ρa medida e da calculada. Trabalhar em log faz cada ponto pesar pelo erro relativo e mantém as resistividades e espessuras positivas.

• Por que log nos dados: a curva de SEV atravessa uma ou mais ordens de grandeza. Em escala linear, os pontos de ρa alta dominariam o ajuste. Em log, um desvio de 5% pesa o mesmo em 10 Ω·m e em 1.000 Ω·m.
• Por que log nos parâmetros: ρ e h só podem ser positivos, e variam por ordens de grandeza. Ajustar log ρ e log h resolve as duas coisas e deixa o problema menos assimétrico.
• Como o mínimo é procurado: a relação entre modelo e curva não é linear. Os métodos de Gauss-Newton, Levenberg-Marquardt (Marquardt, 1963) e de região de confiança linearizam o problema em torno do modelo atual (com o Jacobiano, a sensibilidade de cada ponto a cada parâmetro), dão um passo, recalculam e repetem até o desajuste parar de cair.
• Limites: boa parte dos programas impõe faixas para ρ e h, para o otimizador não fugir para valores sem sentido físico.

O segundo termo da fórmula, com peso β, é a regularização, descrita mais abaixo. Como o problema não é linear, o otimizador pode parar num mínimo local, e por isso o modelo inicial importa.

Erro RMS: o que é e quanto é aceitável

O erro RMS resume, em porcentagem, quanto a curva do modelo se afasta da curva medida. O valor aceitável é o compatível com o ruído das leituras, e não um número fixo.

Como regra prática, uma SEV de boa qualidade costuma ser ajustada com RMS de poucos por cento. No Rhoa, abaixo de 3% o ajuste é marcado como excelente, de 3% a 8% como aceitável e acima de 8% como ruim. Essas faixas são uma convenção do app, não uma norma. Use-as com estas ressalvas:

• Compare com o ruído: se as leituras repetidas variam 3%, um RMS de 0,5% significa que o modelo está ajustando ruído. Se variam 1%, um RMS de 6% indica que algo não está sendo explicado.
• Olhe onde está o desajuste: resíduos aleatórios são normais. Resíduos sistemáticos (todos os pontos finais acima da curva, um trecho inteiro abaixo) apontam para camada faltando, embreagem mal corrigida ou efeito lateral.
• Poucos pontos, RMS fácil: com N pequeno e muitas camadas, quase qualquer curva é ajustada. O RMS só tem sentido junto com o número de pontos e de parâmetros.
• Programas diferentes, RMS diferentes: alguns calculam o RMS sobre os logaritmos, outros sobre valores relativos, e cada um conta os pontos repetidos de um jeito. Números de programas distintos não são diretamente comparáveis.

Quantas camadas usar: parcimônia e overfitting

Use o menor número de camadas que ajusta a curva até o nível do ruído. Camadas a mais quase sempre baixam um pouco o RMS, mas passam a ajustar o erro de medida e criam estruturas que não existem.

A regra de partida é contar os mínimos, os máximos e as mudanças de inclinação da curva e somar um (os tipos de curva ajudam nessa contagem). O teste abaixo mostra por que parar aí. Os “dados” são a curva do modelo H da seção anterior com erros fixos de até ±3% somados a cada ponto. Contra esses dados, o próprio modelo verdadeiro tem RMS de 1,78%. Cada linha é uma inversão feita com o motor local do Rhoa.

Com 2 camadas o RMS fica em 50,9%: falta estrutura. Com 3, cai para 1,79%, praticamente o nível do ruído, e o modelo recupera o verdadeiro. Com 4 e 5 camadas o RMS quase não muda (1,70% com 5), mas aparecem camadas extras finas ou repetidas, que o dado não sustenta. Quando o RMS chega ao nível do ruído, acrescentar camadas é overfitting.

O chute inicial

O modelo inicial é o ponto de partida do otimizador. Num problema não linear e não único, ele decide em qual dos mínimos possíveis a inversão termina.

Um bom chute sai da própria curva: ρ₁ perto da ρa das menores aberturas, o último ρ na tendência do ramo final, as camadas intermediárias com valores mais extremos que os da curva (como mostrou o exemplo do modelo H) e as interfaces numa fração do AB/2 em que a curva muda de comportamento. O passo a passo está em como interpretar uma SEV. Há também métodos automáticos que montam um modelo direto da curva, como o de Zohdy (1989), que usa tantas camadas quantos são os pontos e depois as desloca e ajusta.

Na prática, rode a inversão a partir de chutes diferentes. Se os modelos finais forem parecidos, o resultado é robusto. Se mudarem muito com RMS semelhante, você está diante de equivalência, e a escolha precisa de informação externa.

Regularização e suavidade

Regularizar é acrescentar à função objetivo uma preferência por modelos “bem-comportados”. Ela estabiliza a inversão quando o dado não determina todos os parâmetros, ao custo de um ajuste um pouco pior.

• Amortecimento: o Levenberg-Marquardt e a regressão de cumeeira (ridge regression; Inman, 1975) limitam o tamanho de cada passo, o que evita oscilações quando um parâmetro é mal determinado.
• Suavidade: penaliza diferenças entre camadas vizinhas. No extremo, a inversão de Occam (Constable, Parker e Constable, 1987) usa muitas camadas finas de espessura fixa e procura o modelo mais suave que ainda ajusta os dados. O resultado é uma curva de ρ contínua com a profundidade, em vez de degraus.
• O peso β: pequeno demais, a regularização não ajuda. Grande demais, ela apaga contrastes reais e o RMS sobe. O peso deve ser explícito e o mesmo quando se comparam modelos.

Modelos de poucas camadas (o padrão em SEV) e modelos suaves de muitas camadas são duas leituras do mesmo dado. Os primeiros combinam com a ideia de camadas geológicas. Os segundos mostram melhor onde a resolução acaba.

Equivalência e supressão: por que um RMS baixo não é um modelo único

Pela equivalência, camadas diferentes produzem curvas quase iguais; pela supressão, uma camada pode não deixar marca na curva. Nos dois casos, vários modelos têm RMS baixo, e o ajuste sozinho não escolhe entre eles.

Equivalência em uma curva H

Uma camada condutora entre duas resistivas é determinada quase só pela condutância longitudinal S = h/ρ. Abaixo, quatro modelos com 100 Ω·m × 4 m no topo, 500 Ω·m no substrato e a camada do meio com a mesma S = 0,20 S.

As curvas equivalentes se afastam da referência no máximo 2,8% de RMS, menos que o ruído de muitas SEVs reais, mas o topo do substrato vai de 6 a 10 m. Numa curva K, a grandeza preservada é a resistência transversal T = h · ρ da camada resistiva. A calculadora de Dar Zarrouk ajuda a explorar essas faixas.

Supressão em uma curva A

Um modelo de 20,0 Ω·m × 5 m, 100 Ω·m × 5 m e 1.000 Ω·m tem resistividades crescentes: a camada do meio é intermediária. Invertendo a curva exata dele com só 2 camadas, o Rhoa encontra 20,1 Ω·m × 5,9 m / 849 Ω·m, com RMS de 0,80%. Em dados reais, com alguns por cento de ruído, a camada de 100 Ω·m é indetectável.

Calibração com poços e o SIAGAS

A forma mais eficaz de reduzir a equivalência é fixar com dado direto o que a SEV não resolve: a profundidade de uma interface, a espessura de uma camada ou a resistividade típica de uma litologia.

• Poço próximo: o perfil litológico e o nível estático de um poço do SIAGAS ou de uma sondagem a percussão dizem onde estão o topo rochoso, a argila ou o nível d'água. Fixe essa profundidade e veja se o modelo ainda ajusta a curva com RMS aceitável.
• SEV paramétrica: uma SEV feita ao lado de um poço conhecido dá a resistividade de cada litologia naquela área. Ela vira a chave de interpretação das SEVs vizinhas, melhor que qualquer tabela genérica.
• Coerência entre SEVs: numa seção, prefira entre os modelos equivalentes aqueles que mantêm as camadas contínuas de uma sondagem para a outra.

As faixas de litologia que o Rhoa sugere (cortes em 10, 45, 150, 400 e 900 Ω·m) são uma convenção do app para o padrão de campo PASI Polares/UFSC. Elas ajudam na primeira leitura, mas a calibração local tem prioridade. Veja a tabela de resistividade para as faixas da literatura.

Passo a passo da inversão

A inversão é um ciclo curto de preparar os dados, inverter, criticar o resultado e testar alternativas, repetido até sobrar um modelo simples, coerente e documentado.

1. Revise os dados: plote a curva, marque outliers e confira a geometria. Os erros mais frequentes estão em erros comuns em SEV.
2. Corrija a embreagem: emende os segmentos de MN nas aberturas repetidas. Veja o guia da embreagem de MN.
3. Escolha o número de camadas pela forma da curva.
4. Monte o chute inicial a partir da curva, ou use o automático do programa.
5. Inverta e leia o RMS junto com a distribuição dos resíduos.
6. Teste L − 1 e L + 1 camadas e fique com o menor modelo que chega ao nível do ruído.
7. Teste a equivalência: varie a espessura das camadas mal resolvidas mantendo S ou T e veja quanto o RMS sobe.
8. Calibre com poços, sondagens ou outra SEV paramétrica.
9. Documente: número de camadas, RMS, modelos alternativos, faixas de equivalência e se a embreagem foi corrigida. O modelo de laudo traz essa estrutura.

Para praticar com uma curva real, use os dados de exemplo de SEV.

Como o Rhoa faz a inversão 1D

No Rhoa, a inversão ajusta log₁₀ ρ e log₁₀ h de 2 a 7 camadas por mínimos quadrados com limites, contra a curva corrigida de embreagem, com uma suavização leve entre resistividades vizinhas. O mesmo algoritmo roda no servidor (Python) e no navegador (offline).

Antes do ajuste

• Entram só as leituras com resistência medida, ρa positiva e não marcadas como outlier na caderneta. São precisas pelo menos 4.
• Se a correção de embreagem estiver ligada (padrão), os segmentos de MN são emendados. Em seguida, as leituras de mesmo AB/2 viram um único ponto pela média geométrica. N, no RMS, é o número de aberturas distintas que sobram.

O modelo direto

• Recorrência de Pekeris para a transformada de resistividade e filtro linear para a integral de Hankel (J₁) do arranjo Schlumberger, como descrito acima.
• Em SEV de Wenner Alfa, a curva Wenner sai da Schlumberger pela relação exata ρa,W(a) = 2a·∫ ρa,S(t)/t² dt, de t = a a t = 2a, integrada numericamente. Nesse arranjo não há correção de embreagem, porque o MN muda a cada leitura.
• O filtro é próprio do Rhoa: 61 coeficientes, com abscissas de −3 a +3 em passos de 0,1, ajustados por mínimos quadrados contra a solução analítica exata de 2 camadas. Não é o filtro de Guptasarma: um filtro curto de 19 coeficientes foi testado e descartado por dar ρa negativa em modelos comuns. Os detalhes estão na metodologia.

O ajuste

• Parâmetros: log₁₀ ρ de cada camada e log₁₀ h de cada camada finita. Limites: log₁₀ ρ entre −0,5 e 5 (ρ de cerca de 0,3 a 100.000 Ω·m) e log₁₀ h entre −0,7 e 2,5 (h de cerca de 0,2 a 316 m).
• Resíduo: log₁₀ ρa,obs − log₁₀ ρa,calc em cada AB/2, mais um termo 0,03 · (log₁₀ ρₖ₊₁ − log₁₀ ρₖ) para cada par de camadas vizinhas. A suavização atua só nas resistividades, não nas espessuras. Com esse peso, um contraste de uma década entre duas camadas custa o mesmo que um único ponto desajustado em cerca de 7%: o suficiente para conter oscilações, pouco para apagar contrastes reais.
• Modelo inicial: automático. As resistividades são lidas da própria curva observada em posições igualmente espaçadas ao longo dos pontos (a 1ª camada recebe o primeiro ρa, a última recebe o último) e limitadas a 1–50.000 Ω·m. As espessuras seguem uma progressão geométrica de max(0,5 m; 0,4 × menor AB/2) a max(2 m; 0,25 × maior AB/2).
• Otimizador no servidor: least_squares do SciPy, método TRF (Trust Region Reflective), com até 400 avaliações da função de resíduos.
• Otimizador no navegador: Levenberg-Marquardt com limites, Jacobiano por diferenças finitas e amortecimento adaptativo, com o mesmo teto de 400 avaliações. Filtro, resíduo, limites e modelo inicial são idênticos aos do servidor. Como o caminho até o mínimo não é o mesmo, os resultados podem diferir um pouco.
• RMS: calculado depois do ajuste, pela fórmula relativa da seção sobre RMS, só com os pontos de dados (sem o termo de suavização). Repare que o ajuste minimiza desvios em log, e o RMS informado é relativo em escala linear.
• Saída: camadas com topo, base e litologia sugerida, e uma curva sintética de 50 pontos de 0,8 × o menor a 1,2 × o maior AB/2.

O que você controla na tela de inversão

• Número de camadas (2 a 7, padrão 4), correção de embreagem MN (liga e desliga) e motor de cálculo: automático (servidor, com o navegador como reserva), só servidor ou só local.
• Modelagem direta: depois da inversão, cada camada pode ser editada por digitação, por controles deslizantes em escala log ou pelos botões ÷ 1,1 e × 1,1. Dá para acrescentar e remover camadas e trocar a litologia. A curva sintética, os degraus do modelo, o RMS e a interpretação são recalculados na hora. “Descartar” volta ao resultado da inversão.
• O que não dá para fazer: fornecer o seu próprio modelo inicial ao otimizador ou fixar um parâmetro durante a inversão. Para calibrar com um poço, use a modelagem direta: ponha a interface na profundidade do poço, ajuste as outras camadas e veja se o RMS continua aceitável.
• Comparar e guardar: “Salvar modelo” grava o resultado (inclusive o editado) no histórico de inversões, e o gráfico pode sobrepor o modelo salvo anterior. O app também mostra o tipo de curva, os parâmetros de Dar Zarrouk e o nível d'água estimado do modelo em exibição.

A tabela de parcimônia e os exemplos de equivalência e supressão deste guia foram calculados com esse mesmo motor. Para ver a tela de inversão funcionando, abra a demonstração do Rhoa. Para comparar com o fluxo de outro programa, veja Rhoa vs. IPI2Win.

Referências

• BRAGA, A. C. O. Geofísica aplicada: métodos geoelétricos em hidrogeologia. São Paulo: Oficina de Textos, 2016.
• BRANCH, M. A.; COLEMAN, T. F.; LI, Y. A subspace, interior, and conjugate gradient method for large-scale bound-constrained minimization problems. SIAM Journal on Scientific Computing, v. 21, n. 1, p. 1-23, 1999.
• CONSTABLE, S. C.; PARKER, R. L.; CONSTABLE, C. G. Occam's inversion: a practical algorithm for generating smooth models from electromagnetic sounding data. Geophysics, v. 52, n. 3, p. 289-300, 1987.
• GHOSH, D. P. The application of linear filter theory to the direct interpretation of geoelectrical resistivity sounding measurements. Geophysical Prospecting, v. 19, p. 192-217, 1971.
• GUPTASARMA, D. Optimization of short digital linear filters for increased accuracy. Geophysical Prospecting, v. 30, 1982.
• INMAN, J. R. Resistivity inversion with ridge regression. Geophysics, v. 40, n. 5, p. 798-817, 1975.
• KOEFOED, O. Geosounding principles, 1: resistivity sounding measurements. Amsterdam: Elsevier, 1979.
• MARQUARDT, D. W. An algorithm for least-squares estimation of nonlinear parameters. Journal of the Society for Industrial and Applied Mathematics, v. 11, n. 2, p. 431-441, 1963.
• ORELLANA, E. Prospección geoeléctrica en corriente continua. Madrid: Paraninfo, 1972.
• TELFORD, W. M.; GELDART, L. P.; SHERIFF, R. E. Applied geophysics. 2. ed. Cambridge: Cambridge University Press, 1990.
• ZOHDY, A. A. R. A new method for the automatic interpretation of Schlumberger and Wenner sounding curves. Geophysics, v. 54, n. 2, p. 245-253, 1989.