一种气溶胶光学厚度监测方法与流程

本申请涉及遥感技术领域，尤其涉及一种气溶胶光学厚度监测方法。

背景技术：

气溶胶光学厚度是大气污染程度的重要表征，获取其准确空间分布对于大气环境监测与整治具有重要意义。目前现有的一些研究均基于卫星遥感数据开展气溶胶光学厚度反演，主要基于eos/modis、landsat/oli和gf1/wfv传感器等。

高分四号(gf4)卫星是我国首颗地球同步轨道遥感卫星，开辟了我国地球同步轨道高分辨率对地观测的新领域。其搭载的pms和irs面阵凝视相机具有成像幅宽大、重访周期短、空间分辨率高等优点，具备对大范围气溶胶光学厚度进行高频次监测的潜力。但gf4/pms遥感数据与其它陆地遥感卫星数据在轨道高度、波段设置、成像方式上存在明显差别，不能直接将针对eos/modis、landsat/oli、gf1/wfv传感器等的气溶胶光学厚度反演算法用于gf4卫星遥感数据。

技术实现要素：

本申请提供了一种气溶胶光学厚度监测方法，以解决不能直接将其他传感器的气溶胶光学厚度反演算法用于高分四号卫星的问题。

本申请提供了一种气溶胶光学厚度监测方法，包括以下步骤：

s1：获取高分四号遥感数据，包括高分四号多光谱数据和高分四号红外数据；

s2：将所述数据进行交叉辐射定标，得到各波段表观反射率数据；

s3：根据所述表观反射率数据划分浓密植被区和非浓密植被区；

s4：基于6s辐射传输模型构建所述高分四号多光谱数据蓝波段和红波段查找表；

s5：根据所述表观反射率数据计算红波段地表反射率；

s6：根据modis地表反射率产品计算蓝波段地表反射率；

s7：获取所述查找表中全部气溶胶光学厚度的参数值；

s8：根据所述红波段地表反射率和所述参数值计算出红波段表观反射率，根据所述蓝波段地表反射率和所述参数值计算出蓝波段表观反射率；

s9：将所述红波段地表反射率、所述红波段表观反射率和所述参数值进行插值计算，得到所述浓密植被区的气溶胶光学厚度；

s10：将所述蓝波段地表反射率、所述蓝波段表观反射率和所述参数值进行差值计算，得到所述非浓密植被区的气溶胶光学厚度；

s11：将所述浓密植被区域的气溶胶光学厚度与所述非浓密植被区域的气溶胶光学厚度进行合并，得出完整区域的气溶胶光学厚度分布图。

可选的，s1还包括对所述高分四号遥感数据进行几何校正预处理；所述预处理包括选取一定数量的几何控制点，选用局部三角网方法建立几何校正模型，提高所述高分四号遥感数据的空间定位精度。

可选的，各波段反射率数据由下式得出：

ρtoa,i＝ai·dni+bi；

其中，ρtoa,i为第i波段表观反射率；ai和bi分别为第i波段定标系数增益值和偏移值；dni为第i波段图像中像素观测值。

可选的，s3中还包括根据所述高分四号多光谱数据红波段和近红外波段计算归一化植被指数；所述植被指数由下式得出：

其中，ρtoa,4为高分四号多光谱数据红波段表观反射率；ρtoa,3为高分四号多光谱数据近红外波段表观反射率。

可选的，当所述植被指数大于0.5且所述表观反射率中的短波红外波段表观反射率小于0.1为浓密植被区判别条件。

可选的，所述红波段地表反射率和蓝波段地表反射率由下式得出：

其中，ρs为地表反射率，ρtoa为表观反射率，ρ0为大气程辐射，θs为太阳天顶角，θv为观测天顶角，t(θs)为太阳-地面大气路径透过率，t(θv)为传感器-地面大气路径透过率，s为大气底层向下的半球反射率。

由以上技术方案可知，本申请提供了一种气溶胶光学厚度监测方法，方法为获取高分四号遥感数据，对其进行交叉辐射定标提高遥感辐射亮度数据准确性，计算各个波段的表观反射率数据划分浓密植被地区和非浓密植被地区。基于6s辐射传输模型建立查找表；根据短波红外波段表观反射率推算出红波段地表反射率，根据红波段地表反射率与表观反射率数据以及查找表计算浓密植被地区气溶胶光学厚度。根据modis地表反射率产品计算蓝波段地表反射率，根据蓝波段地表反射率与表观反射率数据以及查找表计算非浓密植被地区气溶胶光学厚度。通过分区建模的方式实现陆地所有区域的气溶胶光学厚度遥感监测。以解决不能直接将其他传感器的气溶胶光学厚度反演算法用于高分四号卫星的问题。

附图说明

为了更清楚地说明本申请的技术方案，下面将对实施例中所需要使用的附图作简单地介绍，显而易见地，对于本领域普通技术人员而言，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。

图1为一种气溶胶光学厚度监测方法流程图；

图2为本申请中交叉辐射定标流程。

具体实施方式

下面将详细地对实施例进行说明，其示例表示在附图中。下面的描述涉及附图时，除非另有表示，不同附图中的相同数字表示相同或相似的要素。以下实施例中描述的实施方式并不代表与本申请相一致的所有实施方式。仅是与权利要求书中所详述的、本申请的一些方面相一致的系统和方法的示例。

参见图1，为一种气溶胶光学厚度监测方法的流程图。本申请提供的一种气溶胶光学厚度监测方法，包括以下步骤：

s1：获取高分四号遥感数据，包括高分四号多光谱数据和高分四号红外数据；

s2：将所述数据进行交叉辐射定标，得到各波段表观反射率数据；

s3：根据所述表观反射率数据划分浓密植被区和非浓密植被区；

s4：基于6s辐射传输模型构建所述高分四号多光谱数据蓝波段和红波段查找表；

s5：根据所述表观反射率数据计算红波段地表反射率；

s6：根据modis地表反射率产品计算蓝波段地表反射率；

s7：获取所述查找表中全部气溶胶光学厚度的参数值；

s8：根据所述红波段地表反射率和所述参数值计算出红波段表观反射率，根据所述蓝波段地表反射率和所述参数值计算出蓝波段表观反射率；

s9：将所述红波段地表反射率、所述红波段表观反射率和所述参数值进行插值计算，得到所述浓密植被区的气溶胶光学厚度；

s10：将所述蓝波段地表反射率、所述蓝波段表观反射率和所述参数值进行差值计算，得到所述非浓密植被区的气溶胶光学厚度；

s11：将所述浓密植被区域的气溶胶光学厚度与所述非浓密植被区域的气溶胶光学厚度进行合并，得出完整区域的气溶胶光学厚度分布图。

进一步地，s1还包括对所述高分四号遥感数据进行几何校正预处理；所述预处理包括选取一定数量的几何控制点，选用局部三角网方法建立几何校正模型，提高所述高分四号遥感数据的空间定位精度。具体地，由于高分四号(gf4)遥感数据数据空间定位存在一定误差以及数据为地球同步轨道卫星进行面阵扫描成像，不同地区变形方式不同。先对数据进行几何校正预处理。对照googlemap等高分辨率影像地图在gf4遥感数据中选取一定数量的几何控制点，针对该数据高轨面阵成像导致不同地区变形方式不同的特点，不选用常用的多项式拟合方法而选用局部三角网方法建立几何校正模型。提高gf4遥感数据的空间定位精度。对gf4/irs短波红外波段运用双线性插值进行空间重采样，使其空间分辨率与gf4/pms多光谱保持一致。

进一步地，参见图2为本申请中交叉辐射定标流程。建立6s模型模拟gf4/pms、gf4/irs、eos/modis表观反射率。分别考虑成像时刻卫星天顶角、相对方位角等角度的影响，以及两个通道之间光谱响应的差异。将eos/modis数据进行卫星天顶角校正，获取gf4/pms、gf4/irs卫星天顶角下对应eos/modis数据。校正公式如下：

其中，modis(θ1,α)是gf4/pms、gf4/irs卫星天顶角下对应eos/modis辐射亮度值，modis(θ,α)是eos/modis卫星天顶角下对应eos/modis辐射亮度值。θ是eos/modis卫星天顶角，α是相对方位角，a1，a2，b1，b2分别为转换系数。

进一步地，对eos/modis数据再进行相对方位角的校正，获取gf4/pms、gf4/irs相对方位角下对应eos/modis数据。之后建立gf4/pms、gf4/irs和eos/modis之间的光谱关系，进行gf4/pms、gf4/irs的光谱差异校正，获得经过角度校正和光谱校正后的gf4/pms、gf4/irs表观反射率。将获得的gf4/pms、gf4/irs表观反射率与同一位置的gf4/pms、gf4/irs原始像素观测值(dn值)进行线性回归分析，获得经过交叉辐射定标后的定标系数。使用反射率定标系数可以消除大气层顶太阳辐照度来源的差异。其中，线性公式如下：

ρtoa,i＝ai·dni+bi；

其中，ρtoa,i为第i波段表观反射率；ai和bi分别为第i波段定标系数增益值和偏移值；dni为第i波段图像中像素观测值。

进一步地，划分浓密植被区和非浓密植被区。浓密植被的可见光波段地表反射率与短波红外波段地表反射率之间具有线性关系，可以依据此关系推算气溶胶光学厚度。但是该关系不适用于非浓密植被地区，因此需要区分浓密植被区和非浓密植被区采用不同针对性算法进行气溶胶光学厚度监测。根据所述高分四号多光谱数据红波段(第3波段)和近红外波段(第4波段)计算归一化植被指数；所述植被指数由下式得出：

其中，ρtoa,4为高分四号多光谱数据红波段表观反射率；ρtoa,3为高分四号多光谱数据近红外波段表观反射率。当所述植被指数大于0.5且所述表观反射率中的短波红外波段表观反射率小于0.1为浓密植被区判别条件。

进一步地，构建gf4/pms查找表，利用大气辐射传输模型构建查找表，其地表反射率表示为：

其中，ρs为地表反射率，ρtoa为表观反射率，ρ0为大气程辐射，θs为太阳天顶角，θv为观测天顶角，t(θs)为太阳-地面大气路径透过率，t(θv)为传感器-地面大气路径透过率，s为大气底层向下的半球反射率。

设置相应的6s参数：目标高度设为零；传感器高度为卫星观测高度；气溶胶光学厚度范围设置为0.0001-1.95，共10个数值；太阳天顶角以10°为步长设9个值；观测天顶角同样以10°为步长设9个值；相对方位角以10°为步长设置19个值(将卫星方位角固定为0°，此时太阳方位角的取值即为相对方位角的值)；气溶胶模式选择大陆型气溶胶模式。建立查找表时，气溶胶光学厚度、太阳天顶角、观测天顶角、相对方位角的设置见表1。

表1查找表参数设置

进一步地，进行浓密植被区气溶胶光学厚度监测。浓密植被短波红外波段地表反射率与红波段地表反射率存在线性关系。此外，在植被覆盖茂密区域，大气气溶胶对短波红外波谱区域的干扰很微小甚至没有，进而可以认为地面目标物的反射率与传感器观测得到的表现反射率基本一致。据此可以由gf4/irs短波红外波段表观反射率(近似为短波红外波段地表反射率)计算得到红波段地表反射率。首先逐像元读取gf4/pms的太阳天顶角、观测天顶角和相对方位角数据，根据查找表计算出10个气溶胶光学厚度值(0.0001，0.1，0.25，0.5，0.75，1.0，1.25，1.5，1.7，1.95)对应的红波段s、ρ0、t(θs)·t(θv)值，然后基于该像元红波段地表反射率计算得到对应于这10个气溶胶光学厚度的表观反射率，根据红波段真实的表观反射率、上面10个表观反射率数据及其对应光学厚度值进行线性插值，得到该像元气溶胶光学厚度。最后对所有浓密植被区像元进行循环计算，得到气溶胶光学厚度值。

进一步地，进行非浓密植被区气溶胶光学厚度监测针对非浓密植被区域，在进行气溶胶光学厚度反演前需对gf4/pms进行角度差异校正。采用半经验核驱动模型进行角度效应校正。模型公式如下：

ρ(θi,θv,φ)＝α0+α1kvol(θi,θv,φ)+α2kgeo(θi,θv,φ)；

式中，kvol为典型辐射传输类型体散射核，kgeo为几何光学核，α1和α2分别为核系数，α0为常数项，θi为太阳天顶角，θv为卫星观测天顶角，为相对方位角，g为相角。核函数与太阳天顶角、卫星天顶角、相对方位角和相角均有关。kvol和kgeo核函数以及相角g如下的表达式分别如下：

在角度校正过程中，根据之前交叉辐射定标时模拟的数据，利用最小二乘法回归拟合roujean模型中的α0、α1和α2。通过求解计算出r1和r0，结合gf4/pms蓝光波段的表观反射率ρs，参见下式计算校正后的蓝光波段反射率ρr。

进一步地，建立gf4/pms蓝光波段地表反射率库。由于gf4/pms卫星与eos/modis在蓝光波段的波段响应函数差异，需要采用以下公式进行通道光谱响应差异校正：

ρs,pms＝a+b*ρs,modis；

式中，ρs,pms、ρs,modis分别为gf4/pms、modis传感器蓝光波段观测的地表反射率；系数a、b利用asd实测数据结合两种传感器的通道光谱响应进行拟合。在系数a、b已知的基础上，将gf4/pms成像日期前后各10天内的eos/modis地表反射率逐日产品进行无云多时相合成，从中提取出蓝波段无云地表反射率数据代入公式中，得到gf4/pms蓝光波段对应的地表反射率。

式中，fλ,gf、fλ,modis分别为gf4/pms、modis传感器在λ波长处的光谱响应；δλgf、δλmodis分别为两种传感器光谱定标采用的波长间隔；ρi,λ为第i中地物在λ处的光谱反射率，利用便携式地物光谱仪asd测量的草地、地砖、水体、水泥、沥青等不同地物的光谱反射率。利用上述两公式拟合出不同地物对应的gf4/pms、modis蓝光波段反射率后，利用光谱响应差异校正公式即可回归得到通道光谱响应差异的校正系数a、b。计算得到gf4/pms蓝光波段对应的地表反射率。逐像元读取gf4/pms的太阳天顶角、观测天顶角和相对方位角数据，根据查找表计算出10个气溶胶光学厚度值(0.0001，0.1，0.25，0.5，0.75，1.0，1.25，1.5，1.7，1.95)对应的蓝波段s、ρ0、t(θs)·t(θv)值，然后基于该像元蓝波段地表反射率计算得到对应于这10个气溶胶光学厚度的表观反射率，根据蓝波段真实表观反射率、上面10个表观反射率数据及其对应光学厚度值进行线性插值，得到该像元气溶胶光学厚度。最后对所有非浓密植被区像元进行循环计算，得到气溶胶光学厚度值。

最终，将浓密植被区域气溶胶光学厚度与非浓密植被区域气溶胶光学厚度进行合并，得出完整区域的气溶胶光学厚度分布图。

由以上技术方案可知，本申请提供了一种气溶胶光学厚度监测方法，方法为获取高分四号遥感数据，对其进行交叉辐射定标提高遥感辐射亮度数据准确性，计算各个波段的表观反射率数据划分浓密植被地区和非浓密植被地区。基于6s辐射传输模型建立查找表；根据短波红外波段表观反射率推算出红波段地表反射率，根据红波段地表反射率与表观反射率数据以及查找表计算浓密植被地区气溶胶光学厚度。根据modis地表反射率产品计算蓝波段地表反射率，根据蓝波段地表反射率与表观反射率数据以及查找表计算非浓密植被地区气溶胶光学厚度。通过分区建模的方式实现陆地所有区域的气溶胶光学厚度遥感监测。以解决不能直接将其他传感器的气溶胶光学厚度反演算法用于高分四号卫星的问题。

本申请提供的实施例之间的相似部分相互参见即可，以上提供的具体实施方式只是本申请总的构思下的几个示例，并不构成本申请保护范围的限定。对于本领域的技术人员而言，在不付出创造性劳动的前提下依据本申请方案所扩展出的任何其他实施方式都属于本申请的保护范围。