CN111413296A - 一种考虑地表非朗伯特性的气溶胶光学厚度遥感反演方法 - Google Patents

一种考虑地表非朗伯特性的气溶胶光学厚度遥感反演方法 Download PDF

Info

Publication number
CN111413296A
CN111413296A CN201910838020.9A CN201910838020A CN111413296A CN 111413296 A CN111413296 A CN 111413296A CN 201910838020 A CN201910838020 A CN 201910838020A CN 111413296 A CN111413296 A CN 111413296A
Authority
CN
China
Prior art keywords
reflectivity
remote sensing
aerosol
kernel
optical thickness
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Pending
Application number
CN201910838020.9A
Other languages
English (en)
Inventor
田信鹏
高志强
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Yantai Institute of Coastal Zone Research of CAS
Original Assignee
Yantai Institute of Coastal Zone Research of CAS
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Yantai Institute of Coastal Zone Research of CAS filed Critical Yantai Institute of Coastal Zone Research of CAS
Priority to CN201910838020.9A priority Critical patent/CN111413296A/zh
Publication of CN111413296A publication Critical patent/CN111413296A/zh
Pending legal-status Critical Current

Links

Images

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01N—INVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N21/00—Investigating or analysing materials by the use of optical means, i.e. using sub-millimetre waves, infrared, visible or ultraviolet light
    • G01N21/17—Systems in which incident light is modified in accordance with the properties of the material investigated
    • G01N21/47—Scattering, i.e. diffuse reflection
    • G01N21/49—Scattering, i.e. diffuse reflection within a body or fluid
    • G01N21/53—Scattering, i.e. diffuse reflection within a body or fluid within a flowing fluid, e.g. smoke
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01B—MEASURING LENGTH, THICKNESS OR SIMILAR LINEAR DIMENSIONS; MEASURING ANGLES; MEASURING AREAS; MEASURING IRREGULARITIES OF SURFACES OR CONTOURS
    • G01B11/00—Measuring arrangements characterised by the use of optical techniques
    • G01B11/02—Measuring arrangements characterised by the use of optical techniques for measuring length, width or thickness
    • G01B11/06—Measuring arrangements characterised by the use of optical techniques for measuring length, width or thickness for measuring thickness ; e.g. of sheet material

Abstract

本发明涉及一种考虑地表非朗伯特性的气溶胶光学厚度遥感反演方法，包括获取遥感数据进行预处理；利用遥感图像地表反射率数据基于先验知识约束完成半经验核驱动模型的每日核系数估算；通过离散余弦变换的惩罚最小二乘估计时间序列平滑滤波算法构建单位时间内核系数地表先验知识数据集；建立气溶胶反演通用查找表；利用半经验核驱动模型结合核系数地表先验知识数据集，确定非朗伯体辐射传输前向模型的地表反射率参数；利用上述结果反演获取气溶胶光学厚度。本发明可获取高精度的空间分辨率为500m的气溶胶光学厚度产品，能提供空间相对更为连续的气溶胶分布，显示了在支持区域污染精细管控、重点城市污染传输通道监测及污染物溯源等领域的应用潜力。

Description

一种考虑地表非朗伯特性的气溶胶光学厚度遥感反演方法
技术领域
本发明涉及大气遥感领域，尤指一种考虑地表非朗伯特性的气溶胶光学厚 度遥感反演方法。
背景技术
气溶胶作为地球气候系统的重要组成部分，是大气中三大活跃成分之一， 其数量与性质在时间及空间尺度上的变化均是巨大的，导致难以用固定模式予 以描述。气溶胶光学厚度(Aerosol Optical Depth,AOD)作为表征气溶胶粒子消 光作用的重要光学参数之一，在地气系统辐射收支平衡、全球气候变化及大气 环境和人类健康等方面具有重要指示作用(张莹等，2013；王中挺等，2016)。 由于大气气溶胶具有较大的时空变异性，卫星遥感成为能够提供可靠气溶胶分 布信息的有效手段。卫星传感器接收到的大气层顶辐射信息是由地表及大气共 同作用、相互耦合的最终结果，并且大气贡献相比地表是十分有限的，尤其是 城市等高反射率区域。遥感定量化过程就是从这单一信息中，区分两个未知量， 其本质是一个病态反演问题。在实际的气溶胶遥感中，需要通过先验知识来确 定气溶胶类型及地表反射率参数，这两者的不确定性也正是气溶胶反演误差的 主要来源。气溶胶反演算法的发展与改进，也正是围绕气溶胶类型的假设及地 表反射率估算这两个方面。
近年来，国内外诸多学者在气溶胶反演中对地表方向反射的估算方面做了 许多努力。Thomas等提出了ORAC(Oxford-RAL Aerosol and Cloud)算法实现了 ATSR-2和AATSR数据的气溶胶参数反演，该算法是使用BRDF/Albedo产品构 建反射率先验知识库实现地气分离。Xue等通过中分辨率成像光谱仪(Moderate Resolution ImagingSpectroradiometer,MODIS)BRDF/反照率产品数据来增加方 程约束条件，基于AATSR星下点及前向两个角度观测数据提出了考虑BRDF效 应的气溶胶协同反演算法。在此基础上，Guang等利用MODIS数据提出了AOD 与地表反射率双星协同反演算法，进一步提高了亮地表上空气溶胶反演精度。 Qin等在假设地表BRDF形状稳定的基础上，基于比值法实现了AATSR数据澳 大利亚地区气溶胶及BRDF的反演。上述诸多算法均需要大气条件不变的时间 序列同步观测数据，导致其普适性较差，且使用的仍是朗伯体前向辐射传输模 型计算表观反射率，未从根本上解决地表非朗伯反射的影响。
发明内容
为解决传统气溶胶算法中地表朗伯反射假定对地表反射率估算带来的不确 定性，本算法发明提出一个耦合地表二向反射特性的气溶胶反演算法，最终构 建一个具有较强普适性、可业务化运行的气溶胶反演算法。
为实现上述目的，本发明所采用的技术方案是：一种考虑地表非朗伯特性 的气溶胶光学厚度遥感反演方法，包括：
获取遥感数据并进行预处理；
地表先验知识数据集的构建：利用长时间序列的遥感图像地表反射率数据， 基于先验知识条件约束完成半经验核驱动模型的每日核系数估算；通过离散余 弦变换的惩罚最小二乘估计DCT-PLS时间序列平滑滤波算法构建单位时间内的 核系数地表先验知识数据集；
查找表构建：构建通用大气参数数据集，再根据MODIS光谱响应函数、对 该大气参数数据集中的数据通过积分运算建立MODIS传感器气溶胶反演查找 表文件；
地表反射率参数的确定：利用RTLS半经验核驱动模型，结合构建的核系数 地表先验知识数据集SANIFD，确定非朗伯体辐射传输前向模型中四个地表反射 率参数；
气溶胶光学厚度反演：根据成像角度信息、地表反射率参数及查找表中大 气参数信息，结合预处理获取的表观反射率数据，基于非朗伯体辐射传输前项 模型，反演获取气溶胶光学厚度。
所述基于先验知识条件约束完成半经验核驱动模型核系数的估算，采用的 代价函数为：J(f)＝(Kf-ρ)TΣT(Kf-ρ)+(f-f0)TMT(f-f0)；
其中，K为三个核矩阵，f为对应的三个核系数向量，ρ为方向反射率向量， Σ为反射率误差的协方差矩阵，f0及M为核系数先验知识的均值向量及协方差 矩阵；
上述代价函数的最小二乘解表达式可写成：
f＝(KTΣ-1K+M-1)-1(KTΣ-1ρ+M-1f0)，通过该式即可完成核系数的估算。
所述半经验核驱动模型的核系数包括：各向同性散射fiso、、体散射fvol及几何 光学散射fgeo三个核系数。
所述单位时间内为每若干天内。
通过DCT-PLS时间序列平滑滤波算法构建核系数地表先验知识数据集，具 体实现公式为：
Figure BDA0002192815110000031
式中，|| ||表示欧几里得范数；
Figure BDA0002192815110000032
表示待平滑的时间序列变量；W表示最佳 平滑滤波权重；s为控制滤波平滑度参数，随着s值的增大，
Figure BDA0002192815110000033
会越平滑；对于 一维数据的计算，惩罚项P通常可用Laplace差分计算。
所述通用大气参数数据集记录如下数据：
不同观测几何数据：18个太阳天顶角、15个观测天顶角、19个相对方位角；
气溶胶光学厚度AOD条件(0.05-3.0范围内)下16个气溶胶光学厚度；
0.3-2.5μm光谱范围的大气参数信息，包括tdd(θs)、tdd(θv)、tdh(θs)、thd(θv)、
Figure BDA0002192815110000034
及ρ0。
所述非朗伯体辐射传输前向模型中四个地表反射率参数包括：rdd、rdh、rhd、 rhh，计算公式如下：
Figure BDA0002192815110000035
其中，Kvol及Kgeo分别为体散射核及几何光学核，fiso、fvol及fgeo分别表示各 向同性散射、体散射及几何光学散射核核系数值，与波段λ有关；R为地表反射 率值，当分别计算时，rdd、rdh、rhd、rhh表示地表不同的反射特性，分别为地表 方向-方向反射率、方向-半球反射率、半球-方向反射率及半球-半球反射率。
非朗伯体辐射传输前项模型为：
Figure BDA0002192815110000041
其中，
Figure BDA0002192815110000042
Figure BDA0002192815110000043
矩阵T(θs)、T(θv)及
Figure BDA0002192815110000044
分别定义为：
T(θs)＝[tdd(θs)tdh(θs)]
Figure BDA0002192815110000045
Figure BDA0002192815110000046
式中，θs、θv、
Figure BDA0002192815110000047
分别为太阳天顶角、传感器天顶角、太阳方位角、 传感器方位角及相对方位角；ρTOA为表观反射率；ρ0为大气路径散射反射率，表 示大气程辐射项；
Figure BDA0002192815110000048
为大气球形反照率；t和r分别表示透过率和反射率，下标 h和d分别表示半球(漫射)和方向(直射)；因此，thd(θv)表征在卫星观测方向 上大气对地表漫射辐射的消光，即大气上行漫射透过率；tdd(θv)表征在卫星观测 方向上大气对地表反射直射辐射的消光，即大气上行直射透过率；同样，tdh(θs) 表征在入射方向上大气对天空漫射辐射的消光，即大气下行漫射透过率；tdd(θs) 则为入射方向上大气对太阳直射辐射的消光，即大气下行直射透过率；rdd、rdh、 rhd、rhh表示地表不同的反射特性，分别为地表方向-方向反射率、方向-半球反射率、半球-方向反射率及半球-半球反射率。
对遥感数据进行预处理包括云检测、冰/雪以及水体识别，采用改进的动态 阈值法，具体表示为：
1)首先使用地表先验知识数据集，计算获取各波段地表反射率值：
Figure BDA0002192815110000051
Figure BDA0002192815110000052
然后，利用不同波段云检测动态阈值模型实现云像元的识别，公式表示为：
Figure BDA0002192815110000053
Figure BDA0002192815110000054
Figure BDA0002192815110000055
Figure BDA0002192815110000056
2)由于采用的半经验核驱动模型中几何光学核及体散射核不适用于冰/雪及 水体地表类型的反射率计算，在气溶胶反演中需要将冰/雪及水体进行掩膜剔除， 具体如下：
使用遥感图像预处理后的第2、5波段表观反射率值
Figure RE-GDA0002232421210000057
按照公式
Figure RE-GDA0002232421210000058
计算，识别冰/雪；使用遥感图像预处理后的第1、2波段 的地表反射率值
Figure RE-GDA0002232421210000059
按照公式
Figure RE-GDA00022324212100000510
计算，识别水体。
本发明的有益效果是：
本发明可以获取高精度的空间分辨率为500m的气溶胶光学厚度产品，且 能提供空间相对更为连续的气溶胶分布信息，显示了在支持区域污染精细管控、 重点城市污染传输通道监测及污染物溯源等领域的应用潜力。为了实现上述以 及相关目的，本发明的一个或多个方面同时包括后面的详细说明及权利要求书 中特别指出的特征。
附图说明
图1为本发明的气溶胶光学厚度反演方法流程示意图。
图2为本发明地表先验知识数据集构建流程示意图。
图3为本发明气溶胶产品各数据层示例图。
图4为本发明获取的气溶胶光学厚度与MODIS气溶胶产品空间分布对比图。
图5为本发明的方法AOD反演结果与MODIS各气溶胶产品(a：DT算法，b：DB算法，c：DTDB算法，d：MAIAC算法)在香河植被AERONET站点的 交叉对比验证图。
图6为本发明的方法AOD反演结果与MODIS各气溶胶产品(a：DT算法， b：DB算法，c：DTDB算法，d：MAIAC算法)在北京城市AERONET站点的 交叉对比验证图。
具体实施方式
下面将结合附图对本发明进行更全面的说明，将阐述本发明实现过程中涉 及到的重要具体细节。然而，本发明可以体现为多种不同形式，提供这些实施 例，是为使本发明更加全面和完整，并将本发明的范围完全地传达给本领域的 普通技术人员。
本发面采用一种考虑地表各向异性反射特性的快速、准确的辐射传输模型， 该模型将辐射过程分为上行和下行分量，并且每个分量中包含直射和漫射两部 分。这种划分可以很好的将入射到地表及地表反射辐射的各个反射率分量清楚 的区分确定，其表观反射率计算公式如下：
Figure BDA0002192815110000061
其中，
Figure BDA0002192815110000062
Figure BDA0002192815110000063
矩阵T(θs)、T(θv)及
Figure BDA0002192815110000064
分别定义为：
T(θs)＝[tdd(θs) tdh(θs)]
Figure BDA0002192815110000065
Figure BDA0002192815110000066
式中，θs、θv、
Figure BDA0002192815110000071
分别为太阳天顶角、传感器天顶角、太阳方位角、 传感器方位角及相对方位角；ρTOA为表观反射率；ρ0为大气路径散射反射率，表 示大气程辐射项；
Figure BDA0002192815110000072
为大气球形反照率；t和r分别表示透过率和反射率，下标 h和d分别表示半球(漫射)和方向(直射)；因此，thd(θv)表征在卫星观测方向 上大气对地表漫射辐射的消光，即大气上行漫射透过率；tdd(θv)表征在卫星观测 方向上大气对地表反射直射辐射的消光，即大气上行直射透过率；同样，tdh(θs) 表征在入射方向上大气对天空漫射辐射的消光，即大气下行漫射透过率；tdd(θs) 则为入射方向上大气对太阳直射辐射的消光，即大气下行直射透过率；rdd、rdh、 rhd、rhh表示地表不同的反射特性，分别为地表方向-方向反射率、方向-半球反射率、半球-方向反射率及半球-半球反射率。
由式(1)可知，在计算表观反射率时需要估算与地表有关的反射率参数及 与大气有关的透过率、程辐射等大气参数。由于地表参数和大气参数的估算是 相互独立的，对于地表各反射率参数的估算可通过半经验核驱动模型实现，大 气参数则可以根据观测几何条件从辐射传输代码预先构建的LUT中获取，关于 这两个参数的获取方法在下面的反演流程中会详细阐述。
如图1所示，为本发明的气溶胶光学厚度反演方法流程示意图。流程如下：
(1)遥感数据的获取：本发明算法基于MODIS实现了高分辨率的气溶胶 光学厚度产品的反演生产，在算法构建及产品生产过程中需要的遥感数据主要 包括：1)MODISLevel-1空间分辨率为500m的表观反射率数据(MOD02 HKM)， 2)MODIS Level-1几何位置数据(MOD03)，3)MODIS Level-2空间分辨率为 500m的地表反射率产品数据(MOD09GA)及BRDF/反照率产品数据 (MCD43A1)，4)MODIS Level-3空间分辨率为500m的植被指数产品数据 (MOD13A1)，5)MODIS Level-2气溶胶产品数据(MOD04、MCD19A2)。以 上所有数据均可从LAADS DAAC免费下载获取 (https://ladsweb.modaps.eosdis.nasa.gov/)。
(2)数据处理：对MOD02HKM数据进行几何校正、云/冰雪检测等预处 理，获取可进行气溶胶反演区域的表观反射率数据。
(3)地表先验知识数据集的构建：在本发明算法中，地表先验知识数据集 为半经验核驱动模型中各向同性散射、体散射及几何光学散射三个核系数值， 构建过程包括两步：1)利用长时间序列的MOD09GA数据，基于先验知识条件 约束算法，实现半经验核驱动模型(公式(8))中fiso、fvol及fgeo三个核系数估 算(每日产品)；2)利用第一步估算的逐日核系数数据，基于离散余弦变换的 惩罚最小二乘估计(Discrete Cosine Transform-Penalized Least Square regression, DCT-PLS)时间序列平滑滤波算法合成构建最终每8天的fiso、fvol及fgeo核系数 地表先验知识数据集(Surface ANIsotropic FactorDatabase,SANIFD)，构建流程 如图2所示。
(4)查找表构建：基于大气辐射传输程序(MODTRAN 5)构建通用大气 参数数据集，并结合MODIS传感器光谱响应函数，构建其气溶胶反演查找表文 件。
(5)地表反射率的确定：基于第(3)步构建的地表先验知识数据集，结 合半经验核驱动模型来确定非朗伯体辐射传输前向模型中各地表反射率参数。
(6)气溶胶光学厚度反演：根据成像角度信息、地表反射率及查找表中大 气参数信息，结合第(2)步获取的表观反射率数据，反演获取气溶胶光学厚度。
具体反演流程：
1)几何校正及表观反射率计算。首先利用MOD03数据中的经纬度数据集 作为控制点，逐像元插值计算MOD02HKM对应的实际经纬度；然后，针对每 一个扫描带数据进行几何校正；最后，对所有扫描带校正结果进行拼接完成几 何校正，并同时实现数据的辐射定标获取真实表观反射率数据图像。
2)云检测。本发明采用改进的动态阈值云检测算法进行MODIS数据的云 检测，其基本过程为：首先使用地表先验知识数据集，计算获取公式(2)-(4) 中各波段地表反射率值
Figure BDA0002192815110000081
然后，利用不同波段云检测动态阈值 模型模拟出云像元对应的表观发射率
Figure BDA0002192815110000082
最后，利用真实表 观反射率数据
Figure BDA0002192815110000083
与模拟表观反射率对比实现云像元的识别， 公式表示为：
Figure BDA0002192815110000091
Figure BDA0002192815110000092
Figure BDA0002192815110000093
Figure BDA0002192815110000094
3)冰/雪及水体识别。由于采用的半经验核驱动模型中几何光学核及体散射 核不适用于冰/雪及水体地表类型的反射率计算，在气溶胶反演中需要将这类地 表类型掩膜剔除。算法中采用的利用MODIS第2、5波段的表观反射率对冰/雪 识别(公式(6))。类似的，对于水体像元的识别则利用第1、2波段的表观反 射率(公式(7))：
Figure BDA0002192815110000095
Figure BDA0002192815110000096
4)查找表构建。为了实现不同传感器气溶胶反演，利用MODTRAN 5大气 辐射传输程序构建了一个通用的大气参数数据集(General Atmospheric Parameters Dataset,GAtmParas)，该数据集中记录有不同观测几何(18个太阳天 顶角、15个观测天顶角、19个相对方位角)及AOD条件下(0.05-3.0范围内 16个气溶胶光学厚度)0.3-2.5μm光谱范围的大气参数信息，包括tdd(θs)、tdd(θv)、 tdh(θs)、thd(θv)、
Figure BDA0002192815110000097
及ρ0六个大气参数。根据MODIS光谱响应函数，利用该大气 参数数据集，通过积分运算即可完成MODIS传感器的气溶胶反演查找表的构建。
5)地表反射率的确定。在实际计算中，可利用半经验核驱动模型(公式(8)) 计算公式(1)中的rdd、rdh、rhd、rhh四个反射率。
Figure BDA0002192815110000098
其中，R为地表反射率值，Kvol及Kgeo分别为体散射核及几何光学核，fiso、 fvol及fgeo分别表示各向同性散射、体散射及几何光学散射核核系数值，与波段λ 有关。
并且，实验证明，在利用RossThick体散射核及LiSparseR几何光学核积分 计算黑空及白空反照率时用式(9)能够很好的拟合。因此，本发明也采用该方 式确定各地表反射率，多项式系数如表1所示：
Figure BDA0002192815110000101
表1公式(9)中多项式hk(θ)＝g0k+g1kθ2+g2kθ3系数取值
Figure BDA0002192815110000102
6)气溶胶产品的生产。首先，根据遥感数据的成像几何条件，通过步骤(4) 构建的查找表查算出对应的大气参数；然后，利用步骤(5)确定的地表反射率 信息，结合公式(1)即可计算获取各气溶胶光学厚度条件下的表观反射率向量； 最后，利用步骤(1)计算的真实表观反射率观测值，求算出对应的气溶胶光学 厚度。本发明生产的气溶胶光学厚度产品示例如图3所示。
为了评估反演方法的性能，交叉对比分析了本发明方法获取的气溶胶光学 厚度与MODIS各气溶胶产品在北京地区的空间分布情况，如图4所示。可以看 出，各算法获取的AOD在空间分布上具有较好的一致性，均表现为西南城市地 区具有较高的AOD载荷，植被地区AOD值较低。本文提出的新算法在植被低 反射率及城市等高反射率地区均能得到AOD结果，空间连续性上较好，且空间 分辨率达到500m，可以提供比MODIS各气溶胶产品更多的细节信息。
此外，为了更加细致地对本发明反演算法精度进行评价，利用2013-2016 年香河(XHS)及北京(BJS)两个AERONET站点地面观测资料，将本发明新 算法获取的气溶胶结果分别与MODIS C6.1 DT、DB、DTDB算法产品及C6 MAIAC算法产品进行交叉对比分析，如图5和6所示。新算法获取气溶胶与 MAIAC在城市及植被均具有较高的精度，在BJS/XHS站点的在期望误差内的 比重(Percentage Within the Expected error,PWE)分别为74.9/82.1％、73.1/81.2％， 能够有效地改善MODIS DT、DB及DTDB算法对AOD的高估，提高AOD反 演精度，尤其是植被地区表现更为明显，高于期望误差内的比重(Percentage Above the Expectederror,PAE)降低了15％左右。与MAIAC算法相比，新算法 与实测数据具有更高的相关系数及更小的RMSE值。此外，由于新算法获取的 AOD空间分辨率为500m，可以提供比MODIS各产品更加详细的细节信息，这 对城市地区的环境监测等更加有利。
综合来看，利用本方法可以获取城市上空高精度、高空间分辨率的气溶胶 光学厚度信息，显示了在支持区域污染精细管控、重点城市污染传输通道监测 及污染物溯源等领域的应用潜力。
以上所述是本发明的优选实施方式，应当指出，对于本技术领域的普通技 术人员来说，在不脱离本发明所述原理的前提下，还可以做出若干改进和润饰， 这些改进和润饰应视为本发明的保护范围。

Claims (9)

1.一种考虑地表非朗伯特性的气溶胶光学厚度遥感反演方法，其特征在于包括：
获取遥感数据并进行预处理；
地表先验知识数据集的构建：利用长时间序列的遥感图像地表反射率数据，基于先验知识条件约束完成半经验核驱动模型的每日核系数估算；通过离散余弦变换的惩罚最小二乘估计DCT-PLS时间序列平滑滤波算法构建单位时间内的核系数地表先验知识数据集；
查找表构建：构建通用大气参数数据集，再根据MODIS光谱响应函数、对该大气参数数据集中的数据通过积分运算建立MODIS传感器气溶胶反演查找表文件；
地表反射率参数的确定：利用RTLS半经验核驱动模型，结合构建的核系数地表先验知识数据集SANIFD，确定非朗伯体辐射传输前向模型中四个地表反射率参数；
气溶胶光学厚度反演：根据成像角度信息、地表反射率参数及查找表中大气参数信息，结合预处理获取的表观反射率数据，基于非朗伯体辐射传输前项模型，反演获取气溶胶光学厚度。
2.如权利1要求所述的一种考虑地表非朗伯特性的气溶胶光学厚度遥感反演方法，其特征在于，所述基于先验知识条件约束完成半经验核驱动模型核系数的估算，采用的代价函数为：J(f)＝(Kf-ρ)TΣT(Kf-ρ)+(f-f0)TMT(f-f0)；
其中，K为三个核矩阵，f为对应的三个核系数向量，ρ为方向反射率向量，Σ为反射率误差的协方差矩阵，f0及M为核系数先验知识的均值向量及协方差矩阵；
上述代价函数的最小二乘解表达式可写成：f＝(KTΣ-1K+M-1)-1(KTΣ-1ρ+M-1f0)，通过该式即可完成核系数的估算。
3.如权利1要求所述的一种考虑地表非朗伯特性的气溶胶光学厚度遥感反演方法，其特征在于所述半经验核驱动模型的核系数包括：各向同性散射fiso、、体散射fvol及几何光学散射fgeo三个核系数。
4.如权利1要求所述的一种考虑地表非朗伯特性的气溶胶光学厚度遥感反演方法，其特征在于所述单位时间内为每若干天内。
5.如权利1要求所述的一种考虑地表非朗伯特性的气溶胶光学厚度遥感反演方法，其特征在于，通过DCT-PLS时间序列平滑滤波算法构建核系数地表先验知识数据集，具体实现公式为：
Figure FDA0002192815100000021
式中，|| ||表示欧几里得范数；
Figure FDA0002192815100000022
表示待平滑的时间序列变量；W表示最佳平滑滤波权重；s为控制滤波平滑度参数，随着s值的增大，
Figure FDA0002192815100000023
会越平滑；对于一维数据的计算，惩罚项P通常可用Laplace差分计算。
6.如权利1要求所述的一种考虑地表非朗伯特性的气溶胶光学厚度遥感反演方法，其特征在于所述通用大气参数数据集记录如下数据：
不同观测几何数据：18个太阳天顶角、15个观测天顶角、19个相对方位角；
气溶胶光学厚度AOD条件(0.05-3.0范围内)下16个气溶胶光学厚度；
0.3-2.5μm光谱范围的大气参数信息，包括tdd(θs)、tdd(θv)、tdh(θs)、thd(θv)、
Figure FDA0002192815100000024
及ρ0。
7.如权利1要求所述的一种考虑地表非朗伯特性的气溶胶光学厚度遥感反演方法，其特征在于所述非朗伯体辐射传输前向模型中四个地表反射率参数包括：rdd、rdh、rhd、rhh，计算公式如下：
Figure FDA0002192815100000025
其中，Kvol及Kgeo分别为体散射核及几何光学核，fiso、fvol及fgeo分别表示各向同性散射、体散射及几何光学散射核核系数值，与波段λ有关；R为地表反射率值，当分别计算时，rdd、rdh、rhd、rhh表示地表不同的反射特性，分别为地表方向-方向反射率、方向-半球反射率、半球-方向反射率及半球-半球反射率。
8.如权利1要求所述的一种考虑地表非朗伯特性的气溶胶光学厚度遥感反演方法，其特征在于，非朗伯体辐射传输前项模型为：
Figure FDA0002192815100000026
其中，
Figure FDA0002192815100000031
Figure FDA0002192815100000032
矩阵T(θs)、T(θv)及
Figure FDA0002192815100000033
分别定义为：
T(θs)＝[tdd(θs) tdh(θs)]
Figure FDA0002192815100000034
Figure FDA0002192815100000035
式中，θs、θv、
Figure FDA0002192815100000036
分别为太阳天顶角、传感器天顶角、太阳方位角、传感器方位角及相对方位角；ρTOA为表观反射率；ρ0为大气路径散射反射率，表示大气程辐射项；
Figure FDA0002192815100000037
为大气球形反照率；t和r分别表示透过率和反射率，下标h和d分别表示半球(漫射)和方向(直射)；因此，thd(θv)表征在卫星观测方向上大气对地表漫射辐射的消光，即大气上行漫射透过率；tdd(θv)表征在卫星观测方向上大气对地表反射直射辐射的消光，即大气上行直射透过率；同样，tdh(θs)表征在入射方向上大气对天空漫射辐射的消光，即大气下行漫射透过率；tdd(θs)则为入射方向上大气对太阳直射辐射的消光，即大气下行直射透过率；rdd、rdh、rhd、rhh表示地表不同的反射特性，分别为地表方向-方向反射率、方向-半球反射率、半球-方向反射率及半球-半球反射率。
9.如权利1要求所述的一种考虑地表非朗伯特性的气溶胶光学厚度遥感反演方法，其特征在于对遥感数据进行预处理包括云检测、冰/雪以及水体识别，采用改进的动态阈值法，具体表示为：
1)首先使用地表先验知识数据集，计算获取各波段地表反射率值：
Figure FDA0002192815100000038
Figure FDA0002192815100000039
然后，利用不同波段云检测动态阈值模型实现云像元的识别，公式表示为：
Figure FDA00021928151000000310
Figure FDA0002192815100000041
Figure FDA0002192815100000042
Figure FDA0002192815100000043
2)由于采用的半经验核驱动模型中几何光学核及体散射核不适用于冰/雪及水体地表类型的反射率计算，在气溶胶反演中需要将冰/雪及水体进行掩膜剔除，具体如下：
使用遥感图像预处理后的第2、5波段表观反射率值
Figure FDA0002192815100000044
按照公式
Figure FDA0002192815100000045
计算，识别冰/雪；使用遥感图像预处理后的第1、2波段的地表反射率值
Figure FDA0002192815100000046
按照公式
Figure FDA0002192815100000047
计算，识别水体。
CN201910838020.9A 2019-09-05 2019-09-05 一种考虑地表非朗伯特性的气溶胶光学厚度遥感反演方法 Pending CN111413296A (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201910838020.9A CN111413296A (zh) 2019-09-05 2019-09-05 一种考虑地表非朗伯特性的气溶胶光学厚度遥感反演方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201910838020.9A CN111413296A (zh) 2019-09-05 2019-09-05 一种考虑地表非朗伯特性的气溶胶光学厚度遥感反演方法

Publications (1)

Publication Number Publication Date
CN111413296A true CN111413296A (zh) 2020-07-14

Family

ID=71490815

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201910838020.9A Pending CN111413296A (zh) 2019-09-05 2019-09-05 一种考虑地表非朗伯特性的气溶胶光学厚度遥感反演方法

Country Status (1)

Country Link
CN (1) CN111413296A (zh)

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN112964666A (zh) * 2021-02-03 2021-06-15 中国电子科技集团公司第二十八研究所 一种基于地表双向反射的大气二氧化碳含量计算方法
CN113066057A (zh) * 2021-03-17 2021-07-02 云南电网有限责任公司电力科学研究院 一种气溶胶光学厚度监测方法

Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102902883A (zh) * 2012-09-24 2013-01-30 北京师范大学 一种基于多角度测量构建二向性反射分布函数(brdf)原型库的方法
CN103267533A (zh) * 2013-05-08 2013-08-28 中国科学院光电研究院 一种实用型高光谱遥感图像大气自校正方法及装置
CN106407656A (zh) * 2016-08-29 2017-02-15 中国科学院遥感与数字地球研究所 一种基于高分辨率卫星影像数据的气溶胶光学厚度反演方法
US20170294011A1 (en) * 2016-04-08 2017-10-12 University Of Electronic Science And Technology Of China Method for Retrieving Atmospheric Aerosol Based on Statistical Segmentation
CN107782700A (zh) * 2017-09-07 2018-03-09 北京师范大学 一种avhrr地表反射率重建方法、系统与装置

Patent Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102902883A (zh) * 2012-09-24 2013-01-30 北京师范大学 一种基于多角度测量构建二向性反射分布函数(brdf)原型库的方法
CN103267533A (zh) * 2013-05-08 2013-08-28 中国科学院光电研究院 一种实用型高光谱遥感图像大气自校正方法及装置
US20170294011A1 (en) * 2016-04-08 2017-10-12 University Of Electronic Science And Technology Of China Method for Retrieving Atmospheric Aerosol Based on Statistical Segmentation
CN106407656A (zh) * 2016-08-29 2017-02-15 中国科学院遥感与数字地球研究所 一种基于高分辨率卫星影像数据的气溶胶光学厚度反演方法
CN107782700A (zh) * 2017-09-07 2018-03-09 北京师范大学 一种avhrr地表反射率重建方法、系统与装置

Non-Patent Citations (4)

* Cited by examiner, † Cited by third party
Title
WENHAN QIN等: ""A fast, accurate algorithm to account for non-lambertian surface effects on TOA radiance"" *
佘璐: ""基于Himawari-8/AHI数据的气溶胶光学厚度反演和沙尘监测研究"", 《中国博士学位论文全文数据库 基础科学辑》 *
张桂欣等: ""AMSR2缺失数据重建及其土壤湿度反演精度评价"", 《农业工程学报》 *
王志慧等著: "《植被生物量与生产力遥感估算技术与应用》", 31 December 2016, 黄河水利出版社 *

Cited By (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN112964666A (zh) * 2021-02-03 2021-06-15 中国电子科技集团公司第二十八研究所 一种基于地表双向反射的大气二氧化碳含量计算方法
CN113066057A (zh) * 2021-03-17 2021-07-02 云南电网有限责任公司电力科学研究院 一种气溶胶光学厚度监测方法
CN113066057B (zh) * 2021-03-17 2023-04-14 云南电网有限责任公司电力科学研究院 一种气溶胶光学厚度监测方法

Similar Documents

Publication Publication Date Title
CN105678085B (zh) 一种pm2.5浓度的估算方法及系统
Li et al. An evaluation of the use of atmospheric and BRDF correction to standardize Landsat data
CN104613944B (zh) 一种基于gwr和bp神经网络的分布式水深预测方法
CN108303044B (zh) 一种叶面积指数获取方法及系统
CN108759820B (zh) 一种基于仿复眼多通道偏振传感器的太阳矢量解算方法
CN112560570B (zh) 一种协同降尺度与数据融合的高分辨率地表温度估算方法
CN100483147C (zh) 高光谱亚像元目标探测方法及装置
CN110174359A (zh) 一种基于高斯过程回归的航空高光谱影像土壤重金属浓度评估方法
Wang et al. Image misregistration error in change measurements
CN109635249B (zh) 水体浊度反演模型建立方法、检测方法及装置
CN112013822A (zh) 基于改进gwr模型的多光谱遥感水深反演方法
CN111413296A (zh) 一种考虑地表非朗伯特性的气溶胶光学厚度遥感反演方法
CN113639716A (zh) 一种基于深度残差收缩网络的水深遥感反演方法
CN113324915A (zh) 一种支持高分辨率气溶胶光学厚度反演的城市复杂地表反射率估算方法
Jiang et al. Synergistic use of optical and InSAR data for urban impervious surface mapping: a case study in Hong Kong
Dong et al. Assessment of the hotspot effect for the PROSAIL model with POLDER hotspot observations based on the hotspot-enhanced kernel-driven BRDF model
CN105841819A (zh) 一种有云条件下的地表温度的估算方法和装置
CN111339959A (zh) 基于sar和光学影像融合的近海浮筏养殖区提取方法
Douxchamps et al. Stereoscopic and velocimetric reconstructions of the free surface topography of antidune flows
CN102073038A (zh) 基于微小地形的遥感影像的地形校正的方法
CN116385871B (zh) 一种森林生物量遥感估算方法、系统、设备及介质
Tian et al. Improving MODIS aerosol estimates over land with the surface BRDF reflectances using the 3-D discrete cosine transform and RossThick-LiSparse models
CN102798851A (zh) 一种基于几何成像的modis lai产品验证方法
Zhen et al. Spatial Resolution Requirements for the Application of Temperature and Emissivity Separation (TES) Algorithm over Urban Areas
CN112052720B (zh) 一种基于直方图聚类的高时空归一化植被指数ndvi的融合模型

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
WD01 Invention patent application deemed withdrawn after publication
WD01 Invention patent application deemed withdrawn after publication

Application publication date: 20200714