Thanks to the blog author for his report
Introduction
This note was written in the wake of the conference and contains a summary of the report with which I spoke. The event took place in June in St. Petersburg. Given that I work a block away from the conference venue, I couldn't help attending this event. In 2016 and 2017, I listened to the reports of the conference, and this year I made a report. Firstly, an interesting (as it seems to me) topic has appeared, which we are developing with , and secondly, after a long study of the legislation of the Russian Federation in terms of the sanctions policy, at the enterprise where I work, there are as many as two licenses .
Before moving on to the topic of my speech, I would like to note the good organization of the event. The visiting page of the conference uses the image of the Kazan Cathedral. The cathedral is one of the main sights of St. Petersburg and it is very clearly visible from the hall where the conference was held.
At the entrance to St. Petersburg State University of Economics, participants were met by assistants from among the students - they did not allow them to get lost. During registration, small souvenirs were given out (a toy - a flashing spike, a pen, stickers with Wolfram symbols). Lunch and coffee breaks were also included in the conference schedule. About delicious coffee and pies, I have already noted on the wall of the group - well done cooks. With this introductory part, I would like to emphasize that the event itself, its format and venue already bring positive emotions.
The report, which was prepared by me and Kirill Belov, is called “Using Wolfram Mathematica to solve problems in applied geophysics. Spectral analysis of seismic data or "where ancient rivers ran". The content of the report covers two parts: firstly, it is the use of algorithms available in for the analysis of geophysical data, and secondly, this is how to place geophysical data in Wolfram Mathematica.
Seismic survey
First you need to make a short digression into geophysics. Geophysics is a science that studies the physical properties of rocks. Well, since the rocks have different properties: electrical, magnetic, elastic, then there are appropriate methods of geophysics: electrical exploration, magnetic exploration, seismic exploration ... In the context of this article, we will only touch on seismic exploration in more detail. Seismic exploration is the main method of searching for oil and gas. The method is based on the excitation of elastic vibrations and the subsequent registration of the response from the rocks that make up the study area. Excitation of vibrations is carried out on land (by dynamite or non-explosive vibration sources of elastic vibrations) or at sea (by air guns). Elastic vibrations propagate through the thickness of rocks, being refracted and reflected at the boundaries of layers with different properties. Reflected waves return to the surface and are recorded by geophones on land (usually electrodynamic devices based on the movement of a magnet suspended in a coil) or hydrophones at sea (based on the piezoelectric effect). By the time of arrival of the waves, one can judge the depths of the geological layers.
Seismic vessel towing equipment
The air gun excites elastic vibrations
Waves pass through the rock mass and are recorded by hydrophones
Research vessel of geophysical exploration "Ivan Gubkin" at the berth near the Blagoveshchensky bridge in St. Petersburg
Seismic signal model
Rocks have different physical properties. For seismic exploration, first of all, elastic properties are important - the velocity of propagation of elastic vibrations and density. If two layers have the same or close properties, then the wave will “not notice” the boundary between them. If the wave velocities in the layers differ, then reflection will occur at the boundary of the layers. The greater the difference in properties, the more intense the reflection. Its intensity will be determined by the reflection coefficient (rc):
where ρ is the rock density, ν is the wave velocity, 1 and 2 denote the upper and lower layers.
One of the simplest and most commonly used seismic signal models is the convolutional model, when the recorded seismic trace is represented as the result of convolution of a sequence of reflection coefficients with a sounding pulse:
where s(t) — seismic trace, i.e. everything recorded by a hydrophone or geophone during a fixed registration time, w(t) - the signal generated by the air gun, n(t) - random noise.
For example, let's calculate a synthetic seismic trace. We will use the Ricker pulse widely used in seismic exploration as the initial signal.
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]Initial seismic impulse
We set two boundaries at depths of 300 ms and 600 ms, and the reflection coefficients will be random numbers
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]]Sequence of reflection coefficients
Calculate and display the seismic trace. Since the reflection coefficients have different signs, we also obtain two sign-alternating reflections on the seismic trace.
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]]simulated track
For this example, you need to make a reservation - in reality, the depth of the layers is determined, of course, in meters, and the calculation of the seismic trace occurs for the time domain. It would be more correct to set the depths in meters and calculate the times of arrival knowing the velocities in the layers. In this case, I immediately set the layers on the time axis.
If we talk about field studies, then as a result of such observations, a huge number of such time series (seismic traces) are recorded. For example, when studying a section 25 km long and 15 km wide, where, as a result of the work, each trace characterizes a cell sized 25x25 meters (such a cell is called a bin), the final data array will contain 600000 traces. With a time sampling step of 1 ms, a recording time of 5 seconds, the final data file will be more than 11 GB, and the volume of the original "raw" material can be hundreds of gigabytes.
How to work with them ?
Plastic bag
The development of the package began on the VK wall of the Russian-speaking support group. Thanks to the community's responses, a solution was found very quickly. And as a result, it turned into a serious development. Corresponding was even marked by the moderators. Currently, the package supports the following types of data that are actively used in the geological industry:
- import of ZMAP and IRAP cartographic data
- import of well measurements in LAS format
- input and output of seismic files of the format
To install the package, you must follow the instructions on the download page of the compiled package, i.e. execute the following code in any :
If[PacletInformation["GeologyIO"] === {}, PacletInstall[URLDownload[
"https://wolfr.am/FiQ5oFih",
FileNameJoin[{CreateDirectory[], "GeologyIO-0.2.2.paclet"}]
]]]After that, the package will be installed to the default folder, the path to which can be obtained as follows:
FileNameJoin[{$UserBasePacletsDirectory, "Repository"}]For example, let's demonstrate the main features of the package. The call is made traditionally for packages in the Wolfram Language:
Get["GeologyIO`"]The package is developed using . This allows you to accompany the main functionality of the package with documentation that does not differ in presentation format from the documentation of Wolfram Mathematica itself and provide the package with test files for a first acquaintance.
Such a file, in particular, is the file "Marmousi.segy" - this is a synthetic model of a geological section, which was developed by the French Institute of Petroleum. Using this model, developers test their own algorithms for wave field modeling, data processing, seismic trace inversion, etc. The Marmousi model itself is stored in the repository where the package was downloaded from. To get the file, run the following code:
If[Not[FileExistsQ["Marmousi.segy"]],
URLDownload["https://wolfr.am/FiQGh7rk", "Marmousi.segy"];]
marmousi = SEGYImport["Marmousi.segy"]Import result - SEGYData object
The SEGY format involves storing various information about observations. First, there are text comments. Information about the place of work, the names of the companies that performed the measurements, etc. are entered here. In our case, this header is called by a request with the TextHeader key. Here is a shortened text title:
Short[marmousi["TextHeader"]]"The Marmousi data set was generated at the Institute …nimum velocity of 1500 m/s and a maximum of 5500 m/s)"
You can display the actual geological model by accessing the seismic traces using the “traces” key (one of the features of the package is case-insensitive keys):
ArrayPlot[Transpose[marmousi["traces"]], PlotTheme -> "Detailed"]Marmousi Model
Currently, the package also allows you to download data in parts from large files, which makes it possible to process files that can reach tens of gigabytes in size. The package also includes functions for exporting data to .segy and partial appending to the end of the file.
Separately, it is worth noting the functionality of the package when working with a complex structure of .segy-files. Since it allows not only accessing individual traces, headers by keys and indexes, but also changing them with subsequent writing to a file. Many of the technical details of the implementation of GeologyIO are beyond the scope of this article and probably deserve a separate description.
The relevance of spectral analysis in seismic exploration
The ability to import seismic data into Wolfram Mathematica allows you to use the built-in signal processing functionality for experimental data. Since each seismic trace is a time series, one of the main tools for studying them is spectral analysis. Among the prerequisites for the analysis of the frequency content of seismic data are, for example, the following:
- Different types of waves are characterized by different frequency composition. This allows you to select useful waves and suppress interference waves.
- Rock properties such as porosity and saturation can affect the frequency content. This allows you to select rocks with better properties.
- Layers with different thickness cause anomalies in different frequency ranges.
The third point is the main one in the context of this article. Below is a code fragment for calculating seismic traces in the case of a layer with varying thickness - a wedge model. This model is traditionally studied in seismic exploration for the analysis of interference effects when waves reflected from many layers are superimposed on each other.
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]]Pinchout Model
The wave speed inside the wedge is 4500 m/s, outside the wedge 4000 m/s, and the density is assumed to be constant 2200 g/cm³. For such a model, we calculate the reflection coefficients and seismic traces.
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]Seismic traces for the wedge model
The sequence of seismic traces depicted in this figure is called a seismic section. As you can see, its interpretation can also be performed on an intuitive level, since the geometry of the reflected waves uniquely corresponds to the model that was specified earlier. If we analyze the traces in more detail, we can see that the traces from the 1st to, approximately, the 30th do not differ - the reflection from the top of the formation and from the sole do not overlap each other. Starting from the 31st trace, the reflections begin to interfere. And, although, in the model, the reflection coefficients do not change horizontally - seismic traces change their intensity when the reservoir thickness changes.
Consider the amplitude of the reflection from the upper boundary of the formation. Starting from the 60th trace, the reflection intensity begins to increase and becomes maximum on the 70th trace. This is how the interference of waves from the top and bottom of the layers manifests itself, leading in some cases to significant anomalies in the seismic record.
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]Graph of the amplitude of the reflected wave from the upper edge of the wedge
It is logical that when the signal is lower-frequency, then the interference begins to appear at large reservoir thicknesses, and in the case of a high-frequency signal, interference occurs at smaller thicknesses. The following code snippet creates a signal at 35Hz, 55Hz, and 85Hz.
waveletSet=Table[(1.0-2.0*(Pi^2)*(f^2)*(t^2))*Exp[-(Pi^2)*(f^2)*(t^2)],
{f,{35,55,85}}];
ListLinePlot[waveletSet,PlotRange->Full,PlotStyle->Black,Frame->True,
PlotLabel->Style["Set of wavelets",Black,20],
LabelStyle->Directive[Black,Italic],
ImageSize->Large,InterpolationOrder->2]Set of initial signals with frequencies 35 Hz, 55 Hz, 85 Hz
After performing the calculation of seismic traces and plotting the amplitudes of the reflected wave, we can see that for different frequencies the anomaly is observed at different reservoir thicknesses.
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]Graphs of the amplitudes of the reflected wave from the upper edge of the wedge for different frequencies
The ability to draw conclusions about the thickness of the reservoir from the results of seismic observations is extremely useful, because one of the main tasks in the exploration of an oil field is to evaluate the most promising points for laying a well (i.e., those areas where the reservoir has a large thickness). In addition, in the geological section there may be such objects that, by their genesis, cause a sharp change in the thickness of the reservoir. This makes spectral analysis an effective tool for studying them. In the next part of the article, we will consider such geological objects in more detail.
Experimental data. Where are they obtained and what to look for in them?
The materials analyzed in the article were obtained on the territory of Western Siberia. The region, as everyone without exception probably knows, is the main oil-producing region of our country. Active development of deposits began in the region in the 60s of the last century. The main method of searching for oil deposits is seismic exploration. It is interesting to consider satellite images of this territory. At a small scale, a huge number of swamps and lakes can be noted, by zooming in on the map, you can see well drilling sites, and by zooming in to the limit, you can also distinguish the clearings of the profiles along which seismic observations were made.
Satellite image of Yandex maps - area of the city of Noyabrsk
Network of well pads at one of the fields
Oil-bearing rocks of Western Siberia occur in a wide range of depths - from 1 km to 5 km. The main volume of rocks containing oil was formed in the Jurassic and Cretaceous. The Jurassic period is probably known to many from the film of the same name. significantly different from today. The Encyclopedia Britannica has a series of paleo-maps that characterize each geological epoch.
Currently,
Jurassic period
Please note that in the Jurassic, the territory of Western Siberia was a sea coast (land, crossed by rivers and a shallow sea). Since the climate was comfortable, it can be assumed that a typical landscape of that time looked like this:
Jurassic Siberia
In this picture, it is not so much animals and birds that are important to us, but the image of the river in the background. The river is the same geological object that we stopped at earlier. The fact is that the activity of the rivers allows the accumulation of well-sorted sandstones, which then become a reservoir for oil. These reservoirs can have a bizarre, complex shape (like a river bed) and they have a variable thickness - the thickness is small near the coast, and increases closer to the center of the channel or in meander areas. So, the rivers formed in the Jurassic are now at a depth of about three kilometers and are the object of a search for oil reservoirs.
Experimental data. Processing and visualization
Let us immediately make a reservation regarding the seismic materials shown in the article - due to the fact that the amount of data used for analysis is significant - only a fragment of the original set of seismic traces is included in the text of the article. This will allow everyone to reproduce the above calculations.
When working with seismic data, a geophysicist usually uses specialized software (there are several industry leaders whose developments are actively used, such as Petrel or Paradigm), which allows you to analyze different types of data and has a convenient graphical interface. Despite all the convenience, such types of software also have their drawbacks - for example, the introduction of modern algorithms into stable versions takes a lot of time, and the possibilities for automating calculations are usually limited. In such a situation, it becomes very convenient to use computer mathematics systems and high-level programming languages that allow the use of a wide algorithmic base and, at the same time, take on a lot of routine. This principle is used to work with seismic data in Wolfram Mathematica. It is not advisable to write rich functionality for interactive work with data - it is more important to provide loading from a generally accepted format, apply the desired algorithms to them and upload them back to an external format.
Following the proposed scheme, we will load the original seismic data and display them in :
Get["GeologyIO`"]
seismic3DZipPath = "seismic3D.zip";
seismic3DSEGYPath = "seismic3D.sgy";
If[FileExistsQ[seismic3DZipPath], DeleteFile[seismic3DZipPath]];
If[FileExistsQ[seismic3DSEGYPath], DeleteFile[seismic3DSEGYPath]];
URLDownload["https://wolfr.am/FiQIuZuH", seismic3DZipPath];
ExtractArchive[seismic3DZipPath];
seismic3DSEGY = SEGYImport[seismic3DSEGYPath]The data loaded and imported in this way is the traces recorded in a section of 10 by 5 kilometers. In the event that the data are obtained using the XNUMXD seismic survey method (waves are recorded not along individual geophysical profiles, but over the entire area simultaneously), it becomes possible to obtain seismic data cubes. These are three-dimensional objects, vertical and horizontal sections of which allow you to study the geological environment in detail. In the considered example, we are dealing with just three-dimensional data. We can get some information from the text header, like this
StringPartition[seismic3DSEGY["textheader"], 80] // TableFormC 1 THIS IS DEMO FILE FOR GEOLOGYIO PACKAGE TEST
C 2
C 3
C 4
C 5 DATE USER NAME: WOLFRAM USER
C6 SURVEY NAME: SOMEWHERE IN SIBERIA
C 7 FILE TYPE 3D SEISMIC VOLUME
C 8
C 9
C10 Z RANGE: FIRST 2200M LAST 2400M
This data set will be enough for us to demonstrate the main stages of data analysis. The traces in the file are recorded sequentially and each of them looks like the following figure - this is the distribution of the amplitudes of the reflected waves along the vertical axis (depth axis).
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]One of the traces of the seismic section
Knowing how many traces are located in each direction of the studied area, you can form a three-dimensional data array and display it using the Image3D[] function
traces=seismic3DSEGY["traces"];
startIL=1050;EndIL=2000;stepIL=2; (*координата Х начала и конца съёмки и шаг трасс*)
startXL=1165;EndXL=1615;stepXL=2; (*координата Y начала и конца съёмки и шаг трасс*)
numIL=(EndIL-startIL)/stepIL+1; (*количество трасс по оис Х*)
numXL=(EndXL-startXL)/stepIL+1; (*количество трасс по оис Y*)
Image3D[ArrayReshape[Abs[traces/Max[Abs[traces[[All,1;;;;4]]]]],{numIL,numXL,101}],ViewPoint->{-1, 0, 0},Background->RGBColor[0,0,0]]XNUMXD image of a seismic data cube. (Vertical axis - depth)
In the event that geological objects of interest create intense seismic anomalies, then visualization tools with transparency can be used. "Unimportant" sections of the recording can be made invisible, leaving only anomalies visible. In Wolfram Mathematica, this can be done with и .
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 of a seismic data cube using Opacity[] and Raster3D[] functions
As in the synthetic example, some geological boundaries (layers) with variable relief can be distinguished on the slices of the original cube.
The main tool for spectral analysis is the Fourier transform. It can be used to estimate the amplitude-frequency spectrum of each trace or group of traces. However, after transferring the data to the frequency domain, information is lost about at what times (read, at what depths) the frequency changes. In order to be able to localize signal changes on the time (depth) axis, the windowed Fourier transform and wavelet decomposition are used. This article uses wavelet decomposition. Wavelet analysis technology began to be actively used in seismic exploration in the 90s. The advantage over the windowed Fourier transform is considered to be better temporal resolution.
Using the following code fragment, one of the seismic traces can be decomposed into individual components:
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]Decomposing a trace into components
To estimate how the reflection energy is distributed at different wave arrival times, scalograms (analogous to a spectrogram) are used. As a rule, in practice it is not necessary to analyze all the components. Usually choose a low-, medium- and high-frequency component.
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]scalogram. Function result
The Wolfram Language uses the wavelet transform function . And the application of this function to the entire set of traces will be carried out using the function . Here it is worth noting one of the strengths of Wolfram Mathematica, the ability to use parallelization . In the above example, there is no need for parallelization - the amount of data is not large, but when working with experimental data sets containing hundreds of thousands of traces, this is a necessity.
tracesCWD=Table[Map[Hilbert[#,0]&,Re[ContinuousWaveletTransform[traces[[i]]][[1]]][[{13,15,18}]]],{i,1,Length@traces}]; After applying the function new data sets appear corresponding to the selected frequencies. In the example above, these are the frequencies: 38Hz, 33Hz, 27Hz. The choice of frequencies is carried out most often on the basis of testing - they receive result maps for different frequency combinations and choose the most informative one from the point of view of a geologist.
If the results need to be shared with colleagues or provided to the customer, then you can use the SEGYExport[] function of the GeologyIO package
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];With three of these cubes available (low, mid, and high), it is common to use RGB blending to visualize the data together. Each component is assigned its own color - red, green, blue. In Wolfram Mathematica, this can be done using the function .
The result is images that can be used for geological interpretation. Meanders, which are fixed on the cut, make it possible to delineate paleochannels, which are more likely to be reservoirs and contain oil reserves. The search and analysis of modern analogues of such a river system makes it possible to determine the most promising parts of the meanders. The channels themselves are characterized by thick layers of well-sorted sandstone and are a good reservoir for oil. Areas outside the "string" anomalies are similar to modern floodplain deposits. Floodplain deposits are mainly represented by clay rocks and drilling into these zones will be inefficient.
RGB slice of the data cube. In the center (slightly to the left of the center) a meandering river can be traced.
RGB slice of the data cube. On the left side you can trace the meandering river.
In some cases, the quality of the seismic data allows for significantly sharper images. It depends on the method of field work, the equipment used by the noise reduction algorithm. In such cases, not only fragments of river systems are visible, but also entire extended paleo-rivers.
RGB blending of the three components of the seismic data cube (horizontal slice). Depth is approximately 2 km.
Satellite image of the Volga River near Saratov
Conclusion
In Wolfram Mathematica, you can analyze seismic data and solve applied problems related to the search for minerals, and the GeologyIO package makes this process more convenient. The structure of seismic data is such that the use of built-in methods for accelerating calculations (, ,…) is very efficient and allows you to process large amounts of data. To a large extent, this is facilitated by the data storage features of the GeologyIO package. By the way, the package can be used not only in the field of applied seismic exploration. Almost the same types of data are used in GPR and seismology. If you have suggestions on how to improve the result, which signal analysis algorithms from the Wolfram Mathematica arsenal are applicable to such data, or if you have critical comments, please leave comments.
Source: habr.com