CN111753439A - 一种国产多角度偏振卫星传感器的气溶胶光学厚度反演方法 - Google Patents

一种国产多角度偏振卫星传感器的气溶胶光学厚度反演方法 Download PDF

Info

Publication number
CN111753439A
CN111753439A CN202010655337.1A CN202010655337A CN111753439A CN 111753439 A CN111753439 A CN 111753439A CN 202010655337 A CN202010655337 A CN 202010655337A CN 111753439 A CN111753439 A CN 111753439A
Authority
CN
China
Prior art keywords
aerosol
angle
optical thickness
data
radiation
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.)
Granted
Application number
CN202010655337.1A
Other languages
English (en)
Other versions
CN111753439B (zh
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.)
Aerospace Information Research Institute of CAS
Original Assignee
Aerospace Information Research Institute 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 Aerospace Information Research Institute of CAS filed Critical Aerospace Information Research Institute of CAS
Priority to CN202010655337.1A priority Critical patent/CN111753439B/zh
Publication of CN111753439A publication Critical patent/CN111753439A/zh
Application granted granted Critical
Publication of CN111753439B publication Critical patent/CN111753439B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00Computer-aided design [CAD]
    • G06F30/20Design optimisation, verification or simulation
    • YGENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02ATECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE
    • Y02A90/00Technologies having an indirect contribution to adaptation to climate change
    • Y02A90/10Information and communication technologies [ICT] supporting adaptation to climate change, e.g. for weather forecasting or climate simulation

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Computer Hardware Design (AREA)
  • Evolutionary Computation (AREA)
  • Geometry (AREA)
  • General Engineering & Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • Image Processing (AREA)
  • Investigating Or Analysing Materials By Optical Means (AREA)

Abstract

本发明公开了一种国产多角度偏振卫星传感器的气溶胶光学厚度反演方法,包括:卫星多角度、多波段数据匹配预处理;多角度偏振卫星影像数据的云像元识别;基于经验正交函数的地表反射贡献估算,构建多角度观测下的散射矩阵进行地表辐射项估算;气溶胶模型构建及反演查找表建立;基于最优模型判识的气溶胶光学厚度反演;本发明首次提出适用于高分五号卫星多角度偏振成像仪(Directional Polarimetric Camera,DPC)的全球陆地上空气溶胶光学厚度反演方法,算法稳定、高效、可验证,并可业务化应用于其他具有多角度探测能力的卫星平台。

Description

一种国产多角度偏振卫星传感器的气溶胶光学厚度反演方法
技术领域
本发明涉及大气科学、遥感科学领域的气溶胶光学厚度反演方法,特别是适用于国产多角度偏振卫星的陆地上空大气气溶胶光学厚度反演方法。
背景技术
大气气溶胶是地球大气中的固体或液体颗粒物质,粒径范围跨越0.001~100 微米5个数量级,形成地球大气中的一个复杂巨系统。气溶胶主要来自于自然和人类活动,自然源气溶胶包括土地扬尘、沙漠地区沙尘传输、火烧爆发、林火燃烧等自然来源排放到大气中的颗粒物,人为源气溶胶则包括工业、农业、交通尾气、建筑施工、城市供暖等活动排放到大气中的部分,既包括直接排放的颗粒物(例如一次气溶胶),也包括气态排放物质经光化学反应、非均相反应等产生的二次气溶胶。气溶胶在气候变化、大气环境、空气质量、人体健康、对地观测、交通运输、生态等众多领域均有重要影响。例如在气候方面,大气气溶胶对气候变化影响复杂且巨大,对太阳辐射和地表辐射同时存在散射和吸收作用,并且由于其构成复杂、时空变化显著,很难定量评估,气溶胶辐射强迫的不确定性超过温室气体等其他因素,是气候变化评估最大的误差来源之一。在大气环境污染方面,气溶胶除了直接降低大气能见度,对交通运输、出行生活产生影响,更重要的是对人体健康有重大危害,不同粒径气溶胶接触人体后附着在呼吸系统的不同部位,这些颗粒物中可能含有酸雾、重金属离子、细菌病毒等,增加了人类罹患呼吸道和心血管疾病的可能性。因此,亟需发展可靠稳定的大气气溶胶空间观测能力,特别是气溶胶的定量化探测方法,是国家和社会的重大需求。
传统上,一般采用地面监测台站对近地面大气进行采样、在位监测,并通过理化、光学等手段分析气溶胶的特性,但这种方法存在一定缺陷:首先是空间覆盖性较差,监测站点主要分布在城市地区,而在广袤的乡村、森林、沙漠、海洋等区域,则很少有分布甚至没有站点,无法形成全球或区域性空间连续的气溶胶监测能力;其次是传统方法一般仅能获取近地面大气的气溶胶样本,不具有中、高层大气的探测能力,更无法获得整层大气的气溶胶信息。遥感是通过非接触方式远距离对目标进行探测的手段,主要利用传感器探测目标物体和电磁波相互作用的辐射特性,从而实现对目标特性的定量或定性探测。利用遥感方法进行大气气溶胶探测,特别是搭载于卫星或航天飞机的遥感传感器,能够有效克服传统采样方法的弊端,并且保持真实大气中气溶胶的自然状态,因此在近几十年得到了快速发展。
气溶胶卫星遥感利用搭载于卫星(极轨卫星、静止卫星等)上的多种先进传感器(主动/被动、光学/微波等)测量和解析多个维度下的气溶胶辐射信号数据(多波段、多角度、辐射、偏振等),实时获取自然状态下整层大气气溶胶的光学特性,实现气溶胶的高精度定量探测。气溶胶卫星传感器历经了从单波段发展到多波段、从单角度发展到多角度、从仅有强度辐射观测到强度和偏振结合观测等阶段,已形成综合利用多波段、多角度、强度、偏振等手段的大气气溶胶探测方案,被认为是最全面的卫星探测手段。2004年12月发射的法国 PARASOL卫星上搭载的POLDER-3就是具备上述探测能力的代表性传感器, POLDER-3在轨稳定运行了约9年的时间,于2013年10月结束探测任务。众多学者基于POLDER-3数据,发展了一系列卫星反演方法进行气溶胶定量化探测,但主要利用的是对细粒子气溶胶敏感的偏振数据,因此POLDER官方仅有细粒子气溶胶光学厚度产品,而缺乏气溶胶总光学厚度产品。2018年5月,我国高分系列中的高分五号卫星在太原卫星发射中心成功发射,其上搭载的多角度偏振成像仪(Directional Polarimetric Camera,DPC)与POLDER-3类似,是具有多波段(可见光到近红外的8个波段)、多角度(目标像元9-12个探测角度)、强度、偏振(STOKES参量I、Q、U)综合探测能力的先进传感器。目前还未有针对DPC/GF-5、能够业务化运行的气溶胶光学厚度反演方法。如何基于我国自主研发的多角度偏振卫星传感器,发展一套稳定可靠的气溶胶光学厚度反演方法,从而实现陆地上空气溶胶含量的定量估计,是提高气候变化评估精度、大气环境污染监测能力等方面所亟需解决的问题。
发明内容
针对现有技术存在的问题,本发明结合卫星遥感的多维度高精度探测能力,利用多波段、多角度、辐射和偏振卫星观测,发展一种适用于国产多角度偏振卫星传感器DPC/GF-5的大气气溶胶光学厚度反演方法,实现陆地上空气溶胶含量的定量获取。
为实现上述目的,本发明的技术方案如下:
一种国产多角度偏振卫星传感器的气溶胶光学厚度反演方法,包括:
步骤一:卫星多角度、多波段数据匹配预处理:首先,根据卫星影像成像时的波段先后顺序,利用观测几何数据判断各角度下多波段影像是否匹配,并通过建立边缘识别算子对观测几何发生跳变的区域进行提取,作为待重组数据;其次,对于待重组数据中的各个像元,根据各观测角度下观测天顶角和方位角随波段和成像序列的变化,进行角度次序定位;最后,将所有待重组数据重新排序,并与其他无需重组的像元整合形成完整的1级辐射和偏振数据,实现多角度、多波段数据匹配,并用于后续的云像元识别、气溶胶反演等步骤;
步骤二:多角度偏振卫星影像数据的云像元识别:首先,提取影像中各像元数据,剔除无效数据并进行多角度、多光谱匹配;其次,通过可见光-近红外通道云敏感波段的表观反射率阈值,逐像元进行光谱维云判识,包括对不同地表类型的像元进行判断;在逐像元识别的基础上,建立影像N×N空间窗口,通过计算所述影像窗口内N×N个像元表观反射率的均一性计算,进行逐窗口空间维云判识;最后,将上述步骤应用在所有有效观测角度上,实现多角度偏振卫星影像数据的云像元识别;
步骤三:基于经验正交函数的地表反射贡献估算:
1)构建多角度观测下的散射矩阵
地-气辐射与卫星传感器入瞳辐射量的交互关系可简化表示为:
Figure BDA0002576565870000031
其中,式中Rtoa表示目标像元的大气层顶表观反射率,Ratm表示目标像元的大气程辐射项,Rsurf表示地表辐射项,μs、μv、φ分别为所述目标像元的太阳天顶角、观测天顶角、相对方位角,λ表示波长;相对于变化情况较复杂的地表像元,大气在一定空间范围内较为稳定,假设在影像N×N窗口对应的范围内大气性质不变,则有:
Figure BDA0002576565870000032
式中,<Rtoa>表示大气层顶表观反射率的偏置量,<Rsurf>为对应的地表辐射项,本发明使用所选窗口内表观反射率的均值作为偏置量,式(2)可进一步变形为:
Figure 1
式中,(x,y)表示窗口内的像元位置下标,Jλ,(x,y)表示像元在波段λ对应的辐射减少量,利用多个观测角度下的辐射观测,可构建散射矩阵:
Figure 2
式中,i、j为角度标号(i,j=1,2,…M),与各角度下的观测几何(μs、μv、φ)相关,Sλ,x,y为像元(x,y)在波长λ处的散射矩阵,可表示为窗口内各像元的地表辐射项及其窗口均值差值的乘积求和形式;
2)地表辐射项估算
通过数学运算对式(4)的散射矩阵求解其特征向量和对应的特征值,实现地表辐射项的估算,假设使用散射矩阵的n个彼此正交的特征向量表征地表辐射项,则有:
Figure BDA0002576565870000041
式中,fk,λ为散射矩阵的特征向量,Ak,λ为对应的特征值,i为角度标识,即特征向量中的元素序号,对应该角度下的观测几何,式(5)经过数学推导可得:
Figure BDA0002576565870000042
式中,Ratm(i,model,τ)表示在气溶胶模型model、气溶胶光学厚度τ下,第i个角度对应的大气程辐射,最终,窗口内地表反射贡献可由式(5)-(6) 联立估算获得;
步骤四:气溶胶模型构建及反演查找表建立:通过6SV辐射传输模式建立查找表,输入参量包括观测几何参数、大气模式、气溶胶模型、气溶胶光学厚度、传感器观测高度、传感器光谱响应函数、地表类型、偏振参数等;其中,构建气溶胶模型的主要参数包括逐波段复折射指数和气溶胶粒子数谱分布,首先计算各气溶胶模型对应的Angstrom波长指数,将6SV默认的550nm波段气溶胶光学厚度,转换到反演所用的各个波段上,与辐射传输前向计算的大气程辐射波段一致,另外,分别输入发射前实验室测量的DPC/GF-5各波段光谱响应函数构建多波段反演查找表;
步骤五:基于最优模型判识的气溶胶光学厚度反演:气溶胶光学厚度反演依据的是构建代价函数,利用步骤(3)建立的查找表,通过正向模拟迭代计算不同气溶胶模型、光学厚度对应的多波段、多角度大气层顶表观反射率,并与实际观测值进行比较得到最优解,本发明所采用的反演代价函数基于最小二乘法建立:
Figure BDA0002576565870000051
式中,Ne为估算地表反射贡献采用的特征向量个数,vλ(j)为权重因子,当第j个角度下的表观反射率数据有效时为1,否则为0,方差项σh,j 2与数据不确定性和噪声相关,代价函数χ2计算所有波段、所有角度的表观反射率总体模拟残差,当代价函数取得最小值时,对应的气溶胶光学厚度即为所求,关于最优气溶胶模型的选择,由步骤(3)可知,在目标像元确定的观测几何下,查找表中多个气溶胶模型对应不同的气溶胶光学厚度,通过代价函数最小化计算出的程辐射项Ratm(j,model,τλ)还存在两个未知量,需要对最优气溶胶模型进行判识;
本发明中采用的最优模型判识及最优光学厚度拟合计算方法如下:首先,对于每个气溶胶模型,根据式(7)计算采用不同特征向量个数表征地表反射贡献时对应的总体模拟残差;然后,根据这些模拟残差建立最优模型判识代价函数:
Figure 100002_3
式中,Nmax为DPC/GF-5数据可使用的特征向量个数最大值;χ2 max为DPC/GF-5数据适用的模型代价函数阈值,满足式(8)的气溶胶模型均为符合反演假设的合理模型;最后,根据模型判识结果对选择出的气溶胶模型对应的光学厚度进行加权平均,得到最终反演结果。
进一步,在步骤五中对所述气溶胶模型对应的光学厚度通过如下公式进行加权平均:
Figure 100002_4
式中,Nselect为根据式(8)选择出的合理模型个数,τNe为各模型中不同特征向量个数对应解算的气溶胶光学厚度,τ为最优拟合的气溶胶光学厚度。
进一步,步骤三中所述<Rtoa>的偏置量采用的是所选窗口内表观反射率的均值作为偏置量。
进一步,步骤三中所述散射矩阵为目标像元地表辐射项的M×M维协方差矩阵,在一定区域内大气均一性假设的条件下,将初始数据集的大气层顶表观反射率Rtoa转换为仅由地表辐射项Rsurf表征的协方差矩阵,通过计算协方差矩阵可以实现地表辐射项的估算。
进一步,步骤三中,所述散射矩阵的正交特征向量的个数小于等于所述散射矩阵维度。
进一步,步骤三中,在式(5)-(6)联立估算前将特征向量按照对应特征值由大到小进行排列,设定阈值滤除较小的特征值及其对应的特征向量。
进一步,步骤三中,分别使用不同特征向量个数表征地表辐射项,其中,可使用的特征向量最大个数通过DPC/GF-5数据的有效观测角度数获得,将计算结果进行加权平均处理估算得到地表反射贡献。
本发明一种国产多角度偏振卫星传感器的气溶胶光学厚度反演方法,采用基于经验正交函数改进方法估算表观反射率中地表反射贡献,解决地气解耦合问题,利用成熟的6SV辐射传输建立气溶胶反演查找表,并针对国产高分五号多角度偏振传感器建立最优气溶胶模型判识方法,最终实现气溶胶光学厚度的高精度反演。
附图说明
图1为实施例中DPC/GF-5获得的东部地区影像;
图2为实施例中DPC/GF-5数据经多角度、多波段匹配后的影像;
图3为实施例中DPC/GF-5数据的云像元识别结果;
图4为6SV辐射传输模式建立气溶胶光学厚度反演查找表方法流程;
图5为实施例中东部地区陆地上空气溶胶光学厚度分布;
图6为实施例中气溶胶光学厚度反演结果与地面验证数据对比。
具体实施方式
为了使本领域的技术人员更好地理解本发明的方案,下面结合本发明示例中的附图对本发明实施例中的技术方案进行清楚、完整地描述,显然,所描述的示例仅仅是本发明的一部分示例,而不是全部的示例。基于本发明的中示例,本领域的普通技术人员在没有做出创造性劳动的前提下,所获得的所有其他实施方式都应当属于本发明保护的范围。
在本实施方式的描述中,术语“内”、“外”、“前”、“后”、“左”、“右”等指示的方位或位置关系均为基于附图所示的方位或位置关系,仅仅是为了便于描述本发明和简化描述,而不是指示或暗示所指的装置或元件必须具有特定的方位、以特定的方位构造和操作,因此不能理解为对本发明的限制。此外,术语“第一”、“第二”仅用于区别类似的对象,而不能理解为特定的顺序或先后次序,应该理解这样的使用在适当情况下可以互换。
为清楚地说明本发明的设计思想,下面结合实施例对本发明进行说明。
实施例
以国产多角度偏振卫星传感器DPC/GF-5实际观测数据为例,描述本发明的具体实施方式。本发明采用的实施例数据是DPC/GF-5观测数据(数据时间:2018 年5月27日,数据空间范围:中国东部地区),数据真彩色影像见图1;根据本发明的方法利用计算机对数据进行处理。
(1)卫星多角度、多波段数据匹配预处理
本发明提出一种基于成像时序排列的多角度、多波段数据匹配方法,适用于国产多角度偏振卫星传感器1级辐射和偏振数据的匹配校正预处理。该方法主要是解决卫星分时成像模式下传感器向前飞行造成的多波段数据难以准确匹配的问题,而现有多角度数据组织形式将这一问题及其带来的后果进一步放大,显著影响大气气溶胶光学厚度等参数的定量反演。本发明提出的这项匹配方法,通过精准定位观测几何/辐射/偏振数据发生跳变的像元,基于逐波段成像时序进行排序和重组,实现多角度、多波段数据的精准匹配。
首先,根据卫星影像成像时的波段先后顺序,利用观测几何数据(观测天顶角和观测方位角)判断各角度下多波段影像是否匹配,并通过建立边缘识别算子对观测几何发生跳变的区域进行提取,作为待重组数据;其次,对于待重组数据中的各个像元,根据各观测角度下观测天顶角和方位角随波段和成像序列的变化,进行角度次序定位;最后,将所有待重组数据重新排序,并与其他无需重组的像元整合形成完整的1级辐射和偏振数据,实现多角度、多波段数据匹配,并用于后续的云像元识别、气溶胶反演等步骤。DPC/GF-5实施例数据多角度、多波段匹配预处理示例结果真彩色影像见图2。
(2)多角度偏振卫星影像数据的云像元识别
本发明采用一种基于光谱维、角度维和空间维联合检测的云像元识别方法,适用于国产多角度偏振卫星传感器影像数据的云判识。与传统卫星遥感影像云识别相比,本发明采用的云像元识别方法的优势主要表现在两个方面:一是多角度数据云识别,由于云的高度有差别,反映在遥感影像上就与观测角度、太阳和传感器的空间位置等因素相关,例如高分辨率影像上可以看到云在地表产生的阴影,多角度观测的优势之一是在短时间内对地表像元进行多次、快速探测,因此可以避免出现在某些特定观测角度下云像元表现为非云像元的特征而被漏识别的情况;二是增加了空间维信息辅助云识别,这是因为大面积厚云通过光谱维阈值方法就可以较好的识别,而薄卷云、高层云以及云边缘等区域像元的光谱特征与晴空像元具有一定的相似性,仅用光谱特征难以有效识别,但这类像元通常具有较大的空间变化特征,因此可利用空间维检测方法进行识别。
首先,提取影像中各像元数据,剔除无效数据(包括填充值和丢失值)并进行多角度、多光谱匹配;其次,通过可见光-近红外通道云敏感波段(443nm、 490nm、670nm、763nm)的表观反射率阈值,逐像元进行光谱维云判识,包括对不同地表类型(植被区域、干旱区域等)的像元进行判断;在逐像元识别的基础上,建立影像N×N空间窗口,通过计算窗口内N×N个像元表观反射率 (443nm、763nm)的均一性(均值和标准差)计算,进行逐窗口空间维云判识;最后,将上述步骤应用在所有有效观测角度上,实现多角度偏振卫星影像数据的云像元识别。DPC/GF-5实施例数据云像元识别示例结果见图3。
(3)基于经验正交函数的地表反射贡献估算
本发明采用的基于经验正交函数的地表反射贡献估算方法源于美国Terra 卫星上的多角度相机MISR。基于经验正交函数的地表反射贡献估算方法的核心思想是:基于一定空间范围内大气具有较好的均一性这一基本假设,将原始数据集的表观反射率转换为地表辐射项,构建基于地表辐射项的协方差矩阵,通过计算经验正交函数(即特征向量)表征该协方差矩阵,从而实现地表辐射项的估计,解决地气解耦合的难题。与传统的单角度探测相比,多角度观测一方面使有效数据量倍增,大幅提高数据等效信噪比,另一方面可以构建散射矩阵,通过散射矩阵的一定数量、彼此正交的特征向量对原始数据集进行重构,在确保数据信息量的同时对数据进行降维,实现快速计算。
由于DPC与MISR的观测机制不同,需要对MISR官方算法进行改进,建立适用于国产多角度偏振传感器的反演算法。与MISR的方法相比,本发明着重在散射矩阵构建、地表辐射项估算等方面进行了改进。
1)构建多角度观测下的散射矩阵
地-气辐射与卫星传感器入瞳辐射量的交互关系可简化表示为:
Figure BDA0002576565870000081
式中Rtoa、Ratm、Rsurf分别表示目标像元的大气层顶表观反射率、大气程辐射项、地表辐射项,μs、μv、φ分别为该像元的太阳天顶角、观测天顶角、相对方位角(由太阳方位角和观测方位角计算),λ表示波长。上述各物理量均随观测角度发生变化。相对于变化情况较复杂的地表像元,大气在一定空间范围内较为稳定,假设在影像N×N窗口对应的范围内大气性质不变,则有:
Figure BDA0002576565870000091
式中,<Rtoa>表示大气层顶表观反射率的偏置量,<Rsurf>为对应的地表辐射项。MISR官方算法中,<Rtoa>偏置量采用的是所选窗口内对应于星下点观测角度的绿光波段最小值的像元表观反射率。而对于DPC/GF-5数据,由于多角度数据匹配形式有所不同,并且影像边缘部分可能出现某些波段或某些观测角度缺失的情况,本发明使用所选窗口内各像元在该波段、该观测角度下的表观反射率均值作为偏置量,则式(2)表达的含义为基于窗口内大气程辐射一致的假设,将各像元的辐射传输转换为整个窗口的辐射传输,物理意义更加明确,并且后续可通过经验正交函数计算得到窗口平均的地表辐射项。式(2)可进一步变形为:
Figure 100002_6
式中,(x,y)表示窗口内的像元位置下标,Jλ,(x,y)表示像元在波段λ对应的辐射减少量。利用多个观测角度下的辐射观测,可构建散射矩阵:
Figure 100002_7
式中,i、j为角度标号(i,j=1,2,…M),与各角度下的观测几何(μs、μv、φ)相关,Sλ,x,y为像元(x,y)在波长λ处的散射矩阵,可表示为窗口内各像元的地表辐射项及其窗口均值差值的乘积求和形式,该散射矩阵实际上是目标像元地表辐射项的M×M维协方差矩阵。上述变换过程通过一定区域内大气均一性假设,将初始数据集的大气层顶表观反射率Rtoa转换为仅由地表辐射项Rsurf表征的协方差矩阵,通过计算协方差矩阵可以实现地表辐射项的估算,即地表反射对大气层顶表观辐射的贡献。
2)地表辐射项估算
通过数学运算对式(4)的散射矩阵求解其特征向量和对应的特征值。假设使用散射矩阵的n个彼此正交的特征向量表征地表辐射项,则有:
Figure BDA0002576565870000094
式中,fk,λ为散射矩阵的特征向量,Ak,λ为对应的特征值,i为角度标识(即特征向量中的元素序号),对应该角度下的观测几何。由于散射矩阵间元素的相关性,正交特征向量的个数一般小于等于散射矩阵维度(即n<=M),因此采用特征值和特征向量表征散射矩阵可以降低计算复杂度,同时最大程度保留有效数据信息。由于正交向量相乘内积为0,式(5)经过数学推导可得:
Figure BDA0002576565870000101
式中,Ratm(i,model,τ)表示在气溶胶模型model、气溶胶光学厚度τ下,第i个角度对应的大气程辐射。最终,窗口内地表反射贡献可由式(5)-(6) 联立估算获得。需要说明,通常情况下较大的特征值对应有效信息,而较小的特征值反映的是噪声信息,因此计算前需首先将特征向量按照对应特征值由大到小进行排列,设定阈值滤除掉一部分特征向量。另外为降低计算不确定性,本发明采用的方案是分别使用不同特征向量个数表征地表辐射项(可使用的特征向量最大个数通过DPC/GF-5数据的有效观测角度数获得),并将计算结果进行加权平均得到地表反射贡献。该步骤可以有效减少特征向量个数取值带来的影响。
(4)气溶胶模型构建及反演查找表建立
气溶胶光学厚度的查找表反演方法的核心思路是基于气溶胶模型,建立各观测几何下不同气溶胶光学厚度与多波段大气程辐射的一一对应关系,即反演查找表,之后结合地表反射估算得到大气层顶表观反射率模拟值,通过设置代价函数,对表观反射率的模拟值和实际观测值进行迭代比较,得到最接近真实状态的大气辐射参数,即查找出待求的气溶胶光学厚度。因此,反演查找表的建立在反演中至关重要,需具备全面性和细致性(理论上所有可能的观测几何、气溶胶模型、光学厚度组合)。
式(6)中不同观测几何下气溶胶模型、气溶胶光学厚度等信息,需要利用本步骤的气溶胶模型构建和查找表建立获得(图4)。查找表采用发展较成熟的 6SV辐射传输模式(6SV-Version 2.1, http://6s.ltdri.org/pages/downloads.html)建立,输入参量包括观测几何参数、大气模式、气溶胶模型、气溶胶光学厚度、传感器观测高度、传感器光谱响应函数、地表类型、偏振参数等。其中,气溶胶模型的构建是关键过程,本发明采用的气溶胶模型来自于代表性文献(基于全球范围、长时间地基观测数据获得),具有较好的典型性和代表性,构建气溶胶模型的主要参数包括逐波段复折射指数(实部、虚部)、气溶胶粒子数谱分布等。另一方面,考虑到反演中各波段的辐射传输计算,需要首先计算各气溶胶模型对应的Angstrom波长指数,将6SV默认的550nm波段气溶胶光学厚度,转换到反演所用的各个波段上,与辐射传输前向计算的大气程辐射波段一致。基于上述设置,通过6SV执行运算后输出的主要有地表辐射参数与大气辐射参数,大气参数包括水汽、臭氧、碳氧化物、氮氧化物、甲烷等吸收性气体,分子瑞利散射,以及气溶胶等主要大气成分的上行、下行、总透过率参数,同时还包括瑞利散射、气溶胶散射辐射参数,如半球反照率、程辐射、偏振反射率、偏振度、相函数、单次散射反照率等。根据上述输入和输出参数,建立逐模型、逐观测几何、逐气溶胶光学厚度对应的查找表,用于后续反演。
(4)基于最优模型判识的气溶胶光学厚度反演
气溶胶光学厚度反演依据的是构建代价函数,利用步骤(3)建立的查找表,通过正向模拟迭代计算不同气溶胶模型、光学厚度对应的多波段、多角度大气层顶表观反射率,并与实际观测值进行比较已得到最优解。本发明所采用的反演代价函数基于最小二乘法建立:
Figure BDA0002576565870000111
式中,Ne为估算地表反射贡献采用的特征向量个数,vλ(j)为权重因子,当第j个角度下的表观反射率数据有效时为1,否则为0,方差项σh,j 2与数据不确定性和噪声相关。代价函数χ2计算所有波段、所有角度的表观反射率总体模拟残差,当代价函数取得最小值时,对应的气溶胶光学厚度即为所求。除地表估算外,气溶胶光学厚度反演的另一个关键问题是最优气溶胶模型的选择,由步骤(3)可知,在目标像元确定的观测几何下,查找表中多个气溶胶模型对应不同的气溶胶光学厚度,即通过代价函数最小化计算出的程辐射项Ratm(j,model, τλ)还存在两个未知量,需要进一步对最优气溶胶模型进行判识。本发明采用的最优模型判识及最优光学厚度拟合计算方法如下:首先,对于每个气溶胶模型,根据式(7)计算采用不同特征向量个数表征地表反射贡献时对应的总体模拟残差;然后,根据这些模拟残差建立最优模型判识代价函数:
Figure 100002_8
式中,Nmax为DPC/GF-5数据可使用的特征向量个数最大值;χ2 max为DPC/GF-5 数据适用的模型代价函数阈值,满足式(8)的气溶胶模型均为符合反演假设的合理模型;最后,根据模型判识结果对上述各选择模型对应的气溶胶光学厚度进行加权平均,得到最终反演结果:
Figure 100002_9
式中,Nselect为根据式(8)选择出的合理模型个数,τNe为各模型中不同特征向量个数对应解算的气溶胶光学厚度,τ为最优拟合的气溶胶光学厚度。 DPC/GF-5实施例数据气溶胶光学厚度反演结果见图5。
实施例数据气溶胶光学厚度反演结果,采用全球地基气溶胶自动观测网(AERONET,https://aeronet.gsfc.nasa.gov)下载的level 1.5级别的气溶胶光学厚度数据(与卫星观测数据同步)进行验证。从验证结果(图6)可以看出,本发明实施例中,DPC/GF-5卫星数据反演的气溶胶光学厚度,与AERONET 数据有较好的一致性和相关性,说明了本发明所提出的适用于国产多角度偏振传感器的陆地气溶胶光学厚度反演方法可靠、稳定。
需要说明的是,除了上述给出的具体示例之外,其中的一些部分可有不同选择。而这些都是本领域技术人员在理解本发明思想的基础上基于其基本技能即可做出的,故在此不再一一例举。
最后,可以理解的是,以上实施方式仅仅是为了说明本发明的原理而采用的示例性实施方式,然而本发明并不局限于此。对于本领域普通技术人员而言,在不脱离本发明的原理和实质的情况下,可以做出各种变型和改进,这些变型和改进也视为本发明的保护范围。

Claims (7)

1.一种国产多角度偏振卫星传感器的气溶胶光学厚度反演方法,其特征在于,包括:
步骤一:卫星多角度、多波段数据匹配预处理:首先,根据卫星影像成像时的波段先后顺序,利用观测几何数据判断各角度下多波段影像是否匹配,并通过建立边缘识别算子对观测几何发生跳变的区域进行提取,作为待重组数据;其次,对于待重组数据中的各个像元,根据各观测角度下观测天顶角和方位角随波段和成像序列的变化,进行角度次序定位;最后,将所有待重组数据重新排序,并与其他无需重组的像元整合形成完整的1级辐射和偏振数据,实现多角度、多波段数据匹配,并用于后续的云像元识别、气溶胶反演等步骤;
步骤二:多角度偏振卫星影像数据的云像元识别:首先,提取影像中各像元数据,剔除无效数据并进行多角度、多光谱匹配;其次,通过可见光-近红外通道云敏感波段的表观反射率阈值,逐像元进行光谱维云判识,包括对不同地表类型的像元进行判断;在逐像元识别的基础上,建立影像N×N空间窗口,通过计算所述影像窗口内N×N个像元表观反射率的均一性计算,进行逐窗口空间维云判识;最后,将上述步骤应用在所有有效观测角度上,实现多角度偏振卫星影像数据的云像元识别;
步骤三:基于经验正交函数的地表反射贡献估算:
1)构建多角度观测下的散射矩阵
地-气辐射与卫星传感器入瞳辐射量的交互关系为:
Figure FDA0002576565860000011
其中,式中Rtoa表示目标像元的大气层顶表观反射率,Ratm表示目标像元的大气程辐射项,Rsurf表示地表辐射项,μs、μv、φ分别为所述目标像元的太阳天顶角、观测天顶角、相对方位角,λ表示波长;相对于变化情况较复杂的地表像元,大气在一定空间范围内较为稳定,假设在影像N×N窗口对应的范围内大气性质不变,则有:
Figure FDA0002576565860000012
式中,<Rtoa>表示大气层顶表观反射率的偏置量,<Rsurf>为对应的地表辐射项,本发明使用所选窗口内表观反射率的均值作为偏置量,式(2)可进一步变形为:
Figure 3
式中,(x,y)表示窗口内的像元位置下标,Jλ,(x,y)表示像元在波段λ对应的辐射减少量,利用多个观测角度下的辐射观测,可构建散射矩阵:
Figure 4
式中,i、j为角度标号(i,j=1,2,…M),与各角度下的观测几何(μs、μv、φ)相关,Sλ,x,y为像元(x,y)在波长λ处的散射矩阵,可表示为窗口内各像元的地表辐射项及其窗口均值差值的乘积求和形式;
2)地表辐射项估算
通过数学运算对式(4)的散射矩阵求解其特征向量和对应的特征值,实现地表辐射项的估算,假设使用散射矩阵的n个彼此正交的特征向量表征地表辐射项,则有:
Figure 5
式中,fk,λ为散射矩阵的特征向量,Ak,λ为对应的特征值,i为角度标识,即特征向量中的元素序号,对应该角度下的观测几何,式(5)经过数学推导可得:
Figure FDA0002576565860000023
式中,Ratm(i,model,τ)表示在气溶胶模型model、气溶胶光学厚度τ下,第i个角度对应的大气程辐射,最终,窗口内地表反射贡献可由式(5)-(6)联立估算获得;
步骤四:气溶胶模型构建及反演查找表建立:通过6SV辐射传输模式建立查找表,输入参量包括观测几何参数、大气模式、气溶胶模型、气溶胶光学厚度、传感器观测高度、传感器光谱响应函数、地表类型、偏振参数等;其中,构建气溶胶模型的主要参数包括逐波段复折射指数和气溶胶粒子数谱分布,首先计算各气溶胶模型对应的Angstrom波长指数,将6SV默认的550nm波段气溶胶光学厚度,转换到反演所用的各个波段上,与辐射传输前向计算的大气程辐射波段一致,另外,分别输入发射前实验室测量的DPC/GF-5各波段光谱响应函数构建多波段反演查找表;
步骤五:基于最优模型判识的气溶胶光学厚度反演:气溶胶光学厚度反演依据的是构建代价函数,利用步骤(3)建立的查找表,通过正向模拟迭代计算不同气溶胶模型、光学厚度对应的多波段、多角度大气层顶表观反射率,并与实际观测值进行比较得到最优解,本发明所采用的反演代价函数基于最小二乘法建立:
Figure 6
式中,Ne为估算地表反射贡献采用的特征向量个数,vλ(j)为权重因子,当第j个角度下的表观反射率数据有效时为1,否则为0,方差项σh,j 2与数据不确定性和噪声相关,代价函数χ2计算所有波段、所有角度的表观反射率总体模拟残差,当代价函数取得最小值时,对应的气溶胶光学厚度即为所求,关于最优气溶胶模型的选择,由步骤(3)可知,在目标像元确定的观测几何下,查找表中多个气溶胶模型对应不同的气溶胶光学厚度,通过代价函数最小化计算出的程辐射项Ratm(j,model,τλ)还存在两个未知量,需要对最优气溶胶模型进行判识;
本发明中采用的最优模型判识及最优光学厚度拟合计算方法如下:首先,对于每个气溶胶模型,根据式(7)计算采用不同特征向量个数表征地表反射贡献时对应的总体模拟残差;然后,根据这些模拟残差建立最优模型判识代价函数:
Figure 7
式中,Nmax为DPC/GF-5数据可使用的特征向量个数最大值;χ2 max为DPC/GF-5数据适用的模型代价函数阈值,满足式(8)的气溶胶模型均为符合反演假设的合理模型;最后,根据模型判识结果对选择出的气溶胶模型对应的光学厚度进行加权平均,得到最终反演结果。
进一步,在步骤五中对所述气溶胶模型对应的光学厚度通过如下公式进行加权平均:
Figure 8
式中,Nselect为根据式(8)选择出的合理模型个数,τNe为各模型中不同特征向量个数对应解算的气溶胶光学厚度,τ为最优拟合的气溶胶光学厚度。
进一步,步骤三中所述<Rtoa>的偏置量采用的是所选窗口内表观反射率的均值作为偏置量。
进一步,步骤三中所述散射矩阵为目标像元地表辐射项的M×M维协方差矩阵,在一定区域内大气均一性假设的条件下,将初始数据集的大气层顶表观反射率Rtoa转换为仅由地表辐射项Rsurf表征的协方差矩阵,通过计算协方差矩阵可以实现地表辐射项的估算。
2.根据权利要求1所述的气溶胶光学厚度反演方法,其特征在于,在步骤五中对所述气溶胶模型对应的光学厚度通过如下公式进行加权平均:
Figure 9
式中,Nselect为根据式(8)选择出的合理模型个数,τNe为各模型中不同特征向量个数对应解算的气溶胶光学厚度,τ为最优拟合的气溶胶光学厚度。
3.根据权利要求1所述的气溶胶光学厚度反演方法,其特征在于,步骤三中所述<Rtoa>的偏置量采用的是所选窗口内表观反射率的均值作为偏置量。
4.根据权利要求1所述的气溶胶光学厚度反演方法,其特征在于,步骤三中所述散射矩阵为目标像元地表辐射项的M×M维协方差矩阵,在一定区域内大气均一性假设的条件下,将初始数据集的大气层顶表观反射率Rtoa转换为仅由地表辐射项Rsurf表征的协方差矩阵,通过计算协方差矩阵可以实现地表辐射项的估算。
5.根据权利要求1所述的气溶胶光学厚度反演方法,其特征在于,步骤三中,所述散射矩阵的正交特征向量的个数小于等于所述散射矩阵维度。
6.根据权利要求1所述的气溶胶光学厚度反演方法,其特征在于,步骤三中,在式(5)-(6)联立估算前将特征向量按照对应特征值由大到小进行排列,设定阈值滤除较小的特征值及其对应的特征向量。
7.根据权利要求1所述的气溶胶光学厚度反演方法,其特征在于,步骤三中,分别使用不同特征向量个数表征地表辐射项,其中,可使用的特征向量最大个数通过DPC/GF-5数据的有效观测角度数获得,将计算结果进行加权平均处理估算得到地表反射贡献。
CN202010655337.1A 2020-07-09 2020-07-09 一种多角度偏振卫星传感器的气溶胶光学厚度反演方法 Active CN111753439B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202010655337.1A CN111753439B (zh) 2020-07-09 2020-07-09 一种多角度偏振卫星传感器的气溶胶光学厚度反演方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202010655337.1A CN111753439B (zh) 2020-07-09 2020-07-09 一种多角度偏振卫星传感器的气溶胶光学厚度反演方法

Publications (2)

Publication Number Publication Date
CN111753439A true CN111753439A (zh) 2020-10-09
CN111753439B CN111753439B (zh) 2023-07-21

Family

ID=72710035

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202010655337.1A Active CN111753439B (zh) 2020-07-09 2020-07-09 一种多角度偏振卫星传感器的气溶胶光学厚度反演方法

Country Status (1)

Country Link
CN (1) CN111753439B (zh)

Cited By (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN112945384A (zh) * 2021-03-01 2021-06-11 中国科学院空天信息创新研究院 一种多角度偏振卫星的数据预处理方法
CN113033025A (zh) * 2021-04-25 2021-06-25 中国人民解放军91977部队 一种海洋大气辐射效应模拟系统和方法
CN113740263A (zh) * 2021-03-10 2021-12-03 仲恺农业工程学院 一种气溶胶光学厚度反演方法及大气颗粒物遥感反演方法
CN113836731A (zh) * 2021-09-28 2021-12-24 中国科学院空天信息创新研究院 陆表稳定目标大气层顶反射率模型的构建方法和装置
CN115630256A (zh) * 2022-12-07 2023-01-20 自然资源部第二海洋研究所 基于暗像元假定的多角度偏振水色卫星大气校正方法
CN116148189A (zh) * 2023-04-14 2023-05-23 自然资源部第二海洋研究所 一种基于被动偏振卫星数据的气溶胶层高获取方法
CN116660106A (zh) * 2023-07-21 2023-08-29 中国科学院空天信息创新研究院 协同星载标量和偏振观测数据的气溶胶参数迭代反演方法
CN117607919A (zh) * 2023-11-17 2024-02-27 中国科学院大气物理研究所 一种基于城市建筑物阴影的气溶胶卫星遥感反演方法

Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103927455A (zh) * 2014-04-24 2014-07-16 中国科学院遥感与数字地球研究所 基于高分一号卫星的陆地气溶胶光学性质反演方法
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
CN110411927A (zh) * 2019-08-02 2019-11-05 中国科学院遥感与数字地球研究所 一种大气细粒子aod和地表偏振反射率协同反演方法

Patent Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103927455A (zh) * 2014-04-24 2014-07-16 中国科学院遥感与数字地球研究所 基于高分一号卫星的陆地气溶胶光学性质反演方法
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 中国科学院遥感与数字地球研究所 一种基于高分辨率卫星影像数据的气溶胶光学厚度反演方法
CN110411927A (zh) * 2019-08-02 2019-11-05 中国科学院遥感与数字地球研究所 一种大气细粒子aod和地表偏振反射率协同反演方法

Cited By (11)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN112945384A (zh) * 2021-03-01 2021-06-11 中国科学院空天信息创新研究院 一种多角度偏振卫星的数据预处理方法
CN113740263A (zh) * 2021-03-10 2021-12-03 仲恺农业工程学院 一种气溶胶光学厚度反演方法及大气颗粒物遥感反演方法
CN113740263B (zh) * 2021-03-10 2023-06-16 仲恺农业工程学院 一种气溶胶光学厚度反演方法及大气颗粒物遥感反演方法
CN113033025A (zh) * 2021-04-25 2021-06-25 中国人民解放军91977部队 一种海洋大气辐射效应模拟系统和方法
CN113836731A (zh) * 2021-09-28 2021-12-24 中国科学院空天信息创新研究院 陆表稳定目标大气层顶反射率模型的构建方法和装置
CN113836731B (zh) * 2021-09-28 2023-06-20 中国科学院空天信息创新研究院 陆表稳定目标大气层顶反射率模型的构建方法和装置
CN115630256A (zh) * 2022-12-07 2023-01-20 自然资源部第二海洋研究所 基于暗像元假定的多角度偏振水色卫星大气校正方法
CN116148189A (zh) * 2023-04-14 2023-05-23 自然资源部第二海洋研究所 一种基于被动偏振卫星数据的气溶胶层高获取方法
CN116660106A (zh) * 2023-07-21 2023-08-29 中国科学院空天信息创新研究院 协同星载标量和偏振观测数据的气溶胶参数迭代反演方法
CN116660106B (zh) * 2023-07-21 2023-10-17 中国科学院空天信息创新研究院 协同星载标量和偏振观测数据的气溶胶参数迭代反演方法
CN117607919A (zh) * 2023-11-17 2024-02-27 中国科学院大气物理研究所 一种基于城市建筑物阴影的气溶胶卫星遥感反演方法

Also Published As

Publication number Publication date
CN111753439B (zh) 2023-07-21

Similar Documents

Publication Publication Date Title
CN111753439B (zh) 一种多角度偏振卫星传感器的气溶胶光学厚度反演方法
US11461994B2 (en) Methods for in-scene shadow compensation using sunlit and skylit illumination factors
US9990705B2 (en) Atmospheric compensation in satellite imagery
Kokhanovsky et al. Space-based remote sensing of atmospheric aerosols: The multi-angle spectro-polarimetric frontier
Mei et al. Retrieval of aerosol optical properties using MERIS observations: Algorithm and some first results
Susaki et al. Validation of MODIS albedo products of paddy fields in Japan
Garestier et al. PolInSAR analysis of X-band data over vegetated and urban areas
Maddy et al. Using Metop-A AVHRR clear-sky measurements to cloud-clear Metop-A IASI column radiances
Sourdeval et al. A variational approach for retrieving ice cloud properties from infrared measurements: application in the context of two IIR validation campaigns
Lisenko Atmospheric correction of multispectral satellite images based on the solar radiation transfer approximation model
Callewaert et al. The Mineral Aerosol Profiling from Infrared Radiances (MAPIR) algorithm: version 4.1 description and evaluation
Wang et al. Aerosol retrieval algorithm based on adaptive land–atmospheric decoupling for polarized remote sensing over land surfaces
Cairns et al. Polarimetric remote sensing of aerosols over land surfaces
Perkins et al. High-speed atmospheric correction for spectral image processing
Pandya et al. Development of a scheme for atmospheric correction of Resourcesat-2 AWiFS data
Liu et al. Vicarious radiometric calibration using a ground radiance-based approach: A case study of sentinel 2A MSI
Kerekes et al. Analytical model of hyperspectral system performance
Chen et al. On the negative polarization's effects at top of atmosphere
Chen et al. High-resolution BRDF and albedo parameters inversion from sentinel-2 multispectral instrument data
Vincent et al. Subpixel target detection in hyperspectral imaging
Wang et al. Study on the Radiation Performance of FY-3D Visible Short Wave Infrared Waveband Based on Reflectivity Basis Method
Shi et al. Multiparameter estimation from landsat observations with topographic consideration
Tarasenkov et al. Assessment of the cloud adjacency effect on retrieval of ground surface reflectance from MODIS satellite data for a surface area near Lake Baikal
Yong et al. TECIS: The First Mission Towards Forest Carbon Mapping By Combination Of Lidar And Multi-Angle Optical Observations
Chen Global estimation of precipitation using opaque microwave bands

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
GR01 Patent grant
GR01 Patent grant