Danke dem Blogautor fĂŒr seinen Vortrag
EinfĂŒhrung
Dieser Beitrag wurde im Anschluss an die Konferenz verfasst und enthÀlt eine Zusammenfassung des Vortrags, den ich gehalten habe. Die Veranstaltung fand im Juni in St. Petersburg statt. Da ich nur einen Block vom Konferenzort arbeite, konnte ich es mir nicht leisten, dieses Ereignis zu verpassen. In den Jahren 2016 und 2017 habe ich den VortrÀgen der Konferenz gelauscht, und in diesem Jahr habe ich selbst einen Vortrag gehalten. Zuerst gab es ein interessantes (wie ich finde) Thema, das wir mit entwickeln, und zweitens, nach einer langen Analyse der Gesetzgebung der Russischen Föderation im Bereich der Sanktionspolitik, hat das Unternehmen, in dem ich arbeite, tatsÀchlich zwei Lizenzen erhalten .
Bevor ich zum Thema meines Vortrags komme, möchte ich die gute Organisation der Veranstaltung hervorheben. Auf der Visitenkarte der Konferenz wird ein Bild der Kasaner Kathedrale verwendet. Die Kathedrale ist eine der Hauptattraktionen Petersburgs und ist vom Saal, in dem die Konferenz stattfand, sehr gut sichtbar.

Am Eingang der SPbGEU wurden die Teilnehmer von Assistenten aus den Reihen der Studenten empfangen â sie sorgten dafĂŒr, dass niemand sich verlief. Bei der Registrierung wurden kleine Souvenirs (Spielzeug â blinkender Spike, Stift, Aufkleber mit Wolfram-Symbolik) verteilt. Mittagessen und Kaffeepause waren ebenfalls im Konferenzzeitplan enthalten. Zum leckeren Kaffee und den Snacks, die ich bereits in der Gruppenwand erwĂ€hnt habe â die Köche waren groĂartig. Mit diesem einleitenden Teil möchte ich betonen, dass die Veranstaltung selbst, ihr Format und der Veranstaltungsort bereits positive Emotionen hervorrufen.
Der Vortrag, den ich zusammen mit Kirill Belov vorbereitet habe, trĂ€gt den Titel âVerwendung von Wolfram Mathematica zur Lösung von Aufgaben der angewandten Geophysik. Spektralanalyse seismischer Daten oder âwo flossen die alten FlĂŒsseâ. Der Inhalt des Vortrags umfasst zwei Teile: erstens die Verwendung von Algorithmen, die in zur Analyse geophysikalischer Daten zur VerfĂŒgung stehen, und zweitens, wie geophysikalische Daten in Wolfram Mathematica integriert werden können.
Seismik
ZunĂ€chst ist es notwendig, einen kleinen Exkurs in die Geophysik zu machen. Geophysik ist die Wissenschaft, die sich mit den physikalischen Eigenschaften von Gesteinen befasst. Da Gesteine unterschiedliche Eigenschaften besitzen â elektrische, magnetische, elastische â existieren entsprechende Methoden der Geophysik: Elektroexploration, Magnetexploration, Seismik... In diesem Artikel werden wir uns jedoch nĂ€her mit der seismischen Exploration befassen. Seismik ist die Hauptmethode zur Suche nach Ăl und Gas. Das Verfahren basiert auf der Erzeugung elastischer Schwingungen und der anschlieĂenden Registrierung der RĂŒckmeldung von den Gesteinen, die das untersuchte Gebiet bilden. Die Erzeugung der Schwingungen erfolgt an Land (mit Dynamit oder nicht detonierenden vibrierenden Quellen elastischer Schwingungen) oder auf dem Meer (mit Pneumpattern). Elastische Schwingungen breiten sich durch die Schichten der Gesteine aus, refragiert und reflektiert an den Grenzen von Schichten mit unterschiedlichen Eigenschaften. Die reflektierten Wellen kehren zur OberflĂ€che zurĂŒck und werden mit Geophonen an Land (ĂŒblicherweise elektrodynamische GerĂ€te, die auf der Bewegung eines Magneten beruhen, der in einer Spule hĂ€ngt) oder Hydrofonen auf dem Meer (die auf dem Fotoeffekt basieren) registriert. Anhand der Ankunftszeiten der Wellen kann auf die Tiefen der geologischen Schichten geschlossen werden.
Das seismische Schiff schleppt AusrĂŒstung

Die Pneumpattern erzeugt elastische Schwingungen

Die Wellen durchdringen die Gesteinsschichten und werden von Hydrofonen registriert

Das Forschungsschiff fĂŒr geophysikalische Erkundung âIwan Gubkinâ am Kai unter der Blagoweschtschenski-BrĂŒcke in St. Petersburg

Modell des seismischen Signals
Gesteine weisen verschiedene physikalische Eigenschaften auf. FĂŒr die seismische Exploration sind vor allem die elastischen Eigenschaften wichtig â die Geschwindigkeit der elastischen Schwingungen und die Dichte. Wenn zwei Schichten identische oder Ă€hnliche Eigenschaften haben, wird die Welle die Grenze zwischen ihnen ânicht bemerkenâ. Wenn sich jedoch die Wellen Geschwindigkeiten in den Schichten unterscheiden, kommt es an der GrenzflĂ€che der Schichten zu einer Reflexion. Je gröĂer der Unterschied in den Eigenschaften, desto intensiver wird die Reflexion sein. Ihre IntensitĂ€t wird durch den Reflexionskoeffizienten (rc) bestimmt:

wobei Ï die Dichte der Gesteine und Μ die Wellen geschwindigkeit ist, 1 und 2 die obere und untere Schicht bezeichnen.
Eines der einfachsten und am hÀufigsten verwendeten Modelle seismischer Signale ist das Faltungsmodell, bei dem die aufgezeichnete seismische Spur als Ergebnis der Faltung einer Folge von Reflexionskoeffizienten mit einem anregenden Impuls dargestellt wird:
![]()
wo s(t) â seismische Spur, d.h. alles, was das Hydrofon oder Geophon wĂ€hrend der festen Registrierungszeit aufgezeichnet hat, w(t) â das Signal, das die Pneumatik erzeugt, n(t) â zufĂ€lliger LĂ€rm.
Berechnen wir zum Beispiel eine synthetische Seismospur. Als Ausgangssignal verwenden wir den in der seismischen Exploration weit verbreiteten Ricker-Impuls.
length=0.050; (*Signal-LĂ€nge*)
dt=0.001;(*Abtastrate des Signals*)
t=Range[-length/2,(length)/2,dt];(*Signal-Zeit*)
f=35;(*Zentralfrequenz*)
welle=(1.0-2.0*(Pi^2)*(f^2)*(t^2))*Exp[-(Pi^2)*(f^2)*(t^2)];
ListLinePlot[welle, 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]UrsprĂŒnglicher seismischer Impuls

Wir setzen zwei Grenzen in Tiefen von 300 ms und 600 ms, und die Reflexionskoeffizienten werden zufÀllige Zahlen sein.
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["Reflexionskoeffizienten",Black,20],
LabelStyle->Directive[Black,Italic]]Folge der Reflexionskoeffizienten

Wir berechnen und visualisieren die seismische Spur. Da die Reflexionskoeffizienten unterschiedliche Vorzeichen haben, erhalten wir auch in der seismischen Spur zwei wechselnde Reflexionen.
traceExamle=ListConvolve[welle[[1;;;;1]],rcExample];
ListPlot[traceExamle,
PlotStyle->Black,Filling->0,Frame->True,Axes->False,
PlotLabel->Style["Seismische Spur",Black,20],
LabelStyle->Directive[Black,Italic]]Simulierte Spur

FĂŒr dieses Beispiel muss ich darauf hinweisen, dass in der RealitĂ€t die Tiefe der Schichten natĂŒrlich in Metern bestimmt wird, wĂ€hrend die Berechnung der seismischen Spur im Zeitbereich erfolgt. Es wĂ€re richtiger, die Tiefen in Metern anzugeben und die Ankunftszeiten zu berechnen, wenn man die Geschwindigkeiten in den Schichten kennt. In diesem Fall habe ich die Schichten sofort auf der Zeitachse festgelegt.
Wenn es um Feldforschung geht, werden durch solche Beobachtungen eine enorme Menge Ă€hnlicher Zeitreihen (seismische Trassen) registriert. Zum Beispiel, bei der Untersuchung eines Gebiets von 25 km LĂ€nge und 15 km Breite, wo jede Trasse aufgrund der Arbeiten eine Zelle von 25x25 Metern charakterisiert (eine solche Zelle wird als Bin bezeichnet), wird die endgĂŒltige Datenmenge 600.000 Trassen enthalten. Bei einer Zeitdiskretisierung von 1 ms und einer Aufzeichnungszeit von 5 Sekunden wird die endgĂŒltige Datendatei ĂŒber 11 GB groĂ sein, wĂ€hrend das Volumen des ursprĂŒnglichen ârohenâ Materials Hunderte von Gigabyte betragen kann.
Wie man damit in ?
Paket
Die Entwicklung des Pakets begann in der VK-Gruppe der russischsprachigen UnterstĂŒtzung. Dank der Antworten der Community wurde die Lösung sehr schnell gefunden. Dies hat sich zu einer ernsthaften Entwicklung entwickelt. Der entsprechende wurde sogar von Moderatoren hervorgehoben. Momentan unterstĂŒtzt das Paket die Arbeit mit folgenden Datentypen, die in der Geologiebranche aktiv verwendet werden:
- Import von kartografischen Daten im ZMAP- und IRAP-Format
- Import von Messungen in Bohrlöchern im LAS-Format
- Ein- und Ausgabe von seismischen Dateien im
Um das Paket zu installieren, befolgen Sie die Anweisungen auf der Seite zum Herunterladen des gebauten Pakets, d.h. fĂŒhren Sie den folgenden Code in einem :
If[PacletInformation["GeologyIO"] === {}, PacletInstall[URLDownload[
"https://wolfr.am/FiQ5oFih",
FileNameJoin[{CreateDirectory[], "GeologyIO-0.2.2.paclet"}]
]]]Nach dem erfolgt die Installation des Pakets in den Standardordner, dessen Pfad folgendermaĂen abgerufen werden kann:
FileNameJoin[{$UserBasePacletsDirectory, "Repository"}]Um die Hauptfunktionen des Pakets zu demonstrieren, wird der Aufruf traditionell fĂŒr Pakete in der Wolfram-Sprache durchgefĂŒhrt:
Get["GeologyIO`"]Das Paket wird mit dem Einsatz von entwickelt. Dies ermöglicht es, die Hauptfunktionen des Pakets mit Dokumentationen zu begleiten, die in der PrĂ€sentationsform nicht von der Dokumentation selbst von Wolfram Mathematica abweichen und dem Paket Testdateien fĂŒr eine erste Bekanntschaft zur VerfĂŒgung zu stellen.


Eine solche Datei ist unter anderem die Datei âMarmousi.segyâ â ein synthetisches Modell eines geologischen Querschnitts, das vom französischen Institut fĂŒr Erdöl entwickelt wurde. Die Entwickler testen mithilfe dieses Modells ihre eigenen Algorithmen zur Modellierung des Wellenfelds, zur Datenverarbeitung, zur Inversion seismischer Spuren usw. Das Modell Marmousi selbst wird in einem Repository gespeichert, von wo aus das Paket heruntergeladen wurde. Um die Datei zu erhalten, fĂŒhren wir den folgenden Code aus:
If[Not[FileExistsQ["Marmousi.segy"]],
URLDownload["https://wolfr.am/FiQGh7rk", "Marmousi.segy"];]
marmousi = SEGYImport["Marmousi.segy"]Das Ergebnis des Imports ist ein SEGYData-Objekt

Das SEGY-Format sieht die Speicherung unterschiedlicher Informationen ĂŒber die Beobachtungen vor. ZunĂ€chst handelt es sich um Textkommentare. Hier werden Informationen ĂŒber den Arbeitsort, die Namen der Firmen, die die Messungen durchgefĂŒhrt haben, usw. eingetragen. In unserem Fall wird dieser Header durch eine Abfrage mit dem SchlĂŒssel TextHeader aufgerufen. Hier ist der verkĂŒrzte Textheader dargestellt:
Short[marmousi["TextHeader"]]âDas Marmousi-Dataset wurde am Institut ⊠Mindestgeschwindigkeit von 1500 m/s und eine Höchstgeschwindigkeit von 5500 m/s erzeugtâ
Die geologische Modell kann angezeigt werden, indem man auf die seismischen Spuren mit dem SchlĂŒssel âtracesâ zugreift (eine der Funktionen des Pakets ist die GroĂ-/Kleinschreibung der SchlĂŒsselunabhĂ€ngigkeit):
ArrayPlot[Transpose[marmousi["traces"]], PlotTheme -> "Detailed"]Modell Marmousi

Derzeit ermöglicht das Paket auch das schrittweise Laden von Daten aus groĂen Dateien, wodurch die Verarbeitung von Dateien mit einer GröĂe von bis zu mehreren zehn Gigabyte möglich wird. Zu den Funktionen des Pakets gehören auch Funktionen zum Export von Daten in .segy und zum teilweisen AnhĂ€ngen an das Ende der Datei.
Besonders hervorzuheben ist die FunktionalitĂ€t des Pakets im Umgang mit komplexen .segy-Dateistrukturen. Es ermöglicht nicht nur den Zugriff auf einzelne Spuren, Header und deren Ănderungen ĂŒber SchlĂŒssel und Indizes, sondern auch deren folgende Speicherung in der Datei. Viele technische Details der GeologyIO-Implementierung gehen ĂŒber den Rahmen dieses Artikels hinaus und verdienen wahrscheinlich eine gesonderte Beschreibung.
Die Relevanz der Spektralanalyse in der seismischen Exploration
Die Möglichkeit, seismische Materialien in Wolfram Mathematica zu importieren, ermöglicht die Nutzung der integrierten Signalverarbeitungsfunktionen fĂŒr experimentelle Daten. Da jede seismische Spur eine Zeitreihe darstellt, ist die spektrale Analyse eines der wichtigsten Werkzeuge zu ihrer Untersuchung. Zu den Faktoren, die die Analyse des Frequenzinhalts seismischer Daten beeinflussen können, gehören beispielsweise folgende:
- Verschiedene Wellentypen weisen unterschiedliche Frequenzzusammensetzungen auf. Dies ermöglicht es, nĂŒtzliche Wellen zu isolieren und Störwellen zu dĂ€mpfen.
- Eigenschaften von Gesteinen wie PorositÀt und SÀttigung können den Frequenzinhalt beeinflussen. Dies ermöglicht die Unterscheidung von Gesteinen mit besseren Eigenschaften.
- Schichten unterschiedlicher Dicke verursachen Anomalien in verschiedenen Frequenzbereichen.
Der dritte Punkt ist im Kontext dieses Artikels der entscheidende. Im Folgenden ist ein Codeausschnitt zur Berechnung seismischer Spuren im Fall einer Schicht mit variierender Dicke â das Keilmodell â aufgefĂŒhrt. Dieses Modell wird traditionell in der Seismik zur Analyse von Interferenzeffekten untersucht, wenn Wellen, die von vielen Schichten reflektiert werden, ĂŒbereinanderliegen.
nx=200;(* Anzahl der Gitterpunkte in X-Richtung *)
ny=200;(* Anzahl der Gitterpunkte in Y-Richtung *)
T=2;(* Gesamte Ausbreitungszeit *)
(* Geschwindigkeit und Dichte *)
modellv=Table[4000,{i,1,ny},{j,1,nx}];(* P-Wellen-Geschwindigkeit in m/s *)
rho=Table[2200,{i,1,ny},{j,1,nx}];(* Dichte in g/cmÂł, konstante Dichte verwendet *)
Table[modellv[[150-Round[i*0.5];;,i]]=4500;,{i,1,200}];
Table[modellv[[;;70,i]]=4500;,{i,1,200}];
(* Modellzeichnung *)
MatrixPlot[modellv,PlotLabel->Style["Modell der Schicht",Schwarz,20],
LabelStyle->Directive[Schwarz,Italic]]Modell eines keilförmigen Schicht

Die Geschwindigkeit der Wellen innerhalb des Keils betrĂ€gt 4500 m/s, auĂerhalb des Keils 4000 m/s, und die Dichte wird konstant mit 2200 g/cmÂł angenommen. FĂŒr dieses Modell berechnen wir die Reflexionskoeffizienten und seismischen Spuren.
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],Schwarz}}},PlotStyle->Directive[Grau,DĂŒnn],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->Rot],PlotRange->All,Frame->True],270Degree]Seismische Spuren fĂŒr das Keilmodell

Die Sequenz der seismischen Linien, die auf diesem Bild dargestellt ist, wird als seismischer Schnitt bezeichnet. Wie zu erkennen ist, kann ihre Interpretation auch auf intuitiver Ebene erfolgen, da die Geometrie der reflektierten Wellen eindeutig dem zuvor festgelegten Modell entspricht. Bei einer detaillierteren Analyse der Linien kann festgestellt werden, dass die Linien von 1 bis etwa 30 sich nicht unterscheiden â die Reflexion von der Oberkante und der Unterkante der Schicht ĂŒberlagern sich nicht. Ab der 31. Linie beginnen die Reflexionen zu interferieren. Und obwohl in dem Modell die Reflexionskoeffizienten horizontal unverĂ€ndert bleiben â Ă€ndern sich die seismischen Linien in ihrer IntensitĂ€t mit der Ănderung der Schichtdicke.
Betrachten wir die Amplitude der Reflexion von der oberen Grenze der Schicht. Ab der 60. Linie beginnt die IntensitĂ€t der Reflexion zu steigen und erreicht bei der 70. Linie ihren Höchstwert. So zeigt sich die Interferenz der Wellen von der Ober- und Unterkante fĂŒr die Schichten, was in einigen FĂ€llen zu erheblichen Anomalien im seismischen Datensatz fĂŒhrt.
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]Diagramm der Amplitude der reflektierten Welle von der oberen Kante des Keils

Es ist logisch, dass bei einem tieferfrequenten Signal die Interferenz bei gröĂeren Schichtdicken auftritt, wĂ€hrend bei einem hochfrequenten Signal die Interferenz bei geringeren Dicken auftritt. Der folgende Codeabschnitt erzeugt ein Signal mit Frequenzen von 35 Hz, 55 Hz und 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]Set von Eingangssignalen mit Frequenzen von 35 Hz, 55 Hz, 85 Hz

Nachdem wir die seismischen Linien berechnet und die Diagramme der Amplituden der reflektierten Welle erstellt haben, können wir sehen, dass die Anomalie bei unterschiedlichen Frequenzen bei unterschiedlichen Schichtdicken auftritt.
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["Amplituden der Reflexion",Black,20],
LabelStyle->Directive[Black,Italic],
Frame->True]Diagramme der Amplituden der reflektierten Welle von der oberen Kante der SchrĂ€ge fĂŒr verschiedene Frequenzen

Die Möglichkeit, aus seismischen Beobachtungen RĂŒckschlĂŒsse auf die Schichtdicke zu ziehen, ist Ă€uĂerst nĂŒtzlich, denn eine der Hauptaufgaben bei der Erkundung von Ălfeldern ist die Bewertung der vielversprechendsten Stellen fĂŒr den Bohrlochbau (d.h. die Abschnitte, in denen die Schicht eine groĂe Dicke hat). DarĂŒber hinaus können im geologischen Profil Objekte vorkommen, die durch ihre Genese einen abrupten Wechsel der Schichtdicken verursachen. Dies macht die spektrale Analyse zu einem effektiven Werkzeug zu ihrer Untersuchung. Im nĂ€chsten Teil des Artikels werden wir solche geologischen Objekte genauer betrachten.
Experimentelle Daten. Wo wurden sie gewonnen und was ist darin zu suchen?
Die in dem Artikel analysierten Materialien stammen aus dem Gebiet der Westsibirischen Region. Diese Region ist, wie wohl jeder weiĂ, das Hauptölproduktionsgebiet unseres Landes. Die aktive ErschlieĂung von LagerstĂ€tten begann in den 60er Jahren des letzten Jahrhunderts. Die Hauptmethode zur Suche nach Ălvorkommen ist die seismische Exploration. Es ist interessant, Satellitenbilder dieses Gebiets zu betrachten. Im kleinen MaĂstab kann man eine riesige Anzahl von SĂŒmpfen und Seen bemerken, vergröĂert man die Karte, so sind Bohrstellen in Form von Kuppeln zu sehen, und wenn man die Karte bis ins kleinste Detail vergröĂert, kann man auch die DurchgĂ€nge der Profile erkennen, entlang derer seismische Beobachtungen durchgefĂŒhrt wurden.
Satellitenbild von Yandex Maps â Gebiet der Stadt Noyabrsk

Netzwerk von Kuppelbohrstellen in einem der LagerstÀtten

ĂlfĂŒhrende Gesteine in der Westsibirischen Region liegen in einem breiten Tiefenbereich â von 1 km bis 5 km. Der Hauptteil der ölhaltigen Gesteine wurde im Jura- und Kreidezeitalter gebildet. Viele wissen wahrscheinlich, dass die Jurazeit durch den gleichnamigen Film bekannt ist. war wesentlich anders als die moderne. In der Encyclopedia Britannica gibt es eine Reihe von PalĂ€okarten, die jede geologische Epoche charakterisieren.
GegenwÀrtige Zeit

Jurazeit

Beachten Sie, dass das Gebiet Westsibirien zur Jurazeit eine KĂŒstenlinie darstellte (LandflĂ€chen durch FlĂŒsse und ein flaches Meer durchzogen). Da das Klima angenehm war, kann man annehmen, dass die typische Landschaft dieser Zeit folgendermaĂen aussah:
Sibirien in der Jurazeit

Auf diesem Bild sind fĂŒr uns nicht so sehr die Tiere und Vögel wichtig, sondern die Darstellung des Flusses im Hintergrund. Der Fluss ist das geologische Objekt, auf dem wir zuvor bereits verweilten. Der Grund dafĂŒr ist, dass die AktivitĂ€t von FlĂŒssen gut sortierte Sandsteine anhĂ€uft, die spĂ€ter zu einem Reservoir fĂŒr Erdöl werden. Diese Reservoirs können eine kurvenreiche, komplexe Form haben (wie auch das Flussbett) und sie haben eine wechselnde Dicke â an den Ufern ist die Dicke gering, wĂ€hrend sie nĂ€her zur Mitte des Flussbetts oder in MĂ€anderabschnitten zunimmt. Die in der Jurazeit gebildeten FlĂŒsse befinden sich nun auf einer Tiefe von etwa drei Kilometern und sind Ziel von Erdölsuchen.
Experimentelle Daten. Verarbeitung und Visualisierung
Wir möchten sofort anmerken, dass bezĂŒglich der seismischen Materialien, die in diesem Artikel gezeigt werden â aufgrund des erheblichen Datenvolumens, das fĂŒr die Analyse verwendet wurde â nur ein Ausschnitt des Originalsatzes seismischer Strecken in den Text des Artikels aufgenommen wurde. Dies ermöglicht allen Interessierten, die dargestellten Berechnungen nachzuvollziehen.
Wenn Geophysiker mit seismischen Daten arbeiten, nutzen sie in der Regel spezialisierte Software (einige MarktfĂŒhrer, deren Entwicklungen aktiv genutzt werden, sind zum Beispiel Petrel oder Paradigm), die es ermöglicht, verschiedene Datentypen zu analysieren und eine benutzerfreundliche grafische OberflĂ€che bietet. Trotz aller Bequemlichkeit haben solche Softwarelösungen auch ihre Nachteile â so dauert die Implementierung moderner Algorithmen in stabile Versionen oft viel Zeit, und die Möglichkeiten zur Automatisierung von Berechnungen sind in der Regel begrenzt. In einer solchen Situation ist die Verwendung von Systemen der Computer-Mathematik und Hochsprache-Programmiersprachen sehr praktisch, die eine breite algorithmische Basis nutzen und zudem viel Routinearbeit ĂŒbernehmen. Nach diesem Prinzip funktioniert auch die Arbeit mit seismischen Daten in Wolfram Mathematica. Es ist nicht sinnvoll, eine umfangreiche FunktionalitĂ€t fĂŒr die interaktive Arbeit mit Daten zu schreiben â wichtiger ist es, den Import aus einem gebrĂ€uchlichen Format zu gewĂ€hrleisten, die gewĂŒnschten Algorithmen anzuwenden und sie anschlieĂend in ein externes Format zurĂŒckzuladen.
GemÀà dem vorgeschlagenen Schema laden wir die ursprĂŒnglichen seismischen Daten und stellen sie dar in :
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]Die auf diese Weise heruntergeladenen und importierten Daten sind Trassen, die auf einem Gebiet von 10 mal 5 Kilometern registriert wurden. Wenn die Daten nach der Methode der dreidimensionalen Seismik gewonnen wurden (die Registrierung der Wellen erfolgt nicht entlang einzelner geophysikalischer Profile, sondern ĂŒber das gesamte Gebiet gleichzeitig), ist es möglich, WĂŒrfel seismischer Daten zu erzeugen. Dies sind dreidimensionale Objekte, deren vertikale und horizontale Schnitte eine detaillierte Untersuchung des geologischen Umfelds ermöglichen. Im betrachteten Beispiel haben wir es mit genau dreidimensionalen Daten zu tun. Einige Informationen können wir aus dem Textheader entnehmen, zum Beispiel so
StringPartition[seismic3DSEGY["textheader"], 80] // TableFormC 1 DIES IST DIE DEMO-DATEI FĂR DEN TEST DES GEOLOGYIO-PAKETS
C 2
C 3
C 4
C 5 DATUM BENUTZERNAME: WOLFRAM BENUTZER
C 6 BEFRAGUNGSNAME: IRGENDWO IN SIBIRIEN
C 7 DATEITYP 3D SEISMISCHER VOLUMEN
C 8
C 9
C10 Z-BEREICH: ERSTER 2200M LETZTER 2400M
Dieser Datensatz reicht aus, um die grundlegenden Schritte der Datenanalyse zu demonstrieren. Die Traces in der Datei sind sequenziell aufgezeichnet, und jeder von ihnen sieht ungefÀhr so aus wie im folgenden Bild - es handelt sich um die Verteilung der Amplituden der reflektierten Wellen entlang der vertikalen Achse (Tiefenachse).
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]Einer der Traces des seismischen Profils

Wenn man weiĂ, wie viele Traces in jeder Richtung des untersuchten Gebiets vorhanden sind, kann man ein dreidimensionales Datenarray erstellen und es mit der Funktion Image3D[] darstellen.
traces=seismic3DSEGY["traces"];
startIL=1050;EndIL=2000;stepIL=2; (*X-Koordinate des Anfangs und Endes der Aufnahme und Schritt der Traces*)
startXL=1165;EndXL=1615;stepXL=2; (*Y-Koordinate des Anfangs und Endes der Aufnahme und Schritt der Traces*)
numIL=(EndIL-startIL)/stepIL+1; (*Anzahl der Traces entlang der X-Achse*)
numXL=(EndXL-startXL)/stepIL+1; (*Anzahl der Traces entlang der Y-Achse*)
Image3D[ArrayReshape[Abs[traces/Max[Abs[traces[[All,1;;;;4]]]]],{numIL,numXL,101}],ViewPoint-> {-1, 0, 0},Background->RGBColor[0,0,0]]Dreidimensionale Darstellung eines WĂŒrfels seismischer Daten. (Vertikale Achse â Tiefe)

Falls geologische Objekte von Interesse intensive seismische Anomalien erzeugen, können Visualisierungstools mit Transparenz verwendet werden. Unwichtige Teile der Aufzeichnung können unsichtbar gemacht werden, wÀhrend nur die Anomalien sichtbar bleiben. In Wolfram Mathematica kann dies erreicht werden mit und .
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]Darstellung des WĂŒrfels seismischer Daten mit den Funktionen Opacity[] und Raster3D[]

Wie im synthetischen Beispiel können an den Schnitten des ursprĂŒnglichen WĂŒrfels einige geologische Grenzen (Schichten) mit variierendem Relief hervorgehoben werden.
Das Hauptwerkzeug der Spektralanalyse ist die Fourier-Transformation. Damit lĂ€sst sich das Amplitudenfrequenzspektrum jeder Spur oder Gruppe von Spuren bewerten. Allerdings geht beim Ăbergang der Daten in den Frequenzbereich die Information darĂŒber verloren, zu welchen Zeiten (d.h. in welchen Tiefen) sich die Frequenz Ă€ndert. Um die Ănderungen des Signals auf der Zeit- (Tiefen-)Achse lokalisieren zu können, werden die Fenster-Fourier-Transformation und die Wavelet-Zerlegung verwendet. In diesem Artikel wird das Wavelet-Zerlegung verwendet. Die Technologie der Wavelet-Analyse begann in den 90er Jahren aktiv in der seismischen Erkundung Anwendung zu finden. Ein Vorteil gegenĂŒber der Fenster-Fourier-Transformation ist die bessere zeitliche Auflösung.
Mit folgendem Codefragment kann die Zerlegung in einzelne Komponenten einer der seismischen Spuren durchgefĂŒhrt werden:
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-Zerlegung",Black,32],
LabelStyle->Directive[Black,Italic],
PlotRange->All,
Frame->True]Zerlegung der Spur in Komponenten

Um zu bewerten, wie die Energie der Reflexionen zu verschiedenen Ankunftszeiten der Wellen verteilt ist, werden Skalenogramme (Àhnlich wie Spektrogramme) verwendet. In der Praxis ist es in der Regel nicht notwendig, alle Komponenten zu analysieren. Normalerweise werden die niedrigen, mittleren und hohen Frequenzkomponenten ausgewÀhlt.
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->{"Zeit","Frequenz(Hz)"},LabelStyle->Directive[Black,Bold,14],
ColorFunction->"RustTones",ImageSize->Large]Skalogramm. Ergebnis der Funktion

In der Wolfram-Sprache wird fĂŒr die Wavelet-Transformation die Funktion . Die Anwendung dieser Funktion auf den gesamten Satz von Spuren erfolgt mit Hilfe der Funktion . Hier sollte man einen der groĂen Vorteile von Wolfram Mathematica, die Möglichkeit der Parallelisierung, erwĂ€hnen. . Im gegebenen Beispiel ist keine Parallelisierung erforderlich â das Datenvolumen ist nicht groĂ, aber beim Arbeiten mit experimentellen DatensĂ€tzen, die Hunderte von Tausenden von Spuren enthalten, ist dies eine Notwendigkeit.
tracesCWD=Table[Map[Hilbert[#,0]&,Re[ContinuousWaveletTransform[traces[[i]]][[1]]][[{13,15,18}]]],{i,1,Length@traces}]; Nach der Anwendung der Funktion Es erscheinen neue DatensĂ€tze, die den gewĂ€hlten Frequenzen entsprechen. Im obigen Beispiel sind dies die Frequenzen: 38 Hz, 33 Hz, 27 Hz. Die Auswahl der Frequenzen erfolgt hĂ€ufig auf Basis von Tests, wobei effektive Karten fĂŒr verschiedene Frequenzkombinationen erstellt und die informativste aus geologischer Sicht ausgewĂ€hlt wird.
Wenn die Ergebnisse mit Kollegen geteilt oder einem Kunden zur VerfĂŒgung gestellt werden mĂŒssen, kann die Funktion SEGYExport[] des GeologyIO-Pakets verwendet werden.
outputdata=seismic3DSEGY;
outputdata["traces",1;;-1]=tracesCWD[[All,3]];
outputdata["textheader"]="Wavelet-Dekompositions-Ergebnis";
outputdata["binaryheader","NumberDataTraces"]=Length[tracesCWD[[All,3]]];
SEGYExport["D:result.segy",outputdata];Mit drei solchen Kuben (niedrigfrequente, mittel- und hochfrequente Komponenten) wird normalerweise eine RGB-Mischung zur gemeinsamen Visualisierung der Daten verwendet. Jeder Komponente wird eine Farbe zugewiesen â rot, grĂŒn, blau. In Wolfram Mathematica kann dies mit der Funktion gemacht werden. .
Das Ergebnis sind Bilder, anhand derer geologische Interpretationen vorgenommen werden können. Die MĂ€ander, die im Schnitt sichtbar werden, ermöglichen es, alte Flussbetten zu konturieren, die mit höherer Wahrscheinlichkeit Reservoirs enthalten und Ăl fĂŒhren können. Die Suche und Analyse moderner Analogien solcher Flusssysteme ermöglicht es, die vielversprechendsten Teile der MĂ€ander zu bestimmen. Eigentliche Flussbetten sind durch mĂ€chtige Schichten gut sortierten Sandsteins charakterisiert und sind gute Reservoirs fĂŒr Ăl. Bereiche auĂerhalb der "Schnur-Anomalien" sind Ă€hnlich den modernen Auenablagerungen. Auenablagerungen bestehen hauptsĂ€chlich aus tonigen Gesteinen, und Bohrungen in diese Zonen werden ineffektiv sein.
RGB-Schnitt des Datenkubus. In der Mitte (etwas links von der Mitte) kann man den mÀandrierenden Fluss verfolgen.

RGB-Schnitt des Datenkubus. In der linken BildhÀlfte kann man den mÀandrierenden Fluss verfolgen.

In einigen FĂ€llen ermöglicht die QualitĂ€t der seismischen Daten wesentlich klarere Bilder. Dies hĂ€ngt von der Methode der Feldarbeiten, der eingesetzten AusrĂŒstung und den angewandten Algorithmen zur GerĂ€uschunterdrĂŒckung ab. In solchen FĂ€llen sind nicht nur Fragmente von Flusssystemen sichtbar, sondern auch ganze ausgedehnte alte FlĂŒsse.
RGB-Mischung der drei Komponenten des seismischen Datenkubus (horizontaler Schnitt). Die Tiefe betrÀgt etwa 2 km.

Satellitenbild des Flusses Wolga in der NĂ€he von Saratow.

Fazit
In Wolfram Mathematica können seismische Daten analysiert und praktische Aufgaben im Zusammenhang mit der Suche nach BodenschĂ€tzen gelöst werden. Das Paket GeologyIO macht diesen Prozess einfacher. Die Struktur der seismischen Daten ist so, dass die Verwendung integrierter Methoden zur Beschleunigung von Berechnungen (, ,âŠ) sehr effizient ist und groĂe Datenmengen verarbeitet werden können. Dies wird erheblich durch die Speichereigenschaften des GeologyIO-Pakets unterstĂŒtzt. Ăbrigens kann das Paket nicht nur im Bereich der angewandten seismischen Untersuchung eingesetzt werden. Praktisch dieselben Datentypen werden in der Georadaruntersuchung und Seismologie verwendet. Wenn Sie VorschlĂ€ge zur Verbesserung der Ergebnisse haben, welche Signalverarbeitungsalgorithmen aus dem Sortiment von Wolfram Mathematica auf solche Daten anwendbar sind oder wenn Sie kritische Anmerkungen haben â hinterlassen Sie bitte Ihre Kommentare.
Quelle: habr.com
