MAGY V1.0 – PROGRAMA DE PROCESSAMENTO DE DADOS MAGNÉTICOS S. NEVES Centro de Geofísica de Évora, Laboratório HERCULES, Universidade de Évora, Rua Romão Ramalho, 59 Évora, Portugal,
[email protected]
O objectivo deste trabalho é a apresentação de um software que permite o processamento e representação de mapas de anomalias magnéticas (intensidade e gradiente do campo magnético total terrestre), em que os dados são importados directamente do magnetómetro. O programa foi desenvolvido através da linguagem Python, sendo utilizadas bibliotecas específicas de processamento numérico, processamento de sinal e de representação gráfica. É apresentado um caso sintético com dois diques verticais magnéticos paralelos entre si, tendo sido adicionado ruído. Este caso foi processado com os métodos de processamento de sinal implementados no programa.
1.
Introdução
A prospeção magnética é uma das técnicas geofísicas mais utilizadas no mapeamento de anomalias em escalas geográficas variáveis e com o objectivo de identificar estruturas de interesse geológico, económico, ambiental e arqueológico. Esta técnica consiste no mapeamento de um ou mais parâmetros do campo magnético terrestre, de modo a analisar as anomalias magnéticas e relacioná-las com as estruturas sub-superficiais. A interpretação dos mapas de anomalias magnéticas consiste na localização geoespacial da fonte, a sua caracterização geométrica, e a determinação das propriedades magnéticas dos materiais que a constituem (sempre que possível, acompanhada da identificação dos materiais geológicos que a compõem). As linguagens de programação de alto nível, como é o caso do Python [2] e Matlab [3], cada vez mais têm-se tornado ferramentas acessíveis aos utilizadores sem grandes conhecimentos de programação. Estas linguagens de programação permitem criar programas de baixa a elevada complexidade, depende unicamente do objetivo do programa, bem como do tempo de cálculo exigido e das capacidades de cálculo do computador. O MAGY foi desenvolvido através da linguagem Python, tendo sido utilizadas bibliotecas numéricas que permitem otimizar os recursos computacionais (numpy [6]) bem como criar as interfaces gráficas funcionais e apelativas (matplotlib [7]). A estrutura do programa é apresentada na Figura 1, sendo constituída por pacotes de funções (input, projecto, processamento, representação, análise, e output). Estes pacotes de funções são descritos no parágrafo seguinte.
Input
Representação Projeto Dados originais Dados
Processamento
Análise Output
Figura 1. Estrutura do programa MAGY.
2.
Input
A função input permite ao utilizador escolher o formato e importar os dados que pretende bem como definir alguns parâmetros de entrada (Figura 2). Nomeadamente, o nome do projeto, o comprimento da linha (L), o afastamento entre linhas (dy), a posição do inicio da prospeção, e definir se o operador do magnetómetro realizou o percurso em zig-zag ou não. Estes parâmetros de entrada permitem ao programa corrigir as posições das amostras, dado que as mesmas são obtidas pelo GPS do magnetómetro e por isso, sujeitas a algumas imprecisões de localização. Após a leitura dos dados e definição de todos os parâmetros, o utilizar pode gravar o projecto ou sair desta opção. O projeto é gravado num ficheiro em formato próprio (.bin). Depois de gravar o projeto, o programa coloca automaticamente os dados em memória, permitindo ao utilizador continuar o seu trabalho. A versão MAGY 1.0 só permite ler dados provenientes de magnetómetros da marca Gem Systems, com dois sensores, e com dados obtidos no campo segundo uma grelha regular (quadrática ou rectangular).
Figura 22. Janela de criaçãão de projecto.
3.
P rojeeto
O projeto é um ficheiiro que contéém uma basee de dados com c toda a informação paramétricaa relativa ao ensaio e magnettométrico (geo ometria, topog grafia, parâmeetros físicos, número de amostras). O projjeto é um ficcheiro (.bin) composto po or dois conten ntores de infoormação. O primeiro coontentor contéém os dados ooriginais, tal co omo são impo ortados do maggnetómetro. Estes dadoos são armazzenados e nuunca são alteerados (dadoss estáticos). O segundo contentor é composto po or um conjunnto de macross (sub-contenttores) que arm mazenam os dados que poderão serr alvo de proocessamento (dados dinâm micos). Estes dados são atualizáveiss em qualquerr momento. 4.
Repreesentação
A represenntação dos dad dos é realizadda com recurso a bibliotecaa Matplotlib [77], aqual se encontra naa base da reprresentação doos diversos grráficos. O programa permitte visualizar os dados em e gráfico de traços ou iimagem. Paraa além das possibilidades p s anteriores, também peermite visualiizar a posiçãoo das amostrras no plano. As Figuras 3 - a) e b) apresentam m exemplos do os gráficos de traços, e de im magem, respeectivamente. O gráfico de imagem uttiliza uma maalha regular dde 0,25 x 0,2 25 m que resu ulta da aplicaação de um polinómio interpolador de quinto grrau. É possív vel definir o intervalo de valores da intensidadee e do gradien nte total do ccampo magnéético terrestre, bem como a paleta de cores do grráfico de imag gem.
a)
b)
Figura 3. a) Grááfico de traços, b) b Gráfico de imag gem.
5.
Proceessamento
O processaamento dos dados d geofísiicos é fundam mental para obter resultaddos de boa qualidade e permitir boaas interpretaçções. Deste modo, m este paccote tem um cconjunto de funções quue permitem ao o utilizador m melhorar a quaalidade do resultado, bem ccomo ajudar a compreennder e interpreetar os mesmoos. A funçção despiking g permite elim minar ondas dee pequeno comprimento dee onda com elevada am mplitude relativ vamente à méédia (em geral, resultado daa existência dde artefactos metálicos no n terreno). Este E tipo de aanomalias im mpossibilita a correta visuaalização dos dados maggnéticos, visto o que por vezzes os dados que se preten ndem analisarr estão num intervalo dee valores muitto inferiores aao intervalo dee valores global. O MAGY tem dispo oníveis dois m métodos para eliminar e ou red duzir o efeito de spiking. O método de d z-score (eliimina os dadoos anómalos), e o método lo ogarítmico (redduz os spikes seguundo um razão o logarítmica) . O método z--score [4] é um ma ferramentaa estatística que permitee avaliar se ass series de daddos têm outlierrs, ou não, baseando na Eq.. (1).
z i xi x / s
Eq. (1)
Onde, xi é a amostrra, x é o valoor médio do conjunto c das amostras, e S é o desvio padrão do mesmo m conjun nto. Se um daada amostra fo or outlier (zi > 3.5) o prograama permite a substituiçção pela média ou pela meddiana da seriee de dados, caaso contrário, as amostras permanecem m com o mesm mo valor. O méttodo logarítmiico [5] realizaa uma suavizzação das gran ndes amplituddes de onda através da Eq. (2), sendo o x o valor daa amostra, e A o fator de su uavização. Quuanto maior for o fator de redução (A) ( maior serrá a suavizaçãão imposta ao os dados. Desste modo, é possível lim mitar a variaação da granddeza represen ntada a um intervalo de vvalores que possibilita a visualização o adequada daas principais anomalias mag gnéticas.
x x mod A. log10 1 A
Eq. (2)
O destripe é um método que reduz o efeito visual de bandeamento do mapa de anomalias magnéticas. Este efeito deve-se ao facto de existirem normalmente variações do valor médio das series individuais contidas no mapa de anomalias. (O mapa de anomalias magnéticas é obtido através de um conjunto de series de dados individuais, paralelas entre si, sendo que as mesmas são obtidas segundo um percurso em linha no campo). Para reduzir este efeito, este método calcula o valor médio do gradiente magnético para cada serie, e realiza uma translação vertical em sinal contrário ao valor médio obtido. Logo, todas as serieis são normalizadas com o valor médio nulo. O desttager é uma função que permite reduzir o efeito de desfasamento (as series impares estão deslocadas de uma determinada distância das series pares, produzindo o efeito de desfasamento). Este efeito acontece devido ao modo de recolha de dados em campo (normalmente em zig-zag), incorreto posicionamento das amostras (obtido por GPS), erros introduzidos pelo operador do magnetómetro ou devido à topografia do terreno. Foi desenvolvido o algoritmo que permite ao utilizador realizar a translação de series par ou impar, mediante uma determinada distância e sentido, sendo o sentido dado pelo sinal do valor da distância. Após a definição de todos os parâmetros, o programa realiza a operação de destagger, ou seja, vai procurar todas as series individuais ímpares ou pares, e realiza a operação de translação, com a distância e sentido definido pelo utilizador. A translação é realizada no mesmo eixo onde estão dispostas as amostras magnéticas. O filtro IIR (Infinite Impulse Response) foi obtido através da biblioteca Scipy [8]. O filtro utiliza o método de filtragem digital butterword, permitindo ao utilizador diversas opções de filtragem: passa baixo, passa alto e passa banda, e atribuir o número de onda a filtrar. É possível escolher o número de polos a utilizar no processo de filtragem. O processo de filtragem é realizado nos dois sentidos de modo a anular o efeito de translação de fase de onda. Este método pode ser utilizado várias vezes sem o incómodo de sair desta opção e voltar a entrar nela. O filtro Median é um método que permite eliminar o ruído magnético de fundo. Este filtro altera o valor de cada amostra magnética pela mediana, sendo esta determinada através de um conjunto de dados contidos num intervalo definido pelo utilizador (janela do filtro). Contudo, este filtro deve ser utilizado com cuidado, pois a sua incorreta utilização pode eliminar anomalias magnéticas importantes do local de estudo [10,11]. O Esta função permite ao utilizador escolher unicamente a dimensão da janela do filtro. Este filtro foi obtido através da biblioteca Scipy [8].
6.
Análise espectral
Esta função permite calcular o espectro de todas as séries individualmente, recorrendo à transformada rápida de Fourier (FFT). O método numérico foi obtido através da biblioteca Scipy [8], sendo depois representada num gráfico próprio pela ferramenta gráfica Matplotlib [7]. A FFT é determinada segundo a Eq. (3); N 1
xk xn e
i 2k
n N
,k 0,..., N 1
Eq. (3)
n 0
Onde, xn é o conjunto de dados originais com uma dimensão N (numero de amostras). Ao conjunto inicial de dados são adicionados zeros com a mesma dimensão do conjunto inicial, de modo a melhorar a visualização das amplitudes com baixo número de onda. A análise espectral permite ao utilizador avaliar as amplitudes do espectro em função do número de onda. Com base na informação apresentada por esta análise o utilizador pode aplicar os diversos métodos de processamento disponíveis e avaliar a sua evolução espectral. 7.
Output
Esta função permite ao utilizar exportar os resultados que pretende. O utilizador pode exportar os dados em formato de imagem, sendo que esta opção está sempre disponível na função representação (parágrafo 4). Também foi desenvolvido o algoritmo que permite ao utilizar exportar o resultado em ficheiro ASCII, sendo posteriormente utilizável noutra aplicação informática. 8.
Aplicação do MAGY a dados sintéticos
Foi cálculado um modelo sintético composto por dois diques magnéticos paralelos (verticalmente polarizados). Ambos os diques têm a mesma geometria, susceptibilidade magnética, e estão posicionados à mesma profundidade, entre os 3,0 e os 6,0 m (Figura 4 - b). Tratando-se de corpos simples, o campo magnético originado pelos diques foi calculado através da Eq. (4) [9]; z z H s 2 It 2 1 2 2 2 2 Eq. (4) z x z 2 x 1 Onde I é a susceptibilidade do corpo magnético, t é a largura do dique, z1 é a distância entre a superfície do terreno e o topo do dique, z2 é a distância entre a superfície do terreno e a base do dique, x é a posição da amostra magnética modelada (Hs). Os corpos magnéticos são constituídos por 12% de magnetite (k=0.5). Normalmente, na prospeção arqueológica é utilizado o gradiente magnético vertical, por isso, foi calculada a componente vertical do campo magnético terrestre para as cotas z=0 e z=0.8 m. Para a cota z = 0 m, os valores utilizados na Equação 8.1 foram; z1=3,0 m, z2= 6,0 m,
t=0,5 m, e I 0,5 0,12 0,6 0,0336 cgs. Os vaalores utilizad dos para a coota z=0,8 m foram; z1=33,8 m, z2= 6,8 m, t=0,5 m, e I 0,5 0,12 0,6 0,03 36 cgs. Posteriiormente foi calculado o valor do grad diente verticaal com base nnos valores anteriores dividido d pela distância entree eles. A quad drícula ao níveel do solo é dee 20 x 20 m (Figura 4-aa), sendo que a disposição das amostrass começou do o canto inferioor esquerdo para o direiito, e as linhas de baixo parra cima (as lin nhas estão afaastadas entre ssi de 1,0 m). As linhas foram f todas ob btidas no messmo sentido, motivo m pelo qu ual o modo ziig-zag não é activado.
a)
b)
Figuraa 4. a) Vista de top po da geometria ddos corpos magnéticos sintéticos; b) Corte transverrsal AA’. [m]
Foi addicionado aos dados sintétiicos originais ruído numériico aleatório, 4 spikes, e variação aleatória a do valor médiio do gradieente vertical magnético (efeito de bandeamennto). A introdu ução deste tip o de artefacto os numéricos permite p simullar situações reais de aquuisição de dad dos em campoo. Esta inform mação foi gravada num ficheeiro texto, e depois foi importada parra o programaa MAGY. Ap pós criar o pro ojecto, os daddos originais foram repreesentados com m a função imaagem (Figura5 5).
Spikes
Figura 55. Imagem dos daados originais.
Verificcou-se que attravés da Figuura 5 os dados contêm muito m ruído. D Deste modo, apresenta-sse um procedimento tipo dee processamen nto de dados magnéticos. m Innicialmente, foi aplicadda a função despiking. d Fooram utilizado os os dois métodos m que o programa permite (z--score e logaritmo). Verificcou-se que o método m z-scorre eliminou ccom sucesso todos os sppikes (Figura 6 – b). O mesm mo não aconteceu com o método m logarítm tmico, tendo os spikes continuando c na n imagem. O motivo pelo qual o último método nãoo funcionou correctameente, deve-se ao a facto que a amplitude do os spikes e dos restantes dad ados não são muito diferrentes. Logo, quando o algooritmo é apliccado, os dadoss são suavizaddos todos de igual formaa, resultando o aparecimentto dos spikes.
a)
b))
Figgura 6. a) Aplicad do o método logaarítmico, b) Apliccado o método z-sscore (valor médiio).
Foi esccolhido o méttodo z-score ppara dar contiinuidade ao processamentoo dos dados. De seguidaa foi aplicado o o método deestripe de mo odo a eliminarr o bandeado (Figura 7). Comparanddo a Figura 6-b com a Figura 7 verrifica-se que os dados esstão menos perturbados.
Figura 7. Apliccação do método destripe aos dado os processados co om o método z-sccore.
Neste exemplo não é utilizado o método desttagger, dado que os dadoss não têm o efeito de traanslação entree series. Para melhorar m o ressultado final ffoi aplicado o filtro IIR, no o domínio doo número de onda. Foi aplicado um filtro passa banda entre 0.5 e 2 ciclo os/m com doiis polos. A utilização desses d valoress deve-se a eliiminação dos grandes e peq quenos númerros de onda, visto que am mbos estão reelacionados coom o ruído num mérico aleatórrio. Atravéés da Figura 8 verifica-see que existe um melhoraamento da quualidade de imagem, seendo as princcipais anomallias magnéticaas reforças. Verificou-se V qque o efeito oscilatório das anomaliias magnéticaas entre linh has (devido ao a adicionado do de ruído numérico) foi reduzido.
Figura 8. 8 Aplicação do ffiltro IIR e mapa de anomalias magnéticas final.
Este programa p estáá em fase de desenvolvim mento, motivo pelo qual teem algumas limitações ao nível do desempenho d coom um númerro elevado dee amostras, quue por vezes poderá bloqquear o compu utador. Pretennde-se melhorrar os actuai s algoritmos bem como introduzir i méétodos mais robustos. O pressente program ma está vocaacionado paraa utilizar o modo: m gradiennte vertical magnético, só sendo po ossível visualiizar os dadoss do valor tottal do campoo magnético terrestre (não ( é possíível aplicar qualquer tip po de proceessamento). Pretende-se disponibilizzar ferramentas que permittam ao utilizaador processarr os dados dee valor total do campo magnético m bem m como a sua interpretação.
9.
Conclusões
Este trabalho tem como objectivo a apresentação do software (MAGY), que permite realizar o processamento e representação de mapas de anomalias magnéticas, intensidade e gradiente magnético do campo terrestre. Foram abordados os métodos de processamento que estão implementados no MAGY bem como a sua aplicação num caso sintético. O caso sintético é constituído por dois diques verticais magnéticos (polarizados verticalmente). Ao campo à superfície produzido por estes corpos foi adicionado ruído numérico aleatório, 4 spikes e oscilação aleatória do valor médio de cada serie (efeito de bandeamento). Verificou-se que o despiking com o método z-score foi bem-sucedido, tendo-se eliminado os spikes que se encontravam nos dados. A aplicação do método logarítmico não produziu não foi tão eficiente, motivo pelo qual não foi escolhido para continuar o processo de tratamento de dados. O método destripe também foi aplicado com sucesso, tendo reduzido o efeito do bandeado, melhorando a qualidade do mapa de anomalias magnéticas. O método destagger não foi aplicado visto que os dados não contêm o efeito de translação entre linhas consecutivas. Foi aplicado o filtro IIR de modo a melhorar a qualidade final do mapa de anomalias magnéticas, tendo-se verificado uma diminuição das oscilações do gradiente magnético e por consequência um melhoramento substancial do mapa de anomalias magnéticas. O programa tem algumas limitações, ao nível do desempenho com elevado número de amostras (que poderá provocar o bloqueio do computador) bem como o número reduzido de funções de processamento de dados (despiking, destripe, desttager, filtro IIR e filtro median), motivo pelo qual o programa será alvo de melhoramento, nomeadamente ao nível dos actuais algoritmos bem como na implementação de métodos de processamento mais robustos. 10. Agradecimentos Este trabalho foi co-financiado pelo Fundo Europeu de Desenvolvimento Regional (FEDER), através do programa INALENTEJO, no âmbito do projecto IMAGOS – Inovative Methodologies in Archaeology, Archaeometry and Geophysics – Optimizing Strategies X APOLLO – Archaeological and Physical On-site Laboratory – Lifting Outputs (ALENT-07-0224-FEDER-001760), do Laboratório HERCULES (Universidade de Évora). 11. Referências 1.
2. 3.
R. E. Chavez, D. L. Argote, M. E. Camara, Geophysical characterization of a rural archaeological site in central mexico. Proc. EAGE 66th conference & exhibition, Paris, France, (2004). Python, https://www.python.org/, acedido a 5-1-2015. Matlab, http://www.mathworks.com/products/matlab, acedido a 5-1-2015
4.
B. Iglewicz and D. Hoaglin, How to Detect and Handle Outliers, The ASQC Basic References in Quality Control: Statistical Techniques. PhD Edward F. Mykytka, Volume 16, (1993). 5. A. Mojica, L. Pastor, C. Camerlynck, N. Florsch, A. Tabbagh, Magnetic prospection of the pre-columbian archaeological site of El Caño in the cultural regrion of Gran Coclé, Panama. Archaeological Prospection, 11 p, 2014. 6. Numpy, http://www.numpy.org, acedido a 8-1-2015. 7. Matplotlib, http://matplotlib.org, acedido a 8-1-205. 8. Scipy, http://www.scipy.org, acedido a 8-1-2015. 9. M. B. Dobrin, Introduction to geophysical prospecting, third edition, McGrawHill book company, pp. 500-501, (1981). 10. M. Ciminale, M. Loddo, Aspects of magnetic data processing, Archaeological Prospectiong, Vol. 8, pp 239-246, 2001 11. A. Eder-Hinterleitner, W. Neubauer, P. Melichar, Restoring magnetic anomalies, Archaeological Prospecting, Vol. 3, pp. 185 – 197, 1996.
