Skip to content

Commit

Permalink
Merge pull request #924 from ebocher/osm_date
Browse files Browse the repository at this point in the history 
OSM and BDTopo data extraction options
  • Loading branch information
@ebocher
ebocher authored Jan 18, 2024
2 parents 23b290f + c4ec7ca commit 720a7ad
Show file tree
Hide file tree
Showing 5 changed files with 113 additions and 8 deletions.
17 changes: 16 additions & 1 deletion bdtopo/src/main/groovy/org/orbisgis/geoclimate/bdtopo/BDTopoV2Workflow.groovy
View file
Original file line number Diff line number Diff line change
Expand Up @@ -85,7 +85,13 @@ 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) {
//def tmp_insee = code.join("_")
if(code.size()==3){
if(code[2]<100){
error("The distance to create a bbox from a point must be greater than 100 meters")
return
}
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
st_setsrid(the_geom, $commune_srid)
Expand Down Expand Up @@ -306,6 +312,13 @@ 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){
error("The distance to create a bbox from a point must be greater than 100 meters")
return
}
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
ST_INTERSECTION(the_geom, ST_MakeEnvelope(${location[1]},${location[0]},${location[3]},${location[2]}, $sourceSRID)) as the_geom, CODE_INSEE from ${inputTables.commune} where the_geom
Expand All @@ -314,6 +327,8 @@ 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{
return
}
def count = h2gis_datasource."$outputTableName".rowCount
if (count > 0) {
Expand Down
15 changes: 15 additions & 0 deletions bdtopo/src/main/groovy/org/orbisgis/geoclimate/bdtopo/BDTopoV3Workflow.groovy
View file
Original file line number Diff line number Diff line change
Expand Up @@ -55,6 +55,14 @@ def filterLinkedShapeFiles(def location, float distance, LinkedHashMap inputTabl
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){
error("The distance to create a bbox from a point must be greater than 100 meters")
return
}
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())
Expand Down Expand Up @@ -261,6 +269,13 @@ Integer loadDataFromPostGIS(Object input_database_properties, Object code, Objec
if (code in Collection) {
def communeColumns = h2gis_datasource.getColumnNames(commune_location)
if(communeColumns.contains("INSEE_COM")) {
if(code.size()==3){
if(code[2]<100){
error("The distance to create a bbox from a point must be greater than 100 meters")
return
}
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
st_setsrid(the_geom, $commune_srid)
Expand Down
18 changes: 12 additions & 6 deletions osm/src/main/groovy/org/orbisgis/geoclimate/osm/WorkflowOSM.groovy
View file
Original file line number Diff line number Diff line change
Expand Up @@ -232,6 +232,9 @@ Map workflow(def input) {
error "The all parameter must be a boolean value"
return null
}

def osm_date= inputParameters.get("date")

def osm_size_area = inputParameters.get("area")

if (!osm_size_area) {
Expand All @@ -241,7 +244,6 @@ Map workflow(def input) {
error "The area of the bounding box to be extracted from OSM must be greater than 0 km²"
return null
}

def overpass_timeout = inputParameters.get("timeout")
if (!overpass_timeout) {
overpass_timeout = 900
Expand Down Expand Up @@ -360,7 +362,7 @@ Map workflow(def input) {
def logTableZones = postfix("log_zones")
Map osmprocessing = osm_processing(h2gis_datasource, processing_parameters, locations.findAll { it }, file_outputFolder, outputFileTables,
outputDatasource, outputTables, outputSRID, downloadAllOSMData, deleteOutputData, deleteOSMFile, logTableZones, osm_size_area,
overpass_timeout, overpass_maxsize)
overpass_timeout, overpass_maxsize, osm_date)
if (!osmprocessing) {
h2gis_datasource.save(logTableZones,"${file_outputFolder.getAbsolutePath() + File.separator}logzones.geojson", true)
return null
Expand Down Expand Up @@ -401,7 +403,7 @@ Map osm_processing(JdbcDataSource h2gis_datasource, def processing_parameters, d
def outputSRID, def downloadAllOSMData,
def deleteOutputData, def deleteOSMFile,
def logTableZones, def bbox_size,
def overpass_timeout, def overpass_maxsize) {
def overpass_timeout, def overpass_maxsize,def overpass_date) {
//Store the zone identifier and the names of the tables
def outputTableNamesResult = [:]
//Create the table to log on the processed zone
Expand All @@ -421,7 +423,7 @@ Map osm_processing(JdbcDataSource h2gis_datasource, def processing_parameters, d
error "Cannot find any geometry to define the zone to extract the OSM data"
return
}
def srid = h2gis_datasource.getSpatialTable(zone).srid
def srid = h2gis_datasource.getSrid(zone)
def reproject = false
if (outputSRID) {
if (outputSRID != srid) {
Expand All @@ -432,7 +434,11 @@ Map osm_processing(JdbcDataSource h2gis_datasource, def processing_parameters, d
}
//Prepare OSM extraction
//TODO set key values ?
def query = "[timeout:$overpass_timeout][maxsize:$overpass_maxsize]" + OSMTools.Utilities.buildOSMQuery(zones.envelope, null, OSMElement.NODE, OSMElement.WAY, OSMElement.RELATION)
def osm_date=""
if(overpass_date){
osm_date = "[date:\"$overpass_date\"]"
}
def query = "[timeout:$overpass_timeout][maxsize:$overpass_maxsize]$osm_date" + OSMTools.Utilities.buildOSMQuery(zones.envelope, null, OSMElement.NODE, OSMElement.WAY, OSMElement.RELATION)

if (downloadAllOSMData) {
//Create a custom OSM query to download all requiered data. It will take more time and resources
Expand All @@ -441,7 +447,7 @@ Map osm_processing(JdbcDataSource h2gis_datasource, def processing_parameters, d
"leisure", "highway", "natural",
"landuse", "landcover",
"vegetation", "waterway", "area", "aeroway", "area:aeroway", "tourism", "sport"]
query = "[timeout:$overpass_timeout][maxsize:$overpass_maxsize]" + OSMTools.Utilities.buildOSMQueryWithAllData(zones.envelope, keysValues, OSMElement.NODE, OSMElement.WAY, OSMElement.RELATION)
query = "[timeout:$overpass_timeout][maxsize:$overpass_maxsize]$osm_date" + OSMTools.Utilities.buildOSMQueryWithAllData(zones.envelope, keysValues, OSMElement.NODE, OSMElement.WAY, OSMElement.RELATION)
}

def extract = OSMTools.Loader.extract(query)
Expand Down
61 changes: 61 additions & 0 deletions osm/src/test/groovy/org/orbisgis/geoclimate/osm/WorflowOSMTest.groovy
View file
Original file line number Diff line number Diff line change
Expand Up @@ -663,6 +663,7 @@ class WorflowOSMTest extends WorkflowAbstractTest {
"input" : [
"locations": [location],//["Pont-de-Veyle"],//[nominatim["bbox"]],//["Lorient"],
"area": 2800,
//"date":"2017-12-31T19:20:00Z",
/*"timeout":182,
"maxsize": 536870918,
"endpoint":"https://lz4.overpass-api.de/api"*/],
Expand Down Expand Up @@ -809,6 +810,66 @@ class WorflowOSMTest extends WorkflowAbstractTest {
assertTrue h2gis.firstRow("select count(*) as count from $building where HEIGHT_WALL>0 and HEIGHT_ROOF>0").count > 0
}

@Test
void testEstimateBuildingWithAllInputHeightFromPoint() {
String directory = folder.absolutePath + File.separator + "test_building_height"
File dirFile = new File(directory)
dirFile.delete()
dirFile.mkdir()
def osm_parmeters = [
"geoclimatedb": [
"folder": dirFile.absolutePath,
"name" : "geoclimate_chain_db",
"delete": false
],
"input" : [
"locations": [[43.726898,7.298452,100]]],
"output" : [
"folder": ["path" : directory,
"tables": ["building", "zone"]]],
"parameters" :
[
rsu_indicators: ["indicatorUse" : ["LCZ"]]
]
]
Map process = OSM.WorkflowOSM.workflow(osm_parmeters)
def tableNames = process.values()
def building = tableNames.building[0]
H2GIS h2gis = H2GIS.open("${directory + File.separator}geoclimate_chain_db;AUTO_SERVER=TRUE")
assertTrue h2gis.firstRow("select count(*) as count from $building where HEIGHT_WALL>0 and HEIGHT_ROOF>0").count > 0
}

@Disabled //Because it takes some time to build the OSM query
@Test
void testEstimateBuildingWithAllInputHeightDate() {
String directory = folder.absolutePath + File.separator + "test_building_height"
File dirFile = new File(directory)
dirFile.delete()
dirFile.mkdir()
def osm_parmeters = [
"geoclimatedb": [
"folder": dirFile.absolutePath,
"name" : "geoclimate_chain_db",
"delete": false
],
"input" : [
"locations": [[43.726898,7.298452,43.727677,7.299632]],
"date":"2015-12-31T19:20:00Z"],
"output" : [
"folder": ["path" : directory,
"tables": ["building", "zone"]]],
"parameters" :
["distance" : 0,
rsu_indicators: ["indicatorUse" : ["LCZ"]]
]
]
Map process = OSM.WorkflowOSM.workflow(osm_parmeters)
def tableNames = process.values()
def building = tableNames.building[0]
H2GIS h2gis = H2GIS.open("${directory + File.separator}geoclimate_chain_db;AUTO_SERVER=TRUE")
assertTrue h2gis.firstRow("select count(*) as count from $building where HEIGHT_WALL>0 and HEIGHT_ROOF>0").count > 0
}

@Test
void testOneSeaLCZ() {
String directory = folder.absolutePath + File.separator + "test_sea_lcz"
Expand Down
10 changes: 9 additions & 1 deletion osmtools/src/main/groovy/org/orbisgis/geoclimate/osmtools/utils/Utilities.groovy
View file
Original file line number Diff line number Diff line change
Expand Up @@ -638,7 +638,8 @@ Geometry geometryFromNominatim(def bbox) {
*
* @author Erwan Bocher (CNRS LAB-STICC)
*
* @param bbox 4 values to define a bbox
* @param bbox 4 values to define a bbox or 3 values to define a point (lat/long)
* with a distance around it, expressed in meters
* @return a JTS polygon
*/
Geometry geometryFromValues(def bbox) {
Expand All @@ -651,6 +652,13 @@ Geometry geometryFromValues(def bbox) {
if (bbox.size() == 4) {
return buildGeometry([bbox[1], bbox[0], bbox[3], bbox[2]]);
}
else if (bbox.size()==3){
if(bbox[2]<100){
error("The distance to create a bbox from a point must be greater than 100 meters")
return null
}
return getAreaFromPoint(bbox[1], bbox[0], bbox[2])
}
}


Expand Down

0 comments on commit 720a7ad

Please sign in to comment.