LoopStructural

Overview

LoopStructural is an open source 3D modelling library providing access to multiple interpolation schemes with a high level and easy to use API for creating geological models. The library has been written for the Loop platform by Lachlan Grose at Monash University.

Loop is an open source 3D probabilistic geological and geophysical modelling platform, initiated by Geoscience Australia and the OneGeology consortium. The project is funded by Australian territory, State and Federal Geological Surveys, the Australian Research Council and the MinEx Collaborative Research Centre.

Examples

LoopStructural can be used to build a 3D geological model using geological relationships between geological objects e.g. faults, folds, unconformities and stratigraphic contacts. The library also provides a high level API to access the fast interpolation schemes that are used by LoopStructural.

Using GeologicalModel

The following example shows how to use the geological model interface to create a geological model from a dataset and evaluate the scalar field and gradient of the interpolator at some random locations.

from LoopStructural import GeologicalModel
from LoopStructural.datatypes import BoundingBox
from LoopStructural.visualisation import Loop3DView
from LoopStructural.datasets import load_claudius

import numpy as np
data, bb = load_claudius()

#bb constaints origin and maximum of axis aligned bounding box
#data is a pandas dataframe with X,Y,Z,val,nx,ny,nz, feature_name

model = GeologicalModel(bb[0,:],bb[1,:])
model.data = data
# nelements specifies the number of discrete interpolation elements
# 'stratí' is the feature name in the data dataframe
model.create_and_add_foliation('strati',nelements=1e5)
model.update()
# get the value of the interpolator at some random locations
locations = np.array(
    [
        np.random.uniform(bb[0, 0], bb[1, 0],5),
        np.random.uniform(bb[0, 1], bb[1, 1],5),
        np.random.uniform(bb[0, 2], bb[1, 2],5),
    ]
).T
val = model.evaluate_feature_value('strati', locations)
# get the gradient of the interpolator
gradient = model.evaluate_feature_gradient('strati',locations)

#Plot the scalar field of the model
model['strati'].scalar_field().plot()

Using InterpolatorBuilder

To access the interpolators directly the InterpolatorBuilder can be used to help assemble an interpolator from a combination of value, gradient and inequality constraints.

from LoopStructural import BoundingBox, InterpolatorBuilder
from LoopStructural.utils import EuclideanTransformation
from LoopStructural.datasets import load_claudius

# load in a dataframe with x,y,z,val,nx,ny,nz to constrain the value and
# gradient normal of the interpolator
data, bb = load_claudius()

# Find the bounding box of the data by finding the extent
bounding_box = BoundingBox().fit(data[["X", "Y", "Z"]])
# assemble an interpolator using the bounding box and the data
interpolator = (
    InterpolatorBuilder(
        interpolatortype="FDI", bounding_box=bounding_box.with_buffer(0.1)
    )
    .add_value_constraints(data.loc[data["val"].notna(), ["X", "Y", "Z", "val"]].values)
    .add_normal_constraints(
        data.loc[data["nx"].notna(), ["X", "Y", "Z", "nx", "ny", "nz"]].values
    )
    .build()
)
# Set the number of elements in the bounding box to 10000 and create a structured grid
bounding_box.nelements = 10000
mesh = bounding_box.structured_grid()
# add the interpolated values to the mesh at the nodes
mesh.properties["v"] = interpolator.evaluate_value(mesh.nodes)

# or the cell centres
# mesh.cell_properties["v"] = interpolator.evaluate_value(mesh.cell_centres)


# visualise the scalar value

mesh.plot()

# We can also add gradient properties and visualise these
mesh.properties["grad"] = interpolator.evaluate_gradient(mesh.nodes)
mesh.vtk().glyph(orient="grad").plot()