0

I'm trying to write a Python script to convert a vector layer using 4326 to EPSG:54002 World Equidistant Cylindrical for road distance calculations, but the script fails to do the conversion, I've tried numerous methods to get it right:

Using the re-project layer algorithm:

reprojectedLayer = processing.run("native:reprojectlayer", {'INPUT': roadLayer ,'TARGET_CRS': 'EPSG:54002', 'OUTPUT':'TEMPORARY_OUTPUT'})['OUTPUT']
reprojectedLayer = processing.run("native:reprojectlayer", {'INPUT': roadLayer ,'TARGET_CRS': QgsCoordinateReferenceSystem('EPSG:54002'), 'OUTPUT':'TEMPORARY_OUTPUT'})['OUTPUT']

Calling the setCrs method on the layer:

roadLayer.setCrs(QgsCoordinateReferenceSystem(54002, QgsCoordinateReferenceSystem.PostgisCrsId))
roadLayer.setCrs(QgsCoordinateReferenceSystem('EPSG:54002'))

Even just writing the layer to a file while specifying a different CRS:

QgsVectorFileWriter.writeAsVectorFormat(roadLayer, sys.argv[2] + "/roads", "utf-8",  QgsCoordinateReferenceSystem('EPSG:54002'), "GPKG")

Unfortunately, none of these have worked. If I save the file after any of the operations above and open it in QGIS it either just opens it fine as 4326 or opens it and displays a message that the CRS is invalid and that it will be defaulting to 4326.

One interesting thing to note is that when I specify the crs while writing the file, this pops up in the console:

proj_create_from_database: crs not found

It doesn't do that with the other methods, only when I write the file. This message doesn't make sense since if I try to reproject the layer to 54002 in QGIS, it works fine. Its only in the script that it fails.

6
  • 54002 isn't a code in the EPSG registry
    – nmtoken
    Feb 12, 2020 at 6:59
  • Hmmm... Its listed in QGIS as a CRS and seems to work fine if I reproject to it using QGIS. the Python code I used to call the reproject layer method is from the processing log window. But I do see that epsg.io has it as ESRI:54002, let me try that.
    – Walter
    Feb 12, 2020 at 7:33
  • Okay no that didn't work, maybe the question I should be asking is, what CRS should I be using for measuring the length of a road (that isn't UTM)
    – Walter
    Feb 12, 2020 at 7:56
  • What if you put the proj definition in your code, rather than using an EPSG /ESRI number.. If that's possible.
    – nmtoken
    Feb 12, 2020 at 9:07
  • 1
    Unfortunately, there's no map projection that's great at road distances. If the data is in a limited area, some large scale (designed for a smallish area) projCRS would probably give okay values, but using a geodesic distance-based function would be better.
    – mkennedy
    Feb 12, 2020 at 18:18

1 Answer 1

2

This script working fine for me :

params = {
    'INPUT': r'path/to/input/file.shp',
    'OUTPUT': r'path/to/the/export/file.shp',
    'TARGET_CRS': QgsCoordinateReferenceSystem('EPSG:54002')
}

processing.run("qgis:reprojectlayer", params)
vlayer = QgsVectorLayer(r'path/to/the/export/file.shp', 'reprojected', 'ogr')
QgsProject.instance().addMapLayer(vlayer)

The memory output work fine too :

params = {
    'INPUT': r'path/to/input/file.shp',
    'OUTPUT': 'memory:',
    'TARGET_CRS': QgsCoordinateReferenceSystem('EPSG:54002')
}

vlayer = processing.run("qgis:reprojectlayer", params)['OUTPUT']
QgsProject.instance().addMapLayer(vlayer)

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge that you have read and understand our privacy policy and code of conduct.

Not the answer you're looking for? Browse other questions tagged or ask your own question.