3

I'm trying to extract data from an OSM.PBF file using the python GDAL/OGR module.

Currently my code looks like this: 

import gdal, ogr

osm = ogr.Open('file.osm.pbf')

## Select multipolygon from the layer
layer = osm.GetLayer(3) 
# Create list to store pubs
pubs = []
for feat in layer:
    if feat.GetField('amenity') == 'pub':
         pubs.append(feat)

While this little bit of code works fine with small.pbf files (15mb). However, when parsing files larger than 50mb I get the following error:

 ERROR 1: Too many features have accumulated in points layer. Use OGR_INTERLEAVED_READING=YES MODE

When I turn this mode on with:

gdal.SetConfigOption('OGR_INTERLEAVED_READING', 'YES')

ogr does not return any features at all anymore, even when parsing small files.

Does anyone know what is going on here?

S Hook
  • 71
  • 1
  • 5
  • According to http://gdal.org/1.11/ogr/drv_osm.html interleaving also requires a special reading pattern. Here is a python example for this special reading pattern: https://lists.osgeo.org/pipermail/gdal-dev/2014-April/038633.html I've never used GDAL/OGR so no idea about the details. If it works feel free to post a more useful answer with example code. – scai Feb 16 '16 at 18:19

1 Answers1

4

Thanks to scai's answer I was able to figure it out.

The special reading pattern required for interleaved reading that is mentioned in gdal.org/1.11/ogr/drv_osm.html is translated into a working python example that can be found below.

This is an example of how to extract all features in an .osm.pbf file that have the 'amenity=pub' tag

import gdal, ogr

gdal.SetConfigOption('OGR_INTERLEAVED_READING', 'YES')
osm = ogr.Open('file.osm.pbf')

# Grab available layers in file
nLayerCount = osm.GetLayerCount()

thereIsDataInLayer = True

pubs = []

while thereIsDataInLayer:

    thereIsDataInLayer = False

    # Cycle through available layers
    for iLayer in xrange(nLayerCount):

        lyr=osm.GetLayer(iLayer)

        # Get first feature from layer
        feat = lyr.GetNextFeature()

        while (feat is not None):

             thereIsDataInLayer = True

             #Do something with feature, in this case store them in a list
             if feat.GetField('amenity') == 'pub':
                 pubs.append(feat)

             #The destroy method is necessary for interleaved reading
             feat.Destroy()

             feat = lyr.GetNextFeature()

As far as I understand it, a while-loop is needed instead of a for-loop because when using the interleaved reading method, it is impossible to obtain the featurecount of a collection.

More clarification on why this piece of code works like it does would be greatly appreciated.

S Hook
  • 71
  • 1
  • 5