如何在 R 中操作 ERA5 数据

Question

这是我第一次使用空间数据，所以我真的需要你的帮助。

我尝试使用 python（格式为 netcdf.zip）从 this 每月下载 ERA5-land 数据。我选择了可以覆盖法国的分区。之后我尝试通过以下代码读取R中的数据：

library(ncdf4)
library(raster)
library(sf)
library(ggplot2)

ncfile <- nc_open(filepath)
print(ncfile)

我得到这样的结果：（例如，我列出变量 t2m，其他有相同的描述）。

33 variables (excluding dimension variables):
        
        short t2m[longitude,latitude,expver,time]   
            scale_factor: 0.000759505975281538
            add_offset: 277.460679817325
            _FillValue: -32767
            missing_value: -32767
            units: K
            long_name: 2 metre temperature
        ..... 

     4 dimensions:
        longitude  Size:151 
            units: degrees_east
            long_name: longitude
        latitude  Size:100 
            units: degrees_north
            long_name: latitude
        expver  Size:2 
            long_name: expver
        time  Size:886 
            units: hours since 1900-01-01 00:00:00.0
            long_name: time
            calendar: gregorian

    2 global attributes:
        Conventions: CF-1.6
        history: 2023-11-21 15:03:18 GMT by grib_to_netcdf-2.24.0: /opt/ecmwf/eccodes/bin/grib_to_netcdf -S param -o /cache/tmp/7a615d9b-759a-4ab0-a842-c978c79c8850-adaptor.mars.internal-1700578884.4900465-17448-6-tmp.nc /cache/tmp/7a615d9b-759a-4ab0-a842-c978c79c8850-adaptor.mars.internal-1700573425.5042555-17448-4-tmp.grib

(1) 我不知道为什么缺失值的数量等于填充值的数量？我检查了一下，结果发现如果我过滤了expver = 2，那么我得到的都是NA，事实上，在expver 1中也存在NA。 当我尝试使用此代码绘制数据时：

t2m_var <- ncvar_get(ncfile, "t2m")
t2m_v1 <- t2m_var[,,1,880]  #longitude, latitude, expver, time
t2m_v2 <- t2m_var[,,2,880]
t2m_slice <- t2m_v1
index <- is.na(t2m_slice)
t2m_slice[index] <- t2m_v2[index]
t2m_raster <- raster(t2m_slice, xmn=min(lon), xmx=max(lon), ymn=min(lat), ymx=max(lat), crs=CRS("+proj=longlat +datum=WGS84"))

(2) 我得到了一张奇怪的图片，其中有许多对应于 NA 的空白区域（它不应该是这样的？因为 ERA5 数据覆盖了所有区域？） 我在想我在处理这些数据时是否遗漏了一些东西。

(3) 你能帮我找出哪一步出了问题以及解决方法吗？ 非常感谢！

Answer 1

您有两个独立的问题，您必须单独解决它们，因为它们不相关或复合。

首先，您有 ERA5 和 ERA5T 数据的混合，后者类型是最近一段时间的“实验”数据。数据最初作为实验发布。只有经过额外的质量控制后，数据才会作为最终产品发布。该过程大约需要 3 个月，因此所有小于 3 个月的数据都是 ERA5T。当您的数据集的时间序列跨越了 3 个月的边界时，您将获得一个额外的维度

expver

，其维度值为 1（最终产品，较旧的数据）和 5（最近的实验数据）。 这里更详细地解释了这种安排。如果您想混合这些数据（这对于许多应用程序来说是完美的），您可以付出一些努力来做到这一点。具有两个

expver

值的 3D 数组是互补的：其中一个有值，另一个

NA

值。文件中没有具体信息指示每个

expver

数组中有多少个有效时间步，但您可以简单地读取数据，然后将所有

NA

值和

 值放在“时间”维度中的切片abind::abind()

剩余数据。

其次，你的图片一点也不奇怪，有趣的地方没有

NA

值。您有 ERA5land 数据，因此有所有陆地区域的数据和海洋的

NA

数据。问题出在你映射数据的方式上。您的图像显然是镜像的，布列塔尼朝上，阿尔卑斯山朝上，比利牛斯山脉朝上，反之亦然。这里的“罪魁祸首”是 R 按列主顺序组织数据（从左上角的值向下，然后按列向右）（阅读 R 手册中的第 5 章将为您提供所有详细信息），而典型的地图软件使用行优先顺序（从左上角开始，然后向下）。如果您传递

raster::raster()

一个 NetCDF 文件，它会正确解释这一点，但如果您给它一个数组，它无法知道数据是否已针对存储中的这种差异进行了纠正。您可以通过使用

expver

排列 3D 数组（因此在合并两个

aperm(., c(2, 1, 3))

切片之后）轻松纠正此问题。