[postgis-tickets] r17808 - Speed up the calculation of cartesian bbox

Raul raul at rmr.ninja
Fri Sep 13 08:40:41 PDT 2019


Author: algunenano
Date: 2019-09-13 08:40:41 -0700 (Fri, 13 Sep 2019)
New Revision: 17808

Modified:
   trunk/NEWS
   trunk/liblwgeom/gbox.c
   trunk/liblwgeom/liblwgeom.h.in
   trunk/liblwgeom/lwgeom_api.c
   trunk/liblwgeom/lwgeom_geos.c
   trunk/liblwgeom/lwinline.h
   trunk/liblwgeom/lwout_geojson.c
   trunk/liblwgeom/lwout_gml.c
Log:
Speed up the calculation of cartesian bbox

Closes https://github.com/postgis/postgis/pull/475
Closes https://trac.osgeo.org/postgis/ticket/4503



Modified: trunk/NEWS
===================================================================
--- trunk/NEWS	2019-09-12 17:47:47 UTC (rev 17807)
+++ trunk/NEWS	2019-09-13 15:40:41 UTC (rev 17808)
@@ -11,6 +11,7 @@
   - #4495, Fix ST_SnapToGrid output having an outdated bbox (Raúl Marín)
   - #4496, Make ST_Simplify(TRIANGLE) collapse if requested (Raúl Marín)
   - #4501, Allow postgis_tiger_geocoder to be installable by non-super users (Regina Obe)
+  - #4503, Speed up the calculation of cartesian bbox (Raúl Marín)
 
 PostGIS 3.0.0alpha4
 2019/08/10

Modified: trunk/liblwgeom/gbox.c
===================================================================
--- trunk/liblwgeom/gbox.c	2019-09-12 17:47:47 UTC (rev 17807)
+++ trunk/liblwgeom/gbox.c	2019-09-13 15:40:41 UTC (rev 17808)
@@ -533,47 +533,115 @@
 	return rv;
 }
 
-int ptarray_calculate_gbox_cartesian(const POINTARRAY *pa, GBOX *gbox )
+static void
+ptarray_calculate_gbox_cartesian_2d(const POINTARRAY *pa, GBOX *gbox)
 {
-	uint32_t i;
-	POINT4D p;
-	int has_z, has_m;
+	const POINT2D *p = getPoint2d_cp(pa, 0);
 
-	if ( ! pa ) return LW_FAILURE;
-	if ( ! gbox ) return LW_FAILURE;
-	if ( pa->npoints < 1 ) return LW_FAILURE;
+	gbox->xmax = gbox->xmin = p->x;
+	gbox->ymax = gbox->ymin = p->y;
 
-	has_z = FLAGS_GET_Z(pa->flags);
-	has_m = FLAGS_GET_M(pa->flags);
+	for (uint32_t i = 1; i < pa->npoints; i++)
+	{
+		p = getPoint2d_cp(pa, i);
+		gbox->xmin = FP_MIN(gbox->xmin, p->x);
+		gbox->xmax = FP_MAX(gbox->xmax, p->x);
+		gbox->ymin = FP_MIN(gbox->ymin, p->y);
+		gbox->ymax = FP_MAX(gbox->ymax, p->y);
+	}
+}
+
+/* Works with X/Y/Z. Needs to be adjusted after if X/Y/M was required */
+static void
+ptarray_calculate_gbox_cartesian_3d(const POINTARRAY *pa, GBOX *gbox)
+{
+	const POINT3D *p = getPoint3d_cp(pa, 0);
+
+	gbox->xmax = gbox->xmin = p->x;
+	gbox->ymax = gbox->ymin = p->y;
+	gbox->zmax = gbox->zmin = p->z;
+
+	for (uint32_t i = 1; i < pa->npoints; i++)
+	{
+		p = getPoint3d_cp(pa, i);
+		gbox->xmin = FP_MIN(gbox->xmin, p->x);
+		gbox->xmax = FP_MAX(gbox->xmax, p->x);
+		gbox->ymin = FP_MIN(gbox->ymin, p->y);
+		gbox->ymax = FP_MAX(gbox->ymax, p->y);
+		gbox->zmin = FP_MIN(gbox->zmin, p->z);
+		gbox->zmax = FP_MAX(gbox->zmax, p->z);
+	}
+}
+
+static void
+ptarray_calculate_gbox_cartesian_4d(const POINTARRAY *pa, GBOX *gbox)
+{
+	const POINT4D *p = getPoint4d_cp(pa, 0);
+
+	gbox->xmax = gbox->xmin = p->x;
+	gbox->ymax = gbox->ymin = p->y;
+	gbox->zmax = gbox->zmin = p->z;
+	gbox->mmax = gbox->mmin = p->m;
+
+	for (uint32_t i = 1; i < pa->npoints; i++)
+	{
+		p = getPoint4d_cp(pa, i);
+		gbox->xmin = FP_MIN(gbox->xmin, p->x);
+		gbox->xmax = FP_MAX(gbox->xmax, p->x);
+		gbox->ymin = FP_MIN(gbox->ymin, p->y);
+		gbox->ymax = FP_MAX(gbox->ymax, p->y);
+		gbox->zmin = FP_MIN(gbox->zmin, p->z);
+		gbox->zmax = FP_MAX(gbox->zmax, p->z);
+		gbox->mmin = FP_MIN(gbox->mmin, p->m);
+		gbox->mmax = FP_MAX(gbox->mmax, p->m);
+	}
+}
+
+int
+ptarray_calculate_gbox_cartesian(const POINTARRAY *pa, GBOX *gbox)
+{
+	if (!pa || pa->npoints == 0)
+		return LW_FAILURE;
+	if (!gbox)
+		return LW_FAILURE;
+
+	int has_z = FLAGS_GET_Z(pa->flags);
+	int has_m = FLAGS_GET_M(pa->flags);
 	gbox->flags = lwflags(has_z, has_m, 0);
 	LWDEBUGF(4, "ptarray_calculate_gbox Z: %d M: %d", has_z, has_m);
+	int coordinates = 2 + has_z + has_m;
 
-	getPoint4d_p(pa, 0, &p);
-	gbox->xmin = gbox->xmax = p.x;
-	gbox->ymin = gbox->ymax = p.y;
-	if ( has_z )
-		gbox->zmin = gbox->zmax = p.z;
-	if ( has_m )
-		gbox->mmin = gbox->mmax = p.m;
-
-	for ( i = 1 ; i < pa->npoints; i++ )
+	switch (coordinates)
 	{
-		getPoint4d_p(pa, i, &p);
-		gbox->xmin = FP_MIN(gbox->xmin, p.x);
-		gbox->xmax = FP_MAX(gbox->xmax, p.x);
-		gbox->ymin = FP_MIN(gbox->ymin, p.y);
-		gbox->ymax = FP_MAX(gbox->ymax, p.y);
-		if ( has_z )
+	case 2:
+	{
+		ptarray_calculate_gbox_cartesian_2d(pa, gbox);
+		break;
+	}
+	case 3:
+	{
+		if (has_z)
 		{
-			gbox->zmin = FP_MIN(gbox->zmin, p.z);
-			gbox->zmax = FP_MAX(gbox->zmax, p.z);
+			ptarray_calculate_gbox_cartesian_3d(pa, gbox);
 		}
-		if ( has_m )
+		else
 		{
-			gbox->mmin = FP_MIN(gbox->mmin, p.m);
-			gbox->mmax = FP_MAX(gbox->mmax, p.m);
+			double zmin = gbox->zmin;
+			double zmax = gbox->zmax;
+			ptarray_calculate_gbox_cartesian_3d(pa, gbox);
+			gbox->mmin = gbox->zmin;
+			gbox->mmax = gbox->zmax;
+			gbox->zmin = zmin;
+			gbox->zmax = zmax;
 		}
+		break;
 	}
+	default:
+	{
+		ptarray_calculate_gbox_cartesian_4d(pa, gbox);
+		break;
+	}
+	}
 	return LW_SUCCESS;
 }
 

Modified: trunk/liblwgeom/liblwgeom.h.in
===================================================================
--- trunk/liblwgeom/liblwgeom.h.in	2019-09-12 17:47:47 UTC (rev 17807)
+++ trunk/liblwgeom/liblwgeom.h.in	2019-09-13 15:40:41 UTC (rev 17808)
@@ -953,22 +953,6 @@
  */
 extern int getPoint2d_p(const POINTARRAY *pa, uint32_t n, POINT2D *point);
 
-/**
-* Returns a POINT3DZ pointer into the POINTARRAY serialized_ptlist,
-* suitable for reading from. This is very high performance
-* and declared const because you aren't allowed to muck with the
-* values, only read them.
-*/
-extern const POINT3DZ* getPoint3dz_cp(const POINTARRAY *pa, uint32_t n);
-
-/**
-* Returns a POINT4D pointer into the POINTARRAY serialized_ptlist,
-* suitable for reading from. This is very high performance
-* and declared const because you aren't allowed to muck with the
-* values, only read them.
-*/
-extern const POINT4D* getPoint4d_cp(const POINTARRAY *pa, uint32_t n);
-
 /*
  * set point N to the given value
  * NOTE that the pointarray can be of any

Modified: trunk/liblwgeom/lwgeom_api.c
===================================================================
--- trunk/liblwgeom/lwgeom_api.c	2019-09-12 17:47:47 UTC (rev 17807)
+++ trunk/liblwgeom/lwgeom_api.c	2019-09-13 15:40:41 UTC (rev 17808)
@@ -362,46 +362,6 @@
 	return 1;
 }
 
-const POINT3DZ*
-getPoint3dz_cp(const POINTARRAY *pa, uint32_t n)
-{
-	if ( ! pa ) return 0;
-
-	if ( ! FLAGS_GET_Z(pa->flags) )
-	{
-		lwerror("getPoint3dz_cp: no Z coordinates in point array");
-		return 0; /*error */
-	}
-
-	if ( n>=pa->npoints )
-	{
-		lwerror("getPoint3dz_cp: point offset out of range");
-		return 0; /*error */
-	}
-
-	return (const POINT3DZ*)getPoint_internal(pa, n);
-}
-
-const POINT4D*
-getPoint4d_cp(const POINTARRAY* pa, uint32_t n)
-{
-	if (!pa) return 0;
-
-	if (!(FLAGS_GET_Z(pa->flags) && FLAGS_GET_M(pa->flags)))
-	{
-		lwerror("getPoint4d_cp: no Z and M coordinates in point array");
-		return 0; /*error */
-	}
-
-	if (n >= pa->npoints)
-	{
-		lwerror("getPoint4d_cp: point offset out of range");
-		return 0; /*error */
-	}
-
-	return (const POINT4D*)getPoint_internal(pa, n);
-}
-
 /*
  * set point N to the given value
  * NOTE that the pointarray can be of any

Modified: trunk/liblwgeom/lwgeom_geos.c
===================================================================
--- trunk/liblwgeom/lwgeom_geos.c	2019-09-12 17:47:47 UTC (rev 17807)
+++ trunk/liblwgeom/lwgeom_geos.c	2019-09-13 15:40:41 UTC (rev 17808)
@@ -243,7 +243,7 @@
 	uint32_t dims = 2;
 	uint32_t i;
 	int append_points = 0;
-	const POINT3DZ* p3d = NULL;
+	const POINT3D *p3d = NULL;
 	const POINT2D* p2d = NULL;
 	GEOSCoordSeq sq;
 
@@ -273,7 +273,7 @@
 	{
 		if (dims == 3)
 		{
-			p3d = getPoint3dz_cp(pa, i);
+			p3d = getPoint3d_cp(pa, i);
 			p2d = (const POINT2D*)p3d;
 			LWDEBUGF(4, "Point: %g,%g,%g", p3d->x, p3d->y, p3d->z);
 		}
@@ -293,7 +293,7 @@
 	{
 		if (dims == 3)
 		{
-			p3d = getPoint3dz_cp(pa, 0);
+			p3d = getPoint3d_cp(pa, 0);
 			p2d = (const POINT2D*)p3d;
 		}
 		else

Modified: trunk/liblwgeom/lwinline.h
===================================================================
--- trunk/liblwgeom/lwinline.h	2019-09-12 17:47:47 UTC (rev 17807)
+++ trunk/liblwgeom/lwinline.h	2019-09-13 15:40:41 UTC (rev 17808)
@@ -90,12 +90,33 @@
 static inline const POINT2D *
 getPoint2d_cp(const POINTARRAY *pa, uint32_t n)
 {
-	if (!pa)
-		return 0;
-
 	return (const POINT2D *)getPoint_internal(pa, n);
 }
 
+/**
+ * Returns a POINT2D pointer into the POINTARRAY serialized_ptlist,
+ * suitable for reading from. This is very high performance
+ * and declared const because you aren't allowed to muck with the
+ * values, only read them.
+ */
+static inline const POINT3D *
+getPoint3d_cp(const POINTARRAY *pa, uint32_t n)
+{
+	return (const POINT3D *)getPoint_internal(pa, n);
+}
+
+/**
+ * Returns a POINT2D pointer into the POINTARRAY serialized_ptlist,
+ * suitable for reading from. This is very high performance
+ * and declared const because you aren't allowed to muck with the
+ * values, only read them.
+ */
+static inline const POINT4D *
+getPoint4d_cp(const POINTARRAY *pa, uint32_t n)
+{
+	return (const POINT4D *)getPoint_internal(pa, n);
+}
+
 static inline LWPOINT *
 lwgeom_as_lwpoint(const LWGEOM *lwgeom)
 {

Modified: trunk/liblwgeom/lwout_geojson.c
===================================================================
--- trunk/liblwgeom/lwout_geojson.c	2019-09-12 17:47:47 UTC (rev 17807)
+++ trunk/liblwgeom/lwout_geojson.c	2019-09-13 15:40:41 UTC (rev 17808)
@@ -752,8 +752,7 @@
 	{
 		for (i=0; i<pa->npoints; i++)
 		{
-			const POINT3DZ *pt;
-			pt = getPoint3dz_cp(pa, i);
+			const POINT3D *pt = getPoint3d_cp(pa, i);
 
 			lwprint_double(
 			    pt->x, precision, x, OUT_DOUBLE_BUFFER_SIZE);

Modified: trunk/liblwgeom/lwout_gml.c
===================================================================
--- trunk/liblwgeom/lwout_gml.c	2019-09-12 17:47:47 UTC (rev 17807)
+++ trunk/liblwgeom/lwout_gml.c	2019-09-13 15:40:41 UTC (rev 17808)
@@ -689,8 +689,7 @@
 	{
 		for (i=0; i<pa->npoints; i++)
 		{
-			const POINT3DZ *pt;
-			pt = getPoint3dz_cp(pa, i);
+			const POINT3D *pt = getPoint3d_cp(pa, i);
 			lwprint_double(
 			    pt->x, precision, x, OUT_DOUBLE_BUFFER_SIZE);
 			lwprint_double(
@@ -1919,8 +1918,7 @@
 	{
 		for (i=0; i<pa->npoints; i++)
 		{
-			const POINT3DZ *pt;
-			pt = getPoint3dz_cp(pa, i);
+			const POINT3D *pt = getPoint3d_cp(pa, i);
 
 			lwprint_double(
 			    pt->x, precision, x, OUT_DOUBLE_BUFFER_SIZE);



More information about the postgis-tickets mailing list