[Bfast-commits] FW: Running bfast on raster bricks

[White, Joanne]

Hello,

I am attempting to run bfast on a raster brick (a time series of AVHRR data).

I have adapted the tutorial code provided in the bfast manual for applying bfastmonitor to raster bricks (code below).

When I try to apply my bfast function to a full image, I get the following error:

        Error in .calcTest(x[1:5], fun, na.rm, forcefun, forceapply) :
                cannot use this function

The test of the single pixel works fine and returns the correct value.

I would appreciate any suggestions on how to fix this.

Many thanks,

Joanne


################################################################
require("raster")
require("bfast")
setwd("T:/ecomonitor/age/prototype/data/fpar")

#read the ENVI format file into a raster brick
data <- brick("test.img")

#read in the list of dates for the time series
date.list <- scan("dates_list.txt", what=character(0))

#assign layer names from date.list
layerNames(data) <- date.list

#helper function for the calc function
xbfast <- function(data, date.list) {
        ndvi <- ts(data, names=date.list, frequency=36, start=c(1985,1), end=c(2007,36))
        result <- bfast(ndvi, h=0.012, season="harmonic", max.iter=1, breaks=2)
    return(cbind(result$Time))
 }

#apply on a single pixel for testing
pixel <- as.vector(data[4])
ndvi <- ts(pixel, names=date.list, frequency=36, start=c(1985,1), end=c(2007,36))
output <- xbfast(ndvi, date.list)
plot(output)

#apply on the full image                                   #This is the part that throws the error
timeofbreak <- calc(data, fun=function(x) {
output <- ts(apply(x, 1, xbfast, date.list))     #Not sure if this is supposed to be "t" as in the example or "ts". Regardless, the error message is the same.
return(output)
})
plot(timeofbreak)
################################################################


Note that the test.img is very small (4 x 4 pixels):
>data
class       : RasterBrick
dimensions  : 2, 2, 4, 828  (nrow, ncol, ncell, nlayers)
resolution  : 1000, 1000  (x, y)
extent      : 787641.7, 789641.7, 5455728, 5457728  (xmin, xmax, ymin, ymax)
coord. ref. : +proj=utm +zone=16 +datum=WGS84 +units=m +no_defs +ellps=WGS84 +towgs84=0,0,0
values      : T:\ecomonitor\age\prototype\data\fpar\test.img
layer names : X1985.01.01, X1985.01.11, X1985.01.21, X1985.02.01, X1985.02.11, X1985.02.21, X1985.03.01,...




-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://lists.r-forge.r-project.org/pipermail/bfast-commits/attachments/20120319/fc8c056f/attachment.html>