UTM coordinates in NAD83 differ from WGS84

[sfroyen]

Steps to reproduce

1. Import Lidar generated contours and slope templates in UTM NAD83(2011) or more recent.
2. Import georeferenced templates in lat/lon WGS84 (e.g., Google satellite images) or use an accurate GPS with OOM on Android.
3. Observe a systematic difference between the Lidar and the GPS positions.

Actual behaviour

NAD83 and WGS84 in recent realizations no longer coincide. However, most lat/lon to UTM coordinate conversion software, including OOM, assumes they do. The error in Colorado is around 2 m, enough to be noticeable (at least to the mapper) on sprint maps.

Expected behaviour

OOM should perform lat/lon to UTM conversion correctly. Although it sounds like PROJ-6 should be able to handle this, it may be simpler to use a hack. For any orienteering map, the shift is a constant that can be implemented by a single delta-x delta-y pair. Input fields for this can added to the map georeferencing screen (in lat/lon degrees, UTM meters or, possibly, even map coordinates mm). There are web sites that will the provide the shift (I'll follow up with one later) and survey precision GPS units usually come with software that does the same.

In principle, there could be data that come in UTM WGS84. I do not know if we need to handle such cases.

Configuration

Mapper Version: latest (Edit: v0.8.4 and/or dev v20190610.2) Operating System: all (Edit: not verified)

[dg0yt]

Can you provide correct pairs of geographic coordinates (WGS84) and projected coordinates (UTM) which are handled incorrectly by Mapper?

---

There is no doubt that "NAD83 and WGS84 in recent realizations no longer coincide." It is for a good reason:
"NAD 83 is defined to remain constant over time for points on the North American Plate, whereas WGS 84 is defined with respect to the average of stations all over the world. Thus the two systems naturally diverge over time." [Wikipedia on NAD].

IIUC your key claims is that:
Mapper incorrectly transforms Lat/Lon from WGS84 (?) to UTM, and an extra x/y shift would be needed for correction.

I would rather expect this is what can be said for NAD83, not UTM. There are "grid shift" files which provide the correction for NAD83. I never tried to make this work with Mapper. But even without these accurate files, it should be possible to to include a local shift adjustment to the NAD83 PROJ.4 specification.

Most parts of the world are not on the North American Plate. As far as I know UTM is generally based on WGS84. But it is hard to find a formal reference, and Wikipedia concludes the history of UTM with: "The WGS84 ellipsoid is now generally used to model the Earth in the UTM coordinate system, which means current UTM northing at a given point can differ up to 200 meters from the old. For different geographic regions, other datum systems (e.g.: ED50, NAD83) can be used." [Wikipedia on UTM]

To add another complication, even WGS84 has undergone revisions.

---

Regarding PROJ version 6 and up, Mapper (and perhaps exisiting maps files) might need extra modifications. One prominent change is the removal of the defaults (file proj_def.dat).

[krticka]

I ran into similar problem when mapping at Yukon, where most of the geodata are using Yukon Albers Projection on GRS80 ellipsoid. The shift with GPS measurements were approximately up to 2 meters (the difference can regionally vary over the continent). However such systematic positional mismatch with precise geodata is not pleasant, I think it is something you can simply live with as shift of 2 m in orienteering mapping is negligible.

Not sure if redefining of +x_0= +y_0= parameters in PROJ4 can help you. I was trying to search for parameters which allows to define local offset, but with no luck. You also need to know approximate direction of the shift. For example OCAD allows to define additional local offset. Personally I used it only once or twice in the past.

[lpechacek]

As to my understanding, which I acquired in u-blox GPS compendium, there are local (national) coordinate systems and global ones. As a matter of principle, the local reference systems, anchored to some place on a continent, will diverge over time from the global reference systems due to continental drift. Moreover, there is reportedly significant tectonic activity in some places which asks for local adjustments. That's at least OCAD's explanation for the Additional local offset field. Local offset setting appears like a useful band aid to me. Just my 0,02 €.

[ghost]

That's at least OCAD's explanation for the Additional local offset field.

Name it "Offset of CRS center" or "Local CRS center" in Mapper (if such feature will be added).

[dg0yt]

Assuming that both the LIDAR data and the GPS data are correct with respect to their CRS and point in time, but there is a mismatch, the questions are:

• Was the CRS correctly specified for the imported data?
• For projected data not using the map's projected CRS: Was the WGS84 transformation given correctly and precisely enough?

If that is too hard, and we want to add a custom offset: We must deal with the manual offset for templates and for data import.

BTW for georeferenced raster templates, it is already possible in dev to switch off georeferencing, and then adjust the position manually.

[sfroyen]

Here's the link to the Web site I was referring to in my initial posting

http://tagis.dep.wv.gov/convert/

As an example, a conversion between lat/lon WGS84 and lan/lon NAD83 produces the output

40.0000000,-105.0000000,LL WGS84(G1150) 39.9999935,-104.9999899,LL NAD83(CORS96)

This shift nicely matches the "error" I observe in my data.

WGS(G1150) and NAD83(CORS96) were used ca. 2000-2010. There are additional changes since then but those are of the order of centimeters.

Note that there is no "simple" transformation that corrects for these smaller shifts as they depend strongly on the local geographic position.

[dg0yt]

Mapper's UTM is linked to WGS84. If you first setup the georefencing for 40N 105E with UTM Zone 13 N, and then switch from "UTM" to "Custom PROJ.4" CRS, the specification will be displayed as:

+proj=utm +datum=WGS84 +zone=13

According to spatialreference.org, "NAD83 / UTM zone 13N" is EPSG:26913, which is defined as

+proj=utm +zone=13 +datum=NAD83 +units=m +no_defs 

in PROJ.4. Spatialreference.org adds +ellps=GRS80.

But even with the NAD83 UTM, Mapper won't display other UTM coordinates than with WGS84 UTM.

With PROJ version 6 in the unstable dev builds of Mapper, there are many ways to do transformations, including simple affine transformation by +xoff=<value> and +yoff=<value>.

What should work with all versions is the +towgs84 parameter. There is a seven parameters form for which I found this old documentation: (http://proj.maptools.org/gen_parms.html)

The seven parameter case uses delta_x, delta_y, delta_z, Rx - rotation X, Ry - rotation Y, Rz - rotation Z, M_BF - Scaling. The three translation parameters are in meters as in the three parameter case. The rotational parameters are in seconds of arc. The scaling is apparently the scale change in parts per million.

So +towgs84=0,0,0,0,0,0,0 should be neutral, and +towgs84=4,-2,0,0,0,0,0 should apply a small x/y shift.

Putting it all together, something like

+proj=utm +zone=13 +datum=NAD83 +units=m +no_defs +ellps=GRS80 +towgs84=4,-2,0,0,0,0,0

might establish the expected local UTM coordinates with the current version of Mapper.

[pkturner]

But even with the NAD83 UTM, Mapper won't display other UTM coordinates than with WGS84 UTM.

Mapper does not accurately convert EPSG:26913 to WGS84, because PROJ conflates WGS84 with NAD83 (when interpreting a spec string). See https://proj.org/development/migration.html. I checked this with some Colorado coordinates. If the map's CRS has an NAD83 datum, then Mapper displays the NAD83 lat-long, slightly at variance with WGS84.

[pkturner]

As an example, a conversion between lat/lon WGS84 and lan/lon NAD83 produces the output

40.0000000,-105.0000000,LL WGS84(G1150) 39.9999935,-104.9999899,LL NAD83(CORS96)

To add detail, the difference between these is 1.1 meter. That's usually negligible, but it can annoy a mapper.

[pkturner]

There are "grid shift" files which provide the correction for NAD83.

Is there a way that Mapper could become aware of which CRSes have NAD83 as datum? I would hope that the needed grid shift correction (whether performed in full or by local approximation) would not depend on a special option in Mapper's Georeferencing dialog.

[dg0yt]

See https://proj.org/development/migration.html

@pkturner Thanks, this is a reference I remembered and was trying to find yesterday, without success. The most interesting point is that "the view that “WGS84 is the true foundation of the world, and everything else can be transformed natively to and from WGS84” is inherently flawed." This is what Mapper is doing at the moment. (At least we are also told that "it is not inherently bad".)

If the map's CRS has an NAD83 datum, then Mapper displays the NAD83 lat-long, slightly at variance with WGS84.

This would explain why I had to find out that "even with the NAD83 UTM, Mapper won't display other UTM coordinates than with WGS84 UTM." But at the moment I fail to see why this is the case - except for that there is no general precise mapping of coordinates between NAD83 and WGS84?

Is there a way that Mapper could become aware of which CRSes have NAD83 as datum?

I think that the PROJ specifications would have a +datum=NAD83. I'm not sure of sure if its possible to get this via PROJ API from a CRS specified with EPSG code, or if there is only a fully resolved specification string.

I would hope that the needed grid shift correction (whether performed in full or by local approximation) would not depend on a special option in Mapper's Georeferencing dialog.

The PROJ specification is actually a generic field, useful for the experts who know how it works. Most of us are not, me included. But before we add additional special fields, we should know if we can inject these inputs into the PROJ projection which links geographic and projected coordinates, or if we really need to add another level of coordinate manipulation. (And note that both GDAL and PROJ support other CRS specification formats, at least WKT. Maybe we want this in Mapper, too?)

Maybe I should ask like this:

1. What else do we need to tell PROJ in order to turn a coordinate received live from GPS in a proper projected coordinated for your desired CRS?
2. How I can present this in the user interface?

[sfroyen]

I've been trying to track down the source of the changes and it appears the both WGS84 and NAD83 coordinates change over time. For instance, https://confluence.qps.nl/qinsy/en/world-geodetic-system-1984-wgs84-29855173.html states that WGS84 updates are such that "Its time evolution in orientation will create no residual global rotation with regards to the crust". It is based on GPS measurements so I assume that an average over the GPS measurements defines the residual. The link above lists (center of the earth?) shifts of almost 1 m from original WGS84 to its 2013 realization. Similarly, NAD83 coordinates have moved w.r.t. the original NAD83 coordinates, presumably using a North American based average.

In order to import data into an existing map, we need to know the coordinate system that was used when the map was created, e.g., UTM NAD83(2011), as well as the coordinate system of the new data. I'm guessing that most new data will be in current or at least recent GPS coordinates. As a default, it would therefore make most sense to shift the GPS coordinates to the original map coordinates. There would also be a need to override the default shift for templates and data that are provided in other (older) coordinates.

It is certainly possible to do this using the appropriate coordinate transformations, which may be supported by the latest versions of PROJ. However, It still seems to me that using a simple (x, y) shift is just as easy. (We would also need a text field to indicate the original map coordinate system.) It is somewhat akin to how the declination is treated.

[dg0yt]

In order to import data into an existing map, we need to know the coordinate system that was used when the map was created, e.g., UTM NAD83(2011), as well as the coordinate system of the new data.

Isn't this what Mapper is doing? We attach the coordinate reference system to each map and template.

to shift the GPS coordinates to the original map coordinates

This is what Mapper is supposed to do for template data, including GPS.

AFAIU the missing piece is providing the details for the desired shift. It does not come magically, by using the correct CRS specifications, but needs extra information, e.g. grid shift files produced from HDTP, or other forms of transformation. We have mentions of a simple X/Y shift here but obviously the general case might require more parameters for some levels of precision.

(We would also need a text field to indicate the original map coordinate system.)

I do not understand the purpose of this proposal. We use formal specifications of the coordinate systems, in order to allow for projecting/transforming coordinates. What is the added value of the text field?

It is somewhat akin to how the declination is treated.

Declination is not a free text field, but another formal parameter which is orthogonal to the coordinate system.

[sfroyen]

You are right. My apprehension is probably caused by my frustration with PROJ. In order to fix the conversion issue in OOM, I must be able to make PROJ reproduce the shifts. I've been struggling with PROJ 6 to do so but I always end up with no shift. Has anyone else succeeded in this? Here are some of the datums that I've tried

NAD83(HARN) UTM zone 13N EPSG:3743 (with transformation: 15930) NAD83(NSRS2007) UTM zone 13N EPSG:3720 (with transformation: 15931) NAD83(2011) UTM zone 13N EPSG:6342

WGS84(TRANSIT) EPSG:7816 (m) EPSG:7815 (ll) WGS84(G730) EPSG:7656 (m) EPSG:7657 (ll) WGS84(G873) EPSG:7658 (m) EPSG:7659 (ll) WGS84(G1150) EPSG:7660 (m) EPSG:7661 (ll) WGS84(G1674) EPSG:7662 (m) EPSG:7663 (ll) WGS84(G1762) EPSG:7664 (m) EPSG:7665 (ll)

WGS84 UTM zone 13N EPSG:32613

I found this comment on https://confluence.qps.nl/qinsy/en/world-geodetic-system-1984-wgs84-29855173.html

Note that the EPSG database contains no specific WGS84 datum realizations.

Which may be part of my issue.

Once I get this resolved, I should be able to provide the necessary information.

[dg0yt]

Mapper (dev) now use the new PROJ 6 programming interface which has support for (time-dependent) "dynamic CRS" and "late binding". However Mapper doesn't use late binding at the moment.

[sfroyen]

A couple of comments:

This link, https://proj.org/operations/transformations/helmert.html , documents PROJ Helmert transformations. The 4-parameter, 2D version would solve my immediate issue.

OOM could, of course, also implement 3D and 4D transformations (7 or 14 parameters). The trouble with these is obtaining the parameter values (in the US there are local variations).

If we use an interface which includes the transformation type as well as the parameter values we can leave the choice to the user.

Note that these transformations are between reference frames. We would need to link templates, map parts(?), GPS input, etc. to their respective reference frames and then provide transformation data between the frames. If, as seems reasonable, we restrict all map objects (but not templates) to one common reference frame, all other items can be tagged with a single transformation.

We need to able to update the parameters as they change over time.

[dg0yt]

1. Mapper doesn't implement a single transformation: All is done via PROJ. In a number of places, e.g. map georeferencing, you may use any PROJ CRS specification you want. However, WKT2 might be a more appropriate general specification language in the future, when EPSG codes are not enough.
2. PROJ just went away from the now-old-school single reference frame (WGS84), and Mapper will probably to do so as well. PROJ does the job of providing the most accurate transformations, depending on region and time. Mapper may need to add the third and fourth dimension in order to make full use of it.
3. With proper test data, I can really help track down errors. For example, with https://github.com/OpenOrienteering/mapper/issues/1570, I helped find another bug in PROJ code and GDAL documentation. It might even be related to your issue, so please try the latest master build...

[dg0yt]

@pkturner Maybe you can find out if NAD83 situation improved with last master/#1570 changes?

[sfroyen]

I was finally able to make proj (really cs2cs) provide sensible results. The command "cs2cs EPSG:7665 EPSG:6342", converting from WGS(G1762) to NAD83(2011) UTM 13N, produced numbers that were reasonable and different from NAD83(HARN). So, I switched the map georeferencing from UTM 13N to EPSG 6342 and reimported the contours. I no longer see differences between the Lidar contours and GPS positions :-)

Lidar based templates are correctly positioned too.

In principle, there is still a problem with the different WGS84 realisations. (It looks like GPS receivers report updated positions when a new realisation is introduced.) Since the first realisation, positions have changed as much as 1 m. However, the adjustments for the last couple of realisations are only 1 cm, so probably not a pressing issue.

This is with OOM 0.9.2. I will try the new build later.

[ghost]

This is with OOM 0.9.2. I will try the new build later.

@sfroyen If you use Linux, you could grab latest master build from here:

• https://software.opensuse.org/download.html?project=home%3Adg0yt&package=openorienteering-mapper-master

You could install it side-by-side on your PC with v0.9.2 release.

[pkturner]

@dg0yt I will look at it on 04-30.

[pkturner]

After a small amount of experimenting, the use of EPSG:7665 found by sfroyen looks promising as a high accuracy replacement for +proj=latlong +datum=WGS84 that Mapper now uses for its geographic CRS.

EPSG:7665 represents a very high accuracy, recently updated version of WGS84, and PROJ 6 knows how to convert between it and NAD83(2011) which is the datum used by UTM 13N as needed by sfroyen. (It also happens to be the basis for the LiDAR I'm using in Tennessee, though 1 m accuracy is not a concern for me.)

My experiments show that when converting to/from EPSG:4326 (equivalent to Mapper's geographic CRS), PROJ does not perform the transformation between NAD83(2011) and WGS84(G1762), but that when converting to/from EPSG:7665, PROJ does insert the transformation when appropriate.

As I understand it, NAD83(2011) is a satellite-based datum for the North American plate, so the transformation between it and recent realizations of WGS84 is represented as a 7 or 14 parameter Helmert transformation rather than require access to a detailed grid transformation file. It's an advantage if Mapper does not need to be packaged with large data files and does not need to access them on the fly over the Internet.

It would be good to try building Mapper with EPSG:7665 as its geographic CRS, and to see what it does with a CRS that is defined using a grid transformation.

[pkturner]

Maybe you can find out if NAD83 situation improved with last master/#1570 changes?

I've begun to think of how to create a test for that. Best would be for the test to start with a map that has a CRS based on NAD83(2011), then import or add a template GPX track and check the projected coordinates of the track.

Based on what I saw with cs2cs, I expect the test to fail with master.

[pkturner]

The key to persuading PROJ to perform the transformation between NAD83 and WGS84 may be to change the handling of GPX tracks, etc. Currently they are assigned a basic WGS84 transformation which PROJ takes to be the rough, or imprecise original WGS84. If they were given a geographic CRS of EPSG:7665 then PROJ should insert the datum transformation.

[dg0yt]

For testing, I would suggest extending GeoreferencingTest::testProjection_data first.

I do not fully understand the effect of generally switching to EPSG:7665. I see that it is a 3D CRS which uses the G1765 realization of the WGS84 datum and the WGS84 ellipsoid. EPSG 4326 is described as 2D CRS. Or do we want to select the WGS84 realization based on record date?

Maybe we should start with using "EPGS:4326" instead of "+proj=latlong +datum=WGS84" if this helps PROJ to find a more accurate transformation? There might be a few more spots than Georeferencing which would need be touched then.

But as noted in #1598, we might want to get away from binding to a single geographic CRS at all.

[pkturner]

I wasn't paying attention to 3D CRS vs. 2D. EPSG 9057 is a 2D CRS that's based on WGS84(G1762).

1598 is good, but it doesn't seem to bring us closer to resolving the present issue.

Here's what cs2cs gives:

$ cs2cs EPSG:4326 EPSG:6342
40 -105
500000.00       4427757.22 0.00
$ cs2cs EPSG:7665 EPSG:6342
40 -105
500000.86       4427756.50 0.89
$ cs2cs EPSG:9057 EPSG:6342
40 -105
500000.86       4427756.50 0.89

Only the latter two include the desired datum transformation. PROJ doesn't handle EPSG:4326 as a datum which merits a careful transformation to NAD83. Actually, record date does come into it, because 30 years ago NAD83 and WGS84 were not that different.

[dg0yt]

You may also list the full pipeline:

projinfo -s EPSG:9057 -t EPSG:6342

which results in:

[...] Conversion from WGS 84 (G1762) (geog2D) to WGS 84 (G1762) (geocentric) + WGS 84 (G1762) to ITRF2008 (1) + ITRF2008 to NAD83(2011) (1) + Conversion from NAD83(2011) (geocentric) to NAD83(2011) (geog2D) + UTM zone 13N, 0.01 m, USA - CONUS and Alaska; PRVI

PROJ string: +proj=pipeline +step +proj=axisswap +order=2,1 +step +proj=unitconvert +xy_in=deg +xy_out=rad +step +proj=cart +ellps=WGS84 +step +proj=helmert +x=0.99343 +y=-1.90331 +z=-0.52655 +rx=0.02591467 +ry=0.00942645 +rz=0.01159935 +s=0.00171504 +dx=0.00079 +dy=-0.0006 +dz=-0.00134 +drx=6.667e-05 +dry=-0.00075744 +drz=-5.133e-05 +ds=-0.00010201 +t_epoch=1997 +convention=coordinate_frame +step +inv +proj=cart +ellps=GRS80 +step +proj=utm +zone=13 +ellps=GRS80

WKT2_2018 string: [...]

[dg0yt]

Edit: source CRS was EPSG:9057

[dg0yt]

IIUC what we need is to attach the exact realization of WGS84 to GNSS data / GPX tracks. We can try that when I finish my current work of adding choice for vector data CRS (needed for GH-1561). And consequently we need to deal with this on coordinate display and vector data export.

[pkturner]

Going a little deeper, cs2cs does

 $ cs2cs EPSG:9057 EPSG:6342
 40 -105
 500000.86       4427756.50 0.89
 40 -105 0 1997.0
 500000.86       4427756.50 0.89 1997.0
 40 -105 0 2000.0
 500001.24       4427756.65 0.86 2020.0

Note that if the epoch is not specified, then the datum transformation is performed for the reference epoch of 1997, yielding an adjustment that is 30 or 40% different from what is correct for 2020.

For testing, I would suggest extending GeoreferencingTest::testProjection_data first.

Enhancements to GeoreferencingTest::testProjection are desirable.

• Accept a parameter indicating the epoch.
• Require greater accuracy. To check that the epoch is handled correctly, accuracy needs to be checked to within about 0.3m (vs. the current 2.2m).

I guess the details of enhancing the test can't be worked out until changes to the Georeferencing class are better known. I think Georeferencing will expose the epoch so that a TemplateTrack can be tied to a CRS combined with a time. Testing this would compare coordinates in two different projected coordinate systems, because what matters is the transformation between the template coordinate system and the map's projected coordinate system.

[dg0yt]

I also did cs2cs tests with the 4th dimension, but it didn't improve the situation with EPSG:4326 ;-)

At the moment, we already HUGE_VAL instead of 0 for the epoch, following feedback in one of our PROJ bug reports.

The test accuracy still provides for Linux builds with PROJ 4.9, but this could be solved by #ifdefs, or by using another data field.

[sfroyen]

Here's a partial list of EPSG numbers related to WGS 84:

1994 G730 EPSG: 7656 (x,y,z), 7657 (lat,lon,h), 9053 2D (lat,lon) 1997 G873 EPSG: 7658 (x,y,z), 7659 (lat,lon,h), 9054 2D (lat,lon) 2002 G1150 EPSG: 7660 (x,y,z), 7661 (lat,lon,h), 9055 2D (lat,lon) 2012 G1674 EPSG: 7662 (x,y,z), 7663 (lat,lon,h), 9056 2D (lat,lon) 2013 G1762 EPSG: 7664 (x,y,z), 7665 (lat,lon,h), 9057 2D (lat,lon)

The first three appear to be identical. Same with the last two. (Checked using cs2cs and projinfo).

I spoke too soon when I said that opening templates worked correctly. Templates using lat/lon are positioned incorrectly. It also makes no difference which EPSG I use. In particular, EPSG 7657 and 7665 produces exactly the same template position even though cs2cs says they should be different:

echo 40N 105W 0 | cs2cs EPSG:7657 EPSG:6342 500000.00 4427757.22 0.00 echo 40N 105W 0 | cs2cs EPSG:7665 EPSG:6342 500000.86 4427756.50 0.89

The first results, which is shifted NW compared to he second, appears to be what I need. I also note that EPSG:4326 produces the "incorrect" result

echo 40N 105W 0 | cs2cs EPSG:4326 EPSG:6342 500000.86 4427756.50 0.00

My georeferenced images are from Pictometry and I do not know what datum their coordinates use.

I will a look at the Mapper code later, but if somebody beats me to it I will not be unhappy :-)

Mapper version v20200502.2

[ghost]

Mapper version v20200502.2

@sfroyen As there no such release builds, it would be better to reference latest commit ID used for self compile Mapper.

[dg0yt]

Mapper version v20200502.2

That's the CI build number. It is a perfect version reference.

However, it only fixes issue 1597, not this one. Here, we need to touch and test several places in the source code. I don't want to make changes in a hurry without understanding the impact - at the moment, every transformation is depending on our default geographic CRS, and there will be people relying on the current transformations.

[ghost]

That's the CI build number

Wow! Where those nightly builds available for download?

[dg0yt]

Wow! Where those nightly builds available for download?

It is off-topic here. Those are not nightly builds. And there is intention in not announcing them as a general service.

[dg0yt]

For GPX files, GDAL reports the following CRS:

GEOGCS["WGS 84",
    DATUM["WGS_1984",
        SPHEROID["WGS 84",6378137,298.257223563,
            AUTHORITY["EPSG","7030"]],
        AUTHORITY["EPSG","6326"]],
    PRIMEM["Greenwich",0,
        AUTHORITY["EPSG","8901"]],
    UNIT["degree",0.0174532925199433,
        AUTHORITY["EPSG","9122"]],
    AXIS["Latitude",NORTH],
    AXIS["Longitude",EAST],
    AUTHORITY["EPSG","4326"]]

which means: EPSG:4326.

[sfroyen]

I've been running some tests with a Ublox F9P based GPS (using their C099-F9P board). It supports RTK corrections over the Internet. I've discovered that GPS coordinates for a given location depends on how you use your GPS:

-Uncorrected GPS reports lat/lon coordinates using WGS84 in the original(?) realization -SBAS corrected GPS uses WGS84 in a recent (less than 6 months old?) realization. -RTK corrected GPS uses whatever your base station is set up to use. In the US, this is frequently NAD83.

Here're coordinates for a point in my neighborhood:

F9P GPS w/RTK corrected measurements (GPS antenna on top of my head) 39.67797133 -105.16596350 WGS84(current) base station 39.67796617 -105.16595033 NAD83(2011) base station/network delta: 5.16E-6 -13.17E-6

Conversion from the WGS84 value using https://tagis.dep.wv.gov/convert/ 39.6779713 -105.1659635 wgs84 (from dropdown) 39.6779649 -105.1659534 nad83 (to dropdown) delta: 6.4E-6 -10.1E-6

Same conversion using PROJ with the command echo 39.6779713N 105.1659635W \ | cs2cs -f "%.8f" EPSG:7665 EPSG:6342 \ | cs2cs -d 7 EPSG:6342 EPSG:7657 produces an identical result. Transformation is WGS84(G1762) -> NAD83(2011) UTM 13N -> WGS84(G730)

Note that skipping the intermediate NAD83(2011), using instead echo 39.6779713N 105.1659635W | cs2cs -d 7 EPSG:7665 EPSG:7657 produces a delta equal to zero!

With OOM set up to use CRS=EPSG:6342 (NAD83 version above), I find that GPS positions agree with contour details only if I use an SBAS or an RTK w/WGS84 base station. If I use RTK w/NAD83 base, I can tell they are off.

Comparing the EPSG and the geographic positions in OOM's georeferencing panel, it looks like OOM is using the latest WGS84 realization (EPSG:7665 or the generic EPSG:4326). SBAS corrected GPS data can therefore be used directly. In principle, uncorrected GPS signals need to be modified but the intrinsic error is probably too large to worry about it.

For RTK corrected GPS there should be an option to transform them.

I'm in dialog with Pictometry about the georeferencing data they supply and I'll add a note here once I have an answer. I'm in good position to test any new code for template coordinate transformations.

[dg0yt]

Thanks for the update. I completed some of the more fundamental changes to our template code (now in the master release), so I should be able to continue with georeferencing options for geospatial vector data.

[lpechacek]

• RTK corrected GPS uses whatever your base station is set up to use. In the US, this is frequently NAD83.

Fun ahead. Appendix B of ZED-F9P Integration manual confirms that:

... The International Terrestrial Reference Frame (ITRF) must be obtained from the reference network and then used to transform the rover position output to match the required reference frame. The rover does not output the position in the local receiver WGS84 (based onITRF2008) datum, it matches the reference receiver (or base) reference frame. The user application needs to do the transformation for use in a mapping application if it does not use the same reference frame. ...

Good news is that it is likely only ZED F9P which behaves this way. The earlier M8P (RTK) product does not exhibit this behavior, according to its documentation.

[pkturner]

There's a change to the behavior of Mapper which is showing up when built with proj-7.0.1. I believe the change arrived in proj-6.3.0.

I'm seeing it when experimenting with georeferencing dialogs on maps that use EPSG:6575 (a projection with an NAD83 datum). The PROJ transformation between the projected and geographic coordinates has begun to include a conversion between the datums. Here in Tennessee the difference is about 0.5 meters.

The difference can be seen in maps that had their georeferencing set with an earlier version of PROJ. In the georeferencing dialog if you change one of the geographic coordinates, and then change the coordinate back to its original value, the result is a small but relevant change to the projected coordinates. I expect there will be a few other effects, perhaps relating to templates.

This change has the potential to improve the issue originally reported by @sfroyen, because this behavior change affects EPSG:6342, which is UTM Zone 13N based on NAD83(2011). In Colorado the difference is more than 1 meter.

[dg0yt]

The key to exploring the transformation pipeline is the projinfo command, e.g.:

projinfo -s EPSG:9057 -t EPSG:6342

And this pipeline is going to improve with new PROJ versions and/or availability of data.

On the other hand, we still need to deal with assigning other CRS than "WGS 84" to geospatial vector data such as GPX, especially when then data doesn't tell us correctly.

[pkturner]

`projinfo -s EPSG:6575 -t EPSG:6342`  

was not helpful at first, because it shows a transformation pipeline that matches the behavior of Mapper prior to proj-6.3. Today I tried

`projinfo -s EPSG:6575 -t "+proj=longlat +datum=WGS84 +type=crs"`  

which yields a pipeline corresponding to the current behavior of Mapper (built with proj-7.0.1). It inserts an inverted "WGS1984(ITRF08)_To_NAD_1983_2011" step that consists of a 7-parameter Helmert transformation, whereas the former has a step called "Ballpark geographic offset from NAD83(2011) to WGS 84" which is a no-op.

My GPS isn't accurate enough for me to bother about 0.85 meters in a GPX track -- I use those tracks for a day and then discard them. The thought that concerns me a little is that my maps have LiDAR-based and orthophoto templates kept for reference. Now that Mapper has changed from showing NAD83 geographic coordinates for my maps, to showing WGS84, I wonder whether that change will affect the alignment of these templates.

[dg0yt]

Looks like there are also others struggling with NAD, GPS, and PROJ: https://lists.osgeo.org/pipermail/proj/2020-October/date.html, in particular

• EPSG v10 update status, https://lists.osgeo.org/pipermail/proj/2020-October/009847.html
• Using latest realization of a datum ensemble ? https://lists.osgeo.org/pipermail/proj/2020-October/009863.html

and following dicussion.

[pkturner]

For testing, I would suggest extending GeoreferencingTest::testProjection_data first.

Well, I went ahead with my original idea of importing a GPX template, and checking its position, in PR #1757. That kind of test is more like the use case raised by this issue.

GeoreferencingTest::testProjection checks projection from/to WGS 84 coordinates. The discussion above demonstrates that WGS 84 is hard to pin down. Also, the conversion to EPSG:6342 depends on the date of GPS operation, which is not a parameter of this test. It seems that testProjection is not suited to situations that involve a change of datums.

I'm hopeful that Mapper can address the present issue and update its use of PROJ, dropping the assumption that every Georeferencing must relate to WGS 84.

[dg0yt]

I'm hopeful that Mapper can address the present issue and update its use of PROJ, dropping the assumption that every Georeferencing must relate to WGS 84.

I'm willing to drop that. See draft #1598, and also taken into consideration by the changes in #1753. #1753 is factored out of my work for making the data CRS configurable for vector data - which would help with aligning GPX data on NAD 83 maps. The data in your PR helps, too.

I just don't want to break too much. PROJ and GDAL still have there own issues. #1755 just showed jow easily we can miss necessary changes - and once more, this also lead to documentation updates in GDAL.