CN106407656A - 一种基于高分辨率卫星影像数据的气溶胶光学厚度反演方法 - Google Patents

Abstract

本发明公开了一种基于高分辨率卫星影像数据的气溶胶光学厚度反演方法，具体为1)根据6S辐射传输模型建立查找表；2)高分辨率数据数据预处理，包括辐射定标，几何校正，云检测，获得原始表观反射率，与观测角度信息，根据观测角度信息与查找表获得大气参数.3)计算每个像元的归一化植被指数，根据植被指数与本发明提出的先验知识，确定每个像元对应的红蓝波段关系。4)根据卫星观测的表观反射率、大气参数与红蓝波段关系反演气溶胶光学厚度。本发明的气溶胶光学厚度遥感反演方法可以有效的对高分辨率卫星进行气溶胶监测，可以为区域及城市的大气环境与污染提供数据源。

Description

一种基于高分辨率卫星影像数据的气溶胶光学厚度反演方法

技术领域

本发明涉及遥感领域，尤指一种基于高分辨率卫星影像数据的气溶胶光学厚度反演方法。

背景技术

大气气溶胶是大气最主要的不确定性的成分之一，可用来估算气溶胶对区域以及全球气候变化的影响。目前各种气溶胶光学厚度反演方法(如暗目标算法、深蓝算法)被广泛的应用于各个卫星影像上，但是由于目前区域尺度的气溶胶光学特性越来越受到重视，区域高分辨率的气溶胶光学特性的监测也越来越得到关注。更加精细的气溶胶分布可以用来监测城市地区的气溶胶空间分布，能够有效快速的查找出污染源，结合其他气象数据可以提供更加精确的大气环境监测与预报。目前，有一些研究者致力于高分辨率的气溶胶光学厚度反演的研究，例如，王中挺等利用基于红蓝波段比值的暗目标方法，对GF-1影像的暗像元进行了反演。Benas等通过合成MERIS(Medium Resolution lmaging Spectrometer)和AATSR(Advanced Along Track Scanning Radiometer)的光谱和角度信息反演了意大利上空1km高分辨率气溶胶光学厚度；Li等通过协同MODIS数据以及环境星(HJ)CCD数据，反演了北京城市上空100m分辨率气溶胶光学特性。尽管如此，目前高分辨率影像气溶胶光学厚度反演还主要集中在暗目标以及多源卫星数据协同反演上，对于单一高分辨率卫星传感器以及城市等高亮地表区域，目前的高分辨率气溶胶光学厚度反演算法并不适用，反演结果的空间覆盖程度以及反演结果精度也有较大的损失。

发明内容

针对目前高分辨率气溶胶反演算法主要集中以暗目标以及多元卫星数据协同反演的问题，提高高分辨率大气气溶胶在环境监测上的运用，提高反演结果空间覆盖率以及精度等问题，提出一种基于高分辨率卫星影像数据的气溶胶光学厚度反演方法。

为实现上述目的，本发明的一种基于高分辨率卫星影像数据的气溶胶光学厚度反演方法，具体为1)根据6S辐射传输模型建立查找表；2) 高分辨率数据数据预处理，包括辐射定标，几何校正，云检测，获得原始表观反射率，与观测角度信息，根据观测角度信息与查找表获得大气参数.3)计算每个像元的归一化植被指数，根据植被指数与本发明提出的先验知识，确定每个像元对应的红蓝波段关系。4)根据卫星观测的表观反射率、大气参数与红蓝波段关系反演气溶胶光学厚度。本发明的气溶胶光学厚度遥感反演方法可以有效的对高分辨率卫星进行气溶胶监测，可以为区域及城市的大气环境与污染提供数据源。

进一步，步骤1)采用6S辐射传输模型进行辐射传输计算得出多组不同的大气参数(F_λ(θ_s)，T_λ(θ)，S_λ和)与观测角度组合而成的查找表。

进一步，步骤2)云检测具体按照公式：

[ρ_blue(D)-ρ_blue(D_r)]＞0.03*[1+(D-D_r)/30]选取像元是否为云像元。其中ρ_blue(D)是D日的像元的蓝波段反射率，D_r是在D日前某个无云的日期。

进一步，步骤3)具体的红蓝波段之间的关系通过MODISBRDF产品计算得到与归一化植被指数(NDVI)的关系，具体可以描述为(图1)：

K＝2.022；b＝-0.002；0.10≤NDVI＜0.50

K＝1.558；b＝0.010；0.50≤NDVI＜0.75

K＝0.974；b＝0.016；0.75≤NDVI＜1.00

不同的NDVI范围对应着不同的红蓝波段关系，其中NDVI＝(NIR-R)/(NIR+R)，通过近红外通道卫星观测的反射率与红波段反射率计算得到。

进一步，步骤4)具体为：根据读取的太阳天顶角、太阳方位角、观测天顶角、观测方位角，在查找表选取相应的数据，进行线性插值，得到红波段与蓝波段、不同气溶胶光学厚度的F_λ(θ_s)，T_λ(θ)，S_λ和大气参数；结合步骤3)得到的纯像元的理论上的红蓝波段的关系，按照方程组

根据卫星观测的太阳天顶角θ_s，传感器天顶角φ_s，是太阳方位角θ_v，是传感器方位角φ_v，对步骤1)得到的查找表进行插值，得到大气参数路径程辐射反射率ρ⁰，太阳——目标、目标——传感器大气路径透射率F、T，下界半球反照率S。再根据卫星观测的表观反射率ρ^TOA计算每个预设气溶胶光学厚度对应的红波段与蓝波段地表反射率ρ^s，与按照步骤3)所确定的红蓝波段地表反射率关系进行比较，最符合理论关系的气溶胶光学厚度值即为所求。

本发明的气溶胶光学厚度遥感反演方法可以有效的利用高分辨率卫星数据对气溶胶光学厚度反演，反演结果可以对城市地区等亮地表地区的气溶胶进行监测，增加了空间覆盖与空间分辨率，提高了气溶胶光学厚度反演精度。

附图说明

图1为不同归一化植被指数下(NDVI)，红蓝波段地表反射率线性关系。

图2为晴朗天气下气溶胶光学厚度反演MODIS官方产品结果(a)与高分辨率卫星(高分一号)气溶胶光学厚度反演结果(b)。

图3为污染天气下气溶胶光学厚度反演MODIS官方产品结果(a)与高分辨率卫星(高分一号)气溶胶光学厚度反演结果(b)。

图4为高分辨率卫星(高分一号)反演结果与AERONET对比验证图。

图5为高分辨率卫星(高分一号)污染源真彩色影像(a)与高分辨率卫星(高分一号)污染源监测结果(b)。

图6为高分辨率卫星(高分一号)华北地区气溶胶光学厚度反演结果。

具体实施方式

与现有的气溶胶反演的传感器相比，高分辨率卫星传感器(例如高分一号卫星)具有极高的空间分辨率，更高的空间分辨率能够反映区域及城市气溶胶光学厚度的空间分布状况，反映城市大气环境状况，能够找出污染源，结合其他气象数据可以对某一区域的污染进行综合分析。

因此，本方法充分发挥高分辨率卫星相机数据高空间分辨率的优势，利用像元NDVI值与红蓝波段地表反射率关系，从观测信号中分离出地表的贡献，得到气溶胶光学厚度。

对于平行大气下的朗伯体地表表面，大气顶层的上行反射率是地表与大气辐射相互作用产生的，是地表与大气耦合的函数。可以表示为

式中：分别是太阳天顶角，太阳方位角，相对方位角，T_g是气体吸收透过率，ρ₀是由分子散射加气溶胶散射所构成的路径程辐射反射率，F_λ(θ_s)、T_λ(θ)分别为太阳——目标、目标——传感器大气路径透射率，S_λ为大气下届半球反照率，是地物目标反射率，ρ_TOA是卫星观测的大气顶层反射率。

对于地表反射率，像元在不同的归一化植被指数下红蓝波段呈现线性关系：

结合公式(1)、(2)以及归一化植被指数与红蓝波段的关系(图2)，气溶胶光学厚度可以从以下方程组计算得到：

反演流程：

1)根据辐射传输模型6S建立查找表；查找表是通过设定不同卫星观测几何参数(8个太阳天顶角、14个观测天顶角、15个太阳与卫星之间的相对方位角)，不同的大气气溶胶参数(Lee和Kim 6种聚类的气溶胶模式、22个0.55微米处大气气溶胶光学厚度值：即0.0-2.0步长为0.1以及2.5)，考虑要反演所需要的波段(红蓝)，并考虑不同气溶胶类型参数，使用6S辐射传输程序进行辐射传输计算得出多组不同的大气参数(F_λ(θ_s)，T_λ(θ)，S_λ和)组合而成的查找表。

2)读取高分辨率影像数据，利用官方提供的定标系数，对原始数据进行定标，得到大气顶层表观反射率。依据云掩模方法对整幅影像进行云检测，提取非云像元进行进一步的反演。云检测方法是多时相的云检测方法(MTCD)，考虑到遥感影像的客观性和连续性，该算法以同一区域相近时相相近日期的卫星遥感地表反射率数据为参考，选取云阴影和典型地表样本点，对样本点进行统计分析，据其动态确定云阴影检测的阈值。

3)确定像元红蓝地表反射率关系：计算每个有效像元的NDVI值，再根据图1所示的关系，依照NDVI选取相应的红蓝波段地表反射率的关系。

4)根据红蓝波段地表反射率关系以及红蓝波段辐射传输方程建立方程组(公式3)。根据读取高分辨率数据的太阳天顶角、太阳方位角、观测天顶角、观测方位角，在查找表选取相应的大气参数，进行线性插值，得到不同波段、不同气溶胶光学厚度的F_λ(θ_s)，T_λ(θ)，S_λ和大气参数。结合方程组(3)，计算得到像元地表反射率与气溶胶光学厚度。

按照上述步骤，以高分一号为例，分别对晴朗天气下与重污染天气下的高分辨率影像进行气溶胶光学厚度反演。晴朗天气反演结果如图2(b)所示。整个华北平原的气溶胶光学厚度水平较低，城中心例如北京、天津的气溶胶光学厚度更高，这与MODIS官方产品结果(图2(a))空间分布一致，但是空间分辨更高，细节更加明显。同样的，在较重的大气污染条件下，反演结果如图3(b)所示整个华北平原都处在气溶胶光学厚度较高的情况下，除了部分气溶胶蔓延至山谷地带之外北部山区的由于地形原因气溶胶光学厚度整体较低，与MODIS官方产品结果(图3(a))空间分布一致。无论区域处在晴空大气下还是污染大气下，反演结果都可以有效的刻画区域气溶胶空间分布及变化情况，空间分辨更高，细节更加突出，有利于对城市地区的气溶胶分布进行研究。

为了更加准确地评估评估反演方法的性能，使用AERONET的地基观测数据进行对比验证。选取的地基数据为卫星过境时间前后30分钟之内，确保卫星过境与地基观测之间大气状况稳定。采用2015年上半年的数据进行验证，共有75的验证点进行验证，图4表示高分辨率气溶胶光学厚度反演结果与地基的相关性较好，相关系数为R＝0.96，据方根误差为0.185，相关性和均方根误差要优于其他高分辨率算法，拟合关系的截距为0.106，表明使用我们改进的基于NDVI红蓝波段地表反射率线性关系进行的气溶胶光学厚度反演的结果，大大提高了高分辨卫星反演精度，无论是在晴朗天气下还是重污染天气下，气溶胶光学厚度反演结果更加接近于地基观测值。

综合来看，利用本方法，高分辨率卫星数据可以有效的用于城区高分辨率气溶胶监测，更高的分辨率能够查找大气污染的污染源(图5)，覆盖更全，能为大气校正以及高分辨率大气环境变化研究提供基础数据(图6)。

Claims (8)

1.一种基于高分辨率卫星影像数据的气溶胶光学厚度反演方法，具体为1)根据6S辐射传输模型建立查找表；2)高分辨率数据数据预处理，包括辐射定标，几何校正，云检测，获得原始表观反射率，与观测角度信息，根据观测角度信息与查找表获得大气参数.3)计算每个像元的归一化植被指数，根据植被指数与本发明提出的先验知识，确定每个像元对应的红蓝波段关系。4)根据卫星观测的表观反射率、大气参数与红蓝波段关系反演气溶胶光学厚度。

2.如权利要求1所述的气溶胶光学厚度遥感反演方法，其特征在于，步骤1)利用6S辐射传输模型建立查找表。

3.如权利要求1所述的气溶胶光学厚度遥感反演方法，其特征在于，步骤2)云检测利用多时相的云检测方法(MTCD)，具体按照公式：

[ρ_blue(D)-ρ_blue(D_r)]＞0.03*[1+(D-D_r)/30]选取像元是否为云像元。其中ρ_blue(D)是D日的像元的蓝波段反射率，D_r是在D日前某个无云的日期。

4.如权利要求1所述的气溶胶光学厚度遥感反演方法，其特征在于，步骤3)根据卫星观测的表观反射率、大气参数与红蓝波段地表反射率关系反演气溶胶光学厚度。

5.如权利要求4所述的红蓝波段地表反射率关系，其特征在于计算每个像元的归一化植被指数，根据像元的植被指数与本发明提出的改进红蓝波段关系先验知识，确定每个像元对应的红蓝波段关系。

6.如权利要求5所述的本发明提出的改进红蓝波段关系先验知识具体表示为：

NDVI＝(NIR-R)/(NIR+R)

K＝2.012；b＝-0.002；0.10≤NDVI＜0.50

K＝1.593；b＝0.009；0.50≤NDVI＜0.75

K＝0.988；b＝0.016；0.75≤NDVI＜1.00

其中NDVI是归一化植被指数，NIR是高分辨率影像近红外波段表观反射率，R是高分辨率影像红外波段表观反射率，K和b是通过NDVI确定的红蓝波段地表反射率关系系数。

7.如权利要求1所述的气溶胶光学厚度遥感反演方法，其特征在于，步骤4)根据红蓝波段地表反射率关系以及红蓝波段辐射传输方程，通过查找表插值出来的大气参数，建立方程组，根据方程组求解得到气溶胶光学厚度。

8.如权利要求6所述的气溶胶光学厚度遥感反演方法，其特征在于，步骤4)查找表插值获得的大气参数、红蓝波段辐射传输方程以及红蓝波段地表反射率关系组成的方程组具体为：

${\mathrm{\ρ}}_{red}^{TOA}={\mathrm{\ρ}}_{red}^0+F_{red}{T}_{red}{\mathrm{\ρ}}_{red}^s/(1-S_{red}{\mathrm{\ρ}}_{red}^s)$

${\mathrm{\ρ}}_{blue}^{TOA}={\mathrm{\ρ}}_{blue}^0+F_{blue}{T}_{blue}{\mathrm{\ρ}}_{blue}^s/(1-S_{blue}{\mathrm{\ρ}}_{blue}^s)$

${\mathrm{\ρ}}_{red}^s=K*{\mathrm{\ρ}}_{blue}^s+b$

根据读取高分辨率遥感影像数据的太阳天顶角、太阳方位角、观测天顶角、观测方位角，在查找表选取相应的数据，进行线性插值，得到不同波段、不同气溶胶光学厚度的F_λ(θ_s)，T_λ(θ)，S_λ和ρ₀ 大气参数。结合方程组(3)，计算得到像元地表反射率与气溶胶光学厚度。其中θ_s，θ_v，分别是太阳天顶角，太阳方位角，相对方位角，T_g是气体吸收透过率，ρ₀是由分子散射加气溶胶散射所构成的路径程辐射反射率，F_λ(θ_s)、T_λ(θ)分别为太阳——目标、目标——传感器大气路径透射率，S_λ为大气下届半球反照率，是地物目标反射率，ρ_TOA是卫星观测的大气顶层反射率。K和b是通过NDVI确定的红蓝波段地表反射率关系系数。