Compute Intersection Between An Object and A GridObject¶

In [1]:

import geolime as geo
from pyproj import CRS

import geolime as geo from pyproj import CRS

In [3]:

dh = geo.datasets.load("rocklea_dome/dh_hyper.geo")

dh = geo.datasets.load("rocklea_dome/dh_hyper.geo")

In [4]:

grid = geo.Voxel(
    name="NewGrid", 
    shape=[75, 150], 
    origin=[545400, 7472700, 430], 
    axis=[[50, 0, 0], [0, 50, 0]]
)

grid = geo.Voxel( name="NewGrid", shape=[75, 150], origin=[545400, 7472700, 430], axis=[[50, 0, 0], [0, 50, 0]] )

In [5]:

gs = geo.GriddedSurface(name="GridSurf", xyz=grid.coords())

gs = geo.GriddedSurface(name="GridSurf", xyz=grid.coords())

In [6]:

gs.set_property_expr(name="XY", expr="sin(X / 2) * sin (Y / 2)")

gs.set_property_expr(name="XY", expr="sin(X / 2) * sin (Y / 2)")

GriddedSurface and other object may not be exaclty superimposed. This means if we want to know the surface coordinates at the XY location of drillholes or blockmodel we need to compute the equivalent of the vertical intersection (or projection) of the object to the surface.

In [7]:

geo.compute_surface_elevation_intersection(
    obj=dh,
    grid_object=gs,
    z_prop="XY",
    new_property_name="intersection",
    ignore_warnings=False
)

geo.compute_surface_elevation_intersection( obj=dh, grid_object=gs, z_prop="XY", new_property_name="intersection", ignore_warnings=False )

In [8]:

geo.plot_2d(
    dh, 
    property="intersection", 
    agg_method="mean", 
    crs=CRS("EPSG:20350"), 
    width=650, 
    height=650, 
    mapbox_zoom=12
)

geo.plot_2d( dh, property="intersection", agg_method="mean", crs=CRS("EPSG:20350"), width=650, height=650, mapbox_zoom=12 )