使用R的netcdf栅格堆栈或栅格砖的时间和地理子集 [英] time and geographical subset of netcdf raster stack or raster brick using R

问题描述

对于以下具有2016年每日全球海表温度的netcdf文件，我正在尝试(i)暂时划分子集，(ii)地理子集，(iii)然后对每个像素采取长期手段并创建一个基本情节.

For the following netcdf file with daily global sea surface temperatures for 2016, I'm trying to (i) subset temporally, (ii) subset geographically, (iii) then take long-term means for each pixel and create a basic plot.

链接到文件:此处

library(raster)
library(ncdf4)

设置工作目录后打开netcdf

open the netcdf after setting my working directory

nc_data <- nc_open('sst.day.mean.2016.v2.nc')

更改时间变量，以便于解释

change the time variable so it's easy to interpret

time <- ncdf4::ncvar_get(nc_data, varid="time")
head(time)

更改为我可以解释的日期

change to dates that I can interpret

time_d <- as.Date(time, format="%j", origin=as.Date("1800-01-01"))

现在我只想将9月1日至10月15日作为子集，但无法弄清楚……

Now I'd like to subset only September 1 to October 15, but can't figure that out...

跟随时间子集，创建栅格图块(或堆栈)和地理子集

Following temporal subset, create raster brick (or stack) and geographical subset

b <- brick('sst.day.mean.2016.v2.nc') # I would change this name to my file with time subest

地理上的子集

b <- crop(b, extent(144, 146, 14, 16))

最后，我想获取我整天数据的每个像素的平均值，将其分配给单个栅格，然后进行简单绘图...

Finally, I'd like to take the average for each pixel across all my days of data, assign this to a single raster, and make a simple plot...

感谢您的帮助和指导.

推荐答案

在b <- brick('sst.day.mean.2016.v2.nc')之后，我们可以键入b来查看光栅图块的信息.

After b <- brick('sst.day.mean.2016.v2.nc'), we can type b to see information of the raster brick.

b
# class       : RasterBrick 
# dimensions  : 720, 1440, 1036800, 366  (nrow, ncol, ncell, nlayers)
# resolution  : 0.25, 0.25  (x, y)
# extent      : 0, 360, -90, 90  (xmin, xmax, ymin, ymax)
# coord. ref. : +proj=longlat +datum=WGS84 +ellps=WGS84 +towgs84=0,0,0 
# data source : C:\Users\basaw\Downloads\sst.day.mean.2016.v2.nc 
# names       : X2016.01.01, X2016.01.02, X2016.01.03, X2016.01.04, X2016.01.05, X2016.01.06, X2016.01.07, X2016.01.08, X2016.01.09, X2016.01.10, X2016.01.11, X2016.01.12, X2016.01.13, X2016.01.14, X2016.01.15, ... 
# Date        : 2016-01-01, 2016-12-31 (min, max)
# varname     : sst 

请注意，Date插槽具有从2016-01-01到2016-12-31的信息，这意味着Z值已经具有日期信息，我们可以使用它来对栅格块进行子集化.

Notice that the Date slot has information from 2016-01-01 to 2016-12-31, which means the Z values already has date information and we can use that to subset the raster brick.

我们可以使用getZ函数访问存储在Z值中的值.输入getZ(b)我们可以看到一系列日期.

We can use the getZ function to access the values stored in the Z values. Type getZ(b) we can see a series of dates.

head(getZ(b))
# [1] "2016-01-01" "2016-01-02" "2016-01-03" "2016-01-04" "2016-01-05" "2016-01-06"

class(getZ(b))
# [1] "Date"

因此，我们可以使用以下代码对栅格砖进行子集化.

We can thus use the following code to subset the raster brick.

b2 <- b[[which(getZ(b) >= as.Date("2016-09-01") & getZ(b) <= as.Date("2016-10-15"))]]

然后我们可以根据您提供的代码裁剪图像.

We can then crop the image based on the code you provided.

b3 <- crop(b2, extent(144, 146, 14, 16))

要计算平均值，只需使用mean函数.

To calculate the average, just use the mean function.

b4 <- mean(b3, na.rm = TRUE)

最后，我们可以绘制平均值.

Finally, we can plot the average.

plot(b4)

这篇关于使用R的netcdf栅格堆栈或栅格砖的时间和地理子集的文章就介绍到这了，希望我们推荐的答案对大家有所帮助，也希望大家多多支持IT屋！