在定量遥感中，影像中的地物光谱信息受大气、地形效应等多种因素影响^[1]，特别是复杂的山区影像，由于地形的起伏，各个像元接收到的有效光照差别较大^[2]，并且影像空间分辨率越高，这种效应就越明显，为了消除这种地形效应的影响，需要进行山区地形辐射校正。对于山区遥感影像，地表接收到的辐射主要有太阳直射辐射、天空的散射辐射、邻近地表的反射和散射4个部分。由于受到坡度、坡向等因素影响，山地阴影区域不能接收到太阳直接辐射，其接收到的辐射主要为天空的散射及邻近地表的反射与散射，因此在对高空间分辨率数据进行地形校正时，以上两项必须加以考虑^[3]。

地形校正算法主要分为3种类型：经验模型、物理模型、半经验模型。其中全色影像的地形校正主要是运用经验模型校正法。经验模型为太阳入射角与卫星传感器接收的辐射之间的经验关系的总结。其优点是所用到的地形因子少，简单易懂且适用性很强，缺点是理论性不完善。波段比模型、Minnaer校正模型都属于该模型。这些模型用于全色影像的地形校正，往往出现过校正或欠校正的效果。物理模型为通过研究光与地表的物理相互作用，基于辐射传输理论建立的模型。其优点是理论基础完善，模型参数物理意义明确，缺点是模型复杂，所需参数较多。这类模型有余弦校正法、SCS校正模型等。半经验模型则是综合了经验模型和物理模型的优点^{[2, 4]}。相比之下，基于辐射传输理论的模型能克服理论性不完善和过校正的缺点，利用定量化信息对影像进行校正^[5]。与此同时，基于辐射传输模型可以在更好地恢复地物像元之间相对关系的基础上对全色遥感影像进行山区地形校正。

利用辐射传输模型校正时，遮蔽因子的精确计算对校正的效果至关重要。通常采用DEM和太阳有效入射角、坡度和坡向角度等地形参量结合的方法计算遮蔽因子，该方法受DEM精度影响较大，且阴影的边缘位置判断不准确，从而会导致地形校正不完整的现象。本文首先用影像头文件的定标系数对全色影像进行定标，然后运用直方图的特殊性计算遮蔽因子并进行地形校正，最后与DEM计算了遮蔽因子得到的地形校正结果作对比，分别从目视效果和统计指标等方面进行对比分析。

1 原理与方法

地形辐射校正的实质是通过各种方法变换，将像元的辐亮度变换到某一相同基准面^[6]。对复杂山区地形校正采用的辐射传输模型需要考虑天空散射辐射和邻近地表反射辐射两个重要影响因素。将地表视为朗伯体，且假定在所有高程下太阳辐射及大气为均匀的情况下，地表反射率公式如下^[7]

式中，d为日地距离；b为遮蔽因子；L_ba为图像程辐射对应的辐亮度，即图像有效像元区域最小值；τ_v、τ_s、E_dir、E_dif分别为大气上行透射率、下行透射率、地面直射辐照度和地面散射辐照度，可直接由6SV计算得到；V_sky(x, y) 为像元对应天空观测因子，由坡度计算得到；E_t为邻近像元的多次反射辐照度；ρ_terrain^(i－1)为环境反射率；β为太阳有效入射角；θ_z为太阳天顶角。太阳有效入射角的计算公式为^[8]

式中，S为坡度；A为坡向, 可通过DEM计算得到；ϕ_S为太阳方位角。以上为本文用来进行地形校正所使用的辐射传输模型。本文中b因子的提取采取两种方法并对其进行分析。

DEM计算遮蔽因子主要是在ENVI的Topographic Modeling-Shaded Relief模块中输入太阳高度角、太阳方位角，通过对仰角的判断遮蔽情况判别阴影区，即可通过波段运算得到值为0-1的遮蔽因子。

直方图计算遮蔽因子主要采用阈值分割法。阈值法是一类被广泛使用的分割方法，它具有操作简便、易于理解并能实现快速分割的优势，其中直方图阈值分割是最为常用的方法。灰度直方图横坐标代表图像的像素灰度值，纵坐标表示具有该灰度级所包含的像素数或出现这个灰度级的概率，直方图曲线可以体现出图像对比度、明亮程度及灰度分布等方面的差异^[9]。当图像中目标物与背景出现明显的谷峰形状时，利用直方图分割可以取得较好的效果^[10]，本文基于该特征通过直方图提取遮蔽因子。

总体来说，图像直方图可以分为单峰直方图、双峰直方图及多峰直方图 3种类型。单峰直方图比较常见，通常是由于目标区域与背景区域面积相差太大，从而像素数量相差悬殊造成的。双峰直方图是选取分割阈值的最理想情况。此类图像明显可以区分两种地物。目标对象与背景区域各自对应一个波峰，两峰之间存在一个波谷，波谷区域的像素量是最少的，同时也代表着两类区域的边界处，可认为谷值为最佳的分割点。多峰直方图则表示该类图像包含的信息非常巨大且结构复杂^[10]。在进行山区遮蔽因子的提取时，对阴影区域的识别尤为重要。直方图阈值法可以运用到山区地形校正中的阴影区域判别与提取。对于裸露的复杂山区遥感影像，地物信息较为单一，山体 (阳坡、阴坡)、阴影边界及阴影是影像的重要组成部分，三者的影像灰度值逐级递减，而且在直方图上表现为明显的双峰直方图特征。两峰之间的波谷对应影像中类别为阴影边界的像元，此类像元在频数上明显少于其余两种地物类别。此时选取双峰之间的谷值作为阈值，通过波段运算将检测出的阴影区域设为0，非阴影区值设为1，便可得到更为精确的遮蔽因子。

2 试验与分析 2.1 试验过程

试验所用数据为SPOT-5 HRG2传感器的全色影像数据，获取日期为2013年6月13日。研究区域位于甘肃省张掖市 (39°38′06″N，100°35′53″E)，海拔1630 m，地形以戈壁及山区为主。SPOT-5的全色数据头文件中获取定标系数C₀为1.320 660(1/W·m²·sr·μm)。试验当天在地面同步测量了大气数据。另外，试验所用DEM数据 (SRTM DEM，分辨率30 m) 下载自地理空间数据云 (http://www.gscloud.cn/)，如图 1所示。

根据上述辐射传输模型，本文选用美国RSI (Research System Inc) 的IDL (interactive data language) 交互式数据语言作为开发平台，同时结合6SV模型计算的辐射分量进行试验。表 1为利用6SV辐射传输模型计算时采用的输入参数。

在计算遮蔽因子时，分别采用DEM和直方图两种方法，并对比其获取遮蔽因子的差异性。利用ENVI的Topographic-Topographic Modeling选取DEM数据，选取shaded relief生成阴影地貌图像。进行波段运算，得到的二值化遮蔽因子如图 2所示。

随后通过直方图的方法计算遮蔽因子，打开全色影像，统计出直方图信息，如图 3所示。

从曲线的形状可以看出该直方图为双峰直方图且有明显的波谷。对波谷进行量取可知影像灰度值小于或大于50的像元均符合正态分布，故选取双峰之间的谷值50作为阈值，通过波段运算计算出遮蔽因子，结果如图 4所示。

2.2 结果分析

分别采用两种方法计算遮蔽因子并利用式 (1) 进行地形校正，得到结果如图 5、图 6所示。

对上述结果分别从目视判读和指标统计两个方面进行分析评价。

2.2.1 目视判读评价

从校正结果可以看出，使用不同方法进行地形校正后，原始影像的阴影均得到了消除，且影像的凹凸感也得到了平滑^[11]。从图 5中原始影像和校正结果影像的对比可以看出，大部分阴影区域并没有得到很好的恢复，且有部分阴影区出现过校正的现象。从图 6中校正前后的对比可以看出，阴影区域均得到了很好的恢复，且阴影区域的山坡、山脊线轮廓清晰。

2.2.2 统计指标评价

本文采用离散分析、回归分析和直方图 3个指标对校正结果进行指标统计评价。

2.2.2.1 离散分析

为了检验影像的整体效果，采用离散度 (方差/均值) 指标对研究区影像进行统计分析。离散度由下式计算^[3]得出

式中，DI表示离散度；M与SD分别表示SPOT-5影像校正前后的辐亮度均值及标准差。

地形校正前，影像离散度相对较大，通过不同方法校正后，离散度均得到了不同程度的减小，说明地形校正减弱了地形对地物信息遮挡的效应，从整体上改善了图像质量^{[3, 7]}。直方图法得到的遮蔽因子进行地形校正的结果离散度小于DEM法，这充分说明当DEM数据精度不高时，直方图在提取遮蔽因子时的准确程度，可以大大提高地形校正的效果。

2.2.2.2 回归分析

影像与太阳有效入射角之间的回归关系反映了地形效应的程度^[12]。一般情况下，地形校正前后拟合的线性回归方程斜率k越大，地形效应越明显^[13-14]。图 7、图 8分别为原始影像和不同校正方法所得结果与太阳入射角的回归曲线示意图。

从图 7、图 8可以看出，经过校正后的回归曲线斜率小于原始影像与太阳入射角的回归曲线斜率，其中直方图法校正结果 (图 8(b)) 的斜率最小，说明经过校正后，地形效应得到明显减小，且直方图法的校正效果优于DEM法。

2.2.2.3 统计直方图

图 9是两种方法校正结果影像的灰度直方图。对比图 3和图 9(a)可以看出，通过DEM计算出的遮蔽因子对地形校正有一定的效果，但是曲线的波谷说明该方法没有很好地去除阴影与非阴影的边界过渡问题。从图 3和图 9(b)的对比可以看出，直方图呈正态分布，说明校正后影像的灰度值过渡平滑自然，阴影区域得到了很好地消除。与DEM计算遮蔽因子不同，直方图计算遮蔽因子的方法是通过判断在全色影像中，不同灰度值出现的频率来分离待校正的阴影区域的值。之所以可以使用直方图分割的方法，是因为在裸露山体的山区地物信息单一，并且在高分辨率的全色影像上能明显区分阴影区域与非阴影区域。将阴影区域看作目标地物，非阴影区域作为背景信息。试验中所选取的影像全部为山区且没有植被等其他地物信息覆盖。阴影区域与非阴影区域的比例均衡，不存在大面积都为阴影区域或只有极小一部分阴影区域的现象。当影像中出现部分城市或部分区域为植被所覆盖，直方图曲线都不出现双峰，也不呈正态分布形状；若影像中绝大多数面积都为阴影区域或非阴影区域，那么直方图曲线只呈正态分布，同样不出现双峰。这样均无法得到分割值，则无法通过该方法计算遮蔽因子。

使用直方图方法计算遮蔽因子进行地形校正，一方面由于计算遮蔽因子的过程中不受DEM的影响；另一方面由于山区这种特殊的地形，能明显区分阴影区域与非阴影区域。特别是在地表裸露地物信息单一的情况下，统计直方图为双峰状，此时可以快速准确地找到分割值，进而计算遮蔽因子，得到较好的校正效果。

3 结语

DEM计算遮蔽因子普遍适用于辐射传输方法的地形校正，但这种方法存在以下问题：

(1)  DEM精度将影响地形阴影判断的精度，如果能获取高精度DEM数据，该方法可以很好地提取阴影区。由于大多数情况无法获取高精度DEM，在选择使用30 m DEM时，DEM精度对阴影判断的精度有很大的影响，导致无法准确地提取阴影区域，造成校正结果不佳。

(2) 在判断阴影区域的过程中一般均使用高程采样点内插方法。由于山区特殊的复杂地形，在高程插值时容易造成高程值计算不准确，并且无法准确确定地形遮蔽角的搜索半径。这些都会造成无法相对准确地提取遮蔽因子，从而导致地形校正效果不理想。虽然延长地形遮蔽角的搜索半径可以使对阴影区域的判断更加精确，但是会增大计算量，降低校正时的计算效率^{[13, 15]}。

本文采用直方图计算遮蔽因子的方法进行地形校正，取得了较好的效果。通过与DEM法校正结果的比较分析可以得出，对于特殊的山区地形影像，使用直方图方法进行地形校正的结果要优于DEM法。