If statement in Python QGIS Raster Calculator

Question

I would like to perform a raster calculation with the if condition and I wrote this code:

import qgis
import gdal 
from qgis.analysis import QgsRasterCalculator, QgsRasterCalculatorEntry
from qgis.core import QgsRasterLayer

inputpath='/Users/macbook/Desktop/TESI/PROGETTO/rasterUfinali/U_1990.tif'
inputpath2='/Users/macbook/Desktop/TESI/PROGETTO/SIGMAmax/SIGMAmax,har.tif'
inputpath3 = '/Users/macbook/Desktop/TESI/PROGETTO/sigm/sigmahar90.tif'
inputpath4 = '/Users/macbook/Desktop/TESI/PROGETTO/lambdamin/lambminhar90.tif'
outputpath = '/Users/macbook/Desktop/TESI/PROGETTO/alfamax/alphamaxhar90.tif'


inputrasterfile = QgsRasterLayer(inputpath)

entries = []
ras = QgsRasterCalculatorEntry()
ras.ref = 'ras@1'
ras.raster = inputrasterfile
ras.bandNumber = 1
entries.append( ras )

inputrasterfile2 = QgsRasterLayer(inputpath2)

ras2 = QgsRasterCalculatorEntry()
ras2.ref = 'ras@2'
ras2.raster = inputrasterfile2
ras2.bandNumber = 1
entries.append( ras2 )


inputrasterfile3 = QgsRasterLayer(inputpath3)

ras3 = QgsRasterCalculatorEntry()
ras3.ref = 'ras@3'
ras3.raster = inputrasterfile3
ras3.bandNumber = 1
entries.append( ras3 )


inputrasterfile4 = QgsRasterLayer(inputpath4)

ras4 = QgsRasterCalculatorEntry()
ras4.ref = 'ras@4'
ras4.raster = inputrasterfile4
ras4.bandNumber = 1
entries.append( ras4 )

if '(ras@2 - ras@1) > 0':
    alpha_max='1'
else: alpha_max = 'ras@3 / (1-2.71828^(ras@4))'

new_path = outputpath
print(alpha_max)

alfamax= QgsRasterCalculator(alpha_max, new_path, 'GTiff', inputrasterfile.extent(), inputrasterfile.width(), inputrasterfile.height(), entries )
alfamax.processCalculation()
iface.addRasterLayer(new_path)

The code works and print a valid output raster, but it doesn't respect the "if condition", because it prints out alpha_max = 1 even where (ras@1 - ras@1)<0.

Can someone tell me where is the mistake?

Answer 1

A mistake I see is that you are not passing your variables into the if statement. You are actually passing a string. Let me explain. Open your python interpreter and past in the following code and run it. You will see that without even defining ras@2 or ras@1 you will always result in alpha_max = 1. This is because you passed your variables in quotation marks.

if '(ras@2 - ras@1) > 0':
    alpha_max='1'
else: alpha_max = 'ras@3 / (1-2.71828^(ras@4))'

print(alpha_max)

You are essentially saying: if "string" > 0, which translates to "if my 'string' is not empty". You are essentially comparing your string to the Boolean True/False. If you're string is not empty, it will always have a value of 1, or True. If you pass an empty string, you would have a value of 0, or False. You need to pass your variables to the if statement not inside a string. Also, it is not clear what aspect of these rasters you are trying to compare. Overall size, height, width, resolution, etc. You will likely need to try something like this, but also substituting in what value you are actually trying to compare for the values of ras1,ras2,ras3,and ras4:

if (ras2 - ras1) > 0:
    alpha_max='1'
else: alpha_max = ras3 / (1-2.71828^(ras4))

print(alpha_max)

Answer 2

QGIS Raster calculator doesn't use If...else statements. It uses this logic instead. So, functional version of your code (with my paths) looks as follows.

import qgis
import gdal 
from qgis.analysis import QgsRasterCalculator, QgsRasterCalculatorEntry
from qgis.core import QgsRasterLayer

#inputpath='/Users/macbook/Desktop/TESI/PROGETTO/rasterUfinali/U_1990.tif'
inputpath='/home/zeito/Desktop/PROGETTO/rasterUfinali/U_1990.tif'

#inputpath2='/Users/macbook/Desktop/TESI/PROGETTO/SIGMAmax/SIGMAmax,har.tif'
inputpath2='/home/zeito/Desktop/PROGETTO/SIGMAmax/SIGMAmax,har.tif'

#inputpath3 = '/Users/macbook/Desktop/TESI/PROGETTO/sigm/sigmahar90.tif'
inputpath3 = '/home/zeito/Desktop/PROGETTO/sigm/sigmahar90.tif'

#inputpath4 = '/Users/macbook/Desktop/TESI/PROGETTO/lambdamin/lambminhar90.tif'
inputpath4 = '/home/zeito/Desktop/PROGETTO/lambdamin/lambminhar90.tif'

#outputpath = '/Users/macbook/Desktop/TESI/PROGETTO/alfamax/alphamaxhar90.tif'
outputpath = '/home/zeito/Desktop/PROGETTO/alfamax/alphamaxhar90.tif'

inputrasterfile = QgsRasterLayer(inputpath)

entries = []
ras1 = QgsRasterCalculatorEntry()
ras1.ref = 'ras1@1'
ras1.raster = inputrasterfile
ras1.bandNumber = 1
entries.append( ras1 )

inputrasterfile2 = QgsRasterLayer(inputpath2)

ras2 = QgsRasterCalculatorEntry()
ras2.ref = 'ras2@1'
ras2.raster = inputrasterfile2
ras2.bandNumber = 1
entries.append( ras2 )

inputrasterfile3 = QgsRasterLayer(inputpath3)

ras3 = QgsRasterCalculatorEntry()
ras3.ref = 'ras3@1'
ras3.raster = inputrasterfile3
ras3.bandNumber = 1
entries.append( ras3 )

inputrasterfile4 = QgsRasterLayer(inputpath4)

ras4 = QgsRasterCalculatorEntry()
ras4.ref = 'ras4@1'
ras4.raster = inputrasterfile4
ras4.bandNumber = 1
entries.append( ras4 )

alpha_max = '(' + ras2.ref + ' - ' + ras1.ref + ' > 0 ' + ')*1 ' + ' + (' + ras2.ref + ' - ' + ras1.ref + ' <= 0 ' + ')* ' + ras3.ref + ' / ( 1 - 2.71828^' + ras4.ref + ')'
print(alpha_max)

new_path = outputpath

alfamax= QgsRasterCalculator(alpha_max, 
                             new_path, 
                             'GTiff', 
                             inputrasterfile.extent(), 
                             inputrasterfile.width(), 
                             inputrasterfile.height(), 
                             entries )

alfamax.processCalculation()
iface.addRasterLayer(new_path)

After run above code in Python Console, I got following result. Observe how is printed alpha_max (it contains used conditions of your formula).