[R-sig-Geo] How to apply this function on each raster stack layer?

[Thiago V. dos Santos]

Frede and Loic,
Thanks a lot for your feedback. I just realized that a simple multiplication 
s1 <- (s * mult) + add
not only works but is much faster than the function I wrote. Raster is a really well-implemented package. Greetings,
--
Thiago V. dos Santos
PhD student
Land and Atmospheric Science
University of Minnesota
http://www.laas.umn.edu/CurrentStudents/MeettheStudents/ThiagodosSantos/index.htm
Phone: (612) 323 9898 

     On Friday, February 6, 2015 2:07 AM, "Dutrieux, Loic" <loic.dutrieux at wur.nl> wrote:
   

 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


   
	[[alternative HTML version deleted]]