fix: world locations for display, export and intrusion classification; deprecate rescale and scale, stage S4 (spec 005)

API.md

@@ -119,7 +119,9 @@ the same thing it used to.
   `evaluate_feature_value`, `evaluate_feature_gradient`,
   `evaluate_fault_displacements`
 - Solve: `update`
-- Geometry: `scale`, `rescale`, `regular_grid`, `bounding_box` (property)
+- Geometry: `scale` (deprecated, use `bounding_box.project`), `rescale`
+  (deprecated, use `bounding_box.reproject`), `regular_grid`,
+  `bounding_box` (property)
 - Stratigraphy: `stratigraphic_column` (property), `stratigraphic_ids`
 - Output: `get_fault_surfaces`, `get_stratigraphic_surfaces`,
   `get_block_model`, `save`, `to_file`, `to_dict`

COMPAT.md

@@ -14,6 +14,8 @@ any module-path move/rename on the 1.x line gets a re-export shim with a
 | `LoopStructural.interpolators._interpolator_builder` | `loop_interpolation._interpolator_builder` | 2026-07-29 | 2026-07-29 | extracted-package integration | 2 minor releases after first shim release | Active | Core maintainers |
 | `LoopStructural.interpolators._interpolator_factory` | `loop_interpolation._interpolator_factory` | 2026-07-29 | 2026-07-29 | extracted-package integration | 2 minor releases after first shim release | Active | Core maintainers |
 | Loading a legacy interpolator payload with no choice of the caller (the loader assumes an identity frame, now with a `FutureWarning`) | Pass `local_origin` and `local_rotation`, or `declare_identity=True`, to `from_dict`, `from_json` or `from_yaml`; at the next major release a call with no choice raises `InterpolatorFormatError` | 2026-11-05 | TBD (next release containing this warning) | spec 005 stage S3 (pull request #397) | The next major release | Active | Core maintainers |
+| `GeologicalModel.rescale(points, inplace=False)` | `GeologicalModel.bounding_box.reproject(points, inplace=False)` | 2026-11-05 | TBD (next release containing this warning) | spec 005 stage S4 | TBD (2 minor releases after the first release that warns) | Active | Core maintainers |
+| `GeologicalModel.scale(points, inplace=False)` | `GeologicalModel.bounding_box.project(points, inplace=False)` | 2026-11-05 | TBD (next release containing this warning) | spec 005 stage S4 | TBD (2 minor releases after the first release that warns) | Active | Core maintainers |
 
 ## Deprecation lifecycle
 
@@ -51,11 +53,12 @@ unchanged, but the new path is preferred going forward.
 | `create_and_add_fault(fault_center=...)` projected the centre into local coordinates; `fault_buffer`/`force_mesh_geometry` were replaced by a default 0.2 buffer; an inferred centre had z = 0; any 3x3 matrix was accepted as a rotation | `fault_center` is a world point; the requested buffer and forced geometry are kept; an inferred centre keeps its elevation; `set_local_transform` and `rotation_xy` raise `RigidRotationError` for a non-rigid matrix. Extra observations added to work around these defects are no longer necessary (spec 005) | 2026-11-04 |
 | Fold `mask_fn`, custom fold providers (`get_deformed_orientation`), `directional_regularisation` direction callables and arrays, and the P2 `steepness_wtfunc` received or returned local (computational-frame) points and directions | They receive world points and give world directions; the interpolator converts to its local frame once. With an identity frame nothing changes (spec 005, stage S2) | 2026-11-04 |
 | `copy()` returned an interpolator with no frame; a transform change after constraints existed gave a stale result | New payloads carry `format_version: 2`, `constraint_space`, `transform` and `bounding_box`. `copy()` keeps the frame. A transform change after constraints exist raises `StaleFrameError` at the next set, evaluate or save (spec 005, stage S3). Loading a legacy payload with no choice is deprecated: see the deprecation table | 2026-11-05 |
+| `GeologicalFeature.get_data()` (used by `Loop2DView.add_data` and `ModelExporter.save`) returned locations moved by the local origin and rotated a second time; the intrusion classification evaluated the model at such points | Both use the stored world locations. With an identity frame nothing changes (spec 005, stage S4) | 2026-11-05 |
 
 The following stable methods remain part of the documented public API surface after the package extraction work and are covered by the snapshot-based contract test: add_onlap_unconformity, add_unconformity, create_and_add_domain_fault, create_and_add_fault, create_and_add_fold_frame, create_and_add_folded_fold_frame, create_and_add_folded_foliation, create_and_add_foliation, evaluate_model, evaluate_model_gradient, from_file, get_fault_surfaces, get_feature_by_name, get_stratigraphic_surfaces, rescale, save, scale.
 
 ## Compatibility debt summary
 
-- Active shims: 8
+- Active shims: 10
 - Oldest active shim added: 2026-07-24
 - Next cleanup milestone: remove the interpolator-path shims after 2 minor releases from 2026-07-29

LoopStructural/export/exporters.py

@@ -89,7 +89,6 @@ def write_feat_surfs(
                 model.bounding_box.bb[1, 2],
             ]
         )
-        model.rescale(verts)
         surf = Surface(
             vertices=verts,
             triangles=faces,
@@ -362,8 +361,6 @@ def _write_cubeface_evtk(model, file_name, data_label, nsteps, real_coords=True)
     # Evaluate model at points
     points, tri = create_box(model.bounding_box, nsteps)
     val = model.evaluate_model(points, scale=False)
-    if real_coords:
-        model.rescale(points)
 
     # Define vertices
     x = np.zeros(points.shape[0])
@@ -427,8 +424,6 @@ def _write_vol_evtk(model, file_name, data_label, nsteps, real_coords=True):
     # project world -> local internally)
     xyz = model.bounding_box.regular_grid(nsteps=nsteps, local=False)
     vals = model.evaluate_model(xyz, scale=False)
-    if real_coords:
-        model.rescale(xyz)
 
     # Define vertices - xyz.shape[0] is length of vector array
     x = np.zeros(xyz.shape[0])
@@ -475,10 +470,6 @@ def _write_vol_gocad(model, file_name, data_label, nsteps, real_coords=True):
     vol_vals = np.reshape(vals, nsteps, order="F")
     bbox = model.bounding_box[:]
 
-    # Convert bounding box to real world scale coords
-    if real_coords:
-        model.rescale(bbox)
-
     # If integer values
     if type(vals[0]) is np.int64:
         d_type = np.int8

LoopStructural/modelling/core/_model_exporter.py

@@ -74,7 +74,7 @@ def get_block_model(
         grid = model.bounding_box.structured_grid(name=name)
 
         grid.cell_properties['stratigraphy'] = model.evaluate_model(
-            model.rescale(model.bounding_box.cell_centres())
+            model.bounding_box.reproject(model.bounding_box.cell_centres())
         )
         return grid, model.stratigraphic_ids()
 

LoopStructural/modelling/core/geological_model.py

@@ -1107,6 +1107,10 @@ 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?
 
+        .. deprecated::
+            Use ``model.bounding_box.reproject`` instead.
+            Deprecated by spec 005, stage S4.
+
         Parameters
         ----------
         points : np.array((N,3),dtype=double)
@@ -1118,7 +1122,12 @@ def rescale(self, points: np.ndarray, *, inplace: bool = False) -> np.ndarray:
         points : np.array((N,3),dtype=double)
 
         """
-
+        # stacklevel 3 skips the public_api wrapper frame
+        warnings.warn(
+            "GeologicalModel.rescale is deprecated and will be removed in a future release; use model.bounding_box.reproject(points) instead",
+            DeprecationWarning,
+            stacklevel=3,
+        )
         return self.bounding_box.reproject(points, inplace=inplace)
 
     # TODO move scale to bounding box/transformer
@@ -1127,6 +1136,10 @@ def scale(self, points: np.ndarray, *, inplace: bool = False) -> np.ndarray:
         """Take points in UTM coordinates and reproject
         into scaled model space
 
+        .. deprecated::
+            Use ``model.bounding_box.project`` instead.
+            Deprecated by spec 005, stage S4.
+
         Parameters
         ----------
         points : np.array((N,3),dtype=float)
@@ -1138,6 +1151,12 @@ def scale(self, points: np.ndarray, *, inplace: bool = False) -> np.ndarray:
         points : np.a::rray((N,3),dtype=double)
 
         """
+        # stacklevel 3 skips the public_api wrapper frame
+        warnings.warn(
+            "GeologicalModel.scale is deprecated and will be removed in a future release; use model.bounding_box.project(points) instead",
+            DeprecationWarning,
+            stacklevel=3,
+        )
         return self.bounding_box.project(np.array(points).astype(float), inplace=inplace)
 
     @public_api(tier="stable")

LoopStructural/modelling/features/_geological_feature.py

@@ -303,6 +303,8 @@ def copy(self, name=None):
     def get_data(self, value_map: dict | None = None) -> list[ValuePoints | VectorPoints]:
         """Return the data associated with this geological feature
 
+        The locations and the vectors are in world coordinates, as the builder stores them.
+
         Parameters
         ----------
         value_map : Optional[dict], optional
@@ -326,15 +328,15 @@ def get_data(self, value_map: dict | None = None) -> list[ValuePoints | VectorPo
 
             data.append(
                 VectorPoints(
-                    locations=self.model.rescale(gradient_constraints[:, :3]),
+                    locations=gradient_constraints[:, :3],
                     vectors=gradient_constraints[:, 3:6],
                     name=f'{self.name}+_gradient',
                 )
             )
         if norm_constraints.shape[0] > 0:
             data.append(
                 VectorPoints(
-                    locations=self.model.rescale(norm_constraints[:, :3]),
+                    locations=norm_constraints[:, :3],
                     vectors=norm_constraints[:, 3:6],
                     name=f'{self.name}_norm',
                 )
@@ -344,17 +346,15 @@ def get_data(self, value_map: dict | None = None) -> list[ValuePoints | VectorPo
                 for name, v in value_map.items():
                     data.append(
                         ValuePoints(
-                            locations=self.model.rescale(
-                                value_constraints[value_constraints == v, :3]
-                            ),
+                            locations=value_constraints[value_constraints == v, :3],
                             values=value_constraints[value_constraints == v, 3],
                             name=f"{name}_value",
                         )
                     )
             else:
                 data.append(
                     ValuePoints(
-                        locations=self.model.rescale(value_constraints[:, :3]),
+                        locations=value_constraints[:, :3],
                         values=value_constraints[:, 3],
                         name=f"{self.name}_value",
                     )

LoopStructural/modelling/intrusions/intrusion_frame_builder.py

@@ -307,6 +307,13 @@ def add_faults_anisotropies(self, fault_list: list | None = None):
                     fault_i_std,
                 ]
 
+    def _network_model_values(self, data: pd.DataFrame) -> np.ndarray:
+        """
+        The id of the unit at each point. The model is evaluated at the world points; no conversion is applied.
+        """
+        xyz = data.loc[:, ["X", "Y", "Z"]].to_numpy(dtype=float)
+        return self.model.evaluate_model(xyz)
+
     def set_intrusion_steps_parameters(self):
         """
         Use intrusion contact data (reference contact for intrusion frame) to compute the parameters related to each step.
@@ -330,9 +337,6 @@ def set_intrusion_steps_parameters(self):
 
         # set data
         intrusion_network_data = self.intrusion_network_data.copy()
-        intrusion_network_data_xyz = (
-            intrusion_network_data.loc[:, ["X", "Y", "Z"]].copy().to_numpy()
-        )
         std_backup = 25
 
         for step_i, step in self.intrusion_steps.items():
@@ -349,8 +353,8 @@ def set_intrusion_steps_parameters(self):
                 "id"
             )
 
-            intrusion_network_data.loc[:, "model_values"] = self.model.evaluate_model(
-                self.model.rescale(intrusion_network_data_xyz, inplace=False)
+            intrusion_network_data.loc[:, "model_values"] = self._network_model_values(
+                intrusion_network_data
             )
 
             # -- check if step is within the same unit. If so, find clusters of data:

ROADMAP.md

@@ -126,6 +126,20 @@ The policy for a transform change, after constraints exist, is "reject at
 use" (decided 2026-11-05). The interpolator raises `StaleFrameError` at the
 next `set`, `evaluate` or `save` call. The interpolator changes no state.
 
+### Deprecation of GeologicalModel.rescale and scale (spec 005, S4)
+Decision date: 2026-11-05 (spec 005, decisions C-013 and C-016).
+
+`rescale` (local to world) and `scale` (world to local) keep their behaviour
+and their signature. Each one emits a `DeprecationWarning`. The warning names
+`bounding_box.reproject` and `bounding_box.project`.
+
+This replaces the status-log entry "2c-3". That entry said that both methods
+stay as signature-stable pass-throughs. The tier stays "stable" until removal.
+
+Removal is not before 2 minor releases after the first release that warns.
+Removal is also not before the known callers (the QGIS plugin) have moved.
+The library no longer calls either method.
+
 ## Stage sequence
 
 - [x] **Stage 0 — Planning infra.** This file, release/versioning/compat
@@ -866,3 +880,4 @@ next `set`, `evaluate` or `save` call. The interpolator changes no state.
   release-please/PyPI/ruff wiring. `tests/unit`: 695 passed/4 skipped, zero
   regressions.
 - **2026-11-05:** Spec 005 stage S3 implemented: saved interpolator frame, legacy-load rule, frame-keeping copies, stale-frame check. Gate result: root pytest 1978 passed, 5 skipped.
+- **2026-11-05:** Spec 005 stage S4 implemented: display, export and intrusion classification use stored world locations; rescale and scale deprecated.

examples/3_fault/plot_1_faulted_intrusion.py

@@ -68,10 +68,7 @@
 viewer.plot_surface(fault, value=0)
 viewer.plot_vector_field(fault)
 viewer.add_points(
-    model.rescale(
-        model.data[model.data["feature_name"] == "strati"][["X", "Y", "Z"]].values,
-        inplace=False,
-    ),
+    model.data[model.data["feature_name"] == "strati"][["X", "Y", "Z"]].values,
     name="prefault",
 )
 viewer.display()
@@ -97,10 +94,7 @@
 viewer.plot_data(strati)
 viewer.plot_surface(fault, value=0.0)
 viewer.add_points(
-    model.rescale(
-        model.data[model.data["feature_name"] == "strati"][["X", "Y", "Z"]].values,
-        inplace=False,
-    ),
+    model.data[model.data["feature_name"] == "strati"][["X", "Y", "Z"]].values,
     name="prefault",
 )
 viewer.display()