Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Improve OSM layers plus some fixes on spatial units creation #958

Merged
merged 4 commits into from
May 4, 2024
Merged

Improve OSM layers plus some fixes on spatial units creation #958

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
4 commits
Select commit Hold shift + click to select a range
ebde173
Some improvements
ebocher Apr 25, 2024
f0b7a38
Fix block construction with invalid geometries
ebocher Apr 26, 2024
8aef2b4
Fix test due to some changes on blocks
ebocher Apr 26, 2024
d633ae8
Update BUILDING.geojson
ebocher Apr 26, 2024
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
4 changes: 2 additions & 2 deletions ...ndicators/src/main/groovy/org/orbisgis/geoclimate/geoindicators/BuildingIndicators.groovy
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -361,14 +361,14 @@ String roadDistance(JdbcDataSource datasource, String building, String inputRoad
// The buffer is created
datasource """DROP TABLE IF EXISTS $build_buffer;
CREATE TABLE $build_buffer AS
SELECT $ID_FIELD_BU, ST_BUFFER($GEOMETRIC_FIELD, $bufferDist) AS $GEOMETRIC_FIELD
SELECT $ID_FIELD_BU, ST_BUFFER($GEOMETRIC_FIELD, $bufferDist, 2) AS $GEOMETRIC_FIELD
FROM $building;
CREATE SPATIAL INDEX IF NOT EXISTS buff_ids ON $build_buffer ($GEOMETRIC_FIELD)""".toString()
// The road surfaces are created
datasource """
DROP TABLE IF EXISTS $road_surf;
CREATE TABLE $road_surf AS
SELECT ST_BUFFER($GEOMETRIC_FIELD, $ROAD_WIDTH::DOUBLE PRECISION/2,'endcap=flat') AS $GEOMETRIC_FIELD
SELECT ST_BUFFER($GEOMETRIC_FIELD, $ROAD_WIDTH::DOUBLE PRECISION/2,'quad_segs=2 endcap=flat') AS $GEOMETRIC_FIELD
FROM $inputRoadTableName;
CREATE SPATIAL INDEX IF NOT EXISTS buff_ids ON $road_surf ($GEOMETRIC_FIELD)""".toString()
// The roads located within the buffer are identified
Expand Down
29 changes: 16 additions & 13 deletions geoindicators/src/main/groovy/org/orbisgis/geoclimate/geoindicators/RsuIndicators.groovy
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -1233,7 +1233,7 @@ String extendedFreeFacadeFraction(JdbcDataSource datasource, String building, St
def notIncBu = postfix "notIncBu"

datasource """DROP TABLE IF EXISTS $extRsuTable; CREATE TABLE $extRsuTable AS SELECT
ST_BUFFER($GEOMETRIC_FIELD, $buffDist) AS $GEOMETRIC_FIELD,
ST_BUFFER($GEOMETRIC_FIELD, $buffDist, 2) AS $GEOMETRIC_FIELD,
$ID_FIELD_RSU FROM $rsu;""".toString()

// The facade area of buildings being entirely included in the RSU buffer is calculated
Expand Down Expand Up @@ -1337,13 +1337,15 @@ String smallestCommunGeometry(JdbcDataSource datasource, String zone, String id_
//Separate road features according the zindex
def roadTable_zindex0_buffer = postfix "road_zindex0_buffer"
def road_tmp = postfix "road_zindex0"

datasource """DROP TABLE IF EXISTS $roadTable_zindex0_buffer, $road_tmp;
CREATE TABLE $roadTable_zindex0_buffer as SELECT st_buffer(the_geom, WIDTH::DOUBLE PRECISION/2)
AS the_geom
FROM $road where ZINDEX=0 ;
CREATE TABLE $roadTable_zindex0_buffer as SELECT ST_CollectionExtract(st_intersection(a.the_geom,b.the_geom),2) AS the_geom,
a.WIDTH, b.${id_zone}
FROM $road as a, $zone AS b WHERE a.the_geom && b.the_geom AND st_intersects(a.the_geom, b.the_geom) and a.ZINDEX=0 ;
CREATE SPATIAL INDEX IF NOT EXISTS ids_$roadTable_zindex0_buffer ON $roadTable_zindex0_buffer(the_geom);
CREATE TABLE $road_tmp AS SELECT ST_CollectionExtract(st_intersection(st_union(st_accum(a.the_geom)),b.the_geom),3) AS the_geom, b.${id_zone} FROM
$roadTable_zindex0_buffer AS a, $zone AS b WHERE a.the_geom && b.the_geom AND st_intersects(a.the_geom, b.the_geom) GROUP BY b.${id_zone};
CREATE TABLE $road_tmp AS SELECT st_union(st_accum(st_buffer(a.the_geom, WIDTH::DOUBLE PRECISION/2,2)))AS the_geom,
${id_zone} FROM
$roadTable_zindex0_buffer AS a group by ${id_zone} ;
DROP TABLE IF EXISTS $roadTable_zindex0_buffer;
""".toString()
tablesToMerge += ["$road_tmp": "select ST_ToMultiLine(the_geom) as the_geom, ${id_zone} from $road_tmp WHERE ST_ISEMPTY(THE_GEOM)=false"]
Expand All @@ -1356,12 +1358,13 @@ String smallestCommunGeometry(JdbcDataSource datasource, String zone, String id_
def railTable_zindex0_buffer = postfix "rail_zindex0_buffer"
def rail_tmp = postfix "rail_zindex0"
datasource """DROP TABLE IF EXISTS $railTable_zindex0_buffer, $rail_tmp;
CREATE TABLE $railTable_zindex0_buffer as SELECT st_buffer(the_geom, WIDTH::DOUBLE PRECISION/2)
AS the_geom
FROM $rail where ZINDEX=0 ;
CREATE TABLE $railTable_zindex0_buffer as SELECT ST_CollectionExtract(st_intersection(a.the_geom,b.the_geom),3) AS the_geom,
a.WIDTH, b.${id_zone}
FROM $rail as a ,$zone AS b WHERE a.the_geom && b.the_geom AND st_intersects(a.the_geom, b.the_geom) and a.ZINDEX=0 ;
CREATE SPATIAL INDEX IF NOT EXISTS ids_$railTable_zindex0_buffer ON $railTable_zindex0_buffer(the_geom);
CREATE TABLE $rail_tmp AS SELECT ST_CollectionExtract(st_intersection(st_union(st_accum(a.the_geom)),b.the_geom),3) AS the_geom, b.${id_zone} FROM
$railTable_zindex0_buffer AS a, $zone AS b WHERE a.the_geom && b.the_geom AND st_intersects(a.the_geom, b.the_geom) GROUP BY b.${id_zone};
CREATE TABLE $rail_tmp AS SELECT st_union(st_accum(st_buffer(a.the_geom, WIDTH::DOUBLE PRECISION/2,2))) AS the_geom,
${id_zone} FROM
$railTable_zindex0_buffer AS a GROUP BY ${id_zone};
DROP TABLE IF EXISTS $railTable_zindex0_buffer;
""".toString()
tablesToMerge += ["$rail_tmp": "select ST_ToMultiLine(the_geom) as the_geom, ${id_zone} from $rail_tmp WHERE ST_ISEMPTY(THE_GEOM)=false"]
Expand Down Expand Up @@ -1435,7 +1438,7 @@ String smallestCommunGeometry(JdbcDataSource datasource, String zone, String id_
def tmp_point_polygonize = postfix "tmp_point_polygonize_zindex0"
datasource """DROP TABLE IF EXISTS $tmp_point_polygonize;
CREATE INDEX ON $tmp_tables($id_zone);
CREATE TABLE $tmp_point_polygonize as select EXPLOD_ID as ${ID_COLUMN_NAME}, st_pointonsurface(st_force2D(the_geom)) as the_geom ,
CREATE TABLE $tmp_point_polygonize as select EXPLOD_ID as ${ID_COLUMN_NAME}, st_pointonsurface(the_geom) as the_geom ,
st_area(the_geom) as area , ${id_zone}
from st_explode ('(select st_polygonize(st_union(st_force2d(
st_precisionreducer(st_node(st_accum(a.the_geom)), 3)))) as the_geom, ${id_zone} from $tmp_tables as a group by ${id_zone})')""".toString()
Expand Down Expand Up @@ -2254,7 +2257,7 @@ String groundLayer(JdbcDataSource datasource, String zone, String id_zone,
def roadTable_zindex0_buffer = postfix "road_zindex0_buffer"
def road_tmp = postfix "road_zindex0"
datasource """DROP TABLE IF EXISTS $roadTable_zindex0_buffer, $road_tmp;
CREATE TABLE $roadTable_zindex0_buffer as SELECT st_buffer(the_geom, WIDTH::DOUBLE PRECISION/2)
CREATE TABLE $roadTable_zindex0_buffer as SELECT st_buffer(the_geom, WIDTH::DOUBLE PRECISION/2, 2)
AS the_geom, surface as type
FROM $road where ZINDEX=0 ;
CREATE SPATIAL INDEX IF NOT EXISTS ids_$roadTable_zindex0_buffer ON $roadTable_zindex0_buffer(the_geom);
Expand Down
73 changes: 44 additions & 29 deletions geoindicators/src/main/groovy/org/orbisgis/geoclimate/geoindicators/SpatialUnits.groovy
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -74,6 +74,7 @@ String createTSU(JdbcDataSource datasource, String zone,
def tsuDataPrepared = prepareTSUData(datasource,
zone, road, rail,
vegetation, water, sea_land_mask, urban_areas, surface_vegetation, surface_hydro, surface_urban_areas, prefixName)

if (!tsuDataPrepared) {
info "Cannot prepare the data for RSU calculation."
return
Expand Down Expand Up @@ -123,30 +124,38 @@ String createTSU(JdbcDataSource datasource, String inputTableName, String inputz
return null
}
if (inputzone) {


//Clip the geometry
datasource.createSpatialIndex(inputTableName)
datasource.createSpatialIndex(inputzone)

String cut_lines = postfix("cut_lines")
datasource.execute("""DROP TABLE IF EXISTS $cut_lines;
CREATE TABLE $cut_lines as SELECT ST_INTERSECTION(a.the_geom, b.the_geom) as the_geom from
$inputTableName as a, $inputzone as b where a.the_geom && b.the_geom and st_intersects(a.the_geom, b.the_geom);
""".toString())

//Create the polygons from the TSU lines
String polygons = postfix("polygons")
// todo: put 'mitre' instead of 'round' whenever JTS has been updated to version 1.20
datasource.execute("""DROP TABLE IF EXISTS $polygons;
CREATE TABLE $polygons as SELECT EXPLOD_ID AS $COLUMN_ID_NAME, ST_SETSRID(ST_BUFFER(ST_BUFFER(the_geom,-0.01, 2), 0.01, 2), $epsg) AS the_geom FROM
datasource.execute(""" DROP TABLE IF EXISTS $polygons;
CREATE TABLE $polygons as
SELECT
ST_SETSRID(st_makevalid(ST_BUFFER(ST_BUFFER(the_geom,-0.01, 2), 0.01, 2)), $epsg) AS the_geom FROM
ST_EXPLODE('(SELECT ST_POLYGONIZE(ST_UNION(ST_NODE(ST_ACCUM(the_geom)))) AS the_geom
FROM $inputTableName)')
FROM $cut_lines)')
""".toString())
datasource.createSpatialIndex(polygons)
datasource.createSpatialIndex(inputTableName)
datasource """
DROP TABLE IF EXISTS $outputTableName;
CREATE TABLE $outputTableName AS
SELECT a.*
FROM $polygons AS a,
$inputzone AS b
WHERE a.the_geom && b.the_geom
AND ST_INTERSECTS(ST_POINTONSURFACE(a.THE_GEOM), b.the_geom) and st_area(a.the_geom) > $area;
DROP TABLE IF EXISTS $polygons;
""".toString()

datasource.execute("""DROP TABLE IF EXISTS $outputTableName;
CREATE TABLE $outputTableName as
select CAST((row_number() over()) as Integer) as $COLUMN_ID_NAME, the_geom from
st_explode('$polygons') where st_area(the_geom) > $area;
DROP TABLE IF EXISTS $cut_lines, $polygons""".toString())
} else {
datasource """
DROP TABLE IF EXISTS $outputTableName;
CREATE TABLE $outputTableName AS
SELECT EXPLOD_ID AS $COLUMN_ID_NAME, ST_SETSRID(ST_BUFFER(ST_BUFFER(the_geom,-0.01, 'join=mitre'), 0.01, 'join=mitre'), $epsg) AS the_geom
SELECT EXPLOD_ID AS $COLUMN_ID_NAME, ST_SETSRID(ST_BUFFER(ST_BUFFER(the_geom,-0.01, 2), 0.01, 2), $epsg) AS the_geom
FROM ST_EXPLODE('(
SELECT ST_POLYGONIZE(ST_UNION(ST_NODE(ST_ACCUM(the_geom)))) AS the_geom
FROM $inputTableName)') where st_area(the_geom) > $area""".toString()
Expand Down Expand Up @@ -339,7 +348,7 @@ String prepareTSUData(JdbcDataSource datasource, String zone, String road, Strin
}
}

if (water && datasource.hasTable(urban_areas)) {
if (urban_areas && datasource.hasTable(urban_areas)) {
if (datasource.getColumnNames(urban_areas).size() > 0) {
debug "Preparing urban areas..."
queryCreateOutputTable += [urban_areas_tmp: "(SELECT ST_ToMultiLine(THE_GEOM) FROM $urban_areas WHERE st_area(the_geom)>=$surface_urban_areas and type not in ('social_building'))"]
Expand Down Expand Up @@ -430,7 +439,7 @@ String createBlocks(JdbcDataSource datasource, String inputTableName, double sna
CREATE INDEX ON $subGraphTableNodes (NODE_ID);
DROP TABLE IF EXISTS $subGraphBlocks;
CREATE TABLE $subGraphBlocks AS
SELECT ST_UNION(ST_ACCUM(ST_buffer(A.THE_GEOM, $snappingTolerance))) AS THE_GEOM
SELECT ST_UNION(ST_ACCUM(ST_buffer(A.THE_GEOM, $snappingTolerance,2))) AS THE_GEOM
FROM $inputTableName A, $subGraphTableNodes B
WHERE A.id_build=B.NODE_ID GROUP BY B.CONNECTED_COMPONENT;
""".toString()
Expand All @@ -446,16 +455,22 @@ String createBlocks(JdbcDataSource datasource, String inputTableName, double sna
//Create the blocks
debug "Creating the block table..."

datasource """DROP TABLE IF EXISTS $outputTableName;
def blocks = postfix("blocks")
datasource.execute("""DROP TABLE IF EXISTS $blocks;
CREATE TABLE $blocks as
SELECT st_force2d(ST_MAKEVALID(THE_GEOM)) as the_geom FROM $subGraphBlocks
UNION ALL (SELECT st_force2d(ST_MAKEVALID(a.the_geom)) as the_geom FROM $inputTableName a
LEFT JOIN $subGraphTableNodes b ON a.id_build = b.NODE_ID WHERE b.NODE_ID IS NULL);""".toString())

//Don't forget to explode the blocks
datasource.execute("""DROP TABLE IF EXISTS $outputTableName;
CREATE TABLE $outputTableName ($columnIdName SERIAL, THE_GEOM GEOMETRY)
AS SELECT CAST((row_number() over()) as Integer), the_geom FROM
((SELECT st_force2d(ST_MAKEVALID(THE_GEOM)) as the_geom FROM $subGraphBlocks)
UNION ALL (SELECT st_force2d(ST_MAKEVALID(a.the_geom)) as the_geom FROM $inputTableName a
LEFT JOIN $subGraphTableNodes b ON a.id_build = b.NODE_ID WHERE b.NODE_ID IS NULL));""".toString()
AS SELECT CAST((row_number() over()) as Integer), the_geom FROM st_explode('$blocks')
where st_area(the_geom) > 0 ;""".toString())

// Temporary tables are deleted
datasource """DROP TABLE IF EXISTS $graphTable,
$subGraphBlocks, $subGraphTableNodes, $subGraphTableEdges;""".toString()
$subGraphBlocks, $subGraphTableNodes, $subGraphTableEdges,$blocks;""".toString()

debug "The blocks have been created"
return outputTableName
Expand Down Expand Up @@ -510,8 +525,8 @@ String spatialJoin(JdbcDataSource datasource, String sourceTable, String targetT
WHERE a.$GEOMETRIC_COLUMN_SOURCE && b.$GEOMETRIC_COLUMN_TARGET AND
ST_INTERSECTS(a.$GEOMETRIC_COLUMN_SOURCE,
b.$GEOMETRIC_COLUMN_TARGET)
ORDER BY ST_AREA(ST_INTERSECTION(ST_BUFFER(ST_PRECISIONREDUCER(a.$GEOMETRIC_COLUMN_SOURCE, 3),0),
ST_BUFFER(ST_PRECISIONREDUCER(b.$GEOMETRIC_COLUMN_TARGET,3),0)))
ORDER BY ST_AREA(ST_INTERSECTION(a.$GEOMETRIC_COLUMN_SOURCE,
b.$GEOMETRIC_COLUMN_TARGET))
DESC LIMIT $nbRelations) AS $idColumnTarget
FROM $sourceTable a""".toString()
} else {
Expand Down Expand Up @@ -662,7 +677,7 @@ String computeSprawlAreas(JdbcDataSource datasource, String grid_indicators, flo
select st_union(st_accum(the_geom)) as the_geom from
$grid_indicators where lcz_warm>=2
and LCZ_PRIMARY NOT IN (101, 102,103,104,106, 107))')
where st_isempty(st_buffer(the_geom, -100)) =false""".toString())
where st_isempty(st_buffer(the_geom, -100,2)) =false""".toString())
datasource.execute("""CREATE TABLE $outputTableName as SELECT CAST((row_number() over()) as Integer) as id,
the_geom
FROM
Expand Down Expand Up @@ -727,7 +742,7 @@ String extractCoolAreas(JdbcDataSource datasource, String grid_indicators,float
CREATE TABLE $outputTableName as SELECT CAST((row_number() over()) as Integer) as id, the_geom FROM ST_EXPLODE('(
SELECT ST_UNION(ST_ACCUM(a.THE_GEOM)) AS THE_GEOM FROM $grid_indicators as a
where
a.LCZ_PRIMARY in (101, 102, 103,104, 106, 107))') ${distance>0?" where st_isempty(st_buffer(the_geom, -$distance)) =false":""};
a.LCZ_PRIMARY in (101, 102, 103,104, 106, 107))') ${distance>0?" where st_isempty(st_buffer(the_geom, -$distance,2)) =false":""};
""".toString())
return outputTableName
}
Expand Down
8 changes: 4 additions & 4 deletions ...cators/src/main/groovy/org/orbisgis/geoclimate/geoindicators/WorkflowGeoIndicators.groovy
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -385,7 +385,7 @@ String computeRSUIndicators(JdbcDataSource datasource, String buildingTable,
"abbey", "condominium", "villa", "dormitory", "sheltered_housing",
"workers_dormitory", "terrace", "residential"]],
floorAreaTypeAndCompositionLcz : [:],
buildingAreaTypeAndCompositionTeb : ["individual_housing": ["house", "detached", "bungalow", "farm", "villa", "terrace"],
buildingAreaTypeAndCompositionTeb : ["individual_housing": ["house", "detached", "bungalow", "farm", "villa", "terrace", "cabin"],
"collective_housing": [ "apartments","barracks","abbey", "dormitory",
"sheltered_housing", "workers_dormitory",
"condominium"],
Expand All @@ -412,8 +412,8 @@ String computeRSUIndicators(JdbcDataSource datasource, String buildingTable,
"ammunition","bunker","casemate","shelter","religious","place_of_worship",
"wayside_shrine","station","stable","sty","greenhouse","kiosk","marketplace",
"marker","warehouse","planetarium","fire_station","water_tower","grandstand",
"transportation","toll_booth"]],
floorAreaTypeAndCompositionTeb : ["individual_housing": ["house", "detached", "bungalow", "farm", "villa", "terrace"],
"transportation","toll_booth", "hut", "shed"]],
floorAreaTypeAndCompositionTeb : ["individual_housing": ["house", "detached", "bungalow", "farm", "villa", "terrace", "cabin"],
"collective_housing": [ "apartments","barracks","abbey", "dormitory",
"sheltered_housing", "workers_dormitory",
"condominium"],
Expand All @@ -440,7 +440,7 @@ String computeRSUIndicators(JdbcDataSource datasource, String buildingTable,
"ammunition","bunker","casemate","shelter","religious","place_of_worship",
"wayside_shrine","station","stable","sty","greenhouse","kiosk","marketplace",
"marker","warehouse","planetarium","fire_station","water_tower","grandstand",
"transportation","toll_booth"]],
"transportation","toll_booth", "hut", "shed"]],
utrfSurfFraction : ["vegetation_fraction_utrf" : ["high_vegetation_fraction",
"low_vegetation_fraction",
"high_vegetation_low_vegetation_fraction",
Expand Down
Loading
Loading