[PROJ] Adjusting priorities in proj.db: grid vs helmert transformations

[Even Rouault]

Craig,

The sorting algorithm of the operations returned by projinfo is given at 
https://proj.org/operations/operations_computation.html#filtering-and-sorting-of-coordinate-operations

I believe you run here in point 8 "consider as more relevant an 
operation whose area of use is larger" that comes before point 9 
"consider as more relevant an operation that has a better accuracy.". 
The rationale is that, until you know exactly which point to transform, 
an operation whose area of use is larger than another one is a better 
default.

But... if you use cs2cs or proj_trans() on the result of 
proj_create_crs_from_crs(), among the candidate operations, they will 
prefer the ones that can be instantiated (that is, in particular, whose 
grids are available) and that have a better accuracy.

Did you actually try with cs2cs ? With your modified record, 
au_icsm_GDA94_GDA2020_conformal_and_distortion.tif is used for me when 
using -35 120 for example.

$ echo -35 120 | PROJ_DEBUG=2 PROJ_NETWORK=ON PROJ_LIB=data src/cs2cs 
EPSG:4283 EPSG:7844
pj_open_lib(proj.db): call fopen(data/proj.db) - succeeded
pj_open_lib(proj.ini): call fopen(data/proj.ini) - succeeded
Using coordinate operation GDA94 to GDA2020 (2)
pj_open_lib(au_icsm_GDA94_GDA2020_conformal_and_distortion.tif): call 
fopen(data/au_icsm_GDA94_GDA2020_conformal_and_distortion.tif) - failed
pj_open_lib(GDA94_GDA2020_conformal_and_distortion.gsb): call 
fopen(data/GDA94_GDA2020_conformal_and_distortion.gsb) - failed
Using 
https://cdn.proj.org/au_icsm_GDA94_GDA2020_conformal_and_distortion.tif
34d59'59.951"S    120d0'0.037"E 0.000

And if you add for example --bbox 120,-40,130,-30 to the projinfo 
invokation to restrict a bit the area of use, you should see the order 
you expect.

Setting an operation accuracy to NULL / unknown makes it less prioritary 
than operations that have a known accuracy (point 6 of the above algorithm)

Even

Le 14/11/2021 à 21:31, Craig de Stigter a écrit :
> Hey folks
>
> We're trying to adjust proj to use a specific transformation between two
> CRSes *by default* (ie. we don't want our calling software to have to be
> explicit about which to use, since there are several places we're calling
> PROJ from and we'd like to make them automatically consistent)
>
> Stock proj 8.1.1 gives us these transform choices:
>
>> $ PROJ_NETWORK=ON projinfo -s EPSG:4283 -t EPSG:7844 --spatial-test
> intersects --summary
>> Candidate operations found: 5
>>
>> EPSG:8048, GDA94 to GDA2020 (1), 0.01 m, Australia including Lord Howe Island, Macquarie Island, Ashmore and Cartier Islands, Christmas Island, Cocos (Keeling) Islands, Norfolk Island. All onshore and offshore.
>> DERIVED_FROM(EPSG):8446, GDA94 to GDA2020 (3), 0.05 m, Australia - Australian Capital Territory; New South Wales; Northern Territory; Queensland; South Australia; Tasmania; Western Australia; Victoria.
>> DERIVED_FROM(EPSG):8447, GDA94 to GDA2020 (2), 0.05 m, Australia - Australian Capital Territory; New South Wales; Northern Territory; Queensland; South Australia; Tasmania; Western Australia; Victoria.
>> ...
>
> The one we want is DERIVED_FROM(EPSG):8447, GDA94 to GDA2020 (2), which is
> third in the list, and thus not chosen by default.
>
> Aiming to bump that to the top of the list, we modified our proj.db build
> to add this SQL:
>
>> -- When transforming from GDA94 to GDA2020,
>> -- use the distortion+conformal grid in preference to the conformal-only
> grid.
>> UPDATE grid_transformation
>> SET accuracy = 0.0099
>> WHERE auth_name = 'EPSG' AND code = 8447;
> We expected this to bump the transform to the top of the list. It did bump
> up by one slot, above the other grid-based transform (EPSG:8446). However,
> it's still below the helmert transformation (EPSG:8048):
>
>> Candidate operations found: 5
>> EPSG:8048, GDA94 to GDA2020 (1), 0.01 m, Australia including Lord Howe
> Island, Macquarie Island, Ashmore and Cartier Islands, Christmas
> Island, Cocos (Keeling) Islands, Norfolk Island. All onshore and
> offshore.
>> DERIVED_FROM(EPSG):8447, GDA94 to GDA2020 (2), 0.0099 m, Australia -
> Australian Capital Territory; New South Wales; Northern Territory;
> Queensland; South Australia; Tasmania; Western Australia; Victoria.
>> DERIVED_FROM(EPSG):8446, GDA94 to GDA2020 (3), 0.05 m, Australia -
> Australian Capital Territory; New South Wales; Northern Territory;
> Queensland; South Australia; Tasmania; Western Australia; Victoria.
>> ...
>
> I couldn't find anything documenting that PROJ prefers Helmert transforms
> over grid-based ones - merely that it prefers operations with defined
> accuracy, and then prefers transforms in order of accuracy. In which case I
> would have expected it to choose EPSG:8447 over EPSG:8048 now that we've
> hacked the accuracy number.
>
> Interestingly, setting the accuracy of the Helmert transform to NULL *did*
> fix this:
>
>> DERIVED_FROM(EPSG):8447, GDA94 to GDA2020 (2), 0.0499 m, Australia -
> Australian Capital Territory; New South Wales; Northern Territory;
> Queensland; South Australia; Tasmania; Western Australia; Victoria.
>> ...
>> EPSG:8048, GDA94 to GDA2020 (1), unknown accuracy, Australia including
> Lord Howe Island, Macquarie Island, Ashmore and Cartier Islands,
> Christmas Island, Cocos (Keeling) Islands, Norfolk Island. All onshore
> and offshore.
>
>
> So my questions are:
>
> 1. Is setting the accuracy of EPSG:8048 to NULL a bad idea? We'll still
> want to use EPSG:8048 as an offshore fallback.
> 2. Any idea why PROJ is still prioritising EPSG:8048 after our first tweak?
> 3. Is there a better way to accomplish this reprioritisation?
>
> Thanks a bunch!
>
-- 
http://www.spatialys.com
My software is free, but my time generally not.