Localised Uniform Conditioning#
Uniform Conditioning in mining geostatistics is used to estimate the amount of recoverable resources (tonnage and metal content) of blocks inside a panel given the panel grade. It works by directly calculating the tonnage and metal content of blocks within a panel that are above a specified cutoff grade. However, this method has some drawbacks:
It doesn’t pinpoint where the recoverable blocks are spatially.
It’s not suitable for estimating small blocks when the data are sparse compared to the block size.
Localized Uniform Conditioning (LUC) is designed to address these limitations. The key idea behind LUC is to leverage the grade ranking obtained from Ordinary Kriging on the block model at SMU scale while ensuring the final grade distribution aligns with the one obtained through Uniform Conditioning.
Library Import#
import numpy as np
import geolime as geo
import pyvista as pv
from pyproj import CRS
from plotly.subplots import make_subplots
import plotly.graph_objects as go
import json
pv.set_jupyter_backend('html')
geo.Project().set_crs(CRS("EPSG:20350"))
Read Surface, Drillholes and Block Model files#
domain_solid = geo.datasets.load("rocklea_dome/domain_mesh.dxf")
domain_solid.name = "channel_hg"
dh = geo.read_file("../data/domained_drillholes.geo")
bm_panel = geo.read_file("../data/block_model_panel.geo")
bm_smu = geo.read_file("../data/block_model_smu.geo")
Here, both the panel and smu estimation grid have been previously created.
Thus, they are only loaded, and not redefined.
Estimation#
domain_solid.contains(dh)
dh.set_region("channel_hg", "in_channel_hg and Fe_pct.notna()")
Change Of Support#
Uniform conditioning (UC) estimates the recoverable resources of blocks inside a panel based on the panel’s estimated grade. This technique requires understanding the relationship between grades at different scales, a concept known as change of support.
The goal is to predict the grade and tonnage of selectively mined units (SMUs).
These SMUs are the smallest mineable units, and their grades exhibit different statistical properties compared to the larger panels or the original sample data points.
A change of support model quantifies how the statistical distribution of grades changes as the volume or area under consideration changes.
The change of support used in Uniform Condition is the Discrete Gaussian Model (DGM) and relates the distribution of grades at different supports through a series of Hermite polynomials. In GeoLime the object which relates to the conversion of raw data to gaussian data and model the gaussian distribution using Hermite polynomials is calld Anamorphosis. It take into account the data from a given property of an object and the number of coefficient to be used for the polynomial fit.
Here we use 50 polynomials, this is an acceptable trade-off between quality and rapidity of execution of the fit. Indeed, more polynomials lead to longer computation time.
anam = geo.Anamorphosis(
data=dh.data("Fe_pct", "channel_hg"),
coefficient_number=50
)
geo.hermite_plot(anam)
Covariance & Neighborhood settings#
Read Covariance Model from json file
with open("../data/cov_model.json") as json_data:
cov_dict = json.load(json_data)
covariance = cov_dict["nugget"] * geo.Nugget()
for i in range(len(cov_dict["structures"])):
covariance += cov_dict["structures"][i]["contribution"] * getattr(
geo,
cov_dict["structures"][i]["structure_type"]
)(
metric=geo.Metric(
scales=list(cov_dict["structures"][0]["anisotropy"]["ranges"].values()),
angles=list(cov_dict["structures"][0]["anisotropy"]["rotation"].values()))
)
neighborhood = geo.MinMaxPointsNeighborhood(
dimension=3,
metric=geo.Metric(
angles=[0, 0, 0],
scales=[3000, 3000, 20]
),
min_n=1,
max_n=20,
)
Localised Uniform Conditioning#
Localised Uniform Conditioning estimator is constructed the same way as the Kriging estimator using a covariance model and a neighborhood model. The estimator performs uniform conditioning on the panel BlockModel and creates property of proportion and mass (generally called metal) content above each given cutoff for each panel. It also creates the property of the kriging result, allowing user to perfom further QA/QC analysis. The estimator also performs the grade localisation on the SMU BlockModel where the property of grade is created.
uc_estimator = geo.LocalisedUniformConditioning(
covariance_model=covariance,
neighborhood_model=neighborhood,
)
For this example, the specified cut-off grades were chosen to cover the entirety of the iron values in the dataset.
cutoff_grades=np.arange(0, 60, 5, dtype=float)
ore = uc_estimator.solve(
obj=dh,
obj_region='channel_hg',
obj_attribute='Fe_pct',
anamorphosis=anam,
support=bm_smu,
support_region='channel_hg',
panel_support=bm_panel,
panel_support_region='channel_hg',
cutoff_grades=cutoff_grades,
support_attribute='Fe_LUC',
discr=np.array([2, 2, 1]),
panel_discr=np.array([2, 2, 1])
)
Based on proportion and metal property above each cutoff the mean grade above each cutoff is then computed for each panel. This is useful for grade tonnage curves and selectivity analyses.
for prop in bm_panel.properties():
if "metal_above_" in prop :
grade_prop = prop.lstrip("metal_above")
bm_panel.set_property(f"mean_grade_above_{grade_prop}", bm_panel[prop] / bm_panel[f"ore_prop_above_{grade_prop}"])
bm_panel.to_file("../data/block_model_panel")
bm_smu.to_file("../data/block_model_smu")
Results Validation#
dh.plot_2d(property="Fe_pct", agg_method="max", interactive_map=True, region='channel_hg', region_agg='any')
Panels#
An estimation of each panel’s grade is done by ordinary kriging during the Uniform Conditioning process.
bm_panel.plot_2d(property="Fe_LUC", agg_method="mean", interactive_map=True, region="channel_hg", region_agg="any")
For each panel it is possible to study the distribution and selectivity at the SMU scale. For example here, we select a random panel in the BlockModel and compute the mass proportion (also called tonnage), mass content (also called metal) and mean grade above the computed cutoff. We use comprehension list to select the property matching the prefix “ore_prop_above_”, “metal_above_” and “mean_grade_above_” in order to only plot these values in the figure.
panel = 5
cutoff_prop = [bm_panel.data(prop, region="channel_hg")[panel] for prop in bm_panel.properties() if "ore_prop_above_" in prop]
cutoff_metal = [bm_panel.data(prop, region="channel_hg")[panel] for prop in bm_panel.properties() if "metal_above_" in prop]
cutoff_mean_grade =[bm_panel.data(prop, region="channel_hg")[panel] for prop in bm_panel.properties() if "mean_grade_above_" in prop]
We then use Plotly to create the 3 different scatter of the tonnage, metal and mean grade above the computed cutoffs.
fig = go.Figure(data=go.Scatter(x=cutoff_grades, y=cutoff_prop, mode='markers+lines', line_color="black", marker_size=5, showlegend=False))
fig.update_layout(xaxis_title="Cutoff grades", yaxis_title="Tonnage", title=f"Tonnage - panel {panel}")
fig.show()
fig = go.Figure(data=go.Scatter(x=cutoff_grades, y=cutoff_metal, mode='markers+lines', line_color="black", marker_size=5, showlegend=False))
fig.update_layout(xaxis_title="Cutoff grade", yaxis_title="Metal quantity", title=f"Metal quantity - panel {panel}")
fig.show()
fig = go.Figure(data=go.Scatter(y=cutoff_mean_grade, x=cutoff_grades, mode='markers+lines', line_color="black", marker_size=5, showlegend=False))
fig.update_layout(yaxis_title="Mean grade", xaxis_title="Cutoff grade", title=f"Mean grades - panel {panel}")
fig.show()
Since the metal quantity and the proportion above different cutoff values are obtained, these can be visualized in maps for specified cutoff values.
cutoff = ["25_0", "40_0"]
properties = [f"ore_prop_above_{cutoff[0]}", f"ore_prop_above_{cutoff[1]}"]
fig = make_subplots(
rows=1,
cols=2,
subplot_titles=[f"Proportion above {cutoff[0]}", f"Proportion above {cutoff[1]}"],
specs=[[{"type": "mapbox"} , {"type": "mapbox"}]]
)
# Add plots
for i in range(2):
obj = bm_panel
prop = properties[i]
trace=geo.plot_2d(
[
{"georef_object": obj, "property": prop, "agg_method": "mean", "region":"channel_hg", "region_agg":"any"},
]
)
fig.add_trace(trace["data"][0],row=1, col=i+1)
mapbox=f"mapbox{i+1}"
fig.update_layout({mapbox:{"zoom":12, "style":'open-street-map', "center":trace["layout"]["mapbox"]["center"]}})
fig.update_traces({"coloraxis":"coloraxis"}, col=i+1)
fig.update_layout({"coloraxis":{"cmin":0,"cmax":1,"showscale":True,"colorscale":"Viridis" }})
fig.update_layout(width=800, height=600)
fig.update_layout(margin=dict(l=0, r=0, t=30, b=10))
#Set global colorscale
fig.update_layout(showlegend=False)
fig.show()
properties = [f"metal_above_{cutoff[0]}",f"metal_above_{cutoff[1]}"]
fig = make_subplots(
rows=1,
cols=2,
subplot_titles=[f"Metal above {cutoff[0]}", f"Metal above {cutoff[1]}"],
specs=[[{"type": "mapbox"} , {"type": "mapbox"}]]
)
# Add plots
for i in range(2):
obj = bm_panel
prop = properties[i]
trace=geo.plot_2d(
[
{"georef_object": obj, "property": prop, "agg_method": "mean", "region":"channel_hg", "region_agg":"any"},
]
)
fig.add_trace(trace["data"][0],row=1, col=i+1)
mapbox=f"mapbox{i+1}"
fig.update_layout({mapbox:{"zoom":12, "style":'open-street-map', "center":trace["layout"]["mapbox"]["center"]}})
fig.update_traces({"coloraxis":"coloraxis"}, col=i+1)
fig.update_layout({"coloraxis":{"showscale":True,"colorscale":"Viridis"}})
fig.update_layout(width=800, height=600)
fig.update_layout(margin=dict(l=0, r=0, t=30, b=10))
#Set global colorscale
fig.update_layout(showlegend=False)
fig.show()
properties = [f"mean_grade_above_{cutoff[0]}",f"mean_grade_above_{cutoff[1]}"]
fig = make_subplots(
rows=1,
cols=2,
subplot_titles=[f"Mean grade above {cutoff[0]}", f"Mean grade above {cutoff[1]}"],
specs=[[{"type": "mapbox"} , {"type": "mapbox"}]]
)
# Add plots
for i in range(2):
obj = bm_panel
prop = properties[i]
trace=geo.plot_2d(
[
{"georef_object": obj, "property": prop, "agg_method": "mean", "region":"channel_hg", "region_agg":"any"},
]
)
fig.add_trace(trace["data"][0],row=1, col=i+1)
mapbox=f"mapbox{i+1}"
fig.update_layout({mapbox:{"zoom":12, "style":'open-street-map', "center":trace["layout"]["mapbox"]["center"]}})
fig.update_traces({"coloraxis":"coloraxis"}, col=i+1)
fig.update_layout({"coloraxis":{"showscale":True,"colorscale":"Viridis" }})
fig.update_layout(width=750, height=600)
fig.update_layout(margin=dict(l=0, r=0, t=30, b=10))
#Set global colorscale
fig.update_layout(showlegend=False)
fig.show()
bm_panel_pv = bm_panel.to_pyvista(properties="Fe_LUC", region="channel_hg")
p = pv.Plotter()
p.add_mesh(bm_panel_pv)
p.set_scale(zscale=20)
p.show()
Selective Mining Units#
To get results at the scale of smus, an ordinary kriging is executed on them. Then, these relative grade values are used when computing more precise grades from panels’ distributions.
geo.plot_2d([{"georef_object": bm_smu, "property": "Fe_LUC", "agg_method": "mean"}])
bm_smu_pv = bm_smu.to_pyvista(properties="Fe_LUC", region="channel_hg")
p = pv.Plotter()
p.add_mesh(bm_smu_pv)
p.set_scale(zscale=20)
p.show()
QA-QC#
Swath Plot on SMU BlockModel can be compared to Drillholes
geo.swath_plot(
[
{"obj":bm_smu, "attribute":"Fe_LUC", "region":"channel_hg", "swath_interval":50, "axis":"Y", "color": "blue"},
{"obj":dh, "attribute":"Fe_pct", "region":"channel_hg", "swath_interval":50, "axis":"Y_M", "color": "red"}
]
)
geo.swath_plot(
[
{"obj":bm_smu, "attribute":"Fe_LUC", "region":"channel_hg", "swath_interval":50, "axis":"X", "color": "blue"},
{"obj":dh, "attribute":"Fe_pct", "region":"channel_hg", "swath_interval":50, "axis":"X_M", "color": "red"}
]
)
Grade Tonnage Curves can also be made and compared between SMU grade and properties at given cutoff for panel. For Grade Tonnage Curves we need to have the rock density in order to be able to compute tonnage values. Density has been extracted for Murchison resource report and is set to 2.4 t/m3.
bm_panel['density'] = 2.4
bm_smu['density'] = 2.4
GTC can be direclty made using grade values and the gtc function on the SMU BlockModel
fig = geo.gtc(bm_smu, 'Fe_LUC', 'density', cutoff_grades, region="channel_hg")
fig
For panel BlockModel, some computation have to made in order to extract tonnage, metal and mean grade at a given cutoff for the entire deposit.
Tonnage above a given cutoff for a block corresponds to :
For the the entire block the tonnage correspond to the sum of this value for each block:
The BlockModel as block of the same size and the density was set to be unique as there is no measurement, this formula can be simplified to
The mean grade, as mentionned above, is computed dividing mass content by ore proportion. So in the same manner, the formula can be simplified to
Cell volume can be extracted using the product of the cell size
cell_volume = np.prod(bm_panel.axis())
We first extract proportion and mass content above cutoff inside list of array using comprehension list
cutoff_prop = [sum(bm_panel.data(prop, region="channel_hg")) for prop in bm_panel.properties() if "ore_prop_above_" in prop]
cutoff_metal = [sum(bm_panel.data(prop, region="channel_hg")) for prop in bm_panel.properties() if "metal_above_" in prop]
Mean grade and tonnage are then computed following determined formulas.
cutoff_mean_grade = [cutoff_metal[i]/cutoff_prop[i] for i in range(len(cutoff_prop))]
tonnage = [cell_volume * bm_panel["density"][0] * ore for ore in cutoff_prop]
We can then update the first GTC obtained on SMU using the computed tonnage and mean grade
fig.add_trace(go.Bar(y=tonnage, x=cutoff_grades, name="LUC Tonnage", marker_color="goldenrod"))
fig.add_trace(go.Scatter(y=cutoff_mean_grade, x = cutoff_grades, name="LUC Mean Grade", mode="lines", line_color="chocolate"), secondary_y=True)
fig.update_layout(
yaxis_title="Tonnage",
yaxis2_title="Mean Grade",
xaxis_title="Cutoff grade",
title = "SMU scale Grade-Tonnage curve",
xaxis_range=[18, 60],
yaxis2_range=[18, 60]
)
fig.show()
Such comparison allows to ensures validity of localisation process as the approximation is quite close the results obtained through the change of support model computation.