I am reading GeoJSON into shapely objects and need help,

I have a massive MultiPolygon with several nested Polygons in it. I believe this is why it is not considered a valid shape?

-> False

-> 'Nested shells[-118.752368927 36.63641068]'

Are the inner polygons actually the problem?

Right now we attempt to fix all invalid poly/multipolygons we get using the first buffer ring, and throw an error if that does not make the shape valid

        first_buffer = as_shape.buffer(0)
        if not first_buffer.is_valid:
            raise Exception('Invalid geo json contains invalid geometries')

        return mapping(first_buffer)

In Many cases this works, but in the case of a multipolygon with several complex holes, the MultiPolygon gets oversimplified (holes go away and we get the outer shell).

What's the right way to make all inner polygons holes while keeping the outer MultiPolygon in tact?

My attempts so far seem overly complicated --

        all_ps = list(as_shape)
        get = {
            'shells': [ g for g in all_ps if g and not len([ g.buffer(0).within(i) for i in all_ps if i!=g ]) > 2],
            'holes': [ g for g in all_ps if g and len([ g.buffer(0).within(i) for i in all_ps if i!=g ]) > 2],

2 Answers 2


You could clip out the nested shells using outer.difference(inner)

from shapely.geometry import Point, MultiPolygon
from shapely.validation import explain_validity

multi = MultiPolygon([Point(0, 0).buffer(5), Point(0, 0).buffer(2)])

multi = MultiPolygon(sorted(multi, reverse=True, key=lambda p: p.area))  # ensure largest is 1st

outer = multi[0]  # 1st part of the multipoly
inner = multi[1:]  # handle overlapping inner polys

for poly in inner:
    outer = outer.difference(poly)


>>> Nested shells[2 0]
>>> Valid Geometry


enter image description here


enter image description here

  • Nice! This is a smart approach. unfortunately the data does not reliably feature the "outer" as the first item. Wish it did, seems like it should be a part of the geojson standard. Commented Oct 12, 2021 at 20:14
  • 1
    Maybe sort by area? Largest == outer? multi = MultiPolygon(sorted(multi, reverse=True, key=lambda p: p.area))
    – user2856
    Commented Oct 12, 2021 at 20:53
  • Thats a perfect solution if we can assume there is only one "shell" - in our case we technically can't do that :( Commented Oct 13, 2021 at 17:49

I am going to answer with what worked for us, but it is also built on a key assumption about your underlying data, so I am not going to accept my own answer since YMMV. If there is some wildly clever way to solve this universally, perhaps using something I don't know about geojson, please share!

What we ended up determining is that we only got this Nested Shells issue on MultiPolygons that had 3 levels of Polygons on top of each other. What we determined was that the "holes" were being represented by the polygons that were "one level" down (covered only by the outer most geometries). So we ended up with something like this

    all_ps: MultiPolygon = shp.geoms

    solid = MultiPolygon(
            poly for poly in all_ps 
            if len([p for p in all_ps if p.covers(poly)]) % 2 != 0
    with open('fire.json', 'w') as fp:

    holes = MultiPolygon(
            poly for poly in all_ps if len([p for p in all_ps if p.covers(poly)]) % 2 == 0 

    return MultiPolygon([
            holes=[nf.exterior.coords for nf in holes if fire.covers(nf)]
        ) for fire in holes

This works for us. It MIGHT have been good to distinguish exterior solids to better capture actual Polygon holes, but we found that this works fine.

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.