geost.models.model_dataset.ModelDataset.get_thickness

ModelDataset.get_thickness(condition: GeotopUnits | DataArray) → DataArray

Calculate the thickness where a voxelmodel or layermodel meets a specified condition. The condition can be a boolean DataArray or Dataset that indicates which voxels or layers to include in the thickness calculation.

Parameters:

condition (GeotopUnits | xr.DataArray) – A boolean DataArray or Dataset indicating which voxels or layers to include in the thickness calculation. The condition should have the same dimensions as the model.

Returns:

xarray.DataArray containing the calculated thickness for each horizontal x,y- location in the model.

Return type:

xr.DataArray

Examples

Determine the thickness in a voxelmodel Dataset where “strat” equals 1100 and “lith” equals 1:

>>> thickness = voxelmodel.gst.get_thickness(
...    (voxelmodel["lithology"] == 1) & (voxelmodel["strat"] == 1100)
... )

Or in a layermodel for a subset of units where the a value is smaller than 10:

>>> thickness = layermodel.gst.get_thickness(
...    (layermodel["layer"].isin(["B", "D"])) & (layermodel["value"] < 10)
... )

If you are working with GeoTOP, you can also use a GeotopUnits object to specify the condition, for example to get the thickness of the “Formatie van Echteld” unit with clay as the lithology:

>>> geotop = geost.read_geotop_from_opendap(bbox=(110_000, 440_000, 120_000, 450_000))
>>> strat_units = geost.bro.geotop_strat_units()
>>> lithok_units = geost.bro.geotop_lithok_units()

>>> thickness = geotop.gst.get_thickness(
...     strat_units.select_description_contains("Formatie van Echteld")
...     & lithok_units.select_unit("k")
... )