Merge pull request #996 from ebocher/incorrect_osm_geom

Fix the issues #994 and #995

docs/CHANGELOG.md

@@ -2,3 +2,5 @@
 
 - Upgrade dependencies H2, H2GIS, JTS and Groovy
 - Fix some tests due to some rounding precision in JTS 1.20
+- Fix bad OSM geometry representation #994
+- Fix height value parsing with OSM #995

osm/src/main/groovy/org/orbisgis/geoclimate/osm/InputDataFormatting.groovy

@@ -846,36 +846,49 @@ static float getHeightRoof(height, heightPattern) {
         def match2_group1 = matcher.group(1)
         def match2_group2 = matcher.group(2)
         if (match1_group1) {
-            new_h = Float.parseFloat(match1_group1) * 12
+            new_h = parseFloat(match1_group1) * 12
         }
         if (match2_group2 == "''") {
-            new_h += Float.parseFloat(match2_group1)
+            new_h += parseFloat(match2_group1)
         }
         return new_h * 0.0254
 
     } else {
         if (match1_group1 && match1_group2 == null) {
-            return Float.parseFloat(match1_group1)
+               return parseFloat(match1_group1)
         }
         //next mach for feet, inch matcher.find();
         else {
             def type = match1_group2.toLowerCase()
             switch (type) {
                 case "m":
-                    return Float.parseFloat(match1_group1)
+                    return parseFloat(match1_group1)
                 case "foot":
-                    return Float.parseFloat(match1_group1) * 0.3048
+                    return parseFloat(match1_group1) * 0.3048
                 case "'":
-                    return Float.parseFloat(match1_group1) * 12 * 0.0254
+                    return parseFloat(match1_group1) * 12 * 0.0254
                 case "''":
-                    return Float.parseFloat(match1_group1) * 0.0254
+                    return parseFloat(match1_group1) * 0.0254
                 default:
                     return 0
             }
         }
     }
 }
 
+/**
+ * Parse to float otherwise return 0
+ * @param value to parse
+ * @return a float value
+ */
+static Float parseFloat(def value){
+    try {
+        return Float.parseFloat(value)
+    }catch (NumberFormatException ex){
+        return 0
+    }
+}
+
 /**
  * This function defines the value of the column nb_lev according to the values of b_lev and r_lev
  * @param row The row of the raw table to examine

osm/src/test/groovy/org/orbisgis/geoclimate/osm/InputDataFormattingTest.groovy

@@ -61,20 +61,20 @@ class InputDataFormattingTest {
         Map extractData = OSM.InputDataLoading.createGISLayers(
                 h2GIS, new File(this.class.getResource("redon.osm").toURI()).getAbsolutePath(), epsg)
 
-        assertEquals 1038, h2GIS.getTable(extractData.building).rowCount
+        assertEquals 1034, h2GIS.getTable(extractData.building).rowCount
         assertEquals 211, h2GIS.getTable(extractData.road).rowCount
         assertEquals 44, h2GIS.getTable(extractData.rail).rowCount
-        assertEquals 136, h2GIS.getTable(extractData.vegetation).rowCount
+        assertEquals 135, h2GIS.getTable(extractData.vegetation).rowCount
         assertEquals 10, h2GIS.getTable(extractData.water).rowCount
         assertEquals 47, h2GIS.getTable(extractData.impervious).rowCount
-        assertEquals 11, h2GIS.getTable(extractData.urban_areas).rowCount
+        assertEquals 10, h2GIS.getTable(extractData.urban_areas).rowCount
         assertEquals 0, h2GIS.getTable(extractData.coastline).rowCount
 
         //Buildings
         Map buildingLayers = OSM.InputDataFormatting.formatBuildingLayer(h2GIS, extractData.building)
         String building = buildingLayers.building
         assertNotNull h2GIS.getTable(building).save(new File(folder, "osm_building_formated.shp").absolutePath, true)
-        assertEquals 1038, h2GIS.getTable(building).rowCount
+        assertEquals 1028, h2GIS.getTable(building).rowCount
         assertTrue h2GIS.firstRow("select count(*) as count from ${building} where NB_LEV is null".toString()).count == 0
         assertTrue h2GIS.firstRow("select count(*) as count from ${building} where NB_LEV<0".toString()).count == 0
         assertTrue h2GIS.firstRow("select count(*) as count from ${building} where HEIGHT_WALL is null".toString()).count == 0
@@ -92,12 +92,12 @@ class InputDataFormattingTest {
         assertNotNull h2GIS.getTable(buiding_imp).save(new File(folder, "osm_building_formated_type.shp").absolutePath, true)
         def rows = h2GIS.rows("select type from ${buiding_imp} where id_build=158 or id_build=982".toString())
         assertEquals(2, rows.size())
-        assertTrue(rows.type == ['residential', 'residential'])
+        assertTrue(rows.type == ['residential', 'slight_construction'])
 
         rows = h2GIS.rows("select type from ${buiding_imp} where id_build=881 or id_build=484 or id_build=610".toString())
 
         assertEquals(3, rows.size())
-        assertTrue(rows.type == ['industrial', 'industrial', 'industrial'])
+        assertTrue(rows.type == ['industrial', 'building', 'slight_construction'])
 
 
         //Roads
@@ -117,7 +117,7 @@ class InputDataFormattingTest {
         //Vegetation
         String vegetation = OSM.InputDataFormatting.formatVegetationLayer(h2GIS, extractData.vegetation)
         assertNotNull h2GIS.getTable(vegetation).save(new File(folder, "osm_vegetation_formated.shp").absolutePath, true)
-        assertEquals 140, h2GIS.getTable(vegetation).rowCount
+        assertEquals 138, h2GIS.getTable(vegetation).rowCount
         assertTrue h2GIS.firstRow("select count(*) as count from ${vegetation} where type is null".toString()).count == 0
         assertTrue h2GIS.firstRow("select count(*) as count from ${vegetation} where HEIGHT_CLASS is null".toString()).count == 0
 
@@ -246,32 +246,32 @@ class InputDataFormattingTest {
         def epsg = 2154
         Map extractData = OSM.InputDataLoading.createGISLayers(h2GIS, new File(this.class.getResource("redon.osm").toURI()).getAbsolutePath(), epsg)
 
-        assertEquals 1038, h2GIS.getTable(extractData.building).rowCount
+        assertEquals 1034, h2GIS.getTable(extractData.building).rowCount
         assertEquals 211, h2GIS.getTable(extractData.road).rowCount
         assertEquals 44, h2GIS.getTable(extractData.rail).rowCount
-        assertEquals 136, h2GIS.getTable(extractData.vegetation).rowCount
+        assertEquals 135, h2GIS.getTable(extractData.vegetation).rowCount
         assertEquals 10, h2GIS.getTable(extractData.water).rowCount
         assertEquals 47, h2GIS.getTable(extractData.impervious).rowCount
 
         //Buildings with estimation state
         Map buildingLayers = OSM.InputDataFormatting.formatBuildingLayer(h2GIS, extractData.building)
         String buildingLayer = buildingLayers.building
         assertNotNull h2GIS.getTable(buildingLayer).save(new File(folder, "osm_building_formated.shp").absolutePath, true)
-        assertEquals 1038, h2GIS.getTable(buildingLayer).rowCount
+        assertEquals 1028, h2GIS.getTable(buildingLayer).rowCount
         assertTrue h2GIS.firstRow("select count(*) as count from ${buildingLayer} where NB_LEV is null").count == 0
         assertTrue h2GIS.firstRow("select count(*) as count from ${buildingLayer} where NB_LEV<0").count == 0
         assertTrue h2GIS.firstRow("select count(*) as count from ${buildingLayer} where NB_LEV=0").count == 0
         assertTrue h2GIS.firstRow("select count(*) as count from ${buildingLayer} where HEIGHT_WALL is null").count == 0
         assertTrue h2GIS.firstRow("select count(*) as count from ${buildingLayer} where HEIGHT_WALL<0").count == 0
         assertTrue h2GIS.firstRow("select count(*) as count from ${buildingLayer} where HEIGHT_ROOF is null").count == 0
         assertTrue h2GIS.firstRow("select count(*) as count from ${buildingLayer} where HEIGHT_ROOF<0").count == 0
-        assertEquals 1033, h2GIS.getTable(buildingLayers.building_estimated).rowCount
-        assertTrue h2GIS.firstRow("select count(*) as count from ${buildingLayers.building} join ${buildingLayers.building_estimated} using (id_build, id_source) where 1=1").count == 1033
+        assertEquals 1023, h2GIS.getTable(buildingLayers.building_estimated).rowCount
+        assertEquals(1023, h2GIS.firstRow("select count(*) as count from ${buildingLayers.building} join ${buildingLayers.building_estimated} using (id_build, id_source) where 1=1").count)
 
         //Buildings without estimation state
         buildingLayers = OSM.InputDataFormatting.formatBuildingLayer(h2GIS, extractData.building)
-        assertEquals 1038, h2GIS.getTable(buildingLayers.building).rowCount
-        assertEquals 1033, h2GIS.getTable(buildingLayers.building_estimated).rowCount
+        assertEquals 1028, h2GIS.getTable(buildingLayers.building).rowCount
+        assertEquals 1023, h2GIS.getTable(buildingLayers.building_estimated).rowCount
     }
 
 
@@ -294,7 +294,7 @@ class InputDataFormattingTest {
 
         //zoneToExtract = [62.2, 28.2, 62.4, 28.4]
 
-        zoneToExtract =[47.0619, -1.8145005, 47.394558, -1.2849174]
+        zoneToExtract =[51.328681,1.195128,51.331121,1.199162]
         Map extractData = OSM.InputDataLoading.extractAndCreateGISLayers(h2GIS, zoneToExtract)
 
         String formatedPlaceName = zoneToExtract.join("_").trim().split("\\s*(,|\\s)\\s*").join("_");
@@ -383,4 +383,10 @@ class InputDataFormattingTest {
         def road = OSM.InputDataFormatting.formatRoadLayer(h2GIS, gISLayers.road)
         h2GIS.getTable(road).save("/tmp/formated_osm_road.shp", true)
     }
+
+    @Test
+    void parseFloat() {
+        def heightPattern = Pattern.compile("((?:\\d+\\/|(?:\\d+|^|\\s)\\.)?\\d+)\\s*([^\\s\\d+\\-.,:;^\\/]+(?:\\^\\d+(?:\$|(?=[\\s:;\\/])))?(?:\\/[^\\s\\d+\\-.,:;^\\/]+(?:\\^\\d+(?:\$|(?=[\\s:;\\/])))?)*)?", Pattern.CASE_INSENSITIVE)
+        assertEquals(0, InputDataFormatting.getHeightRoof("II OSK 1559/12", heightPattern))
+    }
 }

osm/src/test/groovy/org/orbisgis/geoclimate/osm/InputDataLoadingTest.groovy

@@ -55,9 +55,8 @@ class InputDataLoadingTest {
     void createGISLayersTest() {
         def osmfile = new File(this.class.getResource("redon.osm").toURI()).getAbsolutePath()
         Map extract = OSM.InputDataLoading.createGISLayers(h2GIS, osmfile, 2154)
-        assertEquals 1038, h2GIS.getTable(extract.building).rowCount
-
-        assertEquals 136, h2GIS.getTable(extract.vegetation).rowCount
+        assertEquals 1034, h2GIS.getTable(extract.building).rowCount
+        assertEquals 135, h2GIS.getTable(extract.vegetation).rowCount
         assertEquals 211, h2GIS.getTable(extract.road).rowCount
 
         assertEquals 44, h2GIS.getTable(extract.rail).rowCount
@@ -66,7 +65,7 @@ class InputDataLoadingTest {
 
         assertEquals 47, h2GIS.getTable(extract.impervious).rowCount
 
-        assertEquals 11, h2GIS.getTable(extract.urban_areas).rowCount
+        assertEquals 10, h2GIS.getTable(extract.urban_areas).rowCount
     }
 
     //This test is used for debug purpose

osm/src/test/groovy/org/orbisgis/geoclimate/osm/WorflowOSMTest.groovy

@@ -653,7 +653,7 @@ class WorflowOSMTest extends WorkflowAbstractTest {
         def nominatim = org.orbisgis.geoclimate.osmtools.OSMTools.Utilities.getNominatimData(location)
         def grid_size = 250
         location = nominatim.bbox
-        //location=[46.178404,6.095524,46.222959,6.190109]
+        location=[51.2, 1.0, 51.4, 1.2]
         /* location =[ 48.84017284026897,
                     2.3061887733275785,
                     48.878115442982086,
@@ -678,7 +678,7 @@ class WorflowOSMTest extends WorkflowAbstractTest {
                 "parameters"  :
                         ["distance"             : 0,
                          "rsu_indicators"       : [
-                                 "indicatorUse": ["LCZ"]//, "TEB"] //, "UTRF"]
+                                 "indicatorUse": ["LCZ", "TEB"] //, "UTRF"]
 
                          ], "grid_indicators"   : [
                                 "x_size"    : grid_size,

osmtools/src/main/groovy/org/orbisgis/geoclimate/osmtools/Transform.groovy

@@ -330,7 +330,7 @@ String extractWaysAsPolygons(JdbcDataSource datasource, String osmTablesPrefix,
         query += """ a.the_geom ${OSMTools.TransformUtils.createTagList(datasource, columnsSelector, columnsToKeep)}  """
     }
 
-    query += " FROM $waysPolygonTmp AS a, $osmTableTag b WHERE a.id_way=b.id_way and st_isempty(a.the_geom)=false "
+    query += " FROM $waysPolygonTmp AS a, $osmTableTag b WHERE a.id_way=b.id_way and st_isempty(a.the_geom)=false and st_isvalid(a.the_geom) "
 
     if (columnsToKeep) {
         query += " AND b.TAG_KEY IN ('${columnsToKeep.join("','")}') "
@@ -567,7 +567,7 @@ def extractRelationsAsPolygons(JdbcDataSource datasource, String osmTablesPrefix
     } else {
         query += """ st_normalize(a.the_geom) as the_geom ${OSMTools.TransformUtils.createTagList(datasource, columnsSelector, columnsToKeep)} 
         FROM $relationsMpHoles AS a, ${osmTablesPrefix}_relation_tag  b 
-                    WHERE a.id_relation=b.id_relation and st_isempty(a.the_geom)=false  """
+                    WHERE a.id_relation=b.id_relation and st_isempty(a.the_geom)=false and st_isvalid(a.the_geom) """
     }
 
     if (columnsToKeep) {

osmtools/src/test/groovy/org/orbisgis/geoclimate/osmtools/AbstractOSMToolsTest.groovy

@@ -114,8 +114,8 @@ abstract class AbstractOSMToolsTest {
        ('POINT(0 10)', 3),('POINT(10 10)', 4);
 
        CREATE TABLE ${prefix}_way_node(id_way int, id_node int, node_order int);
-       INSERT INTO ${prefix}_way_node VALUES(1, 1, 1),(1, 2, 2),(1, 3, 3),
-       (1, 4, 4),(1, 1, 5);
+       INSERT INTO ${prefix}_way_node VALUES(1, 1, 1),(1, 2, 2),(1, 4, 3),
+       (1, 5, 4),(1, 1, 5);
 
        CREATE TABLE ${prefix}_way(id_way int);
        INSERT INTO ${prefix}_way VALUES(1);

osmtools/src/test/groovy/org/orbisgis/geoclimate/osmtools/TransformTest.groovy

@@ -257,7 +257,7 @@ class TransformTest extends AbstractOSMToolsTest {
     @Test
     void extractWaysAsPolygonsTest() {
         def prefix = "OSM_" + uuid()
-        def epsgCode = 2453
+        def epsgCode = 4326
         def tags = [building: "house"]
         def columnsToKeep = ["water", "building"]
 
@@ -364,7 +364,7 @@ class TransformTest extends AbstractOSMToolsTest {
 
         String result = OSMTools.Transform.extractRelationsAsPolygons(ds, prefix, epsgCode, tags, columnsToKeep)
         assertFalse result.isEmpty()
-        ds.save(result, "/tmp/building.geojson", true)
+        //ds.save(result, "/tmp/building.geojson", true)
         def table = ds.getTable(result)
         assertEquals 0, table.rowCount
 
@@ -562,13 +562,13 @@ class TransformTest extends AbstractOSMToolsTest {
         def tags = ["building"]
         String outputTableName = OSMTools.Transform.toPolygons(ds, prefix, 4326, tags)
         assertEquals 6, ds.firstRow("select count(*) as count from ${outputTableName} where ST_NumInteriorRings(the_geom)  > 0").count as int
-        assertEquals 1032, ds.firstRow("select count(*) as count from ${outputTableName} where ST_NumInteriorRings(the_geom)  = 0").count as int
+        assertEquals 1028, ds.firstRow("select count(*) as count from ${outputTableName} where ST_NumInteriorRings(the_geom)  = 0").count as int
 
         //Create landuse layer
         tags = ["landuse": ["farmland", "forest", "grass", "meadow", "orchard", "vineyard", "village_green", "allotments"],]
         outputTableName = OSMTools.Transform.toPolygons(ds, prefix, 4326, tags)
-        assertEquals 131, ds.firstRow("select count(*) as count from ${outputTableName}").count as int
-        assertEquals 123, ds.firstRow("select count(*) as count from ${outputTableName} where \"landuse\"='grass'").count as int
+        assertEquals 130, ds.firstRow("select count(*) as count from ${outputTableName}").count as int
+        assertEquals 122, ds.firstRow("select count(*) as count from ${outputTableName} where \"landuse\"='grass'").count as int
 
         //Create urban areas layer
         tags = ["landuse": [
@@ -579,8 +579,8 @@ class TransformTest extends AbstractOSMToolsTest {
         ]]
         outputTableName = OSMTools.Transform.toPolygons(ds, prefix, 4326, tags)
 
-        assertEquals 6, ds.firstRow("select count(*) as count from ${outputTableName}").count as int
-        assertEquals 4, ds.firstRow("select count(*) as count from ${outputTableName} where \"landuse\"='residential'").count as int
+        assertEquals 5, ds.firstRow("select count(*) as count from ${outputTableName}").count as int
+        assertEquals 3, ds.firstRow("select count(*) as count from ${outputTableName} where \"landuse\"='residential'").count as int
 
     }