CN108776727B - 一种基于出租车轨迹数据的道路几何特征提取方法 - Google Patents

Abstract

本发明涉及一种基于出租车轨迹数据的道路几何特征提取方法，首先对出租车轨迹数据预处理；接着计算每个轨迹点的行驶方向角；然后遍历原始轨迹点集合，获取轨迹点P _i 的相似轨迹点集合；接着遍历相似轨迹点集合，计算得到轨迹点P _i 的偏移距离，得到P _i 偏移后的坐标；接着完成轨迹点集合的遍历，将P _i 偏移后的坐标更新为轨迹点的当前坐标，计算所有轨迹点偏移距离的平均值，若小于等于阈值则将偏移后的轨迹点转成轨迹线，对每条轨迹线作缓冲区并融合为面要素进行二值化，利用ArcScan工具提取道路网L。本发明通过对轨迹点进行偏移，有效的避免了将邻近不同行驶方向的轨迹聚为一束的情况，提取得到的路网更能反映真实道路结构，具有较高的实用价值。

Description

一种基于出租车轨迹数据的道路几何特征提取方法

技术领域

本发明涉及地理信息系统与智能交通研究领域，特别是一种基于出租车轨迹数据的道路几何特征提取方法。

背景技术

城市路网数据作为基础地理信息数据的关键组成部分，是城市智能交通建设的重要数据源，准确、精细的路网数据能够为城市交通拥堵治理、城市规划提供有效的决策支持，获取现势性强、详细程度高的路网信息是当前交通领域重要的任务。随着移动端定位技术的成熟，安装了GPS定位装置的出租车在短时间内便能获取大量关于城市路网的时空数据，具有实时、低成本、针对性强的优势，从轨迹数据中提取精细化的路网信息，成为交通地理信息领域重要的研究内容。

利用GPS时空轨迹数据提取道路网信息已经有相关探索性工作，主要分为三类，第一类主要采用栅格化方法将轨迹数据栅格化提取道路骨架线；第二类方法采用轨迹聚类的方法；第三类方法基于轨迹增量的方法。2009年，Cao LiLi等提出了模拟轨迹点受到周围轨迹段吸引的路网提取模型，但现有方法提取的路网信息精细化程度并不高，原因主要在于没有考虑双向行驶道路上存在两个行驶方向的车流轨迹，以至于道路交叉口信息提取不完善。

发明内容

有鉴于此，本发明的目的是提出一种基于出租车轨迹数据的道路几何特征提取方法，解决了当前路网信息提取精细化程度低的问题，能够反映路网真实的拓扑关系，并实现道路行驶方向信息的提取。

本发明采用以下方案实现：一种基于出租车轨迹数据的道路几何特征提取方法，具体包括以下步骤：

步骤S1：对GPS轨迹点数据进行预处理，得到轨迹点集合Traj；

步骤S2：计算轨迹点集合Traj中每个轨迹点的行驶方向角；

步骤S3：遍历轨迹点集合Traj，并获取轨迹点集合Traj中的P_i的相似轨迹点集合Traj_s；

步骤S4：遍历P_i的相似轨迹点集合Traj_s，计算Traj_s中每个轨迹点P_s对P_i的欧氏距离d_i和影响权重ω_s，并以此计算轨迹点P_i的偏移距离Δd_i，经坐标分解得到轨迹点P_i偏移后的坐标；重复步骤S3与S4，直至完成轨迹点集合Traj的遍历，得到轨迹点集合Traj中所有轨迹点偏移后的坐标；

步骤S5：将轨迹点集合Traj中所有的轨迹点的坐标更新为偏移后的坐标，计算所有轨迹点偏移距离的平均值Δd，并判断Δd是否大于阈值δ，若大于则返回步骤S2，否则进入步骤S6；

步骤S6：将偏移后的轨迹点转成轨迹线，对每条轨迹线作缓冲区并融合为面要素进行二值化，并提取道路网L。

进一步地，步骤S1中的所述预处理包括以下步骤：

步骤S11：将原始轨迹点数据的经纬度信息转换成墨卡托平面坐标；

步骤S12：以研究区范围的区划图作为底图，将原始轨迹点数据与底图进行叠置，剔除研究区域以外的轨迹点；

步骤S13：根据出租车载客状态剔除非载客轨迹点数据；

步骤S14：将轨迹点数据中的时间转换成以秒为单位的数据，即将轨迹中“yyyy-MM-dd hh:mm:ss”格式的时间转换成以秒为单位的整数。

步骤S15：对直线距离超过阈值σ的轨迹段进行轨迹点加密。

进一步地，步骤S2中，轨迹点P_i的行驶方向角θ定义如下：

轨迹点P_i行驶方向角θ：从轨迹点P_i-1的指北方向线起依顺时针方向至P_i-1和P_i+1连线间的夹角；

如果P_i为轨迹的起点，则θ为从轨迹点P_i的指北方向线起依顺时针方向至P_i和P_i+1连线间的夹角；如果P_i为轨迹的终点，则θ为从轨迹点P_i-1的指北方向线起依顺时针方向至P_i-1和P_i连线间的夹角。

由此，计算Traj中轨迹点P_i的行驶方向角的方法包括以下步骤：

步骤S21：计算轨迹点P_i-1和P_i+1的X、Y坐标差值Δdx、Δdy；

步骤S22：判断Δdx、Δdy的取值范围，行驶方向角θ的计算采用下式：

进一步地，步骤S3中相似轨迹点的定义为：空间位置邻近、行驶方向角相近的轨迹点。

由此，步骤S3中，获取轨迹点集合Traj中的P_i的相似轨迹点集合Traj_s具体包括以下步骤：

步骤S31：以轨迹点P_i的X、Y坐标为圆心，r为半径，获取半径范围内的轨迹点集合T；

步骤S32：遍历轨迹点集合T，剔除与P_i行驶方向角差值的绝对值小于等于角度阈值的轨迹点后，得到P_i的相似轨迹点集合Traj_s。

进一步地，步骤S4中影响权重的定义为：假设同一条道路上的相似轨迹点之间存在吸引作用，用于反映轨迹点受到相似轨迹点影响大小的因子称为影响权重，权重与距离呈负相关，距离越远，影响权重越小。

由此，步骤S4中，计算Traj_s中每个轨迹点P_s对P_i的欧氏距离d_i和影响权重ω_s，并以此计算轨迹点P_i的偏移距离Δd_i，经坐标分解得到轨迹点P_i偏移后的坐标具体包括以下步骤：

步骤S41：过轨迹点P_i作平行于P_i-1-P_i+1线段的直线l，以直线l为y轴,选定车辆行驶方向为y坐标轴正方向，并以P_i为原点O，建立O-XY平面直角坐标系；如果P_i为轨迹的起点，则过P_i作平行于P_i-P_i+1线段的直线l；如果P_i为轨迹的终点，则过P_i作平行于P_i-1-P_i线段的直线l；

步骤S42：遍历P_i的相似轨迹点集合Traj_s,计算中轨迹点P_s在O-XY坐标系中X轴的坐标值d'_s，按d'_s的正负分别累加取平均得到轨迹点P_i的偏移距离Δd_i，其中Δd_i的计算采用下式：

式中，d'_j、d'_k分别为d'_s≥0和d'_s＜0时集合Traj_s中轨迹点的X轴坐标，ω_j、ω_k为Traj_s中轨迹点的影响权重，m、n分别为相似轨迹点集合Traj_s中d'_s≥0和d'_s＜0的轨迹点数量；

其中，轨迹点的影响权重ω_s的计算采用以下公式：

式中，d_s是P_i和P_s的欧式距离，ω_s为相似轨迹点P_s对P_i的影响权重，k为常数；

步骤S43：在墨卡托平面坐标系下对Δd_i进行分解，计算轨迹点P_s偏移后的坐标。

进一步地，所述步骤S5中，所有轨迹点偏移距离的平均值Δd的计算方法采用下式：

式中，M为轨迹点集合Traj中的轨迹点数量。

进一步地，步骤S6具体包括以下步骤：

步骤S61：剔除没有相似轨迹点的离散点并分割轨迹；

步骤S62：将轨迹点集合Traj中的轨迹点按车辆ID和时间顺序连接成轨迹线；

步骤S63：在每条轨迹线两侧作距离为γ的缓冲区，并将缓冲区融合为单个面要素；

步骤S64：将融合后的面要素二值化，并利用ArcScan工具提取道路网L。

与现有技术相比，本发明有以下有益效果：

1、本发明充分考虑并合理解决了当前路网信息提取精细化程度低的问题，与Cao提出的方法相比，仅考虑轨迹点与其相似轨迹点间的相互影响，忽略轨迹点偏移过程中的“弹力”作用，使得轨迹点偏移过程更为简单，避免了复杂的数学运算，有效避免现有道路几何特征提取算法将同一条道路上相反行驶方向的轨迹识别为同一条道路的情况，本发明所提取的路网为道路双向几何特征结构，能够反映路网真实的拓扑关系。

2、本发明对于不同道路结构的城市路网都能准确提取出精细的道路几何特征，得到的实验结果几何特征上贴近于真实的道路网，并能准确提取方向道路的行驶方向信息。

附图说明

图1为本发明实施例的方法流程示意图。

图2为本发明实施例的计算轨迹点行驶方向角的示意图。

图3为本发明实施例的获取相似轨迹点的示意图。

图4为本发明实施例的轨迹点偏移的示意图。

具体实施方式

下面结合附图及实施例对本发明做进一步说明。

应该指出，以下详细说明都是例示性的，旨在对本申请提供进一步的说明。除非另有指明，本文使用的所有技术和科学术语具有与本申请所属技术领域的普通技术人员通常理解的相同含义。

需要注意的是，这里所使用的术语仅是为了描述具体实施方式，而非意图限制根据本申请的示例性实施方式。如在这里所使用的，除非上下文另外明确指出，否则单数形式也意图包括复数形式，此外，还应当理解的是，当在本说明书中使用术语“包含”和/或“包括”时，其指明存在特征、步骤、操作、器件、组件和/或它们的组合。

如图1所示，本实施例提供了一种基于出租车轨迹数据的道路几何特征提取方法，具体包括以下步骤：

步骤S1：对GPS轨迹点数据进行预处理，得到轨迹点集合Traj；

步骤S2：计算轨迹点集合Traj中每个轨迹点的行驶方向角；

步骤S3：遍历轨迹点集合Traj，并获取轨迹点集合Traj中的P_i的相似轨迹点集合Traj_s；

步骤S4：遍历P_i的相似轨迹点集合Traj_s，计算Traj_s中每个轨迹点P_s对P_i的欧氏距离d_i和影响权重ω_s，并以此计算轨迹点P_i的偏移距离Δd_i，经坐标分解得到轨迹点P_i偏移后的坐标；重复步骤S3与S4，直至完成轨迹点集合Traj的遍历，得到轨迹点集合Traj中所有轨迹点偏移后的坐标；

步骤S5：将轨迹点集合Traj中所有的轨迹点的坐标更新为偏移后的坐标，计算所有轨迹点偏移距离的平均值Δd，并判断Δd是否大于阈值δ，若大于则返回步骤S2，否则进入步骤S6；

步骤S6：将偏移后的轨迹点转成轨迹线，对每条轨迹线作缓冲区并融合为面要素进行二值化，并提取道路网L。

在本实施例中，步骤S1中的所述预处理包括以下步骤：

步骤S11：将原始轨迹点数据的经纬度信息转换成墨卡托平面坐标，这一步是为了便于在后续步骤中轨迹点间欧式距离的计算与表示。转换方法如下；

式中，B为维度，L为经度，X表示横坐标，Y表示纵坐标，

R＝6378137m为地球平均半径。

步骤S12：以研究区范围的区划图作为底图，将原始轨迹点数据与底图进行叠置，剔除研究区域以外的轨迹点；

步骤S13：根据出租车载客状态剔除非载客轨迹点数据，这是由于出租车载客状态下通常在路段上行驶，其轨迹点能更好的反映路网结构；

步骤S14：将轨迹点数据中的时间转换成以秒为单位的数据，即将轨迹中“yyyy-MM-dd hh:mm:ss”格式的时间转换成以秒为单位的整数，并将结果乘以10；这一步是为了后续加密的轨迹点提供可用的时间戳范围，转换方法如下：

t＝(dd*86400+hh*3600+mm*60+ss)*10。

步骤S15：对直线距离超过阈值σ的轨迹段进行轨迹点加密，确保轨迹点集合中相同车辆ID的相邻轨迹点的距离间隔在15-25米之间，这一步是为了后续轨迹点偏移之后，道路转弯区域不同车辆ID的轨迹能够汇聚在一起，轨迹点加密的具体方法如下：

x_i＝x_s+i*x_a，y_i＝y_s+i*y_a(1≤i≤m-1)；

式中，x_s,x_e,y_s,y_e分别为前后相邻轨迹点的x、y坐标，m为相邻轨迹点间待插入轨迹点的数量，通过相邻轨迹点间欧式距离除以20向下取整确定，x_i，y_i，t_i为待插入轨迹点的x、y坐标和时间信息，t_s为相邻轨迹点中时间较小的时间戳。

在本实施例中，步骤S2中，轨迹点P_i的行驶方向角θ定义如下：

轨迹点P_i行驶方向角θ：从轨迹点P_i-1的指北方向线起依顺时针方向至P_i-1和P_i+1连线间的夹角；

如果P_i为轨迹的起点，则θ为从轨迹点P_i的指北方向线起依顺时针方向至P_i和P_i+1连线间的夹角；如果P_i为轨迹的终点，则θ为从轨迹点P_i-1的指北方向线起依顺时针方向至P_i-1和P_i连线间的夹角。

由此，计算Traj中轨迹点P_i的行驶方向角的方法包括以下步骤：

步骤S21：计算轨迹点P_i-1和P_i+1的X、Y坐标差值Δdx、Δdy；

步骤S22：判断Δdx、Δdy的取值范围，行驶方向角θ的计算采用下式：

其中，Δdx、Δdy分别为相邻轨迹点X、Y坐标之差；如果轨迹点P_i为一段轨迹的起点，则P_i的行驶方向角由P_i与P_i+1的坐标确定；如果轨迹点P_i为一段轨迹的终点，则P_i的行驶方向角由P_i-1与P_i的坐标确定。

如图2所示，O₁、O₂、O₃、O₄、O₅分别为轨迹点，位于轨迹段中间部分的轨迹点的行驶方向角由前后相邻轨迹点确定，例如轨迹点O₃的行驶方向角等于从轨迹点O₂的指北方向线起依顺时针方向至O₂和O₄连线间的夹角，即θ₃；起始轨迹点O₁的行驶方向角由轨迹点O₁、O₂确定，终止轨迹点O₅的行驶方向角由轨迹点O₄、O₅确定

在本实施例中，步骤S3中相似轨迹点的定义为：空间位置邻近、行驶方向角相近的轨迹点。

由此，步骤S3中，获取轨迹点集合Traj中的P_i的相似轨迹点集合Traj_s具体包括以下步骤：

步骤S31：以轨迹点P_i的X、Y坐标为圆心，r为半径，获取半径范围内的轨迹点集合T；

步骤S32：遍历轨迹点集合T，剔除与P_i行驶方向角差值的绝对值小于等于角度阈值的轨迹点后，得到P_i的相似轨迹点集合Traj_s。

如图3所示，轨迹点O₃的相似轨迹点获取方法为遍历以O₃为圆心，R为20m范围内，与O₃行驶方向角差值的绝对值≤5°的轨迹点即O₃的相似轨迹点。

在本实施例中，步骤S4中影响权重的定义为：假设同一条道路上的相似轨迹点之间存在吸引作用，用于反映轨迹点受到相似轨迹点影响大小的因子称为影响权重，权重与距离呈负相关，距离越远，影响权重越小。

由此，步骤S4中，计算Traj_s中每个轨迹点P_s对P_i的欧氏距离d_i和影响权重ω_s，并以此计算轨迹点P_i的偏移距离Δd_i，经坐标分解得到轨迹点P_i偏移后的坐标具体包括以下步骤：

步骤S41：过轨迹点P_i作平行于P_i-1-P_i+1线段的直线l，以直线l为y轴,选定车辆行驶方向为y坐标轴正方向，并以P_i为原点O，建立O-XY平面直角坐标系；如果P_i为轨迹的起点，则过P_i作平行于P_i-P_i+1线段的直线l；如果P_i为轨迹的终点，则过P_i作平行于P_i-1-P_i线段的直线l；

如图4所示，轨迹点O₃建立O-XY坐标系是以O₃为原点，过轨迹点O₃且平行于O₂-O₄线段的直线为Y轴，与轨迹运动方向呈锐角的方向作为Y轴正方向，通过建立右手直角坐标系确定X轴；轨迹段起点O₁建立O-XY坐标系则以O₁为原点，过O₁且平行于O₁-O₂线段的直线为Y轴，与轨迹运动方向呈锐角的方向作为Y轴正方向，通过建立右手直角坐标系确定X轴；轨迹终点O₅建立O-XY坐标系的方式与轨迹点起点O₁同理。

步骤S42：遍历P_i的相似轨迹点集合Traj_s,计算中轨迹点P_s在O-XY坐标系中X轴的坐标值d'_s，按d'_s的正负分别累加取平均得到轨迹点P_i的偏移距离Δd_i，其中Δd_i的计算采用下式：

式中，d'_j、d'_k分别为d'_s≥0和d'_s＜0时集合Traj_s中轨迹点的X轴坐标，ω_j、ω_k为Traj_s中轨迹点的影响权重，m、n分别为相似轨迹点集合Traj_s中d'_s≥0和d'_s＜0的轨迹点数量；

其中，轨迹点的影响权重ω_s的计算采用以下公式：

式中，d_s是P_i和P_s的欧式距离，ω_s为相似轨迹点P_s对P_i的影响权重，k为常数；由此得到：

步骤S43：在墨卡托平面坐标系下对Δd_i进行分解，计算轨迹点P_s偏移后的坐标。

在本实施例中，所述步骤S5中，所有轨迹点偏移距离的平均值Δd的计算方法采用下式：

式中，M为轨迹点集合Traj中的轨迹点数量。

较佳的，在本实施例提供了一种相似轨迹点集合的快速查询方法：将研究区域虚拟划分为20×20m单元的网格，为每一个网格进行唯一编号，将同一个网格内的所有轨迹点集中存储在Key为网格编号的HashMap中，遍历轨迹点集合中的轨迹点时，根据轨迹点的坐标信息确定其所在网格，提取以此网格为中心的九宫格所包含的轨迹点进行遍历，其他网格轨迹点则忽略，这样能极大的缩短相似轨迹点查询的时间。

在本实施例中，步骤S6具体包括以下步骤：

步骤S61：剔除没有相似轨迹点的离散点并分割轨迹；

步骤S62：将轨迹点集合Traj中的轨迹点按车辆ID和时间顺序连接成轨迹线；

步骤S63：在每条轨迹线两侧作距离为γ的缓冲区，并将缓冲区融合为单个面要素；

步骤S64：将融合后的面要素二值化，并利用ArcScan工具提取道路网L。

以上所述仅为本发明的较佳实施例，凡依本发明申请专利范围所做的均等变化与修饰，皆应属本发明的涵盖范围。

Claims (5)

1.一种基于出租车轨迹数据的道路几何特征提取方法，其特征在于：包括以下步骤：

步骤S1：将原始轨迹点数据的经纬度坐标转换成墨卡托平面坐标；以研究区范围的区划图作为底图，将轨迹点数据与底图进行叠置，剔除研究区域以外的轨迹点；根据出租车载客状态剔除非载客轨迹点数据，并将轨迹点数据中的时间转换成以秒为单位的数据；对相邻2个轨迹之间直线距离超过阈值σ的轨迹段进行轨迹点加密；得到轨迹点集合Traj；

步骤S2：计算轨迹点集合Traj中每个轨迹点的行驶方向角；

步骤S3：遍历轨迹点集合Traj，以轨迹点P_i的X、Y坐标为圆心，r为半径，获取半径范围内的轨迹点集合T；遍历轨迹点集合T，剔除与P_i行驶方向角差值的绝对值小于等于角度阈值的轨迹点后，得到P_i的相似轨迹点集合Traj_s；

步骤S4：遍历P_i的相似轨迹点集合Traj_s，计算Traj_s中每个轨迹点P_s对P_i的欧氏距离d_i和影响权重ω_s，并以此计算轨迹点P_i的偏移距离Δd_i，经坐标分解得到轨迹点P_i偏移后的坐标；重复步骤S3与S4，直至完成轨迹点集合Traj的遍历，得到轨迹点集合Traj中所有轨迹点偏移后的坐标；

步骤S5：将轨迹点集合Traj中所有的轨迹点的坐标更新为偏移后的坐标，计算所有轨迹点偏移距离的平均值Δd，并判断Δd是否大于阈值δ，若大于则返回步骤S2，否则进入步骤S6；

步骤S6：将偏移后的轨迹点转成轨迹线，对每条轨迹线作缓冲区并融合为面要素进行二值化，并提取道路网L。

2.根据权利要求1所述的一种基于出租车轨迹数据的道路几何特征提取方法，其特征在于：步骤S2中，计算Traj中轨迹点P_i的行驶方向角的方法包括以下步骤：

步骤S21：计算轨迹点P_i-1和P_i+1的X、Y坐标差值Δdx、Δdy；

步骤S22：判断Δdx、Δdy的取值范围，行驶方向角θ的计算采用下式：

3.根据权利要求1所述的一种基于出租车轨迹数据的道路几何特征提取方法，其特征在于：步骤S4中，计算Traj_s中每个轨迹点P_s对P_i的欧氏距离d_i和影响权重ω_s，并以此计算轨迹点P_i的偏移距离Δd_i，经坐标分解得到轨迹点P_i偏移后的坐标具体包括以下步骤：

步骤S41：过轨迹点P_i作平行于P_i-1-P_i+1线段的直线l，以直线l为y轴,选定车辆行驶方向为y坐标轴正方向，并以P_i为原点O，建立O-XY平面直角坐标系；如果P_i为轨迹的起点，则过P_i作平行于P_i-P_i+1线段的直线l；如果P_i为轨迹的终点，则过P_i作平行于P_i-1-P_i线段的直线l；

步骤S42：遍历P_i的相似轨迹点集合Traj_s,计算Traj_s中轨迹点P_s在O-XY坐标系中X轴的坐标值d_s'，按d_s'的正负分别累加取平均得到轨迹点P_i的偏移距离Δd_i，其中Δd_i的计算采用下式：

式中，d_j'、d_k'分别为d_s'≥0和d_s'＜0时集合Traj_s中轨迹点的X轴坐标，ω_j、ω_k为Traj_s中轨迹点的影响权重，m、n分别为相似轨迹点集合Traj_s中d_s'≥0和d_s'＜0的轨迹点数量；

其中，轨迹点的影响权重ω_s的计算采用以下公式：

式中，d_s是P_i和P_s的欧式距离，ω_s为相似轨迹点P_s对P_i的影响权重，k为常数；

步骤S43：在墨卡托平面坐标系下对Δd_i进行分解，计算轨迹点P_s偏移后的坐标。

4.根据权利要求1所述的一种基于出租车轨迹数据的道路几何特征提取方法，其特征在于：所述步骤S5中，所有轨迹点偏移距离的平均值Δd的计算方法采用下式：

式中，M为轨迹点集合Traj中的轨迹点数量。

5.根据权利要求1所述的一种基于出租车轨迹数据的道路几何特征提取方法，其特征在于：步骤S6具体包括以下步骤：

步骤S61：剔除没有相似轨迹点的离散点并分割轨迹；

步骤S62：将轨迹点集合Traj中的轨迹点按车辆ID和时间顺序连接成轨迹线；

步骤S63：在每条轨迹线两侧作距离为γ的缓冲区，并将缓冲区融合为单个面要素；

步骤S64：将融合后的面要素二值化，并利用ArcScan工具提取道路网L。

Images