[R-sig-Geo] extracting time series data from a raster brick of AVHRR satellite data

Robert Hijmans r.hijmans at gmail.com
Thu Jul 28 11:53:05 CEST 2011


On Thu, Jul 28, 2011 at 1:53 AM, janverbesselt [via R-sig-geo] <
ml-node+6629146-1537221348-149542 at n2.nabble.com> wrote:

> 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?
>


Ah, bug, happens when values of b are in memory, not when they are on disk.
Fixed in devel version. Thanks, Robert



>
> 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?
>
>
Oh, that not right either (readAll with a RasterStack does not set the
inMemory flag correctly).

modis <- brick("modis.grd")
modis <- readAll(modis)



> 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.
>
>
Precisely by sending these type of reports.
Thanks, Robert.



>
> 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>.
>
>


--
View this message in context: http://r-sig-geo.2731867.n2.nabble.com/extracting-time-series-data-from-a-raster-brick-of-AVHRR-satellite-data-tp6622055p6629299.html
Sent from the R-sig-geo mailing list archive at Nabble.com.



More information about the R-sig-Geo mailing list