Writing a shapefile using OGR, from Shapely Geometry - No feature added Error

Question

I am trying to write a shapefile from a Shapely geometry using PyQGIS in the QGIS 2.6.0 console.

Lines of code before the OGR method to write shapefile:

extnt = poly.extent() # gives a QtRectangle output -> extnt
x_min = extnt.xMinimum()
x_max = extnt.xMaximum()
y_min = extnt.yMinimum()
y_max = extnt.yMaximum()

A_bbox = QgsPoint(x_min, y_min)
B_bbox = QgsPoint(x_max, y_min)
C_bbox = QgsPoint(x_max, y_max)
D_bbox = QgsPoint(x_min, y_max)

where poly is a polygon shapefile.

Source for the Code below :- How to write Shapely geometries to shapefiles? and GitHub post

linelyr = ogr.Geometry(ogr.wkbLineString)
linelyr.AddPoint(x_min, y_min)
linelyr.AddPoint(x_min, y_max)
out_line = ogr.Geometry(ogr.wkbLineString)
out_line.AddGeometry(linelyr)
out_line.ExportToWkb

driver = ogr.GetDriverByName("Esri Shapefile")
ds = driver.CreateDataSource("outputlocation.shp")
layr1 = ds.CreateLayer('', None, ogr.wkbLineString)

layr1.CreateField(ogr.FieldDefn('id', ogr.OFTInteger))
defn = layr1.GetLayerDefn()

feat = ogr.Feature(defn)
feat.SetField('id', 1)

geom = ogr.CreateGeometryFromWkb(out_line.wkb) # or maybe just 'outline'
feat.SetGeometry(geom)

layr1.CreateFeature(feat)
ds = layr1 = feat = geom = None 

The shapefile is created, along with all the necessary files, meaning no errors. But the shapefile, when opened in QGIS has no features. Just a constructed attribute table.

Pls assist...

P. S. Any Other solution besides OGR usage also works, as long as it can be coded in the Python console of QGIS

Answer 1

Some remarks because you mix many things ( if you want yo use your script outside the console of QGIS, install the Python module GDAL (osgeo) in your Python installation: it is installed in the Python version of QGIS if you are on Windows).

1) There are no Shapely geometries in your script, only

• PyQGIS geometries QgsPoint(x_min, y_min), ...
• ogr geometries linelyr = ogr.Geometry(ogr.wkbLineString), ...

2) out_line is unnecessary because linelyr is already an ogr LineString (and your creation of out_line gives invalid geometries (empty LineStrings = LineString(LineString))

linelyr = ogr.Geometry(ogr.wkbLineString)
linelyr.AddPoint(5, 47)
linelyr.AddPoint(5, 55)
print linelyr.ExportToJson()
{ "type": "LineString", "coordinates": [ [ 5.0, 47.0, 0.0 ], [ 5.0, 55.0, 0.0 ] ] }
# and
out_line = ogr.Geometry(ogr.wkbLineString)
out_line.AddGeometry(linelyr)
print out_line.ExportToWkt()
LINESTRING EMPTY

3) In the same way, you don't need geom = ogr.CreateGeometryFromWkb(out_line.wkb) because as said, out_line or linelyr are already ogr LineStrings.

3) At the end of the process, with ogr, if you not close the resulting shapefile, it remains empty.

So:

from osgeo import ogr
# create ogr geometry
linelyr = ogr.Geometry(ogr.wkbLineString)
linelyr.AddPoint(5, 47)
linelyr.AddPoint(5, 55)
# create the shapefile
driver = ogr.GetDriverByName("Esri Shapefile")
ds = driver.CreateDataSource("outputlocation.shp")
layr1 = ds.CreateLayer('',None, ogr.wkbLineString)
# create the field
layr1.CreateField(ogr.FieldDefn('id', ogr.OFTInteger))
# Create the feature and set values
defn = layr1.GetLayerDefn()
feat = ogr.Feature(defn)
feat.SetField('id', 1)
feat.SetGeometry(linelyr)
layr1.CreateFeature(feat)
# close the shapefile
ds.Destroy()

If you want to use Shapely:

from shapely.geometry import Point, LineString
line = LineString([Point(5,47),Point(5,55)])
# conversion to ogr geometry
linelyr = ogr.CreateGeometryFromWkt(line.wkt)
print linelyr.ExportToJson()
{ "type": "LineString", "coordinates": [ [ 5.0, 47.0, 0.0 ], [ 5.0, 55.0, 0.0 ] ] }

If you want to use PyQGIS:

line = QgsGeometry.fromPolyline([QgsPoint(5,47),QgsPoint(5,55)]))
# conversion to ogr geometry
linelyr = ogr.CreateGeometryFromWkt(line.exportToWkt())
print linelyr.ExportToJson()
{ "type": "LineString", "coordinates": [ [ 5.0, 47.0 ], [ 5.0, 55.0 ] ] }