城市建筑物高度信息是城市规划、建设项过程中的重要的数据^[1]，目前利用GPS、全站仪野外测量或利用航空影像立体像对获取建筑物高度的方法，虽然在技术上很成熟，但投资大、成本高、效率低^[2]。随着遥感技术的发展，高分辨率遥感影像中地物的亮度、纹理、形状、空间位置信息更加精细，为城市定量遥感提供了准确的数据源，使建筑物高度信息快速提取成为了可能。通过高分影像提取出阴影斑块，建立建筑物、阴影、太阳、卫星的几何关系模型来估算建筑物高度，成为了当前研究的热点，如何有效、精确地提取建筑物阴影信息成为了建筑物高度反演的关键。

目前，高分辨率遥感影像的阴影检测方法可分为基于模型和基于阴影特征两大类^[3]。其中基于模型的方法将有关场景、目标几何结构、光源信息等辅助数据作为已知条件输入，鉴于光源、场景及目标结构信息不容易获取，且计算量大。因此这类方法局限性较大^{[4, 5]}；而基于阴影特征的方法主要利用待检测图像自身阴影的一些属性，通过一系列的图像增强算法，提高阴影信息的辨识度，然后借助阈值分割技术提取出阴影区域。由于其算法设计简单而有效是目前研究和应用的主流。Ma等^[6]通过对阴影颜色特性分析，在色调、饱和度、明度(hue saturation value，HSV)空间下，将饱和度与明度归一化，建立阴影指数，有效地突出了阴影信息；Polidorio等^[7]基于色调、亮度、饱和度(hue intensity saturation，HIS)颜色模型，将亮度(I)与饱和度(S)进行差值运算，提高了阴影区域与其他地物的对比度，最后通过人工设定阈值，实现阴影区域的分割；赵志明等^[8]利用微分数学形态学建立阴影指数(morphology shadow index，MSI)，并结合尺度分割优化和面向对象分类技术提取了高分辨率遥感影像建筑物及其阴影轮廓；段光耀等^[9]通过分析遥感阴影光谱特性，构建出亮度特征I、彩色不变特征C3、主成分变换第一特征分量PCI等多个分量，并与面向对象分类技术相结合，实现了阴影区域的提取。

本文以郑州市某小区GF-2遥感影像作为数据源，利用微分数学形态学阴影指数能够充分考虑到阴影的亮度、局部对比度、大小和方向性等特征，并结合阴影在HIS色彩空间下，饱和度分量值较大、亮度值较小的特性，提出一种适合于GF-2遥感卫星数据的阴影提取模型。

1 研究区与数据

本研究选用郑州市某小区作为研究区域，该小区有居民楼43栋、两条主干道路以及若干小路，植被茂盛且颜色较深。选用2015-08-14该小区GF-2遥感影像作为数据源，包括多光谱影像数据(4 m)和全色影像数据(1 m)。综合考虑全色波段分辨率高和多光谱图像波谱特征全面的特点，为了提高影像的分辨率又不丢失颜色特征信息，本文对全色数据选择Pansharp方法进行融合，融合后影像分辨率为1 m，如图 1所示。

2 研究方法

本研究选用融合后的GF-2影像数据(分辨率为1 m)作为数据源，提出了一种多特征分量构建的阴影检测方法。

1) 利用阴影在HIS色彩空间下饱和度分量值较大、亮度值较小的特性，构建出阴影指数1(NDSI)。

2) 基于微分思想，建立不同尺度、方向的结构算子，构建出微分形态学阴影指数2(MSI)。

3) 利用OSTU阈值分割技术，检测出类似阴影区域1、2，然后将两个类似阴影区域进行逻辑与运算，提取出精确阴影区域，并计算其在太阳方位角方向上的长度，构建建筑物、阴影、太阳、卫星的几何关系模型，从而反演出研究区内建筑物高度信息，技术流程如图 2所示。

2.1 阴影指数构建 2.1.1 HIS色彩空间下阴影指数构建

在RGB色彩空间下，阴影的亮度较低，容易和其他低亮度地物混淆，不容易将其区分出来。将RGB转换到HIS色彩空间^[10]后，各个分量的独立性增强、数据冗余度较低。在HIS颜色模型中，阴影与其他地物相比，具有以下3个显著特征:①亮度值更低，由于阴影是由不透光物质的遮挡形成的，体现为I值较低；②饱和度更高，由于大气瑞丽散射的影响，阴影区域的光主要来自波长较短的蓝紫光，体现为S值较大；③色调更大，由于影像中阴影区域接近于黑色，体现为H值较大。基于上述3个特征，通过实验对比，发现选取S和I分量构建的阴影指数可以较好地增强阴影信息。阴影指数构建如式(1)所示：

$ {\rm{NDSI = }}\frac{{S - I}}{{S + I}} $

(1)

2.1.2 形态学阴影指数(MSI)构建

形态学阴影指数(MSI)建立在阴影的光谱结构特征和一系列形态学操作的基础上，充分考虑到阴影局部明暗对比度、亮度、形状、大小和方向性等特点^[11]。形态学阴影指数的主要构建步骤如下：

1) 形态学黑帽构建

$ {W_{ - TH}}\left( {d, s} \right) = r_b^{{\rm{se}}}\left( {d, s} \right) - b $

(2)

式中，W_-TH表示形态学黑帽构建后的图像；se为形态学结构算子；b为主成分变换后的影像；d和s分别表示结构算子的方向和尺度。本研究选取的结构算子为矩形，结构算子选取了4个方向，分别为0°、45°、90°、135°；结构算子在尺度上选取从2~20，间隔为2。道路和阴影的主要区别在于方向的延伸，阴影的延伸方向是杂乱无章的，而道路通常是沿着某一两个方向，通过多个方向和多尺度的形态学运算，可以较好地区分两者。

2) 形态学剖面计算

$ \left\{ \begin{array}{l} {M_{PW - TH}}\left( {d, s} \right) = {W_{ - TH}}\left( {d, s} \right)\\ {M_{PW - TH}}\left( {d, 0} \right) = b \end{array} \right. $

(3)

3) 微分形态学计算

$ \begin{array}{l} {D_{{M_{PW - TH}}}}\left( {d, s} \right) = \\ \left| {{M_{PW - TH}}\left( {d, s + \Delta s} \right) - {M_{PW - TH}}\left( {d, s} \right)} \right| \end{array} $

(4)

4) 形态学阴影指数计算

$ {\rm{MSI = }}\frac{{\sum\limits_{d, s} {{D_{{M_{PW - TH}}}}\left( {d, s} \right)} }}{{D \times S}} $

(5)

式中, S=(S_max-S_min)÷ΔS+1；D为阴影计算剖面时的方向数。本研究中D=4，S_min=2，S_max=20，ΔS=2。MSI的构建是基于阴影结构在微分形态学剖面D_MPW-TH大部分方向上具有较大的值，由于阴影周围是亮度值较高的建筑物，两者具有较强的对比度。因此，MSI值较大的区域是阴影^[12]。

2.2 建筑物高度反演 2.2.1 阴影长度计算

阴影长度计算是利用阴影反演建筑物高度的关键，为了提高阴影长度计算的效率和精度，其计算方法如下:

1) 获取影像拍摄时的太阳高度角、太阳方位角、卫星高度角、卫星方位角角信息，分析出太阳光线与阴影的几何关系。

2) 根据推算出的建筑物与阴影之间的几何关系，生成一系列线簇，方向与太阳光线入射方向一致。

3) 将阴影图层与线图层进行求交运算，为保证计算结果的准确性，每个阴影多边形内应不少于8个矢量线。

4) 获取每个阴影图层内部测量线的长度信息，并对其求均值(去除最大和最小值)，作为阴影的长度信息。

2.2.2 高度反演模型

根据太阳、卫星及阴影的成像关系，建筑物高度和阴影的关系，可以分为太阳和卫星在同侧和异侧两种情况。假定如下条件成立：物体垂直于地球表面，物体的影子直接投影在地面上，影子从物体的底部开始^[13]。本文测量阴影长度时，是按照太阳方位角方向进行的，已经考虑了太阳方位角的影响，遥感卫星的方位角一般都接近于90°，且仅决定了阴影显示的是局部还是全部，因此本文不考虑卫星方位角对试验结果造成的影响^{[14, 15]}。如图 3(a)所示，太阳和卫星在建筑物同侧时，影像上只能显示部分阴影。在不考虑卫星方位角和太阳方位角的影响时，可得如下公式：

$ \left\{ \begin{array}{l} A = S - B = H/\tan \theta - H/\tan \omega \\ H = A\tan \omega \tan \theta /\left( {\tan \omega - \tan \omega } \right) \end{array} \right. $

(6)

式中，H为建筑物高度；S为建筑物阴影长度；A为影像上阴影的长度；θ为太阳高度角；ω为卫星高度角。如图 3(b)所示，太阳和卫星在建筑物异侧时，影像上能看到建筑物完整的阴影。根据几何关系得出阴影与建筑物高度之间的关系，在不考虑太阳方位角对影像影响的情况下，可得如下公式：

$ H = A\tan \theta $

(7)

3 实验仿真及误差分析 3.1 阴影提取

本文在MatlabR2015b平台进行实验仿真，完成了实验区内43栋建筑物阴影提取工作。根据阴影在HIS色彩空间下的低亮度、高饱和度的特性，构造出阴影指数1(图 4(a))；采用OSTU阈值分割算法计算出阴影与非阴影之间的阈值，大于该阈值的则认为是阴影区域1(图 4(b))。在0°、45°、90°、135°等4个方向，同时选取了大小为2~20，间隔为2的矩形结构算子，利用微分数学形态学，构造出形态学阴影指数2(图 4(c))，利用阈值分割技术，并结合面向对象的思想去除面积小于50的斑块，从而提取出类似阴影区域(图 4(d))。最后将阴影指数(NDSI)和形态学影指数(MSI)分割得到的两个类似阴影区域做逻辑与运算，计算出精确阴影区域(图 4(e))。为了较为直观地观察提取效果，将提取得到的阴影图斑边界叠加到原始影像(图 4(f))，结果显示提取得到的阴影边界可以较好地与原始影像上的阴影区域套合。

3.2 建筑物高度计算及误差分析

根据提取得到的精确阴影区域，获取每个阴影图斑的长度信息，利用建筑物高度估算模型反演出研究区内建筑物高度信息。为了验证本文反演高度信息的精度，利用全站仪实地测量出研究区内建筑物高度数据，并与反演高度进行对比。建筑物实测高度信息、反演高度信息以及绝对误差统计结果如表 1所示。由表 1的统计结果可知，绝对误差均值为0.56 m，反演高度与实际建筑高度整体相差不大，反演较为准确。

图 5为研究区43栋建筑物反演高度与实测高度的散点图。图 5中虚线表示直线方程Y=X，假设在理想状况下反演高度和实测高度相等，红色实线表示反演高度和实测高度的相关趋势线，相关系数达到95.03%，说明总体反演精度较高。研究区43栋建筑物绝对误差的统计如图 6所示，绝对误差不大于1 m的建筑物有41栋，占到总体95.35%，而绝对误差超过1 m的建筑物只有两栋，约占总体4.65%，证明了本研究方法的稳定性和可靠性。由此可见，本文方法在分辨率为1 m左右的高分影像数据估算建筑物高度是可行的。

4 结束语

本文选用融合后的GF-2影像数据(分辨率为1 m)作为数据源，利用微分数学形态学、阴影在HIS色彩空间下的特点，提出了一种结合特征分量构建与阈值分割的阴影检测方法。

1) 本文算法能够充分考虑到阴影局部明暗对比度、亮度、形状、大小和方向性的特点，可以很好地对阴影进行识别，有效地估算出建筑物高度信息，解决了建筑物与道路光谱特性相互混淆的问题，具有操作简单灵活、快捷高效等优势。

2) 对研究区内的43栋建筑物反演高度与实测高度进行了统计，绝对误差均值为0.56 m，且绝对误差均在1.2 m以内，误差小于1 m的占到总体95.35%，总提取精度达到95.03%，且阴影形状保持较好。

3) 随着高分影像的不断推广和普及，该方法具有大范围检测建筑物高度信息的前景，开发出全自动的阴影提取、长度计算算法，将其用于城市建筑物容积率、密度等城市诸多检测问题。