Index score integration

geest/core/crs_converter.py

@@ -1,6 +1,12 @@
-from qgis.core import QgsCoordinateTransform, QgsCoordinateReferenceSystem, QgsProcessingFeedback, QgsVectorLayer
+from qgis.core import (
+    QgsCoordinateTransform,
+    QgsCoordinateReferenceSystem,
+    QgsProcessingFeedback,
+    QgsVectorLayer,
+)
 from qgis import processing
 
+
 class CRSConverter:
     def __init__(self, layer):
         """
@@ -22,7 +28,9 @@ def convert_to_crs(self, target_crs_epsg):
 
         # Check if the current CRS is the same as the target CRS
         if current_crs != target_crs:
-            print(f"Converting layer from {current_crs.authid()} to {target_crs.authid()}")
+            print(
+                f"Converting layer from {current_crs.authid()} to {target_crs.authid()}"
+            )
 
             layer = processing.run(
                 "native:reprojectlayer",

geest/core/workflow_factory.py

@@ -18,7 +18,7 @@ def create_workflow(self, attributes, feedback: QgsFeedback):
         if analysis_mode == "Spatial Analysis":
             return RasterLayerWorkflow(attributes, feedback)
         elif analysis_mode == "Use Default Index Score":
-            return DefaultIndexScoreWorkflow(attributes, feedback)        
+            return DefaultIndexScoreWorkflow(attributes, feedback)
         elif analysis_mode == "Don’t Use":
             return DontUseWorkflow(attributes, feedback)
         elif analysis_mode == "Temporal Analysis":

geest/core/workflow_job.py

@@ -34,7 +34,9 @@ def run(self) -> bool:
         """
         if not self._workflow:
             QgsMessageLog.logMessage(
-                f"Error: No workflow assigned to {self.description()}", tag="Geest", level=Qgis.Critical,
+                f"Error: No workflow assigned to {self.description()}",
+                tag="Geest",
+                level=Qgis.Critical,
             )
             return False
 
@@ -47,12 +49,16 @@ def run(self) -> bool:
 
             if result:
                 QgsMessageLog.logMessage(
-                    f"Workflow {self.description()} completed.", tag="Geest", level=Qgis.Info
+                    f"Workflow {self.description()} completed.",
+                    tag="Geest",
+                    level=Qgis.Info,
                 )
                 return True
             else:
                 QgsMessageLog.logMessage(
-                    f"Workflow {self.description()} did not complete successfully.", tag="Geest", level=Qgis.Info
+                    f"Workflow {self.description()} did not complete successfully.",
+                    tag="Geest",
+                    level=Qgis.Info,
                 )
                 return False
 

geest/core/workflow_queue_manager.py

@@ -63,7 +63,10 @@ def on_processing_completed(self, success: bool) -> None:
         """
         if success:
             QgsMessageLog.logMessage(
-                "All workflow tasks completed successfully.", tag="Geest", level=Qgis.Info)
+                "All workflow tasks completed successfully.",
+                tag="Geest",
+                level=Qgis.Info,
+            )
         else:
             QgsMessageLog.logMessage(
                 "Workflow processing was canceled.", tag="Geest", level=Qgis.Info

geest/core/workflows/default_index_score_workflow.py

@@ -1,13 +1,17 @@
 import os
+import glob
 from qgis.core import (
-    QgsMessageLog, 
-    Qgis, 
-    QgsFeedback, 
-    QgsFeature, 
-    QgsVectorLayer, 
-    QgsField, 
-    QgsGeometry, 
-    QgsRectangle)
+    QgsMessageLog,
+    Qgis,
+    QgsFeedback,
+    QgsFeature,
+    QgsVectorLayer,
+    QgsField,
+    QgsGeometry,
+    QgsRectangle,
+    QgsRasterLayer,
+    QgsProject,
+)
 from qgis.PyQt.QtCore import QVariant
 import processing  # QGIS processing toolbox
 from .workflow_base import WorkflowBase
@@ -31,24 +35,37 @@ def execute(self):
         Executes the workflow, reporting progress through the feedback object and checking for cancellation.
         """
 
-        QgsMessageLog.logMessage("Executing Use Default Index Score", tag="Geest", level=Qgis.Info)
-        QgsMessageLog.logMessage("----------------------------------", tag="Geest", level=Qgis.Info)
+        QgsMessageLog.logMessage(
+            "Executing Use Default Index Score", tag="Geest", level=Qgis.Info
+        )
+        QgsMessageLog.logMessage(
+            "----------------------------------", tag="Geest", level=Qgis.Info
+        )
         for item in self.attributes.items():
-            QgsMessageLog.logMessage(f"{item[0]}: {item[1]}", tag="Geest", level=Qgis.Info)
-        QgsMessageLog.logMessage("----------------------------------", tag="Geest", level=Qgis.Info)
+            QgsMessageLog.logMessage(
+                f"{item[0]}: {item[1]}", tag="Geest", level=Qgis.Info
+            )
+        QgsMessageLog.logMessage(
+            "----------------------------------", tag="Geest", level=Qgis.Info
+        )
 
         # loop through self.bboxes_layer and the self.areas_layer  and create a raster mask for each feature
         index_score = self.attributes["Default Index Score"]
         for feature in self.bboxes_layer.getFeatures():
-            geom = feature.geometry() # todo this shoudl come from the areas layer
+            geom = feature.geometry()  # todo this shoudl come from the areas layer
             aligned_box = geom
             mask_name = f"bbox_{feature.id()}"
             self.create_raster(
-                geom=geom, 
-                aligned_box=aligned_box, 
+                geom=geom,
+                aligned_box=aligned_box,
                 mask_name=mask_name,
-                index_score=index_score)
+                index_score=index_score,
+            )
         # TODO Jeff copy create_raster_vrt from study_area.py
+        # Create and add the VRT of all generated raster masks if in raster mode
+        self.create_raster_vrt(
+            output_vrt_name=os.path.join(self.workflow_directory, "combined_mask.vrt")
+        )
 
         steps = 10
         for i in range(steps):
@@ -64,31 +81,36 @@ def execute(self):
                 (i + 1) / steps * 100
             )  # Report progress in percentage
             QgsMessageLog.logMessage(
-                f"Assigning index score: {self.attributes['Default Index Score']}", 
-                tag="Geest", level=Qgis.Info)
+                f"Assigning index score: {self.attributes['Default Index Score']}",
+                tag="Geest",
+                level=Qgis.Info,
+            )
 
         self.attributes["result"] = "Use Default Index Score Workflow Completed"
         QgsMessageLog.logMessage(
-            "Use Default Index Score workflow workflow completed", tag="Geest", level=Qgis.Info
+            "Use Default Index Score workflow workflow completed",
+            tag="Geest",
+            level=Qgis.Info,
         )
         return True
 
-
     def create_raster(
-            self, 
-            geom: QgsGeometry, 
-            aligned_box: QgsGeometry, 
-            mask_name: str,
-            index_score: float) -> None:
+        self,
+        geom: QgsGeometry,
+        aligned_box: QgsGeometry,
+        mask_name: str,
+        index_score: float,
+    ) -> None:
         """
         Creates a byte raster mask for a single geometry.
 
         :param geom: Geometry to be rasterized.
         :param aligned_box: Aligned bounding box geometry for the geometry.
         :param mask_name: Name for the output raster file.
-        """    
+        """
         aligned_box = QgsRectangle(aligned_box.boundingBox())
         mask_filepath = os.path.join(self.workflow_directory, f"{mask_name}.tif")
+        index_score = (self.attributes["Default Index Score"] / 100) * 5
 
         # Create a memory layer to hold the geometry
         temp_layer = QgsVectorLayer(
@@ -114,22 +136,81 @@ def create_raster(
         params = {
             "INPUT": temp_layer,
             "FIELD": None,
-            "BURN": 78, # todo Jeff put on likert scale properly
+            "BURN": index_score,  # todo Jeff put on likert scale properly
             "USE_Z": False,
             "UNITS": 1,
             "WIDTH": x_res,
             "HEIGHT": y_res,
             "EXTENT": f"{aligned_box.xMinimum()},{aligned_box.xMaximum()},"
-                    f"{aligned_box.yMinimum()},{aligned_box.yMaximum()}",  # Extent of the aligned bbox
+            f"{aligned_box.yMinimum()},{aligned_box.yMaximum()}",  # Extent of the aligned bbox
             "NODATA": 0,
             "OPTIONS": "",
-            "DATA_TYPE": 0, # byte
+            "DATA_TYPE": 0,  # byte
             "INIT": None,
             "INVERT": False,
             "EXTRA": "",
             "OUTPUT": mask_filepath,
         }
         # Run the rasterize algorithm
         processing.run("gdal:rasterize", params)
-        QgsMessageLog.logMessage(f"Created raster mask: {mask_filepath}", tag="Geest", level=Qgis.Info)
+        QgsMessageLog.logMessage(
+            f"Created raster mask: {mask_filepath}", tag="Geest", level=Qgis.Info
+        )
 
+    def create_raster_vrt(self, output_vrt_name: str = "combined_mask.vrt") -> None:
+        """
+        Creates a VRT file from all generated raster masks and adds it to the QGIS map.
+
+        :param output_vrt_name: The name of the VRT file to create.
+        """
+        QgsMessageLog.logMessage(
+            f"Creating VRT of masks '{output_vrt_name}' layer to the map.",
+            tag="Geest",
+            level=Qgis.Info,
+        )
+        # Directory containing raster masks
+        raster_dir = os.path.dirname(output_vrt_name)
+        raster_files = glob.glob(os.path.join(raster_dir, "*.tif"))
+
+        if not raster_files:
+            QgsMessageLog.logMessage(
+                "No raster masks found to combine into VRT.",
+                tag="Geest",
+                level=Qgis.Warning,
+            )
+            return
+
+        vrt_filepath = os.path.join(raster_dir, output_vrt_name)
+
+        # Define the VRT parameters
+        params = {
+            "INPUT": raster_files,
+            "RESOLUTION": 0,  # Use highest resolution among input files
+            "SEPARATE": False,  # Combine all input rasters as a single band
+            "OUTPUT": vrt_filepath,
+            "PROJ_DIFFERENCE": False,
+            "ADD_ALPHA": False,
+            "ASSIGN_CRS": None,
+            "RESAMPLING": 0,
+            "SRC_NODATA": "0",
+            "EXTRA": "",
+        }
+
+        # Run the gdal:buildvrt processing algorithm to create the VRT
+        processing.run("gdal:buildvirtualraster", params)
+        QgsMessageLog.logMessage(
+            f"Created VRT: {vrt_filepath}", tag="Geest", level=Qgis.Info
+        )
+
+        # Add the VRT to the QGIS map
+        vrt_layer = QgsRasterLayer(vrt_filepath, "Combined Mask VRT")
+
+        if vrt_layer.isValid():
+            QgsProject.instance().addMapLayer(vrt_layer)
+            QgsMessageLog.logMessage(
+                "Added VRT layer to the map.", tag="Geest", level=Qgis.Info
+            )
+        else:
+            QgsMessageLog.logMessage(
+                "Failed to add VRT layer to the map.", tag="Geest", level=Qgis.Critical
+            )

geest/core/workflows/workflow_base.py

@@ -1,6 +1,6 @@
 import os
 from abc import ABC, abstractmethod
-from qgis.core import QgsFeedback, QgsVectorLayer
+from qgis.core import QgsFeedback, QgsVectorLayer, QgsMessageLog, Qgis
 from qgis.PyQt.QtCore import QSettings
 
 
@@ -26,14 +26,21 @@ def __init__(self, attributes: dict, feedback: QgsFeedback):
             raise ValueError("Working directory not set.")
         # This is the lower level directory for this workflow
         self.workflow_directory = self._create_workflow_directory()
-        self.gpkg_path: str = os.path.join(self.working_directory, "study_area", "study_area.gpkg")
+        self.gpkg_path: str = os.path.join(
+            self.working_directory, "study_area", "study_area.gpkg"
+        )
         if not os.path.exists(self.gpkg_path):
             raise ValueError(f"Study area geopackage not found at {self.gpkg_path}.")
-        self.bboxes_layer = QgsVectorLayer(f"{self.gpkg_path}|layername=study_area_bboxes", "study_area_bboxes", "ogr")
-        self.areas_layer = QgsVectorLayer(f"{self.gpkg_path}|layername=study_area_polygons", "study_area_polygons", "ogr")
+        self.bboxes_layer = QgsVectorLayer(
+            f"{self.gpkg_path}|layername=study_area_bboxes", "study_area_bboxes", "ogr"
+        )
+        self.areas_layer = QgsVectorLayer(
+            f"{self.gpkg_path}|layername=study_area_polygons",
+            "study_area_polygons",
+            "ogr",
+        )
         self.output_crs = self.bboxes_layer.crs()
 
-
     @abstractmethod
     def execute(self) -> bool:
         """
@@ -50,12 +57,21 @@ def _create_workflow_directory(self) -> str:
         :return: The path to the workflow directory
         """
         workflow_dir = os.path.join(
-            self.working_directory, "contextual", "workplace_discrimination", "wbl_2024_workplace_index_score")
+            self.working_directory,
+            "contextual",
+            "workplace_discrimination",
+            "wbl_2024_workplace_index_score",
+        )
         if not os.path.exists(workflow_dir):
             try:
                 os.makedirs(workflow_dir)
-                QgsMessageLog.logMessage(f"Created study area directory: {workflow_dir}", tag="Geest", level=Qgis.Info)
+                QgsMessageLog.logMessage(
+                    f"Created study area directory: {workflow_dir}",
+                    tag="Geest",
+                    level=Qgis.Info,
+                )
             except Exception as e:
-                QgsMessageLog.logMessage(f"Error creating directory: {e}", tag="Geest", level=Qgis.Critical)
+                QgsMessageLog.logMessage(
+                    f"Error creating directory: {e}", tag="Geest", level=Qgis.Critical
+                )
         return workflow_dir
-