Wolfram Mathematica geofizikoje

Ačiū tinklaraščio autorei Antonas Ekimenko už jo pranešimą

įvedimas

Ši pastaba buvo parašyta po konferencijos Wolfram Rusijos technologijų konferencija ir yra mano pateiktos ataskaitos santrauka. Renginys vyko birželio mėnesį Sankt Peterburge. Atsižvelgiant į tai, kad dirbu už kvartalo nuo konferencijos vietos, negalėjau nedalyvauti šiame renginyje. 2016 ir 2017 metais klausiausi konferencijos pranešimų, o šiemet skaičiau pranešimą. Pirma, atsirado įdomi (man atrodo) tema, su kuria vystome Kirilas Belovas, ir antra, po ilgo Rusijos Federacijos teisės aktų dėl sankcijų politikos tyrimo, įmonėje, kurioje aš dirbu, atsirado net dvi licencijos. Wolfram Mathematica.

Prieš pereinant prie savo kalbos temos, norėčiau atkreipti dėmesį į gerą renginio organizavimą. Konferencijos apsilankymo puslapyje naudojamas Kazanės katedros vaizdas. Katedra yra viena pagrindinių Sankt Peterburgo lankytinų vietų ir labai aiškiai matoma iš salės, kurioje vyko konferencija.

Wolfram Mathematica geofizikoje

Prie įėjimo į Sankt Peterburgo valstybinį ekonomikos universitetą dalyvius pasitiko asistentai iš studentų tarpo – jie neleido pasiklysti. Registracijos metu buvo išdalinti nedideli suvenyrai (žaislas - mirksintis smaigalys, rašiklis, lipdukai su Wolfram simbolika). Pietūs ir kavos pertraukėlės taip pat buvo įtrauktos į konferencijos tvarkaraštį. Apie skanią kavą ir pyragus jau pažymėjau ant grupės sienos – šefai puikūs. Šia įžangine dalimi noriu pabrėžti, kad pats renginys, jo formatas ir vieta jau kelia teigiamas emocijas.

Mano ir Kirilo Belovo parengtas pranešimas vadinasi „Wolfram Mathematica naudojimas taikomosios geofizikos problemoms spręsti. Seisminių duomenų spektrinė analizė arba „kur tekėjo senovės upės“. Ataskaitos turinys apima dvi dalis: pirma, algoritmų, kuriuos galima rasti, naudojimą Wolfram Mathematica geofiziniams duomenims analizuoti, ir, antra, taip geofizinius duomenis įdėti į Wolfram Mathematica.

Seisminis tyrinėjimas

Pirmiausia turite padaryti trumpą ekskursiją į geofiziką. Geofizika yra mokslas, tiriantis fizines uolienų savybes. Na, o kadangi uolienos turi skirtingas savybes: elektrines, magnetines, elastines, yra atitinkami geofizikos metodai: elektrinis žvalgymas, magnetinis žvalgymas, seisminis žvalgymas... Šio straipsnio kontekste plačiau aptarsime tik seisminę žvalgybą. Seisminiai tyrinėjimai yra pagrindinis naftos ir dujų paieškos būdas. Metodas pagrįstas tampriųjų virpesių sužadinimu ir vėlesniu tiriamą sritį sudarančių uolienų atsako registravimu. Vibracijos sužadinamos sausumoje (su dinamitu arba nesprogiais elastingų virpesių vibracijos šaltiniais) arba jūroje (su pneumatiniais ginklais). Elastiniai virpesiai sklinda per uolienų masę, lūžta ir atsispindi skirtingų savybių sluoksnių ribose. Atsispindinčios bangos grįžta į paviršių ir jas fiksuoja sausumoje esantys geofonai (dažniausiai elektrodinaminiai prietaisai, pagrįsti ritėje pakabinto magneto judėjimu) arba hidrofonai jūroje (remiantis pjezoelektriniu efektu). Pagal bangų atvykimo laiką galima spręsti apie geologinių sluoksnių gylius.

Seisminė laivų vilkimo įranga
Wolfram Mathematica geofizikoje

Pneumatinis pistoletas sužadina elastines vibracijas
Wolfram Mathematica geofizikoje

Bangos praeina per uolienų masę ir yra įrašomos hidrofonais
Wolfram Mathematica geofizikoje

Geofizinių tyrimų laivas „Ivanas Gubkinas“ prieplaukoje prie Blagoveščenskio tilto Sankt Peterburge
Wolfram Mathematica geofizikoje

Seisminio signalo modelis

Uolos turi skirtingas fizines savybes. Seisminiams tyrinėjimams pirmiausia svarbios tamprumo savybės – elastingų virpesių sklidimo greitis ir tankis. Jei du sluoksniai turi tas pačias ar panašias savybes, banga „nepastebės“ ribos tarp jų. Jei bangų greičiai sluoksniuose skiriasi, tada atsispindės ties sluoksnių riba. Kuo didesnis savybių skirtumas, tuo intensyvesnis atspindys. Jo intensyvumas bus nustatomas pagal atspindžio koeficientą (rc):

Wolfram Mathematica geofizikoje

kur ρ – uolienų tankis, ν – bangos greitis, 1 ir 2 nurodo viršutinį ir apatinį sluoksnius.

Vienas iš paprasčiausių ir dažniausiai naudojamų seisminių signalų modelių yra konvoliucijos modelis, kai užfiksuotas seisminis pėdsakas vaizduojamas kaip atspindžio koeficientų sekos konvoliucijos rezultatas su zondavimo impulsu:

Wolfram Mathematica geofizikoje

kur s(t) — seisminis pėdsakas, t.y. viskas, kas buvo įrašyta hidrofonu arba geofonu per fiksuotą įrašymo laiką, w(t) - pneumatinio pistoleto generuojamas signalas, n(t) - atsitiktinis triukšmas.

Kaip pavyzdį apskaičiuokime sintetinį seisminį pėdsaką. Kaip pradinį signalą naudosime Rickerio impulsą, plačiai naudojamą seisminiuose tyrimuose.

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]

Pradinis seisminis impulsas
Wolfram Mathematica geofizikoje

Mes nustatysime dvi ribas 300 ms ir 600 ms gylyje, o atspindžio koeficientai bus atsitiktiniai skaičiai

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

Atspindžio koeficientų seka
Wolfram Mathematica geofizikoje

Apskaičiuokime ir parodykime seisminį pėdsaką. Kadangi atspindžio koeficientai turi skirtingus ženklus, gauname du kintamus seisminio pėdsako atspindžius.

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

Imituojamas takelis
Wolfram Mathematica geofizikoje

Šiam pavyzdžiui būtina padaryti išlygą – realiai sluoksnių gylis nustatomas, žinoma, metrais, o seisminio pėdsako apskaičiavimas vyksta laiko domenui. Teisingiau būtų nustatyti gylius metrais ir skaičiuoti atvykimo laiką žinant greičius sluoksniuose. Šiuo atveju sluoksnius iš karto nustatau laiko ašyje.

Jei mes kalbame apie lauko tyrimus, tada tokių stebėjimų rezultatas yra daugybė panašių laiko eilučių (seisminių pėdsakų). Pavyzdžiui, tiriant 25 km ilgio ir 15 km pločio aikštelę, kur darbo rezultatas kiekvienas pėdsakas apibūdina 25x25 metrų dydžio langelį (toks langelis vadinamas šiukšliadėže), galutiniame duomenų masyve bus 600000 1 pėdsakų. Esant 5 ms atrankos laikui ir 11 sekundžių įrašymo laikui, galutinis duomenų failas bus didesnis nei XNUMX GB, o originalios „žaliavos“ tūris gali siekti šimtus gigabaitų.

Kaip su jais dirbti Wolfram Mathematica?

Pakuotė GeologijaIO

Pradėtas paketo kūrimas klausimas ant rusakalbių paramos grupės VK sienos. Dėl bendruomenės atsakymų sprendimas buvo rastas labai greitai. Ir dėl to tai išaugo į rimtą vystymąsi. Atitinkamas Wolfram bendruomenės sienos įrašas Tai netgi pažymėjo moderatoriai. Šiuo metu paketas palaiko darbą su šiais duomenų tipais, kurie aktyviai naudojami geologijos pramonėje:

  1. žemėlapių duomenų importas ZMAP ir IRAP formatais
  2. matavimų importas LAS formato šuliniuose
  3. seisminių failų formato įvestis ir išvestis SEGY

Norėdami įdiegti paketą, turite vadovautis surinkto paketo atsisiuntimo puslapyje pateiktomis instrukcijomis, t.y. paleiskite šį kodą bet kuriame Mathematica sąsiuvinis:

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

Po to paketas bus įdiegtas numatytajame aplanke, kurio kelią galima gauti taip:

FileNameJoin[{$UserBasePacletsDirectory, "Repository"}]

Kaip pavyzdį parodysime pagrindines paketo galimybes. Tradiciškai skambinama paketams Wolfram kalba:

Get["GeologyIO`"]

Paketas sukurtas naudojant Wolfram Workbench. Tai leidžia prie pagrindinio paketo funkcionalumo pridėti dokumentaciją, kuri savo pateikimo formatu nesiskiria nuo pačios Wolfram Mathematica dokumentacijos, ir pateikti paketą su bandomaisiais failais pirmajai pažinčiai.

Wolfram Mathematica geofizikoje

Wolfram Mathematica geofizikoje

Toks failas visų pirma yra failas „Marmousi.segy“ - tai sintetinis geologinio skyriaus modelis, kurį sukūrė Prancūzijos naftos institutas. Naudodamiesi šiuo modeliu, kūrėjai išbando savo algoritmus bangų lauko modeliavimui, duomenų apdorojimui, seisminių pėdsakų inversijai ir kt. Pats Marmousi modelis yra saugomas saugykloje, iš kurios buvo atsiųstas pats paketas. Norėdami gauti failą, paleiskite šį kodą:

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

Importo rezultatas – SEGYData objektas
Wolfram Mathematica geofizikoje

SEGY formatas apima įvairios informacijos apie stebėjimus saugojimą. Pirma, tai yra tekstiniai komentarai. Tai apima informaciją apie darbų atlikimo vietą, matavimus atlikusių įmonių pavadinimus ir kt. Mūsų atveju ši antraštė iškviečiama užklausa su TextHeader raktu. Štai sutrumpinta teksto antraštė:

Short[marmousi["TextHeader"]]

„Marmousi duomenų rinkinys buvo sukurtas institute ... mažiausias greitis 1500 m/s ir maksimalus 5500 m/s)“

Galite parodyti tikrąjį geologinį modelį, pasiekdami seisminius pėdsakus naudodami klavišą „traces“ (viena iš paketo ypatybių yra ta, kad raktai neskiria didžiųjų ir mažųjų raidžių):

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

Modelis Marmousi
Wolfram Mathematica geofizikoje

Šiuo metu paketas taip pat leidžia įkelti duomenis dalimis iš didelių failų, todėl galima apdoroti failus, kurių dydis gali siekti keliasdešimt gigabaitų. Paketo funkcijos taip pat apima duomenų eksportavimo į .segy ir dalinio pridėjimo prie failo pabaigos funkcijas.

Atskirai verta atkreipti dėmesį į paketo funkcionalumą dirbant su sudėtinga .segy failų struktūra. Kadangi tai leidžia ne tik pasiekti atskirus pėdsakus ir antraštes naudojant raktus ir indeksus, bet ir juos pakeisti bei įrašyti į failą. Daugelis techninių GeologyIO įgyvendinimo detalių nepatenka į šio straipsnio taikymo sritį ir tikriausiai nusipelno atskiro aprašymo.

Spektrinės analizės svarba seisminiuose tyrimuose

Galimybė importuoti seisminius duomenis į Wolfram Mathematica leidžia naudoti integruotą signalų apdorojimo funkciją eksperimentiniams duomenims. Kadangi kiekvienas seisminis pėdsakas reiškia laiko eilutę, viena iš pagrindinių jų tyrimo priemonių yra spektrinė analizė. Tarp būtinųjų seisminių duomenų dažninės sudėties analizės sąlygų galime paminėti, pavyzdžiui:

  1. Skirtingiems bangų tipams būdinga skirtinga dažnių sudėtis. Tai leidžia išryškinti naudingas bangas ir slopinti trukdžių bangas.
  2. Uolienų savybės, tokios kaip poringumas ir prisotinimas, gali paveikti dažnio sudėtį. Tai leidžia identifikuoti geriausių savybių turinčias uolienas.
  3. Skirtingo storio sluoksniai sukelia anomalijas skirtinguose dažnių diapazonuose.

Trečias punktas yra pagrindinis šio straipsnio kontekste. Žemiau pateikiamas kodo fragmentas, skirtas seisminiams pėdsakams apskaičiuoti įvairaus storio sluoksnio atveju – pleišto modelis. Šis modelis tradiciškai tiriamas atliekant seisminius tyrimus, siekiant analizuoti trukdžių efektus, kai iš daugelio sluoksnių atsispindinčios bangos yra viena ant kitos.

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

Išspausto darinio modelis
Wolfram Mathematica geofizikoje

Bangos greitis pleišto viduje yra 4500 m/s, už pleišto 4000 m/s, o tankis laikomas pastoviu 2200 g/cm³. Tokiam modeliui apskaičiuojame atspindžio koeficientus ir seisminius pėdsakus.

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]

Seisminiai pėdsakai pleišto modeliui
Wolfram Mathematica geofizikoje

Šiame paveikslėlyje parodyta seisminių pėdsakų seka vadinama seismine sekcija. Kaip matote, jo interpretacija taip pat gali būti atlikta intuityviu lygmeniu, nes atspindėtų bangų geometrija aiškiai atitinka anksčiau nurodytą modelį. Jei išanalizuosite pėdsakus išsamiau, pastebėsite, kad pėdsakai nuo 1 iki maždaug 30 nesiskiria – atspindys nuo darinio stogo ir iš apačios vienas kito nepersidengia. Pradedant nuo 31 pėdsako, atspindžiai pradeda trukdyti. Ir, nors modelyje atspindžio koeficientai nesikeičia horizontaliai – seisminiai pėdsakai keičia savo intensyvumą, keičiantis darinio storiui.

Panagrinėkime atspindžio amplitudę nuo viršutinės darinio ribos. Pradedant nuo 60-ojo maršruto, atspindžio intensyvumas pradeda didėti, o 70-ajame maršrute tampa didžiausias. Taip pasireiškia bangų trukdžiai iš stogo ir sluoksnių apačios, kai kuriais atvejais sukeliantys reikšmingas seisminio įrašo anomalijas.

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]

Atsispindėjusios bangos amplitudės grafikas nuo viršutinio pleišto krašto
Wolfram Mathematica geofizikoje

Logiška, kad kai signalas yra žemesnio dažnio, trukdžiai pradeda atsirasti esant dideliam formavimo storiui, o aukšto dažnio signalo atveju trikdžiai atsiranda esant mažesniam storiui. Šis kodo fragmentas sukuria signalą, kurio dažniai yra 35 Hz, 55 Hz ir 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]

Šaltinio signalų rinkinys, kurių dažniai yra 35 Hz, 55 Hz, 85 Hz
Wolfram Mathematica geofizikoje

Apskaičiuojant seisminius pėdsakus ir nubraižant atsispindėjusių bangų amplitudės grafikus, matome, kad esant skirtingiems dažniams anomalija stebima esant skirtingam formavimo storiui.

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]

Atsispindėjusios bangos nuo viršutinio pleišto krašto amplitudės grafikai skirtingiems dažniams
Wolfram Mathematica geofizikoje

Galimybė daryti išvadas apie darinio storį iš seisminių stebėjimų rezultatų yra itin naudinga, nes naftos žvalgyboje vienas pagrindinių uždavinių yra įvertinti perspektyviausius gręžinio klojimo taškus (t. storesnis). Be to, geologiniame pjūvyje gali būti objektų, kurių genezė sukelia staigų darinio storio pasikeitimą. Dėl to spektrinė analizė yra veiksminga priemonė joms tirti. Kitoje straipsnio dalyje tokius geologinius objektus nagrinėsime plačiau.

Eksperimentiniai duomenys. Kur juos gavai ir ko juose ieškoti?

Straipsnyje analizuojamos medžiagos gautos Vakarų Sibire. Regionas, kaip tikriausiai žino visi be išimties, yra pagrindinis mūsų šalies naftos gavybos regionas. Aktyvi indėlių plėtra regione prasidėjo praėjusio amžiaus 60-aisiais. Pagrindinis naftos telkinių paieškos būdas yra seisminis tyrinėjimas. Įdomu pažiūrėti šios teritorijos palydovines nuotraukas. Nedideliame mastelyje galite pastebėti daugybę pelkių ir ežerų, padidinus žemėlapį matomos kasetinių gręžinių gręžimo vietos, o padidinus žemėlapį iki ribos, taip pat galima atskirti profilių proskynas, išilgai kurių yra seisminis. buvo atlikti stebėjimai.

Palydovinis „Yandex“ žemėlapių vaizdas – Nojabrsko miesto sritis
Wolfram Mathematica geofizikoje

Šulinių trinkelių tinklas viename iš laukų
Wolfram Mathematica geofizikoje

Naftingos Vakarų Sibiro uolienos yra labai įvairiose gylyje - nuo 1 km iki 5 km. Pagrindinis naftos turinčių uolienų tūris susidarė juros ir kreidos laikais. Juros periodas tikriausiai daugeliui žinomas iš to paties pavadinimo filmo. Juros periodo klimatas gerokai skyrėsi nuo šiuolaikinio. Encyclopedia Britannica turi daugybę paleomapų, apibūdinančių kiekvieną heloginę erą.

Šiandien
Wolfram Mathematica geofizikoje
Juros periodas
Wolfram Mathematica geofizikoje

Atkreipkite dėmesį, kad juros laikais Vakarų Sibiro teritorija buvo jūros pakrantė (žemė, kurią kerta upės ir sekli jūra). Kadangi klimatas buvo patogus, galime manyti, kad tipiškas to meto kraštovaizdis atrodė taip:

Juros periodo Sibiras
Wolfram Mathematica geofizikoje

Šiame paveikslėlyje mums svarbūs ne tiek gyvūnai ir paukščiai, kiek upės vaizdas fone. Upė yra tas pats geologinis objektas, prie kurio sustojome anksčiau. Faktas yra tas, kad upių veikla leidžia kauptis gerai išrūšiuotiems smiltainiams, kurie vėliau taps naftos rezervuaru. Šie rezervuarai gali būti keistos, sudėtingos formos (kaip upės vaga) ir įvairaus storio – prie krantų storis mažas, bet arčiau vagos centro ar vingiuotose vietose didėja. Taigi, juros periode susiformavusios upės dabar yra maždaug trijų kilometrų gylyje ir yra naftos telkinių paieškos objektas.

Eksperimentiniai duomenys. Apdorojimas ir vizualizacija

Iš karto padarykime išlygą dėl straipsnyje pateiktų seisminių medžiagų – dėl to, kad analizei panaudotų duomenų kiekis yra didelis – į straipsnio tekstą įtrauktas tik originalaus seisminių pėdsakų rinkinio fragmentas. Tai leis bet kam atkurti aukščiau pateiktus skaičiavimus.

Dirbdamas su seisminiais duomenimis, geofizikas dažniausiai naudoja specializuotą programinę įrangą (yra keli pramonės lyderiai, kurių plėtra aktyviai naudojama, pavyzdžiui, Petrel ar Paradigm), kuri leidžia analizuoti įvairaus tipo duomenis ir turi patogią grafinę sąsają. Nepaisant visų patogumų, tokio tipo programinė įranga turi ir trūkumų – pavyzdžiui, šiuolaikinių algoritmų įdiegimas stabiliose versijose užima daug laiko, o skaičiavimų automatizavimo galimybės dažniausiai yra ribotos. Esant tokiai situacijai, tampa labai patogu naudotis kompiuterinėmis matematikos sistemomis ir aukšto lygio programavimo kalbomis, kurios leidžia naudotis plačia algoritmine baze ir tuo pačiu įgauna daug rutinos. Tai yra principas, naudojamas dirbant su seisminiais duomenimis Wolfram Mathematica. Interaktyviam darbui su duomenimis netikslinga rašyti turtingą funkcionalumą – svarbiau užtikrinti įkėlimą iš visuotinai priimto formato, pritaikant jiems norimus algoritmus ir įkėlimą atgal į išorinį formatą.

Pagal siūlomą schemą įkelsime originalius seisminius duomenis ir juos parodysime 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]

Tokiu būdu atsisiunčiami ir importuojami duomenys yra maršrutai, užfiksuoti 10 x 5 kilometrų plote. Jei duomenys gauti trimačiu seisminio tyrimo metodu (bangos fiksuojamos ne pagal atskirus geofizinius profilius, o per visą plotą vienu metu), atsiranda galimybė gauti seisminių duomenų kubus. Tai trimačiai objektai, kurių vertikalios ir horizontalios pjūviai leidžia detaliai ištirti geologinę aplinką. Nagrinėjamame pavyzdyje kalbame apie trimačius duomenis. Tam tikros informacijos galime gauti iš teksto antraštės, pavyzdžiui, ši

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

C 1 TAI DEMO FAILAS, skirtas GEOLOGYIO PAKETINIMO TESTYMUI
C 2
C 3
C 4
C 5 DATA NAUDOTOJO VARDAS: WOLFRAM NAUDOTOJAS
C 6 TYRIMO PAVADINIMAS: KAŽkur SIBIRE
C 7 FAILOS TIPAS 3D SEISMINIS TŪRIS
C 8
C 9
C10 Z RANGA: PIRMA 2200M PASKUTINĖ 2400M

Šio duomenų rinkinio mums pakaks, kad parodytume pagrindinius duomenų analizės etapus. Pėdsakai faile įrašomi nuosekliai ir kiekvienas iš jų atrodo panašiai kaip toliau pateiktame paveikslėlyje – tai atspindėtų bangų amplitudės pasiskirstymas išilgai vertikalios ašies (gylio ašies).

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]

Vienas iš seisminių pjūvių pėdsakų
Wolfram Mathematica geofizikoje

Žinodami, kiek pėdsakų yra kiekviena tiriamos srities kryptimi, galite sugeneruoti trimatį duomenų masyvą ir parodyti jį naudodami funkciją 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 seisminių duomenų kubo vaizdas. (Vertikali ašis – gylis)
Wolfram Mathematica geofizikoje

Jei dominančios geologinės ypatybės sukuria intensyvias seismines anomalijas, galima naudoti skaidrias vizualizacijos priemones. „Nesvarbios“ įrašo sritys gali būti nematomos, todėl matomos tik anomalijas. Wolfram Mathematica tai galima padaryti naudojant Nepermatomumas[] и 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]

Seisminių duomenų kubo vaizdas naudojant Opacity[] ir Raster3D[] funkcijas Wolfram Mathematica geofizikoje

Kaip ir sintetiniame pavyzdyje, originalaus kubo atkarpose galima nustatyti kai kurias geologines ribas (sluoksnius) su kintamu reljefu.

Pagrindinis spektrinės analizės įrankis yra Furjė transformacija. Su jo pagalba galite įvertinti kiekvieno pėdsako ar pėdsakų grupės amplitudės-dažnio spektrą. Tačiau perkėlus duomenis į dažnių sritį, prarandama informacija apie tai, kokiu metu (skaitykite, kokiame gylyje) keičiasi dažnis. Kad būtų galima lokalizuoti signalo pokyčius laiko (gylio) ašyje, naudojama langinė Furjė transformacija ir bangelių skaidymas. Šiame straipsnyje naudojamas bangelių skaidymas. Bangų analizės technologija pradėta aktyviai naudoti seisminiuose tyrimuose 90-aisiais. Privalumas prieš languotą Furjė transformaciją laikomas geresne laiko skiriamąja geba.

Naudodami šį kodo fragmentą galite suskaidyti vieną iš seisminių pėdsakų į atskirus komponentus:

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]

Pėdsakų išskaidymas į komponentus
Wolfram Mathematica geofizikoje

Norint įvertinti, kaip atspindžio energija pasiskirsto skirtingais bangų atvykimo laikais, naudojamos skalogramos (analogiškos spektrogramai). Paprastai praktiškai nereikia analizuoti visų komponentų. Paprastai pasirenkami žemo, vidutinio ir aukšto dažnio komponentai.

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]

Skalograma. Funkcijos rezultatas Bangos skalograma[]
Wolfram Mathematica geofizikoje

Wolfram kalba naudoja bangelių transformacijos funkciją ContinuousWaveletTransform[]. Ir šios funkcijos taikymas visam pėdsakų rinkiniui bus atliktas naudojant funkciją Lentelė[]. Čia verta paminėti vieną iš „Wolfram Mathematica“ stipriųjų pusių – galimybę naudoti lygiagretavimą ParallelTable[]. Aukščiau pateiktame pavyzdyje lygiagretinimo nereikia – duomenų kiekis nėra didelis, tačiau dirbant su eksperimentiniais duomenų rinkiniais, kuriuose yra šimtai tūkstančių pėdsakų, tai yra būtinybė.

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

Pritaikius funkciją ContinuousWaveletTransform[] Pasirodo nauji duomenų rinkiniai, atitinkantys pasirinktus dažnius. Aukščiau pateiktame pavyzdyje šie dažniai yra: 38Hz, 33Hz, 27Hz. Dažnių pasirinkimas dažniausiai atliekamas testavimo pagrindu – jie gauna efektyvius skirtingų dažnių derinių žemėlapius ir parenka geologo požiūriu informatyviausią.

Jei jums reikia pasidalyti rezultatais su kolegomis arba pateikti juos klientui, galite naudoti „GeologyIO“ paketo funkciją SEGYExport[]

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

Naudojant tris iš šių kubelių (žemo dažnio, vidutinio dažnio ir aukšto dažnio komponentus), RGB maišymas paprastai naudojamas duomenims kartu vizualizuoti. Kiekvienam komponentui priskiriama savo spalva – raudona, žalia, mėlyna. Wolfram Mathematica tai galima padaryti naudojant funkciją ColorCombine[].

Rezultatas yra vaizdai, iš kurių galima padaryti geologinę interpretaciją. Atkarpoje užfiksuoti vingiai leidžia išskirti paleokanalus, kurie greičiausiai yra rezervuarai ir kuriuose yra naftos atsargų. Šiuolaikinių tokios upių sistemos analogų paieška ir analizė leidžia nustatyti perspektyviausias vingių dalis. Patys kanalai pasižymi storais gerai išrūšiuoto smiltainio sluoksniais ir yra geras naftos rezervuaras. Sritys, esančios už "nėrinių" anomalijų, yra panašios į šiuolaikinius užliejamus telkinius. Salpos telkinius daugiausia sudaro molingos uolienos, todėl gręžimas į šias zonas bus neveiksmingas.

RGB duomenų kubo dalis. Centre (šiek tiek į kairę nuo centro) galite atsekti vingiuojančią upę.
Wolfram Mathematica geofizikoje
RGB duomenų kubo dalis. Kairėje pusėje galite atsekti vingiuojančią upę.
Wolfram Mathematica geofizikoje

Kai kuriais atvejais seisminių duomenų kokybė leidžia gauti žymiai aiškesnius vaizdus. Tai priklauso nuo lauko darbų metodikos, triukšmo mažinimo algoritmo naudojamos įrangos. Tokiais atvejais matomi ne tik upių sistemų fragmentai, bet ir ištisos pratęstos paleoupės.

Trijų seisminių duomenų kubo (horizontalaus pjūvio) komponentų RGB maišymas. Gylis apie 2 km.
Wolfram Mathematica geofizikoje
Palydovinis Volgos upės vaizdas netoli Saratovo
Wolfram Mathematica geofizikoje

išvada

„Wolfram Mathematica“ leidžia analizuoti seisminius duomenis ir spręsti taikomąsias problemas, susijusias su mineralų žvalgymu, o „GeologyIO“ paketas šį procesą daro patogesnį. Seisminių duomenų struktūra yra tokia, kad naudojant integruotus metodus skaičiavimams paspartinti (ParallelTable[], ParallelDo[],…) yra labai efektyvus ir leidžia apdoroti didelius duomenų kiekius. Tai didžiąja dalimi palengvina „GeologyIO“ paketo duomenų saugojimo funkcijos. Beje, paketą galima naudoti ne tik taikomosios seisminės žvalgybos srityje. Beveik tos pačios rūšies duomenys naudojami žemės skverbimosi radare ir seismologijoje. Jei turite pasiūlymų, kaip pagerinti rezultatą, kokie signalų analizės algoritmai iš Wolfram Mathematica arsenalo tinka tokiems duomenims, arba jei turite kokių nors kritinių pastabų, prašome Palikite komentarą.

Šaltinis: www.habr.com

Добавить комментарий