From c832b341608f25eb5e46da61af6d344d15613ed2 Mon Sep 17 00:00:00 2001 From: Eric Fischer Date: Thu, 3 Dec 2015 13:12:34 -0800 Subject: [PATCH 1/3] Better handling of features that cross the antimeridian --- geojson.c | 2 +- projection.c | 24 ++---------------------- projection.h | 4 ++-- tile.cc | 20 ++++++++++++++++++++ 4 files changed, 25 insertions(+), 25 deletions(-) diff --git a/geojson.c b/geojson.c index cc9d499..2d347e8 100644 --- a/geojson.c +++ b/geojson.c @@ -160,7 +160,7 @@ void parse_geometry(int t, json_object *j, unsigned *bbox, long long *fpos, FILE } } else { if (j->length >= 2 && j->array[0]->type == JSON_NUMBER && j->array[1]->type == JSON_NUMBER) { - unsigned x, y; + long long x, y; double lon = j->array[0]->number; double lat = j->array[1]->number; latlon2tile(lat, lon, 32, &x, &y); diff --git a/projection.c b/projection.c index effb347..839370e 100644 --- a/projection.c +++ b/projection.c @@ -2,39 +2,19 @@ #include "projection.h" // http://wiki.openstreetmap.org/wiki/Slippy_map_tilenames -void latlon2tile(double lat, double lon, int zoom, unsigned int *x, unsigned int *y) { +void latlon2tile(double lat, double lon, int zoom, long long *x, long long *y) { double lat_rad = lat * M_PI / 180; unsigned long long n = 1LL << zoom; long long llx = n * ((lon + 180) / 360); long long lly = n * (1 - (log(tan(lat_rad) + 1 / cos(lat_rad)) / M_PI)) / 2; - if (lat >= 85.0511) { - lly = 0; - } - if (lat <= -85.0511) { - lly = n - 1; - } - - if (llx < 0) { - llx = 0; - } - if (lly < 0) { - lly = 0; - } - if (llx >= n) { - llx = n - 1; - } - if (lly >= n) { - lly = n - 1; - } - *x = llx; *y = lly; } // http://wiki.openstreetmap.org/wiki/Slippy_map_tilenames -void tile2latlon(unsigned int x, unsigned int y, int zoom, double *lat, double *lon) { +void tile2latlon(long long x, long long y, int zoom, double *lat, double *lon) { unsigned long long n = 1LL << zoom; *lon = 360.0 * x / n - 180.0; *lat = atan(sinh(M_PI * (1 - 2.0 * y / n))) * 180.0 / M_PI; diff --git a/projection.h b/projection.h index af7002e..20d67ec 100644 --- a/projection.h +++ b/projection.h @@ -1,4 +1,4 @@ -void latlon2tile(double lat, double lon, int zoom, unsigned int *x, unsigned int *y); -void tile2latlon(unsigned int x, unsigned int y, int zoom, double *lat, double *lon); +void latlon2tile(double lat, double lon, int zoom, long long *x, long long *y); +void tile2latlon(long long x, long long y, int zoom, double *lat, double *lon); unsigned long long encode(unsigned int wx, unsigned int wy); void decode(unsigned long long index, unsigned *wx, unsigned *wy); diff --git a/tile.cc b/tile.cc index bdd95a6..6daa00a 100644 --- a/tile.cc +++ b/tile.cc @@ -544,6 +544,26 @@ long long write_tile(char **geoms, char *metabase, char *stringpool, int z, unsi continue; } + if (z == 0) { + if (bbox[0] < 0 || bbox[2] > 1LL << 32) { + // If the geometry extends off the edge of the world, concatenate on another copy + // shifted by 360 degrees, and then make sure both copies get clipped down to size. + + unsigned n = geom.size(); + for (unsigned i = 0; i < n; i++) { + geom.push_back(draw(geom[i].op, geom[i].x - (1LL << 32), geom[i].y)); + } + for (unsigned i = 0; i < n; i++) { + geom.push_back(draw(geom[i].op, geom[i].x + (1LL << 32), geom[i].y)); + } + + bbox[0] = 0; + bbox[2] = 1LL << 32; + + quick = -1; + } + } + if (quick != 1) { if (t == VT_LINE) { geom = clip_lines(geom, z, line_detail, buffer); From 9faa625e753f7f57287bf9551b3a8990a343d1b8 Mon Sep 17 00:00:00 2001 From: Eric Fischer Date: Thu, 3 Dec 2015 14:42:59 -0800 Subject: [PATCH 2/3] Calculate bounding boxes better for features that cross the antimeridian --- geojson.c | 36 ++++++++++++++++++++++++++++++------ 1 file changed, 30 insertions(+), 6 deletions(-) diff --git a/geojson.c b/geojson.c index 2d347e8..041492c 100644 --- a/geojson.c +++ b/geojson.c @@ -138,7 +138,7 @@ void serialize_string(FILE *out, const char *s, long long *fpos, const char *fna *fpos += len + 1; } -void parse_geometry(int t, json_object *j, unsigned *bbox, long long *fpos, FILE *out, int op, const char *fname, json_pull *source, long long *wx, long long *wy, int *initialized) { +void parse_geometry(int t, json_object *j, long long *bbox, long long *fpos, FILE *out, int op, const char *fname, json_pull *source, long long *wx, long long *wy, int *initialized) { if (j == NULL || j->type != JSON_ARRAY) { fprintf(stderr, "%s:%d: expected array for type %d\n", fname, source->line, t); return; @@ -437,7 +437,7 @@ long long addpool(struct memfile *poolfile, struct memfile *treefile, char *s, c return off; } -int serialize_geometry(json_object *geometry, json_object *properties, const char *reading, json_pull *jp, long long *seq, long long *metapos, long long *geompos, long long *indexpos, struct pool *exclude, struct pool *include, int exclude_all, FILE *metafile, FILE *geomfile, FILE *indexfile, struct memfile *poolfile, struct memfile *treefile, const char *fname, int maxzoom, int layer, double droprate, unsigned *file_bbox, json_object *tippecanoe) { +int serialize_geometry(json_object *geometry, json_object *properties, const char *reading, json_pull *jp, long long *seq, long long *metapos, long long *geompos, long long *indexpos, struct pool *exclude, struct pool *include, int exclude_all, FILE *metafile, FILE *geomfile, FILE *indexfile, struct memfile *poolfile, struct memfile *treefile, const char *fname, int maxzoom, int layer, double droprate, long long *file_bbox, json_object *tippecanoe) { json_object *geometry_type = json_hash_get(geometry, "type"); if (geometry_type == NULL) { static int warned = 0; @@ -492,7 +492,7 @@ int serialize_geometry(json_object *geometry, json_object *properties, const cha } } - unsigned bbox[] = {UINT_MAX, UINT_MAX, 0, 0}; + long long bbox[] = {UINT_MAX, UINT_MAX, 0, 0}; int nprop = 0; if (properties != NULL && properties->type == JSON_HASH) { @@ -595,7 +595,13 @@ int serialize_geometry(json_object *geometry, json_object *properties, const cha struct index index; index.start = geomstart; index.end = *geompos; - index.index = encode(bbox[0] / 2 + bbox[2] / 2, bbox[1] / 2 + bbox[3] / 2); + + // Calculate the center even if off the edge of the plane, + // and then mask to bring it back into the addressable area + long long midx = (bbox[0] / 2 + bbox[2] / 2) & ((1LL << 32) - 1); + long long midy = (bbox[1] / 2 + bbox[3] / 2) & ((1LL << 32) - 1); + index.index = encode(midx, midy); + fwrite_check(&index, sizeof(struct index), 1, indexfile, fname); *indexpos += sizeof(struct index); @@ -620,7 +626,7 @@ int serialize_geometry(json_object *geometry, json_object *properties, const cha return 1; } -void parse_json(json_pull *jp, const char *reading, long long *seq, long long *metapos, long long *geompos, long long *indexpos, struct pool *exclude, struct pool *include, int exclude_all, FILE *metafile, FILE *geomfile, FILE *indexfile, struct memfile *poolfile, struct memfile *treefile, char *fname, int maxzoom, int layer, double droprate, unsigned *file_bbox) { +void parse_json(json_pull *jp, const char *reading, long long *seq, long long *metapos, long long *geompos, long long *indexpos, struct pool *exclude, struct pool *include, int exclude_all, FILE *metafile, FILE *geomfile, FILE *indexfile, struct memfile *poolfile, struct memfile *treefile, char *fname, int maxzoom, int layer, double droprate, long long *file_bbox) { long long found_hashes = 0; long long found_features = 0; long long found_geometries = 0; @@ -814,7 +820,7 @@ int read_json(int argc, char **argv, char *fname, const char *layername, int max memfile_write(treefile, &p, sizeof(struct stringpool)); } - unsigned file_bbox[] = {UINT_MAX, UINT_MAX, 0, 0}; + long long file_bbox[] = {UINT_MAX, UINT_MAX, 0, 0}; unsigned midx = 0, midy = 0; long long seq = 0; @@ -1159,6 +1165,24 @@ int read_json(int argc, char **argv, char *fname, const char *layername, int max midlat = (maxlat + minlat) / 2; midlon = (maxlon + minlon) / 2; + // If the bounding box extends off the plane on either side, + // a feature wrapped across the date line, so the width of the + // bounding box is the whole world. + if (file_bbox[0] < 0) { + file_bbox[0] = 0; + file_bbox[2] = (1LL << 32) - 1; + } + if (file_bbox[2] > (1LL << 32) - 1) { + file_bbox[0] = 0; + file_bbox[2] = (1LL << 32) - 1; + } + if (file_bbox[1] < 0) { + file_bbox[1] = 0; + } + if (file_bbox[3] < (1LL << 32) - 1) { + file_bbox[3] = (1LL << 32) - 1; + } + tile2latlon(file_bbox[0], file_bbox[1], 32, &maxlat, &minlon); tile2latlon(file_bbox[2], file_bbox[3], 32, &minlat, &maxlon); From 5506288273e69dfbbc66f37c079cf93876d04c20 Mon Sep 17 00:00:00 2001 From: Eric Fischer Date: Thu, 3 Dec 2015 15:19:54 -0800 Subject: [PATCH 3/3] Bump version number --- CHANGELOG.md | 5 +++++ version.h | 2 +- 2 files changed, 6 insertions(+), 1 deletion(-) diff --git a/CHANGELOG.md b/CHANGELOG.md index 90a4816..3b2e2ec 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -1,3 +1,8 @@ +## 1.4.1 + +* Features that cross the antimeridian are split into two parts instead + of being partially lost off the edge + ## 1.4.0 * More polygon correctness diff --git a/version.h b/version.h index 087d6a1..5c66f11 100644 --- a/version.h +++ b/version.h @@ -1 +1 @@ -#define VERSION "tippecanoe v1.4.0\n" +#define VERSION "tippecanoe v1.4.1\n"