[GIS] Are the functions in OGRPolygon accessible through the python bindings

gdalogrosgeopython

I'm currently trying to use the OGR python bindings to process some shapefiles. I want to obtain the number of interior rings in a polygon feature. Unfortunately I'm unable to use the documented function to get that information.

Here is my code:

from osgeo import ogr
field = ogr.Open("testExtract.shp")
lyr = field.GetLayer("testExtract")
feat = lyr.GetNextFeature()

geom = feat.GetGeometryRef()
featName = geom.GetGeometryName()
print featName
geoCount = geom.GetGeometryCount()
print "geoCount: " + str(geoCount)

print geom.getNumInteriorRings()
print geom.GetNumInteriorRings()
print geom.getNumInteriorRing()
print geom.GetNumInteriorRing()

This is what those print commands output:

POLYGON
geoCount: 6

AttributeError: 'Geometry' object has no attribute 'getNumInteriorRings'
AttributeError: 'Geometry' object has no attribute 'GetNumInteriorRings'
AttributeError: 'Geometry' object has no attribute 'getNumInteriorRing'
AttributeError: 'Geometry' object has no attribute 'GetNumInteriorRing'

Now, the only polygon in that shapefile has interior rings. And the python code is telling me that geom is a polygon with 6 interior and exterior rings, which is correct. And the documentation, http://www.gdal.org/ogr/classOGRPolygon.html#a3b996195adcf9fcd4f33570a37753dde, is saying getNumInteriorRings() is a valid function for an OGRPolygon.

So is the problem that GetNumInteriorRings() is not accessible through the python API? The python documention doesn't mention OGRPolygon as a class, just Geometry. http://gdal.org/python/osgeo.ogr.Geometry-class.html

Edit: I realize that each polygon only has one exterior polygon, so I can get that information. But I'm also interested in a couple of the other functions in OGRPolygon such as getInteriorRing (int) and getExteriorRing () so I am curious if those public member functions for OGRPolygon are accessible in Python.

Best Answer

When you start with a Python module, there are several solutions to find the available functions. One of them is dir:

geom = feat.GetGeometryRef()
print dir(geom)
['AddGeometry', 'AddGeometryDirectly', 'AddPoint', 'AddPoint_2D', 'Area', 'AssignSpatialReference', 'Boundary', 'Buffer', 'Centroid', 'Clone', 'CloseRings', 'Contains', 'ConvexHull', 'Crosses', 'Destroy', 'Difference', 'Disjoint', 'Distance', 'Empty', 'Equal', 'Equals', 'ExportToGML', 'ExportToJson', 'ExportToKML', 'ExportToWkb', 'ExportToWkt', 'FlattenTo2D', 'GetArea', 'GetBoundary', 'GetCoordinateDimension', 'GetDimension', 'GetEnvelope', 'GetEnvelope3D', 'GetGeometryCount', 'GetGeometryName', 'GetGeometryRef', 'GetGeometryType', 'GetPoint', 'GetPointCount', 'GetPoint_2D', 'GetPoints', 'GetSpatialReference', 'GetX', 'GetY', 'GetZ', 'Intersect', 'Intersection', 'Intersects', 'IsEmpty', 'IsRing', 'IsSimple', 'IsValid', 'Length', 'Overlaps', 'PointOnSurface', 'Segmentize', 'SetCoordinateDimension', 'SetPoint', 'SetPoint_2D', 'Simplify', 'SimplifyPreserveTopology', 'SymDifference', 'SymmetricDifference', 'Touches', 'Transform', 'TransformTo', 'Union', 'UnionCascaded', 'Within', 'WkbSize', '__class__', '__del__', '__delattr__', '__dict__', '__doc__', '__format__', '__getattr__', '__getattribute__', '__hash__', '__init__', '__iter__', '__module__', '__new__', '__reduce__', '__reduce_ex__', '__repr__', '__setattr__', '__setstate__', '__sizeof__', '__str__', '__subclasshook__', '__swig_destroy__', '__swig_getmethods__', '__swig_setmethods__', '__weakref__', 'next', 'this']

and no getNumInteriorRings: apparently this method is not exposed in the python bindings.
But if you look at [gdal-dev] Polygon topology:

In fact you have to use the Geometry.GetGeometryCount() that returns 1 (the exterior ring) + the number of interior rings. So polygon.GetGeometryCount() - 1 should return the number of interior rings

So a polygon only:

geom.GetGeometryCount() 
1

For a polygon whith one hole:

enter image description here

geom.GetGeometryCount()
2 # -> one hole

and 

geom.GetGeometryRef(0) # -> exterior ring = LinearRing, in red
geom.GetGeometryRef(1) # -> interior ring = LinearRing

enter image description here

exterior = Geometry(ogr.wkbPolygon)
exterior.AddGeometry(geom.GetGeometryRef(0))

enter image description here

interior = Geometry(ogr.wkbPolygon)
exterior.AddGeometry(geom.GetGeometryRef(1))

enter image description here

But if you look at the geo_interface protocol (a python object property that returns a geojson dictionary, look at Python Geo_interface applications):

 print geom.ExportToJson()
 { "type": "Polygon", "coordinates": [ [ [ 1149490.109727980103344, 691044.609108003089204 ], [ 1191579.109752570046112, 691044.609108003089204 ], [ 1191579.109752570046112, 648030.576115839648992 ], [ 1149490.109727980103344, 648030.576115839648992 ], [ 1149490.109727980103344, 691044.609108003089204 ] ], [ [ 1154115.274565846892074, 686419.444270136067644 ], [ 1154115.274565846892074, 653118.257437493419275 ], [ 1165678.186660513980314, 653118.257437493419275 ], [ 1165678.186660513980314, 686419.444270136067644 ], [ 1154115.274565846892074, 686419.444270136067644 ] ] ] }

it is easier with fiona (also based on GDAL/OGR)

import fiona
layer = fiona.open("myshape.shp")
feat = layer.next()
print feat
{'geometry': {'type': 'Polygon', 'coordinates': [[(1149490.1097279801, 691044.60910800309), (1191579.10975257, 691044.60910800309), (1191579.10975257, 648030.57611583965), (1149490.1097279801, 648030.57611583965), (1149490.1097279801, 691044.60910800309)], [(1154115.2745658469, 686419.44427013607), (1154115.2745658469, 653118.25743749342), (1165678.186660514, 653118.25743749342), (1165678.186660514, 686419.44427013607), (1154115.2745658469, 686419.44427013607)]]}, 'type': 'Feature', 'id': '0', 'properties': OrderedDict([(u'FID', 0.0)])}
print feat['geometry']['type'])
Polygon
print len(feat['geometry']['coordinates']
2
# outer ring
print feat['geometry']['coordinates'][0]
[(1149490.1097279801, 691044.60910800309), (1191579.10975257, 691044.60910800309), (1191579.10975257, 648030.57611583965), (1149490.1097279801, 648030.57611583965), (1149490.1097279801, 691044.60910800309)]
# inner ring
print feat['geometry']['coordinates'][1]
[(1154115.2745658469, 686419.44427013607), (1154115.2745658469, 653118.25743749342), (1165678.186660514, 653118.25743749342), (1165678.186660514, 686419.44427013607), (1154115.2745658469, 686419.44427013607)]

Related Question