Get empty layer from QGIS Raster Calculator PYQGIS

Question

I have a raster file containing different type of land use. For each type of land use I want to apply a coefficient to obtain a new map. Then I am Using QGIS Raster Calculator.

I did it manually with for example this formula :

("lu_4@1"=210)*200+("lu_4@1"=12)*20.897+("lu_4@1"=30)*27.812

Where lu_4 is my land use layer.

I obtain a black layer with nan values. When I apply my precreated style, my layer is displayed correctly.

I wanted to automatized it. When my layer is loaded in the project, I figured out that it is empty. Indeed it has only nan value and when I check with the Value Tool, all the data are no data and I can't see the layer.

Here is my code :

    rast = self.pathabs+'/_Inputs/LandUse/lu_'+ str(self.num_grid )+".tif"

    fileInfo = QFileInfo(rast)
    baseName = fileInfo.baseName()
    bohLayer = QgsRasterLayer(rast, baseName)

    entries = []
    # Define band
    boh1 = QgsRasterCalculatorEntry()
    boh1.ref = lu[:-4]+"@1"
    boh1.raster = bohLayer
    boh1.bandNumber = 1
    entries.append( boh1 )

    folder=self.pathabs+'/Maps/'+self.name_grid
    if not path.exists(folder): # create a folder for the grid if not existent
        os.makedirs(folder)
    path_map = folder+'/presence_'+name_sp+'_'+self.name_grid+'.tif'
    # Process calculation with input extent and resolution
    calc = QgsRasterCalculator( f,
                                path_map, 
                                'GTiff',
                                bohLayer.extent(), 
                                bohLayer.width(), 
                                bohLayer.height(),
                                entries)

    calc.processCalculation()

So rast is my raster path, bohLayer is my raster layer, f is the formula, folder is the output path.

I checked that :

• my formula is correct

• my raster layer is correct too, I added it to the canvas

• extent,width and height correspond to the raster layer, so no problem with these neither

• the output path is correct as i get the layer, even though it is empty

So does anyone know where the problem can come from ? I still wonder if it could come from the QgsRasterCalculator itself even if I don't get a error message. I got a problem with it by trying it directly from the python console. See this post if you want more information.

I am using QGIS 2.18.9

Answer 1

Next code produces raster that you desire. You only need to change for your raster values because, for testing proposal, I used a raster with only three values (1, 2, 3).

from qgis.analysis import QgsRasterCalculator, QgsRasterCalculatorEntry

layer = iface.activeLayer()

entries = []

#("lu_4@1"=210)*200+("lu_4@1"=12)*20.897+("lu_4@1"=30)*27.812

lyr = QgsRasterCalculatorEntry()
lyr.ref = layer.name() +'@1'
lyr.raster = layer
lyr.bandNumber = 1
entries.append( lyr )

expression = '( ( ' + entries[0].ref + ' ) = 1 ) * 200 '  + ' + ( ( '  + entries[0].ref + ' ) = 2 ) * 20.897 '  + ' + ( ( ' + entries[0].ref + ' ) = 3 ) * 27.812 '

print expression

root = '/home/zeito/pyqgis_data/'

calc = QgsRasterCalculator(expression, 
                           root + 'zzz.tif', 
                           'GTiff', 
                           layer.extent(), 
                           layer.width(), 
                           layer.height(), 
                           entries )

calc.processCalculation()

After running the code at Python Console of QGIS, with help of Value Tool plugin, it was corroborated that it worked in each case; as it can be observed at next image: