基于分级栅格化和改进细化算法的轨迹数据路网生成研究

蔡 林，李英冰，邹子昕

(武汉大学 测绘学院，湖北 武汉 430079)

当前在城市中的出租车每天都生成海量轨迹数据，数据中蕴含车辆状态、人类行为等信息，这些信息对城市规划、交通建设等有重要应用价值[1]，从轨迹数据中挖掘居民出行特征等信息[2]、路网等城市建设信息[3]已成为研究的热点。

近年我国城市道路建设迅速发展，但路网信息更新慢，现势性无法满足人们日益增长的出行需求。传统道路更新技术，如测量技术，周期较长且花费成本高，而从遥感影像中提取路网相较于测量技术，周期短，效率高，但相比于从轨迹数据获取路网信息，其过程较复杂、数据成本较高。轨迹数据不仅能反映路网实时变化，且数据成本低，利用轨迹数据进行挖掘与分析，能实现对某一区域路网几何信息的获取，使构建或实时更新道路地图成为可能，如文献[4]利用轨迹地图匹配技术，从路段周边局部范围内的轨迹数据中提取并更新相关道路信息。

目前从轨迹数据中提取路网的方法主要有3种：第一种是将轨迹点连接成轨迹线，进而融合成道路，如文献[5]利用新轨迹与旧轨迹基于Delaunay三角网进行融合获取路网。这类方法是以旧轨迹线为基础获取新道路中心线，其精度会受到原有轨迹线的影响。第二种是聚类算法，将轨迹数据进行聚类，从中提取中心线获取路网，如文献[6]使用一种滚动式聚类算法完成了道路拓扑的生成和交叉口的识别。这类方法在轨迹数据量多的情形下，算法复杂，处理速度较慢。第三种是将轨迹数据转为栅格数据，通过轨迹点落入栅格点中个数确定栅格点的像素值，将栅格图像中骨架线作为道路中心线，如文献[7]在此基础上从轨迹数据中构建了路网，但丢失原始轨迹的连通性，对路网拓扑提取造成困难，面对海量轨迹数据，使用这类方法可实现一群轨迹点的聚类处理，而非单个点，从而减少轨迹数据处理时间，但这类方法提取路网，没有顾及细化后交叉点细节变化及其提取，且忽略栅格点中轨迹点的速度、方向等有效信息，也没有顾及近邻车道、数据缺失对提取骨架线的影响。

笔者基于第三种方法设计了一种分级栅格化方法，不仅提升栅格化速度，且保留利用轨迹点方向、速度等信息；利用图像膨胀算法将近邻车道合并、断开路段连接；使用改进Zhang细化算法对栅格路网图像进行细化，提取道路交叉口、中心线，最后生成矢量路网。

1 轨迹数据的索引与分级栅格化

轨迹数据生成路网栅格图像的传统过程为[8]：①确定栅格图像中栅格大小和研究区域；②建立轨迹点和栅格对应关系；③通过统计落入相应栅格中的轨迹点个数确定栅格点的像素值，从而生成栅格图像。

针对传统方法为了追求计算速度而牺牲方向、速度等大量有用信息，笔者提出一种分级索引的栅格化方法。面对海量轨迹数据，为提高栅格化效率，该方法对轨迹点落入的选取区域进行分级处理，为保留和利用栅格中的信息，建立方向、速度和时间索引，如图1所示。

图1 索引方案

首先根据选取区域大小对其进行分级处理，将选取区域分为3级，分别对应图1中的Div、Grid以及Cell。对于提取道路中心线，Cell大小与车辆一致，约5 m×5 m，然后对Cell建立如图1中的方向索引、时间索引、速度索引，其中方向索引将方向分为8个不同范围，0～7代表不同方向，时间索引是以小时为单位，0～23代表不同时间段，速度索引是根据轨迹点速度大小将速度划分为(0～9)10个不同范围。

其次轨迹数据以天为单位进行处理，对于每个属于选取区域Div的轨迹点，先记录所属于的Grid，然后寻找对应Grid中轨迹点所属的Cell，对应Cell内轨迹点统计数加一，且将轨迹点的方向、时间、速度信息分别落入Cell的方向、时间、速度索引的对应范围内统计数也加一，这样不仅能记录Cell中的轨迹点个数，而且能统计出Cell中不同方向、时段、速度范围内的轨迹点个数。

最后将多天同一Cell内不同轨迹点个数相加，同时将不同索引范围内统计数相加，采用多天的轨迹数据，能够补全部分缺少的道路信息。利用Cell速度索引中轨迹点速度信息和方位索引中轨迹点方位信息，去除Cell统计数中速度属于高低速范围或方向不合理的轨迹点，具有过滤噪声的作用。若Cell中存在轨迹点，则像素值为0，不存在则为255，可得栅格路网图像。

分级栅格化有两个优点:①整个区域被分为不同级别，只对轨迹点落入Cell的前一级做计算，不针对整个选取区域做计算，速度相较于传统栅格化会有所提高。②轨迹点的速度、方向等信息保留于栅格点中，且能利用其判断轨迹点是否合理，能去除噪声点，从而提高栅格路网图像的准确性。但分级栅格化有两个要点：①根据选取区域分级，区域越大，分级次数越多，最小的区域级别一定为Cell，Cell与Div之间的Grid可根据区域分为多级；②Cell大小需控制，太大细节信息易丢失，太小近邻车道明显，会加大中心线提取难度。

2 快速膨胀与改进细化算法在路网提取中的应用

获取的栅格路网图像中存在数据缺失导致路段断开和近邻车道的情形，断开路段需连接，提取道路中心线要进行近邻车道合并，而膨胀算法可以使路段向周围扩张，同时将路段连接和合并车道。一种简单快速的膨胀算法为依次扫描图像中的像素点，寻找灰度值为255的像素点，分析该像素周围8个像素点，仅要8个像素点中存在灰度值为0的像素点，就设置该像素点灰度值为0。

经膨胀算法处理后，路网图像基本显示路网特征，但路段为多像素宽度，提取路段较复杂。而用细化算法处理能保持图像中黑色部分的拓扑结构，处理后表现为单像素宽度，图中骨架线为道路中心线，笔者采用一种改进Zhang细化算法。将膨胀后图像二值化，对于具有道路信息点(像素值为0)，其像素值设为1，称为目标点；对无道路信息点(像素值为255)，将其像素值设为0，称为背景点，可得简单二值图像。

改进Zhang细化算法流程如下：

(1)寻找以目标点为中心的8邻域，记中心点为P1，其邻域的8个点逆时针绕中心点分别记为P2，P3，…，P9，P2在P1上。标记同时满足下列条件的像素点：

①2≤N(P1)≤6

②S(P1)=1或B(P1)∈{65,5,20,80,13,22,52,133,141,54}

③P2×P4×P8=0或S(P2)!=1

④P2×P4×P6=0或S(P4)!=1

(2)去除所有满足标记的点，且将目标点灰度值赋值为0，背景点灰度值赋值为255。

其中N(P1)为P的8邻域内非零邻点的个数，S(P1)为以P为中心，其周围8邻域内任一点开始，依次沿逆时针方向，从0到1的变化次数。B(P1)为P1邻域点的二进制编码值S，计算公式如式(1)所示[9]，Pi的值为其像素值。

(1)

改进Zhang细化算法使用5×5邻域，在加大邻域范围基础上改变原有Zhang算法细化条件，并添加了文献[8]改进条件。某一交叉口3种细化算法效果如图2所示，其中浅色部分为膨胀后交叉口效果，深色部分为细化后效果。图2(a)为Zhang细化算法，细化后非单像素宽度，路段交叉口处效果较差；图2(b)为文献[8]改进细化算法，仅将图2(a)改为单像素宽度；图2(c)为笔者改进Zhang细化算法，细化后为单像素宽度且交叉口较前两种更为直观、符合现实状况，细化效果更好。

图2 细化算法效果图

图2(a)、图2(b)、图2(c)相比，交叉口拓扑结构发生改变，主要是因细化过程不仅取决前次迭代结果，且取决本次迭代中已处理过像素点分布情况，与迭代次数、使用邻域大小、细化条件有很大关系。三者本质都是把多像素宽度路段细化为一条中心线，虽然交叉口等局部细节走向可能不同，但路段整体走向都是一致的。细化后图像中依然存在着一些孤立像素点与较少像素点组成的像素段，是由轨迹数据发生漂移或者轨迹数据发生缺失等噪声造成。对于提取路段来说，这部分像素是无价值的，需过滤。

3 从细化路网图中提取路段构建矢量路网

像素点所在行列为其坐标，提取路段、路段交叉点、端点均用像素坐标表示。提取路段前需进行路段端点、交叉点的获取，其思路为判断某目标像素点(像素值为0)3×3邻域内的其他目标像素点个数，若有2个，则该像素点为端点，3个以上为交叉点。

路段提取分两类，首先是端点至端点或交叉点间的路段，如图3中数字为1部分的路段，其中数字为3部分是交叉点，因细化后图像为单像素宽度，则可从各端点出发，依次在3×3邻域寻找并储存具有灰度值的点，再以储存的点为开始，重复寻找，当遇到其他交叉点或端点停止，可提取端点至端点或交叉点的路段。其次是交叉点至交叉点之间的路段，如图3中数字为2部分，从某一交叉点3×3邻域内的各个像素值为0的点出发，类似于从端点出发提取路段，但当遇到另一交叉点停止，则可提取交叉点至交叉点间所有路段。

图3 路段提取示意图

基于端点和交叉点提取在细化后图像中的全部路段，但由较少像素点组成的路段是非必要的，这类路段主要是由图像中间隔较近交叉点(实际在现实中同为一道路交叉点)之间的路段组成，需将这些路段删除，然后将这部分交叉点合并为一个交叉点，如图4(b)所示，深色部分为交叉点。

经过试验，设置删除长度小于或等于4个像素的路段，其大多都为间隔较近交叉点间的路段，路段删除前后部分交叉口变化情况如图4，图中a、b为两种不同的路口情况。

图4 部分交叉口路段变化

删除较短路段后，首先如果一个交叉点同时与3条路段以上相连接，则该点一定为真正交叉点，按照这种思路对这类交叉点获取，如图4(a)，经过此步，完成了部分真正交叉点的获取；其次若在交叉点周围存在其他交叉点，如图4(b)删除后，只需要从某一交叉点出发，在该交叉点邻域范围内，搜索出其他交叉点，利用这些交叉点插值出一个新交叉点，将在这些路段中交叉点替换为新交叉点，能将一个交叉点与多条路段相连接，获取合并后的真正交叉点，如图4(b)所示，交叉点合并后路段较好连接，一个真交叉点至少与3条以上路段相连接，一条路段的首尾点必然是端点或者真正交叉点。

构建矢量路网前，需用DP(Dougls-Pucker)算法进行路段压缩[10]，DP算法使原始曲线路段更加接近真实的道路和保持光滑形态，但必须选择合适的阈值。本文选取的阈值为像素宽度0.8，代表真实距离为4 m左右，将提取路段进行了DP算法压缩。由于压缩后路段的首尾点一定是端点或者真正交叉点，只需要将每条路段依次生成矢量路段，通过调用GDAL库，生成了shpfile格式的路网矢量图文件。

4 结果分析

笔者利用成都市某5 km×5 km区域20161101～20161108共8天滴滴出租车轨迹数据进行实验，共计66 196 629条，其中最小栅格Cell为5 m×5 m,传统栅格化方法需427 678 ms，分级栅格化只需400 756 ms，两者相比，分级栅格化速度有所提升。

选取区域经分级栅格化、膨胀、细化且去噪后进行矢量化处理的路网如图5所示。

图5 不同类型路网图像

图5(a)与图5(b)相比，图5(b)中道路中心线明显，近邻车道合并，断开路段相连，更有利于从中提取道路中心线，但道路中心线是非单像素宽度。图5(c)为该区域细化且去噪后的路网图像，道路宽度为单像素宽度，与图5(b)中路网拓扑信息基本一致。图5(d)是图5(a)的基础上提取矢量路网，两图中道路中心线基本一致，说明了矢量路网的生成过程是正确的。

为说明生成路网的有效性，笔者选取矢量路网中两个区域与相应百度地图路网相比，如图6所示，矢量路网中的路段与真实地图路网中路段相贴近，路段拓扑结构特征基本一致，但存在部分路段缺失，主要是因为轨迹数据没有覆盖这些区域或区域内轨迹点数量较少造成的。之后将提取真正交叉点的经纬度坐标与真实的道路交叉点坐标进行对比，选取图6区域2的A、B、C 3点，对比结果如表1所示。通过对比，坐标差在30 m左右，这类误差主要由在提取交叉点中将交叉口化为一个点造成的，更进一步证明本文从出租车轨迹数据中提取路网方案是可行的。

图6 百度地图与实验结果对比图

表1 区域2中的交叉点坐标对比

ID(a)中坐标(°)(b)中坐标(°)差距/mA104.078 992,30.713 222104.078 714,30.713 18331B104.078 943,30.715 717104.078 804,30.715 85322C104.076 161,30.714 127104.076 049,30.714 32926

5 结论

针对目前利用轨迹数据实现构建矢量路网的问题，提出利用分级栅格化和改进Zhang细化算法实现路网快速提取；利用细化后路网栅格图完成对路口端点、交叉点、矢量路网的提取；最后以成都市某5 km×5 km区域8天的滴滴出租车数据进行了实验论证。但由于路网复杂性，该方法还存在一定的不足，道路复杂的环形立交路段提取效果较差，路网构造中车道识别也至关重要，研究中将车道合并，提取仅为道路中心线。这些研究有待于今后进一步探索。