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