[GIS] Distance from points to borders of polygons in QGIS big data set

distancepythonqgisvector

In QGIS I have two shapefiles: one point layer representing power plants (654 elements) and one polygon layer representing cities and villages connected to gas grid (1348 elements).

For every power plant I need to find all the areas, which are maximum in 10 km distance from the plant. What I can imagine in QGIS I need to find all the areas to whose border the distance from the point is smaller or equal 10 km.

As a result I want to obtain a table with the list of power plants, selected areas and their distances.

How can I make it in QGIS? I am familiar with Python, but never worked with PostGIS or GRASS.

Best Answer

Next PyQGIS code does the work.

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)

I tested it with shapefiles (power_plants and cities) of next image:

enter image description here

After running the code at the Python Console of QGIS, it was obtained a memory layer (named 'selected_plants') whose attributes table has the list of power plants, their distances and selected areas (expressed as 'city_id').

enter image description here

Editing Note:

Creating in each shapefile its respective field name for each feature (in my example 'name_plant' and 'name_city'), next code put in these fields of memory layer the name of attribute that corresponds.

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 = []
plant_names = []
city_names = []

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)
            p_name = feat_pp.attribute('plant_name')
            plant_names.append(p_name)
            c_name = feat_c.attribute('city_name')
            city_names.append(c_name)

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

uri = "Point?crs=epsg:" + str(epsg) + "&field=id:integer&field=distance:double(20,5)&field=plant_names:string(30)&field=city_names:string(30)""&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], plant_names[i], city_names[i]])

prov.addFeatures(selected_feats)

QgsMapLayerRegistry.instance().addMapLayer(mem_layer)

After running the above code, I got:

enter image description here

Related Question