Fiona fails to write CRS to shapefile

Question

I want to save some geometry to a shapefile with the python27 module fiona and manually specify the coordinate reference system (CRS).

Here's what I tried:

import fiona
import fiona.crs

outpath = "shape.shp"

schema = {'geometry': 'Point', 'properties': {'place': 'str'}}
point = {u"type": u"Point", u"coordinates": [0.0, 51.4769]}
properties = {'place': 'Greenwich'}
wgs84 = fiona.crs.from_epsg(4326)

with fiona.open(outpath, 'w', crs=wgs84, driver='ESRI Shapefile',
                schema=schema) as layer:
    layer.write({'geometry': point, 'properties': properties})

Unfortunately the shape.prj file that's produced is empty, hence the CRS is not specified in the output file.

A little diagnostics: After running the script, the layer has the CRS specified:

In [10]: layer.crs
Out[10]: {'init': 'epsg:4326', 'no_defs': True}

Any idea what's going wrong there?

I can work around this issue by downloading the .PRJ file from spatialreference.org, placing it in the same folder and renaming it shape.prj. I guess I'll figure out how to do that programmatically if there's no direct solution (should be easy).

Answer 1

I had this issue, and fixed it by setting the GDAL_DATA variable correctly, so that fiona can map the projection strings correctly. For reference, I'm using Anaconda, the Spyder IDE, Fiona 1.8.4, Python 3.6.8, and GDAL 2.3.3.

While Anaconda usually sets the GDAL_DATA variable upon entering the virtual environment, using another IDE like Spyder will not preserve it, and thus causes issues where fiona can't export the CRS correctly.

You can test this fix by trying to printing out a EPSG to WKT transformation before & after setting the GDAL_DATA variable explictly.

Without setting GDAL_DATA:

import os
print('GDAL_DATA' in os.environ)
from osgeo import osr
srs = osr.SpatialReference() # Declare a new SpatialReference 
srs.ImportFromEPSG(3413) # Import the EPSG code into the new object srs
print(srs.ExportToWkt()) # Print the result before transformation to ESRI WKT (prints nothing)

Results in:

False

With setting GDAL_DATA:

import os
os.environ['GDAL_DATA'] = 'D:\\ProgramData\\Anaconda3\\envs\\cfm\\Library\\share\\gdal'
print('GDAL_DATA' in os.environ)
from osgeo import , osr
srs = osr.SpatialReference() # Declare a new SpatialReference 
srs.ImportFromEPSG(3413) # Import the EPSG code into the new object srs
print(srs.ExportToWkt()) # Print the result before transformation to ESRI WKT (prints the EPSG converted to a WKT string)

Results in:

True
PROJCS["WGS 84 / NSIDC Sea Ice Polar Stereographic North",GEOGCS["WGS 84",DATUM["WGS_1984",SPHEROID["WGS 84",6378137,298.257223563,AUTHORITY["EPSG","7030"]],AUTHORITY["EPSG","6326"]],PRIMEM["Greenwich",0,AUTHORITY["EPSG","8901"]],UNIT["degree",0.0174532925199433,AUTHORITY["EPSG","9122"]],AUTHORITY["EPSG","4326"]],PROJECTION["Polar_Stereographic"],PARAMETER["latitude_of_origin",70],PARAMETER["central_meridian",-45],PARAMETER["scale_factor",1],PARAMETER["false_easting",0],PARAMETER["false_northing",0],UNIT["metre",1,AUTHORITY["EPSG","9001"]],AXIS["X",EAST],AXIS["Y",NORTH],AUTHORITY["EPSG","3413"]]

Answer 2

You can use GDAL/OGR in your code with fiona:

import fiona

from osgeo import osr              ###
srs = osr.SpatialReference()       ###
srs.SetFromUserInput("EPSG:4326")  ###
wgs84 = srs.ExportToProj4()        ###
# or
# wgs84 = srs.ExportToPrettyWkt()  ###
# or
# wgs84 = srs.ExportToWkt()        ###

outpath  = "shape.shp"

schema = {'geometry': 'Point', 'properties': {'place': 'str'}}
point = {u"type": u"Point", u"coordinates": [0.0, 51.4769]}
properties = {'place': 'Greenwich'}

with fiona.open(outpath, 'w', 
                crs=wgs84,         ###
                # crs_wkt=wgs84,   ### use this line if you use wkt
                driver='ESRI Shapefile',
                schema=schema) as layer:
    layer.write({'geometry': point, 'properties': properties})