PyQGIS QgsPolygon Issue – How to Fix QgsPolygon: CurvePolygon EMPTY Error When Extracting Largest MultiPolygon Part in PyQGIS

areamulti-polygonmultipartpyqgissplitting

There is a shapefile with MultiPolygons, see the image below

input

With this code, I am trying to get the largest part (by area) of each MultiPolygon feature. It utilizes the constGet() method.

layer = iface.activeLayer()

def get_largest_part(feat):
    parts = [part for part in feat.geometry().constGet()]
    largest_part = max(parts, key=lambda x: x.area())
    return largest_part

largest_parts = [get_largest_part(feat) for feat in layer.getFeatures()]

print(largest_parts)

However, the output looks a bit weird. Can somebody explain to me why?

[<QgsPolygon: CurvePolygon EMPTY>, <QgsPolygon: CurvePolygon EMPTY>]

Sometimes when running these both codes my QGIS 3.18.2-Zürich on Windows 8.1 even crashes. The same story happened when using parts() and constGet() on Ubuntu 21.10.

I tried to apply the solution suggested by @KadirŞahbaz in this thread: Getting area of part of MultiPolygon in QGIS, that utilizes the parts() method

layer = iface.activeLayer()

def get_largest_part(feat):
    parts = [part for part in feat.geometry().parts()]
    largest_part = sorted(parts, key=lambda x: x.area())[0]
    return largest_part
    
largest_parts = [get_largest_part(feat) for feat in layer.getFeatures()]
    
print(largest_parts)

however, the output is the same:

[<QgsPolygon: CurvePolygon EMPTY>, <QgsPolygon: CurvePolygon EMPTY>]

Best Answer

I don't know why you get empty polygons, but you can use the following script:

layer = iface.activeLayer()

def get_largest_part(feat):
    geoms = feat.geometry().asGeometryCollection()
    geom = sorted(geoms, key=lambda g: g.area(), reverse=True)[0]
    return QgsGeometry(geom)

largest_parts = [get_largest_part(feat) for feat in layer.getFeatures()]

enter image description here