Results: Add study geotiff in plot xy.

src/View/Results/PlotXY.py

@@ -31,6 +31,20 @@ from PyQt5.QtCore import (
 )
 from PyQt5.QtWidgets import QApplication, QTableView
 
+try:
+    import rasterio
+    import rasterio.control
+    import rasterio.crs
+    import rasterio.sample
+    import rasterio.vrt
+    import rasterio._features
+
+    from rasterio.io import MemoryFile
+    _rasterio_loaded = True
+except Exception as e:
+    print(f"Module 'rasterio' is not available: {e}")
+    _rasterio_loaded = False
+
 _translate = QCoreApplication.translate
 
 logger = logging.getLogger()
@@ -52,6 +66,8 @@ class PlotXY(PamhyrPlot):
         self.line_gl = []
         self.overflow = []
 
+        self._plot_img = {}
+
         self._timestamps = parent._timestamps
         self._current_timestamp = max(self._timestamps)
         self._current_reach_id = reach_id
@@ -153,6 +169,7 @@ class PlotXY(PamhyrPlot):
         reach = results.river.reach(self._current_reach_id)
         reaches = results.river.reachs
 
+        self.draw_geotiff()
         self.draw_profiles(reach, reaches)
         self.draw_water_elevation(reach)
         self.draw_water_elevation_max(reach)
@@ -166,6 +183,7 @@ class PlotXY(PamhyrPlot):
         if reach.geometry.number_profiles == 0:
             self._init = False
             return
+
         self.line_xy = []
         # TODO uncomment to draw all the reaches
         # self.draw_other_profiles(reaches)
@@ -306,6 +324,36 @@ class PlotXY(PamhyrPlot):
             alpha=0.7
         )
 
+    def draw_geotiff(self):
+        if not _rasterio_loaded:
+            return
+
+        lst = self._data[0]._study.river._geotiff.lst
+
+        for img in self._plot_img:
+            self._plot_img[img].remove()
+
+        self._plot_img = {}
+
+        for geotiff in lst:
+            memfile = geotiff.memfile
+            if memfile is None:
+                return
+
+            with memfile.open() as gt:
+                img = gt.read()
+                coords = geotiff.coordinates
+
+                self._plot_img[geotiff] = self.canvas.axes.imshow(
+                    img.transpose((1, 2, 0)),
+                    extent=list(coords.values())
+                )
+
+            if not geotiff.is_enabled():
+                self._plot_img[geotiff].set(alpha=0.5)
+
+        self.idle()
+
     def set_reach(self, reach_id):
         self._current_reach_id = reach_id
         self._current_profile_id = 0