3

I have a vector shape with many different sized polygons each with an other id and values. I want to assign a new table with Information of each polygon containing the name of the majority (higher proportion) of adjacent polygons attribute name. Ergo some kind of polygon neighborhood analysis.

The query should look to the attribute of each bordering polygon and meassure the length of neighborhood polygon and id.

Does anyone have an idea or a tool how to accomplish this task?

I intend to use QGIS, or postgis. I know there is a tool in arcgis called neighborhood statistic.

  • The ArcGIS Neighborhood statistics tools don't work like that - they are raster tools that use a 'neighborhood' to search an underlying raster in which those raster cells intersect the neighborhood - then the summary of the intersecting raster values is assigned to the neighborhood. In your case, you may have a polygon bordered by 5 other polygons - you might want to simply select the surrounding features and take an average (?) of the values of those 'neighbors'? – DPSSpatial Feb 7 '17 at 16:57
  • In ArcGis it is in the analysis toolbox / Proximity toolset / Polygon Neighbors. There you can make some neighborhood analysis. For example the length of coincident edges (edge neighbors). The output is a table with all the relations of each polygon. – Peter Feb 9 '17 at 10:54
3

This blog post has the answer I think you are looking for.

In this case it looks at neighboring polygons; lists their ids, and sums up a defined attribute of them.
You might want to make a copy of the original data and then run the code in the python console of QGIS.

NOTE:
I noticed an error in his code on line 52 here is the corrected if statement:

    # these conditions. So if a feature is not disjoint, it is a neighbor.
    if (f != intersecting_f and not intersecting_f.geometry().disjoint(geom)):
        neighbors.append(intersecting_f[_NAME_FIELD])
        neighbors_sum += intersecting_f[_SUM_FIELD]

Sample Output: Sample output


#UPDATE#

I've updated to reflect your comment.

If you load the layer into QGIS , edit the variables and output file locations, and run the following, you will get:

  1. A layer for each shared edge
    (showing: id, touches, length, sumcolumn)
  2. A table with a row for each polygon (stats summed)
    (showing: id, length of edges summed, sumcolumn summed)

Sample Output:
1.Shared Edges: 1.shared Edges
2.Sums per polygon:

Sums per polygon

'''
FOR EACH POYGON OF INPUT: 
SUM the length of touching sides
SUM a specified column of touching polygons
'''
import itertools
import math
from PyQt4.QtCore import QVariant

#SPECIFY VARIABLES
polycode = 'ID'
polysum = 'SumCol'
#SPECIFY EPSG
epsg = '28355'
#SPECIFY OUTPUTS
outputedges = r"C://temp//DELETEME//t_line_2.shp"
outputedges_sum = r"C://temp//DELETEME//t_line_2_sum.shp"

#USE TJHE LOADED LAYER:
layer = iface.activeLayer()

#CREATE LINE MEMORY LAYER
vlline = QgsVectorLayer("LineString?crs=EPSG:%s" % epsg, "temporary_lines", "memory")
prline = vlline.dataProvider()
# add fields
prline.addAttributes([QgsField('ID', QVariant.String), QgsField('Touches', QVariant.String), QgsField('length', QVariant.Double), QgsField('tosum', QVariant.Double)])
vlline.updateFields() # tell the vector layer to fetch changes from the provid

#LIST OF EXTERIOR RINGS, ALONG WITH PECIFIED COLUMNS
rings = [(QgsGeometry.fromPolyline(elem.geometry().asPolygon()[0]), elem[polycode], elem[polysum])for elem in layer.getFeatures()]


#CREATE TABLE FOR EACH EDGE (SHOW: id, touches, length, sumcolumn)
for i in  itertools.permutations(rings, 2):
    print 'Working on ID: '+ str(i[0][1])
    if i[0][0].intersects(i[1][0]):
        len = i[0][0].intersection(i[1][0]).length()
        #print len
        fetline = QgsFeature()
        fetline.setGeometry(i[0][0].intersection(i[1][0]))
        fetline.setAttributes([str(i[0][1]),str(i[1][1]),len,(i[1][2])])
        prline.addFeatures([fetline])
        vlline.updateExtents()

#WRITE EDGES TABLES
QgsVectorFileWriter.writeAsVectorFormat(vlline,outputedges,"utf-8",None,"ESRI Shapefile")
iface.addVectorLayer(outputedges, "hello_edges", "ogr")

#CREATE Vlayer of summed data (grouped by ID)
vlayer = QgsVectorLayer( "?query=SELECT ID, SUM(length) sum_length, SUM(tosum) sum_tosum  FROM hello_edges GROUP BY ID ", "vlayer", "virtual" )
querylayer = vlayer.dataProvider()
QgsVectorFileWriter.writeAsVectorFormat(vlayer,outputedges_sum,"utf-8",None,"ESRI Shapefile")
iface.addVectorLayer(outputedges_sum, "hello_sums", "ogr")
iface.showAttributeTable(iface.activeLayer());
|improve this answer|||||
  • Thank you, but I need a table with the length of coincident edges (edge neighbors) for each polygon and one value, like population in a table. – Peter Feb 9 '17 at 11:05
  • I've updated the answer to reflect your comment, see above. I think this is what you are after, if not perhaps edit the question to add a some more detail/explanation. – Cushen Feb 10 '17 at 4:56
0

I fixed my problem with the following code in postgis. Select a polygon layer, and be aware of the projection. With this code you get a table with the length of intersections for all pairs of polygons. If you also want some values for calculation like sum() and group by it is also possible.

SELECT a.gid AS gid_a, b.gid AS gid_b, 
ST_Length(ST_CollectionExtract(ST_Intersection(a.geom, b.geom), 2)),
b.value --any column from the table, join or calculation
FROM mypoly a, mypoly b
WHERE a.gid != b.gid AND ST_Touches(a.geom, b.geom);

This code is as good as "Polygon Neighbors" in the ArcGis Toolbox, only the performance could be better, on my iMac 2009 with PostGis, Qgis. Maybe there are faster cpu's or a faster OS :-)

  1. Sample Polygons: all neighborhood lines where calculated seperatly
  2. Sample Output from PostGis:

the calculation was made with district polygons.

|improve this answer|||||

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.