[gdal-dev] Troubles in coord conversions with GDAL+proj6

[Joaquim Manuel Freire Luís]

|>> and the other is that in this later case it outputs lat,long instead of lon,lat
|>
|>That one is expected. SetWellKnownGeogCS( "WGS84" ) now returns a SRS
|>with lat, lon order. You can force "GIS" friendly order by calling
|>oSRS.SetAxisMappingStrategy(OAMS_TRADITIONAL_GIS_ORDER);

Hmm, lat,longs are bad 😊
When did that change? (SetAxisMappingStrategy was introduced in GDAL 2.5 only)

Anyway, I added this to my code

#if (GDAL_VERSION_NUM >= 2500)
	 oSrcSRS.SetAxisMappingStrategy(OAMS_TRADITIONAL_GIS_ORDER);
	 oDstSRS.SetAxisMappingStrategy(OAMS_TRADITIONAL_GIS_ORDER);
#endif

And things changed a bit (now lat is always = 0). And I think I found the issue (Matlab craziness) , which is simply unbelievable. 

OK, starting again. This now gives the correct value

>> prj.SrcProjSRS = '+proj=utm +zone=29 +datum=WGS84 +units=m +no_defs';
>> xy = ogrproj([0 0], prj)

xy =  -13.488743884387199                   0

But on a second run is already wrong

>> xy = ogrproj([0 0], prj)
xy = -13.488864729988810                   0

Now the unbelievable part. Let's use x,y variables instead of [0 0] (note, in Matlab [0 0] is a double precision array) and the result is: !!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

CORRECT
>> x=0;y=0;
>> xy = ogrproj([y y], prj)
xy = -13.488743884387199                   0

WRONG
>> xy = ogrproj([0 0], prj)
xy = -13.488864731071464                   0

CORRECT
>> xy = ogrproj([y y], prj)
xy = -13.488743884387199                   0

So, I'm sorry for the noise but I think it was not all my fault.