TransWikia.com

Reprojection with gdal.Warp doesn't work

Geographic Information Systems Asked by delViento on April 23, 2021

I am trying to reproject a raster. I have tried gdal.Warp and it doesn’t work.

RuntimeError: Cannot find coordinate operations from ‘EPSG:6257’ to ‘EPSG:6686’.

Working in Python 3, Windows 64bit

  • GDAL 3.0.2
  • proj 6.2.1

EPSG source Projection is 6257

PROJCS["MAGNA-SIRGAS / Medellin urban grid",
    GEOGCS["MAGNA-SIRGAS",
        DATUM["Marco_Geocentrico_Nacional_de_Referencia",
            SPHEROID["GRS 1980",6378137,298.257222101,
                AUTHORITY["EPSG","7019"]],
            TOWGS84[0,0,0,0,0,0,0],
            AUTHORITY["EPSG","6686"]],
        PRIMEM["Greenwich",0,
            AUTHORITY["EPSG","8901"]],
        UNIT["degree",0.0174532925199433,
            AUTHORITY["EPSG","9122"]],
        AUTHORITY["EPSG","4686"]],
    PROJECTION["Colombia_Urban"],
    PARAMETER["Latitude of natural origin",6.22920888888889],
    PARAMETER["Longitude of natural origin",-75.5648869444444],
    PARAMETER["False easting",835378.647],
    PARAMETER["False northing",1180816.875],
    PARAMETER["Projection plane origin height",1510],
    UNIT["metre",1,
        AUTHORITY["EPSG","9001"]],
    AXIS["Northing",NORTH],
    AXIS["Easting",EAST],
    AUTHORITY["EPSG","6257"]]

I want to reproject to the following projection (which is the proj of a vector):

EPSG destination Projection is 6686

PROJCS["MAGNA-SIRGAS / Medellin urban grid",
    GEOGCS["MAGNA-SIRGAS",
        DATUM["Marco_Geocentrico_Nacional_de_Referencia",
            SPHEROID["GRS 1980",6378137,298.257222101,
                AUTHORITY["EPSG","7019"]],
            AUTHORITY["EPSG","6686"]],
        PRIMEM["Greenwich",0],
        UNIT["Degree",0.0174532925199433]],
    PROJECTION["IGAC_Plano_Cartesiano"],
    PARAMETER["False_Easting",835378.647],
    PARAMETER["False_Northing",1180816.875],
    PARAMETER["Longitude_Of_Center",-75.5648869444444],
    PARAMETER["Latitude_Of_Center",6.22920888888889],
    PARAMETER["Height",1510],
    UNIT["metre",1,
        AUTHORITY["EPSG","9001"]],
    AXIS["Easting",EAST],
    AXIS["Northing",NORTH]]

I tried:

input_raster = os.path.join(geotiff_path, 'aoi1.tif')
output_raster = os.path.join(geotiff_path, 'aoi1_reproj.tif')

gdal.Warp(output_raster,input_raster,dstSRS='EPSG:6686')

Any ideas? Thanks

EDIT:

when running from the command line:

(rio) C:Usersxxxxxx>projinfo -s epsg:6257 -t epsg:6686

Candidate operations found: 1
-------------------------------------
Operation n┬░1:

unknown id, Inverse of Medellin urban grid + Ballpark geographic offset from MAGNA-SIRGAS to JGD2011 + Japan Plane Rectangular CS zone XVIII, unknown accuracy, World, has ballpark transformation

PROJ string:
Error when exporting to PROJ string: Unsupported conversion method: Colombia Urban

WKT2_2018 string:
CONCATENATEDOPERATION["Inverse of Medellin urban grid + Ballpark geographic offset from MAGNA-SIRGAS to JGD2011 + Japan Plane Rectangular CS zone XVIII",
    SOURCECRS[
        PROJCRS["MAGNA-SIRGAS / Medellin urban grid",
            BASEGEOGCRS["MAGNA-SIRGAS",
                DATUM["Marco Geocentrico Nacional de Referencia",
                    ELLIPSOID["GRS 1980",6378137,298.257222101,
                        LENGTHUNIT["metre",1]]],
                PRIMEM["Greenwich",0,
                    ANGLEUNIT["degree",0.0174532925199433]],
                ID["EPSG",4686]],
            CONVERSION["Medellin urban grid",
                METHOD["Colombia Urban",
                    ID["EPSG",1052]],
                PARAMETER["Latitude of natural origin",6.22920888888889,
                    ANGLEUNIT["degree",0.0174532925199433],
                    ID["EPSG",8801]],
                PARAMETER["Longitude of natural origin",-75.5648869444444,
                    ANGLEUNIT["degree",0.0174532925199433],
                    ID["EPSG",8802]],
                PARAMETER["False easting",835378.647,
                    LENGTHUNIT["metre",1],
                    ID["EPSG",8806]],
                PARAMETER["False northing",1180816.875,
                    LENGTHUNIT["metre",1],
                    ID["EPSG",8807]],
                PARAMETER["Projection plane origin height",1510,
                    LENGTHUNIT["metre",1],
                    ID["EPSG",1039]]],
            CS[Cartesian,2],
                AXIS["northing (N)",north,
                    ORDER[1],
                    LENGTHUNIT["metre",1]],
                AXIS["easting (E)",east,
                    ORDER[2],
                    LENGTHUNIT["metre",1]],
            ID["EPSG",6257]]],
    TARGETCRS[
        PROJCRS["JGD2011 / Japan Plane Rectangular CS XVIII",
            BASEGEOGCRS["JGD2011",
                DATUM["Japanese Geodetic Datum 2011",
                    ELLIPSOID["GRS 1980",6378137,298.257222101,
                        LENGTHUNIT["metre",1]]],
                PRIMEM["Greenwich",0,
                    ANGLEUNIT["degree",0.0174532925199433]],
                ID["EPSG",6668]],
            CONVERSION["Japan Plane Rectangular CS zone XVIII",
                METHOD["Transverse Mercator",
                    ID["EPSG",9807]],
                PARAMETER["Latitude of natural origin",20,
                    ANGLEUNIT["degree",0.0174532925199433],
                    ID["EPSG",8801]],
                PARAMETER["Longitude of natural origin",136,
                    ANGLEUNIT["degree",0.0174532925199433],
                    ID["EPSG",8802]],
                PARAMETER["Scale factor at natural origin",0.9999,
                    SCALEUNIT["unity",1],
                    ID["EPSG",8805]],
                PARAMETER["False easting",0,
                    LENGTHUNIT["metre",1],
                    ID["EPSG",8806]],
                PARAMETER["False northing",0,
                    LENGTHUNIT["metre",1],
                    ID["EPSG",8807]]],
            CS[Cartesian,2],
                AXIS["northing (X)",north,
                    ORDER[1],
                    LENGTHUNIT["metre",1]],
                AXIS["easting (Y)",east,
                    ORDER[2],
                    LENGTHUNIT["metre",1]],
            ID["EPSG",6686]]],
    STEP[
        CONVERSION["Inverse of Medellin urban grid",
            METHOD["Inverse of Colombia Urban",
                ID["INVERSE(EPSG)",1052]],
            PARAMETER["Latitude of natural origin",6.22920888888889,
                ANGLEUNIT["degree",0.0174532925199433],
                ID["EPSG",8801]],
            PARAMETER["Longitude of natural origin",-75.5648869444444,
                ANGLEUNIT["degree",0.0174532925199433],
                ID["EPSG",8802]],
            PARAMETER["False easting",835378.647,
                LENGTHUNIT["metre",1],
                ID["EPSG",8806]],
            PARAMETER["False northing",1180816.875,
                LENGTHUNIT["metre",1],
                ID["EPSG",8807]],
            PARAMETER["Projection plane origin height",1510,
                LENGTHUNIT["metre",1],
                ID["EPSG",1039]],
            ID["INVERSE(EPSG)",6225]]],
    STEP[
        COORDINATEOPERATION["Ballpark geographic offset from MAGNA-SIRGAS to JGD2011",
            SOURCECRS[
                GEOGCRS["MAGNA-SIRGAS",
                    DATUM["Marco Geocentrico Nacional de Referencia",
                        ELLIPSOID["GRS 1980",6378137,298.257222101,
                            LENGTHUNIT["metre",1]]],
                    PRIMEM["Greenwich",0,
                        ANGLEUNIT["degree",0.0174532925199433]],
                    CS[ellipsoidal,2],
                        AXIS["geodetic latitude (Lat)",north,
                            ORDER[1],
                            ANGLEUNIT["degree",0.0174532925199433]],
                        AXIS["geodetic longitude (Lon)",east,
                            ORDER[2],
                            ANGLEUNIT["degree",0.0174532925199433]],
                    ID["EPSG",4686]]],
            TARGETCRS[
                GEOGCRS["JGD2011",
                    DATUM["Japanese Geodetic Datum 2011",
                        ELLIPSOID["GRS 1980",6378137,298.257222101,
                            LENGTHUNIT["metre",1]]],
                    PRIMEM["Greenwich",0,
                        ANGLEUNIT["degree",0.0174532925199433]],
                    CS[ellipsoidal,2],
                        AXIS["geodetic latitude (Lat)",north,
                            ORDER[1],
                            ANGLEUNIT["degree",0.0174532925199433]],
                        AXIS["geodetic longitude (Lon)",east,
                            ORDER[2],
                            ANGLEUNIT["degree",0.0174532925199433]],
                    ID["EPSG",6668]]],
            METHOD["Geographic2D offsets",
                ID["EPSG",9619]],
            PARAMETER["Latitude offset",0,
                ANGLEUNIT["degree",0.0174532925199433],
                ID["EPSG",8601]],
            PARAMETER["Longitude offset",0,
                ANGLEUNIT["degree",0.0174532925199433],
                ID["EPSG",8602]]]],
    STEP[
        CONVERSION["Japan Plane Rectangular CS zone XVIII",
            METHOD["Transverse Mercator",
                ID["EPSG",9807]],
            PARAMETER["Latitude of natural origin",20,
                ANGLEUNIT["degree",0.0174532925199433],
                ID["EPSG",8801]],
            PARAMETER["Longitude of natural origin",136,
                ANGLEUNIT["degree",0.0174532925199433],
                ID["EPSG",8802]],
            PARAMETER["Scale factor at natural origin",0.9999,
                SCALEUNIT["unity",1],
                ID["EPSG",8805]],
            PARAMETER["False easting",0,
                LENGTHUNIT["metre",1],
                ID["EPSG",8806]],
            PARAMETER["False northing",0,
                LENGTHUNIT["metre",1],
                ID["EPSG",8807]],
            ID["EPSG",17818]]],
    USAGE[
        SCOPE["unknown"],
        AREA["World"],
        BBOX[-90,-180,90,180]]]

One Answer

The answer that was given in the gdal-dev mailing list:

EPSG:6257 uses the "Colombia Urban" projection method, whose support has been added only in PROJ 7.2 ( https://proj.org/operations/projections/col_urban.html )

And there must be an error in the target EPSG code in your post: EPSG:6686 is a Japanese CRS (well, you could technically try to do the reprojection, but it is unlikely this is what you want to do)

Correct answer by user30184 on April 23, 2021

Add your own answers!

Ask a Question

Get help from others!

© 2024 TransWikia.com. All rights reserved. Sites we Love: PCI Database, UKBizDB, Menu Kuliner, Sharing RPP