Index score integration

geest/core/crs_converter.py

@@ -1,8 +1,8 @@
 from qgis.core import (
-    QgsCoordinateTransform,
     QgsCoordinateReferenceSystem,
     QgsProcessingFeedback,
-    QgsVectorLayer,
+    QgsMessageLog,
+    Qgis,
 )
 from qgis import processing
 
@@ -28,8 +28,10 @@ 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()}"
+            QgsMessageLog.logMessage(
+                f"Converting layer from {current_crs.authid()} to {target_crs.authid()}",
+                tag="Geest",
+                level=Qgis.Info,
             )
 
             layer = processing.run(
@@ -41,8 +43,18 @@ def convert_to_crs(self, target_crs_epsg):
                 },
                 feedback=QgsProcessingFeedback(),
             )["OUTPUT"]
-            print(f"Layer successfully converted to {target_crs.authid()}")
+            QgsMessageLog.logMessage(
+                f"Layer successfully converted to {target_crs.authid()}",
+                tag="Geest",
+                level=Qgis.Info,
+            )
+
             return layer
         else:
-            print(f"Layer is already in the target CRS: {target_crs.authid()}")
+            QgsMessageLog.logMessage(
+                f"Layer is already in the target CRS: {target_crs.authid()}",
+                tag="Geest",
+                level=Qgis.Info,
+            )
+
             return self.layer

geest/core/study_area.py

@@ -49,7 +49,7 @@ def __init__(
         self.gpkg_path: str = os.path.join(
             self.working_dir, "study_area", "study_area.gpkg"
         )
-        self.counter: int = 0
+
         # Remove the GeoPackage if it already exists to start with a clean state
         if os.path.exists(self.gpkg_path):
             try:
@@ -226,6 +226,9 @@ def process_singlepart_geometry(
         bbox: QgsRectangle = self.grid_aligned_bbox(geom.boundingBox())
 
         # Create a feature for the aligned bounding box
+        study_area_feature: QgsFeature = QgsFeature()
+        study_area_feature.setGeometry(QgsGeometry.fromRect(bbox))
+        study_area_feature.setAttributes([area_name])
         # Always save the study area bounding boxes regardless of mode
         self.save_to_geopackage(
             layer_name="study_area_bboxes",
@@ -241,10 +244,14 @@ def process_singlepart_geometry(
         geom.transform(transform)
 
         # Create a feature for the original part
+        study_area_polygon: QgsFeature = QgsFeature()
+        study_area_polygon.setGeometry(geom)
+        study_area_polygon.setAttributes([area_name])
         # Always save the study area bounding boxes regardless of mode
         self.save_to_geopackage(
             layer_name="study_area_polygons", geom=geom, area_name=normalized_name
         )
+
         # Process the geometry based on the selected mode
         if self.mode == "vector":
             QgsMessageLog.logMessage(
@@ -322,16 +329,15 @@ def grid_aligned_bbox(self, bbox: QgsRectangle) -> QgsRectangle:
             * 100
         )
 
-        y_min -= 100  # Offset by 100m to ensure the grid covers the entire geometry
-        y_max += 100  # Offset by 100m to ensure the grid covers the entire geometry
-        x_min -= 100  # Offset by 100m to ensure the grid covers the entire geometry
-        x_max += 100  # Offset by 100m to ensure the grid covers the entire geometry
-
         # Return the aligned bbox in the output CRS
         return QgsRectangle(x_min, y_min, x_max, y_max)
 
     def save_to_geopackage(
-        self, layer_name: str, geom: QgsGeometry, area_name: str
+        self,
+        features: List[QgsFeature],
+        layer_name: str,
+        fields: List[QgsField],
+        geometry_type: QgsWkbTypes,
     ) -> None:
         """
         Save features to GeoPackage. Create or append the layer as necessary.
@@ -380,7 +386,9 @@ def append_to_layer(
                 level=Qgis.Critical,
             )
 
-    def create_layer_if_not_exists(self, layer_name: str) -> None:
+    def create_layer_if_not_exists(
+        self, layer_name: str, fields: List[QgsField], geometry_type: QgsWkbTypes
+    ) -> None:
         """
         Create a new layer in the GeoPackage if it doesn't already exist.
 
@@ -425,6 +433,11 @@ def create_layer_if_not_exists(self, layer_name: str) -> None:
                     QgsVectorFileWriter.CreateOrOverwriteLayer
                 )
 
+            # Convert list of QgsField objects to QgsFields object
+            qgs_fields = QgsFields()
+            for field in fields:
+                qgs_fields.append(field)
+
             # Create a new GeoPackage layer
             QgsVectorFileWriter.create(
                 fileName=self.gpkg_path,

geest/core/workflows/default_index_score_workflow.py

@@ -1,4 +1,5 @@
 import os
+import glob
 from qgis.core import (
     QgsMessageLog,
     Qgis,
@@ -8,6 +9,8 @@
     QgsField,
     QgsGeometry,
     QgsRectangle,
+    QgsRasterLayer,
+    QgsProject,
 )
 from qgis.PyQt.QtCore import QVariant
 import processing  # QGIS processing toolbox
@@ -46,9 +49,23 @@ def execute(self):
             "----------------------------------", tag="Geest", level=Qgis.Info
         )
 
+        self.workflow_directory = self._create_workflow_directory(
+            "contextual",
+            self.attributes["ID"].lower(),
+        )
+
         # 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():
+            if (
+                self.feedback.isCanceled()
+            ):  # Check for cancellation before each major step
+                QgsMessageLog.logMessage(
+                    "Workflow canceled before processing feature.",
+                    tag="Geest",
+                    level=Qgis.Warning,
+                )
+                return False
             geom = feature.geometry()  # todo this shoudl come from the areas layer
             aligned_box = geom
             mask_name = f"bbox_{feature.id()}"
@@ -59,25 +76,10 @@ def execute(self):
                 index_score=index_score,
             )
         # TODO Jeff copy create_raster_vrt from study_area.py
-
-        steps = 10
-        for i in range(steps):
-            if self.feedback.isCanceled():
-                QgsMessageLog.logMessage(
-                    "Dont use workflow canceled.", tag="Geest", level=Qgis.Warning
-                )
-                return False
-
-            # Simulate progress and work
-            self.attributes["progress"] = f"Dont use workflow Step {i + 1} completed"
-            self.feedback.setProgress(
-                (i + 1) / steps * 100
-            )  # Report progress in percentage
-            QgsMessageLog.logMessage(
-                f"Assigning index score: {self.attributes['Default Index Score']}",
-                tag="Geest",
-                level=Qgis.Info,
-            )
+        # 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")
+        )
 
         self.attributes["result"] = "Use Default Index Score Workflow Completed"
         QgsMessageLog.logMessage(
@@ -101,8 +103,17 @@ def create_raster(
         :param aligned_box: Aligned bounding box geometry for the geometry.
         :param mask_name: Name for the output raster file.
         """
+        if self.feedback.isCanceled():  # Check for cancellation before starting
+            QgsMessageLog.logMessage(
+                "Workflow canceled before creating raster.",
+                tag="Geest",
+                level=Qgis.Warning,
+            )
+            return
+
         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(
@@ -128,7 +139,7 @@ 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,
@@ -148,3 +159,71 @@ def create_raster(
         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.
+        """
+        if self.feedback.isCanceled():  # Check for cancellation before starting
+            QgsMessageLog.logMessage(
+                "Workflow canceled before creating VRT.",
+                tag="Geest",
+                level=Qgis.Warning,
+            )
+            return
+
+        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
+        )
+
+        layer_id = self.attributes["ID"].replace("_", " ")
+
+        # Add the VRT to the QGIS map
+        vrt_layer = QgsRasterLayer(vrt_filepath, f"Combined Mask VRT ({layer_id})")
+
+        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/dont_use_workflow.py

@@ -20,6 +20,12 @@ def execute(self):
         """
         Executes the workflow, reporting progress through the feedback object and checking for cancellation.
         """
+        if self.feedback.isCanceled():
+            QgsMessageLog.logMessage(
+                "Dont use workflow canceled.", tag="Geest", level=Qgis.Warning
+            )
+            return False
+
         QgsMessageLog.logMessage("Executing 'dont use'", tag="Geest", level=Qgis.Info)
 
         steps = 10

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
 
 
@@ -49,7 +49,7 @@ def execute(self) -> bool:
         """
         pass
 
-    def _create_workflow_directory(self) -> str:
+    def _create_workflow_directory(self, *subdirs: str) -> str:
         """
         Creates the directory for this workflow if it doesn't already exist.
         It will be in the scheme of working_dir/dimension/factor/indicator
@@ -58,9 +58,7 @@ def _create_workflow_directory(self) -> str:
         """
         workflow_dir = os.path.join(
             self.working_directory,
-            "contextual",
-            "workplace_discrimination",
-            "wbl_2024_workplace_index_score",
+            *subdirs,
         )
         if not os.path.exists(workflow_dir):
             try: