Good afternoon, <br>I am working on a program to convert a shapefile from ocentric to ographic latitude. The program reads an ESRI shapefile, creates a new shapefile, copying the attribute data, and performs a latitude "shift" for each point in the geometry. OGR is working wonderfully on all geometry types, except for multipolygons with donuts.<br>
<br>Below is the code snippet. Above this the geometry type and feature attributes are read from the input shapefile.<br><br>When printing the point_count (in blue) variable I am able to see that my multipolygon consists of two polygons, each with 2 rings (one interior and one exterior). I iterate through these rings and attempt to get a point count, which I use to iterate over each linearring. The point count for the exterior ring runs without a problem and returns the proper number of points. The point count on the interior ring returns an error. <br>
<br>Stepping up one "level" I can print the geometry (green text) for the multipolygon and the description of 2 polygons (in proper WKT format with an exterior and interior ring) are written. It seems that once I am "into" the linearrings of a multipolygon OGR is only "seeing" the exterior ring.<br>
<br>Thoughts?<br><br>Thanks!<br>Jay<br><br><br>Code Snippet:<br><br>elif geom_type == "MULTIPOLYGON": <br> <br> while feature_in is not None: <br> <br> #Create the output feature from the input<br>
geometry = feature_in.GetGeometryRef()<br> <br> #Dig into the input geometry<br> geom_count = geometry.GetGeometryCount()<br> <br> #Pack the multipolygon lunch for later<br> multipoly = ogr.Geometry(ogr.wkbMultiPolygon)<br>
<br> #Iterate through the polygons<br> for g in range(geom_count):<br> <span style="color: rgb(0, 102, 0);">geom = geometry.GetGeometryRef(g)</span><br> geom_count = geom.GetGeometryCount()<br>
geom_name = geom.GetGeometryName()<br> <br> #For individual polygons in the multipolygon layer<br> if geom_name == "LINEARRING":<br> <br> #Create a linear ring to populate with the points.<br>
ring = ogr.Geometry(ogr.wkbLinearRing)<br> <br> #Get a point count<br> point_count = geom.GetPointCount()<br> <br> #Iterate through the points<br>
for h in range (point_count):<br> pntx = geom.GetX(h)<br> pnty = geom.GetY(h)<br> pntz = geom.GetZ(h)<br> <br> #perform the conversion<br>
if pntz == 0.0:<br> new_y = convert(pnty)<br> ring.AddPoint_2D(pntx, new_y)<br> <br> <br> elif pntz != 0.0:<br>
new_y = convert(pnty)<br> ring.AddPoint(pntx, new_y, pntz)<br> <br> #Add the linearring <br> poly.AddGeometry(ring)<br> <br>
#Assign the spatial reference to the feature<br> poly.AssignSpatialReference(feature_in_spatref)<br> <br> #For multipolygons in the multipolygon layer<br> elif geom_name == "POLYGON":<br>
<br> #Create a polygon to populate with linear ring(s)<br> poly = ogr.Geometry(ogr.wkbPolygon)<br> <br> #Iterate through the polygons<br> for h in range (geom_count):<br>
geom = geom.GetGeometryRef(h)<br><br> <span style="color: rgb(255, 0, 0);"> #Check for multipolyons with donuts</span><br style="color: rgb(255, 0, 0);"><span style="color: rgb(255, 0, 0);"> try:</span><br style="color: rgb(255, 0, 0);">
<span style="color: rgb(255, 0, 0);"> #Get a point count</span><br style="color: rgb(255, 0, 0);"><span style="color: rgb(255, 0, 0);"> <span style="color: rgb(0, 0, 153);">point_count</span> = geom.GetPointCount()</span><br style="color: rgb(255, 0, 0);">
<span style="color: rgb(255, 0, 0);"> </span><br style="color: rgb(255, 0, 0);"><span style="color: rgb(255, 0, 0);"> #Create a linear ring to populate with the points.</span><br style="color: rgb(255, 0, 0);">
<span style="color: rgb(255, 0, 0);"> ring = ogr.Geometry(ogr.wkbLinearRing) </span><br style="color: rgb(255, 0, 0);"><span style="color: rgb(255, 0, 0);"> </span><br style="color: rgb(255, 0, 0);">
<span style="color: rgb(255, 0, 0);"> #Iterate through the points</span><br style="color: rgb(255, 0, 0);"><span style="color: rgb(255, 0, 0);"> for j in range (point_count):</span><br style="color: rgb(255, 0, 0);">
<span style="color: rgb(255, 0, 0);"> pntx = geom.GetX(j)</span><br style="color: rgb(255, 0, 0);"><span style="color: rgb(255, 0, 0);"> pnty = geom.GetY(j)</span><br style="color: rgb(255, 0, 0);">
<span style="color: rgb(255, 0, 0);"> pntz = geom.GetZ(j)</span><br style="color: rgb(255, 0, 0);"><span style="color: rgb(255, 0, 0);"> </span><br style="color: rgb(255, 0, 0);">
<span style="color: rgb(255, 0, 0);"> #perform the conversion</span><br style="color: rgb(255, 0, 0);"><span style="color: rgb(255, 0, 0);"> if pntz == 0.0:</span><br style="color: rgb(255, 0, 0);">
<span style="color: rgb(255, 0, 0);"> new_y = convert(pnty)</span><br style="color: rgb(255, 0, 0);"><span style="color: rgb(255, 0, 0);"> ring.AddPoint_2D(pntx, new_y)</span><br style="color: rgb(255, 0, 0);">
<span style="color: rgb(255, 0, 0);"> </span><br style="color: rgb(255, 0, 0);"><span style="color: rgb(255, 0, 0);"> </span><br style="color: rgb(255, 0, 0);"><span style="color: rgb(255, 0, 0);"> elif pntz != 0.0:</span><br style="color: rgb(255, 0, 0);">
<span style="color: rgb(255, 0, 0);"> new_y = convert(pnty)</span><br style="color: rgb(255, 0, 0);"><span style="color: rgb(255, 0, 0);"> ring.AddPoint(pntx, new_y, pntz)</span><br style="color: rgb(255, 0, 0);">
<span style="color: rgb(255, 0, 0);"> </span><br style="color: rgb(255, 0, 0);"><span style="color: rgb(255, 0, 0);"> #Exit gracefully if multipolygons with donuts are encountered </span><br style="color: rgb(255, 0, 0);">
<span style="color: rgb(255, 0, 0);"> except:</span><br style="color: rgb(255, 0, 0);"><span style="color: rgb(255, 0, 0);"> print "This shapefile is not supported."</span><br>
<br> #Add the linearring<br> poly.AddGeometry(ring)<br> <br> #Add the polygon to the multipolygon<br> multipoly.AddGeometry(poly)<br> <br>
#Set the spatial reference of the feature<br> multipoly.AssignSpatialReference(feature_in_spatref)<div style="visibility: hidden; display: inline;" id="avg_ls_inline_popup"></div><style type="text/css">#avg_ls_inline_popup { position:absolute; z-index:9999; padding: 0px 0px; margin-left: 0px; margin-top: 0px; width: 240px; overflow: hidden; word-wrap: break-word; color: black; font-size: 10px; text-align: left; line-height: 13px;}</style>