3

If have three rasters(as matrix):

r1 <- raster(nrows=10, ncols=10); r1 <- setValues(r1, 1:ncell(r1))
r16 <- raster(nrows=10, ncols=10);r16 <- setValues(r16, 1:ncell(r16))
r30 <- raster(nrows=10, ncols=10);r30 <- setValues(r30, 1:ncell(r30))

I would like to linearly interpolate r1,c16,c30 to find the values in between i.e. r2,r3,r4,......r15 then r17,r18,r19,..........r29.

Is this possible using R?

temor
  • 935
  • 2
  • 10
  • 26

2 Answers2

3

Here is a way to do that

library(raster)
r <- raster(nrows=10, ncols=10); 
values(r) <- NA

x <- sapply(1:30, function(...) r)
x[[1]] <- setValues(r, runif(ncell(r)))
x[[16]] <- setValues(r, runif(ncell(r))) + 10
x[[30]] <- setValues(r, runif(ncell(r))) + 20

s <- stack(x)

z <- approxNA(s)

plot(z)
plot(1:30, z[1])

Here is another way to do it

library(raster)
r <- raster(nrows=10, ncols=10); 
x1 <- setValues(r, runif(ncell(r)))
x16 <- setValues(r, runif(ncell(r))) + 10
x30 <- setValues(r, runif(ncell(r))) + 20

s <- stack(x1, x16, x30)
x <- calc(s, fun=function(y) approx(c(1,16,30), y, 1:30)$y)

But this will fail if there are NA values in the three layers. You would need to adjust the function fun to deal with that (here is an example).

Community
  • 1
  • 1
Robert Hijmans
  • 40,301
  • 4
  • 55
  • 63
  • 1
    Well that seems to make [this monstrosity](https://github.com/johnbaums/things/blob/aae7fe671211dc3a4ca8fe6f3da71c28bad448be/R/interpolate_linear.R) more or less redundant. – jbaums Oct 21 '15 at 19:44
  • 1
    @jbaums the main issue with your function is it does not extrapolate backwards. If for example I have a data for one year that starts on the 20th jan and I want to inter/extra for the whole year,your function will not work. – temor Oct 22 '15 at 14:53
  • @temor - that's true (and I mention it in the details section of the docs). – jbaums Oct 22 '15 at 15:00
  • can you please take a look at this: http://stackoverflow.com/questions/33284379/convert-spline-interpolation-to-linear-interpolation – temor Oct 22 '15 at 15:08
  • I can look at it, but why? – Robert Hijmans Oct 22 '15 at 16:30
0

Is there a way to do this but change xout in approx to the cell values of another raster? I've been trying to get two different rasters into app and it won't seem to work for me.

I'm looking to interpolate values to an depth, where each layer in the rast would be the climate value at a fixed depth. Going from surface to the seafloor. The output would be a single layer raster where the stack cell level values are interpolated to the nearest depth based on depth of input raster z.

I was thinking something like this:

library(raster)
r <- raster(nrows=10, ncols=10); 
x1 <- setValues(r, runif(ncell(r)))
x16 <- setValues(r, runif(ncell(r))) + 10
x30 <- setValues(r, runif(ncell(r))) + 20

## this would be say elevations per cell
z <- setValues(r,sample(1:30,100,replace=TRUE))

s <- stack(x1, x16, x30)
x <- calc(s, fun=function(y) approx(x=c(1,16,30), y=y, xout=z)$y)

I've tried various ways to index the depth raster z with no luck. Any pointers would be great.

Skiptoniam
  • 91
  • 1
  • 7