1.3b: 2-D sections

Importing

import gempy as gp
import gempy_viewer as gpv
import numpy as np

np.random.seed(1234)

Setup the model

Importing the data from CSV-files and setting extent and resolution

data_path = 'https://raw.githubusercontent.com/cgre-aachen/gempy_data/master/'

geo_model: gp.data.GeoModel = gp.create_geomodel(
    project_name='Tutorial_ch1_1_Basics',
    extent=[0, 2000, 0, 2000, 0, 750],
    resolution=[20, 20, 20],  # * Here we define the resolution of the voxels
    refinement=4,  # * Here we define the number of octree levels. If octree levels are defined, the resolution is ignored.
    importer_helper=gp.data.ImporterHelper(
        path_to_orientations=data_path + "/data/input_data/getting_started/simple_fault_model_orientations.csv",
        path_to_surface_points=data_path + "/data/input_data/getting_started/simple_fault_model_points.csv",
        hash_surface_points="4cdd54cd510cf345a583610585f2206a2936a05faaae05595b61febfc0191563",
        hash_orientations="7ba1de060fc8df668d411d0207a326bc94a6cdca9f5fe2ed511fd4db6b3f3526"
    )
)

gp.map_stack_to_surfaces(
    gempy_model=geo_model,
    mapping_object=  # TODO: This mapping I do not like it too much. We should be able to do it passing the data objects directly
    {
        "Fault_Series": 'Main_Fault',
        "Strat_Series": ('Sandstone_2', 'Siltstone', 'Shale', 'Sandstone_1')
    }
)

gp.set_is_fault(
    frame=geo_model.structural_frame,
    fault_groups=['Fault_Series']
)

Surface points hash:  4cdd54cd510cf345a583610585f2206a2936a05faaae05595b61febfc0191563
Orientations hash:  7ba1de060fc8df668d411d0207a326bc94a6cdca9f5fe2ed511fd4db6b3f3526

Structural Groups:	StructuralGroup: Name: Fault_Series Structural Relation: StackRelationType.FAULT Elements: StructuralElement: Name: Main_Fault StructuralGroup: Name: Strat_Series Structural Relation: StackRelationType.ERODE Elements: StructuralElement: Name: Sandstone_2 StructuralElement: Name: Siltstone StructuralElement: Name: Shale StructuralElement: Name: Sandstone_1
Fault Relations:	Fault_Seri... Strat_Seri... Fault_Series Strat_Series
	True False

Add sections

pass section dictionary with startpoint, endpoint and resolution for every section:

gp.set_section_grid(
    grid=geo_model.grid,
    section_dict={
        'section1': ([0, 0], [2000, 2000], [100, 80]),
        'section2': ([800, 0], [800, 2000], [150, 100]),
        'section3': ([0, 200], [1500, 500], [200, 150])
    }  # p1,p2,resolution
)

Active grids: GridTypes.DENSE|SECTIONS|NONE

	start	stop	resolution	dist
section1	(0.0, 0.0)	(2000.0, 2000.0)	(100, 80)	2828.427125
section2	(800.0, 0.0)	(800.0, 2000.0)	(150, 100)	2000.000000
section3	(0.0, 200.0)	(1500.0, 500.0)	(200, 150)	1529.705854

Add topography

gp.set_topography_from_random(
    grid=geo_model.grid,
    fractal_dimension=1.2,
    d_z=np.array([300, 750]),
    topography_resolution=np.array([50, 50])
)

Active grids: GridTypes.DENSE|TOPOGRAPHY|SECTIONS|NONE

Topography(_regular_grid=RegularGrid(resolution=array([20, 20, 20]), extent=array([   0., 2000.,    0., 2000.,    0.,  750.]), values=array([[  50.  ,   50.  ,   18.75],
       [  50.  ,   50.  ,   56.25],
       [  50.  ,   50.  ,   93.75],
       ...,
       [1950.  , 1950.  ,  656.25],
       [1950.  , 1950.  ,  693.75],
       [1950.  , 1950.  ,  731.25]], shape=(8000, 3)), mask_topo=array([], shape=(0, 3), dtype=bool), _transform=None, _base_resolution=array([2, 2, 2])), values_2d=array([[[   0.        ,    0.        ,  525.81364595],
        [   0.        ,   40.81632653,  534.4311714 ],
        [   0.        ,   81.63265306,  545.21138351],
        ...,
        [   0.        , 1918.36734694,  493.45417828],
        [   0.        , 1959.18367347,  505.51902125],
        [   0.        , 2000.        ,  516.62320964]],

       [[  40.81632653,    0.        ,  520.02419799],
        [  40.81632653,   40.81632653,  529.19849936],
        [  40.81632653,   81.63265306,  540.73889549],
        ...,
        [  40.81632653, 1918.36734694,  485.70563878],
        [  40.81632653, 1959.18367347,  498.13591443],
        [  40.81632653, 2000.        ,  509.94870633]],

       [[  81.63265306,    0.        ,  512.29034438],
        [  81.63265306,   40.81632653,  522.38242549],
        [  81.63265306,   81.63265306,  533.74749759],
        ...,
        [  81.63265306, 1918.36734694,  475.96047534],
        [  81.63265306, 1959.18367347,  488.33556325],
        [  81.63265306, 2000.        ,  500.75092545]],

       ...,

       [[1918.36734694,    0.        ,  536.8572301 ],
        [1918.36734694,   40.81632653,  549.32197157],
        [1918.36734694,   81.63265306,  561.72341993],
        ...,
        [1918.36734694, 1918.36734694,  503.24054915],
        [1918.36734694, 1959.18367347,  513.73925164],
        [1918.36734694, 2000.        ,  525.47837607]],

       [[1959.18367347,    0.        ,  534.68908826],
        [1959.18367347,   40.81632653,  545.65518713],
        [1959.18367347,   81.63265306,  558.27313536],
        ...,
        [1959.18367347, 1918.36734694,  501.0502999 ],
        [1959.18367347, 1959.18367347,  512.95973605],
        [1959.18367347, 2000.        ,  524.20361683]],

       [[2000.        ,    0.        ,  530.78028142],
        [2000.        ,   40.81632653,  540.15206274],
        [2000.        ,   81.63265306,  551.79148373],
        ...,
        [2000.        , 1918.36734694,  498.47649577],
        [2000.        , 1959.18367347,  510.05476038],
        [2000.        , 2000.        ,  521.16890329]]], shape=(50, 50, 3)), source=None, values=array([[   0.        ,    0.        ,  525.81364595],
       [   0.        ,   40.81632653,  534.4311714 ],
       [   0.        ,   81.63265306,  545.21138351],
       ...,
       [2000.        , 1918.36734694,  498.47649577],
       [2000.        , 1959.18367347,  510.05476038],
       [2000.        , 2000.        ,  521.16890329]], shape=(2500, 3)), resolution=(50, 50), raster_shape=())

Active grids:

geo_model.grid.active_grids

<GridTypes.DENSE|TOPOGRAPHY|SECTIONS|NONE: 1050>

gpv.plot_section_traces(geo_model)

<function plot_section_traces at 0x7f582f2ef8a0>

geo_model.interpolation_options.mesh_extraction = False
sol = gp.compute_model(geo_model)

Setting Backend To: AvailableBackends.PYTORCH
GPU requested but unavailable; falling back to CPU (GEMPY_GPU_FALLBACK=True)
Setting Backend To: AvailableBackends.PYTORCH
Chunking done: 8 chunks

gpv.plot_2d(geo_model, section_names=['topography'])

/opt/buildAgent/work/3a8738c25f60c3c9/venv/lib/python3.14/site-packages/gempy_viewer/API/_plot_2d_sections_api.py:112: UserWarning: Section contacts not implemented yet. We need to pass scalar field for the sections grid
  warnings.warn(

<gempy_viewer.modules.plot_2d.visualization_2d.Plot2D object at 0x7f5820cff3d0>

gpv.plot_2d(geo_model, section_names=['section1'])

/opt/buildAgent/work/3a8738c25f60c3c9/venv/lib/python3.14/site-packages/gempy_viewer/API/_plot_2d_sections_api.py:112: UserWarning: Section contacts not implemented yet. We need to pass scalar field for the sections grid
  warnings.warn(

<gempy_viewer.modules.plot_2d.visualization_2d.Plot2D object at 0x7f5820c4b450>

gpv.plot_2d(
    model=geo_model,
    section_names=['section1', 'section2', 'section3', 'topography'],
    show_topography=True
)


# sphinx_gallery_thumbnail_number = 4

/opt/buildAgent/work/3a8738c25f60c3c9/venv/lib/python3.14/site-packages/gempy_viewer/API/_plot_2d_sections_api.py:112: UserWarning: Section contacts not implemented yet. We need to pass scalar field for the sections grid
  warnings.warn(

<gempy_viewer.modules.plot_2d.visualization_2d.Plot2D object at 0x7f580ed83e50>