# Distance from points to borders of polygons in QGIS big data set

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.

• Are you saying you want to buffer the power plants and use the buffer geometry to select the cities/villages whose borders touch those boundaries? Then you want the distance each power plant is from the boundary of each city/town within 10km? Mar 15, 2017 at 16:14
• Actually somehow I did not think about buffering but in this case it can work as a really easy solution to pick up interesting cities. But then I still need to calculate certain distance to the border of selected cities, yes. Mar 16, 2017 at 7:06

## 2 Answers

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:

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').

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:

• Thank you very much, looks like really easy solution for my problem! my question is: how is this code creating an id number of feature? In other words, how can I find out which feature in selected_plants table matches my original feature in plants layer? and the same about cities Mar 16, 2017 at 8:04
• I added an 'Editing Note' where I explain how to answer this question. I hope that it helps. By the way, I assumed that you have in attributes table for each shapefile its respective features name (in my case, of course, they are fictional). Mar 16, 2017 at 8:53
• Thank you for explanation, naming like this works. But it looks like the code is giving me only parts of the answers. It is duplicating answers, which you can see in my results here: link and giving one answer where there should be more, which you can see here. (the interested point is selected as yellow) link Mar 16, 2017 at 10:57

I have figured out the solution based on @xunilk answer, with the small change of code (in his version in my case the code was overwriting the attributes of power plants). Now it works as it should!

``````registry = QgsMapLayerRegistry.instance()

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

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

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

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:
distances.append(d)
p_name = feat_pp.attribute('ZE-nr')
plant_names.append(p_name)
c_name = feat_c.attribute('NAMN1')
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()
selected_feats = [QgsFeature() for _ in range(len(distances))]

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)
``````