Merge branch 'main' of https://github.com/EsriDE/urban-heat-risk-index

CONTRIBUTING.md

@@ -0,0 +1,47 @@
+# How to Contribute
+We welcome contributions from anyone and everyone. Please see our [guidelines for contributing](https://github.com/esri/contributing).
+
+If you have any questions, please reach out to our developer community [GeoDev Germany Questions](https://community.esri.com/t5/geodev-germany-questions/bd-p/geodev-germany-questions).
+
+Here’s how you can get involved:
+
+## Fork the Repository
+- Click the “Fork” button at the top right of the repository page to create a copy of the repository in your GitHub account.
+
+## Clone Your Fork
+- Clone your forked repository to your local machine using:
+
+```
+git clone https://github.com/EsriDE/urban-heat-risk-index.git
+```
+
+## Create a Branch
+- Create a new branch for your feature or bug fix:
+
+```
+git checkout -b feature/your-feature-name
+```
+
+## Make Changes
+- Make your changes in the codebase. Ensure your code follows the project’s coding standards.
+
+## Commit Your Changes
+- Commit your changes with a descriptive commit message:
+
+```
+git add .
+git commit -m "Add feature: your feature description"
+```
+
+## Push to Your Fork
+- Push your changes to your forked repository:
+
+```
+git push origin feature/your-feature-name
+```
+
+## Create a Pull Request
+- Go to the original repository on GitHub and click the “New Pull Request” button. Select your branch and submit the pull request for review.
+
+Feedback and Support
+If you have any questions or need help, feel free to open an issue or contact us at [GeoDev Germany Questions](https://community.esri.com/t5/geodev-germany-questions/bd-p/geodev-germany-questions).

README.md

@@ -24,6 +24,9 @@ We gather data from various sources, including satellite imagery, weather statio
 
 Using ArcGIS, we apply spatial analysis techniques to identify patterns and correlations between urban features and temperature variations. This helps us pinpoint specific areas that are most susceptible to the UHI effect.
 
+![Screenshot ModelBuilder Heat Risk Index](https://raw.githubusercontent.com/EsriDE/urban-heat-risk-index/main/doc/img/HRI.svg)
+*Screenshot: ModelBuilder Heat Risk Index*
+
 ### Visualization
 
 The results are visualized through interactive maps and dashboards, making it easy for stakeholders to understand the extent and impact of UHIs in their regions.
@@ -58,12 +61,12 @@ The results are visualized through interactive maps and dashboards, making it ea
 
 ## Contributing
 
-We welcome contributions from the community.
+We welcome contributions from anyone and everyone. Please see our [guidelines for contributing](CONTRIBUTING.md).
 
 ## License
 
 This project is licensed under the Apache V2 License - see the [LICENSE](LICENSE) file for details.
 
 ## Acknowledgments
 
-We would like to thank the ArcGIS community and all the contributors who have helped make this project possible.
+We would like to thank the ArcGIS community and all the contributors who have helped make this project possible.

spatial-data-science/scripts/Heat Risk Index.HRITool.pyt.xml

@@ -0,0 +1,2 @@
+<?xml version="1.0" encoding="UTF-8"?>
+<metadata xml:lang="de"><Esri><CreaDate>20240918</CreaDate><CreaTime>15024700</CreaTime><ArcGISFormat>1.0</ArcGISFormat><SyncOnce>TRUE</SyncOnce></Esri></metadata>

spatial-data-science/scripts/Heat Risk Index.pyt

@@ -0,0 +1,53 @@
+# -*- coding: utf-8 -*-
+
+import arcpy
+from heat_risk_index import initialize_arcpy
+
+
+class Toolbox:
+    def __init__(self):
+        """Define the toolbox (the name of the toolbox is the name of the
+        .pyt file)."""
+        self.label = "HRIToolbox"
+        self.alias = "hritoolbox"
+
+        # List of tool classes associated with this toolbox
+        self.tools = [HRITool]
+
+
+class HRITool:
+    def __init__(self):
+        """Define the tool (tool name is the name of the class)."""
+        self.label = "HRI"
+        self.description = ""
+        initialize_arcpy() 
+
+
+    def getParameterInfo(self):
+        """Define the tool parameters."""
+        params = None
+        return params
+
+    def isLicensed(self):
+        """Set whether the tool is licensed to execute."""
+        return True
+
+    def updateParameters(self, parameters):
+        """Modify the values and properties of parameters before internal
+        validation is performed.  This method is called whenever a parameter
+        has been changed."""
+        return
+
+    def updateMessages(self, parameters):
+        """Modify the messages created by internal validation for each tool
+        parameter. This method is called after internal validation."""
+        return
+
+    def execute(self, parameters, messages):
+        """The source code of the tool."""
+        return
+
+    def postExecute(self, parameters):
+        """This method takes place after outputs are processed and
+        added to the display."""
+        return

spatial-data-science/scripts/Heat Risk Index.pyt.xml

@@ -0,0 +1,2 @@
+<?xml version="1.0" encoding="UTF-8"?>
+<metadata xml:lang="de"><Esri><CreaDate>20240918</CreaDate><CreaTime>14591300</CreaTime><ArcGISFormat>1.0</ArcGISFormat><SyncOnce>TRUE</SyncOnce><ModDate>20240918</ModDate><ModTime>150839</ModTime></Esri><toolbox name="Heat Risk Index" alias="hritoolbox"><arcToolboxHelpPath>c:\program files\arcgis\pro\Resources\Help\gp</arcToolboxHelpPath><toolsets/></toolbox><dataIdInfo><idCitation><resTitle>Heat Risk Index</resTitle></idCitation></dataIdInfo><distInfo><distributor><distorFormat><formatName>ArcToolbox Toolbox</formatName></distorFormat></distributor></distInfo></metadata>

spatial-data-science/scripts/heat_risk_index.py

@@ -0,0 +1,160 @@
+import arcpy
+from arcpy.sa import *
+from sys import argv
+
+def initialize_arcpy():
+    """ Initialize arcpy settings and check necessary extensions. """
+    arcpy.env.overwriteOutput = True
+    arcpy.CheckOutExtension("3D")
+    arcpy.CheckOutExtension("spatial")
+    arcpy.CheckOutExtension("ImageExt")
+    arcpy.CheckOutExtension("ImageAnalyst")
+
+def generate_tessellation(output_feature_class, extent, size, spatial_ref):
+    """ Generate a tessellation grid. """
+    arcpy.management.GenerateTessellation(
+        Output_Feature_Class=output_feature_class,
+        Extent=extent,
+        Size=size,
+        Spatial_Reference=spatial_ref
+    )
+
+def spatial_join(target_features, join_features, out_feature_class):
+    """ Perform a spatial join between two feature sets. """
+    arcpy.analysis.SpatialJoin(
+        target_features=target_features,
+        join_features=join_features,
+        out_feature_class=out_feature_class,
+        match_option="LARGEST_OVERLAP"
+    )
+
+def zonal_statistics(in_zone_data, zone_field, in_value_raster, out_table, statistics_type="MAXIMUM"):
+    """ Calculate zonal statistics as a table. """
+    arcpy.sa.ZonalStatisticsAsTable(
+        in_zone_data,
+        zone_field,
+        in_value_raster,
+        out_table,
+        "DATA",
+        statistics_type
+    )
+
+def copy_raster(in_raster, out_rasterdataset, raster_format="TIFF"):
+    """ Copy a raster to a new dataset. """
+    arcpy.management.CopyRaster(in_raster, out_rasterdataset, format=raster_format)
+
+def reclassify_raster(in_raster, remap, out_raster):
+    """ Reclassify a raster based on value ranges. """
+    reclass_raster = arcpy.sa.Reclassify(in_raster, "Value", remap)
+    reclass_raster.save(out_raster)
+
+def join_field(in_data, in_field, join_table, join_field, fields):
+    """ Join fields from one table to another. """
+    arcpy.management.JoinField(
+        in_data=in_data,
+        in_field=in_field,
+        join_table=join_table,
+        join_field=join_field,
+        fields=fields
+    )
+
+def calculate_field(in_table, field, expression, field_type="FLOAT"):
+    """ Calculate a new field based on an expression. """
+    arcpy.management.CalculateField(
+        in_table=in_table,
+        field=field,
+        expression=expression,
+        field_type=field_type
+    )
+
+def standardize_field(in_table, fields, method="MIN-MAX", min_value=1, max_value=5):
+    """ Standardize fields using the given method. """
+    arcpy.management.StandardizeField(
+        in_table=in_table,
+        fields=fields,
+        method=method,
+        min_value=min_value,
+        max_value=max_value
+    )
+
+def delete_features(in_features):
+    """ Delete features from a feature class or layer. """
+    arcpy.management.DeleteFeatures(in_features=in_features)
+
+def hri_main(land_cover, zensus_2022):
+    """ Main function to execute the HRI analysis. """
+    initialize_arcpy()
+
+    # Generate tessellation
+    tessellation_output = "HRI_Hexagone"
+    generate_tessellation(
+        output_feature_class=tessellation_output, 
+        extent="781745.292120143 6556576.21979931 802689.19726414 6581479.0533047",
+        size="1500 SquareMeters",
+        spatial_ref="PROJCS[\"WGS_1984_Web_Mercator_Auxiliary_Sphere\",GEOGCS[\"GCS_WGS_1984\",...]"
+    )
+
+    # Spatial join
+    spatial_join_output = "HRI_Hexagone_SpatialJoin1"
+    spatial_join(tessellation_output, zensus_2022, spatial_join_output)
+
+    # Zonal statistics for surface temperature
+    surf_temp_max = "surf_temp_max"
+    zonal_statistics(spatial_join_output, "GRID_ID", "Multispectral Landsat.tif", surf_temp_max)
+
+    # Copy raster for tree canopy
+    copy_raster(land_cover, "tree_canopy.tif")
+
+    # Reclassify tree canopy raster
+    reclassify_raster("tree_canopy.tif", "10 1;20 0;30 0;40 0;50 0;60 0;70 0;80 0;90 0;95 0;100 0", "reclass_tree_canopy.tif")
+
+    # Zonal statistics for reclassified tree canopy
+    tree_canopy_count = "tree_canopy_count"
+    zonal_statistics(spatial_join_output, "GRID_ID", "reclass_tree_canopy.tif", tree_canopy_count, statistics_type="SUM")
+
+    # Join fields for tree canopy statistics
+    join_field(spatial_join_output, "GRID_ID", tree_canopy_count, "GRID_ID", ["SUM"])
+
+    # Calculate percentage tree cover and lacking
+    calculate_field(tree_canopy_count, "PCT_Tree_Cover", "(!SUM! / !COUNT!) * 100")
+    calculate_field(tree_canopy_count, "PCT_Lacking", "100 - !PCT_Tree_Cover!")
+
+    # Join the fields for surface temperature and tree canopy data
+    join_field(spatial_join_output, "GRID_ID", surf_temp_max, "GRID_ID", ["MAX"])
+    join_field(spatial_join_output, "GRID_ID", tree_canopy_count, "GRID_ID", ["PCT_Tree_Cover", "PCT_Lacking"])
+
+    # Copy raster for built-up area
+    copy_raster(land_cover, "built_up_area.tif")
+
+    # Reclassify built-up area raster
+    reclassify_raster("built_up_area.tif", "10 0;20 0;30 0;40 0;50 1;60 0;70 0;80 0;90 0;95 0;100 0", "reclass_built_up_area.tif")
+
+    # Zonal statistics for reclassified built-up area
+    built_up_area_count = "built_up_area_count"
+    zonal_statistics(spatial_join_output, "GRID_ID", "reclass_built_up_area.tif", built_up_area_count, statistics_type="SUM")
+
+    # Join the built-up area statistics
+    join_field(spatial_join_output, "GRID_ID", built_up_area_count, "GRID_ID", ["SUM"])
+    calculate_field(built_up_area_count, "PCT_built_up_area", "(!SUM! / !COUNT!) * 100")
+    join_field(spatial_join_output, "GRID_ID", built_up_area_count, "GRID_ID", ["PCT_built_up_area"])
+
+    # Final HRI calculation using standardized values
+    standardize_field(
+        spatial_join_output,
+        [["Einwohner", "Einwohner_MIN_MAX"], ["MAX", "TEMP_MAX_MIN_MAX"], ["PCT_Lacking", "PCT_Lacking_MIN_MAX"]]
+    )
+    calculate_field(spatial_join_output, "HRI", "Sum($feature.TEMP_MAX_MIN_MAX, $feature.PCT_Lacking_MIN_MAX, $feature.Einwohner_MIN_MAX)", field_type="FLOAT")
+
+    # Select layer by location
+    arcpy.management.SelectLayerByLocation(
+        in_layer=[spatial_join_output],
+        overlap_type="WITHIN",
+        select_features="Ortsteile_Bonn",
+        invert_spatial_relationship="INVERT"
+    )
+
+    # Delete features not in Bonn
+    delete_features(spatial_join_output)
+
+if __name__ == '__main__':
+    hri_main(*argv[1:])