diff --git a/docs/development/usability-plan.md b/docs/development/usability-plan.md index 4cbbc62d..3f5cd6ae 100644 --- a/docs/development/usability-plan.md +++ b/docs/development/usability-plan.md @@ -346,7 +346,7 @@ cross-section to files without the 3D view. ### Phase 6: Follow-up changes -Five changes that come from use of phases 3 to 5. Do them in this order. Each +Nine changes that come from use of phases 3 to 5. Do them in this order. Each one is a separate pull request. #### 6.1 Two separate workflows @@ -560,18 +560,18 @@ Rules: Tasks: -- [ ] `names_to_refresh` returns no result that comes from the column when the +- [x] `names_to_refresh` returns no result that comes from the column when the column has no units. `DerivedRefresh` does not raise when there are no units. -- [ ] `update_model` builds the faults when the column has no groups. +- [x] `update_model` builds the faults when the column has no groups. `update_foliation_features` does nothing for an empty column. Check that `model.stratigraphic_column` is safe to leave unset. -- [ ] `check_stratigraphy` does not ask for the geology layer, the structure +- [x] `check_stratigraphy` does not ask for the geology layer, the structure layer or the contacts when the column is empty. The text says that the step is optional if the user models only faults. -- [ ] `model_state`, `valid` and the primary action accept a model that has +- [x] `model_state`, `valid` and the primary action accept a model that has faults and no groups. -- [ ] Step 5 (view and export) works with fault features only. The block model +- [x] Step 5 (view and export) works with fault features only. The block model and the stratigraphic surfaces are not offered when there are no units. Acceptance: a project has a fault trace layer and an empty column. The user @@ -581,11 +581,247 @@ the fault features, and the user can view and export the fault surfaces. Tests: unit tests for `names_to_refresh` and for the checks with an empty column. A QGIS test of `update_model` with faults and no column. +#### 6.6 Fault topology as a button, not a dialog + +Problem: "Calculate topology..." in step 3 opens a dialog +(`FaultTopologyWidget`). The dialog asks again for the fault layer and the fault +ID field. The user already chose these in the "Fault layer" section of the same +page. The dialog also closes by itself and shows a second message box. + +Rules: + +- The button runs the calculation at once. It uses the fault layer and the name + field that are set in the data manager. +- The button is disabled until a fault layer and a name field are set. The + tooltip says why it is disabled. +- No dialog opens. The result shows in the step message bar (for example, + "Calculated fault topology for 12 pairs."), not in a message box. Errors use + the same message bar. +- The label is "Calculate topology", with no ellipsis, because no dialog opens. + +Tasks: + +- [x] Move the calculation from `FaultTopologyWidget._run_topology` to a + function that does not use Qt widgets. It takes the layer, the ID field + and the data manager, and it returns the number of pairs or raises an + error with a clear message. +- [x] `FaultsStep` calls this function from the button. It reads the layer and + the field from `data_manager.get_fault_traces()`. It updates the enabled + state when the fault layer or the field changes. +- [x] Show the result and the errors in the message bar of the step page. +- [x] Remove `FaultTopologyWidget`, `fault_topology_widget.ui` and + `launchers.FAULT_TOPOLOGY`. Keep the toolbar action only if it still has + a use: if it stays, it calls the same function with the fault layer from + the data manager. + +Acceptance: the user sets the fault layer in step 3 and presses "Calculate +topology". The Fault Adjacency tables fill with no dialog and no second choice +of layer. + +Tests: unit tests for the function (a missing layer, a missing field, an empty +layer, and a normal result). A QGIS test that the button is disabled with no +fault layer. + +#### 6.7 Fault topology sets the relationship order wrongly + +Problem: the topology calculation writes the fault pairs in the wrong order. +`FaultTopologyWidget._run_topology` reads `Fault1` and `Fault2` from the +map2loop table and calls `update_fault_relationship(Fault1, Fault2, ABUTTING)`. +In `FaultTopology` the pair `(a, b)` is directional. In +`apply_fault_abutting_relationships`, `(a, b)` ABUTTING means that fault `a` is +cropped by fault `b`. The map2loop pair has no such direction. The code takes +the table order, so the abutting fault and the fault that it abuts can be the +wrong way round. The crop then removes the wrong part of the model. + +Other places that use order, to check at the same time: + +- `new_faults` is `sorted(...)` on strings, so the fault list order is + alphabetical ("10" before "2"), not the order in the layer. +- The code first removes all pairs with `NONE`, and then adds the new pairs. + Pairs that the user set to FAULTED are lost when the user calculates again. + +Status: I read the LoopStructural and plugin code. I did not read map2loop's +table (the package is not installed here), so the cause in the table must be +confirmed first. + +Rules: + +- For each detected pair, decide which fault ends at the other, from the + geometry of the traces (the fault whose end point lies on the other trace is + the abutting fault). Write the pair as `(abutting fault, other fault)`. +- If the geometry does not give a direction, do not guess. Leave the pair as + `NONE` and list it in the message bar, so the user sets it in the Fault + Adjacency tab. +- Keep the layer order for the fault list. +- A new calculation does not delete a relationship that the user set by hand. + +Tasks: + +- [x] Confirm the columns and the order of the map2loop table. Confirmed in + `map2loop/topology.py`: the columns are `Fault1`, `Fault2`, `Type`, + `Angle`. The pairs come from the lower triangle of a buffer adjacency + matrix (buffer 500 map units), so `Fault1` is only the later fault in + the layer. The table has no direction. Pairs are near each other, and + not always touching. +- [x] Add a function that gives the direction of a pair from two traces. It + has no Qt code. +- [x] Use it in the calculation from 6.6. Keep the layer order of the faults. +- [x] Keep user-set relationships when the user calculates again. + +Acceptance: for two faults where one ends at the other, the Fault Adjacency +table shows the right fault as abutting, and the built fault is cropped on the +right side. A second calculation does not change a relationship that the user +set. + +Tests: unit tests for the direction function (a T-junction, a crossing, two +separate traces, and the two input orders giving the same result). A test that +a second calculation keeps a user-set FAULTED pair. + +#### 6.8 Clean up the viewer and add isosurfaces from model features + +Problem: the viewer code is large and has grown in parts. The files in +`gui/visualisation/` have about 4 000 lines. `feature_list_widget.py` (1 450 +lines) mixes the feature tree, the code that builds meshes from the model, and +the code for cross-sections, block models and topography. +`object_properties_widget.py` (1 100 lines) has many places that read +`viewer.meshes` directly and test `current_object_name`. `object_list_widget.py` +(670 lines) holds the object/model view. The user can add only the surfaces that +the model gives (`feature.surfaces()`), so the user cannot choose a value and +make an isosurface of a scalar field. + +Goals: + +- The viewer classes are smaller and simpler, with one clear job for each. +- The user can add many isosurfaces from a model feature, and can choose the + value of each one. + +Rules: + +- One object registry owns the meshes of the viewer (name, mesh, source feature, + source type, isovalue, style). The widgets read and change objects only through + it. No widget reads `viewer.meshes` directly. +- The code that builds a mesh from the model (scalar field, surface, vector + field, isosurface, block model, cross-section) has no Qt code. The widgets + only call it. +- Each mesh object keeps the information that is needed to build it again (the + source feature, the type and the isovalue), so "Update viewer objects" works + for isosurfaces in the same way as for other objects. +- An isosurface is an object like a surface, with its own name, colour, opacity + and visibility. It does not replace the surfaces that the model gives. + +Tasks: + +- [ ] Read the three files and write a short list of what each class does and + which methods are duplicated or not used. Remove dead code first. +- [ ] Add an object registry class (no Qt) for the meshes and their source + information. Move the reads and writes of `viewer.meshes` to it. +- [ ] Split `feature_list_widget.py`: move the mesh builders to a module without + Qt code (`mesh_builders.py`), and move the cross-section, block model and + topography actions out of the tree widget. +- [ ] Simplify `object_properties_widget.py`: one handler for the selected + object, and the controls shown for each source type, not a check for the + object name in each method. +- [ ] Simplify the object/model view (`object_list_widget.py`): group the + objects by source feature, and show the type and the isovalue of each + object. +- [ ] Add "Add isosurface..." to the menu of a model feature. The user enters + one or more values (a list, or a start, an end and a count). The default + values come from the range of the scalar field. A pure function gives the + values from the input. +- [ ] Build the isosurfaces with the scalar field of the feature + (`feature.surfaces(value)`). Give each object a name with its value, for + example `Fault_1_iso_0.50`. A value outside the range of the field gives + a message in the message bar, not an error dialog. +- [ ] Let the user change the value of an existing isosurface in the object + properties. The object is built again. +- [ ] Update the docs for the viewer. + +Acceptance: the user adds five isosurfaces from one feature in one action. They +show in the object list under that feature, each with its value. The user +changes one value and only that object changes. After the model is updated, +"Update viewer objects" builds the isosurfaces again with the same values. + +Tests: unit tests for the registry and for the function that gives the values +(a list, a range, a value outside the range, a repeated value). A QGIS test that +the menu action adds the objects with the right names and isovalues. + +#### 6.9 True clipping relationships with loop_cgal + +Problem: the model has no way to cut one surface with another. A surface that +goes above the DEM stays in the model and in the exports. The user can only hide +it in the viewer. The same is true for a surface that must end at another +surface, for example a unit that an unconformity cuts. `loop_cgal` does mesh +boolean operations (exact geometry), so it can make true cuts. The plugin does +not use `loop_cgal` now. + +Status: the plugin has no `loop_cgal` code. I did not read the `loop_cgal` API. +Confirm which operations it has (clip a mesh with a mesh, keep the part above or +below, and the result type) and how to install it in the QGIS Python +environment before the design is final. + +Rules: + +- `loop_cgal` is optional. If it is not installed, the clip controls are + disabled with a tooltip that says why, and the model builds as before. +- A clipping relationship has a target surface, a clipping surface (the DEM, a + feature surface, or the bounding box) and a side to keep (above or below). +- The cut is made on the output meshes (the surfaces for the viewer and the + export). It does not change the solved interpolator or the model data. +- The cut is a relationship in the model state. It is saved and loaded with the + state, and the user can turn it off. Old state files have no clipping. +- Cutting by the DEM is a one-step action: "Cut all surfaces by the DEM". It + uses the DEM from step 1 and makes one relationship for each surface. +- The viewer and the export (step 5) use the clipped meshes. The unclipped mesh + is kept, so the user can remove the cut. + +Tasks: + +- [ ] Confirm the `loop_cgal` operations and the install method. Add the + dependency check and a clear message when it is missing. +- [ ] Add a module without Qt code (`main/clipping.py`) with a function that + clips a mesh with a mesh and returns the new mesh, and a function that + makes a surface mesh from the DEM in the model bounding box. +- [ ] Add a `ClippingRelationship` record (target, clipping surface, side, on + or off) to the model manager, with save and load. Add the result to the + object registry of 6.8, so each object has a clipped and an unclipped + mesh. +- [ ] Add "Cut all surfaces by the DEM" to the viewer and to step 5. Add a + "Clip by..." action to the menu of one surface for a cut by another + feature surface. +- [ ] Use the clipped meshes in the surface export and in the cross-section and + block model code where they make sense. Say in the docs which outputs the + cut changes (the block model is not cut). +- [ ] Block model filter: the block model is not cut, but the user can choose + to show only the cells below the DEM. Add a "Below DEM" cell array to the + block model (true when the cell centre is below the DEM height at its x, + y), and a "Show only below DEM" check box in the object properties. The + check box hides the other cells with a threshold filter in the viewer. + It does not delete cells, and it does not need `loop_cgal`. The export of + the block model has the same option, off by default. +- [ ] Handle failures: an open or self-crossing mesh, no overlap, and an empty + result. Show the reason in the message bar and keep the unclipped mesh. +- [ ] Docs: say what the cut does, that it changes only the output meshes, and + how to install `loop_cgal`. + +Acceptance: a model has a DEM and five surfaces, some of which go above the +ground. After "Cut all surfaces by the DEM", no surface in the viewer or in the +export is above the DEM. The user turns the cut off and the full surfaces come +back. With "Show only below DEM" on, the block model shows only the cells below +the ground, and with it off all cells show again. The state saves and loads with +the cut. Without `loop_cgal`, the model +builds and the clip controls are disabled. + +Tests: unit tests for the clipping function with simple meshes (a plane cut by a +plane, no overlap, an empty result) that skip when `loop_cgal` is missing. Unit +tests for the relationship save and load. A QGIS test of the DEM cut action. + Order and links between the parts: 6.1 first, because it defines what the model reads in each workflow. 6.2 depends on the same checks, so do it next. 6.3 comes after 6.2, because it must lay out the final content of the pages. 6.4 does not depend on the others. 6.5 comes after 6.2, because it changes the -same checks. +same checks. 6.6 does not depend on the others. 6.7 goes with 6.6, because both change the +same calculation. 6.8 does not depend on the others. 6.9 comes +after 6.8, because it uses the object registry. ### Phase 7: Fold modelling @@ -882,14 +1118,194 @@ viewer. The new viewer shows geological data objects that are similar to those of Geoscience ANALYST (`geoh5` types): points, curves, surfaces, sections, block models, drillholes and orientations. -This phase is large, and most of the work is in a separate repository. The -tasks, the protocol, the risks and the open questions are in the -[viewer plan](viewer-plan.md). +This phase is large, and most of the work is in a separate repository +(`geoviewer`). The tasks, the protocol, the risks and the open questions are +in the viewer plan of that repository (`docs/viewer-plan.md`). Acceptance: from step 5, a user opens the advanced viewer. The viewer shows the model and its input data, and updates after each build. A pick in the viewer selects the feature in the dock and shows the point on the map. +### Phase 9: Demo mode + +Problem: a new user cannot see what the finished workflow looks like before +they use their own data. The dock has five steps and many buttons. The docs +have screenshots, but they go out of date. Trainers and developers also need a +repeatable way to show the plugin, and a way to test the full workflow through +the real UI. + +Aim: a demo mode runs the workflow by itself. It clicks the real buttons, in +the order that the user must use them, with sample data, until a model is +built. The user watches and can pause, step or stop at any time. + +The demo does not use a second code path. It calls the same widgets as a user +does. Thus the demo also shows when a step is broken, and it is a UI test. + +#### Rules + +- **Real controls.** The demo changes a control (a layer picker, a combo, a + check box) and presses a button through the widget. It does not call the + model manager or the data manager directly. If a control is hidden or + disabled, the demo stops with an error. This is a test failure, not a case + to work around. +- **Visible.** Before each action, the demo moves a highlight (an overlay frame + and a one-line caption) to the control. It waits a set time, then it acts. + The caption says what the action does and why ("Build column from the map: + the order of units comes from the geology layer"). +- **Same speed as the work.** The demo waits for the end of each background + task (extraction, thickness, build) with the task signals. It does not wait + a fixed time. If a task fails, the demo stops and shows the error. +- **User control.** A small demo bar shows Pause, Step, Speed (slow, normal, + fast) and Stop. The demo also stops if the user clicks or types in the dock. + The demo never runs by itself at start-up. +- **Safe project.** The demo works in a new, empty QGIS project, and it asks + before it replaces the current project. If the current project has unsaved + changes, the demo asks the user to save them first. Stop removes the demo + layers and restores the plugin state to the state before the demo. +- **Sample data.** The demo uses a small data set that the plugin includes + (see 9.1). It does not need a network. +- **Both entries.** There are two demos: "Build from a geological map" and + "Interpolate surfaces from constraints". The user selects one in the demo + menu. +- **No modal dialogs.** The demo does not open a modal dialog. If a step uses + a dialog (for example "Build column" or "Add Foliation"), the demo drives + the dialog and closes it. A dialog that the demo cannot drive is a problem + in the step (see the general rules, rule 2). + +#### Design + +A demo is a list of steps in data. A step has a target, an action and a +caption. The runner is a small state machine. The list and the runner do not +import QGIS widgets, so unit tests can run them. + +``` +demo script (data) runner dock +------------------ ------ ---- + step: target "step2.build_column" -> find control by id -> highlight + action click wait (speed) click() + wait_for "column_changed" wait for signal <- signal + caption "..." next step +``` + +- **Target ids.** Each control that a demo uses has a stable `objectName` + (for example `step2.build_column`, `step4.primary_button`). The runner finds + the control with `findChild`. A rename of an id fails a test. The ids are + also usable by the other UI tests. +- **Actions.** `click`, `set_text`, `select_layer`, `select_item` (a combo or + a menu entry), `check`, `go_to_step`, `wait_for` (a signal or a condition + with a time limit). +- **Order from the dock.** The script does not hold the order of the workflow + as a second copy. Where the next action is "press Next" or "press the primary + button", the runner reads `choose_primary_action` and the footer text from the + dock. A step that does not match the dock is a failure. +- **Checks.** After a step, the runner can assert a condition (the step status + is "done", the footer shows no problem). A failed check stops the demo and + names the step. This lets the CI run the same script without the delays. +- **Headless mode.** The same runner runs with no highlight and with zero + delay, in the QGIS test job. This is the end-to-end test. + +#### The two scripts + +Map demo (about 20 actions): + +1. Step 1: set the bounding box from the geology layer, set the CRS, and + select the DEM. +2. Step 2: select the geology layer and the unit name field. "Build column" + from the map. Calculate basal contacts and thickness ("Derive from map"). +3. Step 3: select the fault layer and the name field. "Calculate topology". + Show the adjacency tables. +4. Step 4: press the primary button (build and solve). +5. Step 5: "Open 3D view". Show an export action, with the output in a + temporary folder. + +Constraints demo (about 12 actions): + +1. Step 1: set the bounding box. Select "Interpolate surfaces from + constraints". +2. Step 4: "Add Foliation" from a value layer and an orientation layer, then + press the primary button. +3. Step 4: "Add Fold Event" and a folded foliation, if phase 7 is merged. +4. Step 5: "Open 3D view". + +Each demo ends with a message in the message bar and the demo bar shows +"Demo finished. Stop to remove the demo data, or keep it to explore." + +#### Tasks + +Do the parts in this order. Each part is a separate pull request. + +9.1 Sample data + +- [ ] Choose a small open data set that gives a model in less than one + minute (a part of the Hamersley data, or a synthetic area). Check the + licence, and record the source in the docs. +- [ ] Add the layers as one GeoPackage and one small DEM in + `loopstructural/resources/demo/`. Keep the package size small (state the + limit in the pull request, for example 5 MB). +- [ ] A function that adds the layers to a new project, with names and styles, + and a function that removes them. It does not need the widgets. + +9.2 Ids and signals + +- [ ] Give each control that the demos use a stable `objectName`. Add a test + that lists the ids and finds each one in the dock. +- [ ] Check that each action of the demos ends with a signal that the runner + can wait for (column changed, derived data updated, topology + calculated, model built or failed). Add the signal where it is missing. + +9.3 Runner and overlay + +- [ ] `gui/demo/script.py`: the step form, the actions, the speed settings. No + QGIS import. +- [ ] `gui/demo/runner.py`: the state machine (run, pause, step, stop, + error). It takes the time source and the control finder as arguments, so + unit tests can use fakes. +- [ ] `gui/demo/overlay.py`: the highlight frame and the caption. It follows the + control when the dock moves or resizes, and it works for a dock that is + floating, tabbed or in a separate window (`separate_dock_widgets`). +- [ ] `gui/demo/demo_bar.py`: Pause, Step, Speed and Stop. Stop on a user click + in the dock. + +9.4 Scripts and entry + +- [ ] The map demo and the constraints demo as data in `gui/demo/scripts/`. +- [ ] "Demo" in the dock header menu and in the Plugins menu. It asks before it + replaces the project, then loads the sample data and starts the script. +- [ ] Stop removes the demo layers, resets the plugin state and restores the + state before the demo (use the save and load of the application state). + +9.5 Tests and docs + +- [ ] A headless QGIS test that runs each script with zero delay and checks the + result: the model is solved, and the number of features is as expected. + Run it in the QGIS test job, so a change to the UI that breaks the + workflow fails the build. +- [ ] A user guide page in `docs/usage` ("Try the demo"). Say what the demo + does, how to stop it, and what it changes in the project. +- [ ] Use the demo to make the screenshots of the docs, so that they match the + current UI. + +Files: new `gui/demo/` module, `resources/demo/`, `plugin_main.py`, +`gui/modelling/steps/header.py`, and the step pages (ids and signals only). + +Acceptance: a new user opens the demo from the dock menu. The plugin loads the +sample data and clicks through the steps, with a caption for each action, and +it ends with a solved model in the 3D view. The user can pause, step and stop. +After Stop, the project and the plugin state are as they were before the demo. +A script that is out of date (for example, a control was renamed or a step does +not become "done") stops with an error that names the step. + +Tests: unit tests for the runner (run, pause, step, stop; a missing control; a +timeout; a failed check; a failed task) with a fake clock and fake controls. +Unit tests for the script form and for the sample data functions. A QGIS test +that runs both demos end to end with zero delay. + +Order and links: 9.1 and 9.2 do not depend on each other. 9.3 depends on 9.2. +9.4 depends on 9.1 and 9.3. 9.5 comes last. Phase 9 depends on phase 3 (the +steps), phase 4 (the primary button and the two entries) and phase 5 (step 5). +The Fold Event part of the constraints demo depends on phase 7, and the 3D +view action can use the PyVista viewer or the viewer of phase 8. + ## Risks - **Large UI change.** Users of the current version must learn the new layout. @@ -926,6 +1342,15 @@ viewer selects the feature in the dock and shows the point on the map. arguments (`limb_wl`, `axis_wl`, `av_fold_axis`, the profile types). Some of them are keyword arguments, not public API. Pin the LoopStructural version and add a test for each argument. +- **Demo scripts become out of date (phase 9).** A script depends on the ids and + the order of the controls. Run both scripts in the QGIS test job, so that a UI + change that breaks them fails the build. Keep the script as data, with no + copy of the workflow order. +- **Demo changes the user project (phase 9).** The demo adds layers and + changes the plugin state. Use a new project, ask before it replaces the + current one, and restore the state on Stop. +- **Package size (phase 9).** The sample data increases the size of the plugin + package. Keep it small, and decide in the review if it must be a download. ## Open questions @@ -967,3 +1392,9 @@ viewer selects the feature in the dock and shows the point on the map. 11. (7) Must a fault that cuts a fold frame also cut the folded features? LoopStructural calculates the rotation angles in the restored space. Recommendation: yes, use the same faults. Test it before 7.5. +12. (9) Must the sample data be in the plugin package, or a download from the + first demo run? Recommendation: in the package if it is below 5 MB. If it + is larger, download it once and keep it in the QGIS profile folder. +13. (9) Must the demo also work from the Processing toolbox or the Python + console (for example, `loopstructural.run_demo("map")`)? Recommendation: + yes, one function. The test job and the screenshot job use it. diff --git a/loopstructural/gui/map2loop_tools/fault_topology_widget.py b/loopstructural/gui/map2loop_tools/fault_topology_widget.py deleted file mode 100644 index 0e50fdd8..00000000 --- a/loopstructural/gui/map2loop_tools/fault_topology_widget.py +++ /dev/null @@ -1,243 +0,0 @@ -"""Widget for calculating fault topology from a fault layer.""" - -import os -from contextlib import nullcontext - -import geopandas as gpd -from qgis.core import QgsMapLayerProxyModel -from qgis.PyQt import uic -from qgis.PyQt.QtWidgets import QDialog, QMessageBox - -from ..compatibility import configure_layer_combo - - -class FaultTopologyWidget(QDialog): - """Widget for calculating fault topology from a fault layer.""" - - def __init__(self, parent=None, data_manager=None, debug_manager=None): - super().__init__(parent) - self.data_manager = data_manager - # Load the UI file - ui_path = os.path.join(os.path.dirname(__file__), "fault_topology_widget.ui") - uic.loadUi(ui_path, self) - # Set filter for fault layer selection - configure_layer_combo(self.faultLayerComboBox, QgsMapLayerProxyModel.Filter.LineLayer) - self.faultLayerComboBox.layerChanged.connect(self._on_fault_layer_changed) - # react to field changes so we can update the modelling widget via the data manager - try: - # QgsFieldComboBox uses fieldChanged signal - self.faultIdFieldComboBox.fieldChanged.connect(self._on_fault_field_changed) - except Exception: - pass - - self.runButton.clicked.connect(self._run_topology) - # After attempting to guess, synchronise with current data manager state (if any) - self._sync_with_data_manager() - - def _on_fault_layer_changed(self): - layer = self.faultLayerComboBox.currentLayer() - self.faultIdFieldComboBox.setLayer(layer) - # Optionally auto-select a likely ID field - if layer: - fields = [field.name() for field in layer.fields()] - for name in ["id", "ID", "fault_id", "FaultID", "FaultId"]: - if name in fields: - self.faultIdFieldComboBox.setField(name) - break - # Inform the data manager / modelling widgets about the change and preserve other settings - self._update_data_manager_fault_layer() - - def _on_fault_field_changed(self): - # When the selected ID field changes, update the data manager so the modelling widget updates - self._update_data_manager_fault_layer() - - def _sync_with_data_manager(self): - """Set the widget UI to reflect the current fault traces selection in the data manager.""" - if not hasattr(self, 'data_manager') or self.data_manager is None: - print("No data manager to sync with") - return - try: - fault_traces = self.data_manager.get_fault_traces() - except Exception: - fault_traces = None - if not fault_traces: - return - layer = fault_traces.get('layer') - print(f"Syncing fault topology widget with layer: {layer}") - if layer is not None: - try: - self.faultLayerComboBox.setLayer(layer) - except Exception: - pass - # set the name field if available - fault_name_field = fault_traces.get('fault_name_field') - if fault_name_field and layer is not None: - try: - self.faultIdFieldComboBox.setLayer(layer) - self.faultIdFieldComboBox.setField(fault_name_field) - except Exception: - pass - - def _update_data_manager_fault_layer(self): - """Update the ModellingDataManager with the layer/field chosen in this widget. - - Preserve any other fault settings already present in the data manager (dip, displacement, use_z). - """ - if not hasattr(self, 'data_manager') or self.data_manager is None: - return - # Gather current selections from this widget - layer = self.faultLayerComboBox.currentLayer() - name_field = None - try: - name_field = self.faultIdFieldComboBox.currentField() - except Exception: - name_field = None - # Preserve existing settings from data manager if present - existing = self.data_manager.get_fault_traces() or {} - fault_dip_field = existing.get('fault_dip_field') - fault_displacement_field = existing.get('fault_displacement_field') - use_z_coordinate = existing.get('use_z_coordinate', False) - # Call data manager to set the fault trace layer which will notify the modelling UI - try: - self.data_manager.set_fault_trace_layer( - layer, - fault_name_field=name_field, - fault_dip_field=fault_dip_field, - fault_displacement_field=fault_displacement_field, - use_z_coordinate=use_z_coordinate, - ) - except Exception: - # Fail silently to avoid breaking UI if data_manager is not fully initialised - return - - def _run_topology(self): - layer = self.faultLayerComboBox.currentLayer() - if not layer: - QMessageBox.warning(self, "Missing Input", "Please select a fault layer.") - return - id_field = self.faultIdFieldComboBox.currentField() - if not id_field: - QMessageBox.warning(self, "Missing Input", "Please select a fault ID field.") - return - # Convert to GeoDataFrame - gdf = gpd.GeoDataFrame.from_features(layer.getFeatures()) - if gdf.empty: - QMessageBox.warning(self, "No Data", "The selected layer has no features.") - return - # Rename the selected ID field to 'ID' for Topology class compatibility - if id_field != "ID": - gdf = gdf.rename(columns={id_field: "ID"}) - # Use map2loop Topology class - try: - from map2loop.topology import Topology - except ImportError: - QMessageBox.critical(self, "Error", "Could not import map2loop Topology class.") - return - topology = Topology(geology_data=None, fault_data=gdf) - df = topology.fault_fault_relationships - - # Update the modelling FaultTopology (so the Fault Adjacency tab refreshes) - if hasattr(self, 'data_manager') and self.data_manager is not None: - try: - from LoopStructural.modelling.core.fault_topology import FaultRelationshipType - - ft = self.data_manager._fault_topology - model_manager = getattr(self.data_manager, '_model_manager', None) - # Repopulating the whole topology fires one notification per - # add_fault/remove_fault/update_fault_relationship call below. - # Each notification normally triggers a full O(faults^2) - # rescan in the model manager (see - # `batch_fault_topology_updates`), so for a map2loop run with - # many faults that's a lot of redundant, slow work before the - # dialog even closes. Batch it into a single rescan. - batch_cm = ( - model_manager.batch_fault_topology_updates() - if model_manager is not None - else nullcontext() - ) - with batch_cm: - # Remove existing fault-fault relationships (notify observers) - for f1, f2 in list(ft.adjacency.keys()): - try: - ft.update_fault_relationship(f1, f2, FaultRelationshipType.NONE) - except Exception: - pass - - # Remove existing stratigraphy relationships - for unit, fault in list(ft.stratigraphy_fault_relationships.keys()): - try: - ft.update_fault_stratigraphy_relationship(unit, fault, False) - except Exception: - pass - - # Determine faults from the fault layer itself (all IDs present), - # not just the ones map2loop found a relationship for. A fault with - # no detected topological relationship is still a real fault and - # must not be dropped from the fault topology. - new_faults = {str(v) for v in gdf['ID'].unique()} - - # Add new faults; never remove existing ones here, so faults - # without a detected relationship (or ones the user added - # manually) are preserved and only updated, not deleted. - for f in sorted(new_faults): - if f not in ft.faults: - try: - ft.add_fault(f) - except Exception: - pass - - # Add relationships from df. map2loop detects these pairs purely - # from spatial proximity/intersection of fault traces (see - # Topology._calculate_fault_fault_relationships), the same signal - # LoopStructural's own map2loop processor (GeologicalModel.from_processor) - # resolves to a splay or an ABUTTING relationship -- never a FAULTED - # (interpolation-chaining) one. Mark these ABUTTING here too: it's the - # cheap, post-hoc region crop, and matches what this data actually - # represents (near/intersecting faults, not a deliberate cross-cutting - # order). Leave FAULTED for the user to set explicitly in the Fault - # Adjacency tab when they really do want one fault's interpolation to - # depend on another -- auto-marking every detected pair FAULTED made - # every fault chain-build against every nearby fault, which is - # combinatorially expensive and freezes Initialize Model on models - # with more than a handful of faults. - if df is not None and not df.empty: - for _, row in df.iterrows(): - try: - if 'Fault1' in row.index and 'Fault2' in row.index: - f1 = str(row['Fault1']) - f2 = str(row['Fault2']) - else: - f1 = str(row.iloc[0]) - f2 = str(row.iloc[1]) - ft.update_fault_relationship(f1, f2, FaultRelationshipType.ABUTTING) - except Exception: - pass - - # Update unit-fault relationships if available from topology - try: - uf = topology.unit_fault_relationships - if uf is not None and not uf.empty: - for _, r in uf.iterrows(): - try: - unit = r.get('Unit', r.iloc[0]) - fault = r.get('Fault', r.iloc[1]) - ft.update_fault_stratigraphy_relationship( - unit, str(fault), True - ) - except Exception: - pass - except Exception: - # unit-fault relationships not available - pass - - except Exception: - # If anything fails here, still continue to show success of topology run - pass - - QMessageBox.information( - self, - "Success", - f"Calculated fault topology for {len(df) if df is not None else 0} pairs.", - ) - self.close() - return True diff --git a/loopstructural/gui/map2loop_tools/fault_topology_widget.ui b/loopstructural/gui/map2loop_tools/fault_topology_widget.ui deleted file mode 100644 index 79522883..00000000 --- a/loopstructural/gui/map2loop_tools/fault_topology_widget.ui +++ /dev/null @@ -1,64 +0,0 @@ - - - FaultTopologyWidget - - - - 0 - 0 - 400 - 120 - - - - Fault Topology Calculator - - - - - - - - Fault Layer: - - - - - - - - - - Fault ID Field: - - - - - - - - - - - - Calculate Topology - - - - - - - - QgsMapLayerComboBox - QComboBox -
qgsmaplayercombobox.h
- - - QgsFieldComboBox - QComboBox -
qgsfieldcombobox.h
- - - - - diff --git a/loopstructural/gui/map2loop_tools/launchers.py b/loopstructural/gui/map2loop_tools/launchers.py index e14f4aa8..d3a4bf25 100644 --- a/loopstructural/gui/map2loop_tools/launchers.py +++ b/loopstructural/gui/map2loop_tools/launchers.py @@ -6,7 +6,6 @@ BASAL_CONTACTS = 'basal_contacts' THICKNESS = 'thickness' PAINT_STRAT_ORDER = 'paint_strat_order' -FAULT_TOPOLOGY = 'fault_topology' DATA_CONVERSION = 'data_conversion' @@ -16,10 +15,6 @@ def _dialog_class(tool): from loopstructural.gui.data_conversion import AutomaticConversionDialog return AutomaticConversionDialog - if tool == FAULT_TOPOLOGY: - from loopstructural.gui.map2loop_tools.fault_topology_widget import FaultTopologyWidget - - return FaultTopologyWidget from loopstructural.gui import map2loop_tools names = { diff --git a/loopstructural/gui/modelling/steps/pages.py b/loopstructural/gui/modelling/steps/pages.py index d5c6be6b..cc72f27d 100644 --- a/loopstructural/gui/modelling/steps/pages.py +++ b/loopstructural/gui/modelling/steps/pages.py @@ -4,7 +4,8 @@ The existing tabs and dialogs are the contents of the pages. """ -from qgis.core import QgsApplication +from qgis.core import Qgis, QgsApplication +from qgis.gui import QgsMessageBar from qgis.PyQt.QtCore import Qt, pyqtSignal from qgis.PyQt.QtWidgets import ( QComboBox, @@ -28,6 +29,11 @@ ) from loopstructural.main import layer_roles +from loopstructural.main.fault_topology_calc import ( + FaultTopologyError, + calculate_from_data_manager, + result_message, +) from loopstructural.main.workflow_mode import WORKFLOW_MODE_LABELS, WORKFLOW_MODES from . import checks @@ -204,17 +210,59 @@ def __init__(self, parent=None, **kwargs): ) layout.addWidget(self.sections) - topology = QPushButton("Calculate topology...", self) - topology.setToolTip("Find which faults touch each other, from the fault traces.") - topology.clicked.connect(lambda _checked=False: self._show_tool(launchers.FAULT_TOPOLOGY)) + self.message_bar = QgsMessageBar(self) + layout.addWidget(self.message_bar) + + self.topology_button = QPushButton("Calculate topology", self) + self.topology_button.clicked.connect(self._calculate_topology) row = QHBoxLayout() row.addStretch(1) - row.addWidget(topology) + row.addWidget(self.topology_button) layout.addLayout(row) + # The button follows the fault layer and the name field + self.fault_layers.faultTraceLayer.layerChanged.connect(self.update_topology_button) + self.fault_layers.faultNameField.fieldChanged.connect(self.update_topology_button) + if self.data_manager is not None: + self.data_manager.layer_roles.attach(lambda _role, _value: self.update_topology_button()) + self.update_topology_button() self.adjacency = FaultAdjacencyTab(self, data_manager=self.data_manager) layout.addWidget(page_scroll_area(self.adjacency, self), 1) + def _fault_traces(self): + traces = self.data_manager.get_fault_traces() if self.data_manager is not None else None + return traces or {} + + def update_topology_button(self, *_args): + """Enable the button when the fault layer and the name field are set.""" + traces = self._fault_traces() + layer = self.data_manager.get_layer_role(layer_roles.FAULT_TRACES) if self.data_manager else None + if layer is None or not traces.get('fault_name_field'): + self.topology_button.setEnabled(False) + self.topology_button.setToolTip( + "Select a fault layer and the field with the fault name in the Fault layer section." + ) + else: + self.topology_button.setEnabled(True) + self.topology_button.setToolTip( + "Find which faults touch each other, from the fault traces." + ) + + def _calculate_topology(self): + self.message_bar.clearWidgets() + try: + result = calculate_from_data_manager(self.data_manager) + except FaultTopologyError as error: + self.message_bar.pushMessage("Fault topology", str(error), level=Qgis.MessageLevel.Critical) + return + text, warn = result_message(result) + self.message_bar.pushMessage( + "Fault topology", + text, + level=Qgis.MessageLevel.Warning if warn else Qgis.MessageLevel.Success, + duration=0 if warn else 8, + ) + def _fault_layer_summary(self): if self.data_manager is None: return "" diff --git a/loopstructural/main/fault_topology_calc.py b/loopstructural/main/fault_topology_calc.py new file mode 100644 index 00000000..7658a084 --- /dev/null +++ b/loopstructural/main/fault_topology_calc.py @@ -0,0 +1,219 @@ +"""Calculate the fault topology from the fault traces. + +The module has no widgets. The button in step 3 calls `calculate_fault_topology`. +The direction of a pair comes from the geometry of the traces and not from the +order of the map2loop table. +""" + +from collections import namedtuple +from contextlib import nullcontext + +from LoopStructural.modelling.core.fault_topology import FaultRelationshipType + +TopologyResult = namedtuple('TopologyResult', ['pairs', 'undetermined']) +"""``pairs`` is the number of pairs that map2loop found. +``undetermined`` is the list of ``(fault, fault)`` pairs without a direction, +or whose relationship was set by the user and is not changed.""" + + +class FaultTopologyError(Exception): + """The topology cannot be calculated. The text is safe to show to the user.""" + + +def _end_points(geometry): + """Return the end points of a trace as shapely points.""" + from shapely.ops import linemerge + + if geometry is None or geometry.is_empty: + return [] + if geometry.geom_type == 'MultiLineString': + geometry = linemerge(geometry) + boundary = geometry.boundary + if boundary.is_empty: # a closed trace has no end points + return [] + return list(boundary.geoms) if hasattr(boundary, 'geoms') else [boundary] + + +def default_tolerance(geometry_a, geometry_b): + """A distance that is small for the size of the two traces.""" + xmin_a, ymin_a, xmax_a, ymax_a = geometry_a.bounds + xmin_b, ymin_b, xmax_b, ymax_b = geometry_b.bounds + width = max(xmax_a, xmax_b) - min(xmin_a, xmin_b) + height = max(ymax_a, ymax_b) - min(ymin_a, ymin_b) + return max((width**2 + height**2) ** 0.5 * 1e-3, 1e-9) + + +def fault_pair_direction(geometry_a, geometry_b, tolerance=None): + """Find which of two fault traces ends at the other. + + The abutting fault is the fault with an end point on the other trace. + + Parameters + ---------- + geometry_a, geometry_b : shapely geometry + The two traces (a line or a multi-line). + tolerance : float, optional + The largest distance between an end point and the other trace that + counts as "on the trace". The default is 0.1 % of the size of the + two traces. + + Returns + ------- + int or None + 0 if ``geometry_a`` abuts ``geometry_b``, 1 if ``geometry_b`` abuts + ``geometry_a``, and None if the geometry gives no direction (the + traces cross, are separate, or end at each other). + """ + if tolerance is None: + tolerance = default_tolerance(geometry_a, geometry_b) + a_abuts = any(p.distance(geometry_b) <= tolerance for p in _end_points(geometry_a)) + b_abuts = any(p.distance(geometry_a) <= tolerance for p in _end_points(geometry_b)) + if a_abuts and not b_abuts: + return 0 + if b_abuts and not a_abuts: + return 1 + return None + + +def update_topology(topology, traces, pairs, model_manager=None): + """Write the faults and the relationships of the pairs into the topology. + + Parameters + ---------- + topology : FaultTopology + The topology of the modelling data manager. + traces : dict + Fault name to shapely geometry, in the order of the layer. + pairs : iterable of tuple + ``(name, name)`` pairs that are close or touch. The order is not used. + model_manager : GeologicalModelManager, optional + The notifications of the topology are batched when this is given. + + Returns + ------- + list of tuple + The pairs that were not written: the geometry gives no direction, or + the user already set a relationship for the pair. + """ + skipped = [] + batch = ( + model_manager.batch_fault_topology_updates() if model_manager is not None else nullcontext() + ) + with batch: + # Add the faults in the order of the layer. Faults that are in the + # topology already stay, with their relationships. + for name in traces: + if name not in topology.faults: + topology.add_fault(name) + for first, second in pairs: + first, second = str(first), str(second) + if first not in traces or second not in traces or first == second: + continue + # A relationship in either direction was set before: keep it. + if (first, second) in topology.adjacency or (second, first) in topology.adjacency: + skipped.append((first, second)) + continue + direction = fault_pair_direction(traces[first], traces[second]) + if direction is None: + skipped.append((first, second)) + continue + abutting, other = (first, second) if direction == 0 else (second, first) + topology.update_fault_relationship(abutting, other, FaultRelationshipType.ABUTTING) + return skipped + + +def _traces_from_layer(layer, id_field): + """Return ``(geodataframe, traces)`` from a vector layer. + + ``traces`` is a dict of the fault name to one shapely geometry, in the + order of the layer. + """ + import geopandas as gpd + + gdf = gpd.GeoDataFrame.from_features(layer.getFeatures()) + if gdf.empty: + raise FaultTopologyError("The fault layer has no features.") + if id_field not in gdf.columns: + raise FaultTopologyError(f"The fault layer has no field '{id_field}'.") + if id_field != "ID": + gdf = gdf.rename(columns={id_field: "ID"}) + traces = {} + for name, geometry in zip(gdf["ID"], gdf.geometry): + if geometry is None or geometry.is_empty: + continue + name = str(name) + if name in traces: + traces[name] = traces[name].union(geometry) + else: + traces[name] = geometry + return gdf, traces + + +def calculate_fault_topology(layer, id_field, data_manager): + """Calculate the fault topology and write it into the data manager. + + Parameters + ---------- + layer : QgsVectorLayer + The fault traces. + id_field : str + The field with the name of the fault. + data_manager : ModellingDataManager + The data manager that has the fault topology. + + Returns + ------- + TopologyResult + The number of pairs and the pairs without a direction. + + Raises + ------ + FaultTopologyError + If the layer or the field is missing, the layer is empty, or map2loop + is not installed. + """ + if layer is None: + raise FaultTopologyError("Select a fault layer.") + if not id_field: + raise FaultTopologyError("Select the field with the fault name.") + gdf, traces = _traces_from_layer(layer, id_field) + try: + from map2loop.topology import Topology + except ImportError as error: + raise FaultTopologyError("Could not import the map2loop Topology class.") from error + table = Topology(geology_data=None, fault_data=gdf).fault_fault_relationships + pairs = [] + if table is not None and not table.empty: + for _, row in table.iterrows(): + if 'Fault1' in row.index and 'Fault2' in row.index: + pairs.append((row['Fault1'], row['Fault2'])) + else: + pairs.append((row.iloc[0], row.iloc[1])) + skipped = update_topology( + data_manager._fault_topology, + traces, + pairs, + getattr(data_manager, '_model_manager', None), + ) + return TopologyResult(len(pairs), skipped) + + +def calculate_from_data_manager(data_manager): + """Calculate the topology with the fault layer and field of the data manager.""" + traces = data_manager.get_fault_traces() if data_manager is not None else None + traces = traces or {} + return calculate_fault_topology( + traces.get('layer'), traces.get('fault_name_field'), data_manager + ) + + +def result_message(result): + """Return the text for the message bar, and True if it is a warning.""" + text = f"Calculated fault topology for {result.pairs} pairs." + if result.undetermined: + names = ", ".join(f"{a}-{b}" for a, b in result.undetermined) + text += ( + f" These pairs were not set (no direction from the traces, or set before): {names}." + " Set them in the Fault Adjacency table." + ) + return text, bool(result.undetermined) diff --git a/loopstructural/plugin_main.py b/loopstructural/plugin_main.py index cccdaa1c..69efb481 100644 --- a/loopstructural/plugin_main.py +++ b/loopstructural/plugin_main.py @@ -144,7 +144,7 @@ def initGui(self): # -- Actions self.action_fault_topology = QAction( - self.tr("Fault Topology Calculator"), + self.tr("Calculate fault topology"), self.iface.mainWindow(), ) self.action_fault_topology.triggered.connect(self.show_fault_topology_dialog) @@ -423,8 +423,21 @@ def show_paint_strat_order_dialog(self): self._show_tool(launchers.PAINT_STRAT_ORDER) def show_fault_topology_dialog(self): - """Show the fault topology calculator dialog.""" - self._show_tool(launchers.FAULT_TOPOLOGY) + """Calculate the fault topology for the fault layer of the data manager.""" + from loopstructural.gui.messages import push_success, push_warning + from loopstructural.main.fault_topology_calc import ( + FaultTopologyError, + calculate_from_data_manager, + result_message, + ) + + try: + result = calculate_from_data_manager(self.data_manager) + except FaultTopologyError as error: + push_warning("Fault topology", str(error)) + return + text, warn = result_message(result) + (push_warning if warn else push_success)("Fault topology", text) def tr(self, message: str) -> str: """Translate a string using Qt translation API. diff --git a/tests/qgis/test_fault_topology_button.py b/tests/qgis/test_fault_topology_button.py new file mode 100644 index 00000000..b1bc2a5a --- /dev/null +++ b/tests/qgis/test_fault_topology_button.py @@ -0,0 +1,37 @@ +"""The "Calculate topology" button of step 3 follows the fault layer.""" + +from unittest.mock import Mock + +import pytest +from qgis.core import QgsProject, QgsVectorLayer + +from loopstructural.gui.modelling.steps.pages import FaultsStep +from loopstructural.main.data_manager import ModellingDataManager +from loopstructural.main.model_manager import GeologicalModelManager + + +@pytest.fixture +def page(): + project = QgsProject.instance() + model_manager = GeologicalModelManager() + data_manager = ModellingDataManager( + project=project, mapCanvas=Mock(), logger=Mock() + ) + data_manager.set_model_manager(model_manager) + return FaultsStep(None, data_manager=data_manager, model_manager=model_manager) + + +def test_button_is_disabled_without_a_fault_layer(page): + assert not page.topology_button.isEnabled() + assert 'fault layer' in page.topology_button.toolTip() + + +def test_button_is_enabled_with_layer_and_name_field(page): + layer = QgsVectorLayer('LineString?field=name:string', 'faults', 'memory') + QgsProject.instance().addMapLayer(layer) + page.data_manager.set_fault_trace_layer(layer, fault_name_field='name') + page.update_topology_button() + assert page.topology_button.isEnabled() + page.data_manager.set_fault_trace_layer(layer, fault_name_field=None) + page.update_topology_button() + assert not page.topology_button.isEnabled() diff --git a/tests/unit/test_fault_topology_calc.py b/tests/unit/test_fault_topology_calc.py new file mode 100644 index 00000000..21c3d913 --- /dev/null +++ b/tests/unit/test_fault_topology_calc.py @@ -0,0 +1,144 @@ +"""Pytest tests for the fault topology calculation. + +The module has no QGIS code, so the tests use traces from shapely and run in +the fast tests/unit/ job. +""" + +import pytest +from LoopStructural import FaultTopology, StratigraphicColumn +from LoopStructural.modelling.core.fault_topology import FaultRelationshipType +from shapely.geometry import LineString + +from loopstructural.main import fault_topology_calc as calc + +ABUTTING = FaultRelationshipType.ABUTTING +FAULTED = FaultRelationshipType.FAULTED + + +def t_junction(): + """Fault "stem" ends on the middle of fault "bar".""" + bar = LineString([(0, 0), (10, 0)]) + stem = LineString([(5, 0), (5, 10)]) + return bar, stem + + +def test_t_junction_stem_abuts(): + bar, stem = t_junction() + assert calc.fault_pair_direction(stem, bar) == 0 + assert calc.fault_pair_direction(bar, stem) == 1 + + +def test_input_order_gives_the_same_fault(): + bar, stem = t_junction() + first = calc.fault_pair_direction(bar, stem) + second = calc.fault_pair_direction(stem, bar) + assert [bar, stem][first] is stem + assert [stem, bar][second] is stem + + +def test_crossing_has_no_direction(): + a = LineString([(0, 5), (10, 5)]) + b = LineString([(5, 0), (5, 10)]) + assert calc.fault_pair_direction(a, b) is None + + +def test_separate_traces_have_no_direction(): + a = LineString([(0, 0), (1, 0)]) + b = LineString([(0, 5), (1, 5)]) + assert calc.fault_pair_direction(a, b) is None + + +def test_end_to_end_has_no_direction(): + a = LineString([(0, 0), (5, 0)]) + b = LineString([(5, 0), (10, 0)]) + assert calc.fault_pair_direction(a, b) is None + + +def make_topology(): + return FaultTopology(StratigraphicColumn()) + + +def test_update_writes_the_abutting_fault_first(): + bar, stem = t_junction() + topology = make_topology() + # the table order is "bar" then "stem": the wrong way round + skipped = calc.update_topology(topology, {'bar': bar, 'stem': stem}, [('bar', 'stem')]) + assert skipped == [] + assert topology.adjacency == {('stem', 'bar'): ABUTTING} + + +def test_update_keeps_the_layer_order_of_the_faults(): + traces = {name: LineString([(i * 20, 0), (i * 20 + 1, 0)]) for i, name in enumerate(['10', '2', '1'])} + topology = make_topology() + calc.update_topology(topology, traces, []) + assert topology.faults == ['10', '2', '1'] + + +def test_pair_without_direction_is_listed(): + a = LineString([(0, 5), (10, 5)]) + b = LineString([(5, 0), (5, 10)]) + topology = make_topology() + skipped = calc.update_topology(topology, {'a': a, 'b': b}, [('a', 'b')]) + assert skipped == [('a', 'b')] + assert topology.adjacency == {} + + +def test_second_calculation_keeps_user_relationship(): + bar, stem = t_junction() + topology = make_topology() + traces = {'bar': bar, 'stem': stem} + calc.update_topology(topology, traces, [('bar', 'stem')]) + topology.update_fault_relationship('stem', 'bar', FAULTED) + skipped = calc.update_topology(topology, traces, [('bar', 'stem')]) + assert topology.adjacency == {('stem', 'bar'): FAULTED} + assert skipped == [('bar', 'stem')] + + +def test_user_relationship_in_the_other_direction_is_kept(): + bar, stem = t_junction() + topology = make_topology() + traces = {'bar': bar, 'stem': stem} + calc.update_topology(topology, traces, []) + topology.update_fault_relationship('bar', 'stem', FAULTED) + calc.update_topology(topology, traces, [('bar', 'stem')]) + assert topology.adjacency == {('bar', 'stem'): FAULTED} + + +class FakeLayer: + def __init__(self, features=()): + self._features = list(features) + + def getFeatures(self): + return iter(self._features) + + +def test_missing_layer(): + with pytest.raises(calc.FaultTopologyError, match='fault layer'): + calc.calculate_fault_topology(None, 'name', object()) + + +def test_missing_field(): + with pytest.raises(calc.FaultTopologyError, match='field'): + calc.calculate_fault_topology(FakeLayer(), '', object()) + + +def test_empty_layer(): + pytest.importorskip('geopandas') + with pytest.raises(calc.FaultTopologyError, match='no features'): + calc.calculate_fault_topology(FakeLayer(), 'name', object()) + + +def test_no_fault_layer_in_the_data_manager(): + class Manager: + def get_fault_traces(self): + return None + + with pytest.raises(calc.FaultTopologyError): + calc.calculate_from_data_manager(Manager()) + + +def test_result_message(): + text, warn = calc.result_message(calc.TopologyResult(12, [])) + assert text == "Calculated fault topology for 12 pairs." and not warn + text, warn = calc.result_message(calc.TopologyResult(2, [('a', 'b')])) + assert 'a-b' in text and warn