LCOV - code coverage report
Current view: top level - geom/monetdb5 - geom_atoms.c (source / functions) Hit Total Coverage
Test: coverage.info Lines: 426 793 53.7 %
Date: 2024-10-03 20:03:20 Functions: 52 72 72.2 %

          Line data    Source code
       1             : /*
       2             :  * SPDX-License-Identifier: MPL-2.0
       3             :  *
       4             :  * This Source Code Form is subject to the terms of the Mozilla Public
       5             :  * License, v. 2.0.  If a copy of the MPL was not distributed with this
       6             :  * file, You can obtain one at http://mozilla.org/MPL/2.0/.
       7             :  *
       8             :  * Copyright 2024 MonetDB Foundation;
       9             :  * Copyright August 2008 - 2023 MonetDB B.V.;
      10             :  * Copyright 1997 - July 2008 CWI.
      11             :  */
      12             : 
      13             : #include "geom_atoms.h"
      14             : 
      15             : /***********************************************/
      16             : /************* wkb type functions **************/
      17             : /***********************************************/
      18             : 
      19             : /* Creates the string representation (WKT) of a WKB */
      20             : /* return length of resulting string. */
      21             : ssize_t
      22         982 : wkbTOSTR(char **geomWKT, size_t *len, const void *GEOMWKB, bool external)
      23             : {
      24         982 :         const wkb *geomWKB = GEOMWKB;
      25         982 :         char *wkt = NULL;
      26         982 :         size_t dstStrLen = 5;   /* "nil" */
      27             : 
      28             :         /* from WKB to GEOSGeometry */
      29         982 :         GEOSGeom geosGeometry = wkb2geos(geomWKB);
      30             : 
      31         982 :         if (geosGeometry) {
      32         982 :                 size_t l;
      33         982 :                 GEOSWKTWriter *WKT_wr = GEOSWKTWriter_create_r(geoshandle);
      34             :                 //set the number of dimensions in the writer so that it can
      35             :                 //read correctly the geometry coordinates
      36         982 :                 GEOSWKTWriter_setOutputDimension_r(geoshandle, WKT_wr, GEOSGeom_getCoordinateDimension_r(geoshandle, geosGeometry));
      37         982 :                 GEOSWKTWriter_setTrim_r(geoshandle, WKT_wr, 1);
      38         982 :                 wkt = GEOSWKTWriter_write_r(geoshandle, WKT_wr, geosGeometry);
      39         982 :                 if (wkt == NULL) {
      40           0 :                         GDKerror("GEOSWKTWriter_write failed\n");
      41           0 :                         return -1;
      42             :                 }
      43         982 :                 GEOSWKTWriter_destroy_r(geoshandle, WKT_wr);
      44         982 :                 GEOSGeom_destroy_r(geoshandle, geosGeometry);
      45             : 
      46         982 :                 l = strlen(wkt);
      47         982 :                 dstStrLen = l;
      48         982 :                 if (external)
      49         831 :                         dstStrLen += 2; /* add quotes */
      50         982 :                 if (*len < dstStrLen + 1 || *geomWKT == NULL) {
      51         100 :                         *len = dstStrLen + 1;
      52         100 :                         GDKfree(*geomWKT);
      53         100 :                         if ((*geomWKT = GDKmalloc(*len)) == NULL) {
      54           0 :                                 GEOSFree_r(geoshandle, wkt);
      55           0 :                                 return -1;
      56             :                         }
      57             :                 }
      58         982 :                 if (external)
      59         831 :                         snprintf(*geomWKT, *len, "\"%s\"", wkt);
      60             :                 else
      61         151 :                         strcpy(*geomWKT, wkt);
      62         982 :                 GEOSFree_r(geoshandle, wkt);
      63             : 
      64         982 :                 return (ssize_t) dstStrLen;
      65             :         }
      66             : 
      67             :         /* geosGeometry == NULL */
      68           0 :         if (*len < 4 || *geomWKT == NULL) {
      69           0 :                 GDKfree(*geomWKT);
      70           0 :                 if ((*geomWKT = GDKmalloc(*len = 4)) == NULL)
      71             :                         return -1;
      72             :         }
      73           0 :         if (external) {
      74           0 :                 strcpy(*geomWKT, "nil");
      75           0 :                 return 3;
      76             :         }
      77           0 :         strcpy(*geomWKT, str_nil);
      78           0 :         return 1;
      79             : }
      80             : 
      81             : ssize_t
      82           2 : wkbFROMSTR(const char *geomWKT, size_t *len, void **GEOMWKB, bool external)
      83             : {
      84           2 :         wkb **geomWKB = (wkb **) GEOMWKB;
      85           2 :         size_t parsedBytes;
      86           2 :         str err;
      87             : 
      88           2 :         if (external && strncmp(geomWKT, "nil", 3) == 0) {
      89           0 :                 *geomWKB = wkbNULLcopy();
      90           0 :                 if (*geomWKB == NULL)
      91             :                         return -1;
      92             :                 return 3;
      93             :         }
      94           2 :         err = wkbFROMSTR_withSRID(geomWKT, len, geomWKB, 0, &parsedBytes);
      95           2 :         if (err != MAL_SUCCEED) {
      96           1 :                 GDKerror("%s", getExceptionMessageAndState(err));
      97           1 :                 freeException(err);
      98           1 :                 return -1;
      99             :         }
     100           1 :         return (ssize_t) parsedBytes;
     101             : }
     102             : 
     103             : BUN
     104        9001 : wkbHASH(const void *W)
     105             : {
     106        9001 :         const wkb *w = W;
     107        9001 :         int i;
     108        9001 :         BUN h = 0;
     109             : 
     110      158111 :         for (i = 0; i < (w->len - 1); i += 2) {
     111      149110 :                 BUN a = ((unsigned char *) w->data)[i];
     112      149110 :                 BUN b = ((unsigned char *) w->data)[i + 1];
     113             : #if '\377' < 0                                       /* char is signed? */
     114             :                 /* maybe sign extend */
     115      149110 :                 if (a & 0x80)
     116       27436 :                         a |= ~(BUN)0x7f;
     117      149110 :                 if (b & 0x80)
     118       27188 :                         b |= ~(BUN)0x7f;
     119             : #endif
     120      149110 :                 h = (h << 3) ^ (h >> 11) ^ (h >> 17) ^ (b << 8) ^ a;
     121             :         }
     122        9001 :         return h;
     123             : }
     124             : 
     125             : /* returns a pointer to a null wkb */
     126             : const void *
     127         315 : wkbNULL(void)
     128             : {
     129         315 :         return &wkb_nil;
     130             : }
     131             : 
     132             : int
     133        8192 : wkbCOMP(const void *L, const void *R)
     134             : {
     135        8192 :         const wkb *l = L, *r = R;
     136             : 
     137        8192 :         if (l->srid != r->srid)
     138             :                 return -1;
     139             : 
     140        7783 :         int len = l->len;
     141             : 
     142        7783 :         if (len != r->len)
     143        4906 :                 return len - r->len;
     144             : 
     145        2877 :         if (len == ~(int) 0)
     146             :                 return (0);
     147             : 
     148        2552 :         return memcmp(l->data, r->data, len);
     149             : }
     150             : 
     151             : /* read wkb from log */
     152             : void *
     153         469 : wkbREAD(void *A, size_t *dstlen, stream *s, size_t cnt)
     154             : {
     155         469 :         wkb *a = A;
     156         469 :         int len;
     157         469 :         int srid;
     158             : 
     159         469 :         (void) cnt;
     160         469 :         assert(cnt == 1);
     161         469 :         if (mnstr_readInt(s, &len) != 1)
     162             :                 return NULL;
     163         469 :         if (mnstr_readInt(s, &srid) != 1)
     164             :                 return NULL;
     165         469 :         size_t wkblen = (size_t) wkb_size(len);
     166         469 :         if (a == NULL || *dstlen < wkblen) {
     167           0 :                 if ((a = GDKrealloc(a, wkblen)) == NULL)
     168             :                         return NULL;
     169           0 :                 *dstlen = wkblen;
     170             :         }
     171         469 :         a->len = len;
     172         469 :         a->srid = srid;
     173         469 :         if (len > 0 && mnstr_read(s, (char *) a->data, len, 1) != 1) {
     174           0 :                 GDKfree(a);
     175           0 :                 return NULL;
     176             :         }
     177             :         return a;
     178             : }
     179             : 
     180             : /* write wkb to log */
     181             : gdk_return
     182         473 : wkbWRITE(const void *A, stream *s, size_t cnt)
     183             : {
     184         473 :         const wkb *a = A;
     185         473 :         int len = a->len;
     186         473 :         int srid = a->srid;
     187             : 
     188         473 :         (void) cnt;
     189         473 :         assert(cnt == 1);
     190         473 :         if (!mnstr_writeInt(s, len))    /* 64bit: check for overflow */
     191             :                 return GDK_FAIL;
     192         473 :         if (!mnstr_writeInt(s, srid))   /* 64bit: check for overflow */
     193             :                 return GDK_FAIL;
     194         943 :         if (len > 0 &&               /* 64bit: check for overflow */
     195         470 :             mnstr_write(s, (char *) a->data, len, 1) < 0)
     196             :                 return GDK_FAIL;
     197             :         return GDK_SUCCEED;
     198             : }
     199             : 
     200             : var_t
     201        1038 : wkbPUT(BAT *b, var_t *bun, const void *VAL)
     202             : {
     203        1038 :         const wkb *val = VAL;
     204        1038 :         char *base;
     205             : 
     206        1038 :         *bun = HEAP_malloc(b, wkb_size(val->len));
     207        1038 :         base = b->tvheap->base;
     208        1038 :         if (*bun != (var_t) -1) {
     209        1038 :                 memcpy(&base[*bun], val, wkb_size(val->len));
     210        1037 :                 b->tvheap->dirty = true;
     211             :         }
     212        1037 :         return *bun;
     213             : }
     214             : 
     215             : void
     216           0 : wkbDEL(Heap *h, var_t *index)
     217             : {
     218           0 :         HEAP_free(h, *index);
     219           0 : }
     220             : 
     221             : size_t
     222        1370 : wkbLENGTH(const void *P)
     223             : {
     224        1370 :         const wkb *p = P;
     225        1370 :         var_t len = wkb_size(p->len);
     226        1370 :         assert(len <= GDK_int_max);
     227        1370 :         return (size_t) len;
     228             : }
     229             : 
     230             : gdk_return
     231         971 : wkbHEAP(Heap *heap, size_t capacity)
     232             : {
     233         971 :         return HEAP_initialize(heap, capacity, 0, (int) sizeof(var_t));
     234             : }
     235             : 
     236             : /* Non-atom WKB functions */
     237             : wkb *
     238           3 : wkbNULLcopy(void)
     239             : {
     240           3 :         wkb *n = GDKmalloc(sizeof(wkb_nil));
     241           3 :         if (n)
     242           3 :                 *n = wkb_nil;
     243           3 :         return n;
     244             : }
     245             : 
     246             : wkb *
     247           1 : wkbCopy(const wkb* src)
     248             : {
     249           1 :         wkb *n = GDKmalloc(wkb_size(src->len));
     250           1 :         if (n) {
     251           1 :                 n->len = src->len;
     252           1 :                 n->srid = src->srid;
     253           1 :                 memcpy(n->data, src->data, src->len);
     254             :         }
     255           1 :         return n;
     256             : }
     257             : 
     258             : /* returns the size of variable-sized atom wkb */
     259             : var_t
     260        6310 : wkb_size(size_t len)
     261             : {
     262        6310 :         if (len == ~(size_t) 0)
     263         187 :                 len = 0;
     264        6310 :         assert(offsetof(wkb, data) + len <= VAR_MAX);
     265        6310 :         return (var_t) (offsetof(wkb, data) + len);
     266             : }
     267             : 
     268             : /* Creates WKB representation (including srid) from WKT representation */
     269             : /* return number of parsed characters. */
     270             : str
     271        1035 : wkbFROMSTR_withSRID(const char *geomWKT, size_t *len, wkb **geomWKB, int srid, size_t *nread)
     272             : {
     273        1035 :         GEOSGeom geosGeometry = NULL;   /* The geometry object that is parsed from the src string. */
     274        1035 :         GEOSWKTReader *WKT_reader;
     275        1035 :         static const char polyhedralSurface[] = "POLYHEDRALSURFACE";
     276        1035 :         static const char multiPolygon[] = "MULTIPOLYGON";
     277        1035 :         char *geomWKT_new = NULL;
     278        1035 :         size_t parsedCharacters = 0;
     279             : 
     280        1035 :         *nread = 0;
     281             : 
     282             :         /* we always allocate new memory */
     283        1035 :         GDKfree(*geomWKB);
     284        1035 :         *len = 0;
     285        1035 :         *geomWKB = NULL;
     286             : 
     287        1035 :         if (strNil(geomWKT)) {
     288           0 :                 *geomWKB = wkbNULLcopy();
     289           0 :                 if (*geomWKB == NULL)
     290           0 :                         throw(MAL, "wkb.FromText", SQLSTATE(HY013) MAL_MALLOC_FAIL);
     291           0 :                 *len = sizeof(wkb_nil);
     292           0 :                 return MAL_SUCCEED;
     293             :         }
     294             :         //check whether the representation is binary (hex)
     295        1035 :         if (geomWKT[0] == '0') {
     296         289 :                 str ret = wkbFromBinary(geomWKB, &geomWKT);
     297             : 
     298         289 :                 if (ret != MAL_SUCCEED)
     299             :                         return ret;
     300         289 :                 *nread = strlen(geomWKT);
     301         289 :                 *len = (size_t) wkb_size((*geomWKB)->len);
     302         289 :                 return MAL_SUCCEED;
     303             :         }
     304             :         //check whether the geometry type is polyhedral surface
     305             :         //geos cannot handle this type of geometry but since it is
     306             :         //a special type of multipolygon I just change the type before
     307             :         //continuing. Of course this means that isValid for example does
     308             :         //not work correctly.
     309         746 :         if (strncasecmp(geomWKT, polyhedralSurface, strlen(polyhedralSurface)) == 0) {
     310           2 :                 size_t sizeOfInfo = strlen(geomWKT) - strlen(polyhedralSurface) + strlen(multiPolygon) + 1;
     311           2 :                 geomWKT_new = GDKmalloc(sizeOfInfo);
     312           2 :                 if (geomWKT_new == NULL)
     313           0 :                         throw(MAL, "wkb.FromText", SQLSTATE(HY013) MAL_MALLOC_FAIL);
     314           2 :                 snprintf(geomWKT_new, sizeOfInfo, "%s%s", multiPolygon, geomWKT + strlen(polyhedralSurface));
     315           2 :                 geomWKT = geomWKT_new;
     316             :         }
     317             :         ////////////////////////// UP TO HERE ///////////////////////////
     318             : 
     319         746 :         WKT_reader = GEOSWKTReader_create_r(geoshandle);
     320         745 :         if (WKT_reader == NULL) {
     321           0 :                 if (geomWKT_new)
     322           0 :                         GDKfree(geomWKT_new);
     323           0 :                 throw(MAL, "wkb.FromText", SQLSTATE(38000) "Geos operation GEOSWKTReader_create failed");
     324             :         }
     325         745 :         geosGeometry = GEOSWKTReader_read_r(geoshandle, WKT_reader, geomWKT);
     326         746 :         GEOSWKTReader_destroy_r(geoshandle, WKT_reader);
     327             : 
     328         746 :         if (geosGeometry == NULL) {
     329          55 :                 if (geomWKT_new)
     330           0 :                         GDKfree(geomWKT_new);
     331          55 :                 throw(MAL, "wkb.FromText", SQLSTATE(38000) "Geos operation GEOSWKTReader_read failed");
     332             :         }
     333             : 
     334         691 :         if (GEOSGeomTypeId_r(geoshandle, geosGeometry) == -1) {
     335           0 :                 if (geomWKT_new)
     336           0 :                         GDKfree(geomWKT_new);
     337           0 :                 GEOSGeom_destroy_r(geoshandle, geosGeometry);
     338           0 :                 throw(MAL, "wkb.FromText", SQLSTATE(38000) "Geos operation GEOSGeomTypeId failed");
     339             :         }
     340             : 
     341         691 :         GEOSSetSRID_r(geoshandle, geosGeometry, srid);
     342             :         /* the srid was lost with the transformation of the GEOSGeom to wkb
     343             :          * so we decided to store it in the wkb */
     344             : 
     345             :         /* we have a GEOSGeometry with number of coordinates and SRID and we
     346             :          * want to get the wkb out of it */
     347         690 :         *geomWKB = geos2wkb(geosGeometry);
     348         691 :         GEOSGeom_destroy_r(geoshandle, geosGeometry);
     349         691 :         if (*geomWKB == NULL) {
     350           0 :                 if (geomWKT_new)
     351           0 :                         GDKfree(geomWKT_new);
     352           0 :                 throw(MAL, "wkb.FromText", SQLSTATE(38000) "Geos operation geos2wkb failed");
     353             :         }
     354             : 
     355         691 :         *len = (size_t) wkb_size((*geomWKB)->len);
     356         690 :         parsedCharacters = strlen(geomWKT);
     357         690 :         assert(parsedCharacters <= GDK_int_max);
     358             : 
     359         690 :         GDKfree(geomWKT_new);
     360             : 
     361         690 :         *nread = parsedCharacters;
     362         690 :         return MAL_SUCCEED;
     363             : }
     364             : 
     365             : /***********************************************/
     366             : /************* mbr type functions **************/
     367             : /***********************************************/
     368             : 
     369             : #define MBR_WKTLEN 256
     370             : 
     371             : /* TOSTR: print atom in a string. */
     372             : /* return length of resulting string. */
     373             : ssize_t
     374          91 : mbrTOSTR(char **dst, size_t *len, const void *ATOM, bool external)
     375             : {
     376          91 :         const mbr *atom = ATOM;
     377          91 :         char tempWkt[MBR_WKTLEN];
     378          91 :         size_t dstStrLen;
     379             : 
     380          91 :         if (!is_mbr_nil(atom)) {
     381          91 :                 dstStrLen = (size_t) snprintf(tempWkt, MBR_WKTLEN,
     382             :                                               "BOX (%f %f, %f %f)",
     383          91 :                                               atom->xmin, atom->ymin,
     384          91 :                                               atom->xmax, atom->ymax);
     385             :         } else {
     386           0 :                 tempWkt[0] = 0; /* not used */
     387           0 :                 dstStrLen = 0;
     388             :         }
     389             : 
     390          91 :         if (*len < dstStrLen + 4 || *dst == NULL) {
     391           0 :                 GDKfree(*dst);
     392           0 :                 if ((*dst = GDKmalloc(*len = dstStrLen + 4)) == NULL)
     393             :                         return -1;
     394             :         }
     395             : 
     396          91 :         if (dstStrLen > 4) {
     397          91 :                 if (external) {
     398          75 :                         snprintf(*dst, *len, "\"%s\"", tempWkt);
     399          75 :                         dstStrLen += 2;
     400             :                 } else {
     401          16 :                         strcpy(*dst, tempWkt);
     402             :                 }
     403           0 :         } else if (external) {
     404           0 :                 strcpy(*dst, "nil");
     405           0 :                 dstStrLen = 3;
     406             :         } else {
     407           0 :                 strcpy(*dst, str_nil);
     408           0 :                 dstStrLen = 1;
     409             :         }
     410          91 :         return (ssize_t) dstStrLen;
     411             : }
     412             : 
     413             : /* FROMSTR: parse string to mbr. */
     414             : /* return number of parsed characters. */
     415             : ssize_t
     416           2 : mbrFROMSTR(const char *src, size_t *len, void **ATOM, bool external)
     417             : {
     418           2 :         mbr **atom = (mbr **) ATOM;
     419           2 :         size_t nchars = 0;      /* The number of characters parsed; the return value. */
     420           2 :         GEOSGeom geosMbr = NULL;        /* The geometry object that is parsed from the src string. */
     421           2 :         double xmin = 0, ymin = 0, xmax = 0, ymax = 0;
     422           2 :         const char *c;
     423             : 
     424           2 :         if (*len < sizeof(mbr) || *atom == NULL) {
     425           2 :                 GDKfree(*atom);
     426           2 :                 if ((*atom = GDKmalloc(*len = sizeof(mbr))) == NULL)
     427             :                         return -1;
     428             :         }
     429           2 :         if (external && strncmp(src, "nil", 3) == 0) {
     430           0 :                 **atom = mbrNIL;
     431           0 :                 return 3;
     432             :         }
     433           2 :         if (strNil(src)) {
     434           0 :                 **atom = mbrNIL;
     435           0 :                 return 1;
     436             :         }
     437             : 
     438           2 :         if ((strstr(src, "mbr") == src || strstr(src, "MBR") == src
     439           0 :              || strstr(src, "box") == src || strstr(src, "BOX") == src)
     440           2 :             && (c = strstr(src, "(")) != NULL) {
     441             :                 /* Parse the mbr */
     442           2 :                 if ((c - src) != 3 && (c - src) != 4) {
     443           0 :                         GDKerror("ParseException: Expected a string like 'MBR(0 0,1 1)' or 'MBR (0 0,1 1)'\n");
     444           0 :                         return -1;
     445             :                 }
     446             : 
     447           2 :                 if (sscanf(c, "(%lf %lf,%lf %lf)", &xmin, &ymin, &xmax, &ymax) != 4) {
     448           0 :                         GDKerror("ParseException: Not enough coordinates.\n");
     449           0 :                         return -1;
     450             :                 }
     451           0 :         } else if ((geosMbr = GEOSGeomFromWKT(src)) == NULL) {
     452           0 :                 GDKerror("GEOSGeomFromWKT failed\n");
     453           0 :                 return -1;
     454             :         }
     455             : 
     456           2 :         if (geosMbr == NULL) {
     457           2 :                 assert(GDK_flt_min <= xmin && xmin <= GDK_flt_max);
     458           2 :                 assert(GDK_flt_min <= xmax && xmax <= GDK_flt_max);
     459           2 :                 assert(GDK_flt_min <= ymin && ymin <= GDK_flt_max);
     460           2 :                 assert(GDK_flt_min <= ymax && ymax <= GDK_flt_max);
     461           2 :                 (*atom)->xmin = (float) xmin;
     462           2 :                 (*atom)->ymin = (float) ymin;
     463           2 :                 (*atom)->xmax = (float) xmax;
     464           2 :                 (*atom)->ymax = (float) ymax;
     465           2 :                 nchars = strlen(src);
     466             :         }
     467           2 :         if (geosMbr)
     468           0 :                 GEOSGeom_destroy_r(geoshandle, geosMbr);
     469           2 :         assert(nchars <= GDK_int_max);
     470           2 :         return (ssize_t) nchars;
     471             : }
     472             : 
     473             : /* HASH: compute a hash value. */
     474             : /* returns a positive integer hash value */
     475             : BUN
     476          10 : mbrHASH(const void *ATOM)
     477             : {
     478          10 :         const mbr *atom = ATOM;
     479          10 :         return ATOMhash(TYPE_flt, &atom->xmin) ^ ATOMhash(TYPE_flt, &atom->ymin) ^
     480          10 :                 ATOMhash(TYPE_flt, &atom->xmax) ^ ATOMhash(TYPE_flt, &atom->ymax);
     481             : }
     482             : 
     483             : const void *
     484         315 : mbrNULL(void)
     485             : {
     486         315 :         return &mbrNIL;
     487             : }
     488             : 
     489             : /* COMP: compare two mbrs. */
     490             : /* returns int <0 if l<r, 0 if l==r, >0 else */
     491             : int
     492        1824 : mbrCOMP(const void *L, const void *R)
     493             : {
     494             :         /* simple lexicographical ordering on (x,y) */
     495        1824 :         const mbr *l = L, *r = R;
     496        1824 :         int res;
     497        1824 :         if (is_mbr_nil(l))
     498         156 :                 return -!is_mbr_nil(r);
     499        1669 :         if (is_mbr_nil(r))
     500             :                 return 1;
     501        1006 :         if (l->xmin == r->xmin)
     502         536 :                 res = (l->ymin < r->ymin) ? -1 : (l->ymin != r->ymin);
     503             :         else
     504         470 :                 res = (l->xmin < r->xmin) ? -1 : 1;
     505         375 :         if (res == 0) {
     506         360 :                 if (l->xmax == r->xmax)
     507         262 :                         res = (l->ymax < r->ymax) ? -1 : (l->ymax != r->ymax);
     508             :                 else
     509          98 :                         res = (l->xmax < r->xmax) ? -1 : 1;
     510             :         }
     511             :         return res;
     512             : }
     513             : 
     514             : /* read mbr from log */
     515             : void *
     516          10 : mbrREAD(void *A, size_t *dstlen, stream *s, size_t cnt)
     517             : {
     518          10 :         mbr *a = A;
     519          10 :         mbr *c;
     520          10 :         size_t i;
     521          10 :         int v[4];
     522          10 :         flt vals[4];
     523             : 
     524          10 :         if (a == NULL || *dstlen < cnt * sizeof(mbr)) {
     525           0 :                 if ((a = GDKrealloc(a, cnt * sizeof(mbr))) == NULL)
     526             :                         return NULL;
     527           0 :                 *dstlen = cnt * sizeof(mbr);
     528             :         }
     529          24 :         for (i = 0, c = a; i < cnt; i++, c++) {
     530          14 :                 if (!mnstr_readIntArray(s, v, 4)) {
     531           0 :                         if (a != A)
     532           0 :                                 GDKfree(a);
     533           0 :                         return NULL;
     534             :                 }
     535          14 :                 memcpy(vals, v, 4 * sizeof(int));
     536          14 :                 c->xmin = vals[0];
     537          14 :                 c->ymin = vals[1];
     538          14 :                 c->xmax = vals[2];
     539          14 :                 c->ymax = vals[3];
     540             :         }
     541             :         return a;
     542             : }
     543             : 
     544             : /* write mbr to log */
     545             : gdk_return
     546          14 : mbrWRITE(const void *C, stream *s, size_t cnt)
     547             : {
     548          14 :         const mbr *c = C;
     549          14 :         size_t i;
     550          14 :         flt vals[4];
     551          14 :         int v[4];
     552             : 
     553          28 :         for (i = 0; i < cnt; i++, c++) {
     554          14 :                 vals[0] = c->xmin;
     555          14 :                 vals[1] = c->ymin;
     556          14 :                 vals[2] = c->xmax;
     557          14 :                 vals[3] = c->ymax;
     558          14 :                 memcpy(v, vals, 4 * sizeof(int));
     559          14 :                 if (!mnstr_writeIntArray(s, v, 4))
     560             :                         return GDK_FAIL;
     561             :         }
     562             :         return GDK_SUCCEED;
     563             : }
     564             : 
     565             : /* Non-atom mbr functions */
     566             : /* Check if fixed-sized atom mbr is null */
     567             : bool
     568        5369 : is_mbr_nil(const mbr *m)
     569             : {
     570        5369 :         return (m == NULL || is_flt_nil(m->xmin) || is_flt_nil(m->ymin) || is_flt_nil(m->xmax) || is_flt_nil(m->ymax));
     571             : }
     572             : 
     573             : /* MBR FUNCTIONS */
     574             : /* MBR */
     575             : 
     576             : /* Creates the mbr for the given geom_geometry. */
     577             : str
     578         692 : wkbMBR(mbr **geomMBR, wkb **geomWKB)
     579             : {
     580         692 :         GEOSGeom geosGeometry;
     581         692 :         str ret = MAL_SUCCEED;
     582         692 :         bit empty;
     583             : 
     584             :         //check if the geometry is nil
     585         692 :         if (is_wkb_nil(*geomWKB)) {
     586           2 :                 if ((*geomMBR = GDKmalloc(sizeof(mbr))) == NULL)
     587           0 :                         throw(MAL, "geom.MBR", SQLSTATE(HY013) MAL_MALLOC_FAIL);
     588           2 :                 **geomMBR = mbrNIL;
     589           2 :                 return MAL_SUCCEED;
     590             :         }
     591             :         //check if the geometry is empty
     592         689 :         if ((ret = wkbIsEmpty(&empty, geomWKB)) != MAL_SUCCEED) {
     593             :                 return ret;
     594             :         }
     595         690 :         if (empty) {
     596           5 :                 if ((*geomMBR = GDKmalloc(sizeof(mbr))) == NULL)
     597           0 :                         throw(MAL, "geom.MBR", SQLSTATE(HY013) MAL_MALLOC_FAIL);
     598           5 :                 **geomMBR = mbrNIL;
     599           5 :                 return MAL_SUCCEED;
     600             :         }
     601             : 
     602         685 :         geosGeometry = wkb2geos(*geomWKB);
     603         686 :         if (geosGeometry == NULL) {
     604           0 :                 *geomMBR = NULL;
     605           0 :                 throw(MAL, "geom.MBR", SQLSTATE(38000) "Geos problem converting GEOS to WKB");
     606             :         }
     607             : 
     608         686 :         *geomMBR = mbrFromGeos(geosGeometry);
     609             : 
     610         684 :         GEOSGeom_destroy_r(geoshandle, geosGeometry);
     611             : 
     612         686 :         if (*geomMBR == NULL || is_mbr_nil(*geomMBR)) {
     613           0 :                 GDKfree(*geomMBR);
     614           0 :                 *geomMBR = NULL;
     615           0 :                 throw(MAL, "wkb.mbr", SQLSTATE(38000) "Geos failed to create mbr");
     616             :         }
     617             : 
     618             :         return MAL_SUCCEED;
     619             : }
     620             : 
     621             : str
     622          23 : wkbBox2D(mbr **box, wkb **point1, wkb **point2)
     623             : {
     624          23 :         GEOSGeom point1_geom, point2_geom;
     625          23 :         double xmin = 0.0, ymin = 0.0, xmax = 0.0, ymax = 0.0;
     626          23 :         str err = MAL_SUCCEED;
     627             : 
     628             :         //check null input
     629          23 :         if (is_wkb_nil(*point1) || is_wkb_nil(*point2)) {
     630           1 :                 if ((*box = GDKmalloc(sizeof(mbr))) == NULL)
     631           0 :                         throw(MAL, "geom.MakeBox2D", SQLSTATE(HY013) MAL_MALLOC_FAIL);
     632           1 :                 **box = mbrNIL;
     633           1 :                 return MAL_SUCCEED;
     634             :         }
     635             :         //check input not point geometries
     636          22 :         point1_geom = wkb2geos(*point1);
     637          22 :         point2_geom = wkb2geos(*point2);
     638          22 :         if (point1_geom == NULL || point2_geom == NULL) {
     639           0 :                 if (point1_geom)
     640           0 :                         GEOSGeom_destroy_r(geoshandle, point1_geom);
     641           0 :                 if (point2_geom)
     642           0 :                         GEOSGeom_destroy_r(geoshandle, point2_geom);
     643           0 :                 throw(MAL, "geom.MakeBox2D", SQLSTATE(HY013) MAL_MALLOC_FAIL);
     644             :         }
     645          44 :         if (GEOSGeomTypeId_r(geoshandle, point1_geom) + 1 != wkbPoint_mdb ||
     646          22 :             GEOSGeomTypeId_r(geoshandle, point2_geom) + 1 != wkbPoint_mdb) {
     647           1 :                 err = createException(MAL, "geom.MakeBox2D", SQLSTATE(38000) "Geometries should be points");
     648          42 :         } else if (GEOSGeomGetX_r(geoshandle, point1_geom, &xmin) == -1 ||
     649          42 :                    GEOSGeomGetY_r(geoshandle, point1_geom, &ymin) == -1 ||
     650          42 :                    GEOSGeomGetX_r(geoshandle, point2_geom, &xmax) == -1 ||
     651          21 :                    GEOSGeomGetY_r(geoshandle, point2_geom, &ymax) == -1) {
     652             : 
     653           0 :                 err = createException(MAL, "geom.MakeBox2D", SQLSTATE(38000) "Geos error in reading the points' coordinates");
     654             :         } else {
     655             :                 //Assign the coordinates. Ensure that they are in correct order
     656          21 :                 *box = GDKmalloc(sizeof(mbr));
     657          21 :                 if (*box == NULL) {
     658           0 :                         err = createException(MAL, "geom.MakeBox2D", SQLSTATE(HY013) MAL_MALLOC_FAIL);
     659             :                 } else {
     660          21 :                         (*box)->xmin = (float) (xmin < xmax ? xmin : xmax);
     661          21 :                         (*box)->ymin = (float) (ymin < ymax ? ymin : ymax);
     662          21 :                         (*box)->xmax = (float) (xmax > xmin ? xmax : xmin);
     663          31 :                         (*box)->ymax = (float) (ymax > ymin ? ymax : ymin);
     664             :                 }
     665             :         }
     666          22 :         GEOSGeom_destroy_r(geoshandle, point1_geom);
     667          22 :         GEOSGeom_destroy_r(geoshandle, point2_geom);
     668             : 
     669          22 :         return err;
     670             : }
     671             : 
     672             : static str
     673         162 : mbrrelation_wkb(bit *out, wkb **geom1WKB, wkb **geom2WKB, str (*func)(bit *, mbr **, mbr **))
     674             : {
     675         162 :         mbr *geom1MBR = NULL, *geom2MBR = NULL;
     676         162 :         str ret = MAL_SUCCEED;
     677             : 
     678         162 :         if (is_wkb_nil(*geom1WKB) || is_wkb_nil(*geom2WKB)) {
     679          29 :                 *out = bit_nil;
     680          29 :                 return MAL_SUCCEED;
     681             :         }
     682             : 
     683         133 :         ret = wkbMBR(&geom1MBR, geom1WKB);
     684         133 :         if (ret != MAL_SUCCEED) {
     685             :                 return ret;
     686             :         }
     687             : 
     688         133 :         ret = wkbMBR(&geom2MBR, geom2WKB);
     689         133 :         if (ret != MAL_SUCCEED) {
     690           0 :                 GDKfree(geom1MBR);
     691           0 :                 return ret;
     692             :         }
     693             : 
     694         133 :         ret = (*func) (out, &geom1MBR, &geom2MBR);
     695             : 
     696         133 :         GDKfree(geom1MBR);
     697         133 :         GDKfree(geom2MBR);
     698             : 
     699         133 :         return ret;
     700             : }
     701             : 
     702             : /*returns true if the two mbrs overlap */
     703             : str
     704          14 : mbrOverlaps(bit *out, mbr **b1, mbr **b2)
     705             : {
     706          14 :         if (is_mbr_nil(*b1) || is_mbr_nil(*b2))
     707           0 :                 *out = bit_nil;
     708             :         else                    //they cannot overlap if b2 is left, right, above or below b1
     709          28 :                 *out = !((*b2)->ymax < (*b1)->ymin ||
     710          12 :                          (*b2)->ymin > (*b1)->ymax ||
     711          12 :                          (*b2)->xmax < (*b1)->xmin ||
     712          12 :                          (*b2)->xmin > (*b1)->xmax);
     713          14 :         return MAL_SUCCEED;
     714             : }
     715             : 
     716             : /*returns true if the mbrs of the two geometries overlap */
     717             : str
     718           6 : mbrOverlaps_wkb(bit *out, wkb **geom1WKB, wkb **geom2WKB)
     719             : {
     720           6 :         return mbrrelation_wkb(out, geom1WKB, geom2WKB, mbrOverlaps);
     721             : }
     722             : 
     723             : /* returns true if b1 is above b2 */
     724             : str
     725           5 : mbrAbove(bit *out, mbr **b1, mbr **b2)
     726             : {
     727           5 :         if (is_mbr_nil(*b1) || is_mbr_nil(*b2))
     728           0 :                 *out = bit_nil;
     729             :         else
     730           5 :                 *out = ((*b1)->ymin > (*b2)->ymax);
     731           5 :         return MAL_SUCCEED;
     732             : }
     733             : 
     734             : /*returns true if the mbrs of geom1 is above the mbr of geom2 */
     735             : str
     736           6 : mbrAbove_wkb(bit *out, wkb **geom1WKB, wkb **geom2WKB)
     737             : {
     738           6 :         return mbrrelation_wkb(out, geom1WKB, geom2WKB, mbrAbove);
     739             : }
     740             : 
     741             : /* returns true if b1 is below b2 */
     742             : str
     743           5 : mbrBelow(bit *out, mbr **b1, mbr **b2)
     744             : {
     745           5 :         if (is_mbr_nil(*b1) || is_mbr_nil(*b2))
     746           0 :                 *out = bit_nil;
     747             :         else
     748           5 :                 *out = ((*b1)->ymax < (*b2)->ymin);
     749           5 :         return MAL_SUCCEED;
     750             : }
     751             : 
     752             : /*returns true if the mbrs of geom1 is below the mbr of geom2 */
     753             : str
     754           6 : mbrBelow_wkb(bit *out, wkb **geom1WKB, wkb **geom2WKB)
     755             : {
     756           6 :         return mbrrelation_wkb(out, geom1WKB, geom2WKB, mbrBelow);
     757             : }
     758             : 
     759             : /* returns true if box1 is left of box2 */
     760             : str
     761          65 : mbrLeft(bit *out, mbr **b1, mbr **b2)
     762             : {
     763          65 :         if (is_mbr_nil(*b1) || is_mbr_nil(*b2))
     764           0 :                 *out = bit_nil;
     765             :         else
     766          65 :                 *out = ((*b1)->xmax < (*b2)->xmin);
     767          65 :         return MAL_SUCCEED;
     768             : }
     769             : 
     770             : /*returns true if the mbrs of geom1 is on the left of the mbr of geom2 */
     771             : str
     772          84 : mbrLeft_wkb(bit *out, wkb **geom1WKB, wkb **geom2WKB)
     773             : {
     774          84 :         return mbrrelation_wkb(out, geom1WKB, geom2WKB, mbrLeft);
     775             : }
     776             : 
     777             : /* returns true if box1 is right of box2 */
     778             : str
     779          11 : mbrRight(bit *out, mbr **b1, mbr **b2)
     780             : {
     781          11 :         if (is_mbr_nil(*b1) || is_mbr_nil(*b2))
     782           0 :                 *out = bit_nil;
     783             :         else
     784          11 :                 *out = ((*b1)->xmin > (*b2)->xmax);
     785          11 :         return MAL_SUCCEED;
     786             : }
     787             : 
     788             : /*returns true if the mbrs of geom1 is on the right of the mbr of geom2 */
     789             : str
     790          12 : mbrRight_wkb(bit *out, wkb **geom1WKB, wkb **geom2WKB)
     791             : {
     792          12 :         return mbrrelation_wkb(out, geom1WKB, geom2WKB, mbrRight);
     793             : }
     794             : 
     795             : /* returns true if box1 overlaps or is above box2 when only the Y coordinate is considered*/
     796             : str
     797           0 : mbrOverlapOrAbove(bit *out, mbr **b1, mbr **b2)
     798             : {
     799           0 :         if (is_mbr_nil(*b1) || is_mbr_nil(*b2))
     800           0 :                 *out = bit_nil;
     801             :         else
     802           0 :                 *out = ((*b1)->ymin >= (*b2)->ymin);
     803           0 :         return MAL_SUCCEED;
     804             : }
     805             : 
     806             : /*returns true if the mbrs of geom1 overlaps or is above the mbr of geom2 */
     807             : str
     808           0 : mbrOverlapOrAbove_wkb(bit *out, wkb **geom1WKB, wkb **geom2WKB)
     809             : {
     810           0 :         return mbrrelation_wkb(out, geom1WKB, geom2WKB, mbrOverlapOrAbove);
     811             : }
     812             : 
     813             : /* returns true if box1 overlaps or is below box2 when only the Y coordinate is considered*/
     814             : str
     815           0 : mbrOverlapOrBelow(bit *out, mbr **b1, mbr **b2)
     816             : {
     817           0 :         if (is_mbr_nil(*b1) || is_mbr_nil(*b2))
     818           0 :                 *out = bit_nil;
     819             :         else
     820           0 :                 *out = ((*b1)->ymax <= (*b2)->ymax);
     821           0 :         return MAL_SUCCEED;
     822             : }
     823             : 
     824             : /*returns true if the mbrs of geom1 overlaps or is below the mbr of geom2 */
     825             : str
     826           0 : mbrOverlapOrBelow_wkb(bit *out, wkb **geom1WKB, wkb **geom2WKB)
     827             : {
     828           0 :         return mbrrelation_wkb(out, geom1WKB, geom2WKB, mbrOverlapOrBelow);
     829             : }
     830             : 
     831             : /* returns true if box1 overlaps or is left of box2 when only the X coordinate is considered*/
     832             : str
     833          10 : mbrOverlapOrLeft(bit *out, mbr **b1, mbr **b2)
     834             : {
     835          10 :         if (is_mbr_nil(*b1) || is_mbr_nil(*b2))
     836           0 :                 *out = bit_nil;
     837             :         else
     838          10 :                 *out = ((*b1)->xmax <= (*b2)->xmax);
     839          10 :         return MAL_SUCCEED;
     840             : }
     841             : 
     842             : /*returns true if the mbrs of geom1 overlaps or is on the left of the mbr of geom2 */
     843             : str
     844           6 : mbrOverlapOrLeft_wkb(bit *out, wkb **geom1WKB, wkb **geom2WKB)
     845             : {
     846           6 :         return mbrrelation_wkb(out, geom1WKB, geom2WKB, mbrOverlapOrLeft);
     847             : }
     848             : 
     849             : /* returns true if box1 overlaps or is right of box2 when only the X coordinate is considered*/
     850             : str
     851          10 : mbrOverlapOrRight(bit *out, mbr **b1, mbr **b2)
     852             : {
     853          10 :         if (is_mbr_nil(*b1) || is_mbr_nil(*b2))
     854           0 :                 *out = bit_nil;
     855             :         else
     856          10 :                 *out = ((*b1)->xmin >= (*b2)->xmin);
     857          10 :         return MAL_SUCCEED;
     858             : }
     859             : 
     860             : /*returns true if the mbrs of geom1 overlaps or is on the right of the mbr of geom2 */
     861             : str
     862           6 : mbrOverlapOrRight_wkb(bit *out, wkb **geom1WKB, wkb **geom2WKB)
     863             : {
     864           6 :         return mbrrelation_wkb(out, geom1WKB, geom2WKB, mbrOverlapOrRight);
     865             : }
     866             : 
     867             : /* returns true if b1 is contained in b2 */
     868             : str
     869          46 : mbrContained(bit *out, mbr **b1, mbr **b2)
     870             : {
     871          46 :         if (is_mbr_nil(*b1) || is_mbr_nil(*b2))
     872           0 :                 *out = bit_nil;
     873             :         else
     874          75 :                 *out = (((*b1)->xmin >= (*b2)->xmin) && ((*b1)->xmax <= (*b2)->xmax) && ((*b1)->ymin >= (*b2)->ymin) && ((*b1)->ymax <= (*b2)->ymax));
     875          46 :         return MAL_SUCCEED;
     876             : }
     877             : 
     878             : /*returns true if the mbrs of geom1 is contained in the mbr of geom2 */
     879             : str
     880          18 : mbrContained_wkb(bit *out, wkb **geom1WKB, wkb **geom2WKB)
     881             : {
     882          18 :         return mbrrelation_wkb(out, geom1WKB, geom2WKB, mbrContained);
     883             : }
     884             : 
     885             : /*returns true if b1 contains b2 */
     886             : str
     887          26 : mbrContains(bit *out, mbr **b1, mbr **b2)
     888             : {
     889          26 :         return mbrContained(out, b2, b1);
     890             : }
     891             : 
     892             : /*returns true if the mbrs of geom1 contains the mbr of geom2 */
     893             : str
     894           6 : mbrContains_wkb(bit *out, wkb **geom1WKB, wkb **geom2WKB)
     895             : {
     896           6 :         return mbrrelation_wkb(out, geom1WKB, geom2WKB, mbrContains);
     897             : }
     898             : 
     899             : /* returns true if the boxes are the same */
     900             : str
     901          31 : mbrEqual(bit *out, mbr **b1, mbr **b2)
     902             : {
     903          31 :         if (is_mbr_nil(*b1) && is_mbr_nil(*b2))
     904           0 :                 *out = 1;
     905          31 :         else if (is_mbr_nil(*b1) || is_mbr_nil(*b2))
     906           0 :                 *out = 0;
     907             :         else
     908          55 :                 *out = (((*b1)->xmin == (*b2)->xmin) && ((*b1)->xmax == (*b2)->xmax) && ((*b1)->ymin == (*b2)->ymin) && ((*b1)->ymax == (*b2)->ymax));
     909          31 :         return MAL_SUCCEED;
     910             : }
     911             : 
     912             : /*returns true if the mbrs of geom1 and the mbr of geom2 are the same */
     913             : str
     914          12 : mbrEqual_wkb(bit *out, wkb **geom1WKB, wkb **geom2WKB)
     915             : {
     916          12 :         return mbrrelation_wkb(out, geom1WKB, geom2WKB, mbrEqual);
     917             : }
     918             : 
     919             : str
     920          20 : mbrDiagonal(dbl *out, mbr **b)
     921             : {
     922          20 :         double side_a = .0, side_b = .0;
     923             : 
     924          20 :         if (is_mbr_nil(*b)) {
     925           0 :                 *out = dbl_nil;
     926           0 :                 return MAL_SUCCEED;
     927             :         }
     928             : 
     929          20 :         side_a = (*b)->xmax - (*b)->xmin;
     930          20 :         side_b = (*b)->ymax - (*b)->ymin;
     931             : 
     932          20 :         *out = sqrt(pow(side_a, 2.0) + pow(side_b, 2.0));
     933             : 
     934          20 :         return MAL_SUCCEED;
     935             : }
     936             : 
     937             : /* returns the Euclidean distance of the centroids of the boxes */
     938             : str
     939         135 : mbrDistance(dbl *out, mbr **b1, mbr **b2)
     940             : {
     941         135 :         double b1_Cx = 0.0, b1_Cy = 0.0, b2_Cx = 0.0, b2_Cy = 0.0;
     942             : 
     943         135 :         if (is_mbr_nil(*b1) || is_mbr_nil(*b2)) {
     944           0 :                 *out = dbl_nil;
     945           0 :                 return MAL_SUCCEED;
     946             :         }
     947             :         //compute the centroids of the two polygons
     948         135 :         b1_Cx = ((*b1)->xmin + (*b1)->xmax) / 2.0;
     949         135 :         b1_Cy = ((*b1)->ymin + (*b1)->ymax) / 2.0;
     950         135 :         b2_Cx = ((*b2)->xmin + (*b2)->xmax) / 2.0;
     951         135 :         b2_Cy = ((*b2)->ymin + (*b2)->ymax) / 2.0;
     952             : 
     953             :         //compute the euclidean distance
     954         135 :         *out = sqrt(pow(b2_Cx - b1_Cx, 2.0) + pow(b2_Cy - b1_Cy, 2.0));
     955             : 
     956         135 :         return MAL_SUCCEED;
     957             : }
     958             : 
     959             : /*returns the Euclidean distance of the centroids of the mbrs of the two geometries */
     960             : str
     961         170 : mbrDistance_wkb(dbl *out, wkb **geom1WKB, wkb **geom2WKB)
     962             : {
     963         170 :         mbr *geom1MBR = NULL, *geom2MBR = NULL;
     964         170 :         str ret = MAL_SUCCEED;
     965             : 
     966         170 :         if (is_wkb_nil(*geom1WKB) || is_wkb_nil(*geom2WKB)) {
     967          45 :                 *out = dbl_nil;
     968          45 :                 return MAL_SUCCEED;
     969             :         }
     970             : 
     971         125 :         ret = wkbMBR(&geom1MBR, geom1WKB);
     972         125 :         if (ret != MAL_SUCCEED) {
     973             :                 return ret;
     974             :         }
     975             : 
     976         125 :         ret = wkbMBR(&geom2MBR, geom2WKB);
     977         125 :         if (ret != MAL_SUCCEED) {
     978           0 :                 GDKfree(geom1MBR);
     979           0 :                 return ret;
     980             :         }
     981             : 
     982         125 :         ret = mbrDistance(out, &geom1MBR, &geom2MBR);
     983             : 
     984         125 :         GDKfree(geom1MBR);
     985         125 :         GDKfree(geom2MBR);
     986             : 
     987         125 :         return ret;
     988             : }
     989             : 
     990             : /* get Xmin, Ymin, Xmax, Ymax coordinates of mbr */
     991             : str
     992         212 : wkbCoordinateFromMBR(dbl *coordinateValue, mbr **geomMBR, int *coordinateIdx)
     993             : {
     994             :         //check if the MBR is null
     995         212 :         if (is_mbr_nil(*geomMBR) || is_int_nil(*coordinateIdx)) {
     996           8 :                 *coordinateValue = dbl_nil;
     997           8 :                 return MAL_SUCCEED;
     998             :         }
     999             : 
    1000         208 :         switch (*coordinateIdx) {
    1001          52 :         case 1:
    1002          52 :                 *coordinateValue = (*geomMBR)->xmin;
    1003          52 :                 break;
    1004          52 :         case 2:
    1005          52 :                 *coordinateValue = (*geomMBR)->ymin;
    1006          52 :                 break;
    1007          52 :         case 3:
    1008          52 :                 *coordinateValue = (*geomMBR)->xmax;
    1009          52 :                 break;
    1010          52 :         case 4:
    1011          52 :                 *coordinateValue = (*geomMBR)->ymax;
    1012          52 :                 break;
    1013           0 :         default:
    1014           0 :                 throw(MAL, "geom.coordinateFromMBR", SQLSTATE(38000) "Geos unrecognized coordinateIdx: %d\n", *coordinateIdx);
    1015             :         }
    1016             : 
    1017             :         return MAL_SUCCEED;
    1018             : }
    1019             : 
    1020             : str
    1021           0 : wkbCoordinateFromWKB(dbl *coordinateValue, wkb **geomWKB, int *coordinateIdx)
    1022             : {
    1023           0 :         mbr *geomMBR;
    1024           0 :         str ret = MAL_SUCCEED;
    1025           0 :         bit empty;
    1026             : 
    1027           0 :         if (is_wkb_nil(*geomWKB) || is_int_nil(*coordinateIdx)) {
    1028           0 :                 *coordinateValue = dbl_nil;
    1029           0 :                 return MAL_SUCCEED;
    1030             :         }
    1031             : 
    1032             :         //check if the geometry is empty
    1033           0 :         if ((ret = wkbIsEmpty(&empty, geomWKB)) != MAL_SUCCEED) {
    1034             :                 return ret;
    1035             :         }
    1036             : 
    1037           0 :         if (empty) {
    1038           0 :                 *coordinateValue = dbl_nil;
    1039           0 :                 return MAL_SUCCEED;
    1040             :         }
    1041             : 
    1042           0 :         if ((ret = wkbMBR(&geomMBR, geomWKB)) != MAL_SUCCEED)
    1043             :                 return ret;
    1044             : 
    1045           0 :         ret = wkbCoordinateFromMBR(coordinateValue, &geomMBR, coordinateIdx);
    1046             : 
    1047           0 :         GDKfree(geomMBR);
    1048             : 
    1049           0 :         return ret;
    1050             : }
    1051             : 
    1052             : str
    1053           2 : mbrFromString(mbr **w, const char **src)
    1054             : {
    1055           2 :         size_t len = *w ? sizeof(mbr) : 0;
    1056           2 :         char *errbuf;
    1057           2 :         str ex;
    1058             : 
    1059           2 :         if (mbrFROMSTR(*src, &len, (void **) w, false) >= 0)
    1060             :                 return MAL_SUCCEED;
    1061           0 :         GDKfree(*w);
    1062           0 :         *w = NULL;
    1063           0 :         errbuf = GDKerrbuf;
    1064           0 :         if (errbuf) {
    1065           0 :                 if (strncmp(errbuf, "!ERROR: ", 8) == 0)
    1066           0 :                         errbuf += 8;
    1067             :         } else {
    1068             :                 errbuf = "cannot parse string";
    1069             :         }
    1070             : 
    1071           0 :         ex = createException(MAL, "mbr.FromString", SQLSTATE(38000) "Geos %s", errbuf);
    1072             : 
    1073           0 :         GDKclrerr();
    1074             : 
    1075           0 :         return ex;
    1076             : }
    1077             : 
    1078             : /* COMMAND mbr
    1079             :  * Creates the mbr for the given geom_geometry.
    1080             :  */
    1081             : 
    1082             : str
    1083           0 : ordinatesMBR(mbr **res, flt *minX, flt *minY, flt *maxX, flt *maxY)
    1084             : {
    1085           0 :         if ((*res = GDKmalloc(sizeof(mbr))) == NULL)
    1086           0 :                 throw(MAL, "geom.mbr", SQLSTATE(HY013) MAL_MALLOC_FAIL);
    1087           0 :         if (is_flt_nil(*minX) || is_flt_nil(*minY) || is_flt_nil(*maxX) || is_flt_nil(*maxY))
    1088           0 :                 **res = mbrNIL;
    1089             :         else {
    1090           0 :                 (*res)->xmin = *minX;
    1091           0 :                 (*res)->ymin = *minY;
    1092           0 :                 (*res)->xmax = *maxX;
    1093           0 :                 (*res)->ymax = *maxY;
    1094             :         }
    1095             :         return MAL_SUCCEED;
    1096             : }
    1097             : 
    1098             : str
    1099           0 : mbrIntersects(bit *out, mbr** mbr1, mbr** mbr2) {
    1100           0 :         if (((*mbr1)->ymax < (*mbr2)->ymin) || ((*mbr1)->ymin > (*mbr2)->ymax))
    1101           0 :                 (*out) = false;
    1102           0 :         else if (((*mbr1)->xmax < (*mbr2)->xmin) || ((*mbr1)->xmin > (*mbr2)->xmax))
    1103           0 :         (*out) = false;
    1104             :         else
    1105           0 :                 (*out) = true;
    1106           0 :         return MAL_SUCCEED;
    1107             : }
    1108             : 
    1109             : /************************************************/
    1110             : /************* wkba type functions **************/
    1111             : /************************************************/
    1112             : 
    1113             : /* Creates the string representation of a wkb_array */
    1114             : /* return length of resulting string. */
    1115             : ssize_t
    1116           0 : wkbaTOSTR(char **toStr, size_t *len, const void *FROMARRAY, bool external)
    1117             : {
    1118           0 :         const wkba *fromArray = FROMARRAY;
    1119           0 :         int items = fromArray->itemsNum, i;
    1120           0 :         int itemsNumDigits = (int) ceil(log10(items));
    1121           0 :         size_t dataSize;        //, skipBytes=0;
    1122           0 :         char **partialStrs;
    1123           0 :         char *toStrPtr = NULL, *itemsNumStr = GDKmalloc(itemsNumDigits + 1);
    1124             : 
    1125           0 :         if (itemsNumStr == NULL)
    1126             :                 return -1;
    1127             : 
    1128           0 :         dataSize = (size_t) snprintf(itemsNumStr, itemsNumDigits + 1, "%d", items);
    1129             : 
    1130             :         // reserve space for an array with pointers to the partial
    1131             :         // strings, i.e. for each wkbTOSTR
    1132           0 :         partialStrs = GDKzalloc(items * sizeof(char *));
    1133           0 :         if (partialStrs == NULL) {
    1134           0 :                 GDKfree(itemsNumStr);
    1135           0 :                 return -1;
    1136             :         }
    1137             :         //create the string version of each wkb
    1138           0 :         for (i = 0; i < items; i++) {
    1139           0 :                 size_t llen = 0;
    1140           0 :                 ssize_t ds;
    1141           0 :                 ds = wkbTOSTR(&partialStrs[i], &llen, fromArray->data[i], false);
    1142           0 :                 if (ds < 0) {
    1143           0 :                         GDKfree(itemsNumStr);
    1144           0 :                         while (i >= 0)
    1145           0 :                                 GDKfree(partialStrs[i--]);
    1146           0 :                         GDKfree(partialStrs);
    1147           0 :                         return -1;
    1148             :                 }
    1149           0 :                 dataSize += ds;
    1150             : 
    1151           0 :                 if (strNil(partialStrs[i])) {
    1152           0 :                         GDKfree(itemsNumStr);
    1153           0 :                         while (i >= 0)
    1154           0 :                                 GDKfree(partialStrs[i--]);
    1155           0 :                         GDKfree(partialStrs);
    1156           0 :                         if (*len < 4 || *toStr == NULL) {
    1157           0 :                                 GDKfree(*toStr);
    1158           0 :                                 if ((*toStr = GDKmalloc(*len = 4)) == NULL)
    1159             :                                         return -1;
    1160             :                         }
    1161           0 :                         if (external) {
    1162           0 :                                 strcpy(*toStr, "nil");
    1163           0 :                                 return 3;
    1164             :                         }
    1165           0 :                         strcpy(*toStr, str_nil);
    1166           0 :                         return 1;
    1167             :                 }
    1168             :         }
    1169             : 
    1170             :         //add [] around itemsNum
    1171           0 :         dataSize += 2;
    1172             :         //add ", " before each item
    1173           0 :         dataSize += 2 * sizeof(char) * items;
    1174             : 
    1175             :         //copy all partial strings to a single one
    1176           0 :         if (*len < dataSize + 3 || *toStr == NULL) {
    1177           0 :                 GDKfree(*toStr);
    1178           0 :                 *toStr = GDKmalloc(*len = dataSize + 3);        /* plus quotes + termination character */
    1179           0 :                 if (*toStr == NULL) {
    1180           0 :                         for (i = 0; i < items; i++)
    1181           0 :                                 GDKfree(partialStrs[i]);
    1182           0 :                         GDKfree(partialStrs);
    1183           0 :                         GDKfree(itemsNumStr);
    1184           0 :                         return -1;
    1185             :                 }
    1186             :         }
    1187           0 :         toStrPtr = *toStr;
    1188           0 :         if (external)
    1189           0 :                 *toStrPtr++ = '\"';
    1190           0 :         *toStrPtr++ = '[';
    1191           0 :         strcpy(toStrPtr, itemsNumStr);
    1192           0 :         toStrPtr += strlen(itemsNumStr);
    1193           0 :         *toStrPtr++ = ']';
    1194           0 :         for (i = 0; i < items; i++) {
    1195           0 :                 if (i == 0)
    1196           0 :                         *toStrPtr++ = ':';
    1197             :                 else
    1198           0 :                         *toStrPtr++ = ',';
    1199           0 :                 *toStrPtr++ = ' ';
    1200             : 
    1201             :                 //strcpy(toStrPtr, partialStrs[i]);
    1202           0 :                 memcpy(toStrPtr, partialStrs[i], strlen(partialStrs[i]));
    1203           0 :                 toStrPtr += strlen(partialStrs[i]);
    1204           0 :                 GDKfree(partialStrs[i]);
    1205             :         }
    1206             : 
    1207           0 :         if (external)
    1208           0 :                 *toStrPtr++ = '\"';
    1209           0 :         *toStrPtr = '\0';
    1210             : 
    1211           0 :         GDKfree(partialStrs);
    1212           0 :         GDKfree(itemsNumStr);
    1213             : 
    1214           0 :         return (ssize_t) (toStrPtr - *toStr);
    1215             : }
    1216             : 
    1217             : static ssize_t wkbaFROMSTR_withSRID(const char *fromStr, size_t *len, wkba **toArray, int srid);
    1218             : 
    1219             : /* return number of parsed characters. */
    1220             : ssize_t
    1221           0 : wkbaFROMSTR(const char *fromStr, size_t *len, void **TOARRAY, bool external)
    1222             : {
    1223           0 :         wkba **toArray = (wkba **) TOARRAY;
    1224           0 :         if (external && strncmp(fromStr, "nil", 3) == 0) {
    1225           0 :                 size_t sz = wkba_size(~0);
    1226           0 :                 if ((*len < sz || *toArray == NULL)
    1227           0 :                     && (*toArray = GDKmalloc(sz)) == NULL)
    1228             :                         return -1;
    1229           0 :                 **toArray = wkba_nil;
    1230           0 :                 return 3;
    1231             :         }
    1232           0 :         return wkbaFROMSTR_withSRID(fromStr, len, toArray, 0);
    1233             : }
    1234             : 
    1235             : /* returns a pointer to a null wkba */
    1236             : const void *
    1237         315 : wkbaNULL(void)
    1238             : {
    1239         315 :         return &wkba_nil;
    1240             : }
    1241             : 
    1242             : BUN
    1243           0 : wkbaHASH(const void *WARRAY)
    1244             : {
    1245           0 :         const wkba *wArray = WARRAY;
    1246           0 :         int j, i;
    1247           0 :         BUN h = 0;
    1248             : 
    1249           0 :         for (j = 0; j < wArray->itemsNum; j++) {
    1250           0 :                 wkb *w = wArray->data[j];
    1251           0 :                 for (i = 0; i < (w->len - 1); i += 2) {
    1252           0 :                         BUN a = ((unsigned char *) w->data)[i];
    1253           0 :                         BUN b = ((unsigned char *) w->data)[i + 1];
    1254             : #if '\377' < 0                                       /* char is signed? */
    1255             :                         /* maybe sign extend */
    1256           0 :                         if (a & 0x80)
    1257           0 :                                 a |= ~(BUN)0x7f;
    1258           0 :                         if (b & 0x80)
    1259           0 :                                 b |= ~(BUN)0x7f;
    1260             : #endif
    1261           0 :                         h = (h << 3) ^ (h >> 11) ^ (h >> 17) ^ (b << 8) ^ a;
    1262             :                 }
    1263             :         }
    1264           0 :         return h;
    1265             : }
    1266             : 
    1267             : int
    1268           0 : wkbaCOMP(const void *L, const void *R)
    1269             : {
    1270           0 :         const wkba *l = L, *r = R;
    1271           0 :         int i, res = 0;
    1272             : 
    1273             :         //compare the number of items
    1274           0 :         if (l->itemsNum != r->itemsNum)
    1275           0 :                 return l->itemsNum - r->itemsNum;
    1276             : 
    1277           0 :         if (l->itemsNum == ~(int) 0)
    1278             :                 return (0);
    1279             : 
    1280             :         //compare each wkb separately
    1281           0 :         for (i = 0; i < l->itemsNum; i++)
    1282           0 :                 res += wkbCOMP(l->data[i], r->data[i]);
    1283             : 
    1284             :         return res;
    1285             : }
    1286             : 
    1287             : /* read wkb from log */
    1288             : void *
    1289           0 : wkbaREAD(void *A, size_t *dstlen, stream *s, size_t cnt)
    1290             : {
    1291           0 :         wkba *a = A;
    1292           0 :         int items, i;
    1293             : 
    1294           0 :         (void) cnt;
    1295           0 :         assert(cnt == 1);
    1296             : 
    1297           0 :         if (mnstr_readInt(s, &items) != 1)
    1298             :                 return NULL;
    1299             : 
    1300           0 :         size_t wkbalen = (size_t) wkba_size(items);
    1301           0 :         if (a == NULL || *dstlen < wkbalen) {
    1302           0 :                 if ((a = GDKrealloc(a, wkbalen)) == NULL)
    1303             :                         return NULL;
    1304           0 :                 *dstlen = wkbalen;
    1305             :         }
    1306             : 
    1307           0 :         a->itemsNum = items;
    1308             : 
    1309           0 :         for (i = 0; i < items; i++) {
    1310           0 :                 size_t wlen = 0;
    1311           0 :                 a->data[i] = wkbREAD(NULL, &wlen, s, cnt);
    1312             :         }
    1313             : 
    1314             :         return a;
    1315             : }
    1316             : 
    1317             : /* write wkb to log */
    1318             : gdk_return
    1319           0 : wkbaWRITE(const void *A, stream *s, size_t cnt)
    1320             : {
    1321           0 :         const wkba *a = A;
    1322           0 :         int i, items = a->itemsNum;
    1323           0 :         gdk_return ret = GDK_SUCCEED;
    1324             : 
    1325           0 :         (void) cnt;
    1326           0 :         assert(cnt == 1);
    1327             : 
    1328           0 :         if (!mnstr_writeInt(s, items))
    1329             :                 return GDK_FAIL;
    1330           0 :         for (i = 0; i < items; i++) {
    1331           0 :                 ret = wkbWRITE(a->data[i], s, cnt);
    1332             : 
    1333           0 :                 if (ret != GDK_SUCCEED)
    1334           0 :                         return ret;
    1335             :         }
    1336             :         return GDK_SUCCEED;
    1337             : }
    1338             : 
    1339             : var_t
    1340           0 : wkbaPUT(BAT *b, var_t *bun, const void *VAL)
    1341             : {
    1342           0 :         const wkba *val = VAL;
    1343           0 :         char *base;
    1344             : 
    1345           0 :         *bun = HEAP_malloc(b, wkba_size(val->itemsNum));
    1346           0 :         base = b->tvheap->base;
    1347           0 :         if (*bun != (var_t) -1) {
    1348           0 :                 memcpy(&base[*bun], val, wkba_size(val->itemsNum));
    1349           0 :                 b->tvheap->dirty = true;
    1350             :         }
    1351           0 :         return *bun;
    1352             : }
    1353             : 
    1354             : void
    1355           0 : wkbaDEL(Heap *h, var_t *index)
    1356             : {
    1357           0 :         HEAP_free(h, *index);
    1358           0 : }
    1359             : 
    1360             : size_t
    1361           0 : wkbaLENGTH(const void *P)
    1362             : {
    1363           0 :         const wkba *p = P;
    1364           0 :         var_t len = wkba_size(p->itemsNum);
    1365           0 :         assert(len <= GDK_int_max);
    1366           0 :         return (size_t) len;
    1367             : }
    1368             : 
    1369             : gdk_return
    1370         315 : wkbaHEAP(Heap *heap, size_t capacity)
    1371             : {
    1372         315 :         return HEAP_initialize(heap, capacity, 0, (int) sizeof(var_t));
    1373             : }
    1374             : 
    1375             : /* Non-atom WKBA functions */
    1376             : /* returns the size of variable-sized atom wkba */
    1377             : var_t
    1378           0 : wkba_size(int items)
    1379             : {
    1380           0 :         var_t size;
    1381             : 
    1382           0 :         if (items == ~0)
    1383           0 :                 items = 0;
    1384           0 :         size = (var_t) (offsetof(wkba, data) + items * sizeof(wkb *));
    1385           0 :         assert(size <= VAR_MAX);
    1386             : 
    1387           0 :         return size;
    1388             : }
    1389             : 
    1390             : static ssize_t
    1391           0 : wkbaFROMSTR_withSRID(const char *fromStr, size_t *len, wkba **toArray, int srid)
    1392             : {
    1393           0 :         int items, i;
    1394           0 :         size_t skipBytes = 0;
    1395             : 
    1396             : //IS THERE SPACE OR SOME OTHER CHARACTER?
    1397             : 
    1398             :         //read the number of items from the beginning of the string
    1399           0 :         memcpy(&items, fromStr, sizeof(int));
    1400           0 :         skipBytes += sizeof(int);
    1401           0 :         *toArray = GDKmalloc(wkba_size(items));
    1402           0 :         if (*toArray == NULL)
    1403             :                 return -1;
    1404             : 
    1405           0 :         for (i = 0; i < items; i++) {
    1406           0 :                 size_t parsedBytes;
    1407           0 :                 str err = wkbFROMSTR_withSRID(fromStr + skipBytes, len, &(*toArray)->data[i], srid, &parsedBytes);
    1408           0 :                 if (err != MAL_SUCCEED) {
    1409           0 :                         GDKerror("%s", getExceptionMessageAndState(err));
    1410           0 :                         freeException(err);
    1411           0 :                         return -1;
    1412             :                 }
    1413           0 :                 skipBytes += parsedBytes;
    1414             :         }
    1415             : 
    1416           0 :         assert(skipBytes <= GDK_int_max);
    1417           0 :         return (ssize_t) skipBytes;
    1418             : }
    1419             : 
    1420             : /* Only use of WKBA in exported functions */
    1421             : str
    1422           0 : wkbInteriorRings(wkba **geomArray, wkb **geomWKB)
    1423             : {
    1424           0 :         int interiorRingsNum = 0, i = 0;
    1425           0 :         GEOSGeom geosGeometry;
    1426           0 :         str ret = MAL_SUCCEED;
    1427             : 
    1428           0 :         if (is_wkb_nil(*geomWKB)) {
    1429           0 :                 if ((*geomArray = GDKmalloc(wkba_size(~0))) == NULL)
    1430           0 :                         throw(MAL, "geom.InteriorRings", SQLSTATE(HY013) MAL_MALLOC_FAIL);
    1431           0 :                 **geomArray = wkba_nil;
    1432           0 :                 return MAL_SUCCEED;
    1433             :         }
    1434             : 
    1435           0 :         geosGeometry = wkb2geos(*geomWKB);
    1436           0 :         if (geosGeometry == NULL) {
    1437           0 :                 throw(MAL, "geom.InteriorRings", SQLSTATE(38000) "Geos operation  wkb2geos failed");
    1438             :         }
    1439             : 
    1440           0 :         if ((GEOSGeomTypeId_r(geoshandle, geosGeometry) + 1) != wkbPolygon_mdb) {
    1441           0 :                 GEOSGeom_destroy_r(geoshandle, geosGeometry);
    1442           0 :                 throw(MAL, "geom.interiorRings", SQLSTATE(38000) "Geometry not a Polygon");
    1443             : 
    1444             :         }
    1445             : 
    1446           0 :         ret = wkbNumRings(&interiorRingsNum, geomWKB, &i);
    1447             : 
    1448           0 :         if (ret != MAL_SUCCEED) {
    1449           0 :                 GEOSGeom_destroy_r(geoshandle, geosGeometry);
    1450           0 :                 return ret;
    1451             :         }
    1452             : 
    1453           0 :         *geomArray = GDKmalloc(wkba_size(interiorRingsNum));
    1454           0 :         if (*geomArray == NULL) {
    1455           0 :                 GEOSGeom_destroy_r(geoshandle, geosGeometry);
    1456           0 :                 throw(MAL, "geom.InteriorRings", SQLSTATE(HY013) MAL_MALLOC_FAIL);
    1457             :         }
    1458           0 :         (*geomArray)->itemsNum = interiorRingsNum;
    1459             : 
    1460           0 :         for (i = 0; i < interiorRingsNum; i++) {
    1461           0 :                 const GEOSGeometry *interiorRingGeometry;
    1462           0 :                 wkb *interiorRingWKB;
    1463             : 
    1464             :                 // get the interior ring of the geometry
    1465           0 :                 interiorRingGeometry = GEOSGetInteriorRingN_r(geoshandle, geosGeometry, i);
    1466           0 :                 if (interiorRingGeometry == NULL) {
    1467           0 :                         while (--i >= 0)
    1468           0 :                                 GDKfree((*geomArray)->data[i]);
    1469           0 :                         GDKfree(*geomArray);
    1470           0 :                         GEOSGeom_destroy_r(geoshandle, geosGeometry);
    1471           0 :                         *geomArray = NULL;
    1472           0 :                         throw(MAL, "geom.InteriorRings", SQLSTATE(38000) "Geos operation GEOSGetInteriorRingN failed");
    1473             :                 }
    1474             :                 // get the wkb representation of it
    1475           0 :                 interiorRingWKB = geos2wkb(interiorRingGeometry);
    1476           0 :                 if (interiorRingWKB == NULL) {
    1477           0 :                         while (--i >= 0)
    1478           0 :                                 GDKfree((*geomArray)->data[i]);
    1479           0 :                         GDKfree(*geomArray);
    1480           0 :                         GEOSGeom_destroy_r(geoshandle, geosGeometry);
    1481           0 :                         *geomArray = NULL;
    1482           0 :                         throw(MAL, "geom.InteriorRings", SQLSTATE(38000) "Geos operation wkb2geos failed");
    1483             :                 }
    1484             : 
    1485           0 :                 (*geomArray)->data[i] = interiorRingWKB;
    1486             :         }
    1487           0 :         GEOSGeom_destroy_r(geoshandle, geosGeometry);
    1488             : 
    1489           0 :         return MAL_SUCCEED;
    1490             : }

Generated by: LCOV version 1.14