Ticket #210: postgis-nullunion.patch

File postgis-nullunion.patch, 10.6 KB (added by mcayland, 3 years ago)

  • postgis/lwgeom_geos.c

     
    107107        GEOSGeometry * geos_result=NULL; 
    108108        int SRID=-1; 
    109109        size_t offset = 0; 
     110        bits8 *bitmap; 
     111        int bitmask; 
    110112#if POSTGIS_DEBUG_LEVEL > 0 
    111113        static int call=1; 
    112114#endif 
    … …  
    128130 
    129131        nelems = ArrayGetNItems(ARR_NDIM(array), ARR_DIMS(array)); 
    130132 
     133        bitmap = ARR_NULLBITMAP(array); 
     134 
    131135        POSTGIS_DEBUGF(3, "unite_garray: number of elements: %d", nelems); 
    132136 
    133137        if ( nelems == 0 ) PG_RETURN_NULL(); 
    134138 
    135139        /* One-element union is the element itself */ 
    136         if ( nelems == 1 ) PG_RETURN_POINTER((PG_LWGEOM *)(ARR_DATA_PTR(array))); 
     140        if ( nelems == 1 ) 
     141        { 
     142                /* If the element is a NULL then we need to handle it separately */ 
     143                if (bitmap && (*bitmap & 1) == 0) 
     144                        PG_RETURN_NULL(); 
     145                else 
     146                        PG_RETURN_POINTER((PG_LWGEOM *)(ARR_DATA_PTR(array))); 
     147        } 
    137148 
    138149        /* Ok, we really need geos now ;) */ 
    139150        initGEOS(lwnotice, lwnotice); 
    … …  
    145156        ** If they are, we can use UnionCascaded for faster results. 
    146157        */ 
    147158        offset = 0; 
     159        bitmask = 1; 
    148160        for ( i = 0; i < nelems; i++ ) 
    149161        { 
    150                 PG_LWGEOM *pggeom = (PG_LWGEOM *)(ARR_DATA_PTR(array)+offset); 
    151                 int pgtype = TYPE_GETTYPE(pggeom->type); 
    152                 offset += INTALIGN(VARSIZE(pggeom)); 
    153                 if ( ! i ) /* Initialize SRID */ 
     162                /* Don't do anything for NULL values */ 
     163                if ((bitmap && (*bitmap & bitmask) != 0) || !bitmap) 
    154164                { 
    155                         SRID = pglwgeom_getSRID(pggeom); 
    156                         if ( TYPE_HASZ(pggeom->type) ) is3d = 1; 
     165                        PG_LWGEOM *pggeom = (PG_LWGEOM *)(ARR_DATA_PTR(array)+offset); 
     166                        int pgtype = TYPE_GETTYPE(pggeom->type); 
     167                        offset += INTALIGN(VARSIZE(pggeom)); 
     168                        if ( ! i ) /* Initialize SRID */ 
     169                        { 
     170                                SRID = pglwgeom_getSRID(pggeom); 
     171                                if ( TYPE_HASZ(pggeom->type) ) is3d = 1; 
     172                        } 
     173                        if ( pgtype != POLYGONTYPE && pgtype != MULTIPOLYGONTYPE ) 
     174                        { 
     175                                allpolys = 0; 
     176                                break; 
     177                        } 
    157178                } 
    158                 if ( pgtype != POLYGONTYPE && pgtype != MULTIPOLYGONTYPE ) 
    159                 { 
    160                         allpolys = 0; 
    161                         break; 
     179 
     180                /* Advance NULL bitmap */ 
     181                if (bitmap) 
     182                { 
     183                        bitmask <<= 1; 
     184                        if (bitmask == 0x100) 
     185                        { 
     186                                bitmap++; 
     187                                bitmask = 1; 
     188                        } 
    162189                } 
    163190        } 
    164191 
    … …  
    176203                ** First make an array of GEOS Polygons. 
    177204                */ 
    178205                offset = 0; 
     206                bitmap = ARR_NULLBITMAP(array); 
     207                bitmask = 1; 
    179208                for ( i = 0; i < nelems; i++ ) 
    180209                { 
    181                         PG_LWGEOM *pggeom = (PG_LWGEOM *)(ARR_DATA_PTR(array)+offset); 
    182                         int pgtype = TYPE_GETTYPE(pggeom->type); 
    183                         offset += INTALIGN(VARSIZE(pggeom)); 
    184                         if ( pgtype == POLYGONTYPE ) 
     210                        /* Don't do anything for NULL values */ 
     211                        if ((bitmap && (*bitmap & bitmask) != 0) || !bitmap) 
    185212                        { 
    186                                 if ( curgeom == geoms_size ) 
     213                                PG_LWGEOM *pggeom = (PG_LWGEOM *)(ARR_DATA_PTR(array)+offset); 
     214                                int pgtype = TYPE_GETTYPE(pggeom->type); 
     215                                offset += INTALIGN(VARSIZE(pggeom)); 
     216                                if ( pgtype == POLYGONTYPE ) 
    187217                                { 
    188                                         geoms_size *= 2; 
    189                                         geoms = repalloc( geoms, sizeof(GEOSGeom) * geoms_size ); 
    190                                 } 
    191                                 geoms[curgeom] = (GEOSGeometry *)POSTGIS2GEOS(pggeom); 
    192                                 curgeom++; 
    193                         } 
    194                         if ( pgtype == MULTIPOLYGONTYPE ) 
    195                         { 
    196                                 int j = 0; 
    197                                 LWGEOM_INSPECTED *lwgeom = lwgeom_inspect(SERIALIZED_FORM(pggeom)); 
    198                                 for ( j = 0; j < lwgeom->ngeometries; j++ ) 
    199                                 { 
    200                                         LWPOLY *lwpoly = NULL; 
    201                                         int k = 0; 
    202218                                        if ( curgeom == geoms_size ) 
    203219                                        { 
    204220                                                geoms_size *= 2; 
    205221                                                geoms = repalloc( geoms, sizeof(GEOSGeom) * geoms_size ); 
    206222                                        } 
    207                                         /* This builds a LWPOLY on top of the serialized form */ 
    208                                         lwpoly = lwgeom_getpoly_inspected(lwgeom, j); 
    209                                         geoms[curgeom] = LWGEOM2GEOS(lwpoly_as_lwgeom(lwpoly)); 
    210                                         /* We delicately free the LWPOLY and POINTARRAY structs, leaving the serialized form below untouched. */ 
    211                                         for ( k = 0; k < lwpoly->nrings; k++ ) 
     223                                        geoms[curgeom] = (GEOSGeometry *)POSTGIS2GEOS(pggeom); 
     224                                        curgeom++; 
     225                                } 
     226                                if ( pgtype == MULTIPOLYGONTYPE ) 
     227                                { 
     228                                        int j = 0; 
     229                                        LWGEOM_INSPECTED *lwgeom = lwgeom_inspect(SERIALIZED_FORM(pggeom)); 
     230                                        for ( j = 0; j < lwgeom->ngeometries; j++ ) 
    212231                                        { 
    213                                                 lwfree(lwpoly->rings[k]); 
     232                                                LWPOLY *lwpoly = NULL; 
     233                                                int k = 0; 
     234                                                if ( curgeom == geoms_size ) 
     235                                                { 
     236                                                        geoms_size *= 2; 
     237                                                        geoms = repalloc( geoms, sizeof(GEOSGeom) * geoms_size ); 
     238                                                } 
     239                                                /* This builds a LWPOLY on top of the serialized form */ 
     240                                                lwpoly = lwgeom_getpoly_inspected(lwgeom, j); 
     241                                                geoms[curgeom] = LWGEOM2GEOS(lwpoly_as_lwgeom(lwpoly)); 
     242                                                /* We delicately free the LWPOLY and POINTARRAY structs, leaving the serialized form below untouched. */ 
     243                                                for ( k = 0; k < lwpoly->nrings; k++ ) 
     244                                                { 
     245                                                        lwfree(lwpoly->rings[k]); 
     246                                                } 
     247                                                lwpoly_release(lwpoly); 
     248                                                curgeom++; 
    214249                                        } 
    215                                         lwpoly_release(lwpoly); 
    216                                         curgeom++; 
    217250                                } 
    218251                        } 
    219252 
     253                        /* Advance NULL bitmap */ 
     254                        if (bitmap) 
     255                        { 
     256                                bitmask <<= 1; 
     257                                if (bitmask == 0x100) 
     258                                { 
     259                                        bitmap++; 
     260                                        bitmask = 1; 
     261                                } 
     262                        } 
    220263                } 
    221264                /* 
    222265                ** Take our GEOS Polygons and turn them into a GEOS MultiPolygon, 
    223266                ** then pass that into cascaded union. 
    224267                */ 
    225                 g1 = GEOSGeom_createCollection(GEOS_MULTIPOLYGON, geoms, curgeom); 
    226                 if ( g1 ) g2 = GEOSUnionCascaded(g1); 
    227                 if ( g2 ) GEOSSetSRID(g2, SRID); 
    228                 if ( g2 ) result = GEOS2POSTGIS(g2, is3d); 
    229                 /* Clean up the mess. */ 
    230                 if ( g1 ) GEOSGeom_destroy((GEOSGeometry *)g1); 
    231                 if ( g2 ) GEOSGeom_destroy(g2); 
     268                if (curgeom > 0) 
     269                { 
     270                        g1 = GEOSGeom_createCollection(GEOS_MULTIPOLYGON, geoms, curgeom); 
     271                        if ( g1 ) g2 = GEOSUnionCascaded(g1); 
     272                        if ( g2 ) GEOSSetSRID(g2, SRID); 
     273                        if ( g2 ) result = GEOS2POSTGIS(g2, is3d); 
     274                        /* Clean up the mess. */ 
     275                        if ( g1 ) GEOSGeom_destroy((GEOSGeometry *)g1); 
     276                        if ( g2 ) GEOSGeom_destroy(g2); 
     277                } 
     278                else 
     279                { 
     280                        /* All we found were NULLs, so let's return NULL */ 
     281                        PG_RETURN_NULL(); 
     282                } 
    232283        } 
    233284        else 
    234285        { 
    … …  
    237288                ** Heterogeneous result, let's slog through this one union at a time. 
    238289                */ 
    239290                offset = 0; 
     291                bitmap = ARR_NULLBITMAP(array); 
     292                bitmask = 1; 
    240293                for (i=0; i<nelems; i++) 
    241294                { 
    242                         PG_LWGEOM *geom = (PG_LWGEOM *)(ARR_DATA_PTR(array)+offset); 
    243                         offset += INTALIGN(VARSIZE(geom)); 
    244  
    245                         pgis_geom = geom; 
    246  
    247                         POSTGIS_DEBUGF(3, "geom %d @ %p", i, geom); 
    248  
    249                         /* Check is3d flag */ 
    250                         if ( TYPE_HASZ(geom->type) ) is3d = 1; 
    251  
    252                         /* Check SRID homogeneity and initialize geos result */ 
    253                         if ( ! i ) 
     295                        /* Don't do anything for NULL values */ 
     296                        if ((bitmap && (*bitmap & bitmask) != 0) || !bitmap) 
    254297                        { 
    255                                 geos_result = (GEOSGeometry *)POSTGIS2GEOS(geom); 
    256                                 SRID = pglwgeom_getSRID(geom); 
    257                                 POSTGIS_DEBUGF(3, "first geom is a %s", lwgeom_typename(TYPE_GETTYPE(geom->type))); 
    258                                 continue; 
     298                                PG_LWGEOM *geom = (PG_LWGEOM *)(ARR_DATA_PTR(array)+offset); 
     299                                offset += INTALIGN(VARSIZE(geom)); 
     300         
     301                                pgis_geom = geom; 
     302         
     303                                POSTGIS_DEBUGF(3, "geom %d @ %p", i, geom); 
     304         
     305                                /* Check is3d flag */ 
     306                                if ( TYPE_HASZ(geom->type) ) is3d = 1; 
     307         
     308                                /* Check SRID homogeneity and initialize geos result */ 
     309                                if ( ! i ) 
     310                                { 
     311                                        geos_result = (GEOSGeometry *)POSTGIS2GEOS(geom); 
     312                                        SRID = pglwgeom_getSRID(geom); 
     313                                        POSTGIS_DEBUGF(3, "first geom is a %s", lwgeom_typename(TYPE_GETTYPE(geom->type))); 
     314                                        continue; 
     315                                } 
     316                                else 
     317                                { 
     318                                        errorIfSRIDMismatch(SRID, pglwgeom_getSRID(geom)); 
     319                                } 
     320         
     321                                g1 = POSTGIS2GEOS(pgis_geom); 
     322         
     323                                POSTGIS_DEBUGF(3, "unite_garray(%d): adding geom %d to union (%s)", 
     324                                        call, i, lwgeom_typename(TYPE_GETTYPE(geom->type))); 
     325         
     326                                g2 = GEOSUnion(g1, geos_result); 
     327                                if ( g2 == NULL ) 
     328                                { 
     329                                        GEOSGeom_destroy((GEOSGeometry *)g1); 
     330                                        GEOSGeom_destroy((GEOSGeometry *)geos_result); 
     331                                        elog(ERROR,"GEOS union() threw an error!"); 
     332                                } 
     333                                GEOSGeom_destroy((GEOSGeometry *)g1); 
     334                                GEOSGeom_destroy((GEOSGeometry *)geos_result); 
     335                                geos_result = g2; 
    259336                        } 
    260                         else 
    261                         { 
    262                                 errorIfSRIDMismatch(SRID, pglwgeom_getSRID(geom)); 
    263                         } 
    264337 
    265                         g1 = POSTGIS2GEOS(pgis_geom); 
    266  
    267                         POSTGIS_DEBUGF(3, "unite_garray(%d): adding geom %d to union (%s)", 
    268                                        call, i, lwgeom_typename(TYPE_GETTYPE(geom->type))); 
    269  
    270                         g2 = GEOSUnion(g1, geos_result); 
    271                         if ( g2 == NULL ) 
     338                        /* Advance NULL bitmap */ 
     339                        if (bitmap) 
    272340                        { 
    273                                 GEOSGeom_destroy((GEOSGeometry *)g1); 
    274                                 GEOSGeom_destroy((GEOSGeometry *)geos_result); 
    275                                 elog(ERROR,"GEOS union() threw an error!"); 
     341                                bitmask <<= 1; 
     342                                if (bitmask == 0x100) 
     343                                { 
     344                                        bitmap++; 
     345                                        bitmask = 1; 
     346                                } 
    276347                        } 
    277                         GEOSGeom_destroy((GEOSGeometry *)g1); 
    278                         GEOSGeom_destroy((GEOSGeometry *)geos_result); 
    279                         geos_result = g2; 
     348                         
    280349                } 
    281350 
    282                 GEOSSetSRID(geos_result, SRID); 
    283                 result = GEOS2POSTGIS(geos_result, is3d); 
    284                 GEOSGeom_destroy(geos_result); 
     351                /* If geos_result is set then we found at least one non-NULL geometry */ 
     352                if (geos_result) 
     353                { 
     354                        GEOSSetSRID(geos_result, SRID); 
     355                        result = GEOS2POSTGIS(geos_result, is3d); 
     356                        GEOSGeom_destroy(geos_result); 
     357                } 
     358                else 
     359                { 
     360                        /* All we found were NULLs, so let's return NULL */ 
     361                        PG_RETURN_NULL(); 
     362                } 
    285363 
    286364#if POSTGIS_GEOS_VERSION >= 31 
    287365        } 
    … …  
    289367 
    290368        if ( result == NULL ) 
    291369        { 
    292                 elog(ERROR, "Union returned a NULL geometry."); 
    293                 PG_RETURN_NULL(); /* never get here */ 
     370                /* Union returned a NULL geometry */ 
     371                PG_RETURN_NULL(); 
    294372        } 
    295373 
    296374        PG_RETURN_POINTER(result); 

  • postgis/lwgeom_accum.c

     
    2121#include "lwgeom_pg.h" 
    2222 
    2323/* Local prototypes */ 
     24Datum PGISDirectFunctionCall1(PGFunction func, Datum arg1); 
    2425Datum pgis_geometry_accum_transfn(PG_FUNCTION_ARGS); 
    2526Datum pgis_geometry_accum_finalfn(PG_FUNCTION_ARGS); 
    2627Datum pgis_geometry_union_finalfn(PG_FUNCTION_ARGS); 
    … …  
    210211        p = (pgis_abs*) PG_GETARG_POINTER(0); 
    211212 
    212213        geometry_array = pgis_accum_finalfn(p, CurrentMemoryContext, fcinfo); 
    213         result = DirectFunctionCall1( pgis_union_geometry_array, geometry_array ); 
     214        result = PGISDirectFunctionCall1( pgis_union_geometry_array, geometry_array ); 
     215        if (!result) 
     216                PG_RETURN_NULL(); 
    214217 
    215218        PG_RETURN_DATUM(result); 
    216219} 
    … …  
    284287        PG_RETURN_DATUM(result); 
    285288} 
    286289 
     290/** 
     291* A modified version of PostgreSQL's DirectFunctionCall1 which allows NULL results; this 
     292* is required for aggregates that return NULL. 
     293*/ 
     294Datum 
     295PGISDirectFunctionCall1(PGFunction func, Datum arg1) 
     296{ 
     297        FunctionCallInfoData fcinfo; 
     298        Datum           result; 
     299 
     300        InitFunctionCallInfoData(fcinfo, NULL, 1, NULL, NULL); 
     301 
     302        fcinfo.arg[0] = arg1; 
     303        fcinfo.argnull[0] = false; 
     304 
     305        result = (*func) (&fcinfo); 
     306 
     307        /* Check for null result, since caller is clearly not expecting one */ 
     308        if (fcinfo.isnull) 
     309                return (Datum) 0; 
     310 
     311        return result; 
     312}