使用R的raster包,我brick从一个文件中获取,带有以下ncdump标题(我显示一个小的示例文件,实际文件要大得多):
dimensions:
lon = 2 ;
lat = 3 ;
time = UNLIMITED ; // (125000 currently)
variables:
float lon(lon) ;
lon:standard_name = "longitude" ;
lon:long_name = "longitude" ;
lon:units = "degrees_east" ;
lon:axis = "X" ;
float lat(lat) ;
lat:standard_name = "latitude" ;
lat:long_name = "latitude" ;
lat:units = "degrees_north" ;
lat:axis = "Y" ;
double time(time) ;
time:standard_name = "time" ;
time:long_name = "Time" ;
time:units = "seconds since 2001-1-1 00:00:00" ;
time:calendar = "standard" ;
time:axis = "T" ;
short por(time, lat, lon) ;
por:_FillValue = 0s ;
por:missing_value = 0s ;
并在R:
class : RasterBrick
dimensions : 3, 2, 6, 125000 (nrow, ncol, ncell, nlayers)
resolution : 0.008999825, 0.009000778 (x, y)
extent : 6.4955, 6.5135, 44.0955, 44.1225 (xmin, xmax, ymin, ymax)
coord. ref. : +proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0
data source : /home/clima-archive/afantini/chym/chym_output/test.nc
names : X0, X3600, X7200, X10800, X14400, X18000, X21600, X25200, X28800, X32400, X36000, X39600, X43200, X46800, X50400, ...
z-value : 0, 449996400 (min, max)
varname : por
但是,为了更快地访问和更高压缩,已经交换了两个文件维度,因此对于我们需要的那种使用,分块更好.所以文件将是这样的(链接下载1MB文件):
dimensions:
lon = UNLIMITED ; // (2 currently)
lat = 3 ;
time = 125000 ;
variables:
float lon(lon) ;
lon:standard_name = "longitude" ;
lon:long_name = "longitude" ;
lon:units = "degrees_east" ;
lon:axis = "X" ;
float lat(lat) ;
lat:standard_name = "latitude" ;
lat:long_name = "latitude" ;
lat:units = "degrees_north" ;
lat:axis = "Y" ;
double time(time) ;
time:standard_name = "time" ;
time:long_name = "Time" ;
time:units = "seconds since 2001-1-1 00:00:00" ;
time:calendar = "standard" ;
time:axis = "T" ;
short por(lon, lat, time) ;
por:_FillValue = 0s ;
por:missing_value = 0s ;
并在R:
class : RasterBrick
dimensions : 3, 125000, 375000, 2 (nrow, ncol, ncell, nlayers)
resolution : 3600, 0.009000778 (x, y)
extent : -1800, 449998200, 44.0955, 44.1225 (xmin, xmax, ymin, ymax)
coord. ref. : NA
data source : /home/clima-archive/afantini/chym/chym_output/test_swapped.nc
names : X6.5, X6.50899982452393
degrees_east: 6.5, 6.50899982452393
varname : por
如您所见,文件打开就好像列数为125000.我想将列数与层数交换,而不读取所有数据.我想从我应该使用的光栅手册,layer或者lvar,因为:
layer:整数.要在多层文件中使用的图层(变量),或从RasterStack/Brick或SpatialPixelsDataFrame或SpatialGridDataFrame中提取的图层.如果'layer = 0',则返回空的RasterLayer(没有关联的值)
.......
'lvar':整数> 0(默认值= 3).要选择要使用的"级别变量"(第三维变量),如果文件有4个维度(例如深度而不是时间)
但这似乎不起作用,例如layer="time",因为它什么都没改变.
我怎样才能做到这一点？
如果您不介意在打开/读取后进行重塑，我认为您可以使用库读取变量中的数据ncdf4，然后转置它。就像是：
nc <- nc_open(*your_nc_file*)
data <- ncvar_get(nc, por) # "por" is the name of your variable, right ?
data_new <- aperm(data, c(1,3,2)) # "transpose" the matrix
一个可能的问题可能是它data_new不再是一个raster*对象，但您可以很容易地从中重新创建一个对象。
哈特哈，
洛伦佐