CODEWITHSUNDEEP
pythoncoredumpshapefilegdalogr
Jacques MALAPRADE
Jacques MALAPRADE

Reputation: 1013

Segmentation fault (core dumped) on SetField/SetFeature in ogr

I have been trying for the past 2 days to get a python script to complete creating a shapefile of 289995 points with attributes. The points can be created but the script does not complete the attributes. The code runs correctly until the loop: for j, p in enumerate(wCoords):(see code below-2nd loop) when after a while the segmentation fault happens. I have tried to add an if statement which halts the processes at certain intervals to see if I can find the location in the loop cycle when it happens. The loop cycles without fault until 1000 cycles, but not until 10,000 when it it stops without feedback, seeming like an endless loop. The program is to create tree points and then attach tree height attributes to the points. The code is below:

def save_shp(wCoords):
     print 'saving shapefile...'
     driver = ogr.GetDriverByName('ESRI Shapefile')
     if os.path.exists('tree_points.shp'):
         driver.DeleteDataSource('tree_points.shp')
     ds = driver.CreateDataSource('tree_points.shp')
     layer = ds.CreateLayer('trees', geom_type=ogr.wkbPoint)
     layerDefn = layer.GetLayerDefn()
     point = ogr.Geometry(ogr.wkbPoint)

     for i, p in enumerate(wCoords):
         point.AddPoint(p[0],p[1])
         featureIndex = i
         feature = ogr.Feature(layerDefn)
         feature.SetGeometry(point)
         feature.SetFID(featureIndex)
         layer.CreateFeature(feature)

     fieldDefn = ogr.FieldDefn('tree_hts', ogr.OFTReal)
     layer.CreateField(fieldDefn)
     i = feature.GetFieldIndex('tree_hts')#???

     for j, p in enumerate(wCoords):

         feature_n = layer.GetFeature(j)
         feature_n.SetField(i, p[2])#???
         layer.SetFeature(feature_n)

     try:
         ds.Destroy()
     except:
         print 'still core dumping!'

I don't know enough about gdal/ogr to give you any more info than this. Please help. Jacques

Upvotes: 1

Views: 1427

Answers (1)

Mike T
Mike T

Reputation: 43642

A few quick tips:

  • Use the same layer name as the shapefile prefix: ds.CreateLayer('tree_points', ogr.wkbPoint)
  • Do your layer.CreateField(fieldDefn) call before adding any data
  • Loop through your features once, creating geometry and feature objects at the same time

  • Inside your for loop, you need to make a new geometry object, and point index 0:

    point = ogr.Geometry(ogr.wkbPoint)
point.SetPoint_2D(0, p[0], p[1])
  • You don't need ds.Destroy(); save/close using ds = None

Upvotes: 1

Related Questions