|
| 1 | +""" |
| 2 | +Modeling step by step |
| 3 | +^^^^^^^^^^^^^^^^^^^^^ |
| 4 | +
|
| 5 | +This tutorial demonstrates step-by-step geological modeling using the `gempy` and `gempy_viewer` libraries. We will start by importing the necessary packages, loading input data, creating a geological model, and then visualizing the results. |
| 6 | +""" |
| 7 | + |
| 8 | +# %% |
| 9 | +# Import minimal requirements |
| 10 | +# We need to import the `gempy` library for geological modeling and `gempy_viewer` for visualization. |
| 11 | +import gempy as gp |
| 12 | +import gempy_viewer as gpv |
| 13 | + |
| 14 | +# %% |
| 15 | +# Define the path to input data |
| 16 | +# Here, we provide two ways to define the path to the input data: using a URL or a local path. |
| 17 | +# Uncomment the first two lines if you want to use the online data source. |
| 18 | + |
| 19 | +# data_path = 'https://raw.githubusercontent.com/cgre-aachen/gempy_data/master/' |
| 20 | +# path_to_data = data_path + "/data/input_data/jan_models/" |
| 21 | + |
| 22 | +# For this tutorial, we will use the local path: |
| 23 | +data_path = '../../' |
| 24 | +path_to_data = data_path + "/data/input_data/jan_models/" |
| 25 | + |
| 26 | +# %% |
| 27 | +# Create a GeoModel instance |
| 28 | +# We create a GeoModel instance with a specified project name and extent. |
| 29 | +# The ImporterHelper class is used to specify the paths to the orientation and surface points data. |
| 30 | + |
| 31 | +geo_model = gp.create_geomodel( |
| 32 | + project_name='tutorial_model', |
| 33 | + extent=[0, 2500, 0, 1000, 0, 1110], |
| 34 | + refinement=6, |
| 35 | + importer_helper=gp.data.ImporterHelper( |
| 36 | + path_to_orientations=path_to_data + "tutorial_model_orientations.csv", |
| 37 | + path_to_surface_points=path_to_data + "tutorial_model_surface_points.csv" |
| 38 | + ) |
| 39 | +) |
| 40 | + |
| 41 | +# %% |
| 42 | +# Displaying simple data cross section |
| 43 | +# We use the `gempy_viewer` to visualize the initial cross-section of our geological model. |
| 44 | +gpv.plot_2d(geo_model) |
| 45 | + |
| 46 | +# %% |
| 47 | +# Map geological series to surfaces |
| 48 | +# Here, we map the geological series to specific surfaces. This step is crucial for defining the stratigraphic relationships in our model. |
| 49 | +gp.map_stack_to_surfaces( |
| 50 | + gempy_model=geo_model, |
| 51 | + mapping_object={ |
| 52 | + "Strat_Series1": ('rock3'), |
| 53 | + "Strat_Series2": ('rock2', 'rock1'), |
| 54 | + } |
| 55 | +) |
| 56 | + |
| 57 | +# %% |
| 58 | +# Update transformation and interpolation options |
| 59 | +# We update the model with anisotropy settings and specify various interpolation options to refine the model's accuracy. |
| 60 | + |
| 61 | +geo_model.update_transform(auto_anisotropy=gp.data.GlobalAnisotropy.CUBE) |
| 62 | + |
| 63 | +interpolation_options: gp.data.InterpolationOptions = geo_model.interpolation_options |
| 64 | + |
| 65 | +interpolation_options.mesh_extraction = True |
| 66 | +interpolation_options.evaluation_options.compute_scalar_gradient = True |
| 67 | + |
| 68 | +interpolation_options.kernel_options.range = 1 |
| 69 | +interpolation_options.evaluation_options.number_octree_levels_surface = 6 |
| 70 | +interpolation_options.evaluation_options.curvature_threshold = 0.6 |
| 71 | + |
| 72 | +# %% |
| 73 | +# Compute the geological model |
| 74 | +# We use the specified backend (in this case, PyTorch) to compute the model. |
| 75 | +gp.compute_model( |
| 76 | + gempy_model=geo_model, |
| 77 | + engine_config=gp.data.GemPyEngineConfig(backend=gp.data.AvailableBackends.PYTORCH) |
| 78 | +) |
| 79 | + |
| 80 | +# %% |
| 81 | +# Visualize the model: 2D cross-section without scalar field |
| 82 | +# After computing the model, we visualize it again in 2D without the scalar field. |
| 83 | +gpv.plot_2d(geo_model, show_scalar=False, series_n=1) |
| 84 | + |
| 85 | +# %% md |
| 86 | +# Visualize the model: 2D cross-section with scalar field |
| 87 | +# In this cell, we visualize the 2D cross-section with the scalar field enabled. |
| 88 | +# %% |
| 89 | +gpv.plot_2d(geo_model, show_scalar=True, series_n=1) |
| 90 | + |
| 91 | +# %% |
| 92 | +# Visualize the model in 3D |
| 93 | +# Finally, we create a 3D visualization of the geological model without lithological coloring and image. |
| 94 | +gpv.plot_3d(geo_model, show_lith=False, image=False) |
| 95 | +# sphinx_gallery_thumbnail_number = -1 |
| 96 | + |
| 97 | +# %% |
| 98 | +# ### Coming up next |
| 99 | +# Additional: Manually adding data (optional) |
| 100 | +# Here is an example of how you can manually add surface points to the model. Uncomment and modify the code as needed. |
| 101 | + |
| 102 | +''' |
| 103 | +gp.add_surface_points( |
| 104 | + geo_model=geo_model, |
| 105 | + x=[458, 612], |
| 106 | + y=[0, 0], |
| 107 | + z=[-107, -14], |
| 108 | + elements_names=['surface1', 'surface1'] |
| 109 | +) |
| 110 | +''' |
| 111 | + |
| 112 | +# %% |
| 113 | +# Displaying data cross section (optional) |
| 114 | +# You can re-plot the 2D cross-section if needed. |
| 115 | +# gpv.plot_2d(geo_model) |
| 116 | +# %% |
0 commit comments