createOperations(): emulate PROJ < 6 behavior when doing geocentric <--> geographic transformation between datum with unknown transformation (fixes #3360)

[rouault]

We switch from one arbitrary to the other one. The new behavior is perhaps slightly better because, besides matching PROJ <6 results, it matches also the behavior when operating in the geographic coordinate space between datums without known transformations.
@kbevers opinions ? I considered at that point to add some API switch to select which behavior is desired, but that seems overkill and nobody will figure that they need to use such potential switch, right ?

Likely not appropriate for 9.1.x backport

[kbevers]

I don't see why this solution is arbitrary. To me it is obviously correct to do the +proj=cart +inv step using the WGS84 ellipsoid and not the radius of a sphere. This might be easier for human eyes to see than for a computer program that follows strict rules.

In fact, the +proj=longlat +R=6378137 is not really doing anything in the pipeline, as seen below where I've removed it in the last example.

$ echo "4057449.87922115 1236185.10760693 4747287.50008112" | cct -t0 +proj=pipeline +step +inv +proj=cart +ellps=WGS84 +step +proj=longlat +R=6378137 +step +proj=unitconvert +xy_in=rad +xy_out=deg +step +proj=axisswap +order=2,1
 48.4111246681   16.9444298151      -74.6128        0.0000

$ echo "4057449.87922115 1236185.10760693 4747287.50008112" | cct -t0 +proj=pipeline +step +inv +proj=cart +ellps=WGS84 +step +proj=unitconvert +xy_in=rad +xy_out=deg +step +proj=axisswap +order=2,1     
 48.4111246681   16.9444298151      -74.6128        0.0000

Had there been a a projection in play on the output CRS the +R would of course play a role but the +proj=latlong +R=6378137 step would still be unnecessary.

So in conclusion, looks more like a bug fix than a breaking change in behaviour to me.

[rouault]

To me it is obviously correct

well, I'm sill a bit perplexed by that. anyway... Adding the 9.1 backport label then

[kbevers]

Geocentric (ECEF) does not have an ellipsoid... Right...

But WGS84 as a datum does, and that needs to be taken into consideration when converting from cartesian coordinates to geodetic coordinates.