jueves, 16 de marzo de 2017

Selección de features de puntos con distancia mínima a features tipo polígono en PyQGIS

Estos pueden ser seleccionados mediante la función 'distance' de QgsGeometry para las sucesivas interacciones entre rasgos en un double loop. La condición de distancia mínima se incluye en el loop (en este caso 10000 m = 10 km) lo que sirve para poblar las listas con las distancias y los ids (o los nombres) de los rasgos seleccionados.



Los shapefiles a usar en la prueba se visualizan en la imagen a continuación:


El código completo es el siguiente:

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
registry = QgsMapLayerRegistry.instance()

power_plants = registry.mapLayersByName('power_plants')
cities = registry.mapLayersByName('cities')

feats_pp = [ feat for feat in power_plants[0].getFeatures() ]

feats_cities = [ feat for feat in cities[0].getFeatures() ]

selected_feats = []
distances = []
city_id = []

for feat_c in feats_cities:
    for feat_pp in feats_pp:
        d = feat_c.geometry().distance(feat_pp.geometry())
        if d <= 10000:
            selected_feats.append(feat_pp)
            distances.append(d)
            city_id.append(feat_c.id())

epsg = power_plants[0].crs().postgisSrid()

uri = "Point?crs=epsg:" + str(epsg) + "&field=id:integer&field=distance:double(20,5)&field=city_id:integer""&index=yes"

mem_layer = QgsVectorLayer(uri,
                           'selected_plants',
                           'memory')

prov = mem_layer.dataProvider()

for i, feat in enumerate(selected_feats):
    feat.setAttributes([i, distances[i], city_id[i]])

prov.addFeatures(selected_feats)

QgsMapLayerRegistry.instance().addMapLayer(mem_layer)

el cual una vez ejecutado en la Python Console de QGIS permite obtener la memory layer (en verde) de la imagen siguiente:



1 comentario:

Unknown dijo...

Jose Guerrero,agradecido por su valiosa informacion que imparte por este medio, me gustaria encaminar mis conocimientos a traves de sus orientaciones, me dedico al Area de Topografia, Drones Topograficos, recien realice un Curso introductorio a la Geomatica. mi correo es 2008arquitectura@gmail.com