I would like to know if there is a Python library that can be used to convert GeoJSON Multi-polygon to a Polygon.
Note: i want to convert .geojson file Multi-Polygon feature type to Polygon feature type.
I would like to know if there is a Python library that can be used to convert GeoJSON Multi-polygon to a Polygon.
Note: i want to convert .geojson file Multi-Polygon feature type to Polygon feature type.
from osgeo import ogr
str = '''
{ "type": "MultiPolygon",
"coordinates": [
[
[[40, 40], [20, 45], [45, 30], [40, 40]]
],
[
[[20, 35], [10, 30], [10, 10], [30, 5], [45, 20], [20, 35]],
[[30, 20], [20, 15], [20, 25], [30, 20]]
]
]
}
'''
g = ogr.CreateGeometryFromJson(str)
print "Hi! I'm a %s with an Area %s" % (g.GetGeometryName(), g.Area())
print "I have inside me %s feature(s)!\n" % g.GetGeometryCount()
for idx, f in enumerate(g):
print "I'm feature n.%s and I am a %s.\t I have an Area of %s - You can get my json repr with f.ExportToJson()" % (idx, f.GetGeometryName(),f.Area())
>>> Hi! I'm a MULTIPOLYGON with an Area 712.5
>>> I have inside me 2 feature(s)!
>>> I'm feature n.0 and I am a POLYGON. I have an Area of 87.5 - You can get my json repr with f.ExportToJson()
>>> I'm feature n.1 and I am a POLYGON. I have an Area of 625.0 - You can get my json repr with f.ExportToJson()
This is a quick way to get a list of single part geometries from a multipart feature using shapely. The buffer is in there to correct any invalid geometry in both the original feature and in the exploded parts.
import json
from shapely.geometry import shape, mapping
# Sample data
str = '''
{ "type": "MultiPolygon",
"coordinates": [
[[[40, 40], [20, 45], [45, 30], [40, 40]]],
[[[20, 35], [10, 30], [10, 10], [30, 5], [45, 20], [20, 35]],[[30, 20], [20, 15], [20, 25], [30, 20]]]
]
}
'''
g = [x.buffer(0) for x in shape(json.loads(str)).buffer(0).geoms]
What is going on here:
First, starting from the inside, json.loads() loads our geojson string into a python object. shape() then turns this into a shapely feature, a Multipolygon in this case. buffer(0) on this object corrects any invalid geometry (for example, turning a single-part bowtie into a multipolygon). geoms is a property of the resulting geometry from the buffer operation that contains a collection of all the geometries in the multipart geometry. Iterate over this ([x.buffer(0) for x in ...]), again calling buffer on each geometry, to construct a list of valid single part geometries.
What is missing here is attaching any attributes to each of the new geometry parts. If needed, I can add that in a little bit (once I get a better sample put together).
In PyQGIS, you only need to convert from Multi Part to Single Part. I tried out my approach with the Multi Part Polygon of next image:
The following code did the work:
layer = iface.mapCanvas().currentLayer()
features = layer.getFeatures()
for feature in features:
geom = feature.geometry()
print geom.exportToGeoJSON()
# check if feature geometry is multipart
if geom.isMultipart():
new_features = []
temp_feature = QgsFeature(feature)
# set the geometry of each part
for part in geom.asGeometryCollection ():
temp_feature.setGeometry(part)
new_features.append(QgsFeature(temp_feature))
for feature in new_features:
print feature.geometry().exportToWkt()
print feature.geometry().exportToGeoJSON()
After running the code at the Python Console I got:
The first print in the red rectangles is for the GeoJSON Multi-polygon and the other are for single polygons in GeoJSON format.
To test the effective conversion to single parts, I used the WKT notation for the geometries with the QuickWKT plugin (see next image):
The Multi Part geometry was effectively converted to Single Part. I can use these features to get a memory layer or a shapefile.
There are many Python modules to do that: all those which use the geointerface protocol (GeoJSON-like protocol for geo-spatial (GIS) vector data, proposed by Sean Gillies)
Shapely is one of them but we begin with PyGeoIf ("shapely ultralight", easiest to install)
from pygeoif import geometry
import geojson
str = '''
{ "type": "MultiPolygon",
"coordinates": [
[[[40, 40], [20, 45], [45, 30], [40, 40]]],
[[[20, 35], [10, 30], [10, 10], [30, 5], [45, 20], [20, 35]],[[30, 20], [20, 15], [20, 25], [30, 20]]]
]
}
'''
mp = geojson.loads(str)
Multipoly = geometry.as_shape(mp)
for poly in Multipoly.geoms:
print poly
POLYGON((40.0 40.0, 20.0 45.0, 45.0 30.0, 40.0 40.0))
POLYGON((20.0 35.0, 10.0 30.0, 10.0 10.0, 30.0 5.0, 45.0 20.0, 20.0 35.0),(30.0 20.0, 20.0 15.0, 20.0 25.0, 30.0 20.0))
# __geointerface__
for poly in Multipoly.geoms:
print poly.__geo_interface__
{'type': 'Polygon', 'coordinates': (((40.0, 40.0), (20.0, 45.0), (45.0, 30.0), (40.0, 40.0)),)}
{'type': 'Polygon', 'coordinates': (((20.0, 35.0), (10.0, 30.0), (10.0, 10.0), (30.0, 5.0), (45.0, 20.0), (20.0, 35.0)), ((30.0, 20.0), (20.0, 15.0), (20.0, 25.0), (30.0, 20.0)))}
Likewise, you could also use Shapely
from shapely.geometry import shape, mapping
Multipoly = shape(mp)
for poly in Multipoly:
print poly.wkt
POLYGON ((40 40, 20 45, 45 30, 40 40))
POLYGON ((20 35, 10 30, 10 10, 30 5, 45 20, 20 35), (30 20, 20 15, 20 25, 30 20))
for poly in Multipoly:
print mapping(poly)
{'type': 'Polygon', 'coordinates': (((40.0, 40.0), (20.0, 45.0), (45.0, 30.0), (40.0, 40.0)),)}
{'type': 'Polygon', 'coordinates': (((20.0, 35.0), (10.0, 30.0), (10.0, 10.0), (30.0, 5.0), (45.0, 20.0), (20.0, 35.0)), ((30.0, 20.0), (20.0, 15.0), (20.0, 25.0), (30.0, 20.0)))}
And you can use the QGIS QuickWKT plugin to show the results (as Xulnik)
import geojson. The example works with import json as well (at least, the one using shapely). json is included in the standard library. Apart from this (and i might be wrong), great answer (+1)!
str as variable name as it overwrites the str built-in Python function.