Skip to content
Permalink
Browse files

Faster PointsInPolygon(+Unique/Weighted) algorithms

- Avoid use of seperate feature requests for every point
- Use GEOS prepared geometries when testing for point containment

Quick stats:

~1500 point layer

BEFORE: 17 seconds
AFTER: 3 seconds

~900k point layer

BEFORE: 30 mins = canceled at 20%
AFTER: 2.5 mins = 100% complete
  • Loading branch information
nyalldawson committed Nov 13, 2015
1 parent 7c8177e commit 596b56c44130f44a2d3d4e0c43c231d7a92d921d
@@ -76,26 +76,25 @@ def processAlgorithm(self, progress):
geom = QgsGeometry()

current = 0
hasIntersections = False

features = vector.features(polyLayer)
total = 100.0 / float(len(features))
for ftPoly in features:
geom = ftPoly.geometry()
engine = QgsGeometry.createGeometryEngine(geom.geometry())
engine.prepareGeometry()

attrs = ftPoly.attributes()

count = 0
hasIntersections = False
points = spatialIndex.intersects(geom.boundingBox())
if len(points) > 0:
hasIntersections = True

if hasIntersections:
for i in points:
request = QgsFeatureRequest().setFilterFid(i)
ftPoint = pointLayer.getFeatures(request).next()
tmpGeom = QgsGeometry(ftPoint.geometry())
if geom.contains(tmpGeom):
request = QgsFeatureRequest().setFilterFids(points)
fit = pointLayer.getFeatures(request)
ftPoint = QgsFeature()
while fit.nextFeature(ftPoint):
tmpGeom = ftPoint.geometry()
if engine.contains(tmpGeom.geometry()):
count += 1

outFeat.setGeometry(geom)
@@ -80,26 +80,25 @@ def processAlgorithm(self, progress):
geom = QgsGeometry()

current = 0
hasIntersections = False

features = vector.features(polyLayer)
total = 100.0 / float(len(features))
for ftPoly in features:
geom = ftPoly.geometry()
engine = QgsGeometry.createGeometryEngine(geom.geometry())
engine.prepareGeometry()

attrs = ftPoly.attributes()

classes = []
hasIntersections = False
points = spatialIndex.intersects(geom.boundingBox())
if len(points) > 0:
hasIntersections = True

if hasIntersections:
for i in points:
request = QgsFeatureRequest().setFilterFid(i)
ftPoint = pointLayer.getFeatures(request).next()
request = QgsFeatureRequest().setFilterFids(points)
fit = pointLayer.getFeatures(request)
ftPoint = QgsFeature()
while fit.nextFeature(ftPoint):
tmpGeom = QgsGeometry(ftPoint.geometry())
if geom.contains(tmpGeom):
if engine.contains(tmpGeom.geometry()):
clazz = ftPoint.attributes()[classFieldIndex]
if clazz not in classes:
classes.append(clazz)
@@ -86,27 +86,25 @@ def processAlgorithm(self, progress):
geom = QgsGeometry()

current = 0
hasIntersections = False

features = vector.features(polyLayer)
total = 100.0 / float(len(features))
for ftPoly in features:
geom = ftPoly.geometry()
engine = QgsGeometry.createGeometryEngine(geom.geometry())
engine.prepareGeometry()

attrs = ftPoly.attributes()

count = 0
hasIntersections = False
points = spatialIndex.intersects(geom.boundingBox())
if len(points) > 0:
hasIntersections = True

if hasIntersections:
progress.setText(unicode(len(points)))
for i in points:
request = QgsFeatureRequest().setFilterFid(i)
ftPoint = pointLayer.getFeatures(request).next()
request = QgsFeatureRequest().setFilterFids(points)
fit = pointLayer.getFeatures(request)
ftPoint = QgsFeature()
while fit.nextFeature(ftPoint):
tmpGeom = QgsGeometry(ftPoint.geometry())
if geom.contains(tmpGeom):
if engine.contains(tmpGeom.geometry()):
weight = unicode(ftPoint.attributes()[fieldIdx])
try:
count += float(weight)

0 comments on commit 596b56c

Please sign in to comment.
You can’t perform that action at this time.