为QGIS 2空间连接(PyQGIS)创建空间索引



我编写了一些代码来在QGIS 2和2.2中执行一个简单的空间连接(位于缓冲区内的点以获取缓冲区的属性)。但是,我想使用QgsSpatialIndex来加快速度。从这里我可以去哪里:

pointProvider = self.pointLayer.dataProvider()
rotateProvider = self.rotateBUFF.dataProvider()
all_point = pointProvider.getFeatures()
point_spIndex = QgsSpatialIndex()
for feat in all_point:
    point_spIndex.insertFeature(feat)
all_line = rotateProvider.getFeatures()
line_spIndex = QgsSpatialIndex()
for feat in all_line:
    line_spIndex.insertFeature(feat)
rotate_IDX = self.rotateBUFF.fieldNameIndex('bearing')
point_IDX = self.pointLayer.fieldNameIndex('bearing')
self.pointLayer.startEditing()
for rotatefeat in self.rotateBUFF.getFeatures():
    for pointfeat in self.pointLayer.getFeatures():
        if pointfeat.geometry().intersects(rotatefeat.geometry()) == True:
            pointID = pointfeat.id()
bearing = rotatefeat.attributes()[rotate_IDX]
self.pointLayer.changeAttributeValue(pointID, point_IDX, bearing)
self.pointLayer.commitChanges()

要进行这种空间连接,您可以使用QgsSpatialIndex (http://www.qgis.org/api/classQgsSpatialIndex.html) intersects(QgsRectangle)函数来获取候选特征id的列表,或者使用nearestNeighbor (QgsPoint,n)函数来获取n个最近邻居作为特征id的列表。

由于您只想要位于缓冲区内的点,因此intersects函数似乎是最合适的。我没有测试如果一个退化bbox(点)可以使用。如果没有,就在你的点周围画一个非常小的边界框。

intersects函数返回具有与给定矩形相交的边界框的所有特征,因此您必须测试这些候选特征以获得真正的相交。

你的外部循环应该在点上(你想从它们的包含缓冲区中为每个点添加属性值)。

# If degenerate rectangles are allowed, delta could be 0,
# if not, choose a suitable, small value
delta = 0.1
# Loop through the points
for point in all_point:
    # Create a search rectangle
    # Assuming that all_point consist of QgsPoint
    searchRectangle = QgsRectangle(point.x() - delta, point.y()  - delta, point.x() + delta, point.y() + delta)
    # Use the search rectangle to get candidate buffers from the buffer index
    candidateIDs = line_index.intesects(searchRectangle)
    # Loop through the candidate buffers to find the first one that contains the point
    for candidateID in candidateIDs:
        candFeature == rotateProvider.getFeatures(QgsFeatureRequest(candidateID)).next()
        if candFeature.geometry().contains(point):
            # Do something useful with the point - buffer pair
            # No need to look further, so break
            break

最新更新