temporal interpolation with stackApply - seg fault

Michael,

I think that is because of a difference in versions. The error should
go away with a more recent version of 'raster'. You can get one from
R-Forge:  install.packages("raster",
repos="http://R-Forge.R-project.org")

The current version on CRAN, 1.6-22, has some glitches, and I just
uploaded version 1.7-6 to CRAN. This version should be available for
install from there within 24 hrs.

Robert


On Wed, Dec 1, 2010 at 6:04 PM, Michael Loranty <mloranty at whrc.org>  
wrote:

Thanks for the quick response, Robert.

I hadn't realized I was misusing the stackApply. For some reason  
I'm getting
a similar error with calc.  I don't know how to look at the function
definition for calc to diagnose further.

Best,
Mike

r <- raster(ncol=10, nrow=10)
r[]=1:ncell(r)
s <- brick(r,flip(r,direction = "x"),flip(r,direction =
"y"),r,flip(r,direction = "y"),flip(r,direction = "x"))

set.NA <- function(x, na.rm = T){x[x==18] <- NA; return(x)}

s <- stackApply(s,fun = set.NA, indices = 1:nlayers(s))

fill.NA <- function(x,na.rm = F){

+                       fill <- approxfun(x, y = NULL, method =  
'linear',
rule = 2);
+                       new <- fill(1:length(x));
+                       rec <- which(is.na(x));
+                       x[rec] <- new[rec];
+                       return(x)
+                       }

ss < - calc(s,fun=fill.NA)

Error in approxfun(x, y = NULL, method = "linear", rule = 2) :
 need at least two non-NA values to interpolate

s[13]

[1] 13 NA 83 13 83 NA

fill.NA(s[13])

[1] 13 48 83 13 83 83

Michael Loranty Ph.D.
Postdoctoral Fellow
The Woods Hole Research Center
149 Woods Hole Rd
Falmouth, MA 02540
508.444.1560
mloranty at whrc.org








On Dec 1, 2010, at 7:34 PM, Robert J. Hijmans wrote:

Michael,

Thanks for sending an good self-contained example. If you use
stackApply with layers=1:nlayers you are requesting that the  
function
should be applied to each layer separately. That is OK for the first
case (set.NA) although 'calc' would be more direct. But it does not
for fill.NA, because that function needs the values of all layers  
in a
single computations. In this case you must use calc. The below  
works:


r <- raster(ncol=10, nrow=10)
r[] <- 1:ncell(r)

s <- stack(r,flip(r,direction = "x"),flip(r,direction =
"y"),r,flip(r,direction = "y"),flip(r,direction = "x"))
set.NA <- function(x, na.rm = T){x[x==18] <- NA; return(x)}
s <- calc(s, fun=set.NA)

fill.NA <- function(x,na.rm = F){
  fill <- approxfun(x, y = NULL, method = 'linear', rule = 2)
  rec <- which(is.na(x))
  x[rec] <- fill(rec)
  return(x)
}

ss <- calc(s,fun=fill.NA)


Best, Robert

On Wed, Dec 1, 2010 at 1:44 PM, Michael Loranty  
<mloranty at whrc.org> wrote:

Hi All,

I have a raster stack, that represents a time series of temperature
observations, and would like to interpolate temporally to fill  
gaps in
the
data. That is, to apply a gap-filling function across all layers  
in a
stack,
pixel by pixel. My hope was to accomplish this using stackApply,  
however
this doesn't seem to work. The reason seems to be associated with  
how the
data are handled within stackApply which is based upon the type of
function
that is used.

Has anyone else encountered a similar problem or come up with a  
solution?
It
seems to me that it may not be too difficult to make a a function  
to do
this
based upon few minor adjustments to the current raster function.  
A simple
example is below.

Best,
Mike

r <- raster(ncol=10, nrow=10)
r[]=1:ncell(r)
s <- brick(r,flip(r,direction = "x"),flip(r,direction =
"y"),r,flip(r,direction = "y"),flip(r,direction = "x"))

set.NA <- function(x, na.rm = T){x[x==18] <- NA; return(x)}

s <- stackApply(s,fun = set.NA, indices = 1:nlayers(s))

fill.NA <- function(x,na.rm = F){

+                       fill <- approxfun(x, y = NULL, method =  
'linear',
rule = 2);
+                       new <- fill(1:length(x));
+                       rec <- which(is.na(x));
+                       x[rec] <- new[rec];
+                       return(x)
+                       }

s < - stackApply(s,fun=fill.NA,indices = 1:nlayers(s))

Error in approxfun(x, y = NULL, method = "linear", rule = 2) :
 need at least two non-NA values to interpolate

sessionInfo()

R version 2.11.1 (2010-05-31)
x86_64-apple-darwin9.8.0

other attached packages:
[1] rgdal_0.6-27  raster_1.6-22 sp_0.9-72



Michael Loranty Ph.D.
Postdoctoral Fellow
The Woods Hole Research Center
149 Woods Hole Rd
Falmouth, MA 02540
508.444.1560
mloranty at whrc.org