[R-sig-Geo] [raster newbie] split along dimension into RasterBrick layers

Tom Roche Tom_Roche at pobox.com
Fri Feb 22 00:30:03 CET 2013


Suppose one has a netCDF file (name='foo.nc') with a data variable=
'foo(time, lat, lon)' having size(time)=14. I want to convert that to
a RasterBrick with one layer per timestep, each like 'bar(lat, lon)'.
I will then `projectRaster(brick ...)`, and operate on the values from
each layer separately.

However I cannot find the appropriate constructor syntax using either
`brick(x='foo.nc',...)` or `raster(x='foo.nc',...)`. Can this be done?
Or must I instead do something like

# r <- rep(raster(nrow=lat.size, ncol=lon.size), time.size) # fail
# more fail:
# r <- rep(NA, time.size) # =14
# for (i in 1:time.size) {
#   r[i] <- raster(nrow=lat.size, ncol=lon.size, crs=my.crs)
#   values(r[i]) <- foo.timestep(i)
# }
# subquestion: how does one create a vector of S4 classes?

# brute forcing the main question :-(
r.1 <- raster(nrow=lat.size, ncol=lon.size, crs=my.crs)
values(r.1) <- foo.timestep(1)
...
r.14 <- raster(nrow=lat.size, ncol=lon.size, crs=my.crs)
values(r.14) <- foo.timestep(14)
brick( r.1,  r.2,  r.3,  r.4, r.5, r.6, r.7, r.8, r.9, r.10,
      r.11, r.12, r.13, r.14)

Or is there another/better way to do this?

TIA, Tom Roche <Tom_Roche at pobox.com>



More information about the R-sig-Geo mailing list