diff --git a/bdtopo/src/main/groovy/org/orbisgis/geoclimate/bdtopo/AbstractBDTopoWorkflow.groovy b/bdtopo/src/main/groovy/org/orbisgis/geoclimate/bdtopo/AbstractBDTopoWorkflow.groovy index 5fbac7e8f1..2393453556 100644 --- a/bdtopo/src/main/groovy/org/orbisgis/geoclimate/bdtopo/AbstractBDTopoWorkflow.groovy +++ b/bdtopo/src/main/groovy/org/orbisgis/geoclimate/bdtopo/AbstractBDTopoWorkflow.groovy @@ -542,7 +542,9 @@ abstract class AbstractBDTopoWorkflow extends BDTopoUtils { } def allowed_grid_indicators = ["BUILDING_FRACTION", "BUILDING_HEIGHT", "BUILDING_POP", "BUILDING_TYPE_FRACTION", "WATER_FRACTION", "VEGETATION_FRACTION", "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"] + "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"] def allowedOutputIndicators = allowed_grid_indicators.intersect(list_indicators*.toUpperCase()) if (allowedOutputIndicators) { //Update the RSU indicators list according the grid indicators diff --git a/geoindicators/src/main/groovy/org/orbisgis/geoclimate/geoindicators/DataUtils.groovy b/geoindicators/src/main/groovy/org/orbisgis/geoclimate/geoindicators/DataUtils.groovy index 9e046f840f..66b73dcbf8 100644 --- a/geoindicators/src/main/groovy/org/orbisgis/geoclimate/geoindicators/DataUtils.groovy +++ b/geoindicators/src/main/groovy/org/orbisgis/geoclimate/geoindicators/DataUtils.groovy @@ -73,7 +73,7 @@ String joinTables(JdbcDataSource datasource, Map inputTableNamesWithId, String o } leftQuery += " LEFT JOIN $key as $alias ON $alias.$value = $columnKey " } - indexes += "CREATE INDEX IF NOT EXISTS ${key}_ids ON $key USING BTREE($value);" + indexes += "CREATE INDEX IF NOT EXISTS ${key}_ids ON $key ($value);" alias++ } diff --git a/geoindicators/src/main/groovy/org/orbisgis/geoclimate/geoindicators/GenericIndicators.groovy b/geoindicators/src/main/groovy/org/orbisgis/geoclimate/geoindicators/GenericIndicators.groovy index 41c28a5d5e..21b9025751 100644 --- a/geoindicators/src/main/groovy/org/orbisgis/geoclimate/geoindicators/GenericIndicators.groovy +++ b/geoindicators/src/main/groovy/org/orbisgis/geoclimate/geoindicators/GenericIndicators.groovy @@ -447,7 +447,7 @@ String distributionCharacterization(JdbcDataSource datasource, String distribTab def outputTableName = prefix prefixName, BASENAME // Get all columns from the distribution table and remove the geometry column if exists - def allColumns = datasource."$distribTableName".columns + def allColumns = datasource.getColumnNames(distribTableName) if (allColumns.contains(GEOMETRY_FIELD)) { allColumns -= GEOMETRY_FIELD } @@ -455,6 +455,11 @@ String distributionCharacterization(JdbcDataSource datasource, String distribTab def distribColumns = allColumns.minus(inputId.toUpperCase()) def nbDistCol = distribColumns.size() + if(distribColumns.size()==0){ + error("Any columns to compute the distribution from the table $distribTableName".toString()) + return + } + def idxExtrem = nbDistCol - 1 def idxExtrem_1 = nbDistCol - 2 if (extremum.toUpperCase() == "LEAST") { diff --git a/geoindicators/src/main/groovy/org/orbisgis/geoclimate/geoindicators/RsuIndicators.groovy b/geoindicators/src/main/groovy/org/orbisgis/geoclimate/geoindicators/RsuIndicators.groovy index da1c592f7a..e6d8a1497d 100644 --- a/geoindicators/src/main/groovy/org/orbisgis/geoclimate/geoindicators/RsuIndicators.groovy +++ b/geoindicators/src/main/groovy/org/orbisgis/geoclimate/geoindicators/RsuIndicators.groovy @@ -174,9 +174,16 @@ String freeExternalFacadeDensityExact(JdbcDataSource datasource, String building DROP TABLE IF EXISTS $onlyBuildRsu; CREATE TABLE $onlyBuildRsu AS SELECT a.$idRsu, - (a.$FACADE_AREA-b.$FACADE_AREA)/a.$RSU_AREA AS FREE_EXTERNAL_FACADE_DENSITY + a.$FACADE_AREA/a.$RSU_AREA AS FREE_EXTERNAL_FACADE_DENSITY FROM $buildLineRsu AS a LEFT JOIN $sharedLineRsu AS b - ON a.$idRsu = b.$idRsu""".toString() + ON a.$idRsu = b.$idRsu WHERE b.$idRsu IS NULL + union all + SELECT a.$idRsu, + (a.$FACADE_AREA-b.$FACADE_AREA)/a.$RSU_AREA AS FREE_EXTERNAL_FACADE_DENSITY + FROM $buildLineRsu AS a right JOIN $sharedLineRsu AS b + ON a.$idRsu = b.$idRsu + """.toString() + // 5. Join RSU having no buildings and set their value to 0 datasource.createIndex(onlyBuildRsu,idRsu) @@ -207,6 +214,7 @@ String freeExternalFacadeDensityExact(JdbcDataSource datasource, String building * @param datasource A connexion to a database (H2GIS, PostGIS, ...) where are stored the input Table and in which * the resulting database will be stored * @param rsu The name of the input ITable where are stored the RSU + * @param id_rsu Unique identifier column name * @param correlationBuildingTable The name of the input ITable where are stored the buildings and the relationships * between buildings and RSU * @param pointDensity The density of points (nb / free m²) used to calculate the spatial average SVF. Use 0.008f @@ -225,11 +233,11 @@ String freeExternalFacadeDensityExact(JdbcDataSource datasource, String building * @author Jérémy Bernard * @author Erwan Bocher */ -String groundSkyViewFactor(JdbcDataSource datasource, String rsu, String correlationBuildingTable, float pointDensity, +String groundSkyViewFactor(JdbcDataSource datasource, String rsu, String id_rsu, String correlationBuildingTable, float pointDensity, float rayLength, int numberOfDirection, String prefixName) { def GEOMETRIC_COLUMN_RSU = "the_geom" def GEOMETRIC_COLUMN_BU = "the_geom" - def ID_COLUMN_RSU = "id_rsu" + def ID_COLUMN_RSU = id_rsu def HEIGHT_WALL = "height_wall" def BASE_NAME = "ground_sky_view_factor" @@ -318,7 +326,7 @@ String groundSkyViewFactor(JdbcDataSource datasource, String rsu, String correla // The result of the SVF calculation is averaged at RSU scale datasource """ - CREATE TABLE $outputTableName(id_rsu integer, $BASE_NAME double) + CREATE TABLE $outputTableName($ID_COLUMN_RSU integer, $BASE_NAME double) AS (SELECT a.$ID_COLUMN_RSU, CASE WHEN AVG(b.SVF) is not null THEN AVG(b.SVF) ELSE 1 END FROM $rsu a LEFT JOIN $svfPts b @@ -392,13 +400,13 @@ String aspectRatio(JdbcDataSource datasource, String rsuTable, String rsuFreeExt * @return A database table name. * @author Jérémy Bernard */ -String projectedFacadeAreaDistribution(JdbcDataSource datasource, String building, String rsu, +String projectedFacadeAreaDistribution(JdbcDataSource datasource, String building, String rsu,String id_rsu, List listLayersBottom = [0, 10, 20, 30, 40, 50], int numberOfDirection = 12, String prefixName) { def BASE_NAME = "projected_facade_area_distribution" def GEOMETRIC_COLUMN_RSU = "the_geom" def GEOMETRIC_COLUMN_BU = "the_geom" - def ID_COLUMN_RSU = "id_rsu" + def ID_COLUMN_RSU = id_rsu def ID_COLUMN_BU = "id_build" def HEIGHT_WALL = "height_wall" @@ -410,7 +418,7 @@ String projectedFacadeAreaDistribution(JdbcDataSource datasource, String buildin datasource.createSpatialIndex(building,"the_geom") datasource.createSpatialIndex(rsu,"the_geom") datasource.createIndex(building,"id_build") - datasource.createIndex(rsu,"id_rsu") + datasource.createIndex(rsu,ID_COLUMN_RSU) if (360 % numberOfDirection == 0 && numberOfDirection % 2 == 0) { @@ -504,7 +512,7 @@ String projectedFacadeAreaDistribution(JdbcDataSource datasource, String buildin // Intersections between free facades and rsu geometries are calculated datasource """ DROP TABLE IF EXISTS $buildingFreeExpl; - CREATE TABLE $buildingFreeExpl(id_rsu INTEGER, the_geom GEOMETRY, $namesAndType) AS + CREATE TABLE $buildingFreeExpl($ID_COLUMN_RSU INTEGER, the_geom GEOMETRY, $namesAndType) AS (SELECT a.$ID_COLUMN_RSU, ST_INTERSECTION(a.$GEOMETRIC_COLUMN_RSU, ST_TOMULTILINE(b.the_geom)), ${onlyNamesB} FROM $rsu a, $buildingLayer b WHERE a.$GEOMETRIC_COLUMN_RSU && b.the_geom @@ -513,8 +521,8 @@ String projectedFacadeAreaDistribution(JdbcDataSource datasource, String buildin // Intersections facades are exploded to multisegments datasource """DROP TABLE IF EXISTS $rsuInter; - CREATE TABLE $rsuInter(id_rsu INTEGER, the_geom GEOMETRY, $namesAndType) - AS (SELECT id_rsu, the_geom, ${onlyNames} FROM ST_EXPLODE('$buildingFreeExpl'))""".toString() + CREATE TABLE $rsuInter($ID_COLUMN_RSU INTEGER, the_geom GEOMETRY, $namesAndType) + AS (SELECT $ID_COLUMN_RSU, the_geom, ${onlyNames} FROM ST_EXPLODE('$buildingFreeExpl'))""".toString() // The analysis is then performed for each direction ('numberOfDirection' / 2 because calculation @@ -550,8 +558,8 @@ String projectedFacadeAreaDistribution(JdbcDataSource datasource, String buildin sumNamesDir = sumNamesDir.join(",") def query = "DROP TABLE IF EXISTS $finalIndicator; " + - "CREATE TABLE $finalIndicator AS SELECT a.id_rsu," + queryColumns + - " FROM (SELECT id_rsu, CASE WHEN ST_AZIMUTH(ST_STARTPOINT(the_geom), ST_ENDPOINT(the_geom)) >= PI()" + + "CREATE TABLE $finalIndicator AS SELECT a.$ID_COLUMN_RSU," + queryColumns + + " FROM (SELECT $ID_COLUMN_RSU, CASE WHEN ST_AZIMUTH(ST_STARTPOINT(the_geom), ST_ENDPOINT(the_geom)) >= PI()" + "THEN ST_AZIMUTH(ST_STARTPOINT(the_geom), ST_ENDPOINT(the_geom)) - PI() " + "ELSE ST_AZIMUTH(ST_STARTPOINT(the_geom), ST_ENDPOINT(the_geom)) END AS azimuth," + " ST_LENGTH(the_geom) AS length, ${onlyNames} FROM $rsuInter) a" @@ -559,15 +567,15 @@ String projectedFacadeAreaDistribution(JdbcDataSource datasource, String buildin datasource query.toString() - datasource.createIndex(finalIndicator,"id_rsu") + datasource.createIndex(finalIndicator,ID_COLUMN_RSU) // Sum area at RSU scale and fill null values with 0 datasource """ DROP TABLE IF EXISTS $outputTableName; CREATE TABLE ${outputTableName} - AS SELECT a.id_rsu, ${sumNamesDir} + AS SELECT a.$ID_COLUMN_RSU, ${sumNamesDir} FROM $rsu a LEFT JOIN $finalIndicator b - ON a.id_rsu = b.id_rsu - GROUP BY a.id_rsu""".toString() + ON a.$ID_COLUMN_RSU = b.$ID_COLUMN_RSU + GROUP BY a.$ID_COLUMN_RSU""".toString() // Remove all temporary tables created datasource """DROP TABLE IF EXISTS $buildingIntersection, $buildingIntersectionExpl, @@ -820,6 +828,7 @@ String roofAreaDistribution(JdbcDataSource datasource, String rsu, String buildi * the resulting database will be stored * @param rsuTable the name of the input ITable where are stored the RSU geometries, the rsu_geometric_mean_height and * the rsu_projected_facade_area_distribution fields + * @param id_rsu Unique identifier column name * @param projectedFacadeAreaName the name base used for naming the projected facade area field within the * inputA rsuTable (default rsu_projected_facade_area_distribution - note that the field is also constructed * with the direction and vertical height informations) @@ -835,12 +844,13 @@ String roofAreaDistribution(JdbcDataSource datasource, String rsu, String buildi * * @author Jérémy Bernard */ -String effectiveTerrainRoughnessLength(JdbcDataSource datasource, String rsuTable, +String effectiveTerrainRoughnessLength(JdbcDataSource datasource, String rsuTable, String id_rsu, String projectedFacadeAreaName, String geometricMeanBuildingHeightName, - List listLayersBottom = [0, 10, 20, 30, 40, 50], int numberOfDirection = 12, String prefixName) { + List listLayersBottom = [0, 10, 20, 30, 40, 50], + int numberOfDirection = 12, String prefixName) { def GEOMETRIC_COLUMN = "the_geom" - def ID_COLUMN_RSU = "id_rsu" + def ID_COLUMN_RSU = id_rsu def BASE_NAME = "effective_terrain_roughness_length" debug "Executing RSU effective terrain roughness length" @@ -858,14 +868,16 @@ String effectiveTerrainRoughnessLength(JdbcDataSource datasource, String rsuTabl // The name of the outputTableName is constructed def outputTableName = prefix prefixName, "rsu_" + BASE_NAME + def layerSize = listLayersBottom.size() // The lambda_f indicator is first calculated def lambdaQuery = "DROP TABLE IF EXISTS $lambdaTable;" + "CREATE TABLE $lambdaTable AS SELECT $ID_COLUMN_RSU, $geometricMeanBuildingHeightName, (" - for (int i in 1..listLayersBottom.size()) { + for (int i in 1..layerSize) { + //TODO : why an array here and not a variable names[i - 1] = "${projectedFacadeAreaName}_H${listLayersBottom[i - 1]}_${listLayersBottom[i]}" - if (i == listLayersBottom.size()) { - names[listLayersBottom.size() - 1] = - "${projectedFacadeAreaName}_H${listLayersBottom[listLayersBottom.size() - 1]}" + if (i == layerSize) { + names[layerSize - 1] = + "${projectedFacadeAreaName}_H${listLayersBottom[layerSize - 1]}" } for (int d = 0; d < numberOfDirection / 2; d++) { def dirDeg = d * 360 / numberOfDirection @@ -1089,6 +1101,7 @@ String linearRoadOperations(JdbcDataSource datasource, String rsuTable, String r * @param datasource A connexion to a database (H2GIS, PostGIS, ...) where are stored the input Table and in which * the resulting database will be stored * @param rsuTable the name of the input ITable where are stored the effectiveTerrainRoughnessHeight values + * @param id_rsu Unique identifier column name * @param effectiveTerrainRoughnessLength the field name corresponding to the RSU effective terrain roughness class due * to roughness elements (buildings, trees, etc.) (in the rsuTable) * @param prefixName String use as prefix to name the output table @@ -1097,8 +1110,7 @@ String linearRoadOperations(JdbcDataSource datasource, String rsuTable, String r * * @author Jérémy Bernard */ -String effectiveTerrainRoughnessClass(JdbcDataSource datasource, String rsu, String effectiveTerrainRoughnessLength, String prefixName) { - def ID_COLUMN_RSU = "id_rsu" +String effectiveTerrainRoughnessClass(JdbcDataSource datasource, String rsu, String id_rsu, String effectiveTerrainRoughnessLength, String prefixName) { def BASE_NAME = "effective_terrain_roughness_class" debug "Executing RSU effective terrain roughness class" @@ -1108,7 +1120,7 @@ String effectiveTerrainRoughnessClass(JdbcDataSource datasource, String rsu, Str // Based on the lookup Table of Davenport datasource """DROP TABLE IF EXISTS $outputTableName; - CREATE TABLE $outputTableName AS SELECT $ID_COLUMN_RSU, + CREATE TABLE $outputTableName AS SELECT $id_rsu, CASEWHEN($effectiveTerrainRoughnessLength<0.0 OR $effectiveTerrainRoughnessLength IS NULL, null, CASEWHEN($effectiveTerrainRoughnessLength<0.00035, 1, CASEWHEN($effectiveTerrainRoughnessLength<0.01525, 2, @@ -1681,7 +1693,7 @@ String buildingSurfaceDensity(JdbcDataSource datasource, String facadeDensityTab } /** - * Script to compute the distributio of (only) horizontal roof area fraction for each layer of the canopy. If the height + * Script to compute the distribution of (only) horizontal roof area fraction for each layer of the canopy. If the height * roof and height wall differ, take the average value at building height. Note that this process first cut the buildings * according to RSU in order to calculate the exact distribution. * @@ -1694,6 +1706,7 @@ String buildingSurfaceDensity(JdbcDataSource datasource, String facadeDensityTab * @param prefixName String use as prefix to name the output table * @param listLayersBottom the list of height corresponding to the bottom of each vertical layers (default * [0, 5, 10, 15, 20, 25, 30, 35, 40, 45, 50]) + * @param cutBuilding set to true if the building layer must be cuted by the rsu layer * * @return Table name in which the rsu id and their corresponding indicator value are stored * @@ -1701,7 +1714,7 @@ String buildingSurfaceDensity(JdbcDataSource datasource, String facadeDensityTab */ String roofFractionDistributionExact(JdbcDataSource datasource, String rsu, String building, String idRsu, List listLayersBottom = [0, 5, 10, 15, 20, 25, 30, 35, 40, 45, 50], - boolean density = true, String prefixName) { + boolean cutBuilding=true, String prefixName) { def GEOMETRIC_COLUMN_RSU = "the_geom" def GEOMETRIC_COLUMN_BU = "the_geom" def ID_COLUMN_BU = "id_build" @@ -1718,15 +1731,17 @@ String roofFractionDistributionExact(JdbcDataSource datasource, String rsu, Stri // To avoid overwriting the output files of this step, a unique identifier is created // Temporary table names - def buildInter = postfix "build_inter" + def buildInter =building def rsuBuildingArea = postfix "rsu_building_area" def buildFracH = postfix "roof_frac_H" def bufferTable = postfix "buffer_table" - // 1. Create the intersection between buildings and RSU polygons - datasource."$building"."$ID_COLUMN_BU".createIndex() - datasource."$rsu"."$idRsu".createIndex() - datasource """ + if(cutBuilding) { + buildInter = postfix "build_inter" + // 1. Create the intersection between buildings and RSU polygons + datasource."$building"."$ID_COLUMN_BU".createIndex() + datasource."$rsu"."$idRsu".createIndex() + datasource """ DROP TABLE IF EXISTS $buildInter; CREATE TABLE $buildInter AS SELECT a.$ID_COLUMN_BU, a.$idRsu, @@ -1734,9 +1749,10 @@ String roofFractionDistributionExact(JdbcDataSource datasource, String rsu, Stri (a.$HEIGHT_WALL + a.$HEIGHT_ROOF) / 2 AS $BUILD_HEIGHT FROM $building AS a LEFT JOIN $rsu AS b ON a.$idRsu = b.$idRsu""".toString() + } // 2. Calculate the total building roof area within each RSU - datasource."$buildInter"."$idRsu".createIndex() + datasource.createIndex(buildInter, idRsu) datasource """ DROP TABLE IF EXISTS $rsuBuildingArea; CREATE TABLE $rsuBuildingArea @@ -1745,9 +1761,8 @@ String roofFractionDistributionExact(JdbcDataSource datasource, String rsu, Stri GROUP BY $idRsu""".toString() // 3. Calculate the fraction of roof for each level of the canopy (defined by 'listLayersBottom') except the last - datasource."$buildInter"."$BUILD_HEIGHT".createIndex() - datasource."$buildInter"."$idRsu".createIndex() - datasource."$rsuBuildingArea"."$idRsu".createIndex() + datasource.createIndex(buildInter,BUILD_HEIGHT) + datasource.createIndex(rsuBuildingArea,idRsu) def tab_H = [:] def indicToJoin = [:] for (i in 1..(listLayersBottom.size() - 1)) { @@ -1840,14 +1855,17 @@ String roofFractionDistributionExact(JdbcDataSource datasource, String rsu, Stri * [0, 5, 10, 15, 20, 25, 30, 35, 40, 45, 50]) * @param numberOfDirection the number of directions used for the calculation - according to the method used it should * be divisor of 360 AND a multiple of 2 (default 12) + * @param distributionAsIndex set the value to false to avoid normalizing the distribution by the area of the input rsu + * and the range of the bottom layers * @param prefixName String use as prefix to name the output table * * @return A database table name. * @author Jérémy Bernard + * TODO : https://github.com/orbisgis/geoclimate/issues/848 */ String frontalAreaIndexDistribution(JdbcDataSource datasource, String building, String rsu, String idRsu, List listLayersBottom = [0, 5, 10, 15, 20, 25, 30, 35, 40, 45, 50], - int numberOfDirection = 12, String prefixName) { + int numberOfDirection = 12, boolean distributionAsIndex = true, String prefixName) { def GEOMETRIC_FIELD_RSU = "the_geom" def GEOMETRIC_FIELD_BU = "the_geom" def ID_FIELD_BU = "id_build" @@ -1939,7 +1957,8 @@ String frontalAreaIndexDistribution(JdbcDataSource datasource, String building, // Indicator name def indicName = "FRONTAL_AREA_INDEX_H${layer_bottom}_${layer_top}_D${k * angleRangeDeg}_${(k + 1) * angleRangeDeg}" // Define query to sum the projected facade for buildings and shared facades - dirQueryVertFrac[k] = """ + if(distributionAsIndex) { + dirQueryVertFrac[k] = """ CASE WHEN $v > AZIMUTH AND $v-AZIMUTH < PI() THEN CASE WHEN $HEIGHT_WALL >= $layer_top THEN LENGTH*SIN($v-AZIMUTH) @@ -1953,7 +1972,24 @@ String frontalAreaIndexDistribution(JdbcDataSource datasource, String building, ELSE 0 END END AS $indicName""" - dirQueryDiv[k] = """COALESCE(SUM(b.$indicName)/ST_AREA(a.$GEOMETRIC_FIELD_RSU), 0) AS $indicName""" + dirQueryDiv[k] = """COALESCE(SUM(b.$indicName)/ST_AREA(a.$GEOMETRIC_FIELD_RSU), 0) AS $indicName""" + }else{ + dirQueryVertFrac[k] = """ + CASE WHEN $v > AZIMUTH AND $v-AZIMUTH < PI() + THEN CASE WHEN $HEIGHT_WALL >= $layer_top + THEN LENGTH*SIN($v-AZIMUTH) + ELSE LENGTH*SIN($v-AZIMUTH)*($HEIGHT_WALL-$layer_bottom) + END + ELSE CASE WHEN $v - AZIMUTH < -PI() + THEN CASE WHEN $HEIGHT_WALL >= $layer_top + THEN LENGTH*SIN($v+2*PI()-AZIMUTH) + ELSE LENGTH*SIN($v+2*PI()-AZIMUTH)*($HEIGHT_WALL-$layer_bottom) + END + ELSE 0 + END + END AS $indicName""" + dirQueryDiv[k] = """COALESCE(SUM(b.$indicName), 0) AS $indicName""" + } } // Calculates projected surfaces for buildings and shared facades datasource """ @@ -1990,6 +2026,9 @@ String frontalAreaIndexDistribution(JdbcDataSource datasource, String building, dirList.each { k, v -> // Indicator name def indicName = "FRONTAL_AREA_INDEX_H${layer_bottom}_${layer_top}_D${k * angleRangeDeg}_${(k + 1) * angleRangeDeg}" + if(!distributionAsIndex){ + indicName = "FRONTAL_AREA_INDEX_H${layer_bottom}_D${k * angleRangeDeg}_${(k + 1) * angleRangeDeg}" + } // Define query to calculate the vertical fraction of projected facade for buildings and shared facades dirQueryVertFrac[k] = """ CASE WHEN $v-AZIMUTH > 0 AND $v-AZIMUTH < PI() diff --git a/geoindicators/src/main/groovy/org/orbisgis/geoclimate/geoindicators/SpatialUnits.groovy b/geoindicators/src/main/groovy/org/orbisgis/geoclimate/geoindicators/SpatialUnits.groovy index 46107ee6da..48d4fa77fa 100644 --- a/geoindicators/src/main/groovy/org/orbisgis/geoclimate/geoindicators/SpatialUnits.groovy +++ b/geoindicators/src/main/groovy/org/orbisgis/geoclimate/geoindicators/SpatialUnits.groovy @@ -126,7 +126,7 @@ String createTSU(JdbcDataSource datasource, String inputTableName, String inputz CREATE TABLE $outputTableName AS SELECT EXPLOD_ID AS $COLUMN_ID_NAME, ST_SETSRID(a.the_geom, $epsg) AS the_geom FROM ST_EXPLODE('( - SELECT ST_POLYGONIZE(ST_UNION(ST_NODE(ST_ACCUM(the_geom)))) AS the_geom + SELECT ST_BUFFER(ST_POLYGONIZE(ST_UNION(ST_NODE(ST_ACCUM(the_geom)))), -0.01) AS the_geom FROM $inputTableName)') AS a, $inputzone AS b WHERE a.the_geom && b.the_geom @@ -136,9 +136,9 @@ String createTSU(JdbcDataSource datasource, String inputTableName, String inputz datasource """ DROP TABLE IF EXISTS $outputTableName; CREATE TABLE $outputTableName AS - SELECT EXPLOD_ID AS $COLUMN_ID_NAME, ST_SETSRID(ST_FORCE2D(the_geom), $epsg) AS the_geom + SELECT EXPLOD_ID AS $COLUMN_ID_NAME, ST_SETSRID(st_buffer(the_geom, -0.01), $epsg) AS the_geom FROM ST_EXPLODE('( - SELECT ST_POLYGONIZE(ST_UNION(ST_NODE(ST_ACCUM(the_geom)))) AS the_geom + SELECT ST_BUFFER(ST_POLYGONIZE(ST_UNION(ST_NODE(ST_ACCUM(the_geom)))),-0.01) AS the_geom FROM $inputTableName)') where st_area(the_geom) > $area""".toString() } diff --git a/geoindicators/src/main/groovy/org/orbisgis/geoclimate/geoindicators/WorkflowGeoIndicators.groovy b/geoindicators/src/main/groovy/org/orbisgis/geoclimate/geoindicators/WorkflowGeoIndicators.groovy index 53ddb2ddc4..db1146f88e 100644 --- a/geoindicators/src/main/groovy/org/orbisgis/geoclimate/geoindicators/WorkflowGeoIndicators.groovy +++ b/geoindicators/src/main/groovy/org/orbisgis/geoclimate/geoindicators/WorkflowGeoIndicators.groovy @@ -26,6 +26,8 @@ import org.locationtech.jts.geom.Geometry import org.orbisgis.data.jdbc.JdbcDataSource import org.orbisgis.geoclimate.Geoindicators +import java.sql.SQLException + @BaseScript Geoindicators geoindicators /** @@ -357,14 +359,15 @@ String computeRSUIndicators(JdbcDataSource datasource, String buildingTable, [facadeDensListLayersBottom : [0, 10, 20, 30, 40, 50], facadeDensNumberOfDirection : 12, svfPointDensity : 0.008, - svfRayLength : 100, svfNumberOfDirection: 60, + svfRayLength : 100, + svfNumberOfDirection : 60, heightColumnName : "height_roof", inputFields : ["id_build", "the_geom"], levelForRoads : [0], angleRangeSizeBuDirection : 30, svfSimplified : true, indicatorUse : ["LCZ", "UTRF", "TEB"], - surfSuperpositions : ["high_vegetation": ["water", "building", "low_vegetation", "rail","road", "impervious"]], + surfSuperpositions : ["high_vegetation": ["water", "building", "low_vegetation", "rail", "road", "impervious"]], surfPriorities : ["water", "building", "high_vegetation", "low_vegetation", "rail", "road", "impervious"], buildingAreaTypeAndComposition: ["light_industry": ["light_industry"], "heavy_industry": ["heavy_industry"], @@ -419,7 +422,7 @@ String computeRSUIndicators(JdbcDataSource datasource, String buildingTable, def columnIdRsu = "id_rsu" def columnIdBuild = "id_build" def BASE_NAME = "rsu_indicators" - def tablesToDrop=[] + def tablesToDrop = [] // Maps for intermediate or final joins def finalTablesToJoin = [:] def intermediateJoin = [:] @@ -447,9 +450,9 @@ String computeRSUIndicators(JdbcDataSource datasource, String buildingTable, if (!surfaceFractions) { // Calculate all surface fractions indicators // Need to create the smallest geometries used as input of the surface fraction process - def superpositionsTable = Geoindicators.RsuIndicators.smallestCommunGeometry(datasource, + def superpositionsTable = Geoindicators.RsuIndicators.smallestCommunGeometry(datasource, rsu, "id_rsu", buildingTable, road, water, vegetation, - impervious,rail, temporaryPrefName) + impervious, rail, temporaryPrefName) if (!superpositionsTable) { info "Cannot compute the smallest commun geometries" return @@ -464,7 +467,7 @@ String computeRSUIndicators(JdbcDataSource datasource, String buildingTable, return } - tablesToDrop< urbTypoCorrespondenceTabInverted[ini] = fin } - datasource.createIndex(utrfBuild,"I_TYPO") + datasource.createIndex(utrfBuild, "I_TYPO") def queryDistinct = """SELECT DISTINCT I_TYPO AS I_TYPO FROM $utrfBuild""" def mapTypos = datasource.rows(queryDistinct.toString()) def listTypos = [] @@ -909,8 +912,8 @@ Map computeTypologyIndicators(JdbcDataSource datasource, String building_indicat } queryCaseWhenReplace = queryCaseWhenReplace + " 'unknown' " + endCaseWhen utrfBuilding = postfix "UTRF_BUILDING" - datasource.createIndex(utrfBuild,COLUMN_ID_BUILD) - datasource.createIndex(building_indicators,COLUMN_ID_BUILD) + datasource.createIndex(utrfBuild, COLUMN_ID_BUILD) + datasource.createIndex(building_indicators, COLUMN_ID_BUILD) datasource """ DROP TABLE IF EXISTS $utrfBuilding; CREATE TABLE $utrfBuilding AS SELECT a.$COLUMN_ID_BUILD, a.$COLUMN_ID_RSU, a.THE_GEOM, @@ -923,9 +926,9 @@ Map computeTypologyIndicators(JdbcDataSource datasource, String building_indicat def queryCasewhen = [:] queryCasewhen["AREA"] = "" queryCasewhen["FLOOR_AREA"] = "" - datasource.createIndex(rsu_indicators,COLUMN_ID_RSU) - datasource.createIndex(building_indicators,COLUMN_ID_RSU) - datasource.createIndex(utrfBuilding,COLUMN_ID_BUILD) + datasource.createIndex(rsu_indicators, COLUMN_ID_RSU) + datasource.createIndex(building_indicators, COLUMN_ID_RSU) + datasource.createIndex(utrfBuilding, COLUMN_ID_BUILD) queryCasewhen.keySet().each { ind -> def querySum = "" @@ -943,26 +946,26 @@ Map computeTypologyIndicators(JdbcDataSource datasource, String building_indicat WHERE b.$COLUMN_ID_RSU IS NOT NULL GROUP BY b.$COLUMN_ID_RSU """.toString() - tablesToDrop<2.1.0-SNAPSHOT 2.0.9 1.4.11 - 3.0.11 + 3.0.19 4.6.3 2.11.0 1.21