r - 从 netCDF 更快地读取时间序列？

Question

我有一些大型的 netCDF 文件，其中包含 0.5 度分辨率的地球 6 小时数据。

每年有360个纬度点，720个经度点，1420个时间点。我有两个年度文件(12 GB ea)和一个包含 110 年数据(1.3 TB)的文件存储为 netCDF-4(这是 1901 年数据的示例， 1901.nc ，其 use policy 和 original, public files that I started with )。

据我了解，从一个 netCDF 文件中读取应该更快，而不是循环遍历一组文件 separated by year and variable originally provided 。

我想为每个网格点提取一个时间序列，例如距特定纬度和经度的 10 或 30 年。但是，我发现这非常慢。例如，虽然我可以在 0.002 秒内从单个时间点读取 10000 个值的全局切片(维度的顺序是 lat、lon、时间):

## a time series of 10 points from one location:
library(ncdf4)
met.nc <- nc_open("1901.nc")
system.time(a <- ncvar_get(met.nc, "lwdown", start = c(100,100,1), 
                                             count = c(1,1,10)))
   user  system elapsed 
  0.001   0.000   0.090 

## close down session

## a global slice of 10k points from one time
library(ncdf4)
system.time(met.nc <- nc_open("1901.nc"))
system.time(a <- ncvar_get(met.nc, "lwdown", start = c(100,100,1), 
                                             count = c(100,100,1)))
   user  system elapsed 
  0.002   0.000   0.002 

我怀疑这些文件是为了优化空间层的读取而编写的，因为 a) 变量的顺序是 lat, lon, time, b) 这将是生成这些文件的气候模型的逻辑顺序和 c) 因为全局范围是最常见的可视化。

我试图对变量重新排序，以便时间优先:

ncpdq -a time,lon,lat 1901.nc 1901_time.nc

( ncpdq 来自 NCO (netCDF operators) software )

> library(ncdf4)

## first with the original data set:
> system.time(met.nc <- nc_open("test/1901.nc"))
   user  system elapsed 
  0.024   0.045  22.334 
> system.time(a <- ncvar_get(met.nc, "lwdown", start = c(100,100,1), count = c(1, 1, 1000))
+ )
   user  system elapsed 
  0.005   0.027  14.958 

## now with the rearranged dimensions:
> system.time(met_time.nc <- nc_open("test/1901_time.nc"))
   user  system elapsed 
  0.025   0.041  16.704 
> system.time(a <- ncvar_get(met_time.nc, "lwdown", start = c(100,100,1), count = c(1, 1, 1000)))
   user  system elapsed 
  0.001   0.019   9.660 

如何优化一个点的阅读时间序列，而不是一个时间点的大区域层？例如，如果文件的写入方式不同，例如时间、纬度、经度，是否会更快？是否有一种“简单”的方法来转换 netCDF-4 文件中的维度顺序？

更新

(@mdsumner 要求的基准测试)

library(rbenchmark)
library(ncdf4)
nc <- nc_open("1901.nc")
benchmark(timeseries = ncvar_get(nc, "lwdown", 
                                 start = c(1, 1, 50), 
                                 count = c(10, 10, 100)), 
          spacechunk = ncvar_get(nc, "lwdown", 
                                  start = c(1, 1, 50), 
                                  count = c(100, 100, 1)),           
          replications = 1000)

        test replications elapsed relative user.self sys.self user.child
2 spacechunk         1000   0.909    1.000     0.843    0.066          0
1 timeseries         1000   2.211    2.432     1.103    1.105          0
  sys.child
2         0
1         0

更新 2:

我已经开始在这里开发解决方案。这些点点滴滴都在 github.com/ebimodeling/model-drivers/tree/master/met/cruncep 中的一组脚本中

脚本仍然需要一些工作和组织 - 并非所有脚本都有用。但读取速度快如闪电。与上述结果并不完全可比，但在一天结束时，我可以立即从 1.3TB 文件(0.5 度分辨率，2.5 秒)中读取 100 年、每 6 小时的时间序列:

system.time(ts <- ncvar_get(met.nc, "lwdown", start = c(50, 1, 1), count = c(160000, 1, 1)))
   user  system elapsed 
  0.004   0.000   0.004 

(注意:维度顺序已更改，如下所述: How can I specify dimension order when using ncdf4::ncvar_get? )

Answer 1

我认为这个问题的答案不会是重新排序数据，因为它会将数据分块。有关分块 netCDF 文件影响的完整讨论，请参阅 Unidata 的 NetCDF 首席开发人员 Russ Rew 的以下博客文章:

Chunking Data: Why it Matters

Chunking Data: Choosing Shapes

结果是虽然采用不同的分块策略可以大大提高访问速度，但选择正确的策略并非易事。

在较小的样本数据集上，sst.wkmean.1990-present.nc ，我在使用您的基准测试命令时看到以下结果:

1) 不分块:

## test replications elapsed relative user.self sys.self user.child sys.child
## 2 spacechunk         1000   0.841    1.000     0.812    0.029          0         0
## 1 timeseries         1000   1.325    1.576     0.944    0.381          0         0

2)天真分块:

## test replications elapsed relative user.self sys.self user.child sys.child
## 2 spacechunk         1000   0.788    1.000     0.788    0.000          0         0
## 1 timeseries         1000   0.814    1.033     0.814    0.001          0         0

天真的分块只是在黑暗中拍摄；我用了 nccopy效用因此:

$ nccopy -c"lat/100,lon/100,time/100,nbnds/" sst.wkmean.1990-present.nc chunked.nc

nccopy 的 Unidata 文档可以找到实用程序 here .

我希望我可以推荐一种用于分块数据集的特定策略，但它高度依赖于数据。希望上面链接的文章能让您深入了解如何对数据进行分块以实现您正在寻找的结果!

更新

Marcos Hermida 的以下博客文章展示了不同的分块策略如何影响读取特定 netCDF 文件的时间序列时的速度。这应该仅用作可能的起点。

Netcdf-4 Chunking Performance Results on AR-4 3D Data File

关于通过 nccopy 重新分块明显挂；该问题似乎与 4MB 的默认块缓存大小有关。通过将其增加到 4GB(或更多)，您可以将大文件的复制时间从 24 小时以上减少到 11 分钟以下!

Nccopy extremly slow/hangs

Unidata JIRA Trouble Ticket System: NCF-85, Improve use of chunk cache in nccopy utility, making it practical for rechunking large files.

有一点我不确定；在第一个链接中，讨论是关于 chunk cache , 但参数传递给 nccopy, -m , 指定复制缓冲区中的字节数。 -m nccopy 的参数控制块缓存的大小。