Add an option to buffer polygons outward to keep them from collapsing

CHANGELOG.md

@@ -1,3 +1,7 @@
+# 2.36.0
+
+* Add --buffer-polygons-outward option to prevent collapse of polygon spindles
+
 # 2.35.0
 
 * Fix a bug in --detect-longitude-wraparound when there are multiple rings

README.md

@@ -496,6 +496,7 @@ the same layer, enclose them in an `all` expression so they will all be evaluate
  * `-pt` or `--no-tiny-polygon-reduction`: Don't combine the area of very small polygons into small squares that represent their combined area.
  * `-pT` or `--no-tiny-polygon-reduction-at-maximum-zoom`: Combine the area of very small polygons into small squares that represent their combined area only at zoom levels below the maximum.
  * `--tiny-polygon-size=`_size_: Use the specified _size_ for tiny polygons instead of the default 2. Anything above 6 or so will lead to visible artifacts with the default tile detail.
+ * `-aO` or `--buffer-polygons-outward`: Buffer polygons, except those that have been detected to be continous by `--no-simplification-of-shared-nodes`, outward slightly to prevent narrow polygon spindles fom collapsing away.
  * `-av` or `--visvalingam`: Use Visvalingam's simplification algorithm rather than Douglas-Peucker's.
 
 ### Attempts to improve shared polygon boundaries

geometry.cpp

@@ -197,6 +197,27 @@ drawvec reduce_tiny_poly(drawvec &geom, int z, int detail, bool *still_needs_sim
 
 			double area = get_area(geom, i, j);
 
+			// if we are trying to salvage polygons that would otherwise drop out,
+			// also raise the standards for what can qualify as polygon dust
+			if (additional[A_BUFFER_POLYGONS_OUTWARD]) {
+				long long minx = LLONG_MAX;
+				long long maxx = LLONG_MIN;
+				long long miny = LLONG_MAX;
+				long long maxy = LLONG_MIN;
+				for (auto const &d : geom) {
+					minx = std::min(minx, (long long) d.x);
+					maxx = std::max(maxx, (long long) d.x);
+					miny = std::min(miny, (long long) d.y);
+					maxy = std::max(maxy, (long long) d.y);
+				}
+				if (area > 0 && area <= pixel * pixel && area < (maxx - minx) * (maxy - miny) / 3) {
+					// if the polygon doesn't use most of its area,
+					// don't let it be dust, because the shape is
+					// probably something weird and interesting.
+					area = pixel * pixel * 2;
+				}
+			}
+
 			// XXX There is an ambiguity here: If the area of a ring is 0 and it is followed by holes,
 			// we don't know whether the area-0 ring was a hole too or whether it was the outer ring
 			// that these subsequent holes are somehow being subtracted from. I hope that if a polygon
@@ -485,6 +506,26 @@ drawvec impose_tile_boundaries(drawvec &geom, long long extent) {
 	return out;
 }
 
+bool is_shared_node(draw d, int z, int tx, int ty, struct node *shared_nodes_map, size_t nodepos) {
+	if (nodepos > 0) {
+		// offset to global
+		if (z != 0) {
+			d.x += tx * (1LL << (32 - z));
+			d.y += ty * (1LL << (32 - z));
+		}
+
+		// to quadkey
+		struct node n;
+		n.index = encode_quadkey((unsigned) d.x, (unsigned) d.y);
+
+		if (bsearch(&n, shared_nodes_map, nodepos / sizeof(node), sizeof(node), nodecmp) != NULL) {
+			return true;
+		}
+	}
+
+	return false;
+}
+
 drawvec simplify_lines(drawvec &geom, int z, int tx, int ty, int detail, bool mark_tile_bounds, double simplification, size_t retain, drawvec const &shared_nodes, struct node *shared_nodes_map, size_t nodepos) {
 	int res = 1 << (32 - detail - z);
 	long long area = 1LL << (32 - z);
@@ -514,21 +555,8 @@ drawvec simplify_lines(drawvec &geom, int z, int tx, int ty, int detail, bool ma
 				geom[i].necessary = true;
 			}
 
-			if (nodepos > 0) {
-				// offset to global
-				draw d = geom[i];
-				if (z != 0) {
-					d.x += tx * (1LL << (32 - z));
-					d.y += ty * (1LL << (32 - z));
-				}
-
-				// to quadkey
-				struct node n;
-				n.index = encode_quadkey((unsigned) d.x, (unsigned) d.y);
-
-				if (bsearch(&n, shared_nodes_map, nodepos / sizeof(node), sizeof(node), nodecmp) != NULL) {
-					geom[i].necessary = true;
-				}
+			if (is_shared_node(geom[i], z, tx, ty, shared_nodes_map, nodepos)) {
+				geom[i].necessary = true;
 			}
 		}
 	}
@@ -1398,3 +1426,83 @@ drawvec checkerboard_anchors(drawvec const &geom, int tx, int ty, int z, unsigne
 
 	return out;
 }
+
+bool is_continuous_ring(drawvec const &geom, size_t i, size_t j, int z, int tx, int ty, struct node *shared_nodes_map, size_t nodepos) {
+	size_t count = 0;
+
+	// j - 1 to avoid double-counting the duplicate last node
+	for (; i < j - 1; i++) {
+		if (is_shared_node(geom[i], z, tx, ty, shared_nodes_map, nodepos)) {
+			count++;
+
+			// every ring will have three shared nodes; the fourth indicates a shared vertex
+			if (count >= 4) {
+				return true;
+			}
+		}
+	}
+
+	return false;
+}
+
+bool is_continuous_poly(drawvec const &geom, int z, int tx, int ty, struct node *shared_nodes_map, size_t nodepos) {
+	for (size_t i = 0; i < geom.size(); i++) {
+		if (geom[i].op == VT_MOVETO) {
+			size_t j;
+			for (j = i + 1; j < geom.size(); j++) {
+				if (geom[j].op != VT_LINETO) {
+					break;
+				}
+			}
+
+			if (is_continuous_ring(geom, i, j, z, tx, ty, shared_nodes_map, nodepos)) {
+				return true;
+			}
+		}
+	}
+
+	return false;
+}
+
+drawvec buffer_poly(drawvec const &geom, double buffer) {
+	drawvec out = geom;
+
+	if (buffer != 0) {
+		for (size_t i = 0; i < geom.size(); i++) {
+			if (geom[i].op == VT_MOVETO) {
+				size_t j;
+				for (j = i + 1; j < geom.size(); j++) {
+					if (geom[j].op != VT_LINETO) {
+						break;
+					}
+				}
+
+				for (size_t k = i; k < j - 1; k++) {
+					draw p0 = geom[(k + 0 - i) % (j - i - 1) + i];
+					draw p1 = geom[(k + 1 - i) % (j - i - 1) + i];
+					draw p2 = geom[(k + 2 - i) % (j - i - 1) + i];
+
+					double a10 = atan2(p1.y - p0.y, p1.x - p0.x);
+					double a21 = atan2(p2.y - p1.y, p2.x - p1.x);
+
+					double dx = cos(a10 - 90 * M_PI / 180) + cos(a21 - 90 * M_PI / 180);
+					double dy = sin(a10 - 90 * M_PI / 180) + sin(a21 - 90 * M_PI / 180);
+
+					// the angle halfway between the angles
+					// perpendicular to a0->a1 (a10) and a1->a2 (a21)
+					double a2 = atan2(dy, dx);
+
+					out[(k + 1 - i) % (j - i - 1) + i].x = std::round(p1.x + buffer * cos(a2));
+					out[(k + 1 - i) % (j - i - 1) + i].y = std::round(p1.y + buffer * sin(a2));
+				}
+
+				out[j - 1].x = out[i].x;
+				out[j - 1].y = out[i].y;
+
+				i = j - 1;
+			}
+		}
+	}
+
+	return out;
+}

geometry.hpp

@@ -104,4 +104,7 @@ std::string overzoom(mvt_tile tile, int oz, int ox, int oy, int nz, int nx, int
 std::string overzoom(std::string s, int oz, int ox, int oy, int nz, int nx, int ny,
 		     int detail, int buffer, std::set<std::string> const &keep, bool do_compress);
 
+drawvec buffer_poly(drawvec const &geom, double buffer);
+bool is_continuous_poly(drawvec const &geom, int z, int tx, int ty, struct node *shared_nodes_map, size_t nodepos);
+
 #endif

main.cpp

@@ -3181,6 +3181,7 @@ int main(int argc, char **argv) {
 		{"no-tiny-polygon-reduction", no_argument, &prevent[P_TINY_POLYGON_REDUCTION], 1},
 		{"no-tiny-polygon-reduction-at-maximum-zoom", no_argument, &prevent[P_TINY_POLYGON_REDUCTION_AT_MAXZOOM], 1},
 		{"tiny-polygon-size", required_argument, 0, '~'},
+		{"buffer-polygons-outward", no_argument, &additional[A_BUFFER_POLYGONS_OUTWARD], 1},
 		{"no-simplification-of-shared-nodes", no_argument, &prevent[P_SIMPLIFY_SHARED_NODES], 1},
 		{"visvalingam", no_argument, &additional[A_VISVALINGAM], 1},
 

man/tippecanoe.1

@@ -637,6 +637,8 @@ the line or polygon within one tile unit of its proper location. You can probabl
 .IP \(bu 2
 \fB\fC\-\-tiny\-polygon\-size=\fR\fIsize\fP: Use the specified \fIsize\fP for tiny polygons instead of the default 2. Anything above 6 or so will lead to visible artifacts with the default tile detail.
 .IP \(bu 2
+\fB\fC\-aO\fR or \fB\fC\-\-buffer\-polygons\-outward\fR: Buffer polygons, except those that have been detected to be continous by \fB\fC\-\-no\-simplification\-of\-shared\-nodes\fR, outward slightly to prevent narrow polygon spindles fom collapsing away.
+.IP \(bu 2
 \fB\fC\-av\fR or \fB\fC\-\-visvalingam\fR: Use Visvalingam's simplification algorithm rather than Douglas\-Peucker's.
 .RE
 .SS Attempts to improve shared polygon boundaries

options.hpp

@@ -26,6 +26,7 @@
 #define A_HILBERT ((int) 'h')
 #define A_VISVALINGAM ((int) 'v')
 #define A_GENERATE_POLYGON_LABEL_POINTS ((int) 'P')
+#define A_BUFFER_POLYGONS_OUTWARD ((int) 'O')
 
 #define P_SIMPLIFY ((int) 's')
 #define P_SIMPLIFY_LOW ((int) 'S')

serial.cpp

@@ -221,6 +221,7 @@ std::string serialize_feature(serial_feature *sf, long long wx, long long wy) {
 	serialize_int(s, sf->segment);
 
 	write_geometry(sf->geometry, s, wx, wy);
+	serialize_byte(s, sf->continuous);
 
 	if (sf->index != 0) {
 		serialize_ulong_long(s, sf->index);
@@ -276,8 +277,10 @@ serial_feature deserialize_feature(std::string &geoms, unsigned z, unsigned tx,
 	sf.index = 0;
 	sf.label_point = 0;
 	sf.extent = 0;
+	sf.continuous = -1;  // unknown
 
 	sf.geometry = decode_geometry(&cp, z, tx, ty, sf.bbox, initial_x[sf.segment], initial_y[sf.segment]);
+	deserialize_byte(&cp, &sf.continuous);
 
 	if (sf.layer & (1 << FLAG_INDEX)) {
 		deserialize_ulong_long(&cp, &sf.index);

serial.hpp

@@ -43,6 +43,7 @@ struct serial_feature {
 
 	signed char t = 0;
 	signed char feature_minzoom = 0;
+	signed char continuous = -1;  // is this a continuous polygon?
 
 	bool has_id = false;
 	unsigned long long id = 0;

tests/buffer-polygons/4836bc56-0fab-4ffb-81c2-88142a05cec3-buffer.json

@@ -0,0 +1,2 @@
+{ "type": "Feature", "properties": { "felt:color": "#C93535", "felt:id": "d32e3884-e5fa-4ede-9bb0-a00d8094a846", "felt:locked": 0, "felt:ordering": 1695668996295712, "felt:routeMode": "NONE", "felt:showLength": 0, "felt:strokeOpacity": 1, "felt:strokeStyle": "solid", "felt:strokeWidth": 4, "felt:type": "Path" }, "geometry": { "type": "MultiPolygon", "coordinates": [ [ [ [ -122.2605783, 37.8081111 ], [ -122.2368893, 37.7923771 ], [ -122.236884, 37.7923751 ], [ -122.2368781, 37.7923755 ], [ -122.2368733, 37.7923782 ], [ -122.2368708, 37.7923824 ], [ -122.2368713, 37.7923871 ], [ -122.2368747, 37.7923909 ], [ -122.2605599, 37.8081224 ], [ -122.2612464, 37.8182874 ], [ -122.2473437, 37.8342159 ], [ -122.2473417, 37.8342203 ], [ -122.2473427, 37.8342249 ], [ -122.2473465, 37.8342284 ], [ -122.2473521, 37.83423 ], [ -122.2473579, 37.8342291 ], [ -122.2473623, 37.8342261 ], [ -122.2612673, 37.8182951 ], [ -122.2612693, 37.8182895 ], [ -122.2605823, 37.8081175 ], [ -122.2605812, 37.808114 ], [ -122.2605783, 37.8081111 ] ] ] ] } }
+{ "type": "Feature", "properties": { "felt:color": "#C93535", "felt:id": "72fa3875-764b-4785-8ae5-0d10b43acb11", "felt:locked": 0, "felt:ordering": 1695669003828484, "felt:routeMode": "NONE", "felt:showLength": 0, "felt:strokeOpacity": 1, "felt:strokeStyle": "solid", "felt:strokeWidth": 4, "felt:type": "Path" }, "geometry": { "type": "MultiPolygon", "coordinates": [ [ [ [ -122.2766184, 37.8184338 ], [ -122.3006514, 37.8228408 ], [ -122.3006573, 37.8228406 ], [ -122.3006623, 37.8228381 ], [ -122.3006651, 37.822834 ], [ -122.3006649, 37.8228294 ], [ -122.3006618, 37.8228254 ], [ -122.3006566, 37.8228232 ], [ -122.2766285, 37.8184171 ], [ -122.2676189, 37.8056726 ], [ -122.2676155, 37.8056696 ], [ -122.2676108, 37.8056681 ], [ -122.2245393, 37.8002451 ], [ -122.2337976, 37.786019 ], [ -122.2884614, 37.8002603 ], [ -122.3132655, 37.8112459 ], [ -122.3132712, 37.811247 ], [ -122.3132768, 37.8112457 ], [ -122.3132809, 37.8112424 ], [ -122.3132824, 37.8112379 ], [ -122.3132808, 37.8112334 ], [ -122.3132765, 37.8112301 ], [ -122.2884715, 37.8002441 ], [ -122.2884695, 37.8002434 ], [ -122.2337955, 37.7859994 ], [ -122.2337903, 37.7859991 ], [ -122.2337854, 37.7860007 ], [ -122.2337819, 37.7860039 ], [ -122.2245119, 37.8002479 ], [ -122.2245106, 37.800252 ], [ -122.2245119, 37.8002561 ], [ -122.2245153, 37.8002593 ], [ -122.2245202, 37.8002609 ], [ -122.2676018, 37.8056852 ], [ -122.2766111, 37.8184294 ], [ -122.2766142, 37.8184322 ], [ -122.2766184, 37.8184338 ] ] ] ] } }