Fix invalid polygons before snapping to tile coordinates (#566)

pull/571/head
Michael Barry 2023-04-27 08:19:31 -04:00 zatwierdzone przez GitHub
rodzic e0285d2e1c
commit 8d0e06c667
Nie znaleziono w bazie danych klucza dla tego podpisu
ID klucza GPG: 4AEE18F83AFDEB23
12 zmienionych plików z 244 dodań i 19 usunięć

Wyświetl plik

@ -9,6 +9,8 @@ import com.onthegomap.planetiler.geo.GeoUtils;
import com.onthegomap.planetiler.geo.GeometryException;
import com.onthegomap.planetiler.geo.GeometryType;
import com.onthegomap.planetiler.geo.MutableCoordinateSequence;
import com.onthegomap.planetiler.stats.DefaultStats;
import com.onthegomap.planetiler.stats.Stats;
import java.util.ArrayList;
import java.util.BitSet;
import java.util.Collection;
@ -233,12 +235,13 @@ public class FeatureMerge {
* @param minDist the minimum threshold in tile pixels between polygons to combine into a group
* @param buffer the amount (in tile pixels) to expand then contract polygons by in order to combine
* almost-touching polygons
* @param stats for counting data errors
* @return a new list containing all unaltered features in their original order, then each of the merged groups
* ordered by the index of the first element in that group from the input list.
* @throws GeometryException if an error occurs encoding the combined geometry
*/
public static List<VectorTile.Feature> mergeNearbyPolygons(List<VectorTile.Feature> features, double minArea,
double minHoleArea, double minDist, double buffer) throws GeometryException {
double minHoleArea, double minDist, double buffer, Stats stats) throws GeometryException {
List<VectorTile.Feature> result = new ArrayList<>(features.size());
Collection<List<VectorTile.Feature>> groupedByAttrs = groupByAttrs(features, result, GeometryType.POLYGON);
for (List<VectorTile.Feature> groupedFeatures : groupedByAttrs) {
@ -272,7 +275,7 @@ public class FeatureMerge {
if (!(merged instanceof Polygonal) || merged.getEnvelopeInternal().getArea() < minArea) {
continue;
}
merged = GeoUtils.snapAndFixPolygon(merged).reverse();
merged = GeoUtils.snapAndFixPolygon(merged, stats, "merge").reverse();
} else {
merged = polygonGroup.get(0);
if (!(merged instanceof Polygonal) || merged.getEnvelopeInternal().getArea() < minArea) {
@ -289,6 +292,11 @@ public class FeatureMerge {
return result;
}
public static List<VectorTile.Feature> mergeNearbyPolygons(List<VectorTile.Feature> features, double minArea,
double minHoleArea, double minDist, double buffer) throws GeometryException {
return mergeNearbyPolygons(features, minArea, minHoleArea, minDist, buffer, DefaultStats.get());
}
/**
* Returns all the clusters from {@code geometries} where elements in the group are less than {@code minDist} from

Wyświetl plik

@ -1,6 +1,7 @@
package com.onthegomap.planetiler.geo;
import com.onthegomap.planetiler.collection.LongLongMap;
import com.onthegomap.planetiler.stats.Stats;
import java.util.ArrayList;
import java.util.List;
import java.util.stream.Stream;
@ -22,6 +23,7 @@ import org.locationtech.jts.geom.PrecisionModel;
import org.locationtech.jts.geom.TopologyException;
import org.locationtech.jts.geom.impl.PackedCoordinateSequence;
import org.locationtech.jts.geom.impl.PackedCoordinateSequenceFactory;
import org.locationtech.jts.geom.util.GeometryFixer;
import org.locationtech.jts.geom.util.GeometryTransformer;
import org.locationtech.jts.io.WKBReader;
import org.locationtech.jts.precision.GeometryPrecisionReducer;
@ -294,8 +296,8 @@ public class GeoUtils {
* Returns a copy of {@code geom} with coordinates rounded to {@link #TILE_PRECISION} and fixes any polygon
* self-intersections or overlaps that may have caused.
*/
public static Geometry snapAndFixPolygon(Geometry geom) throws GeometryException {
return snapAndFixPolygon(geom, TILE_PRECISION);
public static Geometry snapAndFixPolygon(Geometry geom, Stats stats, String stage) throws GeometryException {
return snapAndFixPolygon(geom, TILE_PRECISION, stats, stage);
}
/**
@ -304,21 +306,29 @@ public class GeoUtils {
*
* @throws GeometryException if an unrecoverable robustness exception prevents us from fixing the geometry
*/
public static Geometry snapAndFixPolygon(Geometry geom, PrecisionModel tilePrecision) throws GeometryException {
public static Geometry snapAndFixPolygon(Geometry geom, PrecisionModel tilePrecision, Stats stats, String stage)
throws GeometryException {
try {
if (!geom.isValid()) {
geom = fixPolygon(geom);
stats.dataError(stage + "_snap_fix_input");
}
return GeometryPrecisionReducer.reduce(geom, tilePrecision);
} catch (IllegalArgumentException e) {
} catch (TopologyException | IllegalArgumentException e) {
// precision reduction fails if geometry is invalid, so attempt
// to fix it then try again
geom = fixPolygon(geom);
geom = GeometryFixer.fix(geom);
stats.dataError(stage + "_snap_fix_input2");
try {
return GeometryPrecisionReducer.reduce(geom, tilePrecision);
} catch (IllegalArgumentException e2) {
} catch (TopologyException | IllegalArgumentException e2) {
// give it one last try but with more aggressive fixing, just in case (see issue #511)
geom = fixPolygon(geom, tilePrecision.gridSize() / 2);
stats.dataError(stage + "_snap_fix_input3");
try {
return GeometryPrecisionReducer.reduce(geom, tilePrecision);
} catch (IllegalArgumentException e3) {
} catch (TopologyException | IllegalArgumentException e3) {
stats.dataError(stage + "_snap_fix_input3_failed");
throw new GeometryException("snap_third_time_failed", "Error reducing precision");
}
}

Wyświetl plik

@ -793,10 +793,11 @@ public class OsmReader implements Closeable, MemoryEstimator.HasEstimate {
if (member.type() == OsmElement.Type.WAY) {
if (poly != null && !poly.isEmpty()) {
rings.add(poly);
} else if (relation.hasTag("type", "multipolygon")) {
} else {
// boundary and land_area relations might not be complete for extracts, but multipolygons should be
LOGGER.warn(
"Missing " + role + " OsmWay[" + member.ref() + "] for " + relation.getTag("type") + " " + this);
stats.dataError("osm_" + relation.getTag("type") + "_missing_way");
LOGGER.trace(
"Missing {} OsmWay[{}] for {} {}", role, member.ref(), relation.getTag("type"), this);
}
}
}

Wyświetl plik

@ -237,7 +237,7 @@ public class FeatureRenderer implements Consumer<FeatureCollector.Feature>, Clos
* See https://docs.mapbox.com/vector-tiles/specification/#simplification for issues that can arise from naive
* coordinate rounding.
*/
geom = GeoUtils.snapAndFixPolygon(geom);
geom = GeoUtils.snapAndFixPolygon(geom, stats, "render");
// JTS utilities "fix" the geometry to be clockwise outer/CCW inner but vector tiles flip Y coordinate,
// so we need outer CCW/inner clockwise
geom = geom.reverse();

Wyświetl plik

@ -0,0 +1,18 @@
package com.onthegomap.planetiler.stats;
/**
* Holder for default {@link Stats} implementation to use for this process.
*/
public class DefaultStats {
private DefaultStats() {}
private static Stats defaultValue = null;
public static Stats get() {
return defaultValue;
}
public static void set(Stats stats) {
defaultValue = stats;
}
}

Wyświetl plik

@ -23,6 +23,7 @@ import java.time.Duration;
import java.util.ArrayList;
import java.util.List;
import java.util.Map;
import java.util.concurrent.ConcurrentHashMap;
import java.util.concurrent.ConcurrentSkipListMap;
import java.util.concurrent.Executors;
import java.util.concurrent.ScheduledExecutorService;
@ -49,6 +50,7 @@ class PrometheusStats implements Stats {
private ScheduledExecutorService executor;
private final String job;
private final Map<String, Path> filesToMonitor = new ConcurrentSkipListMap<>();
private final Map<String, Long> dataErrorCounters = new ConcurrentHashMap<>();
private final Map<String, MemoryEstimator.HasEstimate> heapObjectsToMonitor = new ConcurrentSkipListMap<>();
/** Constructs a new instance but does not start polling (for tests). */
@ -128,6 +130,7 @@ class PrometheusStats implements Stats {
@Override
public void dataError(String errorCode) {
Stats.super.dataError(errorCode);
dataErrors.labels(errorCode).inc();
}
@ -203,6 +206,11 @@ class PrometheusStats implements Stats {
}.register(registry);
}
@Override
public Map<String, Long> dataErrors() {
return dataErrorCounters;
}
@Override
public void close() {
executor.shutdown();

Wyświetl plik

@ -9,6 +9,7 @@ import com.onthegomap.planetiler.util.MemoryEstimator;
import java.nio.file.Path;
import java.time.Duration;
import java.util.Map;
import java.util.concurrent.ConcurrentHashMap;
import java.util.concurrent.ConcurrentSkipListMap;
import java.util.function.LongSupplier;
import java.util.function.Supplier;
@ -26,7 +27,9 @@ public interface Stats extends AutoCloseable {
/** Returns a new stat collector that stores basic stats in-memory to report through {@link #printSummary()}. */
static Stats inMemory() {
return new InMemory();
var stats = new InMemory();
DefaultStats.set(stats);
return stats;
}
/**
@ -34,7 +37,9 @@ public interface Stats extends AutoCloseable {
* <a href="https://github.com/prometheus/pushgateway">prometheus push gateway</a> at {@code destination}.
*/
static Stats prometheusPushGateway(String destination, String job, Duration interval) {
return PrometheusStats.createAndStartPushing(destination, job, interval);
var stats = PrometheusStats.createAndStartPushing(destination, job, interval);
DefaultStats.set(stats);
return stats;
}
/**
@ -47,6 +52,11 @@ public interface Stats extends AutoCloseable {
if (logger.isInfoEnabled()) {
logger.info("");
logger.info("-".repeat(40));
logger.info("data errors:");
dataErrors().entrySet().stream()
.sorted(Map.Entry.<String, Long>comparingByValue().reversed())
.forEachOrdered(entry -> logger.info("\t{}\t{}", entry.getKey(), format.integer(entry.getValue())));
logger.info("-".repeat(40));
timers().printSummary();
logger.info("-".repeat(40));
for (var entry : monitoredFiles().entrySet()) {
@ -156,11 +166,16 @@ public interface Stats extends AutoCloseable {
*/
void counter(String name, String label, Supplier<Map<String, LongSupplier>> values);
/** Returns all the data error counters. */
Map<String, Long> dataErrors();
/**
* Records that an invalid input feature was discarded where {@code errorCode} can be used to identify the kind of
* failure.
*/
void dataError(String errorCode);
default void dataError(String errorCode) {
dataErrors().merge(errorCode, 1L, Long::sum);
}
/**
* A stat collector that stores top-level metrics in-memory to report through {@link #printSummary()}.
@ -174,6 +189,7 @@ public interface Stats extends AutoCloseable {
private final Timers timers = new Timers();
private final Map<String, Path> monitoredFiles = new ConcurrentSkipListMap<>();
private final Map<String, Long> dataErrors = new ConcurrentHashMap<>();
@Override
public void wroteTile(int zoom, int bytes) {}
@ -208,10 +224,12 @@ public interface Stats extends AutoCloseable {
public void counter(String name, String label, Supplier<Map<String, LongSupplier>> values) {}
@Override
public void processedElement(String elemType, String layer) {}
public Map<String, Long> dataErrors() {
return dataErrors;
}
@Override
public void dataError(String errorCode) {}
public void processedElement(String elemType, String layer) {}
@Override
public void gauge(String name, Supplier<Number> value) {}

Wyświetl plik

@ -1710,6 +1710,54 @@ class PlanetilerTests {
), counts);
}
@Test
void testIssue546Terschelling() throws Exception {
Geometry geometry = new WKBReader()
.read(
new InputStreamInStream(Files.newInputStream(TestUtils.pathToResource("issue_546_terschelling.wkb"))));
geometry = GeoUtils.worldToLatLonCoords(geometry);
assertNotNull(geometry);
// automatically checks for self-intersections
var results = runWithReaderFeatures(
Map.of("threads", "1"),
List.of(
newReaderFeature(geometry, Map.of())
),
(in, features) -> features.polygon("ocean")
.setBufferPixels(4.0)
.setMinZoom(0)
);
// this lat/lon is in the middle of an island and should not be covered by
for (int z = 4; z <= 14; z++) {
double lat = 53.391958;
double lon = 5.2438441;
var coord = TileCoord.aroundLngLat(lon, lat, z);
var problematicTile = results.tiles.get(coord);
if (z == 14) {
assertNull(problematicTile);
continue;
}
double scale = Math.pow(2, coord.z());
double tileX = (GeoUtils.getWorldX(lon) * scale - coord.x()) * 256;
double tileY = (GeoUtils.getWorldY(lat) * scale - coord.y()) * 256;
var point = newPoint(tileX, tileY);
assertEquals(1, problematicTile.size());
var geomCompare = problematicTile.get(0).geometry();
geomCompare.validate();
var geom = geomCompare.geom();
assertFalse(geom.covers(point), "z" + z);
}
}
@ParameterizedTest
@ValueSource(strings = {
"",

Wyświetl plik

@ -3,12 +3,16 @@ package com.onthegomap.planetiler.geo;
import static com.onthegomap.planetiler.TestUtils.*;
import static com.onthegomap.planetiler.geo.GeoUtils.*;
import static org.junit.jupiter.api.Assertions.assertEquals;
import static org.junit.jupiter.api.Assertions.assertFalse;
import static org.junit.jupiter.api.Assertions.assertTrue;
import com.onthegomap.planetiler.stats.Stats;
import java.util.List;
import org.geotools.geometry.jts.WKTReader2;
import org.junit.jupiter.api.Test;
import org.junit.jupiter.params.ParameterizedTest;
import org.junit.jupiter.params.provider.CsvSource;
import org.locationtech.jts.algorithm.Orientation;
import org.locationtech.jts.geom.Geometry;
import org.locationtech.jts.geom.LinearRing;
import org.locationtech.jts.geom.Point;
@ -366,7 +370,81 @@ class GeoUtilsTest {
var result = GeoUtils.snapAndFixPolygon(new WKTReader2().read(
"""
MULTIPOLYGON (((198.83750000000003 46.07500000000004, 199.0625 46.375, 199.4375 46.0625, 199.5 46.43750000000001, 199.5625 46, 199.3125 45.5, 198.8912037037037 46.101851851851876, 198.83750000000003 46.07500000000004)), ((198.43750000000003 46.49999999999999, 198.5625 46.43750000000001, 198.6875 46.25, 198.1875 46.25, 198.43750000000003 46.49999999999999)), ((198.6875 46.25, 198.81249999999997 46.062500000000014, 198.6875 46.00000000000002, 198.6875 46.25)), ((196.55199579831933 46.29359243697479, 196.52255639097743 46.941259398496236, 196.5225563909774 46.941259398496236, 196.49999999999997 47.43750000000001, 196.875 47.125, 197 47.5625, 197.47880544905414 46.97729334004497, 197.51505401161464 46.998359569801956, 197.25 47.6875, 198.0625 47.6875, 198.5 46.625, 198.34375 46.546875, 198.34375000000003 46.54687499999999, 197.875 46.3125, 197.875 46.25, 197.875 46.0625, 197.82894736842107 46.20065789473683, 197.25 46.56250000000001, 197.3125 46.125, 196.9375 46.1875, 196.9375 46.21527777777778, 196.73250000000002 46.26083333333334, 196.5625 46.0625, 196.55199579831933 46.29359243697479)), ((196.35213414634146 45.8170731707317, 197.3402027027027 45.93108108108108, 197.875 45.99278846153846, 197.875 45.93750000000002, 197.93749999999997 45.99999999999999, 197.9375 46, 197.90625 45.96874999999999, 197.90625 45.96875, 196.75000000000006 44.81250000000007, 197.1875 45.4375, 196.3125 45.8125, 196.35213414634146 45.8170731707317)), ((195.875 46.124999999999986, 195.8125 46.5625, 196.5 46.31250000000001, 195.9375 46.4375, 195.875 46.124999999999986)), ((196.49999999999997 46.93749999999999, 196.125 46.875, 196.3125 47.125, 196.49999999999997 46.93749999999999)))
"""));
"""),
Stats.inMemory(), "test");
assertEquals(3.146484375, result.getArea(), 1e-5);
}
@Test
void testSnapAndFixIssue546() throws GeometryException, ParseException {
var orig = new WKTReader2().read(
"""
POLYGON(
(
0 0,
2 0,
2 1,
0 1,
0 0
),
(
1.19053596444428 0.9029693332377,
1.06439281777784 0.92184484124482,
1.42787640888855 0.8257737683016,
0.94579086222257 0.5617504680322,
1.19053596444428 0.9029693332377
)
)
""");
var result = GeoUtils.snapAndFixPolygon(orig, Stats.inMemory(), "test");
var point = newPoint(1.14602002962965, 0.76978969252622);
assertTrue(result.isValid());
assertFalse(result.contains(point));
}
@Test
void testSnapAndFixIssue546_2() throws GeometryException, ParseException {
var orig = new WKTReader2().read(
"""
POLYGON(
(
1.19053596444428 0.9029693332377,
1.06439281777784 0.92184484124482,
1.42787640888855 0.8257737683016,
0.94579086222257 0.5617504680322,
1.19053596444428 0.9029693332377
)
)
""");
var result = GeoUtils.snapAndFixPolygon(orig, Stats.inMemory(), "test");
assertTrue(result.isValid());
assertFalse(Orientation.isCCWArea(result.getCoordinates()));
}
@Test
void testSnapAndFixIssue546_3() throws GeometryException, ParseException {
var orig = new WKTReader2().read(
"""
POLYGON(
(
0 0,
2 0,
2 1,
0 1,
0 0
),
(
1.19053596444428 0.9029693332377,
1.06439281777784 0.92184484124482,
1.42787640888855 0.8257737683016,
0.94579086222257 0.5617504680322,
1.19053596444428 0.9029693332377
)
)
""");
var result = GeoUtils.snapAndFixPolygon(orig, Stats.inMemory(), "test");
var point = newPoint(1.14602002962965, 0.76978969252622);
assertTrue(result.isValid());
assertFalse(result.contains(point));
}
}

Wyświetl plik

@ -112,6 +112,16 @@ class PrometheusStatsTest {
assertContainsStat("^planetiler_heap_object_test_size_bytes 10", stats);
}
@Test
void captureDataErrors() {
PrometheusStats stats = new PrometheusStats("job");
assertEquals(Map.of(), stats.dataErrors());
stats.dataError("a");
stats.dataError("a");
stats.dataError("b");
assertEquals(Map.of("a", 2L, "b", 1L), stats.dataErrors());
}
private static Counter.Readable counterAt(int num) {
var result = Counter.newSingleThreadCounter();
result.incBy(num);

Wyświetl plik

@ -0,0 +1,26 @@
package com.onthegomap.planetiler.stats;
import static org.junit.jupiter.api.Assertions.assertEquals;
import static org.junit.jupiter.api.Assertions.assertSame;
import java.util.Map;
import org.junit.jupiter.api.Test;
class StatsTest {
@Test
void testDefaultStats() {
var a = Stats.inMemory();
var b = DefaultStats.get();
assertSame(a, b);
}
@Test
void captureDataErrors() {
var a = Stats.inMemory();
assertEquals(Map.of(), a.dataErrors());
a.dataError("a");
a.dataError("a");
a.dataError("b");
assertEquals(Map.of("a", 2L, "b", 1L), a.dataErrors());
}
}

Plik binarny nie jest wyświetlany.