GDAL读取MAIAC的HDF文件（熟肉）

一个有趣的灵魂W

MAIAC全称Multi-Angle Implementationof Atmospheric Correction algorithm （多角度大气校正算法），主要生产的是气溶胶光学厚度（AOD）产品，在MODIS数据库中的序列号是MCD19A2，目前有c6和c6.1两个版本。空间分辨率1km，时间分辨率为1d。这个算法得到的AOD更为精确，同时获得的AOD范围也更为广泛。有兴趣的同学可以到NASA官网了解更为丰富的MAIAC算法细节（本文不做扩展）。

GDAL(Geospatial Data Abstraction Library)是一个在X/MIT许可协议下的开源栅格空间数据转换库。可以读取主流的栅格格式（以下是部分GDAL支持的栅格格式）。

MAIAC产品在NASA官网中以HDF（Hierarchical Data Format）存储的格式为主。这个格式可以存储不同类型的

图像

和数码数据的文件格式，并且可以在不同类型的机器上传输，同时还有统一处理这种文件格式的函数库。大多数普通计算机都支持这种文件格式。

好了，前言部分介绍完毕，接下来需要获取一个HDF文件进行解析，你可以选择Earth或者Disc下载一个MAIAC的文件。具体可以进NASA的官网注册账号，选择数据并下载（这里也不做展开了）。

1、注册-登录

2、搜索选择数据

3、下载

于是我得到了一个HDF（其实是一堆）数据：

这里有几个注意的：

A2017365代表数据的时间，2017是年份，365代表第365天的数据。h21v03代表数据的位置，后部分不重要了（具体我也忘了，有兴趣的可以参考用户手册）。

一切预备工作齐备，接下来就是利用Python+GDAL进行读取的工作。

首先，导入库：

import gdal

紧接着，重头戏，开始读取：

hdf=r'C:pytempThesisaodrepairmaiachdfMCD19A2.A2017365.h21v03.006.2018119214031.hdf'
ds=gdal.Open(hdf)
sub=ds.GetSubDatasets()

可以发现，HDF文件中包含了13个层。其中第一层是470nm处AOD信息（顺带说一句，Python的计数基本都是从0开始）。

第二层是550nm处的AOD：

数据的位数都是16bit。

其他层都是一些辅助信息，包括不确定性层、质量控制层、角度层等（有兴趣的可以去查一查用户手册，了解更多信息）。

到了这步已经结束了吗？并没有，还需要继续深扒层里的信息和数据，重新利用GDAL打开550nmAOD层的信息，并转为数组，才算正确读取为可以直接运行和计算的部分。

pb=sub[1][0]
b=gdal.Open(pb)
barr=b.ReadAsArray()

这样：

可以看到barr中是一个3维的数组。为什么有3维，这跟它的算法有莫大关系，简而言之，算法里包含了时间维度的检测，所以成了3维的数据。

以上，就是今天要讲解的HDF数据的读取~

彩蛋~

数据打开了，稍加处理（选取每个维度的最大值）：

max函数

最后得到这么一张AOD分布图。

好啦~今天到此为止

往期

PyCharm2019亲测破解方式

中国范围2019年道路网数据分享

分享一套中国区域的矢量图层（到县级）-更新

UTM坐标的遥感图，Python终于把它扶正啦(手动狗头)

机器学习感兴趣么？无监督的遥感图像分类感兴趣吗？来嘛！

微信号：一个有趣的灵魂W

喜欢就关注我们吧