Index: Makefile.in ================================================================== --- Makefile.in +++ Makefile.in @@ -347,11 +347,12 @@ SRC += \ $(TOP)/ext/icu/sqliteicu.h \ $(TOP)/ext/icu/icu.c SRC += \ $(TOP)/ext/rtree/rtree.h \ - $(TOP)/ext/rtree/rtree.c + $(TOP)/ext/rtree/rtree.c \ + $(TOP)/ext/rtree/geopoly.c SRC += \ $(TOP)/ext/session/sqlite3session.c \ $(TOP)/ext/session/sqlite3session.h SRC += \ $(TOP)/ext/rbu/sqlite3rbu.h \ @@ -547,11 +548,12 @@ $(TOP)/ext/fts3/fts3.h \ $(TOP)/ext/fts3/fts3Int.h \ $(TOP)/ext/fts3/fts3_hash.h \ $(TOP)/ext/fts3/fts3_tokenizer.h EXTHDR += \ - $(TOP)/ext/rtree/rtree.h + $(TOP)/ext/rtree/rtree.h \ + $(TOP)/ext/rtree/geopoly.c EXTHDR += \ $(TOP)/ext/icu/sqliteicu.h EXTHDR += \ $(TOP)/ext/rtree/sqlite3rtree.h Index: Makefile.msc ================================================================== --- Makefile.msc +++ Makefile.msc @@ -336,11 +336,11 @@ # the Windows platform. # !IFNDEF OPT_FEATURE_FLAGS !IF $(MINIMAL_AMALGAMATION)==0 OPT_FEATURE_FLAGS = $(OPT_FEATURE_FLAGS) -DSQLITE_ENABLE_FTS3=1 -OPT_FEATURE_FLAGS = $(OPT_FEATURE_FLAGS) -DSQLITE_ENABLE_RTREE=1 +OPT_FEATURE_FLAGS = $(OPT_FEATURE_FLAGS) -DSQLITE_ENABLE_RTREE=1 -DSQLITE_ENABLE_GEOPOLY=1 !ENDIF OPT_FEATURE_FLAGS = $(OPT_FEATURE_FLAGS) -DSQLITE_ENABLE_COLUMN_METADATA=1 !ENDIF # Should the session extension be enabled? If so, add compilation options @@ -1398,10 +1398,11 @@ $(TOP)\ext\fts3\fts3Int.h \ $(TOP)\ext\fts3\fts3_hash.h \ $(TOP)\ext\fts3\fts3_tokenizer.h \ $(TOP)\ext\icu\sqliteicu.h \ $(TOP)\ext\rtree\rtree.h \ + $(TOP)\ext\rtree\geopoly.c \ $(TOP)\ext\rbu\sqlite3rbu.h \ $(TOP)\ext\session\sqlite3session.h # Generated source code files # @@ -1571,11 +1572,12 @@ $(TOP)\ext\fts3\fts3.h \ $(TOP)\ext\fts3\fts3Int.h \ $(TOP)\ext\fts3\fts3_hash.h \ $(TOP)\ext\fts3\fts3_tokenizer.h EXTHDR = $(EXTHDR) \ - $(TOP)\ext\rtree\rtree.h + $(TOP)\ext\rtree\rtree.h \ + $(TOP)\ext\rtree\geopoly.c EXTHDR = $(EXTHDR) \ $(TOP)\ext\icu\sqliteicu.h EXTHDR = $(EXTHDR) \ $(TOP)\ext\rtree\sqlite3rtree.h EXTHDR = $(EXTHDR) \ Index: configure ================================================================== --- configure +++ configure @@ -909,10 +909,11 @@ enable_fts3 enable_fts4 enable_fts5 enable_json1 enable_update_limit +enable_geopoly enable_rtree enable_session enable_gcov ' ac_precious_vars='build_alias @@ -1561,10 +1562,11 @@ --enable-fts3 Enable the FTS3 extension --enable-fts4 Enable the FTS4 extension --enable-fts5 Enable the FTS5 extension --enable-json1 Enable the JSON1 extension --enable-update-limit Enable the UPDATE/DELETE LIMIT clause + --enable-geopoly Enable the GEOPOLY extension --enable-rtree Enable the RTREE extension --enable-session Enable the SESSION extension --enable-gcov Enable coverage testing using gcov Optional Packages: @@ -3930,17 +3932,17 @@ if ${lt_cv_nm_interface+:} false; then : $as_echo_n "(cached) " >&6 else lt_cv_nm_interface="BSD nm" echo "int some_variable = 0;" > conftest.$ac_ext - (eval echo "\"\$as_me:3935: $ac_compile\"" >&5) + (eval echo "\"\$as_me:3937: $ac_compile\"" >&5) (eval "$ac_compile" 2>conftest.err) cat conftest.err >&5 - (eval echo "\"\$as_me:3938: $NM \\\"conftest.$ac_objext\\\"\"" >&5) + (eval echo "\"\$as_me:3940: $NM \\\"conftest.$ac_objext\\\"\"" >&5) (eval "$NM \"conftest.$ac_objext\"" 2>conftest.err > conftest.out) cat conftest.err >&5 - (eval echo "\"\$as_me:3941: output\"" >&5) + (eval echo "\"\$as_me:3943: output\"" >&5) cat conftest.out >&5 if $GREP 'External.*some_variable' conftest.out > /dev/null; then lt_cv_nm_interface="MS dumpbin" fi rm -f conftest* @@ -5142,11 +5144,11 @@ fi rm -rf conftest* ;; *-*-irix6*) # Find out which ABI we are using. - echo '#line 5147 "configure"' > conftest.$ac_ext + echo '#line 5149 "configure"' > conftest.$ac_ext if { { eval echo "\"\$as_me\":${as_lineno-$LINENO}: \"$ac_compile\""; } >&5 (eval $ac_compile) 2>&5 ac_status=$? $as_echo "$as_me:${as_lineno-$LINENO}: \$? = $ac_status" >&5 test $ac_status = 0; }; then @@ -6667,15 +6669,15 @@ # The option is referenced via a variable to avoid confusing sed. lt_compile=`echo "$ac_compile" | $SED \ -e 's:.*FLAGS}\{0,1\} :&$lt_compiler_flag :; t' \ -e 's: [^ ]*conftest\.: $lt_compiler_flag&:; t' \ -e 's:$: $lt_compiler_flag:'` - (eval echo "\"\$as_me:6672: $lt_compile\"" >&5) + (eval echo "\"\$as_me:6674: $lt_compile\"" >&5) (eval "$lt_compile" 2>conftest.err) ac_status=$? cat conftest.err >&5 - echo "$as_me:6676: \$? = $ac_status" >&5 + echo "$as_me:6678: \$? = $ac_status" >&5 if (exit $ac_status) && test -s "$ac_outfile"; then # The compiler can only warn and ignore the option if not recognized # So say no if there are warnings other than the usual output. $ECHO "X$_lt_compiler_boilerplate" | $Xsed -e '/^$/d' >conftest.exp $SED '/^$/d; /^ *+/d' conftest.err >conftest.er2 @@ -7006,15 +7008,15 @@ # The option is referenced via a variable to avoid confusing sed. lt_compile=`echo "$ac_compile" | $SED \ -e 's:.*FLAGS}\{0,1\} :&$lt_compiler_flag :; t' \ -e 's: [^ ]*conftest\.: $lt_compiler_flag&:; t' \ -e 's:$: $lt_compiler_flag:'` - (eval echo "\"\$as_me:7011: $lt_compile\"" >&5) + (eval echo "\"\$as_me:7013: $lt_compile\"" >&5) (eval "$lt_compile" 2>conftest.err) ac_status=$? cat conftest.err >&5 - echo "$as_me:7015: \$? = $ac_status" >&5 + echo "$as_me:7017: \$? = $ac_status" >&5 if (exit $ac_status) && test -s "$ac_outfile"; then # The compiler can only warn and ignore the option if not recognized # So say no if there are warnings other than the usual output. $ECHO "X$_lt_compiler_boilerplate" | $Xsed -e '/^$/d' >conftest.exp $SED '/^$/d; /^ *+/d' conftest.err >conftest.er2 @@ -7111,15 +7113,15 @@ # with a dollar sign (not a hyphen), so the echo should work correctly. lt_compile=`echo "$ac_compile" | $SED \ -e 's:.*FLAGS}\{0,1\} :&$lt_compiler_flag :; t' \ -e 's: [^ ]*conftest\.: $lt_compiler_flag&:; t' \ -e 's:$: $lt_compiler_flag:'` - (eval echo "\"\$as_me:7116: $lt_compile\"" >&5) + (eval echo "\"\$as_me:7118: $lt_compile\"" >&5) (eval "$lt_compile" 2>out/conftest.err) ac_status=$? cat out/conftest.err >&5 - echo "$as_me:7120: \$? = $ac_status" >&5 + echo "$as_me:7122: \$? = $ac_status" >&5 if (exit $ac_status) && test -s out/conftest2.$ac_objext then # The compiler can only warn and ignore the option if not recognized # So say no if there are warnings $ECHO "X$_lt_compiler_boilerplate" | $Xsed -e '/^$/d' > out/conftest.exp @@ -7166,15 +7168,15 @@ # with a dollar sign (not a hyphen), so the echo should work correctly. lt_compile=`echo "$ac_compile" | $SED \ -e 's:.*FLAGS}\{0,1\} :&$lt_compiler_flag :; t' \ -e 's: [^ ]*conftest\.: $lt_compiler_flag&:; t' \ -e 's:$: $lt_compiler_flag:'` - (eval echo "\"\$as_me:7171: $lt_compile\"" >&5) + (eval echo "\"\$as_me:7173: $lt_compile\"" >&5) (eval "$lt_compile" 2>out/conftest.err) ac_status=$? cat out/conftest.err >&5 - echo "$as_me:7175: \$? = $ac_status" >&5 + echo "$as_me:7177: \$? = $ac_status" >&5 if (exit $ac_status) && test -s out/conftest2.$ac_objext then # The compiler can only warn and ignore the option if not recognized # So say no if there are warnings $ECHO "X$_lt_compiler_boilerplate" | $Xsed -e '/^$/d' > out/conftest.exp @@ -9546,11 +9548,11 @@ lt_cv_dlopen_self=cross else lt_dlunknown=0; lt_dlno_uscore=1; lt_dlneed_uscore=2 lt_status=$lt_dlunknown cat > conftest.$ac_ext <<_LT_EOF -#line 9551 "configure" +#line 9553 "configure" #include "confdefs.h" #if HAVE_DLFCN_H #include #endif @@ -9642,11 +9644,11 @@ lt_cv_dlopen_self_static=cross else lt_dlunknown=0; lt_dlno_uscore=1; lt_dlneed_uscore=2 lt_status=$lt_dlunknown cat > conftest.$ac_ext <<_LT_EOF -#line 9647 "configure" +#line 9649 "configure" #include "confdefs.h" #if HAVE_DLFCN_H #include #endif @@ -11627,10 +11629,24 @@ fi if test "${enable_udlimit}" = "yes" ; then OPT_FEATURE_FLAGS="${OPT_FEATURE_FLAGS} -DSQLITE_ENABLE_UPDATE_DELETE_LIMIT" fi + +######### +# See whether we should enable GEOPOLY +# Check whether --enable-geopoly was given. +if test "${enable_geopoly+set}" = set; then : + enableval=$enable_geopoly; enable_geopoly=yes +else + enable_geopoly=no +fi + +if test "${enable_geopoly}" = "yes" ; then + OPT_FEATURE_FLAGS="${OPT_FEATURE_FLAGS} -DSQLITE_ENABLE_GEOPOLY" + enable_rtree=yes +fi ######### # See whether we should enable RTREE # Check whether --enable-rtree was given. if test "${enable_rtree+set}" = set; then : Index: configure.ac ================================================================== --- configure.ac +++ configure.ac @@ -657,10 +657,20 @@ [Enable the UPDATE/DELETE LIMIT clause]), [enable_udlimit=yes],[enable_udlimit=no]) if test "${enable_udlimit}" = "yes" ; then OPT_FEATURE_FLAGS="${OPT_FEATURE_FLAGS} -DSQLITE_ENABLE_UPDATE_DELETE_LIMIT" fi + +######### +# See whether we should enable GEOPOLY +AC_ARG_ENABLE(geopoly, AC_HELP_STRING([--enable-geopoly], + [Enable the GEOPOLY extension]), + [enable_geopoly=yes],[enable_geopoly=no]) +if test "${enable_geopoly}" = "yes" ; then + OPT_FEATURE_FLAGS="${OPT_FEATURE_FLAGS} -DSQLITE_ENABLE_GEOPOLY" + enable_rtree=yes +fi ######### # See whether we should enable RTREE AC_ARG_ENABLE(rtree, AC_HELP_STRING([--enable-rtree], [Enable the RTREE extension]), ADDED ext/rtree/geopoly.c Index: ext/rtree/geopoly.c ================================================================== --- ext/rtree/geopoly.c +++ ext/rtree/geopoly.c @@ -0,0 +1,853 @@ +/* +** 2018-05-25 +** +** The author disclaims copyright to this source code. In place of +** a legal notice, here is a blessing: +** +** May you do good and not evil. +** May you find forgiveness for yourself and forgive others. +** May you share freely, never taking more than you give. +** +****************************************************************************** +** +** This file implements an alternative R-Tree virtual table that +** uses polygons to express the boundaries of 2-dimensional objects. +** +** This file is #include-ed onto the end of "rtree.c" so that it has +** access to all of the R-Tree internals. +*/ +#include + +/* Enable -DGEOPOLY_ENABLE_DEBUG for debugging facilities */ +#ifdef GEOPOLY_ENABLE_DEBUG + static int geo_debug = 0; +# define GEODEBUG(X) if(geo_debug)printf X +#else +# define GEODEBUG(X) +#endif + +#ifndef JSON_NULL /* The following stuff repeats things found in json1 */ +/* +** Versions of isspace(), isalnum() and isdigit() to which it is safe +** to pass signed char values. +*/ +#ifdef sqlite3Isdigit + /* Use the SQLite core versions if this routine is part of the + ** SQLite amalgamation */ +# define safe_isdigit(x) sqlite3Isdigit(x) +# define safe_isalnum(x) sqlite3Isalnum(x) +# define safe_isxdigit(x) sqlite3Isxdigit(x) +#else + /* Use the standard library for separate compilation */ +#include /* amalgamator: keep */ +# define safe_isdigit(x) isdigit((unsigned char)(x)) +# define safe_isalnum(x) isalnum((unsigned char)(x)) +# define safe_isxdigit(x) isxdigit((unsigned char)(x)) +#endif + +/* +** Growing our own isspace() routine this way is twice as fast as +** the library isspace() function. +*/ +static const char geopolyIsSpace[] = { + 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 0, 0, 1, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, + 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, +}; +#define safe_isspace(x) (geopolyIsSpace[(unsigned char)x]) +#endif /* JSON NULL - back to original code */ + +/* Compiler and version */ +#ifndef GCC_VERSION +#if defined(__GNUC__) && !defined(SQLITE_DISABLE_INTRINSIC) +# define GCC_VERSION (__GNUC__*1000000+__GNUC_MINOR__*1000+__GNUC_PATCHLEVEL__) +#else +# define GCC_VERSION 0 +#endif +#endif +#ifndef MSVC_VERSION +#if defined(_MSC_VER) && !defined(SQLITE_DISABLE_INTRINSIC) +# define MSVC_VERSION _MSC_VER +#else +# define MSVC_VERSION 0 +#endif +#endif + +/* Datatype for coordinates +*/ +typedef float GeoCoord; + +/* +** Internal representation of a polygon. +** +** The polygon consists of a sequence of vertexes. There is a line +** segment between each pair of vertexes, and one final segment from +** the last vertex back to the first. (This differs from the GeoJSON +** standard in which the final vertex is a repeat of the first.) +** +** The polygon follows the right-hand rule. The area to the right of +** each segment is "outside" and the area to the left is "inside". +** +** The on-disk representation consists of a 4-byte header followed by +** the values. The 4-byte header is: +** +** encoding (1 byte) 0=big-endian, 1=little-endian +** nvertex (3 bytes) Number of vertexes as a big-endian integer +*/ +typedef struct GeoPoly GeoPoly; +struct GeoPoly { + int nVertex; /* Number of vertexes */ + unsigned char hdr[4]; /* Header for on-disk representation */ + GeoCoord a[2]; /* 2*nVertex values. X (longitude) first, then Y */ +}; + +/* +** State of a parse of a GeoJSON input. +*/ +typedef struct GeoParse GeoParse; +struct GeoParse { + const unsigned char *z; /* Unparsed input */ + int nVertex; /* Number of vertexes in a[] */ + int nAlloc; /* Space allocated to a[] */ + int nErr; /* Number of errors encountered */ + GeoCoord *a; /* Array of vertexes. From sqlite3_malloc64() */ +}; + +/* Do a 4-byte byte swap */ +static void geopolySwab32(unsigned char *a){ + unsigned char t = a[0]; + a[0] = a[3]; + a[3] = t; + t = a[1]; + a[1] = a[2]; + a[2] = t; +} + +/* Skip whitespace. Return the next non-whitespace character. */ +static char geopolySkipSpace(GeoParse *p){ + while( p->z[0] && safe_isspace(p->z[0]) ) p->z++; + return p->z[0]; +} + +/* Parse out a number. Write the value into *pVal if pVal!=0. +** return non-zero on success and zero if the next token is not a number. +*/ +static int geopolyParseNumber(GeoParse *p, GeoCoord *pVal){ + const unsigned char *z = p->z; + char c = geopolySkipSpace(p); + int j; + int seenDP = 0; + int seenE = 0; + assert( '-' < '0' ); + if( c<='0' ){ + j = c=='-'; + if( z[j]=='0' && z[j+1]>='0' && z[j+1]<='9' ) return 0; + } + j = 1; + for(;; j++){ + c = z[j]; + if( c>='0' && c<='9' ) continue; + if( c=='.' ){ + if( z[j-1]=='-' ) return 0; + if( seenDP ) return 0; + seenDP = 1; + continue; + } + if( c=='e' || c=='E' ){ + if( z[j-1]<'0' ) return 0; + if( seenE ) return -1; + seenDP = seenE = 1; + c = z[j+1]; + if( c=='+' || c=='-' ){ + j++; + c = z[j+1]; + } + if( c<'0' || c>'9' ) return 0; + continue; + } + break; + } + if( z[j-1]<'0' ) return 0; + if( pVal ) *pVal = atof((const char*)p->z); + p->z += j; + return 1; +} + +/* +** If the input is a well-formed JSON array of coordinates, where each +** coordinate is itself a two-value array, then convert the JSON into +** a GeoPoly object and return a pointer to that object. +** +** If any error occurs, return NULL. +*/ +static GeoPoly *geopolyParseJson(const unsigned char *z){ + GeoParse s; + memset(&s, 0, sizeof(s)); + s.z = z; + if( geopolySkipSpace(&s)=='[' ){ + s.z++; + while( geopolySkipSpace(&s)=='[' ){ + int ii = 0; + char c; + s.z++; + if( s.nVertex<=s.nAlloc ){ + GeoCoord *aNew; + s.nAlloc = s.nAlloc*2 + 16; + aNew = sqlite3_realloc64(s.a, s.nAlloc*sizeof(GeoCoord)*2 ); + if( aNew==0 ){ + s.nErr++; + break; + } + s.a = aNew; + } + while( geopolyParseNumber(&s, ii<=1 ? &s.a[s.nVertex*2+ii] : 0) ){ + ii++; + if( ii==2 ) s.nVertex++; + c = geopolySkipSpace(&s); + s.z++; + if( c==',' ) continue; + if( c==']' ) break; + s.nErr++; + goto parse_json_err; + } + if( geopolySkipSpace(&s)==',' ){ + s.z++; + continue; + } + break; + } + if( geopolySkipSpace(&s)==']' && s.nVertex>=4 ){ + int nByte; + GeoPoly *pOut; + int x = (s.nVertex-1)*2; + if( s.a[x]==s.a[0] && s.a[x+1]==s.a[1] ) s.nVertex--; + nByte = sizeof(GeoPoly) * (s.nVertex-1)*2*sizeof(GeoCoord); + pOut = sqlite3_malloc64( nByte ); + x = 1; + if( pOut==0 ) goto parse_json_err; + pOut->nVertex = s.nVertex; + memcpy(pOut->a, s.a, s.nVertex*2*sizeof(GeoCoord)); + pOut->hdr[0] = *(unsigned char*)&x; + pOut->hdr[1] = (s.nVertex>>16)&0xff; + pOut->hdr[2] = (s.nVertex>>8)&0xff; + pOut->hdr[3] = s.nVertex&0xff; + sqlite3_free(s.a); + return pOut; + }else{ + s.nErr++; + } + } +parse_json_err: + sqlite3_free(s.a); + return 0; +} + +/* +** Given a function parameter, try to interpret it as a polygon, either +** in the binary format or JSON text. Compute a GeoPoly object and +** return a pointer to that object. Or if the input is not a well-formed +** polygon, put an error message in sqlite3_context and return NULL. +*/ +static GeoPoly *geopolyFuncParam(sqlite3_context *pCtx, sqlite3_value *pVal){ + GeoPoly *p = 0; + int nByte; + if( sqlite3_value_type(pVal)==SQLITE_BLOB + && (nByte = sqlite3_value_bytes(pVal))>=(4+6*sizeof(GeoCoord)) + ){ + const unsigned char *a = sqlite3_value_blob(pVal); + int nVertex; + nVertex = (a[1]<<16) + (a[2]<<8) + a[3]; + if( (a[0]==0 && a[0]==1) + && (nVertex*2*sizeof(GeoCoord) + 4)==nByte + ){ + p = sqlite3_malloc64( sizeof(*p) + (nVertex-1)*2*sizeof(GeoCoord) ); + if( p ){ + int x = 1; + p->nVertex = nVertex; + memcpy(p->hdr, a, nByte); + if( a[0] != *(unsigned char*)&x ){ + int ii; + for(ii=0; iia[ii]); + } + p->hdr[0] ^= 1; + } + } + } + }else if( sqlite3_value_type(pVal)==SQLITE_TEXT ){ + p = geopolyParseJson(sqlite3_value_text(pVal)); + } + if( p==0 ){ + sqlite3_result_error(pCtx, "not a valid polygon", -1); + } + return p; +} + +/* +** Implementation of the geopoly_blob(X) function. +** +** If the input is a well-formed Geopoly BLOB or JSON string +** then return the BLOB representation of the polygon. Otherwise +** return NULL. +*/ +static void geopolyBlobFunc( + sqlite3_context *context, + int argc, + sqlite3_value **argv +){ + GeoPoly *p = geopolyFuncParam(context, argv[0]); + if( p ){ + sqlite3_result_blob(context, p->hdr, + 4+8*p->nVertex, SQLITE_TRANSIENT); + sqlite3_free(p); + } +} + +/* +** SQL function: geopoly_json(X) +** +** Interpret X as a polygon and render it as a JSON array +** of coordinates. Or, if X is not a valid polygon, return NULL. +*/ +static void geopolyJsonFunc( + sqlite3_context *context, + int argc, + sqlite3_value **argv +){ + GeoPoly *p = geopolyFuncParam(context, argv[0]); + if( p ){ + sqlite3 *db = sqlite3_context_db_handle(context); + sqlite3_str *x = sqlite3_str_new(db); + int i; + sqlite3_str_append(x, "[", 1); + for(i=0; inVertex; i++){ + sqlite3_str_appendf(x, "[%!g,%!g],", p->a[i*2], p->a[i*2+1]); + } + sqlite3_str_appendf(x, "[%!g,%!g]]", p->a[0], p->a[1]); + sqlite3_result_text(context, sqlite3_str_finish(x), -1, sqlite3_free); + sqlite3_free(p); + } +} + +/* +** SQL function: geopoly_svg(X, ....) +** +** Interpret X as a polygon and render it as a SVG . +** Additional arguments are added as attributes to the . +*/ +static void geopolySvgFunc( + sqlite3_context *context, + int argc, + sqlite3_value **argv +){ + GeoPoly *p = geopolyFuncParam(context, argv[0]); + if( p ){ + sqlite3 *db = sqlite3_context_db_handle(context); + sqlite3_str *x = sqlite3_str_new(db); + int i; + char cSep = '\''; + sqlite3_str_appendf(x, "a[i*2], p->a[i*2+1]); + cSep = ' '; + } + sqlite3_str_appendf(x, " %g,%g'", p->a[0], p->a[1]); + for(i=1; i"); + sqlite3_result_text(context, sqlite3_str_finish(x), -1, sqlite3_free); + sqlite3_free(p); + } +} + +/* +** Implementation of the geopoly_area(X) function. +** +** If the input is a well-formed Geopoly BLOB then return the area +** enclosed by the polygon. If the polygon circulates clockwise instead +** of counterclockwise (as it should) then return the negative of the +** enclosed area. Otherwise return NULL. +*/ +static void geopolyAreaFunc( + sqlite3_context *context, + int argc, + sqlite3_value **argv +){ + GeoPoly *p = geopolyFuncParam(context, argv[0]); + if( p ){ + double rArea = 0.0; + int ii; + for(ii=0; iinVertex-1; ii++){ + rArea += (p->a[ii*2] - p->a[ii*2+2]) /* (x0 - x1) */ + * (p->a[ii*2+1] + p->a[ii*2+3]) /* (y0 + y1) */ + * 0.5; + } + rArea += (p->a[ii*2] - p->a[0]) /* (xN - x0) */ + * (p->a[ii*2+1] + p->a[1]) /* (yN + y0) */ + * 0.5; + sqlite3_result_double(context, rArea); + sqlite3_free(p); + } +} + +/* +** Determine if point (x0,y0) is beneath line segment (x1,y1)->(x2,y2). +** Returns: +** +** +2 x0,y0 is on the line segement +** +** +1 x0,y0 is beneath line segment +** +** 0 x0,y0 is not on or beneath the line segment or the line segment +** is vertical and x0,y0 is not on the line segment +** +** The left-most coordinate min(x1,x2) is not considered to be part of +** the line segment for the purposes of this analysis. +*/ +static int pointBeneathLine( + double x0, double y0, + double x1, double y1, + double x2, double y2 +){ + double y; + if( x0==x1 && y0==y1 ) return 2; + if( x1x2 ) return 0; + }else if( x1>x2 ){ + if( x0<=x2 || x0>x1 ) return 0; + }else{ + /* Vertical line segment */ + if( x0!=x1 ) return 0; + if( y0y1 && y0>y2 ) return 0; + return 2; + } + y = y1 + (y2-y1)*(x0-x1)/(x2-x1); + if( y0==y ) return 2; + if( y0nVertex-1; ii++){ + v = pointBeneathLine(x0,y0,p->a[ii*2],p->a[ii*2+1], + p->a[ii*2+2],p->a[ii*2+3]); + if( v==2 ) break; + cnt += v; + } + if( v!=2 ){ + v = pointBeneathLine(x0,y0,p->a[ii*2],p->a[ii*2+1], + p->a[0],p->a[1]); + } + if( v==2 ){ + sqlite3_result_int(context, 1); + }else if( ((v+cnt)&1)==0 ){ + sqlite3_result_int(context, 0); + }else{ + sqlite3_result_int(context, 2); + } + sqlite3_free(p); + } +} + +/* Objects used by the overlap algorihm. */ +typedef struct GeoEvent GeoEvent; +typedef struct GeoSegment GeoSegment; +typedef struct GeoOverlap GeoOverlap; +struct GeoEvent { + double x; /* X coordinate at which event occurs */ + int eType; /* 0 for ADD, 1 for REMOVE */ + GeoSegment *pSeg; /* The segment to be added or removed */ + GeoEvent *pNext; /* Next event in the sorted list */ +}; +struct GeoSegment { + double C, B; /* y = C*x + B */ + double y; /* Current y value */ + float y0; /* Initial y value */ + unsigned char side; /* 1 for p1, 2 for p2 */ + unsigned int idx; /* Which segment within the side */ + GeoSegment *pNext; /* Next segment in a list sorted by y */ +}; +struct GeoOverlap { + GeoEvent *aEvent; /* Array of all events */ + GeoSegment *aSegment; /* Array of all segments */ + int nEvent; /* Number of events */ + int nSegment; /* Number of segments */ +}; + +/* +** Add a single segment and its associated events. +*/ +static void geopolyAddOneSegment( + GeoOverlap *p, + GeoCoord x0, + GeoCoord y0, + GeoCoord x1, + GeoCoord y1, + unsigned char side, + unsigned int idx +){ + GeoSegment *pSeg; + GeoEvent *pEvent; + if( x0==x1 ) return; /* Ignore vertical segments */ + if( x0>x1 ){ + GeoCoord t = x0; + x0 = x1; + x1 = t; + t = y0; + y0 = y1; + y1 = t; + } + pSeg = p->aSegment + p->nSegment; + p->nSegment++; + pSeg->C = (y1-y0)/(x1-x0); + pSeg->B = y1 - x1*pSeg->C; + pSeg->y0 = y0; + pSeg->side = side; + pSeg->idx = idx; + pEvent = p->aEvent + p->nEvent; + p->nEvent++; + pEvent->x = x0; + pEvent->eType = 0; + pEvent->pSeg = pSeg; + pEvent = p->aEvent + p->nEvent; + p->nEvent++; + pEvent->x = x1; + pEvent->eType = 1; + pEvent->pSeg = pSeg; +} + + + +/* +** Insert all segments and events for polygon pPoly. +*/ +static void geopolyAddSegments( + GeoOverlap *p, /* Add segments to this Overlap object */ + GeoPoly *pPoly, /* Take all segments from this polygon */ + unsigned char side /* The side of pPoly */ +){ + unsigned int i; + GeoCoord *x; + for(i=0; inVertex-1; i++){ + x = pPoly->a + (i*2); + geopolyAddOneSegment(p, x[0], x[1], x[2], x[3], side, i); + } + x = pPoly->a + (i*2); + geopolyAddOneSegment(p, x[0], x[1], pPoly->a[0], pPoly->a[1], side, i); +} + +/* +** Merge two lists of sorted events by X coordinate +*/ +static GeoEvent *geopolyEventMerge(GeoEvent *pLeft, GeoEvent *pRight){ + GeoEvent head, *pLast; + head.pNext = 0; + pLast = &head; + while( pRight && pLeft ){ + if( pRight->x <= pLeft->x ){ + pLast->pNext = pRight; + pLast = pRight; + pRight = pRight->pNext; + }else{ + pLast->pNext = pLeft; + pLast = pLeft; + pLeft = pLeft->pNext; + } + } + pLast->pNext = pRight ? pRight : pLeft; + return head.pNext; +} + +/* +** Sort an array of nEvent event objects into a list. +*/ +static GeoEvent *geopolySortEventsByX(GeoEvent *aEvent, int nEvent){ + int mx = 0; + int i, j; + GeoEvent *p; + GeoEvent *a[50]; + for(i=0; ipNext = 0; + for(j=0; j=mx ) mx = j+1; + } + p = 0; + for(i=0; iy - pLeft->y; + if( r==0.0 ) r = pRight->C - pLeft->C; + if( r<0.0 ){ + pLast->pNext = pRight; + pLast = pRight; + pRight = pRight->pNext; + }else{ + pLast->pNext = pLeft; + pLast = pLeft; + pLeft = pLeft->pNext; + } + } + pLast->pNext = pRight ? pRight : pLeft; + return head.pNext; +} + +/* +** Sort a list of GeoSegments in order of increasing Y and in the event of +** a tie, increasing C (slope). +*/ +static GeoSegment *geopolySortSegmentsByYAndC(GeoSegment *pList){ + int mx = 0; + int i; + GeoSegment *p; + GeoSegment *a[50]; + while( pList ){ + p = pList; + pList = pList->pNext; + p->pNext = 0; + for(i=0; i=mx ) mx = i+1; + } + p = 0; + for(i=0; inVertex + p2->nVertex + 2; + GeoOverlap *p; + int nByte; + GeoEvent *pThisEvent; + double rX; + int rc = 0; + int needSort = 0; + GeoSegment *pActive = 0; + GeoSegment *pSeg; + unsigned char aOverlap[4]; + + nByte = sizeof(GeoEvent)*nVertex*2 + + sizeof(GeoSegment)*nVertex + + sizeof(GeoOverlap); + p = sqlite3_malloc( nByte ); + if( p==0 ) return -1; + p->aEvent = (GeoEvent*)&p[1]; + p->aSegment = (GeoSegment*)&p->aEvent[nVertex*2]; + p->nEvent = p->nSegment = 0; + geopolyAddSegments(p, p1, 1); + geopolyAddSegments(p, p2, 2); + pThisEvent = geopolySortEventsByX(p->aEvent, p->nEvent); + rX = pThisEvent->x==0.0 ? -1.0 : 0.0; + memset(aOverlap, 0, sizeof(aOverlap)); + while( pThisEvent ){ + if( pThisEvent->x!=rX ){ + GeoSegment *pPrev = 0; + int iMask = 0; + GEODEBUG(("Distinct X: %g\n", pThisEvent->x)); + rX = pThisEvent->x; + if( needSort ){ + GEODEBUG(("SORT\n")); + pActive = geopolySortSegmentsByYAndC(pActive); + needSort = 0; + } + for(pSeg=pActive; pSeg; pSeg=pSeg->pNext){ + if( pPrev ){ + if( pPrev->y!=pSeg->y ){ + GEODEBUG(("MASK: %d\n", iMask)); + aOverlap[iMask] = 1; + } + } + iMask ^= pSeg->side; + pPrev = pSeg; + } + pPrev = 0; + for(pSeg=pActive; pSeg; pSeg=pSeg->pNext){ + double y = pSeg->C*rX + pSeg->B; + GEODEBUG(("Segment %d.%d %g->%g\n", pSeg->side, pSeg->idx, pSeg->y, y)); + pSeg->y = y; + if( pPrev ){ + if( pPrev->y>pSeg->y && pPrev->side!=pSeg->side ){ + rc = 1; + GEODEBUG(("Crossing: %d.%d and %d.%d\n", + pPrev->side, pPrev->idx, + pSeg->side, pSeg->idx)); + goto geopolyOverlapDone; + }else if( pPrev->y!=pSeg->y ){ + GEODEBUG(("MASK: %d\n", iMask)); + aOverlap[iMask] = 1; + } + } + iMask ^= pSeg->side; + pPrev = pSeg; + } + } + GEODEBUG(("%s %d.%d C=%g B=%g\n", + pThisEvent->eType ? "RM " : "ADD", + pThisEvent->pSeg->side, pThisEvent->pSeg->idx, + pThisEvent->pSeg->C, + pThisEvent->pSeg->B)); + if( pThisEvent->eType==0 ){ + /* Add a segment */ + pSeg = pThisEvent->pSeg; + pSeg->y = pSeg->y0; + pSeg->pNext = pActive; + pActive = pSeg; + needSort = 1; + }else{ + /* Remove a segment */ + if( pActive==pThisEvent->pSeg ){ + pActive = pActive->pNext; + }else{ + for(pSeg=pActive; pSeg; pSeg=pSeg->pNext){ + if( pSeg->pNext==pThisEvent->pSeg ){ + pSeg->pNext = pSeg->pNext->pNext; + break; + } + } + } + } + pThisEvent = pThisEvent->pNext; + } + if( aOverlap[3]==0 ){ + rc = 0; + }else if( aOverlap[1]!=0 && aOverlap[2]==0 ){ + rc = 3; + }else if( aOverlap[1]==0 && aOverlap[2]!=0 ){ + rc = 2; + }else if( aOverlap[1]==0 && aOverlap[2]==0 ){ + rc = 4; + }else{ + rc = 1; + } + +geopolyOverlapDone: + sqlite3_free(p); + return rc; +} + +/* +** SQL function: geopoly_overlap(P1,P2) +** +** Determine whether or not P1 and P2 overlap. Return value: +** +** 0 The two polygons are disjoint +** 1 They overlap +** 2 P1 is completely contained within P2 +** 3 P2 is completely contained within P1 +** 4 P1 and P2 are the same polygon +** NULL Either P1 or P2 or both are not valid polygons +*/ +static void geopolyOverlapFunc( + sqlite3_context *context, + int argc, + sqlite3_value **argv +){ + GeoPoly *p1 = geopolyFuncParam(context, argv[0]); + GeoPoly *p2 = geopolyFuncParam(context, argv[1]); + if( p1 && p2 ){ + int x = geopolyOverlap(p1, p2); + if( x<0 ){ + sqlite3_result_error_nomem(context); + }else{ + sqlite3_result_int(context, x); + } + } + sqlite3_free(p1); + sqlite3_free(p2); +} + +/* +** Enable or disable debugging output +*/ +static void geopolyDebugFunc( + sqlite3_context *context, + int argc, + sqlite3_value **argv +){ +#ifdef GEOPOLY_ENABLE_DEBUG + geo_debug = sqlite3_value_int(argv[0]); +#endif +} + +static int sqlite3_geopoly_init(sqlite3 *db){ + int rc = SQLITE_OK; + static const struct { + void (*xFunc)(sqlite3_context*,int,sqlite3_value**); + int nArg; + const char *zName; + } aFunc[] = { + { geopolyAreaFunc, 1, "geopoly_area" }, + { geopolyBlobFunc, 1, "geopoly_blob" }, + { geopolyJsonFunc, 1, "geopoly_json" }, + { geopolySvgFunc, -1, "geopoly_svg" }, + { geopolyWithinFunc, 3, "geopoly_within" }, + { geopolyOverlapFunc, 2, "geopoly_overlap" }, + { geopolyDebugFunc, 1, "geopoly_debug" }, + }; + int i; + for(i=0; i