Monday, 15 February 2016

QGIS python script loop for raster calculator


How can I write a QGIS script for making a loop of raster calculator through a set of rasters? I know I can read multiple raster by applying:


lddLrs = qgis.utils.iface.legendInterface().layers()
for lyr in lddLrs:

And I want to apply an easy operation such as “Raster/Raster*Raster” and repeat it through multiple rasters and get the result with the same name of the input raster with a suffix.


Any Idea?


For using the script to make a new tool within the "processing Toolbox - Scripts - Create New Script" I made the following script but it doesn't work. Any idea to fix it?


##NoDATA Background=name  

##lyr=multiple raster
##OUT=output raster

import qgis
from qgis.analysis import QgsRasterCalculatorEntry, QgsRasterCalculator

# Get layer object
layer = processing.getObject(lyr)
lddLrs = qgis.utils.iface.legendInterface().layers()
for lyr in lddLrs:

entries = []
ras = QgsRasterCalculatorEntry()
ras.ref = 'lyr@1'
ras.raster = lyr
ras.bandNumber = 1
entries.append( ras )
calc = QgsRasterCalculator( '("lyr@1" / "lyr@1") * "lyr@1"', +lyr.name() + "_suffix.tif", 'GTiff', lyr.extent(), lyr.width(), lyr.height(), entries )
calc.processCalculation()

resulting tool




Answer



You could use something like the following:


import processing

lddLrs = qgis.utils.iface.legendInterface().layers()
path = "path/to/results//"
for lyr in lddLrs:
processing.runalg("saga:rastercalculator", lyr, None, "(a/a)*a", True, 7, path + lyr.name() + "_suffix.tif")




EDIT :


If you want to use the QgsRasterCalculator, you could try using the following (based on this post):


from qgis.analysis import QgsRasterCalculatorEntry, QgsRasterCalculator

lddLrs = qgis.utils.iface.legendInterface().layers()
path = "path/to/results"
for lyr in lddLrs:
entries = []
ras = QgsRasterCalculatorEntry()
ras.ref = 'ras@1'

ras.raster = lyr
ras.bandNumber = 1
entries.append( ras )
calc = QgsRasterCalculator( '(ras@1 / ras@1) * ras@1', path + lyr.name() + "_suffix.tif", 'GTiff', lyr.extent(), lyr.width(), lyr.height(), entries )
calc.processCalculation()



UPDATE:


You can use the following code in a script which allows you to select individual rasters and run the calculator on them. I've changed the output to a folder so that multiple (or single) rasters can be output to the same directory but you can change this:


##NoDATA Background=name

##Raster=multiple raster
##OUT=folder

import glob, qgis
from PyQt4.QtCore import QFileInfo
from qgis.analysis import QgsRasterCalculatorEntry, QgsRasterCalculator
from qgis.core import QgsMapLayerRegistry, QgsRasterLayer

# Split rasters
layers = Raster.split(';')

output_path = OUT + "/"
suffix = "_suffix.tif"

for ras in layers:
# Get layer object
lyr = processing.getObjectFromUri(ras)
entries = []
ras = QgsRasterCalculatorEntry()
ras.ref = 'lyr@1'
ras.raster = lyr

ras.bandNumber = 1
entries.append( ras )
calc = QgsRasterCalculator( '(lyr@1 / lyr@1) * lyr@1', output_path + lyr.name() + suffix, 'GTiff', lyr.extent(), lyr.width(), lyr.height(), entries )
calc.processCalculation()

for results in glob.glob(output_path + "*" + suffix):
fileInfo = QFileInfo(results)
path = fileInfo.filePath()
baseName = fileInfo.baseName()
layer = QgsRasterLayer(path, baseName)

QgsMapLayerRegistry.instance().addMapLayer(layer)

No comments:

Post a Comment

arcpy - Changing output name when exporting data driven pages to JPG?

Is there a way to save the output JPG, changing the output file name to the page name, instead of page number? I mean changing the script fo...