[R-sig-Geo] R+SAGA formula to subtract two raster

[Alexander Brenning]

Hi Alessandro,

SAGA has a grid_calculus library. RSAGA provides a wrapper function 
rsaga.grid.calculus, and even rsaga.linear.combination as a more 
convenient version when you want to apply linear models to stacks of 
grids. See ?rsaga.grid.calculus, there's also some example code.

Here some code for your problem of subtracting grids. I am playing 
around with two small identical grids temp1.sgrd and temp2.sgrd and want 
to show that their difference is zero.


library(RSAGA)

rsaga.grid.calculus(in.grids = c("temp1","temp2"),
     out.grid = "result", formula = ~a-b)
# formula could also be  formula = "a+b"
# 'a' refers to the first grid from in.grids,
# 'b' to the second one, etc.

# Same result:
#rsaga.linear.combination( c("temp1", "temp2"), "result",
#    coef = c(1,-1))

# Let's have a look at the data in R
# (if the grid is not too large...):
rsaga.sgrd.to.esri("result", prec = 1)
res = read.ascii.grid("result", return.header = FALSE)
summary(as.vector(res)) # in my case, all ==0 because temp1=temp2


Note that there's also a wrapper function rsaga.close.gaps for the 
module that you are using in the following piece of code:

 > # filter the missing values using neighbours:
 >
 > rsaga.geoprocessor(lib="grid_tools", module=7,
 > param=list(INPUT="DEM_1.sgrd", RESULT="DEM_1_f.sgrd", THRESHOLD=0.1))


I hope this helps.
  Alex




Alessandro wrote:
> Hi All,
> 
>  
> 
> I need help to write a code in R to subtract two raster in SAGA.
> 
>  
> 
> I have "DEM_1_f.sgrd" and "DSM_1_f.sgrd" in SAGA format and I wish to
> subtract (DSM_1_f.sgrd - DEM_1_f.sgrd) to obtain a new layer 
> 
>  
> 
> A part of the code is this:
> 
>  
> 
> #A# DEM - 1) Shapes to Grid
> 
> rsaga.get.usage("grid_gridding", 3)
> 
> rsaga.geoprocessor(lib="grid_gridding", module=3,
> param=list(GRID="DEM_1.sgrd", INPUT="ground.shp", FIELD=0, LINE_TYPE=0,
> USER_CELL_SIZE=dem.pixelsize, USER_X_EXTENT_MIN=ground at bbox[1,1],
> USER_X_EXTENT_MAX=ground at bbox[1,2], USER_Y_EXTENT_MIN=ground at bbox[2,1],
> USER_Y_EXTENT_MAX=ground at bbox[2,2]))
> 
> # filter the missing values using neighbours:
> 
> rsaga.geoprocessor(lib="grid_tools", module=7,
> param=list(INPUT="DEM_1.sgrd", RESULT="DEM_1_f.sgrd", THRESHOLD=0.1))
> 
>  
> 
> #A# DSM - 2) Shapes to Grid
> 
> rsaga.get.usage("grid_gridding", 3)
> 
> rsaga.geoprocessor(lib="grid_gridding", module=3,
> param=list(GRID="DSM_1.sgrd", INPUT="vegetation.shp", FIELD=0, LINE_TYPE=0,
> USER_CELL_SIZE=dem.pixelsize, USER_X_EXTENT_MIN=vegetation at bbox[1,1],
> USER_X_EXTENT_MAX=vegetation at bbox[1,2],
> USER_Y_EXTENT_MIN=vegetation at bbox[2,1],
> USER_Y_EXTENT_MAX=vegetation at bbox[2,2]))
> 
> # filter the missing values using neighbours:
> 
> rsaga.geoprocessor(lib="grid_tools", module=7,
> param=list(INPUT="DSM_1.sgrd", RESULT="DSM_1_f.sgrd", THRESHOLD=0.1))
> 
>  
> 
> #DCM = DSM - DEM
> 
>  
> 
> Rsaga.geoprocessor(???)
> 
>  
> 
> Thanks for help 
> 
>  
> 
> Ale
> 
>  
> 
> PS: is there a link with a SAGA and R codes banks
> 
>  
> 
>  
> 
> 
> 	[[alternative HTML version deleted]]
> 
> _______________________________________________
> R-sig-Geo mailing list
> R-sig-Geo at stat.math.ethz.ch
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo
> 
> 

-- 
Alexander Brenning
brenning at uwaterloo.ca - T +1-519-888-4567 ext 35783
Department of Geography and Environmental Management
University of Waterloo
200 University Ave. W - Waterloo, ON - Canada N2L 3G1
http://www.fes.uwaterloo.ca/geography/faculty/brenning/