CN115311854B - 一种基于数据融合的车辆时空轨迹重构方法 - Google Patents

一种基于数据融合的车辆时空轨迹重构方法 Download PDF

Info

Publication number
CN115311854B
CN115311854B CN202210868085.XA CN202210868085A CN115311854B CN 115311854 B CN115311854 B CN 115311854B CN 202210868085 A CN202210868085 A CN 202210868085A CN 115311854 B CN115311854 B CN 115311854B
Authority
CN
China
Prior art keywords
track
vehicle
time
space
road section
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN202210868085.XA
Other languages
English (en)
Other versions
CN115311854A (zh
Inventor
任刚
诸赛
邓玥
马景峰
曹奇
李大韦
吴辰旸
晁晓宁
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Southeast University
Original Assignee
Southeast University
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Southeast University filed Critical Southeast University
Priority to CN202210868085.XA priority Critical patent/CN115311854B/zh
Publication of CN115311854A publication Critical patent/CN115311854A/zh
Application granted granted Critical
Publication of CN115311854B publication Critical patent/CN115311854B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • G—PHYSICS
    • G08—SIGNALLING
    • G08G—TRAFFIC CONTROL SYSTEMS
    • G08G1/00—Traffic control systems for road vehicles
    • G08G1/01—Detecting movement of traffic to be counted or controlled
    • G08G1/0104—Measuring and analyzing of parameters relative to traffic conditions
    • G08G1/0125—Traffic data processing
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S19/00—Satellite radio beacon positioning systems; Determining position, velocity or attitude using signals transmitted by such systems
    • G01S19/38—Determining a navigation solution using signals transmitted by a satellite radio beacon positioning system
    • G01S19/39—Determining a navigation solution using signals transmitted by a satellite radio beacon positioning system the satellite radio beacon positioning system transmitting time-stamped messages, e.g. GPS [Global Positioning System], GLONASS [Global Orbiting Navigation Satellite System] or GALILEO
    • G01S19/393—Trajectory determination or predictive tracking, e.g. Kalman filtering
    • G—PHYSICS
    • G08—SIGNALLING
    • G08G—TRAFFIC CONTROL SYSTEMS
    • G08G1/00—Traffic control systems for road vehicles
    • G08G1/01—Detecting movement of traffic to be counted or controlled
    • G08G1/017—Detecting movement of traffic to be counted or controlled identifying vehicles
    • G—PHYSICS
    • G08—SIGNALLING
    • G08G—TRAFFIC CONTROL SYSTEMS
    • G08G1/00—Traffic control systems for road vehicles
    • G08G1/01—Detecting movement of traffic to be counted or controlled
    • G08G1/052—Detecting movement of traffic to be counted or controlled with provision for determining speed or overspeed
    • Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
    • Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
    • Y02T—CLIMATE CHANGE MITIGATION TECHNOLOGIES RELATED TO TRANSPORTATION
    • Y02T10/00—Road transport of goods or passengers
    • Y02T10/10—Internal combustion engine [ICE] based vehicles
    • Y02T10/40—Engine management systems

Landscapes

  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Engineering & Computer Science (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Remote Sensing (AREA)
  • Chemical & Material Sciences (AREA)
  • Analytical Chemistry (AREA)
  • Computer Networks & Wireless Communication (AREA)
  • Traffic Control Systems (AREA)

Abstract

本发明公开了一种基于数据融合的车辆时空轨迹重构方法，方法包括：(1)浮动车时空轨迹重构算法，利用卡口AVI数据和GPS数据的关联性，简化地图匹配流程，基于分段三次Hermit插值法，对已知轨迹的缺失或未知的轨迹段进行修复，得到完整的浮动车时空轨迹曲线，并从中提取交通流基本图参数，作为全时空车辆轨迹重构算法的输入。(2)全时空车辆轨迹重构算法，在Newell冲击波模型中引入车辆顺序变量，用车辆顺序的变化表征车辆的超车行为，推导其与空间、时间和累计车辆数的数学关系，将轨迹重构问题转化为车辆顺序的求解问题，实现目标时段内路段上所有车辆的时空轨迹重构，有助于准确把握道路交通运行状态，实时而精准的进行交通管理与控制。

Description

一种基于数据融合的车辆时空轨迹重构方法
技术领域
本发明涉及交通规划与管理领域，具体涉及一种基于数据融合的城市快速路车辆时空轨迹重构方法。
背景技术
根据来源划分，当前主流的交通数据可以分为传统的固定位置检测数据和移动检测数据两大类。固定位置检测数据以视频监测卡口的自动车辆识别数据(AVI)为主，移动检测数据包含浮动车GPS数据、手机信令数据等。然而，由于资金限制，卡口视频设备安装位置十分有限，大多数安装在交叉口上下游位置，并且由于环境和设备本身的影响，采集数据会存在较多的错误和缺失。浮动车GPS数据也受限于设备和环境遮挡的影响，检测数据存在采样率低、采样时间间隔不统一、轨迹漂移和缺失问题。此外，由于交通数据行业分散，不同源交通数据监管缺乏共享机制与统一标准，导致了“交通信息孤岛”现象，单一数据源的挖掘分析效率低下，结果的实际应用价值不大。现有研究中大多数针对单一交通数据源进行轨迹进行挖掘分析，采用数据融合方法的研究较少，且结果有待商榷。因此，如何将异源数据有效融合，提高数据质量和使用效率，进一步深入挖掘潜在的交通信息，是一个亟待解决的问题。
车辆轨迹重构问题可以分为两类。第一类是针对已知部分轨迹的车辆轨迹重构，可以描述为：基于离散或不完整的车辆轨迹，利用有效的插值或递推算法对缺失或未知的轨迹段进行修复，从而恢复出完整平滑的车辆时空轨迹曲线。第二类是针对未知轨迹车辆的轨迹重构，可以描述为：基于已知的完整车辆轨迹，利用交通模型原理或数据驱动方法推断出未知轨迹车辆的时空轨迹。尽管移动检测数据可以在短时间间隔内提供空间上连续的信息，但更提供该类数据的车辆到道路交通流中渗透率较低，只占整个交通流的一小部分，重构出全部车辆的时空轨迹仍然是一个巨大挑战。典型的轨迹重构方法是利用卡口AVI数据或固定线圈检测数据对车辆轨迹进行估计，重构出的车辆轨迹的准确性与卡口或线圈检测器的布设密度存在较大的关联性。而随着交通数据采集技术的多元化发展，近年来许多学者尝试采用多源数据融合的方法对车辆轨迹进行重构，有效克服了单一数据源低采样率和估计精度不高等缺陷，这也成为轨迹重构的重要研究方向之一。
针对上述问题，面向我国城市道路交通状态评估与优化的实际需求，本发明以城市快速路为研究对象，提出一种基于数据融合的城市快速路车辆时空轨迹重构方法。
发明内容
本发明目的：在于提供一种基于数据融合的车辆时空轨迹重构方法，利用路网GIS数据、浮动车GPS数据和视频卡口AVI数据，以城市快速路为研究范围，在路段层面重建车辆的时空轨迹图谱，从而把握准确的交通流状态，为交通管控和决策提供参考。
为实现以上功能，本发明设计一种基于数据融合的车辆时空轨迹重构方法，针对目标路段上各浮动车，执行以下步骤1-步骤4，针对目标路段上所有车辆，执行步骤5，获得目标路段上所有车辆的行驶轨迹，进一步构建考虑超车行为的车辆时空轨迹重构模型，完成目标路段的交通管控：
步骤1：采集目标路段各卡口的AVI数据，以及目标路段上各浮动车的GPS数据，分别构建AVI数据集、GPS数据集，以浮动车的车牌号为标识，分别从GPS数据集中提取各浮动车的车牌号所对应的轨迹点，并按照每辆浮动车的车牌号所对应的轨迹点的采样时间排序，获得每辆浮动车离散的轨迹点序列；
步骤2：基于地图匹配方法中的几何匹配法，利用AVI数据和GPS数据中的关联特征，确定各浮动车的轨迹点所对应的匹配道路，完成浮动车轨迹的修正；
根据球面距离计算公式计算各浮动车的轨迹点之间的距离，完成浮动车轨迹的时空坐标转换，获得离散的浮动车时空轨迹；
步骤3：基于分段三次hermit插值法，对步骤2所获得的浮动车时空轨迹的缺失段进行插值，获得连续且平滑的浮动车时空轨迹；
步骤4：基于三角型基本图模型，并基于步骤3所获得的浮动车时空轨迹，提取交通流参数，交通流参数包括自由流速度、交通流激波波速、阻塞密度；
步骤5：在LWR模型和Newell模型中引入车辆顺序变量，用车辆顺序的变化表示车辆的超车行为，通过推导车辆顺序变量与空间、时间、累计车辆数的数学关系，将车辆轨迹重构问题转化为车辆顺序的求解问题，并以各交通流参数为输入、以所有车辆在目标路段上的重构轨迹为输出，构建考虑超车行为的车辆时空轨迹重构模型，并应用车辆时空轨迹重构模型，完成目标路段的交通管控。
有益效果：相对于现有技术，本发明的优点包括：
本发明基于插值算法的浮动车时空轨迹重构方法。研究面向数据融合的GPS轨迹修正方法，利用卡口AVI数据和GPS数据在时间、空间和身份标识等纬度上的关联性简化地图匹配的候选路径集，提高了GPS轨迹修正的效率。通过时空坐标转换得到离散的浮动车时空轨迹。基于分段三次Hermit插值等四种插值方法，结合实际时空轨迹特征进行修正，对浮动车时空轨迹缺失段进行插值补全，可得到所有浮动车的时空轨迹图谱。
本发明提出考虑超车行为的车辆时空轨迹重构算法。基于LWR模型和Newell模型基本原理，研究累计车辆数、空间和时间三者的数学关系，引入车辆顺序变量，用车辆顺序的变化表征车辆的超车行为，推导了车辆空间位置与车辆顺序相对于时间的数学关系式，建立了考虑超车行为的车辆时空轨迹重构模型。提出的轨迹重构算法可以对研究路段和研究时段内所有车辆的时空轨迹进行重建，并进行可视化展示，进一步可提取交通流特征参数，实现交通流状态的精准把控，为交通控制与管理提供参考价值。
附图说明
图1是基于数据融合的车辆时空轨迹重构方法的流程图；
图2是根据本发明实施例提供的车辆轨迹数据提取的流程图；
图3(a)是根据本发明实施例提供的GPS轨迹修正流程图；
图3(b)是根据本发明实施例提供的时空坐标转换图；
图4是根据本发明实施例提供的交通流激波波速提取过程图
图5是根据本发明实施例提供的全时空轨迹重构场景图；
图6是根据本发明实施例提供的基于分段三次Hermit插值得到的浮动车时空轨迹图；
图7是根据本发明实施例提供的上下游卡口累计流量函数分段线性化结果图；
图8是根据本发明实施例提供的浮动车观测轨迹与重构轨迹对比图；
图9是根据本发明实施例提供的全时空车辆轨迹重构结果图；
图10是根据本发明实施例提供的浮动车轨迹重构三种误差结果分布图；
图11是根据本发明实施例提供的从时空轨迹图提取交通流特征参数示例图。
具体实施方式
下面结合附图对本发明作进一步描述。以下实施例仅用于更加清楚地说明本发明的技术方案，而不能以此来限制本发明的保护范围。
本发明实施例提供的一种基于数据融合的车辆时空轨迹重构方法，参照图1，针对目标路段上各浮动车，执行以下步骤1-步骤4，针对目标路段上所有车辆，执行步骤5，获得目标路段上所有车辆的行驶轨迹，进一步构建考虑超车行为的车辆时空轨迹重构模型，完成目标路段的交通管控：
步骤1：采集目标路段各卡口的AVI数据，以及目标路段上各浮动车的GPS数据，分别构建AVI数据集、GPS数据集，以浮动车的车牌号为标识，分别从GPS数据集中提取各浮动车的车牌号所对应的轨迹点，并按照每辆浮动车的车牌号所对应的轨迹点的采样时间排序，获得每辆浮动车离散的轨迹点序列；
GPS数据集中包括每辆浮动车的轨迹点，且GPS数据集中每辆浮动车的轨迹点按时间戳排序，以浮动车的车牌号为标识，从GPS数据中提取每辆浮动车的车牌号所对应的轨迹点，并按照每辆浮动车的车牌号所对应的轨迹点的采样时间排序，参照图2，获得每辆浮动车离散的轨迹点序列Pprobe＝(p1，p2，...，pi，...，pn)，其中p1，p2，...，pi，...，pn表示车牌号所对应的各轨迹点，n为轨迹点总数。
步骤2：基于地图匹配方法中的几何匹配法，利用AVI数据和GPS数据中的关联特征，确定各浮动车的轨迹点所对应的匹配道路，完成浮动车轨迹的修正；
根据球面距离计算公式计算各浮动车的轨迹点之间的距离，完成浮动车轨迹的时空坐标转换，获得离散的浮动车时空轨迹；
步骤2的具体步骤如下：
步骤2.1：分别取AVI数据集中上游卡口和下游卡口AVI数据的车牌和时间字段，判断浮动车j的车牌号是否在其中且连续通过上游卡口和下游卡口，若是，则认为浮动车j通过目标路段；否则认为浮动车j未通过目标路段，取下一辆浮动车的车牌号重复上述判断过程，直至所有浮动车判断完毕；
步骤2.2：参照图3(a)，通过卡口AVI数据约束确定匹配路段后，针对轨迹点pi，基于几何匹配法获得轨迹点pi在目标路段的道路中心线上的投影点ei，取下一个轨迹点重复上述匹配过程，直至所有轨迹点匹配完毕，则输出修正后的轨迹点序列(e1，e2，...，ei，...，en)，其中e1，e2，...，ei，...，en表示修正后的各轨迹点，n为轨迹点总数；
步骤2.3：参照图3(b)，根据球面距离计算公式，计算轨迹点pi、pj之间的距离dij如下式：
dij≈cos-1(sinβi·sinβj+cosβi·cosβj·cosΔα)·R
式中，Δα＝αi-αj，其中αi、αj分别为轨迹点pi、pj的经度，βi、βj分别为轨迹点pi、pj的纬度，R为地球半径；
构建以目标路段的上游卡口为起点的离散浮动车时空轨迹序列(t1，x1)，(t2，x2)，...，(tn，xn)，其中t1，t2，...，tn、x1，x2，...，xn分别为轨迹点p1，p2，...，pn的时间、空间位置。
步骤3：基于分段三次hermit插值法，对步骤2所获得的浮动车时空轨迹的缺失段进行插值，获得连续且平滑的浮动车时空轨迹；
步骤3的具体步骤如下：
步骤3.1：根据分段三次Hermit插值法，在每个分段区间上，插值函数都是三次多项式并且要求在插值点处插值函数连续，插值函数的一阶导连续。对于浮动车时空轨迹序列(t1，x1)，(t2，x2)，...，(ti，xi)，...，(tn，xn)，在区间[ti-1，ti]上的插值条件如下式：
其中，f(t)为插值函数，其表达式以分段三次多项式表示为如下形式：
步骤3.2：构建Hermit插值法的基函数如下式：
式中，Φ′i-1(t)、Φ′i(t)均为三次多项式，称为三次Hermit插值多项式的基函数；
步骤3.3：对步骤3.2所构建的Hermit插值法的基函数两边求导：
结合步骤3.1所述插值条件，可得下式：
记Δi＝ti-ti-1，基函数的表达式如下式：
针对上式两边求导可得下式：
其中c和d表达式如下式：
求解获得的表达式如下式：
分别求解获得Φi-1(t)、Φi(t)的表达式如下式：
多项式xi(t)的表达式如下式：
通过求解每个插值区间上的多项式表达式，获得分段三次Hermit多项式表达式X(t)；
步骤3.4：分段三次Hermit插值只有在被插值函数在插值节点处的函数值和导数值已知时才可以使用，对于本发明的场景，每辆浮动车的每个轨迹点都对应记录一个瞬时速度值，该值可以近似认为是当前轨迹点对应时空函数的导数值。并且车辆的瞬时速度均大于等于零，两个轨迹点之间插值曲线的形状基本被限定，大多数情况下不会出现三次样条插值类似的震荡现象。而对于少部分的曲线震荡，添加与分段三次样条曲线类似的约束条件，t＝ti时刻插值点空间位置x(ti)与x(ti-1)进行比较，并且x(ti)要满足小于右边已知点对应的空间位置；
针对两轨迹点之间时空轨迹曲线的曲线震荡，引入约束条件如下式，进行插值结果修正：
xi＝max{x(ti)，x(ti-1)}
步骤3.5：针对浮动车j的时空轨迹序列(t1，x1)，(t2，x2)，...，(ti，xi)，...，(tn，xn)，取一个轨迹点(ti，xi)作为误差计算对象，用时空轨迹序列中剩余的n-1个轨迹点进行插值补齐，插值的时间步长取Δt＝0.1s，误差指标采用平均绝对误差MAE和均方根误差RMSE，浮动车j的平均绝对误差MAEj和均方根误差RMSEj如下式：
式中，为第i个轨迹点的插值结果；
目标路段上所有浮动车的平均绝对误差MAE和均方根误差RMSE如下式：
式中，m为目标路段上的浮动车总数。
步骤4：基于三角型基本图模型，并基于步骤3所获得的浮动车时空轨迹，提取交通流参数，交通流参数包括自由流速度、交通流激波波速、阻塞密度；
步骤4的具体步骤如下：
步骤4.1：自由流速度vf指的是在目标道路上不受其它车辆干扰，并且由驾驶员主观意愿自由选择的行驶速度。基于Hermit插值法，获得目标路段上每辆浮动车的时空轨迹曲线，从浮动车i进入目标路段的起始时刻开始，以步长Δt为窗口向由平滑，对于第j个窗口，对应的自由流速度vj如下式：
根据自由流速度vj获得浮动车i在目标路段上速度的分布情况，基于帕累托法则，取80％分位的速度作为浮动车i自由流速度的估计值，记为vi-80th；重复上述过程获得每辆浮动车自由流速度的估计值，目标路段的自由流速度的估计值如下式：
式中，m为目标路段上的浮动车总数；
步骤4.2：交通流激波指的是在拥挤乃至阻塞状态下与车流前进方向相反的波，类似于管道内的水流突然受阻时的后涌现象。本发明根据目标路段上各浮动车的时空轨迹曲线中拥堵发生处拟合获得交通流激波波速w；
参照图4，从各条浮动车时空轨迹图上找到进入拥挤状态的拐点，用这些拐点来拟合出交通流激波，对应拟合直线的斜率代表波速大小。由于交通发生拥挤状态时并不一定会产生阻塞，因此本文以阻塞状态下浮动车时空轨迹的拐点来拟合交通流激波。
步骤4.3：在已知自由流速度vf和交通激波波速w的情况下，要计算阻塞密度kj，还需要得到道路的饱和流率qm，借鉴《通行能力手册》(HCM 2010)关于快速路基本通行能力描述，基本通行能力描述的是道路和交通都处于理想条件下，标准车以最小的车头间距连续行驶的理想交通流，在单位时间内通过道路断面的最大车辆数。本发明用饱和流率qm代替基本通行能力。那么根据三角型基本图的特点，阻塞密度kj的表达式如下式：
式中，km为车流达到饱和流率对应最佳密度，qm为目标路段的饱和流率。
步骤5：在LWR模型和Newell模型中引入车辆顺序变量，用车辆顺序的变化表示车辆的超车行为，通过推导车辆顺序变量与空间、时间、累计车辆数的数学关系，将车辆轨迹重构问题转化为车辆顺序的求解问题，并以各交通流参数为输入、以所有车辆在目标路段上的重构轨迹为输出，构建考虑超车行为的车辆时空轨迹重构模型，并应用车辆时空轨迹重构模型，完成目标路段的交通管控。
具体分析方法如下：
步骤5.1：Newell运动波模型的特性推导。完整的LWR模型表达式如下：
式中，q(x，t)，k(x，t)，u(x，t)分别是流率q、密度k和速度u关于时间和空间位置的函数；u＝ue(k)为平衡状态下交通流速度-密度函数关系。
Newell模型则在LWR模型的基础上基于交通流三角型基本图进行简化，以时间t和空间位置x对应的累计车辆数N(x，t)为状态变量，交通流守恒方程可以转化为：
式中，Q(k)＝min{vfk，w(kj-k)}，是流率q关于密度k的函数。
参照图5，针对长度为l的目标路段，目标路段为均质路段，且目标路段的上游卡口和下游卡口各设置了AVI设备，用于采集通过车辆的数据，目标路段交通流特征符合三角型基本图。
以t＝0时刻，x＝l处的车辆为参考，将[0，t]时段路段上游卡口通过的累计车辆数记为U(t)＝N(0，t)-N0，其中N0为t＝0时路段内部初始的车辆数，路段下游卡口通过的累计车辆数记为D(t)＝N(l，t)。在Newell的理论中，在t＞0且x∈[0，l]时的任一位置的累计车辆数N(x，t)可由Daganzo提出的变分理论来求解，累计车辆数N(x，t)公式如下：
上式表明累计车辆数N(x，t)由上游条件(方程中的第一项)或者下游条件(方程中的第二项)决定。用Nf(x，t)表示非拥挤状态下的累计车辆数，用Nc(x，t)表示拥挤状态下的累计车辆数，可以得到：
则累计车辆数N(x，t)公式转化为：
N(x，t)＝min{Nf(x，t)，Nc(x，t)}
为方便计算，假设N(x，t)在几何上是光滑曲面，即N(x，t)关于x和t可微。显然，对任意时刻t，N(x，t)、U(t)和D(t)都是单调不减的。用Nf(x，t)表示非拥挤状态下的累计车辆数，用Nc(x，t)表示拥挤状态下的累计车辆数，N(x，t)可以表示为：N(x，t)＝min{Nf(x，t)，Nc(x，t)}；
结合交通流基本图，得到下式：
对于任意时刻t，在几何上累计车辆数N(x，t)关于x的图像应当为Nf(x)和Nc(x)曲线叠加后取下方的部分，两曲线交界于位置[x1，x2]处(x1＝x2时两曲线交于一点)，且两交点连线对应的斜率为临界值km。整个目标路段被划分为三段：
(1)当x∈[0，x1]时，N(x，t)＝Nf(x，t)＜Nc(x，t)，交通流处于非拥挤状态；
(2)当x∈[x1，x2]时，N(x，t)＝Nf(x，t)＝Nc(x，t)，交通流处于临界状态，此时流率达到最大值；
(3)当x∈[x2，l]时，N(x，t)＝Nc(x，t)＜Nf(x，t)，交通流处于拥挤状态。
考虑极端情况下，当x1＝x2＝0时，N(x，t)＝Nc(x，t)，整个目标道路交通流均为拥挤状态；当x1＝x2＝l时，N(x，t)＝Nf(x，t)，整个目标道路交通流处于非拥挤状态。
将xf(n，t)和xc(n，t)记为Nf(x，t)和Nc(x，t)的反函数，则Nf(x，t)、Nc(x，t)表达式转换为：
将目标路段分成三段来分析：
(1)当x∈[0，x1]时，N(x，t)＝Nf(x，t)＜Nc(x，t)且xf(n，t)＜xc(n，t)，则有n＝N(x(n，t)，t)＝Nf(xf(n，t)，t)，故x(n，t)＝xf(n，t)，即对于这部分路段，满足x(n，t)＝min{xf(n，t)，xc(n，t)}；
(2)当x∈[x1，x2]时，N(x，t)＝Nf(x，t)＝Nc(x，t)且xf(n，t)＝xc(n，t)，则有n＝N(x(n，t)，t)＝Nf(xf(n，t)，t)＝Nc(xc(n，t)，t)，故x(n，t)＝xf(n，t)＝xc(n，t)，该部分路段同样满足x(n，t)＝min{xf(n，t)，xc(n，t)}；
(3)当x∈[x2，l]时，N(x，t)＝Nc(x，t)＜Nf(x，t)且xf(n，t)＞xc(n，t)，则有n＝N(x(n，t)，t)＝Nc(xc(n，t)，t)，故x(n，t)＝xc(n，t)，即对于这部分路段，满足x(n，t)＝min{xf(n，t)，xc(n，t)}。
综上所述，对于N(x，t)的反函数x(n，t)，N(x，t)＝min{Nf(x，t)，Nc(x，t)}的形式同样适用，即：
x(n，t)＝min{xf(n，t)，xc(n，t)}
步骤5.2：车辆顺序变量引入。针对t＝0时，x＝l处的车辆，用δi(t)表示t时刻车辆i的顺序，对应时刻的车辆空间位置为xi(t)，则事实上δi(t)等同于[0，t]时段内，x＝xi(t)处的累计车辆数，即：
δi(t)＝N(xi(t)，t)
考虑x(n，t)为N(x，t)的反函数，则有下式：
xi(t)＝x(δi(t)，t)
用和/>表示Nf(x，t)和Nc(x，t)的反函数，则有下式：
则x(n，t)＝min{xf(n，t)，xc(n，t)}可以写为下式形式：
则可得下式：
由上式可知，在已知交通基本图参数的情况下，要求得交通流中任一车辆i在任一时刻t的空间位置xi(t)，等同于计算任一时刻t车辆i的顺序δi(t)，通过比较拥挤和非拥挤情况下的车辆位置即可得到最终结果。
在实际的道路交通流中，车辆的换道、超车行为是普遍存在的。由于本发明只研究沿道路前进方向的位置变化，故对车辆的换道行为不展开研究，只考虑车辆在车流前进方向上顺序变化。用ei和di分别表示车辆进入目标道路上游卡口和离开下游卡口的时间，则当δi(t)＝δi(ei)＝δi(dt)时，表明车辆i在通过目标路段时没有被超越或者超越其它车辆，也可能超越该车的车辆数和该车辆超越的车辆数相等，但是这种情况对于路段较短的情况比较特殊，不作考虑。
用变量Ci来表示车辆i通过目标路段时顺序的变化量，即：
Ci＝δi(dt)-δi(ei)
当Ci＞0时，车辆i通过目标路段时被Ci辆车超越；当Ci＝0时，车辆i通过目标路段时顺序没有变化；当Ci＜0时，车辆i通过目标路段时超越了Ci辆车。
考虑较一般的情况，对于车辆i和车辆j，在时段[t1，t2]内，若δi(t1)＜δj(t1)，δi(t2)＞δj(t2)，则表明在该时段内车辆i被车辆j超越了，在空间位置上表现为xi(t1)＞xj(t1)，xi(t2)＜xj(t2)。因此，Ci和δi(t)可以用于刻画车辆的超车行为。
当交通流运行稳定时，可以对δi(t)随时间的变化进行简化，本发明用线性模型来表示δi(t)，即：
δi(t)＝αit+βi
当t＝ei时
δi(ei)＝N(xi(ei)，ei)＝N(0，ei)＝U(ei)+N0
当t＝di时
δi(di)＝N(xi(di)，di)＝N(l，di)＝D(di)
根据线性函数斜率计算公式，可以得到：
则化简后δi(t)的表达式为：
步骤5.3：考虑超车行为的车辆时空轨迹重构算法。根据步骤S52，重构车辆时空轨迹相当于求得任一时刻t交通流中任一车辆i空间位置xi(t)。因为时间t是连续变量，则xi(t)也是连续的，当计算xi(t)的时间步长Δt足够小时，xi(t)可以近似认为是车辆i的时空轨迹。在计算xi(t)之前，需要先求得以下参数：
(1)交通流参数vf，kj和w；
(2)t＝0时刻路段内部初始车辆数N0；
(3)路段上下游卡口累计车辆数关于时间的函数U(t)和D(t)；
(4)任一时刻t车辆i的顺序δi(t)。
首先交通流参数vf，kj和w可以由步骤4的方法利用重构的浮动车时空轨迹得到。本发明用每辆车通过上下游卡口时，对应累计车辆数差的平均值n0来估计t＝0时刻路段内部初始车辆数N0，即
式中，h为通过上下游卡口且均被记录数据的车辆数。
路段上下游卡口累计流量函数U(t)和D(t)实际上是分段函数，因为函数值对应的是车辆数，均为整数。为了方便计算，将U(t)和D(t)用分段线性差值方法进行连续化处理，将其转换成分段线性函数，使得每个时刻t都能对应一个函数值，D(t)与U(t)的处理方法相同。
计算任一时刻t车辆i的顺序δi(t)，求解出交通流非拥挤和拥挤状态下对应的位置和/>比较后可得最终的位置xi(t)。
综合上述分析，在得到交通流参数vf，kj，w和函数U(t)，D(t)后，车辆时空轨迹重构的算法步骤如下：
(1)计算车辆i通过目标路段时顺序的变化量Ci；
(2)将车辆i进入和离开目标路段的时间区间[ei，di]以步长Δt＝0.1s划分；
(3)当t＝μΔt时，计算此时顺序δi(t)：
其中，μ为整数，t＝μΔt∈[ei，di]；
(4)当t＝μΔt时，计算在非拥挤情况下，车辆i的位置/>
已知则
其中，U-1(t)为U(t)的逆函数；
(5)当t＝μΔt时，计算在拥挤情况下，车辆i的位置
已知 的表达式不能直接求得，可以将该式转化为函数方程：
其中，唯一未知数求解方法为二分法、牛顿法、梯度法中的一种；
(6)比较和/>大小，取小者作为t＝μΔt时车辆i的预测空间位置。因为不确定/>和/>的单调性，所以增加约束，将xi(μΔt)与上一个预测空间位置比较，取较大者，以满足xi(t)单调不减：
(7)返回(3)，计算t＝μΔt+Δt时刻的δi(t)、和/>直到t+Δt此时，可以得到车辆i以Δt＝0.1s为时间间隔的位置序列：xi(ei)，xi(ei+Δt)，xi(ei+2Δt)，...，xi(di)；
步骤5.4：轨迹重构误差验证。利用浮动车GPS数据集和卡口AVI数据集的关联性，考虑用重构得到的浮动车轨迹来进行车辆轨迹重构精度的评价。评价指标为平均绝对误差(MAE)和均方根误差(RMSE)，即比较每个时刻t对应的浮动车i的实际空间位置和重构轨迹的预测空间位置之差。
式中，xi(t)为t时刻浮动车i的实际空间位置，/>为t时刻浮动车i的预测空间位置。
此外，车辆重构轨迹误差还可用相对误差来表示，即浮动车i的实际时空轨迹和重构轨迹围成的面积与实际轨迹围成的面积之比，计算公式为：
式中，Si1和Si2分别表示浮动车i实际轨迹围成的面积和实际轨迹与重构轨迹围成的面积。因为Δt＝0.1s足够小，所以近似认为xi(t)和所围成的面积是以Δt为宽的多个矩形面积之和，可以作为积分的近似计算。
具体实施例：
1.数据准备
研究路段：深圳机场南路位于深圳宝安区宝安航城街道，是深圳宝安国际机场、深圳机场客运码头连接城市道路网的重要快速主干道，西起领航路，东至机场宝安立交，全长约4.552公里。道路主车道为双向六车道，辅道双向四车道，中间设有高架，长约2公里，道路全程最高限速80km/h。研究路段选取了机场南路由西向东方向的高架路段。路段在机场南路机场人行天桥与下角山人行天桥之间，两天桥均设有车辆身份识别监控摄像头，记录由西向东的过车数据。快速路段分别有一出口和入口，出口距离上游卡口约600米，入口距离上游卡口约1700米。
GPS数据：选取深圳市2016年9月1日晚高峰时段出租车GPS数据，平均采样间隔为15s。
AVI数据：选取深圳市2016年9月1日晚高峰时段研究路段上下游卡口AVI数据，卡口的具体信息如表1所示。
表1
卡口名称 经度 纬度 ID
机场人行天桥西往东 113.824647 22.614894 20507302
下角山人行天桥东侧西往东 113.844703 22.622033 20507303
GIS数据：从开源平台Open Street Map上下载对应于研究路段的GIS数据，并通过ArcGIS软件进行预处理。根据一般的通勤时间，将晚高峰时段定为17:30-19:30。获取在该时段内通过目标路段的浮动车GPS和AVI数据方法如下：
(1)以卡口ID为条件，筛选17:30-19:30内连续通过上下游卡口的AVI数据；
(2)对于浮动车数据，以车牌号遍历，判断车牌号是否存在于筛选后的AVI数据中，保留判断为“是”的浮动车数据，并以车辆到达上下游卡口的时间为条件进一步筛选；
(3)从筛选结果中保留路段内部轨迹点数大于5的浮动车数据，最终得到晚高峰时段内通过研究路段的浮动车GPS数据集。
数据筛选后进一步除噪，最后共得到浮动车GPS数据162条，AVI数据3042条，浮动车的渗透率约为5.3％。以15min为间隔统计晚高峰各时段内AVI记录的车流量，如表2所示。各时段内车辆数相差不大，17:30:00-18:30:00时段内车流量略高于18:30:00-19:30:00时段，在整个晚高峰时段内研究路段较为畅通。
表2
时段 车辆数 时段 车辆数
17:30:00-17:45:00 408 18:30:00-18:45:00 328
17:45:00-18:00:00 447 18:45:00-19:00:00 362
18:00:00-18:15:00 337 19:00:00-19:15:00 372
18:15:00-18:30:00 406 19:15:00-19:30:00 342
2.浮动车轨迹重构结果
经过数据处理后，得到晚高峰时段通过研究路段的浮动车GPS轨迹数据。根据步骤2的方法对GPS轨迹点进行经纬度的修正，将GPS轨迹点与研究路段相匹配。以上下游卡口的经纬度坐标为起终点，对每辆车的GPS轨迹进行时空坐标转换，得到转换后的浮动车时空轨迹数据集。基于分段三次Hermit插值法的浮动车轨迹重构误差如表3所示。根据表3中的误差结果，分段三次Hermit插值法的MAE和RMSE分别为24.546和41.010，误差较小。
表3
MAE RMSE
24.546 41.010
以分段三次Hermit插值法进行浮动车时空轨迹重构，得到162辆浮动车完整的时空轨迹曲线。为了更清楚地展示车辆的运行状态变化，以0.1s为窗口滑动，近似计算车辆的瞬时速度，并给对应轨迹点根据速度大小附上不同深浅的颜色，用热力图方法呈现，如图6所示。
从轨迹时空图上看，在研究时段内浮动车在整个交通流中分布较均匀，不存在轨迹数据在某一个时段内集中的情况，因此能够大致看出研究路段整体交通流在时空图上的分布情况。研究路段在整个晚高峰期间内交通较为畅通，大部分浮动车均能以接近自由流的速度在该路段行驶。其中，在18:00-18:25时段内，距离上游卡口约1750米处出现拥挤状态，车流以低速缓慢前进，并且在该时段内拥挤状态向上游传播至距离上游卡口约1500米处。路段距离上游1700米左右有一匝道接入，并且1700米后部分路段存在交织区，产生拥堵的可能原因是在该时段内时段，路段在此处入口匝道的车流量较大，匝道接入的交通流对主路上车流产生影响，导致主路车辆减速，产生了与车流前进方向相反的交通流激波。
3.基本图参数估计结果
(1)自由流速度估计
选取凌晨2016年9月1日0:00-2:00通过该路段的浮动车数据作为自由流速度估计的数据集。经过轨迹插值重构后，共得到105辆浮动车的轨迹数据。剔除离群数据后得到的各浮动车的80％分位速度估计值分布，平均值为18.70，即该路段自由流速度vf估计值为18.70m/s，约为67.32km/h，比道路限速略低。
(2)交通流激波波速估计
从图6来看，研究时段内没有出现明显的阻塞流，拥挤区域有少部分浮动车有停车行为，因此本发明从有明显停车行为的浮动车轨迹来拟合交通流激波波速，经计算，w≈4.15m/s。
(3)阻塞密度估计
研究路段西向东方向为3车道，限速为80km/h，假定该路段高峰时段内无中、大型车通行，根据HCM2010取单车道饱和流率为2100pcu/h，那么该路段的饱和流率qm＝6300pcu/h。已估计得到vf＝18.70m/s，w＝4.15m/s，则阻塞密度估计值为：
即每条车道平均阻塞密度约为171.7pcu/km。假设平均车长为5m，那么阻塞情况下平均车头间距约为5.79m。
4.全时空轨迹重构与误差分析结果
(1)数据集构建及输入参数设置
在轨迹冲重构之前，先基于筛选后的AVI数据构造模型的输入数据集，包含字段及对应示意如表4，数据集示例见表5。起始时刻设置为17:30:00，参考车为该时刻通过下游卡口的第一辆车。因此在17:30:00时，路段内部已经存在一部分初始车辆，故表5内有部分数据U(t)和ei字段数据有缺失。
表4
字段 含义
ID 车牌号
U(t) [0，t]时段路段上游卡口通过的累计车辆数
D(t) [0，t]时段路段下游卡口通过的累计车辆数
ei 车辆i经过上游卡口的时间
di 车辆i经过下游卡口的时间
表5
ID U(t) D(t) ei di
粤S****9 1 2016/9/1 17:30:00
粤B****8 2 2016/9/1 17:30:00
粤S****P 3 2016/9/1 17:30:03
… … … … …
粤S****6 22 87 2016/9/1 17:30:38 2016/9/1 17:32:40
… … … … …
粤B****2 2948 3042 2016/9/1 19:28:09 2016/9/1 19:29:58
基于输入数据集，可以得到路段上下游卡口累计流量函数U(t)和D(t)，以17：30：00为基准，用分段线性差值方法对函数U(t)和D(t)进行连续化处理，结果如图7所示。对于研究时段内任一时刻，均能求得U(t)和D(t)函数值，为车辆轨迹重构提供数据支撑。路段内部初始车辆数N0估计值为79辆。路段长度在GIS中求出为1990米。此外，研究路段的交通流基本图参数如表6所示。
表6
参数 vf(m/s) w(m/s) kj(pcu/m)
估计值 18.70 4.15 0.515
(2)全时空轨迹重构与误差结果
基于(1)中构造的数据集和输入参数，根据步骤5的轨迹重构算法对晚高峰时段内通过路段的所有车辆(不包括路段中间从出口离开和入口进入的车辆)进行轨迹重构。其中，计算在拥挤情况下车辆i的位置采用牛顿法求解，收敛误差设置为0.01m。
用重构的162辆完整的浮动车时空轨迹作为观测数据集对算法精度进行验证，所有浮动车时空轨迹的观测值和预测值如图8所示，整体上浮动车时空轨迹的观测值和预测值较为比较一致，但是每个时段的误差分布不够均匀。首先在17:30:00-17:40:00、17:50:00-18:00:00这两个时段内，浮动车在同一时刻对应空间位置的预测值小于观测值，并且在上下游卡口所记录的行程时间内，预测轨迹还没有达到下游卡口的位置。产生这一误差的可能原因在于，在非拥挤状态下，车辆的位置预测值与自由流速度有很大关系，而该时段内车辆的平均速度超过了估计的自由流速度。在18:00:00-18:20:00时段内，对于距离上游卡口1500-1750处的拥挤状态部分的轨迹重构，预测的车辆时空轨迹相对观测值在空间上有约100米的滞后，而在时间上较为一致。
因为研究时段内车辆数较多，时空轨迹较密集，为了更清楚地展示交通流运行状态变化，同样以0.1s为窗口滑动，近似计算车辆的瞬时速度，并给对应轨迹点根据速度大小附上不同深浅的颜色，用热力图方法呈现，如图9所示。
从结果来看，整个研究路段在晚高峰时段内较为畅通。其中，18:30-19:30时段内总体交通流比17:30-18:30时段内畅通。在18:00:00-18:25:00时段内，距离上游卡口约1700-2000米处出现了拥挤流，平均车辆运行速度较小，车流以小于5m/s的速度行进，主要原因在于附近入口匝道接入车流对主路车流造成影响。除此之外，从预测效果来看，在18:00之前和18:38之后预测的车辆轨迹绝大部分都滞后于观测轨迹，因为在对应的旅行时间内，预测车辆位置还没有达到下游卡口处。原因可能是路段的自由流速度估计值偏小，因为估计的自由流速度为18.7m/s(67.32km/h)，小于路段限速(80km/h)，在道路畅通状态下，按估计的自由流速度预测的位置会滞后于实际车辆位置。
图10为所有浮动车时空轨迹重构的相对误差RE、MAE和RMSE的计算结果分布情况，误差的平均值分别为REmean＝7.10％，MAEmean＝71.37，RMSEmean＝84.24，从各浮动车轨迹预测的误差分布来看，整体处于较稳定的水平，没有出现异常的起伏，有一半左右的浮动车轨迹预测值的误差处在较低水平，另一半的车辆的误差在平均值附近有较大起伏。其中MAE和RMSE的误差相对较大，意味着对任一时刻平均每辆车的空间位置估计值的平均误差有100多米。而RE误差接近10％，说明预测的轨迹精度达到了合理水平，轨迹重构算法具有较高的可信度和精度。
(3)交通流特征参数提取示例
结合时空轨迹图，除了从宏观上能看出交通流运行状态变化情况，从微观角度也能提取图中任意一块时空区域的交通流特征参数，如图11所示。
选取任一辆车的时空轨迹，通过判断它从上方超越其它车辆轨迹的交叉点个数，可以得知该车在路段内行驶时超车数以及产生超车行为的大致时间和位置。例如，图中红线标识出的是车牌为“粤B****W”车辆的时空轨迹，在路段运行期间，它一共发生4次超车行为，大致的时空坐标分别为(19:26:29，260m)、(19:27:10，900m)、(19:27:20，1150m)和(19:27:55，1700m)。
选取时空轨迹图中任意一块区域，通过统计某时刻某空间范围内或某空间位置出在某一时间范围内的轨迹通过数，可以估计对应时刻的车流密度和流率。例如图中红色框线标出位置，在19:25:30时，距离上游卡口1250-1750米范围内共有12辆车在运行，可以估算出此时该范围的车流密度大约为32pcu/km。距离上游卡口1500米处，在19:25:30-19:26:00时段内，通过了11辆车，则平均车头时距为2.73s/pcu，对应的流率估计值为1318pcu/h。
本申请中所述的“轨迹重构”的含义指的是车辆时空轨迹的重构，而非实际的运行轨迹重构。上面结合附图对本发明的实施方式作了详细说明，但是本发明并不限于上述实施方式，在本领域普通技术人员所具备的知识范围内，还可以在不脱离本发明宗旨的前提下做出各种变化。

Claims (1)

1.一种基于数据融合的车辆时空轨迹重构方法，其特征在于，针对目标路段上各浮动车，执行以下步骤1-步骤4，针对目标路段上所有车辆，执行步骤5，获得目标路段上所有车辆的行驶轨迹，进一步构建考虑超车行为的车辆时空轨迹重构模型，完成目标路段的交通管控：
步骤1：采集目标路段各卡口的AVI数据，以及目标路段上各浮动车的GPS数据，分别构建AVI数据集、GPS数据集，以浮动车的车牌号为标识，分别从GPS数据集中提取各浮动车的车牌号所对应的轨迹点，并按照每辆浮动车的车牌号所对应的轨迹点的采样时间排序，获得每辆浮动车离散的轨迹点序列；
步骤1中，GPS数据集中包括每辆浮动车的轨迹点，且GPS数据集中每辆浮动车的轨迹点按时间戳排序，以浮动车的车牌号为标识，从GPS数据中提取每辆浮动车的车牌号所对应的轨迹点，并按照每辆浮动车的车牌号所对应的轨迹点的采样时间排序，获得每辆浮动车离散的轨迹点序列Pprobe＝(p1，p2，...，pi，...，pn)，其中p1，p2，...，pi，...，pn表示车牌号所对应的各轨迹点，n为轨迹点总数；
步骤2：基于地图匹配方法中的几何匹配法，利用AVI数据和GPS数据中的关联特征，确定各浮动车的轨迹点所对应的匹配道路，完成浮动车轨迹的修正；
根据球面距离计算公式计算各浮动车的轨迹点之间的距离，完成浮动车轨迹的时空坐标转换，获得离散的浮动车时空轨迹；
步骤2的具体步骤如下：
步骤2.1：分别取AVI数据集中上游卡口和下游卡口AVI数据的车牌和时间字段，判断浮动车j的车牌号是否在其中且连续通过上游卡口和下游卡口，若是，则认为浮动车j通过目标路段；否则认为浮动车j未通过目标路段，取下一辆浮动车的车牌号重复上述判断过程，直至所有浮动车判断完毕；
步骤2.2：通过卡口AVI数据约束确定匹配路段后，针对轨迹点pi，基于几何匹配法获得轨迹点pi在目标路段的道路中心线上的投影点ei，取下一个轨迹点重复上述匹配过程，直至所有轨迹点匹配完毕，则输出修正后的轨迹点序列(e1，e2，...，ei，...，en)，其中e1，e2，...，ei，...，en表示修正后的各轨迹点，n为轨迹点总数；
步骤2.3：根据球面距离计算公式，计算轨迹点pi、pj之间的距离dij如下式：
dij≈cos-1(sinβi·sinβj+cosβi·cosβj·cosΔα)·R
式中，Δα＝αi-αj，其中αi、αj分别为轨迹点pi、pj的经度，βi、βj分别为轨迹点pi、pj的纬度，R为地球半径；
构建以目标路段的上游卡口为起点的离散浮动车时空轨迹序列(t1，x1)，(t2，x2)，...，(tn，xn)，其中t1，t2，...，tn、x1，x2，...，xn分别为轨迹点p1，p2，...，pn的时间、空间位置；
步骤3：基于分段三次hermit插值法，对步骤2所获得的浮动车时空轨迹的缺失段进行插值，获得连续且平滑的浮动车时空轨迹；
步骤3的具体步骤如下：
步骤3.1：根据分段三次Hermit插值法，对于浮动车时空轨迹序列(t1，x1)，(t2，x2)，...，(ti，xi)，...，(tn，xn)，在区间[ti-1，ti]上的插值条件如下式：
其中，f(t)为插值函数，其表达式以分段三次多项式表示为如下形式：
步骤3.2：构建Hermit插值法的基函数如下式：
式中，均为三次多项式，称为三次Hermit插值多项式的基函数；
步骤3.3：对步骤3.2所构建的Hermit插值法的基函数两边求导：
结合步骤3.1所述插值条件，可得下式：
记Δi＝ti-ti-1，基函数的表达式如下式：
针对上式两边求导可得下式：
其中c和d表达式如下式：
求解获得的表达式如下式：
分别求解获得的表达式如下式：
多项式xi(t)的表达式如下式：
通过求解每个插值区间上的多项式表达式，获得分段三次Hermit多项式表达式X(t)；
步骤3.4：针对两轨迹点之间时空轨迹曲线的曲线震荡，引入约束条件如下式，进行插值结果修正：
xi＝max{x(ti)，x(ti-1)}
步骤3.5：针对浮动车j的时空轨迹序列(t1，x1)，(t2，x2)，...，(ti，xi)，...，(tn，xn)，取一个轨迹点(ti，xi)作为误差计算对象，用时空轨迹序列中剩余的n-1个轨迹点进行插值补齐，插值的时间步长取Δt＝0.1s，误差指标采用平均绝对误差MAE和均方根误差RMSE，浮动车j的平均绝对误差MAEj和均方根误差RMSEj如下式：
式中，为第i个轨迹点的插值结果；
目标路段上所有浮动车的平均绝对误差MAE和均方根误差RMSE如下式：
式中，m为目标路段上的浮动车总数；
步骤4：基于三角型基本图模型，并基于步骤3所获得的浮动车时空轨迹，提取交通流参数，交通流参数包括自由流速度、交通流激波波速、阻塞密度；
步骤4的具体步骤如下：
步骤4.1：基于Hermit插值法，获得目标路段上每辆浮动车的时空轨迹曲线，从浮动车i进入目标路段的起始时刻开始，以步长Δt为窗口向由平滑，对于第j个窗口，对应的自由流速度vj如下式：
根据自由流速度vj获得浮动车i在目标路段上速度的分布情况，基于帕累托法则，取80％分位的速度作为浮动车i自由流速度的估计值，记为vi-80th；重复上述过程获得每辆浮动车自由流速度的估计值，目标路段的自由流速度的估计值如下式：
式中，m为目标路段上的浮动车总数；
步骤4.2：根据目标路段上各浮动车的时空轨迹曲线中拥堵发生处拟合获得交通流激波波速w；
步骤S4.3：阻塞密度kj的表达式如下式：
式中，km为车流达到饱和流率对应最佳密度，qm为目标路段的饱和流率；
步骤5：在LWR模型和Newell模型中引入车辆顺序变量，用车辆顺序的变化表示车辆的超车行为，通过推导车辆顺序变量与空间、时间、累计车辆数的数学关系，将车辆轨迹重构问题转化为车辆顺序的求解问题，并以各交通流参数为输入、以所有车辆在目标路段上的重构轨迹为输出，构建考虑超车行为的车辆时空轨迹重构模型，并应用车辆时空轨迹重构模型，完成目标路段的交通管控；
步骤5的具体步骤如下：
步骤5.1：以时间t和空间位置x所对应的累计车辆数N(x，t)为状态变量，将[0，t]时段目标路段上游卡口通过的累计车辆数记为U(t)＝N(0，t)-N0，其中N0为t＝0时目标路段初始的车辆数，目标路段下游卡口通过的累计车辆数记为D(t)＝N(l，t)，l为目标路段长度；
步骤5.2：针对车辆i，用ei和di分别表示车辆i进入目标路段和离开目标路段的时间，计算车辆i通过目标路段时车辆顺序的变化量Ci如下式：
Ci＝δi(di)-δi(ei)
式中，δi(di)为车辆i离开目标路段时的车辆顺序，δi(ei)为车辆i进入目标路段时的车辆顺序；
当Ci＞0时，表示车辆i通过目标路段时被Ci辆车超越；当Ci＝0时，表示车辆i通过目标路段时车辆顺序没有变化；当Ci＜0时，表示车辆i通过目标路段时超越了Ci辆车；
步骤5.2：计算t＝0时刻目标路段内部初始车辆数N0；以每辆车辆通过上下游卡口时，对应累计车辆数差的平均值n0来估计t＝0时刻目标路段内部初始车辆数N0，N0表达式如下式：
式中，h为通过上下游卡口且均被记录数据的车辆数；
步骤5.3：将车辆i进入和离开目标路段的时间区间[ei，di]以步长Δt＝0.1s划分；
步骤5.4：当t＝μΔt∈[ei，di]时，计算车辆i的车辆顺序δi(μΔt)如下式：
式中，μ为整数；
当t＝μΔt时，在非拥挤情况下，车辆i的空间位置如下式：
其中，U-1(t)为U(t)的逆函数；
当t＝μΔt时，在拥挤情况下，车辆i的空间位置如下式：
其中，的求解方法为二分法、牛顿法、梯度法中的一种；
步骤5.5：比较和/>大小，取小者作为t＝μΔt时车辆i的预测空间位置，并增加约束条件，具体如下式：
步骤5.6：重复步骤5.4，计算t＝μΔt+Δt时刻的δi(t)、和/>直到t+Δt≥di，获得车辆i以Δt＝0.1s为时间间隔的空间位置序列：xi(ei)，xi(ei+Δt)，xi(ei+2Δt)，...，xi(di)，作为车辆i在目标路段上的重构轨迹；
步骤5.7：基于浮动车GPS数据集和卡口AVI数据集的关联性，针对重构得到的浮动车重构轨迹进行车辆轨迹重构精度的评价；
步骤5.7中评价车辆轨迹重构精度的具体方法如下：
基于平均绝对误差MAE和均方根误差RMSE，进行车辆轨迹重构精度的评价，即比较每个时刻t对应的浮动车i的实际空间位置和重构轨迹的预测空间位置之差：
式中，xi(t)为t时刻浮动车i的实际空间位置，/>为t时刻浮动车i重构轨迹的预测空间位置；
以相对误差来评价重构轨迹误差，即浮动车i的实际时空轨迹和重构轨迹围成的面积与实际时空轨迹围成的面积之比，浮动车i重构轨迹的相对误差errori如下式：
式中，Si1和Si2分别表示浮动车i实际时空轨迹围成的面积和实际时空轨迹与重构轨迹围成的面积。
CN202210868085.XA 2022-07-22 2022-07-22 一种基于数据融合的车辆时空轨迹重构方法 Active CN115311854B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202210868085.XA CN115311854B (zh) 2022-07-22 2022-07-22 一种基于数据融合的车辆时空轨迹重构方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202210868085.XA CN115311854B (zh) 2022-07-22 2022-07-22 一种基于数据融合的车辆时空轨迹重构方法

Publications (2)

Publication Number Publication Date
CN115311854A CN115311854A (zh) 2022-11-08
CN115311854B true CN115311854B (zh) 2023-08-25

Family

ID=83856728

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202210868085.XA Active CN115311854B (zh) 2022-07-22 2022-07-22 一种基于数据融合的车辆时空轨迹重构方法

Country Status (1)

Country Link
CN (1) CN115311854B (zh)

Families Citing this family (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN115879294B (zh) * 2022-11-30 2024-04-26 西部科学城智能网联汽车创新中心(重庆)有限公司 一种基于多车环境感知的全样车流轨迹生成方法及系统
CN116415199B (zh) * 2023-04-13 2023-10-20 广东铭太信息科技有限公司 基于审计中间表的业务数据离群分析方法
CN116543559A (zh) * 2023-05-22 2023-08-04 东南大学 一种考虑超车行为的车辆轨迹重构方法
CN116597397B (zh) * 2023-07-17 2023-10-24 腾讯科技（深圳）有限公司 预测车辆轨迹的模型训练方法和装置及存储介质
CN117334051B (zh) * 2023-10-26 2024-05-10 江苏中路交通发展有限公司 一种高速路车辆轨迹重构方法及系统

Citations (9)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN105788252A (zh) * 2016-03-22 2016-07-20 连云港杰瑞电子有限公司 基于定点检测器和信号配时数据融合的城市干道车辆轨迹重构方法
CN106652458A (zh) * 2017-02-20 2017-05-10 东南大学 基于虚拟车辆轨迹重构的在线城市道路路径行程时间估计方法
CN108447256A (zh) * 2018-03-22 2018-08-24 连云港杰瑞电子有限公司 基于电警和定点检测器数据融合的干道车辆轨迹重构方法
CN108492562A (zh) * 2018-04-12 2018-09-04 连云港杰瑞电子有限公司 基于定点检测与电警数据融合的交叉口车辆轨迹重构方法
CN109064741A (zh) * 2018-08-01 2018-12-21 北京航空航天大学 基于多源数据融合的干道车辆运行轨迹重构的方法
CN111243277A (zh) * 2020-03-09 2020-06-05 山东大学 基于车牌识别数据的通勤车辆时空轨迹重构方法及系统
WO2021073524A1 (zh) * 2019-10-15 2021-04-22 同济大学 一种拥堵交通流溯源分析方法
CN113129604A (zh) * 2021-03-19 2021-07-16 同济大学 一种基于网联车辆轨迹数据的信控交叉口运行评估方法
CN113724489A (zh) * 2021-07-22 2021-11-30 东南大学 基于多源数据的交通拥堵溯源方法

Patent Citations (9)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN105788252A (zh) * 2016-03-22 2016-07-20 连云港杰瑞电子有限公司 基于定点检测器和信号配时数据融合的城市干道车辆轨迹重构方法
CN106652458A (zh) * 2017-02-20 2017-05-10 东南大学 基于虚拟车辆轨迹重构的在线城市道路路径行程时间估计方法
CN108447256A (zh) * 2018-03-22 2018-08-24 连云港杰瑞电子有限公司 基于电警和定点检测器数据融合的干道车辆轨迹重构方法
CN108492562A (zh) * 2018-04-12 2018-09-04 连云港杰瑞电子有限公司 基于定点检测与电警数据融合的交叉口车辆轨迹重构方法
CN109064741A (zh) * 2018-08-01 2018-12-21 北京航空航天大学 基于多源数据融合的干道车辆运行轨迹重构的方法
WO2021073524A1 (zh) * 2019-10-15 2021-04-22 同济大学 一种拥堵交通流溯源分析方法
CN111243277A (zh) * 2020-03-09 2020-06-05 山东大学 基于车牌识别数据的通勤车辆时空轨迹重构方法及系统
CN113129604A (zh) * 2021-03-19 2021-07-16 同济大学 一种基于网联车辆轨迹数据的信控交叉口运行评估方法
CN113724489A (zh) * 2021-07-22 2021-11-30 东南大学 基于多源数据的交通拥堵溯源方法

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
基于改进三次Hermite插值的车辆时空轨迹重构研究;陈志军;吴超仲;吕能超;马杰;;交通信息与安全(06);全文 *

Also Published As

Publication number Publication date
CN115311854A (zh) 2022-11-08

Similar Documents

Publication Publication Date Title
CN115311854B (zh) 一种基于数据融合的车辆时空轨迹重构方法
CN109670277B (zh) 一种基于多模态数据融合与多模型集成的旅行时间预测方法
Chen et al. Assessing right-turning vehicle-pedestrian conflicts at intersections using an integrated microscopic simulation model
Geroliminis et al. Identification and analysis of queue spillovers in city street networks
Ma et al. Estimating maximum queue length for traffic lane groups using travel times from video-imaging data
Anand et al. Traffic density estimation under heterogeneous traffic conditions using data fusion
Lu et al. Effects of winter weather on traffic operations and optimization of signalized intersections
CN103295414A (zh) 一种基于海量历史gps轨迹数据的公交车到站时间预测方法
CN110299011A (zh) 一种基于收费数据的高速公路任意断面的交通流预测方法
CN111915887A (zh) 一种基于多源异构交通数据的集成与处理系统及其方法
Aljamal et al. Real-time estimation of vehicle counts on signalized intersection approaches using probe vehicle data
CN111881557B (zh) 基于道路平均速度的车流仿真方法
CN102855755B (zh) 一种基于运行车速预测的城市干道车队离散模型的建模方法
Wang et al. Forecasting traffic volume at a designated cross-section location on a freeway from large-regional toll collection data
Oskarbski et al. Estimating the average speed of public transport vehicles based on traffic control system data
Bassani et al. Experimental analysis of operational data for roundabouts through advanced image processing
Anusha et al. Dynamical systems approach for queue and delay estimation at signalized intersections under mixed traffic conditions
Hawas et al. Optimized multistage fuzzy-based model for incident detection and management on urban streets
CN107886192B (zh) 基于固定式与移动式车辆检测数据的数据与信息融合方法
CN113538902B (zh) 基于交通状态的交叉口车辆轨迹数据修复方法
Derevitskiy et al. Traffic estimation on full graph of transport network using GPS data of bus movements
Hegyi et al. Freeway Traffic Management and Control.
KR20230063080A (ko) 교통량 증대를 위한 신호교차로 노면조건 및 교통패턴 분석 결과에 따른 신호등 운영 시스템
Wang et al. Trajectory-based vehicle emission evaluation for signalized intersection using roadside LiDAR data
Guilbert et al. Re-identification by Inductive Loop Detector: Experimentation on target origin—Destination matrix

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant