前言
这两天帮一个朋友处理了些 nc 数据，本以为很简单的事情，没想到里面涉及到了很多的细节和坑，无论是“知难行易”还是“知易行难”都不能充分的说明问题，还是“知行合一”来的更靠谱些，既要知道理论又要知道如何实现，于是经过不太充分的研究后总结成此文，以记录如何使用 python 处理 nc 数据。
一、nc 数据介绍
nc 全称 netCDF（The Network Common Data Form），可以用来存储一系列的数组，就是这么简单（参考https://www.unidata.ucar.edu/software/netcdf/docs/netcdf_introduction.html）。
既然 nc 可以用来一系列的数组，所以经常被用来存储科学观测数据，最好还是长时间序列的。
试想一下一个科学家每隔一分钟采集一次实验数据并存储了下来，如果不用这种格式存储，时间长了可能就需要创建一系列的 csv 或者 txt 等，而采用 nc 一个文件就可以搞定，是不是很方便。
更方便的是如果这个科学实验与气象、水文、温度等地理信息稍微沾点边的，完全也可以用 nc 进行存储， GeoTiff 顶多能多存几个波段（此处波段可以认为是气象、水文等不同信号），而 nc 可以存储不同波段的长时间观测结果，是不是非常方便。
可以使用 gdal 查看数据信息，执行：
gdalinfo name.nc
即可得到如下信息：
每一个 SUBDATASET 表示记录的是一种格式的数据（气象、水文等等），如果要想查看此 SUBDATASET 的具体信息，可以执行：
gdalinfo NETCDF:name.nc:SUBDATASET_NAME
此处的 SUBDATASET_NAME 为上面的 T2、PSFC 等等，可以得到如下信息：
此处只有一个 Band ，每一个 Band 记录了一个时间点（或者其他区分形式）的一条记录，这个记录是一个数组。
所以看到这里，各位应该已经明白了，可以直接使用 GDAL 处理 nc 数据，比如直接使用 gdalwarp 将某个 SUBDATASET 转成 GeoTiff 等等，此处暂且不表，各位只需要查阅一下 gdalwarp 手册即可知道如何处理。
明白了以上信息基本也就清楚了如何处理此数据。
二、数据处理
python 是运用非常广泛，自然其下各种类库非常丰富，专业一点的说法就叫生态丰富。
2.1 netCDF4
此框架可以直接将 nc 读取成数组（详细信息参考https://github.com/Unidata/netcdf4-python）。读取方式如下：
dataset = netCDF4.Dataset('name.nc')# open the dataset
这样即可读出整个 nc 中的数据信息，如果需要获取某个 SUBDATASET 只需要使用 dataset[SUBDATASET_NAME] 即可，返回的是一个三维数组，表示不同时间段（或其他区分方式下）的数据信息。
我们可以对此数组做各种操作，如求平均值、方差等等，又让我想起了大学里的那一堆枯燥但又让人很有兴趣的实验课程。当然，此处如果使用 numpy 框架进行处理，会起到事半功倍的效果，如求长时间序列下的平均值：
到这里跟地信有关的同志都会看出一个问题，此框架只能对数据进行处理，而不能进行与位置有关的操作，这就导致数据无法变成直白的地图可视化效果。其实任何数据都是相通的，我们可以采用此种方式处理完后转为 GeoTiff 等，当然我们也可以直接采用 GeoTiff 的处理流程来进行处理。
2.2 rasterio
rasterio 是 Mapbox 开源的空间数据处理框架，功能非常强大，此处不细说，只表如何处理我们的 nc 数据。
当然第一种方式就是使用 netCDF4 处理完之后，使用此框架写入 GeoTiff，但是这样不太优雅，而且使用了两个框架，明显过于麻烦，我们直接使用此框架从读数据开始处理。
此处读的时候就有技巧了，要像采用 gdalinfo 读取 SUBDATASET 一样来直接读取此 SUBDATASET 数据，如下：
即给 open 函数传入 NETCDF:name.nc:SUBDATASET_NAME，采用 src.read(range(1, dim + 1)) 可以直接读出此范围内所有 Band （时间点）的信息，范围可以自己设定，注意从 0 开始，当然也可以仅读取某个 Band 的信息。
src.meta 记录了此 SUBDATASET 的元数据信息，与 gdalinfo 看到的基本相同。
这样我们就可以继续将此数据使用 numpy 等框架进行处理，处理完之后更重要的是要写入 GeoTiff 中（直白的说就是添加空间信息）。
也很简单，如下即可：
newfile 为存储路径，res_arr 为计算结果数组，注意尺寸不要发生变化（width*height），out_meta 为目标文件的元数据描述信息，可以直接将上面 src.meta 进行简单处理即可。
crs 表示目标数据空间投影信息，transform 表示目标文件 空间范围信息，可以通过经纬度信息和图像尺寸等计算得到。
dst.write_band 将数据写入对应波段，当然此处也可以写入多个波段，根据计算结果而定，同样从 1 开始。
三、总结
本文简单介绍了 nc 数据的特点及如何使用 python 处理 nc 数据。每个目标都有多条路可以达到，重要的是找到那条自己喜欢的和适合自己的路，然而话又说回来，即使走的不是想要的那条路，不是一样可以达到目标嘛！所以关键是要找到自己的目标。
作者：魏守峰
网友评论