Merge pull request #977 from ebocher/save_urban_cool_areas

Improve urban_sprawl_areas and urban_cool_areas computation
orbisgis · Jun 10, 2024 · 7910c41 · 7910c41
2 parents d5d7c5b + 1e61e2a
commit 7910c41
Show file tree

Hide file tree

Showing 62 changed files with 1,239 additions and 1,213 deletions.
diff --git a/bdtopo/pom.xml b/bdtopo/pom.xml
@@ -1,4 +1,5 @@
-<project xmlns="http://maven.apache.org/POM/4.0.0" xmlns:xsi="http://www.w3.org/2001/XMLSchema-instance" xsi:schemaLocation="http://maven.apache.org/POM/4.0.0 http://maven.apache.org/xsd/maven-4.0.0.xsd">
+<project xmlns="http://maven.apache.org/POM/4.0.0" xmlns:xsi="http://www.w3.org/2001/XMLSchema-instance"
+         xsi:schemaLocation="http://maven.apache.org/POM/4.0.0 http://maven.apache.org/xsd/maven-4.0.0.xsd">
     <parent>
         <groupId>org.orbisgis.geoclimate</groupId>
         <artifactId>geoclimate-parent</artifactId>

diff --git a/bdtopo/src/main/groovy/org/orbisgis/geoclimate/bdtopo/AbstractBDTopoWorkflow.groovy b/bdtopo/src/main/groovy/org/orbisgis/geoclimate/bdtopo/AbstractBDTopoWorkflow.groovy
@@ -28,7 +28,6 @@ import org.orbisgis.data.H2GIS
 import org.orbisgis.data.api.dataset.ITable
 import org.orbisgis.data.jdbc.JdbcDataSource
 import org.orbisgis.geoclimate.Geoindicators
-
 import org.orbisgis.geoclimate.worldpoptools.WorldPopTools
 
 import java.sql.Connection
@@ -291,7 +290,7 @@ abstract class AbstractBDTopoWorkflow extends BDTopoUtils {
      * @param message
      * @throws Exception
      */
-    void saveLogZoneTable(JdbcDataSource dataSource,String databaseFolder, String location, String message) throws Exception {
+    void saveLogZoneTable(JdbcDataSource dataSource, String databaseFolder, String location, String message) throws Exception {
         def logTableZones = postfix("log_zones")
         //Create the table to log on the processed zone
         dataSource.execute("""DROP TABLE IF EXISTS $logTableZones;
@@ -310,7 +309,7 @@ abstract class AbstractBDTopoWorkflow extends BDTopoUtils {
                             '${Geoindicators.version()}',
                             '${Geoindicators.buildNumber()}')""")
         }
-        dataSource.save(logTableZones, databaseFolder+File.separator+"log_zones_"+id_zone+".fgb", true )
+        dataSource.save(logTableZones, databaseFolder + File.separator + "log_zones_" + id_zone + ".fgb", true)
     }
 
     /**
@@ -349,7 +348,7 @@ abstract class AbstractBDTopoWorkflow extends BDTopoUtils {
     def linkDataFromFolder(def inputFolder, def inputWorkflowTableNames,
                            H2GIS h2gis_datasource, def inputSRID) throws Exception {
         def folder = new File(inputFolder)
-        if(!folder.exists()){
+        if (!folder.exists()) {
             throw new Exception("The input folder doesn't exist")
         }
         if (folder.isDirectory()) {
@@ -558,8 +557,8 @@ abstract class AbstractBDTopoWorkflow extends BDTopoUtils {
                                                    "ROAD_FRACTION", "IMPERVIOUS_FRACTION", "UTRF_AREA_FRACTION", "UTRF_FLOOR_AREA_FRACTION", "LCZ_FRACTION", "LCZ_PRIMARY", "FREE_EXTERNAL_FACADE_DENSITY",
                                                    "BUILDING_HEIGHT_WEIGHTED", "BUILDING_SURFACE_DENSITY",
                                                    "BUILDING_HEIGHT_DIST", "FRONTAL_AREA_INDEX", "SEA_LAND_FRACTION", "ASPECT_RATIO",
-                                                   "SVF", "HEIGHT_OF_ROUGHNESS_ELEMENTS", "TERRAIN_ROUGHNESS_CLASS", "SPRAWL_AREAS",
-                                                   "SPRAWL_DISTANCES", "SPRAWL_COOL_DISTANCE"]
+                                                   "SVF", "HEIGHT_OF_ROUGHNESS_ELEMENTS", "TERRAIN_ROUGHNESS_CLASS", "URBAN_SPRAWL_AREAS",
+                                                   "URBAN_SPRAWL_DISTANCES", "URBAN_SPRAWL_COOL_DISTANCE"]
                     def allowedOutputIndicators = allowed_grid_indicators.intersect(list_indicators*.toUpperCase())
                     if (allowedOutputIndicators) {
                         //Update the RSU indicators list according the grid indicators
@@ -650,7 +649,8 @@ abstract class AbstractBDTopoWorkflow extends BDTopoUtils {
                 "road_traffic",
                 "population",
                 "ground_acoustic",
-                "sprawl_areas"]
+                "urban_sprawl_areas",
+                "urban_cool_areas"]
     }
 
 
@@ -860,7 +860,7 @@ abstract class AbstractBDTopoWorkflow extends BDTopoUtils {
                     x_size, y_size, srid, grid_indicators_params.rowCol)
             if (gridTableName) {
                 String rasterizedIndicators = Geoindicators.WorkflowGeoIndicators.rasterizeIndicators(h2gis_datasource, gridTableName,
-                        grid_indicators_params.indicators, grid_indicators_params.lcz_lod,
+                        grid_indicators_params.indicators,
                         results.building, results.road, results.vegetation,
                         results.water, results.impervious,
                         results.rsu_lcz, results.rsu_utrf_area, "", "",
@@ -871,7 +871,10 @@ abstract class AbstractBDTopoWorkflow extends BDTopoUtils {
                     Map sprawl_indic = Geoindicators.WorkflowGeoIndicators.sprawlIndicators(h2gis_datasource, rasterizedIndicators, "id_grid", grid_indicators_params.indicators,
                             Math.max(x_size, y_size).floatValue())
                     if (sprawl_indic) {
-                        results.put("sprawl_areas", sprawl_indic.sprawl_areas)
+                        results.put("urban_sprawl_areas", sprawl_indic.urban_sprawl_areas)
+                        if (sprawl_indic.urban_cool_areas) {
+                            results.put("urban_cool_areas", sprawl_indic.urban_cool_areas)
+                        }
                         results.put("grid_indicators", sprawl_indic.grid_indicators)
                     }
                     info("End computing grid_indicators")

diff --git a/bdtopo/src/main/groovy/org/orbisgis/geoclimate/bdtopo/BDTopo.groovy b/bdtopo/src/main/groovy/org/orbisgis/geoclimate/bdtopo/BDTopo.groovy
@@ -18,9 +18,6 @@
  *
  */
 package org.orbisgis.geoclimate.bdtopo
-
-import org.orbisgis.geoclimate.utils.LoggerUtils
-
 /**
  * Main module for BDTopo processing tasks
  */
@@ -143,7 +140,7 @@ abstract class BDTopo extends BDTopoUtils {
      * Meteorological Society 93, no. 12 (2012): 1879-1900.
      *
      */
-    static Map workflow(def input, int version) throws Exception{
+    static Map workflow(def input, int version) throws Exception {
         if (version == 2) {
             BDTopoV2Workflow bdtopo_v2_workflow = new BDTopoV2Workflow()
             return bdtopo_v2_workflow.execute(input)

diff --git a/bdtopo/src/main/groovy/org/orbisgis/geoclimate/bdtopo/BDTopoV2Workflow.groovy b/bdtopo/src/main/groovy/org/orbisgis/geoclimate/bdtopo/BDTopoV2Workflow.groovy
@@ -39,7 +39,7 @@ import java.sql.SQLException
 
 @Override
 Integer loadDataFromPostGIS(Object input_database_properties, Object code, Object distance, Object inputTables,
-                            Object inputSRID, H2GIS h2gis_datasource) throws Exception{
+                            Object inputSRID, H2GIS h2gis_datasource) throws Exception {
     def commune_location = inputTables.commune
     if (!commune_location) {
         throw new Exception("The commune table must be specified to run Geoclimate")
@@ -71,7 +71,7 @@ Integer loadDataFromPostGIS(Object input_database_properties, Object code, Objec
     if (commune_srid == 0 && inputSRID) {
         commune_srid = inputSRID
     } else if (commune_srid <= 0) {
-        throw new Exception( "Cannot find a SRID value for the layer commune.\n" +
+        throw new Exception("Cannot find a SRID value for the layer commune.\n" +
                 "Please set a valid OGC prj or use the parameter srid to force it.")
     }
 
@@ -80,11 +80,11 @@ Integer loadDataFromPostGIS(Object input_database_properties, Object code, Objec
     //Check if code is a string or a bbox
     //The zone is a osm bounding box represented by ymin,xmin , ymax,xmax,
     if (code in Collection) {
-        if(code.size()==3){
-            if(code[2]<100){
+        if (code.size() == 3) {
+            if (code[2] < 100) {
                 throw new Exception("The distance to create a bbox from a point must be greater than 100 meters")
             }
-            code = BDTopoUtils.bbox(code[0], code[1],code[2])
+            code = BDTopoUtils.bbox(code[0], code[1], code[2])
         }
         String inputTableName = """(SELECT
                     ST_INTERSECTION(st_setsrid(the_geom, $commune_srid), ST_MakeEnvelope(${code[1]},${code[0]},${code[3]},${code[2]}, $commune_srid)) as the_geom, CODE_INSEE  from $commune_location where 
@@ -139,7 +139,7 @@ Integer loadDataFromPostGIS(Object input_database_properties, Object code, Objec
 
         if (inputTables.troncon_voie_ferree) {
             //Extract troncon_voie_ferree
-                def inputTableName = "(SELECT ID, st_setsrid(the_geom, $commune_srid) as the_geom, NATURE, LARGEUR, NB_VOIES, POS_SOL, FRANCHISST FROM ${inputTables.troncon_voie_ferree}  WHERE st_setsrid(the_geom, $commune_srid) && 'SRID=$commune_srid;$geomToExtract'::GEOMETRY AND ST_INTERSECTS(st_setsrid(the_geom, $commune_srid), 'SRID=$commune_srid;$geomToExtract'::GEOMETRY))"
+            def inputTableName = "(SELECT ID, st_setsrid(the_geom, $commune_srid) as the_geom, NATURE, LARGEUR, NB_VOIES, POS_SOL, FRANCHISST FROM ${inputTables.troncon_voie_ferree}  WHERE st_setsrid(the_geom, $commune_srid) && 'SRID=$commune_srid;$geomToExtract'::GEOMETRY AND ST_INTERSECTS(st_setsrid(the_geom, $commune_srid), 'SRID=$commune_srid;$geomToExtract'::GEOMETRY))"
             outputTableName = "TRONCON_VOIE_FERREE"
             debug "Loading in the H2GIS database $outputTableName"
             IOMethods.exportToDataBase(sourceConnection, inputTableName, h2gis_datasource.getConnection(), outputTableName, -1, 1000)
@@ -233,7 +233,7 @@ int getVersion() {
 }
 
 @Override
-Map formatLayers(JdbcDataSource datasource, Map layers, float distance, float hLevMin) throws Exception{
+Map formatLayers(JdbcDataSource datasource, Map layers, float distance, float hLevMin) throws Exception {
     if (!hLevMin) {
         hLevMin = 3
     }
@@ -290,7 +290,7 @@ Map formatLayers(JdbcDataSource datasource, Map layers, float distance, float hL
 
 @Override
 def filterLinkedShapeFiles(def location, float distance, LinkedHashMap inputTables,
-                           int sourceSRID, int inputSRID, H2GIS h2gis_datasource) throws Exception{
+                           int sourceSRID, int inputSRID, H2GIS h2gis_datasource) throws Exception {
     def formatting_geom = "the_geom"
     if (sourceSRID == 0 && sourceSRID != inputSRID) {
         formatting_geom = "st_setsrid(the_geom, $inputSRID) as the_geom"
@@ -303,11 +303,11 @@ def filterLinkedShapeFiles(def location, float distance, LinkedHashMap inputTabl
     //Check if code is a string or a bbox
     //The zone is a osm bounding box represented by ymin,xmin , ymax,xmax,
     if (location in Collection) {
-        if(location.size()==3){
-            if(location[2]<100){
+        if (location.size() == 3) {
+            if (location[2] < 100) {
                 throw new Exception("The distance to create a bbox from a point must be greater than 100 meters")
             }
-            location = BDTopoUtils.bbox(location[0], location[1],location[2])
+            location = BDTopoUtils.bbox(location[0], location[1], location[2])
         }
         debug "Loading in the H2GIS database $outputTableName"
         h2gis_datasource.execute("""DROP TABLE IF EXISTS $outputTableName ; CREATE TABLE $outputTableName as  SELECT
@@ -317,7 +317,7 @@ def filterLinkedShapeFiles(def location, float distance, LinkedHashMap inputTabl
         debug "Loading in the H2GIS database $outputTableName"
         h2gis_datasource.execute("""DROP TABLE IF EXISTS $outputTableName ; CREATE TABLE $outputTableName as SELECT $formatting_geom, 
             CODE_INSEE FROM ${inputTables.commune} WHERE CODE_INSEE='$location' or lower(nom)='${location.toLowerCase()}'""".toString())
-    }else{
+    } else {
         throw new Exception("Invalid location data type. Please set a text value or a collection of coordinates to specify a bbox")
     }
     def count = h2gis_datasource.getRowCount(outputTableName)

diff --git a/bdtopo/src/main/groovy/org/orbisgis/geoclimate/bdtopo/BDTopoV3Workflow.groovy b/bdtopo/src/main/groovy/org/orbisgis/geoclimate/bdtopo/BDTopoV3Workflow.groovy
@@ -39,7 +39,7 @@ import java.sql.SQLException
 
 @Override
 def filterLinkedShapeFiles(def location, float distance, LinkedHashMap inputTables,
-                           int sourceSRID, int inputSRID, H2GIS h2gis_datasource) throws Exception{
+                           int sourceSRID, int inputSRID, H2GIS h2gis_datasource) throws Exception {
     def formatting_geom = "the_geom"
     if (sourceSRID == 0 && sourceSRID != inputSRID) {
         formatting_geom = "st_setsrid(the_geom, $inputSRID) as the_geom"
@@ -54,34 +54,31 @@ def filterLinkedShapeFiles(def location, float distance, LinkedHashMap inputTabl
     if (location in Collection) {
         debug "Loading in the H2GIS database $outputTableName"
         def communeColumns = h2gis_datasource.getColumnNames(inputTables.commune)
-        if(communeColumns.contains("INSEE_COM")) {
-            if(location.size()==3){
-                if(location[2]<100){
+        if (communeColumns.contains("INSEE_COM")) {
+            if (location.size() == 3) {
+                if (location[2] < 100) {
                     throw new IllegalArgumentException("The distance to create a bbox from a point must be greater than 100 meters")
                 }
-                location = BDTopoUtils.bbox(location[0], location[1],location[2])
+                location = BDTopoUtils.bbox(location[0], location[1], location[2])
             }
 
             h2gis_datasource.execute("""DROP TABLE IF EXISTS $outputTableName ; CREATE TABLE $outputTableName as  SELECT
                     ST_INTERSECTION(the_geom, ST_MakeEnvelope(${location[1]},${location[0]},${location[3]},${location[2]}, $sourceSRID)) as the_geom, INSEE_COM AS CODE_INSEE  from ${inputTables.commune} where the_geom 
                     && ST_MakeEnvelope(${location[1]},${location[0]},${location[3]},${location[2]}, $sourceSRID) """.toString())
-        }
-        else {
+        } else {
             throw new Exception("Cannot find a column insee_com or code_insee to filter the commune")
         }
     } else if (location instanceof String) {
         debug "Loading in the H2GIS database $outputTableName"
         def communeColumns = h2gis_datasource.getColumnNames(inputTables.commune)
-        if(communeColumns.contains("INSEE_COM")){
-        h2gis_datasource.execute("""DROP TABLE IF EXISTS $outputTableName ; CREATE TABLE $outputTableName as SELECT $formatting_geom, 
+        if (communeColumns.contains("INSEE_COM")) {
+            h2gis_datasource.execute("""DROP TABLE IF EXISTS $outputTableName ; CREATE TABLE $outputTableName as SELECT $formatting_geom, 
             INSEE_COM AS CODE_INSEE FROM ${inputTables.commune} WHERE INSEE_COM='$location' or lower(nom)='${location.toLowerCase()}'""".toString())
-        }
-        else {
+        } else {
             throw new Exception("Cannot find a column insee_com or code_insee to filter the commune")
         }
 
-    }
-    else {
+    } else {
         throw new IllegalArgumentException("Invalid location")
     }
     def count = h2gis_datasource.getRowCount(outputTableName)
@@ -215,7 +212,7 @@ def filterLinkedShapeFiles(def location, float distance, LinkedHashMap inputTabl
 
 @Override
 Integer loadDataFromPostGIS(Object input_database_properties, Object code, Object distance, Object inputTables, Object inputSRID,
-                            H2GIS h2gis_datasource) throws Exception{
+                            H2GIS h2gis_datasource) throws Exception {
     def commune_location = inputTables.commune
     if (!commune_location) {
         throw new Exception("The commune table must be specified to run Geoclimate")
@@ -259,12 +256,12 @@ Integer loadDataFromPostGIS(Object input_database_properties, Object code, Objec
     //The zone is a osm bounding box represented by ymin,xmin , ymax,xmax,
     if (code in Collection) {
         def communeColumns = h2gis_datasource.getColumnNames(commune_location)
-        if(communeColumns.contains("INSEE_COM")) {
-            if(code.size()==3){
-                if(code[2]<100){
+        if (communeColumns.contains("INSEE_COM")) {
+            if (code.size() == 3) {
+                if (code[2] < 100) {
                     throw new Exception("The distance to create a bbox from a point must be greater than 100 meters")
                 }
-                code = BDTopoUtils.bbox(code[0], code[1],code[2])
+                code = BDTopoUtils.bbox(code[0], code[1], code[2])
             }
             String inputTableName = """(SELECT
                     ST_INTERSECTION(st_setsrid(the_geom, $commune_srid), ST_MakeEnvelope(${code[1]},${code[0]},${code[3]},${code[2]}, $commune_srid)) as the_geom, INSEE_COM as CODE_INSEE   from $commune_location where 
@@ -273,20 +270,19 @@ Integer loadDataFromPostGIS(Object input_database_properties, Object code, Objec
                     st_intersects(st_setsrid(the_geom, $commune_srid), ST_MakeEnvelope(${code[1]},${code[0]},${code[3]},${code[2]}, $commune_srid)))""".toString()
             debug "Loading in the H2GIS database $outputTableName"
             IOMethods.exportToDataBase(sourceConnection, inputTableName, h2gis_datasource.getConnection(), outputTableName, -1, 100)
-        }else {
+        } else {
             throw new Exception("Cannot find a column insee_com or code_insee to filter the commune")
         }
     } else if (code instanceof String) {
         def communeColumns = h2gis_datasource.getColumnNames(commune_location)
-        if(communeColumns.contains("insee_com")){
+        if (communeColumns.contains("insee_com")) {
             String inputTableName = "(SELECT st_setsrid(the_geom, $commune_srid) as the_geom, INSEE_COM as CODE_INSEE FROM $commune_location WHERE INSEE_COM='$code' or lower(nom)='${code.toLowerCase()}')"
             debug "Loading in the H2GIS database $outputTableName"
             IOMethods.exportToDataBase(sourceConnection, inputTableName, h2gis_datasource.getConnection(), outputTableName, -1, 1000)
-        }
-        else {
+        } else {
             throw new Exception("Cannot find a column insee_com to filter the commune")
         }
-    }else{
+    } else {
         throw new Exception("Invalid location data type. Please set a text value or a collection of coordinates to specify a bbox")
     }
     def count = h2gis_datasource.getRowCount(outputTableName)
@@ -427,7 +423,7 @@ int getVersion() {
 }
 
 @Override
-Map formatLayers(JdbcDataSource datasource, Map layers, float distance, float hLevMin) throws Exception{
+Map formatLayers(JdbcDataSource datasource, Map layers, float distance, float hLevMin) throws Exception {
     if (!hLevMin) {
         hLevMin = 3
     }