1

I have a workflow to eliminate small polygons within a layer and merge them with adjacent polygons. My input is the result of a polygonization of a raster (has no overlaps) and it has one field "class".

Input:

enter image description here

In the first step I loop over all features and if a polygon is of class 0 (blue), and is too small, i change its class to a suitable neighbour.

enter image description here

Next I dissolve the features, according to their changed class values.

enter image description here

As you can see, the blue areas are dissolved successfully, however the island in the middle appears to have changed value. Actually there are now two polygons laying on top of each other there, a light green one, and the original one. The light green basically filled a hole during the dissolve process. This happens only for a couple of polygons, and not all of them.

Here's the code to change the values of the polygons:

def KFB_for_1_class (input, field_, c, area_threshold):
    feats = []
    dis_temp = r"z:\dis.shp"

    layer = iface.addVectorLayer(input,"","ogr")
    layer.startEditing()

    layer_cpy = iface.addVectorLayer(input,"","ogr")

    print "cleaning class "+str(c)
    for f in layer.getFeatures():
        if f[field_] == c and f["area"] < area_threshold:            
            org_int = f[field_]
            candidates = []
            for f2 in layer_cpy.getFeatures():
                if f.geometry().touches(f2.geometry()):
                    if f2[field_] not in candidates and f2[field_] != org_int:
                        candidates.append(f2[field_])

            candidates.sort()
            if len(candidates) >0:
                f[field_] = candidates[0]
            for i in (candidates):
                if abs(org_int - i) <= abs(f[field_] - org_int):
                    f[field_] = i
                    f["area"] = area_threshold + 1

            layer.updateFeature(f)
            feats.append(f["id"])

            print "feature " +str(f["id"]) + " changed from " + str(org_int) + " to " + str(f["int"])

    print feats


    layer.commitChanges()

To dissolve I use this:

processing.runalg("qgis:dissolve", input, False , "class", dis_temp)

What I tried:

Both did not work. The different tools all ended up creating the same error.

How can I resolve this?

PolyGeo
  • 65,136
  • 29
  • 109
  • 338
Ben the bear
  • 455
  • 4
  • 10
  • You should paste more of your code snippet. It's difficult to tell anything or try to reproduce your workflow there... – gisnside Sep 29 '17 at 13:33
  • Are you able to confirm there was only a single polygon present the source data within the problem area? Shot in the dark, sometimes I try save as new shapefile and then run the process again on a new file. This method has cleared up some interesting challenges, the odds are similar to roulette. – user92055 Sep 29 '17 at 13:49
  • @gisnside I added the code that alters the values in the class field. – Ben the bear Sep 29 '17 at 14:02
  • @user92055 I am able to confirm, that both in the original data (polygonization of raster) and in the altered dataset (see code above), only 1 polygon is present. Only after the dissolve there appear to be 2. – Ben the bear Sep 29 '17 at 14:05
  • I'm unable to get this to run because feats is not defined – user92055 Sep 29 '17 at 15:38
  • @user92055 feats is just an empty list, I added it in the code. You need your layer to have a field "area" and "id" besides the class field for the code to run. However none of those take part in the problem with the dissolve tool. – Ben the bear Sep 29 '17 at 18:43
  • well i'm not much help on this because it works for me. I ran it on some random landuse data and am unable to replicate the error you are getting. Any way you can post the data from your screenshots? – user92055 Sep 29 '17 at 19:21
  • 1
    Is it possible you have multipart feature involved somewhere ? Have a try single-partying your data (Multipart to singleparts in QGIs processing window) before working your data out – gisnside Oct 01 '17 at 21:50
  • To check for multipart features is a good idea. The workflow that leads up to the problem does not create any though and I know that. I still tried it and when I ran the multipart to singlepart tool, I get the exact same amount of records. Interstingly enough it seems, that the dissolve works afterwards though. Maybe that tool repairs something funky within the data. – Ben the bear Oct 02 '17 at 07:48
  • Another thing I tried to get a bit closer to the problem, is extract the classes individually and dissolve them seperately. That way the issue does not appear to happen either. So it must be some sort of topology handling the dissolve tools cannot do with the data set. – Ben the bear Oct 02 '17 at 07:50

0 Answers0