[gdal-dev] GDAL, ESRI Shapefiles and nested polygons

[Even Rouault]

Le vendredi 08 août 2014 15:19:38, Graeme Wilkie a écrit :
> Hi,
>  
> You were correct, it is returning multipolygon (9sorry, I should have
> included that information). 
> I've made the change you suggested so the code looks like
>  Int index = 0;
> OGRGeometryH hGeometry = OGR_F_GetGeometryRef(hFeature);
> int shape_type = wkbFlatten(ORG_G_GetGeometry(hGeometry);
> int num_geom = OGR_G_GetGeometryCount(hGeometry) - 1; // Returns a valid
> result int innerIndex = 0;
>  
> do
> {
>                 OGRGeometryH hRing = OGR_G_GetGeometryRef(hGeometry,
> index); OGRGeometryH tmpGeo = OGR_G_GetGeometryRef(hRing, innerIndex); Int
> num_vertices = OGR_G_GetPointCount(tmpGeo);
>                 //.
> //. Do something with the data
> //.
> Index++;
> }
> // Destroy the hfeature here
> )
>  
> If innerIndex is 0 then I can read the data for that polygon. If innerIndex
> in 1 then tmpGeo is NULL even though there are 2 inner polygons in the
> test data. 
> Have I do something stupid ?
>  
> The test data is a simple square with 2 squares inside it created using
> shapefile.py 

You still didn't include the WKT... Well I will assume that you have square T3 
included inside square T2 included inside square T1.

In that case the structure of the WKT should be ((T1,T2),T3)
i.e. a multipolygon with 2 polygons : (T1,T2) and T3
T1 is an outer ring
T2 an inner ring of T1
T3 an outer ring
That's the way polygons are modelised in SFS specification.

> If anybody has a small piece of example could it would be a great help.
>  
> Regards
>  
> Graeme
> 
> ________________________________
>  From: Even Rouault <even.rouault at spatialys.com>
> To: gdal-dev at lists.osgeo.org
> Cc: Graeme Wilkie <graeme.wilkie at btinternet.com>
> Sent: Thursday, August 7, 2014 9:26 PM
> Subject: Re: [gdal-dev] GDAL, ESRI Shapefiles and nested polygons
> 
> Le jeudi 07 août 2014 22:20:34, Graeme Wilkie a écrit :
> > Hi,
> > 
> > I'm trying to use GDAL to read ESRI shapefiles. I've got it working for
> > points, lines, polygons, polylines and multipoint but I'm having problems
> > with polygons that contain polygons in the same feature.
> > 
> > 
> > 
> > The following is a cut down code sample of the code I'm using. The
> > problem I have is that I can read the number of polygons in the feature
> > (i.e. 3 polygons, nested) using OGR_G_GetGemoetryCount. When I use
> > OGR_R_GetGeometryRef to get the handle/pointer for each of the polygons
> > in turn I get what appears to be a valid handle but when I use it to get
> > the number of points/vertices in the polygon it always returns 0.
> 
> Graeme,
> 
> You didn't display the WKT that corresponds to the geometry, so I will just
> make the guess that the geometries returned is a multipolygon and not a
> polygon.
> 
> So you might need one more level of OGR_G_GetGeometryCount and
> OGR_G_GetGeometryRef to go to the ring level.
> 
> Best regards,
> 
> Even

-- 
Spatialys - Geospatial professional services
http://www.spatialys.com