[PROJ] question about EPSG:3003

classic Classic list List threaded Threaded
2 messages Options
Reply | Threaded
Open this post in threaded view
|

[PROJ] question about EPSG:3003

a.furieri
I notice an unpleasant regression in the proj-string returned by
proj_as_proj_string() for EPSG:3003 Monte Mario / Italy zone 1
in all previous versions, accordingly to the string returned by
calling the GDAL function OSRExportToProj4(), this was defined
as follows:

+proj=tmerc +lat_0=0 +lon_0=9 +k=0.9996 +x_0=1500000 +y_0=0 \
+ellps=intl +towgs84=-104.1,-49.1,-9.9,0.971,-2.917,0.714,-11.68 \
+units=m +no_defs

but in PROJ.6 the string returned by calling proj_as_proj_string()
is now:

+proj=tmerc +lat_0=0 +lon_0=9 +k=0.9996 +x_0=1500000 +y_0=0 \
+ellps=intl +units=m +no_defs +type=crs

as you can notice, the +towgs84 member has been suppressed,
and the sad practical consequence is that any transformation
based on the above proj-string (as e.g. from EPSG:4326 to
EPSG:3003) is now shifted by about 75 metres (at least, this
is true on Tuscany).

if I remember well after so many years, PROJ.4/GDAL initially
supported a definition for 3003 lacking the +towgs84 member,
and this caused endless troubles to all italian users.
but since 2012 (more or less) PROJ.4/GDAL started supporting
a revised proj-string including the +towgs84 member, and it
was an highly appreciated improvement.
now PROJ.6 seems to mark a bad regression to the dark past.

note: this is not at all an issue for new SpatiaLite versions
fully based on PROJ.6, because in this case ST_Transform()
will use by default the CRS definitions retrieved from the
PROJ.6 own SQLite database.
but nevertheless a very dangerous scenario exists.
imagine some SpatiaLite database created by a recent
version based on PROJ.6
in this case the "spatial_ref_sys" table will be
populated by inserting the proj-strings returned
by proj_as_proj_string().
if such a database will be connected in a second
time to some earlier version of SpatiaLite still
based on old PROJ.4 (a not at all uncommon practice
in the  peculiar world of spatialite) any call to
ST_Transform() will then be based on the proj-strings
stored into "spatial_ref_sys", and any transformation
to/from srid=3003 will be inexorably affected by
severe shifts.

any idea/comment/suggestion about this ?

bye Sandro

p.s. the same identical issue affects srids
3004, 23032 and 23033
_______________________________________________
PROJ mailing list
[hidden email]
https://lists.osgeo.org/mailman/listinfo/proj
Reply | Threaded
Open this post in threaded view
|

Re: question about EPSG:3003

Even Rouault-2
Hi Sandro,

First, if you instanciate with proj_create() / proj_create_from_database() and
export with proj_as_proj_string(), you'll never get a towgs84 clause.
Since the CRS by itself has no transformation clause.
However if you look at the output of the projinfo utility or
OSRImportFromEPSG(), you may, sometimes, get a towgs84 clause when exporting
to PROJ string (e.g with EPSG:4322 / WGS72), because they use the
proj_crs_create_bound_crs_to_WGS84() function which tries to find the most
precise transformation from the CRS to EPSG:4326, and, that's the important
point, which is valid for the whole area of use of the CRS.

For EPSG:3003, there are however 3 possibilities (you can ignore the last one
in the below output), none of which cover the whole area of use of the
projected CRS.

$ src/projinfo -s EPSG:3003 -t EPSG:4326 --spatial-test intersects --summary
Candidate operations found: 4
unknown id, Inverse of Italy zone 1 + Monte Mario to WGS 84 (4), 4 m, Italy -
mainland
unknown id, Inverse of Italy zone 1 + Monte Mario to WGS 84 (2), 4 m, Italy -
Sardinia onshore
unknown id, Inverse of Italy zone 1 + Monte Mario to WGS 84 (11), 10 m, Italy
- Sicily Strait west of 13°E
unknown id, Inverse of Italy zone 1 + Ballpark geographic offset from Monte
Mario to WGS 84, unknown accuracy, World, has ballpark transformation

Previous PROJ version hard-coded the transformation that has the largest area
of use, but that's not necessarily ideal in all use cases (Sardinia or Sicily)

So, this is intended behaviour. I'm campaigning here and on all relevants
channels: "With PROJ 6, PROJ strings to express CRS are very suboptimal. They
should be restricted to express a transformation pipeline"

If you want Spatialite to have the same behaviour as before regarding
generation of PROJ strings, you could possibly use the
proj_create_operations() API to list the possible transformations, by setting
proj_operation_factory_context_set_spatial_criterion(PROJ_SPATIAL_CRITERION_PARTIAL_INTERSECTION)
on the factory_ctx, pickup the one that has the largest area of use, create
the final bound CRS with that transformation with proj_crs_create_bound_crs(),
and export it to a PROJ string.
Admitedly, we could have a releaxed proj_crs_create_bound_crs_to_WGS84()
implementation in PROJ itself to do that, but I'm not sure if we really need
to replicate/encourage past suboptimal behaviour.

Even

--
Spatialys - Geospatial professional services
http://www.spatialys.com
_______________________________________________
PROJ mailing list
[hidden email]
https://lists.osgeo.org/mailman/listinfo/proj