[PROJ] Different transformation results between PROJ 9.0.0 and PROJ 9.0.1

[Even Rouault]

Alexander,

in 9.0.0 the operation that gets applied is:

unknown id, Inverse of WGS 84 to EGM96 height (1) + Inverse of NAD83 to 
WGS 84 (43) + NAD83 to NAD83(HARN) (3) + NAD83(HARN) to NAVD88 height 
(5) + Inverse of NAD83 to NAD83(HARN) (3) + SPCS83 California zone 1 (US 
Survey feet), 3.15 m, unknown domain of validity, at least one grid missing

PROJ string:
+proj=pipeline
   +step +proj=axisswap +order=2,1
   +step +proj=unitconvert +xy_in=deg +xy_out=rad
   +step +proj=vgridshift +grids=us_nga_egm96_15.tif +multiplier=1
   +step +inv +proj=vgridshift +grids=us_noaa_g1999u05.tif +multiplier=1
   +step +inv +proj=hgridshift +grids=us_noaa_cnhpgn.tif
   +step +proj=lcc +lat_0=39.3333333333333 +lon_0=-122 
+lat_1=41.6666666666667
         +lat_2=40 +x_0=2000000.0001016 +y_0=500000.0001016 +ellps=GRS80
   +step +proj=unitconvert +xy_in=m +xy_out=us-ft


and in 9.0.1/master:

unknown id, Inverse of WGS 84 to EGM96 height (1) + Inverse of 
NAD83(HARN) to WGS 84 (3) + NAD83(HARN) to NAVD88 height (5) + Inverse 
of NAD83 to NAD83(HARN) (3) + SPCS83 California zone 1 (US Survey feet), 
2.1 m, unknown domain of validity, at least one grid missing

PROJ string:
+proj=pipeline
   +step +proj=axisswap +order=2,1
   +step +proj=unitconvert +xy_in=deg +xy_out=rad
   +step +proj=vgridshift +grids=us_nga_egm96_15.tif +multiplier=1
   +step +proj=cart +ellps=WGS84
   +step +inv +proj=helmert +x=-0.991 +y=1.9072 +z=0.5129 
+rx=-0.0257899075194932
         +ry=-0.0096500989602704 +rz=-0.0116599432323421 +s=0
         +convention=coordinate_frame
   +step +inv +proj=cart +ellps=GRS80
   +step +inv +proj=vgridshift +grids=us_noaa_g1999u05.tif +multiplier=1
   +step +inv +proj=hgridshift +grids=us_noaa_cnhpgn.tif
   +step +proj=lcc +lat_0=39.3333333333333 +lon_0=-122 
+lat_1=41.6666666666667
         +lat_2=40 +x_0=2000000.0001016 +y_0=500000.0001016 +ellps=GRS80
   +step +proj=unitconvert +xy_in=m +xy_out=us-ft

So the difference is that in 9.0.1/master NAD83(HARN) to WGS 84 
transformation is modeled by a Helmert transformation, whereas in 9.0 
its is a null transformation. This is likely due to a new record in the 
EPSG database incorporated in >= 9.0.1

(I got this by looking at the pipeline name emitted by PROJ_DEBUG=2 on 
cs2cs, and then getting the list of transformations from projinfo)


Regarding the differences between specifying or not --bbox, this was a 
bug for which I've queued a fix in 
https://github.com/OSGeo/PROJ/pull/3282. With that fix I get the same 
result when specifying --bbox -122.1,40.9,-121.9,41.1 or not


Even

Le 04/08/2022 à 11:34, Alexander Nehrbass via PROJ a écrit :
>
> Hi everyone,
>
> I’m observing a significant difference between PROJ 9.0.0 and PROJ 
> 9.0.1/master[c482ea8] when doing some coordinate transformations:
>
> $ echo 41.0 -122.0 0.0 | PROJ_NETWORK=ON ./cs2cs EPSG:4326+5773 
> EPSG:2225+5703 -d 8
>
> 6561666.44282083    2247576.19215388 0.21377754 (in PROJ 9.0.0)
> 6561670.17486934    2247574.26932194 0.70338464 (PROJ 9.0.1/master)
>
> When running these commands with PROJ_DEBUG=3, I can see that 
> different transformation pipelines are being used. PROJ 9.0.0 uses an 
> 11 step pipeline whereas PROJ 9.0.1 uses a 7 step pipeline.
>
> While trying to understand where this difference might come from, I 
> started to play around with the –bbox option. Doing the above 
> transformation with a bounding box surprisingly does NOT lead to 
> different transformation results between PROJ 9.0.0 and PROJ 
> 9.0.1/master:
>
> $ echo 41.0 -122.0 0.0 | PROJ_NETWORK=ON ./cs2cs EPSG:4326+5773 
> EPSG:2225+5703 -d 8 --bbox -122.1,40.9,-121.9,41.1
>
> 6561666.44282083    2247576.19215388 0.31567764
>
> As you can see the horizontal coordinates are the same as when doing 
> the transformation without the bounding box in PROJ 9.0.0, however 
> there is a difference in the vertical coordinate. I also tried 
> changing the size of the bounding box which didn’t seem to have an 
> influence.
>
> My questions at the moment are:
> 1. Where does the difference between PROJ 9.0.0 and 9.0.1 when doing 
> the transformation without the bounding box come from?
> 2. How does using the bbox option affect the transformation results? 
> In particular:
> 2a. Why are the differences between PROJ 9.0.0 and 9.0.1 gone when 
> using a bounding box?
> 2b. Why does the vertical coordinate in PROJ 9.0.0 differ when doing 
> the transformation with and without the bounding box?
> 3. Related bonus question: I was under the impression that passing the 
> -bbox option to projinfo would result in projinfo returning exactly 
> one transformation. However, this seems in correct. Was I under the 
> wrong impression? E.g.:
>
> $ ./projinfo -s EPSG:4326+5773 -t EPSG:2225+5703 --spatial-test 
> intersects --bbox -122.1,40.9,-121.9,41.1
>
>
> Any help with finding answers to these question is highly appreciated.
>
> Thanks in advance and best regards,
> Alex
> _______________________________________________
> PROJ mailing list
> PROJ at lists.osgeo.org
> https://lists.osgeo.org/mailman/listinfo/proj

-- 
http://www.spatialys.com
My software is free, but my time generally not.