[gdal-dev] Difference between ogr2ogr and Python reprojection
Even Rouault
even.rouault at spatialys.com
Tue Mar 14 11:25:54 PDT 2023
Nicolas,
ogr2ogr uses the OGRGeometryFactory::transformWithOptions() method that
has logic to deal with poles, antimerdian, etc. whereas
OGRGeometry::transform() does basic point-to-point reprojection
You can achieve what you want by using the ogr.GeomTransformer(ct)
method that returns a ogr.GeomTransformer object that has a
Transform(ogr.Geometry) method. See below test from GDAL's autotest suite:
def test_ogr_geomtransfomer_default():
src = osr.SpatialReference()
src.ImportFromEPSG(32660)
dst = osr.SpatialReference()
dst.SetAxisMappingStrategy(osr.OAMS_TRADITIONAL_GIS_ORDER)
dst.ImportFromEPSG(4326)
ct = osr.CoordinateTransformation(src, dst)
geom = ogr.CreateGeometryFromWkt(
"LINESTRING(832864.275023695 0,835092.849076364 0)"
)
transformer = ogr.GeomTransformer(ct)
geom_dst = transformer.Transform(geom)
assert (
geom_dst.ExportToWkt()
== "MULTILINESTRING ((179.99 0.0,180 0),(-180 0,-179.99 0.0))"
)
Even
Le 14/03/2023 à 19:15, Dion-Degodez,Nicolas (ECCC) a écrit :
>
> _Context_
>
> We are trying to reproject a GeoJSON polygon to EPSG:4326.
>
> In some cases, depending on whether we use the ogr2ogr command line or
> the equivalent in python, the resulting polygon is different.
>
> In those situations, the command line gives us the appropriate
> polygon, whereas the python code gives us a distorted polygon.
>
> _Question_
>
> Any idea why there is a difference between the command line and the
> python results? Is there a mistake in our interpretation of the
> code/how we converted it in python?
>
> _Sample GeoJSON_: https://pastebin.com/hUiVVjFv
>
> Proj4: "+proj=stere +lat_0=90 +lat_ts=60 +lon_0=264 +x_0=0 +y_0=0
> +R=6371229 +units=m +no_defs"
>
> _Command line that works_: ogr2ogr -s_srs "+proj=stere +lat_0=90
> +lat_ts=60 +lon_0=264 +x_0=0 +y_0=0 +R=6371229 +units=m +no_defs"
> -t_srs EPSG:4326 result_command_line.json geojson.json
>
> Result: https://pasteboard.co/cTpZBRamEYv8.jpg
>
> _Python code_ (we would expect the same output as the command line):
> https://pastebin.com/QgWDBL1W
>
> Result: https://pasteboard.co/8iSv0vbnQlHw.jpg
>
> Thank you!
>
> Nicolas D.
>
>
> _______________________________________________
> gdal-dev mailing list
> gdal-dev at lists.osgeo.org
> https://lists.osgeo.org/mailman/listinfo/gdal-dev
--
http://www.spatialys.com
My software is free, but my time generally not.
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.osgeo.org/pipermail/gdal-dev/attachments/20230314/826099ff/attachment.htm>
More information about the gdal-dev
mailing list