Obrigado ao autor do blog Anton Ekimenko por sua palestra1. Introdução
Esta nota foi escrita na sequência da
Wolfram Russian Technology Conference e contém uma sinopse do relatório que apresentei. O evento ocorreu em junho na cidade de São Petersburgo. Dado o fato de trabalhar a uma quadra do local da conferência, não pude deixar de participar deste evento. Em 2016 e 2017, ouvi os relatórios da conferência e, neste ano, fiz uma apresentação. Primeiro, apareceu um tópico interessante (como me parece) que estamos desenvolvendo com
Kirill Belov e, depois, após um longo estudo da legislação da Federação Russa em relação à política de sanções, duas licenças da
Wolfram Mathematica apareceram na empresa em que trabalho.
Antes de passar para o tópico da minha apresentação, gostaria de observar a boa organização do evento. A página da conferência usa a imagem da Catedral de Kazan. A catedral é uma das principais atrações de São Petersburgo e é claramente visível do salão em que a conferência foi realizada.
Na entrada da Universidade Estadual de Economia de São Petersburgo, os participantes foram recebidos por assistentes estudantis - eles não permitiram que se perdessem. Durante o registro, foram entregues pequenas lembranças (um adesivo de brinquedo, uma caneta, adesivos com os símbolos Wolfram). Os intervalos para almoço e café também foram incluídos na programação da conferência. Sobre deliciosos cafés e tortas, eu já notei na parede do grupo - chefs bem feitos. Com esta parte introdutória, gostaria de enfatizar que o próprio evento, seu formato e local já trazem emoções positivas.
O relatório, preparado por mim e Kirill Belov, é chamado de “Usando o Wolfram Mathematica para resolver problemas de geofísica aplicada. Análise espectral de dados sísmicos ou "onde os rios antigos corriam". O conteúdo do relatório abrange duas partes: primeiro, o uso de algoritmos disponíveis no
Wolfram Mathematica para a análise de dados geofísicos e, em segundo lugar, é assim que você coloca os dados geofísicos no Wolfram Mathematica.
Exploração sísmica
Primeiro você precisa fazer uma pequena excursão na geofísica. Geofísica é uma ciência que estuda as propriedades físicas das rochas. Bem, como as rochas têm propriedades diferentes: elétrica, magnética, elástica, existem métodos correspondentes de geofísica: exploração elétrica, exploração magnética, exploração sísmica ... No contexto deste artigo, tocaremos apenas na exploração sísmica com mais detalhes. A exploração sísmica é o principal método para encontrar petróleo e gás. O método baseia-se na excitação de vibrações elásticas e no registro subsequente da resposta das rochas que compõem a área de estudo. A excitação de vibrações é realizada em terra (dinamite ou fontes de vibração não explosivas de vibrações elásticas) ou no mar (pistolas pneumáticas). As vibrações elásticas se propagam através da espessura das rochas, sendo refratadas e refletidas nos limites das camadas com propriedades diferentes. As ondas refletidas retornam à superfície e são registradas por geofones terrestres (geralmente dispositivos eletrodinâmicos com base no movimento de um ímã suspenso em uma bobina) ou hidrofones no mar (com base no efeito piezoelétrico). No momento da chegada das ondas, é possível julgar as profundezas das camadas geológicas.
O navio sísmico reboca o equipamento:
Pistola de ar excita vibrações elásticas:
As ondas passam através da massa rochosa e são registradas pelos hidrofones:
Navio de pesquisa geofísica "Ivan Gubkin" no cais na ponte Blagoveshchensky em São Petersburgo:
Modelo de sinal sísmico
Rochas têm propriedades físicas diferentes. Para a exploração sísmica, as propriedades elásticas são mais importantes - a velocidade de propagação das vibrações e da densidade elásticas. Se duas camadas tiverem propriedades iguais ou semelhantes, a onda não "notará" o limite entre elas. Se as velocidades das ondas nas camadas diferirem, a reflexão ocorrerá no limite das camadas. Quanto maior a diferença nas propriedades, mais intensa a reflexão. Sua intensidade será determinada pelo coeficiente de reflexão (rc):
onde ρ é a densidade da rocha, ν é a velocidade da onda, 1 e 2 denotam as camadas superior e inferior.Um dos modelos de sinal sísmico mais simples e mais comumente usado é o modelo de convolução, quando um traço sísmico registrado é apresentado como resultado da convolução de uma sequência de coeficientes de reflexão com um pulso de sonda:
onde s (t) é a faixa sísmica, ou seja, tudo o que um hidrofone ou um geofone gravou durante um tempo fixo de gravação,
w (t) é o sinal que a pistola de ar gera,
n (t) é ruído aleatório.
Calculamos, por exemplo, uma pesquisa sísmica sintética. Usaremos o pulso Ricker amplamente utilizado na exploração sísmica como o sinal inicial.
length=0.050; (*Signal lenght*) dt=0.001;(*Sample rate of signal*) t=Range[-length/2,(length)/2,dt];(*Signal time*) f=35;(*Central frequency*) wavelet=(1.0-2.0*(Pi^2)*(f^2)*(t^2))*Exp[-(Pi^2)*(f^2)*(t^2)]; ListLinePlot[wavelet, Frame->True,PlotRange->Full,Filling->Axis,PlotStyle->Black, PlotLabel->Style["Initial wavelet",Black,20], LabelStyle->Directive[Black,Italic], FillingStyle->{White,Black},ImageSize->Large,InterpolationOrder->2]
Fonte pulso sísmico
Definimos dois limites em profundidades de 300 ms e 600 ms, e os coeficientes de reflexão serão números aleatórios.
rcExample=ConstantArray[0,1000]; rcExample[[300]]=RandomReal[{-1,0}]; rcExample[[600]]=RandomReal[{0,1}]; ListPlot[rcExample,Filling->0,Frame->True,Axes->False,PlotStyle->Black, PlotLabel->Style["Reflection Coefficients",Black,20], LabelStyle->Directive[Black,Italic]]
A sequência dos coeficientes de reflexão:
Calcule e exiba a trilha sísmica. Como os coeficientes de reflexão possuem sinais diferentes, obtemos duas reflexões alternadas no caminho sísmico.
traceExamle=ListConvolve[wavelet[[1;;;;1]],rcExample]; ListPlot[traceExamle, PlotStyle->Black,Filling->0,Frame->True,Axes->False, PlotLabel->Style["Seismic trace",Black,20], LabelStyle->Directive[Black,Italic]]
Trilha simulada:
Para este exemplo, uma reserva deve ser feita - na realidade, a profundidade das camadas é determinada, obviamente, em metros, e o cálculo do caminho sísmico ocorre para o domínio do tempo. Seria mais correto definir as profundidades em metros e calcular os tempos de chegada sabendo a velocidade nos estratos. Nesse caso, eu imediatamente defino as camadas no eixo do tempo.Se falamos de estudos de campo, como resultado de tais observações, um grande número de séries temporais (trilhas sísmicas) é registrado. Por exemplo, ao explorar um site com 25 km de comprimento e 15 km de largura, onde, como resultado do trabalho, cada trilha caracteriza uma célula de 25x25 metros (essa célula é chamada bin), o conjunto final de dados conterá 600.000 trilhas. Com uma etapa de amostragem de tempo de 1 ms, um tempo de gravação de 5 segundos, o arquivo final de dados será superior a 11 GB e a quantidade de material "bruto" original pode ser de centenas de gigabytes.
Como trabalhar com eles no
Wolfram Mathematica ?
O desenvolvimento do pacote começou com uma
pergunta no muro VK do grupo de suporte de língua russa. Graças às respostas da comunidade, uma solução foi encontrada muito rapidamente. E, como resultado, tornou-se um desenvolvimento sério. O
post correspondente
no mural da Wolfram Comunity foi anotado pelos moderadores. Atualmente, o pacote oferece suporte ao trabalho com os seguintes tipos de dados usados ativamente no setor geológico:
- importar dados de mapa do formato ZMAP e IRAP
- importação de medições em poços no formato LAS
- entrada e saída de arquivos sísmicos no formato SEGY
Para instalar o pacote, você deve seguir as instruções na página de download do pacote montado, ou seja, execute o seguinte código em qualquer
notebook Mathematica :
If[PacletInformation["GeologyIO"] === {}, PacletInstall[URLDownload[ "https://wolfr.am/FiQ5oFih", FileNameJoin[{CreateDirectory[], "GeologyIO-0.2.2.paclet"}] ]]]
Depois disso, o pacote será instalado na pasta padrão, cujo caminho pode ser obtido da seguinte maneira:
FileNameJoin[{$UserBasePacletsDirectory, "Repository"}]
Como exemplo, demonstramos os principais recursos do pacote. A chamada é tradicionalmente para pacotes no Wolfram Language:
Get["GeologyIO`"]
O pacote está sendo desenvolvido usando o
Wolfram Workbench . Isso permite acompanhar a principal funcionalidade do pacote com a documentação, que no formato de apresentação não difere da documentação do próprio Wolfram Mathematica e fornece ao pacote arquivos de teste para o primeiro conhecido.
Esse arquivo, em particular, é o arquivo Marmousi.segy, um modelo sintético de uma seção geológica desenvolvida pelo Instituto Francês de Petróleo. Usando esse modelo, os desenvolvedores testam seus próprios algoritmos para modelar o campo de onda, processamento de dados, inversão de traços sísmicos etc. O próprio modelo Marmousi é armazenado no repositório, de onde o próprio pacote foi baixado. Para obter o arquivo, execute o seguinte código:
If[Not[FileExistsQ["Marmousi.segy"]], URLDownload["https://wolfr.am/FiQGh7rk", "Marmousi.segy"];] marmousi = SEGYImport["Marmousi.segy"]
O resultado da importação é um objeto SEGYData:
O formato SEGY envolve o armazenamento de várias informações sobre as observações. Em primeiro lugar, estes são comentários de texto. Informações sobre o local de trabalho, os nomes das empresas que realizaram as medições etc. são inseridas aqui. No nosso caso, esse cabeçalho é chamado por uma solicitação com a chave TextHeader. Aqui está um cabeçalho de texto abreviado:
Short[marmousi["TextHeader"]]
“O conjunto de dados de Marmousi foi gerado no Instituto ... velocidade máxima de 1500 m / se máxima de 5500 m / s)”
O modelo geológico real pode ser exibido entrando em contato com traços sísmicos usando a tecla "traços" (um dos recursos do pacote é a distinção entre maiúsculas e minúsculas das chaves):
ArrayPlot[Transpose[marmousi["traces"]], PlotTheme -> "Detailed"]
Modelo Marmousi:
No momento, o pacote também permite baixar dados em partes de arquivos grandes, para que seja possível processar arquivos cujo tamanho pode atingir dezenas de gigabytes. O pacote também inclui funções para exportar dados para .segy e reescrever parcialmente para o final do arquivo.
Separadamente, vale a pena notar a funcionalidade do pacote ao trabalhar com a estrutura complexa dos arquivos .segy. Uma vez que permite não apenas acessar chaves e índices para rastreamentos individuais, cabeçalhos, mas também alterá-los com a gravação subsequente em um arquivo. Muitos dos detalhes técnicos da implementação do GeologyIO estão fora do escopo deste artigo e provavelmente merecem uma descrição separada.
A relevância da análise espectral na exploração sísmica
A capacidade de importar materiais sísmicos para o Wolfram Mathematica permite usar a funcionalidade de processamento de sinal incorporada para dados experimentais. Como cada trilha sísmica é uma série temporal, uma das principais ferramentas para estudá-las é a análise espectral. Entre os pré-requisitos para a análise da composição de frequência dos dados sísmicos estão, por exemplo, os seguintes:
- Diferentes tipos de ondas são caracterizados por diferentes composições de frequência. Isso nos permite selecionar ondas úteis e suprimir ondas de interferência.
- Propriedades de rochas como porosidade e saturação podem afetar a composição da frequência. Isso permite que você selecione rochas com melhores propriedades.
- Camadas com diferentes espessuras causam anomalias em diferentes faixas de frequência.
O terceiro parágrafo é central no contexto deste artigo. Abaixo está um trecho de código para calcular traços sísmicos no caso de uma camada com espessura variável - um modelo de cunha. Este modelo é tradicionalmente estudado em exploração sísmica para a análise de efeitos de interferência, quando ondas refletidas de várias camadas se sobrepõem.
nx=200;(* Number of grid points in X direction*) ny=200;(* Number of grid points in Y direction*) T=2;(*Total propagation time*) (*Velocity and density*) modellv=Table[4000,{i,1,ny},{j,1,nx}];(* P-wave velocity in m/s*) rho=Table[2200,{i,1,ny},{j,1,nx}];(* Density in g/cm^3, used constant density*) Table[modellv[[150-Round[i*0.5];;,i]]=4500;,{i,1,200}]; Table[modellv[[;;70,i]]=4500;,{i,1,200}]; (*Plotting model*) MatrixPlot[modellv,PlotLabel->Style["Model of layer",Black,20], LabelStyle->Directive[Black,Italic]]
Modelo de formação de explosão:
A velocidade da onda dentro da cunha é de 4500 m / s, fora da cunha é de 4000 m / s, e a densidade é assumida como sendo constante de 2200 g / cm³. Para esse modelo, calculamos os coeficientes de reflexão e traços sísmicos.
rc=Table[N[(modellv[[All,i]]-PadLeft[modellv[[All,i]],201,4000][[1;;200]])/(modellv[[All,i]]+PadLeft[modellv[[All,i]],201,4500][[1;;200]])],{i,1,200}]; traces=Table[ListConvolve[wavelet[[1;;;;1]],rc[[i]]],{i,1,200}]; starttrace=10; endtrace=200; steptrace=10; trasenum=Range[starttrace,endtrace,steptrace]; traserenum=Range[Length@trasenum]; tracedist=0.5; Rotate[Show[ Reverse[Table[ ListLinePlot[traces[[trasenum[[i]]]]*50+trasenum[[i]]*tracedist,Filling->{1->{trasenum[[i]]*tracedist,{RGBColor[0.97,0.93,0.68],Black}}},PlotStyle->Directive[Gray,Thin],PlotRange->Full,InterpolationOrder->2,Axes->False,Background->RGBColor[0.97,0.93,0.68]], {i,1,Length@trasenum}]],ListLinePlot[Transpose[{ConstantArray[45,80],Range[80]}],PlotStyle->Red],PlotRange->All,Frame->True],270Degree]
Trilhas sísmicas para o modelo de cunha:
A sequência de traços sísmicos mostrada nesta figura é chamada de seção sísmica. Como você pode ver, sua interpretação pode ser realizada em um nível intuitivo, pois a geometria das ondas refletidas corresponde exclusivamente ao modelo que foi definido anteriormente. Se você analisar as faixas com mais detalhes, poderá ver que as faixas do 1º ao 30º não diferem - o reflexo do topo da formação e da sola não se sobrepõe. A partir da rota 31, as reflexões começam a interferir. E, embora, no modelo, os coeficientes de reflexão não sejam alterados horizontalmente - as trilhas sísmicas alteram sua intensidade com uma alteração na espessura do reservatório.
Considere a amplitude de reflexão a partir do limite superior do reservatório. A partir da 60ª rota, a intensidade da reflexão começa a aumentar e, na 70ª rota, torna-se máxima. É assim que a interferência das ondas do telhado e da sola nas formações se manifesta, levando em alguns casos a anomalias significativas do registro sísmico.
ListLinePlot[GaussianFilter[Abs[traces[[All,46]]],3][[;;;;2]], InterpolationOrder->2,Frame->True,PlotStyle->Black, PlotLabel->Style["Amplitude of reflection",Black,20], LabelStyle->Directive[Black,Italic], PlotRange->All]
Gráfico da amplitude da onda refletida a partir da borda superior da cunha
É lógico que quando o sinal é de frequência mais baixa, a interferência começa a aparecer em grandes espessuras da formação e, no caso de um sinal de alta frequência, a interferência ocorre em espessuras menores. O fragmento de código a seguir cria um sinal com frequências de 35 Hz, 55 Hz e 85 Hz.
waveletSet=Table[(1.0-2.0*(Pi^2)*(f^2)*(t^2))*Exp[-(Pi^2)*(f^2)*(t^2)], {f,{35,55,85}}]; ListLinePlot[waveletSet,PlotRange->Full,PlotStyle->Black,Frame->True, PlotLabel->Style["Set of wavelets",Black,20], LabelStyle->Directive[Black,Italic], ImageSize->Large,InterpolationOrder->2]
Um conjunto de sinais de fonte com frequências de 35 Hz, 55 Hz, 85 Hz
Após calcular os traços sísmicos e traçar as amplitudes da onda refletida, podemos ver que, para diferentes frequências, uma anomalia é observada em diferentes espessuras da formação.
tracesSet=Table[ListConvolve[waveletSet[[j]][[1;;;;1]],rc[[i]]],{j,1,3},{i,1,200}]; lowFreq=ListLinePlot[GaussianFilter[Abs[tracesSet[[1]][[All,46]]],3][[;;;;2]],InterpolationOrder->2,PlotStyle->Black,PlotRange->All]; medFreq=ListLinePlot[GaussianFilter[Abs[tracesSet[[2]][[All,46]]],3][[;;;;2]],InterpolationOrder->2,PlotStyle->Black,PlotRange->All]; highFreq=ListLinePlot[GaussianFilter[Abs[tracesSet[[3]][[All,46]]],3][[;;;;2]],InterpolationOrder->2,PlotStyle->Black,PlotRange->All]; Show[lowFreq,medFreq,highFreq,PlotRange->{{0,100},All}, PlotLabel->Style["Amplitudes of reflection",Black,20], LabelStyle->Directive[Black,Italic], Frame->True]
Gráficos das amplitudes da onda refletida a partir da borda superior da cunha para diferentes frequências
A capacidade de tirar conclusões sobre a espessura do reservatório com base nos resultados de observações sísmicas é extremamente útil, porque uma das principais tarefas na prospecção de um campo de petróleo é avaliar os pontos mais promissores para a colocação de um poço (ou seja, as áreas em que o reservatório tem uma grande espessura). Além disso, na seção geológica, pode haver objetos que, devido à sua gênese, causem uma mudança acentuada na espessura da formação. Isso faz da análise espectral uma ferramenta eficaz para estudá-los. Na próxima parte do artigo, consideraremos esses objetos geológicos com mais detalhes.
Dados experimentais. Onde eles são recebidos e o que procurar neles?
Os materiais analisados no artigo foram obtidos no território da Sibéria Ocidental. A região, como todos, sem exceção, provavelmente sabem, é a principal região produtora de petróleo do nosso país. O desenvolvimento ativo dos nascimentos metropolitanos começou na região nos anos 60 do século passado. O principal método para encontrar campos de petróleo é a exploração sísmica. É interessante considerar imagens de satélite deste território. Em pequena escala, pode-se observar um grande número de pântanos e lagos, aumentando o mapa, você pode ver poços de cluster para perfurar poços e aumentando o mapa até o limite, também é possível distinguir glades de perfis pelas quais foram feitas observações sísmicas.
Imagem de satélite dos mapas Yandex - a área da cidade Noyabrsk:
Uma rede de sites de cluster em um dos campos:
As rochas contendo óleo da Sibéria Ocidental ocorrem em uma ampla gama de profundidades - de 1 a 5 km. O principal volume de rochas contendo óleo é formado nos tempos Jurássico e Cretáceo. O período jurássico é provavelmente conhecido por muitos pelo filme de mesmo nome.
O clima do período jurássico foi significativamente diferente do período moderno. Na Enciclopédia Britânica, há uma série de cartões de memória que caracterizam todas as épocas gelógicas.
Hora atual:
Período Jurássico:
Observe que, no Jurássico, o território da Sibéria Ocidental era uma costa marítima (terra atravessada por rios e mar raso). Como o clima era confortável, pode-se supor que uma paisagem típica da época era a seguinte:
Sibéria do Jurássico:
Nesta foto, tanto animais e pássaros são importantes para nós como a imagem do rio ao fundo. Um rio é o próprio objeto geológico sobre o qual paramos antes. O fato é que a atividade dos rios permite o acúmulo de arenitos bem classificados, que depois se tornam um reservatório de petróleo. Esses reservatórios podem ter uma forma bizarra e complexa (como um leito do rio) e têm uma espessura variável - ao longo da costa, a espessura é pequena e aumenta mais perto do centro do canal ou em seções dos meandros. Portanto, os rios formados no período jurássico estão agora a uma profundidade de cerca de três quilômetros e são objeto da busca por reservatórios de petróleo.
Dados experimentais. Processamento e visualização
Nós imediatamente fazemos uma reserva em relação aos materiais sísmicos mostrados no artigo - devido ao fato de a quantidade de dados usada para análise ser significativa - apenas um fragmento do conjunto original de traços sísmicos é colocado no texto do artigo. Isso permitirá que todos reproduzam os cálculos acima.
Ao trabalhar com dados sísmicos, os geofísicos geralmente usam software especializado (existem vários líderes do setor cujo desenvolvimento é usado ativamente, por exemplo, Petrel ou Paradigm), que permite analisar diferentes tipos de dados e possui uma interface gráfica conveniente. Apesar de toda a conveniência, esses tipos de software também têm suas desvantagens - por exemplo, a implementação de algoritmos modernos em versões estáveis leva muito tempo e as possibilidades de automação de cálculos são geralmente limitadas. Em tal situação, torna-se muito conveniente usar sistemas matemáticos de computador e linguagens de programação de alto nível que permitem usar uma ampla base algorítmica e, ao mesmo tempo, adotar muita rotina. Usando esse princípio, o trabalho com dados sísmicos no Wolfram Mathematica é construído. É impraticável escrever uma rica funcionalidade do trabalho interativo com dados - é mais importante garantir o carregamento de um formato geralmente aceito, aplicar os algoritmos desejados a eles e enviá-los novamente para um formato externo.
Seguindo o esquema proposto, carregue os dados sísmicos originais e exiba-os no
Wolfram Mathematica :
Get["GeologyIO`"] seismic3DZipPath = "seismic3D.zip"; seismic3DSEGYPath = "seismic3D.sgy"; If[FileExistsQ[seismic3DZipPath], DeleteFile[seismic3DZipPath]]; If[FileExistsQ[seismic3DSEGYPath], DeleteFile[seismic3DSEGYPath]]; URLDownload["https://wolfr.am/FiQIuZuH", seismic3DZipPath]; ExtractArchive[seismic3DZipPath]; seismic3DSEGY = SEGYImport[seismic3DSEGYPath]
Os dados baixados e importados dessa maneira são as faixas registradas em um site de 10 a 5 quilômetros. No caso de os dados serem obtidos pelo método de exploração sísmica tridimensional (as ondas não são registradas ao longo de perfis geofísicos individuais, mas por toda a área ao mesmo tempo), torna-se possível obter cubos de dados sísmicos. São objetos tridimensionais, seções verticais e horizontais, que permitem um estudo detalhado do ambiente geológico. No exemplo considerado, estamos lidando apenas com dados tridimensionais. Podemos obter algumas informações do cabeçalho do texto, por exemplo, como este
StringPartition[seismic3DSEGY["textheader"], 80] // TableForm
C 1 ESTE É UM ARQUIVO DE DEMONSTRAÇÃO PARA O ENSAIO DE PACOTE GEOLOGYIO
C 2
C 3
C 4
C 5 DATA NOME DO USUÁRIO: USUÁRIO WOLFRAM
C 6 NOME DA PESQUISA: EM ALGUM LUGAR DA SIBÉRIA
C 7 VOLUME SÍSMICO DE TIPO DE ARQUIVO 3D
C 8
C 9
C10 Z RANGE: PRIMEIROS 2200M ÚLTIMOS 2400M
Esse conjunto de dados será suficiente para demonstrarmos os principais estágios da análise de dados. As trilhas no arquivo são gravadas seqüencialmente e cada uma delas se parece com a figura a seguir - esta é a distribuição das amplitudes das ondas refletidas ao longo do eixo vertical (eixo de profundidade).
ListLinePlot[seismic3DSEGY["traces"][[100]], InterpolationOrder -> 2, PlotStyle -> Black, PlotLabel -> Style["Seismic trace", Black, 20], LabelStyle -> Directive[Black, Italic], PlotRange -> All, Frame -> True, ImageSize -> 1200, AspectRatio -> 1/5]
Uma das seções sísmicas:
Sabendo quantos traços estão localizados em cada direção da área estudada, você pode formar um array de dados tridimensional e exibi-lo usando a função Image3D []
traces=seismic3DSEGY["traces"]; startIL=1050;EndIL=2000;stepIL=2; (* *) startXL=1165;EndXL=1615;stepXL=2; (* Y *) numIL=(EndIL-startIL)/stepIL+1; (* *) numXL=(EndXL-startXL)/stepIL+1; (* Y*) Image3D[ArrayReshape[Abs[traces/Max[Abs[traces[[All,1;;;;4]]]]],{numIL,numXL,101}],ViewPoint->{-1, 0, 0},Background->RGBColor[0,0,0]]
Imagem tridimensional de um cubo de dados sísmicos (eixo vertical - profundidade):
Caso objetos geológicos de interesse criem intensas anomalias sísmicas, podem ser usadas ferramentas de visualização com transparência. As seções “sem importância” da gravação podem ficar invisíveis, deixando apenas as anomalias visíveis. No Wolfram Mathematica, isso pode ser feito usando o
Opacity [] e o
Raster3D [] .
data = ArrayReshape[Abs[traces/Max[Abs[traces[[All,1;;;;4]]]]],{numIL,numXL,101}]; Graphics3D[{Opacity[0.1], Raster3D[data, ColorFunction->"RainbowOpacity"]}, Boxed->False, SphericalRegion->True, ImageSize->840, Background->None]
Imagem de um cubo de dados sísmicos usando as funções Opacity [] e Raster3D []
Como no exemplo sintético, nas fatias do cubo original, alguns limites geológicos (camadas) com relevo variável podem ser distinguidos.
A principal ferramenta de análise espectral é a transformada de Fourier. Usando-o, você pode avaliar o espectro de amplitude-frequência de cada caminho ou grupo de caminhos. No entanto, após a transferência de dados para o domínio da frequência, as informações são perdidas sobre a que horas (leia a que profundidade) a frequência muda. Para poder localizar alterações de sinal no eixo do tempo (profundo), são utilizadas a transformação de Fourier da janela e a decomposição da wavelet. Este artigo usa decomposição de wavelet. A tecnologia de análise Wavelet começou a ser usada ativamente na exploração sísmica nos anos 90. A vantagem sobre a janela de transformação de Fourier é considerada a melhor resolução de tempo.
Usando o seguinte fragmento de código, você pode se decompor em componentes individuais de um dos rastreamentos sísmicos:
cwd=ContinuousWaveletTransform[seismicSection["traces"][[100]]] Show[ ListLinePlot[Re[cwd[[1]]],PlotRange->All], ListLinePlot[seismicSection["traces"][[100]], PlotStyle->Black,PlotRange->All],ImageSize->{1500,500},AspectRatio->Full, PlotLabel->Style["Wavelet decomposition",Black,32], LabelStyle->Directive[Black,Italic], PlotRange->All, Frame->True]
Decomposição de componentes
Para avaliar como a energia de reflexão é distribuída em diferentes momentos da chegada das ondas, são utilizados escalogramas (análogos do espectrograma). Como regra, na prática, não há necessidade de analisar todos os componentes. Geralmente escolha o componente de baixa, média e alta frequência.
freq=(500/(#*contWD["Wavelet"]["FourierFactor"]))&/@(Thread[{Range[contWD["Octaves"]],1}]/.contWD["Scales"])//Round; ticks=Transpose[{Range[Length[freq]],freq}]; WaveletScalogram[contWD,Frame->True,FrameTicks->{{ticks,Automatic},Automatic},FrameTicksStyle->Directive[Orange,12], FrameLabel->{"Time","Frequency(Hz)"},LabelStyle->Directive[Black,Bold,14], ColorFunction->"RustTones",ImageSize->Large]
Scalogram. Resultado da função
WaveletScalogram []
O Wolfram Language usa a função ContinuousWaveletTransform [] para transformações de wavelet . E a aplicação dessa função a todo o conjunto de rastreios é realizada com o uso da função Table [] . Aqui vale a pena notar que um dos pontos fortes do Wolfram Mathematica é a capacidade de usar a paralelização ParallelTable [] . No exemplo dado, a paralelização não é necessária - o volume de dados não é grande, mas ao trabalhar com conjuntos de dados experimentais contendo centenas de milhares de rastreios, isso é uma necessidade. tracesCWD=Table[Map[Hilbert[#,0]&,Re[ContinuousWaveletTransform[traces[[i]]][[1]]][[{13,15,18}]]],{i,1,Length@traces}];
Após aplicar a função ContinuousWaveletTransform [] , novas matrizes de dados correspondentes às frequências selecionadas são exibidas. No exemplo acima, estas são frequências: 38Hz, 33Hz, 27Hz. A escolha das frequências é feita com mais frequência com base em testes - eles obtêm mapas eficazes para diferentes combinações de frequências e escolhem os mais informativos do ponto de vista de um geólogo.Se você precisar compartilhar os resultados com colegas ou fornecê-los ao cliente, poderá usar a função SEGYExport [] do GeologyIO outputdata=seismic3DSEGY; outputdata["traces",1;;-1]=tracesCWD[[All,3]]; outputdata["textheader"]="Wavelet Decomposition Result"; outputdata["binaryheader","NumberDataTraces"]=Length[tracesCWD[[All,3]]]; SEGYExport["D:\\result.segy",outputdata];
Tendo três desses cubos disponíveis (componentes de baixa frequência, média e alta frequência), como regra, a mistura RGB é usada para visualização conjunta de dados. Cada componente recebe sua própria cor - vermelho, verde, azul. No Wolfram Mathematica, isso pode ser feito usando a função ColorCombine [] .Como resultado, são obtidas imagens que podem ser usadas para executar uma interpretação geológica. Os meandros pendurados na seção permitem delinear o paleorus, que provavelmente é um reservatório e contém reservas de petróleo. A busca e análise de análogos modernos de um sistema desse tipo de rio nos permite determinar as partes mais promissoras dos meandros. As próprias camas são caracterizadas por grossas camadas de arenito bem classificado e são um bom reservatório de óleo. Áreas fora das anomalias de "renda" são semelhantes aos depósitos modernos de várzea. Os sedimentos das planícies de inundação são representados principalmente por rochas argilosas e a perfuração nessas zonas será ineficaz.Fatia RGB do cubo de dados. No centro (um pouco à esquerda do centro), você pode traçar o rio sinuoso.
Fatia RGB do cubo de dados. No lado esquerdo, você pode traçar o rio sinuoso.
Em alguns casos, a qualidade dos dados sísmicos permite imagens significativamente mais nítidas. Depende da metodologia do trabalho de campo, equipamento adotado pelo algoritmo de redução de ruído. Nesses casos, não apenas fragmentos do sistema fluvial são visíveis, mas também paleorecs estendidos inteiros.Mistura RGB dos três componentes do cubo de dados sísmicos (fatia horizontal). A profundidade é de aproximadamente 2 km.
Imagem de satélite do rio Volga na região de Saratov
Conclusão
No Wolfram Mathematica, você pode analisar dados sísmicos e resolver problemas de aplicativos associados à pesquisa de minerais, e o pacote GeologyIO permite tornar esse processo mais conveniente. A estrutura dos dados sísmicos é tal que o uso de métodos internos para acelerar cálculos ( ParallelTable [] , ParallelDo [], ...) é muito eficiente e permite processar grandes quantidades de dados. Em grande parte, isso é facilitado pelos recursos de armazenamento de dados do GeologyIO. A propósito, o pacote pode ser usado não apenas no campo da exploração sísmica aplicada. Quase os mesmos tipos de dados são usados em georadar e sismologia.Se você tiver sugestões sobre como melhorar o resultado, quais algoritmos de análise de sinal do Wolfram Mathematica são aplicáveis a esses dados, ou se você tiver comentários críticos, deixe comentários.