Skip to content
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
4 changes: 2 additions & 2 deletions LoopStructural/export/exporters.py
Original file line number Diff line number Diff line change
Expand Up @@ -423,7 +423,7 @@ def _write_vol_evtk(model, file_name, data_label, nsteps, real_coords=True):

"""
# Define grid spacing
xyz = model.bounding_box.regular_grid(nsteps)
xyz = model.bounding_box.regular_grid(nsteps=nsteps)
vals = model.evaluate_model(xyz, scale=False)
if real_coords:
model.rescale(xyz)
Expand Down Expand Up @@ -465,7 +465,7 @@ def _write_vol_gocad(model, file_name, data_label, nsteps, real_coords=True):

"""
# Define grid spacing in model scale coords
xyz = model.bounding_box.regular_grid(nsteps)
xyz = model.bounding_box.regular_grid(nsteps=nsteps)

vals = model.evaluate_model(xyz, scale=False)
# Use FORTRAN style indexing for GOCAD VOXET
Expand Down
3 changes: 3 additions & 0 deletions LoopStructural/interpolators/_surfe_wrapper.py
Original file line number Diff line number Diff line change
Expand Up @@ -205,3 +205,6 @@ def evaluate_gradient(self, evaluation_points):
@property
def nx(self):
return self.get_data_locations().shape[0]
@property
def n_elements(self)->int:
return self.get_data_locations().shape[0]
147 changes: 103 additions & 44 deletions LoopStructural/modelling/core/geological_model.py
Original file line number Diff line number Diff line change
Expand Up @@ -6,7 +6,7 @@

import numpy as np
import pandas as pd
from typing import List
from typing import List, Optional
import pathlib
from ...modelling.features.fault import FaultSegment

Expand Down Expand Up @@ -654,7 +654,9 @@ def set_stratigraphic_column(self, stratigraphic_column, cmap="tab20"):

def create_and_add_foliation(
self,
series_surface_data: str,
series_surface_name: str,
*,
series_surface_data: pd.DataFrame = None,
interpolatortype: str = "FDI",
nelements: int = 1000,
tol=None,
Expand All @@ -664,8 +666,18 @@ def create_and_add_foliation(
"""
Parameters
----------
series_surface_data : string
series_surface_name : string
corresponding to the feature_name in the data
series_surface_data : pd.DataFrame, optional
data frame containing the surface data
interpolatortype : str
the type of interpolator to use, default is 'FDI'
nelements : int
the number of elements to use in the series surface
tol : float, optional
tolerance for the solver, if not specified uses the model default
faults : list, optional
list of faults to be used in the series surface, if not specified uses the model faults
kwargs

Returns
Expand Down Expand Up @@ -696,16 +708,18 @@ def create_and_add_foliation(
bounding_box=self.bounding_box,
interpolatortype=interpolatortype,
nelements=nelements,
name=series_surface_data,
name=series_surface_name,
model=self,
**kwargs,
)
# add data
series_data = self.data[self.data["feature_name"] == series_surface_data]
if series_data.shape[0] == 0:
if series_surface_data is None:
series_surface_data = self.data.loc[self.data["feature_name"] == series_surface_name]

if series_surface_data.shape[0] == 0:
logger.warning("No data for {series_surface_data}, skipping")
return
series_builder.add_data_from_data_frame(series_data)
series_builder.add_data_from_data_frame(series_surface_data)
self._add_faults(series_builder, features=faults)

# build feature
Expand All @@ -721,7 +735,9 @@ def create_and_add_foliation(

def create_and_add_fold_frame(
self,
foldframe_data,
fold_frame_name:str,
*,
fold_frame_data=None,
interpolatortype="FDI",
nelements=1000,
tol=None,
Expand All @@ -731,10 +747,24 @@ def create_and_add_fold_frame(
"""
Parameters
----------
foldframe_data : string
fold_frame_name : string
unique string in feature_name column

kwargs
fold_frame_data : pandas data frame
if not specified uses the model data
interpolatortype : str
the type of interpolator to use, default is 'FDI'
nelements : int
the number of elements to use in the fold frame
tol : float, optional
tolerance for the solver
buffer : float
buffer to add to the bounding box of the fold frame
**kwargs : dict
additional parameters to be passed to the
:class:`LoopStructural.modelling.features.builders.StructuralFrameBuilder`
and :meth:`LoopStructural.modelling.features.builders.StructuralFrameBuilder.setup`
and the interpolator, such as `domain` or `tol`


Returns
-------
Expand All @@ -751,14 +781,18 @@ def create_and_add_fold_frame(
fold_frame_builder = StructuralFrameBuilder(
interpolatortype=interpolatortype,
bounding_box=self.bounding_box.with_buffer(buffer),
name=foldframe_data,
name=fold_frame_name,
frame=FoldFrame,
nelements=nelements,
model=self,
**kwargs,
)
# add data
fold_frame_data = self.data[self.data["feature_name"] == foldframe_data]
if fold_frame_data is None:
fold_frame_data = self.data.loc[self.data["feature_name"] == fold_frame_name]
if fold_frame_data.shape[0] == 0:
logger.warning(f"No data for {fold_frame_name}, skipping")
return
fold_frame_builder.add_data_from_data_frame(fold_frame_data)
self._add_faults(fold_frame_builder[0])
self._add_faults(fold_frame_builder[1])
Expand All @@ -775,7 +809,9 @@ def create_and_add_fold_frame(

def create_and_add_folded_foliation(
self,
foliation_data,
foliation_name,
*,
foliation_data=None,
interpolatortype="DFI",
nelements=10000,
buffer=0.1,
Expand Down Expand Up @@ -832,15 +868,18 @@ def create_and_add_folded_foliation(
bounding_box=self.bounding_box.with_buffer(buffer),
nelements=nelements,
fold=fold,
name=foliation_data,
name=foliation_name,
svario=svario,
model=self,
**kwargs,
)

if foliation_data is None:
foliation_data = self.data.loc[self.data["feature_name"] == foliation_name]
if foliation_data.shape[0] == 0:
logger.warning(f"No data for {foliation_name}, skipping")
return
series_builder.add_data_from_data_frame(
self.data[self.data["feature_name"] == foliation_data]
)
foliation_data )
self._add_faults(series_builder)
# series_builder.add_data_to_interpolator(True)
# build feature
Expand All @@ -858,7 +897,9 @@ def create_and_add_folded_foliation(

def create_and_add_folded_fold_frame(
self,
fold_frame_data,
fold_frame_name: str,
*,
fold_frame_data: Optional[pd.DataFrame] = None,
interpolatortype="FDI",
nelements=10000,
fold_frame=None,
Expand All @@ -869,14 +910,22 @@ def create_and_add_folded_fold_frame(

Parameters
----------
fold_frame_data : string
fold_frame_name : string
name of the feature to be added

fold_frame_data : pandas data frame, optional
data frame containing the fold frame data, if not specified uses the model data
interpolatortype : str
the type of interpolator to use, default is 'FDI' (unused) 5/6/2025
fold_frame : StructuralFrame, optional
the fold frame for the fold if not specified uses last feature added

kwargs : dict
parameters passed to child functions
nelements : int
the number of elements to use in the fold frame
tol : float, optional
tolerance for the solver, if not specified uses the model default
**kwargs : dict
additional parameters to be passed to the
:class:`LoopStructural.modelling.features.builders.StructuralFrameBuilder`
and :meth:`LoopStructural.modelling.features.builders.StructuralFrameBuilder.setup`

Returns
-------
Expand Down Expand Up @@ -917,15 +966,15 @@ def create_and_add_folded_fold_frame(
interpolatortype=interpolatortypes,
bounding_box=self.bounding_box.with_buffer(kwargs.get("buffer", 0.1)),
nelements=[nelements, nelements, nelements],
name=fold_frame_data,
name=fold_frame_name,
fold=fold,
frame=FoldFrame,
model=self,
**kwargs,
)
fold_frame_builder.add_data_from_data_frame(
self.data[self.data["feature_name"] == fold_frame_data]
)
if fold_frame_data is None:
fold_frame_data = self.data[self.data["feature_name"] == fold_frame_name]
fold_frame_builder.add_data_from_data_frame(fold_frame_data)

for i in range(3):
self._add_faults(fold_frame_builder[i])
Expand All @@ -947,6 +996,7 @@ def create_and_add_intrusion(
self,
intrusion_name,
intrusion_frame_name,
*,
intrusion_frame_parameters={},
intrusion_lateral_extent_model=None,
intrusion_vertical_extent_model=None,
Expand Down Expand Up @@ -1224,7 +1274,7 @@ def add_onlap_unconformity(self, feature: GeologicalFeature, value: float) -> Ge
return uc_feature

def create_and_add_domain_fault(
self, fault_surface_data, nelements=10000, interpolatortype="FDI", **kwargs
self, fault_surface_data,*, nelements=10000, interpolatortype="FDI", **kwargs
):
"""
Parameters
Expand Down Expand Up @@ -1252,7 +1302,7 @@ def create_and_add_domain_fault(
)

# add data
unconformity_data = self.data[self.data["feature_name"] == fault_surface_data]
unconformity_data = self.data.loc[self.data["feature_name"] == fault_surface_data]

domain_fault_feature_builder.add_data_from_data_frame(unconformity_data)
# look through existing features if there is a fault before an
Expand All @@ -1275,8 +1325,10 @@ def create_and_add_domain_fault(

def create_and_add_fault(
self,
fault_surface_data,
displacement,
fault_name: str,
displacement: float,
*,
fault_data:Optional[pd.DataFrame] = None,
interpolatortype="FDI",
tol=None,
fault_slip_vector=None,
Expand All @@ -1299,9 +1351,12 @@ def create_and_add_fault(
"""
Parameters
----------
fault_surface_data : string
fault_name : string
name of the fault surface data in the dataframe
displacement : displacement magnitude
displacement magnitude of the fault, in model units
fault_data : pd.DataFrame, optional
data frame containing the fault data, if not specified uses the model data
major_axis : [type], optional
[description], by default None
minor_axis : [type], optional
Expand All @@ -1328,7 +1383,7 @@ def create_and_add_fault(
if "fault_vectical_radius" in kwargs and intermediate_axis is None:
intermediate_axis = kwargs["fault_vectical_radius"]

logger.info(f'Creating fault "{fault_surface_data}"')
logger.info(f'Creating fault "{fault_name}"')
logger.info(f"Displacement: {displacement}")
logger.info(f"Tolerance: {tol}")
logger.info(f"Fault function: {faultfunction}")
Expand All @@ -1353,7 +1408,7 @@ def create_and_add_fault(
# tol *= 0.1*minor_axis

if displacement == 0:
logger.warning(f"{fault_surface_data} displacement is 0")
logger.warning(f"{fault_name} displacement is 0")

if "data_region" in kwargs:
kwargs.pop("data_region")
Expand All @@ -1363,14 +1418,18 @@ def create_and_add_fault(
interpolatortype,
bounding_box=self.bounding_box,
nelements=kwargs.pop("nelements", 1e4),
name=fault_surface_data,
name=fault_name,
model=self,
**kwargs,
)
fault_frame_data = self.data.loc[self.data["feature_name"] == fault_surface_data].copy()
if fault_data is None:
fault_data = self.data.loc[self.data["feature_name"] == fault_name]
if fault_data.shape[0] == 0:
logger.warning(f"No data for {fault_name}, skipping")
return

self._add_faults(fault_frame_builder, features=faults)
# add data
fault_frame_data = self.data.loc[self.data["feature_name"] == fault_surface_data].copy()

if fault_center is not None and ~np.isnan(fault_center).any():
fault_center = self.scale(fault_center, inplace=False)
Expand All @@ -1381,7 +1440,7 @@ def create_and_add_fault(
if intermediate_axis:
intermediate_axis = intermediate_axis / self.scale_factor
fault_frame_builder.create_data_from_geometry(
fault_frame_data=fault_frame_data,
fault_frame_data=fault_data,
fault_center=fault_center,
fault_normal_vector=fault_normal_vector,
fault_slip_vector=fault_slip_vector,
Expand Down Expand Up @@ -1418,7 +1477,7 @@ def create_and_add_fault(
return fault

# TODO move rescale to bounding box/transformer
def rescale(self, points: np.ndarray, inplace: bool = False) -> np.ndarray:
def rescale(self, points: np.ndarray, *, inplace: bool = False) -> np.ndarray:
"""
Convert from model scale to real world scale - in the future this
should also do transformations?
Expand All @@ -1441,7 +1500,7 @@ def rescale(self, points: np.ndarray, inplace: bool = False) -> np.ndarray:
return points

# TODO move scale to bounding box/transformer
def scale(self, points: np.ndarray, inplace: bool = False) -> np.ndarray:
def scale(self, points: np.ndarray, *, inplace: bool = False) -> np.ndarray:
"""Take points in UTM coordinates and reproject
into scaled model space

Expand All @@ -1467,7 +1526,7 @@ def scale(self, points: np.ndarray, inplace: bool = False) -> np.ndarray:
points /= self.scale_factor
return points

def regular_grid(self, nsteps=None, shuffle=True, rescale=False, order="C"):
def regular_grid(self, *, nsteps=None, shuffle=True, rescale=False, order="C"):
"""
Return a regular grid within the model bounding box

Expand All @@ -1483,7 +1542,7 @@ def regular_grid(self, nsteps=None, shuffle=True, rescale=False, order="C"):
"""
return self.bounding_box.regular_grid(nsteps=nsteps, shuffle=shuffle, order=order)

def evaluate_model(self, xyz: np.ndarray, scale: bool = True) -> np.ndarray:
def evaluate_model(self, xyz: np.ndarray, *, scale: bool = True) -> np.ndarray:
"""Evaluate the stratigraphic id at each location

Parameters
Expand Down Expand Up @@ -1559,7 +1618,7 @@ def evaluate_model(self, xyz: np.ndarray, scale: bool = True) -> np.ndarray:
logger.error(f"Model does not contain {group}")
return strat_id

def evaluate_model_gradient(self, points: np.ndarray, scale: bool = True) -> np.ndarray:
def evaluate_model_gradient(self, points: np.ndarray, *, scale: bool = True) -> np.ndarray:
"""Evaluate the gradient of the stratigraphic column at each location

Parameters
Expand Down
Loading
Loading