We are using pyproj to transform between coordinate reference systems. What is a good way to check if the transformation is simply scale/translate/units? Basically, I want to know if a grid in the first crs will still be in a grid in the target crs.
Edit: I want to add some clarity about what I am trying to detect. Suppose I have a 3x4 grid of 12 points represented by 3 x
values and 4 y
values. When I transform these points, sometimes the grid is distorted such that 12 x,y
pairs are required. But sometimes, the transformation yields a 3x4 grid that can still be represented using 3 x
values and 4 y
values. It is these latter cases that we would like to detect.
An example is EPSG:4236
to EPSG:3395
. Other transforms distort the grid, such as EPSG:4236
to EPSG:2193
.
I've listed out a lot of detail below, but I think this is going to be a pretty easy question. Ideally, there is something simple in the pyproj.Transformer
object that I can inspect, or an easy way to compare projection definitions in some format.
Empirical Test
Obviously, I can transform a few points and then check if they are still gridded. I'd like to find a better way.
Compare DATUM
Using
pyproj.CRS('EPSG:4236').datum
, the datums are- EPSG:4236:
DATUM["World Geodetic System 1984", ELLIPSOID["WGS 84",6378137,298.257223563, LENGTHUNIT["metre",1]], ID["EPSG",6326]]
- EPSG:3395:
DATUM["World Geodetic System 1984", ELLIPSOID["WGS 84",6378137,298.257223563, LENGTHUNIT["metre",1]], ID["EPSG",6326]]
- EPSG:2193:
DATUM["New Zealand Geodetic Datum 2000", ELLIPSOID["GRS 1980",6378137,298.257222101, LENGTHUNIT["metre",1]], ID["EPSG",6167]]
So the DATUM for 4326 and 3395 are the exact same. Is it the case that two crs with the same exact datum string is always a scale/translate? Is it enough to compare just the DATUM portion (
World Geodetic System 1984
vsNew Zealand Geodetic Datum 2000
) or just the ID portion (EPSG 6326
vsEPSG 6167
)? Or, is this not going to work at all (can two CRS share the same DATUM but still require distortion when transforming)?Furthermore, even if that does work, does every crs have a "datum"? It would be nice to detect this even with custom crs.
- EPSG:4236:
Compare PROJ4
The proj4 strings:
- EPSG:4326:
+proj=longlat +datum=WGS84 +no_defs +type=crs
- EPSG:3395:
+proj=merc +lon_0=0 +k=1 +x_0=0 +y_0=0 +datum=WGS84 +units=m +no_defs +type=crs
- EPSG:2193:
+proj=tmerc +lat_0=0 +lon_0=173 +k=0.9996 +x_0=1600000 +y_0=10000000 +ellps=GRS80 +units=m +no_defs +type=crs
Anything I can do with this? Seems like there is not quite sufficient info, here, but maybe some values have "defaults"? Which parts of these strings indicate scaling/transformation as opposed to other computations?
- EPSG:4326:
Compare WKT
The WKT strings or JSON representations have the same basic parameters, but they are more detailed and may be simpler to use. I won't post them here entirely. But, for example, EPSG:3395 has EPSG:4326 as the
base_crs
, but of course EPSG:2193 does not. Would comparingbase_crs
be sufficient?Inspect the pyproj
Transformer
objectI'm using a pyproj
Transformer
, so maybe the best thing to do is simply inspect it.- Using the transformer
definition
directly might be better than trying to compare proj4 strings; which part(s) should I use? - Alternatively, the
operations
is empty for the 4236 to 3395 transformation, but not empty for 4326 to 2193. Is it as simple as that?
>>> t = pyproj.Transformer.from_crs('EPSG:4236', 'EPSG:3395') >>> t.definition 'proj=pipeline step proj=axisswap order=2,1 step proj=unitconvert xy_in=deg xy_out=rad step proj=merc lon_0=0 k=1 x_0=0 y_0=0 ellps=WGS84' >>> t.operations () >>> t = pyproj.Transformer.from_crs('EPSG:4236', 'EPSG:2193') >>> t.definition 'proj=pipeline step proj=axisswap order=2,1 step proj=unitconvert xy_in=deg xy_out=rad step proj=tmerc lat_0=0 lon_0=173 k=0.9996 x_0=1600000 y_0=10000000 ellps=GRS80 step proj=axisswap order=2,1' >>> t.operations (<Coordinate Operation: Transformation> Inverse of NZGD2000 to WGS 84 (1) Area of Use: - name: New Zealand - bounds: (160.6, -55.95, -171.2, -25.88), <Coordinate Operation: Conversion> New Zealand Transverse Mercator 2000 Area of Use: - name: New Zealand - onshore - bounds: (166.37, -47.33, 178.63, -34.1))
- Using the transformer