Different Projection Results in AWS Managed Services
Greg Troxel
gdt at lexort.com
Mon Feb 17 15:34:02 PST 2025
Pat Blair <pblair at geocomm.com> writes:
> Below I am pasting some SQL that I hope will illustrate what we are
> seeing. When we transform a coordinate from EPSG:26860 (NAD83 /
> Minnesota Central (ftUS)) to WGS84 in two different environments (“AWS
> RDS” and “AWS Aurora”) we get a different result with a displacement
> of about 0.12 meters.
Multiple things going on here. You need two things to start:
Geodetic clarity on what you are doing. You are calling 0.12m a
discrepancy which means the allowable fuzz in understanding is zero.
And, unless you are working unusually closely with the US DoD, you
can't possibly be obtaining positions in WGS84 to anywhere near the
0.1 m level, because there are no mechanisms for accurate access to
that datum. (You'd need carrier phase data from the official
stations, which as far as I know they don't publish.)
The ability to compare test vectors to authoritative sources such as
NGS.
> The proj4 definitions for the projected
> coordinate system look to be the same on both systems, though since
> these are managed database services we don’t have direct access to any
> of the projection databases, grids, etc.
Hence, these questions need to be should be directed to your hosted
service provider. But once you achieve clarity, you may not need to.
A key point is whether proj has access to and is using downloadable
grids. You can't tell from the info string if they are present already.
Note you are using 8.0.1 and 9.1.0. Install both locally and see if you
can reproduce the differences. This may have nothing to do with hosted
services.
> Many thanks in advance if you can provide some advice! We would
> certainly be happy to provide any more information that might be
> helpful in describing what we are seeing.
So....
EPSG:26860 that you wrote about above is "NAD83(HARN) / Nebraska
(ftUS)". But it seems you mean EPSG:26850. This just says NAD83,
which strikes me as surprising today.
https://epsg.org/crs_4269/NAD83.html
You are concerned about 0.12 meter, but it's hard to believe you have
coordinates in the original NAD83 datum that are better than 1 m.
It's likely there is a code for your SPCS with NAD83(2011) epoch
2010.0 instead. But first the question is what frame are your data in
really. I have found that data from a GIS environment often lacks
clarity about these fine points -- doesn't have meter-level accuracy
anyway. (With MassGIS some data is sub-meter, but it, perhaps not
surprising, is recent and is labeled more carefully. An example is
15cm imagery which has been ground-controlled to NAD83(2011) epoch
2010.0).
WGS84 is an ensemble. There are a number of members, and
unfortunately one of them is the original, WGS84(TRANSIT), about 2m
different than modern members. When you use an ensemble you mean that
**your data is in one of the members but you don't know which**. In
my view, transforming **to** an ensemble is non-sensical and I
recommend that you come up with a new plan where you don't do that.
It does not make sense to transform to WGS84 and consider any
discrepancy under 2m any kind of error. Generally people want WGS84
because some standard refers to it, and my guess at what the standard
authors would have wanted had they really understood is the most
recent realization. Either that or they do not believe their standard
will be used with data that is accurate to better than a few meters.
(Which you might or might not be doing :-)
Because NAD83 is an ensemble and particularly because WGS84 is an
ensemble, and because of how proj currently works (which I think is
suboptimal but have not offered a patch to address), proj will use a
null transform from NAD83 to WGS84. If you expect your coordinates to
be in WGS84(G2159), that's wrong. You can instead transform to some
ITRF, or perhaps now you can transform to a specific WGS84. But you
have to run proj to look at the transforms and understand what ought
to be.
WGS84 realizations and ITRF realizations (that you might use as WGS84
proxies) are dynamic datums. Points have velocities (2-3 cm/y in NA),
and at the 0.12m level this cannot be neglected. Coordinates in WGS84
realizations and ITRF realizations are typically expressed at some
reference epoch, or a specific epoch. That's hard to start with, and
the reference epoch of an ensemble puts you through the thin ice, not
just on it.
> /**
> AWS RDS
> */
> -- POINT(-95.1331380928918 45.773407072930254)
FWIW, running your transform on my system (NetBSD 10 amd64, proj 9.5.1,
postgis 3.5.3), I got the same result as RDS. But, here I ask to
transform to NAD83(1986), to NAD83(2011) epoch 2010.0 (3d, but that's
ignored I think), to WGS84(ensemble) as you did, and to ITRF2014.
First, let's transform from the CRS you said the data is in to the
modern (but still ftUS :) flavor with NAD83(2011) epoch 2010.0:
gdt=> SELECT ST_AsText(ST_Transform(ST_SetSRID('POINT(2399319.7374562174 611365.0841263086)'::geometry, 26850), 6501)) as crs6501;
crs6501
---------------------------------------------
POINT(2399319.3357411437 611365.1332752574)
(1 row)
Just eyeballing, that's 45 cm right there.
Let's transform to NAD83(1986), NAD83(2011) epoch 2010.0,
WGS84(ensemble) and ITF2014.
gdt=> SELECT ST_AsText(ST_Transform(ST_SetSRID('POINT(2399319.7374562174 611365.0841263086)'::geometry, 26850), 4269)) as crs4269;
crs4269
---------------------------------------------
POINT(-95.13313652618976 45.77340681030725)
(1 row)
gdt=> SELECT ST_AsText(ST_Transform(ST_SetSRID('POINT(2399319.7374562174 611365.0841263086)'::geometry, 26850), 6319)) as crs6319;
crs6319
---------------------------------------------
POINT(-95.13313810253163 45.77340693280185)
(1 row)
gdt=> SELECT ST_AsText(ST_Transform(ST_SetSRID('POINT(2399319.7374562174 611365.0841263086)'::geometry, 26850), 4326)) as crs4326;
crs4326
---------------------------------------------
POINT(-95.1331380928918 45.773407072930254)
(1 row)
gdt=> SELECT ST_AsText(ST_Transform(ST_SetSRID('POINT(2399319.7374562174 611365.0841263086)'::geometry, 26850), 7912)) as crs7912;
crs7912
----------------------------------------------
POINT(-95.13267275410031 45.773464527082645)
(1 row)
Now the same, but relabeling the source as NAD83(2011) epoch 2010.0:
gdt=> SELECT ST_AsText(ST_Transform(ST_SetSRID('POINT(2399319.7374562174 611365.0841263086)'::geometry, 6501), 4269)) as crs4269;
crs4269
----------------------------------------------
POINT(-95.13313494984804 45.773406687812546)
(1 row)
gdt=> SELECT ST_AsText(ST_Transform(ST_SetSRID('POINT(2399319.7374562174 611365.0841263086)'::geometry, 6501), 6319)) as crs6319;
crs6319
---------------------------------------------
POINT(-95.13313652618976 45.77340681030725)
(1 row)
gdt=> SELECT ST_AsText(ST_Transform(ST_SetSRID('POINT(2399319.7374562174 611365.0841263086)'::geometry, 6501), 6319)) as crs4326;
crs4326
---------------------------------------------
POINT(-95.13313652618976 45.77340681030725)
(1 row)
gdt=> SELECT ST_AsText(ST_Transform(ST_SetSRID('POINT(2399319.7374562174 611365.0841263086)'::geometry, 6501), 7912)) as crs7912;
crs7912
--------------------------------------------
POINT(-95.13267117775811 45.7734644045764)
(1 row)
Use NCAT or other NGS tools to transform your point to NAD83(2011) epoch
2010.0. Also do it in proj. You might need to use proj to transform to
the base CRS of your SPC, EPSG:4269. That's just a projection, and not
a datum shift. You might see how different "NAD83(1986)" is from
"NAD83(2011) epoch 2010.0".
https://www.ngs.noaa.gov/NCAT/
Then, you can use various tools to transform that to WGS84(G2296) or
WGS84(G2139). NGS doesn't deal as much in WGS84. As I understand it
G2296 is equivalent to ITRF2020 and G2139 is equivalent to ITRF2014.
This is probably useful
https://www.ngs.noaa.gov/cgi-bin/HTDP/htdp.prl?f1=4&f2=1
This is old, but it gives you some sense of the difficulties.
https://www.ngs.noaa.gov/CORS/Articles/WGS84NAD83.pdf
> /**
> AWS Aurora.
> */
> -- POINT(-95.13313652618976 45.77340681030725)
This matches my transform to NAD83(1986), which is a projection and not
a datum shift.
proj 9 had major improvements in US national datum transforms, and seems
to be transforming between NAD83(1986) and NAD83(2011) epoch 2010.0. My
impression is that it is close to what NCAT produces. Such transforms
are almost certainly gridded as they are really about distortions in the
early frames from errors, greatly reduced later.
What's surprising me and I need to dig in more is that it appears the
transformation from 6319 to 4326 is not null, but it's also quite
different from 7912 (ITRF2014).
This may raise more questions than answers, but you really need to drive
definitional uncertainty to near zero if you want to talk about 12 cm.
Please report back with a summary of your understanding after
investigating, as I suspect it is helpful to others, as well as in
understanding the current state of proj. Dealing with ensembles and
with older datums is hard.
Greg
More information about the postgis-users
mailing list