Merci à l'auteur du blog Anton Ekimenko pour son discoursPrésentation
Cette note a été rédigée à la suite de la
Conférence technologique russe de Wolfram et contient un résumé du rapport que j'ai présenté. L'événement a eu lieu en juin dans la ville de Saint-Pétersbourg. Étant donné que je travaille à un pâté de maisons du lieu de la conférence, je n'ai pas pu m'empêcher d'assister à cet événement. En 2016 et 2017, j'ai écouté les rapports de la conférence, et cette année j'ai fait une présentation. Premièrement, un sujet intéressant (il me semble) est apparu que nous développons avec
Kirill Belov , et deuxièmement, après une longue étude de la législation de la Fédération de Russie concernant la politique de sanctions, deux licences
Wolfram Mathematica sont apparues dans l'entreprise où je travaille.
Avant de passer au sujet de ma présentation, je voudrais noter la bonne organisation de l'événement. La page de la conférence utilise l'image de la cathédrale de Kazan. La cathédrale est l'une des principales attractions de Saint-Pétersbourg et elle est très clairement visible depuis la salle dans laquelle s'est tenue la conférence.
À l'entrée de l'Université d'économie d'État de Saint-Pétersbourg, les participants ont été accueillis par des assistants étudiants - ils ne leur ont pas permis de se perdre. Lors de l'inscription, de petits souvenirs ont été distribués (un jouet - des adhésions clignotantes, un stylo, des autocollants avec des symboles Wolfram). Le déjeuner et les pauses-café étaient également inclus dans le programme de la conférence. À propos de délicieux cafés et tartes, j'ai déjà noté sur le mur du groupe - des chefs bien faits. Avec cette partie introductive, je voudrais souligner que l'événement lui-même, son format et son lieu apportent déjà des émotions positives.
Le rapport, préparé par moi et Kirill Belov, s'intitule «Utilisation de Wolfram Mathematica pour résoudre des problèmes de géophysique appliquée. Analyse spectrale des données sismiques ou "où coulaient les anciennes rivières". Le contenu du rapport comprend deux parties: premièrement, l'utilisation des algorithmes disponibles dans
Wolfram Mathematica pour l'analyse des données géophysiques, et deuxièmement, comment placer les données géophysiques dans Wolfram Mathematica.
Exploration sismique
Vous devez d'abord faire une petite excursion en géophysique. La géophysique est une science qui étudie les propriétés physiques des roches. Eh bien, puisque les roches ont des propriétés différentes: électriques, magnétiques, élastiques, il existe des méthodes de géophysique correspondantes: l'exploration électrique, l'exploration magnétique, l'exploration sismique ... Dans le cadre de cet article, nous n'aborderons que l'exploration sismique plus en détail. L'exploration sismique est la principale méthode de recherche de pétrole et de gaz. La méthode est basée sur l'excitation des vibrations élastiques et l'enregistrement subséquent de la réponse des roches composant la zone d'étude. L'excitation des vibrations s'effectue sur terre (dynamite ou sources de vibrations non explosives de vibrations élastiques) ou en mer (canons à air). Les vibrations élastiques se propagent à travers l'épaisseur des roches, étant réfractées et réfléchies aux limites de couches aux propriétés différentes. Les ondes réfléchies reviennent à la surface et sont enregistrées par des géophones terrestres (généralement des appareils électrodynamiques basés sur le mouvement d'un aimant suspendu dans une bobine) ou des hydrophones dans la mer (basés sur l'effet piézoélectrique). Au moment de l'arrivée des vagues, on peut juger de la profondeur des couches géologiques.
Un navire sismique remorque du matériel:
Le pistolet à air excite les vibrations élastiques:
Les vagues traversent la masse rocheuse et sont enregistrées par des hydrophones:
Navire de recherche en exploration géophysique «Ivan Gubkin» sur le quai du pont Blagoveshchensky à Saint-Pétersbourg:
Modèle de signal sismique
Les roches ont des propriétés physiques différentes. Pour l'exploration sismique, les propriétés élastiques sont les plus importantes - la vitesse de propagation des vibrations élastiques et la densité. Si deux couches ont des propriétés identiques ou similaires, l'onde ne «remarquera» pas la frontière entre elles. Si les vitesses des vagues dans les couches diffèrent, une réflexion se produira à la limite des couches. Plus la différence de propriétés est grande, plus la réflexion est intense. Son intensité sera déterminée par le coefficient de réflexion (rc):
où ρ est la densité de la roche, ν est la vitesse des vagues, 1 et 2 désignent les couches supérieure et inférieure.L'un des modèles de signaux sismiques les plus simples et les plus couramment utilisés est le modèle de convolution, lorsqu'une trace sismique enregistrée est présentée à la suite de la convolution d'une séquence de coefficients de réflexion avec une impulsion de sonde:
où s (t) est la piste sismique, c'est-à-dire tout ce qu'un hydrophone ou un géophone enregistré sur une durée d'enregistrement fixe,
w (t) est le signal généré par le pistolet à air,
n (t) est un bruit aléatoire.
Nous calculons par exemple un levé sismique synthétique. Nous utiliserons l'impulsion Ricker largement utilisée en exploration sismique comme signal initial.
length=0.050; (*Signal lenght*) dt=0.001;(*Sample rate of signal*) t=Range[-length/2,(length)/2,dt];(*Signal time*) f=35;(*Central frequency*) wavelet=(1.0-2.0*(Pi^2)*(f^2)*(t^2))*Exp[-(Pi^2)*(f^2)*(t^2)]; ListLinePlot[wavelet, Frame->True,PlotRange->Full,Filling->Axis,PlotStyle->Black, PlotLabel->Style["Initial wavelet",Black,20], LabelStyle->Directive[Black,Italic], FillingStyle->{White,Black},ImageSize->Large,InterpolationOrder->2]
Impulsion sismique source
Nous avons fixé deux limites à des profondeurs de 300 ms et 600 ms, et les coefficients de réflexion seront des nombres aléatoires.
rcExample=ConstantArray[0,1000]; rcExample[[300]]=RandomReal[{-1,0}]; rcExample[[600]]=RandomReal[{0,1}]; ListPlot[rcExample,Filling->0,Frame->True,Axes->False,PlotStyle->Black, PlotLabel->Style["Reflection Coefficients",Black,20], LabelStyle->Directive[Black,Italic]]
La séquence des coefficients de réflexion:
Calculez et affichez la piste sismique. Les coefficients de réflexion ayant des signes différents, nous obtenons deux réflexions alternées sur le trajet sismique.
traceExamle=ListConvolve[wavelet[[1;;;;1]],rcExample]; ListPlot[traceExamle, PlotStyle->Black,Filling->0,Frame->True,Axes->False, PlotLabel->Style["Seismic trace",Black,20], LabelStyle->Directive[Black,Italic]]
Piste simulée:
Pour cet exemple, une réservation doit être faite - en réalité, la profondeur des couches est bien sûr déterminée en mètres, et le calcul du chemin sismique se fait pour le domaine temporel. Il serait plus correct de régler les profondeurs en mètres et de calculer les temps d'arrivée en connaissant la vitesse dans les strates. Dans ce cas, j'ai immédiatement mis les calques sur l'axe du temps.Si nous parlons d'études sur le terrain, à la suite de ces observations, un grand nombre de ces séries chronologiques (traces sismiques) sont enregistrées. Par exemple, lors de l'exploration d'un site de 25 km de long et 15 km de large, où, à la suite du travail, chaque piste caractérise une cellule de 25 x 25 mètres (une telle cellule est appelée bac), le tableau de données final contiendra 600000 pistes. Avec une étape d'échantillonnage temporel de 1 ms, une durée d'enregistrement de 5 secondes, le fichier de données final sera supérieur à 11 Go et la quantité de matière «brute» d'origine peut atteindre des centaines de gigaoctets.
Comment travailler avec eux dans
Wolfram Mathematica ?
Le développement du paquet a commencé par une
question sur le mur VK du groupe de soutien russophone. Grâce aux réponses de la communauté, une solution a été trouvée très rapidement. Et en conséquence, il est devenu un développement sérieux. Le
message correspondant
sur le mur de Wolfram Comunity a même
été noté par les modérateurs. Actuellement, le package prend en charge l'utilisation des types de données suivants qui sont activement utilisés dans l'industrie géologique:
- importer des données cartographiques au format ZMAP et IRAP
- importation de mesures dans des puits au format LAS
- entrée et sortie de fichiers sismiques au format SEGY
Pour installer le package, vous devez suivre les instructions sur la page de téléchargement du package assemblé, c'est-à-dire exécutez le code suivant sur n'importe quel
bloc -
notes Mathematica :
If[PacletInformation["GeologyIO"] === {}, PacletInstall[URLDownload[ "https://wolfr.am/FiQ5oFih", FileNameJoin[{CreateDirectory[], "GeologyIO-0.2.2.paclet"}] ]]]
Après cela, le package sera installé dans le dossier par défaut, dont le chemin peut être obtenu comme suit:
FileNameJoin[{$UserBasePacletsDirectory, "Repository"}]
À titre d'exemple, nous montrons les principales fonctionnalités du package. L'appel est traditionnellement pour les packages en Wolfram Language:
Get["GeologyIO`"]
Le package est développé à l'aide de
Wolfram Workbench . Cela vous permet d'accompagner la fonctionnalité principale du package de documentation, qui ne diffère pas dans le format de présentation de la documentation de Wolfram Mathematica lui-même et de fournir au package des fichiers de test pour la première connaissance.
Un tel fichier, en particulier, est le fichier Marmousi.segy, un modèle synthétique d'une section géologique développé par l'Institut Français du Pétrole. En utilisant ce modèle, les développeurs testent leurs propres algorithmes pour modéliser le champ d'onde, le traitement des données, l'inversion des traces sismiques, etc. Le modèle Marmousi lui-même est stocké dans le référentiel, d'où le package lui-même a été téléchargé. Afin d'obtenir le fichier, exécutez le code suivant:
If[Not[FileExistsQ["Marmousi.segy"]], URLDownload["https://wolfr.am/FiQGh7rk", "Marmousi.segy"];] marmousi = SEGYImport["Marmousi.segy"]
Le résultat de l'importation est un objet SEGYData:
Le format SEGY implique le stockage de diverses informations sur les observations. Tout d'abord, ce sont des commentaires de texte. Les informations sur le lieu de travail, les noms des entreprises qui ont effectué les mesures, etc. sont saisies ici. Dans notre cas, cet en-tête est appelé par une requête avec la clé TextHeader. Voici un titre de texte abrégé:
Short[marmousi["TextHeader"]]
«L'ensemble de données Marmousi a été généré à l'Institut ... vitesse minimale de 1500 m / s et maximum de 5500 m / s)»
Le modèle géologique réel peut être affiché en contactant des traces sismiques à l'aide de la clé «traces» (l'une des caractéristiques du package est la non-casse des clés):
ArrayPlot[Transpose[marmousi["traces"]], PlotTheme -> "Detailed"]
Modèle Marmousi:
À l'heure actuelle, le package vous permet également de télécharger des données en plusieurs parties à partir de gros fichiers, de sorte qu'il devient possible de traiter des fichiers dont la taille peut atteindre des dizaines de gigaoctets. Le package comprend également des fonctions d'exportation des données vers .segy et une réécriture partielle à la fin du fichier.
Séparément, il convient de noter la fonctionnalité du package lorsque vous travaillez avec la structure complexe des fichiers .segy. Puisqu'il permet non seulement d'accéder aux clés et aux index des traces individuelles, des en-têtes, mais aussi de les modifier avec l'écriture ultérieure dans un fichier. De nombreux détails techniques de la mise en œuvre de GeologyIO dépassent le cadre de cet article et méritent probablement une description distincte.
La pertinence de l'analyse spectrale dans l'exploration sismique
La possibilité d'importer des matériaux sismiques dans Wolfram Mathematica vous permet d'utiliser la fonctionnalité de traitement du signal intégrée pour les données expérimentales. Chaque piste sismique étant une série temporelle, l'un des principaux outils pour les étudier est l'analyse spectrale. Par exemple, les conditions préalables à l'analyse de la composition en fréquence des données sismiques sont les suivantes:
- Différents types d'ondes sont caractérisés par une composition de fréquences différente. Cela nous permet de sélectionner les ondes utiles et de supprimer les ondes d'interférence.
- Les propriétés des roches telles que la porosité et la saturation peuvent affecter la composition en fréquence. Cela vous permet de sélectionner des roches avec de meilleures propriétés.
- Des couches d'épaisseurs différentes provoquent des anomalies dans différentes gammes de fréquences.
Le troisième paragraphe est central dans le contexte de cet article. Ci-dessous, un extrait de code pour calculer les traces sismiques dans le cas d'une couche d'épaisseur variable - un modèle en coin. Ce modèle est traditionnellement étudié en exploration sismique pour l'analyse des effets d'interférence, lorsque les ondes réfléchies par de nombreuses couches se chevauchent.
nx=200;(* Number of grid points in X direction*) ny=200;(* Number of grid points in Y direction*) T=2;(*Total propagation time*) (*Velocity and density*) modellv=Table[4000,{i,1,ny},{j,1,nx}];(* P-wave velocity in m/s*) rho=Table[2200,{i,1,ny},{j,1,nx}];(* Density in g/cm^3, used constant density*) Table[modellv[[150-Round[i*0.5];;,i]]=4500;,{i,1,200}]; Table[modellv[[;;70,i]]=4500;,{i,1,200}]; (*Plotting model*) MatrixPlot[modellv,PlotLabel->Style["Model of layer",Black,20], LabelStyle->Directive[Black,Italic]]
Modèle de formation des explosions:
La vitesse des vagues à l'intérieur du coin est de 4500 m / s, à l'extérieur du coin de 4000 m / s, et la densité est supposée constante à 2200 g / cm³. Pour un tel modèle, nous calculons les coefficients de réflexion et les traces sismiques.
rc=Table[N[(modellv[[All,i]]-PadLeft[modellv[[All,i]],201,4000][[1;;200]])/(modellv[[All,i]]+PadLeft[modellv[[All,i]],201,4500][[1;;200]])],{i,1,200}]; traces=Table[ListConvolve[wavelet[[1;;;;1]],rc[[i]]],{i,1,200}]; starttrace=10; endtrace=200; steptrace=10; trasenum=Range[starttrace,endtrace,steptrace]; traserenum=Range[Length@trasenum]; tracedist=0.5; Rotate[Show[ Reverse[Table[ ListLinePlot[traces[[trasenum[[i]]]]*50+trasenum[[i]]*tracedist,Filling->{1->{trasenum[[i]]*tracedist,{RGBColor[0.97,0.93,0.68],Black}}},PlotStyle->Directive[Gray,Thin],PlotRange->Full,InterpolationOrder->2,Axes->False,Background->RGBColor[0.97,0.93,0.68]], {i,1,Length@trasenum}]],ListLinePlot[Transpose[{ConstantArray[45,80],Range[80]}],PlotStyle->Red],PlotRange->All,Frame->True],270Degree]
Pistes sismiques pour le modèle à coin:
La séquence de traces sismiques représentée sur cette figure est appelée une section sismique. Comme vous pouvez le voir, son interprétation peut être effectuée à un niveau intuitif, car la géométrie des ondes réfléchies correspond uniquement au modèle défini précédemment. Si vous analysez les pistes plus en détail, vous pouvez voir que les pistes du 1er au 30 environ ne diffèrent pas - la réflexion du haut de la formation et de la semelle ne se chevauchent pas. À partir de la route 31, les réflexions commencent à interférer. Et, bien que, dans le modèle, les coefficients de réflexion ne changent pas horizontalement - les pistes sismiques changent leur intensité avec un changement dans l'épaisseur du réservoir.
Considérez l'amplitude de réflexion de la limite supérieure du réservoir. À partir du 60e itinéraire, l'intensité de la réflexion commence à augmenter et au 70e itinéraire devient maximum. C'est ainsi que se manifeste l'interférence des ondes du toit et de la semelle pour les formations, conduisant dans certains cas à des anomalies significatives d'enregistrement sismique.
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]
Graphique de l'amplitude de l'onde réfléchie depuis le bord supérieur du coin
Il est logique que lorsque le signal est à basse fréquence, alors des interférences commencent à apparaître à de grandes épaisseurs de la formation, et dans le cas d'un signal à haute fréquence, des interférences se produisent à des épaisseurs plus petites. Le fragment de code suivant crée un signal avec des fréquences de 35 Hz, 55 Hz et 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]
Un ensemble de signaux source avec des fréquences de 35 Hz, 55 Hz, 85 Hz
Après avoir calculé les traces sismiques et tracé les amplitudes de l'onde réfléchie, nous pouvons voir que pour différentes fréquences une anomalie est observée à différentes épaisseurs de la formation.
tracesSet=Table[ListConvolve[waveletSet[[j]][[1;;;;1]],rc[[i]]],{j,1,3},{i,1,200}]; lowFreq=ListLinePlot[GaussianFilter[Abs[tracesSet[[1]][[All,46]]],3][[;;;;2]],InterpolationOrder->2,PlotStyle->Black,PlotRange->All]; medFreq=ListLinePlot[GaussianFilter[Abs[tracesSet[[2]][[All,46]]],3][[;;;;2]],InterpolationOrder->2,PlotStyle->Black,PlotRange->All]; highFreq=ListLinePlot[GaussianFilter[Abs[tracesSet[[3]][[All,46]]],3][[;;;;2]],InterpolationOrder->2,PlotStyle->Black,PlotRange->All]; Show[lowFreq,medFreq,highFreq,PlotRange->{{0,100},All}, PlotLabel->Style["Amplitudes of reflection",Black,20], LabelStyle->Directive[Black,Italic], Frame->True]
Graphiques des amplitudes de l'onde réfléchie depuis le bord supérieur du coin pour différentes fréquences
La capacité de tirer des conclusions sur l'épaisseur du réservoir sur la base des résultats des observations sismiques est extrêmement utile, car l'une des principales tâches de la prospection d'un gisement de pétrole est d'évaluer les points les plus prometteurs pour la pose d'un puits (c'est-à-dire les sections où le réservoir a une grande épaisseur). De plus, dans la section géologique, il peut y avoir de tels objets qui, en raison de leur genèse, provoquent un changement brusque de l'épaisseur de la formation. Cela fait de l'analyse spectrale un outil efficace pour les étudier. Dans la prochaine partie de l'article, nous examinerons ces objets géologiques plus en détail.
Données expérimentales. Où sont-ils reçus et que rechercher en eux?
Les matériaux analysés dans l'article ont été obtenus sur le territoire de la Sibérie occidentale. La région, comme chacun le sait sans exception, est probablement la principale région productrice de pétrole de notre pays. Le développement actif des naissances en métro a commencé dans la région dans les années 60 du siècle dernier. La principale méthode de recherche de gisements de pétrole est l'exploration sismique. Il est intéressant de considérer les images satellites de ce territoire. À petite échelle, un grand nombre de marécages et de lacs peuvent être notés, en augmentant la carte, vous pouvez voir des puits de grappe pour le forage de puits, et en augmentant la carte à la limite, vous pouvez également distinguer des clairières de profils le long desquels des observations sismiques ont été effectuées.
Image satellite des cartes Yandex - la région de la ville de Noyabrsk:
Un réseau de sites de cluster sur l'un des domaines:
Les roches pétrolifères de la Sibérie occidentale se trouvent dans une large gamme de profondeurs - de 1 km à 5 km. Le volume principal des roches contenant du pétrole se forme au Jurassique et au Crétacé. La période jurassique est probablement connue de beaucoup par le film du même nom.
Le climat de la période jurassique était significativement différent de celui moderne. Dans l'Encyclopedia Britannica, il existe une série de paléocartes qui caractérisent chaque époque gélogique.
Heure actuelle:
Période Jurassique:
Veuillez noter que dans le Jurassique, le territoire de la Sibérie occidentale était une côte maritime (terre traversée par des rivières et une mer peu profonde). Comme le climat était confortable, on peut supposer qu'un paysage typique de l'époque ressemblait à ceci:
Sibérie du Jurassique:
Sur cette photo, pas tant les animaux et les oiseaux sont importants pour nous, que l'image de la rivière en arrière-plan. Une rivière est l'objet géologique même sur lequel nous nous sommes arrêtés plus tôt. Le fait est que l'activité des rivières permet l'accumulation de grès bien triés, qui deviennent alors un réservoir de pétrole. Ces réservoirs peuvent avoir une forme bizarre et complexe (comme un lit de rivière) et ils ont une épaisseur variable - le long de la côte, l'épaisseur est petite et augmente plus près du centre du canal ou dans des sections des méandres. Ainsi, les rivières formées au Jurassique sont maintenant à une profondeur d'environ trois kilomètres et font l'objet de la recherche de gisements de pétrole.
Données expérimentales. Traitement et visualisation
Nous émettons immédiatement une réserve concernant les matériaux sismiques indiqués dans l'article - en raison du fait que la quantité de données utilisées pour l'analyse est importante - seul un fragment de l'ensemble original de traces sismiques est placé dans le texte de l'article. Cela permettra à chacun de reproduire les calculs ci-dessus.
Lorsqu'ils travaillent avec des données sismiques, les géophysiciens utilisent généralement des logiciels spécialisés (il existe plusieurs leaders de l'industrie dont le développement est activement utilisé, par exemple Petrel ou Paradigm), ce qui vous permet d'analyser différents types de données et dispose d'une interface graphique pratique. Malgré toute la commodité, ces types de logiciels ont également leurs inconvénients - par exemple, la mise en œuvre d'algorithmes modernes dans des versions stables prend beaucoup de temps et les possibilités d'automatisation des calculs sont généralement limitées. Dans une telle situation, il devient très pratique d'utiliser des systèmes de mathématiques informatiques et des langages de programmation de haut niveau qui vous permettent d'utiliser une large base algorithmique et, en même temps, de prendre beaucoup de routine. En utilisant ce principe, le travail avec les données sismiques dans Wolfram Mathematica est construit. Il n'est pas pratique d'écrire une riche fonctionnalité de travail interactif avec des données - il est plus important d'assurer le chargement à partir d'un format généralement accepté, de leur appliquer les algorithmes souhaités et de les télécharger à nouveau dans un format externe.
En suivant le schéma proposé, chargez les données sismiques d'origine et affichez-les dans
Wolfram Mathematica :
Get["GeologyIO`"] seismic3DZipPath = "seismic3D.zip"; seismic3DSEGYPath = "seismic3D.sgy"; If[FileExistsQ[seismic3DZipPath], DeleteFile[seismic3DZipPath]]; If[FileExistsQ[seismic3DSEGYPath], DeleteFile[seismic3DSEGYPath]]; URLDownload["https://wolfr.am/FiQIuZuH", seismic3DZipPath]; ExtractArchive[seismic3DZipPath]; seismic3DSEGY = SEGYImport[seismic3DSEGYPath]
Les données téléchargées et importées de cette manière sont les traces enregistrées sur un site de 10 kilomètres sur 5. Dans le cas où les données sont obtenues par la méthode d'exploration sismique tridimensionnelle (l'enregistrement des vagues n'est pas effectué le long de profils géophysiques séparés, mais sur toute la zone en même temps), il devient possible d'obtenir des cubes de données sismiques. Ce sont des objets tridimensionnels dont les coupes verticales et horizontales permettent une étude détaillée de l'environnement géologique. Dans l'exemple considéré, nous ne traitons que de données tridimensionnelles. Nous pouvons obtenir des informations à partir du titre du texte, par exemple comme ceci
StringPartition[seismic3DSEGY["textheader"], 80] // TableForm
C 1 CECI EST UN FICHIER DEMO POUR LE TEST DU PAQUET GEOLOGYIO
C 2
C 3
C 4
C 5 DATE NOM D'UTILISATEUR: UTILISATEUR WOLFRAM
C 6 NOM DE L'ENQUÊTE: QUELQUE PART EN SIBÉRIA
C 7 TYPE DE FICHIER 3D VOLUME SISMIQUE
C 8
C 9
GAMME C10 Z: PREMIER 2200M DERNIER 2400M
Cet ensemble de données sera suffisant pour nous montrer les principales étapes de l'analyse des données. Les pistes du fichier sont enregistrées séquentiellement et chacune d'elles ressemble à la figure suivante - il s'agit de la distribution des amplitudes des ondes réfléchies le long de l'axe vertical (axe de profondeur).
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]
L'une des sections sismiques:
Sachant combien de traces se trouvent dans chaque direction de la zone étudiée, vous pouvez former un tableau de données en trois dimensions et l'afficher à l'aide de la fonction Image3D []
traces=seismic3DSEGY["traces"]; startIL=1050;EndIL=2000;stepIL=2; (* *) startXL=1165;EndXL=1615;stepXL=2; (* Y *) numIL=(EndIL-startIL)/stepIL+1; (* *) numXL=(EndXL-startXL)/stepIL+1; (* Y*) Image3D[ArrayReshape[Abs[traces/Max[Abs[traces[[All,1;;;;4]]]]],{numIL,numXL,101}],ViewPoint->{-1, 0, 0},Background->RGBColor[0,0,0]]
Image tridimensionnelle d'un cube de données sismiques (axe vertical - profondeur):
Dans le cas où des objets géologiques d'intérêt créent des anomalies sismiques intenses, des outils de visualisation avec transparence peuvent être utilisés. Les sections «sans importance» de l'enregistrement peuvent être rendues invisibles, ne laissant apparaître que des anomalies. Dans Wolfram Mathematica, cela peut être fait en utilisant
Opacity [] et
Raster3D [] .
data = ArrayReshape[Abs[traces/Max[Abs[traces[[All,1;;;;4]]]]],{numIL,numXL,101}]; Graphics3D[{Opacity[0.1], Raster3D[data, ColorFunction->"RainbowOpacity"]}, Boxed->False, SphericalRegion->True, ImageSize->840, Background->None]
Image d'un cube de données sismiques utilisant les fonctions Opacity [] et Raster3D []
Comme dans l'exemple synthétique, sur les tranches du cube d'origine, certaines limites géologiques (couches) à relief variable peuvent être distinguées.
Le principal outil d'analyse spectrale est la transformée de Fourier. En l'utilisant, vous pouvez évaluer le spectre amplitude-fréquence de chaque chemin ou groupe de chemins. Cependant, après avoir transféré des données dans le domaine fréquentiel, des informations sont perdues sur les moments (lus à quelles profondeurs) la fréquence change. Afin de pouvoir localiser les changements de signal sur l'axe du temps (profond), la transformée de Fourier de fenêtre et la décomposition en ondelettes sont utilisées. Cet article utilise la décomposition en ondelettes. La technologie d'analyse des ondelettes a commencé à être activement utilisée dans l'exploration sismique dans les années 90. L'avantage par rapport à la transformée de Fourier de la fenêtre est considéré comme la meilleure résolution temporelle.
À l'aide du fragment de code suivant, vous pouvez vous décomposer en composants individuels de l'une des traces sismiques:
cwd=ContinuousWaveletTransform[seismicSection["traces"][[100]]] Show[ ListLinePlot[Re[cwd[[1]]],PlotRange->All], ListLinePlot[seismicSection["traces"][[100]], PlotStyle->Black,PlotRange->All],ImageSize->{1500,500},AspectRatio->Full, PlotLabel->Style["Wavelet decomposition",Black,32], LabelStyle->Directive[Black,Italic], PlotRange->All, Frame->True]
Décomposition des composants
Pour évaluer la répartition de l'énergie de réflexion à différents moments d'arrivée des ondes, des scalogrammes (analogiques du spectrogramme) sont utilisés. En règle générale, dans la pratique, il n'est pas nécessaire d'analyser tous les composants. Choisissez généralement la composante basse, moyenne et haute fréquence.
freq=(500/(#*contWD["Wavelet"]["FourierFactor"]))&/@(Thread[{Range[contWD["Octaves"]],1}]/.contWD["Scales"])//Round; ticks=Transpose[{Range[Length[freq]],freq}]; WaveletScalogram[contWD,Frame->True,FrameTicks->{{ticks,Automatic},Automatic},FrameTicksStyle->Directive[Orange,12], FrameLabel->{"Time","Frequency(Hz)"},LabelStyle->Directive[Black,Bold,14], ColorFunction->"RustTones",ImageSize->Large]
Scalogramme. Résultat de la fonction
WaveletScalogram []
Wolfram Language utilise la fonction ContinuousWaveletTransform [] pour les transformées en ondelettes . Et l'application de cette fonction à l'ensemble des traces s'effectue à l'aide de la fonction Table [] . Ici, il convient de noter que l'une des forces de Wolfram Mathematica est la possibilité d'utiliser la parallélisation ParallelTable [] . Dans l'exemple donné, la parallélisation n'est pas nécessaire - le volume de données n'est pas important, mais lorsque vous travaillez avec des ensembles de données expérimentales contenant des centaines de milliers de traces, c'est une nécessité. tracesCWD=Table[Map[Hilbert[#,0]&,Re[ContinuousWaveletTransform[traces[[i]]][[1]]][[{13,15,18}]]],{i,1,Length@traces}];
Après avoir appliqué la fonction ContinuousWaveletTransform [] , de nouveaux tableaux de données correspondant aux fréquences sélectionnées apparaissent. Dans l'exemple ci-dessus, ce sont des fréquences: 38 Hz, 33 Hz, 27 Hz. Le choix des fréquences s'effectue le plus souvent sur la base de tests - ils obtiennent des cartes efficaces pour différentes combinaisons de fréquences et choisissent la plus informative du point de vue d'un géologue.Si vous avez besoin de partager les résultats avec des collègues ou de les fournir au client, vous pouvez utiliser la fonction SEGYExport [] de GeologyIO outputdata=seismic3DSEGY; outputdata["traces",1;;-1]=tracesCWD[[All,3]]; outputdata["textheader"]="Wavelet Decomposition Result"; outputdata["binaryheader","NumberDataTraces"]=Length[tracesCWD[[All,3]]]; SEGYExport["D:\\result.segy",outputdata];
Disposant de trois cubes de ce type (composants basse fréquence, moyenne fréquence et haute fréquence), le mixage RVB est généralement utilisé pour la visualisation conjointe des données. Chaque composant a sa propre couleur - rouge, vert, bleu. Dans Wolfram Mathematica, cela peut être fait en utilisant la fonction ColorCombine [] .En conséquence, des images sont obtenues qui peuvent être utilisées pour effectuer une interprétation géologique. Les méandres qui pendent sur la section vous permettent de tracer le paléore, qui est plus susceptible d'être des réservoirs et de contenir des réserves de pétrole. La recherche et l'analyse d'analogues modernes d'un tel système fluvial nous permettent de déterminer les parties les plus prometteuses des méandres. Les lits eux-mêmes sont caractérisés par d'épaisses couches de grès bien triés et constituent un bon réservoir de pétrole. Les zones en dehors des anomalies «dentelles» sont similaires aux dépôts modernes des plaines inondables. Les sédiments des plaines inondables sont principalement représentés par des roches argileuses et le forage dans ces zones sera inefficace.Tranche RVB du cube de données. Au centre (un peu à gauche du centre), vous pouvez tracer la rivière sinueuse.
Tranche RVB du cube de données. Sur le côté gauche, vous pouvez tracer la rivière sinueuse.
Dans certains cas, la qualité des données sismiques permet des images nettement plus nettes. Cela dépend de la méthodologie du travail sur le terrain, des équipements adoptés par l'algorithme de réduction du bruit. Dans de tels cas, non seulement des fragments du système fluvial sont visibles, mais aussi des paléorecs étendus entiers.Mélange RVB des trois composants du cube de données sismiques (tranche horizontale). La profondeur est d'environ 2 km.
Image satellite de la Volga dans la région de Saratov
Conclusion
Dans Wolfram Mathematica, vous pouvez analyser les données sismiques et résoudre les problèmes d'application associés à la recherche de minéraux, et le package GeologyIO vous permet de rendre ce processus plus pratique. La structure des données sismiques est telle que l'utilisation de méthodes intégrées pour accélérer les calculs ( ParallelTable [] , ParallelDo [], ...) est très efficace et vous permet de traiter de grandes quantités de données. Dans une large mesure, cela est facilité par les fonctionnalités de stockage du package GeologyIO. Soit dit en passant, le paquet peut être utilisé non seulement dans le domaine de l'exploration sismique appliquée. Presque les mêmes types de données sont utilisés en géoradar et en sismologie. Si vous avez des suggestions sur la façon d'améliorer le résultat, quels algorithmes d'analyse de signal Wolfram Mathematica sont applicables à ces données, ou si vous avez des commentaires critiques, laissez des commentaires.