基于非线性地表反射率模型的气溶胶遥感定量反演方法与流程

1.本发明涉及物联网技术领域，具体而言，涉及一种基于非线性地表反射率模型的气溶胶遥感定量反演方法。


背景技术：

2.在如今对地观测时代，从大气粒子的辐射传输理论，到卫星探测技术、气溶胶反演原理，再到产品生产投入业务化运行，整个流程已相对成熟。基于卫星遥感观测手段开发的全球覆盖的气溶胶光学厚度(aerosol optical depth,aod)产品已非常丰富，为研究任意尺度的空气质量提供数据支撑。气溶胶反演算法主要以经典的暗像元算法(dark target,dt)和深蓝算法(deep blue,db)为主(hsu等，2013；levy等，2013)。前者用于暗像元，后者用于亮地表像元。
3.dt算法假定气溶胶光学厚度一般随波长的增大而减小，认为在短波近红外通道几乎不受大气信号的影响，可见光通道0.47和0.65地表反射率与短波近红外通道的表观反射率或地表反射率之间存在某种线性统计关系，在已知地表反射率的前提下，进而实现气溶胶反演。该思想被modis团队应用在terra/modis和aqua/modis气溶胶业务化算法中。dt算法也被移植到许多具备相似观测能力的传感器上，例如，suomi-npp/viirs、fy3d/mersi-ii。也有一些学者将dt算法应用至静止卫星上h8。
4.深蓝算法(dark blue algorithm,db)主要针对沙漠、干旱、半干旱的亮地表而设计，利用深蓝波段(或近紫外)地表反射率比较低、而气溶胶信号贡献较大的特性，从长时间的表观反射率数据中挑选出最小值，进行大气订正后得到地表反射率库。在过去几十年，基于dt算法开展的研究比比皆是，针对载荷实际状态建立地表反射率模型是其中比较重要的一个改进。总体而言，新改进的地表反射率模型依然是随着观测几何和植被指数变化的线性关系，当临近像元地表类型或者观测角度变化较大时，就会出现地表反射率跳变的问题，进而导致定量反演的aod出现异常。


技术实现要素：

5.为了解决上述问题，本发明的目的在于提供一种基于非线性地表反射率模型的气溶胶遥感定量反演方法，通过筛选辅助aod数据集中小于0.2的样本作为目标载荷大气订正的背景值，并循环太阳、卫星观测角和方位角变化，借助大气辐射传输模型获得不同通道的地表反射率，通过分析不同季节下通道间地表反射率随观测几何、地表类型、地表植被参数等关键参数的变化规律和显著性，建立适用于每组地表反射率关系的数学模型，尽可能精确地表征地表反射率随季节、角度的变化规律，从而实现aod高精度反演。
6.为实现上述目的，本发明提供了以下技术方案：
7.根据本发明的第一方面，提供了一种基于非线性地表反射率模型的气溶胶遥感定量反演方法，该气溶胶遥感定量反演方法包括以下步骤：
8.1.读取目标卫星-载荷逐轨道一级数据的观测时间和辅助卫星-载荷aod产品的观
测时间，若二者时间差小于5分钟，执行下述过程，否则终止；
9.2.从辅助卫星-载荷aod数据集中筛选出aod数值小于0.2的样本，作为目标载荷大气订正的背景值；
10.3.读取目标卫星-载荷观测逐轨道一级数据每个像元的计数值并换算为大气层顶表观反射率，订正后得到气体吸收订正后的表观反射率；
11.4.读取目标卫星-载荷观测二级云检测产品，若云检测产品未标记为“冰雪”，则读取目标卫星-载荷观测二级云检测产品，若云检测产品标记为“一定是晴空”，进行下述过程，否则终止；
12.5.读取目标卫星-载荷多个通道反射率，若不符合内陆水体的判识条件，进行下述过程，否则终止；
13.6.读取辅助卫星-载荷aod产品中经度、纬度、海陆标识和aod 4个数据集，筛选出预设范围内目标卫星-载荷逐轨道一级数据的有效像元，对于辅助卫星-载荷每个像元，对筛选到的目标卫星-载荷的有效像元进行平均，得到时空匹配的目标卫星-载荷与辅助卫星-载荷的观测记录；
14.7.根据观测记录，循环太阳、卫星观测角和方位角变化，获得目标卫星-载荷的地表反射率，通过分析确定观测角分类设置的断点以及地表参数分类设置的断点，建立地表反射率模型。
15.作为本发明的进一步方案，从辅助卫星-载荷aod数据集中筛选出aod数值小于0.2的样本，包括：
16.读取辅助卫星-载荷aod逐轨道产品，对于每个像元，执行下述判断条件，若满足，继续进行下述过程，否则终止；
17.aod
auxiliary
＜0.2
18.其中，aod
auxiliary
为辅助卫星-载荷aod数据，0.2为筛选样本的阈值。
19.作为本发明的进一步方案，读取目标卫星-载荷观测逐轨道一级数据每个像元的计数值并换算为大气层顶表观反射率，包括：
20.读取目标卫星-载荷观测逐轨道一级数据每个像元于多个可见光通道、近红外通道、短波近红外通道的计数值，通过定标系数和物理量转换公式，将计数值换算为大气层顶表观反射率；
21.其中，所述可见光通道包含蓝光波段、绿光波段和红光波段，所述近红外通道包含两个近红外波段。
22.作为本发明的进一步方案，订正后得到气体吸收订正后的表观反射率时，利用含水汽、臭氧、二氧化碳的再分析资料对目标卫星-载荷相关波段的表观反射率进行气体吸收订正，得到气体吸收订正后的表观反射率，其中，所述表观反射率的波段包括可见光通道、近红外通道、短波近红外通道。
23.作为本发明的进一步方案，读取目标卫星-载荷观测二级云检测产品，对于每个像元，若云检测产品标记为“冰雪”，不进行下述过程，过程终止。
24.作为本发明的进一步方案，读取目标卫星-载荷观测二级云检测产品，若云检测产品标记为“一定是晴空”，对于每个像元，执行下述内陆水体判断条件，若不满足，继续执行下述过程，否则终止；
25.判断是否为内体水体，若满足下述条件，过程终止，否则进行下述过程；
[0026][0027]
其中，ρ
nir1
和ρ
red
分别为近红外通道1和红光表观反射率，t1为判断为内陆水体的阈值；ρ
swir
为短波近红外通道表观反射率，t2为短波近红外通道反射率判断阈值。
[0028]
作为本发明的进一步方案，对筛选到的目标卫星-载荷的有效像元进行平均，得到时空匹配的{时间、经度、纬度、θs、θv、θ、海陆标识、aod、ρ
blue
、ρ
green
、ρ
red
、ρ
nir1
、ρ
nir2
、ρ
swir
、ndvi、ndvi
swir
、afri}观测记录，其中，θs为太阳天顶角、θv为卫星天顶角、为相对方位角、θ为散射角、ρ
blue
、ρ
green
、ρ
red
分别为目标卫星-载荷的蓝光、绿光和红光通道表观反射率，ρ
nir1
和ρ
nir2
分别为目标卫星-载荷的近红外通道1和近红外通道2表观反射率，ρ
swir
为目标卫星-载荷的短波近红外通道表观反射率，ndvi为近红外通道1表观反射率和红光通道表观反射率计算得到，ndvi
swir
为近红外通道1表观反射率和短波近红外通道表观反射率计算得到，afri为蓝光通道表观反射率和短波近红外通道表观反射率计算得到。
[0029]
作为本发明的进一步方案，建立地表反射率模型之前，根据目标卫星-载荷地表反射率关系显著变化的观测角，联合与与两组关系分析拟合参数和相关性变化的拐点，确定观测角分类设置的断点以及地表参数分类设置的断点，其中，根据选择相关性的强弱将参数划分为主变量和次变量，在主变量每个类别下，按照次变量分类标准，分别统计每类样本下与与与相关性最强的拟合关系，确定地表反射率模型；为计算得到蓝光通道地表反射率，为计算得到红光通道地表反射率，为计算得到短波近红外通道地表反射率。
[0030]
作为本发明的进一步方案，所述地表反射率模型用于代入到气溶胶反演算法中进行气溶胶光学厚度反演，所述气溶胶反演算法为基于多角度反射率分析的气溶胶反演算法。
[0031]
作为本发明的进一步方案，所述气溶胶反演算法为基于可见光和近红外反射率的气溶胶反演算法。
[0032]
根据本发明的第二方面，还提供了一种计算机设备，包括存储器和处理器，该存储器中存储有计算机程序，该计算机程序被处理器执行时执行上述任一项根据本发明的基于非线性地表反射率模型的气溶胶遥感定量反演方法。
[0033]
根据本发明的第三方面，还提供了一种计算机可读存储介质，存储有计算机程序指令，该计算机程序指令被执行时实现上述任一项根据本发明的基于非线性地表反射率模型的气溶胶遥感定量反演方法。
[0034]
与现有技术相比较而言，本发明提供的一种基于非线性地表反射率模型的气溶胶遥感定量反演方法具有以下有益效果：
[0035]
本发明通过筛选小于0.2的样本作为背景值，建立非线性地表反射率模型，考虑了地表反射率随观测几何和地表参数的非线性变化规律，提高了aod反演精度。同时，本发明利用ndvi、ndvi
swir
和afri等地表参数，对地表反射率随地表参数的变化规律进行了分析，更好地刻画了地表反射率的变化规律。本发明方法简单易行，适用性广，可用于各种卫星遥感数据的气溶胶遥感定量反演。
[0036]
本技术的这些方面或其他方面在以下实施例的描述中会更加简明易懂。应当理解的是，以上的一般描述和后文的细节描述仅是的和解释性的，并不能限制本技术。
附图说明
[0037]
为了更清楚地说明本技术实施例中的技术方案，下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本技术的一些实施例。在附图中：
[0038]
图1为本发明的实施例中一种基于非线性地表反射率模型的气溶胶遥感定量反演方法的流程图；
[0039]
图2为本发明的实施例中一种基于非线性地表反射率模型的气溶胶遥感定量反演方法中fy3d/mersi-ii与随θs和θ拟合关系的示意图；
[0040]
图3为本发明的实施例中一种基于非线性地表反射率模型的气溶胶遥感定量反演方法中fy3d/mersi-ii与随θs和θ拟合关系的示意图；
[0041]
图4为本发明的实施例中一种基于非线性地表反射率模型的气溶胶遥感定量反演方法中fy3d/mersi-ii与随ndvi、ndvi
swir
、afri拟合关系的示意图；
[0042]
图5为本发明的实施例中一种基于非线性地表反射率模型的气溶胶遥感定量反演方法中fy3d/mersi-ii与随ndvi、ndvi
swir
、afri拟合关系的示意图；
[0043]
图6为本发明的实施例中一种基于非线性地表反射率模型的气溶胶遥感定量反演方法中fy3d/mersi-ii在不同θs下与随ndvi
swir
拟合关系的示意图；
[0044]
图7为本发明的实施例中一种基于非线性地表反射率模型的气溶胶遥感定量反演方法中fy3d/mersi-ii在不同θs下与随ndvi
swir
拟合关系的示意图。
[0045]
图8为本发明中执行基于非线性地表反射率模型的气溶胶遥感定量反演方法的计算机设备的一个实施例的硬件结构示意图。
[0046]
图9为本发明实施例提供的基于非线性地表反射率模型的气溶胶遥感定量反演方法的计算机可读存储介质的示意图。
[0047]
本技术目的实现、功能特点及优点将结合实施例，参照附图做进一步说明。
具体实施方式
[0048]
下面，结合附图以及具体实施方式，对本技术做进一步描述，需要说明的是，在不相冲突的前提下，以下描述的各实施例之间或各技术特征之间可以任意组合形成新的实施例。
[0049]
应当理解，此处所描述的具体实施例仅仅用以解释本技术，并不用于限定本技术。
[0050]
下面将结合本技术实施例中的附图，对本技术实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例是本技术一部分实施例，而不是全部的实施例。基于本技术中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本技术保护的范围。
[0051]
附图中所示的流程图仅是示例说明，不是必须包括所有的内容和操作/步骤，也不
是必须按所描述的顺序执行。例如，有的操作/步骤还可以分解、组合或部分合并，因此实际执行的顺序有可能根据实际情况改变。
[0052]
下面结合附图，对本技术的一些实施方式作详细说明。在不冲突的情况下，下述的实施例及实施例中的特征可以相互组合。
[0053]
由于新改进的地表反射率模型依然是随着观测几何和植被指数变化的线性关系，当临近像元地表类型或者观测角度变化较大时，就会出现地表反射率跳变的问题，进而导致定量反演的aod出现异常。为了解决上述问题，本发明提供了一种基于非线性地表反射率模型的气溶胶遥感定量反演方法。
[0054]
在一些实施方式中，基于非线性地表反射率模型的气溶胶遥感定量反演方法可以应用在计算机设备，该计算机设备可以是pc、便携计算机、移动终端等具有显示和处理功能的设备，当然也不限于此。
[0055]
本技术实施例提供了一种基于非线性地表反射率模型的气溶胶遥感定量反演方法，该基于非线性地表反射率模型的气溶胶遥感定量反演方法包括以下步骤：
[0056]
步骤s10、读取目标卫星-载荷逐轨道一级数据的观测时间和辅助卫星-载荷aod产品的观测时间，若二者时间差小于5分钟，执行下述过程，否则终止；
[0057]
步骤s20、从辅助卫星-载荷aod数据集中筛选出aod数值小于0.2的样本，作为目标载荷大气订正的背景值；
[0058]
步骤s30、读取目标卫星-载荷观测逐轨道一级数据每个像元的计数值并换算为大气层顶表观反射率，订正后得到气体吸收订正后的表观反射率；
[0059]
步骤s40、读取目标卫星-载荷观测二级云检测产品，若云检测产品未标记为“冰雪”，则读取目标卫星-载荷观测二级云检测产品，若云检测产品标记为“一定是晴空”，进行下述过程，否则终止；
[0060]
步骤s50、读取目标卫星-载荷多个通道反射率，若不符合内陆水体的判识条件，进行下述过程，否则终止；
[0061]
步骤s60、读取辅助卫星-载荷aod产品中经度、纬度、海陆标识和aod 4个数据集，筛选出预设范围内目标卫星-载荷逐轨道一级数据的有效像元，对于辅助卫星-载荷每个像元，对筛选到的目标卫星-载荷的有效像元进行平均，得到时空匹配的目标卫星-载荷与辅助卫星-载荷的观测记录；
[0062]
步骤s70、根据观测记录，循环太阳、卫星观测角和方位角变化，获得目标卫星-载荷的地表反射率，通过分析确定观测角分类设置的断点以及地表参数分类设置的断点，建立地表反射率模型。
[0063]
本发明提供了一种基于非线性地表反射率模型的气溶胶遥感定量反演方法，通过筛选辅助aod数据集中小于0.2的样本作为目标载荷大气订正的背景值，并循环太阳、卫星观测角和方位角变化，借助大气辐射传输模型获得不同通道的地表反射率，通过分析不同季节下通道间地表反射率随观测几何、地表类型、地表植被参数等关键参数的变化规律和显著性，建立适用于每组地表反射率关系的数学模型，尽可能精确地表征地表反射率随季节、角度的变化规律，从而实现aod高精度反演。
[0064]
在进行气溶胶遥感定量反演操作时，包括以下步骤：
[0065]
1)读取目标卫星-载荷逐轨道一级数据的观测时间，读取辅助卫星-载荷aod产品
的观测时间，若前者与后者观测时间差小于5分钟，执行下述过程，否则终止；
[0066]
2)读取辅助卫星-载荷aod逐轨道产品，对于每个像元，执行下述判断条件，若满足，继续进行下述过程，否则终止；
[0067]
aod
auxiliary
＜0.2
[0068]
其中，aod
auxiliary
为辅助卫星-载荷aod数据，0.2为筛选样本的阈值。
[0069]
3)读取目标卫星-载荷观测逐轨道一级数据每个像元于多个可见光、近红外、短波近红外的计数值，通过定标系数和物理量转换公式，将计数值换算为大气层顶表观反射率。其中所述可见光通道包含蓝光波段、绿光波段和红光波段，分别简写为blue，green和red，近红外通道包含两个，分别简写为nir1和nir2，短波近红外通道简写为swir；
[0070]
4)利用水汽、臭氧、二氧化碳等辅助资料，对目标卫星-载荷相关波段的表观反射率进行气体吸收订正，得到气体吸收订正后的表观反射率。其中所述相关波段包括blue，green，red，nir1，nir2和swir；
[0071]
5)读取目标卫星-载荷观测二级云检测产品，对于每个像元，若云检测产品标记为“冰雪”，不进行下述过程，否则继续执行下述过程；
[0072]
6)读取目标卫星-载荷观测二级云检测产品，对于每个像元，若云检测产品标记为“一定是晴空”，进行下述过程，否则终止；
[0073]
7)对于每个像元，执行下述内陆水体判断条件，若不满足，继续执行下述过程，否则终止；判断是否为内体水体，若满足下述条件，过程终止，否则进行下述过程；
[0074][0075]
其中，ρ
nir1
和ρ
red
分别为近红外通道1和红光表观反射率，t1为判断为内陆水体的阈值；ρ
swir
为短波近红外通道表观反射率，t2为短波近红外通道反射率判断阈值。
[0076]
8)读取辅助卫星-载荷aod产品中经度、纬度、海陆标识和aod 4个数据集，对于每个像元，筛选出半径5km范围内目标卫星-载荷逐轨道一级数据的有效像元，若像元数大于2个，执行下述判断条件，否则终止；
[0077]
9)对于辅助卫星-载荷每个像元，对筛选到的目标卫星-载荷的有效像元进行平均，得到二者时空匹配的{时间、经度、纬度、θs、θv、θ、海陆标识、aod、ρ
blue
、ρ
green
、ρ
red
、ρ
nir1
、ρ
nir2
、ρ
swir
、ndvi、ndvi
swir
、afri}观测记录；其中，θs为太阳天顶角、θv为卫星天顶角、为相对方位角、θ为散射角、ρ
blue
、ρ
green
、ρ
red
分别为目标卫星-载荷的蓝光、绿光和红光通道表观反射率，ρ
nir1
和ρ
nir2
分别为目标卫星-载荷的近红外通道1和近红外通道2表观反射率，ρ
swir
为目标卫星-载荷的短波近红外通道表观反射率，ndvi为近红外通道1表观反射率和红光通道表观反射率计算得到，ndvi
swir
为近红外通道1表观反射率和短波近红外通道表观反射率计算得到，afri为蓝光通道表观反射率和短波近红外通道表观反射率计算得到；
[0078]
10)分别统计θs、θv、θ、ndvi、ndvi
swir
、afri 6个参数的最值、平均值分布特征，选择其中一个参数作为自变量，对每个参数取值设置19个断点分为20类，统计每个参数在每类样本下与与的线性拟合关系(斜率、截距)和相关性；依次变化自变量，直至将6个参数遍历完毕；
[0079]
11)对于观测几何θs、θv、θ，对比与与的拟合关系和相关
性在这三个变量下的变化情况，根据相关性高低，确定目标卫星-载荷地表反射率关系随哪个观测角显著变化；
[0080]
12)对于ndvi、ndvi
swir
、afri，对比与与的拟合关系和相关性在这三个变量下的变化情况，根据相关性高低，确定目标卫星-载荷地表反射率关系随哪个地表参数显著变化；
[0081]
13)对于目标卫星-载荷地表反射率关系显著变化的观测角而言，联合与与两组关系分析拟合参数和相关性变化的拐点，确定观测角分类设置的断点；
[0082]
14)对于目标卫星-载荷地表反射率关系显著变化的地表参数而言，联合与与两组关系分析拟合参数和相关性变化的拐点，确定地表参数分类设置的断点；
[0083]
15)对比13)和14)两组关系中相关性，选择相关性最强的参数作为主变量，相关性较弱的参数作为次变量；在主变量每个类别下，按照次变量分类标准，分别统计每类样本下与与相关性最强的拟合关系(可以是幂指数、对数、多项式、线性等关系中的一种)，进而确定整个地表反射率模型；
[0084]
16)基于上述地表反射率模型，代入到气溶胶反演算法中，实现气溶胶光学厚度反演。
[0085]
在本实施例中，所述地表反射率模型用于代入到气溶胶反演算法中进行气溶胶光学厚度反演，所述气溶胶反演算法为基于多角度反射率分析的气溶胶反演算法，所述气溶胶反演算法为基于可见光和近红外反射率的气溶胶反演算法。
[0086]
为了便于本领域普通技术人员理解和实施本发明，下面以本发明的方法应用在风云三号卫星mersi-ii传感器(中分辨率光谱成像仪)上，选择aqua/modis传感器作为辅助数据源，对本发明实施例进行详细说明。
[0087]
1)从mersi逐轨道一级数据文件名上读取观测时间，读取aqua/modis l2 aod产品逐轨道数据的观测时间，若前者与后者观测时间差小于5分钟，执行下述过程，否则终止；
[0088]
2)读取aqua/modis l2 aod逐轨道产品，对于每个像元，执行下述判断条件，若满足，继续进行下述过程，否则终止；
[0089]
aod
modis
＜0.2
[0090]
其中，aod
modis
为辅助卫星-载荷aod数据，0.2为筛选样本的阈值。
[0091]
3)读取fy-3d/mersi-ii l1数据每个像元在三个可见光通道(0.47μm、0.55μm和0.65μm)、两个近红外通道(0.865μm、1.03μm)、一个短波近红外通道(2.13μm)的计数值，通过定标系数和物理量转换公式，将以上通道计算值换算为大气层顶表观反射率；
[0092]
4)利用水汽、臭氧、二氧化碳等辅助资料，对0.47μm、0.55μm和0.65μm、0.865μm、1.03μm和2.13μm分别进行水汽吸收、臭氧吸收和二氧化碳吸收订正，得到气体吸收订正后的表观反射率；
[0093]
5)读取mersi-ii二级云检测产品，对于每个像元，若云检测产品标记为“冰雪”，不进行下述过程，否则继续执行下述过程；
[0094]
6)读取mersi-ii二级云检测产品，对于每个像元，若云检测产品标记为“一定是晴
空”，进行下述过程，否则终止；
[0095]
7)对于每个像元，执行下述内陆水体判断条件，若不满足，继续执行下述过程，否则终止；判断是否为内体水体，若满足下述条件，过程终止，否则进行下述过程；
[0096][0097]
其中，ρ
nir1
和ρ
red
分别为0.865μm和0.65μm表观反射率，t1为0.1，t2为0.08。
[0098]
8)读取aqua/aod产品中经度、纬度、海陆标识和aod 4个数据集，对于每个像元，筛选出半径5km范围内mersi-ii l1有效像元，若像元数大于2个，执行下述判断条件，否则终止；
[0099]
9)对于aqua/aod每个像元，对筛选到的mersi-ii有效像元进行平均，得到二者时空匹配的{时间、经度、纬度、θs、θv、θ、海陆标识、aod、ρ
blue
、ρ
green
、ρ
red
、ρ
nir1
、ρ
nir2
、ρ
swir
、ndvi、ndvi
swir
、afri}观测记录；其中，θs为太阳天顶角、θv为卫星天顶角、为相对方位角、θ为散射角、ρ
blue
、ρ
green
、ρ
red
分别为mersi-ii蓝光、绿光和红光通道表观反射率，ρ
nir1
和ρ
nir2
分别为mersi-ii近红外通道1和近红外通道2表观反射率，ρ
swir
为mersi-ii短波近红外通道表观反射率，ndvi为近红外通道1表观反射率和红光通道表观反射率计算得到，ndvi
swir
为近红外通道1表观反射率和短波近红外通道表观反射率计算得到，afri为蓝光通道表观反射率和短波近红外通道表观反射率计算得到；
[0100]
10)分别统计θs、θv、θ、ndvi、ndvi
swir
、afri 6个参数的最值、平均值分布特征，选择其中一个参数作为自变量，对每个参数取值设置19个断点分为20类，统计每个参数在每类样本下与与的线性拟合关系(斜率、截距)和相关性；依次变化自变量，直至将6个参数遍历完毕；
[0101]
11)对于观测几何θs、θv、θ，对比与与的拟合关系和相关性在这三个变量下的变化情况，根据相关性高低，确定mersi-ii地表反射率关系随哪个观测角显著变化；以θs和θ为例，图2和图3分别展示了与与与两组拟合关系随之变化趋势，从相关性r2来看，两组关系随θs显著变化。
[0102]
12)对于ndvi、ndvi
swir
、afri，对比与与的拟合关系和相关性在这三个变量下的变化情况，根据相关性高低，确定mersi-ii地表反射率关系随哪个地表参数显著变化；图4分别展示了与随着ndvi、ndvi
swir
、afri拟合参数的变化趋势，图5分别展示了与随着ndvi、ndvi
swir
、afri拟合参数的变化趋势；结合图4和图5的表现，以及ndvi、ndvi
swir
、afri值域在每类样本中的分布情况，选择mersi-ii地表反射率关系随着ndvi
swir
参数变化。
[0103]
13)对于mersi-ii地表反射率关系显著变化的观测角而言，联合与与两组关系分析拟合参数和相关性变化的拐点，确定观测角分类设置的断点；结合图2和图3，可确定θs在30
°
和55
°
作为断点，图2为2022年7月fy3d/mersi-ii与随θs(上图)和θ(下图)拟合关系的示意图，图3为2022年7月fy3d/mersi-ii与随θs(上图)和θ(下图)拟合关系的示意图；
[0104]
14)对于mersi-ii地表反射率关系显著变化的地表参数而言，联合与
与两组关系分析拟合参数和相关性变化的拐点，确定地表参数分类设置的断点；结合图4可确定ndvi
swir
的断点分别为-0.05，0.25和0.75，图4为2022年7月fy3d/mersi-ii与随ndvi(上图)、ndvi
swir
(中图)、afri(下图)拟合关系的示意图；结合图5可确定ndvi
swir
的断点分别为-0.05和0.75，图5为2022年7月fy3d/mersi-ii与随ndvi(上图)、ndvi
swir
(中图)、afri(下图)拟合关系示意图；
[0105]
15)对比13)和14)两组关系中相关性，选择相关性最强的参数作为主变量，相关性较弱的参数作为次变量；在主变量每个类别下，按照次变量分类标准，分别统计每类样本下与与相关性最强的拟合关系(可以是幂指数、对数、多项式、线性等关系中的一种)，进而确定整个地表反射率模型；图6分别展示了θs小于30
°
、30～55
°
、大于55
°
三组样本下与拟合关系(第一列斜率，第二列截距，第三列相关性)随着ndvi
swir
变化情况，图7和图6相似，但针对与图6为图6fy3d/mersi-ii在不同θs下与随ndvi
swir
拟合关系的示意图，图6中第一行为θs小于30
°
，第二行θs为30～55
°
，第三行为θs大于55
°
。图7为fy3d/mersi-ii在不同θs下与随ndvi
swir
拟合关系的示意图，图7中第一行为θs小于30
°
，第二行θs为30～55
°
，第三行为θs大于55
°
。
[0106]
基于上述地表反射率模型，代入到气溶胶反演算法中，实现气溶胶光学厚度反演。
[0107]
本发明通过筛选小于0.2的样本作为背景值，建立非线性地表反射率模型，考虑了地表反射率随观测几何和地表参数的非线性变化规律，提高了aod反演精度。同时，本发明利用ndvi、ndvi
swir
和afri等地表参数，对地表反射率随地表参数的变化规律进行了分析，更好地刻画了地表反射率的变化规律。本发明方法简单易行，适用性广，可用于各种卫星遥感数据的气溶胶遥感定量反演。
[0108]
本发明实施例的第三个方面，还提供了一种计算机设备1000，包括存储器1001和处理器1002，该存储器中存储有计算机程序，该计算机程序被该处理器执行时实现上述任意一项实施例的方法。
[0109]
如图8所示，为本发明提供的执行基于非线性地表反射率模型的气溶胶遥感定量反演方法的计算机设备的一个实施例的硬件结构示意图。以如图8所示的计算机设备1000为例，在该计算机设备中包括一个处理器1002以及一个存储器1001，并还可以包括：输入装置和输出装置。处理器1002、存储器1001、输入装置和输出装置可以通过总线或者其他方式连接，图8中以通过总线连接为例。输入装置可接收输入的数字或字符信息，以及产生与基于非线性地表反射率模型的气溶胶遥感定量反演有关的信号输入。输出装置可包括显示屏等显示设备。
[0110]
存储器1001作为一种非易失性计算机可读存储介质，可用于存储非易失性软件程序、非易失性计算机可执行程序以及模块，如本技术实施例中的基于非线性地表反射率模型的气溶胶遥感定量反演方法对应的程序指令/模块。存储器1001可以包括存储程序区和存储数据区，其中，存储程序区可存储操作系统、至少一个功能所需要的应用程序；存储数据区可存储基于非线性地表反射率模型的气溶胶遥感定量反演方法的使用所创建的数据等。此外，存储器1001可以包括高速随机存取存储器，还可以包括非易失性存储器，例如至少一个磁盘存储器件、闪存器件、或其他非易失性固态存储器件。在一些实施例中，存储器
1001可选包括相对于处理器1002远程设置的存储器，这些远程存储器可以通过网络连接至本地模块。上述网络的实例包括但不限于互联网、企业内部网、局域网、移动通信网及其组合。
[0111]
处理器1002在一些实施例中可以是中央处理器(central processing unit，cpu)、控制器、微控制器、微处理器、或其他数据处理芯片。该处理器1002通常用于控制计算机设备1000的总体操作。本实施例中，处理器1002用于运行存储器1001中存储的程序代码或者处理数据。本实施例计算机设备的多个计算机设备1000的处理器1002通过运行存储在存储器1001中的非易失性软件程序、指令以及模块，从而执行服务器的各种功能应用以及数据处理，即实现上述方法实施例的基于非线性地表反射率模型的气溶胶遥感定量反演方法的步骤。
[0112]
通过以上的实施方式的描述，本领域普通技术人员可以清楚地了解到各实施方式可借助软件加通用硬件平台的方式来实现，当然也可以通过硬件。本领域普通技术人员可以理解实现上述实施例方法中的全部或部分流程是可以通过计算机程序来指令相关的硬件来完成。
[0113]
本发明实施例的第四个方面，还提供了一种计算机可读存储介质，图9示出了根据本发明实施例提供的基于非线性地表反射率模型的气溶胶遥感定量反演方法的计算机可读存储介质的示意图。如图9所示，计算机可读存储介质2000存储有计算机程序指令2001，该计算机程序指令2001可以被处理器执行。该计算机程序指令2001被执行时实现上述任意一项实施例的方法，即实现上述方法实施例的基于非线性地表反射率模型的气溶胶遥感定量反演方法的步骤。
[0114]
应当理解，在相互不冲突的情况下，以上针对根据本发明的基于非线性地表反射率模型的气溶胶遥感定量反演方法阐述的所有实施方式、特征和优势同样地适用于根据本发明的基于非线性地表反射率模型的气溶胶遥感定量反演系统和存储介质。
[0115]
本领域技术人员还将明白的是，结合这里的公开所描述的各种的逻辑块、模块、电路和算法步骤可以被实现为电子硬件、计算机软件或两者的组合。为了清楚地说明硬件和软件的这种可互换性，已经就各种示意性组件、方块、模块、电路和步骤的功能对其进行了一般性的描述。这种功能是被实现为软件还是被实现为硬件取决于具体应用以及施加给整个系统的设计约束。本领域技术人员可以针对每种具体应用以各种方式来实现的功能，但是这种实现决定不应被解释为导致脱离本发明实施例公开的范围。
[0116]
最后需要说明的是，本文的计算机可读存储介质(例如，存储器)可以是易失性存储器或非易失性存储器，或者可以包括易失性存储器和非易失性存储器两者。作为例子而非限制性的，非易失性存储器可以包括只读存储器(rom)、可编程rom(prom)、电可编程rom(eprom)、电可擦写可编程rom(eeprom)或快闪存储器。易失性存储器可以包括随机存取存储器(ram)，该ram可以充当外部高速缓存存储器。作为例子而非限制性的，ram可以以多种形式获得，比如同步ram(dram)、动态ram(dram)、同步dram(sdram)、双数据速率sdram(ddr sdram)、增强sdram(esdram)、同步链路dram(sldram)、以及直接rambus ram(drram)。所公开的方面的存储设备意在包括但不限于这些和其它合适类型的存储器。
[0117]
结合这里的公开所描述的各种的逻辑块、模块和电路可以利用被设计成用于执行这里功能的下列部件来实现或执行：通用处理器、数字信号处理器(dsp)、专用集成电路
(asic)、现场可编程门阵列(fpga)或其它可编程逻辑器件、分立门或晶体管逻辑、分立的硬件组件或者这些部件的任何组合。通用处理器可以是微处理器，但是可替换地，处理器可以是任何传统处理器、控制器、微控制器或状态机。处理器也可以被实现为计算设备的组合，例如，dsp和微处理器的组合、多个微处理器、一个或多个微处理器结合dsp和/或任何其它这种配置。
[0118]
以上仅为本技术的优选实施例，并非因此限制本技术的专利范围，凡是利用本技术说明书及附图内容所作的等效结构或等效流程变换，或直接或间接运用在其他相关的技术领域，均同理包括在本技术的专利保护范围内。

技术特征：
1.一种基于非线性地表反射率模型的气溶胶遥感定量反演方法，其特征在于，该气溶胶遥感定量反演方法包括以下步骤：读取目标卫星-载荷逐轨道一级数据的观测时间和辅助卫星-载荷aod产品的观测时间，若二者时间差小于5分钟，执行下述过程，否则终止；从辅助卫星-载荷aod数据集中筛选出aod数值小于0.2的样本，作为目标载荷大气订正的背景值；读取目标卫星-载荷观测逐轨道一级数据每个像元的计数值并换算为大气层顶表观反射率，订正后得到气体吸收订正后的表观反射率；读取目标卫星-载荷观测二级云检测产品，若云检测产品未标记为“冰雪”，则读取目标卫星-载荷观测二级云检测产品，若云检测产品标记为“一定是晴空”，进行下述过程，否则终止；读取目标卫星-载荷多个通道反射率，若不符合内陆水体的判识条件，进行下述过程，否则终止；读取辅助卫星-载荷aod产品中经度、纬度、海陆标识和aod4个数据集，筛选出预设范围内目标卫星-载荷逐轨道一级数据的有效像元，对于辅助卫星-载荷每个像元，对筛选到的目标卫星-载荷的有效像元进行平均，得到时空匹配的目标卫星-载荷与辅助卫星-载荷的观测记录；根据观测记录，循环太阳、卫星观测角和方位角变化，获得目标卫星-载荷的地表反射率，通过分析确定观测角分类设置的断点以及地表参数分类设置的断点，建立地表反射率模型。2.根据权利要求1所述的基于非线性地表反射率模型的气溶胶遥感定量反演方法，其特征在于，从辅助卫星-载荷aod数据集中筛选出aod数值小于0.2的样本，包括：读取辅助卫星-载荷aod逐轨道产品，对于每个像元，执行下述判断条件，若满足，继续执行，否则终止；aod
auxiliary
<0.2其中，aod
auxiliary
为辅助卫星-载荷aod数据，0.2为筛选样本的阈值。3.根据权利要求2所述的基于非线性地表反射率模型的气溶胶遥感定量反演方法，其特征在于，读取目标卫星-载荷观测逐轨道一级数据每个像元的计数值并换算为大气层顶表观反射率，包括：读取目标卫星-载荷观测逐轨道一级数据每个像元于多个可见光通道、近红外通道、短波近红外通道的计数值，通过定标系数和物理量转换公式，将计数值换算为大气层顶表观反射率；其中，所述可见光通道包含蓝光波段、绿光波段和红光波段，所述近红外通道包含两个近红外波段。4.根据权利要求3所述的基于非线性地表反射率模型的气溶胶遥感定量反演方法，其特征在于，订正后得到气体吸收订正后的表观反射率时，利用含水汽、臭氧、二氧化碳的再分析资料对目标卫星-载荷相关波段的表观反射率进行气体吸收订正，得到气体吸收订正后的表观反射率，其中，所述表观反射率的波段包括可见光通道、近红外通道、短波近红外通道。
5.根据权利要求1所述的基于非线性地表反射率模型的气溶胶遥感定量反演方法，其特征在于，读取目标卫星-载荷观测二级云检测产品，对于每个像元，若云检测产品标记为“冰雪”，过程终止。6.根据权利要求5所述的基于非线性地表反射率模型的气溶胶遥感定量反演方法，其特征在于，读取目标卫星-载荷观测二级云检测产品，若云检测产品标记为“一定是晴空”，对于每个像元，执行下述内陆水体判断条件，若不满足，继续执行下述过程，否则终止；判断是否为内体水体，若满足下述条件，过程终止，否则进行下述过程：其中，ρ
nir1
和ρ
red
分别为近红外通道1和红光表观反射率，t1为判断为内陆水体的阈值；ρ
swir
为短波近红外通道表观反射率，t2为短波近红外通道反射率判断阈值。7.根据权利要求6所述的基于非线性地表反射率模型的气溶胶遥感定量反演方法，其特征在于，对筛选到的目标卫星-载荷的有效像元进行平均，得到时空匹配的{时间、经度、纬度、θ
s
、θ
v
、θ、海陆标识、aod、ρ
blue
、ρ
green
、ρ
red
、ρ
nir1
、ρ
nir2
、ρ
swir
、ndvi、ndvi
swir
、afri}观测记录；其中，θ
s
为太阳天顶角、θ
v
为卫星天顶角、为相对方位角、θ为散射角、ρ
blue
、ρ
green
、ρ
red
分别为目标卫星-载荷的蓝光、绿光和红光通道表观反射率，ρ
nir1
和ρ
nir2
分别为目标卫星-载荷的近红外通道1和近红外通道2表观反射率，ρ
swir
为目标卫星-载荷的短波近红外通道表观反射率，ndvi为近红外通道1表观反射率和红光通道表观反射率计算得到，ndvi
swir
为近红外通道1表观反射率和短波近红外通道表观反射率计算得到，afri为蓝光通道表观反射率和短波近红外通道表观反射率计算得到。8.根据权利要求7所述的基于非线性地表反射率模型的气溶胶遥感定量反演方法，其特征在于，还包括：分别统计θ
s
、θ
v
、θ、ndvi、ndvi
swir
、afri参数的最值、平均值分布特征，选择其中一个参数作为自变量，对每个参数取值设置19个断点分为20类，统计每个参数在每类样本下与与的线性拟合关系和相关性，依次变化自变量，直至将6个参数遍历完毕；其中，线性拟合关系包括斜率和截距；为计算得到蓝光通道地表反射率，为计算得到红光通道地表反射率，为计算得到短波近红外通道地表反射率。9.一种计算机设备，其特征在于，所述计算机设备包括多个计算机设备，各计算机设备包括存储器、处理器以及存储在存储器上并可在处理器上运行的计算机程序，其特征在于，所述多个计算机设备的处理器执行所述计算机程序时共同实现权利要求1至8中任一项所述基于非线性地表反射率模型的气溶胶遥感定量反演方法的步骤。10.一种计算机可读存储介质，其上存储有计算机程序，其特征在于，所述存储介质存储的所述计算机程序被处理器执行时实现权利要求1至8中任一项所述基于非线性地表反射率模型的气溶胶遥感定量反演方法的步骤。

技术总结
本发明涉及物联网技术领域，具体涉及一种基于非线性地表反射率模型的气溶胶遥感定量反演方法，从辅助卫星-载荷AOD数据集中筛选出AOD数值小于0.2的样本，作为目标载荷大气订正的背景值，循环太阳、卫星观测角和方位角变化，借助大气辐射传输模型获得不同通道的地表反射率，通过分析不同季节下通道间地表反射率随观测几何、地表类型、地表植被参数等关键参数的变化规律和显著性，建立适用于每组地表反射率关系的数学模型，考虑了地表反射率随观测几何和地表参数的非线性变化规律，提高了AOD反演精度。同时，对地表反射率随地表参数的变化规律进行了分析，更好地刻画了地表反射率的变化规律，可用于各种卫星遥感数据的气溶胶遥感定量反演。定量反演。定量反演。


技术研发人员：司一丹 陈林 徐娜 张兴赢
受保护的技术使用者：国家卫星气象中心（国家空间天气监测预警中心）
技术研发日：2023.05.19
技术公布日：2023/9/23