How to speed up PyQGIS code for finding intersection of features in the same layer

Question

I an deleting one of the intersecting polygons from a vector layer using PyQGIS code mentioned below. It is too slow as I have many features in the layer. Is there a way to speed it up?

index= QgsSpatialIndex(layer.getFeatures())
request = QgsFeatureRequest()
clause = QgsFeatureRequest.OrderByClause('id', ascending=True)
orderby = QgsFeatureRequest.OrderBy([clause])
request.setOrderBy(orderby)
for feat in layer.getFeatures(request):
    feats=[f for f in layer.getFeatures(request)]
    for f in feats:
        if index.intersects(f.geometry().boundingBox()):
            layer.select(feat.id())
 selIDs=layer.selectedFeatureIds()
 layer.startEditing()
 layer.deleteFeatures(selIDs1)
 layer.commitChanges()

Answer 1

@GermánCarrillo already provided a pretty useful link to a great answer. Here is something more specific about your script, where there are a few issues:

1. Since you are looking for intersecting features, there is no need to sort the layer by ids (or is there something you did not mention?)
2. You are iterating over the whole layer for each feature and checking the index (incorrectly) afterwards, so there is pretty no performance gain in the index usage. Instead, for each feature, iterate over the index only.
3. You can store the features geometries directly inside the index by using flags=QgsSpatialIndex.FlagStoreFeatureGeometries so you can avoid an inperformant feature request and instead get the geometry from the index much faster
4. Dont forget to exclude self-intersections; otherwise your entire layer will be deleted
5. If the index matches this does not mean that the features actually intersect, that only means the bounding box intersects, so you need to do an actual intersects-test on the matches.
6. There is no need to do a selection if you only want to delete the features

Try this:

layer = iface.activeLayer()
index = QgsSpatialIndex(layer.getFeatures(), flags=QgsSpatialIndex.FlagStoreFeatureGeometries)
featids = []
for feat in layer.getFeatures():
    for fid in index.intersects(feat.geometry().boundingBox()): # iterate over the index-matches. The index returns the IDs of the features where the boundingbox intersects
        if fid == feat.id(): # ignore self-intersections
            continue
        fgeom = index.geometry(fid) # get the feature geometry by the id from the index
        if fgeom.intersects(feat.geometry()): # now check if not only the bounding box intersects, but if the actual features geometries intersect
            featids.append(feat.id()) # if so append the id to a list for selection afterwards
layer.startEditing()
layer.deleteFeatures(featids)
layer.commitChanges()