CN104463970B - 一种基于遥感影像确定城市三维重心的方法及其应用 - Google Patents

Abstract

本发明提供了一种基于遥感影像确定城市三维重心的方法，所述方法利用城市的遥感影像，结合城市地籍、地形数据，提取目标城市的各个建筑的轮廓和相应高度。然后，按预定的正方形格网将城市建筑高度和轮廓数据划分成一个个标准单元，确定各个格网(单元)的建筑高度，进而构建建筑高度等值线。最后，基于目标区域的建筑高度等值线构建目标区域的三维模型，并基于三维模型计算目标区域的三维重心。本发明能够准确快速地确定出城市或城市某一区域的三维重心。因为三维重心可从平面尺度、立体高度等角度对城市的区位中心进行测度，所以利用三维重心来进行商业中心的布局、基础设施建设项目的选址比利用二维重心来辅助决策要更加科学。

Description

一种基于遥感影像确定城市三维重心的方法及其应用

技术领域

本发明涉及图像处理领域，具体涉及一种基于城市影像确定城市三维重心的方法及该方法的应用。

背景技术

城市作为人类活动最为活跃的空间区域，其形态在城市形成与发展的时空过程中不断发生着变化，同时也反映着这一过程的特征和规律。在城市建设、土地管理等领域，都需要掌握城市的发展趋势状况，以便基于城市的当前发展状况进一步对新建建筑进行规划、设计、选址等工作。在对城市的发展趋势进行研究过程中，获得城市重心及其变化过程是至关重要的。

很多学者都对城市形态的演变和重心转移规律进行了深入研究。在相关研究中，张树清等在矢量景观理论中对于单个地理对象的空间重心给出了理论定义和公式推导，匡文慧等人在此基础上提出空间重心转移距离和重心转移角度两项指标，用以进一步分析重心转移的距离和方向性，并在近几年的城市土地扩张和形态演化等相关研究中得到了广泛应用。同时，重心的研究方法也被拓展应用于人口重心、经济重心、景观空间重心、生态服务重心和土地利用重心等的研究中，它是将重心作为一种加权平均中心，以人口、经济总量、各类土地资源面积等作为权重参考，以经纬度坐标或者直角坐标为基础计算得到的。

当前对城市重心的研究多基于城市平面形态获得城市二维重心。在城市垂向伸展日益显著的背景下，非常有必要对城市的三维重心进行变化分析。因为三维重心更能够反映出城市人口的集聚和经济发展的态势，以及城市区位中心的移动方向。

但是，现有的城市重心研究方法只能够获得城市二维重心的发展和变化，还无法得到城市三维重心。

发明内容

针对上述问题，本发明提出了一种基于城市的遥感影像确定城市三维重心的方法，并且，本发明将该方法应用于商业中心的布局、基础设施建设项目的选址中，以便提高商业中心和基础设施的利用效率。

具体而言，本发明提供了一种基于遥感影像确定城市三维重心的方法，其特征在于，所述方法包括如下步骤：

步骤1，获取目标城市的一个或多个目标区域的影像以及所述目标区域的城市地籍和/或地形数据；

步骤2，对所获取的影像进行预处理；

步骤3，基于所获取的影像以及所述目标区域的城市地籍和/或地形数据提取所述目标区域的各个建筑的轮廓；

步骤4，基于所获取的影像以及所述目标区域的城市地籍和/或地形数据提取所述目标区域的各个建筑的高度；

步骤5，以预定尺度绘制覆盖所述目标区域的正方形格网，所述正方形格网具有若干正方形标准单元；

步骤6，基于所获得的所述目标城市的各个建筑的轮廓和高度数据，采用所绘制的所述目标区域的正方形格网将各建筑划分在所述正方形格网的一个个的标准单元内；

步骤7，计算所述正方形格网的每个标准单元内的建筑总体积，从而计算出每个标准单元内的平均建筑高度；

步骤8，基于所述正方形格网的每个标准单元内的平均建筑高度，生成所述目标区域的建筑高度等值线；

步骤9，基于所生成的所述目标区域的建筑高度等值线获得所述目标区域的三维立体模型；

步骤10，沿水平面方向对所述三维立体模型进行分割，以获得与水平面平行的若干建筑薄层；

步骤11，分别计算每个建筑薄层的重心的平面坐标(x,y)和高度坐标(z)；

步骤12，基于各个建筑薄层重心的三维坐标及各个薄层的体积计算各薄层的组合体的重心的三维坐标。

进一步地，所述步骤4中，提取所述目标区域的各个建筑的高度的步骤包括：

步骤4-1，获取所述目标区域内若干已知建筑物的高度数据；

步骤4-2，基于所述影像，对于所述目标区域中的每个建筑，评估该建筑的多条特征线的提取难度；

步骤4-3，对于每个建筑，基于对多条特征线的提取难度的评估，选取最低难度的特征线；

步骤4-4，从所述影像中提取选定的特征线；

步骤4-5，基于所述目标区域内已知的建筑物高度数据反算所述影像的摄取参数；

步骤4-6，基于所选定特征线的长度、所述摄取参数，计算所述目标区域内，每个未知建筑的高度。

进一步地，所述方法还包括判断所述影像是单景影像还是多景影像的镶嵌，如果所述影像为多景影像的镶嵌，则在所述高度提取的步骤中，分别针对各景影像进行参数计算和高度提取。

进一步地，所述步骤11中，每个建筑薄层的重心的平面坐标等于该建筑薄层的形心的坐标，每个建筑薄层的重心的高度坐标为：

其中，i为“建筑薄层”所在的层数；h为“建筑薄层”的厚度，即两条建筑高度等值线数值之差。

进一步地，所述目标区域的三维重心的坐标为：

式中，为城市三维重心的坐标，x_i、y_i、z_i为第i个“建筑薄层”的三维重心的坐标，v_i为第i个“建筑薄层”的体积，等于“建筑薄层”面积与厚度h的乘积。

进一步地，所述方法还包括计算所述目标城市的三维重心的转移距离，计算公式如下：

其中，D_t+1为三维空间重心转移距离，为基于在第一时间获得的所述目标城市的影像计算所得的所述目标城市的三维重心的坐标， 为基于在第二时间获得的所述目标城市的影像计算所得的所述目标城市的三维重心的坐标。

进一步地，所述方法还包括计算三维重心转移距离在xoy平面上的投影距离D’_t+1：

另一方面，本发明提供一种基于上述方法获得的城市三维重心的应用，其特征在于，所述应用包括：基于所获得的城市三维重心为商业中心进行选址，使得所述商业中心位于所述三维重心的预定距离之内。

另一方面，本发明提供一种基于上述方法获得的城市三维重心的应用，其特征在于，所述应用包括：基于所获得的三维重心为基础设施建设项目选址服务。

另一方面，本发明提供一种基于上述方法获得的城市三维重心的应用，其特征在于，所述应用包括：基于所获得的三维重心进行城市发展方向和重点发展区位的选择。

在所述高度提取的过程中，所采用的特征线包括：第一特征线——目标建筑的屋顶位移点与其阴影点连线、第二特征线——目标建筑的阴影全长、第三特征线——目标建筑的高差引起的屋顶像点位移，以及第四特征线——目标建筑被遮挡后的阴影长。

进一步地，本发明中的高分辨率影像数据的来源较为广泛，可以是卫星影像、航拍影像以及Google Earth、百度地图等网络上可共享的卫星数据。特别地，Google Earth和百度地图所提供的影像是比较典型的无参数高分辨率影像，可以利用专门的影像截取软件获取待提取区域范围内的相应影像数据，具有速度快、质量高、拼接误差小的特点，利于之后的影像配准。所述影像的摄取参数包括所述影像摄取时，摄取装置的摄取方向相对于地面的高度角和方位角，以及所述影像摄取时太阳的高度角和方位角。

由于三维重心是城市空间的真正重心，它反映出了真正的城市区位中心，根据城市内部地域结构相关理论，城市区位重心可用来指导中心商务区的布局、城市基础设施建设、城市发展方向和重点发展区位的选择等。

技术效果

本发明基于对二维平面重心的求取方法的深入分析，推导了一种由城市建筑高度等值线计算城市实体三维重心的方法，将城市三维实体看作由一层层“薄板”搭建成的组合体，利用组合体的重心计算公式求取城市三维重心。本发明以南京市为实例进行了验证，将所求得的城市二维重心和三维重心的转移结果进行了比较分析。结果表明，该计算方法可方便计算城市实体的三维重心。城市三维重心充分考虑到城市垂直方向上的伸展，可为更全面地反映城市空间变化的规律和特点提供一条有效途径。

另外，本发明在对高度进行提取时，能够实现在无影像摄取参数的情况下，计算所获影像内的建筑物的高度信息。本发明可采用网上易于获取的无参数的影像进行建筑高度的提取，而无需花经费购买带完整元数据(参数)的影像，适用性广、成本低。本发明可以应用于城市规划和管理、土地管理、数字城市构建等领域，能够显著节约社会成本。

下面结合说明书附图1、2、3和具体实施案例对本发明做进一步阐述。

附图说明

图1为本发明的三维重心提取方法的示意性流程图；

图2为以卫星影像为例所示出的卫星成像原理示意图；

图3为对建筑薄层进行二维重心计算的示意性视图；

图4为通过对三维重心的计算原理示意图；

图5为以南京市为例所计算的城市的二维重心和三维重心的平面投影随时间的变化。

具体实施方式

本发明对三维重心求解的基本原理。本文在研究二维重心提取方法的基础上，提出城市三维重心提取的基本思路(图1)：首先，基于城市高分辨率遥感影像和城镇地籍图、地形图等辅助资料，提取城市中所有建筑轮廓和建筑高度；其次，将整个城市建成区划分成合适边长的若干个规则格网，本文选择500m边长的规则格网，对跨格网建筑进行分割，根据各格网内的建筑轮廓和建筑高度计算格网内的建筑总体积，依据格网面积计算各格网的平均建筑高度；再次，用反距离权重法(本领域常用的一种空间数据插值方法)对格网数据进行插值处理，并采用等值线追踪法(用于将高度值相等的点连接起来)生成建筑高度等值线，将建筑高度等值线三维立体化，生成三维城市简化模型；最后，将三维城市简化模型看成是由若干“薄板”构成的组合体，利用组合体重心求解原理计算三维城市的重心。

下面结合图1为例，具体说明本发明的三维重心确定方法的具体流程。

如图1所示，首先，需要获取南京市建成区的遥感影像以及相关的城市地籍、地形数据。这两种数据可以相互补充也可以彼此结合使用。

然后，对所获取的遥感影像进行预处理，例如，调整所述遥感影像的亮度、对比度和色彩平衡。完成增强处理后，利用地形图或正射影像图进行影像的配准工作，也可采用地面实测点坐标进行配准。配准后可生成直接用于特征点、线提取的底图，也可以方便轮廓的确定。本发明对图像的增强和配准，采用的是现有技术，这里不再累述。

接下来，基于所获取的遥感影像并结合地籍、地形数据(部分不易获取高分辨率影像的区域可利用地籍和地形数据中的建筑楼层信息进行建筑高度的换算)提取南京市建成区的各个建筑的轮廓，并且提取该区域的各个建筑的高度。建筑轮廓可以采用现有的建筑轮廓提取方法进行，并且进行适当的人工干预。但是，由于遥感影像中没有直接的建筑高度数据可供使用，所以，本发明采用了一种新的建筑高度提取方法来获取各个建筑的高度数据，建筑高度的提取方法，在后文中详细描述。

然后，在计算机中，以500m为边长构建正方形格网，格网的数量要能够覆盖目标区域。

接下来，按正方形格网将目标区域的所有建筑划分进规则格网单元中。

然后，计算各个格网内的建筑总体积，以及每个格网内的平均建筑高度，并追踪目标区域的建筑高度等值线。

然后，基于所追踪的目标区域的建筑高度等值线获得目标区域的三维立体模型。然后，如图4所示，沿水平面方向对所述三维立体模型进行分割，以获得与水平面平行的若干建筑薄层。

接下来，分别计算每个建筑薄层的重心的平面坐标和高度坐标。对于均质物体(密度、比重为定值)、在相同的重力场下，其重心、质心和形心是重合的，在本研究中，将城市建成区的三维形态看成是均质的，因此求重心的问题可转化成求形心的问题。矢量景观理论对于单个地理对象的空间重心给出了理论定义和公式推导。

如图3所示，S为城市建成区，其面积为A，(x，y)是S内部某点的坐标，C点为S的重心，其坐标为()，dxdy是(x，y)点的微分面积。

的公式推导如下：

式中，Q_x、Q_y可理解为绕x轴和y轴的力矩。

然后，如图4所示，可对各个建筑薄层的重心进行平面坐标和高度坐标(三维坐标)的求取，进而结合各建筑薄层的体积计算各薄层的组合体的三维坐标。

各“薄板”重心的平面坐标x、y通过求平面重心的方法获得，各“薄板”重心的z坐标为：

式(4)中，i为“薄板”所在的层数；h为“薄板”的厚度(高度)，即两条建筑高度等值线数值之差。基于各部件(薄板)的组合体的三维重心的坐标为：

式(5)中，为城市三维重心的坐标，x_i、y_i、z_i为第i个“薄板”的三维重心的坐标，v_i为第i个“薄板”的体积，等于“薄板”面积与厚度h的乘积。

为进一步分析不同时期城市三维空间重心转移的距离和方向性，引入如下变量：

三维空间重心转移距离(D_t+1)：

三维重心转移距离在xoy平面上的投影距离(D′_t+1)：

三维重心高差变化(Z_t+1)：

三维重心在三维空间的转移角度可分解为xoy平面上的转移角度和在垂直方向上(z轴)的转移角度。其中，xoy平面上的转移角度：

当且

当且

当

垂直方向上的转移角度：

上式中，β_t+1为从t到t+1时期三维空间重心转移在xoy平面上投影后从正东方向逆时针旋转的角度，α_t+1为三维空间重心转移在垂直方向上的变化角度。

下面详细介绍在本实施例中，各个建筑的高度提取的具体过程。

在高度提取中，假定如下条件是成立的：①建筑物垂直于地球表面；②建筑物的影子直接投影在地面上；③影子从建筑物的底部开始。当太阳和卫星处于建筑物的同侧时，卫星成像如图2所示，图中MO为建筑的实际高度H，O点为建筑物顶面角点M点的垂直投影点，建筑立面MNPO在地面的投影为BEPO，其地面阴影为ADPO。A点是M点的阴影点，B点是M点在影像上的成像点位置，C点为MN在地面的投影BE和AO相交的点。卫星高度角ω＝∠MBO，太阳高度角θ＝∠MAO，卫星和太阳的方位角分别为α和γ，卫星和太阳的方位角交角是∠BOA＝α－γ。在图2中，A、B、O、C点为特征点，AO、BO、AB、AC线为特征线，其中AO为阴影全长(l₁)，BO为建筑高差引起的建筑顶面角点像点位移(l₂)，AB为屋顶角点成像点与其阴影点的连线长(l₃)，AC为建筑遮挡后的可见阴影长(l₄)分别根据四类特征线段推导建筑高度的计算模型。以下公式中，ω为卫星高度角，θ为太阳高度角，α为卫星方位角，γ为太阳方位角，δ为阴影顺时针到建筑物的夹角。

(1)利用阴影全长l₁计算建筑物高度：

H＝AO×tanθ (13)

(2)利用建筑高差引起的建筑顶面角点像点位移l₂计算建筑物高度：

H＝BO×tanω (14)

(3)利用屋顶角点成像点与其阴影点的连线长l₃计算建筑物高度：

∠ABO＝α-γ (15)

由余弦定理可得：

AB²＝AO²+BO²-2·AO·BO·cos(α-γ)

＝H²(ctanθ²+ctanω²-2ctanθ·ctanω·cos(α-γ)) (16)

(4)利用建筑遮挡后的可见阴影长l4计算建筑物高度：

∠BCO＝180°-δ (18)

δ＝∠CBO+∠BOC＝∠CBO+α-γ (19)

由式(6)、(7)可得：

∠CBO＝δ-α+γ (20)

由正弦定理可得：

求得建筑高度H为：

从以上公式可以得出，已知l₁，建筑高度仅和θ有关；已知l₂，建筑高度仅和ω有关；已知l₃，建筑高度和ω、θ、α、γ四个角度有关；已知l₄，建筑高度和ω、θ、α、γ、δ五个角度有关。在高分辨率影像成像过程中，每景影像的范围很小，成像时间很短，而我国又处于中低纬度，所以可以认为整景影像上每个像元的ω、θ、α、γ角度大小都相等，因此前式(13)、(14)和(17)式可以简写成：

H＝l₁×k₁ (23)

H＝l₂×k₂ (24)

H＝l₃×k₃ (25)

由于Google Earth影像无法直接获取ω、θ、α、γ等参数，可通过影像上的部分已知建筑高度反算θ和ω角度，α和γ角可通过在影像上量算图2中的BO和AO角获得。即同一景影像上各个建筑的高度和其阴影全长、建筑高差引起的建筑顶面角点像点位移、以及屋顶角点成像点与其阴影点的连线长均成正比。式(23)、(24)、(25)中的比例系数可运用区域中已知的某建筑的高度反推，即本发明先利用公式(23)-(25)中的一个计算出相应系数k，然后再利用求得的系数k，计算目标建筑的高度。对于一幢建筑来说，只要四类特征线中的一类能够准确判读(第四类并不常用)，即可运用相应的转换模型换算成建筑高度。

实际量算中应注意区分Google Earth上每景高分辨率影像的边界，在每景影像上分别计算ω、θ、α、γ等角度参数和k₁，k₂，k₃系数。也就是说，本发明的方法还包括，如果采用的影像为Google Earth图像，则提取所述影像的每景边界，并且，如果当影像边界不易或无法区分时，将影像划分成一定边长的正方形网格，每个网格求一组拍摄参数以及相应系数。例如，网格的边长设定为5km为宜，当相邻网格间的系数值产生突变时，说明两个网格处于不同的景上，这时应仔细分辨两景影像间的镶嵌线，并严格区分镶嵌线两侧系数值。

在高分辨率遥感影像上，受影像成像参数和地面复杂程度的影响，四类特征线段量取的难易和准确程度不一，为了快速准确获取建筑实际高度，在影像上选择最易准确判读的特征点和特征线，然后基于相应计算模型进行建筑高度的计算。

图5示出了基于本实施例中的方法所获得的南京城市三维重心的变化。从图中可以看出，实际上三维重心的平面投影与二维重心并不是重合的，而是存在一定的偏差。因此，在进行城市规划和土地利用规划时，利用三维重心来为商业中心和基础设施建设项目选址比利用二维重心来选址要更加科学，因为三维重心更能体现出人口的聚集度。

虽然上面结合本发明的优选实施例对本发明的原理进行了详细的描述，本领域技术人员应该理解，上述实施例仅仅是对本发明的示意性实现方式的解释，并非对本发明包含范围的限定。实施例中的细节并不构成对本发明范围的限制，在不背离本发明的精神和范围的情况下，任何基于本发明技术方案的等效变换、简单替换等显而易见的改变，均落在本发明保护范围之内。

Claims (10)

1.一种基于影像确定城市三维重心的方法，其特征在于，所述方法包括如下步骤：

步骤1，获取目标城市的一个或多个目标区域的影像以及所述目标区域的城市地籍和/或地形数据；

步骤2，对所获取的影像进行预处理；

步骤3，基于所获取的影像以及所述目标区域的城市地籍和/或地形数据提取所述目标区域的各个建筑的轮廓；

步骤4，基于所获取的影像以及所述目标区域的城市地籍和/或地形数据提取所述目标区域的各个建筑的高度；

步骤5，以预定尺度绘制覆盖所述目标区域的正方形格网，所述正方形格网具有若干正方形标准单元；

步骤6，基于所获得的所述目标城市的各个建筑的轮廓和高度数据，采用所绘制的所述目标区域的正方形格网将各建筑划分在所述正方形格网的一个个的标准单元内；

步骤7，计算所述正方形格网的每个标准单元内的建筑总体积，从而计算出每个标准单元内的平均建筑高度；

步骤8，基于所述正方形格网的每个标准单元内的平均建筑高度，生成所述目标区域的建筑高度等值线；

步骤9，基于所生成的所述目标区域的建筑高度等值线获得所述目标区域的三维立体模型；

步骤10，沿水平面方向对所述三维立体模型进行分割，以获得与水平面平行的若干建筑薄层；

步骤11，分别计算每个建筑薄层的重心的平面坐标(x,y)和高度坐标(z)；

步骤12，基于各个建筑薄层重心的三维坐标及各个薄层的体积计算各薄层的组合体的重心的三维坐标。

2.根据权利要求1所述的确定城市三维重心的方法，其特征在于，所述步骤4中，提取所述目标区域的各个建筑的高度的步骤包括：

步骤4-1，获取所述目标区域内若干已知建筑物的高度数据；

步骤4-2，基于所述影像，对于所述目标区域中的每个建筑，评估该建筑的多条特征线的提取难度；

步骤4-3，对于每个建筑，基于对多条特征线的提取难度的评估，选取最低难度的特征线；

步骤4-4，从所述影像中提取选定的特征线；

步骤4-5，基于所述目标区域内已知的建筑物高度数据反算所述影像的摄取参数；

步骤4-6，基于所选定特征线的长度、所述摄取参数，计算所述目标区域内，每个未知建筑的高度。

3.根据权利要求2所述的确定城市三维重心的方法，其特征在于，所述方法还包括判断所述影像是单景影像还是多景影像的镶嵌，如果所述影像为多景影像的镶嵌，则在所述高度提取的步骤中，分别针对各景影像进行参数计算和高度提取。

4.根据权利要求1所述的确定城市三维重心的方法，其特征在于，

所述步骤11中，每个建筑薄层的重心的平面坐标等于该建筑薄层的形心的坐标，每个建筑薄层的重心的高度坐标为：

其中，i为“建筑薄层”所在的层数；h为“建筑薄层”的厚度，即两条建筑高度等值线数值之差。

5.根据权利要求4所述的确定城市三维重心的方法，其特征在于，所述目标区域的三维重心的坐标为：

式中，为城市三维重心的坐标，x_i、y_i、z_i为第i个“建筑薄层”的三维重心的坐标，v_i为第i个“建筑薄层”的体积，等于“建筑薄层”面积与厚度h的乘积。

6.根据权利要求1所述的确定城市三维重心的方法，其特征在于，所述方法还包括计算所述目标城市的三维重心的转移距离，计算公式如下：

其中，D_t+1为三维空间重心转移距离，为基于在第一时间获得的所述目标城市的影像计算所得的所述目标城市的三维重心的坐标， 为基于在第二时间获得的所述目标城市的影像计算所得的所述目标城市的三维重心的坐标。

7.根据权利要求6所述的确定城市三维重心的方法，其特征在于，所述方法还包括计算三维重心转移距离在xoy平面上的投影距离D′_t+1：

。

8.一种基于权利要求1-7中任意一项所述方法获得的城市三维重心的应用，其特征在于，所述应用包括：基于所获得的城市三维重心为商业中心进行选址，使得所述商业中心位于所述三维重心的预定距离之内。

9.一种基于权利要求1-7中任意一项所述方法获得的城市三维重心的应用，其特征在于，所述应用包括：基于所获得的三维重心为基础设施建设项目选址服务。

10.一种基于权利要求1-7中任意一项所述方法获得的城市三维重心的应用，其特征在于，所述应用包括：基于所获得的三维重心进行城市发展方向和重点发展区位的选择。