CN101034477A - 遥感数字图像上阴影消除和阴影中像元遥感值恢复方法 - Google Patents

Abstract

本发明涉及一种遥感数字图像上阴影消除和阴影中像元遥感值恢复的方法。首先，判断和消除遥感数字图像的山体与云层阴影。再次，在此基础上最后完成遥感数字图像上阴影消除和阴影中像元遥感值恢复。其具体步骤是：输入遥感数字图像和对应的数字地形图，寻找阴影并计算参数，在遥感数字图像上寻找分别位于阳光下和阴影(云层与山体)中的同类像元P1、P2并分别计算其归一化天空光地形照射系数：G1和G2以及太阳直射光地形照射系数(归一化)F1；继而计算该阴影区水平地面自然光的散射/直射比，然后计算全部阴影区水平地面自然光的散射/直射比，插值计算所有无阴影地区各点水平地面上自然光的散射/直射比值，最终输出对应水平地面上自然光散射/直射比和对应水平地面上遥感数字图像上阴影消除后的图像和阴影中像元遥感值恢复数值，大大地提高了散射/直射比获取的效率和精度。

Description

遥感数字图像上阴影消除和阴影中像元遥感值恢复方法

技术领域

本发明涉及一种遥感与地理信息系统技术，特别是一种遥感数字图像上阴影消除和阴影中像元遥感值恢复方法。

背景技术

散射/直射比是表征大气状况的重要参数，在大气环境、气象、环境遥感和遥感图像处理等方面有特殊意义。由于地形对天空的遮挡和屏蔽，自然地表上的天空散射光照度分布受地形影响将发生变化。因地形的影响，不仅改变了地面上自然光的总照度并且改变了地表自然光中天空散射光和太阳直射光照度的比例，即：散射/直射比。散射/直射比的常规获取方法是人工逐点实测，费时费力，效率低。判断和消除遥感数字图像的山体与云层阴影和阴影中像元遥感值的恢复一直是遥感图像处理的难题，常规、传统处理方法，尚未解决该难题。

发明内容

本发明的目的在于客服现有技术中存在的难题，提供一种遥感数字图像上阴影消除和阴影中像元遥感值恢复方法。

为了达到上述发明目的，本发明采用下述技术方案：一种遥感数字图像上阴影消除和阴影中像元遥感值恢复方法，其特征在于，首先，判断和消除遥感数字图像的山体与云层阴影；再次，在此基础上最后完成遥感数字图像上阴影消除和阴影中像元遥感值恢复，其具体步骤是：

(1).输入遥感数字图像和对应的数字地形图；

(2).寻找阴影并计算参数：在遥感数字图像上寻找分别位于阳光下和阴影(云层与山体)中的同类像元P1、P2并分别计算其归一化天空光地形照射系数：G1和G2以及太阳直射光地形照射系数(归一化)F1；

(3).计算该阴影区水平地面自然光的散射/直射比：

L＝DN2*F1/(DN1*G2-DN2*G1)

(4).计算全部阴影区水平地面自然光的散射/直射比；

(5).插值计算所有无阴影地区各点水平地面上自然光的散射/直射比值；

(6).输出：对应水平地面上自然光散射/直射比；

(7).输出：对应水平地面上遥感数字图像上阴影消除后的图像和阴影中像元遥感值恢复数值

上述的判断和消除遥感数字图像的山体与云层阴影的方法如下：

设山体和云层阴影中一像元P_kl，其遥感值DN_kl，大气程辐射遥感值D_Akl，像元点水平地面上的散射/直射比L_kl、太阳直射光水平地形变换系数：F_kl、天空散射光地形变换系数：G_kl为已知。

则：阴影中太阳直射光照度为0：DN_skl＝0

像元坡面天空散射光遥感值：DN_dkl＝(DN_kl-D_Akl)                            (1)

像元水平地面天空散射光遥感值：DN`_dkl＝DN_dkl/G_kl＝(DN_kl-D_Akl)/G_kl     (2)

则：当像元P_kl不在阴影中时：

像元水平地面太阳直射光遥感值：DN`<sub>skl</sub>＝DN`_dkl/L_kl＝(DN_kl-D_Akl)/L_klG_kl(3)

像元坡面太阳直射光遥感值：

DN_skl＝DN`_skl；F_kl＝DN_dkl；F_kl/L_klG_kl＝(DN_kl-D_Akl)F_kl/L_klG_kl    (4)

坡面像元遥感值：

DN_kl＝DN_dkl+DN_skll＝DN_dkl+DN_dklF_kl/L_kl G_kl

    ＝DN_dkl(1+F_kl/L_kl G_kl)＝(DN_kl-D_Akl)(1+F_kl/L_kl G_kl)             (5)

像元水平地面遥感值：

DN`<sub>k</sub>＝DN`_skl+DN`_dkl＝(DN_kl-D_Akl)/L_klG_kl+(DN_kl-D_Akl)/G_kl

    ＝(1+L_kl)(DN_kl-D_Akl)/G_klL_kl                                        (6)

公式中：F_kl＝1-tgα_kl·ctgθ_kl·cosω_kl，其中α_kl、θ_kl、ω_kl、分别为：像元P_kl的坡角、太阳高度角和该像元坡向、太阳入射方向之间的夹角。

上述的步骤(2)中，阴影中任意像元点P_ij的归一化天空光地形照射系数G_ij和阳光下像元点P1(P_ij)的归一化太阳直射光地形照射系数F的计算步骤如下：

(1)像元P_ij归一化天空光地形照射系数的计算：

(a)分别计算像元点P_ij在各个方向上的最大地形高度角为β_k(k＝1，2，……n)，

n为方位间隔：n＝2π/Δt，Δt为方位角步长，为可整除360的数；

β_k为点P_ij在第k方向上的最大地形高度(屏蔽)角，

β_k＝MAX(h_l，L＝1、2、3…ML)                                            (7)

h_l为点P_ij在第k个方向上的第L个点的地形高度(屏蔽)角，

h_l＝tg(Z_PL/S_PL)^-1                                                     (8)

Z_PL，S_PL：分别为点P_ij与第k个方向上的第L个点的高差和水平距离。

(b)求和以计算该点的归一化天空光地形屏蔽系数：

${\mathrm{\&Omega;}}_{\mathrm{ij}}=\mathrm{\&Delta;t}\underset{k=1}{\overset{n}{\mathrm{\&Sigma;}}}\sin {\mathrm{\&beta;}}_k/{2\mathrm{\&pi;}}_k---\left(9\right)$

(c)计算像元点的归一化天空光地形照射系数：G_ij＝1-Ω_ij                (10)

(2)阳光下像元点P1(P_ij)的归一化太阳直射光地形照射系数的计算：

F1＝F_ij＝1-tgα_ij·ctgθ_ij·cosω_ij  ω_ij＝AL_ij-A_ij            (11)

像元地面太阳高度角θ_ij、方位角AL_ij由卫片注记给出，像元地面坡度α_ij、坡向A_ij以及像元经、纬度由卫片对应数字地形图生成的数字地面模型DTM给出，ω_ij代表太阳方位角与地面坡向之夹角。

本项发明与现有技术相比所具有的优点与积极效果

水平地面自然光的散射/直射比是表达太阳光受到环境(大气、地形、地理位置)影响程度的重要参数。随着大气(气象)状况、地理位置和地形条件的不同水平地面上各地的自然光散射/直射比一般是不同的。水平地面自然光的散射/直射比的常规获取方法是人工逐点实测。本项发明是利用遥感数字图像和数字地形图，自动获取区域(大范围)水平地面自然光的散射/直射比。大大地提高了散射/直射比自动化智能化获取的效率和精度。

本项发明在大气环境质量评价、遥感图像的地形影响修正、地表光辐射的地形再分配及其区域分布研究等方面均有重要意义。

附图说明

图1为本发明遥感数字图像上阴影消除和阴影中像元遥感值恢复的原理方法的计算机流程图。

图2为福建龙岩地区1∶10万数字地形图(1024×1024)；

图3为中国福建龙岩地区太阳直射光地形分布图像(1024×1024，太阳方向东南)。

图4为福建龙岩地区天空散射光地形分布图像(1024×1024)。

具体实施方式

本发明的一个优选实施例结合附图详细说明如下：

本例为福建龙岩地区1∶10万数字地形图(1024×1024)以说明遥感数字图像上阴影消除和阴影中像元遥感值恢复的原理方法。

参见图2，其出示了本实施所要研究地区，即中国福建龙岩地区1∶10万数字地形图(1024×1024)。

参见图3，其出示了中国福建龙岩地区太阳直射光地形分布图像(1024×1024，太阳方向东南)。

请参阅图1，它是本发明遥感数字图像上阴影消除和阴影中像元遥感值恢复的原理方法的计算机流程图，在GIS中，运行步骤如下：

(1).输入遥感数字图像和对应的数字地形图。

(2).寻找阴影并计算参数：在遥感数字图像上寻找分别位于阳光下和阴影(云层与山体)中的同类像元P1、P2并分别计算其归一化天空光地形照射系数：G1和G2以及太阳直射光地形照射系数(归一化)F1。

(3).计算该阴影区水平地面自然光的散射/直射比：

L＝DN2*F1/(DN1*G2-DN2*G1)

(4).计算全部阴影区水平地面自然光的散射/直射比。

(5).插值计算所有无阴影地区各点水平地面上自然光的散射/直射比值。

(6).输出：对应水平地面上自然光散射/直射比。

(7).输出：对应水平地面上遥感数字图像上阴影消除后的图像和阴影中像元遥感值恢复数值。见图4。

上述的判断和消除遥感数字图像的山体与云层阴影的方法如下：

设山体和云层阴影中一像元P_kl，其遥感值DN_kl，大气程辐射遥感值D_Akl，像元点水平地面上的散射/直射比L_kl、太阳直射光水平地形变换系数：F_kl、天空散射光地形变换系数：G_kl为已知。

则：阴影中太阳直射光照度为0：DN_skl＝0

像元坡面天空散射光遥感值：DN_dkl＝(DN_kl-D_Akl)                            (1)

像元水平地面天空散射光遥感值：DN`_dkl＝DN_dkl/G_kl＝(DN_kl-D_Akl)/G_kl     (2)

则：当像元P_kl不在阴影中时：

像元水平地面太阳直射光遥感值：DN`<sub>skl</sub>＝DN`_dkl/L_kl＝(DN_kl-D_Akl)/L_klG_kl(3)

像元坡面太阳直射光遥感值：

DN_skl＝DN`_skl；F_kl＝DN_dkl；F_kl/L_klG_kl＝(DN_kl-D_Akl)F_kl/L_klG_kl    (4)

坡面像元遥感值：

DN_kl＝DN_dkl+DN_skll＝DN_dkl+DN_dklF_kl/L_kl G_kl

    ＝DN_dkl(1+F_kl/L_kl G_kl)＝(DN_kl-D_Akl)(1+F_kl/L_kl G_kl)           (5)

像元水平地面遥感值：

DN`<sub>k</sub>＝DN`_skl+DN`_dkl＝(DN_kl-D_Akl)/L_klG_kl+(DN_kl-D_Akl)/G_kl

    ＝(1+L_kl)(DN_kl-D_Akl)/G_kl L_kl                                     (6)

公式中：F_kl＝1-tgα_kl·ctgθ_kl·cosω_kl，其中α_kl、θ_kl、ω_kl、分别为：像元P_kl的坡角、太阳高度角和该像元坡向、太阳入射方向之间的夹角。

上述的步骤(2)中，阴影中任意像元点P_ij的归一化天空光地形照射系数G_ij和阳光下像元点P1(P_ij)的归一化太阳直射光地形照射系数F的计算步骤如下：

(1)像元P_ij归一化天空光地形照射系数的计算：

(a)分别计算像元点P_ij在各个方向上的最大地形高度角为β_k(k＝1，2，……n)；

n为方位间隔：n＝2π/Δt，Δt为方位角步长，为可整除360的数；

β_k为点P_ij在第k方向上的最大地形高度(屏蔽)角，

β_k＝MAX(h_l，L＝1、2、3…ML)                                           (7)

h_l为点P_ij在第k个方向上的第L个点的地形高度(屏蔽)角，

h_l＝tg(Z_PL/S_PL)^-1                                                    (8)

Z_PL，S_PL：分别为点P_ij与第k个方向上的第L个点的高差和水平距离。

(b)求和以计算该点的归一化天空光地形屏蔽系数：

${\mathrm{\&Omega;}}_{\mathrm{ij}}=\mathrm{\&Delta;t}\underset{k=1}{\overset{n}{\mathrm{\&Sigma;}}}\sin {\mathrm{\&beta;}}_k/{2\mathrm{\&pi;}}_k---\left(9\right)$

(c)计算像元点的归一化天空光地形照射系数：G_ij＝1-Ω_ij                 (10)

(2)阳光下像元点P1(P_ij)的归一化太阳直射光地形照射系数的计算：

F1＝F_ij＝1-tgα_ij·ctgθ_ij·cosω_ij ω_ij＝AL_ij-A_ij              (11)

像元地面太阳高度角θ_ij、方位角AL_ij由卫片注记给出，像元地面坡度α_ij、坡向A_ij以及像元经、纬度由卫片对应DTM(数字地形图生成的数字地面模型)给出，ω_ij代表太阳方位角与地面坡向之夹角。

实例与讨论

山体和云层阴影中某像元：其遥感值DN_kl＝11，大气程辐射遥感值D_Akl＝4，像元点水平地面上的散射/直射比L_kl＝0。08、太阳直射光水平地形变换系数：F_kl＝0.883、天空散射光地形变换系数：G_kl＝0..735

则：像元坡面天空散射光遥感值DN_dkl＝11-4＝7

像元水平地面天空散射光遥感值：DN`_dkl＝DN_dkl/G_kl＝7/0.735＝9

则：当像元P_kl不在阴影中时：

像元水平地面太阳直射光遥感值：DN`<sub>skl</sub>＝DN`_dkl/L_kl＝9/0.08＝113

像元坡面太阳直射光遥感值：DN_skl＝DN`_skl F_kl＝113×0.881＝100

坡面像元遥感值：DN_kl＝DN_dkl+DN_skll＝7+100＝107

像元水平地面遥感值：DN`<sub>k</sub>＝DN`_skl+DN`_dk l＝9+113＝122

遥感数字图像的水平地形变换(TM：1波段)：

行	列	像元遥感值	坡度(度)	坡向(度)	F值	G值	水平地面像元遥感值
								31	148	83	23.7	83	0.69297	0.82279	116
39	239	165	39.1	190.0	1.55536	0.85220	114
								23	61	113	31.7	153.8	1.56395	0.71841	79
26	96	113	29.4	348.2	0.51337	0.83611	199

从上面可以看出：分别位于阴披和阳坡的同类像元(前2)遥感值经水平地形变换后接近一致；而分别位于阴披和阳坡的非同类像元(后2)遥感值经水平地形变换后显著分离。

Claims (3)

1.遥感数字图像上阴影消除和阴影中像元遥感值恢复方法，其特征在于，首先，判断和消除遥感数字图像的山体与云层阴影；再次，在此基础上最后完成遥感数字图像上阴影消除和阴影中像元遥感值恢复；其具体步骤是：

(1).输入遥感数字图像和对应的数字地形图；

(2).寻找阴影并计算参数：在遥感数字图像上寻找分别位于阳光下和阴影中的同类像元P1、P2并分别计算其归一化天空光地形照射系数：G1和G2以及太阳直射光地形照射系数F1；

(3).计算该阴影区水平地面自然光的散射/直射比：

L＝DN2*F1/(DN1*G2-DN2*G1)

(4).计算全部阴影区水平地面自然光的散射/直射比；

(5).插值计算所有无阴影地区各点水平地面上自然光的散射/直射比值；

(6).输出：对应水平地面上自然光散射/直射比；

(7).输出：对应水平地面上遥感数字图像上阴影消除后的图像和阴影中像元遥感值恢复数值。

2.根据权利要求1所述的遥感数字图像上阴影消除和阴影中像元遥感值恢复方法，其特征在于所述的判断和消除遥感数字图像的山体与云层阴影的方法如下：

设山体和云层阴影中一像元P_kl，其遥感值DN_kl，大气程辐射遥感值D_Akl，像元点水平地面上的散射/直射比L_kl、太阳直射光水平地形变换系数：F_kl、天空散射光地形变换系数：G_kl为已知。

则：阴影中太阳直射光照度为0：DN_skl＝0，

像元坡面天空散射光遥感值：DN_dkl＝(DN_kl-D_Akl)，

像元水平地面天空散射光遥感值：DN`_dkl＝DN_dkl/G_kl＝(DN_kl-D_Akl)/G_kl，

则：当像元P_kl不在阴影中时：

像元水平地面太阳直射光遥感值：DN`<sub>skl</sub>＝DN`_dkl/L_kl＝(DN_kl-D_Akl)/L_klG_kl，

像元坡面太阳直射光遥感值：

DN_skl＝DN`_skl；F_kl＝DN_dkl；F_kl/L_klG_kl＝(DN_kl-D_Akl)F_kl/L_klG_kl，

坡面像元遥感值：

DN_kl＝DN_dkl+DN_skll＝DN_dkl+DN_dklF_kl/L_klG_kl

    ＝DN_dkl(1+F_kl/L_klG_kl)＝(DN_kl-D_Akl)(1+F_kl/L_klG_kl)，

像元水平地面遥感值：

DN`<sub>k</sub>＝DN`_skl+DN`_dkl＝(DN_kl-D_Akl)/L_klG_kl+(DN_kl-D_Akl)/G_kl

        ＝(1+L_kl)(DN_kl-D_Akl)/G_klL_kl，

公式中：F_kl＝1-tgα_kl·ctgθ_kl·cosω_kl，其中α_kl、θ_kl、ω_kl、分别为：像元P_kl的坡角、太阳高度角和该像元坡向、太阳入射方向之间的夹角。

3.根据权利要求1所述的遥感数字图像上阴影消除和阴影中像元遥感值恢复方法，其特征在于所述的步骤(2)中，阴影中任意像元点P_ij的归一化天空光地形照射系数G_ij和阳光下像元点P1(P_ij)的归一化太阳直射光地形照射系数F的计算步骤如下：(1)像元P_ij归一化天空光地形照射系数的计算：

(a)分别计算像元点P_ij在各个方向上的最大地形高度角为β_k(k＝1，2，……n)；

n为方位间隔：n＝2π/Δt，Δt为方位角步长，为可整除360的数；

β_k为点P_ij在第k方向上的最大地形高度(屏蔽)角，

β_k＝MAX(h_l，L＝1、2、3…ML)                        (7)

h_l为点P_ij在第k个方向上的第L个点的地形高度(屏蔽)角，

h_l＝tg(Z_PL/S_PL)^-1                                   (8)

Z_PL，S_PL：分别为点P_ij与第k个方向上的第L个点的高差和水平距离。

(b)求和以计算该点的归一化天空光地形屏蔽系数：

${\mathrm{\&Omega;}}_{\mathrm{ij}}=\mathrm{\&Delta;t}\underset{k=1}{\overset{n}{\mathrm{\&Sigma;}}}\sin {\mathrm{\&beta;}}_k/2{\mathrm{\&pi;}}_k---\left(9\right)$

(c)计算像元点的归一化天空光地形照射系数：G_ij＝1-Ω_ij    (10)

(2)阳光下像元点P1(P_ij)的归一化太阳直射光地形照射系数的计算：

F1＝F_ij＝1-tgα_ij·ctgθ_ij·cosω_ij         ω_ij＝AL_ij-A_ij    (11)

像元地面太阳高度角θ_ij、方位角AL_ij由卫片注记给出，像元地面坡度α_ij、坡向A_ij以及像元经、纬度由卫片对应DTM(数字地形图生成的数字地面模型)给出，ω_ij代表太阳方位角与地面坡向之夹角。

Images