2

我有一个具有多边形特征的图层。每个特征都有属性和值。我还有一个坐标列表,我想知道坐标所在的特征(或多边形)。

有人可以指导我如何去做吗?API 中是否有可以帮助我实现目标的函数,或者我应该使用一些计算几何算法自己来完成它?我知道如何做后者,但如果已经有内置函数,它会节省我一些时间。

谢谢。

4

2 回答 2

2
while provider.nextFeature(feature):
    if (feature.geometry().contains(QgsGeometry.fromPoint(QgsPoint(lon, lat)))):
        print 'Contained in feature %d' % feature.id()
于 2013-08-25T21:06:15.287 回答
1

我最终设法自己做到了。

import sys
import os
from qgis.core import *
import matplotlib.pyplot as plt
from matplotlib.path import Path
import matplotlib.patches as patches

LATITUDE = 1.29306
LONGITUDE = 103.856

QgsApplication.setPrefixPath("/usr", True)
QgsApplication.initQgis()

layer=QgsVectorLayer("/home/shubham/SMART/zones/mtz1092p.shp", "mtz1092p", "ogr")
if not layer.isValid():
    print "Layer failed to load!"
provider = layer.dataProvider()

def findFeatureId(point):
    feat = QgsFeature()
    allAttrs = provider.attributeIndexes()
    provider.select(allAttrs)
    while provider.nextFeature(feat):
        geom = feat.geometry()
        x = geom.asPolygon()
        if len(x) == 0:
            print "Feature ID %d has no ring" % feat.id()
        else:
            codes = []
            codes.append(Path.MOVETO)
            for i in range (0, len(x[0]) - 2):
                codes.append(Path.LINETO)
            codes.append(Path.CLOSEPOLY)
            path = Path(x[0], codes)
            if (path.contains_point(point, None, 0.0)):
                print "Point contained in feature ID %d" %feat.id()

if __name__ == "__main__":
    crsSrc = QgsCoordinateReferenceSystem(4326) # WGS84
    crsDest = QgsCoordinateReferenceSystem(3414)# SVY21
    xform = QgsCoordinateTransform(crsSrc, crsDest)
    pt = xform.transform(QgsPoint(LONGITUDE, LATITUDE))
    findFeatureId(pt)
于 2013-05-30T08:59:54.563 回答