How to apply this function on each raster stack layer?

Fri, Feb 6, 2015 12:06 AM

Hi Thiago,

The approach you took seems correct. See example below.

library(raster)
fn <- system.file("external/test.grd", package="raster")
s  <- stack(fn, fn)

# And these are sample coefficients:
mult <- c(0.0003342, 0.0005534) 
add  <- c(0.1, 0.2) 

# define calc function
fun <- function(x) {
    y <- x * mult + add
    return(y)
}

out <- calc(x = s, fun = fun)

# For each band individually
s2 <- s
s2[[1]] <- (s2[[1]] * mult[1]) + add[1]
s2[[2]] <- (s2[[2]] * mult[2]) + add[2]

all.equal(out, s2)


Best regards,
---
Lo?c Dutrieux
PhD student
Laboratory of Geo Information Science and Remote Sensing
Wageningen University
The Netherlands

From: R-sig-Geo <r-sig-geo-bounces at r-project.org> on behalf of Frede Aakmann T?gersen <frtog at vestas.com>
Sent: Friday, February 6, 2015 7:25 AM
To: Thiago V. dos Santos; R-SIG list
Subject: Re: [R-sig-Geo] How to apply this function on each raster stack layer?

Hi Thiago

I think that the forceapply argument of calc() is what you need; see ?calc. Here is the relevant part from the man:

     The intent of some functions can be ambiguous. Consider:

     'library(raster)'

     'r <- raster(volcano)'

     'calc(r, function(x) x * 1:10)'

     In this case, the cell values are multiplied in a vectorized
     manner and a single layer is returned where the first cell has
     been multiplied with one, the second cell with two, the 11th cell
     with one again, and so on. But perhaps the intent was to create 10
     new layers ('x*1, x*2, ...')? This can be achieved by using
     argument 'forceapply=TRUE'

     'calc(r, function(x) x * 1:10), forceapply=TRUE'


Yours sincerely / Med venlig hilsen


Frede Aakmann T?gersen
Specialist, M.Sc., Ph.D.
Plant Performance & Modeling

Technology & Service Solutions
T +45 9730 5135
M +45 2547 6050
frtog at vestas.com
http://www.vestas.com

Company reg. name: Vestas Wind Systems A/S
This e-mail is subject to our e-mail disclaimer statement.
Please refer to www.vestas.com/legal/notice
If you have received this e-mail in error please contact the sender.


> -----Original Message-----
> From: R-sig-Geo [mailto:r-sig-geo-bounces at r-project.org] On Behalf Of
> Thiago V. dos Santos
> Sent: 5. februar 2015 19:22
> To: R-SIG list
> Subject: [R-sig-Geo] How to apply this function on each raster stack layer?
>
> Hi all,
>
> I am trying to write an efficient script to calibrate hundreds of Landsat 8
> images. At a certain point of the calibration steps, I need to apply some
> coefficients in each layer of a raster stack.
>
>
> This is one sample stack:fn <- system.file("external/test.grd",
> package="raster")
>
> s  <- stack(fn, fn)
>
> And these are sample coefficients:
> mult <- c(0.0003342, 0.0005534)
> add  <- c(0.1, 0.2)
>
>
> What I need to is to apply each index of the coefficients to the
> correspondent index of the stack layer, like in this example:
>
> s[[1]] <- (s[[1]] * mult[1]) + add[1]
> s[[2]] <- (s[[2]] * mult[2]) + add[2]
>
> This is my poor attempt, which obviously does not work:
>
> cal.fun <- function(x) {
> x <- (x * mult) + add
> }
>
> s.cal <- calc(s, cal.fun, progress='text')
>
>
> Any ideas on how to do that?
>
> Many thanks,
> --
> Thiago V. dos Santos
> PhD student
> Land and Atmospheric Science
> University of Minnesota
> http://www.laas.umn.edu/CurrentStudents/MeettheStudents/ThiagodosSa
> ntos/index.htm
> Phone: (612) 323 9898
>
> _______________________________________________
> R-sig-Geo mailing list
> R-sig-Geo at r-project.org
> https://stat.ethz.ch/mailman/listinfo/r-sig-geo

_______________________________________________
R-sig-Geo mailing list
R-sig-Geo at r-project.org
https://stat.ethz.ch/mailman/listinfo/r-sig-geo

How to apply this function on each raster stack layer?

Thread (4 messages)