6

I am relatively new using QGIS.

Unfortunately I can´t figure out the following yet. I have two polygon shapefiles, one is a grid I have created using QGIS 2.12.3-Lyon (Creat graticule) and the second represents native communities´ lands. I want to assign to each graticule´s cell the name and attributes of the corresponding overlapping community. Some cells overlap with more than one communal land. In these cases I want to assign the name of the community that overlaps with the cell centroid.

How can I do this?

|improve this question
  • 2
    How do you want a solution handle cases with many communities overlapping a cell, but none of them contains the cell centroid? In case of a single community the cell gets a match, in case others somehow touching the cell there would be no match at all. Assigning the community with the largest overlap may be a further criterion. – Detlev Feb 2 '16 at 18:33
  • True! Actually, it might be more useful, as you mention, to assign a community to a cell based on the largest overlap. How can I do this? – Renzo Giudice Feb 4 '16 at 14:47
1

Given 2 polygon layers. The features of the first one should receive the attributes of this feature of the second one, which has the largest overlap area.

The following Python script loops over both feature sets (brute force). If there is an overlap, it assignes the attributes if the overlap area is larger than the previous.

Spatial indexing and preselection with spatial filtering may enhance the execution speed.

# get layers
grid_lyr = QgsMapLayerRegistry.instance().mapLayersByName('grid')[0]
poly_lyr = QgsMapLayerRegistry.instance().mapLayersByName('community')[0]

# get field definitions
fields = poly_lyr.fields()

# create result layer
res_lyr = QgsVectorLayer('Polygon?crs=EPSG:4326', 'result', 'memory')
res_prov = res_lyr.dataProvider()
res_prov.addAttributes(fields)
res_lyr.updateFields()

# holds the attributes to update with
attributes = []
feats = []

# loop over grid features
for g in grid_lyr.getFeatures():
    overlap_area = -9999
    # loop over the features to assign to, if they overlap
    for poly in poly_lyr.getFeatures():
        # build overlap geometry
        overlap = g.geometry().intersection(poly.geometry())
        # if geometries overlap
        if not overlap.isEmpty():
            # and if overlap is larger than the previous one
            if overlap_area < overlap.area():
                overlap_area = overlap.area()
                # remember attributes
                attributes = poly.attributes()

    # create new feature
    feat = QgsFeature(fields)
    feat.setGeometry(g.geometry())
    feat.setAttributes(attributes)
    feats.append(feat)

# finally add features to result layer and display it 
res_prov.addFeatures(feats)
QgsMapLayerRegistry.instance().addMapLayer(res_lyr)

The test case:

Test case

|improve this answer
  • Thanks Detlev! This really seems to be a way to go, however when trying to implement this with my data I got this error at almost the beginning of the code: >>> for g in grid_lyr.getFeatures(): overlap_area = -9999 File "<input>", line 3 overlap_area = -9999 ^ IndentationError: expected an indented block – Renzo Giudice Mar 19 '16 at 14:50
  • @RenzoGiudice maybe that the leading spaces got lost while copy & pasting. The lines should be indented by multiples of 4 blanks. See the code on this page, and copy appropriate blocks of 4 blanks, if they get lost by pasting. Some editors delete leading blanks on lines, or set tabs when breaking a line. Both can cause IndentationError, – Detlev Mar 19 '16 at 15:14
  • Ok @Detlev, it worked, no problems with indentation. However, I´m getting a couple of error messages. (1): ...res_lyr.updatesFields() Traceback (most recent call last): File "<input>", line 1, in <module> AttributeError: 'QgsVectorLayer' object has no attribute 'updatesFields' and (2): res_prov.addFeatures(feats) File "<input>", line 13 res_prov.addFeatures(feats) ^ SyntaxError: invalid syntax. Everything else worked fine. I even produced a result but it is empty, as expected, no visual result, nor attributes in the attribute table. – Renzo Giudice Mar 21 '16 at 13:23
  • The first error is 'updatesFields', correct is updateFields(). I don't see what should be wrong with the second line, perhaps it's a follow up error. Can you upload your code somewhere into the cloud? – Detlev Mar 21 '16 at 14:04

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.