CN113324915A - 一种支持高分辨率气溶胶光学厚度反演的城市复杂地表反射率估算方法 - Google Patents

一种支持高分辨率气溶胶光学厚度反演的城市复杂地表反射率估算方法 Download PDF

Info

Publication number
CN113324915A
CN113324915A CN202110562764.XA CN202110562764A CN113324915A CN 113324915 A CN113324915 A CN 113324915A CN 202110562764 A CN202110562764 A CN 202110562764A CN 113324915 A CN113324915 A CN 113324915A
Authority
CN
China
Prior art keywords
reflectivity
wave band
blue
red
earth surface
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
CN202110562764.XA
Other languages
English (en)
Other versions
CN113324915B (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.)
Wuhan University WHU
Original Assignee
Wuhan University WHU
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 Wuhan University WHU filed Critical Wuhan University WHU
Priority to CN202110562764.XA priority Critical patent/CN113324915B/zh
Publication of CN113324915A publication Critical patent/CN113324915A/zh
Application granted granted Critical
Publication of CN113324915B publication Critical patent/CN113324915B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

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
    • 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/25—Colour; Spectral properties, i.e. comparison of effect of material on the light at two or more different wavelengths or wavelength bands
    • G01N21/31—Investigating relative effect of material at wavelengths characteristic of specific elements or molecules, e.g. atomic absorption spectrometry
    • 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/25—Colour; Spectral properties, i.e. comparison of effect of material on the light at two or more different wavelengths or wavelength bands
    • G01N21/31—Investigating relative effect of material at wavelengths characteristic of specific elements or molecules, e.g. atomic absorption spectrometry
    • G01N21/35—Investigating relative effect of material at wavelengths characteristic of specific elements or molecules, e.g. atomic absorption spectrometry using infrared light
    • 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/25—Colour; Spectral properties, i.e. comparison of effect of material on the light at two or more different wavelengths or wavelength bands
    • G01N21/31—Investigating relative effect of material at wavelengths characteristic of specific elements or molecules, e.g. atomic absorption spectrometry
    • G01N21/35—Investigating relative effect of material at wavelengths characteristic of specific elements or molecules, e.g. atomic absorption spectrometry using infrared light
    • G01N21/3504—Investigating relative effect of material at wavelengths characteristic of specific elements or molecules, e.g. atomic absorption spectrometry using infrared light for analysing gases, e.g. multi-gas analysis
    • 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/25—Colour; Spectral properties, i.e. comparison of effect of material on the light at two or more different wavelengths or wavelength bands
    • G01N21/31—Investigating relative effect of material at wavelengths characteristic of specific elements or molecules, e.g. atomic absorption spectrometry
    • G01N21/35—Investigating relative effect of material at wavelengths characteristic of specific elements or molecules, e.g. atomic absorption spectrometry using infrared light
    • G01N21/3563—Investigating relative effect of material at wavelengths characteristic of specific elements or molecules, e.g. atomic absorption spectrometry using infrared light for analysing solids; Preparation of samples therefor
    • 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/25—Colour; Spectral properties, i.e. comparison of effect of material on the light at two or more different wavelengths or wavelength bands
    • G01N21/31—Investigating relative effect of material at wavelengths characteristic of specific elements or molecules, e.g. atomic absorption spectrometry
    • G01N21/35—Investigating relative effect of material at wavelengths characteristic of specific elements or molecules, e.g. atomic absorption spectrometry using infrared light
    • G01N21/359—Investigating relative effect of material at wavelengths characteristic of specific elements or molecules, e.g. atomic absorption spectrometry using infrared light using near infrared light
    • 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/55—Specular reflectivity
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/20—Design optimisation, verification or simulation
    • 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
    • G01N2021/1765—Method using an image detector and processing of image signal
    • 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
    • G01N2021/178—Methods for obtaining spatial resolution of the property being measured
    • 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
    • G01N2021/1793—Remote sensing
    • G01N2021/1795—Atmospheric mapping of gases
    • 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
    • G01N2021/1793—Remote sensing
    • G01N2021/1797—Remote sensing in landscape, e.g. crops

Landscapes

  • Physics & Mathematics (AREA)
  • Spectroscopy & Molecular Physics (AREA)
  • General Physics & Mathematics (AREA)
  • Chemical & Material Sciences (AREA)
  • Pathology (AREA)
  • Analytical Chemistry (AREA)
  • Biochemistry (AREA)
  • General Health & Medical Sciences (AREA)
  • Health & Medical Sciences (AREA)
  • Immunology (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Engineering & Computer Science (AREA)
  • Theoretical Computer Science (AREA)
  • Computer Hardware Design (AREA)
  • Evolutionary Computation (AREA)
  • Geometry (AREA)
  • General Engineering & Computer Science (AREA)
  • Image Processing (AREA)
  • Investigating Or Analysing Materials By Optical Means (AREA)

Abstract

本发明涉及一种支持高分辨率气溶胶光学厚度反演的城市复杂地表反射率估算方法。首先为每个像元构建像元数据库，根据归一化植被指数(NDVI)将城市复杂地表分为三类，然后对每类地表构建精确的地表反射率估算方案。本发明实现的高空间分辨率地表反射率确定方案，可根据城市地区不同地表的反射特性和变化特性，精确的获取不同地表的地表反射率，且在伪不变地表上考虑了地表的双向反射特性。基于该发明获取的地表反射率可用于支撑城市复杂地表的高分辨率卫星AOD反演。此外，该方法或该方法获取的地表反射率可以推广到具有同类型传感器的卫星上，进一步提高高空间分辨率卫星在AOD反演中的应用。

Description

一种支持高分辨率气溶胶光学厚度反演的城市复杂地表反射 率估算方法
技术领域
本发明属于卫星被动遥感领域，特别是涉及一种支持高分辨率气溶胶光学厚度反演的城市复杂地表反射率估算方法。
背景技术
气溶胶是悬浮在大气中直径小于100um的固态或液态颗粒物的总称，已对全球气候变化、城市生态环境和人类健康造成重要影响。现有的粗分辨率(>1km)气溶胶光学厚度(AOD)产品，无法反映大气气溶胶在小空间尺度上的变化，尤其是不能够满足城市地区的精细化和精准化大气污染监测的需求。高空间分辨率卫星观测能提供高分辨率气溶胶产品，但其反演算法需要将地表贡献从卫星观测中剥离，以获取微弱的气溶胶散射信号，因此获取准确的高分辨率地表反射率是高分辨率卫星气溶胶反演的关键问题之一。
地表反射率并不是固定的，它与地表覆盖类型、视图几何(即太阳天顶角、卫星天顶角、相对方位角)等相关。Wei等利用landsat逐场景的地表反射率产品，基于第二最小值的原则构建了逐月的地表反射率先验知识库。Chen等利用一定时段内的sentinel-2的地表反射率库和地基AERONET站点的AOD确定最干净的像元，并用其来代替一个月内的地表反射率。上述方法都假设地表为朗伯面，而城市的复杂地表(如人造地表)是非朗伯面，存在明显的双向反射特性。Lyapustin等基于中分辨率卫星MODIS重访周期短的特点，在短期内(8-16天)获取多个(至少4个)不同角度观测的地表反射率数据，基于半经验RTLS模型建立了地表双向反射率分布函数(BRDF)，例如MODIS的MOD19A3产品。然而，高分辨率卫星影像由于其空间分辨率较高而卫星视场较窄，所以重访周期较长(如Landsat-8卫星的重访周期为16天)，很难在短期内获得足够的有效观测，这为准确估算城市复杂地表的双向反射特性带来了较大的困难。
因此，急需研究一种考虑地表双向反射特性，适用于高分辨率卫星影像的地表反射率估算方法，为复杂城市地表的气溶胶光学厚度反演提供精确的地表反射率，进而提高气溶胶反演的精度。
发明内容
本发明针对现有技术的不足，提供一种支持高分辨率气溶胶光学厚度反演的城市复杂地表反射率估算方法。本发明首先为每个像元构建像元数据库，根据归一化植被指数(NDVI)将城市复杂地表分为三类，然后分别对每类地表构建精确的地表反射率估算方案，可广泛应用于高分辨率卫星传感器的地表反射率准确估算。
为了达到上述目的，本发明提供的技术方案是一种支持高分辨率气溶胶光学厚度反演的城市复杂地表反射率估算方法，包括以下步骤：
步骤1，构建像元数据库，具体包括以下几个子步骤：
步骤1.1，获取影像；
步骤1.2，剔除无效数据；
步骤1.3，进行大气校正；
步骤1.4，基于红光波段和近红外波段的地表反射率计算NDVI；
步骤1.5，构建像元数据库；
步骤2，确定城市复杂地表的类型，具体包括以下几个子步骤：
步骤2.1，标记浓密植被地表；
步骤2.2，标记伪不变地表；
步骤2.3，标记稀疏植被地表；
步骤3，根据地表类型分别估算反射率，具体包括以下几个子步骤：
步骤3.1，估算浓密植被地表的反射率；
步骤3.2，估算伪不变地表的反射率；
步骤3.3，估算稀疏植被地表的反射率。
而且，所述步骤1.1是获取多年的高分辨率多光谱卫星影像，包含蓝光波段、红光波段、近红外波段、短波红外波段的表观反射率。
而且，所述步骤1.2剔除无效数据是逐影像进行云层、云阴影、水体、冰雪的检测，并将检测结果标记为无效数据。
而且，所述步骤1.3是对有效的蓝光波段、红光波段、近红外波段的表观反射率数据进行大气校正，得到蓝光波段、红光波段、近红外波段的地表反射率，此处的地表反射率是基于一次卫星观测及大气校正的结果，其准确度需要后续处理以进一步提高，在短波红外波段，大气(瑞利、气溶胶)散射对表观反射率的贡献很小，可以被忽略，所以短波红外波段的表观反射率即地表反射率。
而且，所述步骤1.5中由于多景影像是相互重叠的，且其中经纬度相同的像元在不同影像上的地表反射率是变化的，因此收集经纬度相同像元的地表反射率、NDVI、观测的月份和观测时刻的视图几何数据(包含太阳天顶角、太阳方位角、卫星天顶角、卫星方位角)，为每个像元分别建立数据库。
而且，所述步骤2.1中若某像元在某一景影像中的NDVI大于0.55，则标记该影像上的该像元为浓密植被地表，此标记仅对当前影像有效。
而且，所述步骤2.2中若某像元在所有影像中NDVI都小于0.2，则标记该像元为伪不变地表，此标记对所有影像均有效。
而且，所述步骤2.3是将影像中未被标记的部分定为稀疏植被地表。
而且，所述步骤3.1是采用暗目标法估算浓密植被地表的反射率，浓密植被的蓝、红光波段与短波红外波段的地表反射率存在稳定的线性关系，因此根据当前影像提供的短波红外波段的表观反射率，基于线性关系可以估算当前影像中浓密植被地表的的蓝光波段和红光波段的地表反射率。
而且，所述步骤3.2中伪不变地表主要包含城市区域内的建筑物、道路等城市人造地表，表现出显著的双向反射率特性，由于伪不变地表物理和光学特性可在长时间内保持不变，因此采用半经验核驱动(RTLS)模型进行估算，具体步骤如下：
步骤3.2.1，针对每一个伪不变地表像元，从像元数据库中获取蓝光波段和红光波段的地表反射率，以及对应的视图几何数据(包含太阳天顶角、太阳方位角、卫星天顶角、卫星方位角)。
步骤3.2.2，将步骤3.2.1获取的视图几何数据代入RossThick核和LiSparse核，计算每个伪不变地表像元的两个内核参数。
步骤3.2.3，将步骤3.2.2获取的两个内核参数和步骤3.2.1获取的蓝光波段、红光波段的地表反射率，代入半经验核驱动模型，通过最小二乘的方法，分别确定该像元的蓝光波段和红光波段在半经验核驱动模型中的三个核系数，在每一个伪不变地表像元上，都应建立蓝光波段和红光波段两个模型。
步骤3.2.4，基于步骤3.2.2计算得到的两个内核参数和步骤3.2.3得到的半经验核驱动模型中的三个核系数，计算当前影像中伪不变地表像元的地表反射率，或将步骤3.2.3得到的这三个核系数代入辐射传输模型，模拟在地表各向异性状态下的地表-大气-天顶反射率辐射传输特性。
而且，所述步骤3.3是基于中值算法合成稀疏植被地表的逐月地表反射率。稀疏植被地表包含农田、草地等地表覆盖类型，其地表反射率随时间变化较大，但在多年同月份的比较中变化较小，因此可以假设多年同月份的地表反射率保持不变，用同月份的中值来估算稀疏植被地表的反射率，具体步骤如下：
步骤3.3.1，针对每一个稀疏植被地表像元，从像元数据库中获取蓝光波段和红光波段的地表反射率，以及观测时的月份。
步骤3.3.2，将获取的蓝光波段和红光波段的地表反射率按月份分组，并统计每个分组中地表反射率的中值，作为该像元在各月份的地表反射率。在每一个稀疏植被地表像元上，都应得到蓝光波段和红光波段在12个月份上的中值。
步骤3.3.3，基于每个稀疏植被像元在各月份上的中值与观测时的月份，确定当前影像中稀疏植被地表像元的地表反射率。
与现有技术相比，本发明具有如下优点：本发明实现的高空间分辨率地表反射率确定方案，可根据城市地区不同地表的反射特性和变化特性，精确地获取不同地表的地表反射率，且在伪不变地表上考虑了地表的双向反射特性。基于该发明获取的地表反射率可用于支撑城市复杂地表的高分辨率卫星AOD反演。此外，该方法或该方法获取的地表反射率可以推广到具有同类型传感器的卫星上，进一步提高高空间分辨率卫星在AOD反演中的应用。
附图说明
图1为本发明实施例的流程图。
图2为本发明实施例估算的地表反射率精度评价图，其中图2(a)为蓝色波段的地表反射率精度评价图，图2(b)为红色波段的地表反射率精度评价图。
具体实施方式
本发明提供一种支持高分辨率气溶胶光学厚度反演的城市复杂地表反射率估算方法，首先为每个像元构建像元数据库，根据归一化植被指数(NDVI)将城市复杂地表分为三类，然后分别对每类地表构建精确的地表反射率估算方案，可广泛应用于高分辨率卫星传感器的地表反射率准确估算。
下面结合附图和实施例对本发明的技术方案作进一步说明。
如图1所示，本发明实施例的流程包括以下步骤：
步骤1，构建像元数据库，具体包括以下几个子步骤：
步骤1.1，获取影像，获取多年的高分辨率多光谱卫星影像，包含蓝光波段、红光波段、近红外波段、短波红外波段的表观反射率。
步骤1.2，剔除无效数据，逐影像进行云层、云阴影、水体、冰雪的检测，并将检测结果标记为无效数据。
步骤1.3，大气校正，对有效的蓝光波段、红光波段、近红外波段的表观反射率数据进行大气校正，得到蓝光波段、红光波段、近红外波段的地表反射率，此处的地表反射率是基于一次卫星观测及大气校正的结果，其准确度需要后续处理以进一步提高，在短波红外波段，大气(瑞利、气溶胶)散射对表观反射率的贡献很小，可以被忽略，所以短波红外波段的表观反射率即地表反射率。
步骤1.4，基于红光波段和近红外波段的地表反射率计算NDVI。
步骤1.5，构建像元数据库，多景影像是相互重叠的，且其中经纬度相同的像元在不同影像上的地表反射率是变化的，收集经纬度相同像元的地表反射率、NDVI、观测的月份和观测时刻的视图几何数据(包含太阳天顶角、太阳方位角、卫星天顶角、卫星方位角)，为每个像元分别建立数据库。
步骤2，确定城市复杂地表的类型，具体包括以下几个子步骤：
步骤2.1，标记浓密植被地表，若某像元在某一景影像中的NDVI大于0.55，则标记该影像上的该像元为浓密植被地表，此标记仅对当前影像有效。
步骤2.2，标记伪不变地表，若某像元在所有影像中NDVI都小于0.2，则标记该像元为伪不变地表，此标记对所有影像均有效。
步骤2.3，标记稀疏植被地表，将影像中未被标记的部分定为稀疏植被地表。
步骤3，根据地表类型分别估算反射率，具体包括以下几个子步骤：
步骤3.1，估算浓密植被地表，采用暗目标法估算浓密植被地表的反射率，浓密植被的蓝、红光波段与短波红外波段的地表反射率存在稳定的线性关系，因此根据当前影像提供的短波红外波段的表观反射率，基于线性关系可以估算当前影像中浓密植被地表的的蓝光波段和红光波段的地表反射率。
步骤3.2，估算伪不变地表，基于半经验核驱动模型估算伪不变地表的反射率，伪不变地表主要包含城市区域内的建筑物、道路等城市人造地表，表现出显著的双向反射率特性，但伪不变地表物理和光学特性可在长时间内保持不变，可采用半经验核驱动(RTLS)模型进行估算，具体步骤如下：
步骤3.2.1，针对每一个伪不变地表像元，从像元数据库中获取蓝光波段和红光波段的地表反射率，以及对应的视图几何数据(包含太阳天顶角、太阳方位角、卫星天顶角、卫星方位角)。
步骤3.2.2，将步骤3.2.1获取的视图几何数据代入RossThick核和LiSparse核，计算每个伪不变地表像元的两个内核参数。
步骤3.2.3，将步骤3.2.2获取的两个内核和步骤3.2.1获取的蓝光波段、红光波段的地表反射率，代入半经验核驱动模型，通过最小二乘的方法，分别确定该像元的蓝光波段和红光波段在半经验核驱动模型中的三个核系数，在每一个伪不变地表像元上，都应建立蓝光波段和红光波段两个模型。
步骤3.2.4，基于步骤3.2.2计算得到的两个内核参数和步骤3.2.3得到的半经验核驱动模型中的三个核系数，计算当前影像中伪不变地表像元的地表反射率，或将步骤3.2.3得到的这三个核系数代入辐射传输模型，模拟在地表各向异性状态下的地表-大气-天顶反射率辐射传输特性。
步骤3.3，估算稀疏植被地表，基于中值算法合成稀疏植被地表的逐月地表反射率。稀疏植被地表包含农田、草地等地表覆盖类型，其地表反射率随时间变化较大，但在多年同月份的比较中变化较小，因此可以假设多年同月份的地表反射率保持不变，用同月份的中值来估算稀疏植被地表的反射率，具体步骤如下：
步骤3.3.1，针对每一个稀疏植被地表像元，从像元数据库中获取蓝光波段和红光波段的地表反射率，以及观测时的月份。
步骤3.3.2，将获取的蓝光波段和红光波段的地表反射率按月份分组，并统计每个分组中地表反射率的中值，作为该像元在各月份的地表反射率。在每一个稀疏植被地表像元上，都应得到蓝光波段和红光波段在12个月份上的中值。
步骤3.3.3，基于每个稀疏植被像元在各月份上的中值与观测时的月份，确定当前影像中稀疏植被地表像元的地表反射率。
下面以高分辨卫星Sentinel-2为例，进一步阐述本发明的技术方案。
1.实施目标
基于Sentinel-2卫星的影像，估算2017-2020年北京地区蓝光波段和红光波段的30米分辨率的地表反射率，可用于高分辨率气溶胶光学厚度反演。
2.数据选择
样例数据选取由Sentinel-2A/B两颗卫星携带的多光谱成像仪(MSI)获取的L1C级数据，选择研究区域是北京市。MSI传感器共设置有13个光谱波段，有4个可见光波段、4个红边波段、1个水汽波段、1个近红外波段和3个短波红外波段，各波段的空间分辨率为10-60m。为了减小传感器噪声的影响，进一步提高信噪比，所有的波段将被重采样为30m。
3.实施过程
(1)构建像元数据库
①收集2017-2020年内北京地区所有云量小于30％的Sentinel-2A/B的L1C级数据，其中包含蓝光波段表观反射率
Figure BDA0003079648810000071
红光波段表观反射率
Figure BDA0003079648810000072
近红外波段表观反射率
Figure BDA0003079648810000073
短波红外波段表观反射率
Figure BDA0003079648810000074
以及观测的月份m和观测时刻的太阳天顶角θ0、太阳方位角
Figure BDA0003079648810000075
卫星天顶角θ、卫星方位角
Figure BDA0003079648810000076
②利用Fmask4.2算法，对收集的每一景影像分别进行云层、云阴影、内陆水体和冰雪的检测，并将检测结果标记为无效数据。
③采用大气校正代码(LaSRC)，对有效的蓝光波段、红光波段、近红外波段的表观反射率数据进行大气校正，获取每一景大气校正后影像的蓝光波段地表反射率ρ′BLUE、红光波段地表反射率ρ′RED、近红外波段地表反射率ρ′NIR，将短波红外的表观反射率当做其地表反射率ρ′SWIR。
④根据数据库中的红光波段和近红外波段的地表反射率，计算归一化植被指数(NDVI)：
Figure BDA0003079648810000077
⑤为经纬度相同的像元构建数据库，其中包含每一景的影像上对应位置的经过大气校正的蓝光波段地表反射率ρ′BLUE、经过大气校正的红光波段地表反射率ρ′RED、短波红外波段地表反射率ρ′SWIR、归一化植被指数NDVI，以及观测月份m、观测时刻的太阳天顶角θ0、太阳方位角
Figure BDA0003079648810000078
卫星天顶角θ、卫星方位角
Figure BDA0003079648810000079
(2)确定城市复杂地表的类型
①浓密植被地表：逐影像筛选NDVI大于0.55的像元，标记对应影像上的对应像元为浓密植被地表。
②伪不变地表：筛选在所有影像中NDVI都小于0.2的像元，标记所有影像上的该像元为伪不变地表。
③稀疏植被地表：逐影像筛选未标记的像元，标记对应影像上的对应像元为稀疏植被地表。
(3)不同地表类型下地表反射率估算方案的确定
①确定浓密植被地表的反射率
基于蓝光波段、红光波段与短波红外波段ρ′SWIR地表反射率间稳定的线性关系，可通过公式(2)和公式(3)估算浓密植被地表像元的准确的蓝光波段地表反射率ρBLUE和准确的红光波段地表反射率ρRED：
ρBLUE＝0.25×ρ′SWIR (2)
ρRED＝0.50×ρ′SWIR (3)
②确定伪不变地表的反射率
采用半经验核驱动(RTLS)模型估算伪不变地表像元的准确的蓝光波段地表反射率ρBLUE和准确的红光波段地表反射率ρRED，具体步骤如下：
a)针对每个伪不变像元，收集经大气校正后的蓝光波段地表反射率ρ′BLUE、红光波段地表反射率ρ′RED，以及观测时刻的太阳天顶角θ0、太阳方位角
Figure BDA0003079648810000081
卫星天顶角θ、卫星方位角
Figure BDA0003079648810000082
b)基于Ross-LiBRDF模型，计算模型中RossThick内核参数KVOL和LiSparse内核参数KGEO：
Figure BDA0003079648810000083
Figure BDA0003079648810000084
其中：
Figure BDA0003079648810000085
Figure BDA0003079648810000086
cosξ＝cosθ0cosθ+sinθ0sinθcosφ (8)
Figure BDA0003079648810000087
每一个伪不变地表像元上都能得到一组内核参数KVOL和KGEO。
c)将公式(4)和公式(5)计算的内核参数KVOL、KGEO和经大气校正后的的蓝光波段地表反射率ρ′BLUE、红光波段地表反射率ρ′RED，分别代入公式(10)和(11)，并基于最小二乘原则求得蓝光波段的系数fISO,BLUE、fVOL,BLUE、fGEO,BLUE和红光波段的系数fISO,RED、fVOL,RED、fGEO,RED：
ρ′BLUE＝fISO,BLUE+fVOL,BLUE×KVOL+fGEO,BLUE×KGEO (10)
ρ′RED＝fISO,RED+fVOL,RED×KVOL+fGEO,RED×KGEO (11)
其中，系数fISO,BLUE和fISO,RED代表地表的各项同性散射，系数fVOL,BLUE和fVOL,RED代表均一地表的体积散射，系数fGEO,BLUE和fGEO,RED代表几何-光学表面散射。
d)针对某个伪不变像元，将内核参数KVOL和KGEO，以及蓝光波段的系数fISO,BLUE、fVOL,BLUE、fGEO,BLUE和红光波段的系数fISO,RED、fVOL,RED、fGEO,RED带入公式(12)和(13)估算伪不变地表像元的准确的蓝光波段地表反射率ρBLUE和红光波段地表反射率ρRED：
ρBLUE＝fISO,BLUE+fVOL,BLUE×KVOL+fGEO,BLUE×KGEO (12)
ρRED＝fISO,RED+fVOL,RED×KVOL+fGEO,RED×KGEO (13)
此处也可以基于辐射传输模型，将红、蓝波段的三个核系数，代入辐射传输模型，模拟在地表各向异性状态下的地表-大气-天顶反射率辐射传输特性。
③确定稀疏植被地表的反射率
取稀疏植被地表在同月份中地表反射率的中值为准确的蓝光波段地表反射率ρBLUE和准确的红光波段地表反射率ρRED，具体步骤如下：
a)针对每个稀疏植被像元，收集经大气校正后的蓝光波段地表反射率ρ′BLUE、红光波段地表反射率ρ′RED，以及观测的月份m；
b)针对某个稀疏植被像元，依次统计其所有年份中第M(M∈(1,2,3，…12))月的经大气校正后的蓝光波段地表反射率的中值ρBLUE,MED和红光波段地表反射率的中值ρRED,MED：
ρBLUE,MED(M)＝Median{ρ′BLUE(m＝M)} (14)
ρRED,MED(M)＝Median{ρ′RED(m＝M)} (15)
c)针对某个稀疏植被像元，根据其观测月份m估算稀疏植被地表像元的准确的蓝光波段地表反射率ρBLUE和红光波段地表反射率ρRED：
ρBLUE＝ρBLUE,MED(m) (16)
ρRED＝ρRED,MED(m) (17)
4.精度评价
北京市区设有4个AERONET站点，可用于监测气溶胶的光学和物理特性。基于获取卫星影像时间、太阳-卫星视图几何和天顶反射率，结合对应地面AERONET站点提供的气溶胶光学厚度、臭氧浓度、水汽柱浓度、地表高程参数等，输入6S辐射传输模型中，进行精确的大气校正，可获取各站点上对应卫星影像在不同波段上的精确的地表反射率。与本发明获取的地表反射率进行对比，结果展示如图2。结果表明，本方法获取的地表反射率具有较高的精度，在红色波段的可决系数(R2)高达0.982，由于蓝光波段更容易受气溶胶的影响，其可决系数达到了0.823，均方根误差(RMSE)为0.007，平均绝对误差(MAE)为0.005，均能满足气溶胶监测的需要。
本文中所描述的具体实施例仅仅是对本发明精神作举例说明。本发明所属技术领域的技术人员可以对所描述的具体实施例做各种各样的修改或补充或采用类似的方式替代，但并不会偏离本发明的精神或者超越所附权利要求书所定义的范围。

Claims (10)

1.一种支持高分辨率气溶胶光学厚度反演的城市复杂地表反射率估算方法，其特征在于，包括如下步骤：
步骤1，构建像元数据库，具体包括以下几个子步骤：
步骤1.1，获取影像；
步骤1.2，剔除无效数据；
步骤1.3，进行大气校正；
步骤1.4，基于红光波段和近红外波段的地表反射率计算NDVI；
步骤1.5，构建像元数据库；
步骤2，确定城市复杂地表的类型，具体包括以下几个子步骤：
步骤2.1，标记浓密植被地表；
步骤2.2，标记伪不变地表；
步骤2.3，将影像中未被标记的部分标记为稀疏植被地表；
步骤3，根据地表类型分别估算反射率，具体包括以下几个子步骤：
步骤3.1，估算浓密植被地表的反射率；
步骤3.2，估算伪不变地表的反射率；
步骤3.3，估算稀疏植被地表的反射率。
2.如权利要求1所述的一种支持高分辨率气溶胶光学厚度反演的城市复杂地表反射率估算方法，其特征在于：所述步骤1.1是获取多年的高分辨率多光谱卫星影像，包含蓝光波段、红光波段、近红外波段、短波红外波段的表观反射率。
3.如权利要求1所述的一种支持高分辨率气溶胶光学厚度反演的城市复杂地表反射率估算方法，其特征在于：所述步骤1.2剔除无效数据是逐影像进行云层、云阴影、水体、冰雪的检测，并将检测结果标记为无效数据。
4.如权利要求1所述的一种支持高分辨率气溶胶光学厚度反演的城市复杂地表反射率估算方法，其特征在于：所述步骤1.3是对有效的蓝光波段、红光波段、近红外波段的表观反射率数据进行大气校正，得到蓝光波段、红光波段、近红外波段的地表反射率，此处的地表反射率是基于一次卫星观测及大气校正的结果，其准确度需要后续处理以进一步提高，在短波红外波段，大气散射对表观反射率的贡献很小，可以被忽略，所以短波红外波段的表观反射率即地表反射率。
5.如权利要求1所述的一种支持高分辨率气溶胶光学厚度反演的城市复杂地表反射率估算方法，其特征在于：所述步骤1.5中由于多景影像是相互重叠的，且其中经纬度相同的像元在不同影像上的地表反射率是变化的，因此收集经纬度相同像元的地表反射率、NDVI、观测的月份、观测时刻的视图几何数据，包括太阳天顶角、太阳方位角、卫星天顶角和卫星方位角，为每个像元分别建立数据库。
6.如权利要求1所述的一种支持高分辨率气溶胶光学厚度反演的城市复杂地表反射率估算方法，其特征在于：所述步骤2.1中若某像元在某一景影像中的NDVI大于α，则标记该影像上的该像元为浓密植被地表，此标记仅对当前影像有效。
7.如权利要求1所述的一种支持高分辨率气溶胶光学厚度反演的城市复杂地表反射率估算方法，其特征在于：所述步骤2.2中若某像元在所有影像中NDVI都小于β，则标记该像元为伪不变地表，此标记对所有影像均有效。
8.如权利要求1所述的一种支持高分辨率气溶胶光学厚度反演的城市复杂地表反射率估算方法，其特征在于：所述步骤3.1是采用暗目标法估算浓密植被地表的反射率，浓密植被的蓝、红光波段与短波红外波段的地表反射率存在稳定的线性关系，因此根据当前影像提供的短波红外波段的表观反射率，基于线性关系可以估算当前影像中浓密植被地表的的蓝光波段和红光波段的地表反射率，计算公式如下：
ρBLUE＝λ×ρ′SWIR+ν (2)
ρRED＝μ×ρ′SWIR+τ (3)
其中，ρ′SWIR为短波红外波段，ρBLUE为浓密植被地表像元的准确的蓝光波段地表反射率，ρRED为浓密植被地表像元的准确的红光波段地表反射率。
9.如权利要求1所述的一种支持高分辨率气溶胶光学厚度反演的城市复杂地表反射率估算方法，其特征在于：所述步骤3.2中伪不变地表主要包含城市区域内的建筑物、道路等城市人造地表，表现出显著的双向反射率特性，由于伪不变地表物理和光学特性可在长时间内保持不变，因此采用半经验核驱动模型进行估算，具体步骤如下：
步骤3.2.1，针对每一个伪不变地表像元，从像元数据库中获取蓝光波段的地表反射率和红光波段的地表反射率，以及对应的视图几何数据，包含太阳天顶角、太阳方位角、卫星天顶角和卫星方位角；
步骤3.2.2，将步骤3.2.1获取的视图几何数据代入RossThick核和LiSparse核，计算每个伪不变地表像元的两个内核参数KVOL和KGEO：
Figure FDA0003079648800000031
Figure FDA0003079648800000032
其中：
Figure FDA0003079648800000033
Figure FDA0003079648800000034
cosξ＝cosθ0cosθ+sinθ0sinθcosφ (8)
Figure FDA0003079648800000035
式中，θ0为太阳天顶角，
Figure FDA0003079648800000036
为太阳方位角，θ为卫星天顶角，
Figure FDA0003079648800000037
为卫星方位角，每一个伪不变地表像元上都能得到一组内核参数KVOL和KGEO；
步骤3.2.3，将步骤3.2.2获取的两个内核参数和步骤3.2.1获取的蓝光波段、红光波段的地表反射率，代入半经验核驱动模型，通过最小二乘的方法，分别确定该像元的蓝光波段和红光波段在半经验核驱动模型中的三个核系数，在每一个伪不变地表像元上，都应建立蓝光波段和红光波段两个模型；
ρ′BLUE＝fISO,BLUE+fVOL,BLUE×KVOL+fGEO,BLUE×KGEO (10)
ρ′RED＝fISO,RED+fVOL,RED×KVOL+fGEO,RED×KGEO (11)
式中，ρ′BLUE是从像元数据库中获取的蓝光波段的地表反射率，ρ′RED是从像元数据库中获取的红光波段的地表反射率，fISO,BLUE、fVOL,BLUE、fGEO,BLUE是像元的蓝光波段在半经验核驱动模型中的三个核系数，fISO,RED、fVOL,RED、fGEO,RED是像元的红光波段在半经验核驱动模型中的三个核系数，系数fISO,BLUE和fISO,RED代表地表的各项同性散射，系数fVOL,BLUE和fVOL,RED代表均一地表的体积散射，系数fGEO,BLUE和fGEO,RED代表几何-光学表面散射；
步骤3.2.4，基于步骤3.2.2计算得到的两个内核参数和步骤3.2.3得到的半经验核驱动模型中的三个核系数，计算当前影像中伪不变地表像元的地表反射率，或将步骤3.2.3得到的这三个核系数代入辐射传输模型，模拟在地表各向异性状态下的地表-大气-天顶反射率辐射传输特性；
ρBLUE＝fISO,BLUE+fVOL,BLUE×KVOL+fGEO,BLUE×KGEO (12)
ρRED＝fISO,RED+fVOL,RED×KVOL+fGEO,RED×KGEO (13)
式中，ρBLUE为伪不变地表像元的准确的蓝光波段地表反射率，ρRED为伪不变地表像元的准确的红光波段地表反射率。
10.如权利要求1所述的一种支持高分辨率气溶胶光学厚度反演的城市复杂地表反射率估算方法，其特征在于：所述步骤3.3是基于中值算法合成稀疏植被地表的逐月地表反射率，稀疏植被地表包含农田、草地等地表覆盖类型，具体步骤如下：
步骤3.3.1，针对每一个稀疏植被地表像元，从像元数据库中获取蓝光波段和红光波段的地表反射率，以及观测时的月份；
步骤3.3.2，将获取的蓝光波段和红光波段的地表反射率按月份分组，并统计每个分组中地表反射率的中值，作为该像元在各月份的地表反射率，在每一个稀疏植被地表像元上，都应得到蓝光波段和红光波段在12个月份上的中值，
ρBLUE,MED(M)＝Median{ρ′BLUE(m＝M)} (14)
ρRED,MED(M)＝Median{ρ′RED(m＝M)} (15)
式中，m为观测的月份，M∈(1,2,3，…12)，ρ′BLUE为从像元数据库中获取的蓝光波段的地表反射率，ρ′RED为从像元数据库中获取的红光波段的地表反射率，ρBLUE,MED为经大气校正后的蓝光波段地表反射率的中值，ρRED,MED为经大气校正后的红光波段地表反射率的中值；
步骤3.3.3，基于每个稀疏植被像元在各月份上的中值与观测时的月份，确定当前影像中稀疏植被地表像元的地表反射率，
ρBLUE＝ρBLUE,MED(m) (16)
ρRED＝ρRED,MED(m) (17)
式中，ρBLUE为稀疏植被地表像元的准确的蓝光波段地表反射率，ρRED为稀疏植被地表像元的准确的红光波段地表反射率。
CN202110562764.XA 2021-05-24 2021-05-24 一种支持高分辨率气溶胶光学厚度反演的城市复杂地表反射率估算方法 Active CN113324915B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202110562764.XA CN113324915B (zh) 2021-05-24 2021-05-24 一种支持高分辨率气溶胶光学厚度反演的城市复杂地表反射率估算方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202110562764.XA CN113324915B (zh) 2021-05-24 2021-05-24 一种支持高分辨率气溶胶光学厚度反演的城市复杂地表反射率估算方法

Publications (2)

Publication Number Publication Date
CN113324915A true CN113324915A (zh) 2021-08-31
CN113324915B CN113324915B (zh) 2022-10-21

Family

ID=77416428

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202110562764.XA Active CN113324915B (zh) 2021-05-24 2021-05-24 一种支持高分辨率气溶胶光学厚度反演的城市复杂地表反射率估算方法

Country Status (1)

Country Link
CN (1) CN113324915B (zh)

Cited By (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN114266968A (zh) * 2021-12-16 2022-04-01 河南大学 一种城市不同土地覆盖类型遥感自动解译方法
CN114970214A (zh) * 2022-07-28 2022-08-30 南京航天宏图信息技术有限公司 一种气溶胶光学厚度反演方法及装置
CN117347282A (zh) * 2023-08-22 2024-01-05 中南大学 星基气溶胶光学厚度反演方法、装置及系统和存储介质
CN117607919A (zh) * 2023-11-17 2024-02-27 中国科学院大气物理研究所 一种基于城市建筑物阴影的气溶胶卫星遥感反演方法

Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20030150992A1 (en) * 2002-02-13 2003-08-14 Chavez Pat S. Field based spectral radiometer
CN106407656A (zh) * 2016-08-29 2017-02-15 中国科学院遥感与数字地球研究所 一种基于高分辨率卫星影像数据的气溶胶光学厚度反演方法
CN106680273A (zh) * 2016-08-29 2017-05-17 中国科学院遥感与数字地球研究所 一种高空间分辨率卫星地表反射率反演方法
CN110501716A (zh) * 2019-07-29 2019-11-26 武汉大学 基于单光子激光雷达背景噪声率的地表分类方法
CN111537510A (zh) * 2020-05-09 2020-08-14 东北林业大学 一种基于空间信息技术的农田防护林防风效应计量方法

Patent Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20030150992A1 (en) * 2002-02-13 2003-08-14 Chavez Pat S. Field based spectral radiometer
CN106407656A (zh) * 2016-08-29 2017-02-15 中国科学院遥感与数字地球研究所 一种基于高分辨率卫星影像数据的气溶胶光学厚度反演方法
CN106680273A (zh) * 2016-08-29 2017-05-17 中国科学院遥感与数字地球研究所 一种高空间分辨率卫星地表反射率反演方法
CN110501716A (zh) * 2019-07-29 2019-11-26 武汉大学 基于单光子激光雷达背景噪声率的地表分类方法
CN111537510A (zh) * 2020-05-09 2020-08-14 东北林业大学 一种基于空间信息技术的农田防护林防风效应计量方法

Non-Patent Citations (6)

* Cited by examiner, † Cited by third party
Title
N. C. HSU ET AL.: "Enhanced Deep Blue aerosol retrieval algorithm: The second generation", 《JOURNAL OF GEOPHYSICAL RESEARCH: ATMOSPHERES》 *
N. C. HSU ET AL.: "Enhanced Deep Blue aerosol retrieval algorithm: The second generation", 《JOURNAL OF GEOPHYSICAL RESEARCH: ATMOSPHERES》, 12 August 2013 (2013-08-12), pages 9296 - 9315 *
ROUJEAN ET AL.: "A bidirectional reflectance model of the Earth"s surface for the correction of remote sensing data", 《JOURNAL OF GEOPHYSICAL RESEARCH》 *
ROUJEAN ET AL.: "A bidirectional reflectance model of the Earth"s surface for the correction of remote sensing data", 《JOURNAL OF GEOPHYSICAL RESEARCH》, 31 December 1992 (1992-12-31), pages 466 - 468 *
王利民等: ""基于暗目标法和 GF-1 的农作物光合有效辐射反演"", 《农业工程学报》 *
王利民等: ""基于暗目标法和 GF-1 的农作物光合有效辐射反演"", 《农业工程学报》, 23 November 2016 (2016-11-23), pages 184 - 191 *

Cited By (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN114266968A (zh) * 2021-12-16 2022-04-01 河南大学 一种城市不同土地覆盖类型遥感自动解译方法
CN114266968B (zh) * 2021-12-16 2023-01-31 河南大学 一种城市不同土地覆盖类型遥感自动解译方法
CN114970214A (zh) * 2022-07-28 2022-08-30 南京航天宏图信息技术有限公司 一种气溶胶光学厚度反演方法及装置
CN117347282A (zh) * 2023-08-22 2024-01-05 中南大学 星基气溶胶光学厚度反演方法、装置及系统和存储介质
CN117347282B (zh) * 2023-08-22 2024-05-28 中南大学 星基气溶胶光学厚度反演方法、装置及系统和存储介质
CN117607919A (zh) * 2023-11-17 2024-02-27 中国科学院大气物理研究所 一种基于城市建筑物阴影的气溶胶卫星遥感反演方法

Also Published As

Publication number Publication date
CN113324915B (zh) 2022-10-21

Similar Documents

Publication Publication Date Title
CN113324915B (zh) 一种支持高分辨率气溶胶光学厚度反演的城市复杂地表反射率估算方法
Mumby et al. Mapping marine environments with IKONOS imagery: enhanced spatial resolution can deliver greater thematic accuracy
Kuc et al. Sentinel-2 imagery for mapping and monitoring imperviousness in urban areas
Bryson et al. Kite aerial photography for low-cost, ultra-high spatial resolution multi-spectral mapping of intertidal landscapes
Li et al. An evaluation of the use of atmospheric and BRDF correction to standardize Landsat data
Richards et al. Error correction and registration of image data
CN108303044A (zh) 一种叶面积指数获取方法及系统
Kant et al. Satellite-based analysis of the role of land use/land cover and vegetation density on surface temperature regime of Delhi, India
Dai et al. Estimating river surface elevation from ArcticDEM
Chu et al. Technical framework for shallow-water bathymetry with high reliability and no missing data based on time-series sentinel-2 images
Baldi et al. Validation and comparison of different techniques for the derivation of digital elevation models and volcanic monitoring (Vulcano Island, Italy)
Hsu et al. Cross-estimation of Soil Moisture Using Thermal Infrared Images with Different Resolutions.
Unger et al. Modeling of the urban heat island pattern based on the relationship between surface and air temperatures
Tao et al. A Landsat-derived annual inland water clarity dataset of China between 1984 and 2018
Jiang et al. Synergistic use of optical and InSAR data for urban impervious surface mapping: a case study in Hong Kong
CN110009584B (zh) 基于参考光谱匹配的多光谱遥感影像大气校正系统及方法
CN111413296A (zh) 一种考虑地表非朗伯特性的气溶胶光学厚度遥感反演方法
Dehm et al. SUAS based multispectral imagery for monitoring wetland inundation and vegetation
Scambos et al. Improving digital elevation models over ice sheets using AVHRR-based photoclinometry
Fraser et al. 3D building reconstruction from high-resolution Ikonos stereo-imagery
CN113034519A (zh) 滨海滩涂湿地植被恢复建设区域确定方法和装置
Zhang et al. Prediction of Soil Organic Carbon Content Using Sentinel-1/2 and Machine Learning Algorithms in Swamp Wetlands in Northeast China
Mozafari et al. Air pollution estimation using aerosol optical thickness by oli images in Tehran
Zhou et al. Comparison between in situ and MODIS-derived spectral reflectances of snow and sea ice in the Amundsen Sea, Antarctica
CN110705089B (zh) 一种细模态气溶胶参数反演方法

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