Wolfram Mathematica во геофизика

Благодарност до авторот на блогот Антон Екименко за неговиот извештај

Вовед

Оваа белешка беше напишана во пресрет на конференцијата Волфрам руска технолошка конференција и содржи резиме на извештајот што го дадов. Настанот се случи во јуни во Санкт Петербург. Имајќи предвид дека работам на еден блок од местото на конференцијата, не можев а да не присуствувам на овој настан. Во 2016 и 2017 година слушав извештаи од конференции, а оваа година одржав презентација. Прво, се појави една интересна (ми се чини) тема со која ја развиваме Кирил Белов, и второ, по долго проучување на законодавството на Руската Федерација во врска со политиката на санкции, во претпријатието каде што работам се појавија дури две лиценци Волфрам Математика.

Пред да преминам на темата на мојот говор, би сакал да ја забележам добрата организација на настанот. Посетената страница на конференцијата користи слика на катедралата Казан. Катедралата е една од главните атракции на Санкт Петербург и е многу јасно видлива од салата во која се одржа конференцијата.

На влезот на Државниот економски универзитет во Санкт Петербург, учесниците ги сретнаа асистенти од студентите - не им дозволија да се изгубат. При регистрацијата беа поделени мали сувенири (играчка - трепкачки шил, пенкало, налепници со симболи на Волфрам). Во распоредот на конференцијата беа вклучени и паузите за ручек и кафе. Веќе забележав за вкусно кафе и пити на ѕидот на групата - готвачите се одлични. Со овој воведен дел сакам да истакнам дека самиот настан, неговиот формат и локација веќе носат позитивни емоции.

Извештајот што го подготвивме јас и Кирил Белов се вика „Користење на Wolfram Mathematica за решавање проблеми во применетата геофизика. Спектрална анализа на сеизмички податоци или „каде течеле античките реки“. Содржината на извештајот опфаќа два дела: прво, употребата на алгоритми достапни во Волфрам Математика за анализа на геофизички податоци, и второ, вака се ставаат геофизички податоци во Wolfram Mathematica.

Сеизмички истражувања

Прво треба да направите кратка екскурзија во геофизиката. Геофизиката е наука која ги проучува физичките својства на карпите. Па, бидејќи карпите имаат различни својства: електрични, магнетни, еластични, постојат соодветни методи на геофизика: електрично истражување, магнетно истражување, сеизмичко истражување... Во контекст на оваа статија, ќе разговараме само за сеизмичкото истражување подетално. Сеизмичкото истражување е главниот метод за барање нафта и гас. Методот се заснова на возбудување на еластичните вибрации и последователно снимање на одговорот од карпите што ја сочинуваат областа на проучување. Вибрациите се возбудуваат на копно (со динамит или неексплозивни извори на вибрации на еластични вибрации) или на море (со воздушни пиштоли). Еластичните вибрации се шират низ карпестата маса, се прекршуваат и се рефлектираат на границите на слоевите со различни својства. Рефлектираните бранови се враќаат на површината и се снимаат со геофони на копно (обично електродинамички уреди базирани на движење на магнет суспендиран во калем) или хидрофони во морето (врз основа на пиезоелектричниот ефект). До моментот на пристигнување на брановите, може да се суди за длабочините на геолошките слоеви.

Опрема за влечење на сеизмички бродови


Воздушниот пиштол возбудува еластични вибрации


Брановите минуваат низ карпестата маса и се снимаат со хидрофони


Геофизички истражувачки брод „Иван Губкин“ на пристаништето кај мостот Благовешченски во Санкт Петербург

Модел на сеизмички сигнал

Карпите имаат различни физички својства. За сеизмичко истражување примарно се важни еластичните својства - брзината на ширење на еластичните вибрации и густината. Ако два слоја имаат исти или слични својства, тогаш бранот „нема да ја забележи“ границата меѓу нив. Ако брзините на брановите во слоевите се разликуваат, тогаш рефлексијата ќе се појави на границата на слоевите. Колку е поголема разликата во својствата, толку е поинтензивен одразот. Неговиот интензитет ќе се определи со коефициентот на рефлексија (rc):

каде ρ е густината на карпата, ν е брзината на бранот, 1 и 2 ги означуваат горните и долните слоеви.

Еден од наједноставните и најчесто користените модели на сеизмички сигнал е моделот на конволуција, кога снимената сеизмичка трага е претставена како резултат на конволуција на низа од коефициенти на рефлексија со пулс за сондирање:

каде s(t) — сеизмичка трага, т.е. сè што е снимено со хидрофон или геофон за време на фиксно време за снимање, w(t) - сигналот генериран од воздушниот пиштол, n(t) - случаен шум.

Дозволете ни да пресметаме синтетичка сеизмичка трага како пример. Како почетен сигнал ќе го користиме пулсот Рикер, широко користен во сеизмичките истражувања.

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]

Почетен сеизмички импулс


Ќе поставиме две граници на длабочини од 300 ms и 600 ms, а коефициентите на рефлексија ќе бидат случајни броеви

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

Низа на коефициенти на рефлексија


Да ја пресметаме и прикажеме сеизмичката трага. Бидејќи коефициентите на рефлексија имаат различни знаци, добиваме две наизменични рефлексии на сеизмичката трага.

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

Симулирана песна


За овој пример, потребно е да се направи резервација - во реалноста, длабочината на слоевите се одредува, се разбира, во метри, а пресметката на сеизмичката трага се случува за временскиот домен. Би било поправилно да се постават длабочините во метри и да се пресметаат времињата на пристигнување знаејќи ги брзините во слоевите. Во овој случај, веднаш ги поставив слоевите на временската оска.

Ако зборуваме за теренски истражувања, тогаш како резултат на ваквите набљудувања се евидентирани огромен број слични временски серии (сеизмички траги). На пример, кога се проучува локација долга 25 километри и широка 15 километри, каде што, како резултат на работата, секоја трага карактеризира ќелија со димензии 25x25 метри (таквата ќелија се нарекува канта), финалната низа на податоци ќе содржи 600000 траги. Со време на земање примероци од 1 ms и време на снимање од 5 секунди, конечната датотека со податоци ќе биде повеќе од 11 GB, а обемот на оригиналниот „суровина“ може да биде стотици гигабајти.

Како да се работи со нив Волфрам Математика?

Пакет ГеологијаИО

Започна развојот на пакетот прашање на ѕидот VK на групата за поддршка што зборува руски. Благодарение на одговорите на заедницата, решението беше пронајдено многу брзо. И како резултат на тоа, прерасна во сериозен развој. Соодветни Објава на ѕидот на заедницата на Волфрам Дури беше одбележан од модераторите. Во моментов, пакетот поддржува работа со следниве типови податоци кои активно се користат во геолошката индустрија:

1. увоз на податоци од картата во формати ZMAP и IRAP
2. увоз на мерења во бунари со формат LAS
3. формат на влез и излез на сеизмички датотеки СЕГИ

За да го инсталирате пакетот, мора да ги следите упатствата на страницата за преземање на склопениот пакет, т.е. извршете го следниов код во кој било Тетратка за математика:

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

После тоа, пакетот ќе се инсталира во стандардната папка, патеката до која може да се добие на следниов начин:

FileNameJoin[{$UserBasePacletsDirectory, "Repository"}]

Како пример, ќе ги покажеме главните способности на пакетот. Повикот се врши традиционално за пакети на јазикот Волфрам:

Get["GeologyIO`"]

Пакетот е развиен со користење Волфрам работна маса. Ова ви овозможува да ја придружите главната функционалност на пакетот со документација, која во однос на форматот на презентација не се разликува од документацијата на самата Wolfram Mathematica и да го обезбедите пакетот со тест-датотеки за првото запознавање.

Таква датотека, особено, е датотеката „Marmousi.segy“ - ова е синтетички модел на геолошки дел, кој го разви Францускиот институт за нафта. Користејќи го овој модел, програмерите ги тестираат сопствените алгоритми за моделирање на бранови полиња, обработка на податоци, инверзија на сеизмички траги итн. Самиот модел Marmousi е зачуван во складиштето од каде што е преземен самиот пакет. За да ја добиете датотеката, извршете го следниов код:

If[Not[FileExistsQ["Marmousi.segy"]], 
URLDownload["https://wolfr.am/FiQGh7rk", "Marmousi.segy"];]
marmousi = SEGYImport["Marmousi.segy"]

Увези резултат - објект SEGYData


Форматот SEGY вклучува складирање на различни информации за набљудувањата. Прво, ова се текстуални коментари. Тука спаѓаат информации за локацијата на работата, имињата на компаниите што ги извршиле мерењата итн. Во нашиот случај, ова заглавие се повикува со барање со клучот TextHeader. Еве скратено заглавие на текстот:

Short[marmousi["TextHeader"]]

„Збирката на податоци Marmousi беше генерирана во Институтот ... минимална брзина од 1500 m/s и максимум 5500 m/s)“

Можете да го прикажете вистинскиот геолошки модел со пристап до сеизмичките траги користејќи го копчето „траги“ (една од карактеристиките на пакетот е тоа што клучевите не се чувствителни на букви):

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

Моделот Мармузи


Во моментов, пакетот исто така ви овозможува да вчитате податоци во делови од големи датотеки, што овозможува обработка на датотеки чија големина може да достигне десетици гигабајти. Функциите на пакетот вклучуваат и функции за извоз на податоци во .segy и делумно додавање на крајот на датотеката.

Одделно, вреди да се забележи функционалноста на пакетот при работа со сложената структура на датотеките .segy. Бидејќи ви овозможува не само пристап до поединечни траги и заглавија користејќи клучеви и индекси, туку и да ги менувате и потоа да ги запишете во датотека. Многу од техничките детали за имплементацијата на GeologyIO се надвор од опсегот на овој напис и веројатно заслужуваат посебен опис.

Релевантност на спектралната анализа во сеизмичките истражувања

Способноста за увоз на сеизмички податоци во Wolfram Mathematica ви овозможува да користите вградена функционалност за обработка на сигналот за експериментални податоци. Бидејќи секоја сеизмичка трага претставува временска серија, една од главните алатки за нивно проучување е спектралната анализа. Меѓу предусловите за анализа на фреквентниот состав на сеизмичките податоци, можеме да го наведеме, на пример, следново:

1. Различни типови бранови се карактеризираат со различен фреквентен состав. Ова ви овозможува да нагласите корисни бранови и да ги потиснете брановите на пречки.
2. Карпестите својства како што се порозноста и заситеноста може да влијаат на составот на фреквенцијата. Ова овозможува да се идентификуваат карпите со најдобри својства.
3. Слоевите со различни дебелини предизвикуваат аномалии во различни фреквентни опсези.

Третата точка е главната во контекст на оваа статија. Подолу е фрагмент од код за пресметување на сеизмички траги во случај на слој со различна дебелина - модел на клин. Овој модел традиционално се изучува во сеизмичкото истражување за да се анализираат ефектите на пречки кога брановите рефлектирани од многу слоеви се надредени еден на друг.

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

Модел на формација со штипкање


Брзината на бранот внатре во клинот е 4500 m/s, надвор од клинот 4000 m/s, а густината се претпоставува дека е константна 2200 g/cm³. За таков модел пресметуваме коефициенти на рефлексија и сеизмички траги.

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]

Сеизмички траги за клин модел


Редоследот на сеизмички траги прикажани на оваа слика се нарекува сеизмички пресек. Како што можете да видите, неговото толкување може да се изврши и на интуитивно ниво, бидејќи геометријата на рефлектираните бранови јасно одговара на моделот што беше наведен претходно. Ако подетално ги анализирате трагите, ќе забележите дека трагите од 1 до приближно 30 не се разликуваат - отсјајот од покривот на формацијата и од дното не се преклопуваат. Почнувајќи од 31-та трага, рефлексиите почнуваат да се мешаат. И, иако во моделот, коефициентите на рефлексија не се менуваат хоризонтално - сеизмичките траги го менуваат својот интензитет како што се менува дебелината на формацијата.

Да ја разгледаме амплитудата на рефлексијата од горната граница на формацијата. Почнувајќи од 60-та рута, интензитетот на рефлексијата почнува да се зголемува и на 70-та рута станува максимален. Така се манифестира мешањето на брановите од покривот и дното на слоевите, што во некои случаи доведува до значителни аномалии во сеизмичкиот запис.

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]

График на амплитудата на рефлектираниот бран од горниот раб на клинот


Логично е дека кога сигналот е со пониска фреквенција, пречки почнуваат да се појавуваат на големи дебелини на формацијата, а во случај на сигнал со висока фреквенција, пречки се јавуваат на помали дебелини. Следниот фрагмент од код создава сигнал со фреквенции од 35 Hz, 55 Hz и 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]

Збир на изворни сигнали со фреквенции од 35 Hz, 55Hz, 85Hz


Со пресметување на сеизмички траги и исцртување графикони на амплитуди на рефлектираните бранови, можеме да видиме дека за различни фреквенции се забележува аномалија на различни дебелини на формацијата.

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]

Графикони на амплитудите на рефлектираниот бран од горниот раб на клинот за различни фреквенции


Способноста да се извлечат заклучоци за дебелината на формацијата од резултатите од сеизмичките набљудувања е исклучително корисна, бидејќи една од главните задачи во истражувањето на нафтата е да се проценат најперспективните точки за поставување бунар (т.е. оние области каде што е формацијата подебели). Покрај тоа, во геолошкиот дел може да има објекти чија генеза предизвикува остра промена во дебелината на формацијата. Ова ја прави спектралната анализа ефективна алатка за нивно проучување. Во следниот дел од статијата подетално ќе ги разгледаме ваквите геолошки објекти.

Експериментални податоци. Каде ги добивте и што да барате во нив?

Материјалите анализирани во статијата се добиени во Западен Сибир. Регионот, како што веројатно сите без исклучок знаат, е главниот регион за производство на нафта во нашата земја. Активниот развој на наоѓалиштата започна во регионот во 60-тите години на минатиот век. Главниот метод за пребарување на наоѓалишта на нафта е сеизмичко истражување. Интересно е да се погледнат сателитски снимки од оваа територија. Во мал обем, можете да забележите огромен број на мочуришта и езера; со зголемување на картата, можете да ги видите местата за дупчење на бунари во кластери, а со зголемување на картата до граница, можете да ги разликувате и расчистувањата на профилите по кои сеизмички беа извршени набљудувања.

Сателитска слика на мапи на Yandex - градска област Нојабрск


Мрежа од подлоги за бунари на едно од полињата


Нафтаните карпи од Западен Сибир се појавуваат во широк опсег на длабочини - од 1 км до 5 км. Главниот волумен на карпите што содржат нафта бил формиран во времето на Јура и Креда. Периодот Јура на многумина веројатно им е познат од истоимениот филм. Јура клима значително се разликуваше од модерниот. Encyclopedia Britannica има серија палеомапи кои ја карактеризираат секоја хелогична ера.

Денес

Јура период


Ве молиме имајте предвид дека во времето на Јура, територијата на Западен Сибир беше морски брег (земја што ја преминуваа реки и плитко море). Бидејќи климата беше удобна, можеме да претпоставиме дека типичен пејзаж од тоа време изгледал вака:

Јура Сибир


На оваа слика, она што ни е важно не се толку животните и птиците, туку сликата на реката во позадина. Реката е истиот геолошки објект на кој застанавме порано. Факт е дека активноста на реките овозможува да се акумулираат добро сортирани песочник, кои потоа ќе станат резервоар за нафта. Овие акумулации можат да имаат бизарна, сложена форма (како речно корито) и имаат променлива дебелина - во близина на бреговите дебелината е мала, но поблиску до центарот на каналот или во меандерските области се зголемува. Значи, реките формирани во Јура сега се на длабочина од околу три километри и се предмет на потрага по резервоари за нафта.

Експериментални податоци. Обработка и визуелизација

Веднаш да направиме резерва во однос на сеизмичките материјали прикажани во статијата - поради фактот што количината на податоци што се користат за анализата е значителна - само фрагмент од оригиналниот сет на сеизмички траги е вклучен во текстот на статијата. Ова ќе му овозможи на секој да ги репродуцира горенаведените пресметки.

Кога работи со сеизмички податоци, геофизичарот обично користи специјализиран софтвер (има неколку индустриски лидери чии случувања активно се користат, на пример Петрел или Парадигма), што ви овозможува да анализирате различни видови податоци и има удобен графички интерфејс. И покрај сета погодност, овие типови софтвери имаат и свои недостатоци - на пример, имплементацијата на современи алгоритми во стабилни верзии трае многу време, а можностите за автоматизирање на пресметките обично се ограничени. Во таква ситуација, станува многу погодно да се користат компјутерски математички системи и програмски јазици на високо ниво, кои овозможуваат користење на широка алгоритамска база и, во исто време, преземаат многу рутина. Ова е принципот што се користи за работа со сеизмички податоци во 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]

Податоците што се преземаат и увезуваат на овој начин се маршрутите снимени на површина со димензии 10 на 5 километри. Доколку податоците се добијат со метод на тридимензионален сеизмички преглед (брановите се снимаат не по одделни геофизички профили, туку на целата област истовремено), станува возможно да се добијат коцки за сеизмички податоци. Станува збор за тродимензионални објекти, чии вертикални и хоризонтални пресеци овозможуваат детално проучување на геолошката средина. Во разгледаниот пример, се работи со тродимензионални податоци. Можеме да добиеме некои информации од заглавието на текстот, вака

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

C 1 ОВА Е ДЕМО-ДАТОЈКА ЗА ТЕСТ НА ПАКЕТ НА GEOLOGYIO
Ц 2
Ц 3
Ц 4
C 5 ДАТУМ ИМЕ НА КОРИСНИКОТ: КОРИСНИК НА WOLFRAM
C 6 ИМЕ НА АНКЕТА: НЕКАДЕ ВО СИБИР
C 7 ТИП НА ДАТОТЕКА 3D СЕИЗМИЧЕН ВОЛУМЕН
Ц 8
Ц 9
ОПЕГА C10 Z: ПРВИ 2200 М ПОСЛЕДНИ 2400 м

Овој сет на податоци ќе биде доволен за да ги покажеме главните фази на анализа на податоците. Трагите во датотеката се снимаат последователно и секоја од нив изгледа нешто како следната слика - ова е распределбата на амплитудите на рефлектираните бранови долж вертикалната оска (оска на длабочина).

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]

Една од трагите на сеизмичкиот дел


Знаејќи колку траги се наоѓаат во секоја насока на проучуваната област, можете да генерирате тродимензионална низа со податоци и да ја прикажете со помош на функцијата 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]]

XNUMXD слика на коцка за сеизмички податоци (Вертикална оска - длабочина)


Доколку геолошките карактеристики од интерес создаваат интензивни сеизмички аномалии, тогаш може да се користат алатки за визуелизација со транспарентност. „Неважните“ области на снимката може да се направат невидливи, оставајќи ги видливи само аномалии. Во Wolfram Mathematica ова може да се направи со користење Непроѕирност[] и Растер 3Д[].

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]

Слика на коцка за сеизмички податоци со користење на функциите Opacity[] и Raster3D[]

Како и во синтетичкиот пример, некои геолошки граници (слоеви) со променлив релјеф може да се идентификуваат на пресеците на оригиналната коцка.

Главната алатка за спектрална анализа е Фуриеовата трансформација. Со негова помош, можете да го оцените амплитудно-фреквентниот спектар на секоја трага или група на траги. Меѓутоа, по пренесувањето на податоците во доменот на фреквенција, информациите се губат за тоа во кои моменти (читај во кои длабочини) се менува фреквенцијата. За да може да се локализираат промените на сигналот на оската на времето (длабочината), се користат прозорците на Фуриеовата трансформација и разложувањето на брановидовите. Овој напис користи разградување на брановидни бранови. Технологијата за анализа на брановидници почна активно да се користи во сеизмичките истражувања во 90-тите. Предноста во однос на прозорската Фуриеова трансформација се смета за подобра временска резолуција.

Користејќи го следниот фрагмент од код, можете да разложите една од сеизмичките траги на поединечни компоненти:

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]

Распаѓање на трага на компоненти


За да се процени како се распределува енергијата на рефлексијата во различни времиња на пристигнување на брановите, се користат скалограми (аналогни на спектрограм). Како по правило, во пракса нема потреба да се анализираат сите компоненти. Вообичаено, се избираат компоненти со ниска, средна и висока фреквенција.

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]

Скалограм. Резултат на функцијата WaveletScalogram[]


Јазикот на Волфрам ја користи функцијата за трансформација на брановидни линии ContinuousWaveletTransform[]. И примената на оваа функција на целиот сет на траги ќе се изврши со помош на функцијата Табела[]. Тука вреди да се забележи една од силните страни на Wolfram Mathematica - способноста да се користи паралелизација Паралелна табела[]. Во горниот пример, нема потреба од паралелизација - обемот на податоци не е голем, но кога се работи со експериментални множества на податоци што содржат стотици илјади траги, ова е неопходност.

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

По примената на функцијата ContinuousWaveletTransform[] Се појавуваат нови збирки податоци што одговараат на избраните фреквенции. Во горниот пример, овие фреквенции се: 38Hz, 33Hz, 27Hz. Изборот на фреквенции најчесто се врши врз основа на тестирање - тие добиваат ефективни карти за различни комбинации на фреквенции и ја избираат најинформативната од гледна точка на геолог.

Ако треба да ги споделите резултатите со колегите или да ги доставите на клиентот, можете да ја користите функцијата SEGYExport[] од пакетот 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];

Со три од овие коцки (компоненти со ниска фреквенција, средна фреквенција и висока фреквенција), RGB мешањето обично се користи за да се визуелизираат податоците заедно. Секоја компонента е доделена своја боја - црвена, зелена, сина. Во Wolfram Mathematica ова може да се направи со помош на функцијата ColorCombine[].

Резултатот се слики од кои може да се направи геолошка интерпретација. Меандрите што се забележани на делот овозможуваат да се разграничат палеоканали, кои поверојатно се акумулации и содржат резерви на нафта. Пребарувањето и анализата на современите аналози на таков речен систем ни овозможува да ги одредиме најперспективните делови од меандрите. Самите канали се карактеризираат со дебели слоеви на добро сортиран песочник и се добар резервоар за нафта. Областите надвор од аномалиите на „тантела“ се слични на современите наноси на поплавни рамнини. Наслагите во поплавните рамнини главно се претставени со глинести карпи и дупчењето во овие зони ќе биде неефикасно.

RGB парче од коцката за податоци. Во центарот (малку лево од центарот) можете да ја следите меандрираната река.

RGB парче од коцката за податоци. На левата страна можете да ја проследите меандрираната река.


Во некои случаи, квалитетот на сеизмичките податоци овозможуваат значително појасни слики. Ова зависи од методологијата на работа на терен, опремата што ја користи алгоритмот за намалување на бучавата. Во такви случаи, не се видливи само фрагменти од речните системи, туку и цели продолжени палеореки.

RGB мешање на три компоненти на коцка за сеизмички податоци (хоризонтално парче). Длабочина околу 2 км.

Сателитска слика на реката Волга во близина на Саратов

Заклучок

Wolfram Mathematica ви овозможува да ги анализирате сеизмичките податоци и да решавате применети проблеми поврзани со истражување на минерали, а пакетот GeologyIO го прави овој процес поудобен. Структурата на сеизмичките податоци е таква што со користење на вградени методи за забрзување на пресметките (Паралелна табела[], ParallelDo[],…) е многу ефикасен и ви овозможува да обработувате големи количини на податоци. Во голема мера, ова е олеснето со карактеристиките за складирање податоци на пакетот GeologyIO. Патем, пакетот може да се користи не само во областа на применетите сеизмички истражувања. Скоро исти типови на податоци се користат во радарите кои продираат по земја и сеизмологијата. Остави коментар.

Извор: www.habr.com