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
4

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!

  • 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.

protected by Community Apr 14 '17 at 21:04

Thank you for your interest in this question. Because it has attracted low-quality or spam answers that had to be removed, posting an answer now requires 10 reputation on this site (the association bonus does not count).

Would you like to answer one of these unanswered questions instead?

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