[postgis-tickets] r17541 - PROJ6: Speed improvements

Raul raul at rmr.ninja
Wed Jun 19 08:12:21 PDT 2019


Author: algunenano
Date: 2019-06-19 08:12:21 -0700 (Wed, 19 Jun 2019)
New Revision: 17541

Modified:
   trunk/NEWS
   trunk/liblwgeom/liblwgeom.h.in
   trunk/liblwgeom/lwgeom_transform.c
   trunk/libpgcommon/lwgeom_cache.h
   trunk/libpgcommon/lwgeom_transform.c
   trunk/libpgcommon/lwgeom_transform.h
   trunk/postgis/lwgeom_in_gml.c
   trunk/postgis/lwgeom_transform.c
Log:
PROJ6: Speed improvements

Caches more elements to avoid recalculation per row

References #4372


Modified: trunk/NEWS
===================================================================
--- trunk/NEWS	2019-06-19 15:09:25 UTC (rev 17540)
+++ trunk/NEWS	2019-06-19 15:12:21 UTC (rev 17541)
@@ -13,7 +13,8 @@
   - #4388, AddRasterConstraints: Ignore NULLs when generating constraints (Raúl Marín)
   - #4327, Avoid pfree'ing the result of getenv (Raúl Marín)
   - #4406, Throw on invalid characters when decoding geohash (Raúl Marín)
-  - #4372, Avoid resource leaks with PROJ6 (Raúl Marín)
+  - #4429, Avoid resource leaks with PROJ6 (Raúl Marín)
+  - #4372, PROJ6: Speed improvements (Raúl Marín)
 
 PostGIS 3.0.0alpha2
 2019/06/02

Modified: trunk/liblwgeom/liblwgeom.h.in
===================================================================
--- trunk/liblwgeom/liblwgeom.h.in	2019-06-19 15:09:25 UTC (rev 17540)
+++ trunk/liblwgeom/liblwgeom.h.in	2019-06-19 15:12:21 UTC (rev 17541)
@@ -30,8 +30,8 @@
 #define _LIBLWGEOM_H 1
 
 #include <stdarg.h>
+#include <stdint.h>
 #include <stdio.h>
-#include <stdint.h>
 
 #include "../postgis_config.h"
 
@@ -42,8 +42,28 @@
     projPJ pj_from;
     projPJ pj_to;
 } PJ;
+
+typedef PJ LWPROJ;
+
 #else
 #include "proj.h"
+
+/* For PROJ6 we cache several extra values to avoid calls to proj_get_source_crs
+ * or proj_get_target_crs since those are very costly
+ */
+typedef struct LWPROJ
+{
+	PJ* pj;
+        /* CRSs are swapped: Used in transformation calls */
+	uint8_t source_swapped;
+	uint8_t target_swapped;
+        /* Source crs is geographic: Used in geography calls (source srid == dst srid) */
+        uint8_t source_is_latlong;
+
+        /* Source ellipsoid parameters */
+        double source_semi_major_metre;
+        double source_semi_minor_metre;
+} LWPROJ;
 #endif
 
 #if POSTGIS_PROJ_VERSION < 49
@@ -2358,7 +2378,7 @@
  *
  * Eg: "+proj=longlat +ellps=WGS84 +datum=WGS84 +no_defs"
  */
-projPJ lwproj_from_string(const char* txt);
+projPJ projpj_from_string(const char* txt);
 
 #endif
 /**
@@ -2366,10 +2386,21 @@
  * @param geom the geometry to transform
  * @param PJ the input and output
  */
-int lwgeom_transform(LWGEOM *geom, PJ* pj);
-int ptarray_transform(POINTARRAY *pa, PJ* pj);
+int lwgeom_transform(LWGEOM *geom, LWPROJ* pj);
+int ptarray_transform(POINTARRAY *pa, LWPROJ* pj);
 
+#if POSTGIS_PROJ_VERSION >= 60
 
+/**
+ * Allocate a new LWPROJ containing the reference to the PROJ's PJ
+ * If extra_geography_data is true, it will generate the following values for
+ * the source srs: is_latlong (geometric or not) and spheroid values
+ */
+LWPROJ *lwproj_from_PJ(PJ *pj, int8_t extra_geography_data);
+
+#endif
+
+
 /*******************************************************************************
  * GEOS-dependent extra functions on LWGEOM
  ******************************************************************************/

Modified: trunk/liblwgeom/lwgeom_transform.c
===================================================================
--- trunk/liblwgeom/lwgeom_transform.c	2019-06-19 15:09:25 UTC (rev 17540)
+++ trunk/liblwgeom/lwgeom_transform.c	2019-06-19 15:12:21 UTC (rev 17541)
@@ -48,9 +48,8 @@
 
 #if POSTGIS_PROJ_VERSION < 60
 
-
 static int
-point4d_transform(POINT4D *pt, PJ* pj)
+point4d_transform(POINT4D *pt, LWPROJ *pj)
 {
 	POINT3D orig_pt = {pt->x, pt->y, pt->z}; /* Copy for error report*/
 
@@ -90,7 +89,7 @@
  * from inpj projection to outpj projection
  */
 int
-ptarray_transform(POINTARRAY *pa, PJ* pj)
+ptarray_transform(POINTARRAY *pa, LWPROJ *pj)
 {
 	uint32_t i;
 	POINT4D p;
@@ -110,9 +109,9 @@
 {
 	char *pj_errstr;
 	int rv;
-	PJ pj;
+	LWPROJ pj;
 
-	pj.pj_from = lwproj_from_string(instr);
+	pj.pj_from = projpj_from_string(instr);
 	if (!pj.pj_from)
 	{
 		pj_errstr = pj_strerrno(*pj_get_errno_ref());
@@ -121,7 +120,7 @@
 		return LW_FAILURE;
 	}
 
-	pj.pj_to = lwproj_from_string(outstr);
+	pj.pj_to = projpj_from_string(outstr);
 	if (!pj.pj_to)
 	{
 		pj_free(pj.pj_from);
@@ -142,7 +141,7 @@
  * from inpj projection to outpj projection
  */
 int
-lwgeom_transform(LWGEOM *geom, PJ* pj)
+lwgeom_transform(LWGEOM *geom, LWPROJ *pj)
 {
 	uint32_t i;
 
@@ -199,7 +198,7 @@
 }
 
 projPJ
-lwproj_from_string(const char *str1)
+projpj_from_string(const char *str1)
 {
 	if (!str1 || str1[0] == '\0')
 	{
@@ -212,6 +211,132 @@
 
 #else /* POSTGIS_PROJ_VERION >= 60 */
 
+static uint8_t
+proj_crs_is_swapped(const PJ *pj_crs)
+{
+	PJ *pj_cs;
+	uint8_t rv = LW_FALSE;
+
+	if (proj_get_type(pj_crs) == PJ_TYPE_COMPOUND_CRS)
+	{
+		PJ *pj_horiz_crs = proj_crs_get_sub_crs(NULL, pj_crs, 0);
+		if (!pj_horiz_crs)
+			lwerror("%s: proj_crs_get_sub_crs returned NULL", __func__);
+		pj_cs = proj_crs_get_coordinate_system(NULL, pj_horiz_crs);
+		proj_destroy(pj_horiz_crs);
+	}
+	else if (proj_get_type(pj_crs) == PJ_TYPE_BOUND_CRS)
+	{
+		PJ *pj_src_crs = proj_get_source_crs(NULL, pj_crs);
+		if (!pj_src_crs)
+			lwerror("%s: proj_get_source_crs returned NULL", __func__);
+		pj_cs = proj_crs_get_coordinate_system(NULL, pj_src_crs);
+		proj_destroy(pj_src_crs);
+	}
+	else
+	{
+		pj_cs = proj_crs_get_coordinate_system(NULL, pj_crs);
+	}
+	if (!pj_cs)
+		lwerror("%s: proj_crs_get_coordinate_system returned NULL", __func__);
+	int axis_count = proj_cs_get_axis_count(NULL, pj_cs);
+	if (axis_count > 0)
+	{
+		const char *out_name, *out_abbrev, *out_direction;
+		double out_unit_conv_factor;
+		const char *out_unit_name, *out_unit_auth_name, *out_unit_code;
+		/* Read only first axis, see if it is degrees / north */
+		proj_cs_get_axis_info(NULL,
+				      pj_cs,
+				      0,
+				      &out_name,
+				      &out_abbrev,
+				      &out_direction,
+				      &out_unit_conv_factor,
+				      &out_unit_name,
+				      &out_unit_auth_name,
+				      &out_unit_code);
+		rv = (strcasecmp(out_direction, "north") == 0);
+	}
+	proj_destroy(pj_cs);
+	return rv;
+}
+
+LWPROJ *
+lwproj_from_PJ(PJ *pj, int8_t extra_geography_data)
+{
+	PJ *pj_source_crs = proj_get_source_crs(NULL, pj);
+	uint8_t source_is_latlong = LW_FALSE;
+	double out_semi_major_metre = DBL_MAX, out_semi_minor_metre = DBL_MAX;
+
+	if (!pj_source_crs)
+	{
+		lwerror("%s: unable to access source crs", __func__);
+		return NULL;
+	}
+	uint8_t source_swapped = proj_crs_is_swapped(pj_source_crs);
+
+	/* We only care about the extra values if there is no transformation */
+	if (!extra_geography_data)
+	{
+		proj_destroy(pj_source_crs);
+	}
+	else
+	{
+		PJ *pj_ellps;
+		double out_inv_flattening;
+		int out_is_semi_minor_computed;
+
+		PJ_TYPE pj_type = proj_get_type(pj_source_crs);
+		if (pj_type == PJ_TYPE_UNKNOWN)
+		{
+			proj_destroy(pj_source_crs);
+			lwerror("%s: unable to access source crs type", __func__);
+			return NULL;
+		}
+		source_is_latlong = (pj_type == PJ_TYPE_GEOGRAPHIC_2D_CRS) || (pj_type == PJ_TYPE_GEOGRAPHIC_3D_CRS);
+
+		pj_ellps = proj_get_ellipsoid(NULL, pj_source_crs);
+		proj_destroy(pj_source_crs);
+		if (!pj_ellps)
+		{
+			lwerror("%s: unable to access source crs ellipsoid", __func__);
+			return NULL;
+		}
+		if (!proj_ellipsoid_get_parameters(NULL,
+						   pj_ellps,
+						   &out_semi_major_metre,
+						   &out_semi_minor_metre,
+						   &out_is_semi_minor_computed,
+						   &out_inv_flattening))
+		{
+			proj_destroy(pj_ellps);
+			lwerror("%s: unable to access source crs ellipsoid parameters", __func__);
+			return NULL;
+		}
+		proj_destroy(pj_ellps);
+	}
+
+	PJ *pj_target_crs = proj_get_target_crs(NULL, pj);
+	if (!pj_target_crs)
+	{
+		lwerror("%s: unable to access target crs", __func__);
+		return NULL;
+	}
+	uint8_t target_swapped = proj_crs_is_swapped(pj_target_crs);
+	proj_destroy(pj_target_crs);
+
+	LWPROJ *lp = malloc(sizeof(LWPROJ));
+	lp->pj = pj;
+	lp->source_swapped = source_swapped;
+	lp->target_swapped = target_swapped;
+	lp->source_is_latlong = source_is_latlong;
+	lp->source_semi_major_metre = out_semi_major_metre;
+	lp->source_semi_minor_metre = out_semi_minor_metre;
+
+	return lp;
+}
+
 int
 lwgeom_transform_from_str(LWGEOM *geom, const char* instr, const char* outstr)
 {
@@ -235,14 +360,18 @@
 		return LW_FAILURE;
 	}
 
-	int ret = lwgeom_transform(geom, pj);
+	LWPROJ *lp = lwproj_from_PJ(pj, LW_FALSE);
+
+	int ret = lwgeom_transform(geom, lp);
+
 	proj_destroy(pj);
+	free(lp);
 
 	return ret;
 }
 
 int
-lwgeom_transform(LWGEOM* geom, PJ* pj)
+lwgeom_transform(LWGEOM *geom, LWPROJ *pj)
 {
 	uint32_t i;
 
@@ -301,56 +430,8 @@
 	return LW_SUCCESS;
 }
 
-static int
-proj_crs_is_swapped(const PJ* pj_crs)
-{
-	PJ *pj_cs;
-	int rv = LW_FALSE;
-
-	if (proj_get_type(pj_crs) == PJ_TYPE_COMPOUND_CRS)
-	{
-		PJ *pj_horiz_crs = proj_crs_get_sub_crs(NULL, pj_crs, 0);
-		if (!pj_horiz_crs) lwerror("%s: proj_crs_get_sub_crs returned NULL", __func__);
-		pj_cs = proj_crs_get_coordinate_system(NULL, pj_horiz_crs);
-		proj_destroy(pj_horiz_crs);
-	}
-	else if (proj_get_type(pj_crs) == PJ_TYPE_BOUND_CRS)
-	{
-		PJ *pj_src_crs = proj_get_source_crs(NULL, pj_crs);
-		if (!pj_src_crs) lwerror("%s: proj_get_source_crs returned NULL", __func__);
-		pj_cs = proj_crs_get_coordinate_system(NULL, pj_src_crs);
-		proj_destroy(pj_src_crs);
-	}
-	else
-	{
-		pj_cs = proj_crs_get_coordinate_system(NULL, pj_crs);
-	}
-	if (!pj_cs) lwerror("%s: proj_crs_get_coordinate_system returned NULL", __func__);
-	int axis_count = proj_cs_get_axis_count(NULL, pj_cs);
-	if (axis_count > 0)
-	{
-		const char *out_name, *out_abbrev, *out_direction;
-		double out_unit_conv_factor;
-		const char *out_unit_name, *out_unit_auth_name, *out_unit_code;
-		/* Read only first axis, see if it is degrees / north */
-		proj_cs_get_axis_info(NULL, pj_cs, 0,
-			&out_name,
-			&out_abbrev,
-			&out_direction,
-			&out_unit_conv_factor,
-			&out_unit_name,
-			&out_unit_auth_name,
-			&out_unit_code
-			);
-		rv = (strcasecmp(out_direction, "north") == 0);
-	}
-	proj_destroy(pj_cs);
-	return rv;
-}
-
-
 int
-ptarray_transform(POINTARRAY* pa, PJ* pj)
+ptarray_transform(POINTARRAY *pa, LWPROJ *pj)
 {
 	uint32_t i;
 	POINT4D p;
@@ -359,24 +440,9 @@
 	size_t point_size = ptarray_point_size(pa);
 	int has_z = ptarray_has_z(pa);
 	double *pa_double = (double*)(pa->serialized_pointlist);
-	int input_swapped, output_swapped;
 
-	PJ* pj_source_crs = proj_get_source_crs(NULL, pj);
-	PJ* pj_target_crs = proj_get_target_crs(NULL, pj);
-
-	if (!(pj_source_crs && pj_target_crs))
-	{
-		lwerror("ptarray_transform: unable to access source and target crs");
-		return LW_FAILURE;
-	}
-
-	input_swapped = proj_crs_is_swapped(pj_source_crs);
-	output_swapped = proj_crs_is_swapped(pj_target_crs);
-	proj_destroy(pj_source_crs);
-	proj_destroy(pj_target_crs);
-
 	/* Convert to radians if necessary */
-	if (proj_angular_input(pj, PJ_FWD))
+	if (proj_angular_input(pj->pj, PJ_FWD)) /* CACHE THIS TOO */
 	{
 		for (i = 0; i < pa->npoints; i++)
 		{
@@ -385,7 +451,7 @@
 		}
 	}
 
-	if (input_swapped)
+	if (pj->source_swapped)
 		ptarray_swap_ordinates(pa, LWORD_X, LWORD_Y);
 
 	/*
@@ -396,15 +462,21 @@
 	* double *t, size_t st, size_t nt)
 	*/
 
-	n_converted = proj_trans_generic(
-		pj, PJ_FWD,
-		pa_double,     point_size, n_points, /* X */
-		pa_double + 1, point_size, n_points, /* Y */
-		has_z ? pa_double + 2 : NULL,
-		has_z ? point_size : 0,
-		has_z ? n_points : 0, /* Z */
-		NULL, 0, 0 /* M */
-		);
+	n_converted = proj_trans_generic(pj->pj,
+					 PJ_FWD,
+					 pa_double,
+					 point_size,
+					 n_points, /* X */
+					 pa_double + 1,
+					 point_size,
+					 n_points, /* Y */
+					 has_z ? pa_double + 2 : NULL,
+					 has_z ? point_size : 0,
+					 has_z ? n_points : 0, /* Z */
+					 NULL,
+					 0,
+					 0 /* M */
+	);
 
 	if (n_converted != n_points)
 	{
@@ -413,7 +485,7 @@
 		return LW_FAILURE;
 	}
 
-	int pj_errno_val = proj_errno(pj);
+	int pj_errno_val = proj_errno(pj->pj);
 	if (pj_errno_val)
 	{
 		lwerror("transform: %s (%d)",
@@ -421,11 +493,11 @@
 		return LW_FAILURE;
 	}
 
-	if (output_swapped)
+	if (pj->target_swapped)
 		ptarray_swap_ordinates(pa, LWORD_X, LWORD_Y);
 
 	/* Convert radians to degrees if necessary */
-	if (proj_angular_output(pj, PJ_FWD))
+	if (proj_angular_output(pj->pj, PJ_FWD))
 	{
 		for (i = 0; i < pa->npoints; i++)
 		{

Modified: trunk/libpgcommon/lwgeom_cache.h
===================================================================
--- trunk/libpgcommon/lwgeom_cache.h	2019-06-19 15:09:25 UTC (rev 17540)
+++ trunk/libpgcommon/lwgeom_cache.h	2019-06-19 15:12:21 UTC (rev 17541)
@@ -67,7 +67,7 @@
 {
 	int32_t srid_from;
 	int32_t srid_to;
-	PJ* projection;
+	LWPROJ *projection;
 	MemoryContext projection_mcxt;
 }
 PROJSRSCacheItem;

Modified: trunk/libpgcommon/lwgeom_transform.c
===================================================================
--- trunk/libpgcommon/lwgeom_transform.c	2019-06-19 15:09:25 UTC (rev 17540)
+++ trunk/libpgcommon/lwgeom_transform.c	2019-06-19 15:12:21 UTC (rev 17541)
@@ -77,7 +77,7 @@
 typedef struct struct_PJHashEntry
 {
 	MemoryContext ProjectionContext;
-	PJ* projection;
+	LWPROJ *projection;
 }
 PJHashEntry;
 
@@ -87,8 +87,8 @@
 
 static HTAB *CreatePJHash(void);
 static void DeletePJHashEntry(MemoryContext mcxt);
-static PJ* GetPJHashEntry(MemoryContext mcxt);
-static void AddPJHashEntry(MemoryContext mcxt, PJ* projection);
+static LWPROJ *GetPJHashEntry(MemoryContext mcxt);
+static void AddPJHashEntry(MemoryContext mcxt, LWPROJ *projection);
 
 /* Internal Cache API */
 /* static PROJPortalCache *GetPROJSRSCache(FunctionCallInfo fcinfo) ; */
@@ -127,7 +127,7 @@
 }
 
 static void
-PROJSRSDestroyPJ(PJ* pj)
+PROJSRSDestroyPJ(LWPROJ *pj)
 {
 #if POSTGIS_PROJ_VERSION < 60
 /* Ape the Proj 6+ API for versions < 6 */
@@ -137,7 +137,8 @@
 		pj_free(pj->pj_to);
 	free(pj);
 #else
-	proj_destroy(pj);
+	proj_destroy(pj->pj);
+	free(pj);
 #endif
 }
 
@@ -163,7 +164,7 @@
 #endif
 
 	/* Lookup the PJ pointer in the global hash table so we can free it */
-	PJ* projection = GetPJHashEntry(context);
+	LWPROJ *projection = GetPJHashEntry(context);
 
 	if (!projection)
 		elog(ERROR, "PROJSRSCacheDelete: Trying to delete non-existant projection object with MemoryContext key (%p)", (void *)context);
@@ -278,7 +279,8 @@
 	return hash_create("PostGIS PROJ Backend MemoryContext Hash", PROJ_BACKEND_HASH_SIZE, &ctl, (HASH_ELEM | HASH_FUNCTION));
 }
 
-static void AddPJHashEntry(MemoryContext mcxt, PJ* projection)
+static void
+AddPJHashEntry(MemoryContext mcxt, LWPROJ *projection)
 {
 	bool found;
 	void **key;
@@ -301,7 +303,8 @@
 	}
 }
 
-static PJ* GetPJHashEntry(MemoryContext mcxt)
+static LWPROJ *
+GetPJHashEntry(MemoryContext mcxt)
 {
 	void **key;
 	PJHashEntry *he;
@@ -356,7 +359,7 @@
 	return false;
 }
 
-static PJ *
+static LWPROJ *
 GetProjectionFromPROJCache(PROJPortalCache *cache, int32_t srid_from, int32_t srid_to)
 {
 	uint32_t i;
@@ -644,8 +647,8 @@
 	PJ* projection = malloc(sizeof(PJ));
 	pj_from_str = from_strs.proj4text;
 	pj_to_str = to_strs.proj4text;
-	projection->pj_from = lwproj_from_string(pj_from_str);
-	projection->pj_to = lwproj_from_string(pj_to_str);
+	projection->pj_from = projpj_from_string(pj_from_str);
+	projection->pj_to = projpj_from_string(pj_to_str);
 
 	if (!projection->pj_from)
 		elog(ERROR,
@@ -657,7 +660,7 @@
 		    "could not form projection from 'srid=%d' to 'srid=%d'",
 		    srid_from, srid_to);
 #else
-	PJ* projection = NULL;
+	PJ *projpj = NULL;
 	/* Try combinations of ESPG/SRTEXT/PROJ4TEXT until we find */
 	/* one that gives us a usable transform. Note that we prefer */
 	/* EPSG numbers over SRTEXT and SRTEXT over PROJ4TEXT */
@@ -669,15 +672,20 @@
 		pj_to_str   = pgstrs_get_entry(&to_strs,   i % 3);
 		if (!(pj_from_str && pj_to_str))
 			continue;
-		projection = proj_create_crs_to_crs(NULL, pj_from_str, pj_to_str, NULL);
-		if (projection && !proj_errno(projection))
+		projpj = proj_create_crs_to_crs(NULL, pj_from_str, pj_to_str, NULL);
+		if (projpj && !proj_errno(projpj))
 			break;
 	}
+	if (!projpj)
+	{
+		elog(ERROR, "could not form projection (PJ) from 'srid=%d' to 'srid=%d'", srid_from, srid_to);
+		return;
+	}
+	LWPROJ *projection = lwproj_from_PJ(projpj, srid_from == srid_to);
 	if (!projection)
 	{
-		elog(ERROR,
-		    "could not form projection from 'srid=%d' to 'srid=%d'",
-		    srid_from, srid_to);
+		elog(ERROR, "could not form projection (LWPROJ) from 'srid=%d' to 'srid=%d'", srid_from, srid_to);
+		return;
 	}
 #endif
 
@@ -836,7 +844,7 @@
 }
 
 int
-GetPJUsingFCInfo(FunctionCallInfo fcinfo, int32_t srid_from, int32_t srid_to, PJ **pj)
+GetPJUsingFCInfo(FunctionCallInfo fcinfo, int32_t srid_from, int32_t srid_to, LWPROJ **pj)
 {
 	PROJPortalCache *proj_cache = NULL;
 
@@ -862,19 +870,12 @@
 }
 
 static int
-proj_pj_is_latlong(const PJ* pj)
+proj_pj_is_latlong(const LWPROJ *pj)
 {
 #if POSTGIS_PROJ_VERSION < 60
 	return pj_is_latlong(pj->pj_from);
 #else
-	PJ_TYPE pj_type;
-	PJ *pj_src_crs = proj_get_source_crs(NULL, pj);
-	if (!pj_src_crs)
-		elog(ERROR, "%s: proj_get_source_crs returned NULL", __func__);
-	pj_type = proj_get_type(pj_src_crs);
-	proj_destroy(pj_src_crs);
-	return (pj_type == PJ_TYPE_GEOGRAPHIC_2D_CRS) ||
-	       (pj_type == PJ_TYPE_GEOGRAPHIC_3D_CRS);
+	return pj->source_is_latlong;
 #endif
 }
 
@@ -881,7 +882,7 @@
 static int
 srid_is_latlong(FunctionCallInfo fcinfo, int32_t srid)
 {
-	PJ* pj;
+	LWPROJ *pj;
 	if ( GetPJUsingFCInfo(fcinfo, srid, srid, &pj) == LW_FAILURE)
 		return LW_FALSE;
 	return proj_pj_is_latlong(pj);
@@ -924,12 +925,8 @@
 int
 spheroid_init_from_srid(FunctionCallInfo fcinfo, int32_t srid, SPHEROID *s)
 {
-	PJ* pj;
-#if POSTGIS_PROJ_VERSION >= 60
-	double out_semi_major_metre, out_semi_minor_metre, out_inv_flattening;
-	int out_is_semi_minor_computed;
-	PJ *pj_ellps, *pj_crs;
-#elif POSTGIS_PROJ_VERSION >= 48
+	LWPROJ *pj;
+#if POSTGIS_PROJ_VERSION >= 48 && POSTGIS_PROJ_VERSION < 60
 	double major_axis, minor_axis, eccentricity_squared;
 #endif
 
@@ -937,25 +934,9 @@
 		return LW_FAILURE;
 
 #if POSTGIS_PROJ_VERSION >= 60
-	if (!proj_pj_is_latlong(pj))
+	if (!pj->source_is_latlong)
 		return LW_FAILURE;
-	pj_crs = proj_get_source_crs(NULL, pj);
-	if (!pj_crs)
-	{
-		lwerror("%s: proj_get_source_crs returned NULL", __func__);
-	}
-	pj_ellps = proj_get_ellipsoid(NULL, pj_crs);
-	if (!pj_ellps)
-	{
-		proj_destroy(pj_crs);
-		lwerror("%s: proj_get_ellipsoid returned NULL", __func__);
-	}
-	proj_ellipsoid_get_parameters(NULL, pj_ellps,
-		&out_semi_major_metre, &out_semi_minor_metre,
-		&out_is_semi_minor_computed, &out_inv_flattening);
-	proj_destroy(pj_ellps);
-	proj_destroy(pj_crs);
-	spheroid_init(s, out_semi_major_metre, out_semi_minor_metre);
+	spheroid_init(s, pj->source_semi_major_metre, pj->source_semi_minor_metre);
 
 #elif POSTGIS_PROJ_VERSION >= 48
 	if (!pj_is_latlong(pj->pj_from))

Modified: trunk/libpgcommon/lwgeom_transform.h
===================================================================
--- trunk/libpgcommon/lwgeom_transform.h	2019-06-19 15:09:25 UTC (rev 17540)
+++ trunk/libpgcommon/lwgeom_transform.h	2019-06-19 15:12:21 UTC (rev 17541)
@@ -32,7 +32,7 @@
 void SetPROJLibPath(void);
 bool IsInPROJCache(ProjCache cache, int32_t srid_from, int32_t srid_to);
 PJ *GetPJFromPROJCache(ProjCache cache, int32_t srid_from, int32_t srid_to);
-int GetPJUsingFCInfo(FunctionCallInfo fcinfo, int32_t srid_from, int32_t srid_to, PJ **pj);
+int GetPJUsingFCInfo(FunctionCallInfo fcinfo, int32_t srid_from, int32_t srid_to, LWPROJ **pj);
 int spheroid_init_from_srid(FunctionCallInfo fcinfo, int32_t srid, SPHEROID *s);
 void srid_check_latlong(FunctionCallInfo fcinfo, int32_t srid);
 srs_precision srid_axis_precision(FunctionCallInfo fcinfo, int32_t srid, int precision);

Modified: trunk/postgis/lwgeom_in_gml.c
===================================================================
--- trunk/postgis/lwgeom_in_gml.c	2019-06-19 15:09:25 UTC (rev 17540)
+++ trunk/postgis/lwgeom_in_gml.c	2019-06-19 15:12:21 UTC (rev 17541)
@@ -311,8 +311,8 @@
 	text_in = GetProj4String(srid_in);
 	text_out = GetProj4String(srid_out);
 
-	pj.pj_from = lwproj_from_string(text_in);
-	pj.pj_to = lwproj_from_string(text_out);
+	pj.pj_from = projpj_from_string(text_in);
+	pj.pj_to = projpj_from_string(text_out);
 
 	lwfree(text_in);
 	lwfree(text_out);
@@ -337,21 +337,39 @@
 gml_reproject_pa(POINTARRAY *pa, int32_t srid_in, int32_t srid_out)
 {
 	PJ *pj;
+	LWPROJ *lwp;
 	char text_in[32];
 	char text_out[32];
 
-	if (srid_in == SRID_UNKNOWN) return pa; /* nothing to do */
-	if (srid_out == SRID_UNKNOWN) gml_lwpgerror("invalid GML representation", 3);
+	if (srid_in == SRID_UNKNOWN)
+		return pa; /* nothing to do */
 
+	if (srid_out == SRID_UNKNOWN)
+	{
+		gml_lwpgerror("invalid GML representation", 3);
+		return NULL;
+	}
+
 	snprintf(text_in, 32, "EPSG:%d", srid_in);
 	snprintf(text_out, 32, "EPSG:%d", srid_out);
 	pj = proj_create_crs_to_crs(NULL, text_in, text_out, NULL);
 
-	if (ptarray_transform(pa, pj) == LW_FAILURE)
+	lwp = lwproj_from_PJ(pj, LW_FALSE);
+	if (!lwp)
 	{
+		proj_destroy(pj);
+		gml_lwpgerror("Could not create LWPROJ*", 57);
+		return NULL;
+	}
+
+	if (ptarray_transform(pa, lwp) == LW_FAILURE)
+	{
+		proj_destroy(pj);
 		elog(ERROR, "gml_reproject_pa: reprojection failed");
+		return NULL;
 	}
 	proj_destroy(pj);
+	free(lwp);
 
 	return pa;
 }

Modified: trunk/postgis/lwgeom_transform.c
===================================================================
--- trunk/postgis/lwgeom_transform.c	2019-06-19 15:09:25 UTC (rev 17540)
+++ trunk/postgis/lwgeom_transform.c	2019-06-19 15:12:21 UTC (rev 17541)
@@ -50,7 +50,7 @@
 	GSERIALIZED* geom;
 	GSERIALIZED* result=NULL;
 	LWGEOM* lwgeom;
-	PJ* pj;
+	LWPROJ *pj;
 	int32 srid_to, srid_from;
 
 	srid_to = PG_GETARG_INT32(1);
@@ -216,7 +216,7 @@
 
 	if (srid_from != srid_to)
 	{
-		PJ* pj;
+		LWPROJ *pj;
 		if (GetPJUsingFCInfo(fcinfo, srid_from, srid_to, &pj) == LW_FAILURE)
 		{
 			PG_FREE_IF_COPY(geom, 0);



More information about the postgis-tickets mailing list