ვოლფრამ მათემატიკა გეოფიზიკაში

მადლობა ბლოგის ავტორს ანტონ ეკიმენკო მისი მოხსენებისთვის

შესავალი

ეს შენიშვნა დაიწერა კონფერენციის ფონზე ვოლფრამის რუსული ტექნოლოგიების კონფერენცია და შეიცავს ჩემს მიერ მიცემული მოხსენების შეჯამებას. ღონისძიება ივნისში სანქტ-პეტერბურგში გაიმართა. იმის გათვალისწინებით, რომ ვმუშაობ კონფერენციის საიტიდან ერთი ბლოკით, არ შემეძლო არ დავესწრო ამ ღონისძიებას. 2016 და 2017 წლებში მოვისმინე კონფერენციის მოხსენებები, წელს კი გავაკეთე პრეზენტაცია. ჯერ ერთი საინტერესო (მეჩვენება) თემა გაჩნდა, რომელსაც ჩვენ ვამუშავებთ კირილ ბელოვიდა მეორეც, რუსეთის ფედერაციის კანონმდებლობის ხანგრძლივი შესწავლის შემდეგ სანქციების პოლიტიკასთან დაკავშირებით, იმ საწარმოში, სადაც მე ვმუშაობ, გამოჩნდა ორი ლიცენზია. ვოლფრამ მათემატიკა.

სანამ ჩემი გამოსვლის თემაზე გადავიდოდე, მინდა აღვნიშნო ღონისძიების კარგი ორგანიზება. კონფერენციის სავიზიტო გვერდზე გამოყენებულია ყაზანის ტაძრის გამოსახულება. საკათედრო ტაძარი პეტერბურგის ერთ-ერთი მთავარი ღირსშესანიშნაობაა და ძალიან კარგად ჩანს დარბაზიდან, რომელშიც კონფერენცია გაიმართა.

პეტერბურგის სახელმწიფო ეკონომიკური უნივერსიტეტის შესასვლელთან მონაწილეებს სტუდენტების ასისტენტები დახვდნენ - დაკარგვის უფლებას არ აძლევდნენ. რეგისტრაციისას დაურიგდათ მცირე ზომის სუვენირები (სათამაშო - მოციმციმე წვეტი, კალამი, სტიკერები ვოლფრამის სიმბოლოებით). სადილისა და ყავის შესვენებები კონფერენციის განრიგშიც იყო შეტანილი. ჯგუფის კედელზე უგემრიელესი ყავისა და ღვეზელების შესახებ უკვე აღვნიშნე - მზარეულები მშვენიერი არიან. ამ შესავალი ნაწილით მინდა ხაზგასმით აღვნიშნო, რომ თავად ღონისძიება, მისი ფორმატი და ლოკაცია უკვე დადებით ემოციებს მოაქვს.

მოხსენებას, რომელიც მოვამზადეთ ჩემი და კირილ ბელოვის მიერ, ჰქვია „ვოლფრამ მათემატიკის გამოყენება გამოყენებითი გეოფიზიკის ამოცანების გადასაჭრელად. სეისმური მონაცემების სპექტრული ანალიზი ან „სად მიედინებოდა უძველესი მდინარეები“. მოხსენების შინაარსი მოიცავს ორ ნაწილს: პირველი, ალგორითმების გამოყენებას, რომლებიც ხელმისაწვდომია ვოლფრამ მათემატიკა გეოფიზიკური მონაცემების გასაანალიზებლად და მეორეც, ეს არის ის, თუ როგორ უნდა ჩავდოთ გეოფიზიკური მონაცემები 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 გბ-ზე მეტი, ხოლო ორიგინალური "ნედლეულის" მოცულობა შეიძლება იყოს ასობით გიგაბაიტი.

როგორ ვიმუშაოთ მათთან ვოლფრამ მათემატიკა?

პაკეტი გეოლოგიაIO

პაკეტის შემუშავება დაიწყო კითხვა რუსულენოვანი დამხმარე ჯგუფის VK კედელზე. საზოგადოების პასუხების წყალობით, გამოსავალი ძალიან სწრაფად მოიძებნა. და შედეგად, ის გადაიზარდა სერიოზულ განვითარებაში. Შესაბამისი Wolfram Community კედლის პოსტი მოდერატორებმაც კი აღნიშნეს. ამჟამად, პაკეტი მხარს უჭერს მუშაობას მონაცემთა შემდეგ ტიპებთან, რომლებიც აქტიურად გამოიყენება გეოლოგიურ ინდუსტრიაში:

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

”მარმუსის მონაცემთა ნაკრები შეიქმნა ინსტიტუტში ... მინიმალური სიჩქარე 1500 მ/წმ და მაქსიმუმ 5500 მ/წმ)”

თქვენ შეგიძლიათ აჩვენოთ ფაქტობრივი გეოლოგიური მოდელი სეისმურ კვალზე წვდომით "კვალი" ღილაკის გამოყენებით (პაკეტის ერთ-ერთი მახასიათებელია ის, რომ გასაღებები არ არის რეგისტრირებული):

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 მ/წმ, სოლის გარეთ 4000 მ/წმ, ხოლო სიმკვრივე მიჩნეულია მუდმივად 2200 გ/სმ³. ასეთი მოდელისთვის ჩვენ ვიანგარიშებთ ასახვის კოეფიციენტებს და სეისმურ კვალს.

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 ჰც, 55 ჰც და 85 ჰც სიხშირით.

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 ჰც, 55 ჰც, 85 ჰც სიხშირით


სეისმური კვალის გამოთვლით და არეკლილი ტალღის ამპლიტუდების გრაფიკების გამოსახვით, ჩვენ ვხედავთ, რომ სხვადასხვა სიხშირეზე ანომალია შეინიშნება ფორმირების სხვადასხვა სისქეზე.

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

C 1 ეს არის დემო ფაილი GEOLOGYIO პაკეტის ტესტისთვის
C ჰიტები
C ჰიტები
C ჰიტები
C 5 თარიღი მომხმარებლის სახელი: WOLFRAM USER
C 6 კვლევის დასახელება: სადღაც ციმბირში
C 7 ფაილის ტიპი 3D სეისმური მოცულობა
C ჰიტები
C ჰიტები
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-ში ეს შეიძლება გაკეთდეს გამოყენებით გამჭვირვალობა[] и რასტერი 3D[].

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


ვოლფრამის ენა იყენებს ფუნქციას ვაილეტების ტრანსფორმაციისთვის უწყვეტი ტალღის ტრანსფორმაცია[]. და ამ ფუნქციის გამოყენება კვალის მთელ კომპლექტზე განხორციელდება ფუნქციის გამოყენებით მაგიდა[]. აქვე აღსანიშნავია Wolfram Mathematica-ს ერთ-ერთი ძლიერი მხარე - პარალელიზაციის გამოყენების უნარი. პარალელური ცხრილი[]. ზემოხსენებულ მაგალითში არ არის საჭირო პარალელიზება - მონაცემთა მოცულობა არ არის დიდი, მაგრამ ასიათასობით კვალის შემცველ ექსპერიმენტულ მონაცემთა ნაკრებებთან მუშაობისას ეს აუცილებლობაა.

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

ფუნქციის გამოყენების შემდეგ უწყვეტი ტალღის ტრანსფორმაცია[] არჩეული სიხშირეების შესაბამისი ახალი მონაცემთა ნაკრები გამოჩნდება. ზემოთ მოცემულ მაგალითში ეს სიხშირეებია: 38Hz, 33Hz, 27Hz. სიხშირეების არჩევა ყველაზე ხშირად ტესტირების საფუძველზე ხდება - ისინი იღებენ ეფექტურ რუკებს სიხშირის სხვადასხვა კომბინაციისთვის და ირჩევენ ყველაზე ინფორმაციულს გეოლოგის თვალსაზრისით.

თუ გჭირდებათ კოლეგებთან შედეგების გაზიარება ან მომხმარებლისთვის მიწოდება, შეგიძლიათ გამოიყენოთ GeologyIO პაკეტის 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];

ამ სამი კუბით (დაბალი სიხშირის, საშუალო სიხშირის და მაღალი სიხშირის კომპონენტებით), RGB შერევა, როგორც წესი, გამოიყენება მონაცემთა ერთად ვიზუალიზაციისთვის. თითოეულ კომპონენტს ენიჭება საკუთარი ფერი - წითელი, მწვანე, ლურჯი. Wolfram Mathematica-ში ეს შეიძლება გაკეთდეს ფუნქციის გამოყენებით ფერების კომბინაცია[].

შედეგი არის სურათები, საიდანაც შესაძლებელია გეოლოგიური ინტერპრეტაციის გაკეთება. მონაკვეთზე დაფიქსირებული მეანდრები შესაძლებელს ხდის პალეოარხების გამოკვეთას, რომლებიც, სავარაუდოდ, რეზერვუარებია და ნავთობის მარაგს შეიცავს. ასეთი მდინარის სისტემის თანამედროვე ანალოგების ძიება და ანალიზი საშუალებას გვაძლევს განვსაზღვროთ მეანდრების ყველაზე პერსპექტიული ნაწილები. თავად არხები ხასიათდება კარგად დალაგებული ქვიშაქვის სქელი ფენებით და ნავთობის კარგი რეზერვუარია. „მაქმანის“ ანომალიების მიღმა არსებული ტერიტორიები თანამედროვე ჭალის საბადოების მსგავსია. ჭალის საბადოები ძირითადად წარმოდგენილია თიხიანი ქანებით და ამ ზონებში ბურღვა არაეფექტური იქნება.

მონაცემთა კუბის RGB ნაჭერი. ცენტრში (ცენტრიდან ოდნავ მარცხნივ) შეგიძლიათ თვალყური ადევნოთ მეანდერ მდინარეს.

მონაცემთა კუბის RGB ნაჭერი. მარცხენა მხარეს შეგიძლიათ თვალყური ადევნოთ მეანდერ მდინარეს.


ზოგიერთ შემთხვევაში, სეისმური მონაცემების ხარისხი იძლევა მნიშვნელოვნად უფრო მკაფიო გამოსახულებებს. ეს დამოკიდებულია საველე სამუშაოების მეთოდოლოგიაზე, ხმაურის შემცირების ალგორითმის მიერ გამოყენებულ აღჭურვილობაზე. ასეთ შემთხვევებში ჩანს არა მხოლოდ მდინარის სისტემების ფრაგმენტები, არამედ მთელი გაშლილი პალეო-მდინარეები.

სეისმური მონაცემების კუბის სამი კომპონენტის RGB შერევა (ჰორიზონტალური ნაჭერი). სიღრმე დაახლოებით 2 კმ.

მდინარე ვოლგის სატელიტური სურათი სარატოვთან

დასკვნა

Wolfram Mathematica-ს შეუძლია სეისმური მონაცემების ანალიზი და მინერალების ძიებასთან დაკავშირებული გამოყენებითი პრობლემების გადაჭრა, ხოლო GeologyIO პაკეტი ამ პროცესს უფრო მოსახერხებელს ხდის. სეისმური მონაცემების სტრუქტურა ისეთია, რომ გამოთვლების დასაჩქარებლად ჩაშენებული მეთოდების გამოყენებით (პარალელური ცხრილი[], ParallelDo[],…) არის ძალიან ეფექტური და გაძლევთ საშუალებას დაამუშაოთ დიდი რაოდენობით მონაცემები. დიდწილად, ამას ხელს უწყობს GeologyIO პაკეტის მონაცემთა შენახვის მახასიათებლები. სხვათა შორის, პაკეტის გამოყენება შესაძლებელია არა მხოლოდ გამოყენებითი სეისმური ძიების სფეროში. თითქმის იგივე ტიპის მონაცემები გამოიყენება სახმელეთო რადარებსა და სეისმოლოგიაში. თუ გაქვთ წინადადებები, როგორ გააუმჯობესოთ შედეგი, Wolfram Mathematica-ს არსენალიდან რომელი სიგნალის ანალიზის ალგორითმები გამოიყენება ასეთ მონაცემებზე, ან თუ გაქვთ რაიმე კრიტიკული კომენტარი, გთხოვთ. დატოვე კომენტარი.

წყარო: www.habr.com