基于高分辨率遥感影像的城市地表信息提取

Extracting Information of Urban Land Surface with High Resolution Remote Sensing Image

摘要: 高分辨率遥感影像具有高空间分辨率、高清晰度、信息量丰富等优点，为城市地表土地覆盖和土地利用信息的提取提供了可能。本文以高分辨率遥感影像数据worldview为主要数据源，以基于立体像对提取DEM、基于分类的面向对象影像分析方法和基于阴影的建筑物高度获取为主要技术，自动提取武汉市十里铺区域的地表信息，得到了该区域的DEM、土地分类情况以及建筑物高度情况等信息，精度高，速度快，结果表明利用高分辨率的遥感影像对复杂的城市地表信息获取是可行的。

Abstract: With the advantages of high spatial resolution and definition, and rich information, land classification and land utilization of urban surface can be reached by using the high-resolution remote sensing data. Based on the high-resolution remote sensing image data worldview as the main data source, the methods of stereo images, object-oriented technique and architecture shadow are used to extract the information of DEM, land classification and building heights in Wuhan Shilipu district. The results show the high precision and rapid speed when extracting the surface information, and high-resolution remote sensing image is feasible for complex urban surface information acquisition.

1. 引言

下垫面的变化是城市化进程的一个主要方面，下垫面的变化主要表现为天然地表被建筑物、道路等所覆盖，且建筑物分布复杂，高低各异，道路交错，往往一个小区域内就包含了多种地表信息，因此中低分辨率的遥感影像资料难以用于城市小区域的地表信息获取 [1] 。高分辨率卫星遥感影像特别是 IKONOS、QUICKBIRD、worldview等的出现，为用户提供了地物景观的结构、纹理、几何和位置等信息，使利用遥感手段提取地表各种景观信息的潜力大大增强 [2] 。通常认为地面分辨率高于5 m的遥感影像即可认为是高分辨率遥感影像 [3] ，而遥感技术的不断发展，使目前部分商业遥感卫星的多波段影像空间分辨率可达米级，全色波段影像空间分辨率甚至可达亚米级。在这样高的空间分辨率下，影像的空间结构、地表纹理等细节，以及地物景观的边缘信息可以清晰地呈现出来，即在获取到地物光谱信息的同时也可获取地物的位置、形状、纹理、大小、阴影等信息 [4] ，使在较小的空间尺度上观察地表的细节变化成为可能，从而为精确提取城市地表信息提供了基础及条件。Yang X. [5] 充分分析了高分辨率遥感影像在城市监测、分析、调控中发挥的重要作用。Zhou等 [6] 用高分辨率LIDAR数据对城市下垫面景观进行分类，在地块尺度上的分类结果精度可达92.3%；Pu等 [7] 通过对IKONOS影像进行地表分类的研究表明，基于面向对象的方法能够使高分辨率遥感影像提取地表信息的精度大大提高；杜培军 [8] 用SPOT7影像对城市地表结构类型进行分类，提取结果较为满意；杨朝斌等 [9] 充分论证了高分辨率遥感影像在城市LUCC应用中的适用性及面向对象分类方法在影像信息提取中的精确性；刘东升等 [10] 详细阐述了国产高分遥感数据的处理方法及在资源、环境监测，应急救援等方面的实际应用情况。大量研究表明，高分辨率遥感影像以其高效便捷及高实效性的特点，推动了城市地表信息提取技术的发展；如何使提取的信息精度更高则成为了目前研究的热点之一。本文针对快速城市化导致城市地表信息快速改变的现状，利用worldview高分辨率影像资料，以基于立体像对提取DEM、基于分类的面向对象影像分析方法和基于阴影的建筑物高度获取为主要技术，以期快速、高效、精准地获取地表信息。

2. 研究方法

2.1. 立体像对提取DEM数据

立体像对通过传感器来模拟人眼从而构造出立体视觉感。基于立体像对的DEM提取原理可描述为 [11] ：在地表上空两个点拍摄地面某一点会形成一个角度，这两个点确定后，所形成的角度越大，地物就越高；反之，角度越小，地物就越低。卫星在成像过程中，会受到多种因素的影响，使所成的影像里面个别点发生不同程度的几何变形和失真 [12] 。传感器校正是恢复传感器状态和地面相对位置关系的重要环节，其精度对DEM提取影响很大，建立几何模型来消除影像中的畸变，是立体像对提取DEM中的必要程序 [13] 。有理函数模型RFM (Rational Function Model)是一种与传感器无关的通用几何模型 [14] ，该模型通过一组有理函数来描述目标与像点之间的映射关系。本文利用有理函数模型进行几何定位，借助ENVI 5.0软件的DEM Extraction模块交互完成DEM提取。

首先，选取需要处理的立体像对数据，读取影像附带的RPB文件信息，自动解算得到有理函数模型的参数RPCs、最大及最小高程值。其次，在地表选取控制点，选取容易识别的点为宜，如在道路的交叉处等，尽量均匀分布于研究范围。然后，匹配选取的连接点，为了消除两幅图像在Y方向的视觉差，将匹配问题从二维降为一维，Y方向的视差越小则连接点选取的精度越高。最后，生成核线影像，得到DEM数据。

2.2. 基于面向对象的土地分类利用提取方法

面向对象的分类方法充分利用高分辨率遥感影像的全色和多光谱数据的空间、纹理和光谱信息，对图像进行分割和分类 [15] ，可建立不同尺度的分类层次，具体更多的分类依据，使分类方式也更加灵活，操作性更强。基于面向对象的土地利用分类主要内容包括影像分割、特征提取、分类图导出等步骤 [16] 。

1) 影像分割

在利用遥感影像进行土地信息分类时，先要对得到的遥感影像资料进行数据的预处理，预处理的准确性直接影响之后信息提取的结果，本次预处理采用几何校正、图像融合和大气校正三个步骤。预处理后采用面向对象的分割算法——多尺度分割 [17] [18] 方法对影像数据进行分割。

多尺度分割方法根据不同的空间尺度参数对影像进行分割，从而得到多个有意义的多边形区域，同时考虑尺度、色调、形状、紧致度及平滑度等，采用区域合并算法，使分割后影像对象的权重异质性最小。权重异质性是光谱异质性和形状异质性的组合，可用公式表示为：

$f=w\ast h_{\text{color}}+\left(1-w\right)\ast h_{\text{shape}}$ (1)

式中，f为异质性大小，w为光谱权值( $0<w<\text{1}$ )， $h_{\text{color}}$ 为光谱异质性， $h_{\text{shape}}$ 为形状异质性，包括光滑度(smooth)和紧致度(compact)两个方面，计算公式如下：

$h_{\text{shape}}=w\ast h_{\text{compact}}+\left(1-w\right)\ast h_{\text{smooth}}$ (2)

式中， $h_{\text{compact}}$ 为紧致度， $h_{\text{smooth}}$ 为光滑度，决定因素包括对象的边长、总体像元数、外接矩形的最短边界等。若紧致度的比重较大，光滑度的比重较小，则分割后对象的形状比较接近矩形；若光滑度的比重较大，紧致度的比重较小，则分割后的对象边界较为平滑。

2) 特征提取

遥感影像多尺度分割完成后进行影像信息的提取。采用面向对象分类方法中的监督分类法对地物信息进行提取。利用ENVI EX模块目视解译法选取并定义数据。在分类过程中需要选择具体的监督分类法、设置相应的参数值。EX模块提供两种分类方法：① K邻近法，相对传统的最邻近方法产生更小的敏感异常和噪声数据集，相应的K参数值一般设置为3、5、7较好；② 支持向量法，是建立在统计学习理论基础上的一种方法，需要定义一系列参数，计算过程比K邻近法麻烦些。利用分类方法可得到监督分类结果，可输出为矢量和栅格两种类型数据。

2.3. 基于阴影的建筑物高度获取方法

当卫星高度角大于太阳高度角时，在遥感影像上比较容易识别出建筑物阴影，因其在影像上表现为低亮度值。图1中，H为所研究建筑物的高度，M为该建筑物房顶角点， $\beta$ 和 $\theta$ 分别为太阳高度角、太阳方位角， $\alpha$ 和 $\phi$ 分别为卫星高度角、卫星方位角， $\delta$ 为卫星和太阳的方位角交角， $M^{\prime }$ 为建筑物房顶角点阴影投影在影像上的位置， $M^{″}$ 为建筑物房顶角点在影像上成像的位置， $M^{\prime }{M}^{″}$ 为影像上的建筑物房顶角点与阴影之间的距离。根据卫星成像的几何关系可得如下公式：

$\{\begin{array}{l}{M}^{\prime }N=H\cot \beta \\ M^{″}N=H\cot \alpha \\ \delta =\theta -\phi \end{array}$ (3)

Figure 1. Relationship among sun, satellite altitude and azimuth angle, and building shadow

图1. 太阳和卫星高度角、方位角与建筑物阴影之间的关系

在 $\Delta M^{\prime }{M}^{″}N$ 中，应用余弦定理可得：

$M^{\prime }{M^{″}}^2=M^{\prime }{N}^2+M^{″}{N}^2-2M^{\prime }N\cdot M^{″}N\cdot \cos \delta$ (4)

即

$L^2=H^2{\cot }^2\beta +H^2{\cot }^2\alpha -2H^2\cot \beta \cot \alpha \cos \delta$ (5)

因此，建筑物的高度H为

$H=L\cdot {\left({\cot }^2\beta +{\cot }^2\alpha -2\cot \beta \cot \alpha \cos \delta \right)}^{-\frac{1}{2}}$ (6)

令

$k={\left({\cot }^2\beta +{\cot }^2\alpha -2\cot \beta \cot \alpha \cos \delta \right)}^{-\frac{1}{2}}$ (7)

则有

$H=L\cdot k$ (8)

为消除太阳方位角和卫星方位角的影响，对遥感影像作正射校正，使太阳方位角和卫星方位角相等，则有

$H=L/\left(\cot \beta -\cot \alpha \right)$ (9)

从理论上讲，知道太阳以及卫星的高度角、方位角，影像上的建筑物房顶角点与其阴影之间的距离，即可计算建筑物高度。

3. 研究区域及提取结果

3.1. 研究区域

武汉市作为湖北省的龙头，长江中游的中枢，城市化发展迅速。武汉市在建国初期，城市化呈现出比较合理的带状组团式发展格局，20世纪80年代，城市开始以长江为轴线发展，20世纪90年代以来，开始了规模巨大的城市化进程。从“九五”到“十五”期间，城市建成面积达355.65平方公里，超过了之前城市规划中2010年343.30平方公里的预测数据，城区面积扩大了1.48倍。本文选择了武汉市十里铺区域，该区域是城市化发展较典型的一个区域，地理位置在东经114.21˚~114.23˚，北纬30.55˚~30.57˚之间，面积为2.51平方公里，东与桃花岛小区相邻，北邻琴台大道，西于龙阳大道相距约1 km，地理位置如图2所示。

Figure 2. Geographical position of the study area

图2. 研究区域地理位置图

3.2. 高空间分辨率DEM提取结果

城市区域高精度DEM数据难以获取，野外实测难度大，所以采用worldview-1高空间分辨率卫星立体像对数据来提取5米分辨率的DEM数据。worldview-1卫星是Digital Globe公司在2007年9月18日于SLC-2W美国加州范登堡空军基地成功发射的一颗高分辨率商业遥感卫星，这是Digital Globe公司首次在民用领域将传感器的空间分辨率提高到0.5米，每天可采集数据超过100平方公里。本文获取的影像采集时间为2012年3月14日，采用的坐标系统为WGS84地理坐标系统，具体信息见表1。利用手持GPS设备(麦哲伦Explorist 610)在研究区域范围内实测了11个位置点的高程数据，尽量分散均匀，如图3所示，选取6个点作为DEM提取的地面控制点，剩下5个点作为验证数据。然后选取10个容易分辨的地物点，如道路交叉处等，作为连接点。利用ENVI的DEM Extration模块提取5米空间分辨率的DEM数据，如图4(b)。

Table 1. The worldview-1 stereopair information of the study area

表1. 研究区域worldview-1立体像对信息

Figure 3. The worldview-1 stereopair data of the study area

图3. worldview-1立体像对研究区域数据

对提取的DEM数据进行精度分析，选取5个验证点所在栅格的高程值与实测值进行对比，如表2所示，大部分点的相对误差都小于5%，平均相对误差为3.71%。并且与从测绘部门获取的研究区域内10 m精度的实测DEM数据如图4(a)进行对比，可以看到提取的5米分辨率DEM数据的高低趋势基本和实测DEM数据保持一致，说明提取得到的高空间分辨率DEM数据较为可靠，可以为后续研究使用。

Table 2. Comparison between DEM and measured value of verification point

表2. DEM与验证点实测值对比

(a) (b)

Figure 4. Comparison of 10m measured DEM (a) and extracted 5 m DEM (b)

图4. (a) 10 m实测DEM和(b) 提取的5 m DEM对比

3.3. 高精度土地利用分类图

本次研究利用获取的worldview-2高分卫星影像进行土地利用分类，worldview-2高分卫星于2009年10月6日发射升空，能够提供0.5米分辨率的全色影像和1.8米分辨率的多光谱影像，本次获取的研究区域影像如图5所示。

首先对worldview-2全色影像和多光谱影像进行融合，得到0.5米分辨率的多光谱影像，随后利用ENVI的Feature Extration工具根据临近像素亮度、纹理、颜色等对融合影像进行多尺度分割，经过多次试算，分割阈值设为50，然后再利用Full Lambda-Schedule算法在结合光谱和空间信息的基础上迭代合并临近的小斑块，经试算，阈值设为94，最后利用监督分类方法，通过观察和识别图像上的地物特征，选取建筑物、树木、草地、裸土、不透水面五类地物各10个样本作为训练样本数据，选取K邻近法进行分类，K值设为3。得到的矢量分类图如图6所示。通过在Google Earth上随机选取20个地物清晰、易辨别的位置点对分类结果进行对比验证，结果有18个点的地物一致，准确率达到90%。

Figure 5. worldview-2 satellite image of the study area

图5. 研究区域worldview-2卫星影像

Figure 6. Land use classification of the study area

图6. 研究区域土地利用分类图

3.4. 建筑物高度提取

从获取的worldview-2卫星影像的头文件可以得到卫星过境时与太阳的高度角为70.8˚，方位角为148.4˚。然后根据worldview-2全色影像中建筑物顶角点与其阴影的像素个数计算距离，利用公式(9)可以计算得到建筑物高度，从而可以得到整个研究区域内的建筑物高度图层，如图7。

Figure 7. Height of buildings in the study area

图7. 研究区域建筑物高度图

为了验证计算得到的建筑物高度的精度，在研究区域内选取识别度较高的2座高楼和2座矮楼的高度，并利用全站仪测量其高度，并与计算结果相对比，结果见表3。从表3可知，4个建筑物高度的相对误差都小于20%，说明计算结果较为可靠。

Table 3. Comparison between extracted values and measured values of building heights

表3. 建筑物高度遥感提取值与实测值对比情况

4. 结论

本文通过高分辨率worldview遥感数据利用立体像对提取了城市地表的DEM信息、基于面向对象的分类方法提取了土地利用信息、结合建筑物阴影获取了建筑物高度信息，计算精度较高，成果质量可靠。由结果可知，利用上述方法提取高分辨率遥感影像获取城市地表信息具有较高的实用性，可为今后城市地表信息提取提供借鉴和参考。