Bedload: Dynamic display for the riverbed over time in results window + correction of import of .ST files in geometry tab when we import a bedload file

src/Model/Geometry/Reach.py

@@ -752,8 +752,11 @@ class Reach(SQLSubModel):
                                 list_profile.append(list_point_profile)
                                 list_point_profile = []
                             else:
-                                line.append("")
-                                list_point_profile.append(line[:3])
+                                # MAGE-8 ST files can append hydraulic data
+                                # after the regular x, y, z and point-name
+                                # fields.  ProfileXYZ still expects exactly
+                                # those first four fields.
+                                list_point_profile.append(line[:4])
                         else:
                             pass
 

src/View/Results/PlotAC.py

@@ -95,7 +95,7 @@ class PlotAC(PamhyrPlot):
 
     def draw_profile(self, reach, profile):
         x = profile.geometry.get_station()
-        z = profile.geometry.z()
+        z = self.get_bottom_elevation(profile)
 
         self.line_rk, = self.canvas.axes.plot(
             x, z,
@@ -110,7 +110,7 @@ class PlotAC(PamhyrPlot):
         result = self.results[self._current_res_id]
 
         x = profile.geometry.get_station()
-        z = profile.geometry.z()
+        z = self.get_bottom_elevation(profile)
         rk = reach.geometry.get_rk()
         water_z = result.get("table")["Z"][
             result.get_timestamp_id(
@@ -139,7 +139,7 @@ class PlotAC(PamhyrPlot):
         result = self.results[self._current_res_id]
 
         x = profile.geometry.get_station()
-        z = profile.geometry.z()
+        z = self.get_bottom_elevation(profile)
         rk = reach.geometry.get_rk()
 
         water_z = np.max(
@@ -162,7 +162,7 @@ class PlotAC(PamhyrPlot):
             return
 
         station = profile.geometry.get_station()
-        elevation = profile.geometry.z()
+        elevation = self.get_bottom_elevation(profile)
         gl = profile.geometry.names()
 
         self.annotation = []
@@ -239,9 +239,11 @@ class PlotAC(PamhyrPlot):
             self._current_reach_id)
         profile = reach.profile(self._current_profile_id)
         x = profile.geometry.get_station()
-        z = profile.geometry.z()
+        z = self.get_bottom_elevation(profile)
 
+        self.update_river_bottom(reach, profile, x, z)
         self.update_water(reach, profile, x, z)
+        self.update_gl()
         self.update_idle()
 
     def update(self):
@@ -252,7 +254,7 @@ class PlotAC(PamhyrPlot):
                     .river.reach(self._current_reach_id)
         profile = reach.profile(self._current_profile_id)
         x = profile.geometry.get_station()
-        z = profile.geometry.z()
+        z = self.get_bottom_elevation(profile)
 
         self.update_river_bottom(reach, profile, x, z)
         self.update_water(reach, profile, x, z)
@@ -279,6 +281,69 @@ class PlotAC(PamhyrPlot):
     def update_river_bottom(self, reach, profile, x, z):
         self.line_rk.set_data(x, z)
 
+    def get_bottom_elevation(self, profile):
+        """Return the cross-section elevation at the selected timestamp."""
+        initial_z = np.asarray(profile.geometry.z(), dtype=float)
+        result = self.results[self._current_res_id]
+
+        try:
+            sediment_timestamps = sorted(result.get("sediment_timestamps"))
+        except (KeyError, TypeError):
+            return initial_z
+
+        if not sediment_timestamps:
+            return initial_z
+
+        initial_timestamp = sediment_timestamps[0]
+        current_timestamp = min(
+            sediment_timestamps,
+            key=lambda timestamp: abs(timestamp - self._current_timestamp)
+        )
+        initial_layers = profile.get_ts_key(initial_timestamp, "sl")
+        current_layers = profile.get_ts_key(current_timestamp, "sl")
+
+        if not initial_layers or not current_layers:
+            return initial_z
+
+        initial_thickness = self._get_sediment_thickness(initial_layers)
+        current_thickness = self._get_sediment_thickness(current_layers)
+        point_count = len(initial_z)
+        wet_mask = self._get_wet_mask(profile, initial_z)
+
+        if wet_mask is None:
+            return initial_z
+
+        if len(initial_thickness) == len(current_thickness) == 1:
+            elevation_delta = current_thickness[0] - initial_thickness[0]
+            return initial_z + np.where(wet_mask, elevation_delta, 0.0)
+
+        if (len(initial_thickness) == len(current_thickness) == point_count):
+            elevation_delta = np.asarray(current_thickness) \
+                - np.asarray(initial_thickness)
+            return initial_z + np.where(wet_mask, elevation_delta, 0.0)
+
+        return initial_z
+
+    def _get_wet_mask(self, profile, elevation):
+        result = self.results[self._current_res_id]
+
+        try:
+            water_elevation = result.get("table")["Z"][
+                result.get_timestamp_id(self._current_timestamp),
+                profile.global_index
+            ]
+        except (KeyError, IndexError, TypeError):
+            return None
+
+        return elevation <= water_elevation
+
+    @staticmethod
+    def _get_sediment_thickness(section_layers):
+        return [
+            sum(float(layer[0]) for layer in point_layers)
+            for point_layers in section_layers
+        ]
+
     def update_water(self, reach, profile, x, z):
         result = self.results[self._current_res_id]
         water_z = result.get("table")["Z"][

src/View/ui/MainWindow.ui

@@ -188,7 +188,7 @@
    </widget>
    <widget class="QMenu" name="menuSediment">
     <property name="title">
-     <string>&amp;Sediment</string>
+     <string>&amp;Bedload</string>
     </property>
     <addaction name="action_menu_boundary_conditions_sediment"/>
     <addaction name="action_menu_sediment_layers"/>