From 75462f6a1bb7d07a65b17f4b4be39bd6ad844a46 Mon Sep 17 00:00:00 2001 From: Bertold Van den Bergh Date: Wed, 14 Aug 2019 14:51:57 +0200 Subject: [PATCH] Split up point reading and polygon checking. Add function to get list of points --- database/builder/makedb.sh | 1 + demo.c | 14 +- library/zonedetect.c | 372 ++++++++++++++++++++++++++++--------- library/zonedetect.h | 3 + 4 files changed, 297 insertions(+), 93 deletions(-) diff --git a/database/builder/makedb.sh b/database/builder/makedb.sh index 88421a6..8cfe749 100755 --- a/database/builder/makedb.sh +++ b/database/builder/makedb.sh @@ -10,6 +10,7 @@ mkdir -p naturalearth; cd naturalearth cd .. ./builder C naturalearth/ne_10m_admin_0_countries_lakes ./out/country16.bin 16 "Made with Natural Earth, placed in the Public Domain." ./builder C naturalearth/ne_10m_admin_0_countries_lakes ./out/country21.bin 21 "Made with Natural Earth, placed in the Public Domain." +./builder C naturalearth/ne_10m_admin_0_countries_lakes ./out/country9.bin 9 "Made with Natural Earth, placed in the Public Domain." mkdir timezone; cd timezone #wget https://github.com/evansiroky/timezone-boundary-builder/releases/download/2018i/timezones.shapefile.zip diff --git a/demo.c b/demo.c index bed9513..6c6ed59 100644 --- a/demo.c +++ b/demo.c @@ -29,7 +29,7 @@ #include #include "zonedetect.h" -void printResults(ZoneDetectResult *results, float safezone) +void printResults(ZoneDetect *cd, ZoneDetectResult *results, float safezone) { if(!results) { printf("No results\n"); @@ -40,6 +40,7 @@ void printResults(ZoneDetectResult *results, float safezone) while(results[index].lookupResult != ZD_LOOKUP_END) { printf("%s:\n", ZDLookupResultToString(results[index].lookupResult)); printf(" meta: %u\n", results[index].metaId); + printf(" polygon: %u\n", results[index].polygonId); if(results[index].data) { for(unsigned int i = 0; i < results[index].numFields; i++) { if(results[index].fieldNames[i] && results[index].data[i]) { @@ -47,6 +48,15 @@ void printResults(ZoneDetectResult *results, float safezone) } } } + size_t length = 0; + float* list = ZDPolygonToList(cd, results[index].polygonId, &length); + if(list){ + printf("Coordinate list length: %lu\n", length/2); + for(size_t i=0; i= (uint32_t)library->length) { @@ -342,105 +347,292 @@ static void ZDDecodePoint(uint64_t point, int32_t* lat, int32_t* lon){ *lat = (int32_t)ZDDecodeUnsignedToSigned(latu); *lon = (int32_t)ZDDecodeUnsignedToSigned(lonu); } - -static ZDLookupResult ZDPointInPolygon(const ZoneDetect *library, uint32_t polygonIndex, int32_t latFixedPoint, int32_t lonFixedPoint, uint64_t *distanceSqrMin) -{ - int32_t pointLat = 0, pointLon = 0, diffLat = 0, diffLon = 0, firstLat = 0, firstLon = 0, prevLat = 0, prevLon = 0; - lonFixedPoint -= 3; - int prevQuadrant = 0, winding = 0; - uint8_t done = 0, first = 1; +struct Reader{ + const ZoneDetect *library; + uint32_t polygonIndex; - uint32_t referenceStart=0, referenceEnd=0; - int32_t referenceDirection = 0; - - uint64_t numVertices = 1; + uint64_t numVertices; + + uint8_t done, first; + uint32_t referenceStart, referenceEnd; + int32_t referenceDirection; - if(library->version == 0){ - if(!ZDDecodeVariableLengthUnsigned(library, &polygonIndex, &numVertices)) return ZD_LOOKUP_PARSE_ERROR; - if(numVertices > 1000000) return ZD_LOOKUP_PARSE_ERROR; + int32_t pointLat, pointLon; + int32_t firstLat, firstLon; +}; + +static void ZDReaderInit(struct Reader *reader, const ZoneDetect *library, uint32_t polygonIndex){ + memset(reader, 0, sizeof(*reader)); + + reader->library = library; + reader->polygonIndex = polygonIndex; + + reader->first = 1; +} + +static int ZDReaderGetPoint(struct Reader *reader, int32_t *pointLat, int32_t *pointLon){ + int32_t diffLat = 0, diffLon = 0; + +readNewPoint: + if(reader->done){ + return 0; } - do{ - uint64_t point; - uint8_t referenceDone = 0; + if(reader->first && reader->library->version == 0){ + if(!ZDDecodeVariableLengthUnsigned(reader->library, &reader->polygonIndex, &reader->numVertices)) return -1; + if(!reader->numVertices) return -1; + } + + uint64_t point; + uint8_t referenceDone = 0; - if(library->version == 1){ - if(!referenceDirection){ - if(!ZDDecodeVariableLengthUnsigned(library, &polygonIndex, &point)) return ZD_LOOKUP_PARSE_ERROR; - }else{ - if(referenceDirection > 0){ - /* Read reference forward */ - if(!ZDDecodeVariableLengthUnsigned(library, &referenceStart, &point)) return ZD_LOOKUP_PARSE_ERROR; - if(referenceStart >= referenceEnd){ - referenceDone = 1; - } - }else if(referenceDirection < 0){ - /* Read reference backwards */ - if(!ZDDecodeVariableLengthUnsignedReverse(library, &referenceStart, &point)) return ZD_LOOKUP_PARSE_ERROR; - if(referenceStart <= referenceEnd){ - referenceDone = 1; - } + if(reader->library->version == 1){ + if(!reader->referenceDirection){ + if(!ZDDecodeVariableLengthUnsigned(reader->library, &reader->polygonIndex, &point)) return -1; + }else{ + if(reader->referenceDirection > 0){ + /* Read reference forward */ + if(!ZDDecodeVariableLengthUnsigned(reader->library, &reader->referenceStart, &point)) return -1; + if(reader->referenceStart >= reader->referenceEnd){ + referenceDone = 1; + } + }else if(reader->referenceDirection < 0){ + /* Read reference backwards */ + if(!ZDDecodeVariableLengthUnsignedReverse(reader->library, &reader->referenceStart, &point)) return -1; + if(reader->referenceStart <= reader->referenceEnd){ + referenceDone = 1; } } + } - if(!point){ - /* This is a special marker, it is not allowed in reference mode */ - if(referenceDirection){ - printf("BUG, marker in reference mode?\n"); - return ZD_LOOKUP_PARSE_ERROR; - } + if(!point){ + /* This is a special marker, it is not allowed in reference mode */ + if(reader->referenceDirection){ + printf("BUG, marker in reference mode?\n"); + return -1; + } - uint64_t value; - if(!ZDDecodeVariableLengthUnsigned(library, &polygonIndex, &value)) return ZD_LOOKUP_PARSE_ERROR; + uint64_t value; + if(!ZDDecodeVariableLengthUnsigned(reader->library, &reader->polygonIndex, &value)) return -1; - if(value == 0){ - done = 1; - }else if(value == 1){ - int32_t diff; - int64_t start; - if(!ZDDecodeVariableLengthUnsigned(library, &polygonIndex, (uint64_t*)&start)) return ZD_LOOKUP_PARSE_ERROR; - if(!ZDDecodeVariableLengthSigned(library, &polygonIndex, &diff)) return ZD_LOOKUP_PARSE_ERROR; - - referenceStart = library->dataOffset+(uint32_t)start; - referenceEnd = library->dataOffset+(uint32_t)(start + diff); - referenceDirection = diff; - if(diff < 0){ - referenceStart--; - referenceEnd--; - } - continue; + if(value == 0){ + reader->done = 1; + }else if(value == 1){ + int32_t diff; + int64_t start; + if(!ZDDecodeVariableLengthUnsigned(reader->library, &reader->polygonIndex, (uint64_t*)&start)) return -1; + if(!ZDDecodeVariableLengthSigned(reader->library, &reader->polygonIndex, &diff)) return -1; + + reader->referenceStart = reader->library->dataOffset+(uint32_t)start; + reader->referenceEnd = reader->library->dataOffset+(uint32_t)(start + diff); + reader->referenceDirection = diff; + if(diff < 0){ + reader->referenceStart--; + reader->referenceEnd--; } - }else{ - ZDDecodePoint(point, &diffLat, &diffLon); - if(referenceDirection < 0){ - diffLat = -diffLat; - diffLon = -diffLon; - } + goto readNewPoint; + } + }else{ + ZDDecodePoint(point, &diffLat, &diffLon); + if(reader->referenceDirection < 0){ + diffLat = -diffLat; + diffLon = -diffLon; + } + } + } + + if(reader->library->version == 0){ + if(!ZDDecodeVariableLengthSigned(reader->library, &reader->polygonIndex, &diffLat)) return -1; + if(!ZDDecodeVariableLengthSigned(reader->library, &reader->polygonIndex, &diffLon)) return -1; + } + + if(!reader->done){ + reader->pointLat += diffLat; + reader->pointLon += diffLon; + if(reader->first) { + reader->firstLat = reader->pointLat; + reader->firstLon = reader->pointLon; + } + } else { + /* Close the polygon (the closing point is not encoded) */ + reader->pointLat = reader->firstLat; + reader->pointLon = reader->firstLon; + } + + reader->first = 0; + + if(reader->library->version == 0){ + reader->numVertices--; + if(!reader->numVertices){ + reader->done = 1; + } + } + + if(referenceDone){ + reader->referenceDirection = 0; + } + + if(pointLat){ + *pointLat = reader->pointLat; + } + + if(pointLon){ + *pointLon = reader->pointLon; + } + + return 1; +} + +static int ZDFindPolygon(const ZoneDetect *library, uint32_t wantedId, uint32_t* metadataIndexPtr, uint32_t* polygonIndexPtr){ + uint32_t polygonId = 0; + uint32_t bboxIndex = library->bboxOffset; + + uint32_t metadataIndex = 0, polygonIndex = 0; + + while(bboxIndex < library->metadataOffset) { + uint64_t polygonIndexDelta; + int32_t metadataIndexDelta; + int32_t tmp; + if(!ZDDecodeVariableLengthSigned(library, &bboxIndex, &tmp)) break; + if(!ZDDecodeVariableLengthSigned(library, &bboxIndex, &tmp)) break; + if(!ZDDecodeVariableLengthSigned(library, &bboxIndex, &tmp)) break; + if(!ZDDecodeVariableLengthSigned(library, &bboxIndex, &tmp)) break; + if(!ZDDecodeVariableLengthSigned(library, &bboxIndex, &metadataIndexDelta)) break; + if(!ZDDecodeVariableLengthUnsigned(library, &bboxIndex, &polygonIndexDelta)) break; + + metadataIndex += (uint32_t)metadataIndexDelta; + polygonIndex += (uint32_t)polygonIndexDelta; + + if(polygonId == wantedId){ + if(metadataIndexPtr){ + metadataIndex += library->metadataOffset; + *metadataIndexPtr = metadataIndex; + } + if(polygonIndexPtr){ + polygonIndex += library->dataOffset; + *polygonIndexPtr = polygonIndex; + } + return 1; + } + + polygonId ++; + } + + return 0; +} + +static int32_t* ZDPolygonToListInternal(const ZoneDetect *library, uint32_t polygonIndex, size_t* length){ + struct Reader reader; + ZDReaderInit(&reader, library, polygonIndex); + + size_t listLength = 2 * 100; + size_t listIndex = 0; + + int32_t* list = malloc(sizeof(int32_t) * listLength); + if(!list){ + goto fail; + } + + while(1){ + int32_t pointLat, pointLon; + int result = ZDReaderGetPoint(&reader, &pointLat, &pointLon); + if(result < 0){ + goto fail; + }else if(result == 0){ + break; + } + + if(listIndex >= listLength){ + listLength *= 2; + if(listLength >= 1048576){ + goto fail; + } + + list = realloc(list, sizeof(int32_t) * listLength); + if(!list){ + goto fail; } } - if(library->version == 0){ - numVertices--; - if(!numVertices){ - done = 1; - } + list[listIndex++] = pointLat; + list[listIndex++] = pointLon; + } - if(!ZDDecodeVariableLengthSigned(library, &polygonIndex, &diffLat)) return ZD_LOOKUP_PARSE_ERROR; - if(!ZDDecodeVariableLengthSigned(library, &polygonIndex, &diffLon)) return ZD_LOOKUP_PARSE_ERROR; - } + if(length){ + *length = listIndex; + } - if(!done){ - pointLat += diffLat; - pointLon += diffLon; - if(first) { - firstLat = pointLat; - firstLon = pointLon; - } - } else { - /* Close the polygon (the closing point is not encoded) */ - pointLat = firstLat; - pointLon = firstLon; + return list; + +fail: + if(list){ + free(list); + } + return NULL; +} + +float* ZDPolygonToList(const ZoneDetect *library, uint32_t polygonId, size_t* lengthPtr){ + uint32_t polygonIndex; + int32_t* data = NULL; + float* flData = NULL; + + if(!ZDFindPolygon(library, polygonId, NULL, &polygonIndex)){ + goto fail; + } + + size_t length = 0; + data = ZDPolygonToListInternal(library, polygonIndex, &length); + + if(!data){ + goto fail; + } + + flData = malloc(sizeof(float) * length); + if(!flData){ + goto fail; + } + + for(size_t i = 0; iprecision); + flData[i+1] = ZDFixedPointToFloat(lon, 180, library->precision); + } + + if(lengthPtr){ + *lengthPtr = length; + } + + return flData; + +fail: + if(data){ + free(data); + } + if(flData){ + free(flData); + } + return NULL; +} + +static ZDLookupResult ZDPointInPolygon(const ZoneDetect *library, uint32_t polygonIndex, int32_t latFixedPoint, int32_t lonFixedPoint, uint64_t *distanceSqrMin) +{ + int32_t pointLat, pointLon, prevLat = 0, prevLon = 0; + int prevQuadrant = 0, winding = 0; + + uint8_t first = 1; + + struct Reader reader; + ZDReaderInit(&reader, library, polygonIndex); + + while(1){ + int result = ZDReaderGetPoint(&reader, &pointLat, &pointLon); + if(result < 0){ + return ZD_LOOKUP_PARSE_ERROR; + }else if(result == 0){ + break; } /* Check if point is ON the border */ @@ -555,15 +747,8 @@ static ZDLookupResult ZDPointInPolygon(const ZoneDetect *library, uint32_t polyg prevQuadrant = quadrant; prevLat = pointLat; prevLon = pointLon; - - if(first){ - first = 0; - } - - if(referenceDone){ - referenceDirection = 0; - } - }while(!done); + first = 0; + }; if(winding == -4) { return ZD_LOOKUP_IN_ZONE; @@ -690,6 +875,8 @@ ZoneDetectResult *ZDLookup(const ZoneDetect *library, float lat, float lon, floa return NULL; } + uint32_t polygonId = 0; + while(bboxIndex < library->metadataOffset) { int32_t minLat, minLon, maxLat, maxLon, metadataIndexDelta; uint64_t polygonIndexDelta; @@ -720,6 +907,7 @@ ZoneDetectResult *ZDLookup(const ZoneDetect *library, float lat, float lon, floa if(newResults) { results = newResults; + results[numResults].polygonId = polygonId; results[numResults].metaId = metadataIndex; results[numResults].numFields = library->numFields; results[numResults].fieldNames = library->fieldNames; @@ -735,6 +923,8 @@ ZoneDetectResult *ZDLookup(const ZoneDetect *library, float lat, float lon, floa /* The data is sorted along minLat */ break; } + + polygonId++; } /* Clean up results */ diff --git a/library/zonedetect.h b/library/zonedetect.h index d23bbdf..3b91dff 100644 --- a/library/zonedetect.h +++ b/library/zonedetect.h @@ -52,6 +52,7 @@ typedef enum { typedef struct { ZDLookupResult lookupResult; + uint32_t polygonId; uint32_t metaId; uint8_t numFields; char **fieldNames; @@ -78,6 +79,8 @@ ZD_EXPORT const char *ZDLookupResultToString(ZDLookupResult result); ZD_EXPORT int ZDSetErrorHandler(void (*handler)(int, int)); ZD_EXPORT const char *ZDGetErrorString(int errZD); +ZD_EXPORT float* ZDPolygonToList(const ZoneDetect *library, uint32_t polygonId, size_t* length); + #ifdef __cplusplus } #endif