1

I'm running Python 2.7 on Ubuntu 17.10, with osgeo v2.2.1 installed via apt.

My code loads osgeo and 1) tries to create a shapefile with 1 layer and 1 field, 2) create a polygon (4 points geometry), then 3) add the polygon to the shapefile's layer. Everything runs without trouble until the third function when I use layer.CreateFeature(feature) :

driver = ogr.GetDriverByName('ESRI Shapefile')
ds = driver.Open(shapefile_name,-1)
srs = osr.SpatialReference()
srs.ImportFromEPSG(4326)
layer = ds.GetLayer()
print "layer",layer,", number of features :",layer.GetFeatureCount()

This prints :

layer <osgeo.ogr.Layer; proxy of <Swig Object of type 'OGRLayerShadow *' at 0x7f79aa499c90> > , number of features : 0

I then create a geometry from the polygon I already created and passed to the function beforehand :

defn = layer.GetLayerDefn()

geom = ogr.CreateGeometryFromWkb(poly)
print "geom",geom

This prints :

geom POLYGON ((-106 24 0,-100 26 0,-103 20 0,-106 20 0))

Code then goes on to create the feature :

feat = ogr.Feature(defn)
feat.SetField('polygon_id', polygon_name)
feat.SetGeometry(geom)
print feat
print "Created feature",feat.GetField('polygon_id')

This feature is created :

<osgeo.ogr.Feature; proxy of <Swig Object of type 'OGRFeatureShadow *' at 0x7f35cd988d50> >
Created feature polygon_1

But when I add it to the layer, nothing happens :

layer.CreateFeature(feat)
print "number of features : ",layer.GetFeatureCount()

Prints :

number of features :  0

What did I miss ?

Tehem
  • 21
  • 3
  • Does it help if you flush first (`ds.FlushCache()`), or close the file and reopen? – Rutger Kassies Mar 09 '18 at 10:52
  • Thank you for your interest in my question. The function I'm creating the shapefile and its layer in ends with `driver = ds = layer = None` And the shapefile and layer are attributed right before creating the geometry (as shown) ; I tried to FlushCache the driver anyway at several places, including : `feat.SetGeometry(geom) ds.FlushCache() layer.CreateFeature(feat)` I also tried ResetReading the layer after the fist GetFeatureCount `print "layer",layer,", number of features :",layer.GetFeatureCount() layer.ResetReading() defn = layer.GetLayerDefn()` Nothing improves the situation. – Tehem Mar 10 '18 at 14:11
  • Does it changes anything if you try to change the -1 to 1 like this in the `driver.Open(shape_filename, 1)`? – Gabriella Giordano Apr 18 '18 at 17:40
  • my suggestion is to include ALL of your code...it will prevent us from guessing what is there and what is not. – user1269942 Sep 10 '18 at 23:50
  • @user1269942 This is the code, it's reproductible. – Tehem Oct 22 '18 at 20:29

1 Answers1

0

The vector tutorial in the gdal documentation provides a good example for writing a wkbPoint feature, but writing a wkbPolygon feature is a little trickier. Below is (c++) code works for me. (Note: TMultline is just an externally defined linked list of vertices.)

bool TTdhOGR_target::WriteLinePts (OGRFeature *poFeature, TMultiLine *lineParam) {
  OGRLinearRing lineString;
  OGRPolygon currPoly;
  TVertex *currPt = lineParam->get_firstPt();
  while (currPt) {
    lineString.addPoint(currPt->getX(), currPt->getY());
    currPt = (TVertex*)currPt->next();
    }
  currPoly.addRing(&lineString);
  return (poFeature->SetGeometry(&currPoly) == OGRERR_NONE);
  }
timh
  • 3
  • 1
tdhgis
  • 1
  • 1
  • Your answer could be improved with additional supporting information. Please [edit] to add further details, such as citations or documentation, so that others can confirm that your answer is correct. You can find more information on how to write good answers [in the help center](/help/how-to-answer). – Community Jan 23 '22 at 03:09