extracting time series data from a raster brick of AVHRR satellite data

row=40
col=58

test = mybrick[i,j]

s=ts(test, start=2001, end=c(2003, 36), frequency=36)

row=40 
col=58 

test = mybrick[row,col] 

## read in Brick
modis <- brick("modis.grd")
nlayers(modis)

dim(modis)

(system.time(test <- extract(modis,1)))

(system.time(test <- modis[1,1]))

Any advice? Should I do this via Python or can this be done via R but more
efficient? 

library(raster)

b <- brick(nc=10, nr=10, nl=648)
b <- setValues(b, matrix(rep(1:100, 648), nc=648))
system.time(v <- b[1])

d <- writeRaster(b, 'test.grd', overwrite=TRUE)
inMemory(d)

system.time(v <- d[1])  # too slow

d <- readAll(d)
inMemory(d)

system.time(v <- d[1])  # OK

e <- writeRaster(b, 'test.nc', overwrite=TRUE)

inMemory(e)

system.time(v <- e[1]) # OK

b

nc <- writeRaster(b,'modis.nc',overwrite=TRUE) # ,datatype='INT2U'

inMemory(modis)

system.time(m <- getValues(modis))

Jan,

Any advice? Should I do this via Python or can this be done via R but

more efficient?

Interesting problem, that I had not noticed before. I think this is because
a loop that can perhaps be omitted or else perhaps needs to move to a C
routine. I need to investigate that. For now, for speed, either read the
values into memory (if possible; perhaps a spatial subset, via 'crop'), or
perhaps use ncdf files (see example below).

library(raster)

Loading required package: sp
raster version 1.9-3 (27-July-2011)

b <- brick(nc=10, nr=10, nl=648)
b <- setValues(b, matrix(rep(1:100, 648), nc=648))
system.time(v <- b[1])

   user  system elapsed
      0       0       0

d <- writeRaster(b, 'test.grd', overwrite=TRUE)
inMemory(d)

[1] FALSE

system.time(v <- d[1])  # too slow

   user  system elapsed
  61.54    0.14   61.98

d <- readAll(d)
inMemory(d)

[1] TRUE

system.time(v <- d[1])  # OK

   user  system elapsed
      0       0       0

e <- writeRaster(b, 'test.nc', overwrite=TRUE)

Loading required package: ncdf

inMemory(e)

[1] FALSE

system.time(v <- e[1]) # OK

   user  system elapsed
   0.01    0.00    0.02

Robert


------------------------------
 If you reply to this email, your message will be added to the discussion
below:

http://r-sig-geo.2731867.n2.nabble.com/extracting-time-series-data-from-a-raster-brick-of-AVHRR-satellite-data-tp6622055p6627415.html
 To unsubscribe from extracting time series data from a raster brick of
AVHRR satellite data, click here<http://r-sig-geo.2731867.n2.nabble.com/template/NamlServlet.jtp?macro=unsubscribe_by_code&node=6622055&code=amFudmVyYmVzc2VsdEBnbWFpbC5jb218NjYyMjA1NXwtMTAxMTIxMDUxNQ==>.

Robert,

I have tried the different options - and based on your example below with
the simulated raster brick everything works fine.
However when applying it to my dataset a few things/errors occur.

This is my dataset:

b

class       : RasterStack
dimensions  : 17, 13, 262  (nrow, ncol, nlayers)
resolution  : 231.6564, 231.6564  (x, y)
extent      : 305554.7, 308566.3, 5713341, 5717279  (xmin, xmax, ymin,
ymax)
coord. ref. : +proj=sinu +lon_0=0 +x_0=0 +y_0=0 +a=6371007.181
+b=6371007.181 +units=m +no_defs
min values  : 2874  535 2526 3779 3626 4520 5747 4820  502 4962 ...
max values  : 7671 9126 7660 7181 8723 8217 8513 8581 8979 8854 ...

1/ when trying:

nc <- writeRaster(b,'modis.nc',overwrite=TRUE) # ,datatype='INT2U'

Error in `dataType<-`(`*tmp*`, value = "FLT4S") :
  Cannot set datatype of a RasterStack

I tried setting the datatype to ,datatype='INT2U' but this gave the same
error. How could this be solved?

2/ When doing the following in a new R session:
modis <- stack("modis.grd")
modis <- readAll(modis)

inMemory(modis)

[1] FALSE

somehow when loading the file from disk in a new R session and then doing
readAll() is does not load it into memory and when trying modis[1] is
remains slow. How can I make sure that the file loads into memory?

3/ therefore I currently opted for the following:

system.time(m <- getValues(modis))

   user  system elapsed
  0.347   0.030   0.453

## create a time series
d <- as.vector(m[70,]) # extract cell 70 from the matrix m

Thanks for your help. I really like the 'raster' package so let me know how
I can help testing the package.

On Wed, Jul 27, 2011 at 9:23 PM, Robert Hijmans [via R-sig-geo] <[hidden
email] <http://user/SendEmail.jtp?type=node&node=6629146&i=0>> wrote:

Jan,

Any advice? Should I do this via Python or can this be done via R but

more efficient?

Interesting problem, that I had not noticed before. I think this is
because a loop that can perhaps be omitted or else perhaps needs to move to
a C routine. I need to investigate that. For now, for speed, either read the
values into memory (if possible; perhaps a spatial subset, via 'crop'), or
perhaps use ncdf files (see example below).

library(raster)

Loading required package: sp
raster version 1.9-3 (27-July-2011)

b <- brick(nc=10, nr=10, nl=648)
b <- setValues(b, matrix(rep(1:100, 648), nc=648))
system.time(v <- b[1])

   user  system elapsed
      0       0       0

d <- writeRaster(b, 'test.grd', overwrite=TRUE)
inMemory(d)

[1] FALSE

system.time(v <- d[1])  # too slow

   user  system elapsed
  61.54    0.14   61.98

d <- readAll(d)
inMemory(d)

[1] TRUE

system.time(v <- d[1])  # OK

   user  system elapsed
      0       0       0

e <- writeRaster(b, 'test.nc', overwrite=TRUE)

Loading required package: ncdf

inMemory(e)

[1] FALSE

system.time(v <- e[1]) # OK

   user  system elapsed
   0.01    0.00    0.02

Robert


------------------------------
 If you reply to this email, your message will be added to the discussion
below:

http://r-sig-geo.2731867.n2.nabble.com/extracting-time-series-data-from-a-raster-brick-of-AVHRR-satellite-data-tp6622055p6627415.html
 To unsubscribe from extracting time series data from a raster brick of
AVHRR satellite data, click here.

------------------------------
 If you reply to this email, your message will be added to the discussion
below:

http://r-sig-geo.2731867.n2.nabble.com/extracting-time-series-data-from-a-raster-brick-of-AVHRR-satellite-data-tp6622055p6629146.html
 To unsubscribe from extracting time series data from a raster brick of
AVHRR satellite data, click here<http://r-sig-geo.2731867.n2.nabble.com/template/NamlServlet.jtp?macro=unsubscribe_by_code&node=6622055&code=ci5oaWptYW5zQGdtYWlsLmNvbXw2NjIyMDU1fC0xMDIyNzY3NTQ3>.