readshapefile returns a multipolygon as several polygons

[Behrouz-Babaki]

I wanted to draw a choropleth map using basemap. My plan was to read the polygons from the shapefile vis readshapefile and then further process them using shapely, descartes, etc.

The problem that I encountered was that readshapefile returned multiple polygons for each MultiPolygon. This causes confusion, as for the enclosed polygon we will have more than one record.

I think the fix will be to return a MultiPolygon instead of multiple polygons.

Here is an example:

from collections import OrderedDict
import fiona
from mpl_toolkits.basemap import Basemap


rec1 = dict()
rec2 = dict()

rec1['id'] = 0
rec1['type'] = 'Feature'
rec1['properties'] = OrderedDict([(u'NAME', 'outer'),(u'VALUE', 100)])
rec1['geometry'] = dict()
rec1['geometry']['type'] = 'MultiPolygon'
rec1['geometry']['coordinates'] = [[[(-50,40),(50,40),
                                    (50,-40),(-50,-40),
                                    (-50,40)
                                   ]],
                                  [[(-45,35),(45,35),
                                    (45,-35),(-45,-35),
                                    (-45,35)
                                   ]]]

rec2['id'] = 1
rec2['type'] = 'Feature'
rec2['properties'] = OrderedDict([(u'NAME', 'inner'),(u'VALUE', 1)])
rec2['geometry'] = dict()
rec2['geometry']['type'] = 'Polygon'
rec2['geometry']['coordinates'] = [[(-45,35),(45,35),
                                    (45,-35),(-45,-35),
                                    (-45,35)
                                   ]]
driver = 'ESRI Shapefile'
crs = {'no_defs': True, 'ellps': 'WGS84', 'datum': 'WGS84', 'proj': 'longlat'}
schema = {'geometry': 'Polygon',
         'properties': {'NAME': 'str',
                'VALUE': 'float:15.2'}}

with fiona.open(
     '/tmp/foo.shp',
     'w',
     driver=driver,
     crs=crs,
     schema=schema) as c:
    c.write(rec1)
    c.write(rec2)

m = Basemap()
m.readshapefile('/tmp/foo','foo')
for info,region  in zip(m.foo_info, m.foo):
    print info
    print region

for which the output is:

{'RINGNUM': 1, 'NAME': 'outer', 'VALUE': 100.0, 'SHAPENUM': 1}
[(-50.0, 40.0), (50.0, 40.0), (50.0, -40.0), (-50.0, -40.0), (-50.0, 40.0)]
{'RINGNUM': 2, 'NAME': 'outer', 'VALUE': 100.0, 'SHAPENUM': 1}
[(-45.0, 35.0), (-45.0, -35.0), (45.0, -35.0), (45.0, 35.0), (-45.0, 35.0)]
{'RINGNUM': 1, 'NAME': 'inner', 'VALUE': 1.0, 'SHAPENUM': 2}
[(-45.0, 35.0), (45.0, 35.0), (45.0, -35.0), (-45.0, -35.0), (-45.0, 35.0)]