1 /*
2 ** 2018-05-25
3 **
4 ** The author disclaims copyright to this source code. In place of
5 ** a legal notice, here is a blessing:
6 **
7 ** May you do good and not evil.
8 ** May you find forgiveness for yourself and forgive others.
9 ** May you share freely, never taking more than you give.
10 **
11 ******************************************************************************
12 **
13 ** This file implements an alternative R-Tree virtual table that
14 ** uses polygons to express the boundaries of 2-dimensional objects.
15 **
16 ** This file is #include-ed onto the end of "rtree.c" so that it has
17 ** access to all of the R-Tree internals.
18 */
19 #include <stdlib.h>
20
21 /* Enable -DGEOPOLY_ENABLE_DEBUG for debugging facilities */
22 #ifdef GEOPOLY_ENABLE_DEBUG
23 static int geo_debug = 0;
24 # define GEODEBUG(X) if(geo_debug)printf X
25 #else
26 # define GEODEBUG(X)
27 #endif
28
29 #ifndef JSON_NULL /* The following stuff repeats things found in json1 */
30 /*
31 ** Versions of isspace(), isalnum() and isdigit() to which it is safe
32 ** to pass signed char values.
33 */
34 #ifdef sqlite3Isdigit
35 /* Use the SQLite core versions if this routine is part of the
36 ** SQLite amalgamation */
37 # define safe_isdigit(x) sqlite3Isdigit(x)
38 # define safe_isalnum(x) sqlite3Isalnum(x)
39 # define safe_isxdigit(x) sqlite3Isxdigit(x)
40 #else
41 /* Use the standard library for separate compilation */
42 #include <ctype.h> /* amalgamator: keep */
43 # define safe_isdigit(x) isdigit((unsigned char)(x))
44 # define safe_isalnum(x) isalnum((unsigned char)(x))
45 # define safe_isxdigit(x) isxdigit((unsigned char)(x))
46 #endif
47
48 /*
49 ** Growing our own isspace() routine this way is twice as fast as
50 ** the library isspace() function.
51 */
52 static const char geopolyIsSpace[] = {
53 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 0, 0, 1, 0, 0,
54 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
55 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
56 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
57 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
58 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
59 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
60 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
61 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
62 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
63 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
64 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
65 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
66 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
67 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
68 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0,
69 };
70 #define safe_isspace(x) (geopolyIsSpace[(unsigned char)x])
71 #endif /* JSON NULL - back to original code */
72
73 /* Compiler and version */
74 #ifndef GCC_VERSION
75 #if defined(__GNUC__) && !defined(SQLITE_DISABLE_INTRINSIC)
76 # define GCC_VERSION (__GNUC__*1000000+__GNUC_MINOR__*1000+__GNUC_PATCHLEVEL__)
77 #else
78 # define GCC_VERSION 0
79 #endif
80 #endif
81 #ifndef MSVC_VERSION
82 #if defined(_MSC_VER) && !defined(SQLITE_DISABLE_INTRINSIC)
83 # define MSVC_VERSION _MSC_VER
84 #else
85 # define MSVC_VERSION 0
86 #endif
87 #endif
88
89 /* Datatype for coordinates
90 */
91 typedef float GeoCoord;
92
93 /*
94 ** Internal representation of a polygon.
95 **
96 ** The polygon consists of a sequence of vertexes. There is a line
97 ** segment between each pair of vertexes, and one final segment from
98 ** the last vertex back to the first. (This differs from the GeoJSON
99 ** standard in which the final vertex is a repeat of the first.)
100 **
101 ** The polygon follows the right-hand rule. The area to the right of
102 ** each segment is "outside" and the area to the left is "inside".
103 **
104 ** The on-disk representation consists of a 4-byte header followed by
105 ** the values. The 4-byte header is:
106 **
107 ** encoding (1 byte) 0=big-endian, 1=little-endian
108 ** nvertex (3 bytes) Number of vertexes as a big-endian integer
109 **
110 ** Enough space is allocated for 4 coordinates, to work around over-zealous
111 ** warnings coming from some compiler (notably, clang). In reality, the size
112 ** of each GeoPoly memory allocate is adjusted as necessary so that the
113 ** GeoPoly.a[] array at the end is the appropriate size.
114 */
115 typedef struct GeoPoly GeoPoly;
116 struct GeoPoly {
117 int nVertex; /* Number of vertexes */
118 unsigned char hdr[4]; /* Header for on-disk representation */
119 GeoCoord a[8]; /* 2*nVertex values. X (longitude) first, then Y */
120 };
121
122 /* The size of a memory allocation needed for a GeoPoly object sufficient
123 ** to hold N coordinate pairs.
124 */
125 #define GEOPOLY_SZ(N) (sizeof(GeoPoly) + sizeof(GeoCoord)*2*((N)-4))
126
127 /* Macros to access coordinates of a GeoPoly.
128 ** We have to use these macros, rather than just say p->a[i] in order
129 ** to silence (incorrect) UBSAN warnings if the array index is too large.
130 */
131 #define GeoX(P,I) (((GeoCoord*)(P)->a)[(I)*2])
132 #define GeoY(P,I) (((GeoCoord*)(P)->a)[(I)*2+1])
133
134
135 /*
136 ** State of a parse of a GeoJSON input.
137 */
138 typedef struct GeoParse GeoParse;
139 struct GeoParse {
140 const unsigned char *z; /* Unparsed input */
141 int nVertex; /* Number of vertexes in a[] */
142 int nAlloc; /* Space allocated to a[] */
143 int nErr; /* Number of errors encountered */
144 GeoCoord *a; /* Array of vertexes. From sqlite3_malloc64() */
145 };
146
147 /* Do a 4-byte byte swap */
geopolySwab32(unsigned char * a)148 static void geopolySwab32(unsigned char *a){
149 unsigned char t = a[0];
150 a[0] = a[3];
151 a[3] = t;
152 t = a[1];
153 a[1] = a[2];
154 a[2] = t;
155 }
156
157 /* Skip whitespace. Return the next non-whitespace character. */
geopolySkipSpace(GeoParse * p)158 static char geopolySkipSpace(GeoParse *p){
159 while( safe_isspace(p->z[0]) ) p->z++;
160 return p->z[0];
161 }
162
163 /* Parse out a number. Write the value into *pVal if pVal!=0.
164 ** return non-zero on success and zero if the next token is not a number.
165 */
geopolyParseNumber(GeoParse * p,GeoCoord * pVal)166 static int geopolyParseNumber(GeoParse *p, GeoCoord *pVal){
167 char c = geopolySkipSpace(p);
168 const unsigned char *z = p->z;
169 int j = 0;
170 int seenDP = 0;
171 int seenE = 0;
172 if( c=='-' ){
173 j = 1;
174 c = z[j];
175 }
176 if( c=='0' && z[j+1]>='0' && z[j+1]<='9' ) return 0;
177 for(;; j++){
178 c = z[j];
179 if( safe_isdigit(c) ) continue;
180 if( c=='.' ){
181 if( z[j-1]=='-' ) return 0;
182 if( seenDP ) return 0;
183 seenDP = 1;
184 continue;
185 }
186 if( c=='e' || c=='E' ){
187 if( z[j-1]<'0' ) return 0;
188 if( seenE ) return -1;
189 seenDP = seenE = 1;
190 c = z[j+1];
191 if( c=='+' || c=='-' ){
192 j++;
193 c = z[j+1];
194 }
195 if( c<'0' || c>'9' ) return 0;
196 continue;
197 }
198 break;
199 }
200 if( z[j-1]<'0' ) return 0;
201 if( pVal ){
202 #ifdef SQLITE_AMALGAMATION
203 /* The sqlite3AtoF() routine is much much faster than atof(), if it
204 ** is available */
205 double r;
206 (void)sqlite3AtoF((const char*)p->z, &r, j, SQLITE_UTF8);
207 *pVal = r;
208 #else
209 *pVal = (GeoCoord)atof((const char*)p->z);
210 #endif
211 }
212 p->z += j;
213 return 1;
214 }
215
216 /*
217 ** If the input is a well-formed JSON array of coordinates with at least
218 ** four coordinates and where each coordinate is itself a two-value array,
219 ** then convert the JSON into a GeoPoly object and return a pointer to
220 ** that object.
221 **
222 ** If any error occurs, return NULL.
223 */
geopolyParseJson(const unsigned char * z,int * pRc)224 static GeoPoly *geopolyParseJson(const unsigned char *z, int *pRc){
225 GeoParse s;
226 int rc = SQLITE_OK;
227 memset(&s, 0, sizeof(s));
228 s.z = z;
229 if( geopolySkipSpace(&s)=='[' ){
230 s.z++;
231 while( geopolySkipSpace(&s)=='[' ){
232 int ii = 0;
233 char c;
234 s.z++;
235 if( s.nVertex>=s.nAlloc ){
236 GeoCoord *aNew;
237 s.nAlloc = s.nAlloc*2 + 16;
238 aNew = sqlite3_realloc64(s.a, s.nAlloc*sizeof(GeoCoord)*2 );
239 if( aNew==0 ){
240 rc = SQLITE_NOMEM;
241 s.nErr++;
242 break;
243 }
244 s.a = aNew;
245 }
246 while( geopolyParseNumber(&s, ii<=1 ? &s.a[s.nVertex*2+ii] : 0) ){
247 ii++;
248 if( ii==2 ) s.nVertex++;
249 c = geopolySkipSpace(&s);
250 s.z++;
251 if( c==',' ) continue;
252 if( c==']' && ii>=2 ) break;
253 s.nErr++;
254 rc = SQLITE_ERROR;
255 goto parse_json_err;
256 }
257 if( geopolySkipSpace(&s)==',' ){
258 s.z++;
259 continue;
260 }
261 break;
262 }
263 if( geopolySkipSpace(&s)==']'
264 && s.nVertex>=4
265 && s.a[0]==s.a[s.nVertex*2-2]
266 && s.a[1]==s.a[s.nVertex*2-1]
267 && (s.z++, geopolySkipSpace(&s)==0)
268 ){
269 GeoPoly *pOut;
270 int x = 1;
271 s.nVertex--; /* Remove the redundant vertex at the end */
272 pOut = sqlite3_malloc64( GEOPOLY_SZ((sqlite3_int64)s.nVertex) );
273 x = 1;
274 if( pOut==0 ) goto parse_json_err;
275 pOut->nVertex = s.nVertex;
276 memcpy(pOut->a, s.a, s.nVertex*2*sizeof(GeoCoord));
277 pOut->hdr[0] = *(unsigned char*)&x;
278 pOut->hdr[1] = (s.nVertex>>16)&0xff;
279 pOut->hdr[2] = (s.nVertex>>8)&0xff;
280 pOut->hdr[3] = s.nVertex&0xff;
281 sqlite3_free(s.a);
282 if( pRc ) *pRc = SQLITE_OK;
283 return pOut;
284 }else{
285 s.nErr++;
286 rc = SQLITE_ERROR;
287 }
288 }
289 parse_json_err:
290 if( pRc ) *pRc = rc;
291 sqlite3_free(s.a);
292 return 0;
293 }
294
295 /*
296 ** Given a function parameter, try to interpret it as a polygon, either
297 ** in the binary format or JSON text. Compute a GeoPoly object and
298 ** return a pointer to that object. Or if the input is not a well-formed
299 ** polygon, put an error message in sqlite3_context and return NULL.
300 */
geopolyFuncParam(sqlite3_context * pCtx,sqlite3_value * pVal,int * pRc)301 static GeoPoly *geopolyFuncParam(
302 sqlite3_context *pCtx, /* Context for error messages */
303 sqlite3_value *pVal, /* The value to decode */
304 int *pRc /* Write error here */
305 ){
306 GeoPoly *p = 0;
307 int nByte;
308 if( sqlite3_value_type(pVal)==SQLITE_BLOB
309 && (nByte = sqlite3_value_bytes(pVal))>=(4+6*sizeof(GeoCoord))
310 ){
311 const unsigned char *a = sqlite3_value_blob(pVal);
312 int nVertex;
313 nVertex = (a[1]<<16) + (a[2]<<8) + a[3];
314 if( (a[0]==0 || a[0]==1)
315 && (nVertex*2*sizeof(GeoCoord) + 4)==(unsigned int)nByte
316 ){
317 p = sqlite3_malloc64( sizeof(*p) + (nVertex-1)*2*sizeof(GeoCoord) );
318 if( p==0 ){
319 if( pRc ) *pRc = SQLITE_NOMEM;
320 if( pCtx ) sqlite3_result_error_nomem(pCtx);
321 }else{
322 int x = 1;
323 p->nVertex = nVertex;
324 memcpy(p->hdr, a, nByte);
325 if( a[0] != *(unsigned char*)&x ){
326 int ii;
327 for(ii=0; ii<nVertex; ii++){
328 geopolySwab32((unsigned char*)&GeoX(p,ii));
329 geopolySwab32((unsigned char*)&GeoY(p,ii));
330 }
331 p->hdr[0] ^= 1;
332 }
333 }
334 }
335 if( pRc ) *pRc = SQLITE_OK;
336 return p;
337 }else if( sqlite3_value_type(pVal)==SQLITE_TEXT ){
338 const unsigned char *zJson = sqlite3_value_text(pVal);
339 if( zJson==0 ){
340 if( pRc ) *pRc = SQLITE_NOMEM;
341 return 0;
342 }
343 return geopolyParseJson(zJson, pRc);
344 }else{
345 if( pRc ) *pRc = SQLITE_ERROR;
346 return 0;
347 }
348 }
349
350 /*
351 ** Implementation of the geopoly_blob(X) function.
352 **
353 ** If the input is a well-formed Geopoly BLOB or JSON string
354 ** then return the BLOB representation of the polygon. Otherwise
355 ** return NULL.
356 */
geopolyBlobFunc(sqlite3_context * context,int argc,sqlite3_value ** argv)357 static void geopolyBlobFunc(
358 sqlite3_context *context,
359 int argc,
360 sqlite3_value **argv
361 ){
362 GeoPoly *p = geopolyFuncParam(context, argv[0], 0);
363 if( p ){
364 sqlite3_result_blob(context, p->hdr,
365 4+8*p->nVertex, SQLITE_TRANSIENT);
366 sqlite3_free(p);
367 }
368 }
369
370 /*
371 ** SQL function: geopoly_json(X)
372 **
373 ** Interpret X as a polygon and render it as a JSON array
374 ** of coordinates. Or, if X is not a valid polygon, return NULL.
375 */
geopolyJsonFunc(sqlite3_context * context,int argc,sqlite3_value ** argv)376 static void geopolyJsonFunc(
377 sqlite3_context *context,
378 int argc,
379 sqlite3_value **argv
380 ){
381 GeoPoly *p = geopolyFuncParam(context, argv[0], 0);
382 if( p ){
383 sqlite3 *db = sqlite3_context_db_handle(context);
384 sqlite3_str *x = sqlite3_str_new(db);
385 int i;
386 sqlite3_str_append(x, "[", 1);
387 for(i=0; i<p->nVertex; i++){
388 sqlite3_str_appendf(x, "[%!g,%!g],", GeoX(p,i), GeoY(p,i));
389 }
390 sqlite3_str_appendf(x, "[%!g,%!g]]", GeoX(p,0), GeoY(p,0));
391 sqlite3_result_text(context, sqlite3_str_finish(x), -1, sqlite3_free);
392 sqlite3_free(p);
393 }
394 }
395
396 /*
397 ** SQL function: geopoly_svg(X, ....)
398 **
399 ** Interpret X as a polygon and render it as a SVG <polyline>.
400 ** Additional arguments are added as attributes to the <polyline>.
401 */
geopolySvgFunc(sqlite3_context * context,int argc,sqlite3_value ** argv)402 static void geopolySvgFunc(
403 sqlite3_context *context,
404 int argc,
405 sqlite3_value **argv
406 ){
407 GeoPoly *p;
408 if( argc<1 ) return;
409 p = geopolyFuncParam(context, argv[0], 0);
410 if( p ){
411 sqlite3 *db = sqlite3_context_db_handle(context);
412 sqlite3_str *x = sqlite3_str_new(db);
413 int i;
414 char cSep = '\'';
415 sqlite3_str_appendf(x, "<polyline points=");
416 for(i=0; i<p->nVertex; i++){
417 sqlite3_str_appendf(x, "%c%g,%g", cSep, GeoX(p,i), GeoY(p,i));
418 cSep = ' ';
419 }
420 sqlite3_str_appendf(x, " %g,%g'", GeoX(p,0), GeoY(p,0));
421 for(i=1; i<argc; i++){
422 const char *z = (const char*)sqlite3_value_text(argv[i]);
423 if( z && z[0] ){
424 sqlite3_str_appendf(x, " %s", z);
425 }
426 }
427 sqlite3_str_appendf(x, "></polyline>");
428 sqlite3_result_text(context, sqlite3_str_finish(x), -1, sqlite3_free);
429 sqlite3_free(p);
430 }
431 }
432
433 /*
434 ** SQL Function: geopoly_xform(poly, A, B, C, D, E, F)
435 **
436 ** Transform and/or translate a polygon as follows:
437 **
438 ** x1 = A*x0 + B*y0 + E
439 ** y1 = C*x0 + D*y0 + F
440 **
441 ** For a translation:
442 **
443 ** geopoly_xform(poly, 1, 0, 0, 1, x-offset, y-offset)
444 **
445 ** Rotate by R around the point (0,0):
446 **
447 ** geopoly_xform(poly, cos(R), sin(R), -sin(R), cos(R), 0, 0)
448 */
geopolyXformFunc(sqlite3_context * context,int argc,sqlite3_value ** argv)449 static void geopolyXformFunc(
450 sqlite3_context *context,
451 int argc,
452 sqlite3_value **argv
453 ){
454 GeoPoly *p = geopolyFuncParam(context, argv[0], 0);
455 double A = sqlite3_value_double(argv[1]);
456 double B = sqlite3_value_double(argv[2]);
457 double C = sqlite3_value_double(argv[3]);
458 double D = sqlite3_value_double(argv[4]);
459 double E = sqlite3_value_double(argv[5]);
460 double F = sqlite3_value_double(argv[6]);
461 GeoCoord x1, y1, x0, y0;
462 int ii;
463 if( p ){
464 for(ii=0; ii<p->nVertex; ii++){
465 x0 = GeoX(p,ii);
466 y0 = GeoY(p,ii);
467 x1 = (GeoCoord)(A*x0 + B*y0 + E);
468 y1 = (GeoCoord)(C*x0 + D*y0 + F);
469 GeoX(p,ii) = x1;
470 GeoY(p,ii) = y1;
471 }
472 sqlite3_result_blob(context, p->hdr,
473 4+8*p->nVertex, SQLITE_TRANSIENT);
474 sqlite3_free(p);
475 }
476 }
477
478 /*
479 ** Compute the area enclosed by the polygon.
480 **
481 ** This routine can also be used to detect polygons that rotate in
482 ** the wrong direction. Polygons are suppose to be counter-clockwise (CCW).
483 ** This routine returns a negative value for clockwise (CW) polygons.
484 */
geopolyArea(GeoPoly * p)485 static double geopolyArea(GeoPoly *p){
486 double rArea = 0.0;
487 int ii;
488 for(ii=0; ii<p->nVertex-1; ii++){
489 rArea += (GeoX(p,ii) - GeoX(p,ii+1)) /* (x0 - x1) */
490 * (GeoY(p,ii) + GeoY(p,ii+1)) /* (y0 + y1) */
491 * 0.5;
492 }
493 rArea += (GeoX(p,ii) - GeoX(p,0)) /* (xN - x0) */
494 * (GeoY(p,ii) + GeoY(p,0)) /* (yN + y0) */
495 * 0.5;
496 return rArea;
497 }
498
499 /*
500 ** Implementation of the geopoly_area(X) function.
501 **
502 ** If the input is a well-formed Geopoly BLOB then return the area
503 ** enclosed by the polygon. If the polygon circulates clockwise instead
504 ** of counterclockwise (as it should) then return the negative of the
505 ** enclosed area. Otherwise return NULL.
506 */
geopolyAreaFunc(sqlite3_context * context,int argc,sqlite3_value ** argv)507 static void geopolyAreaFunc(
508 sqlite3_context *context,
509 int argc,
510 sqlite3_value **argv
511 ){
512 GeoPoly *p = geopolyFuncParam(context, argv[0], 0);
513 if( p ){
514 sqlite3_result_double(context, geopolyArea(p));
515 sqlite3_free(p);
516 }
517 }
518
519 /*
520 ** Implementation of the geopoly_ccw(X) function.
521 **
522 ** If the rotation of polygon X is clockwise (incorrect) instead of
523 ** counter-clockwise (the correct winding order according to RFC7946)
524 ** then reverse the order of the vertexes in polygon X.
525 **
526 ** In other words, this routine returns a CCW polygon regardless of the
527 ** winding order of its input.
528 **
529 ** Use this routine to sanitize historical inputs that that sometimes
530 ** contain polygons that wind in the wrong direction.
531 */
geopolyCcwFunc(sqlite3_context * context,int argc,sqlite3_value ** argv)532 static void geopolyCcwFunc(
533 sqlite3_context *context,
534 int argc,
535 sqlite3_value **argv
536 ){
537 GeoPoly *p = geopolyFuncParam(context, argv[0], 0);
538 if( p ){
539 if( geopolyArea(p)<0.0 ){
540 int ii, jj;
541 for(ii=1, jj=p->nVertex-1; ii<jj; ii++, jj--){
542 GeoCoord t = GeoX(p,ii);
543 GeoX(p,ii) = GeoX(p,jj);
544 GeoX(p,jj) = t;
545 t = GeoY(p,ii);
546 GeoY(p,ii) = GeoY(p,jj);
547 GeoY(p,jj) = t;
548 }
549 }
550 sqlite3_result_blob(context, p->hdr,
551 4+8*p->nVertex, SQLITE_TRANSIENT);
552 sqlite3_free(p);
553 }
554 }
555
556 #define GEOPOLY_PI 3.1415926535897932385
557
558 /* Fast approximation for sine(X) for X between -0.5*pi and 2*pi
559 */
geopolySine(double r)560 static double geopolySine(double r){
561 assert( r>=-0.5*GEOPOLY_PI && r<=2.0*GEOPOLY_PI );
562 if( r>=1.5*GEOPOLY_PI ){
563 r -= 2.0*GEOPOLY_PI;
564 }
565 if( r>=0.5*GEOPOLY_PI ){
566 return -geopolySine(r-GEOPOLY_PI);
567 }else{
568 double r2 = r*r;
569 double r3 = r2*r;
570 double r5 = r3*r2;
571 return 0.9996949*r - 0.1656700*r3 + 0.0075134*r5;
572 }
573 }
574
575 /*
576 ** Function: geopoly_regular(X,Y,R,N)
577 **
578 ** Construct a simple, convex, regular polygon centered at X, Y
579 ** with circumradius R and with N sides.
580 */
geopolyRegularFunc(sqlite3_context * context,int argc,sqlite3_value ** argv)581 static void geopolyRegularFunc(
582 sqlite3_context *context,
583 int argc,
584 sqlite3_value **argv
585 ){
586 double x = sqlite3_value_double(argv[0]);
587 double y = sqlite3_value_double(argv[1]);
588 double r = sqlite3_value_double(argv[2]);
589 int n = sqlite3_value_int(argv[3]);
590 int i;
591 GeoPoly *p;
592
593 if( n<3 || r<=0.0 ) return;
594 if( n>1000 ) n = 1000;
595 p = sqlite3_malloc64( sizeof(*p) + (n-1)*2*sizeof(GeoCoord) );
596 if( p==0 ){
597 sqlite3_result_error_nomem(context);
598 return;
599 }
600 i = 1;
601 p->hdr[0] = *(unsigned char*)&i;
602 p->hdr[1] = 0;
603 p->hdr[2] = (n>>8)&0xff;
604 p->hdr[3] = n&0xff;
605 for(i=0; i<n; i++){
606 double rAngle = 2.0*GEOPOLY_PI*i/n;
607 GeoX(p,i) = x - r*geopolySine(rAngle-0.5*GEOPOLY_PI);
608 GeoY(p,i) = y + r*geopolySine(rAngle);
609 }
610 sqlite3_result_blob(context, p->hdr, 4+8*n, SQLITE_TRANSIENT);
611 sqlite3_free(p);
612 }
613
614 /*
615 ** If pPoly is a polygon, compute its bounding box. Then:
616 **
617 ** (1) if aCoord!=0 store the bounding box in aCoord, returning NULL
618 ** (2) otherwise, compute a GeoPoly for the bounding box and return the
619 ** new GeoPoly
620 **
621 ** If pPoly is NULL but aCoord is not NULL, then compute a new GeoPoly from
622 ** the bounding box in aCoord and return a pointer to that GeoPoly.
623 */
geopolyBBox(sqlite3_context * context,sqlite3_value * pPoly,RtreeCoord * aCoord,int * pRc)624 static GeoPoly *geopolyBBox(
625 sqlite3_context *context, /* For recording the error */
626 sqlite3_value *pPoly, /* The polygon */
627 RtreeCoord *aCoord, /* Results here */
628 int *pRc /* Error code here */
629 ){
630 GeoPoly *pOut = 0;
631 GeoPoly *p;
632 float mnX, mxX, mnY, mxY;
633 if( pPoly==0 && aCoord!=0 ){
634 p = 0;
635 mnX = aCoord[0].f;
636 mxX = aCoord[1].f;
637 mnY = aCoord[2].f;
638 mxY = aCoord[3].f;
639 goto geopolyBboxFill;
640 }else{
641 p = geopolyFuncParam(context, pPoly, pRc);
642 }
643 if( p ){
644 int ii;
645 mnX = mxX = GeoX(p,0);
646 mnY = mxY = GeoY(p,0);
647 for(ii=1; ii<p->nVertex; ii++){
648 double r = GeoX(p,ii);
649 if( r<mnX ) mnX = (float)r;
650 else if( r>mxX ) mxX = (float)r;
651 r = GeoY(p,ii);
652 if( r<mnY ) mnY = (float)r;
653 else if( r>mxY ) mxY = (float)r;
654 }
655 if( pRc ) *pRc = SQLITE_OK;
656 if( aCoord==0 ){
657 geopolyBboxFill:
658 pOut = sqlite3_realloc64(p, GEOPOLY_SZ(4));
659 if( pOut==0 ){
660 sqlite3_free(p);
661 if( context ) sqlite3_result_error_nomem(context);
662 if( pRc ) *pRc = SQLITE_NOMEM;
663 return 0;
664 }
665 pOut->nVertex = 4;
666 ii = 1;
667 pOut->hdr[0] = *(unsigned char*)ⅈ
668 pOut->hdr[1] = 0;
669 pOut->hdr[2] = 0;
670 pOut->hdr[3] = 4;
671 GeoX(pOut,0) = mnX;
672 GeoY(pOut,0) = mnY;
673 GeoX(pOut,1) = mxX;
674 GeoY(pOut,1) = mnY;
675 GeoX(pOut,2) = mxX;
676 GeoY(pOut,2) = mxY;
677 GeoX(pOut,3) = mnX;
678 GeoY(pOut,3) = mxY;
679 }else{
680 sqlite3_free(p);
681 aCoord[0].f = mnX;
682 aCoord[1].f = mxX;
683 aCoord[2].f = mnY;
684 aCoord[3].f = mxY;
685 }
686 }
687 return pOut;
688 }
689
690 /*
691 ** Implementation of the geopoly_bbox(X) SQL function.
692 */
geopolyBBoxFunc(sqlite3_context * context,int argc,sqlite3_value ** argv)693 static void geopolyBBoxFunc(
694 sqlite3_context *context,
695 int argc,
696 sqlite3_value **argv
697 ){
698 GeoPoly *p = geopolyBBox(context, argv[0], 0, 0);
699 if( p ){
700 sqlite3_result_blob(context, p->hdr,
701 4+8*p->nVertex, SQLITE_TRANSIENT);
702 sqlite3_free(p);
703 }
704 }
705
706 /*
707 ** State vector for the geopoly_group_bbox() aggregate function.
708 */
709 typedef struct GeoBBox GeoBBox;
710 struct GeoBBox {
711 int isInit;
712 RtreeCoord a[4];
713 };
714
715
716 /*
717 ** Implementation of the geopoly_group_bbox(X) aggregate SQL function.
718 */
geopolyBBoxStep(sqlite3_context * context,int argc,sqlite3_value ** argv)719 static void geopolyBBoxStep(
720 sqlite3_context *context,
721 int argc,
722 sqlite3_value **argv
723 ){
724 RtreeCoord a[4];
725 int rc = SQLITE_OK;
726 (void)geopolyBBox(context, argv[0], a, &rc);
727 if( rc==SQLITE_OK ){
728 GeoBBox *pBBox;
729 pBBox = (GeoBBox*)sqlite3_aggregate_context(context, sizeof(*pBBox));
730 if( pBBox==0 ) return;
731 if( pBBox->isInit==0 ){
732 pBBox->isInit = 1;
733 memcpy(pBBox->a, a, sizeof(RtreeCoord)*4);
734 }else{
735 if( a[0].f < pBBox->a[0].f ) pBBox->a[0] = a[0];
736 if( a[1].f > pBBox->a[1].f ) pBBox->a[1] = a[1];
737 if( a[2].f < pBBox->a[2].f ) pBBox->a[2] = a[2];
738 if( a[3].f > pBBox->a[3].f ) pBBox->a[3] = a[3];
739 }
740 }
741 }
geopolyBBoxFinal(sqlite3_context * context)742 static void geopolyBBoxFinal(
743 sqlite3_context *context
744 ){
745 GeoPoly *p;
746 GeoBBox *pBBox;
747 pBBox = (GeoBBox*)sqlite3_aggregate_context(context, 0);
748 if( pBBox==0 ) return;
749 p = geopolyBBox(context, 0, pBBox->a, 0);
750 if( p ){
751 sqlite3_result_blob(context, p->hdr,
752 4+8*p->nVertex, SQLITE_TRANSIENT);
753 sqlite3_free(p);
754 }
755 }
756
757
758 /*
759 ** Determine if point (x0,y0) is beneath line segment (x1,y1)->(x2,y2).
760 ** Returns:
761 **
762 ** +2 x0,y0 is on the line segement
763 **
764 ** +1 x0,y0 is beneath line segment
765 **
766 ** 0 x0,y0 is not on or beneath the line segment or the line segment
767 ** is vertical and x0,y0 is not on the line segment
768 **
769 ** The left-most coordinate min(x1,x2) is not considered to be part of
770 ** the line segment for the purposes of this analysis.
771 */
pointBeneathLine(double x0,double y0,double x1,double y1,double x2,double y2)772 static int pointBeneathLine(
773 double x0, double y0,
774 double x1, double y1,
775 double x2, double y2
776 ){
777 double y;
778 if( x0==x1 && y0==y1 ) return 2;
779 if( x1<x2 ){
780 if( x0<=x1 || x0>x2 ) return 0;
781 }else if( x1>x2 ){
782 if( x0<=x2 || x0>x1 ) return 0;
783 }else{
784 /* Vertical line segment */
785 if( x0!=x1 ) return 0;
786 if( y0<y1 && y0<y2 ) return 0;
787 if( y0>y1 && y0>y2 ) return 0;
788 return 2;
789 }
790 y = y1 + (y2-y1)*(x0-x1)/(x2-x1);
791 if( y0==y ) return 2;
792 if( y0<y ) return 1;
793 return 0;
794 }
795
796 /*
797 ** SQL function: geopoly_contains_point(P,X,Y)
798 **
799 ** Return +2 if point X,Y is within polygon P.
800 ** Return +1 if point X,Y is on the polygon boundary.
801 ** Return 0 if point X,Y is outside the polygon
802 */
geopolyContainsPointFunc(sqlite3_context * context,int argc,sqlite3_value ** argv)803 static void geopolyContainsPointFunc(
804 sqlite3_context *context,
805 int argc,
806 sqlite3_value **argv
807 ){
808 GeoPoly *p1 = geopolyFuncParam(context, argv[0], 0);
809 double x0 = sqlite3_value_double(argv[1]);
810 double y0 = sqlite3_value_double(argv[2]);
811 int v = 0;
812 int cnt = 0;
813 int ii;
814 if( p1==0 ) return;
815 for(ii=0; ii<p1->nVertex-1; ii++){
816 v = pointBeneathLine(x0,y0,GeoX(p1,ii), GeoY(p1,ii),
817 GeoX(p1,ii+1),GeoY(p1,ii+1));
818 if( v==2 ) break;
819 cnt += v;
820 }
821 if( v!=2 ){
822 v = pointBeneathLine(x0,y0,GeoX(p1,ii), GeoY(p1,ii),
823 GeoX(p1,0), GeoY(p1,0));
824 }
825 if( v==2 ){
826 sqlite3_result_int(context, 1);
827 }else if( ((v+cnt)&1)==0 ){
828 sqlite3_result_int(context, 0);
829 }else{
830 sqlite3_result_int(context, 2);
831 }
832 sqlite3_free(p1);
833 }
834
835 /* Forward declaration */
836 static int geopolyOverlap(GeoPoly *p1, GeoPoly *p2);
837
838 /*
839 ** SQL function: geopoly_within(P1,P2)
840 **
841 ** Return +2 if P1 and P2 are the same polygon
842 ** Return +1 if P2 is contained within P1
843 ** Return 0 if any part of P2 is on the outside of P1
844 **
845 */
geopolyWithinFunc(sqlite3_context * context,int argc,sqlite3_value ** argv)846 static void geopolyWithinFunc(
847 sqlite3_context *context,
848 int argc,
849 sqlite3_value **argv
850 ){
851 GeoPoly *p1 = geopolyFuncParam(context, argv[0], 0);
852 GeoPoly *p2 = geopolyFuncParam(context, argv[1], 0);
853 if( p1 && p2 ){
854 int x = geopolyOverlap(p1, p2);
855 if( x<0 ){
856 sqlite3_result_error_nomem(context);
857 }else{
858 sqlite3_result_int(context, x==2 ? 1 : x==4 ? 2 : 0);
859 }
860 }
861 sqlite3_free(p1);
862 sqlite3_free(p2);
863 }
864
865 /* Objects used by the overlap algorihm. */
866 typedef struct GeoEvent GeoEvent;
867 typedef struct GeoSegment GeoSegment;
868 typedef struct GeoOverlap GeoOverlap;
869 struct GeoEvent {
870 double x; /* X coordinate at which event occurs */
871 int eType; /* 0 for ADD, 1 for REMOVE */
872 GeoSegment *pSeg; /* The segment to be added or removed */
873 GeoEvent *pNext; /* Next event in the sorted list */
874 };
875 struct GeoSegment {
876 double C, B; /* y = C*x + B */
877 double y; /* Current y value */
878 float y0; /* Initial y value */
879 unsigned char side; /* 1 for p1, 2 for p2 */
880 unsigned int idx; /* Which segment within the side */
881 GeoSegment *pNext; /* Next segment in a list sorted by y */
882 };
883 struct GeoOverlap {
884 GeoEvent *aEvent; /* Array of all events */
885 GeoSegment *aSegment; /* Array of all segments */
886 int nEvent; /* Number of events */
887 int nSegment; /* Number of segments */
888 };
889
890 /*
891 ** Add a single segment and its associated events.
892 */
geopolyAddOneSegment(GeoOverlap * p,GeoCoord x0,GeoCoord y0,GeoCoord x1,GeoCoord y1,unsigned char side,unsigned int idx)893 static void geopolyAddOneSegment(
894 GeoOverlap *p,
895 GeoCoord x0,
896 GeoCoord y0,
897 GeoCoord x1,
898 GeoCoord y1,
899 unsigned char side,
900 unsigned int idx
901 ){
902 GeoSegment *pSeg;
903 GeoEvent *pEvent;
904 if( x0==x1 ) return; /* Ignore vertical segments */
905 if( x0>x1 ){
906 GeoCoord t = x0;
907 x0 = x1;
908 x1 = t;
909 t = y0;
910 y0 = y1;
911 y1 = t;
912 }
913 pSeg = p->aSegment + p->nSegment;
914 p->nSegment++;
915 pSeg->C = (y1-y0)/(x1-x0);
916 pSeg->B = y1 - x1*pSeg->C;
917 pSeg->y0 = y0;
918 pSeg->side = side;
919 pSeg->idx = idx;
920 pEvent = p->aEvent + p->nEvent;
921 p->nEvent++;
922 pEvent->x = x0;
923 pEvent->eType = 0;
924 pEvent->pSeg = pSeg;
925 pEvent = p->aEvent + p->nEvent;
926 p->nEvent++;
927 pEvent->x = x1;
928 pEvent->eType = 1;
929 pEvent->pSeg = pSeg;
930 }
931
932
933
934 /*
935 ** Insert all segments and events for polygon pPoly.
936 */
geopolyAddSegments(GeoOverlap * p,GeoPoly * pPoly,unsigned char side)937 static void geopolyAddSegments(
938 GeoOverlap *p, /* Add segments to this Overlap object */
939 GeoPoly *pPoly, /* Take all segments from this polygon */
940 unsigned char side /* The side of pPoly */
941 ){
942 unsigned int i;
943 GeoCoord *x;
944 for(i=0; i<(unsigned)pPoly->nVertex-1; i++){
945 x = &GeoX(pPoly,i);
946 geopolyAddOneSegment(p, x[0], x[1], x[2], x[3], side, i);
947 }
948 x = &GeoX(pPoly,i);
949 geopolyAddOneSegment(p, x[0], x[1], pPoly->a[0], pPoly->a[1], side, i);
950 }
951
952 /*
953 ** Merge two lists of sorted events by X coordinate
954 */
geopolyEventMerge(GeoEvent * pLeft,GeoEvent * pRight)955 static GeoEvent *geopolyEventMerge(GeoEvent *pLeft, GeoEvent *pRight){
956 GeoEvent head, *pLast;
957 head.pNext = 0;
958 pLast = &head;
959 while( pRight && pLeft ){
960 if( pRight->x <= pLeft->x ){
961 pLast->pNext = pRight;
962 pLast = pRight;
963 pRight = pRight->pNext;
964 }else{
965 pLast->pNext = pLeft;
966 pLast = pLeft;
967 pLeft = pLeft->pNext;
968 }
969 }
970 pLast->pNext = pRight ? pRight : pLeft;
971 return head.pNext;
972 }
973
974 /*
975 ** Sort an array of nEvent event objects into a list.
976 */
geopolySortEventsByX(GeoEvent * aEvent,int nEvent)977 static GeoEvent *geopolySortEventsByX(GeoEvent *aEvent, int nEvent){
978 int mx = 0;
979 int i, j;
980 GeoEvent *p;
981 GeoEvent *a[50];
982 for(i=0; i<nEvent; i++){
983 p = &aEvent[i];
984 p->pNext = 0;
985 for(j=0; j<mx && a[j]; j++){
986 p = geopolyEventMerge(a[j], p);
987 a[j] = 0;
988 }
989 a[j] = p;
990 if( j>=mx ) mx = j+1;
991 }
992 p = 0;
993 for(i=0; i<mx; i++){
994 p = geopolyEventMerge(a[i], p);
995 }
996 return p;
997 }
998
999 /*
1000 ** Merge two lists of sorted segments by Y, and then by C.
1001 */
geopolySegmentMerge(GeoSegment * pLeft,GeoSegment * pRight)1002 static GeoSegment *geopolySegmentMerge(GeoSegment *pLeft, GeoSegment *pRight){
1003 GeoSegment head, *pLast;
1004 head.pNext = 0;
1005 pLast = &head;
1006 while( pRight && pLeft ){
1007 double r = pRight->y - pLeft->y;
1008 if( r==0.0 ) r = pRight->C - pLeft->C;
1009 if( r<0.0 ){
1010 pLast->pNext = pRight;
1011 pLast = pRight;
1012 pRight = pRight->pNext;
1013 }else{
1014 pLast->pNext = pLeft;
1015 pLast = pLeft;
1016 pLeft = pLeft->pNext;
1017 }
1018 }
1019 pLast->pNext = pRight ? pRight : pLeft;
1020 return head.pNext;
1021 }
1022
1023 /*
1024 ** Sort a list of GeoSegments in order of increasing Y and in the event of
1025 ** a tie, increasing C (slope).
1026 */
geopolySortSegmentsByYAndC(GeoSegment * pList)1027 static GeoSegment *geopolySortSegmentsByYAndC(GeoSegment *pList){
1028 int mx = 0;
1029 int i;
1030 GeoSegment *p;
1031 GeoSegment *a[50];
1032 while( pList ){
1033 p = pList;
1034 pList = pList->pNext;
1035 p->pNext = 0;
1036 for(i=0; i<mx && a[i]; i++){
1037 p = geopolySegmentMerge(a[i], p);
1038 a[i] = 0;
1039 }
1040 a[i] = p;
1041 if( i>=mx ) mx = i+1;
1042 }
1043 p = 0;
1044 for(i=0; i<mx; i++){
1045 p = geopolySegmentMerge(a[i], p);
1046 }
1047 return p;
1048 }
1049
1050 /*
1051 ** Determine the overlap between two polygons
1052 */
geopolyOverlap(GeoPoly * p1,GeoPoly * p2)1053 static int geopolyOverlap(GeoPoly *p1, GeoPoly *p2){
1054 sqlite3_int64 nVertex = p1->nVertex + p2->nVertex + 2;
1055 GeoOverlap *p;
1056 sqlite3_int64 nByte;
1057 GeoEvent *pThisEvent;
1058 double rX;
1059 int rc = 0;
1060 int needSort = 0;
1061 GeoSegment *pActive = 0;
1062 GeoSegment *pSeg;
1063 unsigned char aOverlap[4];
1064
1065 nByte = sizeof(GeoEvent)*nVertex*2
1066 + sizeof(GeoSegment)*nVertex
1067 + sizeof(GeoOverlap);
1068 p = sqlite3_malloc64( nByte );
1069 if( p==0 ) return -1;
1070 p->aEvent = (GeoEvent*)&p[1];
1071 p->aSegment = (GeoSegment*)&p->aEvent[nVertex*2];
1072 p->nEvent = p->nSegment = 0;
1073 geopolyAddSegments(p, p1, 1);
1074 geopolyAddSegments(p, p2, 2);
1075 pThisEvent = geopolySortEventsByX(p->aEvent, p->nEvent);
1076 rX = pThisEvent->x==0.0 ? -1.0 : 0.0;
1077 memset(aOverlap, 0, sizeof(aOverlap));
1078 while( pThisEvent ){
1079 if( pThisEvent->x!=rX ){
1080 GeoSegment *pPrev = 0;
1081 int iMask = 0;
1082 GEODEBUG(("Distinct X: %g\n", pThisEvent->x));
1083 rX = pThisEvent->x;
1084 if( needSort ){
1085 GEODEBUG(("SORT\n"));
1086 pActive = geopolySortSegmentsByYAndC(pActive);
1087 needSort = 0;
1088 }
1089 for(pSeg=pActive; pSeg; pSeg=pSeg->pNext){
1090 if( pPrev ){
1091 if( pPrev->y!=pSeg->y ){
1092 GEODEBUG(("MASK: %d\n", iMask));
1093 aOverlap[iMask] = 1;
1094 }
1095 }
1096 iMask ^= pSeg->side;
1097 pPrev = pSeg;
1098 }
1099 pPrev = 0;
1100 for(pSeg=pActive; pSeg; pSeg=pSeg->pNext){
1101 double y = pSeg->C*rX + pSeg->B;
1102 GEODEBUG(("Segment %d.%d %g->%g\n", pSeg->side, pSeg->idx, pSeg->y, y));
1103 pSeg->y = y;
1104 if( pPrev ){
1105 if( pPrev->y>pSeg->y && pPrev->side!=pSeg->side ){
1106 rc = 1;
1107 GEODEBUG(("Crossing: %d.%d and %d.%d\n",
1108 pPrev->side, pPrev->idx,
1109 pSeg->side, pSeg->idx));
1110 goto geopolyOverlapDone;
1111 }else if( pPrev->y!=pSeg->y ){
1112 GEODEBUG(("MASK: %d\n", iMask));
1113 aOverlap[iMask] = 1;
1114 }
1115 }
1116 iMask ^= pSeg->side;
1117 pPrev = pSeg;
1118 }
1119 }
1120 GEODEBUG(("%s %d.%d C=%g B=%g\n",
1121 pThisEvent->eType ? "RM " : "ADD",
1122 pThisEvent->pSeg->side, pThisEvent->pSeg->idx,
1123 pThisEvent->pSeg->C,
1124 pThisEvent->pSeg->B));
1125 if( pThisEvent->eType==0 ){
1126 /* Add a segment */
1127 pSeg = pThisEvent->pSeg;
1128 pSeg->y = pSeg->y0;
1129 pSeg->pNext = pActive;
1130 pActive = pSeg;
1131 needSort = 1;
1132 }else{
1133 /* Remove a segment */
1134 if( pActive==pThisEvent->pSeg ){
1135 pActive = pActive->pNext;
1136 }else{
1137 for(pSeg=pActive; pSeg; pSeg=pSeg->pNext){
1138 if( pSeg->pNext==pThisEvent->pSeg ){
1139 pSeg->pNext = pSeg->pNext->pNext;
1140 break;
1141 }
1142 }
1143 }
1144 }
1145 pThisEvent = pThisEvent->pNext;
1146 }
1147 if( aOverlap[3]==0 ){
1148 rc = 0;
1149 }else if( aOverlap[1]!=0 && aOverlap[2]==0 ){
1150 rc = 3;
1151 }else if( aOverlap[1]==0 && aOverlap[2]!=0 ){
1152 rc = 2;
1153 }else if( aOverlap[1]==0 && aOverlap[2]==0 ){
1154 rc = 4;
1155 }else{
1156 rc = 1;
1157 }
1158
1159 geopolyOverlapDone:
1160 sqlite3_free(p);
1161 return rc;
1162 }
1163
1164 /*
1165 ** SQL function: geopoly_overlap(P1,P2)
1166 **
1167 ** Determine whether or not P1 and P2 overlap. Return value:
1168 **
1169 ** 0 The two polygons are disjoint
1170 ** 1 They overlap
1171 ** 2 P1 is completely contained within P2
1172 ** 3 P2 is completely contained within P1
1173 ** 4 P1 and P2 are the same polygon
1174 ** NULL Either P1 or P2 or both are not valid polygons
1175 */
geopolyOverlapFunc(sqlite3_context * context,int argc,sqlite3_value ** argv)1176 static void geopolyOverlapFunc(
1177 sqlite3_context *context,
1178 int argc,
1179 sqlite3_value **argv
1180 ){
1181 GeoPoly *p1 = geopolyFuncParam(context, argv[0], 0);
1182 GeoPoly *p2 = geopolyFuncParam(context, argv[1], 0);
1183 if( p1 && p2 ){
1184 int x = geopolyOverlap(p1, p2);
1185 if( x<0 ){
1186 sqlite3_result_error_nomem(context);
1187 }else{
1188 sqlite3_result_int(context, x);
1189 }
1190 }
1191 sqlite3_free(p1);
1192 sqlite3_free(p2);
1193 }
1194
1195 /*
1196 ** Enable or disable debugging output
1197 */
geopolyDebugFunc(sqlite3_context * context,int argc,sqlite3_value ** argv)1198 static void geopolyDebugFunc(
1199 sqlite3_context *context,
1200 int argc,
1201 sqlite3_value **argv
1202 ){
1203 #ifdef GEOPOLY_ENABLE_DEBUG
1204 geo_debug = sqlite3_value_int(argv[0]);
1205 #endif
1206 }
1207
1208 /*
1209 ** This function is the implementation of both the xConnect and xCreate
1210 ** methods of the geopoly virtual table.
1211 **
1212 ** argv[0] -> module name
1213 ** argv[1] -> database name
1214 ** argv[2] -> table name
1215 ** argv[...] -> column names...
1216 */
geopolyInit(sqlite3 * db,void * pAux,int argc,const char * const * argv,sqlite3_vtab ** ppVtab,char ** pzErr,int isCreate)1217 static int geopolyInit(
1218 sqlite3 *db, /* Database connection */
1219 void *pAux, /* One of the RTREE_COORD_* constants */
1220 int argc, const char *const*argv, /* Parameters to CREATE TABLE statement */
1221 sqlite3_vtab **ppVtab, /* OUT: New virtual table */
1222 char **pzErr, /* OUT: Error message, if any */
1223 int isCreate /* True for xCreate, false for xConnect */
1224 ){
1225 int rc = SQLITE_OK;
1226 Rtree *pRtree;
1227 sqlite3_int64 nDb; /* Length of string argv[1] */
1228 sqlite3_int64 nName; /* Length of string argv[2] */
1229 sqlite3_str *pSql;
1230 char *zSql;
1231 int ii;
1232
1233 sqlite3_vtab_config(db, SQLITE_VTAB_CONSTRAINT_SUPPORT, 1);
1234
1235 /* Allocate the sqlite3_vtab structure */
1236 nDb = strlen(argv[1]);
1237 nName = strlen(argv[2]);
1238 pRtree = (Rtree *)sqlite3_malloc64(sizeof(Rtree)+nDb+nName+2);
1239 if( !pRtree ){
1240 return SQLITE_NOMEM;
1241 }
1242 memset(pRtree, 0, sizeof(Rtree)+nDb+nName+2);
1243 pRtree->nBusy = 1;
1244 pRtree->base.pModule = &rtreeModule;
1245 pRtree->zDb = (char *)&pRtree[1];
1246 pRtree->zName = &pRtree->zDb[nDb+1];
1247 pRtree->eCoordType = RTREE_COORD_REAL32;
1248 pRtree->nDim = 2;
1249 pRtree->nDim2 = 4;
1250 memcpy(pRtree->zDb, argv[1], nDb);
1251 memcpy(pRtree->zName, argv[2], nName);
1252
1253
1254 /* Create/Connect to the underlying relational database schema. If
1255 ** that is successful, call sqlite3_declare_vtab() to configure
1256 ** the r-tree table schema.
1257 */
1258 pSql = sqlite3_str_new(db);
1259 sqlite3_str_appendf(pSql, "CREATE TABLE x(_shape");
1260 pRtree->nAux = 1; /* Add one for _shape */
1261 pRtree->nAuxNotNull = 1; /* The _shape column is always not-null */
1262 for(ii=3; ii<argc; ii++){
1263 pRtree->nAux++;
1264 sqlite3_str_appendf(pSql, ",%s", argv[ii]);
1265 }
1266 sqlite3_str_appendf(pSql, ");");
1267 zSql = sqlite3_str_finish(pSql);
1268 if( !zSql ){
1269 rc = SQLITE_NOMEM;
1270 }else if( SQLITE_OK!=(rc = sqlite3_declare_vtab(db, zSql)) ){
1271 *pzErr = sqlite3_mprintf("%s", sqlite3_errmsg(db));
1272 }
1273 sqlite3_free(zSql);
1274 if( rc ) goto geopolyInit_fail;
1275 pRtree->nBytesPerCell = 8 + pRtree->nDim2*4;
1276
1277 /* Figure out the node size to use. */
1278 rc = getNodeSize(db, pRtree, isCreate, pzErr);
1279 if( rc ) goto geopolyInit_fail;
1280 rc = rtreeSqlInit(pRtree, db, argv[1], argv[2], isCreate);
1281 if( rc ){
1282 *pzErr = sqlite3_mprintf("%s", sqlite3_errmsg(db));
1283 goto geopolyInit_fail;
1284 }
1285
1286 *ppVtab = (sqlite3_vtab *)pRtree;
1287 return SQLITE_OK;
1288
1289 geopolyInit_fail:
1290 if( rc==SQLITE_OK ) rc = SQLITE_ERROR;
1291 assert( *ppVtab==0 );
1292 assert( pRtree->nBusy==1 );
1293 rtreeRelease(pRtree);
1294 return rc;
1295 }
1296
1297
1298 /*
1299 ** GEOPOLY virtual table module xCreate method.
1300 */
geopolyCreate(sqlite3 * db,void * pAux,int argc,const char * const * argv,sqlite3_vtab ** ppVtab,char ** pzErr)1301 static int geopolyCreate(
1302 sqlite3 *db,
1303 void *pAux,
1304 int argc, const char *const*argv,
1305 sqlite3_vtab **ppVtab,
1306 char **pzErr
1307 ){
1308 return geopolyInit(db, pAux, argc, argv, ppVtab, pzErr, 1);
1309 }
1310
1311 /*
1312 ** GEOPOLY virtual table module xConnect method.
1313 */
geopolyConnect(sqlite3 * db,void * pAux,int argc,const char * const * argv,sqlite3_vtab ** ppVtab,char ** pzErr)1314 static int geopolyConnect(
1315 sqlite3 *db,
1316 void *pAux,
1317 int argc, const char *const*argv,
1318 sqlite3_vtab **ppVtab,
1319 char **pzErr
1320 ){
1321 return geopolyInit(db, pAux, argc, argv, ppVtab, pzErr, 0);
1322 }
1323
1324
1325 /*
1326 ** GEOPOLY virtual table module xFilter method.
1327 **
1328 ** Query plans:
1329 **
1330 ** 1 rowid lookup
1331 ** 2 search for objects overlapping the same bounding box
1332 ** that contains polygon argv[0]
1333 ** 3 search for objects overlapping the same bounding box
1334 ** that contains polygon argv[0]
1335 ** 4 full table scan
1336 */
geopolyFilter(sqlite3_vtab_cursor * pVtabCursor,int idxNum,const char * idxStr,int argc,sqlite3_value ** argv)1337 static int geopolyFilter(
1338 sqlite3_vtab_cursor *pVtabCursor, /* The cursor to initialize */
1339 int idxNum, /* Query plan */
1340 const char *idxStr, /* Not Used */
1341 int argc, sqlite3_value **argv /* Parameters to the query plan */
1342 ){
1343 Rtree *pRtree = (Rtree *)pVtabCursor->pVtab;
1344 RtreeCursor *pCsr = (RtreeCursor *)pVtabCursor;
1345 RtreeNode *pRoot = 0;
1346 int rc = SQLITE_OK;
1347 int iCell = 0;
1348
1349 rtreeReference(pRtree);
1350
1351 /* Reset the cursor to the same state as rtreeOpen() leaves it in. */
1352 resetCursor(pCsr);
1353
1354 pCsr->iStrategy = idxNum;
1355 if( idxNum==1 ){
1356 /* Special case - lookup by rowid. */
1357 RtreeNode *pLeaf; /* Leaf on which the required cell resides */
1358 RtreeSearchPoint *p; /* Search point for the leaf */
1359 i64 iRowid = sqlite3_value_int64(argv[0]);
1360 i64 iNode = 0;
1361 rc = findLeafNode(pRtree, iRowid, &pLeaf, &iNode);
1362 if( rc==SQLITE_OK && pLeaf!=0 ){
1363 p = rtreeSearchPointNew(pCsr, RTREE_ZERO, 0);
1364 assert( p!=0 ); /* Always returns pCsr->sPoint */
1365 pCsr->aNode[0] = pLeaf;
1366 p->id = iNode;
1367 p->eWithin = PARTLY_WITHIN;
1368 rc = nodeRowidIndex(pRtree, pLeaf, iRowid, &iCell);
1369 p->iCell = (u8)iCell;
1370 RTREE_QUEUE_TRACE(pCsr, "PUSH-F1:");
1371 }else{
1372 pCsr->atEOF = 1;
1373 }
1374 }else{
1375 /* Normal case - r-tree scan. Set up the RtreeCursor.aConstraint array
1376 ** with the configured constraints.
1377 */
1378 rc = nodeAcquire(pRtree, 1, 0, &pRoot);
1379 if( rc==SQLITE_OK && idxNum<=3 ){
1380 RtreeCoord bbox[4];
1381 RtreeConstraint *p;
1382 assert( argc==1 );
1383 geopolyBBox(0, argv[0], bbox, &rc);
1384 if( rc ){
1385 goto geopoly_filter_end;
1386 }
1387 pCsr->aConstraint = p = sqlite3_malloc(sizeof(RtreeConstraint)*4);
1388 pCsr->nConstraint = 4;
1389 if( p==0 ){
1390 rc = SQLITE_NOMEM;
1391 }else{
1392 memset(pCsr->aConstraint, 0, sizeof(RtreeConstraint)*4);
1393 memset(pCsr->anQueue, 0, sizeof(u32)*(pRtree->iDepth + 1));
1394 if( idxNum==2 ){
1395 /* Overlap query */
1396 p->op = 'B';
1397 p->iCoord = 0;
1398 p->u.rValue = bbox[1].f;
1399 p++;
1400 p->op = 'D';
1401 p->iCoord = 1;
1402 p->u.rValue = bbox[0].f;
1403 p++;
1404 p->op = 'B';
1405 p->iCoord = 2;
1406 p->u.rValue = bbox[3].f;
1407 p++;
1408 p->op = 'D';
1409 p->iCoord = 3;
1410 p->u.rValue = bbox[2].f;
1411 }else{
1412 /* Within query */
1413 p->op = 'D';
1414 p->iCoord = 0;
1415 p->u.rValue = bbox[0].f;
1416 p++;
1417 p->op = 'B';
1418 p->iCoord = 1;
1419 p->u.rValue = bbox[1].f;
1420 p++;
1421 p->op = 'D';
1422 p->iCoord = 2;
1423 p->u.rValue = bbox[2].f;
1424 p++;
1425 p->op = 'B';
1426 p->iCoord = 3;
1427 p->u.rValue = bbox[3].f;
1428 }
1429 }
1430 }
1431 if( rc==SQLITE_OK ){
1432 RtreeSearchPoint *pNew;
1433 pNew = rtreeSearchPointNew(pCsr, RTREE_ZERO, (u8)(pRtree->iDepth+1));
1434 if( pNew==0 ){
1435 rc = SQLITE_NOMEM;
1436 goto geopoly_filter_end;
1437 }
1438 pNew->id = 1;
1439 pNew->iCell = 0;
1440 pNew->eWithin = PARTLY_WITHIN;
1441 assert( pCsr->bPoint==1 );
1442 pCsr->aNode[0] = pRoot;
1443 pRoot = 0;
1444 RTREE_QUEUE_TRACE(pCsr, "PUSH-Fm:");
1445 rc = rtreeStepToLeaf(pCsr);
1446 }
1447 }
1448
1449 geopoly_filter_end:
1450 nodeRelease(pRtree, pRoot);
1451 rtreeRelease(pRtree);
1452 return rc;
1453 }
1454
1455 /*
1456 ** Rtree virtual table module xBestIndex method. There are three
1457 ** table scan strategies to choose from (in order from most to
1458 ** least desirable):
1459 **
1460 ** idxNum idxStr Strategy
1461 ** ------------------------------------------------
1462 ** 1 "rowid" Direct lookup by rowid.
1463 ** 2 "rtree" R-tree overlap query using geopoly_overlap()
1464 ** 3 "rtree" R-tree within query using geopoly_within()
1465 ** 4 "fullscan" full-table scan.
1466 ** ------------------------------------------------
1467 */
geopolyBestIndex(sqlite3_vtab * tab,sqlite3_index_info * pIdxInfo)1468 static int geopolyBestIndex(sqlite3_vtab *tab, sqlite3_index_info *pIdxInfo){
1469 int ii;
1470 int iRowidTerm = -1;
1471 int iFuncTerm = -1;
1472 int idxNum = 0;
1473
1474 for(ii=0; ii<pIdxInfo->nConstraint; ii++){
1475 struct sqlite3_index_constraint *p = &pIdxInfo->aConstraint[ii];
1476 if( !p->usable ) continue;
1477 if( p->iColumn<0 && p->op==SQLITE_INDEX_CONSTRAINT_EQ ){
1478 iRowidTerm = ii;
1479 break;
1480 }
1481 if( p->iColumn==0 && p->op>=SQLITE_INDEX_CONSTRAINT_FUNCTION ){
1482 /* p->op==SQLITE_INDEX_CONSTRAINT_FUNCTION for geopoly_overlap()
1483 ** p->op==(SQLITE_INDEX_CONTRAINT_FUNCTION+1) for geopoly_within().
1484 ** See geopolyFindFunction() */
1485 iFuncTerm = ii;
1486 idxNum = p->op - SQLITE_INDEX_CONSTRAINT_FUNCTION + 2;
1487 }
1488 }
1489
1490 if( iRowidTerm>=0 ){
1491 pIdxInfo->idxNum = 1;
1492 pIdxInfo->idxStr = "rowid";
1493 pIdxInfo->aConstraintUsage[iRowidTerm].argvIndex = 1;
1494 pIdxInfo->aConstraintUsage[iRowidTerm].omit = 1;
1495 pIdxInfo->estimatedCost = 30.0;
1496 pIdxInfo->estimatedRows = 1;
1497 pIdxInfo->idxFlags = SQLITE_INDEX_SCAN_UNIQUE;
1498 return SQLITE_OK;
1499 }
1500 if( iFuncTerm>=0 ){
1501 pIdxInfo->idxNum = idxNum;
1502 pIdxInfo->idxStr = "rtree";
1503 pIdxInfo->aConstraintUsage[iFuncTerm].argvIndex = 1;
1504 pIdxInfo->aConstraintUsage[iFuncTerm].omit = 0;
1505 pIdxInfo->estimatedCost = 300.0;
1506 pIdxInfo->estimatedRows = 10;
1507 return SQLITE_OK;
1508 }
1509 pIdxInfo->idxNum = 4;
1510 pIdxInfo->idxStr = "fullscan";
1511 pIdxInfo->estimatedCost = 3000000.0;
1512 pIdxInfo->estimatedRows = 100000;
1513 return SQLITE_OK;
1514 }
1515
1516
1517 /*
1518 ** GEOPOLY virtual table module xColumn method.
1519 */
geopolyColumn(sqlite3_vtab_cursor * cur,sqlite3_context * ctx,int i)1520 static int geopolyColumn(sqlite3_vtab_cursor *cur, sqlite3_context *ctx, int i){
1521 Rtree *pRtree = (Rtree *)cur->pVtab;
1522 RtreeCursor *pCsr = (RtreeCursor *)cur;
1523 RtreeSearchPoint *p = rtreeSearchPointFirst(pCsr);
1524 int rc = SQLITE_OK;
1525 RtreeNode *pNode = rtreeNodeOfFirstSearchPoint(pCsr, &rc);
1526
1527 if( rc ) return rc;
1528 if( p==0 ) return SQLITE_OK;
1529 if( i==0 && sqlite3_vtab_nochange(ctx) ) return SQLITE_OK;
1530 if( i<=pRtree->nAux ){
1531 if( !pCsr->bAuxValid ){
1532 if( pCsr->pReadAux==0 ){
1533 rc = sqlite3_prepare_v3(pRtree->db, pRtree->zReadAuxSql, -1, 0,
1534 &pCsr->pReadAux, 0);
1535 if( rc ) return rc;
1536 }
1537 sqlite3_bind_int64(pCsr->pReadAux, 1,
1538 nodeGetRowid(pRtree, pNode, p->iCell));
1539 rc = sqlite3_step(pCsr->pReadAux);
1540 if( rc==SQLITE_ROW ){
1541 pCsr->bAuxValid = 1;
1542 }else{
1543 sqlite3_reset(pCsr->pReadAux);
1544 if( rc==SQLITE_DONE ) rc = SQLITE_OK;
1545 return rc;
1546 }
1547 }
1548 sqlite3_result_value(ctx, sqlite3_column_value(pCsr->pReadAux, i+2));
1549 }
1550 return SQLITE_OK;
1551 }
1552
1553
1554 /*
1555 ** The xUpdate method for GEOPOLY module virtual tables.
1556 **
1557 ** For DELETE:
1558 **
1559 ** argv[0] = the rowid to be deleted
1560 **
1561 ** For INSERT:
1562 **
1563 ** argv[0] = SQL NULL
1564 ** argv[1] = rowid to insert, or an SQL NULL to select automatically
1565 ** argv[2] = _shape column
1566 ** argv[3] = first application-defined column....
1567 **
1568 ** For UPDATE:
1569 **
1570 ** argv[0] = rowid to modify. Never NULL
1571 ** argv[1] = rowid after the change. Never NULL
1572 ** argv[2] = new value for _shape
1573 ** argv[3] = new value for first application-defined column....
1574 */
geopolyUpdate(sqlite3_vtab * pVtab,int nData,sqlite3_value ** aData,sqlite_int64 * pRowid)1575 static int geopolyUpdate(
1576 sqlite3_vtab *pVtab,
1577 int nData,
1578 sqlite3_value **aData,
1579 sqlite_int64 *pRowid
1580 ){
1581 Rtree *pRtree = (Rtree *)pVtab;
1582 int rc = SQLITE_OK;
1583 RtreeCell cell; /* New cell to insert if nData>1 */
1584 i64 oldRowid; /* The old rowid */
1585 int oldRowidValid; /* True if oldRowid is valid */
1586 i64 newRowid; /* The new rowid */
1587 int newRowidValid; /* True if newRowid is valid */
1588 int coordChange = 0; /* Change in coordinates */
1589
1590 if( pRtree->nNodeRef ){
1591 /* Unable to write to the btree while another cursor is reading from it,
1592 ** since the write might do a rebalance which would disrupt the read
1593 ** cursor. */
1594 return SQLITE_LOCKED_VTAB;
1595 }
1596 rtreeReference(pRtree);
1597 assert(nData>=1);
1598
1599 oldRowidValid = sqlite3_value_type(aData[0])!=SQLITE_NULL;;
1600 oldRowid = oldRowidValid ? sqlite3_value_int64(aData[0]) : 0;
1601 newRowidValid = nData>1 && sqlite3_value_type(aData[1])!=SQLITE_NULL;
1602 newRowid = newRowidValid ? sqlite3_value_int64(aData[1]) : 0;
1603 cell.iRowid = newRowid;
1604
1605 if( nData>1 /* not a DELETE */
1606 && (!oldRowidValid /* INSERT */
1607 || !sqlite3_value_nochange(aData[2]) /* UPDATE _shape */
1608 || oldRowid!=newRowid) /* Rowid change */
1609 ){
1610 geopolyBBox(0, aData[2], cell.aCoord, &rc);
1611 if( rc ){
1612 if( rc==SQLITE_ERROR ){
1613 pVtab->zErrMsg =
1614 sqlite3_mprintf("_shape does not contain a valid polygon");
1615 }
1616 goto geopoly_update_end;
1617 }
1618 coordChange = 1;
1619
1620 /* If a rowid value was supplied, check if it is already present in
1621 ** the table. If so, the constraint has failed. */
1622 if( newRowidValid && (!oldRowidValid || oldRowid!=newRowid) ){
1623 int steprc;
1624 sqlite3_bind_int64(pRtree->pReadRowid, 1, cell.iRowid);
1625 steprc = sqlite3_step(pRtree->pReadRowid);
1626 rc = sqlite3_reset(pRtree->pReadRowid);
1627 if( SQLITE_ROW==steprc ){
1628 if( sqlite3_vtab_on_conflict(pRtree->db)==SQLITE_REPLACE ){
1629 rc = rtreeDeleteRowid(pRtree, cell.iRowid);
1630 }else{
1631 rc = rtreeConstraintError(pRtree, 0);
1632 }
1633 }
1634 }
1635 }
1636
1637 /* If aData[0] is not an SQL NULL value, it is the rowid of a
1638 ** record to delete from the r-tree table. The following block does
1639 ** just that.
1640 */
1641 if( rc==SQLITE_OK && (nData==1 || (coordChange && oldRowidValid)) ){
1642 rc = rtreeDeleteRowid(pRtree, oldRowid);
1643 }
1644
1645 /* If the aData[] array contains more than one element, elements
1646 ** (aData[2]..aData[argc-1]) contain a new record to insert into
1647 ** the r-tree structure.
1648 */
1649 if( rc==SQLITE_OK && nData>1 && coordChange ){
1650 /* Insert the new record into the r-tree */
1651 RtreeNode *pLeaf = 0;
1652 if( !newRowidValid ){
1653 rc = rtreeNewRowid(pRtree, &cell.iRowid);
1654 }
1655 *pRowid = cell.iRowid;
1656 if( rc==SQLITE_OK ){
1657 rc = ChooseLeaf(pRtree, &cell, 0, &pLeaf);
1658 }
1659 if( rc==SQLITE_OK ){
1660 int rc2;
1661 pRtree->iReinsertHeight = -1;
1662 rc = rtreeInsertCell(pRtree, pLeaf, &cell, 0);
1663 rc2 = nodeRelease(pRtree, pLeaf);
1664 if( rc==SQLITE_OK ){
1665 rc = rc2;
1666 }
1667 }
1668 }
1669
1670 /* Change the data */
1671 if( rc==SQLITE_OK && nData>1 ){
1672 sqlite3_stmt *pUp = pRtree->pWriteAux;
1673 int jj;
1674 int nChange = 0;
1675 sqlite3_bind_int64(pUp, 1, cell.iRowid);
1676 assert( pRtree->nAux>=1 );
1677 if( sqlite3_value_nochange(aData[2]) ){
1678 sqlite3_bind_null(pUp, 2);
1679 }else{
1680 GeoPoly *p = 0;
1681 if( sqlite3_value_type(aData[2])==SQLITE_TEXT
1682 && (p = geopolyFuncParam(0, aData[2], &rc))!=0
1683 && rc==SQLITE_OK
1684 ){
1685 sqlite3_bind_blob(pUp, 2, p->hdr, 4+8*p->nVertex, SQLITE_TRANSIENT);
1686 }else{
1687 sqlite3_bind_value(pUp, 2, aData[2]);
1688 }
1689 sqlite3_free(p);
1690 nChange = 1;
1691 }
1692 for(jj=1; jj<pRtree->nAux; jj++){
1693 nChange++;
1694 sqlite3_bind_value(pUp, jj+2, aData[jj+2]);
1695 }
1696 if( nChange ){
1697 sqlite3_step(pUp);
1698 rc = sqlite3_reset(pUp);
1699 }
1700 }
1701
1702 geopoly_update_end:
1703 rtreeRelease(pRtree);
1704 return rc;
1705 }
1706
1707 /*
1708 ** Report that geopoly_overlap() is an overloaded function suitable
1709 ** for use in xBestIndex.
1710 */
geopolyFindFunction(sqlite3_vtab * pVtab,int nArg,const char * zName,void (** pxFunc)(sqlite3_context *,int,sqlite3_value **),void ** ppArg)1711 static int geopolyFindFunction(
1712 sqlite3_vtab *pVtab,
1713 int nArg,
1714 const char *zName,
1715 void (**pxFunc)(sqlite3_context*,int,sqlite3_value**),
1716 void **ppArg
1717 ){
1718 if( sqlite3_stricmp(zName, "geopoly_overlap")==0 ){
1719 *pxFunc = geopolyOverlapFunc;
1720 *ppArg = 0;
1721 return SQLITE_INDEX_CONSTRAINT_FUNCTION;
1722 }
1723 if( sqlite3_stricmp(zName, "geopoly_within")==0 ){
1724 *pxFunc = geopolyWithinFunc;
1725 *ppArg = 0;
1726 return SQLITE_INDEX_CONSTRAINT_FUNCTION+1;
1727 }
1728 return 0;
1729 }
1730
1731
1732 static sqlite3_module geopolyModule = {
1733 3, /* iVersion */
1734 geopolyCreate, /* xCreate - create a table */
1735 geopolyConnect, /* xConnect - connect to an existing table */
1736 geopolyBestIndex, /* xBestIndex - Determine search strategy */
1737 rtreeDisconnect, /* xDisconnect - Disconnect from a table */
1738 rtreeDestroy, /* xDestroy - Drop a table */
1739 rtreeOpen, /* xOpen - open a cursor */
1740 rtreeClose, /* xClose - close a cursor */
1741 geopolyFilter, /* xFilter - configure scan constraints */
1742 rtreeNext, /* xNext - advance a cursor */
1743 rtreeEof, /* xEof */
1744 geopolyColumn, /* xColumn - read data */
1745 rtreeRowid, /* xRowid - read data */
1746 geopolyUpdate, /* xUpdate - write data */
1747 rtreeBeginTransaction, /* xBegin - begin transaction */
1748 rtreeEndTransaction, /* xSync - sync transaction */
1749 rtreeEndTransaction, /* xCommit - commit transaction */
1750 rtreeEndTransaction, /* xRollback - rollback transaction */
1751 geopolyFindFunction, /* xFindFunction - function overloading */
1752 rtreeRename, /* xRename - rename the table */
1753 rtreeSavepoint, /* xSavepoint */
1754 0, /* xRelease */
1755 0, /* xRollbackTo */
1756 rtreeShadowName /* xShadowName */
1757 };
1758
sqlite3_geopoly_init(sqlite3 * db)1759 static int sqlite3_geopoly_init(sqlite3 *db){
1760 int rc = SQLITE_OK;
1761 static const struct {
1762 void (*xFunc)(sqlite3_context*,int,sqlite3_value**);
1763 signed char nArg;
1764 unsigned char bPure;
1765 const char *zName;
1766 } aFunc[] = {
1767 { geopolyAreaFunc, 1, 1, "geopoly_area" },
1768 { geopolyBlobFunc, 1, 1, "geopoly_blob" },
1769 { geopolyJsonFunc, 1, 1, "geopoly_json" },
1770 { geopolySvgFunc, -1, 1, "geopoly_svg" },
1771 { geopolyWithinFunc, 2, 1, "geopoly_within" },
1772 { geopolyContainsPointFunc, 3, 1, "geopoly_contains_point" },
1773 { geopolyOverlapFunc, 2, 1, "geopoly_overlap" },
1774 { geopolyDebugFunc, 1, 0, "geopoly_debug" },
1775 { geopolyBBoxFunc, 1, 1, "geopoly_bbox" },
1776 { geopolyXformFunc, 7, 1, "geopoly_xform" },
1777 { geopolyRegularFunc, 4, 1, "geopoly_regular" },
1778 { geopolyCcwFunc, 1, 1, "geopoly_ccw" },
1779 };
1780 static const struct {
1781 void (*xStep)(sqlite3_context*,int,sqlite3_value**);
1782 void (*xFinal)(sqlite3_context*);
1783 const char *zName;
1784 } aAgg[] = {
1785 { geopolyBBoxStep, geopolyBBoxFinal, "geopoly_group_bbox" },
1786 };
1787 int i;
1788 for(i=0; i<sizeof(aFunc)/sizeof(aFunc[0]) && rc==SQLITE_OK; i++){
1789 int enc;
1790 if( aFunc[i].bPure ){
1791 enc = SQLITE_UTF8|SQLITE_DETERMINISTIC|SQLITE_INNOCUOUS;
1792 }else{
1793 enc = SQLITE_UTF8|SQLITE_DIRECTONLY;
1794 }
1795 rc = sqlite3_create_function(db, aFunc[i].zName, aFunc[i].nArg,
1796 enc, 0,
1797 aFunc[i].xFunc, 0, 0);
1798 }
1799 for(i=0; i<sizeof(aAgg)/sizeof(aAgg[0]) && rc==SQLITE_OK; i++){
1800 rc = sqlite3_create_function(db, aAgg[i].zName, 1,
1801 SQLITE_UTF8|SQLITE_DETERMINISTIC|SQLITE_INNOCUOUS, 0,
1802 0, aAgg[i].xStep, aAgg[i].xFinal);
1803 }
1804 if( rc==SQLITE_OK ){
1805 rc = sqlite3_create_module_v2(db, "geopoly", &geopolyModule, 0, 0);
1806 }
1807 return rc;
1808 }
1809