Automatisches stochastisches 3D-Tonfraktionsmodell aus tTEM-Vermessungs- und Bohrlochdaten

Wissenschaftliche Berichte Band 12, Artikelnummer: 17112 (2022) Diesen Artikel zitieren

752 Zugriffe

1 Zitate

Details zu den Metriken

In den meisten urbanisierten und landwirtschaftlich genutzten Gebieten Mitteleuropas besteht der flache Untergrund aus quartären Ablagerungen, die oft die am intensivsten genutzten Schichten sind (Wasserpumpen, flache Geothermie, Materialaushub). Alle diese Lagerstätten sind oft komplex miteinander verflochten, was zu einer hohen räumlichen Variabilität und hohen Komplexität führt. Geophysikalische Daten können eine schnelle und zuverlässige Informationsquelle über den Untergrund sein. Dennoch kann die Integration dieser Daten zeitaufwändig sein, es fehlt eine realistische Interpolation im vollständigen 3D-Raum und die endgültige Unsicherheit wird oft nicht dargestellt. In dieser Studie schlagen wir eine neue Methodik vor, um Bohrlöcher und geophysikalische Daten mit Unsicherheit in einem automatischen Rahmen zu kombinieren. Eine räumlich variierende Übersetzerfunktion, die den Tonanteil anhand des spezifischen Widerstands vorhersagt, wird unter Verwendung der Bohrlochbeschreibung als Kontrollpunkte invertiert. Es wird mit einem 3D-stochastischen Interpolationsrahmen kombiniert, der auf einem Multiple-Point-Statistics-Algorithmus und einer Gaußschen Zufallsfunktion basiert. Dieser neuartige Arbeitsablauf ermöglicht eine robuste Einbeziehung der Daten und ihrer Unsicherheit und erfordert weniger Benutzereingriffe als die bereits vorhandenen Arbeitsabläufe. Die Methodik wird anhand von bodengestützten Schleppdaten (tTEM) und Bohrlochdaten aus dem oberen Aaretal in der Schweiz veranschaulicht. An diesem Standort wurde ein realistisches 3D-Modell mit hoher räumlicher Auflösung der Tonfraktion über das gesamte Tal erstellt. Der sehr dichte Datensatz ermöglichte den Nachweis der Qualität der vorhergesagten Werte und ihrer entsprechenden Unsicherheiten mittels Kreuzvalidierung.

Flache alluviale quartäre Grundwasserleiter werden häufig zur Wasserversorgung oder zur oberflächennahen Nutzung geothermischer Energie genutzt. In diesem Zusammenhang müssen häufig vielfältige damit verbundene Fragen geklärt werden, beispielsweise die Bewertung der Grundwasserressourcen, die Untersuchung der Schadstoffmigration, die Bewertung von Wechselwirkungen mit Oberflächenwasser oder die Vermeidung einer Überlappung von Einflusszonen zwischen benachbarten Geothermiebohrungen. All diese Fragen können nur dann richtig beantwortet werden, wenn die Struktur und die internen Heterogenitäten dieser Grundwasserleiter modelliert werden.

Der Aufbau dieser Modelle erfolgt häufig in mehreren Schritten1,2. Bei kleinräumigen Modellen ist die Verwendung von Bohrlöchern als einzige Datenquelle üblich. Allerdings ist ein solcher Ansatz oft blind für den Großteil der räumlichen Heterogenität und kann daher zu unzureichenden Modellen und falschen Schlussfolgerungen führen. Bohrlöcher sind nur eine Informationsquelle, um Rückschlüsse auf die lokale und vertikale Verteilung der Fazies zu ziehen. Sie sind oft nur begrenzt hilfreich, um komplexe 3D-Strukturen abzuschätzen. Wenn das Interessengebiet groß ist, ist es oft schwierig, zeitaufwändig und teuer, die Anzahl der Bohrlöcher zu erhöhen, um die Unsicherheit auf ein akzeptables Maß zu reduzieren. Eine Lösung besteht darin, kostengünstigere geophysikalische Daten und Bohrlochdaten zu kombinieren. Geophysikalische Daten werden häufig manuell interpretiert und zu einem Strukturmodell kombiniert, das dann mit lithologischen oder Faziessimulationen und schließlich mit physikalischen Eigenschaften gefüllt wird. Ein solcher Arbeitsablauf hat sich beispielsweise bei der Erstellung geologischer Modelle aus elektromagnetischen oder seismischen Luftdaten als effizient erwiesen3,4. Doch diese Modellierungsschritte sind komplex und erfordern oft unterschiedliche Software. Auch wenn einige stochastische Methoden häufig verwendet werden4,5, breiten sich die Unsicherheiten nicht immer über den gesamten Arbeitsablauf aus. Oftmals werden einige Schritte als deterministisch betrachtet und das resultierende geologische Modell ist dasjenige, das den größten Teil des verfügbaren umfassenden Wissens abdeckt6,7,8,9. Wenn schließlich mit einem großen Datensatz gearbeitet wird, kann die manuelle Erstellung des Strukturmodells unter Verwendung sowohl von Bohrlöchern als auch geophysikalischer Daten zeitaufwändig sein.

Daher besteht Bedarf an einem automatischen Ansatz, der mehrere Datentypen integrieren und strukturelle oder parametrische Modelle erstellen kann. Beispielsweise könnte die schnelle Generierung von 3D-Tonmodellen mit einem automatischen Algorithmus von großem Nutzen für lokale Behörden sein, die oft nicht über die Kapazitäten verfügen, die Daten manuell zu integrieren.

In den meisten Fällen handelt es sich bei den entlang von Bohrlöchern verfügbaren Daten jedoch hauptsächlich um lithologische Beschreibungen und nicht um Widerstands- oder Dichtewerte. Aber diese physikalischen Parameter sind es, die in die geophysikalischen Daten einfließen. Aufgrund der Vielzahl von Faktoren, die den Widerstand beeinflussen10, und der oft unvollständigen Beschreibung der lithologischen Fazies ist es schwierig, spezifische Widerstände und Lithologien direkt zu verknüpfen. Der grundlegendste Zusammenhang ist das Archie-Gesetz11, das den spezifischen Widerstand mit der Sättigung, der Wasserleitfähigkeit, der Tortuosität und der Porosität verknüpft. Diese empirische Beziehung basiert auf der Annahme, dass die Matrix nicht leitend ist. Diese Annahme gilt nicht mehr, sobald an ihrer Oberfläche stromleitende Tonminerale vorhanden sind. Diese Porenoberflächenleitfähigkeit hängt von der Oberfläche, der Korngröße, der Tonart und dem Tongehalt ab. Die Schätzung all dieser Parameter, die eine hohe räumliche Variabilität aufweisen können, ist schwierig. Eine aktuelle Rezension12 weist ebenfalls auf Skalierungsprobleme hin. Die meisten empirischen Laborgesetze werden im Kernmaßstab gemessen, wobei die Probe im Bereich von einem Dutzend Kubikzentimetern liegt. Und die Hochskalierung solcher Gesetze auf die Feldskala ist nicht einfach. Schließlich sind die meisten lithologischen Beschreibungen im Zusammenhang mit Bohrlöchern qualitativ und nicht quantitativ12. Außerdem ist es viel zu einfach, eine Funktion anzuwenden, die die lithologische Beschreibung mit einem einzelnen Widerstand verknüpft. Jede lithologische Beschreibung kann mit einem breiten Spektrum an spezifischen Widerständen in Verbindung gebracht werden, wobei es einige Überschneidungen zwischen verschiedenen Lithologien gibt13.

Eine Antwort auf dieses Problem besteht darin, die Wahrscheinlichkeit zu definieren, dass eine bestimmte lithologische Fazies vom Widerstandswert abhängig ist. Diese Wahrscheinlichkeitsverteilungen können aus einer Stichprobe von Bohrlöchern und übereinstimmenden Widerstandsmodellen geschätzt werden13. Allerdings berücksichtigt ein solcher Ansatz nicht die räumliche Verteilung der Bohrlöcher und ignoriert die mögliche räumliche Abhängigkeit zwischen der Art des Sediments und seinem spezifischen Widerstand. Die Wahrscheinlichkeitsverteilungsfunktion (PDF) wird aus Bohrlöchern und Widerstandsmodellen im gesamten Gebiet berechnet. Wenn eine bestimmte lithologische Fazies in einem Teilbereich der Domäne immer einen höheren Widerstand aufweist, wird die über die gesamte Domäne berechnete PDF dies nicht widerspiegeln.

Um dieses Problem zu lösen, schlugen Foged et al.14 eine Methode vor, die auf der Umkehrung einer räumlich variierenden Übersetzerfunktion zwischen spezifischem Widerstand und CF basiert, anstatt zu versuchen, eine einzelne PDF zu schätzen. Die Funktion liefert die beste Übereinstimmung zwischen dem beobachteten CF in Bohrlöchern und dem aus den Widerstandsmodellen berechneten CF. Diese Methode hat den Vorteil, dass sie nicht auf einer vorherigen Parameterschätzung beruht und nur aus beobachteten Daten abgeleitet wird. Es hat auch den Vorteil, dass die lithologische Beschreibung in eine kontinuierliche Variable umgewandelt wird, was die Hochskalierung ermöglicht und gleichzeitig die Co-Lokalisierung der Funktion berücksichtigt. Auch wenn mit dieser Methode der CF an der Position geschätzt werden kann, an der geophysikalische Widerstandsmodelle vorhanden sind, muss er dennoch interpoliert werden, um ein vollständiges kontinuierliches 3D-Modell zu erhalten.

Insofern verwendeten Vilhelmsen et al.15 nach der Vorhersage der CF einen geostatistischen Ansatz. Es gibt verschiedene geostatistische Methoden16,17,18, mit denen Daten interpoliert und realistische Modelle und Simulationen erstellt werden können. Ihre Verwendung ist in Disziplinen wie Risikobewertung, Ressourcenmanagement, Bergbau, Erdöltechnik oder geologischer Modellierung weit verbreitet. Multiple Point Statistics19,20 (MPS) ist eine dieser geostatistischen Techniken. Es handelt sich um eine nichtparametrische Methode, die auf der Verwendung eines Trainingsbildes (TI) beruht, um die räumliche Variabilität einer oder mehrerer Variablen abzuleiten. MPS hat sich in einer Vielzahl von Situationen als in der Lage erwiesen, realistische komplexe räumliche Variabilität zu erzeugen21,22,23,24,25. Vilhelmsen et al.15 gruppierten die CF in Einheiten und verwendeten sie als TI in einem MPS-Verfahren. Die harten Daten wurden als die Zonen definiert, in denen der Cluster am sichersten ist, hauptsächlich bei sehr niedrigem und hohem spezifischem Widerstand. Die anderen Bereiche werden nur durch weiche Daten eingeschränkt. Dies ermöglicht es, eine Art von Unsicherheit in den CF-Daten widerzuspiegeln und Variationen in unsicheren Modellzonen anzuzeigen. Diese Methode erfordert jedoch die Wahl des diskreten Schwellenwerts, bei dem die Zugehörigkeit zu einem Cluster sicher wird. Darüber hinaus ist die Ausbreitung von Unsicherheit aus den Daten nicht möglich. Die Unsicherheit über den spezifischen Widerstand oder CF geht nach der Clusterbildung verloren. Daher wird ein Wert, der nahe an den Clustergrenzen liegt, aber mit einer großen Unsicherheit verbunden ist, als genauso sicher angesehen wie ein Wert, der zweifelsfrei in einen Cluster eintreten würde. Darüber hinaus können aufgrund der Form der Übersetzerfunktion kleine Änderungen des Widerstandswerts in der Übergangszone zwischen Ton und Nicht-Ton große Auswirkungen auf den geschätzten CF haben. Da schließlich die Unsicherheit der Daten als Normalisierung in der Zielfunktion der CF-Inversion verwendet wird, argumentieren wir, dass sie auch bei der Anwendung der Funktion und der Interpolation der Ergebnisse berücksichtigt werden sollte. Dies zeigt die Grenzen der Verwendung eines deterministischen TI, wenn die Unsicherheit berücksichtigt werden sollte.

In diesem Artikel schlagen wir einen erweiterten Arbeitsablauf zur automatischen Generierung eines 3D-Tonfraktionsmodells mit einer robusten Unsicherheitsübertragung von den Daten auf das endgültige Modell vor. Die Arbeit ist wie folgt aufgebaut. Wir stellen zunächst die drei Hauptschritte der Methodik vor: die deterministische CF-Inversion, das stochastische Interpolationsrahmenwerk und schließlich die Kreuzvalidierungsimplementierung. Anschließend präsentieren wir die Anwendung der Methodik für einen bodengestützten, geschleppten transienten elektromagnetischen (tTEM) geophysikalischen Datensatz26, der im oberen Aaretal in der Schweiz erfasst wurde.

Die vorgeschlagene Methodik ermöglicht die automatische Generierung eines 3D-CF-Modells eines vollständigen quartären Tals. Ein wichtiger Aspekt dieses Ansatzes besteht darin, sich auf die Daten selbst zu verlassen, um die meisten Parameter automatisch abzuleiten. In den folgenden Unterabschnitten dieses Dokuments werden alle Schritte im Detail beschrieben. Aber lassen Sie uns zunächst einen kurzen Überblick über den gesamten Arbeitsablauf geben, wie in Abb. 1 dargestellt.

Überblick über die Hauptschritte der vorgeschlagenen Methodik. Der schattierte Bereich hebt den stochastischen Teil des Workflows hervor.

Die Eingaben sind einerseits ein dichter Datensatz von Widerstandsprotokollen, die durch Umkehrung der geophysikalischen Messungen (hier tTEM-Daten) erhalten wurden, und andererseits ein viel spärlicherer Datensatz von geologischen Protokollen, die zur Schätzung der Referenz verwendet werden können CF entlang der Bohrlöcher. Um ein 3D-CF-Modell über die gesamte Domäne zu erstellen, invertieren wir zunächst eine CF-Übersetzerfunktion (Schritt 1 in Abb. 1) unter Verwendung der Methode von Foged et al.14 und beider Datensätze. Die Übersetzerfunktion wird dann zur Schätzung des CF aus dem spezifischen Widerstand für alle verfügbaren geophysikalischen Sondierungen (einschließlich Standorte, an denen keine Bohrlochdaten verfügbar sind) verwendet. Die wichtigste Neuerung des vorgeschlagenen Arbeitsablaufs ist die Methodik zur Interpolation der resultierenden CF-Werte in 3D (Schritt 2 in Abb. 1). Im Vergleich zu früheren Arbeiten gibt es zwei wesentliche Verbesserungen. Zum einen soll vermieden werden, die Daten klassifizieren zu müssen. die CF-Werte werden direkt verwendet. Der zweite und vielleicht wichtigste Punkt besteht darin, dass wir die verschiedenen Quellen von Unsicherheiten berücksichtigen und sie im gesamten Arbeitsablauf verbreiten. Hierzu erfolgt die Interpolation mit einem Multiple Point Statistics-Algorithmus, wobei die CF-Daten selbst in 3D als Trainingsbild (TI) verwendet werden. Die dichte geophysikalische Abdeckung ermöglicht die direkte Nutzung der Daten, ohne dass externe Informationen hinzugefügt werden müssen. Um Unsicherheiten sowohl in den Bohrlöchern als auch in den tTEM-Daten zu berücksichtigen, verwenden wir Gaußsche Zufallsfunktionen zur Simulation von Fehlerkarten, die vor den MPS-Simulationen in die Konditionierungsdaten und TI-Bilder einbezogen werden. Abschließend (Schritt 3 in Abb. 1) wird eine Kreuzvalidierung durchgeführt, um die Qualität der Ergebnisse zu überprüfen. Ein Teil der Konditionierungsdaten wird bei der Interpolation nicht verwendet und mit den Ergebnissen verglichen. Anschließend bewerten wir die Leistung der Methode sowohl im Hinblick auf den vorhergesagten Wert als auch auf die vorhergesagte Unsicherheit.

Dieser Teil unserer Methodik (Schritt 1 in Abb. 1) folgt eng der Arbeit von Christiansen et al.27, erweitert durch Foged et al.14. Wir wenden ihre Idee an, eine Funktion umzukehren, die die Lücke zwischen spezifischem Widerstand und CF füllt. Da die Methode in diesen früheren Veröffentlichungen ausführlich beschrieben wurde, wird hier nur eine kurze Zusammenfassung gegeben. Dieser Ansatz wurde ursprünglich zur Bewertung des Nitratkontaminationsrisikos27 entwickelt und als Accumulated Clay Thickness (ACT) bezeichnet. Später wurde eine Erweiterung vorgeschlagen, um den CF mit derselben Methodik zu schätzen14. Die Hauptannahme ist, dass die größten Änderungen des spezifischen Widerstands in gesättigten und nicht verfestigten Lagerstätten durch Schwankungen in der Tonmenge verursacht werden. Vor diesem Hintergrund kann eine Übersetzerfunktion konstruiert werden, die den spezifischen Widerstand und den CF verknüpft.

Allerdings müssen einige Herausforderungen besprochen werden. Erstens sind die aus den geophysikalischen Daten abgeleiteten Widerstandsmodelle aufgrund der Grundfläche des Instruments, der Auflösungsfähigkeit und der Inversionsverfahren immer eine geglättete Version der Realität. Zweitens kann sich die Zusammensetzung des Tonsediments selbst ändern und an verschiedenen Orten mit ähnlichen lithologischen Beschreibungen unterschiedliche Widerstandswerte verursachen. Dies bedeutet, dass eine Übersetzerfunktion, die nur auf dem spezifischen Widerstand als Eingabe basiert, den Zusammenhang zwischen spezifischem Widerstand und hydraulischer Leitfähigkeit nicht richtig darstellen kann.

(a) Darstellung der automatischen Anpassung der Übersetzerfunktion an einen Datenpunkt. (b) Beispiel eines lithologischen Protokolls der USCS-Hauptkomponente. (c) Zugehöriger geschätzter CF in einem regelmäßigen Schichtgitter (\(CF_{log}\)), abgeleitet aus der lithologischen Beschreibung des Bohrlochs. Die Linie ist der vorhergesagte CF unter Verwendung der Widerstandsmodelle und der invertierten Übersetzerfunktionen. (\(CF_{rho}\)) (d) Widerstandsmodell aus der glatten Inversion der tTEM-Daten an der Bohrlochposition. Für vier Ebenen wird die verwendete invertierte Übersetzerfunktion angezeigt.

Foged et al.14 schlugen daher ein Verfahren vor, das auf der Invertierung zweier Parameter einer Übersetzerfunktion basiert, die einen CF aus einem Widerstandswert vorhersagt und es diesen Parametern ermöglicht, im Raum zu variieren. Mit anderen Worten: Zwei identische Widerstandswerte, die weit voneinander entfernt sind, entsprechen nicht unbedingt demselben CF-Wert. Die grundlegenden Eingabedaten in diesem Verfahren sind geologische Bohrlochprotokolle, wie in Abb. 2b dargestellt, und geophysikalische Widerstandsprofile an denselben Orten, wie sie in Abb. 2d als rote Kurve dargestellt sind. Die beiden Datensätze wurden mit unabhängigen Methoden erfasst und liegen nebeneinander.

Basierend auf der geologischen Beschreibung wird der Tonanteil \(CF_{log}\) entlang des Bohrlochs geschätzt, indem das geologische Protokoll in regelmäßige Tiefenintervalle unterteilt und der Anteil an lithologischem Ton für jedes Intervall berechnet wird. Da die geologische Beschreibung oft qualitativ ist, ist der resultierende CF unsicher und das Ergebnis ist ein Mittelwert und eine Spanne von CFs für jedes Intervall. Dieses Ergebnis ist in Abb. 2c als graue Kästchen dargestellt.

Der nächste Schritt besteht darin, eine parametrisierte Übersetzerfunktion zu definieren:

Dabei ist erfc die komplementäre Fehlerfunktion, \(\rho\) der spezifische Widerstand und \(m_{up}\) und \(m_{low}\) die Widerstandswerte, bei denen die Funktion ein Gewicht von 0,975 zurückgibt und 0,025 bzw. Die Parameter \(m_{up}\) und \(m_{low}\) werden während der Inversion identifiziert und können als Widerstandsgrenzen nur für Ton/nur Sand betrachtet werden. Abbildung 2a zeigt die Invertierung von \(m_{up}\) und \(m_{low}\) für einen Datenpunkt. Vor der Inversion stimmt der vorhergesagte \(CF_{rho}\) nicht mit der beobachteten Tonmenge \(CF_{log}\) überein, die im zusammengefügten Bohrloch beschrieben wurde, was zu einer Fehlanpassung führt. Durch die Inversion werden \(m_{up}\) und \(m_{low}\) angepasst, um die Fehlanpassung zu reduzieren.

Das vorherige Beispiel zeigt das Prinzip für einen einzelnen Datenpunkt, das Problem ist jedoch in 3D korreliert und zielt darauf ab, eine globale Fehlanpassung zu minimieren. Zwei benachbarte Übersetzerfunktionen dürfen keine drastisch unterschiedlichen Parameter haben. Abbildung 2d zeigt beispielsweise, wie sich die Übersetzerfunktion für ein einzelnes Bohrloch mit der Tiefe ändert. Die Funktion ist in Abb. 2d für vier verschiedene Tiefen dargestellt, jedoch für alle Schichten definiert.

Um ein solches Ergebnis zu erhalten, folgt das Inversionsverfahren den folgenden Schritten. Zunächst wird der globale Datenrest als durchschnittlicher quadratischer Fehler zwischen der vorhergesagten CF (\(CF_{rho}\)) und der gemessenen CF in den Bohrlöchern (\(CF_{log}\)) definiert:

Der Datenrest wird durch \(\sigma ^2\), die kombinierte Varianz von \(CF_{rho}\) und \(CF_{log}\) und N die Anzahl der Konditionierungsdaten normalisiert. Anschließend werden Regularisierungsparameter hinzugefügt, um sicherzustellen, dass benachbarte Punkte keine scharfen und unrealistischen Variationen aufweisen. Die Regularisierung begrenzt die räumliche Variation von \(m_{up}\) und \(m_{low}\). Die Regularisierung wird ausgedrückt als:

Dabei ist \(N_{con}\) die Anzahl der Parameterpaare, A die Differenz im logarithmischen Raum zwischen den beiden untersuchten Parameterpaaren und \(e_i\) ein abstandsabhängiger Faktor. Je weiter zwei Punkte voneinander entfernt sind, desto größer kann ihre Differenz sein. Schließlich lautet die vollständige Zielfunktion:

Die optimalen Parameter (\(m_{up}\) und \(m_{low}\)) werden durch Minimierung von Q mithilfe eines iterativen Gauß-Newton-Schemas mit einer Marquardt-Modifikation27 erhalten. Die Übereinstimmung zwischen dem vorhergesagten \(CF_{rho}\) und \(CF_{log}\) ist nicht perfekt, da die Parameter der Übersetzerfunktion nicht nur von den am selben Ort befindlichen Daten (Bohrlochwiderstandsmodelle), sondern auch von beeinflusst werden alle benachbarten (über die Regularisierung). Die Inversion sucht nach einem globalen Minimum.

In dem in 2 gezeigten Beispiel können wir darstellen, wie sich die Übersetzerfunktion anpasst, um Strukturen zu identifizieren, die im Widerstandsmodell nicht gut dargestellt werden, wie beispielsweise die flache Tonschicht.

Am Ende des Schritts zur Schätzung des Tonanteils wird jedem Widerstandsprotokoll ein geschätztes \(CF_{rho}\)-Protokoll zugeordnet, das durch Anwendung der Übersetzerfunktion mit den optimalen Parametern \(m_{up}\) und \(m_ {niedrig}\).

Der letzte Schritt besteht darin, die Unsicherheit von den Widerstandswerten auf den CF zu übertragen. Tatsächlich wurden die spezifischen Widerstände durch eine geophysikalische Inversion ermittelt, die es ermöglicht, die Unsicherheit der spezifischen Widerstände abzuschätzen. Die Übersetzerfunktionen können dann nicht nur auf die Widerstandsdaten, sondern auch auf deren Unsicherheiten angewendet werden, was zu einem \(CF_{rho}\) an den tTEM-Erfassungspunkten mit Unsicherheit führt. Für den nächsten Simulationsschritt wird \(CF_{log}\) als Konditionierungsdaten für die Zellen verwendet, in denen wir Bohrlochinformationen haben. \(CF_{rho}\) wird an anderer Stelle verwendet.

Um ein vollständiges 3D-Modell von CF und Widerstand zu erstellen, müssen die im vorherigen Schritt erhaltenen Werte interpoliert werden, um den Raum abzudecken, in dem keine geophysikalischen oder Bohrlochdaten erfasst wurden. Dieser Teil des Arbeitsablaufs entspricht Schritt 2 in Abb. 1 und ist der wichtigste neuartige Teil des vorgeschlagenen Arbeitsablaufs. Wir verwenden dabei die Direct Sampling MPS-Methode28. Der Hauptvorteil der Verwendung eines MPS-Ansatzes besteht darin, dass er automatisch die räumlichen Muster der regionalen Strukturen aus dem im vorherigen Schritt erstellten sehr dichten Datensatz lernen und diese zur Darstellung der räumlichen Variabilität und Unsicherheit in den interpolierten Bereichen verwenden kann. Bei der Standardanwendung von MPS-Techniken wird jedoch davon ausgegangen, dass TI und Hard Data (HD) deterministisch sind. Hier schlagen wir eine Methode vor, um noch einen Schritt weiter zu gehen und die Unsicherheit in diesen Eingabedaten zu berücksichtigen. Wir beziehen daher nichtdeterministische Hard Data (HD) und Trainingsbilder (TI) in den MPS-Algorithmus ein. Die Gesamtmethode, wie in Schritt 2 in Abb. 1 dargestellt, besteht darin, den MPS-Algorithmus viele Male mit unterschiedlichem TI und unterschiedlichem HD anzuwenden, um ein Ensemble interpolierter 3D-Modelle von CF zu erzeugen. Aus diesem Simulationsensemble werden wir Wahrscheinlichkeitsverteilungen für CF und Widerstand an jedem Ort im 3D-Bereich ableiten.

Das allgemeine Prinzip von MPS-Algorithmen besteht darin, ein Simulationsgitter iterativ zu füllen und gleichzeitig die Muster des TI zu reproduzieren. Der Direct-Sampling-Algorithmus28 ist ein vielseitiger MPS-Algorithmus, der auf einer randomisierten und bedingten Neuabtastung des TI basiert. In diesem Artikel verwenden wir die Deesse-Implementierung29 des Direct Sampling. Um eine Simulation zu generieren, wird das 3D-Gitter zunächst mit den verfügbaren Hard Data gefüllt. Der Algorithmus besucht dann zufällig alle verbleibenden Orte im Raster. Um einen Wert an einem bestimmten Ort zu simulieren, werden die n nächsten HD- und bereits simulierten Werte extrahiert, um ein Datenmuster zu definieren. Das Datenmuster wird dann verwendet, um nach Standorten mit ähnlichen Mustern im TI zu suchen. Während dieser Suche wird der Abstand zwischen dem Datenmuster und den im TI gefundenen Mustern berechnet. Wenn der Abstand unter einem Schwellenwert (t) liegt, werden die beiden Muster als ähnlich betrachtet und der Pixelwert an der fehlenden Stelle wird aus dem TI übernommen und in das Simulationsgitter kopiert. Um den Algorithmus zu beschleunigen und das direkte Kopieren und Einfügen der TI zu vermeiden, wird nur ein Bruchteil (f) der TI gescannt. Die drei Parameter n, t und f werden vorab vom Benutzer ausgewählt. Meerschman et al.30 bieten einige praktische Empfehlungen für die Auswahl dieser Parameter. Eine Stärke des Direct-Sampling-Algorithmus besteht darin, dass er kontinuierliche und mehrere Variablen gleichzeitig verarbeiten kann. Dies bedeutet, dass man einen TI bereitstellen kann, der mehrere Variablen enthält, und der Algorithmus eine oder mehrere Variablen simulieren kann, die von den anderen verfügbaren Daten abhängig sind. Dies wird in mehreren Artikeln ausführlich beschrieben28,31. Es kann verwendet werden, um das Vorhandensein von Trends im Trainingsbild und im Interpolationsgitter mithilfe von Hilfsvariablen zu beschreiben22,24.

Um CF und Widerstand zu interpolieren, haben wir DeeSse verwendet und ein Simulationsproblem mit vier Variablen entworfen: zwei Hilfsvariablen und zwei Hauptvariablen. Jeder von ihnen hat eine bestimmte Anzahl von Nachbarn n und einen Schwellenwert t, während der gescannte Bruchteil f allen Variablen gemeinsam ist. Die beiden sekundären Variablen sind (1) die Tiefe der Zelle und (2) der Nordwert der Zelle. Ihr Zweck besteht darin, unser Vorwissen über räumliche Trends zu beschreiben. Diese Wahl spiegelt unsere Erwartung wider, bestimmte Muster in einer bestimmten Tiefe oder einem bestimmten Bereich in unserem interpolierten Bereich besser zu finden. Beispielsweise werden die tiefsten Strukturen im TI tendenziell am tiefsten in der Simulation reproduziert. Die von uns ausgewählten Schwellenwerte sind jedoch hoch und n ist klein, was den Algorithmus nur in der Mitte einschränkt. Wir haben uns auch dafür entschieden, nur die Nordrichtung einzubeziehen, da wir davon ausgehen, dass die meisten Variationen des Musters entlang des Tals auftreten (ungefähr von Nord nach Süd ausgerichtet). Diese Annahme basiert auf einer visuellen Inspektion der Widerstandsmodelle, geologischen Kenntnissen und Bohrlöchern. Die wichtigsten Simulationsvariablen sind (3) die logarithmische Transformation des spezifischen Widerstands und (4) der Tonanteil. Für diese Variablen haben wir eine höhere Anzahl von Nachbarn n und einen kleineren Schwellenwert t gewählt, da dies die interessierenden Variablen sind. Alle vier dieser Variablen werden simuliert. Wir aktivieren auch die Gaußsche Pyramide-Option des DeeSse-Codes29. Auf diese Weise werden die räumlichen Muster in den Daten auf mehreren Skalen analysiert und modelliert. Der Algorithmus verwendet Gaußsche Filter, um eine Pyramide aus nebeneinander angeordneten Bildern mit gröberem Maßstab zu erstellen, die gemeinsam für Training und Simulation verwendet werden32. Diese Methode verbessert die Qualität der Simulation zwischen dicht bedeckten und spärlicheren Gebieten und kann gleichzeitig auf unterschiedliche Variationsmaßstäbe reagieren.

Um den MPS-Algorithmus auszuführen und die CF zu interpolieren, müssen wir einen 3D-TI und HD bereitstellen. Bei den HD handelt es sich lediglich um die punktuellen 3D-Daten, die aus Schritt 1 der Bohrlöcher und geophysikalischen Sondierungen abgeleitet wurden. Da die räumliche Dichte der tTEM-Daten sehr hoch ist (siehe Beispielanwendung unten), ist dieser Datensatz oft sehr dicht und kann direkt als TI verwendet werden. Diese Situation entspricht dem sogenannten Lückenfüllproblem, bei dem wir nur einige Teile eines bereits sehr dichten Datensatzes interpolieren müssen. MPS und insbesondere die direkte Probenahme haben sich bei diesen Problemen als sehr effizient erwiesen31,33. In diesen Fällen werden die gleichen Daten wie bei HD und TI verwendet. Dieser Modellierungsentscheidung liegt die Annahme zugrunde, dass die Datenabdeckung ausreicht, um die räumlichen Statistiken der Variablen, die interpoliert werden müssen, ordnungsgemäß darzustellen.

In der Praxis werden die 3D-HD-Punkte vorab im 3D-Simulationsgitter platziert. Bei diesen Punkten handelt es sich um Schätzungen der CF- und Widerstandswerte, die entweder aus den Bohrlöchern oder aus der Anwendung der Übersetzerfunktion auf die tTEM-Daten stammen. Die Hilfsvariablen (Tiefe und Hochwert) werden ebenfalls für jeden Standort im Raster berechnet und gespeichert. An diesem Punkt der Methodik könnten wir einfach den MPS-Algorithmus anwenden und ein Ensemble stochastischer Simulationen erhalten, die die Unsicherheit aufgrund der Interpolation darstellen. Dies wird normalerweise bei der Verwendung von MPS durchgeführt. Der TI und der HD sind deterministisch.

Mit den CF- und Widerstandswerten ist jedoch eine Unsicherheit verbunden, die bereits im vorherigen Schritt des Arbeitsablaufs geschätzt wurde. Um auch dieser Unsicherheitsquelle Rechnung zu tragen, bestünde eine Möglichkeit darin, nur die sichersten Punkte als HD zu betrachten, zum Beispiel die extrem niedrigen oder hohen spezifischen Widerstände, die mit ziemlicher Sicherheit jeweils mit völlig tonhaltigen bzw. nicht tonhaltigen Punkten verbunden sind. Diese Idee hat jedoch zwei Hauptnachteile. Zunächst muss eine Auswahl getroffen werden, um die oberen und unteren Widerstandsgrenzen zu bestimmen. Dies widerspräche der Idee, ein möglichst automatisiertes Verfahren zu implementieren. Zweitens muss die Wahrscheinlichkeitsverteilung der Werte im TI der HD-Verteilung ähnlich sein. Wenn sie unterschiedlich sind, neigen die Simulationen dazu, die Ton- und Nicht-Ton-Punkte in der Simulation zu überrepräsentieren. Im Allgemeinen würden wir nicht alle verfügbaren HDs nutzen und könnten einer Zelle einen spezifischen Widerstand oder einen CF zuweisen, der völlig außerhalb des aus den Feldmessungen abgeleiteten Unsicherheitsbereichs liegt. Wir brauchen daher eine bessere Möglichkeit, die Unsicherheit dieser Daten zu berücksichtigen.

Wir überwinden diese Herausforderungen durch eine kombinierte Verwendung von MPS und der Gaußschen Zufallsfunktion (GRF), die es uns ermöglicht, MPS mit einem Trainingsbild (TI) und einem harten Datensatz (HD) durchzuführen, die nicht deterministisch sind. GRF-Modelle sind bekannt18,34. Die räumliche Variabilität wird mithilfe parametrischer Multi-Gauß-Verteilungen modelliert. Diese Modelle werden mit einem Kovarianz- oder Variogrammmodell definiert, das die räumliche Variabilität darstellt. Es können mehrere Realisierungen mit oder ohne Konditionierungsdaten generiert werden. In der vorgeschlagenen Methodik wird das GRF-Modell verwendet, um die Unsicherheit der CF- und Widerstandsdatenmessungen darzustellen. In der Praxis werden für jede MPS-Simulation von CF und spezifischem Widerstand unterschiedliche TI und HD generiert.

Prinzip für die Generierung eines Ensembles aus harten Daten und Trainingsbildern, um die Unsicherheit auf dem HD zu berücksichtigen.

Abbildung 3 skizziert das allgemeine Prinzip der Generierung der TIs. Der ursprüngliche HD (und TI) wird durch das Hinzufügen von Rauschen im Bereich der geschätzten Unsicherheit für diese Datenpunkte leicht gestört. Die Simulation des Lärms erfolgt mittels bedingungsloser 3D-GRF-Simulationen. Die GRF-Simulationen haben einen Mittelwert von 0 und eine Varianz von 1. Um die mögliche räumliche Korrelation des Rauschens zu berücksichtigen, wird automatisch ein anisotropes 3D-Variogramm an die TI-Daten angepasst. Wir verwenden einen Trust-Region-Reflective-Algorithmus35, um die Schwelle und die Mehrrichtungsbereiche zweier Beiträge (eines Gaußschen und eines exponentiellen) zu optimieren. Die Schwelle wird dann normalisiert, um eine Varianz von 1 zu erhalten. Die Werte von CF und Widerstand für den 3D-TI (und 3D-HD) bei jeder Iteration werden dann durch Addition des HD mit einem neu skalierten Rauschen definiert:

wobei \(\text {TI}_{i}(x,y,z)\) der TI-Wert für die Simulation i an der Position (x, y, z) in 3D ist, \(\text {data}(x Dabei ist ,y,z)\) der deterministische Schätzwert der interessierenden Variablen an diesem Ort, \(\sigma (x,y,z)\) die Standardabweichung der geschätzten Unsicherheit und \(GRF_{i}( x,y,z)\) ist der unbedingte simulierte GRF-Wert an Position (x, y, z) für Iteration i. Der endgültige simulierte Wert (z. B. CF) variiert innerhalb seines Unsicherheitsbereichs. In dem in Abb. 3 dargestellten Beispiel wird der aus der geschätzten CF abgeleitete anfängliche TI leicht modifiziert, indem ein korreliertes Zufallsrauschen hinzugefügt wird, das von den lokalen Unsicherheiten abhängt. Der gleiche Vorgang wird für die geophysikalischen und Bohrlochdatensätze unter Verwendung desselben GRF-Modells durchgeführt.

Zusammenfassend lässt sich sagen, dass für jede Simulation ein eindeutiger 3D-TI und der entsprechende 3D-HD-Satz generiert und als Eingabe an den Direct-Sampling-Algorithmus übergeben werden, der die fehlenden Teile im 3D-Gitter simuliert. Dadurch entsteht ein Ensemble von 3D-Realisierungen. Anschließend können wir alle Simulationen stapeln und den Mittelwert und die Standardabweichung für jeden Standort berechnen. Die Zellen, in denen wir Konditionierungsdaten haben, haben den gleichen Mittelwert und die gleiche Standardabweichung wie die ursprünglichen Konditionierungsdaten.

Um die Leistung der vorgeschlagenen Methodik zu überprüfen, haben wir einen Kreuzvalidierungsschritt implementiert. Die Kreuzvalidierung ermöglicht die Quantifizierung der mit einem Modell verbundenen Fehler. Das Prinzip ist einfach: Es wird eine Teilmenge der Daten erstellt; es enthält einen zufälligen Teil des Originaldatensatzes. Anschließend wird die stochastische Interpolation nur unter Verwendung der Teilmenge der HD-Daten angewendet und die resultierende Simulation mit den ausgeschlossenen Daten verglichen. Es können verschiedene Fehlerindikatoren verwendet werden36. In unserem Fall erzeugt eine zufällige Stichprobe der dichten CF-Daten keine ausreichend großen Lücken, um eine repräsentative Fehlerschätzung zu erstellen. Die fehlenden Punkte werden durch die benachbarten zu stark eingeschränkt. Um eine größere Störung in den Daten zu erzeugen, haben wir jedem Punkt eine Gruppe zugewiesen, basierend auf einem 3D-k-Mittelwert-Clustering37. Die Clusterbildung erfolgt anhand der räumlichen Koordinaten (x,y,z) des Modells. Der Zweck dieses Schritts besteht darin, mehrere räumliche Gruppen ähnlicher Größe zu erstellen, die zufällig von den Simulationen ausgeschlossen werden. In unserem Fall haben wir 28 Gruppen definiert. Der simulierte Wert und die Standardabweichung dieser Zonen werden dann mit dem realen Wert verglichen. Der Fehler und der normalisierte Fehler für jeden Punkt im Modell sind definiert als

wobei \(\varepsilon\) und \(\varepsilon _{norm}\) jeweils der Fehler und der normalisierte Fehler sind, n die Gesamtzahl der Simulationen ist, true der Datenpunkt ist, der nicht im interpolierten Datensatz enthalten ist, \(sim_i\ ) ist der simulierte Wert und Sigma ist die Standardabweichung von \(sim_i\) über die n Simulationen, die diese Datenpunkte nicht umfassten. Diese beiden Indikatoren werden punktuell berechnet. Der normalisierte Fehler ist ein wichtiger Indikator, da er zeigt, wie gut wir die Unsicherheit einschätzen, wobei das Ziel ein Verhältnis nahe 1 ist.

(a) Verteilung der Schichtqualität und -tiefe nach der Normalisierung. (b) Übersicht über die Abdeckung der tTEM-Daten und allgemeine Lagekarte (erstellt mit QGIS V3.22.9, qgis.org. Basiskarte ist kostenlos beim Bundesamt für Landestopographie swisstopo erhältlich).

Das Untersuchungsgebiet liegt in der Zentralschweiz und umfasst einen Abschnitt von rund 20 km des Oberen Aaretals zwischen den Städten Thun und Bern. Das Obere Aaretal weist eine für Alpentäler typische quartäre Geologie auf: Das Grundgebirge ist meist einige hundert Meter tief und mit einem komplexen Geflecht aus glazialen, fluvio-glazialen und fluvialen Ablagerungen bedeckt38,39. Im Schweizer Becken wurden mehrere Zyklen von Gletschervorstößen und -rückgängen identifiziert, die zu mehrfachen Veränderungen in den Ablagerungsprozessen führten40. Im besonderen Fall des Oberen Aaretals konnten in Bohrlöchern mindestens vier Glazialzyklen identifiziert werden41. Eine vollständige Beschreibung der Lithostratigraphie allein anhand von Bohrlochbeschreibungen ist jedoch nahezu unmöglich, da ähnliche Ablagerungen unterschiedlichen Alters übereinander liegen und miteinander verflochten sein können. Überall im Tal ist ein Oberflächengrundwasserleiter vorhanden, der aus den jüngsten alluvialen Ablagerungen besteht. Der Grundwasserspiegel liegt zwischen einigen Dutzend Zentimetern und bis zu 2 Metern unter der Oberfläche. Allgemeine Kenntnisse über den Grundwasserleiter legen nahe, dass die Mächtigkeit des oberen Grundwasserleiters je nach Gebiet zwischen 4 und 20 m liegt. Darüber hinaus wurde in einigen tiefen Bohrlöchern ein zweiter tieferer Grundwasserleiter identifiziert, der vom flachen Grundwasserleiter durch eine Tonschicht getrennt ist. Die genaue Ausdehnung, Konnektivität und Mächtigkeit beider Grundwasserleiter ist nicht genau bekannt. Allerdings sind in der Gegend 6.350 Pumpbrunnen (flache Erdwärme oder Trinkwasser) und 5.300 Injektionsbrunnen (Wiederinjektion von Wasser nach der Erdwärmepumpe) im Einsatz. In diesem Zusammenhang wäre die Erstellung eines 3D-Modells für die lokalen Behörden von großem Nutzen, um die Gefährdung des oberen Grundwasserleiters einzuschätzen. Aufgrund der Grundwasserspiegelhöhe kann nahezu in der gesamten Höhe des Modells von gesättigten Verhältnissen ausgegangen werden.

Im Untersuchungsgebiet sind 1542 Bohrlöcher lithologisch beschrieben. Das obere Aaretal war auch eines der Testgebiete, die für das Geoquat-Projekt des Schweizerischen Landestopografischen Instituts39 ausgewählt wurden. In diesem Zusammenhang führten sie eine Digitalisierung und Standardisierung der Bohrlochdaten durch. Alle Bohrlochbeschreibungen wurden in Standard-USCS-Beschreibungen umgewandelt42. Darüber hinaus wurde für jede Schicht ein QC-Wert hinzugefügt, der die Zuverlässigkeit der geologischen Informationen beurteilt: Jeder Schicht wurde eine Note von eins bis fünf zugewiesen, abhängig von der Genauigkeit der Beschreibung. Einer entspricht einer Schicht mit nur einer grundlegenden Beschreibung, fünf entspricht einer Schicht, in der Labormessungen durchgeführt wurden und eine vollständige mehrphasige lithologische Beschreibung verfügbar ist. Eine Schätzung des Tongehalts in den Bohrlöchern wurde anhand der USCS-Richtlinien vorgenommen und die Unsicherheit entsprechend dem Gehalt skaliert. Eine schlecht beschriebene Tonschicht (Qualität 1) hat beispielsweise einen CF-Wert von \(1 \pm 0,5\). Wir gehen davon aus, dass bis zu 50 % der Tonschicht aus Nicht-Ton-Material bestehen könnte. Andererseits hat ein gut beschriebenes (Qualität 5) einen CF-Wert von \(1 \pm 0,08\). Fehlende, künstliche oder unbeschriebene Schichten werden nicht berücksichtigt.

Abbildung 4a zeigt die Verteilung der Qualität der Schichten und deren zugehörige Tiefe. Da der Großteil der Bohrlochexploration entweder zur oberflächennahen geothermischen Exploration oder zu geotechnischen Zwecken durchgeführt wird, stellen wir fest, dass etwa 60 % der Daten über 10 m Tiefe liegen und dass der Anteil ungefähr einer mit der Tiefe abnehmenden Potenzgesetzverteilung folgt.

Im Januar 2020 wurde im oberen Aaretal eine bodengestützte transiente elektromagnetische Schleppuntersuchung (tTEM) durchgeführt. Die Daten wurden unter Verwendung verschiedener Regularisierungen (scharf und glatt) invertiert und die resultierenden Widerstandsmodelle werden in dieser Studie verwendet. Alle Verarbeitungsabläufe, Details zur Regularisierung und Daten werden in Neven et al.26 beschrieben. Der Datensatz umfasst etwa 1500 Hektar mit einem Linienabstand von 25 Metern. Die Probenahmefrequenz nach der Verarbeitung beträgt etwa 1 Widerstandsmodell alle 10 Meter entlang der Leitungen. Die Invertierung der Daten erfolgte mit dem Inversionscode AarhusInv (https://hgg.au.dk/software/aarhusinv)43. Das durchschnittliche Residuum der Inversion beträgt 0,52 für die scharfe Regularisierung, was bedeutet, dass wir tendenziell eine hervorragende Übereinstimmung zwischen den vorhergesagten Daten aus unseren Widerstandsmodellen und den Feldmessungen haben. Um die Untersuchungstiefe abzuschätzen, haben wir die Jacobi-Sensitivitätsmatrix der letzten Iteration44 verwendet. Auf diese Weise können wir die genaue Tiefe ermitteln, in der jedes Widerstandsmodell in den aufgezeichneten Daten nur unzureichend dargestellt wird. Unterhalb dieser Tiefe sind die Widerstandsmodelle ausgeblendet. Für eine vollständige Beschreibung der Invertierung des Datensatzes und der Qualitätsprüfung verweisen wir die Leser wiederum auf Neven et al.26.

Die Umkehrung der Übersetzerfunktion wurde am scharfen Inversionsdatensatz tTEM20AAR26 durchgeführt. 57'862 Widerstandsmodelle aus 30 Schichten wurden berücksichtigt und mit der Standarduntersuchungstiefe geblendet. Abbildung 4b zeigt die geophysikalische Datenabdeckung im Tal nach der Verarbeitung. Die meisten Felder außerhalb der Städte sind kartiert. Städte sowie die Gebiete rund um Straßen und Bahnlinien sind deutlich sichtbar, da die elektromagnetische Kopplung in solchen Umgebungen den Einsatz einer induktiven Methode verbietet und unbedeckt bleibt. Die Unsicherheit des spezifischen Widerstands wurde mithilfe der Kovarianzmatrix der letzten Iteration45 geschätzt. Auch wenn diese Methode die Unsicherheit, die aus einer stochastischen Inversion abgeleitet werden kann, nicht vollständig ersetzt, hat sie den Vorteil, dass sie die relative Unsicherheit auf den Zellen des Modells widerspiegeln kann und gleichzeitig relativ schnell zu berechnen ist.

Ergebnisse der invertierten Übersetzerfunktion an 50'000 TEM-Punkten. Darüber hinaus wurde die gemittelte Unsicherheit mit einem gleitenden Mittelwert von 0,05 CF Breite berechnet.

Der durchschnittliche Datenfehler [siehe Gl. (2)] oder Residuum der Anwendung der Übersetzerfunktion auf die geophysikalischen Modelle beträgt 0,38. Dies bedeutet, dass die Differenz zwischen dem vorhergesagten \(CF_{rho}\) und \(CF_{log}\) mehr als das Zweifache innerhalb des Unsicherheitsbereichs beträgt. Die Hauptursache für Fehlanpassungen ist das Auflösungsproblem. Aufgrund des Fußabdrucks der geophysikalischen Ausrüstung und der gedämpften Umkehrung der kleinsten Quadrate, die an den TEM-Daten durchgeführt wurde, tendieren wir dazu, einen sanfteren Übergang in den geophysikalischen Daten zu haben als in der Realität. Da wir die Ergebnisse mit Bohrlochdaten vergleichen, die tendenziell diese scharfen Übergänge widerspiegeln, können sie eine Ursache für eine Fehlanpassung zwischen dem vorhergesagten \(CF_{rho}\) und dem tatsächlichen \(CF_{log}\) sein. Einige dünne Schichten aus Widerstandsmaterialien in einer großen leitfähigen Schicht werden von der Auflösung des Geräts nicht erfasst, erhöhen jedoch den Restwert. Der Hauptvorteil der14-Methodik besteht jedoch darin, dass sich die Übersetzerfunktion anpasst, um diese Effekte auszugleichen. Abbildung 5 zeigt die Ergebnisse der invertierten Übersetzerfunktion und die gemittelte Unsicherheit an 50.000 TEM-Punkten. Der größte Teil der Unsicherheit konzentriert sich auf die zentralen Punkte, die einen mittleren spezifischen Widerstand aufweisen, da eine kleine Unsicherheit über den spezifischen Widerstand in diesem Bereich die vorhergesagte CF stark beeinflussen wird. Der Wertebereich, in dem der Übergang zwischen Ton und Nicht-Ton mit der höchsten Wahrscheinlichkeit stattfindet, liegt zwischen 90 und 40 Ohm. Solche Werte stimmen mit der erwarteten manuellen Interpretation überein, da der quartäre Ton von Aare einen spezifischen Widerstand zwischen 10 und 50 Ohm aufweist. Die Modelle werden auf die aus den Widerstandsmodellen ermittelte Standarduntersuchungstiefe geblendet und dann an den geostatistischen Simulationsschritt weitergeleitet.

Die hier vorgestellten Ergebnisse basieren auf 200 Simulationsschleifen mit jeweils 10 Realisierungen. Für jede Schleife wird mit einer GRF-Simulation ein neuer TI- und HD-Satz generiert und ein anderer Zufallsteilsatz für die Kreuzvalidierung gezogen. Die Gruppen wurden mithilfe eines K-Means-Algorithmus aus den räumlichen Koordinaten der CF-Daten gebildet. Wie bereits erwähnt, wird der Deesse-Algorithmus durch drei Hauptparameter gesteuert: den gescannten Anteil, den Schwellenwert und die Anzahl der Nachbarn. Der TI-Scan-Anteil ist allen Variablen gemeinsam und wurde auf 5 % festgelegt. Das Muster wird mit dem TI verglichen, bis der gescannte Bruchteil erreicht ist oder bis der Schwellenwert für alle Variablen eingehalten wird. Der Schwellenwert ist für die Hilfsvariablen auf 10 % auf drei benachbarten Knoten und für die Hauptvariablen auf 1 % auf 24 Knoten festgelegt. Wenn der Schwellenwert nie erreicht wird, wird der beste Kandidat ausgewählt und der Punkt wird am Ende der Simulation zur erneuten Simulation markiert. Der Scanpfad im TI ist zufällig. Das Ergebnis ist ein Satz von 2000 Simulationen, basierend auf 200 verschiedenen TI. Das resultierende 3D-Modell hat eine Zellengröße von 50 m x 50 m und eine vertikale Auflösung von 2 m. Die abgedeckte Fläche beträgt 35 km2, bei einer Rechenzeit von 10 h auf einem CPU-Cluster. Vertikale Querschnitte des gemittelten Modells sind in Abb. 6 dargestellt. Ein klarer Trend im Modell ist zwischen der Nord- und der Südseite des Modells zu erkennen. Eine solche Variation war aufgrund der unterschiedlichen Form des Tals zu erwarten und wird durch die Bohrlöcher bestätigt. Die Unsicherheit der Daten ist in den tieferen Zellen des Modells höher, wo keine Bohrlöcher und keine geophysikalischen Daten vorhanden sind, die die Simulationen einschränken. Sie bezeichnen auch die Übergangszonen zwischen Ton- und Nicht-Ton-Bereichen und spiegeln die Unsicherheit über die genaue Tiefe des Übergangs wider. Die Form der unterirdischen Strukturen steht im Einklang mit bestehenden quartären Ablagerungsprozessen. Allerdings reicht die visuelle Konsistenz nicht aus, um einem Modell zu vertrauen.

3D-Ansicht von Schnitten durch das Tonfraktionsmodell mit der damit verbundenen Unsicherheit. Die Z-Skala ist übertrieben.

Histogramm der Fehler und des normalisierten Fehlers für die beiden Hauptvariablen: Widerstand (Rho) und Tonanteil (CF).

Der Fehler und die normalisierten Fehlerverteilungen [siehe Gl. (6) und (7)], die mit Kreuzvalidierung berechnet wurden, sind in Abb. 7 dargestellt. Während der Simulation sind die beiden Hauptvariablen der spezifische Widerstand und der Tonanteil. Die Fehlerschätzung während des Kreuzvalidierungsschritts wird für beide durchgeführt. Der Fehler liegt für die beiden Variablen bei Null, was darauf hindeutet, dass wir im Durchschnitt dazu neigen, einen korrekten Wert ohne Verzerrung vorherzusagen. Darüber hinaus gelang es uns in 44 % bzw. 37 % der Fälle, einen Wert vorherzusagen, der im 10 %-Bereich des tatsächlichen Werts lag. Bezogen auf den normalisierten Fehler beträgt der Mittelwert des normalisierten Fehlers für den spezifischen Widerstand 0,98 (\(\sigma = 0,56\)). Eine solche Verteilung legt nahe, dass die Standardabweichungen der Simulationen in der gleichen Größenordnung liegen wie der tatsächliche Fehler. Die Unsicherheit ist daher gut vorhersehbar. Andererseits beträgt der Mittelwert des CF-normalisierten Fehlers 0,89 (\(\sigma = 0,5\)), was bedeutet, dass wir dazu neigen, den Fehler in den simulierten Daten im Vergleich zu den realen Daten leicht zu überschätzen. Aber insgesamt zeigt die Kreuzvalidierung, dass im Durchschnitt der korrekte Wert mit einer Standardabweichung simuliert wird, die die mögliche Unsicherheit gut widerspiegelt. Ein überzeugendes Ergebnis ist, dass der spezifische Widerstand genau vorhergesagt wird, was zeigt, dass diese Methode sogar für die reine Interpolation von Widerstandskarten von großem Nutzen sein kann.

Um alle geologischen Daten aus den Quartärformationen der Schweiz zu homogenisieren und zu digitalisieren, erstellte das GeoQuat-Projekt einen Prototyp und eine Demonstrationsstudie zum Aaretal. In diesem Zusammenhang wurde ein deterministisches geologisches Modell des Gebiets erstellt, das eine manuelle Interpretation von Bohrlöchern, geophysikalischen Daten und geologischem Wissen durch die Verwendung geologischer Querschnitte39 nutzte. Beim Vergleich der beiden Modelle ist zu berücksichtigen, dass in das stochastische Modell neue Daten einfließen, die zum Zeitpunkt der Erstellung des deterministischen Modells nicht verfügbar waren: Selbst wenn die Bohrlochdatenbank dieselbe ist, wurden die tTEM-Daten nur erfasst im Jahr 2020. Nach der Einbeziehung der Bohrlöcher korrelierten sie manuell Einheiten und Fazies und nutzten die Interpolation des nächsten Nachbarn, um das Modell auf ein vollständiges 3D-Volumen zu erweitern39.

Vergleich zwischen dem bestehenden geologischen Modell und diesem Studienmodell. Die Lage der Schnittpunkte der Abschnitte 3 und 4 ist skizziert.

Abbildung 8 zeigt den Vergleich einiger Querschnitte im deterministischen Modell und im mit MPS generierten Modell. Das deterministische geologische Modell und die Bohrlöcher werden anhand ihrer USCS-Primärkomponenten angezeigt. Die Transparenz des MPS-Modells spiegelt die Unsicherheit wider. Je höher die Transparenz, desto höher die damit verbundene Unsicherheit. Zumindest für großräumige Strukturen stimmen beide Modelle mit den Bohrlochdaten überein. Wir können sehen, dass das MPS-generierte geologische Modell Strukturen anzeigt, die geologisch gesehen realistischer sind. Zu den Querschnitten können einige Anmerkungen gemacht werden:

Abschnitt 1 ist übereinstimmend. Das Modell wird durch zwei tiefe Bohrlöcher, die sich durch das gesamte Modell erstrecken, gut eingegrenzt. Auch wenn die Interpolation des nächsten Nachbarn des deterministischen Modells einige unrealistische Blockformen aufweist, sind die globalen Strukturen zwischen den beiden Modellen ähnlich.

Abschnitt 2 unterscheidet sich drastisch zwischen dem deterministischen und dem stochastischen Modell. Das Fehlen von Bohrlöchern zur Einschränkung des deterministischen Modells führt zur Ausbreitung einer Kiesschicht aus entfernten Bohrlöchern, die an diesem Standort höchstwahrscheinlich nicht vorhanden ist. Der aus den tTEM-Daten berechnete invertierte Widerstand ist gering und hängt höchstwahrscheinlich mit dicken Tonkörpern zusammen. Ein solcher Abschnitt zeigt aufgrund der zusätzlichen Daten eine Verbesserung der Modelle.

Abschnitt 3 verdeutlicht die Bedeutung der Unsicherheit. In dem Gebiet wurden keine tTEM-Daten erfasst und beide Modelle basieren auf denselben Daten. Die oberen paar Meter der Modelle sind bei beiden Modellen eindeutig als oberer Kieskörper zu erkennen. Die Verwendung der vorgeschlagenen Methode ergab jedoch einen möglichen Widerstandskörper (Sand oder Kies) in tieferen Schichten. Die Unsicherheit ist hoch, wie die Farbpalette zeigt. Das Vorhandensein dieser Schicht kann nicht sicher festgestellt werden. Obwohl es sich wahrscheinlich um eine Möglichkeit handelt, die bei der Beurteilung des Gebiets berücksichtigt werden muss.

Abschnitt 4 zeigt, wie eine lokale Quartärstruktur das deterministische Modell drastisch beeinflussen kann. Das 100 m tiefe Bohrloch in der Mitte des Querschnitts weist auf eine Kiesschicht hin, die im Vergleich zu den anderen Bohrlöchern in der Nähe des Abschnitts ungewöhnlich dick ist. Diese übermäßige Verdickung ist äußerst lokal und ist laut persönlicher Mitteilung eines Quartärgeologen höchstwahrscheinlich auf eine Massenbewegung durch Schwerkraft vor einigen Dutzendtausend Jahren zurückzuführen. Aufgrund dieser speziellen Bohrung überschätzte das deterministische geologische Modell die Mächtigkeit des oberen Grundwasserleiters in einem sehr weiten Bereich.

Schließlich konnten wir die Qualität des Modells testen, indem wir die beobachtete Geologie mehrerer neu erstellter Bohrlöcher verglichen, die in keinem der beiden Modelle berücksichtigt wurden, da sie erst kürzlich gebohrt wurden. Die neuen flachen Bohrlöcher (maximal 10 m) stimmen gut mit beiden Modellen überein: 85 % bzw. 94 % für das deterministische und stochastische Modell. Diese Ergebnisse sind nicht überraschend, da die obere Kiesschicht im gesamten Gebiet mehr oder weniger vorhanden ist. Die meisten Fehler sind auf eine geringfügige Über- oder Unterschätzung der Übergangstiefe zurückzuführen. Es wurden jedoch auch vereinzelt Tiefbohrungen gebohrt. In Abb. 9 wird ein neues Bohrlochprotokoll mit den beiden Modellen verglichen. Das Bohrloch liegt mitten im Dorf Wichtrach, wo keine tTEM-Daten erfasst werden können. Das stochastische Modell spiegelt die im Bohrloch beobachtete Geologie gut wider. Das stochastische Modell sagt das Vorhandensein des Kiesgrundwasserleiters an der Oberfläche und hohe Widerstände wie das deterministische Modell bis zu einer Tiefe von etwa 25 m voraus. Doch dann sagt das stochastische Modell das Vorhandensein einer massiven Tonschicht bis zu einer Tiefe von 80 m voraus, bevor es in tieferen Tiefen auf eine große Unsicherheit hinweist. Diese Vorhersage stimmt gut mit den geologischen Beobachtungen überein, die ebenfalls eine massive Tonschicht zeigen. Im Gegenteil, das deterministische Modell sagte eine massive Kiesschicht direkt unter einer dünnen Tonschicht voraus, gefolgt von Schluff. Dies wurde im Bohrloch nicht beobachtet. Beide Modelle werden durch die gleichen Bohrlochdaten in der Umgebung eingeschränkt, aber während das stochastische Modell den regionalen Trend durch die Integration aller Daten des Gebiets ableitet, ist das deterministische Modell nur durch nahegelegene Punkte und Querschnitte eingeschränkt.

Vergleich zwischen einer neuen Bohrung im Dorf Wichtrach (in keinem der beiden Modelle enthalten) und der entsprechenden Zelle im stochastischen Modell und im deterministischen geologischen Modell.

In dieser Studie haben wir eine Methodik vorgestellt, die es ermöglicht, die Unsicherheit von den geophysikalischen Daten und Bohrlochprotokollen bis zum endgültigen Modell der Tonfraktion zu übertragen. Darüber hinaus werden alle Eingabedaten innerhalb ihres Unsicherheitsbereichs als harte Daten verwendet. Ein verrauschter Datenpunkt im Feld führt zu einer unsichereren Zelle im Modell. Diese Funktion ist ein wichtiger Punkt, um eine robuste Fehlerschätzung für das endgültige Modell zu erhalten. Darüber hinaus ist die Integration unsicherer Bohrlöcher aufgrund schlechter oder fehlender Beschreibungen, unsicherer Positionierung oder fehlender Daten ein großer Vorteil. Der Einsatz von MPS ermöglicht die Generierung von Modellen mit komplexen Mustern, die mit Ansätzen wie Kriging nicht reproduziert werden können. Diese Methode hat auch den Vorteil, dass sie sich auf die Daten selbst verlässt, um ohne Vorkenntnisse räumliche Muster abzuleiten. Dies gilt natürlich nur für Regionen, in denen die Datenabdeckung ausreichend ist. In einer Region mit geringerer Datenabdeckung oder zur Steigerung der Modellqualität kann jedoch auch Vorwissen durch den Einsatz von Hilfsvariablen hinzugefügt werden. Nach der Datenaufbereitung kann die Erstellung hochauflösender Valley-Scale-Modelle in wenigen Tagen erreicht werden, verglichen mit Monaten bei deterministischen Modellen.

In Bezug auf die geophysikalische Technik hat die Anwendung von tTEM im oberen Aaretal bewiesen, dass die Technik zuverlässig ist und relevante Informationen über die 3D-Struktur des Untergrunds liefern kann. Dies bestätigt die vorherige Schlussfolgerung von Sandersen et al.46. Darüber hinaus argumentieren wir jedoch, dass eine ordnungsgemäße und effiziente Möglichkeit zur Integration dieser Art großer 3D-Datensätze mit Bohrlochdaten nur durch automatische Methoden erreicht werden kann. Wir gehen auch davon aus, dass die vorgeschlagene Methodik für die Integration dichter luftgestützter Datensätze sehr nützlich sein könnte.

Allerdings hat die Abschätzung der Lithologie anhand des spezifischen Widerstands auch ihre Grenzen. Erstens: Wenn das Gebiet erhebliche Schwankungen im Salzgehalt aufweist, kann dies zu Nichteinzigartigkeit führen. Der gleiche spezifische Widerstand kann gesättigtem Widerstandsmaterial oder Ton entsprechen. Um die Beschaffenheit des Untergrunds zu ermitteln, würden wir uns dann nur noch auf Bohrlöcher stützen. Die vorgeschlagene Methode müsste um eine 3D-Schätzung des Salzgehalts als zusätzliche Hilfsvariable erweitert werden. Mithilfe zusätzlicher hydrogeologischer Informationen kann dies möglich sein. Darüber hinaus ist die vorgeschlagene Methode effizient, um zwischen Widerstandsschichten (Sand oder Kies) und leitfähigen Schichten (Ton) zu unterscheiden, aber der Widerstandsunterschied zwischen Sand und Kies ermöglicht es der Methode nicht, diese klar zu unterscheiden. Das ist der Hauptgrund, warum wir uns entschieden haben, uns ausschließlich auf die Identifizierung der Tonschichten zu konzentrieren.

Auch ohne Simulation des CF-Werts sollte unser Ansatz bei der Interpolation von 1D- oder 2D-Widerstandsmodellen in 3D berücksichtigt werden, insbesondere wenn es um scharfe Inversionen geht. Traditionell erfolgt dies mittels Kriging. Da Kriging jedoch der beste unverzerrte lineare Schätzer ist, führt es zu einer reibungslosen Interpolation zwischen scharfen Modellen und fügt zwangsläufig Artefakte in den 3D-Widerstandsmodellen hinzu. Im Zusammenhang mit quartären Ablagerungen kann der Bereich des Variogramms manchmal kleiner sein oder dem Abstand zwischen zwei Erfassungsgebieten entsprechen, was zu einer falschen oder unvollständigen Interpolation führt. Die vorgeschlagene Methode hat den Vorteil, dass sie sich direkt auf die 1D-Sondierungen verlässt und nicht auf ein durch Kriging erstelltes Raster und interpoliert zwischen den Modellen mit kohärenter Schärfe. Eine ähnliche von Vilhemlmsen et al.15 vorgeschlagene Methodik basiert auf einem Indikator-Kryping auf dem Clusterwert mit einem kurzen Bereich, um die Lücken zwischen den Erfassungslinien im TI zu füllen, was zu einer möglichen Überglättung führen kann. Die hier vorgestellte Methodik hat den Vorteil, dass sie vor der Simulation nicht auf eine Vorinterpolation der Daten angewiesen ist.

Natürlich ist die Inversion selbst bereits eine Interpretation der Daten, und darüber hinaus ist die mit den geophysikalischen Invertierungsmodellen verbundene Unsicherheit oft schwer genau abzuschätzen. Eine weitere Verbesserung könnte die Entwicklung eines Algorithmus sein, der die geophysikalischen Messungen direkt und ohne Inversion bearbeitet, anstatt invertierte Modelle zu verwenden. Auf diese Weise könnten wir die Fehler zwischen den realen Felddaten und dem endgültigen CF-Modell direkt vergleichen. Es würde auch den Arbeitsablauf vereinfachen, indem eine doppelte Inversion vermieden würde, einmal für das Widerstandsmodell selbst und einmal für den CF-Algorithmus.

Der Vergleich zwischen dem bestehenden Modell und dem automatischen Modell (Abb. 8) zeigt das Risiko der Verwendung eines deterministischen Ansatzes. Für Querschnitt 2 sagt das deterministische Modell einen dicken Kieskörper in etwa 40 m Tiefe voraus. Dieser Körper wurde dort nur auf der Grundlage von Bohrlöchern platziert, die Hunderte von Metern entfernt sind (in der Nähe von Querschnitt 1), und nach Einbeziehung geophysikalischer Daten sind wir nun zuversichtlich, dass es dort keinen solchen Körper gibt. Natürlich handelte es sich um eine Interpretation, die auf den damals verfügbaren Daten basierte, es gab jedoch keine Hinweise auf Unsicherheit. In diesem Zusammenhang ist es riskant, Entscheidungen auf der Grundlage eines solchen Modells zu treffen und sollte vermieden werden.

Immer mehr Länder entwickeln zentrale Datenbanken für geologische Messungen. In diesem Zusammenhang ist der Einsatz eines agilen und zuverlässigen Datenaggregationsalgorithmus ein vielversprechender Ansatz. Die Möglichkeit, räumliche Muster aus Daten ohne oder mit geringem Vorwissen abzuleiten, verhindert die Einfügung von Strukturen aus willkürlichen Entscheidungen, die normalerweise bei einer manuellen Interpretation geophysikalischer Modelle erfolgen. Die Standardisierung von Beschreibungsmethoden durch unterschiedliche Bohrunternehmen und Open-Access-Daten sind die Schlüsselpunkte, um die Anwendung dieser Methode schneller und anpassungsfähiger für neue Feldgebiete zu machen.

Eine Verbesserung der Methode könnte darin bestehen, Permeabilitätswerte vorherzusagen. Wenn Permeabilitätsmessungen durchgeführt werden, sei es durch Packungs- oder Pumptests, könnte der Simulation eine neue Variable hinzugefügt werden, um Permeabilitätsfelder direkt zu erzeugen. Eine andere Möglichkeit wäre die Definition einer gekoppelten Wahrscheinlichkeitsdichtefunktion zwischen spezifischem Widerstand, Tonanteil und Permeabilität aus einer Reihe von Vorkenntnissen oder Feldkenntnissen. Dann wären wir in der Lage, direkt parametrische Felder zu erzeugen. Schließlich glauben wir, wie bereits erwähnt, dass die hier vorgestellte Interpolationsmethode auf andere geophysikalische 3D-Modelle angewendet werden könnte, um Teilkarten zu füllen. Aufgrund von Feldeinschränkungen (unzugängliche Bereiche, beschädigte Daten, unterschiedliche Erfassungen usw.) sind einige Bereiche häufig weniger dicht abgedeckt als andere. Unter Verwendung nichtdeterministischer MPS mit Gaußschen Pyramiden könnten vollständige und homogene parametrische Modelle mit geeigneter Unsicherheitsquantifizierung generiert werden.

In diesem Artikel haben wir gezeigt, dass ein neuartiger Arbeitsablauf, der den Tonfraktionsschätzungsalgorithmus14 mit unserem modifizierten Multiple-Points-Statistic-Algorithmus kombiniert, eine robuste Methode für die automatische Generierung eines 3D-Tonfraktionsmodells ist. Die resultierenden 3D-Modelle können von lokalen Behörden oder Projektmanagern genutzt werden, um die Entwicklung und unterirdische Nutzung der Region besser zu planen. Mithilfe dieser Modelle können beispielsweise potenziell hochdurchlässige Zonen für die Grundwassergewinnung oder geothermische Entwicklungen lokalisiert oder das mögliche Vorhandensein von Baumaterialien im Untergrund beurteilt werden. Da die Methode Unsicherheitsschätzungen liefert, kann sie auch bei der Gestaltung der Erfassung weiterer Daten hilfreich sein.

Unsere Methode hat den Vorteil, dass sie datengesteuert ist und nicht auf eine manuelle Interpretation der Strukturen angewiesen ist. Der Arbeitsablauf ist automatisch und umfasst: (1) die Umkehrung der Übersetzerfunktion für die Widerstandsmodelle und die Bohrlochdaten, (2) die automatisierte Anpassung des für das GRF-Modell verwendeten Variogramms und (3) die Generierung verschiedener TI und HD für jede MPS-Simulation. Der Benutzer muss nur wenige Parameter wie die Anzahl n der Nachbarn oder den Schwellenwert t für die MPS-Simulationen auswählen, es können jedoch Standardwerte verwendet werden. Darüber hinaus werden Strukturen reproduziert, die mit klassischer Zweipunktinterpolation wie Kriging nicht modelliert werden können. Durch den Vergleich der an neu erstellten Bohrlöchern beobachteten Geologie mit unserem Modell kann diese Studie den allgemeinen Trend stets gut vorhersagen, wie beispielsweise in Abb. 9 dargestellt. Der Arbeitsablauf kann problemlos in öffentliche Datenbanken integriert werden, sodass die Behörden das 3D-Modell regelmäßig automatisch aktualisieren können wenn neue Felddaten verfügbar sind.

Schließlich beweist die Arbeit auch die Effizienz der Towed Transient ElectroMagnetic (TEM)-Methode. Während in den meisten früheren Veröffentlichungen luftgestützte TEM-Daten verwendet wurden, haben wir Schlepp-TEM-Daten integriert. Dadurch konnten wir eine Auflösung im Raum und in der Tiefe erreichen, die derzeit mit keiner anderen geophysikalischen Methode erreichbar ist. Eine solche Datenerfassung sollte verstärkt in Betracht gezogen werden, auch bei Problemen mittlerer oder kleinerer Größenordnung. Eine einfache, schnelle und kostengünstige Integration mehrerer, oft vorhandener Datentypen kann unabhängig vom Zweck nur von Vorteil sein.

Die in diesem Artikel verwendeten tTEM-Daten sind im Open-Access-Archiv Zenodo unter der folgenden URL verfügbar: https://doi.org/10.5281/ZENODO.4269887. Der Datensatz ist in Neven et al.26 dokumentiert. Der Bohrlochdatensatz ist beim Schweizer Bundesamt für Landestopographie (Swisstopo) erhältlich. Es gelten jedoch Einschränkungen hinsichtlich der Verfügbarkeit dieser Daten, die für die aktuelle Studie unter Lizenz verwendet wurden und daher nicht öffentlich verfügbar sind. Daten sind jedoch auf begründete Anfrage bei Swistopo erhältlich.

Die in dieser Studie verwendeten geostatistischen Codes sind im Python-Geone-Paket implementiert, das unter https://github.com/randlab/geone frei verfügbar ist. Um alle Funktionen von DeeSse nutzen zu können, kann eine akademische Lizenz kostenlos erworben werden. Das CF-Modell wird frei verfügbar sein.

Ringrose, P. & Bentley, M. Reservoir Model Design (Springer, 2015).

Google Scholar

Pyrcz, MJ & Deutsch, CV Geostatistical Reservoir Modeling (Oxford University Press, 2014).

Google Scholar

Jørgensen, F. et al. Eine Methode zur kognitiven 3D-geologischen Voxelmodellierung von AEM-Daten. Stier. Ing. Geol. Umgebung. 72, 421–432. https://doi.org/10.1007/s10064-013-0487-2 (2013).

Artikel Google Scholar

Jørgensen, F., Høyer, A.-S., Sandersen, PB, He, X. & Foged, N. Kombination geologischer 3D-Modellierungstechniken, um Variationen in Geologie, Datentyp und Dichte zu berücksichtigen: Ein Beispiel aus Süddänemark. Berechnen. Geosci. 81, 53–63. https://doi.org/10.1016/j.cageo.2015.04.010 (2015).

Artikel Google Scholar

Wellmann, JF, Varga, MDI, Murdie, RE & Gessner, K. Unsicherheitsschätzung für ein geologisches Modell des Sandstein-Grünsteingürtels, Westaustralien: Erkenntnisse aus der integrierten geologischen und geophysikalischen Inversion in einem Bayes'schen Inferenzrahmen. Geol. Soc. 453, 41–56. https://doi.org/10.1144/sp453.12 (2017).

Artikel Google Scholar

Sophocleous, M. et al. Integrierte numerische Modellierung für das beckenweite Wassermanagement: Der Fall des Rattlesnake Creek-Beckens im Süden von Kansas. J. Hydrol. 214, 179–196. https://doi.org/10.1016/S0022-1694(98)00289-3 (1999).

Artikel ADS Google Scholar

Henriksen, HJ et al. Methodik zur Konstruktion, Kalibrierung und Validierung eines nationalen hydrologischen Modells für Dänemark. J. Hydrol. 280, 52–71. https://doi.org/10.1016/S0022-1694(03)00186-0 (2003).

Artikel ADS Google Scholar

Kollet, SJ & Maxwell, RM Integrierte Oberflächen-Grundwasserströmungsmodellierung: Eine Grenzbedingung für die Überlandströmung an der freien Oberfläche in einem parallelen Grundwasserströmungsmodell. Adv. Wasserressource. 29, 945–958. https://doi.org/10.1016/j.advwatres.2005.08.006 (2006).

Artikel ADS Google Scholar

Lemieux, J.-M., Sudicky, EA, Peltier, WR & Tarasov, L. Dynamik der Grundwasserneubildung und -versickerung über der kanadischen Landschaft während der Wisconsin-Eiszeit. J. Geophys. Res. 113, F01011. https://doi.org/10.1029/2007JF000838 (2008).

Artikel ADS Google Scholar

Linde, N., Binley, A., Tryggvason, A., Pedersen, LB & Revil, A. Verbesserte hydrogeophysikalische Charakterisierung durch gemeinsame Umkehrung des elektrischen Widerstands zwischen Löchern und Daten zur Laufzeit des Bodenradars. Wasserressource. Res. 42, 1–10. https://doi.org/10.1029/2006WR005131 (2006).

Artikel Google Scholar

Archie, G. Das elektrische Widerstandsprotokoll als Hilfe bei der Bestimmung einiger Lagerstätteneigenschaften. Trans. AIME 146, 54–62. https://doi.org/10.2118/942054-G (1942).

Artikel Google Scholar

Knight, R., Gottschalk, I. & Dewar, N. Gesteinsphysik im Feldmaßstab für oberflächennahe Anwendungen. Enzykl. Geol. 1, 884–899. https://doi.org/10.1016/B978-0-12-409548-9.12514-X (2021).

Artikel Google Scholar

Knight, R. et al. Kartierung von Grundwasserleitersystemen mit luftgestützter Elektromagnetik im zentralen Tal von Kalifornien. Grundwasser 56, 893–908. https://doi.org/10.1111/gwat.12656 (2018).

Artikel CAS Google Scholar

Foged, N., Marker, PA, Christansen, AV, Bauer-Gottwein, P. & Jørgensen, F. Groß angelegte 3D-Modellierung durch Integration von Widerstandsmodellen und Bohrlochdaten durch Inversion. Hydrol. Erdsystem. Wissenschaft. 18, 4349–4362. https://doi.org/10.5194/hess-18-4349-2014 (2014).

Artikel ADS Google Scholar

Vilhelmsen, TN et al. Kombination von Clustering-Methoden mit MPS zur Schätzung der strukturellen Unsicherheit für hydrologische Modelle. Vorderseite. Erdwissenschaft. 7, 181. https://doi.org/10.3389/feart.2019.00181 (2019).

Artikel ADS Google Scholar

Matheron, G. Prinzipien der Geostatistik. Wirtschaft. Geol. 58, 1246–1266. https://doi.org/10.2113/gsecongeo.58.8.1246 (1963).

Artikel CAS Google Scholar

De Marsily, G. et al. Umgang mit räumlicher Heterogenität. Hydrogeol. J. 13, 161–183. https://doi.org/10.1007/s10040-004-0432-3 (2005).

Artikel ADS Google Scholar

Chiles, J.-P. & Delfiner, P. Geostatistik: Modellierung räumlicher Unsicherheit Vol. 497 (Wiley, 2009).

MATH Google Scholar

Strebelle, S. Bedingte Simulation komplexer geologischer Strukturen mithilfe von Mehrpunktstatistiken. Mathematik. Geol. 34, 1–21. https://doi.org/10.1023/A:1014009426274 (2002).

Artikel MathSciNet MATH Google Scholar

Mariethoz, G. & Caers, J. Mehrpunkt-Geostatistik: Stochastische Modellierung mit Trainingsbildern (Wiley, 2014).

Buchen Sie Google Scholar

Strebelle, S., Payrazyan, K. & Caers, J. Modellierung eines Tiefsee-Trübungsreservoirs unter Berücksichtigung seismischer Daten unter Verwendung von Mehrpunkt-Geostatistik. SPE Annu. Technik. Konf. Ausstellung: https://doi.org/10.2118/77425-MS (2002).

Artikel Google Scholar

Pirot, G., Straubhaar, J. & Renard, P. Simulation von Zeitreihen für geflochtene Flusshöhenmodelle mit Mehrpunktstatistiken. Geomorphologie 214, 148–156 (2014).

Artikel ADS Google Scholar

de Carvalho, PRM, da Costa, JFCL, Rasera, LG & Varella, LES Geostatistische Faziessimulation mit geometrischen Mustern eines Erdölreservoirs. Stoch. Umgebung. Res. Risikobewertung. 31, 1805–1822. https://doi.org/10.1007/s00477-016-1243-5 (2017).

Artikel Google Scholar

Dall'Alba, V. et al. 3D-Mehrpunkt-Statistiksimulationen des kontinentalen Pliozän-Grundwasserleiters im Roussillon mit DeeSse. Hydrol. Erdsystem. Wissenschaft. 24, 4997–5013. https://doi.org/10.5194/hess-24-4997-2020 (2020).

Artikel ADS Google Scholar

Neven, A., DallAlba, V., Juda, P., Straubhaar, J. & Renard, P. Schätzung des Eisvolumens und der Basaltopographie mithilfe geostatistischer Methoden und bodendurchdringender Radarmessungen: Anwendung auf die Gletscher Tsanfleuron und Scex Rouge, Schweiz Alpen. Kryosphäre 1, 5169–5186. https://doi.org/10.5194/tc-15-5169-2021 (2021).

Artikel ADS Google Scholar

Neven, A., Maurya, PK, Christiansen, AV & Renard, P. ttem20aar: Ein geophysikalischer Benchmark-Datensatz für nicht konsolidierte fluvioglaziale Sedimente. Erdsystem. Wissenschaft. Daten 13, 2743–2752 (2021).

Artikel ADS Google Scholar

Christiansen, AV, Foged, N. & Auken, E. Ein Konzept zur Berechnung der akkumulierten Tondicke aus lithologischen Bohrlochprotokollen und Widerstandsmodellen zur Bewertung der Nitratanfälligkeit. J. Appl. Geophys. 108, 69–77. https://doi.org/10.1016/j.jappgeo.2014.06.010 (2014).

Artikel ADS Google Scholar

Mariethoz, G., Renard, P. & Straubhaar, J. Die direkte Stichprobenmethode zur Durchführung geostatistischer Mehrpunktsimulationen. Wasserressource. Res. 46, 1–14. https://doi.org/10.1029/2008WR007621 (2010).

Artikel Google Scholar

Straubhaar, J. DeeSse-Benutzerhandbuch. Technik. Rep. (Zentrum für Hydrogeologie und Geothermie (CHYN), Universität Neuenburg: Neuchâtel, 2019).

Meerschman, E. et al. Ein praktischer Leitfaden zur Durchführung statistischer Mehrpunktsimulationen mit dem Direktstichprobenalgorithmus. Berechnen. Geosci. 52, 307–324. https://doi.org/10.1016/j.cageo.2012.09.019 (2013).

Artikel ADS Google Scholar

Mariethoz, G., McCabe, MF & Renard, P. Raumzeitliche Rekonstruktion von Lücken in multivariaten Feldern unter Verwendung des direkten Stichprobenansatzes. Wasserressource. Res. 48, 12115. https://doi.org/10.1029/2012WR012115 (2012).

Artikel Google Scholar

Straubhaar, J., Renard, P. & Chugunova, T. Mehrpunktstatistiken unter Verwendung von Bildern mit mehreren Auflösungen. Stoch. Umgebung. Res. Risikobewertung. 1, 1–23. https://doi.org/10.1007/s00477-020-01770-8 (2020).

Artikel Google Scholar

Oriani, F., Borghi, A., Straubhaar, J., Mariethoz, G. & Renard, P. Simulation fehlender Daten innerhalb von Durchflusszeitreihen unter Verwendung von Mehrpunktstatistiken. Umgebung. Modell. Softw. 86, 264–276. https://doi.org/10.1016/j.envsoft.2016.10.002 (2016).

Artikel Google Scholar

Dietrich, CR & Newsam, GN Eine schnelle und exakte Methode für mehrdimensionale Gaußsche stochastische Simulationen. Wasserressource. Res. 29, 2861–2869. https://doi.org/10.1029/93WR01070 (1993).

Artikel ADS Google Scholar

Branch, MA, Coleman, TF & Li, Y. Eine Unterraum-, Innenraum- und konjugierte Gradientenmethode für großräumige Minimierungsprobleme mit beschränkten Einschränkungen. SIAM J. Sci. Berechnen. 21, 1–23. https://doi.org/10.1137/S1064827595289108 (1999).

Artikel ADS MathSciNet MATH Google Scholar

Juda, P., Renard, P. & Straubhaar, J. Ein Rahmen für die Kreuzvalidierung kategorialer geostatistischer Simulationen. Erde-Weltraum-Wissenschaft. 7, 1152. https://doi.org/10.1029/2020ea001152 (2020).

Artikel Google Scholar

Hartigan, JA & Wong, MA Algorithmus AS 136: Ein K-Means-Clustering-Algorithmus. Appl. Stat. 28, 100. https://doi.org/10.2307/2346830 (1979).

Artikel MATH Google Scholar

Kellerhals, P., Haefeli, C. & Tröhler, B. Hydrogeologie Aaretal, zwischen Thun und Bern. Wasseru. Energiewirtschaftsamt des Kt. Bern (WEA) (Schweizerische Geologische Dokumentationsstele, Landeshydrologie und-geologie, 1981).

Google Scholar

Volken, S., Preisig, G. & Gaehwiler, M. GeoQuat: Entwicklung eines Systems für nachhaltiges Management, 3D-Modellierung und Anwendung quartärer Lagerstättendaten. Schweizer Bulletin für Angewandte Geologie. https://doi.org/10.5169/SEALS-658182 (2016).

Graf, HR & Burkhalter, R. Quartäre Ablagerungen: Konzept für eine stratigraphische Klassifikation und Nomenklatur – ein Beispiel aus der Nordschweiz. Schweizer J. Geosci. 109, 137–147. https://doi.org/10.1007/s00015-016-0222-7 (2016).

Artikel CAS Google Scholar

Schlüchter, C. Die vollständigste quartäre Aufzeichnung des Schweizer Alpenvorlandes. Paläogeogr. Paläoklimatol. Paläoökol. 72, 141–146. https://doi.org/10.1016/0031-0182(89)90138-7 (1989).

Artikel Google Scholar

Casagrande, A. Klassifizierung und Identifizierung von Böden. Trans. Bin. Soc. Bauingenieur. 113, 901–930. https://doi.org/10.1061/TACEAT.0006109 (1948).

Artikel Google Scholar

Auken, E. et al. Ein Überblick über einen äußerst vielseitigen Vorwärts- und stabilen Umkehralgorithmus für elektromagnetische und elektrische Daten in der Luft, am Boden und im Bohrloch. Entdecken. Geophys. 46, 223–235. https://doi.org/10.1071/eg13097 (2015).

Artikel ADS Google Scholar

Christiansen, AV & Auken, E. Ein globales Maß für die Untersuchungstiefe. Geophysik 77, 171–177. https://doi.org/10.1190/geo2011-0393.1 (2012).

Artikel Google Scholar

Alumbaugh, DL & Newman, GA Bildbewertung für 2D- und 3D-elektromagnetische Inversion. Geophysik 65, 1455–1467. https://doi.org/10.1190/1.1444834 (2000).

Artikel ADS Google Scholar

Sandersen, PBE et al. Nutzung der Towed Transient ElectroMagnetic-Methode (tTEM) zur Erzielung beispielloser oberflächennaher Details bei der geologischen Kartierung. Ing. Geol. 288, 106125. https://doi.org/10.1016/j.enggeo.2021.106125 (2021).

Artikel Google Scholar

Referenzen herunterladen

Diese Forschung wurde vom Schweizerischen Nationalfonds durch das Projekt Phenix (Fördernummer 182600) unterstützt. Die Autoren danken Nikolaj Foged für seine Expertise bei der Anwendung des ACT-Algorithmus.

Zentrum für Hydrogeologie und Geothermie, Universität Neuenburg, Neuenburg, Schweiz

Alexis Neven & Philippe Renard

Fachbereich Geowissenschaften, Universität Aarhus, Aarhus C, Dänemark

Anders Vest Christiansen

Fachbereich Geowissenschaften, Universität Oslo, Oslo, Norwegen

Philippe Renard

Sie können diesen Autor auch in PubMed Google Scholar suchen

AN führte die Datenanalyse durch, entwarf und implementierte die geostatistische Methodik. AVC war an der Anwendung des Tonfraktionsalgorithmus auf die Daten beteiligt und überwachte die Arbeit. PR erhielt die Finanzierung der Studie. Er überwachte die Arbeit, beteiligte sich an der Gestaltung des Arbeitsablaufs und war an der Datenaufbereitung, dem Verfassen und Bearbeiten der Arbeit beteiligt. Alle Autoren haben das Manuskript überprüft.

Korrespondenz mit Alexis Neven.

Die Autoren geben an, dass keine Interessenkonflikte bestehen.

Springer Nature bleibt neutral hinsichtlich der Zuständigkeitsansprüche in veröffentlichten Karten und institutionellen Zugehörigkeiten.

Open Access Dieser Artikel ist unter einer Creative Commons Attribution 4.0 International License lizenziert, die die Nutzung, Weitergabe, Anpassung, Verbreitung und Reproduktion in jedem Medium oder Format erlaubt, sofern Sie den/die ursprünglichen Autor(en) und die Quelle angemessen angeben. Geben Sie einen Link zur Creative Commons-Lizenz an und geben Sie an, ob Änderungen vorgenommen wurden. Die Bilder oder anderes Material Dritter in diesem Artikel sind in der Creative-Commons-Lizenz des Artikels enthalten, sofern in der Quellenangabe für das Material nichts anderes angegeben ist. Wenn Material nicht in der Creative-Commons-Lizenz des Artikels enthalten ist und Ihre beabsichtigte Nutzung nicht durch gesetzliche Vorschriften zulässig ist oder über die zulässige Nutzung hinausgeht, müssen Sie die Genehmigung direkt vom Urheberrechtsinhaber einholen. Um eine Kopie dieser Lizenz anzuzeigen, besuchen Sie http://creativecommons.org/licenses/by/4.0/.

Nachdrucke und Genehmigungen

Neven, A., Christiansen, AV & Renard, P. Automatisches stochastisches 3D-Tonfraktionsmodell aus tTEM-Vermessungs- und Bohrlochdaten. Sci Rep 12, 17112 (2022). https://doi.org/10.1038/s41598-022-21555-z

Zitat herunterladen

Eingegangen: 03. November 2021

Angenommen: 28. September 2022

Veröffentlicht: 12. Oktober 2022

DOI: https://doi.org/10.1038/s41598-022-21555-z

Jeder, mit dem Sie den folgenden Link teilen, kann diesen Inhalt lesen:

Leider ist für diesen Artikel derzeit kein Link zum Teilen verfügbar.

Bereitgestellt von der Content-Sharing-Initiative Springer Nature SharedIt

Durch das Absenden eines Kommentars erklären Sie sich damit einverstanden, unsere Nutzungsbedingungen und Community-Richtlinien einzuhalten. Wenn Sie etwas als missbräuchlich empfinden oder etwas nicht unseren Bedingungen oder Richtlinien entspricht, kennzeichnen Sie es bitte als unangemessen.