Compute zonal statistics for portion of polygons

Question

I am working with a shapefile containing the district polygons of the political map of India and a raster file containing the elevation data of the same country. My goal is to compute the average and variance altitude within 10 km from the borders of each district. Using Zonal Statistics I was able to compute the statistics I am looking for at the district level, the further step I need is disregard the pixels lying in the central part.

I have tried to get rid of the internal parts of the districts with add rings, but the user should draw the internal polygons by himself. This is not what I am looking for: I want to have an internal polygon with exactly the same angles as the original one but with smaller sides, computed automatically for all India.

How can I create, alternatively, a mask in the shapefile for the area that I don't want (or that I want, the one that most fits the case) to take into account and then use it to compute mean and variance with the Zonal Statistics plugin?

The ideal answer would provide me a PyQGIS code but it's fine also using icons in QGIS.

Answer 1

I tried out the next code with a polygon shapefile (as active layer) where, first, it is produced a "negative" buffer (-10000 m) as memory layer. Afterward, this memory layer is used as difference layer for the active layer for producing a polygon without the internal part (it is obtained a central ring).

layer = iface.activeLayer()

feat = layer.getFeatures().next()

#determining buffer -10000 m
buffer = feat.geometry().buffer(-10000,-1)

#Extract CRS from layer
CRS = layer.crs().postgisSrid()

URI = "Polygon?crs=epsg:"+str(CRS)+"&field=id:integer""&index=yes"

#Create polygon layer for buffer
mem_layer = QgsVectorLayer(URI,
                           "buffer",
                           "memory")

#Prepare mem_layer for editing
mem_layer.startEditing()

#Set feature for buffer
feat2 = QgsFeature()

#Set geometry for buffer
feat2.setGeometry(buffer)

#set attributes values for buffer
feat2.setAttributes([1])

mem_layer.addFeature(feat2, True)

#stop editing and save changes
mem_layer.commitChanges()

#preparing layer for difference
feat3 = mem_layer.getFeatures().next()

difference = feat.geometry().difference(feat3.geometry())

#Create polygon layer for difference
mem_layer2 = QgsVectorLayer(URI,
                           "diference",
                           "memory")

#Prepare mem_layer2 for editing
mem_layer2.startEditing()

#Set feature for difference
feat2 = QgsFeature()

#Set geometry for difference
feat2.setGeometry(difference)

#set attributes values for difference
feat2.setAttributes([1])

mem_layer2.addFeature(feat2, True)

#stop editing and save changes
mem_layer2.commitChanges()

#add Map Layer to Registry
QgsMapLayerRegistry.instance().addMapLayer(mem_layer2)

Next image was obtained after running the code at Python Console of QGIS:

It works as expected.

Editing Note:

Next version of the script works for all features of a polygon vector layer.

layer = iface.activeLayer()

feats = [ feat for feat in layer.getFeatures() ]

#determine buffer -10000 m
buffer = [ feat.geometry().buffer(-10000,-1) for feat in feats ] 

#Extract CRS from layer
CRS = layer.crs().postgisSrid()

URI = "Polygon?crs=epsg:"+str(CRS)+"&field=id:integer""&index=yes"

#Create polygon layer for buffer
mem_layer = QgsVectorLayer(URI,
                           "buffer",
                           "memory")

#Prepare mem_layer for editing
mem_layer.startEditing()

n = len(feats)

#Set feature for buffer
feats2 = [ QgsFeature() for i in range(n) ]

for i in range(n): 
    #Set geometry for buffer
    feats2[i].setGeometry(buffer[i])

    #set attributes values for buffer
    feats2[i].setAttributes([i])

    mem_layer.addFeature(feats2[i], True)

#stop editing and save changes
mem_layer.commitChanges()

feats3 = [ feat3 for feat3 in mem_layer.getFeatures() ]

difference = []

for i in range(n):
    difference.append(feats[i].geometry().difference(feats3[i].geometry()))

#Create polygon layer for difference
mem_layer2 = QgsVectorLayer(URI,
                           "difference",
                           "memory")

#Prepare mem_layer for editing
mem_layer2.startEditing()

#Set feature for difference
feats2 = [ QgsFeature() for i in range(n) ] 

for i in range(n):
    #Set geometry for difference
    feats2[i].setGeometry(difference[i])

    #set attributes values for difference
    feats2[i].setAttributes([i])

    mem_layer2.addFeature(feats2[i], True)

#stop editing and save changes
mem_layer2.commitChanges()

#add Map Layer to Registry
QgsMapLayerRegistry.instance().addMapLayer(mem_layer2)

It was run with a vector layer where, in some features, the distance from the borders is lesser than 10 km; as it can be observed at the next image: