13

I would like to call the QGIS interpolation plugin function (TIN method) (Raster->Interpolate) from python console.

I can't find the corresponding function within QGIS API or within the processing algorithm list. I found the SAGA Triangulation algorithm, which works fine but is 5-10 x slower and speed is important in my case.

Any idea how to execute it?

  • 2
    Although I do not require this, it would be a useful thing to know. I followed this link: (gis.stackexchange.com/questions/11216/…). I got as far as from rasterinterpolation import rasterinterpolation but not sure which module to call (or how to even call). – Joseph Oct 9 '14 at 13:32
  • Could you clarify your requirements a little more? Are you simply looking for a way to compute a new interpolated raster layer from an input raster layer? – underdark Oct 11 '14 at 7:58
  • I have a similar problem: I want to create a countour model that starts with Raster\interpolation followed by Saga\contours from grid. Question is - how to add the rasterinrepolator in the "processing modeler" window? – H.Wiener Jan 9 '15 at 10:25
5

I was able to provide a full solution in the following question:

How to compute an interpolation raster from the python console in QGIS?

I will repost the answer here as well, because of the large interest it seems to attract:

Answer:

The documentation on pyqgis is not very self-explanatory, but i figured out how to properly call the associated interpolation classes (QgsInterpolator, QgsTINInterpolator, QgsIDWInterpolator, QgsGridFileWriter) from python. I am going to describe every step of the script in great detail:

Step 1:

Import the core and analysis module and get the desired vector layer for interpolation by selecting it with a mouseclick in the layer tab.

import qgis.core
import qgis.analysis

layer = qgis.utils.iface.activeLayer()

Step 2:

Prepare the interpolation classes with the necessary Parameters. The exact parameters for initialization of the LayerData struct can be found in the QGIS API docs (searchterm: QgsInterpolator).

layer_data = QgsInterpolator.LayerData()
layer_data.vectorLayer = layer
layer_data.zCoordInterpolation=False
layer_data.InterpolationAttribute =0
layer_data.mInputType = 1

Please notice that I don't use the z Coordinate, I get the first available field (index = 0) as interpolation attribute, and use POINTS as input type.

Step 3:

Choose your interpolation engine. Here you can choose between the TIN-Interpolation method (QgsTINInterpolator) and IDW-Interpolation (QgsIDWInterpolator). I took the QgsTINInterpolator in my code.

tin_interpolator = QgsTINInterpolator([layer_data])

Keep in mind that you have to pass a python list of layer_data to the interpolation engine! This also allows you to add multiple layer_data scenarios.

Step 4:

Setup the parameters that are needed for the export of the interpolation-output (see documentation of QgsGridFileWriter). Those include similar information as the interpolation gui (filepath, extent, resolution, number of colums and rows).

export_path ="C:/SomeFolder/output.asc"
rect = layer.extent()
res = 10
ncol = int( ( rect.xMaximum() - rect.xMinimum() ) / res )
nrows = int( (rect.yMaximum() - rect.yMinimum() ) / res)

output = QgsGridFileWriter(tin_interpolator,export_path,rect,ncol, nrows,res,res)
output.writeFile(True)  

iface.addRasterLayer(export_path, "interpolation_output") 

Be aware of the file extension of your output-raster as QgsGridFileWriter only writes ASCII-grids (.asc). The data gets written to disk by calling the writeFile() method. After export you can add the grid-file as raster to the canvas.

Full script for reference:

import qgis.analysis
import qgis.core

layer = qgis.utils.iface.activeLayer() 
layer_data = QgsInterpolator.LayerData()
layer_data.vectorLayer = layer
layer_data.zCoordInterpolation=False
layer_data.InterpolationAttribute =0
layer_data.mInputType = 1


tin_interpolator = QgsTINInterpolator([layer_data])

export_path = "E:/GIS_Workbench/script_output/test.asc"

rect = layer.extent()
res = 10
ncol = int( ( rect.xMaximum() - rect.xMinimum() ) / res )
nrows = int( (rect.yMaximum() - rect.yMinimum() ) / res)
output = QgsGridFileWriter(tin_interpolator,export_path,rect,ncol,nrows,res,res)
output.writeFile(True)

Keep in mind that the QGIS-API is currently rewritten to version 3.0 and the used interpolation-classes are moved from qgis.analysis to qgis.core! This will have a huge impact on the functionality of this script so that it must be rewritten for version 3.0!

|improve this answer|||||
  • 1
    I try your example code, but only working by layer_data.InterpolationAttribute =0, I try with another field index, but and only get 0. – Leonard Apr 14 '17 at 20:29
  • That is right - I encountered this problem too, but I did not have enough time to investigate into the cause. My solution was to feed the script a layer that had just the single desired field. You might try the QGIS API docs for an improved solution. – root676 Apr 14 '17 at 21:11
3
+50

You could do this if you have the Raster Interpolation plugin installed using Plugin Manager.

from rasterinterpolation.core.rasterinterpolator import RasterInterpolator
rastLayer = iface.activeLayer()
interpolator = RasterInterpolator(rastLayer,0,1)
a= interpolator.linear(QgsPoint(10.662629, 76.225421))
print a

Note: i don't really know what the above code does other than the fact that it printed a value. But it would probably help you to understand the usage.

|improve this answer|||||

Not the answer you're looking for? Browse other questions tagged or ask your own question.