Switch to EPSG:7912 for GeoDataFrame

[jpswinski]

The SlideRule Python client currently uses EPSG:4326 for the GeoDataFrames it returns.

Consider using EPSG:7912. David: "They are effectively the same, but the 7912 tag is more descriptive of the datum/realization/epoch."

[scottyhq]

Just wanted to add a little more context here. EPSG:4326 is an "Ellipsoidal 2D" coordinate reference with respect to WGS84 datum ensemble, and a lot of software today anticipates this for geolocated data (lat,lon,value). In particular, geodataframes are often serialized to geojson format, which requires EPSG:4326 https://datatracker.ietf.org/doc/html/rfc7946 and consequently visualization libraries assume this is the coordinate reference. So there may be an advantage in EPSG:4326 as default when it comes to quick visualizations and exploration of IS2 data.

For geodetic applications of IS2 where it's desirable to strive for sub-meter accuracy and compare to other datasets (e.g. digital elevation models like copernicus dem, lidar, etc), CRS conversions should be as accurate as possible. So it's possible to use "Ellipsoidal 3D" CRS definitions. Instead of changing the default value, maybe we add documentation and an examples of various CRS choices?

The data provider, NSIDC states "This data set (ATL03) contains height above the WGS 84 ellipsoid (ITRF2014 reference frame)". A great reference on ITRF is here https://qpssoftware.scrollhelp.site/qinsy/9.5/international-terrestrial-reference-frame-2014-itr

as @tsutterley pointed out EPSG:7912 actually uses the GRS80 ellipsoid. This is likely within millimeters of using the most recent WGS84 ellipsoid realization EPSG:7665

[tsutterley]

@scottyhq and I were talking about this issue this week. If anything we should use EPSG:7661 as ICESat-2 is in reference to WGS84 G1150 (at least according to the ATL03 ATBD)

[jpswinski]

Still need to decide which transformation is best for ArcticDEM and REMA since it uses the ensemble WGS84 datum.

We are basically choosing between the following 2 options:

(2D)
PROJ string:
+proj=pipeline
  +step +proj=axisswap +order=2,1
  +step +proj=unitconvert +xy_in=deg +xy_out=rad
  +step +proj=stere +lat_0=90 +lat_ts=70 +lon_0=-45 +x_0=0 +y_0=0 +ellps=WGS8

and

(3D)
PROJ string:
+proj=pipeline
  +step +proj=axisswap +order=2,1
  +step +proj=unitconvert +xy_in=deg +z_in=m +xy_out=rad +z_out=m
  +step +proj=cart +ellps=GRS80
  +step +proj=helmert +x=1.0053 +y=-1.9092 +z=-0.5416 +rx=0.0267814
        +ry=-0.0004203 +rz=0.0109321 +s=0.00037 +dx=0.0008 +dy=-0.0006
        +dz=-0.0014 +drx=6.67e-05 +dry=-0.0007574 +drz=-5.13e-05 +ds=-7e-05
        +t_epoch=2010 +convention=coordinate_frame
  +step +inv +proj=cart +ellps=GRS80
  +step +proj=stere +lat_0=90 +lat_ts=70 +lon_0=-45 +x_0=0 +y_0=0 +ellps=WGS84

[jpswinski]

Code is currently using EPSG:7912

[scottyhq]

@scottyhq and I were talking about this issue this week. If anything we should use EPSG:7661 as ICESat-2 is in reference to WGS84 G1150 (at least according to the ATL03 ATBD)

@tsutterley can you link to the doc / excerpt of where this is documented?

NSIDC docs still specify ensemble epsg:4326 unfortunately :( (https://nsidc.org/data/atl03/versions/6)

[tsutterley]

from the ATBD for ATL03:

The geolocation of ATLAS photons provides the starting point for the ATL03 science level product. Each photon event will have been placed within a geodetic coordinate system; the elevations are given in the ITRF2014 reference frame and the geographic coordinates (latitude, longitude, and height) are referenced to the WGS-84 ellipsoid based on the G1150 model (ae = 6378137m, 1/f = 298.257223563).

and Table 4.1 in the altimetry comparison guide provided by ICESat-2

Mission	Ellipsoid	Semi-major axis (a)	Inverse flattening (1/f)
ICESat-2	WGS-84 (G1150)	6378137.0 m	298.257223563

[jpswinski]

Does merging #272 resolve this issue?

[dshean]

I think it may be good enough for now. As @tsutterley points out, the EPSG:7912 CRS isn't exactly what ICESat-2 uses. I don't know of an existing EPSG code for ITRF 2014 and WGS84 ellipsoid (G1150) and 2010.0 epoch. I think we may need to define all of this with WKT and use that to set the CRS of our returned geodataframe. We can also apply for a new code to use for this, which I think is justifiable as the CRS for the ICESat-2 mission.