Благодаря на автора на блога за доклада му
въведение
Тази бележка е написана след конференцията и съдържа резюме на доклада, който изнесох. Събитието се състоя през юни в Санкт Петербург. Като се има предвид, че работя на една пресечка от сайта на конференцията, нямаше как да не присъствам на това събитие. През 2016 и 2017 г. слушах доклади от конференции, а тази година изнесох презентация. Първо, появи се интересна (струва ми се) тема, която развиваме , и второ, след дълго проучване на законодателството на Руската федерация по отношение на санкционната политика, в предприятието, където работя, се появиха цели два лиценза .
Преди да премина към темата на моето изказване, искам да отбележа добрата организация на събитието. Страницата за посещение на конференцията използва изображение на Казанската катедрала. Катедралата е една от основните забележителности на Санкт Петербург и се вижда много добре от залата, в която се проведе конференцията.
На входа на Санкт Петербургския държавен икономически университет участниците бяха посрещнати от асистенти сред студентите - те не им позволиха да се изгубят. По време на регистрацията бяха раздадени малки сувенири (играчка - мигащ шип, химикал, стикери със символи на Wolfram). Обядът и кафе паузите също бяха включени в програмата на конференцията. Вече отбелязах за вкусно кафе и пайове на стената на групата - готвачите са страхотни. С тази уводна част искам да подчертая, че самото събитие, неговият формат и място вече носят положителни емоции.
Докладът, който подготвихме аз и Кирил Белов, се казва „Използване на 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 следи. С време за вземане на проби от 5 ms и време за запис от 11 секунди, крайният файл с данни ще бъде повече от XNUMX GB, а обемът на оригиналния „суров“ материал може да бъде стотици гигабайти.
Как да работим с тях ?
пакет
Разработката на пакета започна на стената на VK на рускоговорящата група за поддръжка. Благодарение на отговорите на общността решение беше намерено много бързо. И в резултат прерасна в сериозно развитие. Кореспондентски Дори беше отбелязано от модераторите. В момента пакетът поддържа работа със следните типове данни, които се използват активно в геоложката индустрия:
- импортиране на картографски данни във формати ZMAP и IRAP
- импортиране на измервания в кладенци във формат LAS
- вход и изход на формат на сеизмични файлове
За да инсталирате пакета, трябва да следвате инструкциите на страницата за изтегляне на сглобения пакет, т.е. изпълнете следния код във всеки :
If[PacletInformation["GeologyIO"] === {}, PacletInstall[URLDownload[
"https://wolfr.am/FiQ5oFih",
FileNameJoin[{CreateDirectory[], "GeologyIO-0.2.2.paclet"}]
]]]След това пакетът ще бъде инсталиран в папката по подразбиране, пътят до който може да бъде получен, както следва:
FileNameJoin[{$UserBasePacletsDirectory, "Repository"}]Като пример ще демонстрираме основните възможности на пакета. Обаждането се извършва традиционно за пакети на езика Wolfram:
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 ви позволява да използвате вградена функционалност за обработка на сигнали за експериментални данни. Тъй като всяка сеизмична следа представлява времева серия, един от основните инструменти за тяхното изследване е спектралният анализ. Сред предпоставките за анализ на честотния състав на сеизмичните данни можем да назовем например следното:
- Различните видове вълни се характеризират с различен честотен състав. Това ви позволява да подчертавате полезни вълни и да потискате вълните на смущения.
- Свойствата на скалата като порьозност и наситеност могат да повлияят на честотния състав. Това прави възможно идентифицирането на скали с най-добри свойства.
- Слоевете с различна дебелина причиняват аномалии в различни честотни диапазони.
Третата точка е основната в контекста на тази статия. По-долу е даден кодов фрагмент за изчисляване на сеизмични следи в случай на слой с различна дебелина - клиновиден модел. Този модел традиционно се изучава при сеизмични проучвания, за да се анализират ефектите на смущение, когато вълните, отразени от много слоеве, се наслагват една върху друга.
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, 55 Hz, 85 Hz
Чрез изчисляване на сеизмични следи и начертаване на графики на амплитудите на отразените вълни можем да видим, че за различни честоти се наблюдава аномалия при различни дебелини на пласта.
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 км. Основният обем от скали, съдържащи нефт, се е образувал през периода юра и креда. Джурасик периодът вероятно е познат на мнозина от едноименния филм. се различава значително от съвременния. Енциклопедия Британика има поредица от палеомапи, които характеризират всяка хелогична ера.
настояще
юрски период
Моля, имайте предвид, че през юрските времена територията на Западен Сибир е била морско крайбрежие (земя, пресечена от реки и плитко море). Тъй като климатът беше удобен, можем да предположим, че типичният пейзаж от онова време изглеждаше така:
Джурасик Сибир
В тази картина за нас са важни не толкова животните и птиците, а изображението на реката на заден план. Реката е същият геоложки обект, на който спряхме по-рано. Факт е, че дейността на реките позволява да се натрупват добре сортирани пясъчници, които след това ще се превърнат в резервоар за нефт. Тези резервоари могат да имат странна, сложна форма (като речно корито) и имат променлива дебелина - близо до бреговете дебелината е малка, но по-близо до центъра на канала или в зоните на меандра се увеличава. И така, реките, образувани през юрския период, сега са на дълбочина около три километра и са обект на търсене на нефтени резервоари.
Експериментални данни. Обработка и визуализация
Нека веднага направим уговорка по отношение на показаните в статията сеизмични материали - поради факта, че количеството данни, използвани за анализа, е значително - в текста на статията е включен само фрагмент от оригиналния набор от сеизмични следи. Това ще позволи на всеки да възпроизведе горните изчисления.
Когато работи със сеизмични данни, геофизикът обикновено използва специализиран софтуер (има няколко лидери в индустрията, чиито разработки се използват активно, например Petrel или Paradigm), който ви позволява да анализирате различни видове данни и има удобен графичен интерфейс. Въпреки всички удобства, този тип софтуер има и своите недостатъци - например внедряването на съвременни алгоритми в стабилни версии отнема много време, а възможностите за автоматизиране на изчисленията обикновено са ограничени. В такава ситуация става много удобно да се използват системи за компютърна математика и езици за програмиране на високо ниво, които позволяват използването на широка алгоритмична база и в същото време поемат много рутина. Това е принципът, използван за работа със сеизмични данни в 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] // TableFormC 1 ТОВА Е ДЕМО ФАЙЛ ЗА ТЕСТ НА ПАКЕТА GEOLOGYIO
C 2
C 3
C 4
C 5 ДАТА ПОТРЕБИТЕЛСКО ИМЕ: WOLFRAM ПОТРЕБИТЕЛ
C 6 ИМЕ НА АНКЕТАТА: НЯКЪДЕ В СИБИР
C 7 ТИП НА ФАЙЛ 3D СЕИЗМИЧЕН ОБЕМ
C 8
C 9
C10 Z ДИАПАЗОН: ПЪРВИ 2200M ПОСЛЕДНИ 2400M
Този набор от данни ще ни бъде достатъчен, за да демонстрираме основните етапи на анализа на данните. Следите във файла се записват последователно и всяка от тях изглежда като следната фигура - това е разпределението на амплитудите на отразените вълни по вертикалната ос (оста на дълбочината).
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 това може да стане с помощта на и .
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]Скалограма. Резултат от функцията
Езикът Wolfram използва функцията за вълнова трансформация . И прилагането на тази функция към целия набор от следи ще се извърши с помощта на функцията . Тук си струва да се отбележи една от силните страни на Wolfram Mathematica - възможността за използване на паралелизация . В горния пример няма нужда от паралелизиране - обемът на данните не е голям, но при работа с експериментални масиви от данни, съдържащи стотици хиляди следи, това е необходимост.
tracesCWD=Table[Map[Hilbert[#,0]&,Re[ContinuousWaveletTransform[traces[[i]]][[1]]][[{13,15,18}]]],{i,1,Length@traces}]; След прилагане на функцията Появяват се нови набори от данни, съответстващи на избраните честоти. В примера по-горе тези честоти са: 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 това може да стане с помощта на функцията .
Резултатът са изображения, от които може да се направи геоложка интерпретация. Меандрите, които са записани на разреза, позволяват да се очертаят палеоканали, които е по-вероятно да бъдат резервоари и да съдържат запаси от нефт. Търсенето и анализът на съвременни аналози на такава речна система ни позволява да определим най-обещаващите части на меандрите. Самите канали се характеризират с дебели слоеве от добре сортирани пясъчници и са добър резервоар за нефт. Зоните извън "дантелените" аномалии са подобни на съвременните наводнени отлагания. Отлаганията в заливните равнини са представени главно от глинести скали и сондажите в тези зони ще бъдат неефективни.
RGB срез на куба с данни. В центъра (малко вляво от центъра) можете да проследите криволичещата река.
RGB срез на куба с данни. От лявата страна можете да проследите криволичещата река.
В някои случаи качеството на сеизмичните данни позволява значително по-ясни изображения. Това зависи от методологията на теренната работа, оборудването, използвано от алгоритъма за намаляване на шума. В такива случаи се виждат не само фрагменти от речни системи, но и цели разширени палеореки.
RGB смесване на три компонента на куб със сеизмични данни (хоризонтален срез). Дълбочина около 2 км.
Сателитна снимка на река Волга близо до Саратов
Заключение
Wolfram Mathematica ви позволява да анализирате сеизмични данни и да решавате приложни проблеми, свързани с проучването на минерали, а пакетът GeologyIO прави този процес по-удобен. Структурата на сеизмичните данни е такава, че използването на вградени методи за ускоряване на изчисленията (, ,…) е много ефективен и ви позволява да обработвате големи количества данни. До голяма степен това се улеснява от функциите за съхранение на данни на пакета GeologyIO. Между другото, пакетът може да се използва не само в областта на приложното сеизмично проучване. Почти същите типове данни се използват в наземните радари и сеизмологията. Ако имате предложения как да подобрите резултата, кои алгоритми за анализ на сигнала от арсенала на Wolfram Mathematica са приложими за такива данни или ако имате критични коментари, моля Оставете коментар.
Източник: www.habr.com