Wolfram Mathematica em Geofísica

Obrigado ao autor do blog Anton Ekimenko pelo seu relatório

Introdução

Esta nota foi escrita após a conferência Conferência de Tecnologia Russa Wolfram e contém um resumo do relatório que apresentei. O evento aconteceu em junho em São Petersburgo. Considerando que trabalho a um quarteirão do local da conferência, não poderia deixar de comparecer a este evento. Em 2016 e 2017, ouvi relatórios de conferências e este ano fiz uma apresentação. Em primeiro lugar, apareceu um tema interessante (me parece), que estamos desenvolvendo com Kirill Belove, em segundo lugar, após um longo estudo da legislação da Federação Russa em relação à política de sanções, na empresa onde trabalho, apareceram até duas licenças Wolfram Mathematica.

Antes de passar ao tema do meu discurso, gostaria de destacar a boa organização do evento. A página de visita da conferência utiliza uma imagem da Catedral de Kazan. A catedral é uma das principais atrações de São Petersburgo e é claramente visível do salão onde ocorreu a conferência.

Wolfram Mathematica em Geofísica

Na entrada da Universidade Estadual de Economia de São Petersburgo, os participantes foram recebidos por assistentes dentre os estudantes - eles não permitiram que se perdessem. Durante a inscrição, foram distribuídos pequenos souvenirs (um brinquedo - uma ponta piscante, uma caneta, adesivos com símbolos Wolfram). Os intervalos para almoço e café também foram incluídos na programação da conferência. Já comentei sobre deliciosos cafés e tortas na parede do grupo - os chefs são ótimos. Com esta parte introdutória, gostaria de ressaltar que o evento em si, seu formato e localização já trazem emoções positivas.

O relatório que foi preparado por mim e Kirill Belov se chama “Usando o Wolfram Mathematica para resolver problemas em geofísica aplicada. Análise espectral de dados sísmicos ou “onde corriam rios antigos”. O conteúdo do relatório abrange duas partes: em primeiro lugar, a utilização de algoritmos disponíveis em Wolfram Mathematica para analisar dados geofísicos e, em segundo lugar, é assim que colocar dados geofísicos no Wolfram Mathematica.

Exploração sísmica

Primeiro você precisa fazer uma pequena excursão pela geofísica. Geofísica é a ciência que estuda as propriedades físicas das rochas. Pois bem, como as rochas têm propriedades diferentes: eléctricas, magnéticas, elásticas, existem métodos geofísicos correspondentes: prospecção eléctrica, prospecção magnética, prospecção sísmica... No contexto deste artigo, discutiremos apenas a prospecção sísmica com mais detalhes. A exploração sísmica é o principal método de busca de petróleo e gás. O método baseia-se na excitação de vibrações elásticas e posterior registro da resposta das rochas que compõem a área de estudo. As vibrações são excitadas em terra (com dinamite ou fontes de vibração não explosivas de vibrações elásticas) ou no mar (com armas de ar comprimido). As vibrações elásticas se propagam pelo maciço rochoso, sendo refratadas e refletidas nos limites de camadas com diferentes propriedades. As ondas refletidas retornam à superfície e são registradas por geofones em terra (geralmente dispositivos eletrodinâmicos baseados no movimento de um ímã suspenso em uma bobina) ou hidrofones no mar (baseados no efeito piezoelétrico). No momento da chegada das ondas, pode-se avaliar a profundidade das camadas geológicas.

Equipamento de reboque de embarcações sísmicas
Wolfram Mathematica em Geofísica

A pistola de ar excita vibrações elásticas
Wolfram Mathematica em Geofísica

As ondas passam pelo maciço rochoso e são registradas por hidrofones
Wolfram Mathematica em Geofísica

Navio de pesquisa geofísica "Ivan Gubkin" no cais próximo à ponte Blagoveshchensky em São Petersburgo
Wolfram Mathematica em Geofísica

Modelo de sinal sísmico

As rochas têm propriedades físicas diferentes. Para a exploração sísmica, as propriedades elásticas são principalmente importantes - a velocidade de propagação das vibrações elásticas e a densidade. Se duas camadas tiverem propriedades iguais ou semelhantes, a onda “não notará” a fronteira entre elas. Se as velocidades das ondas nas camadas forem diferentes, a reflexão ocorrerá nos limites das camadas. Quanto maior a diferença nas propriedades, mais intensa será a reflexão. Sua intensidade será determinada pelo coeficiente de refletância (rc):

Wolfram Mathematica em Geofísica

onde ρ é a densidade da rocha, ν é a velocidade da onda, 1 e 2 indicam as camadas superior e inferior.

Um dos modelos de sinais sísmicos mais simples e utilizados é o modelo de convolução, quando o traço sísmico registrado é representado como o resultado da convolução de uma sequência de coeficientes de reflexão com um pulso de sondagem:

Wolfram Mathematica em Geofísica

onde está(t) — traço sísmico, ou seja, tudo o que foi gravado por um hidrofone ou geofone durante um tempo fixo de gravação, o(a) - o sinal gerado pela pistola de ar comprimido, n(t) - ruído aleatório.

Vamos calcular um traço sísmico sintético como exemplo. Usaremos o pulso de Ricker, amplamente utilizado na exploração sísmica, como 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]

Impulso sísmico inicial
Wolfram Mathematica em Geofísica

Definiremos 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]]

Sequência de coeficientes de reflexão
Wolfram Mathematica em Geofísica

Vamos calcular e exibir o traço sísmico. Como os coeficientes de reflexão têm sinais diferentes, obtemos duas reflexões alternadas no traço 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]]

Pista simulada
Wolfram Mathematica em Geofísica

Para este exemplo é necessário fazer uma ressalva - na realidade, a profundidade das camadas é determinada, claro, em metros, e o cálculo do traço sísmico ocorre para o domínio do tempo. Seria mais correto definir as profundidades em metros e calcular os tempos de chegada conhecendo as velocidades nas camadas. Nesse caso, defino imediatamente as camadas no eixo do tempo.

Se falamos de pesquisa de campo, então, como resultado de tais observações, um grande número de séries temporais semelhantes (traços sísmicos) são registrados. Por exemplo, ao estudar um local com 25 km de comprimento e 15 km de largura, onde, como resultado do trabalho, cada traço caracteriza uma célula de 25x25 metros (tal célula é chamada de bin), a matriz de dados final conterá 600000 traços. Com um tempo de amostragem de 1 ms e um tempo de gravação de 5 segundos, o arquivo de dados final terá mais de 11 GB e o volume do material “bruto” original pode ser de centenas de gigabytes.

Como trabalhar com eles Wolfram Mathematica?

Pacote GeologiaIO

O desenvolvimento do pacote começou pergunta na parede VK do grupo de apoio de língua russa. Graças às respostas da comunidade, foi encontrada uma solução muito rapidamente. E, como resultado, tornou-se um desenvolvimento sério. Correspondente Postagem no mural da Comunidade Wolfram Foi até marcado pelos moderadores. Atualmente, o pacote suporta trabalhar com os seguintes tipos de dados que são usados ​​ativamente na indústria geológica:

  1. importação de dados de mapas nos formatos ZMAP e IRAP
  2. importação de medições em poços no formato LAS
  3. entrada e saída de formato de arquivos sísmicos 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 Caderno matemático:

If[PacletInformation["GeologyIO"] === {}, PacletInstall[URLDownload[
    "https://wolfr.am/FiQ5oFih", 
    FileNameJoin[{CreateDirectory[], "GeologyIO-0.2.2.paclet"}]
]]]

Após o qual o pacote será instalado na pasta padrão, cujo caminho pode ser obtido da seguinte forma:

FileNameJoin[{$UserBasePacletsDirectory, "Repository"}]

Como exemplo, demonstraremos as principais capacidades do pacote. A chamada é feita tradicionalmente para pacotes na Wolfram Language:

Get["GeologyIO`"]

O pacote é desenvolvido usando Bancada Wolfram. Isso permite acompanhar as principais funcionalidades do pacote com documentação, que em termos de formato de apresentação não difere da documentação do próprio Wolfram Mathematica, e fornecer ao pacote arquivos de teste para o primeiro contato.

Wolfram Mathematica em Geofísica

Wolfram Mathematica em Geofísica

Tal arquivo, em particular, é o arquivo “Marmousi.segy” - este é um modelo sintético de uma seção geológica, que foi desenvolvido pelo Instituto Francês de Petróleo. Usando este modelo, os desenvolvedores testam seus próprios algoritmos para modelagem de campos de ondas, processamento de dados, inversão de traços sísmicos, etc. O próprio modelo Marmousi é armazenado no repositório de onde o 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"]

Resultado da importação - objeto SEGYData
Wolfram Mathematica em Geofísica

O formato SEGY envolve o armazenamento de diversas informações sobre as observações. Em primeiro lugar, estes são comentários de texto. Isso inclui informações sobre o local da obra, os nomes das empresas que realizaram as medições, etc. No nosso caso, este 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 Marmousi foi gerado no Instituto...velocidade mínima de 1500 m/s e máxima de 5500 m/s)”

Você pode exibir o modelo geológico real acessando os traços sísmicos usando a chave “traços” (uma das características do pacote é que as chaves não diferenciam maiúsculas de minúsculas):

ArrayPlot[Transpose[marmousi["traces"]], PlotTheme -> "Detailed"]

Modelo Marmousi
Wolfram Mathematica em Geofísica

Atualmente, o pacote também permite carregar dados em partes de arquivos grandes, possibilitando o processamento de arquivos cujo tamanho pode chegar a dezenas de gigabytes. As funções do pacote também incluem funções para exportar dados para .segy e anexá-los parcialmente ao final do arquivo.

Separadamente, vale destacar a funcionalidade do pacote ao trabalhar com a estrutura complexa de arquivos .segy. Uma vez que permite não apenas acessar rastreamentos e cabeçalhos individuais usando chaves e índices, mas também alterá-los e gravá-los em um arquivo. Muitos dos detalhes técnicos da implementação do GeologyIO estão além do escopo deste artigo e provavelmente merecem uma descrição separada.

Relevância da análise espectral na exploração sísmica

A capacidade de importar dados sísmicos para o Wolfram Mathematica permite usar a funcionalidade integrada de processamento de sinais para dados experimentais. Como cada traço sísmico representa uma série temporal, uma das principais ferramentas para estudá-los é a análise espectral. Dentre os pré-requisitos para análise da composição frequencial dos dados sísmicos, podemos citar, por exemplo, os seguintes:

  1. Diferentes tipos de ondas são caracterizados por diferentes composições de frequência. Isso permite destacar ondas úteis e suprimir ondas de interferência.
  2. Propriedades da rocha, como porosidade e saturação, podem afetar a composição de frequência. Isso permite identificar rochas com as melhores propriedades.
  3. Camadas com diferentes espessuras causam anomalias em diferentes faixas de frequência.

O terceiro ponto é o principal no contexto deste artigo. Abaixo está um fragmento de código para cálculo de traços sísmicos no caso de uma camada com espessura variável - um modelo de cunha. Este modelo é tradicionalmente estudado na exploração sísmica para analisar os efeitos de interferência quando as ondas refletidas de muitas camadas são sobrepostas umas às outras.

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 pinçada
Wolfram Mathematica em Geofísica

A velocidade da onda dentro da cunha é 4500 m/s, fora da cunha 4000 m/s, e a densidade é considerada constante 2200 g/cm³. Para tal modelo, calculamos 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]

Traços sísmicos para o modelo em cunha
Wolfram Mathematica em Geofísica

A sequência de traços sísmicos mostrada nesta figura é chamada de seção sísmica. Como você pode ver, sua interpretação também pode ser feita de forma intuitiva, pois a geometria das ondas refletidas corresponde claramente ao modelo especificado anteriormente. Se você analisar os traços com mais detalhes, notará que os traços de 1 a aproximadamente 30 não diferem - os reflexos do telhado da formação e do fundo não se sobrepõem. A partir do 31º traço, os reflexos começam a interferir. E, embora no modelo os coeficientes de reflexão não mudem horizontalmente - os traços sísmicos mudam de intensidade à medida que muda a espessura da formação.

Consideremos a amplitude de reflexão do limite superior da formação. A partir do 60º percurso a intensidade da reflexão começa a aumentar e no 70º percurso torna-se máxima. É assim que se manifesta a interferência das ondas do teto e do fundo das camadas, levando em alguns casos a anomalias significativas no 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 na borda superior da cunha
Wolfram Mathematica em Geofísica

É lógico que quando o sinal é de frequência mais baixa, a interferência começa a aparecer em grandes espessuras de formação e, no caso de um sinal de alta frequência, a interferência ocorre em espessuras menores. O trecho 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 origem com frequências de 35 Hz, 55 Hz, 85 Hz
Wolfram Mathematica em Geofísica

Ao calcular os traços sísmicos e traçar gráficos das amplitudes das ondas refletidas, podemos ver que para diferentes frequências uma anomalia é observada em diferentes espessuras de 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 da borda superior da cunha para diferentes frequências
Wolfram Mathematica em Geofísica

A capacidade de tirar conclusões sobre a espessura da formação a partir dos resultados das observações sísmicas é extremamente útil, porque uma das principais tarefas na exploração de petróleo é avaliar os pontos mais promissores para a colocação de um poço (ou seja, aquelas áreas onde a formação é mais grosso). Além disso, na seção geológica podem existir objetos cuja gênese provoca uma mudança brusca na espessura da formação. Isso torna a 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 você os conseguiu e o que procurar neles?

Os materiais analisados ​​no artigo foram obtidos na Sibéria Ocidental. A região, como provavelmente todos, sem exceção, sabem, é a principal região produtora de petróleo do nosso país. O desenvolvimento ativo de jazidas começou na região na década de 60 do século passado. O principal método de busca de depósitos de petróleo é a exploração sísmica. É interessante observar imagens de satélite deste território. Em pequena escala, você pode notar um grande número de pântanos e lagos; ampliando o mapa, você pode ver locais de perfuração de poços agrupados, e ampliando o mapa até o limite, você também pode distinguir as clareiras dos perfis ao longo dos quais a sísmica observações foram realizadas.

Imagem de satélite de mapas Yandex - área da cidade de Noyabrsk
Wolfram Mathematica em Geofísica

Uma rede de poços em um dos campos
Wolfram Mathematica em Geofísica

As rochas petrolíferas da Sibéria Ocidental ocorrem em uma ampla variedade de profundidades - de 1 km a 5 km. O principal volume de rochas contendo petróleo foi formado nos tempos Jurássico e Cretáceo. O período Jurássico é provavelmente conhecido por muitos pelo filme de mesmo nome. Clima jurássico era significativamente diferente do moderno. A Enciclopédia Britânica possui uma série de paleomapas que caracterizam cada era helógica.

Actualmente,
Wolfram Mathematica em Geofísica
Período jurássico
Wolfram Mathematica em Geofísica

Observe que na época jurássica, o território da Sibéria Ocidental era uma costa marítima (terra atravessada por rios e um mar raso). Como o clima era confortável, podemos supor que uma paisagem típica da época era assim:

Sibéria Jurássica
Wolfram Mathematica em Geofísica

Nesta fotografia, o que importa para nós não são tanto os animais e os pássaros, mas a imagem do rio ao fundo. O rio é o mesmo objeto geológico em que paramos anteriormente. O fato é que a atividade dos rios permite o acúmulo de arenitos bem selecionados, que se tornarão reservatório de petróleo. Esses reservatórios podem ter uma forma bizarra e complexa (como o leito de um rio) e têm espessura variável - perto das margens a espessura é pequena, mas mais perto do centro do canal ou em áreas de meandros ela aumenta. Assim, os rios formados no Jurássico estão hoje com cerca de três quilômetros de profundidade e são objeto de busca por reservatórios de petróleo.

Dados experimentais. Processamento e visualização

Façamos imediatamente uma ressalva relativamente aos materiais sísmicos apresentados no artigo - devido ao facto de a quantidade de dados utilizados para a análise ser significativa - apenas um fragmento do conjunto original de vestígios sísmicos está incluído no texto do artigo. Isso permitirá que qualquer pessoa reproduza os cálculos acima.

Ao trabalhar com dados sísmicos, um geofísico costuma utilizar software especializado (existem vários líderes do setor cujos desenvolvimentos são ativamente utilizados, por exemplo Petrel ou Paradigm), que permite analisar diferentes tipos de dados e possui uma interface gráfica conveniente. Apesar de toda a comodidade, 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 automatização de cálculos costumam ser limitadas. Em tal situação, torna-se muito conveniente utilizar sistemas matemáticos computacionais e linguagens de programação de alto nível, que permitem a utilização de uma ampla base algorítmica e, ao mesmo tempo, exigem muita rotina. Este é o princípio usado para trabalhar com dados sísmicos no Wolfram Mathematica. Não é apropriado escrever funcionalidades ricas para trabalho interativo com dados - é mais importante garantir o carregamento de um formato geralmente aceito, aplicando os algoritmos desejados a eles e enviando-os de volta para um formato externo.

Seguindo o esquema proposto, carregaremos os dados sísmicos originais e os exibiremos em 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 desta forma são as rotas registradas em uma área de 10 por 5 quilômetros. Se os dados forem obtidos através de um método de levantamento sísmico tridimensional (as ondas são registradas não ao longo de perfis geofísicos individuais, mas em toda a área simultaneamente), torna-se possível obter cubos de dados sísmicos. São objetos tridimensionais, cujas seções verticais e horizontais permitem um estudo detalhado do ambiente geológico. No exemplo considerado, estamos lidando com dados tridimensionais. Podemos obter algumas informações do cabeçalho do texto, como esta

StringPartition[seismic3DSEGY["textheader"], 80] // TableForm

C 1 ESTE É UM ARQUIVO DE DEMO PARA TESTE DE PACOTE GEOLOGYIO
C 2
C 3
C 4
C 5 DATA NOME DE USUÁRIO: USUÁRIO WOLFRAM
NOME DA PESQUISA C 6: EM ALGUM LUGAR DA SIBÉRIA
C 7 TIPO DE ARQUIVO VOLUME SÍSMICO 3D
C 8
C 9
GAMA C10 Z: PRIMEIROS 2200M ÚLTIMOS 2400M

Este conjunto de dados será suficiente para demonstrarmos as principais etapas da análise dos dados. Os traços no arquivo são registrados sequencialmente e cada um deles 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]

Um dos traços da seção sísmica
Wolfram Mathematica em Geofísica

Sabendo quantos traços estão localizados em cada direção da área estudada, você pode gerar uma matriz de dados tridimensional e exibi-la 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 XNUMXD de um cubo de dados sísmicos (eixo vertical - profundidade)
Wolfram Mathematica em Geofísica

Se as características geológicas de interesse criarem anomalias sísmicas intensas, então ferramentas de visualização com transparência poderão ser utilizadas. Áreas “sem importância” da gravação podem ficar invisíveis, deixando apenas anomalias visíveis. No Wolfram Mathematica isso pode ser feito usando Opacidade[] и 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 cubo de dados sísmicos usando funções Opacity[] e Raster3D[] Wolfram Mathematica em Geofísica

Tal como no exemplo sintético, em secções do cubo original podem-se identificar alguns limites geológicos (camadas) com relevo variável.

A principal ferramenta para análise espectral é a transformada de Fourier. Com sua ajuda, você pode avaliar o espectro amplitude-frequência de cada traço ou grupo de traços. No entanto, após a transferência dos dados para o domínio da frequência, perde-se informação sobre em que momentos (leia-se em que profundidades) a frequência muda. Para poder localizar mudanças de sinal no eixo do tempo (profundidade), são utilizadas a transformada de Fourier em janela e a decomposição wavelet. Este artigo usa decomposição wavelet. A tecnologia de análise wavelet começou a ser usada ativamente na exploração sísmica na década de 90. A vantagem sobre a transformada de Fourier em janela é considerada uma melhor resolução de tempo.

Usando o seguinte fragmento de código, você pode decompor um dos traços sísmicos em componentes individuais:

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 um traço em componentes
Wolfram Mathematica em Geofísica

Para avaliar como a energia de reflexão é distribuída em diferentes tempos de chegada das ondas, são utilizados escalogramas (análogos a um espectrograma). Via de regra, na prática não há necessidade de analisar todos os componentes. Normalmente, os componentes de baixa, média e alta frequência são selecionados.

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]

Escalograma. Resultado da função Escalograma Wavelet[]
Wolfram Mathematica em Geofísica

A Wolfram Language usa a função para transformação wavelet Transformação Wavelet Contínua[]. E a aplicação desta função a todo o conjunto de traços será realizada através da função Mesa[]. Aqui vale destacar um dos pontos fortes do Wolfram Mathematica - a capacidade de usar paralelização Tabela Paralela[]. No exemplo acima, não há necessidade de paralelização - o volume de dados não é grande, mas ao trabalhar com conjuntos de dados experimentais contendo centenas de milhares de rastreamentos, isso é uma necessidade.

tracesCWD=Table[Map[Hilbert[#,0]&,Re[ContinuousWaveletTransform[traces[[i]]][[1]]][[{13,15,18}]]],{i,1,Length@traces}]; 

Depois de aplicar a função Transformação Wavelet Contínua[] Novos conjuntos de dados aparecem correspondentes às frequências selecionadas. No exemplo acima, essas frequências são: 38Hz, 33Hz, 27Hz. A escolha das frequências é mais frequentemente realizada com base em testes - eles obtêm mapas eficazes para diferentes combinações de frequências e selecionam o mais informativo do ponto de vista de um geólogo.

Se precisar compartilhar os resultados com colegas ou fornecê-los ao cliente, você pode usar a função SEGYExport[] do pacote 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];

Com três desses cubos (componentes de baixa frequência, média frequência e alta frequência), a combinação RGB é normalmente usada para visualizar os dados juntos. Cada componente recebe sua própria cor - vermelho, verde, azul. No Wolfram Mathematica isso pode ser feito usando a função Combinar cores[].

O resultado são imagens a partir das quais a interpretação geológica pode ser feita. Os meandros registrados no trecho permitem delinear paleocanais, mais provavelmente reservatórios e contendo reservas de petróleo. A busca e análise de análogos modernos desse sistema fluvial nos permite determinar as partes mais promissoras dos meandros. Os próprios canais são caracterizados por espessas camadas de arenito bem selecionado e são um bom reservatório de petróleo. As áreas fora das anomalias "rendadas" são semelhantes aos depósitos modernos de várzea. Os depósitos de várzea são representados principalmente por rochas argilosas e a perfuração nestas zonas será ineficaz.

Fatia RGB do cubo de dados. No centro (um pouco à esquerda do centro) você pode traçar o rio sinuoso.
Wolfram Mathematica em Geofísica
Fatia RGB do cubo de dados. No lado esquerdo você pode traçar o rio sinuoso.
Wolfram Mathematica em Geofísica

Em alguns casos, a qualidade dos dados sísmicos permite imagens significativamente mais nítidas. Isso depende da metodologia de trabalho de campo, do equipamento utilizado pelo algoritmo de redução de ruído. Nesses casos, não apenas fragmentos de sistemas fluviais são visíveis, mas também paleo-rios inteiros e extensos.

Mistura RGB de três componentes de um cubo de dados sísmicos (fatia horizontal). Profundidade aproximadamente 2 km.
Wolfram Mathematica em Geofísica
Imagem de satélite do rio Volga perto de Saratov
Wolfram Mathematica em Geofísica

Conclusão

O Wolfram Mathematica permite analisar dados sísmicos e resolver problemas aplicados relacionados à exploração mineral, e o pacote GeologyIO torna esse processo mais conveniente. A estrutura dos dados sísmicos é tal que o uso de métodos integrados para acelerar os cálculos (Tabela Paralela[], ParaleloDo[],…) é muito eficiente e permite processar grandes quantidades de dados. Em grande medida, isso é facilitado pelos recursos de armazenamento de dados do pacote GeologyIO. Aliás, o pacote pode ser utilizado não apenas na área de exploração sísmica aplicada. Quase os mesmos tipos de dados são usados ​​em radar de penetração no solo e sismologia. Se você tiver sugestões sobre como melhorar o resultado, quais algoritmos de análise de sinal do arsenal do Wolfram Mathematica são aplicáveis ​​a esses dados, ou se você tiver algum comentário crítico, por favor Deixe um comentário.

Fonte: habr.com

Adicionar um comentário