Plotlyの3Dサーフェスプロットを使用して、地質表面を視覚化する

Visualize geological surfaces using Plotly's 3D surface plot.

Pythonのデータ可視化ライブラリを使用した地下の可視化

地球科学において、地下に存在する地質層の完全な理解が不可欠です。層の正確な位置とジオメトリを知ることにより、油田・ガス田の新しい潜在的なターゲットや、二酸化炭素のキャプチャ・ストレージの潜在的な場所を特定することが容易になります。これらの表面を精製するために、地震データからウェルログ派生の層の天井までの範囲を組み合わせて使用するなど、様々な方法があります。最も頻繁に、これらの技術はお互いを補完して最終表面を精製するために組み合わされます。

この記事は、地理的空間的変動を理解し、可視化するためのウェルログ測定値の外挿プロセスに焦点を当てた私の以前の記事から続きます。本記事では、対話型のPlotlyチャートを使用して、3Dサーフェスを作成する方法について説明します。

地質学的表面をモデリングすることは複雑なプロセスであり、しばしば複数の反復と改良を必要とするため、本記事ではPythonを使用してこのデータを可視化する非常に単純な例を示します。

地質層の形成を可視化する方法を見るために、エクイナーのボルブデータセットから2つのデータセットを使用します。

• エンジニアリングリーダーは何を気にしているのか？
• H1Bビザはデータ分析の洞察に基づいて承認されますか？
• データアナリストの仕事内容はどのように見えますか？

最初のデータセットは、井戸ピックから派生した層との28の交点から派生し、低解像度の表面を生成するためにクリギングの入力として使用されます。一方、2番目のデータセットは、解釈済みの地震データから派生し、はるかに高い解像度の表面を提供します。

両方のデータセットは、Equinor Volveデータセットから取得されています。詳細は、本記事の下部に記載されています。

また、以下のミニシリーズ内の記事も確認できます:

• PlotlyとPython：ペトロフィジカル＆地質データの対話型ヒートマップの作成
• pykrigeとmatplotlibを使用して地質的変化の空間的可視化を活用する
• Plotly Expressを使用した3Dラインプロットでのウェルパスの可視化

ライブラリとデータのインポート

データを処理する前に、必要なライブラリをインポートする必要があります。これらは以下の通りです:

• pandas — csv形式のデータを読み込むために使用する
• matplotlib — 可視化を作成するために使用する
• pykrige — クリギングを実行するために使用する
• numpy — 数値計算に使用する
• plotly.graph_objects — 3D表面を可視化するために使用する

import pandas as pdimport matplotlib.pyplot as pltimport plotly.graph_objects as gofrom pykrige import OrdinaryKrigingimport numpy as np

次に、pd.read_csv()を使用してデータを読み込むことができます。

このデータには、Volveフィールド内のすべてのウェルからの地質表面に関する情報が含まれているため、query()を使用して必要な形成物のデータを抽出することができます。今回は、Hugin Formationを見ていきます。

df = pd.read_csv('Data/Volve/Well_picks_Volve_v1 copy.csv')df_hugin = df.query('SURFACE == "Hugin Fm. VOLVE Top"')df_hugin

上記のコードを実行すると、以下の表が返されます。一部の井戸がHugin Formationに複数回遭遇していることに注意してください。これは、ウェルボアまたは形成ジオメトリによる複数回の遭遇のためである可能性があります。

TVDSSを外挿して地質学的な表面を生成する

以前の記事では、krigingと呼ばれるプロセスを使用して、測定ポイントの間の「隙間を埋める」方法に焦点を当ててきました。このプロセスの詳細については、この記事を参照してください。

データをロードしたら、pykrigeのOrdinaryKrigingメソッドを呼び出すことでkrigingプロセスを実行できます。

この呼び出しでは、井戸が地下層内で遭遇した位置の東西位置を表すxデータと南北位置を表すyデータを渡します。

Hugin層の表面を生成したいため、TVDSS（True Vertical Depth Subsea）の測定値を使用する必要があります。これにより、選択した基準面よりも表面がどの程度深いかが正確に反映されます。

OK = OrdinaryKriging(x=df_hugin['Easting'],                       y=df_hugin['Northing'],                       z=df_hugin['TVDSS'],                      variogram_model='linear',                      verbose=True, enable_plotting=True)

モデルが生成されたら、井戸/侵入ポイントの範囲全体をカバーする2つの配列に適用できます。

これにより、値のグリッドが生成され、上記で生成したOrdinaryKrigingオブジェクトに渡すことができます。

gridx = np.arange(433986, 438607, 50, dtype='float64')gridy = np.arange(6477539, 6479393, 50,dtype='float64')zstar, ss = OK.execute('grid', gridx, gridy)

最後に、matplotlibのimshowメソッドを使用して、表面の簡単な2Dマップビューを生成できます。

fig, ax = plt.subplots(figsize=(15,5))# Create a 2D image plot of the data in 'zstar'# The 'extent' parameter sets the bounds of the image in data coordinates# 'origin' parameter sets the part of the image that corresponds to the origin of the axesimage = ax.imshow(zstar, extent=(433986, 438607, 6477539, 6479393),                   origin='lower')# Set the labels for the x-axis and y-axisax.set_xlabel('X Location (m)', fontsize=14, fontweight='bold')ax.set_ylabel('Y Location (m)', fontsize=14, fontweight='bold')# Add contourscontours = ax.contour(gridx, gridy, zstar, colors='black')colorbar = fig.colorbar(image)colorbar.set_label('DTC (us/ft)', fontsize=14, fontweight='bold')# Display the plotplt.show()

Plotlyを使用したシンプルな3D表面プロットの作成

2D表面を3Dに変換するには、Plotlyを使用する必要があります。これを行うために、matplotlibを使用することもできますが、私の経験から、Plotlyを使用して3Dビジュアリゼーションを作成する方が簡単でスムーズでインタラクティブです。

以下のコードでは、まず座標のグリッドを作成します。これには、numpyのlinspace関数を使用します。この関数は、指定された範囲内の均等に間隔を開けた数値セットを作成します。私たちのデータセットと例では、範囲はxgrid_extentとygrid_extentの最小値から最大値まで広がります。

この範囲で使用される値の総数は、上記で作成したzstarグリッド内に存在するx要素とy要素の数に相当します。

グリッドが形成されたら、Plotlyを呼び出します。

最初に、フィギュアオブジェクトを作成し、fig.add_traceを使用して3Dサーフェスプロットをフィギュアに追加します。

これが追加されたら、プロットのレイアウトを微調整して、軸ラベル、定義された幅と高さ、およびパディングを持つようにします。

xgrid_extent = [433986、438607]ygrid_extent = [6477539、6479393]x = np.linspace（xgrid_extent [0]、xgrid_extent [1]、zstar.shape [1]）y = np.linspace（ygrid_extent [0]、ygrid_extent [1]、zstar.shape [0]）fig = go.Figure（）fig.add_trace（go.Surface（z = zstar、x = x、y = y））fig.update_layout（scene = dict（                    xaxis_title ='Xロケーション'、                    yaxis_title ='Yロケーション'、                    zaxis_title ='深さ'）、                  width = 1000、                  height = 800、                  margin = dict（r = 20、l = 10、b = 10、t = 10））fig.show（）

上記のコードを実行すると、掘削ウェルボアからの複数のエンカウントに基づくヒュギン層の地質的表面を示す以下のインタラクティブプロットが返されます。

Plotlyで完全に解釈された表面を表示する

ボルブデータセットには、地質解釈、地震データを含む、完全に解釈された表面がいくつかあります。

このデータには、フィールド全体のデータポイントのxおよびyの場所、およびTVDSSデータ（z）が含まれています。

Volveデータポータルで提供されるデータは、.datファイルの形式ですが、テキストエディタ内での少しの操作により、CSVファイルに簡単に変換してpandasを使用して読み込むことができます。

hugin_formation_surface = pd.read_csv('Data/Volve/Hugin_Fm_Top+ST10010ZC11_Near_190314_adj2_2760_EasyDC+STAT+DEPTH.csv')

データを読み込んだら、x、y、およびzデータを変数に抽出することで、自分自身を簡単にすることができます。

x = hugin_formation_surface ['x'] y = hugin_formation_surface ['y'] z = hugin_formation_surface ['z']

次に、xおよびyデータの最小値と最大値の間に正規に間隔を開けたグリッドを作成する必要があります。これは、numpyのmeshgridを使用して行うことができます。

xi = np.linspace(x.min()、x.max()、100)yi = np.linspace(y.min()、y.max()、100)xi、yi = np.meshgrid（xi、yi）

データポイント間を補間するには、いくつかの方法があります。選択する方法は、データがどのように形成されているか（定期的にサンプリングされたデータポイント対非定期的にサンプリングされたデータポイント）、データサイズ、および計算能力に依存します。

ここでのように、大規模なデータセットがある場合、Radial Basis Functionなどのいくつかの方法を使用すると、計算のコストが高くなります。

この例では、scipy内のLinearNDInterpolatorメソッドを使用してモデルを構築し、z（TVDSS）変数に適用します。

点間を補間するためには、LinearNDInterpolatorが1次元配列を期待するため、xi、yiを1次元配列に変形する必要があります。

xir = xi.ravel()yir = yi.ravel()interp = LinearNDInterpolator((x, y), z)zi = interp(xir, yir)

これが計算されたら、Plotlyグラフオブジェクトを使って3Dサーフェスを作成することができます。

fig = go.Figure()fig.add_trace(go.Surface(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,                  height=800,                  margin=dict(r=20, l=10, b=10, t=10))fig.update_traces(contours_z=dict(show=True,                                   usecolormap=True,                                   project_z=True,                                  highlightcolor="white"))fig.show()

上記のコードを実行すると、次のHugin Formationの3Dサーフェスプロットが表示されます。

このプロットを、ウェルボア測定から生成されたプロットと比較すると、中央の谷を含めた全体的な形状に類似性があることが確認できます。しかし、地震によるサーフェスは、ウェルから得られた層の形状のものよりもはるかに詳細な情報を提供します。

概要

この短いチュートリアルでは、Plotlyの3Dサーフェスプロットを使用して、地質サーフェスのインタラクティブな3D可視化を生成する方法を紹介しました。ウェルログから得られた層の形状からは、非常に基本的な3Dサーフェスを生成することができます。これは、Hugin Formationに交差したウェルボアが制限されているため、低解像度のサーフェスになるためです。

対照的に、地震による地平線など、より詳細な測定点がある場合は、より現実的な印象を形成することができます。

両方の方法は同じくらい有効ですが、ウェルログから派生した層の形状から外挿する場合、その地域のその地層の完全なイメージを生成できなかった可能性があることに注意する必要があります。

使用されたデータセット

このチュートリアルで使用されたデータは、Equinorが2018年にリリースしたVolveデータセットのサブセットです。データセットの詳細、ライセンスなどは、以下のリンクで確認できます。

Volveフィールドデータセット

Equinorは、Volveフィールドからの完全なデータセットをリリースしました（2008-2016）。研究や調査のためにダウンロードするには、ここをクリックしてください。

www.equinor.com

Volveデータライセンスは、CC BY 4.0ライセンスに基づいています。ライセンス契約の詳細については、以下のリンクを参照してください。

https://cdn.sanity.io/files/h61q9gi9/global/de6532f6134b9a953f6c41bac47a0c055a3712d3.pdf?equinor-hrs-terms-and-conditions-for-licence-to-data-volve.pdf

読んでくれてありがとうございます。この記事を購読して、私の記事をメールで受け取ることをお勧めします。 ここで購読できます！

第二に、メンバーシップに登録することで、完全なVoAGI体験を得て、他の何千人ものライターと私をサポートすることができます。月額$5の費用がかかりますが、素晴らしいVoAGI記事に完全にアクセスできるだけでなく、自分の文章でお金を稼ぐチャンスもあります。

私のリンクを使用して登録すると、費用の一部が直接私に支払われますが、追加料金はかかりません。もしそうしていただける場合、サポートいただきありがとうございます。

We will continue to update VoAGI; if you have any questions or suggestions, please contact us!