/* ** 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 ** ** Enough space is allocated for 4 coordinates, to work around over-zealous ** warnings coming from some compiler (notably, clang). In reality, the size ** of each GeoPoly memory allocate is adjusted as necessary so that the ** GeoPoly.a[] array at the end is the appropriate size. */ typedef struct GeoPoly GeoPoly; struct GeoPoly { int nVertex; /* Number of vertexes */ unsigned char hdr[4]; /* Header for on-disk representation */ GeoCoord a[8]; /* 2*nVertex values. X (longitude) first, then Y */ }; /* The size of a memory allocation needed for a GeoPoly object sufficient ** to hold N coordinate pairs. */ #define GEOPOLY_SZ(N) (sizeof(GeoPoly) + sizeof(GeoCoord)*2*((N)-4)) /* Macros to access coordinates of a GeoPoly. ** We have to use these macros, rather than just say p->a[i] in order ** to silence (incorrect) UBSAN warnings if the array index is too large. */ #define GeoX(P,I) (((GeoCoord*)(P)->a)[(I)*2]) #define GeoY(P,I) (((GeoCoord*)(P)->a)[(I)*2+1]) /* ** 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( 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){ char c = geopolySkipSpace(p); const unsigned char *z = p->z; int j = 0; int seenDP = 0; int seenE = 0; if( c=='-' ){ j = 1; c = z[j]; } if( c=='0' && z[j+1]>='0' && z[j+1]<='9' ) return 0; for(;; j++){ c = z[j]; if( safe_isdigit(c) ) 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 ){ #ifdef SQLITE_AMALGAMATION /* The sqlite3AtoF() routine is much much faster than atof(), if it ** is available */ double r; (void)sqlite3AtoF((const char*)p->z, &r, j, SQLITE_UTF8); *pVal = r; #else *pVal = (GeoCoord)atof((const char*)p->z); #endif } p->z += j; return 1; } /* ** If the input is a well-formed JSON array of coordinates with at least ** four coordinates and 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, int *pRc){ GeoParse s; int rc = SQLITE_OK; 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 ){ rc = SQLITE_NOMEM; 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==']' && ii>=2 ) break; s.nErr++; rc = SQLITE_ERROR; goto parse_json_err; } if( geopolySkipSpace(&s)==',' ){ s.z++; continue; } break; } if( geopolySkipSpace(&s)==']' && s.nVertex>=4 && s.a[0]==s.a[s.nVertex*2-2] && s.a[1]==s.a[s.nVertex*2-1] && (s.z++, geopolySkipSpace(&s)==0) ){ GeoPoly *pOut; int x = 1; s.nVertex--; /* Remove the redundant vertex at the end */ pOut = sqlite3_malloc64( GEOPOLY_SZ((sqlite3_int64)s.nVertex) ); 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); if( pRc ) *pRc = SQLITE_OK; return pOut; }else{ s.nErr++; rc = SQLITE_ERROR; } } parse_json_err: if( pRc ) *pRc = rc; 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, /* Context for error messages */ sqlite3_value *pVal, /* The value to decode */ int *pRc /* Write error here */ ){ 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)==(unsigned int)nByte ){ p = sqlite3_malloc64( sizeof(*p) + (nVertex-1)*2*sizeof(GeoCoord) ); if( p==0 ){ if( pRc ) *pRc = SQLITE_NOMEM; if( pCtx ) sqlite3_result_error_nomem(pCtx); }else{ int x = 1; p->nVertex = nVertex; memcpy(p->hdr, a, nByte); if( a[0] != *(unsigned char*)&x ){ int ii; for(ii=0; iihdr[0] ^= 1; } } } if( pRc ) *pRc = SQLITE_OK; return p; }else if( sqlite3_value_type(pVal)==SQLITE_TEXT ){ const unsigned char *zJson = sqlite3_value_text(pVal); if( zJson==0 ){ if( pRc ) *pRc = SQLITE_NOMEM; return 0; } return geopolyParseJson(zJson, pRc); }else{ if( pRc ) *pRc = SQLITE_ERROR; return 0; } } /* ** 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], 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], 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],", GeoX(p,i), GeoY(p,i)); } sqlite3_str_appendf(x, "[%!g,%!g]]", GeoX(p,0), GeoY(p,0)); 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; if( argc<1 ) return; p = geopolyFuncParam(context, argv[0], 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, ""); sqlite3_result_text(context, sqlite3_str_finish(x), -1, sqlite3_free); sqlite3_free(p); } } /* ** SQL Function: geopoly_xform(poly, A, B, C, D, E, F) ** ** Transform and/or translate a polygon as follows: ** ** x1 = A*x0 + B*y0 + E ** y1 = C*x0 + D*y0 + F ** ** For a translation: ** ** geopoly_xform(poly, 1, 0, 0, 1, x-offset, y-offset) ** ** Rotate by R around the point (0,0): ** ** geopoly_xform(poly, cos(R), sin(R), -sin(R), cos(R), 0, 0) */ static void geopolyXformFunc( sqlite3_context *context, int argc, sqlite3_value **argv ){ GeoPoly *p = geopolyFuncParam(context, argv[0], 0); double A = sqlite3_value_double(argv[1]); double B = sqlite3_value_double(argv[2]); double C = sqlite3_value_double(argv[3]); double D = sqlite3_value_double(argv[4]); double E = sqlite3_value_double(argv[5]); double F = sqlite3_value_double(argv[6]); GeoCoord x1, y1, x0, y0; int ii; if( p ){ for(ii=0; iinVertex; ii++){ x0 = GeoX(p,ii); y0 = GeoY(p,ii); x1 = (GeoCoord)(A*x0 + B*y0 + E); y1 = (GeoCoord)(C*x0 + D*y0 + F); GeoX(p,ii) = x1; GeoY(p,ii) = y1; } sqlite3_result_blob(context, p->hdr, 4+8*p->nVertex, SQLITE_TRANSIENT); sqlite3_free(p); } } /* ** Compute the area enclosed by the polygon. ** ** This routine can also be used to detect polygons that rotate in ** the wrong direction. Polygons are suppose to be counter-clockwise (CCW). ** This routine returns a negative value for clockwise (CW) polygons. */ static double geopolyArea(GeoPoly *p){ double rArea = 0.0; int ii; for(ii=0; iinVertex-1; ii++){ rArea += (GeoX(p,ii) - GeoX(p,ii+1)) /* (x0 - x1) */ * (GeoY(p,ii) + GeoY(p,ii+1)) /* (y0 + y1) */ * 0.5; } rArea += (GeoX(p,ii) - GeoX(p,0)) /* (xN - x0) */ * (GeoY(p,ii) + GeoY(p,0)) /* (yN + y0) */ * 0.5; return rArea; } /* ** 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], 0); if( p ){ sqlite3_result_double(context, geopolyArea(p)); sqlite3_free(p); } } /* ** Implementation of the geopoly_ccw(X) function. ** ** If the rotation of polygon X is clockwise (incorrect) instead of ** counter-clockwise (the correct winding order according to RFC7946) ** then reverse the order of the vertexes in polygon X. ** ** In other words, this routine returns a CCW polygon regardless of the ** winding order of its input. ** ** Use this routine to sanitize historical inputs that that sometimes ** contain polygons that wind in the wrong direction. */ static void geopolyCcwFunc( sqlite3_context *context, int argc, sqlite3_value **argv ){ GeoPoly *p = geopolyFuncParam(context, argv[0], 0); if( p ){ if( geopolyArea(p)<0.0 ){ int ii, jj; for(ii=1, jj=p->nVertex-1; iihdr, 4+8*p->nVertex, SQLITE_TRANSIENT); sqlite3_free(p); } } #define GEOPOLY_PI 3.1415926535897932385 /* Fast approximation for sine(X) for X between -0.5*pi and 2*pi */ static double geopolySine(double r){ assert( r>=-0.5*GEOPOLY_PI && r<=2.0*GEOPOLY_PI ); if( r>=1.5*GEOPOLY_PI ){ r -= 2.0*GEOPOLY_PI; } if( r>=0.5*GEOPOLY_PI ){ return -geopolySine(r-GEOPOLY_PI); }else{ double r2 = r*r; double r3 = r2*r; double r5 = r3*r2; return 0.9996949*r - 0.1656700*r3 + 0.0075134*r5; } } /* ** Function: geopoly_regular(X,Y,R,N) ** ** Construct a simple, convex, regular polygon centered at X, Y ** with circumradius R and with N sides. */ static void geopolyRegularFunc( sqlite3_context *context, int argc, sqlite3_value **argv ){ double x = sqlite3_value_double(argv[0]); double y = sqlite3_value_double(argv[1]); double r = sqlite3_value_double(argv[2]); int n = sqlite3_value_int(argv[3]); int i; GeoPoly *p; if( n<3 || r<=0.0 ) return; if( n>1000 ) n = 1000; p = sqlite3_malloc64( sizeof(*p) + (n-1)*2*sizeof(GeoCoord) ); if( p==0 ){ sqlite3_result_error_nomem(context); return; } i = 1; p->hdr[0] = *(unsigned char*)&i; p->hdr[1] = 0; p->hdr[2] = (n>>8)&0xff; p->hdr[3] = n&0xff; for(i=0; ihdr, 4+8*n, SQLITE_TRANSIENT); sqlite3_free(p); } /* ** If pPoly is a polygon, compute its bounding box. Then: ** ** (1) if aCoord!=0 store the bounding box in aCoord, returning NULL ** (2) otherwise, compute a GeoPoly for the bounding box and return the ** new GeoPoly ** ** If pPoly is NULL but aCoord is not NULL, then compute a new GeoPoly from ** the bounding box in aCoord and return a pointer to that GeoPoly. */ static GeoPoly *geopolyBBox( sqlite3_context *context, /* For recording the error */ sqlite3_value *pPoly, /* The polygon */ RtreeCoord *aCoord, /* Results here */ int *pRc /* Error code here */ ){ GeoPoly *pOut = 0; GeoPoly *p; float mnX, mxX, mnY, mxY; if( pPoly==0 && aCoord!=0 ){ p = 0; mnX = aCoord[0].f; mxX = aCoord[1].f; mnY = aCoord[2].f; mxY = aCoord[3].f; goto geopolyBboxFill; }else{ p = geopolyFuncParam(context, pPoly, pRc); } if( p ){ int ii; mnX = mxX = GeoX(p,0); mnY = mxY = GeoY(p,0); for(ii=1; iinVertex; ii++){ double r = GeoX(p,ii); if( rmxX ) mxX = (float)r; r = GeoY(p,ii); if( rmxY ) mxY = (float)r; } if( pRc ) *pRc = SQLITE_OK; if( aCoord==0 ){ geopolyBboxFill: pOut = sqlite3_realloc64(p, GEOPOLY_SZ(4)); if( pOut==0 ){ sqlite3_free(p); if( context ) sqlite3_result_error_nomem(context); if( pRc ) *pRc = SQLITE_NOMEM; return 0; } pOut->nVertex = 4; ii = 1; pOut->hdr[0] = *(unsigned char*)ⅈ pOut->hdr[1] = 0; pOut->hdr[2] = 0; pOut->hdr[3] = 4; GeoX(pOut,0) = mnX; GeoY(pOut,0) = mnY; GeoX(pOut,1) = mxX; GeoY(pOut,1) = mnY; GeoX(pOut,2) = mxX; GeoY(pOut,2) = mxY; GeoX(pOut,3) = mnX; GeoY(pOut,3) = mxY; }else{ sqlite3_free(p); aCoord[0].f = mnX; aCoord[1].f = mxX; aCoord[2].f = mnY; aCoord[3].f = mxY; } }else{ memset(aCoord, 0, sizeof(RtreeCoord)*4); } return pOut; } /* ** Implementation of the geopoly_bbox(X) SQL function. */ static void geopolyBBoxFunc( sqlite3_context *context, int argc, sqlite3_value **argv ){ GeoPoly *p = geopolyBBox(context, argv[0], 0, 0); if( p ){ sqlite3_result_blob(context, p->hdr, 4+8*p->nVertex, SQLITE_TRANSIENT); sqlite3_free(p); } } /* ** State vector for the geopoly_group_bbox() aggregate function. */ typedef struct GeoBBox GeoBBox; struct GeoBBox { int isInit; RtreeCoord a[4]; }; /* ** Implementation of the geopoly_group_bbox(X) aggregate SQL function. */ static void geopolyBBoxStep( sqlite3_context *context, int argc, sqlite3_value **argv ){ RtreeCoord a[4]; int rc = SQLITE_OK; (void)geopolyBBox(context, argv[0], a, &rc); if( rc==SQLITE_OK ){ GeoBBox *pBBox; pBBox = (GeoBBox*)sqlite3_aggregate_context(context, sizeof(*pBBox)); if( pBBox==0 ) return; if( pBBox->isInit==0 ){ pBBox->isInit = 1; memcpy(pBBox->a, a, sizeof(RtreeCoord)*4); }else{ if( a[0].f < pBBox->a[0].f ) pBBox->a[0] = a[0]; if( a[1].f > pBBox->a[1].f ) pBBox->a[1] = a[1]; if( a[2].f < pBBox->a[2].f ) pBBox->a[2] = a[2]; if( a[3].f > pBBox->a[3].f ) pBBox->a[3] = a[3]; } } } static void geopolyBBoxFinal( sqlite3_context *context ){ GeoPoly *p; GeoBBox *pBBox; pBBox = (GeoBBox*)sqlite3_aggregate_context(context, 0); if( pBBox==0 ) return; p = geopolyBBox(context, 0, pBBox->a, 0); if( p ){ sqlite3_result_blob(context, p->hdr, 4+8*p->nVertex, SQLITE_TRANSIENT); 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,GeoX(p1,ii), GeoY(p1,ii), GeoX(p1,ii+1),GeoY(p1,ii+1)); if( v==2 ) break; cnt += v; } if( v!=2 ){ v = pointBeneathLine(x0,y0,GeoX(p1,ii), GeoY(p1,ii), GeoX(p1,0), GeoY(p1,0)); } 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(p1); } /* Forward declaration */ static int geopolyOverlap(GeoPoly *p1, GeoPoly *p2); /* ** SQL function: geopoly_within(P1,P2) ** ** Return +2 if P1 and P2 are the same polygon ** Return +1 if P2 is contained within P1 ** Return 0 if any part of P2 is on the outside of P1 ** */ static void geopolyWithinFunc( sqlite3_context *context, int argc, sqlite3_value **argv ){ GeoPoly *p1 = geopolyFuncParam(context, argv[0], 0); GeoPoly *p2 = geopolyFuncParam(context, argv[1], 0); if( p1 && p2 ){ int x = geopolyOverlap(p1, p2); if( x<0 ){ sqlite3_result_error_nomem(context); }else{ sqlite3_result_int(context, x==2 ? 1 : x==4 ? 2 : 0); } } sqlite3_free(p1); sqlite3_free(p2); } /* 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; i<(unsigned)pPoly->nVertex-1; i++){ x = &GeoX(pPoly,i); geopolyAddOneSegment(p, x[0], x[1], x[2], x[3], side, i); } x = &GeoX(pPoly,i); 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; sqlite3_int64 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_malloc64( 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], 0); GeoPoly *p2 = geopolyFuncParam(context, argv[1], 0); 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 } /* ** This function is the implementation of both the xConnect and xCreate ** methods of the geopoly virtual table. ** ** argv[0] -> module name ** argv[1] -> database name ** argv[2] -> table name ** argv[...] -> column names... */ static int geopolyInit( sqlite3 *db, /* Database connection */ void *pAux, /* One of the RTREE_COORD_* constants */ int argc, const char *const*argv, /* Parameters to CREATE TABLE statement */ sqlite3_vtab **ppVtab, /* OUT: New virtual table */ char **pzErr, /* OUT: Error message, if any */ int isCreate /* True for xCreate, false for xConnect */ ){ int rc = SQLITE_OK; Rtree *pRtree; sqlite3_int64 nDb; /* Length of string argv[1] */ sqlite3_int64 nName; /* Length of string argv[2] */ sqlite3_str *pSql; char *zSql; int ii; sqlite3_vtab_config(db, SQLITE_VTAB_CONSTRAINT_SUPPORT, 1); /* Allocate the sqlite3_vtab structure */ nDb = strlen(argv[1]); nName = strlen(argv[2]); pRtree = (Rtree *)sqlite3_malloc64(sizeof(Rtree)+nDb+nName+2); if( !pRtree ){ return SQLITE_NOMEM; } memset(pRtree, 0, sizeof(Rtree)+nDb+nName+2); pRtree->nBusy = 1; pRtree->base.pModule = &rtreeModule; pRtree->zDb = (char *)&pRtree[1]; pRtree->zName = &pRtree->zDb[nDb+1]; pRtree->eCoordType = RTREE_COORD_REAL32; pRtree->nDim = 2; pRtree->nDim2 = 4; memcpy(pRtree->zDb, argv[1], nDb); memcpy(pRtree->zName, argv[2], nName); /* Create/Connect to the underlying relational database schema. If ** that is successful, call sqlite3_declare_vtab() to configure ** the r-tree table schema. */ pSql = sqlite3_str_new(db); sqlite3_str_appendf(pSql, "CREATE TABLE x(_shape"); pRtree->nAux = 1; /* Add one for _shape */ pRtree->nAuxNotNull = 1; /* The _shape column is always not-null */ for(ii=3; iinAux++; sqlite3_str_appendf(pSql, ",%s", argv[ii]); } sqlite3_str_appendf(pSql, ");"); zSql = sqlite3_str_finish(pSql); if( !zSql ){ rc = SQLITE_NOMEM; }else if( SQLITE_OK!=(rc = sqlite3_declare_vtab(db, zSql)) ){ *pzErr = sqlite3_mprintf("%s", sqlite3_errmsg(db)); } sqlite3_free(zSql); if( rc ) goto geopolyInit_fail; pRtree->nBytesPerCell = 8 + pRtree->nDim2*4; /* Figure out the node size to use. */ rc = getNodeSize(db, pRtree, isCreate, pzErr); if( rc ) goto geopolyInit_fail; rc = rtreeSqlInit(pRtree, db, argv[1], argv[2], isCreate); if( rc ){ *pzErr = sqlite3_mprintf("%s", sqlite3_errmsg(db)); goto geopolyInit_fail; } *ppVtab = (sqlite3_vtab *)pRtree; return SQLITE_OK; geopolyInit_fail: if( rc==SQLITE_OK ) rc = SQLITE_ERROR; assert( *ppVtab==0 ); assert( pRtree->nBusy==1 ); rtreeRelease(pRtree); return rc; } /* ** GEOPOLY virtual table module xCreate method. */ static int geopolyCreate( sqlite3 *db, void *pAux, int argc, const char *const*argv, sqlite3_vtab **ppVtab, char **pzErr ){ return geopolyInit(db, pAux, argc, argv, ppVtab, pzErr, 1); } /* ** GEOPOLY virtual table module xConnect method. */ static int geopolyConnect( sqlite3 *db, void *pAux, int argc, const char *const*argv, sqlite3_vtab **ppVtab, char **pzErr ){ return geopolyInit(db, pAux, argc, argv, ppVtab, pzErr, 0); } /* ** GEOPOLY virtual table module xFilter method. ** ** Query plans: ** ** 1 rowid lookup ** 2 search for objects overlapping the same bounding box ** that contains polygon argv[0] ** 3 search for objects overlapping the same bounding box ** that contains polygon argv[0] ** 4 full table scan */ static int geopolyFilter( sqlite3_vtab_cursor *pVtabCursor, /* The cursor to initialize */ int idxNum, /* Query plan */ const char *idxStr, /* Not Used */ int argc, sqlite3_value **argv /* Parameters to the query plan */ ){ Rtree *pRtree = (Rtree *)pVtabCursor->pVtab; RtreeCursor *pCsr = (RtreeCursor *)pVtabCursor; RtreeNode *pRoot = 0; int rc = SQLITE_OK; int iCell = 0; rtreeReference(pRtree); /* Reset the cursor to the same state as rtreeOpen() leaves it in. */ resetCursor(pCsr); pCsr->iStrategy = idxNum; if( idxNum==1 ){ /* Special case - lookup by rowid. */ RtreeNode *pLeaf; /* Leaf on which the required cell resides */ RtreeSearchPoint *p; /* Search point for the leaf */ i64 iRowid = sqlite3_value_int64(argv[0]); i64 iNode = 0; rc = findLeafNode(pRtree, iRowid, &pLeaf, &iNode); if( rc==SQLITE_OK && pLeaf!=0 ){ p = rtreeSearchPointNew(pCsr, RTREE_ZERO, 0); assert( p!=0 ); /* Always returns pCsr->sPoint */ pCsr->aNode[0] = pLeaf; p->id = iNode; p->eWithin = PARTLY_WITHIN; rc = nodeRowidIndex(pRtree, pLeaf, iRowid, &iCell); p->iCell = (u8)iCell; RTREE_QUEUE_TRACE(pCsr, "PUSH-F1:"); }else{ pCsr->atEOF = 1; } }else{ /* Normal case - r-tree scan. Set up the RtreeCursor.aConstraint array ** with the configured constraints. */ rc = nodeAcquire(pRtree, 1, 0, &pRoot); if( rc==SQLITE_OK && idxNum<=3 ){ RtreeCoord bbox[4]; RtreeConstraint *p; assert( argc==1 ); geopolyBBox(0, argv[0], bbox, &rc); if( rc ){ goto geopoly_filter_end; } pCsr->aConstraint = p = sqlite3_malloc(sizeof(RtreeConstraint)*4); pCsr->nConstraint = 4; if( p==0 ){ rc = SQLITE_NOMEM; }else{ memset(pCsr->aConstraint, 0, sizeof(RtreeConstraint)*4); memset(pCsr->anQueue, 0, sizeof(u32)*(pRtree->iDepth + 1)); if( idxNum==2 ){ /* Overlap query */ p->op = 'B'; p->iCoord = 0; p->u.rValue = bbox[1].f; p++; p->op = 'D'; p->iCoord = 1; p->u.rValue = bbox[0].f; p++; p->op = 'B'; p->iCoord = 2; p->u.rValue = bbox[3].f; p++; p->op = 'D'; p->iCoord = 3; p->u.rValue = bbox[2].f; }else{ /* Within query */ p->op = 'D'; p->iCoord = 0; p->u.rValue = bbox[0].f; p++; p->op = 'B'; p->iCoord = 1; p->u.rValue = bbox[1].f; p++; p->op = 'D'; p->iCoord = 2; p->u.rValue = bbox[2].f; p++; p->op = 'B'; p->iCoord = 3; p->u.rValue = bbox[3].f; } } } if( rc==SQLITE_OK ){ RtreeSearchPoint *pNew; pNew = rtreeSearchPointNew(pCsr, RTREE_ZERO, (u8)(pRtree->iDepth+1)); if( pNew==0 ){ rc = SQLITE_NOMEM; goto geopoly_filter_end; } pNew->id = 1; pNew->iCell = 0; pNew->eWithin = PARTLY_WITHIN; assert( pCsr->bPoint==1 ); pCsr->aNode[0] = pRoot; pRoot = 0; RTREE_QUEUE_TRACE(pCsr, "PUSH-Fm:"); rc = rtreeStepToLeaf(pCsr); } } geopoly_filter_end: nodeRelease(pRtree, pRoot); rtreeRelease(pRtree); return rc; } /* ** Rtree virtual table module xBestIndex method. There are three ** table scan strategies to choose from (in order from most to ** least desirable): ** ** idxNum idxStr Strategy ** ------------------------------------------------ ** 1 "rowid" Direct lookup by rowid. ** 2 "rtree" R-tree overlap query using geopoly_overlap() ** 3 "rtree" R-tree within query using geopoly_within() ** 4 "fullscan" full-table scan. ** ------------------------------------------------ */ static int geopolyBestIndex(sqlite3_vtab *tab, sqlite3_index_info *pIdxInfo){ int ii; int iRowidTerm = -1; int iFuncTerm = -1; int idxNum = 0; for(ii=0; iinConstraint; ii++){ struct sqlite3_index_constraint *p = &pIdxInfo->aConstraint[ii]; if( !p->usable ) continue; if( p->iColumn<0 && p->op==SQLITE_INDEX_CONSTRAINT_EQ ){ iRowidTerm = ii; break; } if( p->iColumn==0 && p->op>=SQLITE_INDEX_CONSTRAINT_FUNCTION ){ /* p->op==SQLITE_INDEX_CONSTRAINT_FUNCTION for geopoly_overlap() ** p->op==(SQLITE_INDEX_CONTRAINT_FUNCTION+1) for geopoly_within(). ** See geopolyFindFunction() */ iFuncTerm = ii; idxNum = p->op - SQLITE_INDEX_CONSTRAINT_FUNCTION + 2; } } if( iRowidTerm>=0 ){ pIdxInfo->idxNum = 1; pIdxInfo->idxStr = "rowid"; pIdxInfo->aConstraintUsage[iRowidTerm].argvIndex = 1; pIdxInfo->aConstraintUsage[iRowidTerm].omit = 1; pIdxInfo->estimatedCost = 30.0; pIdxInfo->estimatedRows = 1; pIdxInfo->idxFlags = SQLITE_INDEX_SCAN_UNIQUE; return SQLITE_OK; } if( iFuncTerm>=0 ){ pIdxInfo->idxNum = idxNum; pIdxInfo->idxStr = "rtree"; pIdxInfo->aConstraintUsage[iFuncTerm].argvIndex = 1; pIdxInfo->aConstraintUsage[iFuncTerm].omit = 0; pIdxInfo->estimatedCost = 300.0; pIdxInfo->estimatedRows = 10; return SQLITE_OK; } pIdxInfo->idxNum = 4; pIdxInfo->idxStr = "fullscan"; pIdxInfo->estimatedCost = 3000000.0; pIdxInfo->estimatedRows = 100000; return SQLITE_OK; } /* ** GEOPOLY virtual table module xColumn method. */ static int geopolyColumn(sqlite3_vtab_cursor *cur, sqlite3_context *ctx, int i){ Rtree *pRtree = (Rtree *)cur->pVtab; RtreeCursor *pCsr = (RtreeCursor *)cur; RtreeSearchPoint *p = rtreeSearchPointFirst(pCsr); int rc = SQLITE_OK; RtreeNode *pNode = rtreeNodeOfFirstSearchPoint(pCsr, &rc); if( rc ) return rc; if( p==0 ) return SQLITE_OK; if( i==0 && sqlite3_vtab_nochange(ctx) ) return SQLITE_OK; if( i<=pRtree->nAux ){ if( !pCsr->bAuxValid ){ if( pCsr->pReadAux==0 ){ rc = sqlite3_prepare_v3(pRtree->db, pRtree->zReadAuxSql, -1, 0, &pCsr->pReadAux, 0); if( rc ) return rc; } sqlite3_bind_int64(pCsr->pReadAux, 1, nodeGetRowid(pRtree, pNode, p->iCell)); rc = sqlite3_step(pCsr->pReadAux); if( rc==SQLITE_ROW ){ pCsr->bAuxValid = 1; }else{ sqlite3_reset(pCsr->pReadAux); if( rc==SQLITE_DONE ) rc = SQLITE_OK; return rc; } } sqlite3_result_value(ctx, sqlite3_column_value(pCsr->pReadAux, i+2)); } return SQLITE_OK; } /* ** The xUpdate method for GEOPOLY module virtual tables. ** ** For DELETE: ** ** argv[0] = the rowid to be deleted ** ** For INSERT: ** ** argv[0] = SQL NULL ** argv[1] = rowid to insert, or an SQL NULL to select automatically ** argv[2] = _shape column ** argv[3] = first application-defined column.... ** ** For UPDATE: ** ** argv[0] = rowid to modify. Never NULL ** argv[1] = rowid after the change. Never NULL ** argv[2] = new value for _shape ** argv[3] = new value for first application-defined column.... */ static int geopolyUpdate( sqlite3_vtab *pVtab, int nData, sqlite3_value **aData, sqlite_int64 *pRowid ){ Rtree *pRtree = (Rtree *)pVtab; int rc = SQLITE_OK; RtreeCell cell; /* New cell to insert if nData>1 */ i64 oldRowid; /* The old rowid */ int oldRowidValid; /* True if oldRowid is valid */ i64 newRowid; /* The new rowid */ int newRowidValid; /* True if newRowid is valid */ int coordChange = 0; /* Change in coordinates */ if( pRtree->nNodeRef ){ /* Unable to write to the btree while another cursor is reading from it, ** since the write might do a rebalance which would disrupt the read ** cursor. */ return SQLITE_LOCKED_VTAB; } rtreeReference(pRtree); assert(nData>=1); oldRowidValid = sqlite3_value_type(aData[0])!=SQLITE_NULL;; oldRowid = oldRowidValid ? sqlite3_value_int64(aData[0]) : 0; newRowidValid = nData>1 && sqlite3_value_type(aData[1])!=SQLITE_NULL; newRowid = newRowidValid ? sqlite3_value_int64(aData[1]) : 0; cell.iRowid = newRowid; if( nData>1 /* not a DELETE */ && (!oldRowidValid /* INSERT */ || !sqlite3_value_nochange(aData[2]) /* UPDATE _shape */ || oldRowid!=newRowid) /* Rowid change */ ){ geopolyBBox(0, aData[2], cell.aCoord, &rc); if( rc ){ if( rc==SQLITE_ERROR ){ pVtab->zErrMsg = sqlite3_mprintf("_shape does not contain a valid polygon"); } goto geopoly_update_end; } coordChange = 1; /* If a rowid value was supplied, check if it is already present in ** the table. If so, the constraint has failed. */ if( newRowidValid && (!oldRowidValid || oldRowid!=newRowid) ){ int steprc; sqlite3_bind_int64(pRtree->pReadRowid, 1, cell.iRowid); steprc = sqlite3_step(pRtree->pReadRowid); rc = sqlite3_reset(pRtree->pReadRowid); if( SQLITE_ROW==steprc ){ if( sqlite3_vtab_on_conflict(pRtree->db)==SQLITE_REPLACE ){ rc = rtreeDeleteRowid(pRtree, cell.iRowid); }else{ rc = rtreeConstraintError(pRtree, 0); } } } } /* If aData[0] is not an SQL NULL value, it is the rowid of a ** record to delete from the r-tree table. The following block does ** just that. */ if( rc==SQLITE_OK && (nData==1 || (coordChange && oldRowidValid)) ){ rc = rtreeDeleteRowid(pRtree, oldRowid); } /* If the aData[] array contains more than one element, elements ** (aData[2]..aData[argc-1]) contain a new record to insert into ** the r-tree structure. */ if( rc==SQLITE_OK && nData>1 && coordChange ){ /* Insert the new record into the r-tree */ RtreeNode *pLeaf = 0; if( !newRowidValid ){ rc = rtreeNewRowid(pRtree, &cell.iRowid); } *pRowid = cell.iRowid; if( rc==SQLITE_OK ){ rc = ChooseLeaf(pRtree, &cell, 0, &pLeaf); } if( rc==SQLITE_OK ){ int rc2; pRtree->iReinsertHeight = -1; rc = rtreeInsertCell(pRtree, pLeaf, &cell, 0); rc2 = nodeRelease(pRtree, pLeaf); if( rc==SQLITE_OK ){ rc = rc2; } } } /* Change the data */ if( rc==SQLITE_OK && nData>1 ){ sqlite3_stmt *pUp = pRtree->pWriteAux; int jj; int nChange = 0; sqlite3_bind_int64(pUp, 1, cell.iRowid); assert( pRtree->nAux>=1 ); if( sqlite3_value_nochange(aData[2]) ){ sqlite3_bind_null(pUp, 2); }else{ GeoPoly *p = 0; if( sqlite3_value_type(aData[2])==SQLITE_TEXT && (p = geopolyFuncParam(0, aData[2], &rc))!=0 && rc==SQLITE_OK ){ sqlite3_bind_blob(pUp, 2, p->hdr, 4+8*p->nVertex, SQLITE_TRANSIENT); }else{ sqlite3_bind_value(pUp, 2, aData[2]); } sqlite3_free(p); nChange = 1; } for(jj=1; jjnAux; jj++){ nChange++; sqlite3_bind_value(pUp, jj+2, aData[jj+2]); } if( nChange ){ sqlite3_step(pUp); rc = sqlite3_reset(pUp); } } geopoly_update_end: rtreeRelease(pRtree); return rc; } /* ** Report that geopoly_overlap() is an overloaded function suitable ** for use in xBestIndex. */ static int geopolyFindFunction( sqlite3_vtab *pVtab, int nArg, const char *zName, void (**pxFunc)(sqlite3_context*,int,sqlite3_value**), void **ppArg ){ if( sqlite3_stricmp(zName, "geopoly_overlap")==0 ){ *pxFunc = geopolyOverlapFunc; *ppArg = 0; return SQLITE_INDEX_CONSTRAINT_FUNCTION; } if( sqlite3_stricmp(zName, "geopoly_within")==0 ){ *pxFunc = geopolyWithinFunc; *ppArg = 0; return SQLITE_INDEX_CONSTRAINT_FUNCTION+1; } return 0; } static sqlite3_module geopolyModule = { 3, /* iVersion */ geopolyCreate, /* xCreate - create a table */ geopolyConnect, /* xConnect - connect to an existing table */ geopolyBestIndex, /* xBestIndex - Determine search strategy */ rtreeDisconnect, /* xDisconnect - Disconnect from a table */ rtreeDestroy, /* xDestroy - Drop a table */ rtreeOpen, /* xOpen - open a cursor */ rtreeClose, /* xClose - close a cursor */ geopolyFilter, /* xFilter - configure scan constraints */ rtreeNext, /* xNext - advance a cursor */ rtreeEof, /* xEof */ geopolyColumn, /* xColumn - read data */ rtreeRowid, /* xRowid - read data */ geopolyUpdate, /* xUpdate - write data */ rtreeBeginTransaction, /* xBegin - begin transaction */ rtreeEndTransaction, /* xSync - sync transaction */ rtreeEndTransaction, /* xCommit - commit transaction */ rtreeEndTransaction, /* xRollback - rollback transaction */ geopolyFindFunction, /* xFindFunction - function overloading */ rtreeRename, /* xRename - rename the table */ rtreeSavepoint, /* xSavepoint */ 0, /* xRelease */ 0, /* xRollbackTo */ rtreeShadowName /* xShadowName */ }; static int sqlite3_geopoly_init(sqlite3 *db){ int rc = SQLITE_OK; static const struct { void (*xFunc)(sqlite3_context*,int,sqlite3_value**); signed char nArg; unsigned char bPure; const char *zName; } aFunc[] = { { geopolyAreaFunc, 1, 1, "geopoly_area" }, { geopolyBlobFunc, 1, 1, "geopoly_blob" }, { geopolyJsonFunc, 1, 1, "geopoly_json" }, { geopolySvgFunc, -1, 1, "geopoly_svg" }, { geopolyWithinFunc, 2, 1, "geopoly_within" }, { geopolyContainsPointFunc, 3, 1, "geopoly_contains_point" }, { geopolyOverlapFunc, 2, 1, "geopoly_overlap" }, { geopolyDebugFunc, 1, 0, "geopoly_debug" }, { geopolyBBoxFunc, 1, 1, "geopoly_bbox" }, { geopolyXformFunc, 7, 1, "geopoly_xform" }, { geopolyRegularFunc, 4, 1, "geopoly_regular" }, { geopolyCcwFunc, 1, 1, "geopoly_ccw" }, }; static const struct { void (*xStep)(sqlite3_context*,int,sqlite3_value**); void (*xFinal)(sqlite3_context*); const char *zName; } aAgg[] = { { geopolyBBoxStep, geopolyBBoxFinal, "geopoly_group_bbox" }, }; int i; for(i=0; i