Hi Even, Thank you for your detailed reply. So if I want to use 1169 or also 1660, that is the first candidate listed by
projinfo -s EPSG:3003 -t EPSG:4326 --spatial-test intersects I can transform the coordinates from 3003 to 4265, that is the geodetic SRS used by 3003 sourceProj = osr.SpatialReference() sourceProj.ImportFromEPSG(3003) sourceGeo = osr.SpatialReference() sourceGeo.ImportFromEPSG(4265) proj2geo = osr.CoordinateTransformation(sourcePro,sourceGeo) point.Transform(proj2geo) then I can use EPSG:1169. I tried the code snippet and it seems to be working properly. Do you think this is a correct approach? All the best Stefano Il giorno ven 5 mar 2021 alle ore 12:11 Even Rouault < [email protected]> ha scritto: > Stefano, > > urn:ogc:def:coordinateOperation:EPSG::1169 is a transformation between > EPSG:4265 (the geographic CRS of Monte Mario) and EPSG:4326. You cannot use > it directly when the source CRS is EPSG:3003, which is a projected CRS > based on EPSG:4265. > > You could instead use > > options.SetOperation("""+proj=pipeline > +step +inv +proj=tmerc +lat_0=0 +lon_0=9 +k=0.9996 +x_0=1500000 +y_0=0 > +ellps=intl > +step +proj=push +v_3 > +step +proj=cart +ellps=intl > +step +proj=helmert +x=-225 +y=-65 +z=9 > +step +inv +proj=cart +ellps=WGS84 > +step +proj=pop +v_3 > +step +proj=unitconvert +xy_in=rad +xy_out=deg > +step +proj=axisswap +order=2,1""") > > as suggested by the output of projinfo -s EPSG:3003 -t EPSG:4326 > --spatial-test intersects > > Even > Le 05/03/2021 à 09:42, Stefano Iacovella a écrit : > > Hi all, > > I am trying to use GDAL, via Python bindings, to convert coordinates. > I would like to select a specific transformation when converting among CRS > with different Datums. Reading the documentation my understanding is that I > need to create a osr.CoordinateTransformation with options parameter > valorized properly. > I am using binaries from OSGEO4WIN, the release are 3.1.4 for GDAL and > 6.3.2 for PROJ > This is a code snippet I am using to debug: > > from osgeo import ogr > from osgeo import osr > > source = osr.SpatialReference() > source.ImportFromEPSG(3003) > > target = osr.SpatialReference() > target.ImportFromEPSG(4326) > > options = osr.CoordinateTransformationOptions() > options.SetOperation("urn:ogc:def:coordinateOperation:EPSG::1169") > transform = osr.CoordinateTransformation(source, target, options) > > point = ogr.CreateGeometryFromWkt("POINT (1622987 4756952)") > point.Transform(transform) > > print point.ExportToWkt() > > I have no errors running it but the coordinates of the point, after the > transformation, are infinite. If I create the CoordinateTransformation > without options like this: > > transform = osr.CoordinateTransformation(source, target) > > The point coordinates after transformation are ok but of course the > conversion is performed using the Coordinate Operation automatically > selected by the system. > Am I doing something wrong? > > Thank you in advance for your help. > > All the best, > > Stefano > > > _______________________________________________ > gdal-dev mailing > [email protected]https://lists.osgeo.org/mailman/listinfo/gdal-dev > > -- http://www.spatialys.com > >
_______________________________________________ gdal-dev mailing list [email protected] https://lists.osgeo.org/mailman/listinfo/gdal-dev
