diff --git a/LoopStructural/export/exporters.py b/LoopStructural/export/exporters.py index bf0c30a0..a8eb3093 100644 --- a/LoopStructural/export/exporters.py +++ b/LoopStructural/export/exporters.py @@ -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) @@ -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 diff --git a/LoopStructural/interpolators/_surfe_wrapper.py b/LoopStructural/interpolators/_surfe_wrapper.py index 6ebc1bcf..4784434a 100644 --- a/LoopStructural/interpolators/_surfe_wrapper.py +++ b/LoopStructural/interpolators/_surfe_wrapper.py @@ -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] \ No newline at end of file diff --git a/LoopStructural/modelling/core/geological_model.py b/LoopStructural/modelling/core/geological_model.py index 72233b09..bb6d425a 100644 --- a/LoopStructural/modelling/core/geological_model.py +++ b/LoopStructural/modelling/core/geological_model.py @@ -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 @@ -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, @@ -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 @@ -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 @@ -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, @@ -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 ------- @@ -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]) @@ -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, @@ -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 @@ -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, @@ -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 ------- @@ -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]) @@ -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, @@ -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 @@ -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 @@ -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, @@ -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 @@ -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}") @@ -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") @@ -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) @@ -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, @@ -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? @@ -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 @@ -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 @@ -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 @@ -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 diff --git a/LoopStructural/modelling/features/_base_geological_feature.py b/LoopStructural/modelling/features/_base_geological_feature.py index f1217f0b..e49291a1 100644 --- a/LoopStructural/modelling/features/_base_geological_feature.py +++ b/LoopStructural/modelling/features/_base_geological_feature.py @@ -244,7 +244,7 @@ def min(self): if self.model is None: return 0 - return np.nanmin(self.evaluate_value(self.model.regular_grid((10, 10, 10)))) + return np.nanmin(self.evaluate_value(self.model.regular_grid(nsteps=(10, 10, 10)))) def max(self): """Calculate the maximum value of the geological feature @@ -257,7 +257,7 @@ def max(self): """ if self.model is None: return 0 - return np.nanmax(self.evaluate_value(self.model.regular_grid((10, 10, 10)))) + return np.nanmax(self.evaluate_value(self.model.regular_grid(nsteps=(10, 10, 10)))) def __tojson__(self): regions = [r.name for r in self.regions] diff --git a/LoopStructural/modelling/intrusions/intrusion_builder.py b/LoopStructural/modelling/intrusions/intrusion_builder.py index 745048a5..774f798f 100644 --- a/LoopStructural/modelling/intrusions/intrusion_builder.py +++ b/LoopStructural/modelling/intrusions/intrusion_builder.py @@ -82,7 +82,7 @@ def create_grid_for_evaluation(self, spacing=None): if spacing is None: spacing = self.model.nsteps - grid_points = self.model.regular_grid(spacing, shuffle=False) + grid_points = self.model.regular_grid(nsteps=spacing, shuffle=False) grid_points_coord0 = self.intrusion_frame[0].evaluate_value(grid_points) diff --git a/LoopStructural/modelling/intrusions/intrusion_frame_builder.py b/LoopStructural/modelling/intrusions/intrusion_frame_builder.py index 9d8ae8ff..3ce29068 100644 --- a/LoopStructural/modelling/intrusions/intrusion_frame_builder.py +++ b/LoopStructural/modelling/intrusions/intrusion_frame_builder.py @@ -168,7 +168,7 @@ def create_grid_for_indicator_fxs(self, spacing=None): if spacing is None: spacing = self.model.nsteps - grid_points = self.model.regular_grid(spacing, shuffle=False) + grid_points = self.model.regular_grid(nsteps=spacing, shuffle=False) self.grid_to_evaluate_ifx = grid_points diff --git a/tests/integration/test_fold_models.py b/tests/integration/test_fold_models.py index ece4d922..178944b4 100644 --- a/tests/integration/test_fold_models.py +++ b/tests/integration/test_fold_models.py @@ -16,7 +16,7 @@ def test_average_fold_axis(): fold_frame = model.create_and_add_fold_frame("s1", nelements=10000) stratigraphy = model.create_and_add_folded_foliation( "s0", - fold_frame, + fold_frame=fold_frame, nelements=10000, av_fold_axis=True, # fold_axis=[-6.51626577e-06, -5.00013645e-01, -8.66017526e-01], @@ -39,7 +39,7 @@ def test_fixed_fold_axis(): fold_frame = model.create_and_add_fold_frame("s1", nelements=10000) stratigraphy = model.create_and_add_folded_foliation( "s0", - fold_frame, + fold_frame=fold_frame, nelements=10000, # av_fold_axis=True fold_axis=[-6.51626577e-06, -5.00013645e-01, -8.66017526e-01], @@ -58,7 +58,7 @@ def test_fixed_wavelength(): fold_frame = model.create_and_add_fold_frame("s1", nelements=10000) stratigraphy = model.create_and_add_folded_foliation( "s0", - fold_frame, + fold_frame=fold_frame, nelements=10000, # av_fold_axis=True fold_axis=[-6.51626577e-06, -5.00013645e-01, -8.66017526e-01],