一种可见光遥感卫星影像的大气订正方法和系统与流程

本发明属于遥感卫星技术领域，尤其涉及一种可见光遥感卫星影像的大气订正方法和系统。

背景技术：

大气是遥感成像对地面物体辐射探测、收集和传输的必经介质，地面物体辐射在大气传输过程中会与大气发生一系列的相互作用，导致辐射传输特性发生改变，使遥感影像所记录的信息不能真实反应地面物体辐射信息，从而对遥感影像的质量和遥感影像的解译及定量分析产生很大的影响，因此，消除大气对地面物体辐射的影响，恢复地面物体的本来面目，就成了定量遥感必不可少的一步。

目前，可见光遥感卫星影像的大气订正方法主要有基于地面线性回归模型法、基于大气辐射传输模型法和基于遥感影像特征的订正方法。其中，基于地面线性回归模型法需要以不同地物的大量的野外光谱测量为前提，测量成本较高且无法满足大批量定量化的需求；基于大气辐射传输模型法需要获取不同影像成像时的大气参数，但很多区域无法获取大气参数或即使获取了大气参数也多为观测点的大气参数，很难实现大面积区域遥感的定量化需求；基于遥感影像特征的订正方法可以实现大面积大批量的定量化需要，但其通常需要多个遥感影像波段的配合(如暗目标法需要蓝光波段、红光波段和2.1μm的中波红外波段，多角度法需要每个波段成像时的多个角度信息)，无法实现可见光遥感卫星单波段影像或波段不全影像的大气订正，严重影响可见光遥感卫星影像的定量化应用。

技术实现要素：

本发明的技术解决问题：克服现有技术的不足，提供一种可见光遥感卫星影像的大气订正方法和系统，适用于可见光遥感卫星影像任意波段的大气订正。

为了解决上述技术问题，本发明公开了一种可见光遥感卫星影像的大气订正方法，包括：

将可见光遥感卫星的灰度值影像转换为表观辐亮度影像，得到表观反射率影像；

对所述表观反射率影像进行分类，得到分类影像；

根据大气辐射影响范围距离，确定大气稳定区域内的影像像元范围，将确定的大气稳定区域内的影像像元范围作为大气格网间隔；并，按大气格网间隔对所述表观反射率影像进行大气格网划分，得到空间格网划分影像；

对所述分类影像和空间格网划分影像进行空间叠置分析，获得影像中每个像元的分类属性、空间格网属性、像元表观反射率值及相应位置关联；

根据太阳-大气-地表-大气-遥感卫星相互作用的大气辐射传输方程，建立影像大气订正模型；

根据获得的影像中每个像元的分类属性、空间格网属性、像元表观反射率值及相应位置关联，通过影像大气订正模型，进行影像大气订正反演，确定任意地物分类内所有大气格网内的地物目标反射率。

在上述可见光遥感卫星影像的大气订正方法中，根据可见光遥感卫星影像提供的绝对辐射定标系数和绝对辐射定标计算公式，将可见光遥感卫星的灰度值影像转换为表观辐亮度影像，得到表观反射率影像。

在上述可见光遥感卫星影像的大气订正方法中，根据获得的影像中每个像元的分类属性、空间格网属性、像元表观反射率值及相应位置关联，通过影像大气订正模型，进行影像大气订正反演，确定任意地物分类内所有大气格网内的地物目标反射率，包括：

根据获得的影像中每个像元的分类属性、空间格网属性、像元表观反射率值及相应位置关联，采用带有约束函数的粒子群优化算法进行解算，得到大气格网内同类地物的整体平均反射率全局最优解；

根据整体平均反射率全局最优解和相应大气格网内同类地物的个体差异表观反射率，确定任意地物分类内所有大气格网内的地物目标反射率。

在上述可见光遥感卫星影像的大气订正方法中，根据可见光遥感卫星影像提供的绝对辐射定标系数和绝对辐射定标计算公式，将可见光遥感卫星的灰度值影像转换为表观辐亮度影像，得到表观反射率影像，包括：

根据可见光遥感卫星影像提供的绝对辐射定标系数和绝对辐射定标计算公式，将可见光遥感卫星的灰度值影像dn转换为表观辐亮度影像l；

通过如下公式(1)，得到表观反射率影像ρtoa：

其中，表示日-地距离校正因子；

j表示儒略历天数；es表示影像处理波段的大气外太阳光谱辐照度；μs＝cos(θs)，表示影像成像时太阳天顶角θs的余弦。

在上述可见光遥感卫星影像的大气订正方法中，对所述表观反射率影像进行分类，得到分类影像，包括：

采用监督分类法或非监督分类法，对所述表观反射率影像进行分类，得到分类影像，记为class(k)；其中，k∈[1,k]，表示分类数。

在上述可见光遥感卫星影像的大气订正方法中，根据大气辐射影响范围距离，确定大气稳定区域内的影像像元范围，将确定的大气稳定区域内的影像像元范围作为大气格网间隔；并，按大气格网间隔对所述表观反射率影像进行大气格网划分，得到空间格网划分影像，包括：

根据大气辐射影响范围距离，确定大气稳定区域内的影像像元范围dgrid；

通过如下公式(2)，将确定的大气稳定区域内的影像像元范围作为大气格网间隔ngrid：

其中，gsd表示可见光遥感卫星影像的地面采样分辨率；

按大气格网间隔ngrid，对所述表观反射率影像进行大气格网划分，并依顺序将大气格网记为grid(l)；其中，l∈[1,l]，表示第l个格网。

在上述可见光遥感卫星影像的大气订正方法中，根据太阳-大气-地表-大气-遥感卫星相互作用的大气辐射传输方程，建立影像大气订正模型，包括：

根据太阳-大气-地表-大气-遥感卫星相互作用的大气辐射传输方程，根据如下公式(3)，确定每个像元(i,j)的表观反射率计算模型：·

其中，θv为遥感卫星观测天顶角，μv＝cos(θv)；φs为太阳方位角；ρp为大气路径反射率；ρ为目标反射率；ρe为目标周围区域的平均反射率；s为大气球面反照率；τ为大气的光学厚度；td(θs)为太阳-大气-地表即向下方向的大气漫射散射因子；td(θv)为地表-大气-卫星即向上方向的大气漫射散射因子；为向下方向的大气散射透过率；

将同一大气格网内同一类地物的目标反射率ρ表示为大气格网内同类地物的平均反射率与微小量之和ρ＝ρe+δρ，确定：

确定大气格网内同类地物的整体平均表观反射率中，ρp、t1、t2、s是太阳和卫星观测几何状态参数θs、θv、φs、φv及大气光学厚度τ的函数；

对公式(4)进行简化后，得到：

相应的，本发明还公开了一种可见光遥感卫星影像的大气订正系统，包括：

预处理模块，用于将可见光遥感卫星的灰度值影像转换为表观辐亮度影像，得到表观反射率影像；

分类模块，用于对所述表观反射率影像进行分类，得到分类影像；

格网划分模块，用于根据大气辐射影响范围距离，确定大气稳定区域内的影像像元范围，将确定的大气稳定区域内的影像像元范围作为大气格网间隔；并，按大气格网间隔对所述表观反射率影像进行大气格网划分，得到空间格网划分影像；

分析模块，用于对所述分类影像和空间格网划分影像进行空间叠置分析，获得影像中每个像元的分类属性、空间格网属性、像元表观反射率值及相应位置关联；

建模模块，用于根据太阳-大气-地表-大气-遥感卫星相互作用的大气辐射传输方程，建立影像大气订正模型；

反演模块，用于根据获得的影像中每个像元的分类属性、空间格网属性、像元表观反射率值及相应位置关联，通过影像大气订正模型，进行影像大气订正反演，确定任意地物分类内所有大气格网内的地物目标反射率。

本发明具有以下优点：

(1)本发明利用可见光近遥感卫星影像同类地物的光谱特征和大气尺度的空间特征，将同一地物的光谱特征和大气尺度特征相结合，通过同一大气格网内同一类地物中每个像元的个体属性和同一类地物所有像元的整体属性所具有的粒子群体特性，解决了可见光遥感卫星影像任意波段的大气订正问题。

(2)本发明既考虑了可见光遥感卫星影像地面目标所属类别的差异也考虑了所属类别内大气辐射的差异，通过同类地物的光谱相近性与大气辐射区域的相似性，有效地解决了基于遥感影像特征的订正方法的不足。

(3)本发明为可见光遥感卫星影像大气订正提出的方法，兼顾了批量数据处理问题和可实行性，对可见光近遥感卫星影像大气订正具有很强的工程应用意义。

附图说明

图1是本发明实施例中一种可见光遥感卫星影像的大气订正方法的步骤流程图。

具体实施方式

为使本发明的目的、技术方案和优点更加清楚，下面将结合附图对本发明公开的实施方式作进一步详细描述。

参照图1，示出了本发明实施例中一种可见光遥感卫星影像的大气订正方法的步骤流程图。在本实施例中，所述可见光遥感卫星影像的大气订正方法，包括：

步骤101，将可见光遥感卫星的灰度值影像转换为表观辐亮度影像，得到表观反射率影像。

在本实施例中，可以根据可见光遥感卫星影像提供的绝对辐射定标系数和绝对辐射定标计算公式，将可见光遥感卫星的灰度值影像转换为表观辐亮度影像，得到表观反射率影像。

具体的，根据可见光遥感卫星影像提供的绝对辐射定标系数和绝对辐射定标计算公式，将可见光遥感卫星的灰度值影像dn转换为表观辐亮度影像l；然后，通过如下公式(1)，得到表观反射率影像ρtoa：

其中，表示日-地距离校正因子；

j表示儒略历天数；es表示影像处理波段的大气外太阳光谱辐照度；μs＝cos(θs)，表示影像成像时太阳天顶角θs的余弦。

步骤102，对所述表观反射率影像进行分类，得到分类影像。

在本实施例中，可以采用监督分类法或非监督分类法等任一适当的分类方法，对所述表观反射率影像进行分类，得到分类影像，记为class(k)；其中，k∈[1,k]，表示分类数。

步骤103，根据大气辐射影响范围距离，确定大气稳定区域内的影像像元范围，将确定的大气稳定区域内的影像像元范围作为大气格网间隔；并，按大气格网间隔对所述表观反射率影像进行大气格网划分，得到空间格网划分影像。

在本实施例中，在小尺度大气运动中，可认为一定范围的大气相对稳定、大气各种辐射影响基本相同；为此，设置一定的大气辐射影响范围距离，确定大气稳定区域内的影像像元范围即大气格网间隔，按大气格网间隔将步骤101获得的表观反射率影像进行大气格网划分。也即，可以根据大气辐射影响范围距离，确定大气稳定区域内的影像像元范围dgrid；然后，通过如下公式(2)，将确定的大气稳定区域内的影像像元范围作为大气格网间隔ngrid：

其中，gsd表示可见光遥感卫星影像的地面采样分辨率；

按大气格网间隔ngrid，对所述表观反射率影像进行大气格网划分，并依顺序将大气格网记为grid(l)；其中，l∈[1,l]，表示第l个格网。

步骤104，对所述分类影像和空间格网划分影像进行空间叠置分析，获得影像中每个像元的分类属性、空间格网属性、像元表观反射率值及相应位置关联。

在本实施例中，可以将步骤102得到的分类影像和步骤103得到的空间格网划分影像，进行空间叠置分析，获得影像中每个像元(i,j)的分类属性、空间格网属性、像元表观反射率值及相应位置关联(i,j)∈{class(k),grid(l),ρtoa}。

步骤105，根据太阳-大气-地表-大气-遥感卫星相互作用的大气辐射传输方程，建立影像大气订正模型。

在本实施例中，对于地表非均匀、朗伯体的表面，依据太阳-大气-地表-大气-遥感卫星相互作用的大气辐射传输方程，建立每个像元的表观反射率计算模型；因同一类地物的目标反射率在同一景影像中差异很小，再对经步骤103影像空间格网划分后同一大气格网内同一类地物的表观反射率进行建模。

具体的：

对于地表非均匀、朗伯体的表面，根据太阳-大气-地表-大气-遥感卫星相互作用的大气辐射传输方程，按照如下公式(3)确定每个像元(i,j)的表观反射率计算模型：

其中，θv为遥感卫星观测天顶角，μv＝cos(θv)；φs为太阳方位角；ρp为大气路径反射率；ρ为目标反射率；ρe为目标周围区域的平均反射率；s为大气球面反照率；τ为大气的光学厚度；td(θs)为太阳-大气-地表即向下方向的大气漫射散射因子；td(θv)为地表-大气-卫星即向上方向的大气漫射散射因子；为向下方向的大气散射透过率。

因同一类地物的目标反射率ρ在同一景影像中差异很小，且经步骤103影像空间格网划分后，同一大气格网内同一类地物的目标反射率差异将更小，故可将同一大气格网内同一类地物的目标反射率ρ表示为大气格网内同类地物的平均反射率与微小量之和ρ＝ρe+δρ，为此表观反射率可进一步表示为：

由上式可知，同一大气格网内同一类地物的表观反射率ρtoa可以表示为大气格网内同类地物的整体平均表观反射率与大气格网内同类地物的微小个体差异表观反射率δρtoa的线性之和。

大气格网内同类地物的整体平均表观反射率是ρe、ρp、t1、t2、s的函数其中ρp、t1、t2、s是太阳和卫星观测几何状态参数θs、θv、φs、φv及大气光学厚度τ的函数；在可见光遥感卫星成像时，可获得太阳和卫星观测几何状态参数θs、θv、φs、φv，为此在几何状态参数确定的情况下，大气格网内同类地物的整体平均表观反射率可进一步简化为ρe、τ的二维函数

大气格网内同类地物的个体差异表观反射率δρtoa在整体平均表观反射率确定后，也即二维变量ρe、τ确定后，其是目标反射率微小量δρ的一维线性函数δρtoa＝f(δρ)。

步骤106，根据获得的影像中每个像元的分类属性、空间格网属性、像元表观反射率值及相应位置关联，通过影像大气订正模型，进行影像大气订正反演，确定任意地物分类内所有大气格网内的地物目标反射率。

在本实施例中，可以根据获得的影像中每个像元的分类属性、空间格网属性、像元表观反射率值及相应位置关联，采用带有约束函数的粒子群优化算法进行解算，得到大气格网内同类地物的整体平均反射率全局最优解；根据整体平均反射率全局最优解和相应大气格网内同类地物的个体差异表观反射率，确定任意地物分类内所有大气格网内的地物目标反射率。

具体的：

利用同一大气格网内同一类地物中每个像元的个体属性和同一类地物所有像元的整体属性所具有的粒子群体特性，通过带有约束函数的粒子群优化算法，获得大气格网内同类地物的整体平均反射率全局最优解ρe，然后再利用整体平均反射率全局最优解位置(ρe,τ)和相应大气格网内同类地物的个体差异表观反射率δρtoa，即可确定任意地物分类class(k)内所有大气格网grid(l)内的地物目标反射率ρ＝ρe+δρ。

a)大气订正最优化问题

由步骤105分析知，大气格网l∈[1,l]内ngrid×ngrid个同类地物的整体平均表观反射率为大气格网内同类地物平均反射率ρe和大气光学厚度τ的二维函数通过6s、modtrean等大气辐射传输模型可以获得(其中ρe＝0表示表观反射率只由大气路径散射确定、ρe＝1.0表示目标反射率为强反射地物几乎为100％反射即地物反射率的最大值)和(其中τ＝0为纯分子大气、τ＝2.5为沙尘爆发期间或由生物燃烧释放烟尘造成非常浑浊的大气)情况下整体平均表观反射率的预测估计值为此可将大气格网内同类地物的目标反射率ρ的求解转化求解为已知整体平均表观反射率和预测估计整体平均表观反射率的约束最优化问题：

其中，约束条件考虑了大气格网grid(l)内的大气辐射影响与其相邻8-连接领域grid(l±4)内的大气辐射影响相近、大气格网grid(l)内的大气辐射影响与其相邻24-连接领域grid(l±12)内的大气辐射影响相似以及同一类地物目标反射率ρ在大气格网grid(l)相邻8-连接领域grid(l±4)内相近和整体平均反射率ρe在大气格网grid(l)相邻24-连接领域grid(l±4)内相近等情况，其相应约束参数ετ、ερ、εe依据具体可见光遥感卫星的传感器响应函数、太阳和卫星观测几何状态参数以及先验地物分类光谱特性进行设置。

b)大气订正最优解

同一大气格网内同一类地物的表观反射率ρtoa为大气格网内同类地物的整体平均表观反射率与大气格网内同类地物的微小个体差异表观反射率δρtoa之和，也即同一大气格网内同一类地物中每个像元的个体属性与同一类地物所有像元的整体属性具有ngrid×ngrid个粒子组成的粒子群特性，为此可通过带约束的粒子群优化算法，初始化为一群随机粒子(随机解)，通过群体中个体之间的协作和信息共享使整个群体的运动问题在问题求解空间中利用公式(6)迭代产生从无序到有序的演化过程从而获得最优解找到最优解，在每一次迭代中，粒子通过两个极值来更新自己：一个是粒子本身所找到的最优解，这个解即为个体极值另一个极值是整个粒子群目前找到的最优解，这个极值即为全局最优解ρe,best(t)。

xq(t+1)＝xq(t)+vq(t+1)

其中：xq＝(ρe,q,τq)为第q个粒子的当前位置；

为第q个粒子的当前飞行速度；

c1∈[0,2.05]为调节粒子向自身最优位置飞行的步长位置；

c2∈[0,2.05]为调节粒子向全局最优位置飞行的步长位置；

为收缩因子，其中c＝c1+c2＞4。

c)整景影像大气订正处理

对所有地物分类class(k)、k∈[1,k]依此采用步骤b)获取所有大气格网grid(l)、l∈[1,l]内同类地物的整体平均反射率全局最优解ρe；然后再利用整体平均反射率全局最优解位置(ρe,τ)和相应大气格网内同类地物的个体差异表观反射率即可确定任意分类class(k)内所有大气格网grid(l)内的地物目标反射率ρ＝ρe+δρ。

相应的，本发明还公开了一种可见光遥感卫星影像的大气订正系统，包括：预处理模块，用于将可见光遥感卫星的灰度值影像转换为表观辐亮度影像，得到表观反射率影像；分类模块，用于对所述表观反射率影像进行分类，得到分类影像；格网划分模块，用于根据大气辐射影响范围距离，确定大气稳定区域内的影像像元范围，将确定的大气稳定区域内的影像像元范围作为大气格网间隔；并，按大气格网间隔对所述表观反射率影像进行大气格网划分，得到空间格网划分影像；分析模块，用于对所述分类影像和空间格网划分影像进行空间叠置分析，获得影像中每个像元的分类属性、空间格网属性、像元表观反射率值及相应位置关联；建模模块，用于根据太阳-大气-地表-大气-遥感卫星相互作用的大气辐射传输方程，建立影像大气订正模型；反演模块，用于根据获得的影像中每个像元的分类属性、空间格网属性、像元表观反射率值及相应位置关联，通过影像大气订正模型，进行影像大气订正反演，确定任意地物分类内所有大气格网内的地物目标反射率。

对于系统实施例而言，由于其与方法实施例相对应，所以描述的比较简单，相关之处参见方法实施例部分的说明即可。

本说明中的各个实施例均采用递进的方式描述，每个实施例重点说明的都是与其他实施例的不同之处，各个实施例之间相同相似的部分互相参见即可。

以上所述，仅为本发明最佳的具体实施方式，但本发明的保护范围并不局限于此，任何熟悉本技术领域的技术人员在本发明揭露的技术范围内，可轻易想到的变化或替换，都应涵盖在本发明的保护范围之内。

本发明说明书中未作详细描述的内容属于本领域专业技术人员的公知技术。