最近在做一些数据方面的研究,需要处理很多 netCDF 格式(例如 nc 格式、nc4 格式)的栅格数据。
处理这类数据,使用 Python 应该比较简单,虽然 Python 有性能的问题,但是凭借着众多的数据处理库,还是使用 Python 处理 netCDF 数据较为方便。
其中 netCDF4 库是处理 nc 和 nc4 格式数据的基本库,而 Xarray 库,则是整合了 netCDF4、NumPy 等常用库,称为了处理 netCDF 格式最方便的库了。接下来介绍下 Xarray 库的安装、使用以及踩坑的经验。
1.Xarray 库
1.1 Xarray 简介
xarray 库应该是目前用于处理 netCDF 数据最为方便的库了。xarray 可以看作是对 numpy 和 pandas 的扩展,同时也借鉴了许多 pandas 的优点,具有良好的易用性。xarray 与 netCDF4 在使用上存在许多相似之处,但是从输出的结果我们可以发现,同样是 Dataset,但是 xarray 的结果显然有更好的可读性。
1.2 netCDF 数据结构
netCDF 数据是一种多维数据,包含经度、维度、时间轴、气温、太阳辐射等各类数据。
比较好理解的方式是,将经度、维度和时间想象为自变量;而气温、太阳辐射等各类数据想象为因变量。
加入我们想要知道某一天北京市的温度,那么就可以:
f(经度, 维度, 时间) = 北京市某一天的气温
确定了相关的自变量,就等得出因变量的值。
1.3 Xarray 数据结构
xarray 有两个核心数据结构:DataArray 和 Dataset。
- DataArray(一个带有标签的多维数组)
values
获取数组的具体数值dims
获取维度的名字,如('x', 'y', 'z')
coords
获取一个类似于字典的结果,里面包含各个坐标attrs
获取原始数据的属性,比如变量的名字、单位等
- Dataset(可以简单的理解为由多个 DataArray 组成的集合)
dims
获取维度的名字,结果类似于字典,如{'x': 6, 'y': 6, 'time': 8}
data_vars
获取物理量的名字coords
获取一个类似于字典的结果,里面包含各个坐标attrs
获取原始数据的属性,比如变量的名字、单位等
Data variables 变量,对应着真实的物理数据。比如在气象学中要做出一个气压图,也就是获取”东经 xx 度、北纬 yy 度的大气压为多少帕”,显然这就是一个二维单值函数,两个维度分别是经度和纬度,函数值为大气压
1.4 Xarray 安装
根据 Xarray 的官方文档,需要以下必要的依赖:
根据我这段时间使用的经验来看,建议看装以下的依赖,这些依赖较为常用,剩下的就根据自己的数据处理情况自行选择安装:
pip install numpy
pip install pandas
pip install packaging
pip install netcdf4
pip install xarray
pip install dask
pip install scipy
1.5 读取数据
按需引入头部以及它们的常用缩写:
import numpy as np
import pandas as pd
import netCDF4 as nc
from netCDF4 import Dataset
import xarray as xr
加载 nc 文件(需要注意一些问题):
# 打开 nc 文件
ds = xr.open_dataset("filename.nc")
print(ds)
# 读取多文件数据集(此函数自动将多个文件连接并合并到一个 xarray 数据集中。推荐使用 xarray 打开多个文件的方法)
# 需要安装 dask 库
xr.open_mfdataset('myfiles/*.nc', parallel=True)
💡 中文路径和文件名不能出现中文,否则 xarray 库会报错。(如果要支持中文的话,就要转而使用 netCDF4 库)。 填写路径的时候, \ 为了防止转义,要使用双 \\ 才可以。
其中,读取数据使用的方法有两种,也就是上面将的两种数据结构方法。这两种 Xarray 数据结构的区别是:
- xarray.open_dataset() 读取 Dataset 类型数据,即能读取多个物理量。
- xarray.open_dataarray() 读取 DataArray 类型数据,即只能读取单个物理量。
💡 如果 nc 文件中含有多个物理量,用 open_dataarray() 读取会报错,因此建议统一都用 open_dataset() 来读取文件。
从文件中读取数据,假如数据中含有一个名为 longitude 的物理量:
ds = xarray.open_dataset()
ds.longitude 或 ds[longitude]
1.6 索引核心方法
在 xarray 的官方文档中给出了如下几种索引方式:
1.6.1 根据位置索引
位置索引是最直接也是最简单的索引方式,但是位置索引只对 DataArray
有效,对 DataSet
无效。下面用两种不同方法获取相同的值。
1.6.2 根据维度名字索引
通过维度的名字就可以不必按照指定的维度顺序对数据进行切片了。对 DataArray 和 DataSet 都有效,且方法一致。
# 通过数字索引
ds.isel(longitude=1, time=0)
# 通过标签索引
ds.sel(longitude=0.75, time='2018-01-01')
💡 如果遇到 KeyError(key) from err 错误。可能是存在小浮点差异,请尝试通过 method=”nearest” 解决。
1.7 插值
Xarray 库还带有多种插值功能,这对于我处理栅格数据并进行插值,就非常的方便。xarray 中对 scipy 的插值函数进行了进一步的封装,可以让我们方便的调用。只需要对 DataArray
或 DataSet
使用 interp()
函数就可以实现插值了,就像索引一样简单。不管是一维数据还是多维数据都可以轻松搞定。
In [1]: da = xr.DataArray(np.sin(0.3 * np.arange(12).reshape(4, 3)),
...: [('time', np.arange(4)),
...: ('space', [0.1, 0.2, 0.3])])
...:
# label lookup
In [2]: da.sel(time=3)
Out[2]:
<xarray.DataArray (space: 3)>
array([ 0.42738 , 0.14112 , -0.157746])
Coordinates:
time int64 3
* space (space) float64 0.1 0.2 0.3
# interpolation
In [3]: da.interp(time=2.5)
Out[3]:
<xarray.DataArray (space: 3)>
array([0.700614, 0.502165, 0.258859])
Coordinates:
* space (space) float64 0.1 0.2 0.3
time float64 2.5
da.sel(time=3)
,索引时间维的值为 3(12行)。da.interp(time=2.5)
,将时间维原本不存在的 2.5 插值了出来(22行)。
数据合并也支持维度的拼接和变量的合并,具体参阅 Xarray 库的说明即可。
2.Xarray 库使用例子
最后给大家一个使用用例,方便大家入门理解。
本例采用 “全球高分辨率(3小时,10公里)地表太阳辐射数据集(1983-2018)” 的部分数据作为演示,大家可以自行下载部分数据进行测试。
老哥真棒
\( ̄︶ ̄*\))