[postgis-users] ST_Transform weirdness in 3.1.5

[Jeff Hoffmann]

Thanks at least for looking into this.  One other thing I just
noticed: If I transform from EPSG:4269 -> EPSG:4326 first, it seems to
work as expected.

select st_asewkt(st_transform('SRID=4326;POLYGON((-81.11425169432798
25.139365690162226,-79.88375846050846 25.97,-79.8759485474324
25.97511821217247,-81.11505127181294
25.97943803781699,-81.11425169432798 25.139365690162226))'::geometry,
3087));

                          st_asewkt
-----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
 SRID=3087;POLYGON((690648.4490696276
129676.1163013646,811467.6067641322 225287.6132890527,812229.019377012
225881.33726899035,688387.9040682416
222840.47302733874,690648.4490696276 129676.1163013646))
(1 row)

I have another server configured as follows (also v6.2.1 of PROJ) that
works fine, one of the differences is the postgis version, which is
why I asked about whether it could be a postgis problem.

                                             version
-------------------------------------------------------------------------------------------------
 PostgreSQL 12.4 on x86_64-pc-linux-gnu, compiled by
x86_64-pc-linux-gnu-gcc (GCC) 7.4.0, 64-bit
(1 row)


 postgis_full_version
---------------------------------------------------------------------------------------------------------------------------------------------------------------------
 POSTGIS="3.0.2 2fb2a18" [EXTENSION] PGSQL="120"
GEOS="3.8.1-CAPI-1.13.3" PROJ="6.2.1" LIBXML="2.9.9" LIBJSON="0.12.99"
LIBPROTOBUF="1.3.0" WAGYU="0.4.3 (Internal)"
(1 row)

Could this be something with the proj configuration? Does calling
ST_Transform using an SRID use some sort of internal proj
configuration that may be weird or does it pull the proj4text from the
spatial_ref_sys table? There are other differences that make this
weird (big one is the weird server is on ARM, everything else is x86)
so as much as anything I'm just trying to figure out where to start
looking if nothing else.

On Mon, Apr 4, 2022 at 12:01 PM Paul Ramsey <pramsey at cleverelephant.ca> wrote:
>
> Not super obvious as I'm not seeing a reproduction here on a more recent proj. Anyone with proj6 able to reproduce?
> P
>
> > On Apr 4, 2022, at 9:18 AM, Jeff Hoffmann <jeff at propertykey.com> wrote:
> >
> > I am running PostgreSQL via Amazon RDS. Here is the version
> > information for the instance that exhibits the problem:
> >
> >                                                   version
> > -------------------------------------------------------------------------------------------------------------
> > PostgreSQL 12.8 on aarch64-unknown-linux-gnu, compiled by
> > aarch64-unknown-linux-gnu-gcc (GCC) 7.4.0, 64-bit
> > (1 row)
> >
> > postgis_full_version
> > ---------------------------------------------------------------------------------------------------------------------------------------------------------------
> > POSTGIS="3.1.5 0" [EXTENSION] PGSQL="120" GEOS="3.8.2-CAPI-1.13.4"
> > PROJ="6.2.1" LIBXML="2.9.9" LIBJSON="0.12.99" LIBPROTOBUF="1.3.0"
> > WAGYU="0.5.0 (Internal)"
> > (1 row)
> >
> > Here is an example of a query that produces weird results:
> >
> >
> > select st_asewkt(st_transform('SRID=4269;POLYGON((-81.11425169432798
> > 25.139365690162226,-79.88375846050846
> > 25.135204935817377,-79.8759485474324
> > 25.97511821217247,-81.11505127181294
> > 25.97943803781699,-81.11425169432798 25.139365690162226))'::geometry,
> > 3087));
> >
> >        st_asewkt
> > -----------------------------------------------------------------------------------------------------------------------------------------------------------------------
> > SRID=3087;POLYGON((690648.3378382007 129676.3544158864,Infinity
> > Infinity,Infinity Infinity,688388.1492655604
> > 222840.53520809134,690648.3378382007 129676.3544158864))
> > (1 row)
> >
> > Two of the points are coming up as (Infinty,Infinity).  When I run the
> > same query using the proj4text found in the spatial_ref_sys table
> > instead of the srid I get the following output:
> >
> > select st_asewkt(st_transform('SRID=4269;POLYGON((-81.11425169432798
> > 25.139365690162226,-79.88375846050846
> > 25.135204935817377,-79.8759485474324
> > 25.97511821217247,-81.11505127181294
> > 25.97943803781699,-81.11425169432798 25.139365690162226))'::geometry,
> > '+proj=aea +lat_1=24 +lat_2=31.5 +lat_0=24 +lon_0=-84 +x_0=400000
> > +y_0=0 +ellps=GRS80 +towgs84=0,0,0,0,0,0,0 +units=m +no_defs'));
> >
> >                     st_asewkt
> > -------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------
> > POLYGON((690648.3378382007 129676.3544158864,814558.219082916
> > 132733.37665360732,812229.019377012
> > 225881.33726899035,688388.1492655604
> > 222840.53520809134,690648.3378382007 129676.3544158864))
> > (1 row)
> >
> > This looks pretty close to right.  What is the disconnect between
> > using the srid to call the transform vs the proj4text version of
> > st_transform?
> >
> > --
> > Jeff Hoffmann
> > PropertyKey
> > jeff at propertykey.com
> > _______________________________________________
> > postgis-users mailing list
> > postgis-users at lists.osgeo.org
> > https://lists.osgeo.org/mailman/listinfo/postgis-users
>
> _______________________________________________
> postgis-users mailing list
> postgis-users at lists.osgeo.org
> https://lists.osgeo.org/mailman/listinfo/postgis-users



-- 
Jeff Hoffmann
Head Plate Spinner
PropertyKey
Office: 612-800-6444
jeff at propertykey.com