<div dir="ltr"><div dir="ltr">Hi,<br></div><div dir="ltr"><br></div><div>(resending, as apparently my yesterday message wasn't delivered)<br></div><div dir="ltr"><br>Now that the RFC2 implementation has been merged<br>( <a href="https://github.com/OSGeo/proj.4/pull/1175">https://github.com/OSGeo/proj.4/pull/1175</a> ), come the next steps. One of <br>them will be interation in GDAL, and then libgeotiff, so there will be fixes <br>and complementary work that will be done as a result of those integrations.<br><br>For PROJ itself, one of the main thing I'd like to address is the fact that we <br>have now two CRS databases:<br>- the old one in the data/epsg and data/IGNF files<br>- the new one in the SQLite3 proj.db database<br><br>The new one contains the content of the old one (with some differences, like <br>the EPSG database in proj.db is fresher than in the old one). So the plan is <br>to kill the old 'epsg' and 'IGNF' file that are used by the +init=file:code <br>syntax and go instead into the database to fetch their definitions. If no <br>match is found, so for other authorities than EPSG, IGNF or ESRI currently, <br>then the existing text file based mechanism would still be used.<br><br>All good ? Well, not quite. One of the main issue I foresee before doing any <br>implementation is axis order. In the old epsg file, geographic CRS use the <br>longitude,latitude order, whereas proj.db and the new code I've developped <br>respect the axis order declared by the authority, so in the case of EPSG, <br>latitude,longitude for geographic CRS. This also applies for a few projected <br>CRS: 'epsg' file uses generally easting,northing order (except for a few <br>Transverse Mercator South Orientated projected CRS where it uses the expected <br>westing,southing order), whereas EPSG might have northing,easting order for a <br>few ones, or southing,westing for a few Krovak-based CRS.<br><br>What to do ? I can see different options:<br><br>- just use the 'correct' axis order, and mark this as a breaking change. That <br>is if you use +init=epsg:4326, you'd have to provide / will get latitude first<br><br>- emulate the previous axis order in an attempt of being backward compatible. <br>Not necessarily easy to do since that means we have to remember somewhere the <br>old axis order and add the inversion when needed in some glue code between old <br>and new code. And that would undermine one of the interest of the new work to <br>be more axis aware than previously.<br><br>- have some switch in proj utilies ('proj', 'cs2cs', 'cct') to specify the <br>input and output order (like -input_axis=long,lat -output_axis=east,north) so <br>that they can adapt user input/output to the expected axis of the CRS used.<br>That would still be a breaking change since existing scripts would have to add <br>the new switches to preserve existing behaviour.<br><br>The same question also arises at the API level. For example,<br>proj_create_crs_to_crs(ctxt, "epsg:4326", "epsg:32632") currently expects the <br>input to be in longitude, latitude order.<br><br>There is also the question of units. Traditionnaly, in code paths that would <br>expect angular values, PROJ expects radians (frontend utilities do the <br>conversion into/from degrees). The new code assume sthe unit to be the one <br>declared by the coordinate system in the geographic CRS definition, ie degrees <br>in 99% of the case (or grad for a couple odd & old systems like "NTF (Paris)" <br>/ EPSG:4807)<br><br>Let's consider a practical example. "projinfo -s EPSG:4326 -t EPSG:4807" <br>outputs a transformation pipeline which is<br>"<br>+proj=pipeline +step +proj=axisswap +order=2,1 +step +proj=unitconvert <br>+xy_in=deg +xy_out=rad +step +inv +proj=hgridshift +grids=ntf_r93.gsb +step <br>+proj=longlat +ellps=clrk80ign +pm=paris +step +proj=unitconvert +xy_in=rad <br>+xy_out=grad +step +proj=axisswap +order=2,1"<br><br>The traduction in plain English of this is:<br>- take the input coordinate and reverse the first 2 compoments<br>- convert from degree to radian<br>- do a horizontal grid shift transformation<br>- apply a prime meridian shift<br>- convert from radian to grad<br>- take the ouput coordinate and reverse the first 2 compoments<br><br>That is it expects the input to be lat,long in degrees and output lat,long in <br>grads, as in the offical definition of those CRS.<br><br>Should the PJ* object returned by proj_create_crs_to_crs(ctxt, "epsg:4326", <br>"epsg:4807") just use that pipeline, and consequently expect the coordinates <br>passed into and received from proj_trans() to also use those orders and units.<br><br>Thoughts ?<br><br>Even<br><br>PS: the database might also contain geographic CRS definitions that are <br>longitude,latitude by design. For example I've added an entry for the <br>OGC:CRS84, and if you do "projinfo -s OGC:CRS84 -t EPSG:4807", then the <br>pipeline doesn't contain the initial axisswap operations.<br><br></div></div>