Visualising the Subsurface utilizing Python Information Visualisation Libraries
Inside geoscience, it’s important to have a full understanding of the geological surfaces current throughout the subsurface. By realizing the precise location of the formation and its geometry, it turns into simpler to establish potential new targets for oil and gasoline exploration in addition to potential carbon seize and storage areas. There are a number of strategies we will use to refine these surfaces starting from seismic information to nicely log derived formation tops. Most frequently, these methods are utilized in mixture with one another to refine the ultimate floor.
This text continues from my earlier articles, which targeted on the method of extrapolating nicely log measurements throughout a area to know and visualise geospatial variations. On this article, we’ll take a look at how we will create 3D surfaces utilizing interactive Plotly charts.
As modelling geological surfaces is a fancy course of and sometimes includes a number of iterations and refinement, this text demonstrates a quite simple instance of how we will visualise this information with Python.
To see how we will use Plotly to visualise our geological formation high throughout an space we will probably be utilizing two units of information.
The primary information set is derived from 28 intersections with the formation derived from wellbore picks, that are used as enter for kriging to supply a low-resolution floor. In distinction, the second set of information is derived from interpreted seismic information, offering a a lot greater decision floor.
Each units of information are from the Equinor Volve dataset. Particulars of that are on the backside of this text.
You may also try the next articles inside this mini-series on the hyperlinks under:
Earlier than making an attempt to do any work with the information, we first must import the libraries we’ll want. These are:
- pandas — to learn our information, which is in
csvformat - matplotlib to create our visualisation
- pykrige to hold out the kriging
- numpy for some numerical calculations
- plotly.graph_objects to visualise floor in 3D
import pandas as pd
import matplotlib.pyplot as plt
import plotly.graph_objects as gofrom pykrige import OrdinaryKriging
import numpy as np
Subsequent, we will load the information utilizing pd.read_csv().
As this information incorporates details about geological surfaces from all the wells throughout the Volve discipline, we will use question() to extract the information for the formation we want. On this case, we will probably be wanting on the Hugin Formation.
df = pd.read_csv('Information/Volve/Well_picks_Volve_v1 copy.csv')
df_hugin = df.question('SURFACE == "Hugin Fm. VOLVE High"')
df_hugin
Once we run the above code, we get again the next desk. It’s possible you’ll discover that some wells have encountered the Hugin Formation a number of occasions. That is seemingly because of the wells penetrating the formation a number of occasions both resulting from wellbore or formation geometry.
In my earlier articles, I’ve targeted on how we will use a course of referred to as kriging to “fill within the gaps” between the measurement factors. We gained’t be masking the main points of this course of on this article; nevertheless, you’ll be able to test this article for more information.
As soon as our information has been loaded, we will run the kriging course of by calling upon pykrige’s OrdinaryKriging methodology.
Inside this name, we cross in our x and y information which represents the Easting and Northing place of the place the wellbore encountered the formation throughout the subsurface.
As we wish to generate a floor of the Hugin Formation, we have to use the TVDSS — True Vertical Depth Subsea — measurement. This offers a real reflection of how deep down the floor is under the chosen datum.
OK = OrdinaryKriging(x=df_hugin['Easting'],
y=df_hugin['Northing'],
z=df_hugin['TVDSS'],
variogram_model='linear',
verbose=True, enable_plotting=True)
As soon as the mannequin has been generated, we will apply it to 2 arrays that cowl the complete vary of the nicely/penetration factors.
This enables us to generate a grid of values that are then handed into the OrdinaryKriging object we generated above.
gridx = np.arange(433986, 438607, 50, dtype='float64')
gridy = np.arange(6477539, 6479393, 50,dtype='float64')zstar, ss = OK.execute('grid', gridx, gridy)
To complete issues off, we will generate a easy 2D map view of our floor utilizing matplotlib’s imshow methodology.
fig, ax = plt.subplots(figsize=(15,5))# Create a 2D picture plot of the information in 'zstar'
# The 'extent' parameter units the bounds of the picture in information coordinates
# 'origin' parameter units the a part of the picture that corresponds to the origin of the axes
picture = ax.imshow(zstar, extent=(433986, 438607, 6477539, 6479393),
origin='decrease')
# Set the labels for the x-axis and y-axis
ax.set_xlabel('X Location (m)', fontsize=14, fontweight='daring')
ax.set_ylabel('Y Location (m)', fontsize=14, fontweight='daring')
# Add contours
contours = ax.contour(gridx, gridy, zstar, colours='black')
colorbar = fig.colorbar(picture)
colorbar.set_label('DTC (us/ft)', fontsize=14, fontweight='daring')
# Show the plot
plt.present()
To transform our 2D floor to 3D, we have to use Plotly. We might use matplotlib to do that; nevertheless, from my expertise, it’s simpler, smoother and extra interactive to generate the 3D visualisations with Plotly.
Within the code under, we first create our grid of coordinates. To do that, we use numpy’s linspace operate. This operate will create a set of evenly spaced numbers over a specified vary. For our dataset and instance, the vary extends from the minimal to the utmost of xgrid_extent and ygrid_extent.
The full variety of values used inside this vary will probably be equal to the variety of x and y parts current within the zstar grid we created above.
As soon as our grid is shaped, we then name upon the Plotly.
First, we create our determine object after which use fig.add_trace so as to add our 3D floor plot to the determine.
As soon as this has been added, we have to tweak the structure of our plot in order that we’ve got axis labels, an outlined width and top, and a few padding.
xgrid_extent = [433986, 438607]
ygrid_extent = [6477539, 6479393]x = np.linspace(xgrid_extent[0], xgrid_extent[1], zstar.form[1])
y = np.linspace(ygrid_extent[0], ygrid_extent[1], zstar.form[0])
fig = go.Determine()
fig.add_trace(go.Floor(z=zstar, x=x, y=y))
fig.update_layout(scene = dict(
xaxis_title='X Location',
yaxis_title='Y Location',
zaxis_title='Depth'),
width=1000,
top=800,
margin=dict(r=20, l=10, b=10, t=10))
fig.present()
Once we run the code above, we get again the next interactive plot displaying the geological floor of the Hugin formation primarily based on the a number of encounters from the drilled wellbores.
The Volve dataset has numerous totally interpreted surfaces which have been generated from geological interpretations, together with seismic information.
This information incorporates the x and y areas of information factors throughout the sphere, in addition to our TVDSS information (z).
The info equipped on the Volve information portal is within the type of a .dat file, nevertheless, with a little bit of manipulation inside a textual content editor, it might probably simply be transformed to a CSV file and loaded utilizing pandas.
hugin_formation_surface = pd.read_csv('Information/Volve/Hugin_Fm_Top+ST10010ZC11_Near_190314_adj2_2760_EasyDC+STAT+DEPTH.csv')
As soon as we’ve got the information loaded, we will make issues simpler for ourselves by extracting the x, y and z information to variables.
x = hugin_formation_surface['x']
y = hugin_formation_surface['y']
z = hugin_formation_surface['z']
We then must create a commonly spaced grid between our min and max positions throughout the x and y information areas. This may be completed utilizing numpy’s meshgrid.
xi = np.linspace(x.min(), x.max(), 100)
yi = np.linspace(y.min(), y.max(), 100)
xi, yi = np.meshgrid(xi, yi)
There are a number of methods to interpolate between a sequence of information factors. The strategy chosen will rely on the shape your information is in (commonly sampled information factors vs irregularly sampled), information measurement and computing energy.
If we’ve got a big dataset just like the one right here, it is going to be far more computationally costly with a few of the strategies such because the Radial Foundation Perform.
For this instance, we’ll use the LinearNDInterpolator methodology inside scipy to construct our mannequin, after which apply it to our z (TVDSS) variable.
To ensure that us to interpolate between factors, we have to reshape xi, yi to 1D arrays for interpolation, as LinearNDInterpolator expects 1D array.
xir = xi.ravel()
yir = yi.ravel()interp = LinearNDInterpolator((x, y), z)
zi = interp(xir, yir)
As soon as that has been computed, we will transfer on to creating our 3D Floor with Plotly Graph Objects.
fig = go.Determine()fig.add_trace(go.Floor(z=zi, x=xi, y=yi, colorscale='Viridis'))
fig.update_layout(scene = dict(
xaxis_title='Easting (m)',
yaxis_title='Northing (m)',
zaxis_title='Depth',
zaxis=dict(autorange='reversed')),
width=1000,
top=800,
margin=dict(r=20, l=10, b=10, t=10))
fig.update_traces(contours_z=dict(present=True,
usecolormap=True,
project_z=True,
highlightcolor="white"))
fig.present()
Once we run the above code, we get again the next 3D floor plot of the Hugin Formation.
Once we examine this plot to the one generated from wellbore measurements, we will undoubtedly see similarities within the general form with the valley within the center. Nonetheless, the seismic-derived floor gives a lot higher element than the well-derived formation tops one.
