PostGIS  2.2.8dev-r@@SVN_REVISION@@
gserialized_gist_2d.c
Go to the documentation of this file.
1 /**********************************************************************
2  *
3  * PostGIS - Spatial Types for PostgreSQL
4  * Copyright 2009 Paul Ramsey <pramsey@cleverelephant.ca>
5  *
6  * This is free software; you can redistribute and/or modify it under
7  * the terms of the GNU General Public Licence. See the COPYING file.
8  *
9  **********************************************************************/
10 
11 /*
12 ** R-Tree Bibliography
13 **
14 ** [1] A. Guttman. R-tree: a dynamic index structure for spatial searching.
15 ** Proceedings of the ACM SIGMOD Conference, pp 47-57, June 1984.
16 ** [2] C.-H. Ang and T. C. Tan. New linear node splitting algorithm for
17 ** R-Trees. Advances in Spatial Databases - 5th International Symposium,
18 ** 1997
19 ** [3] N. Beckmann, H.-P. Kriegel, R. Schneider, B. Seeger. The R*tree: an
20 ** efficient and robust access method for points and rectangles.
21 ** Proceedings of the ACM SIGMOD Conference. June 1990.
22 ** [4] A. Korotkov, "A new double sorting-based node splitting algorithm for R-tree",
23 ** http://syrcose.ispras.ru/2011/files/SYRCoSE2011_Proceedings.pdf#page=36
24 */
25 
26 #include "postgres.h"
27 #include "access/gist.h" /* For GiST */
28 #include "access/itup.h"
29 #include "access/skey.h"
30 
31 #include "../postgis_config.h"
32 
33 #include "liblwgeom.h" /* For standard geometry types. */
34 #include "lwgeom_pg.h" /* For debugging macros. */
35 #include "gserialized_gist.h" /* For utility functions. */
36 
37 /* Fall back to older finite() if necessary */
38 #ifndef HAVE_ISFINITE
39 # ifdef HAVE_GNU_ISFINITE
40 # define _GNU_SOURCE
41 # else
42 # define isfinite finite
43 # endif
44 #endif
45 
46 #include <float.h> /* For FLT_MAX */
47 
48 /*
49 ** When is a node split not so good? If more than 90% of the entries
50 ** end up in one of the children.
51 */
52 #define LIMIT_RATIO 0.1
53 
54 /*
55 ** 0 == don't use it
56 ** 1 == use it
57 */
58 #define KOROTKOV_SPLIT 1
59 
60 /*
61 ** For debugging
62 */
63 #if POSTGIS_DEBUG_LEVEL > 0
64 static int g2d_counter_leaf = 0;
65 static int g2d_counter_internal = 0;
66 #endif
67 
68 /*
69 ** GiST 2D key stubs
70 */
71 Datum box2df_out(PG_FUNCTION_ARGS);
72 Datum box2df_in(PG_FUNCTION_ARGS);
73 
74 /*
75 ** GiST 2D index function prototypes
76 */
77 Datum gserialized_gist_consistent_2d(PG_FUNCTION_ARGS);
78 Datum gserialized_gist_compress_2d(PG_FUNCTION_ARGS);
79 Datum gserialized_gist_decompress_2d(PG_FUNCTION_ARGS);
80 Datum gserialized_gist_penalty_2d(PG_FUNCTION_ARGS);
81 Datum gserialized_gist_picksplit_2d(PG_FUNCTION_ARGS);
82 Datum gserialized_gist_union_2d(PG_FUNCTION_ARGS);
83 Datum gserialized_gist_same_2d(PG_FUNCTION_ARGS);
84 Datum gserialized_gist_distance_2d(PG_FUNCTION_ARGS);
85 
86 /*
87 ** GiST 2D operator prototypes
88 */
89 Datum gserialized_same_2d(PG_FUNCTION_ARGS);
90 Datum gserialized_within_2d(PG_FUNCTION_ARGS);
91 Datum gserialized_contains_2d(PG_FUNCTION_ARGS);
92 Datum gserialized_overlaps_2d(PG_FUNCTION_ARGS);
93 Datum gserialized_left_2d(PG_FUNCTION_ARGS);
94 Datum gserialized_right_2d(PG_FUNCTION_ARGS);
95 Datum gserialized_above_2d(PG_FUNCTION_ARGS);
96 Datum gserialized_below_2d(PG_FUNCTION_ARGS);
97 Datum gserialized_overleft_2d(PG_FUNCTION_ARGS);
98 Datum gserialized_overright_2d(PG_FUNCTION_ARGS);
99 Datum gserialized_overabove_2d(PG_FUNCTION_ARGS);
100 Datum gserialized_overbelow_2d(PG_FUNCTION_ARGS);
101 Datum gserialized_distance_box_2d(PG_FUNCTION_ARGS);
102 Datum gserialized_distance_centroid_2d(PG_FUNCTION_ARGS);
103 
104 /*
105 ** true/false test function type
106 */
107 typedef bool (*box2df_predicate)(const BOX2DF *a, const BOX2DF *b);
108 
109 
110 static char* box2df_to_string(const BOX2DF *a)
111 {
112  char *rv = NULL;
113 
114  if ( a == NULL )
115  return pstrdup("<NULLPTR>");
116 
117  rv = palloc(128);
118  sprintf(rv, "BOX2DF(%.12g %.12g, %.12g %.12g)", a->xmin, a->ymin, a->xmax, a->ymax);
119  return rv;
120 }
121 
122 /* Allocate a new copy of BOX2DF */
123 static BOX2DF* box2df_copy(BOX2DF *b)
124 {
125  BOX2DF *c = (BOX2DF*)palloc(sizeof(BOX2DF));
126  memcpy((void*)c, (void*)b, sizeof(BOX2DF));
127  POSTGIS_DEBUGF(5, "copied box2df (%p) to box2df (%p)", b, c);
128  return c;
129 }
130 
131 
132 
133 /* Enlarge b_union to contain b_new. If b_new contains more
134  dimensions than b_union, expand b_union to contain those dimensions. */
135 static void box2df_merge(BOX2DF *b_union, BOX2DF *b_new)
136 {
137 
138  POSTGIS_DEBUGF(5, "merging %s with %s", box2df_to_string(b_union), box2df_to_string(b_new));
139  /* Adjust minimums */
140  if (b_union->xmin > b_new->xmin || isnan(b_union->xmin))
141  b_union->xmin = b_new->xmin;
142  if (b_union->ymin > b_new->ymin || isnan(b_union->ymin))
143  b_union->ymin = b_new->ymin;
144  /* Adjust maximums */
145  if (b_union->xmax < b_new->xmax || isnan(b_union->xmax))
146  b_union->xmax = b_new->xmax;
147  if (b_union->ymax < b_new->ymax || isnan(b_union->ymax))
148  b_union->ymax = b_new->ymax;
149 
150  POSTGIS_DEBUGF(5, "merge complete %s", box2df_to_string(b_union));
151  return;
152 }
153 
154 #if KOROTKOV_SPLIT < 1
155 static bool box2df_intersection(const BOX2DF *a, const BOX2DF *b, BOX2DF *n)
156 {
157  POSTGIS_DEBUGF(5, "calculating intersection of %s with %s", box2df_to_string(a), box2df_to_string(b));
158 
159  if( a == NULL || b == NULL || n == NULL )
160  return FALSE;
161 
162  n->xmax = Min(a->xmax, b->xmax);
163  n->ymax = Min(a->ymax, b->ymax);
164  n->xmin = Max(a->xmin, b->xmin);
165  n->ymin = Max(a->ymin, b->ymin);
166 
167  POSTGIS_DEBUGF(5, "intersection is %s", box2df_to_string(n));
168 
169  if ( (n->xmax < n->xmin) || (n->ymax < n->ymin) )
170  return FALSE;
171 
172  return TRUE;
173 }
174 #endif
175 
176 static float box2df_size(const BOX2DF *a)
177 {
178  float result;
179 
180  if ( a == NULL )
181  return (float)0.0;
182 
183  if ( (a->xmax <= a->xmin) || (a->ymax <= a->ymin) )
184  {
185  result = (float) 0.0;
186  }
187  else
188  {
189  result = (((double) a->xmax)-((double) a->xmin)) * (((double) a->ymax)-((double) a->ymin));
190  }
191 
192  return result;
193 }
194 
195 static float box2df_union_size(const BOX2DF *a, const BOX2DF *b)
196 {
197  float result;
198 
199  POSTGIS_DEBUG(5,"entered function");
200 
201  if ( a == NULL && b == NULL )
202  {
203  elog(ERROR, "box2df_union_size received two null arguments");
204  return 0.0;
205  }
206 
207  if ( a == NULL )
208  return box2df_size(b);
209 
210  if ( b == NULL )
211  return box2df_size(a);
212 
213  result = ((double)Max(a->xmax,b->xmax) - (double)Min(a->xmin,b->xmin)) *
214  ((double)Max(a->ymax,b->ymax) - (double)Min(a->ymin,b->ymin));
215 
216  POSTGIS_DEBUGF(5, "union size of %s and %s is %.8g", box2df_to_string(a), box2df_to_string(b), result);
217 
218  return result;
219 }
220 
221 
222 /* Convert a double-based GBOX into a float-based BOX2DF,
223  ensuring the float box is larger than the double box */
224 static inline int box2df_from_gbox_p(GBOX *box, BOX2DF *a)
225 {
226  a->xmin = next_float_down(box->xmin);
227  a->xmax = next_float_up(box->xmax);
228  a->ymin = next_float_down(box->ymin);
229  a->ymax = next_float_up(box->ymax);
230  return LW_SUCCESS;
231 }
232 
233 /***********************************************************************
234 ** BOX3DF tests for 2D index operators.
235 */
236 
237 /* Ensure all minimums are below maximums. */
238 static inline void box2df_validate(BOX2DF *b)
239 {
240  float tmp;
241  POSTGIS_DEBUGF(5,"validating box2df (%s)", box2df_to_string(b));
242  if ( b->xmax < b->xmin )
243  {
244  tmp = b->xmin;
245  b->xmin = b->xmax;
246  b->xmax = tmp;
247  }
248  if ( b->ymax < b->ymin )
249  {
250  tmp = b->ymin;
251  b->ymin = b->ymax;
252  b->ymax = tmp;
253  }
254  return;
255 }
256 
257 static bool box2df_overlaps(const BOX2DF *a, const BOX2DF *b)
258 {
259  if ( ! a || ! b ) return FALSE; /* TODO: might be smarter for EMPTY */
260 
261  if ( (a->xmin > b->xmax) || (b->xmin > a->xmax) ||
262  (a->ymin > b->ymax) || (b->ymin > a->ymax) )
263  {
264  return FALSE;
265  }
266 
267  return TRUE;
268 }
269 
270 static bool box2df_contains(const BOX2DF *a, const BOX2DF *b)
271 {
272  if ( ! a || ! b ) return FALSE; /* TODO: might be smarter for EMPTY */
273 
274  if ( (a->xmin > b->xmin) || (a->xmax < b->xmax) ||
275  (a->ymin > b->ymin) || (a->ymax < b->ymax) )
276  {
277  return FALSE;
278  }
279 
280  return TRUE;
281 }
282 
283 static bool box2df_within(const BOX2DF *a, const BOX2DF *b)
284 {
285  if ( ! a || ! b ) return FALSE; /* TODO: might be smarter for EMPTY */
286 
287  POSTGIS_DEBUG(5, "entered function");
288  return box2df_contains(b,a);
289 }
290 
291 static bool box2df_equals(const BOX2DF *a, const BOX2DF *b)
292 {
293  if ( a && b ) {
294  if ( (a->xmin != b->xmin) || (a->xmax != b->xmax) ||
295  (a->ymin != b->ymin) || (a->ymax != b->ymax) )
296  {
297  return FALSE;
298  }
299  return TRUE;
300  } else if ( a || b ) {
301  /* one empty, one not */
302  return FALSE;
303  } else {
304  /* both empty */
305  return TRUE;
306  }
307 }
308 
309 static bool box2df_overleft(const BOX2DF *a, const BOX2DF *b)
310 {
311  if ( ! a || ! b ) return FALSE; /* TODO: might be smarter for EMPTY */
312 
313  /* a.xmax <= b.xmax */
314  return a->xmax <= b->xmax;
315 }
316 
317 static bool box2df_left(const BOX2DF *a, const BOX2DF *b)
318 {
319  if ( ! a || ! b ) return FALSE; /* TODO: might be smarter for EMPTY */
320 
321  /* a.xmax < b.xmin */
322  return a->xmax < b->xmin;
323 }
324 
325 static bool box2df_right(const BOX2DF *a, const BOX2DF *b)
326 {
327  if ( ! a || ! b ) return FALSE; /* TODO: might be smarter for EMPTY */
328 
329  /* a.xmin > b.xmax */
330  return a->xmin > b->xmax;
331 }
332 
333 static bool box2df_overright(const BOX2DF *a, const BOX2DF *b)
334 {
335  if ( ! a || ! b ) return FALSE; /* TODO: might be smarter for EMPTY */
336 
337  /* a.xmin >= b.xmin */
338  return a->xmin >= b->xmin;
339 }
340 
341 static bool box2df_overbelow(const BOX2DF *a, const BOX2DF *b)
342 {
343  if ( ! a || ! b ) return FALSE; /* TODO: might be smarter for EMPTY */
344 
345  /* a.ymax <= b.ymax */
346  return a->ymax <= b->ymax;
347 }
348 
349 static bool box2df_below(const BOX2DF *a, const BOX2DF *b)
350 {
351  if ( ! a || ! b ) return FALSE; /* TODO: might be smarter for EMPTY */
352 
353  /* a.ymax < b.ymin */
354  return a->ymax < b->ymin;
355 }
356 
357 static bool box2df_above(const BOX2DF *a, const BOX2DF *b)
358 {
359  if ( ! a || ! b ) return FALSE; /* TODO: might be smarter for EMPTY */
360 
361  /* a.ymin > b.ymax */
362  return a->ymin > b->ymax;
363 }
364 
365 static bool box2df_overabove(const BOX2DF *a, const BOX2DF *b)
366 {
367  if ( ! a || ! b ) return FALSE; /* TODO: might be smarter for EMPTY */
368 
369  /* a.ymin >= b.ymin */
370  return a->ymin >= b->ymin;
371 }
372 
376 static double box2df_distance_leaf_centroid(const BOX2DF *a, const BOX2DF *b)
377 {
378  /* The centroid->centroid distance between the boxes */
379  double a_x = (a->xmax + a->xmin) / 2.0;
380  double a_y = (a->ymax + a->ymin) / 2.0;
381  double b_x = (b->xmax + b->xmin) / 2.0;
382  double b_y = (b->ymax + b->ymin) / 2.0;
383 
384  /* This "distance" is only used for comparisons, */
385  return sqrt((a_x - b_x) * (a_x - b_x) + (a_y - b_y) * (a_y - b_y));
386 }
387 
388 #if POSTGIS_PGSQL_VERSION < 95
389 
393 static double box2df_distance_node_centroid(const BOX2DF *node, const BOX2DF *query)
394 {
395  BOX2DF q;
396  double qx, qy;
397  double d = 0.0;
398 
399  /* Turn query into point */
400  q.xmin = q.xmax = (query->xmin + query->xmax) / 2.0;
401  q.ymin = q.ymax = (query->ymin + query->ymax) / 2.0;
402  qx = q.xmin;
403  qy = q.ymin;
404 
405  /* Check for overlap */
406  if ( box2df_overlaps(node, &q) == LW_TRUE )
407  return 0.0;
408 
409  /* Above or below */
410  if ( qx >= node->xmin && qx <= node->xmax )
411  {
412  if( qy > node->ymax )
413  d = qy - node->ymax;
414  else if ( qy < node->ymin )
415  d = node->ymin - qy;
416  return d;
417  }
418  /* Left or right */
419  else if ( qy >= node->ymin && qy <= node->ymax )
420  {
421  if ( qx > node->xmax )
422  d = qx - node->xmax;
423  else if ( qx < node->xmin )
424  d = node->xmin - qx;
425  return d;
426  }
427  /* Corner quadrants */
428  else
429  {
430  /* below/left of xmin/ymin */
431  if ( qx < node->xmin && qy < node->ymin )
432  {
433  d = (node->xmin - qx) * (node->xmin - qx) +
434  (node->ymin - qy) * (node->ymin - qy);
435  }
436  /* above/left of xmin/ymax */
437  else if ( qx < node->xmin && qy > node->ymax )
438  {
439  d = (node->xmin - qx) * (node->xmin - qx) +
440  (node->ymax - qy) * (node->ymax - qy);
441  }
442  /* above/right of xmax/ymax */
443  else if ( qx > node->xmax && qy > node->ymax )
444  {
445  d = (node->xmax - qx) * (node->xmax - qx) +
446  (node->ymax - qy) * (node->ymax - qy);
447  }
448  /* below/right of xmax/ymin */
449  else if ( qx > node->xmin && qy < node->ymin )
450  {
451  d = (node->xmax - qx) * (node->xmax - qx) +
452  (node->ymin - qy) * (node->ymin - qy);
453  }
454  else
455  {
456  /*ERROR*/
457  elog(ERROR, "%s: reached unreachable code", __func__);
458  }
459  }
460 
461  return sqrt(d);
462 }
463 #endif
464 
465 /* Quick distance function */
466 static inline double pt_distance(double ax, double ay, double bx, double by)
467 {
468  return sqrt((ax - bx) * (ax - bx) + (ay - by) * (ay - by));
469 }
470 
474 static double box2df_distance(const BOX2DF *a, const BOX2DF *b)
475 {
476  /* Check for overlap */
477  if ( box2df_overlaps(a, b) )
478  return 0.0;
479 
480  if ( box2df_left(a, b) )
481  {
482  if ( box2df_above(a, b) )
483  return pt_distance(a->xmax, a->ymin, b->xmin, b->ymax);
484  if ( box2df_below(a, b) )
485  return pt_distance(a->xmax, a->ymax, b->xmin, b->ymin);
486  else
487  return (double)b->xmin - (double)a->xmax;
488  }
489  if ( box2df_right(a, b) )
490  {
491  if ( box2df_above(a, b) )
492  return pt_distance(a->xmin, a->ymin, b->xmax, b->ymax);
493  if ( box2df_below(a, b) )
494  return pt_distance(a->xmin, a->ymax, b->xmax, b->ymin);
495  else
496  return (double)a->xmin - (double)b->xmax;
497  }
498  if ( box2df_above(a, b) )
499  {
500  if ( box2df_left(a, b) )
501  return pt_distance(a->xmax, a->ymin, b->xmin, b->ymax);
502  if ( box2df_right(a, b) )
503  return pt_distance(a->xmin, a->ymin, b->xmax, b->ymax);
504  else
505  return (double)a->ymin - (double)b->ymax;
506  }
507  if ( box2df_below(a, b) )
508  {
509  if ( box2df_left(a, b) )
510  return pt_distance(a->xmax, a->ymax, b->xmin, b->ymin);
511  if ( box2df_right(a, b) )
512  return pt_distance(a->xmin, a->ymax, b->xmax, b->ymin);
513  else
514  return (double)b->ymin - (double)a->ymax;
515  }
516 
517  return FLT_MAX;
518 }
519 
520 
527 static int
528 gserialized_datum_get_box2df_p(Datum gsdatum, BOX2DF *box2df)
529 {
530  GSERIALIZED *gpart;
531  uint8_t flags;
532  int result = LW_SUCCESS;
533 
534  POSTGIS_DEBUG(4, "entered function");
535 
536  /*
537  ** The most info we need is the 8 bytes of serialized header plus the
538  ** of floats necessary to hold the bounding box.
539  */
540  if (VARATT_IS_EXTENDED(gsdatum))
541  {
542  gpart = (GSERIALIZED*)PG_DETOAST_DATUM_SLICE(gsdatum, 0, 8 + sizeof(BOX2DF));
543  }
544  else
545  {
546  gpart = (GSERIALIZED*)PG_DETOAST_DATUM(gsdatum);
547  }
548 
549  flags = gpart->flags;
550 
551  POSTGIS_DEBUGF(4, "got flags %d", gpart->flags);
552 
553  /* Do we even have a serialized bounding box? */
554  if ( FLAGS_GET_BBOX(flags) )
555  {
556  /* Yes! Copy it out into the box! */
557  POSTGIS_DEBUG(4, "copying box out of serialization");
558  memcpy(box2df, gpart->data, sizeof(BOX2DF));
559  result = LW_SUCCESS;
560  }
561  else
562  {
563  /* No, we need to calculate it from the full object. */
564  GBOX gbox;
565  GSERIALIZED *g = (GSERIALIZED*)PG_DETOAST_DATUM(gsdatum);
566  LWGEOM *lwgeom = lwgeom_from_gserialized(g);
567  if ( lwgeom_calculate_gbox(lwgeom, &gbox) == LW_FAILURE )
568  {
569  POSTGIS_DEBUG(4, "could not calculate bbox, returning failure");
570  lwgeom_free(lwgeom);
571  return LW_FAILURE;
572  }
573  lwgeom_free(lwgeom);
574  result = box2df_from_gbox_p(&gbox, box2df);
575  }
576 
577  if ( result == LW_SUCCESS )
578  {
579  POSTGIS_DEBUGF(4, "got box2df %s", box2df_to_string(box2df));
580  }
581 
582  return result;
583 }
584 
585 
590 static int
591 gserialized_datum_predicate_2d(Datum gs1, Datum gs2, box2df_predicate predicate)
592 {
593  BOX2DF b1, b2, *br1=NULL, *br2=NULL;
594  POSTGIS_DEBUG(3, "entered function");
595 
596  if (gserialized_datum_get_box2df_p(gs1, &b1) == LW_SUCCESS) br1 = &b1;
597  if (gserialized_datum_get_box2df_p(gs2, &b2) == LW_SUCCESS) br2 = &b2;
598 
599  if ( predicate(br1, br2) )
600  {
601  POSTGIS_DEBUGF(3, "got boxes %s and %s", br1 ? box2df_to_string(&b1) : "(null)", br2 ? box2df_to_string(&b2) : "(null)");
602  return LW_TRUE;
603  }
604  return LW_FALSE;
605 }
606 
607 
608 
609 
610 /***********************************************************************
611 * GiST 2-D Index Operator Functions
612 */
613 
615 Datum gserialized_distance_centroid_2d(PG_FUNCTION_ARGS)
616 {
617  BOX2DF b1, b2;
618  Datum gs1 = PG_GETARG_DATUM(0);
619  Datum gs2 = PG_GETARG_DATUM(1);
620 
621  POSTGIS_DEBUG(3, "entered function");
622 
623  /* Must be able to build box for each argument (ie, not empty geometry). */
624  if ( (gserialized_datum_get_box2df_p(gs1, &b1) == LW_SUCCESS) &&
626  {
627  double distance = box2df_distance_leaf_centroid(&b1, &b2);
628  POSTGIS_DEBUGF(3, "got boxes %s and %s", box2df_to_string(&b1), box2df_to_string(&b2));
629  PG_RETURN_FLOAT8(distance);
630  }
631  PG_RETURN_FLOAT8(FLT_MAX);
632 }
633 
635 Datum gserialized_distance_box_2d(PG_FUNCTION_ARGS)
636 {
637  BOX2DF b1, b2;
638  Datum gs1 = PG_GETARG_DATUM(0);
639  Datum gs2 = PG_GETARG_DATUM(1);
640 
641  POSTGIS_DEBUG(3, "entered function");
642 
643  /* Must be able to build box for each argument (ie, not empty geometry). */
644  if ( (gserialized_datum_get_box2df_p(gs1, &b1) == LW_SUCCESS) &&
646  {
647  double distance = box2df_distance(&b1, &b2);
648  POSTGIS_DEBUGF(3, "got boxes %s and %s", box2df_to_string(&b1), box2df_to_string(&b2));
649  PG_RETURN_FLOAT8(distance);
650  }
651  PG_RETURN_FLOAT8(FLT_MAX);
652 }
653 
655 Datum gserialized_same_2d(PG_FUNCTION_ARGS)
656 {
657  if ( gserialized_datum_predicate_2d(PG_GETARG_DATUM(0), PG_GETARG_DATUM(1), box2df_equals) == LW_TRUE )
658  PG_RETURN_BOOL(TRUE);
659 
660  PG_RETURN_BOOL(FALSE);
661 }
662 
664 Datum gserialized_within_2d(PG_FUNCTION_ARGS)
665 {
666  POSTGIS_DEBUG(3, "entered function");
667  if ( gserialized_datum_predicate_2d(PG_GETARG_DATUM(0), PG_GETARG_DATUM(1), box2df_within) == LW_TRUE )
668  PG_RETURN_BOOL(TRUE);
669 
670  PG_RETURN_BOOL(FALSE);
671 }
672 
674 Datum gserialized_contains_2d(PG_FUNCTION_ARGS)
675 {
676  POSTGIS_DEBUG(3, "entered function");
677  if ( gserialized_datum_predicate_2d(PG_GETARG_DATUM(0), PG_GETARG_DATUM(1), box2df_contains) == LW_TRUE )
678  PG_RETURN_BOOL(TRUE);
679 
680  PG_RETURN_BOOL(FALSE);
681 }
682 
684 Datum gserialized_overlaps_2d(PG_FUNCTION_ARGS)
685 {
686  if ( gserialized_datum_predicate_2d(PG_GETARG_DATUM(0), PG_GETARG_DATUM(1), box2df_overlaps) == LW_TRUE )
687  PG_RETURN_BOOL(TRUE);
688 
689  PG_RETURN_BOOL(FALSE);
690 }
691 
693 Datum gserialized_left_2d(PG_FUNCTION_ARGS)
694 {
695  if ( gserialized_datum_predicate_2d(PG_GETARG_DATUM(0), PG_GETARG_DATUM(1), box2df_left) == LW_TRUE )
696  PG_RETURN_BOOL(TRUE);
697 
698  PG_RETURN_BOOL(FALSE);
699 }
700 
702 Datum gserialized_right_2d(PG_FUNCTION_ARGS)
703 {
704  if ( gserialized_datum_predicate_2d(PG_GETARG_DATUM(0), PG_GETARG_DATUM(1), box2df_right) == LW_TRUE )
705  PG_RETURN_BOOL(TRUE);
706 
707  PG_RETURN_BOOL(FALSE);
708 }
709 
711 Datum gserialized_above_2d(PG_FUNCTION_ARGS)
712 {
713  if ( gserialized_datum_predicate_2d(PG_GETARG_DATUM(0), PG_GETARG_DATUM(1), box2df_above) == LW_TRUE )
714  PG_RETURN_BOOL(TRUE);
715 
716  PG_RETURN_BOOL(FALSE);
717 }
718 
720 Datum gserialized_below_2d(PG_FUNCTION_ARGS)
721 {
722  if ( gserialized_datum_predicate_2d(PG_GETARG_DATUM(0), PG_GETARG_DATUM(1), box2df_below) == LW_TRUE )
723  PG_RETURN_BOOL(TRUE);
724 
725  PG_RETURN_BOOL(FALSE);
726 }
727 
729 Datum gserialized_overleft_2d(PG_FUNCTION_ARGS)
730 {
731  if ( gserialized_datum_predicate_2d(PG_GETARG_DATUM(0), PG_GETARG_DATUM(1), box2df_overleft) == LW_TRUE )
732  PG_RETURN_BOOL(TRUE);
733 
734  PG_RETURN_BOOL(FALSE);
735 }
736 
738 Datum gserialized_overright_2d(PG_FUNCTION_ARGS)
739 {
740  if ( gserialized_datum_predicate_2d(PG_GETARG_DATUM(0), PG_GETARG_DATUM(1), box2df_overright) == LW_TRUE )
741  PG_RETURN_BOOL(TRUE);
742 
743  PG_RETURN_BOOL(FALSE);
744 }
745 
747 Datum gserialized_overabove_2d(PG_FUNCTION_ARGS)
748 {
749  if ( gserialized_datum_predicate_2d(PG_GETARG_DATUM(0), PG_GETARG_DATUM(1), box2df_overabove) == LW_TRUE )
750  PG_RETURN_BOOL(TRUE);
751 
752  PG_RETURN_BOOL(FALSE);
753 }
754 
756 Datum gserialized_overbelow_2d(PG_FUNCTION_ARGS)
757 {
758  if ( gserialized_datum_predicate_2d(PG_GETARG_DATUM(0), PG_GETARG_DATUM(1), box2df_overbelow) == LW_TRUE )
759  PG_RETURN_BOOL(TRUE);
760 
761  PG_RETURN_BOOL(FALSE);
762 }
763 
764 
765 /***********************************************************************
766 * GiST Index Support Functions
767 */
768 
769 /*
770 ** GiST support function. Given a geography, return a "compressed"
771 ** version. In this case, we convert the geography into a geocentric
772 ** bounding box. If the geography already has the box embedded in it
773 ** we pull that out and hand it back.
774 */
776 Datum gserialized_gist_compress_2d(PG_FUNCTION_ARGS)
777 {
778  GISTENTRY *entry_in = (GISTENTRY*)PG_GETARG_POINTER(0);
779  GISTENTRY *entry_out = NULL;
780  BOX2DF bbox_out;
781  int result = LW_SUCCESS;
782 
783  POSTGIS_DEBUG(4, "[GIST] 'compress' function called");
784 
785  /*
786  ** Not a leaf key? There's nothing to do.
787  ** Return the input unchanged.
788  */
789  if ( ! entry_in->leafkey )
790  {
791  POSTGIS_DEBUG(4, "[GIST] non-leafkey entry, returning input unaltered");
792  PG_RETURN_POINTER(entry_in);
793  }
794 
795  POSTGIS_DEBUG(4, "[GIST] processing leafkey input");
796  entry_out = palloc(sizeof(GISTENTRY));
797 
798  /*
799  ** Null key? Make a copy of the input entry and
800  ** return.
801  */
802  if ( DatumGetPointer(entry_in->key) == NULL )
803  {
804  POSTGIS_DEBUG(4, "[GIST] leafkey is null");
805  gistentryinit(*entry_out, (Datum) 0, entry_in->rel,
806  entry_in->page, entry_in->offset, FALSE);
807  POSTGIS_DEBUG(4, "[GIST] returning copy of input");
808  PG_RETURN_POINTER(entry_out);
809  }
810 
811  /* Extract our index key from the GiST entry. */
812  result = gserialized_datum_get_box2df_p(entry_in->key, &bbox_out);
813 
814  /* Is the bounding box valid (non-empty, non-infinite)? If not, return input uncompressed. */
815  if ( result == LW_FAILURE )
816  {
817  POSTGIS_DEBUG(4, "[GIST] empty geometry!");
818  PG_RETURN_POINTER(entry_in);
819  }
820 
821  POSTGIS_DEBUGF(4, "[GIST] got entry_in->key: %s", box2df_to_string(&bbox_out));
822 
823  /* Check all the dimensions for finite values */
824  if ( ! isfinite(bbox_out.xmax) || ! isfinite(bbox_out.xmin) ||
825  ! isfinite(bbox_out.ymax) || ! isfinite(bbox_out.ymin) )
826  {
827  POSTGIS_DEBUG(4, "[GIST] infinite geometry!");
828  PG_RETURN_POINTER(entry_in);
829  }
830 
831  /* Enure bounding box has minimums below maximums. */
832  box2df_validate(&bbox_out);
833 
834  /* Prepare GISTENTRY for return. */
835  gistentryinit(*entry_out, PointerGetDatum(box2df_copy(&bbox_out)),
836  entry_in->rel, entry_in->page, entry_in->offset, FALSE);
837 
838  /* Return GISTENTRY. */
839  POSTGIS_DEBUG(4, "[GIST] 'compress' function complete");
840  PG_RETURN_POINTER(entry_out);
841 }
842 
843 
844 /*
845 ** GiST support function.
846 ** Decompress an entry. Unused for geography, so we return.
847 */
849 Datum gserialized_gist_decompress_2d(PG_FUNCTION_ARGS)
850 {
851  POSTGIS_DEBUG(5, "[GIST] 'decompress' function called");
852  /* We don't decompress. Just return the input. */
853  PG_RETURN_POINTER(PG_GETARG_POINTER(0));
854 }
855 
856 
857 
858 /*
859 ** GiST support function. Called from gserialized_gist_consistent below.
860 */
861 static inline bool gserialized_gist_consistent_leaf_2d(BOX2DF *key, BOX2DF *query, StrategyNumber strategy)
862 {
863  bool retval;
864 
865  POSTGIS_DEBUGF(4, "[GIST] leaf consistent, strategy [%d], count[%i]",
866  strategy, g2d_counter_leaf++);
867 
868  switch (strategy)
869  {
870 
871  /* Basic overlaps */
872  case RTOverlapStrategyNumber:
873  retval = (bool) box2df_overlaps(key, query);
874  break;
875  case RTSameStrategyNumber:
876  retval = (bool) box2df_equals(key, query);
877  break;
878  case RTContainsStrategyNumber:
879  case RTOldContainsStrategyNumber:
880  retval = (bool) box2df_contains(key, query);
881  break;
882  case RTContainedByStrategyNumber:
883  case RTOldContainedByStrategyNumber:
884  retval = (bool) box2df_contains(query, key);
885  break;
886 
887  /* To one side */
888  case RTAboveStrategyNumber:
889  retval = (bool) box2df_above(key, query);
890  break;
891  case RTBelowStrategyNumber:
892  retval = (bool) box2df_below(key, query);
893  break;
894  case RTRightStrategyNumber:
895  retval = (bool) box2df_right(key, query);
896  break;
897  case RTLeftStrategyNumber:
898  retval = (bool) box2df_left(key, query);
899  break;
900 
901  /* Overlapping to one side */
902  case RTOverAboveStrategyNumber:
903  retval = (bool) box2df_overabove(key, query);
904  break;
905  case RTOverBelowStrategyNumber:
906  retval = (bool) box2df_overbelow(key, query);
907  break;
908  case RTOverRightStrategyNumber:
909  retval = (bool) box2df_overright(key, query);
910  break;
911  case RTOverLeftStrategyNumber:
912  retval = (bool) box2df_overleft(key, query);
913  break;
914 
915  default:
916  retval = FALSE;
917  }
918 
919  return (retval);
920 }
921 
922 /*
923 ** GiST support function. Called from gserialized_gist_consistent below.
924 */
925 static inline bool gserialized_gist_consistent_internal_2d(BOX2DF *key, BOX2DF *query, StrategyNumber strategy)
926 {
927  bool retval;
928 
929  POSTGIS_DEBUGF(4, "[GIST] internal consistent, strategy [%d], count[%i], query[%s], key[%s]",
930  strategy, g2d_counter_internal++, box2df_to_string(query), box2df_to_string(key) );
931 
932  switch (strategy)
933  {
934 
935  /* Basic overlaps */
936  case RTOverlapStrategyNumber:
937  retval = (bool) box2df_overlaps(key, query);
938  break;
939  case RTSameStrategyNumber:
940  case RTContainsStrategyNumber:
941  case RTOldContainsStrategyNumber:
942  retval = (bool) box2df_contains(key, query);
943  break;
944  case RTContainedByStrategyNumber:
945  case RTOldContainedByStrategyNumber:
946  retval = (bool) box2df_overlaps(key, query);
947  break;
948 
949  /* To one side */
950  case RTAboveStrategyNumber:
951  retval = (bool)(!box2df_overbelow(key, query));
952  break;
953  case RTBelowStrategyNumber:
954  retval = (bool)(!box2df_overabove(key, query));
955  break;
956  case RTRightStrategyNumber:
957  retval = (bool)(!box2df_overleft(key, query));
958  break;
959  case RTLeftStrategyNumber:
960  retval = (bool)(!box2df_overright(key, query));
961  break;
962 
963  /* Overlapping to one side */
964  case RTOverAboveStrategyNumber:
965  retval = (bool)(!box2df_below(key, query));
966  break;
967  case RTOverBelowStrategyNumber:
968  retval = (bool)(!box2df_above(key, query));
969  break;
970  case RTOverRightStrategyNumber:
971  retval = (bool)(!box2df_left(key, query));
972  break;
973  case RTOverLeftStrategyNumber:
974  retval = (bool)(!box2df_right(key, query));
975  break;
976 
977  default:
978  retval = FALSE;
979  }
980 
981  return (retval);
982 }
983 
984 /*
985 ** GiST support function. Take in a query and an entry and see what the
986 ** relationship is, based on the query strategy.
987 */
989 Datum gserialized_gist_consistent_2d(PG_FUNCTION_ARGS)
990 {
991  GISTENTRY *entry = (GISTENTRY*) PG_GETARG_POINTER(0);
992  StrategyNumber strategy = (StrategyNumber) PG_GETARG_UINT16(2);
993  bool result;
994  BOX2DF query_gbox_index;
995 
996  /* PostgreSQL 8.4 and later require the RECHECK flag to be set here,
997  rather than being supplied as part of the operator class definition */
998  bool *recheck = (bool *) PG_GETARG_POINTER(4);
999 
1000  /* We set recheck to false to avoid repeatedly pulling every "possibly matched" geometry
1001  out during index scans. For cases when the geometries are large, rechecking
1002  can make things twice as slow. */
1003  *recheck = false;
1004 
1005  POSTGIS_DEBUG(4, "[GIST] 'consistent' function called");
1006 
1007  /* Quick sanity check on query argument. */
1008  if ( DatumGetPointer(PG_GETARG_DATUM(1)) == NULL )
1009  {
1010  POSTGIS_DEBUG(4, "[GIST] null query pointer (!?!), returning false");
1011  PG_RETURN_BOOL(FALSE); /* NULL query! This is screwy! */
1012  }
1013 
1014  /* Quick sanity check on entry key. */
1015  if ( DatumGetPointer(entry->key) == NULL )
1016  {
1017  POSTGIS_DEBUG(4, "[GIST] null index entry, returning false");
1018  PG_RETURN_BOOL(FALSE); /* NULL entry! */
1019  }
1020 
1021  /* Null box should never make this far. */
1022  if ( gserialized_datum_get_box2df_p(PG_GETARG_DATUM(1), &query_gbox_index) == LW_FAILURE )
1023  {
1024  POSTGIS_DEBUG(4, "[GIST] null query_gbox_index!");
1025  PG_RETURN_BOOL(FALSE);
1026  }
1027 
1028  /* Treat leaf node tests different from internal nodes */
1029  if (GIST_LEAF(entry))
1030  {
1032  (BOX2DF*)DatumGetPointer(entry->key),
1033  &query_gbox_index, strategy);
1034  }
1035  else
1036  {
1038  (BOX2DF*)DatumGetPointer(entry->key),
1039  &query_gbox_index, strategy);
1040  }
1041 
1042  PG_RETURN_BOOL(result);
1043 }
1044 
1045 
1046 /*
1047 ** GiST support function. Take in a query and an entry and return the "distance"
1048 ** between them.
1049 **
1050 ** Given an index entry p and a query value q, this function determines the
1051 ** index entry's "distance" from the query value. This function must be
1052 ** supplied if the operator class contains any ordering operators. A query
1053 ** using the ordering operator will be implemented by returning index entries
1054 ** with the smallest "distance" values first, so the results must be consistent
1055 ** with the operator's semantics. For a leaf index entry the result just
1056 ** represents the distance to the index entry; for an internal tree node, the
1057 ** result must be the smallest distance that any child entry could have.
1058 **
1059 ** Strategy 13 = true distance tests <->
1060 ** Strategy 14 = box-based distance tests <#>
1061 */
1063 Datum gserialized_gist_distance_2d(PG_FUNCTION_ARGS)
1064 {
1065  GISTENTRY *entry = (GISTENTRY*) PG_GETARG_POINTER(0);
1066  BOX2DF query_box;
1067  BOX2DF *entry_box;
1068  StrategyNumber strategy = (StrategyNumber) PG_GETARG_UINT16(2);
1069  double distance;
1070 #if POSTGIS_PGSQL_VERSION >= 95
1071  bool *recheck = (bool *) PG_GETARG_POINTER(4);
1072 #endif
1073 
1074  POSTGIS_DEBUG(4, "[GIST] 'distance' function called");
1075 
1076  /* We are using '13' as the gist true-distance <-> strategy number
1077  * and '14' as the gist distance-between-boxes <#> strategy number */
1078  if ( strategy != 13 && strategy != 14 ) {
1079  elog(ERROR, "unrecognized strategy number: %d", strategy);
1080  PG_RETURN_FLOAT8(FLT_MAX);
1081  }
1082 
1083  /* Null box should never make this far. */
1084  if ( gserialized_datum_get_box2df_p(PG_GETARG_DATUM(1), &query_box) == LW_FAILURE )
1085  {
1086  POSTGIS_DEBUG(4, "[GIST] null query_gbox_index!");
1087  PG_RETURN_FLOAT8(FLT_MAX);
1088  }
1089 
1090  /* Get the entry box */
1091  entry_box = (BOX2DF*)DatumGetPointer(entry->key);
1092 
1093 #if POSTGIS_PGSQL_VERSION >= 95
1094 
1095  /* Box-style distance test */
1096  if ( strategy == 14 ) /* operator <#> */
1097  {
1098  distance = box2df_distance(entry_box, &query_box);
1099  }
1100  /* True distance test (formerly centroid distance) */
1101  else if ( strategy == 13 ) /* operator <-> */
1102  {
1103  /* In all cases, since we only have keys (boxes) we'll return */
1104  /* the minimum possible distance, which is the box2df_distance */
1105  /* and let the recheck sort things out in the case of leaves */
1106  distance = box2df_distance(entry_box, &query_box);
1107 
1108  if (GIST_LEAF(entry))
1109  *recheck = true;
1110  }
1111  else
1112  {
1113  elog(ERROR, "%s: reached unreachable code", __func__);
1114  PG_RETURN_NULL();
1115  }
1116 #else
1117  /* Box-style distance test */
1118  if ( strategy == 14 )
1119  {
1120  distance = (double)box2df_distance(entry_box, &query_box);
1121  PG_RETURN_FLOAT8(distance);
1122  }
1123 
1124  /* Treat leaf node tests different from internal nodes */
1125  if (GIST_LEAF(entry))
1126  {
1127  /* Calculate distance to leaves */
1128  distance = (double)box2df_distance_leaf_centroid(entry_box, &query_box);
1129  }
1130  else
1131  {
1132  /* Calculate distance for internal nodes */
1133  distance = (double)box2df_distance_node_centroid(entry_box, &query_box);
1134  }
1135 #endif
1136 
1137  PG_RETURN_FLOAT8(distance);
1138 }
1139 
1140 /*
1141 ** GiST support function. Calculate the "penalty" cost of adding this entry into an existing entry.
1142 ** Calculate the change in volume of the old entry once the new entry is added.
1143 ** TODO: Re-evaluate this in light of R*Tree penalty approaches.
1144 */
1146 Datum gserialized_gist_penalty_2d(PG_FUNCTION_ARGS)
1147 {
1148  GISTENTRY *origentry = (GISTENTRY*) PG_GETARG_POINTER(0);
1149  GISTENTRY *newentry = (GISTENTRY*) PG_GETARG_POINTER(1);
1150  float *result = (float*) PG_GETARG_POINTER(2);
1151  BOX2DF *gbox_index_orig, *gbox_index_new;
1152  float size_union, size_orig;
1153 
1154  POSTGIS_DEBUG(4, "[GIST] 'penalty' function called");
1155 
1156  gbox_index_orig = (BOX2DF*)DatumGetPointer(origentry->key);
1157  gbox_index_new = (BOX2DF*)DatumGetPointer(newentry->key);
1158 
1159  /* Drop out if we're dealing with null inputs. Shouldn't happen. */
1160  if ( (gbox_index_orig == NULL) && (gbox_index_new == NULL) )
1161  {
1162  POSTGIS_DEBUG(4, "[GIST] both inputs NULL! returning penalty of zero");
1163  *result = 0.0;
1164  PG_RETURN_FLOAT8(*result);
1165  }
1166 
1167  /* Calculate the size difference of the boxes. */
1168  size_union = box2df_union_size(gbox_index_orig, gbox_index_new);
1169  size_orig = box2df_size(gbox_index_orig);
1170  *result = size_union - size_orig;
1171 
1172  POSTGIS_DEBUGF(4, "[GIST] 'penalty', union size (%.12f), original size (%.12f), penalty (%.12f)", size_union, size_orig, *result);
1173 
1174  PG_RETURN_POINTER(result);
1175 }
1176 
1177 /*
1178 ** GiST support function. Merge all the boxes in a page into one master box.
1179 */
1181 Datum gserialized_gist_union_2d(PG_FUNCTION_ARGS)
1182 {
1183  GistEntryVector *entryvec = (GistEntryVector *) PG_GETARG_POINTER(0);
1184  int *sizep = (int *) PG_GETARG_POINTER(1); /* Size of the return value */
1185  int numranges, i;
1186  BOX2DF *box_cur, *box_union;
1187 
1188  POSTGIS_DEBUG(4, "[GIST] 'union' function called");
1189 
1190  numranges = entryvec->n;
1191 
1192  box_cur = (BOX2DF*) DatumGetPointer(entryvec->vector[0].key);
1193 
1194  box_union = box2df_copy(box_cur);
1195 
1196  for ( i = 1; i < numranges; i++ )
1197  {
1198  box_cur = (BOX2DF*) DatumGetPointer(entryvec->vector[i].key);
1199  box2df_merge(box_union, box_cur);
1200  }
1201 
1202  *sizep = sizeof(BOX2DF);
1203 
1204  POSTGIS_DEBUGF(4, "[GIST] 'union', numranges(%i), pageunion %s", numranges, box2df_to_string(box_union));
1205 
1206  PG_RETURN_POINTER(box_union);
1207 
1208 }
1209 
1210 /*
1211 ** GiST support function. Test equality of keys.
1212 */
1214 Datum gserialized_gist_same_2d(PG_FUNCTION_ARGS)
1215 {
1216  BOX2DF *b1 = (BOX2DF*)PG_GETARG_POINTER(0);
1217  BOX2DF *b2 = (BOX2DF*)PG_GETARG_POINTER(1);
1218  bool *result = (bool*)PG_GETARG_POINTER(2);
1219 
1220  POSTGIS_DEBUG(4, "[GIST] 'same' function called");
1221 
1222  *result = box2df_equals(b1, b2);
1223 
1224  PG_RETURN_POINTER(result);
1225 }
1226 
1227 #if KOROTKOV_SPLIT > 0
1228 /*
1229  * Adjust BOX2DF b boundaries with insertion of addon.
1230  */
1231 static void
1232 adjustBox(BOX2DF *b, BOX2DF *addon)
1233 {
1234  if (b->xmax < addon->xmax || isnan(b->xmax))
1235  b->xmax = addon->xmax;
1236  if (b->xmin > addon->xmin || isnan(b->xmin))
1237  b->xmin = addon->xmin;
1238  if (b->ymax < addon->ymax || isnan(b->ymax))
1239  b->ymax = addon->ymax;
1240  if (b->ymin > addon->ymin || isnan(b->ymin))
1241  b->ymin = addon->ymin;
1242 }
1243 
1244 /*
1245  * Trivial split: half of entries will be placed on one page
1246  * and another half - to another
1247  */
1248 static void
1249 fallbackSplit(GistEntryVector *entryvec, GIST_SPLITVEC *v)
1250 {
1251  OffsetNumber i,
1252  maxoff;
1253  BOX2DF *unionL = NULL,
1254  *unionR = NULL;
1255  int nbytes;
1256 
1257  maxoff = entryvec->n - 1;
1258 
1259  nbytes = (maxoff + 2) * sizeof(OffsetNumber);
1260  v->spl_left = (OffsetNumber *) palloc(nbytes);
1261  v->spl_right = (OffsetNumber *) palloc(nbytes);
1262  v->spl_nleft = v->spl_nright = 0;
1263 
1264  for (i = FirstOffsetNumber; i <= maxoff; i = OffsetNumberNext(i))
1265  {
1266  BOX2DF *cur = (BOX2DF *) DatumGetPointer(entryvec->vector[i].key);
1267 
1268  if (i <= (maxoff - FirstOffsetNumber + 1) / 2)
1269  {
1270  v->spl_left[v->spl_nleft] = i;
1271  if (unionL == NULL)
1272  {
1273  unionL = (BOX2DF *) palloc(sizeof(BOX2DF));
1274  *unionL = *cur;
1275  }
1276  else
1277  adjustBox(unionL, cur);
1278 
1279  v->spl_nleft++;
1280  }
1281  else
1282  {
1283  v->spl_right[v->spl_nright] = i;
1284  if (unionR == NULL)
1285  {
1286  unionR = (BOX2DF *) palloc(sizeof(BOX2DF));
1287  *unionR = *cur;
1288  }
1289  else
1290  adjustBox(unionR, cur);
1291 
1292  v->spl_nright++;
1293  }
1294  }
1295 
1296  if (v->spl_ldatum_exists)
1297  adjustBox(unionL, (BOX2DF *) DatumGetPointer(v->spl_ldatum));
1298  v->spl_ldatum = BoxPGetDatum(unionL);
1299 
1300  if (v->spl_rdatum_exists)
1301  adjustBox(unionR, (BOX2DF *) DatumGetPointer(v->spl_rdatum));
1302  v->spl_rdatum = BoxPGetDatum(unionR);
1303 
1304  v->spl_ldatum_exists = v->spl_rdatum_exists = false;
1305 }
1306 
1307 /*
1308  * Represents information about an entry that can be placed to either group
1309  * without affecting overlap over selected axis ("common entry").
1310  */
1311 typedef struct
1312 {
1313  /* Index of entry in the initial array */
1314  int index;
1315  /* Delta between penalties of entry insertion into different groups */
1316  float delta;
1317 } CommonEntry;
1318 
1319 /*
1320  * Context for g_box_consider_split. Contains information about currently
1321  * selected split and some general information.
1322  */
1323 typedef struct
1324 {
1325  int entriesCount; /* total number of entries being split */
1326  BOX2DF boundingBox; /* minimum bounding box across all entries */
1327 
1328  /* Information about currently selected split follows */
1329 
1330  bool first; /* true if no split was selected yet */
1331 
1332  float leftUpper; /* upper bound of left interval */
1333  float rightLower; /* lower bound of right interval */
1334 
1335  float4 ratio;
1336  float4 overlap;
1337  int dim; /* axis of this split */
1338  float range; /* width of general MBR projection to the
1339  * selected axis */
1341 
1342 /*
1343  * Interval represents projection of box to axis.
1344  */
1345 typedef struct
1346 {
1347  float lower,
1348  upper;
1349 } SplitInterval;
1350 
1351 /*
1352  * Interval comparison function by lower bound of the interval;
1353  */
1354 static int
1355 interval_cmp_lower(const void *i1, const void *i2)
1356 {
1357  float lower1 = ((const SplitInterval *) i1)->lower,
1358  lower2 = ((const SplitInterval *) i2)->lower;
1359 
1360  if (isnan(lower1))
1361  {
1362  if (isnan(lower2))
1363  return 0;
1364  else
1365  return 1;
1366  }
1367  else if (isnan(lower2))
1368  {
1369  return -1;
1370  }
1371  else
1372  {
1373  if (lower1 < lower2)
1374  return -1;
1375  else if (lower1 > lower2)
1376  return 1;
1377  else
1378  return 0;
1379  }
1380 }
1381 
1382 
1383 /*
1384  * Interval comparison function by upper bound of the interval;
1385  */
1386 static int
1387 interval_cmp_upper(const void *i1, const void *i2)
1388 {
1389  float upper1 = ((const SplitInterval *) i1)->upper,
1390  upper2 = ((const SplitInterval *) i2)->upper;
1391 
1392  if (isnan(upper1))
1393  {
1394  if (isnan(upper2))
1395  return 0;
1396  else
1397  return -1;
1398  }
1399  else if (isnan(upper2))
1400  {
1401  return 1;
1402  }
1403  else
1404  {
1405  if (upper1 < upper2)
1406  return -1;
1407  else if (upper1 > upper2)
1408  return 1;
1409  else
1410  return 0;
1411  }
1412 }
1413 
1414 /*
1415  * Replace negative value with zero.
1416  */
1417 static inline float
1418 non_negative(float val)
1419 {
1420  if (val >= 0.0f)
1421  return val;
1422  else
1423  return 0.0f;
1424 }
1425 
1426 /*
1427  * Consider replacement of currently selected split with the better one.
1428  */
1429 static inline void
1431  float rightLower, int minLeftCount,
1432  float leftUpper, int maxLeftCount)
1433 {
1434  int leftCount,
1435  rightCount;
1436  float4 ratio,
1437  overlap;
1438  float range;
1439 
1440  POSTGIS_DEBUGF(5, "consider split: dimNum = %d, rightLower = %f, "
1441  "minLeftCount = %d, leftUpper = %f, maxLeftCount = %d ",
1442  dimNum, rightLower, minLeftCount, leftUpper, maxLeftCount);
1443 
1444  /*
1445  * Calculate entries distribution ratio assuming most uniform distribution
1446  * of common entries.
1447  */
1448  if (minLeftCount >= (context->entriesCount + 1) / 2)
1449  {
1450  leftCount = minLeftCount;
1451  }
1452  else
1453  {
1454  if (maxLeftCount <= context->entriesCount / 2)
1455  leftCount = maxLeftCount;
1456  else
1457  leftCount = context->entriesCount / 2;
1458  }
1459  rightCount = context->entriesCount - leftCount;
1460 
1461  /*
1462  * Ratio of split - quotient between size of lesser group and total
1463  * entries count.
1464  */
1465  ratio = ((float4) Min(leftCount, rightCount)) /
1466  ((float4) context->entriesCount);
1467 
1468  if (ratio > LIMIT_RATIO)
1469  {
1470  bool selectthis = false;
1471 
1472  /*
1473  * The ratio is acceptable, so compare current split with previously
1474  * selected one. Between splits of one dimension we search for minimal
1475  * overlap (allowing negative values) and minimal ration (between same
1476  * overlaps. We switch dimension if find less overlap (non-negative)
1477  * or less range with same overlap.
1478  */
1479  if (dimNum == 0)
1480  range = context->boundingBox.xmax - context->boundingBox.xmin;
1481  else
1482  range = context->boundingBox.ymax - context->boundingBox.ymin;
1483 
1484  overlap = (leftUpper - rightLower) / range;
1485 
1486  /* If there is no previous selection, select this */
1487  if (context->first)
1488  selectthis = true;
1489  else if (context->dim == dimNum)
1490  {
1491  /*
1492  * Within the same dimension, choose the new split if it has a
1493  * smaller overlap, or same overlap but better ratio.
1494  */
1495  if (overlap < context->overlap ||
1496  (overlap == context->overlap && ratio > context->ratio))
1497  selectthis = true;
1498  }
1499  else
1500  {
1501  /*
1502  * Across dimensions, choose the new split if it has a smaller
1503  * *non-negative* overlap, or same *non-negative* overlap but
1504  * bigger range. This condition differs from the one described in
1505  * the article. On the datasets where leaf MBRs don't overlap
1506  * themselves, non-overlapping splits (i.e. splits which have zero
1507  * *non-negative* overlap) are frequently possible. In this case
1508  * splits tends to be along one dimension, because most distant
1509  * non-overlapping splits (i.e. having lowest negative overlap)
1510  * appears to be in the same dimension as in the previous split.
1511  * Therefore MBRs appear to be very prolonged along another
1512  * dimension, which leads to bad search performance. Using range
1513  * as the second split criteria makes MBRs more quadratic. Using
1514  * *non-negative* overlap instead of overlap as the first split
1515  * criteria gives to range criteria a chance to matter, because
1516  * non-overlapping splits are equivalent in this criteria.
1517  */
1518  if (non_negative(overlap) < non_negative(context->overlap) ||
1519  (range > context->range &&
1520  non_negative(overlap) <= non_negative(context->overlap)))
1521  selectthis = true;
1522  }
1523 
1524  if (selectthis)
1525  {
1526  /* save information about selected split */
1527  context->first = false;
1528  context->ratio = ratio;
1529  context->range = range;
1530  context->overlap = overlap;
1531  context->rightLower = rightLower;
1532  context->leftUpper = leftUpper;
1533  context->dim = dimNum;
1534  POSTGIS_DEBUG(5, "split selected");
1535  }
1536  }
1537 }
1538 
1539 /*
1540  * Return increase of original BOX2DF area by new BOX2DF area insertion.
1541  */
1542 static float
1543 box_penalty(BOX2DF *original, BOX2DF *new)
1544 {
1545  float union_width,
1546  union_height;
1547 
1548  union_width = Max(original->xmax, new->xmax) -
1549  Min(original->xmin, new->xmin);
1550  union_height = Max(original->ymax, new->ymax) -
1551  Min(original->ymin, new->ymin);
1552  return union_width * union_height - (original->xmax - original->xmin) *
1553  (original->ymax - original->ymin);
1554 }
1555 
1556 /*
1557  * Compare common entries by their deltas.
1558  */
1559 static int
1560 common_entry_cmp(const void *i1, const void *i2)
1561 {
1562  float delta1 = ((const CommonEntry *) i1)->delta,
1563  delta2 = ((const CommonEntry *) i2)->delta;
1564 
1565  if (delta1 < delta2)
1566  return -1;
1567  else if (delta1 > delta2)
1568  return 1;
1569  else
1570  return 0;
1571 }
1572 
1573 /*
1574  * --------------------------------------------------------------------------
1575  * Double sorting split algorithm. This is used for both boxes and points.
1576  *
1577  * The algorithm finds split of boxes by considering splits along each axis.
1578  * Each entry is first projected as an interval on the X-axis, and different
1579  * ways to split the intervals into two groups are considered, trying to
1580  * minimize the overlap of the groups. Then the same is repeated for the
1581  * Y-axis, and the overall best split is chosen. The quality of a split is
1582  * determined by overlap along that axis and some other criteria (see
1583  * g_box_consider_split).
1584  *
1585  * After that, all the entries are divided into three groups:
1586  *
1587  * 1) Entries which should be placed to the left group
1588  * 2) Entries which should be placed to the right group
1589  * 3) "Common entries" which can be placed to any of groups without affecting
1590  * of overlap along selected axis.
1591  *
1592  * The common entries are distributed by minimizing penalty.
1593  *
1594  * For details see:
1595  * "A new double sorting-based node splitting algorithm for R-tree", A. Korotkov
1596  * http://syrcose.ispras.ru/2011/files/SYRCoSE2011_Proceedings.pdf#page=36
1597  * --------------------------------------------------------------------------
1598  */
1600 Datum gserialized_gist_picksplit_2d(PG_FUNCTION_ARGS)
1601 {
1602  GistEntryVector *entryvec = (GistEntryVector *) PG_GETARG_POINTER(0);
1603  GIST_SPLITVEC *v = (GIST_SPLITVEC *) PG_GETARG_POINTER(1);
1604  OffsetNumber i,
1605  maxoff;
1606  ConsiderSplitContext context;
1607  BOX2DF *box,
1608  *leftBox,
1609  *rightBox;
1610  int dim,
1611  commonEntriesCount;
1612  SplitInterval *intervalsLower,
1613  *intervalsUpper;
1614  CommonEntry *commonEntries;
1615  int nentries;
1616 
1617  POSTGIS_DEBUG(3, "[GIST] 'picksplit' entered");
1618 
1619  memset(&context, 0, sizeof(ConsiderSplitContext));
1620 
1621  maxoff = entryvec->n - 1;
1622  nentries = context.entriesCount = maxoff - FirstOffsetNumber + 1;
1623 
1624  /* Allocate arrays for intervals along axes */
1625  intervalsLower = (SplitInterval *) palloc(nentries * sizeof(SplitInterval));
1626  intervalsUpper = (SplitInterval *) palloc(nentries * sizeof(SplitInterval));
1627 
1628  /*
1629  * Calculate the overall minimum bounding box over all the entries.
1630  */
1631  for (i = FirstOffsetNumber; i <= maxoff; i = OffsetNumberNext(i))
1632  {
1633  box = (BOX2DF *) DatumGetPointer(entryvec->vector[i].key);
1634  if (i == FirstOffsetNumber)
1635  context.boundingBox = *box;
1636  else
1637  adjustBox(&context.boundingBox, box);
1638  }
1639 
1640  POSTGIS_DEBUGF(4, "boundingBox is %s", box2df_to_string(
1641  &context.boundingBox));
1642 
1643  /*
1644  * Iterate over axes for optimal split searching.
1645  */
1646  context.first = true; /* nothing selected yet */
1647  for (dim = 0; dim < 2; dim++)
1648  {
1649  float leftUpper,
1650  rightLower;
1651  int i1,
1652  i2;
1653 
1654  /* Project each entry as an interval on the selected axis. */
1655  for (i = FirstOffsetNumber; i <= maxoff; i = OffsetNumberNext(i))
1656  {
1657  box = (BOX2DF *) DatumGetPointer(entryvec->vector[i].key);
1658  if (dim == 0)
1659  {
1660  intervalsLower[i - FirstOffsetNumber].lower = box->xmin;
1661  intervalsLower[i - FirstOffsetNumber].upper = box->xmax;
1662  }
1663  else
1664  {
1665  intervalsLower[i - FirstOffsetNumber].lower = box->ymin;
1666  intervalsLower[i - FirstOffsetNumber].upper = box->ymax;
1667  }
1668  }
1669 
1670  /*
1671  * Make two arrays of intervals: one sorted by lower bound and another
1672  * sorted by upper bound.
1673  */
1674  memcpy(intervalsUpper, intervalsLower,
1675  sizeof(SplitInterval) * nentries);
1676  qsort(intervalsLower, nentries, sizeof(SplitInterval),
1678  qsort(intervalsUpper, nentries, sizeof(SplitInterval),
1680 
1681  /*----
1682  * The goal is to form a left and right interval, so that every entry
1683  * interval is contained by either left or right interval (or both).
1684  *
1685  * For example, with the intervals (0,1), (1,3), (2,3), (2,4):
1686  *
1687  * 0 1 2 3 4
1688  * +-+
1689  * +---+
1690  * +-+
1691  * +---+
1692  *
1693  * The left and right intervals are of the form (0,a) and (b,4).
1694  * We first consider splits where b is the lower bound of an entry.
1695  * We iterate through all entries, and for each b, calculate the
1696  * smallest possible a. Then we consider splits where a is the
1697  * uppper bound of an entry, and for each a, calculate the greatest
1698  * possible b.
1699  *
1700  * In the above example, the first loop would consider splits:
1701  * b=0: (0,1)-(0,4)
1702  * b=1: (0,1)-(1,4)
1703  * b=2: (0,3)-(2,4)
1704  *
1705  * And the second loop:
1706  * a=1: (0,1)-(1,4)
1707  * a=3: (0,3)-(2,4)
1708  * a=4: (0,4)-(2,4)
1709  */
1710 
1711  /*
1712  * Iterate over lower bound of right group, finding smallest possible
1713  * upper bound of left group.
1714  */
1715  i1 = 0;
1716  i2 = 0;
1717  rightLower = intervalsLower[i1].lower;
1718  leftUpper = intervalsUpper[i2].lower;
1719  while (true)
1720  {
1721  /*
1722  * Find next lower bound of right group.
1723  */
1724  while (i1 < nentries && (rightLower == intervalsLower[i1].lower ||
1725  isnan(intervalsLower[i1].lower)))
1726  {
1727  leftUpper = Max(leftUpper, intervalsLower[i1].upper);
1728  i1++;
1729  }
1730  if (i1 >= nentries)
1731  break;
1732  rightLower = intervalsLower[i1].lower;
1733 
1734  /*
1735  * Find count of intervals which anyway should be placed to the
1736  * left group.
1737  */
1738  while (i2 < nentries && intervalsUpper[i2].upper <= leftUpper)
1739  i2++;
1740 
1741  /*
1742  * Consider found split.
1743  */
1744  g_box_consider_split(&context, dim, rightLower, i1, leftUpper, i2);
1745  }
1746 
1747  /*
1748  * Iterate over upper bound of left group finding greates possible
1749  * lower bound of right group.
1750  */
1751  i1 = nentries - 1;
1752  i2 = nentries - 1;
1753  rightLower = intervalsLower[i1].upper;
1754  leftUpper = intervalsUpper[i2].upper;
1755  while (true)
1756  {
1757  /*
1758  * Find next upper bound of left group.
1759  */
1760  while (i2 >= 0 && (leftUpper == intervalsUpper[i2].upper ||
1761  isnan(intervalsUpper[i2].upper)))
1762  {
1763  rightLower = Min(rightLower, intervalsUpper[i2].lower);
1764  i2--;
1765  }
1766  if (i2 < 0)
1767  break;
1768  leftUpper = intervalsUpper[i2].upper;
1769 
1770  /*
1771  * Find count of intervals which anyway should be placed to the
1772  * right group.
1773  */
1774  while (i1 >= 0 && intervalsLower[i1].lower >= rightLower)
1775  i1--;
1776 
1777  /*
1778  * Consider found split.
1779  */
1780  g_box_consider_split(&context, dim,
1781  rightLower, i1 + 1, leftUpper, i2 + 1);
1782  }
1783  }
1784 
1785  /*
1786  * If we failed to find any acceptable splits, use trivial split.
1787  */
1788  if (context.first)
1789  {
1790  POSTGIS_DEBUG(4, "no acceptable splits, trivial split");
1791  fallbackSplit(entryvec, v);
1792  PG_RETURN_POINTER(v);
1793  }
1794 
1795  /*
1796  * Ok, we have now selected the split across one axis.
1797  *
1798  * While considering the splits, we already determined that there will be
1799  * enough entries in both groups to reach the desired ratio, but we did
1800  * not memorize which entries go to which group. So determine that now.
1801  */
1802 
1803  POSTGIS_DEBUGF(4, "split direction: %d", context.dim);
1804 
1805  /* Allocate vectors for results */
1806  v->spl_left = (OffsetNumber *) palloc(nentries * sizeof(OffsetNumber));
1807  v->spl_right = (OffsetNumber *) palloc(nentries * sizeof(OffsetNumber));
1808  v->spl_nleft = 0;
1809  v->spl_nright = 0;
1810 
1811  /* Allocate bounding boxes of left and right groups */
1812  leftBox = palloc0(sizeof(BOX2DF));
1813  rightBox = palloc0(sizeof(BOX2DF));
1814 
1815  /*
1816  * Allocate an array for "common entries" - entries which can be placed to
1817  * either group without affecting overlap along selected axis.
1818  */
1819  commonEntriesCount = 0;
1820  commonEntries = (CommonEntry *) palloc(nentries * sizeof(CommonEntry));
1821 
1822  /* Helper macros to place an entry in the left or right group */
1823 #define PLACE_LEFT(box, off) \
1824  do { \
1825  if (v->spl_nleft > 0) \
1826  adjustBox(leftBox, box); \
1827  else \
1828  *leftBox = *(box); \
1829  v->spl_left[v->spl_nleft++] = off; \
1830  } while(0)
1831 
1832 #define PLACE_RIGHT(box, off) \
1833  do { \
1834  if (v->spl_nright > 0) \
1835  adjustBox(rightBox, box); \
1836  else \
1837  *rightBox = *(box); \
1838  v->spl_right[v->spl_nright++] = off; \
1839  } while(0)
1840 
1841  /*
1842  * Distribute entries which can be distributed unambiguously, and collect
1843  * common entries.
1844  */
1845  for (i = FirstOffsetNumber; i <= maxoff; i = OffsetNumberNext(i))
1846  {
1847  float lower,
1848  upper;
1849 
1850  /*
1851  * Get upper and lower bounds along selected axis.
1852  */
1853  box = (BOX2DF *) DatumGetPointer(entryvec->vector[i].key);
1854  if (context.dim == 0)
1855  {
1856  lower = box->xmin;
1857  upper = box->xmax;
1858  }
1859  else
1860  {
1861  lower = box->ymin;
1862  upper = box->ymax;
1863  }
1864 
1865  if (upper <= context.leftUpper || isnan(upper))
1866  {
1867  /* Fits to the left group */
1868  if (lower >= context.rightLower || isnan(lower))
1869  {
1870  /* Fits also to the right group, so "common entry" */
1871  commonEntries[commonEntriesCount++].index = i;
1872  }
1873  else
1874  {
1875  /* Doesn't fit to the right group, so join to the left group */
1876  PLACE_LEFT(box, i);
1877  }
1878  }
1879  else
1880  {
1881  /*
1882  * Each entry should fit on either left or right group. Since this
1883  * entry didn't fit on the left group, it better fit in the right
1884  * group.
1885  */
1886  Assert(lower >= context.rightLower);
1887 
1888  /* Doesn't fit to the left group, so join to the right group */
1889  PLACE_RIGHT(box, i);
1890  }
1891  }
1892 
1893  POSTGIS_DEBUGF(4, "leftBox is %s", box2df_to_string(leftBox));
1894  POSTGIS_DEBUGF(4, "rightBox is %s", box2df_to_string(rightBox));
1895 
1896  /*
1897  * Distribute "common entries", if any.
1898  */
1899  if (commonEntriesCount > 0)
1900  {
1901  /*
1902  * Calculate minimum number of entries that must be placed in both
1903  * groups, to reach LIMIT_RATIO.
1904  */
1905  int m = ceil(LIMIT_RATIO * (double) nentries);
1906 
1907  /*
1908  * Calculate delta between penalties of join "common entries" to
1909  * different groups.
1910  */
1911  for (i = 0; i < commonEntriesCount; i++)
1912  {
1913  box = (BOX2DF *) DatumGetPointer(entryvec->vector[
1914  commonEntries[i].index].key);
1915  commonEntries[i].delta = Abs(box_penalty(leftBox, box) -
1916  box_penalty(rightBox, box));
1917  }
1918 
1919  /*
1920  * Sort "common entries" by calculated deltas in order to distribute
1921  * the most ambiguous entries first.
1922  */
1923  qsort(commonEntries, commonEntriesCount, sizeof(CommonEntry), common_entry_cmp);
1924 
1925  /*
1926  * Distribute "common entries" between groups.
1927  */
1928  for (i = 0; i < commonEntriesCount; i++)
1929  {
1930  box = (BOX2DF *) DatumGetPointer(entryvec->vector[
1931  commonEntries[i].index].key);
1932 
1933  /*
1934  * Check if we have to place this entry in either group to achieve
1935  * LIMIT_RATIO.
1936  */
1937  if (v->spl_nleft + (commonEntriesCount - i) <= m)
1938  PLACE_LEFT(box, commonEntries[i].index);
1939  else if (v->spl_nright + (commonEntriesCount - i) <= m)
1940  PLACE_RIGHT(box, commonEntries[i].index);
1941  else
1942  {
1943  /* Otherwise select the group by minimal penalty */
1944  if (box_penalty(leftBox, box) < box_penalty(rightBox, box))
1945  PLACE_LEFT(box, commonEntries[i].index);
1946  else
1947  PLACE_RIGHT(box, commonEntries[i].index);
1948  }
1949  }
1950  }
1951  v->spl_ldatum = PointerGetDatum(leftBox);
1952  v->spl_rdatum = PointerGetDatum(rightBox);
1953 
1954  POSTGIS_DEBUG(4, "[GIST] 'picksplit' completed");
1955 
1956  PG_RETURN_POINTER(v);
1957 }
1958 
1959 #else /* !KOROTOV_SPLIT */
1960 
1961 typedef struct
1962 {
1963  BOX2DF *key;
1964  int pos;
1965 }
1966 KBsort;
1967 
1968 static int
1969 compare_KB(const void* a, const void* b)
1970 {
1971  BOX2DF *abox = ((KBsort*)a)->key;
1972  BOX2DF *bbox = ((KBsort*)b)->key;
1973  float sa = (abox->xmax - abox->xmin) * (abox->ymax - abox->ymin);
1974  float sb = (bbox->xmax - bbox->xmin) * (bbox->ymax - bbox->ymin);
1975 
1976  if ( sa==sb ) return 0;
1977  return ( sa>sb ) ? 1 : -1;
1978 }
1979 
1986 Datum gserialized_gist_picksplit_2d(PG_FUNCTION_ARGS)
1987 {
1988  GistEntryVector *entryvec = (GistEntryVector *) PG_GETARG_POINTER(0);
1989 
1990  GIST_SPLITVEC *v = (GIST_SPLITVEC *) PG_GETARG_POINTER(1);
1991  OffsetNumber i;
1992  OffsetNumber *listL, *listR, *listB, *listT;
1993  BOX2DF *unionL, *unionR, *unionB, *unionT;
1994  int posL, posR, posB, posT;
1995  BOX2DF pageunion;
1996  BOX2DF *cur;
1997  char direction = ' ';
1998  bool allisequal = true;
1999  OffsetNumber maxoff;
2000  int nbytes;
2001 
2002  POSTGIS_DEBUG(3, "[GIST] 'picksplit' entered");
2003 
2004  posL = posR = posB = posT = 0;
2005 
2006  maxoff = entryvec->n - 1;
2007  cur = (BOX2DF*) DatumGetPointer(entryvec->vector[FirstOffsetNumber].key);
2008 
2009  memcpy((void *) &pageunion, (void *) cur, sizeof(BOX2DF));
2010 
2011  /* find MBR */
2012  for (i = OffsetNumberNext(FirstOffsetNumber); i <= maxoff; i = OffsetNumberNext(i))
2013  {
2014  cur = (BOX2DF *) DatumGetPointer(entryvec->vector[i].key);
2015 
2016  if ( allisequal == true && (
2017  pageunion.xmax != cur->xmax ||
2018  pageunion.ymax != cur->ymax ||
2019  pageunion.xmin != cur->xmin ||
2020  pageunion.ymin != cur->ymin
2021  ) )
2022  allisequal = false;
2023 
2024  if (pageunion.xmax < cur->xmax)
2025  pageunion.xmax = cur->xmax;
2026  if (pageunion.xmin > cur->xmin)
2027  pageunion.xmin = cur->xmin;
2028  if (pageunion.ymax < cur->ymax)
2029  pageunion.ymax = cur->ymax;
2030  if (pageunion.ymin > cur->ymin)
2031  pageunion.ymin = cur->ymin;
2032  }
2033 
2034  POSTGIS_DEBUGF(4, "pageunion is %s", box2df_to_string(&pageunion));
2035 
2036  nbytes = (maxoff + 2) * sizeof(OffsetNumber);
2037  listL = (OffsetNumber *) palloc(nbytes);
2038  listR = (OffsetNumber *) palloc(nbytes);
2039  unionL = (BOX2DF *) palloc(sizeof(BOX2DF));
2040  unionR = (BOX2DF *) palloc(sizeof(BOX2DF));
2041 
2042  if (allisequal)
2043  {
2044  POSTGIS_DEBUG(4, " AllIsEqual!");
2045 
2046  cur = (BOX2DF*) DatumGetPointer(entryvec->vector[OffsetNumberNext(FirstOffsetNumber)].key);
2047 
2048 
2049  if (memcmp((void *) cur, (void *) &pageunion, sizeof(BOX2DF)) == 0)
2050  {
2051  v->spl_left = listL;
2052  v->spl_right = listR;
2053  v->spl_nleft = v->spl_nright = 0;
2054  memcpy((void *) unionL, (void *) &pageunion, sizeof(BOX2DF));
2055  memcpy((void *) unionR, (void *) &pageunion, sizeof(BOX2DF));
2056 
2057  for (i = FirstOffsetNumber; i <= maxoff; i = OffsetNumberNext(i))
2058  {
2059  if (i <= (maxoff - FirstOffsetNumber + 1) / 2)
2060  {
2061  v->spl_left[v->spl_nleft] = i;
2062  v->spl_nleft++;
2063  }
2064  else
2065  {
2066  v->spl_right[v->spl_nright] = i;
2067  v->spl_nright++;
2068  }
2069  }
2070  v->spl_ldatum = PointerGetDatum(unionL);
2071  v->spl_rdatum = PointerGetDatum(unionR);
2072 
2073  PG_RETURN_POINTER(v);
2074  }
2075  }
2076 
2077  listB = (OffsetNumber *) palloc(nbytes);
2078  listT = (OffsetNumber *) palloc(nbytes);
2079  unionB = (BOX2DF *) palloc(sizeof(BOX2DF));
2080  unionT = (BOX2DF *) palloc(sizeof(BOX2DF));
2081 
2082 #define ADDLIST( list, unionD, pos, num ) do { \
2083  if ( pos ) { \
2084  if ( unionD->xmax < cur->xmax ) unionD->xmax = cur->xmax; \
2085  if ( unionD->xmin > cur->xmin ) unionD->xmin = cur->xmin; \
2086  if ( unionD->ymax < cur->ymax ) unionD->ymax = cur->ymax; \
2087  if ( unionD->ymin > cur->ymin ) unionD->ymin = cur->ymin; \
2088  } else { \
2089  memcpy( (void*)unionD, (void*) cur, sizeof( BOX2DF ) ); \
2090  } \
2091  list[pos] = num; \
2092  (pos)++; \
2093 } while(0)
2094 
2095  for (i = FirstOffsetNumber; i <= maxoff; i = OffsetNumberNext(i))
2096  {
2097  cur = (BOX2DF*) DatumGetPointer(entryvec->vector[i].key);
2098 
2099  if (cur->xmin - pageunion.xmin < pageunion.xmax - cur->xmax)
2100  ADDLIST(listL, unionL, posL,i);
2101  else
2102  ADDLIST(listR, unionR, posR,i);
2103  if (cur->ymin - pageunion.ymin < pageunion.ymax - cur->ymax)
2104  ADDLIST(listB, unionB, posB,i);
2105  else
2106  ADDLIST(listT, unionT, posT,i);
2107  }
2108 
2109  POSTGIS_DEBUGF(4, "unionL is %s", box2df_to_string(unionL));
2110  POSTGIS_DEBUGF(4, "unionR is %s", box2df_to_string(unionR));
2111  POSTGIS_DEBUGF(4, "unionT is %s", box2df_to_string(unionT));
2112  POSTGIS_DEBUGF(4, "unionB is %s", box2df_to_string(unionB));
2113 
2114  /* bad disposition, sort by ascending and resplit */
2115  if ( (posR==0 || posL==0) && (posT==0 || posB==0) )
2116  {
2117  KBsort *arr = (KBsort*)palloc( sizeof(KBsort) * maxoff );
2118  posL = posR = posB = posT = 0;
2119  for (i = FirstOffsetNumber; i <= maxoff; i = OffsetNumberNext(i))
2120  {
2121  arr[i-1].key = (BOX2DF*) DatumGetPointer(entryvec->vector[i].key);
2122  arr[i-1].pos = i;
2123  }
2124  qsort( arr, maxoff, sizeof(KBsort), compare_KB );
2125  for (i = FirstOffsetNumber; i <= maxoff; i = OffsetNumberNext(i))
2126  {
2127  cur = arr[i-1].key;
2128  if (cur->xmin - pageunion.xmin < pageunion.xmax - cur->xmax)
2129  ADDLIST(listL, unionL, posL,arr[i-1].pos);
2130  else if ( cur->xmin - pageunion.xmin == pageunion.xmax - cur->xmax )
2131  {
2132  if ( posL>posR )
2133  ADDLIST(listR, unionR, posR,arr[i-1].pos);
2134  else
2135  ADDLIST(listL, unionL, posL,arr[i-1].pos);
2136  }
2137  else
2138  ADDLIST(listR, unionR, posR,arr[i-1].pos);
2139 
2140  if (cur->ymin - pageunion.ymin < pageunion.ymax - cur->ymax)
2141  ADDLIST(listB, unionB, posB,arr[i-1].pos);
2142  else if ( cur->ymin - pageunion.ymin == pageunion.ymax - cur->ymax )
2143  {
2144  if ( posB>posT )
2145  ADDLIST(listT, unionT, posT,arr[i-1].pos);
2146  else
2147  ADDLIST(listB, unionB, posB,arr[i-1].pos);
2148  }
2149  else
2150  ADDLIST(listT, unionT, posT,arr[i-1].pos);
2151  }
2152  pfree(arr);
2153  }
2154 
2155  /* which split more optimal? */
2156  if (Max(posL, posR) < Max(posB, posT))
2157  direction = 'x';
2158  else if (Max(posL, posR) > Max(posB, posT))
2159  direction = 'y';
2160  else
2161  {
2162  float sizeLR, sizeBT;
2163  BOX2DF interLR, interBT;
2164 
2165  if ( box2df_intersection(unionL, unionR, &interLR) == FALSE )
2166  sizeLR = 0.0;
2167  else
2168  sizeLR = box2df_size(&interLR);
2169 
2170  if ( box2df_intersection(unionB, unionT, &interBT) == FALSE )
2171  sizeBT = 0.0;
2172  else
2173  sizeBT = box2df_size(&interBT);
2174 
2175  if (sizeLR < sizeBT)
2176  direction = 'x';
2177  else
2178  direction = 'y';
2179  }
2180 
2181  POSTGIS_DEBUGF(4, "split direction '%c'", direction);
2182 
2183  if (direction == 'x')
2184  {
2185  pfree(unionB);
2186  pfree(listB);
2187  pfree(unionT);
2188  pfree(listT);
2189 
2190  v->spl_left = listL;
2191  v->spl_right = listR;
2192  v->spl_nleft = posL;
2193  v->spl_nright = posR;
2194  v->spl_ldatum = PointerGetDatum(unionL);
2195  v->spl_rdatum = PointerGetDatum(unionR);
2196  }
2197  else
2198  {
2199  pfree(unionR);
2200  pfree(listR);
2201  pfree(unionL);
2202  pfree(listL);
2203 
2204  v->spl_left = listB;
2205  v->spl_right = listT;
2206  v->spl_nleft = posB;
2207  v->spl_nright = posT;
2208  v->spl_ldatum = PointerGetDatum(unionB);
2209  v->spl_rdatum = PointerGetDatum(unionT);
2210  }
2211 
2212  POSTGIS_DEBUG(4, "[GIST] 'picksplit' completed");
2213 
2214  PG_RETURN_POINTER(v);
2215 }
2216 
2217 #endif
2218 
2219 
2220 /*
2221 ** The BOX32DF key must be defined as a PostgreSQL type, even though it is only
2222 ** ever used internally. These no-op stubs are used to bind the type.
2223 */
2225 Datum box2df_in(PG_FUNCTION_ARGS)
2226 {
2227  ereport(ERROR,(errcode(ERRCODE_FEATURE_NOT_SUPPORTED),
2228  errmsg("function box2df_in not implemented")));
2229  PG_RETURN_POINTER(NULL);
2230 }
2231 
2233 Datum box2df_out(PG_FUNCTION_ARGS)
2234 {
2235  BOX2DF *box = (BOX2DF *) PG_GETARG_POINTER(0);
2236  char *result = box2df_to_string(box);
2237  PG_RETURN_CSTRING(result);
2238 }
static float non_negative(float val)
Datum gserialized_above_2d(PG_FUNCTION_ARGS)
Datum gserialized_overabove_2d(PG_FUNCTION_ARGS)
static void box2df_validate(BOX2DF *b)
uint8_t data[1]
Definition: liblwgeom.h:368
#define LIMIT_RATIO
Datum gserialized_gist_penalty_2d(PG_FUNCTION_ARGS)
static BOX2DF * box2df_copy(BOX2DF *b)
static bool box2df_overlaps(const BOX2DF *a, const BOX2DF *b)
static double box2df_distance_node_centroid(const BOX2DF *node, const BOX2DF *query)
Calculate the The node_box_edge->query_centroid distance between the boxes.
Datum gserialized_gist_compress_2d(PG_FUNCTION_ARGS)
LWGEOM * lwgeom_from_gserialized(const GSERIALIZED *g)
Allocate a new LWGEOM from a GSERIALIZED.
static bool box2df_left(const BOX2DF *a, const BOX2DF *b)
double xmax
Definition: liblwgeom.h:277
static void fallbackSplit(GistEntryVector *entryvec, GIST_SPLITVEC *v)
void lwgeom_free(LWGEOM *geom)
Definition: lwgeom.c:1050
#define LW_SUCCESS
Definition: liblwgeom.h:65
static bool box2df_below(const BOX2DF *a, const BOX2DF *b)
Datum gserialized_left_2d(PG_FUNCTION_ARGS)
static bool gserialized_gist_consistent_leaf_2d(BOX2DF *key, BOX2DF *query, StrategyNumber strategy)
static void g_box_consider_split(ConsiderSplitContext *context, int dimNum, float rightLower, int minLeftCount, float leftUpper, int maxLeftCount)
static int interval_cmp_lower(const void *i1, const void *i2)
static void adjustBox(BOX2DF *b, BOX2DF *addon)
static bool box2df_right(const BOX2DF *a, const BOX2DF *b)
float next_float_down(double d)
Definition: lwgeom_api.c:143
static bool box2df_above(const BOX2DF *a, const BOX2DF *b)
static float box2df_size(const BOX2DF *a)
bool(* box2df_predicate)(const BOX2DF *a, const BOX2DF *b)
static bool box2df_within(const BOX2DF *a, const BOX2DF *b)
#define LW_FAILURE
Definition: liblwgeom.h:64
float next_float_up(double d)
Definition: lwgeom_api.c:159
int lwgeom_calculate_gbox(const LWGEOM *lwgeom, GBOX *gbox)
Calculate bounding box of a geometry, automatically taking into account whether it is cartesian or ge...
Definition: lwgeom.c:652
static int box2df_from_gbox_p(GBOX *box, BOX2DF *a)
double ymin
Definition: liblwgeom.h:278
Datum gserialized_overlaps_2d(PG_FUNCTION_ARGS)
Datum gserialized_distance_centroid_2d(PG_FUNCTION_ARGS)
Datum gserialized_overbelow_2d(PG_FUNCTION_ARGS)
double xmin
Definition: liblwgeom.h:276
static double box2df_distance(const BOX2DF *a, const BOX2DF *b)
Calculate the box->box distance.
static double box2df_distance_leaf_centroid(const BOX2DF *a, const BOX2DF *b)
Calculate the centroid->centroid distance between the boxes.
#define FLAGS_GET_BBOX(flags)
Definition: liblwgeom.h:126
static char * box2df_to_string(const BOX2DF *a)
#define LW_FALSE
Definition: liblwgeom.h:62
#define LW_TRUE
Return types for functions with status returns.
Definition: liblwgeom.h:61
Datum gserialized_below_2d(PG_FUNCTION_ARGS)
Datum gserialized_gist_same_2d(PG_FUNCTION_ARGS)
static void box2df_merge(BOX2DF *b_union, BOX2DF *b_new)
Datum gserialized_within_2d(PG_FUNCTION_ARGS)
static int gserialized_datum_get_box2df_p(Datum gsdatum, BOX2DF *box2df)
Peak into a GSERIALIZED datum to find the bounding box.
Datum gserialized_overleft_2d(PG_FUNCTION_ARGS)
Datum gserialized_gist_consistent_2d(PG_FUNCTION_ARGS)
Datum box2df_out(PG_FUNCTION_ARGS)
PG_FUNCTION_INFO_V1(gserialized_distance_centroid_2d)
double ymax
Definition: liblwgeom.h:279
Datum gserialized_same_2d(PG_FUNCTION_ARGS)
static int interval_cmp_upper(const void *i1, const void *i2)
static double pt_distance(double ax, double ay, double bx, double by)
Datum distance(PG_FUNCTION_ARGS)
static bool box2df_overright(const BOX2DF *a, const BOX2DF *b)
static bool box2df_contains(const BOX2DF *a, const BOX2DF *b)
static int common_entry_cmp(const void *i1, const void *i2)
Datum gserialized_gist_picksplit_2d(PG_FUNCTION_ARGS)
static bool box2df_overabove(const BOX2DF *a, const BOX2DF *b)
static bool box2df_overbelow(const BOX2DF *a, const BOX2DF *b)
Datum gserialized_gist_decompress_2d(PG_FUNCTION_ARGS)
static bool box2df_overleft(const BOX2DF *a, const BOX2DF *b)
#define FALSE
Definition: dbfopen.c:168
Datum gserialized_overright_2d(PG_FUNCTION_ARGS)
#define PLACE_RIGHT(box, off)
static float box2df_union_size(const BOX2DF *a, const BOX2DF *b)
static int gserialized_datum_predicate_2d(Datum gs1, Datum gs2, box2df_predicate predicate)
Support function.
static bool box2df_equals(const BOX2DF *a, const BOX2DF *b)
Datum gserialized_right_2d(PG_FUNCTION_ARGS)
Datum box2df_in(PG_FUNCTION_ARGS)
#define TRUE
Definition: dbfopen.c:169
Datum gserialized_contains_2d(PG_FUNCTION_ARGS)
Datum gserialized_gist_union_2d(PG_FUNCTION_ARGS)
Datum gserialized_gist_distance_2d(PG_FUNCTION_ARGS)
#define PLACE_LEFT(box, off)
static bool gserialized_gist_consistent_internal_2d(BOX2DF *key, BOX2DF *query, StrategyNumber strategy)
static float box_penalty(BOX2DF *original, BOX2DF *new)
Datum gserialized_distance_box_2d(PG_FUNCTION_ARGS)
This library is the generic geometry handling section of PostGIS.
uint8_t flags
Definition: liblwgeom.h:367