2

I am trying to get output from r.resample in a stand-alone PYQGIS application. The code runs and does not show any errors, but it does not save an output file (i.e. 'c:/path/to/output.tif'). As far as I can tell, I have set GRASS up ok. I can run the same code through the QGIS console and get a result.

I found a possible answer here, saying that it might be because grass commands need an extent. However, I tried that (see code below) and the command still does not run.

I am using Windows 10, Python 3.7 and QGIS 3.6.

What is the reason it is not producing an output raster?

from sys import argv
import os, sys, time, csv, ogr,subprocess
from osgeo import ogr

# GRASS7 set-up ################################
gisdb = os.path.join(os.path.expanduser("~"), "grassdata")
if sys.platform.startswith('win'):
    # MS Windows
    grass7bin = r'C:\OSGeo4W64\bin\grass76.bat'
# query GRASS GIS itself for its GISBASE
startcmd = [grass7bin, '--config', 'path']
try:
    p = subprocess.Popen(startcmd, shell=False,
                         stdout=subprocess.PIPE, stderr=subprocess.PIPE)
    out, err = p.communicate()
except OSError as error:
    sys.exit("ERROR: Cannot find GRASS GIS start script"
             " {cmd}: {error}".format(cmd=startcmd[0], error=error))
if p.returncode != 0:
    sys.exit("ERROR: Issues running GRASS GIS start script"
             " {cmd}: {error}"
             .format(cmd=' '.join(startcmd), error=err))
out2 = out.decode("utf-8")

gisbase = out2.strip(os.linesep)
print("GIS base: ",gisbase)
os.environ['GISBASE'] = gisbase
# define GRASS-Python environment
grass_pydir = os.path.join(gisbase, "etc", "python")
sys.path.append(grass_pydir)
# GRASS7 set-up complete #########################

app = QgsApplication([], False)   # second argument to False disables the GUI
app.initQgis()
# start processing
import processing
from processing import Processing
Processing.initialize()

#### r.resample starts here ####################
raster_source = path/to/raster.tif
rlay = QgsRasterLayer(raster_source, "raster_name", "gdal")
if not rlay.isValid():
    print("Layer failed to load!") # layer is valid

ext = rlay.extent()
xmin = ext.xMinimum()
xmax = ext.xMaximum()
ymin = ext.yMinimum()
ymax = ext.yMaximum()
coords = "%f,%f,%f,%f" %(xmin, xmax, ymin, ymax)

# You can see what parameters are needed by the algorithm
# using: processing.algorithmHelp("grass7:r.resample")
 params = {
     'input' : rlay,
     'output' : 'c:/path/to/output.tif',
     'GRASS_REGION_PARAMETER':coords,
     'GRASS_REGION_CELLSIZE_PARAMETER':200,
     'GRASS_RASTER_FORMAT_OPT':'',
     'GRASS_RASTER_FORMAT_META':''
                    }

res = processing.run("grass7:r.resample", params)
## Process runs but no output is saved

Edit: some commands output the final file, but not all do. For example, gdal:proximity also does not provide an output.

  • Before app.initQgis() QgsApplication.setPrefixPath("/path/to/QGIS", True) should be called, where /pat/to/QGIS is os.path.join(os.getenv("OSGEO4W_ROOT"), "apps\\qgis") under windows. Don't know if that solves your problem. – Andreas Müller Jul 9 '19 at 8:34

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Browse other questions tagged or ask your own question.