From c553a946d20b911426f46f0d6553326556bde6c4 Mon Sep 17 00:00:00 2001 From: Peter Hanecak Date: Sat, 23 Mar 2024 19:15:20 +0100 Subject: [PATCH 1/5] OSM Lake IDs at low zooms --- .../java/org/openmaptiles/layers/Water.java | 199 +++++++++++++++++- .../org/openmaptiles/layers/WaterTest.java | 103 +++++++-- 2 files changed, 271 insertions(+), 31 deletions(-) diff --git a/src/main/java/org/openmaptiles/layers/Water.java b/src/main/java/org/openmaptiles/layers/Water.java index 8a3d9a9c..6655a9df 100644 --- a/src/main/java/org/openmaptiles/layers/Water.java +++ b/src/main/java/org/openmaptiles/layers/Water.java @@ -42,10 +42,20 @@ OR TORT (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE import com.onthegomap.planetiler.config.PlanetilerConfig; import com.onthegomap.planetiler.expression.MultiExpression; import com.onthegomap.planetiler.geo.GeometryException; +import com.onthegomap.planetiler.reader.SimpleFeature; import com.onthegomap.planetiler.reader.SourceFeature; import com.onthegomap.planetiler.stats.Stats; import com.onthegomap.planetiler.util.Translations; +import java.util.ArrayList; import java.util.List; +import java.util.Map; +import java.util.concurrent.ConcurrentHashMap; +import java.util.function.Consumer; +import org.locationtech.jts.geom.Envelope; +import org.locationtech.jts.geom.Geometry; +import org.locationtech.jts.geom.TopologyException; +import org.locationtech.jts.geom.util.GeometryFixer; +import org.locationtech.jts.index.strtree.STRtree; import org.openmaptiles.OpenMapTilesProfile; import org.openmaptiles.generated.OpenMapTilesSchema; import org.openmaptiles.generated.Tables; @@ -62,7 +72,8 @@ public class Water implements Tables.OsmWaterPolygon.Handler, OpenMapTilesProfile.NaturalEarthProcessor, OpenMapTilesProfile.OsmWaterPolygonProcessor, - ForwardingProfile.FeaturePostProcessor { + ForwardingProfile.FeaturePostProcessor, + OpenMapTilesProfile.FinishHandler { /* * At low zoom levels, use natural earth for oceans and major lakes, and at high zoom levels @@ -71,12 +82,20 @@ public class Water implements * which infers ocean polygons by preprocessing all coastline elements. */ + // smallest NE lake is around 4.42E-13, this is roughly 1/2 of that and approx 17% of OSM features are smaller than this: + private static final double SMALLEST_OSM_LAKE_AREA = Math.pow(4, -21); + private final MultiExpression.Index classMapping; private final PlanetilerConfig config; + private final Stats stats; + private final Map neLakeIndexes = new ConcurrentHashMap<>(); + private final Map> neLakeNameMaps = new ConcurrentHashMap<>(); + private final List neAllLakeInfos = new ArrayList<>(); public Water(Translations translations, PlanetilerConfig config, Stats stats) { this.classMapping = FieldMappings.Class.index(); this.config = config; + this.stats = stats; } @Override @@ -86,21 +105,57 @@ record WaterInfo(int minZoom, int maxZoom, String clazz) {} case "ne_110m_ocean" -> new WaterInfo(0, 1, FieldValues.CLASS_OCEAN); case "ne_50m_ocean" -> new WaterInfo(2, 4, FieldValues.CLASS_OCEAN); case "ne_10m_ocean" -> new WaterInfo(5, 5, FieldValues.CLASS_OCEAN); - - // TODO: get OSM ID from low-zoom natural earth lakes - case "ne_110m_lakes" -> new WaterInfo(0, 1, FieldValues.CLASS_LAKE); - case "ne_50m_lakes" -> new WaterInfo(2, 3, FieldValues.CLASS_LAKE); - case "ne_10m_lakes" -> new WaterInfo(4, 5, FieldValues.CLASS_LAKE); default -> null; }; if (info != null) { - features.polygon(LAYER_NAME) - .setBufferPixels(BUFFER_SIZE) - .setZoomRange(info.minZoom, info.maxZoom) - .setAttr(Fields.CLASS, info.clazz); + setupNeWaterFeature(features, info.minZoom, info.maxZoom, info.clazz, null); + return; + } + + LakeInfo lakeInfo = switch (table) { + case "ne_110m_lakes" -> new LakeInfo(0, 1, FieldValues.CLASS_LAKE); + case "ne_50m_lakes" -> new LakeInfo(2, 3, FieldValues.CLASS_LAKE); + case "ne_10m_lakes" -> new LakeInfo(4, 5, FieldValues.CLASS_LAKE); + default -> null; + }; + if (lakeInfo != null) { + try { + var geom = feature.worldGeometry(); + lakeInfo.geom = geom; + lakeInfo.name = feature.getString("name"); + + var index = neLakeIndexes.computeIfAbsent(table, t -> new STRtree()); + var neLakeNameMap = neLakeNameMaps.computeIfAbsent(table, t -> new ConcurrentHashMap<>()); + + // need to externally synchronize inserts into the STRTree and ArrayList + synchronized (this) { + index.insert(geom.getEnvelopeInternal(), lakeInfo); + neAllLakeInfos.add(lakeInfo); + } + if (lakeInfo.name != null) { + if (!neLakeNameMap.containsKey(lakeInfo.name) || + lakeInfo.geom.getArea() > neLakeNameMap.get(lakeInfo.name).geom.getArea()) { + // on name collision, bigger lake gets on the name list + neLakeNameMap.put(lakeInfo.name, lakeInfo); + } + } + } catch (GeometryException e) { + e.log(stats, "omt_water_ne", + "Error getting geometry for natural earth feature " + table + " " + feature.getTag("ogc_fid")); + // make sure we have this NE lake even if without OSM ID + setupNeWaterFeature(features, lakeInfo.minZoom, lakeInfo.maxZoom, lakeInfo.clazz, null); + } } } + private void setupNeWaterFeature(FeatureCollector features, int minZoom, int maxZoom, String clazz, Long osmId) { + features.polygon(LAYER_NAME) + .setBufferPixels(BUFFER_SIZE) + .setZoomRange(minZoom, maxZoom) + .setAttr(Fields.CLASS, clazz) + .setAttr(Fields.ID, osmId); + } + @Override public void processOsmWater(SourceFeature feature, FeatureCollector features) { features.polygon(LAYER_NAME) @@ -121,6 +176,108 @@ public void process(Tables.OsmWaterPolygon element, FeatureCollector features) { .setAttr(Fields.INTERMITTENT, element.isIntermittent() ? 1 : 0) .setAttrWithMinzoom(Fields.BRUNNEL, Utils.brunnel(element.isBridge(), element.isTunnel()), 12) .setAttr(Fields.CLASS, clazz); + + fillOsmIdIntoNeLake(element); + } + } + + void fillOsmIdIntoNeLake(Tables.OsmWaterPolygon element) { + try { + // match by name: + boolean match = false; + if (element.name() != null) { + for (var map : neLakeNameMaps.values()) { + var lakeInfo = map.get(element.name()); + if (lakeInfo != null) { + match = true; + fillOsmIdIntoNeLake(element, element.source().worldGeometry(), lakeInfo); + } + } + } + if (match) { + return; + } + + // if OSM lake is too small for Z6 (e.g. area bellow ~4px) we assume there is no matching NE lake + Geometry geom = element.source().worldGeometry(); + if (geom.getArea() < SMALLEST_OSM_LAKE_AREA) { + return; + } + + // match by intersection: + Envelope envelope = geom.getEnvelopeInternal(); + + for (var index : neLakeIndexes.values()) { + var items = index.query(envelope); + if (items.size() == 1) { + if (items.getFirst()instanceof LakeInfo lakeInfo) { + fillOsmIdIntoNeLake(element, geom, lakeInfo); + } + } else if (!items.isEmpty()) { + for (var item : items) { + if (item instanceof LakeInfo lakeInfo) { + fillOsmIdIntoNeLake(element, geom, lakeInfo); + } + } + } + } + } catch (GeometryException e) { + e.log(stats, "omt_water", + "Error getting geometry for OSM feature " + element.source().id()); + } + } + + void fillOsmIdIntoNeLake(Tables.OsmWaterPolygon element, Geometry geom, LakeInfo lakeInfo) throws GeometryException { + Geometry neGeom = lakeInfo.geom; + Geometry intersection; + try { + if (!neGeom.intersects(geom)) { + return; + } + intersection = neGeom.intersection(geom); + } catch (TopologyException e) { + try { + Geometry fixedGeom = GeometryFixer.fix(geom); + if (!neGeom.intersects(fixedGeom)) { + return; + } + intersection = neGeom.intersection(fixedGeom); + } catch (TopologyException e2) { + throw new GeometryException("fix_omt_water_topology_error", + "error fixing polygon: " + e2 + "; original error: " + e); + } + } + + // should match following in OpenMapTiles: Distinct on keeps just the first occurence -> order by 'area_ratio DESC' + double areaRatio = intersection.getArea() / neGeom.getArea(); + if (areaRatio > lakeInfo.areaRatio) { + lakeInfo.osmId = element.source().id(); + lakeInfo.areaRatio = areaRatio; + } + } + + @Override + public void finish(String sourceName, FeatureCollector.Factory featureCollectors, + Consumer emit) { + if (OpenMapTilesProfile.NATURAL_EARTH_SOURCE.equals(sourceName)) { + var timer = stats.startStage("ne_lake_index"); + for (var index : neLakeIndexes.values()) { + index.build(); + } + timer.stop(); + } else if (OpenMapTilesProfile.OSM_SOURCE.equals(sourceName)) { + var timer = stats.startStage("ne_lakes"); + for (var item : neAllLakeInfos) { + var features = featureCollectors.get(SimpleFeature.fromWorldGeometry(item.geom)); + setupNeWaterFeature(features, item.minZoom, item.maxZoom, item.clazz, item.osmId); + for (var feature : features) { + emit.accept(feature); + } + } + neLakeNameMaps.clear(); + neLakeIndexes.clear(); + neAllLakeInfos.clear(); + timer.stop(); } } @@ -128,4 +285,26 @@ public void process(Tables.OsmWaterPolygon element, FeatureCollector features) { public List postProcess(int zoom, List items) throws GeometryException { return items.size() > 1 ? FeatureMerge.mergeOverlappingPolygons(items, config.minFeatureSize(zoom)) : items; } + + /** + * Information to hold onto from processing an NE lake to determine OSM ID later. + */ + private static class LakeInfo { + String name; + int minZoom; + int maxZoom; + String clazz; + Geometry geom; + Long osmId; + double areaRatio; + + public LakeInfo(int minZoom, int maxZoom, String clazz) { + this.name = null; + this.minZoom = minZoom; + this.maxZoom = maxZoom; + this.clazz = clazz; + this.osmId = null; + this.areaRatio = 0; + } + } } diff --git a/src/test/java/org/openmaptiles/layers/WaterTest.java b/src/test/java/org/openmaptiles/layers/WaterTest.java index ee900dc5..7dc852d3 100644 --- a/src/test/java/org/openmaptiles/layers/WaterTest.java +++ b/src/test/java/org/openmaptiles/layers/WaterTest.java @@ -2,8 +2,10 @@ import static com.onthegomap.planetiler.TestUtils.rectangle; +import com.onthegomap.planetiler.FeatureCollector; import com.onthegomap.planetiler.geo.GeoUtils; import com.onthegomap.planetiler.reader.SimpleFeature; +import java.util.ArrayList; import java.util.HashMap; import java.util.List; import java.util.Map; @@ -16,7 +18,7 @@ class WaterTest extends AbstractLayerTest { @Test void testWaterNaturalEarth() { assertFeatures(0, List.of(Map.of( - "class", "lake", + "class", "ocean", "intermittent", "", "_layer", "water", "_type", "polygon", @@ -25,49 +27,108 @@ void testWaterNaturalEarth() { rectangle(0, 10), Map.of(), OpenMapTilesProfile.NATURAL_EARTH_SOURCE, - "ne_110m_lakes", + "ne_110m_ocean", 0 ))); - assertFeatures(0, List.of(Map.of( + assertFeatures(6, List.of(Map.of( "class", "ocean", - "intermittent", "", "_layer", "water", "_type", "polygon", - "_minzoom", 0 + "_maxzoom", 5 )), process(SimpleFeature.create( rectangle(0, 10), Map.of(), OpenMapTilesProfile.NATURAL_EARTH_SOURCE, - "ne_110m_ocean", + "ne_10m_ocean", 0 ))); + } - assertFeatures(6, List.of(Map.of( - "class", "lake", - "_layer", "water", - "_type", "polygon", - "_maxzoom", 5 - )), process(SimpleFeature.create( - rectangle(0, 10), + @Test + void testLakeNaturalEarthByIntersection() { + final var polygon = rectangle(0, 0.1); + // NE lakes: + process(SimpleFeature.create( + polygon, + Map.of(), + OpenMapTilesProfile.NATURAL_EARTH_SOURCE, + "ne_110m_lakes", + 0 + )); + process(SimpleFeature.create( + polygon, Map.of(), OpenMapTilesProfile.NATURAL_EARTH_SOURCE, "ne_10m_lakes", 0 - ))); + )); + // OSM lake to take the ID from: + process(SimpleFeature.create( + polygon, + new HashMap<>(Map.of( + "natural", "water", + "water", "reservoir" + )), + OpenMapTilesProfile.OSM_SOURCE, + null, + 123 + )); - assertFeatures(6, List.of(Map.of( - "class", "ocean", + List features = new ArrayList<>(); + profile.finish(OpenMapTilesProfile.OSM_SOURCE, new FeatureCollector.Factory(params, stats), features::add); + assertFeatures(0, List.of(Map.of( + "class", "lake", + "intermittent", "", + "id", 123L, "_layer", "water", "_type", "polygon", + "_minzoom", 0, + "_maxzoom", 1 + ), Map.of( + "class", "lake", + "intermittent", "", + "id", 123L, + "_layer", "water", + "_type", "polygon", + "_minzoom", 4, "_maxzoom", 5 - )), process(SimpleFeature.create( - rectangle(0, 10), - Map.of(), + )), features); + } + + @Test + void testLakeNaturalEarthByName() { + final var polygon = rectangle(0, 0.1); + // NE lake: + process(SimpleFeature.create( + polygon, + Map.of("name", "Test Lake"), OpenMapTilesProfile.NATURAL_EARTH_SOURCE, - "ne_10m_ocean", + "ne_50m_lakes", 0 - ))); + )); + // OSM lake to take the ID from: + process(SimpleFeature.create( + polygon, + new HashMap<>(Map.of( + "name", "Test Lake", + "natural", "water", + "water", "reservoir" + )), + OpenMapTilesProfile.OSM_SOURCE, + null, + 123 + )); + + List features = new ArrayList<>(); + profile.finish(OpenMapTilesProfile.OSM_SOURCE, new FeatureCollector.Factory(params, stats), features::add); + assertFeatures(0, List.of(Map.of( + "class", "lake", + "id", 123L, + "_layer", "water", + "_minzoom", 2, + "_maxzoom", 3 + )), features); } @Test From 57828cc99e5aae4ec394c327acd633a2e7dc4786 Mon Sep 17 00:00:00 2001 From: Peter Hanecak Date: Sat, 23 Mar 2024 19:18:59 +0100 Subject: [PATCH 2/5] testLakeZoomLevels(): NE lakes removed, checked now in testLakeNaturalEarth*() --- .../org/openmaptiles/layers/WaterTest.java | 23 +------------------ 1 file changed, 1 insertion(+), 22 deletions(-) diff --git a/src/test/java/org/openmaptiles/layers/WaterTest.java b/src/test/java/org/openmaptiles/layers/WaterTest.java index 7dc852d3..3264e55e 100644 --- a/src/test/java/org/openmaptiles/layers/WaterTest.java +++ b/src/test/java/org/openmaptiles/layers/WaterTest.java @@ -316,28 +316,7 @@ void testOceanZoomLevels() { @Test void testLakeZoomLevels() { - assertCoversZoomRange(0, 14, "water", - process(SimpleFeature.create( - rectangle(0, 10), - Map.of(), - OpenMapTilesProfile.NATURAL_EARTH_SOURCE, - "ne_110m_lakes", - 0 - )), - process(SimpleFeature.create( - rectangle(0, 10), - Map.of(), - OpenMapTilesProfile.NATURAL_EARTH_SOURCE, - "ne_50m_lakes", - 0 - )), - process(SimpleFeature.create( - rectangle(0, 10), - Map.of(), - OpenMapTilesProfile.NATURAL_EARTH_SOURCE, - "ne_10m_lakes", - 0 - )), + assertCoversZoomRange(6, 14, "water", process(SimpleFeature.create( rectangle(0, 10), Map.of( From 58c86c66b9ad7bfda7a04023bf9b8c9367c8850b Mon Sep 17 00:00:00 2001 From: Peter Hanecak Date: Wed, 27 Mar 2024 21:53:54 +0100 Subject: [PATCH 3/5] OSM lake area limit increased 2x to match smallest NE lake --- src/main/java/org/openmaptiles/layers/Water.java | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/src/main/java/org/openmaptiles/layers/Water.java b/src/main/java/org/openmaptiles/layers/Water.java index 6655a9df..adb4ec5c 100644 --- a/src/main/java/org/openmaptiles/layers/Water.java +++ b/src/main/java/org/openmaptiles/layers/Water.java @@ -82,8 +82,8 @@ public class Water implements * which infers ocean polygons by preprocessing all coastline elements. */ - // smallest NE lake is around 4.42E-13, this is roughly 1/2 of that and approx 17% of OSM features are smaller than this: - private static final double SMALLEST_OSM_LAKE_AREA = Math.pow(4, -21); + // smallest NE lake is around 4.42E-13, this is roughly that and approx. 23% of OSM features are smaller than this: + private static final double SMALLEST_OSM_LAKE_AREA = Math.pow(4, -20.5); private final MultiExpression.Index classMapping; private final PlanetilerConfig config; @@ -210,7 +210,7 @@ void fillOsmIdIntoNeLake(Tables.OsmWaterPolygon element) { for (var index : neLakeIndexes.values()) { var items = index.query(envelope); if (items.size() == 1) { - if (items.getFirst()instanceof LakeInfo lakeInfo) { + if (items.getFirst() instanceof LakeInfo lakeInfo) { fillOsmIdIntoNeLake(element, geom, lakeInfo); } } else if (!items.isEmpty()) { From acbaa22dbac1dabbd50f09e411e9d74f2200fd28 Mon Sep 17 00:00:00 2001 From: Peter Hanecak Date: Wed, 27 Mar 2024 21:56:16 +0100 Subject: [PATCH 4/5] OSM lake area limit increased 2x even further to match smallest observed during Planet processing --- src/main/java/org/openmaptiles/layers/Water.java | 5 +++-- 1 file changed, 3 insertions(+), 2 deletions(-) diff --git a/src/main/java/org/openmaptiles/layers/Water.java b/src/main/java/org/openmaptiles/layers/Water.java index adb4ec5c..56fa91db 100644 --- a/src/main/java/org/openmaptiles/layers/Water.java +++ b/src/main/java/org/openmaptiles/layers/Water.java @@ -82,8 +82,9 @@ public class Water implements * which infers ocean polygons by preprocessing all coastline elements. */ - // smallest NE lake is around 4.42E-13, this is roughly that and approx. 23% of OSM features are smaller than this: - private static final double SMALLEST_OSM_LAKE_AREA = Math.pow(4, -20.5); + // smallest NE lake is around 4.42E-13, smallest matching OSM lake is 9.34E-13, this is slightly bellow that + // and approx. 33% of OSM features are smaller than this, hence to save some CPU cycles: + private static final double SMALLEST_OSM_LAKE_AREA = Math.pow(4, -20); private final MultiExpression.Index classMapping; private final PlanetilerConfig config; From af13750adfa3e91ef3f8b658d7ab7e74d61552ad Mon Sep 17 00:00:00 2001 From: Peter Hanecak Date: Wed, 27 Mar 2024 21:58:13 +0100 Subject: [PATCH 5/5] clean-up: SMALLEST_OSM_LAKE_AREA renamed to OSM_ID_MATCH_AREA_LIMIT to better match what it is --- src/main/java/org/openmaptiles/layers/Water.java | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/src/main/java/org/openmaptiles/layers/Water.java b/src/main/java/org/openmaptiles/layers/Water.java index 56fa91db..c0262415 100644 --- a/src/main/java/org/openmaptiles/layers/Water.java +++ b/src/main/java/org/openmaptiles/layers/Water.java @@ -84,7 +84,7 @@ public class Water implements // smallest NE lake is around 4.42E-13, smallest matching OSM lake is 9.34E-13, this is slightly bellow that // and approx. 33% of OSM features are smaller than this, hence to save some CPU cycles: - private static final double SMALLEST_OSM_LAKE_AREA = Math.pow(4, -20); + private static final double OSM_ID_MATCH_AREA_LIMIT = Math.pow(4, -20); private final MultiExpression.Index classMapping; private final PlanetilerConfig config; @@ -201,7 +201,7 @@ void fillOsmIdIntoNeLake(Tables.OsmWaterPolygon element) { // if OSM lake is too small for Z6 (e.g. area bellow ~4px) we assume there is no matching NE lake Geometry geom = element.source().worldGeometry(); - if (geom.getArea() < SMALLEST_OSM_LAKE_AREA) { + if (geom.getArea() < OSM_ID_MATCH_AREA_LIMIT) { return; }