2

The estimation of the Topographic Position Index (TPI) is scale dependent, but when estimating this index in QGIS (Raster > Analyses > DEM (Terrain models)), I don't see where I can change this scale, nor its standard value.

Do you know how do I obtain/change that(if possible)? Is it related with the pixel value used in the input DEM raster?

  • 1
    I agree that scale (or radius) is important. Just tried to read QGIS documentation and gdaldem tpi. I felt the same as you probably did, that gdal was simply looking at its neighbor cells. You might have to change the size of DEM pixels... probably not a preferred option. Instead, SAGA will give you more control. Processing | SAGA | Terrain Analysis | Morphometry | Topographic position index (tpi). – Kazuhito Dec 1 '16 at 16:18
  • @Kazuhito, yes that's what I thought it was considering, unfortunately... I didn't try with SAGA though. And probably GRASS will have more hypotheses as well! Thanks for your suggestion! – mtao Dec 1 '16 at 22:20
  • 1
    I would just use raster algebra to obtain the desired result. The TIP is just the focal cell subtracted from the focal mean window. You can calculate a focal mean, using any window desired, then subtract the DEM from the resulting focal raster and you have the TPI. – Jeffrey Evans Dec 2 '16 at 17:03
  • 1
    Raster algebra solutions are detailed at gis.stackexchange.com/questions/6056. – whuber Dec 3 '16 at 18:38
6

I do not know much about GRASS, but given there are variety of tools such as r.neighbors (focal statistics), I'm sure you already have found the solution with GRASS.

BTW, this is just to follow up on the other lead - comments by Jeffrey and whuber, both recommended raster algebra.

I tried to put it in R script, which can be run by copy/paste with Script Editor - via Processing Toolbox | R scripts | Tools | Create new R script.

TPI.rsx -- please copy this to Script Editor

##Raster_Tools= group
##Input_Raster_Layer= raster
##Band_number= number 1
##Outer_Window_Size= number 3
##Inner_Window_Size= optional number 1
##Compute_edges= boolean TRUE
##Topographic_Position_Index= output raster

library(raster)
ow <- Outer_Window_Size
iw  <- Inner_Window_Size
edge <- Compute_edges
Layer <- Input_Raster_Layer[[Band_number]]
Outer_window <- focal(Layer, w=matrix(1, nrow=ow, ncol=ow), 
                      fun=mean, na.rm=FALSE, NAonly=TRUE, pad=edge)
if(iw>=3){Inner_window = focal(Layer, w=matrix(1, nrow=iw, ncol=iw), 
                      fun=mean, na.rm=FALSE, NAonly=TRUE, pad=edge)}
if(iw==1){Inner_window = Layer}
Topographic_Position_Index <- Inner_window - Outer_window

If you save this file into ~/.qgis2/processing/rscripts/ as TPI.rsxthere would appear a new geoalgorithm TPI under Raster Tools group.

enter image description here

Note

  1. Outer window size: default= 3 (cells)...odd number 3, 5, 7, ...
  2. Inner window size: default= 1 (cell)...odd number 1, 3, 5, ...

If we run it with Outer=3 and Inner= 1, it produces a raster layer which is equivalent to the gdaldem tpi. Please modify Outer window size as necessary.

Usually Inner window size would always be 1, unless you need to work with very fine grid (e.g. LiDAR) which can make your TPI analysis difficult.

| improve this answer | |
  • Very helpful, thanks Kazuhito. Still there is a problem and I'm not sure where: Following your instructions I get an error when running the tool: "invalid syntax (, line 17) See log for more details" (well actually this is the log...). I would appreciate your further help. – esther Jul 24 '17 at 14:33
  • @esther Apologies... I cannot test it myself right now. It can be from various causes; most likely my faulty script, or recent change in raster package, or configuration in the QGIS processing. Will you post it as new question? If the focus is on the error message, I guess it will not be marked as dupulicate. – Kazuhito Jul 25 '17 at 8:08

Your Answer

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

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