CN110261873A - 一种基于分段统计的大气气溶胶反演方法 - Google Patents

一种基于分段统计的大气气溶胶反演方法 Download PDF

Info

Publication number
CN110261873A
CN110261873A CN201910454155.5A CN201910454155A CN110261873A CN 110261873 A CN110261873 A CN 110261873A CN 201910454155 A CN201910454155 A CN 201910454155A CN 110261873 A CN110261873 A CN 110261873A
Authority
CN
China
Prior art keywords
pixel
aerosol
section
reflectivity
pixel set
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
CN201910454155.5A
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.)
University of Electronic Science and Technology of China
Original Assignee
University of Electronic Science and Technology of China
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 University of Electronic Science and Technology of China filed Critical University of Electronic Science and Technology of China
Priority to CN201910454155.5A priority Critical patent/CN110261873A/zh
Publication of CN110261873A publication Critical patent/CN110261873A/zh
Priority to US16/885,021 priority patent/US20210318253A1/en
Priority to US17/489,300 priority patent/US11733175B2/en
Pending legal-status Critical Current

Links

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S17/00—Systems using the reflection or reradiation of electromagnetic waves other than radio waves, e.g. lidar systems
    • G01S17/88—Lidar systems specially adapted for specific applications
    • G01S17/95—Lidar systems specially adapted for specific applications for meteorological use
    • 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/84—Systems specially adapted for particular applications
    • G01N21/88—Investigating the presence of flaws or contamination
    • G01N21/94—Investigating contamination, e.g. dust
    • 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
    • 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
    • G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V8/00—Prospecting or detecting by optical means
    • G01V8/02—Prospecting
    • 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
    • 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
    • G01N2021/3531—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 without instrumental source, i.e. radiometric
    • 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/84—Systems specially adapted for particular applications
    • G01N21/85—Investigating moving fluids or granular solids
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01N—INVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N2201/00—Features of devices classified in G01N21/00
    • G01N2201/02—Mechanical
    • G01N2201/021—Special mounting in general
    • G01N2201/0214—Airborne
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
    • G01V3/00—Electric or magnetic prospecting or detecting; Measuring magnetic field characteristics of the earth, e.g. declination, deviation
    • G01V3/15—Electric or magnetic prospecting or detecting; Measuring magnetic field characteristics of the earth, e.g. declination, deviation specially adapted for use during transport, e.g. by a person, vehicle or boat
    • G01V3/16—Electric or magnetic prospecting or detecting; Measuring magnetic field characteristics of the earth, e.g. declination, deviation specially adapted for use during transport, e.g. by a person, vehicle or boat specially adapted for use from aircraft

Landscapes

  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Immunology (AREA)
  • Pathology (AREA)
  • Health & Medical Sciences (AREA)
  • Chemical & Material Sciences (AREA)
  • Analytical Chemistry (AREA)
  • Biochemistry (AREA)
  • General Health & Medical Sciences (AREA)
  • Spectroscopy & Molecular Physics (AREA)
  • Engineering & Computer Science (AREA)
  • General Life Sciences & Earth Sciences (AREA)
  • Geophysics (AREA)
  • Computer Networks & Wireless Communication (AREA)
  • Electromagnetism (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Remote Sensing (AREA)
  • Image Analysis (AREA)
  • Image Processing (AREA)

Abstract

该发明公开了一种基于分段统计的大气气溶胶反演方法,属于遥感信息技术领域。本发明首先得到多波段遥感影像和反演波段对应的气溶胶光学厚度查找表，然后按照中红外2.1微米波段的表观反射率区间对像元进行划分和筛选，然后对得到的像元集合按照中红外1.6微米波段的表观反射率区间进行进一步划分和筛选，将最终得到的像元集合按照像元数量划分为两类，像元较多的划为一类，其他为另一类；将像元较多的一类作为基准部分进行反演，以清洁区段的地表反射率作为整个像元集合的地表反射率，反演得到气溶胶厚度值，然后以这些像元作为基准，对另一类进行反演。较传统的暗像元算法，本发明对亮地表地区也具有良好的反演效果。具有更广泛的适用性。

Description

一种基于分段统计的大气气溶胶反演方法
技术领域
本发明属于遥感信息技术领域，更为具体地讲，涉及一种基于分段统计的大气气溶胶反演 法。
背景技术
大气气溶胶是指大气与悬浮在其中的固体和液体微粒共同组成的多相体系，是指悬浮在地 球大气中的具有一定稳定性的，沉降速度小的，尺度范围在微米到几十微米之间的分子团、液 态或固态粒子所组成的混合物。气溶胶粒子主要来源于工业活动、生物燃烧等人为源以及沙尘、 近海海洋粒子等自然源。
研究表明，大气气溶胶不仅影响地球表层系统的辐射收支平衡，而且通过对太阳短波辐射 的散射调节地气系统的反射率从而导致地气系统的降温过程，通过对太阳辐射的吸收作用导致 升温过程，通过对云的形成微物理过程的影响，可以改变云的微物理性质，对大气化学过程以 及生物地球化学循环也起着重要的作用。同时，气溶胶粒子中包含一些对人体有害的粒子，对 人类的健康状况造成严重的影响，尤其是对于人口密集、工业集中的城市地区。所以，大气气 溶胶的反演对于全球气候变化的研究和大气污染的监测和治理工作具有重大意义。
目前国内外学者在气溶胶遥感反演方面取得了一定的成果，并且在新型传感器的设计中也 越来越考虑到对气溶胶的监测。但是目前大多数的反演算法都是针对大尺度平均性质的气溶胶， 而且具有很多局限性，例如，只适用于暗地表。目前还没有一种有效的适用于城市等亮地表地 区的气溶胶反演算法。
发明内容
本发明的目的在于克服现有技术的不足，提供一种基于分段统计的大气气溶胶反演法，采 用分段统计的方式划分为两个部分，采用不同方式进行气溶胶厚度值反演，可以提高亮地表地 区反演结果的精度和分辨率。
为实现上述发明目的，一种基于分段统计的大气气溶胶反演方法，该方法包括以下步骤：
S1：通过卫星获取所分析区域的多波段遥感影像，并进行辐射定标得到含有表观反射率的 多波段遥感影像，所述多波段遥感影像中包括中红外1.6微米波段、中红外2.1微米波段，从 所述多波段遥感影像中选取一个波段作为反演波段；
S2：确定大气辐射传输模型；根据反演波段构建对应的气溶胶光学厚度查找表，查找表中 包含每个气溶胶厚度值AOT所对应的用于计算表观反射率的参数；
S3：对于多波段遥感影像中的所有像素，采用分段统计方法进行分段，分段方法包括以下 步骤：
S3.1：将中红外2.1微米波段的表观反射率范围按照预设间隔λ1划分为N个区间
S3.2：根据多波段遥感影像中每个像元在中红外2.1微米波段的表观反射率，将各个像元 划分到对应的区间，得到每个区间对应的像元集合xn；
S3.3：对于步骤S3.2得到的每个像元集合xn，如果像元集合xn中的像元数量|xn|≥T1，T1表示预设阈值，则保留该像元集合，否则删除；记筛选得到的像元集合为xn′，n′＝1,2,…,N′， N′表示筛选得到的像元集合数量；
S3.4：将中红外1.6微米波段的表观反射率范围按照预设间隔λ2划分为M个区间γm， m＝1,2,…,M；
S3.5：对每个像元集合xn′，根据其每个像元在中红外1.6微米波段的表观反射率，将各 个像元划分到对应的区间，得到该像元集合下每个区间γm对应的像元集合yn′,m；
S3.6：如果像元集合yn′,m中的像元数量|yn′,m|≥T2，T2表示预设阈值，则保留该像元集合， 否则删除；记筛选得到的像元集合为yk，k＝1,2,…,K，K表示筛选得到的像元集合数量；
S3.7：对于步骤S3.6中筛选得到的K个像元集合yk，根据预设阈值T3进行分类，如果像元集合yk中的像元数量|yk|≥T3，T3表示预设阈值，则将该像元集合yk划入集合Ys中，否则划入集合Yt中；
S4：记集合Ys中第p个像元集合为ys,p，p＝1,2,…,|Ys|，|Ys|表示集合Ys中的像元集合数 量，对每个像元集合ys,p依次进行反演，具体方法包括以下步骤：
S4.1：将反演波段的表观反射率范围按照预设间隔λ3划分为D个区间ωd，d＝1,2,…,D； 对于像元集合ys,p，根据其每个像元在反演波段的表观反射率，将各个像元划分到对应的区间， 得到该像元集合下每个区间ωd对应的像元集合zs,p,d；
S4.2：从D个像元集合zs,p,d中搜索得到像元数量|zs,p,d|≥T4的像元集合，T4表示预设阈 值，从搜索得到的像元集合中选择区间序号d最小的像元集合d*作为像元集合ys,p中的清洁 区段；
S4.3：令清洁区段在反演波段的表观反射率其对应的气溶胶厚度值 AOT0为预设的清洁区段气溶胶厚度值；在气溶胶光学厚度查找表中查找 得到对应的参数，根据表观反射率计算清洁区段对应的地表反射率；
S4.4：将清洁区段对应的地表反射率作为整个像元集合ys,p的地表反射率，在步骤 2确定的大气辐射传输模型框架下反演得到像元集合ys,p中每个像元的气溶胶厚度值；
S5：对集合Yt进行反演，其具体方法为：
S5.1：遍历集合Ys中的每个像元，搜索预设半径内的无值像元，令无值像元的气溶胶厚 度值设置为该像元的气溶胶厚度值；记所有被赋值的像元集合为Ps；
S5.2：对于集合Ps中每个像元，在气溶胶光学厚度查找表中查找得到所填充气溶胶厚度 值对应的参数，然后根据表观反射率计算得到该像元的地表反射率；
S5.3：按照预设的边长将多波段遥感影像划分为网格；
S5.4：记集合Yt中第q个像元集合为yt,q，q＝1,2,…,|Yt|，|Yt|表示集合Yt中的像元集合 数量；对于每个像元集合yt,q，在步骤S5.3划分出的每个网格中求取与Ps的交集，令像元集 合yt,q在当前网格中所有像元的地表反射率等于该交集中所有像元的地表反射率均值；
S5.5：在步骤2确定的大气辐射传输模型框架下对步骤S5.4被赋予地表反射率的每个像 元进行气溶胶厚度值反演，得到该像元的气溶胶厚度值。
进一步的，所述步骤2中采用6S辐射传输模型；
表观反射率ρTOA的计算公式为：
其中，ρ0表示大气层辐射反射率，ρs为地表反射率，S表示大气底层向下的半球反射率，Ts表示入射光线从大气顶层到达地表面的总的透过率，Tv表示向上进入卫星传感器视场方向 的总透过率；
构建的气溶胶光学厚度查找表包括：AOT、ρ0、S、Ts·Tv之间的关系。
进一步的，所述步骤S3.1中间隔λ1的取值范围为0.002≤λ1≤0.01。
进一步的，所述步骤S3.4中间隔λ2的取值范围为0.002≤λ2≤0.01。
进一步的，所述步骤S4.2中阈值T4的取值范围为0.02|Ys|≤T4≤0.1|Ys|。
进一步的，所述步骤S4.3中间隔λ3的取值范围为0.002≤λ3≤0.01。
进一步的，所述步骤S4.4和步骤S5.5中气溶胶厚度值的反演方法为：
从查找表中查找得到各个气溶胶厚度值AOT所对应的参数，根据得到的地表反射率计算 每个气溶胶厚度值AOT对应的表观反射率ρi，其中i＝1,2,…,L，L表示查找表中气溶胶厚 度值的个数；对于每个像元，从L个表观反射率ρi中搜索与其自身表观反射率最接近的表观 反射率，其对应的气溶胶厚度值AOT即为该像元的气溶胶厚度值。
本发明基于分段统计的大气气溶胶反演法，首先得到多波段遥感影像和反演波段对应的气 溶胶光学厚度查找表，然后按照中红外2.1微米波段的表观反射率区间对像元进行划分和筛选， 然后对得到的像元集合按照中红外1.6微米波段的表观反射率区间进行进一步划分和筛选，将 最终得到的像元集合按照像元数量划分为两类，像元较多的划为一类，其他为另一类；将像元 较多的一类作为基准部分进行反演，采用方法是先从像元集合中搜索出清洁区段，以清洁区段 的地表反射率作为整个像元集合的地表反射率，反演得到气溶胶厚度值，然后以这些像元作为 基准，对另一类进行反演。
本发明在进行地表反射率的计算时采用了一种新方法，该方法在计算地表反射率时，是依 据清洁像元和基准区段来确定地表反射率的，并不依赖于暗像元，因此只要满足完成分段统计 即可进行反演，所以较传统的暗像元算法，本发明对亮地表地区也具有良好的反演效果。具有 更广泛的适用性。
附图说明
图1是本发明基于分段统计的大气气溶胶反演法的具体实施方式流程图；
图2是本实施例所采用的影像图；
图3是本发明中像素分段的流程图；
图4是本实施例中波段7的像元划分示例图；
图5是地表反射率波谱曲线示例图；
图6是对集合Ys中每个像元集合ys,p进行反演的流程图；
图7是集合Yt进行反演的流程图；
图8是图2所示影像的网格化结果；
图9是采用本发明对图2所示多波段遥感图像的反演结果。
具体实施方式
下面结合附图对本发明的具体实施方式进行描述，以便本领域的技术人员更好地理解本发 明。需要特别提醒注意的是，在以下的描述中，当已知功能和设计的详细描述也许会淡化本发 明的主要内容时，这些描述在这里将被忽略。
实施例
图1是本发明基于分段统计的大气气溶胶反演法的具体实施方式流程图。如图1所示，本 发明基于分段统计的大气气溶胶反演法包括以下步骤：
S101：获取多波段遥感影像：
通过卫星获取所分析区域的多波段遥感影像，并进行辐射定标得到含有表观反射率的多波 段遥感影像，多波段遥感影像中包括中红外1.6微米波段、中红外2.1微米波段，从多波段遥 感影像中选取一个波段作为反演波段。
目前获取多波段遥感影像的卫星为LANDSAT卫星，在LANDSAT 8OLI影像中，包含了9 个波段，分别是Band1海蓝波段、Band2蓝绿波段、Band3绿波段、Band4红波段、Band5近红外、Band6中红外1.6微米波段、Band7中红外2.1微米波段、Band8全色、Band9卷云、Band10热红外、Band11热红外，本发明中需要使用的是Band6中红外和Band7中红外。本实施例中所选择的反演波段为Band1蓝色波段，这是因为LANDSAT 8OLI影像中的Band1对气溶胶更敏感，可以提高反演的精度。
本实施例中所采用的是北京地区2014年第103天的影像。图2是本实施例所采用的影像 图。如图2所示，该影像图中包含了亮地表区域(图中颜色较浅区域)。
S102：获取气溶胶光学厚度查找表：
将步骤S101所获得的多波段遥感影像输入大气辐射传输模型，根据反演波段得到对应的 气溶胶光学厚度查找表，查找表中包含每个气溶胶厚度值AOT所对应的用于计算表观反射率 的参数。本实施例中，假设下垫面是朗伯体表面，其表观反射率ρTOA的计算公式为：
其中，ρ0表示大气层辐射反射率，ρs为地表反射率，S表示大气底层向下的半球反射率， Ts表示入射光线从大气顶层到达地表面的总的透过率，Tv表示向上进入卫星传感器视场方向 的总透过率。因此可知，本实施例中查找表中需要包含的参数包括大气程辐射反射率、大气半 球反射率、入射光线从大气顶层到达地表面的总的透过率和向上进入卫星传感器视场方向的总 透过率的乘积。
目前行业内已经开发出多种大气辐射传输模型，本实施例中采用6S辐射传输模型。除了 多波段遥感影像外，其他需要输入6S辐射传输模型的参数包括：多波段遥感影像所对应的日 期、天顶角、方位角，气溶胶类型和大气模式，气溶胶厚度值的步长等等。本实施例中选取的 大气模式是中纬度夏季，气溶胶类型是大陆型，所选择的反演波段为Band1海蓝波段。表1 是本实施例中获取的查找表。
AOT ρ<sub>0</sub> S T<sub>s</sub>·T<sub>v</sub>
0.01 0.09094 0.75823 0.174
0.02 0.09178 0.75323 0.17572
… … … …
2.99 0.24386 0.07449 0.306720
S103：像元分段统计：
对于多波段遥感影像中的所有像素，采用分段统计方法进行分段。因为中红外1.6微米波 段、Band7中红外2.1微米波段具有比较长的波长，所以它们的表观反射率可以近似等于地表 反射率，因此本发明采用中红外1.6微米波段、中红外2.1微米波段的表观反射率来对像素进 行分段。由于本实施例中采用的LANDSAT 8 OLI影像，因此中红外1.6微米波段为波段6，中 红外2.1微米波段为波段7。
图3是本发明中像素分段的流程图。如图3所示，本发明中像素分段的具体步骤为：
S301：划分波段7表观反射率区间：
将波段7的表观反射率范围按照预设间隔λ1划分为N个区间划分间 隔λ1可以根据实际需要来确定，一般来说，间隔λ1的取值范围为0.002≤λ1≤0.01。
S302：像元划分：
根据多波段遥感影像中每个像元在波段7的表观反射率，将各个像元划分到对应的区间， 得到每个区间对应的像元集合xn。记某个像元在波段7的表观反射率为ρ，那么其对应的 区间序号即为 表示向上取整。
图4是本实施例中波段7的部分像元以0.005的间隔划分示例图。图4中，图(a)是间隔 划分的整体图，图(b)是间隔划分的局部图。如图4所示，通过对波段7的表观反射率范围划 分区间，从而可以对多波段遥感影像中的像元进行区段划分。
S303：像元集合筛选：
依次对步骤S302得到的N个像元集合xn进行判定，如果像元集合xn中的像元数量 |xn|≥T1，T1表示预设阈值，则保留该像元集合，否则删除。记筛选得到的像元集合为xn′， n′＝1,2,…,N′，N′表示筛选得到的像元集合数量。阈值T1的数值也可以根据实际需要来设置，显然，当多波段遥感影像的分辨率越高，间隔λ1越大，每个区间内所包含的像元数量就越多，相应地阈值T1也可以设置得略高一些，反之阈值T1可以设置得略低一些。本实施例中T1＝100000。
S304：划分波段6表观反射率区间：
将波段6的表观反射率范围按照预设间隔λ2划分为M个区间γm，m＝1,2,…,M。
同样地，划分间隔λ2可以根据实际需要来确定，一般来说，预设间隔λ1的取值范围为：0.002≤λ1≤0.01。
S305：像元集合细分：
对步骤S303筛选得到的每个像元集合xn′分别进行细分，细分方法为：对每个像元集合xn′， 根据其每个像元在波段6的表观反射率，将各个像元划分到对应的区间，得到该像元集合下每 个区间γm对应的像元集合yn′,m。显然，对像元集合进行细分后，共计可以得到N′×M个像 元集合。
S306：像元集合再筛选：
对于步骤S305得到的N′×M个像元集合进行筛选。如果像元集合yn′,m中的像元数量|yn′,m|≥T2，T2表示预设阈值，则保留该像元集合，否则删除。记筛选得到的像元集合为yk， k＝1,2,…,K，K表示筛选得到的像元集合数量。与阈值T1一样地，阈值T2也是根据实际情 况来设置的，由于像元集合yn′,m是细分后的像元集合，显然像元集合yn′,m中的像元数量要更 少，因此T2也应当小于T1。本实施例中T2＝10000
S307：像元集合分类：
对于步骤S306中筛选得到的K个像元集合yk，根据预设阈值T3进行分类，如果像元集 合yk中的像元数量|yk|≥T3，T3表示预设阈值，则将该像元集合yk划入集合Ys中，否则划入 集合Yt中。相应地，阈值T3也是根据实际需要进行设置的。本实施例中，T3的取值范围为100000＜T3＜200000，本实例选择T3＝200000，显然，集合Ys和集合Yt包含了许多个小 的像元集合。
根据以上步骤就将多波段遥感影像中的像元按照波段6和波段7的表观反射率进行了分段。 图5是地表反射率波谱曲线示例图。图5中展示了ENVI标准波谱库中一些地物的地表反射率 波谱曲线，通过研究大量地物的波谱曲线可以发现同种地物或不同地物在不同的条件下具有不 同的地表反射率曲线，它们可能在中红外1.6微米波段和中红外2.1微米波段同时有重合，但 是这种情况很少，尤其是对于常见的地物。而且可以通过后续的两种对分段结果不同的处理方 法来减小误差，同时也可以通过规定每个区段的最小像元数量来减小误差。因此，尽管本发明 中所采用的分段方法存在误差，但是仍旧可以利用中红外1.6微米波段和中红外2.1微米波段 将多波段遥感影像中的像元分割成许多区段。如果这些区段没有误差的存在，那么属于同一区 段的像元应该具有近似相同的地表反射率曲线，它们在反演波段的地表反射率应该近似相等。
S104：对集合Ys进行反演：
本发明将集合Ys作为基准部分，采用一种基于清洁像元反演污染像元的方法进行处理。
在气溶胶的反演中，地表反射率的获得至关重要。实验证明基准集合Ys中的每个像元区 段存在很少的误差，这主要因为基准集合Ys中的每个像元集合都是同一种常见地物，例如， 植被和裸土等，同时还可以通过设置一个大的阈值可以进一步排除一些误差。而且Ys中每个 像元集合含有清洁像元的概率也比较大。记Ys中第p个像元集合为ys,p，p＝1,2,…,|Ys|，|Ys| 表示集合Ys中的像元集合数量。那么根据分段方法的特点，可以得到其中，ρs,p表 示每个像元集合ys,p在波段1的近似地表反射率，表示像元集合ys,p中清洁像元在反演波 段(本实施例为波段1)的地表反射率。因此，只要求得每个像元集合中清洁像元的地表反射 率即可得到像元集合中每个像元的地表反射率。因此，首先要得到每个像元集合中清洁像元的 地表反射率。
图6是对集合Ys中每个像元集合ys,p进行反演的流程图。如图6所示，像元集合ys,p进行 反演的具体方法为：
S601：像元集合细分：
将波段1的表观反射率范围按照预设间隔λ3划分为D个区间ωd，d＝1,2,…,D。同样地，划分间隔λ3可以根据实际需要来确定，间隔λ3的取值范围为0.005≤λ3≤0.02。对于像元集合ys,p，根据其每个像元在波段1的表观反射率，将各个像元划分到对应的区间，得到该像元集合下每个区间ωd对应的像元集合zs,p,d。
S602：搜索清洁区段：
从D个像元集合zs,p,d中搜索得到像元数量|zs,p,d|≥T4的像元集合，T4表示预设阈值，从 搜索得到的像元集合中选择区间序号d最小的像元集合d*作为像元集合ys,p中的清洁区段。
阈值T4也是根据实际情况来设置的，一般来说，阈值T4的取值范围为 0.02|Ys|≤T4≤0.1|Ys|，可以根据多波段遥感影像的污染范围进行选择，污染范围越大该值可 以设置的越小，本实例选择T4＝0.03|Ys|。
S603：计算地表反射率：
令清洁区段在波段1的表观反射率其对应的气溶胶厚度值AOT0为预设的清洁区段气溶胶厚度值。显然，AOT0应当是一个极小 值，是预先根据需要来设置的。在气溶胶光学厚度查找表中查找得到对应的参数，计算清洁区段对应的地表反射率。显然，本实施例中地表反射率是通过表观反射率计算 公式反向求解得到的。
S604：对集合ys,p进行反演：
根据之前的分析可知，清洁区段对应的地表反射率就是整个像元集合ys,p的地表反 射率，因此就能以此来反演得到像元集合ys,p中每个像元的气溶胶厚度值。反演方法可以根据 需要来进行选择，本实施例中采用的反演方法为：
从查找表中查找得到各个气溶胶厚度值AOT所对应的参数，计算每个气溶胶厚度值 AOT对应的表观反射率ρi，其中i＝1,2,…,L，L表示查找表中气溶胶厚度值的个数。
对于每个像元，从L个表观反射率ρi中搜索与其自身表观反射率最接近的表观反射率， 其对应的气溶胶厚度值AOT即为该像元的气溶胶厚度值。
S105：对集合Yt进行反演：
由于步骤S104中对集合Ys中的像元进行了反演，得到了各个像元的气溶胶厚度值，对于 集合Yt中的像元，以集合Ys中像元为基准，进行反演。图7是集合Yt进行反演的流程图。如 图7所示，集合Yt进行反演的具体步骤为：
S701：基准像元扩展：
遍历集合Ys中的每个像元，搜索预设半径内的无值像元，令无值像元的气溶胶厚度值设 置为该像元的气溶胶厚度值；记所有被填充赋值的像元集合为Ps。搜索半径r一般根据需要进 行设置，搜索半径越小，所填充的气溶胶厚度值越准确，但是后续剩余的无值像元数量会增多， 一般其取值范围为1≤r≤10。
S702：计算集合Ps中像元的地表反射率：
对于集合Ps中每个像元，在气溶胶光学厚度查找表中查找得到所填充气溶胶厚度值对应 的参数，然后根据表观反射率计算得到该像元的地表反射率。
S703：多波段遥感图像网格化：
考虑到地物在地理空间上具有一定的相关性，按照预设的边长将多波段遥感影像划分为网 格。预设边长可以根据多波段遥感图像的分辨率来设置。图8是图2所示影像的网格化结果。
S704：获取地表反射率：
记集合Yt中第q个像元集合为yt,q，q＝1,2,…,|Yt|，|Yt|表示集合Yt中的像元集合数量。
对于每个像元集合yt,q，在步骤S703划分出的每个网格中求取与Ps的交集，令像元集合yt,q在 当前网格中所有像元的地表反射率等于该交集中所有像元的地表反射率均值。由此，就可以得 到集合Yt中像元的地表反射率。
S705：气溶胶厚度值反演：
对步骤S705被赋予地表反射率的每个像元进行气溶胶厚度值反演，得到该像元的气溶胶 厚度值。
S106：无值像元插值。
由于本发明在步骤S103中对像元集合进行了筛选，因此集合Ys和Yt中并未包含多波段遥 感图像的所有像元，在步骤S105中也难以实现对所有像元的反演。对于剩余的无值像元，采 用插值方法补充其气溶胶厚度值。插值方法是一种图像处理领域的常用方法，其具体过程在此 再赘述。图9是采用本发明对图2所示多波段遥感图像的反演结果。如图9所示的结果可知， 采用本发明对亮地表地区也具有良好的反演效果，较传统的暗像元算法具有更广泛的适用性。
为了验证其准确性，在http://aeronet.gsfc.nasa.gov/这个网站上下载了AERONET (Aerosol Robotic Network)地面观测网络的数据。将本发明的反演结果和地面观测数据以及 暗像元算法的反演结果进行比较，所选择的两个站点均位于城市中，即属于亮地表区域。表2 是图2所示影像图本发明反演结果、暗像元算法反演结果与地面观测数据对比。
站点 Beijing Beijing-RADI
观测数据 1.351 1.389
暗像元法反演结果 1.66 1.73
本发明反演结果 1.38 1.42
暗像元法误差 0.309 0.341
本发明误差 0.029 0.031
表2
如表2所示，本发明的反演结果与实际站点的观测数据更为接近，其误差更小。可见，对 于亮地表地区，本发明具有良好的反演效果。因此与暗像元法相比，本发明的适用范围更广。
尽管上面对本发明说明性的具体实施方式进行了描述，以便于本技术领域的技术人员理解 本发明，但应该清楚，本发明不限于具体实施方式的范围，对本技术领域的普通技术人员来讲， 只要各种变化在所附的权利要求限定和确定的本发明的精神和范围内，这些变化是显而易见的， 一切利用本发明构思的发明创造均在保护之列。

Claims (7)

1.一种基于分段统计的大气气溶胶反演方法，该方法包括以下步骤：
S1：通过卫星获取所分析区域的多波段遥感影像，并进行辐射定标得到含有表观反射率的多波段遥感影像，所述多波段遥感影像中包括中红外1.6微米波段、中红外2.1微米波段，从所述多波段遥感影像中选取一个波段作为反演波段；
S2：确定大气辐射传输模型；根据反演波段构建对应的气溶胶光学厚度查找表，查找表中包含每个气溶胶厚度值AOT所对应的用于计算表观反射率的参数；
S3：对于多波段遥感影像中的所有像素，采用分段统计方法进行分段，分段方法包括以下步骤：
S3.1：将中红外2.1微米波段的表观反射率范围按照预设间隔λ1划分为N个区间n＝1,2,…,N；
S3.2：根据多波段遥感影像中每个像元在中红外2.1微米波段的表观反射率，将各个像元划分到对应的区间，得到每个区间对应的像元集合xn；
S3.3：对于步骤S3.2得到的每个像元集合xn，如果像元集合xn中的像元数量|xn|≥T1，T1表示预设阈值，则保留该像元集合，否则删除；记筛选得到的像元集合为xn′，n′＝1,2,…,N′，N′表示筛选得到的像元集合数量；
S3.4：将中红外1.6微米波段的表观反射率范围按照预设间隔λ2划分为M个区间γm，m＝1,2,…,M；
S3.5：对每个像元集合xn′，根据其每个像元在中红外1.6微米波段的表观反射率，将各个像元划分到对应的区间，得到该像元集合下每个区间γm对应的像元集合yn′,m；
S3.6：如果像元集合yn′,m中的像元数量|yn′,m|≥T2，T2表示预设阈值，则保留该像元集合，否则删除；记筛选得到的像元集合为yk，k＝1,2,…,K，K表示筛选得到的像元集合数量；
S3.7：对于步骤S3.6中筛选得到的K个像元集合yk，根据预设阈值T3进行分类，如果像元集合yk中的像元数量|yk|≥T3，T3表示预设阈值，则将该像元集合yk划入集合Ys中，否则划入集合Yt中；
S4：记集合Ys中第p个像元集合为ys,p，p＝1,2,…,|Ys|，|Ys|表示集合Ys中的像元集合数量，对每个像元集合ys,p依次进行反演，具体方法包括以下步骤：
S4.1：将反演波段的表观反射率范围按照预设间隔λ3划分为D个区间ωd，d＝1,2,…,D；对于像元集合ys,p，根据其每个像元在反演波段的表观反射率，将各个像元划分到对应的区间，得到该像元集合下每个区间ωd对应的像元集合zs,p,d；
S4.2：从D个像元集合zs,p,d中搜索得到像元数量|zs,p,d|≥T4的像元集合，T4表示预设阈值，从搜索得到的像元集合中选择区间序号d最小的像元集合d*作为像元集合ys,p中的清洁区段；
S4.3：令清洁区段在反演波段的表观反射率其对应的气溶胶厚度值AOT0为预设的清洁区段气溶胶厚度值；在气溶胶光学厚度查找表中查找得到对应的参数，根据表观反射率计算清洁区段对应的地表反射率；
S4.4：将清洁区段对应的地表反射率作为整个像元集合ys,p的地表反射率，在步骤2确定的大气辐射传输模型框架下反演得到像元集合ys,p中每个像元的气溶胶厚度值；
S5：对集合Yt进行反演，其具体方法为：
S5.1：遍历集合Ys中的每个像元，搜索预设半径内的无值像元，令无值像元的气溶胶厚度值设置为该像元的气溶胶厚度值；记所有被赋值的像元集合为Ps；
S5.2：对于集合Ps中每个像元，在气溶胶光学厚度查找表中查找得到所填充气溶胶厚度值对应的参数，然后根据表观反射率计算得到该像元的地表反射率；
S5.3：按照预设的边长将多波段遥感影像划分为网格；
S5.4：记集合Yt中第q个像元集合为yt,q，q＝1,2,…,|Yt|，|Yt|表示集合Yt中的像元集合数量；对于每个像元集合yt,q，在步骤S5.3划分出的每个网格中求取与Ps的交集，令像元集合yt,q在当前网格中所有像元的地表反射率等于该交集中所有像元的地表反射率均值；
S5.5：在步骤2确定的大气辐射传输模型框架下对步骤S5.4被赋予地表反射率的每个像元进行气溶胶厚度值反演，得到该像元的气溶胶厚度值。
2.如权利要求1所述的一种基于分段统计的大气气溶胶反演方法,其特征在于所述步骤2中采用6S辐射传输模型；
表观反射率ρTOA的计算公式为：
其中，ρ0表示大气层辐射反射率，ρs为地表反射率，S表示大气底层向下的半球反射率，Ts表示入射光线从大气顶层到达地表面的总的透过率，Tv表示向上进入卫星传感器视场方向的总透过率；
构建的气溶胶光学厚度查找表包括：AOT、ρ0、S、Ts·Tv之间的关系。
3.如权利要求1所述的一种基于分段统计的大气气溶胶反演方法,其特征在于所述步骤S3.1中间隔λ1的取值范围为0.002≤λ1≤0.01。
4.如权利要求1所述的一种基于分段统计的大气气溶胶反演方法,其特征在于所述步骤S3.4中间隔λ2的取值范围为0.002≤λ2≤0.01。
5.如权利要求1所述的一种基于分段统计的大气气溶胶反演方法,其特征在于所述步骤S4.2中阈值T4的取值范围为0.02|Ys|≤T4≤0.1|Ys|。
6.如权利要求1所述的一种基于分段统计的大气气溶胶反演方法,其特征在于所述步骤S4.3中间隔λ3的取值范围为0.002≤λ3≤0.01。
7.如权利要求1所述的一种基于分段统计的大气气溶胶反演方法,其特征在于所述步骤S4.4和步骤S5.5中气溶胶厚度值的反演方法为：
从查找表中查找得到各个气溶胶厚度值AOT所对应的参数，根据得到的地表反射率计算每个气溶胶厚度值AOT对应的表观反射率ρi，其中i＝1,2,…,L，L表示查找表中气溶胶厚度值的个数；对于每个像元，从L个表观反射率ρi中搜索与其自身表观反射率最接近的表观反射率，其对应的气溶胶厚度值AOT即为该像元的气溶胶厚度值。
CN201910454155.5A 2019-05-29 2019-05-29 一种基于分段统计的大气气溶胶反演方法 Pending CN110261873A (zh)

Priority Applications (3)

Application Number Priority Date Filing Date Title
CN201910454155.5A CN110261873A (zh) 2019-05-29 2019-05-29 一种基于分段统计的大气气溶胶反演方法
US16/885,021 US20210318253A1 (en) 2019-05-29 2020-05-27 Method for retrieving atmospheric aerosol based on statistical segmentation
US17/489,300 US11733175B2 (en) 2019-05-29 2021-09-29 Atmospheric aerosol inversion based on statistical and OTSU segmentation

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201910454155.5A CN110261873A (zh) 2019-05-29 2019-05-29 一种基于分段统计的大气气溶胶反演方法

Publications (1)

Publication Number Publication Date
CN110261873A true CN110261873A (zh) 2019-09-20

Family

ID=67915784

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201910454155.5A Pending CN110261873A (zh) 2019-05-29 2019-05-29 一种基于分段统计的大气气溶胶反演方法

Country Status (2)

Country Link
US (2) US20210318253A1 (zh)
CN (1) CN110261873A (zh)

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111323352A (zh) * 2020-04-09 2020-06-23 中南大学 融合细颗粒物浓度数据的区域pm2.5遥感反演模型
CN113740263A (zh) * 2021-03-10 2021-12-03 仲恺农业工程学院 一种气溶胶光学厚度反演方法及大气颗粒物遥感反演方法

Families Citing this family (9)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN115824971A (zh) * 2022-11-17 2023-03-21 中国科学院空天信息创新研究院 实现陆表遥感地表反射率产品质量可追溯的方法及装置
CN116067857B (zh) * 2023-03-10 2023-07-14 南京航天宏图信息技术有限公司 Fy-4a气溶胶光学厚度的反演方法、装置、设备和介质
CN116519557B (zh) * 2023-07-05 2023-09-19 北京数慧时空信息技术有限公司 一种气溶胶光学厚度反演方法
CN116561487B (zh) * 2023-07-10 2023-09-19 中国科学院空天信息创新研究院 气溶胶瞬时短波直接辐射效应反演方法
CN116660106B (zh) * 2023-07-21 2023-10-17 中国科学院空天信息创新研究院 协同星载标量和偏振观测数据的气溶胶参数迭代反演方法
CN117408142B (zh) * 2023-10-11 2024-03-22 中国科学院空天信息创新研究院 一种基于归一化气溶胶指数的气溶胶光学厚度反演方法
CN117110216B (zh) * 2023-10-19 2024-01-30 航天宏图信息技术股份有限公司 气溶胶光学厚度遥感反演方法、装置及电子设备
CN117152636B (zh) * 2023-10-29 2024-03-15 自然资源部第二海洋研究所 一种基于双波段关系的浅海底质反射率遥感监测方法
CN117935063B (zh) * 2024-03-22 2024-06-11 航天宏图信息技术股份有限公司 基于区域相似标准差法rssdm夜间气溶胶光学厚度反演方法、装置及设备

Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102636143A (zh) * 2012-04-20 2012-08-15 中国科学院遥感应用研究所 一种气溶胶光学厚度遥感反演方法
CN104019753A (zh) * 2014-06-17 2014-09-03 杭州电子科技大学 一种基于modis数据反演城市大气气溶胶光学厚度的方法
CN105674901A (zh) * 2016-04-08 2016-06-15 电子科技大学 基于统计学分段的大气气溶胶反演方法

Patent Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102636143A (zh) * 2012-04-20 2012-08-15 中国科学院遥感应用研究所 一种气溶胶光学厚度遥感反演方法
CN104019753A (zh) * 2014-06-17 2014-09-03 杭州电子科技大学 一种基于modis数据反演城市大气气溶胶光学厚度的方法
CN105674901A (zh) * 2016-04-08 2016-06-15 电子科技大学 基于统计学分段的大气气溶胶反演方法

Cited By (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111323352A (zh) * 2020-04-09 2020-06-23 中南大学 融合细颗粒物浓度数据的区域pm2.5遥感反演模型
CN111323352B (zh) * 2020-04-09 2021-08-24 中南大学 融合细颗粒物浓度数据的区域pm2.5遥感反演模型
CN113740263A (zh) * 2021-03-10 2021-12-03 仲恺农业工程学院 一种气溶胶光学厚度反演方法及大气颗粒物遥感反演方法
CN113740263B (zh) * 2021-03-10 2023-06-16 仲恺农业工程学院 一种气溶胶光学厚度反演方法及大气颗粒物遥感反演方法

Also Published As

Publication number Publication date
US20220018787A1 (en) 2022-01-20
US11733175B2 (en) 2023-08-22
US20210318253A1 (en) 2021-10-14

Similar Documents

Publication Publication Date Title
CN105674901B (zh) 基于统计学分段的大气气溶胶反演方法
CN110261873A (zh) 一种基于分段统计的大气气溶胶反演方法
Lu et al. Spatiotemporal variation of surface urban heat islands in relation to land cover composition and configuration: A multi-scale case study of Xi’an, China
Soegaard et al. Towards a spatial CO2 budget of a metropolitan region based on textural image classification and flux measurements
Liu et al. Analysis of the urban heat island effect in Shijiazhuang, China using satellite and airborne data
Saida et al. Assessment of groundwater vulnerability to nitrate based on the optimised DRASTIC models in the GIS environment (case of Sidi Rached Basin, Algeria)
Guo et al. Impacts of climate and land use/cover change on streamflow using SWAT and a separation method for the Xiying River Basin in northwestern China
CN104063718B (zh) 在作物识别和面积估算中选择遥感数据和分类算法的方法
Yao et al. Artificial neural network (ANN) for multi-source PM2. 5 estimation using surface, MODIS, and meteorological data
Yan et al. The challenge of the urban compact form: Three-dimensional index construction and urban land surface temperature impacts
Taheri et al. Energy-based approaches in estimating actual evapotranspiration focusing on land surface temperature: A review of methods, concepts, and challenges
Dimitrov et al. An application of the LCZ approach in surface urban heat island mapping in Sofia, Bulgaria
CN112784395B (zh) 一种预测和模拟河流水体总磷浓度的方法
Wei et al. Spatial-temporal variation in paddy evapotranspiration in subtropical climate regions based on the SEBAL model: A case study of the Ganfu Plain irrigation system, southern China
Li et al. Estimating subpixel surface heat fluxes through applying temperature-sharpening methods to MODIS data
Zhang et al. Estimation of shortwave solar radiation on clear-sky days for a valley glacier with Sentinel-2 time series
Choi et al. Retrieval of hourly PM2. 5 using top-of-atmosphere reflectance from geostationary ocean color imagers I and II
Han et al. CALIOP-based quantification of Central Asian dust transport
Costanzini et al. Identification of SUHI in Urban Areas by Remote Sensing Data and Mitigation Hypothesis through Solar Reflective Materials
Bessagnet et al. An evaluation of the CHIMERE Chemistry Transport Model to simulate dust outbreaks across the Northern Hemisphere in March 2014
Mahdawi et al. Water recharges suitability in Kabul aquifer system within the upper Indus Basin
Ye et al. Extensive evaluation of four satellite precipitation products and their hydrologic applications over the Yarlung Zangbo River
Tang et al. Estimation of ultrahigh resolution pm2. 5 mass concentrations based on mie scattering theory by using landsat8 oli images over pearl river delta
Zhang et al. Construction of urban thermal environment network based on land surface temperature downscaling and local climate zones
Gohil et al. Validation and improvement of the WRF building environment parametrization (BEP) Urban Scheme

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
RJ01 Rejection of invention patent application after publication

Application publication date: 20190920

RJ01 Rejection of invention patent application after publication