From c837d6c1ea8ae232600a1f737265db9122ef3a83 Mon Sep 17 00:00:00 2001 From: Lachlan Grose Date: Tue, 17 Jun 2025 09:19:34 +1000 Subject: [PATCH] fix: allow data to be specified in create and add function. --- .../modelling/core/geological_model.py | 92 +++++++++---------- .../modelling/intrusions/test_intrusions.py | 6 +- tests/unit/modelling/test_geological_model.py | 2 +- 3 files changed, 50 insertions(+), 50 deletions(-) diff --git a/LoopStructural/modelling/core/geological_model.py b/LoopStructural/modelling/core/geological_model.py index e8a10a8f..d098fa70 100644 --- a/LoopStructural/modelling/core/geological_model.py +++ b/LoopStructural/modelling/core/geological_model.py @@ -125,9 +125,7 @@ def __init__( self.features = [] self.feature_name_index = {} self._data = pd.DataFrame() # None - - self.stratigraphic_column = None self.tol = 1e-10 * np.max(self.bounding_box.maximum - self.bounding_box.origin) @@ -179,8 +177,40 @@ def __str__(self): def _ipython_key_completions_(self): return self.feature_name_index.keys() - + def prepare_data(self, data: pd.DataFrame) -> pd.DataFrame: + data = data.copy() + data[['X', 'Y', 'Z']] = self.bounding_box.project(data[['X', 'Y', 'Z']].to_numpy()) + if "type" in data: + logger.warning("'type' is deprecated replace with 'feature_name' \n") + data.rename(columns={"type": "feature_name"}, inplace=True) + if "feature_name" not in data: + logger.error("Data does not contain 'feature_name' column") + raise BaseException("Cannot load data") + for h in all_heading(): + if h not in data: + data[h] = np.nan + if h == "w": + data[h] = 1.0 + if h == "coord": + data[h] = 0 + if h == "polarity": + data[h] = 1.0 + # LS wants polarity as -1 or 1, change 0 to -1 + data.loc[data["polarity"] == 0, "polarity"] = -1.0 + data.loc[np.isnan(data["w"]), "w"] = 1.0 + if "strike" in data and "dip" in data: + logger.info("Converting strike and dip to vectors") + mask = np.all(~np.isnan(data.loc[:, ["strike", "dip"]]), axis=1) + data.loc[mask, gradient_vec_names()] = ( + strikedip2vector(data.loc[mask, "strike"], data.loc[mask, "dip"]) + * data.loc[mask, "polarity"].to_numpy()[:, None] + ) + data.drop(["strike", "dip"], axis=1, inplace=True) + data[['X', 'Y', 'Z', 'val', 'nx', 'ny', 'nz', 'gx', 'gy', 'gz', 'tx', 'ty', 'tz']] = data[ + ['X', 'Y', 'Z', 'val', 'nx', 'ny', 'nz', 'gx', 'gy', 'gz', 'tx', 'ty', 'tz'] + ].astype(float) + return data @classmethod def from_processor(cls, processor): """Builds a model from a :class:`LoopStructural.modelling.input.ProcessInputData` object @@ -473,40 +503,9 @@ def data(self, data: pd.DataFrame): raise BaseException("Cannot load data") logger.info(f"Adding data to GeologicalModel with {len(data)} data points") self._data = data.copy() - self._data[['X','Y','Z']] = self.bounding_box.project(self._data[['X','Y','Z']].to_numpy()) - - - if "type" in self._data: - logger.warning("'type' is deprecated replace with 'feature_name' \n") - self._data.rename(columns={"type": "feature_name"}, inplace=True) - if "feature_name" not in self._data: - logger.error("Data does not contain 'feature_name' column") - raise BaseException("Cannot load data") - for h in all_heading(): - if h not in self._data: - self._data[h] = np.nan - if h == "w": - self._data[h] = 1.0 - if h == "coord": - self._data[h] = 0 - if h == "polarity": - self._data[h] = 1.0 - # LS wants polarity as -1 or 1, change 0 to -1 - self._data.loc[self._data["polarity"] == 0, "polarity"] = -1.0 - self._data.loc[np.isnan(self._data["w"]), "w"] = 1.0 - if "strike" in self._data and "dip" in self._data: - logger.info("Converting strike and dip to vectors") - mask = np.all(~np.isnan(self._data.loc[:, ["strike", "dip"]]), axis=1) - self._data.loc[mask, gradient_vec_names()] = ( - strikedip2vector(self._data.loc[mask, "strike"], self._data.loc[mask, "dip"]) - * self._data.loc[mask, "polarity"].to_numpy()[:, None] - ) - self._data.drop(["strike", "dip"], axis=1, inplace=True) - self._data[['X', 'Y', 'Z', 'val', 'nx', 'ny', 'nz', 'gx', 'gy', 'gz', 'tx', 'ty', 'tz']] = ( - self._data[ - ['X', 'Y', 'Z', 'val', 'nx', 'ny', 'nz', 'gx', 'gy', 'gz', 'tx', 'ty', 'tz'] - ].astype(float) - ) + # self._data[['X','Y','Z']] = self.bounding_box.project(self._data[['X','Y','Z']].to_numpy()) + + def set_model_data(self, data): logger.warning("deprecated method. Model data can now be set using the data attribute") @@ -623,7 +622,7 @@ def create_and_add_foliation( 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_surface_data) + series_builder.add_data_from_data_frame(self.prepare_data(series_surface_data)) self._add_faults(series_builder, features=faults) # build feature @@ -697,7 +696,7 @@ def create_and_add_fold_frame( 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) + fold_frame_builder.add_data_from_data_frame(self.prepare_data(fold_frame_data)) self._add_faults(fold_frame_builder[0]) self._add_faults(fold_frame_builder[1]) self._add_faults(fold_frame_builder[2]) @@ -783,7 +782,10 @@ def create_and_add_folded_foliation( logger.warning(f"No data for {foliation_name}, skipping") return series_builder.add_data_from_data_frame( -foliation_data ) + self.prepare_data( + foliation_data + ) + ) self._add_faults(series_builder) # series_builder.add_data_to_interpolator(True) # build feature @@ -878,7 +880,7 @@ def create_and_add_folded_fold_frame( ) 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) + fold_frame_builder.add_data_from_data_frame(self.prepare_data(fold_frame_data)) for i in range(3): self._add_faults(fold_frame_builder[i]) @@ -1331,7 +1333,7 @@ def create_and_add_fault( 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 @@ -1344,7 +1346,7 @@ def create_and_add_fault( if intermediate_axis: intermediate_axis = intermediate_axis fault_frame_builder.create_data_from_geometry( - fault_frame_data=fault_data, + fault_frame_data=self.prepare_data(fault_data), fault_center=fault_center, fault_normal_vector=fault_normal_vector, fault_slip_vector=fault_slip_vector, @@ -1397,9 +1399,8 @@ def rescale(self, points: np.ndarray, *, inplace: bool = False) -> np.ndarray: points : np.array((N,3),dtype=double) """ - + return self.bounding_box.reproject(points,inplace=inplace) - # TODO move scale to bounding box/transformer def scale(self, points: np.ndarray, *, inplace: bool = False) -> np.ndarray: @@ -1418,7 +1419,6 @@ def scale(self, points: np.ndarray, *, inplace: bool = False) -> np.ndarray: """ return self.bounding_box.project(np.array(points).astype(float),inplace=inplace) - def regular_grid(self, *, nsteps=None, shuffle=True, rescale=False, order="C"): """ diff --git a/tests/unit/modelling/intrusions/test_intrusions.py b/tests/unit/modelling/intrusions/test_intrusions.py index 377ccc2e..0f075125 100644 --- a/tests/unit/modelling/intrusions/test_intrusions.py +++ b/tests/unit/modelling/intrusions/test_intrusions.py @@ -65,9 +65,9 @@ def test_intrusion_builder(): model.data = data model.nsteps = [10, 10, 10] - intrusion_data = data[data["feature_name"] == "tabular_intrusion"] - intrusion_frame_data = model.data[model.data["feature_name"] == "tabular_intrusion_frame"] - + intrusion_data = model.prepare_data(data[data["feature_name"] == "tabular_intrusion"]) + intrusion_frame_data = model.prepare_data(model.data[model.data["feature_name"] == "tabular_intrusion_frame"]) + conformable_feature = model.create_and_add_foliation("stratigraphy") intrusion_frame_parameters = { diff --git a/tests/unit/modelling/test_geological_model.py b/tests/unit/modelling/test_geological_model.py index c09405e9..09e4cf97 100644 --- a/tests/unit/modelling/test_geological_model.py +++ b/tests/unit/modelling/test_geological_model.py @@ -17,7 +17,7 @@ def test_rescale_model_data(): model.set_model_data(data) # Check that the model data is rescaled to local coordinates expected = data[['X', 'Y', 'Z']].values - bb[None, 0, :] - actual = model.data[['X', 'Y', 'Z']].values + actual = model.prepare_data(model.data)[['X', 'Y', 'Z']].values assert np.allclose(actual, expected, atol=1e-6) def test_access_feature_model(): data, bb = load_claudius()