Estimation#
import geolime as geo
import numpy as np
from pyproj import CRS
import pyvista as pv
pv.set_jupyter_backend('panel')
geo.Project().set_crs(CRS("EPSG:20350"))
/tmp/ipykernel_3212/3172071959.py:6: PyVistaDeprecationWarning: `panel` backend is deprecated and is planned for future removal.
pv.set_jupyter_backend('panel')
Read Surface and Drillholes files#
dh = geo.read_file("../data/dh_pop_classif.geo")
domain_solid = geo.datasets.load("rocklea_dome/domain_mesh.dxf")
gs = geo.read_file("../data/surfaces.geo")
nu, nv = gs.shape
size_x, size_y = gs.axis()
size_z = 1
gs.properties()
['X',
'Y',
'Z',
'U',
'V',
'hard_top_Upper',
'hard_thick_Upper',
'hard_top_OreZone',
'hard_thick_OreZone',
'hard_top_Lower']
z_max = np.max(dh["Z_B"])
z_min = np.min(dh["Z_E"])
z_max, z_min
(478.6, 389.2)
Create Estimation Grid#
vx = geo.Voxel(
name="Grid3D",
shape=[nu - 2, nv - 2, (z_max - z_min) // size_z],
origin=[
gs.origin[0].round() + size_x,
gs.origin[1].round() + size_y,
z_min.round()
],
axis=[
[size_x, 0, 0],
[0, size_y, 0],
[0, 0, size_z],
]
)
bm = geo.BlockModel(name="Block3D", xyz=vx.coords())
bm.element_count()
761484
Flag grid cells position compared to surfaces#
geo.compute_surface_elevation_intersection(
obj=bm,
grid_object=gs,
z_prop="hard_top_OreZone",
new_property_name="Z_top_OreZone"
)
geo.compute_surface_elevation_intersection(
obj=bm,
grid_object=gs,
z_prop="hard_top_Lower",
new_property_name="Z_top_Lower"
)
geo.compute_surface_elevation_intersection(
obj=bm,
grid_object=gs,
z_prop="hard_top_Upper",
new_property_name="Z_top_Upper"
)
geo.compute_block_proportion_inside_solid(
grid_object=bm,
solid=domain_solid,
block_proportion_name='domain_prop', # name of the new property
block_discretization=(3,3,3) # proportion will be calculated on 81ths
)
bm.set_property_expr(name="domain_code", expr="1")
bm.set_property_expr(name="domain_code", expr="3", region="(Z > Z_top_Lower) & (Z < Z_top_OreZone)")
bm.set_property_expr(name="domain_code", expr="2", region="(Z >= Z_top_OreZone)")
bm.set_region_condition(
name="OreZone_Domain",
condition="(Z > Z_top_Lower) & (Z < Z_top_OreZone)"
)
bm.set_region_condition(
name="OVB",
condition="(Z >= Z_top_OreZone)"
)
bm.set_region_condition(
name="BDR",
condition="(Z <= Z_top_Lower) and (domain_prop < 0.1)"
)
bm.set_region_condition(
name="Ore",
condition="domain_prop >= 0.1"
)
dh_pv = dh.to_pyvista('Fe_pct')
bm_pv = bm.to_pyvista('domain_code')
domain_solid_pv = domain_solid.to_pyvista()
/opt/hostedtoolcache/Python/3.9.17/x64/lib/python3.9/site-packages/pyvista/core/grid.py:873: PyVistaDeprecationWarning: `UniformGrid` is deprecated. Use `ImageData` instead.
warnings.warn(
p = pv.Plotter()
p.add_mesh(dh_pv.tube(radius=10))
p.add_mesh(bm_pv, opacity=0.3)
p.add_mesh(domain_solid_pv, style='wireframe')
p.set_scale(zscale=10)
p.remove_scalar_bar(title='domain_code')
p.show()