CN103903019A - 一种多车道车辆轨迹时空图的自动生成方法 - Google Patents

Abstract

本发明涉及一种多车道车辆轨迹时空图的自动生成方法。首先通过摄像机离线标定获取摄像机内外参数、根据实际监控需求设置感兴趣区域，并在车流上游起始位置设置三条虚拟检测带；然后基于这三级虚拟检测带进行车辆特征点的提取并按车辆分组，进行特征点跟踪与分组修正并完成车辆轨迹的高精度提取，进而为每条车道自动生成驶过车辆的纵向时空轨迹图并生成一幅含所有车道的横向时空轨迹图，最后自动判断是否处于车道线修正阶段，如果是则基于轨迹聚类法自动修正车道线位置。与现有技术相比，本发明所述方法提高了目标亮度变化的鲁棒性，提高了车辆检测与轨迹跟踪的精度，能够很好地实现车道线的在线自动修正，避免车辆轨迹点所属车道号的错误判断。

Description

一种多车道车辆轨迹时空图的自动生成方法

技术领域

本发明属交通运输监控技术领域，利用计算机视频处理技术及模式识别技术自动提取道路运动车辆轨迹，在此基础上自动判别车辆所属车道并生成各车道车辆轨迹纵向时空图和路段横向时空图。

背景技术

车辆轨迹时空图是车辆行驶过程中的“时间-位置”坐标表示，其中通常取时间为横坐标轴，车辆在车道中的空间位置（纵向或横向）为纵坐标轴。车辆轨迹时空图在交通技术领域中具有重要应用价值。不仅可为交通流量、交通密度、平均速度以及车辆旅行时间等多种交通参数的获取提供统一的支持平台，也能够为交通波检测及交通流特性分析提供高效支持。

车辆轨迹时空图的生成以车辆的定位与跟踪技术为基础。现有的车辆轨迹定位与跟踪技术主要有GPS技术以及视频图像处理技术等。但基于GPS的车辆轨迹定位与跟踪技术需要车载装置的支持，并且其定位精度往往不能很好地精确到车道级别，因此其应用受到很大限制；基于视频的车辆检测与跟踪技术能够很好地解决这些问题，并且具有获取信息量丰富、设备成本低廉、安装维护方便等优点，因而在目前交通监控领域中应用最为广泛。

传统的视频车辆检测跟踪方法主要有基于运动前景提取与跟踪的方法和基于特征点提取与跟踪的方法等。基于运动前景提取与跟踪的方法能够在车流畅通时很好地获取车辆运行轨迹，然而当出现交通拥堵、车流缓慢以及相互遮挡情况时，该方法的效果会急剧下降甚至失效。基于特征点提取与跟踪的方法能够较好地解决这一问题，但其计算量大，在没有特殊硬件支持的条件下算法实时性较差，同时来自于相邻车辆运动特征高度相似的干扰以及车辆活动阴影的干扰也在很大程度上影响了特征点所属车辆的准确判断并降低了车辆跟踪的精度。

发明内容

本发明的目的在于克服现有技术在精度和速度上的不足，提供一种适用于交通监控与分析、能够快速自动生成高精度多车道车辆轨迹时空图的方法，不但可为交通领域研究人员提供一种切实有效的基础研究工具，同时也为交通监控性能和交通分析能力的进一步提高提供有力支持。

为实现上述目的，本发明采取以下技术方案：以来自道路交通监控固定摄像机的视频序列为输入，所述道路交通监控固定摄像机的水平朝向与车流方向允许存在一定角度θ，0≤θ＜90。首先通过摄像机离线标定获取摄像机内外参数、根据实际监控需求设置感兴趣区域（亦可为整幅图像），并在车流上游起始位置设置三条虚拟检测带；工作过程中基于这三级虚拟检测带进行车辆特征点的提取并结合Adaboost智能分类器实现特征点的按车辆分组，然后进行特征点跟踪与分组修正并完成车辆轨迹的高精度提取，进而为每条车道自动生成驶过车辆的纵向（沿车道方向）时空轨迹图并生成一幅含所有车道的横向（垂直于车道方向）时空轨迹图，最后自动判断是否处于车道线修正阶段，如果是则基于轨迹聚类法自动修正车道线位置。

一种多车道交通时空轨迹图的自动生成方法，其特征在于包括以下步骤：

步骤1，设置初始参数。

步骤1.1，根据实际监控需求设置感兴趣区域并手工指定各车道线的初始位置。

步骤1.2，设置垂直于车道的虚拟检测带。

在感兴趣区域内车辆上游位置设置三条垂直于车道的虚拟检测带，各条检测带在垂直车道线方向均贯穿路面的感兴趣区域，各条检测带的高度，即沿车道线方向的尺寸为h_VDL，h_VDL＜＜h_ROI，h_ROI为感兴趣区域的高度，即沿车道线方向的尺寸。

步骤1.3，对摄像机进行标定以获取内外参数。

步骤2，在帧图像中基于三级虚拟检测带进行车辆特征点的自动提取和按车辆分组。

读入道路交通监控摄像机的视频序列图像，并在帧图像中基于三级虚拟检测带进行车辆特征点的自动提取和按车辆分组。首先对三级虚拟检测带进行局部背景建模，然后以三级递增的形式提取驶入车辆的前景特征点，并在车辆驶离第三级检测带时利用Adaboost智能分类器进行车辆识别并根据识别结果对特征点进行分组。Adaboost智能分类器的施用位置由三级虚拟检测带生产的三级PVIs（Panoramic View Images，全景视图）确定。具体方法如下：

步骤2.1，基于局部变更新率背景建模提取局部前景并消除活动阴影。

（1）对三级虚拟检测带进行变更新率局部动态背景建模。

本发明设计了更新率自适应的动态混合高斯背景建模方法，对各虚拟检测带内的每个像素，根据其在最近历史帧中的采样值{X¹,X²,…,X^t}采用K（通常K∈{3,4,5}）个高斯分布的叠加来进行建模。使用历史帧进行背景建模后，对于新像素值X^t+1按ω_i,t/σ_i,t降序与该点各高斯分布进行匹配。对匹配成功的高斯分布更新如下：

$\left\{\begin{array}{c}{\mathrm{\ω}}_{i,t+1}=(1-\mathrm{\α}){\mathrm{\ω}}_{i,t}+\mathrm{\α}\\ {\mathrm{\μ}}_{i,t+1}=(1-\mathrm{\ρ}){\mathrm{\μ}}_{i,t}+\mathrm{\ρ}{X}^{t+1}\\ {\mathrm{\σ}}_{i,t+1}^2=(1-\mathrm{\ρ}){\mathrm{\σ}}_{i,t}^2+\mathrm{\ρ}{(X^{t+1}-{\mathrm{\μ}}_{i,t})}^2\\ \mathrm{\ρ}=\mathrm{\α\η}\left(X^{t+1}\right|{\mathrm{\μ}}_{i,t},{\mathrm{\σ}}_{i,t})\end{array}\right.$

对于其它未匹配的高斯分布，均值和方差保持不变，权值更新为ω_i,t+1＝(1-α)ω_i,t，α为自适应参数更新速率，0≤α≤1。实验表明，交通畅通时的更新速率取0.001～0.005能获得好的更新效果。然而考虑到出现交通拥堵情况时，可能会有排队车辆恰好停在（或非常缓慢的移动在）虚拟检测带上，此时如果仍按畅通时的变化率进行背景更新，则会将车辆误判为背景。为避免这种情况的发生，本发明提出了Log-Sigmoid更新率动态自适应方法，其计算公式如下：

$\mathrm{\α}=\frac{{\mathrm{\α}}_0}{1+e^{-\mathrm{\β}(l-\mathrm{\ξ})}}$

式中，α₀为交通畅通时的更新速率；参数β决定了Log-Sigmoid型函数的陡度；l为由步骤3跟踪程序反馈而来的车体区域特征点组位置到虚拟检测带3的最小像素距离；考虑到车辆之间会存在一定的间隙距离，所以在l大于0但小于等于间隙距离时已很可能有车压线，此时应停止更新背景。此间隙一般不会超过1辆车长，本发明以场景中车辆的平均长度统计值ξ来作为间距估计值，ξ基于前2级虚拟检测带中心线生成的当前时刻之前的PVIs计算得到，公式如下：

$\mathrm{\ξ}=\frac{1}{m}\underset{i=1}{\overset{m}{\mathrm{\Σ}}}\frac{\mathrm{\Δd}\·l_i}{F_s\·(t_i^2-t_i^1)}$

式中，m为两级PVIs中成功匹配出来的车辆Blob（前景块）个数，l_i为第i个车辆Blob的长度，

分别为第i辆车在两级PVIs中的Blob中心t坐标，Δd为两级虚拟带的距离，F_s为帧频。

对于与任何分布匹配均不成功的新像素值X^t+1，以X^t+1为均值、以一个较大值为方差构造一个新的高斯分布并代替ω_i,t/σ_i,t最小的分布来完成更新。对各高斯分布按ω_i,t/σ_i,t进行降序排列并选择最前面的B个分布作为当前点的背景模型，其中：

$B=\arg \underset{b}{\min }(\underset{k=1}{\overset{b}{\mathrm{\Σ}}}{\mathrm{\ω}}_{i,t}>T)$

式中，T为计入背景的最小阈值。

（2）提取局部前景并消除活动阴影。

基于背景差分法分割各虚拟检测带上的活动前景，并采用前景像素值与对应位置处的背景像素值的局部纹理特征相似性原理有效去除前景中包含的车辆活动阴影，提高稀疏特征点的提取精度。

步骤2.2，提取特征点并按车辆进行分组。

（1）提取特征点。

利用文献Tomasi方法，在第一级虚拟检测带的局部前景图像序列中对驶过的车辆进行易跟踪特征点的初步提取。鉴于检测带沿车道方向高度较窄，一级特征点选取不够充分，因此，在第二、三级检测带进行逐级递增式补选特征点，很好地保证车辆前景特征点的提取效果。

（2）对上述由三级检测带得到的特征点按车辆进行分组。

鉴于实际场景中不可避免会出现车辆的相互遮挡问题，本发明设计了基于三级前景PVIs和Adaboost智能分类器相结合的特征点分组方法较好地解决这一问题。具体方法如下：

首先，按下式方法分别将视频序列在三级虚拟检测带中心线上的前景像素信息以时间顺序依次排列生成三级前景PVIs：

$\begin{array}{c}{P}_t^k=X_k^t,k\∈\left\{\mathrm{1,2,3}\right\}\\ P_t^k={[P_{t,0}^k,P_{t,1}^k,\·\·\·,P_{t,N}^k]}^T\\ X_k^t=[X_{0,y_k}^t,X_{1,y_k}^t,\·\·\·,X_{N,y_k}^t]\end{array}$

其中，

为第k级PVI的t时刻列像素向量，

为第k级虚拟检测带中心线位置的t时刻像素向量，N为虚拟检测带中心线的像素级长度，y_k为第k级虚拟检测带中心线y坐标。

然后，在目标刚刚驶离第三级检测带时，利用Adaboost智能分类器在原始图像中进行车辆识别并根据识别结果对特征点进行分组。Adaboost智能分类器在原始图像中的施用区域由该目标在三级PVIs中的Blob（下面记为B^k,k∈{1,2,3}）及相互匹配情况来确定。施用区域左下顶点图像坐标

及区域宽度w^Adaboost计算如下：

$\begin{array}{c}{x}_{\mathrm{Lb}}^{\mathrm{Adaboost}}=x_{\mathrm{Lb}}^{B^3}\\ y_{\mathrm{Lb}}^{\mathrm{Adaboost}}=y^3\\ w^{\mathrm{Adaboost}}=w^{B^3}\end{array}$

其中，

为B³外接矩形的左下角顶点图像坐标，y³为第三级虚拟检测带中心线y坐标，

为B³外接矩形的宽度。

施用区域高度h^Adaboost的计算，分为如下三种情况：

第1种情况：目标在三级PVIs中形成的前景块B³、B²和B¹均相互匹配成功，则：

$h^{\mathrm{Adaboost}}=\frac{l^{B^3}\·|y^3-y^2|}{(t^3-t^2)\·F_s}$

其中，

为B³外接矩形的长度，y²为第二级虚拟检测带中心线y坐标，t³和t²分别是目标到达第三、二级虚拟检测带中心线的时刻值，F_s为帧频；

第2种情况：B³只与B²、B¹中的一个匹配成功，记为B'，则：

$h^{\mathrm{Adaboost}}=\frac{l^{B^3}\·|y^3-y^{\′}|}{(t^3-t^{\′})\·F_s}$

其中，y'为B'所在虚拟检测带中心线y坐标，t'为目标到达B'所在虚拟检测带中心线的时刻值；

第3种情况：B²与B¹匹配成功，但均不能与B³匹配，则：

$h^{\mathrm{Adaboost}}=\frac{|y^2-y^1|}{(t^2-t^1)\·F_s}\·(\max (l^{B^3},l^{B^2},l^{B^1})+\mathrm{\Δh})$

其中，Δh为分离距离修正量，构造为能够有效防止原本粘连的目标在第三级虚拟检测带处发生分离而引起的目标漏检。

（3）同车特征点组的表示。

步骤（2）能将绝大多数的粘连车辆分离并确定各车辆的车体区域。对落在同一车体区域的特征点采用边加权动态图G＝{V,E,W}表述。其中，V为同车特征点集；E为V内各特征点间的优化边集，采用Delaunay三角剖分法构造得到；W为边权集，各边权值取为该边所连特征点对的运动特征相关度，初始化为车体区域长度l_vehicle，并在后继跟踪中动态修改。

步骤3，特征点跟踪与分组修正。

对分组之后的特征点采用在后继帧的图像金字塔中计算Kanade-Lucas-Tomasi稀疏光流并结合Uniform Local Binary Pattern局部二元纹理特征进行特征点匹配的方法进行跟踪。同时，在跟踪过程中基于运动特征相关度对少数的误分组特征点进行必要的修正和干扰点剔除，以充分保证跟踪精度。具体方法如下：

步骤3.1，基于金字塔Kanade-Lucas-Tomasi稀疏光流与Uniform Local Binary Pattern纹理特征匹配的特征点跟踪。

在后继帧的图像金字塔中计算Kanade-Lucas-Tomasi稀疏光流进行特征点跟踪。但Kanade-Lucas-Tomasi稀疏光流法要求被跟踪目标的亮度恒定，且被跟踪目标的运动相对于帧频来说为缓慢运动，然而在实际交通场景中，特征点很容易受到车速过快或者车灯光照变化等因素的影响而导致跟踪漂移。大量的实验表明在复杂的车辆跟踪场景中，局部二元图（Uniform Local Binary Pattern）纹理特征对目标光照变化具有鲁棒性强及目标表达更准确的优势。因此，本发明将其与Kanade-Lucas-Tomasi稀疏光流法进行了有机结合，有效解决了跟踪漂移问题。方法如下：在稀疏光流跟踪过程中，如果发现第i帧的某个特征点在后继帧中的位置发生了变动，则分别以变动前后的这两个点为中心，在各自的帧图像中分别计算其n×n邻域窗口的Uniform Local Binary Pattern局部二元图纹理特征值并进行匹配。若匹配成功则表明两点处的局部纹理特征一致，位置变动后的特征点有效；否则，表明该特征点在跟踪过程中发生了漂移，此时将其剔除。

步骤3.2，基于运动特征相关度对误分组特征点进行必要的修正和干扰点剔除。

在跟踪过程中，采用状态空间转移模型来表述特征点对的运动关系并结合Kalman滤波对步骤2结果中存在的误分组特征点进行必要的分组修正，同时对一些可能出现的干扰点进行剔除。正常行驶的车辆一般可视作刚体，同一车辆上的特征点在运动方向、速度、加速度等运动特征上具有较高相似度，但不同车辆上的特征点的运动特性一般差异相对较大。因此，可基于运动特征的相关程度动态修改边加权动态图中各特征点间的边权值，从而实现分组修正和干扰剔除。分析和实验表明，特征点的上述运动特性差异可通过它们在跟踪过程中的彼此距离及其距离变化率等参数来描述。对于本发明所述交通场景，车辆向远离摄像机的一方行驶，在二维视频帧图中同一车辆上的两个特征点即使车辆以恒定的速度移动，两点间距也会因摄像机角度和透视效果的影响而随时间变化。本发明采用状态空间转移模型来表述特征点对的运动关系并结合Kalman滤波来完成特征点的分组修正和筛选剔除。仅以x轴方向为例（y轴方向同理），对于具有边连接关系的特征点对（i和j），以

表示其在t时刻的运动关系状态向量

，其中和

分别为两点间的x轴方向距离、距离变化速率和变化加速度，则此特征点对系统可表述为：

其中，

为t-1时刻特征点对的状态向量，

为t时刻的状态测量向量，

和η分别为系统噪声向量和测量噪声向量，且～N(0,Q),P(η)～N(0,R)，Q和R分别为系统噪声和测量噪声的协方差矩阵。A为转移矩阵，H为观测矩阵。因此，在t时刻状态预测值及状态修正值

可由点对在上一时刻t-1的历史信息递推得到，具体如下：

$\begin{array}{c}{\hat{X}}_{i,j}^{\left(t\right)}=A{X}_{i,j}^{(t-1)}\\ X_{i,j}^{\left(t\right)}={\hat{X}}_{i,j}^{\left(t\right)}+K(Y_{i,j}^{\left(t\right)}-H{\hat{X}}_{i,j}^{\left(t\right)})\end{array}$

式中，K为Kalman修正矩阵；

为误差向量，其第一分量表征了t时刻的点对距离测量值与预测值之间的差异，如果该差异值为正数，则表明点对可能出现了分离的趋势。为避免误判，本发明不在此时直接断开该点对的边连接，而是动态修改该边的权值，使其减1。然后判断边权，如果其值已由初值减为0，则表明该点对不属于同一车辆，此时断开该点对的连接，实现特征点的分组修正。进而判断这两个特征点，如果某点已无任何边与其相连，则将其作为干扰点剔除。

步骤4，车辆轨迹提取与多个车道的时空轨迹图自动生成。

步骤4.1，对每辆车在其驶出感兴趣区域时由其特征点轨迹提取该车辆轨迹。

对于此时该车的边加权动态图G＝{V,E,W}，V中保存的是该车的全部稳定特征点，E中保存着这些特征点之间的最优边。事实表明，如果一个特征点越靠近车辆中心位置，那么它在E中边的数量就越多；反之，越靠近车辆边缘则边数越少。据此，本发明提出了一种加权平均法来计算V中各特征点的分布中心

$\begin{array}{c}\stackrel{\‾}{x}=\frac{1}{2n}\underset{1\≤i\≤n}{\mathrm{\Σ}}{m}_i{x}_i\\ \stackrel{\‾}{y}=\frac{1}{2n}\underset{1\≤i\≤n}{\mathrm{\Σ}}{m}_i{y}_i\end{array}$

其中，n为V的尺度，x_i和y_i为第i个特征点的图像坐标，m_i为第i个特征点的权重。

进而，计算V中各个特征点i与分布中心

的欧氏距离按照

升序重排各特征点并选出前K个点构成分布中心的K-最近邻点集C。最终按下式选取第j个特征点的轨迹作为当前车辆轨迹：

$j=\underset{p\∈C}{\arg \max }{l}_p$

其中，l_p为K-最近邻点集C中特征点p的轨迹长度。

步骤4.2，根据各辆车的轨迹数据自动生成多个车道的时空轨迹图。

根据实际车道数为每条车道自动生成驶过车辆的纵向时空轨迹图，并生成一幅含所有车道的横向时空轨迹图。

（1）生成纵向时空轨迹图。

以第i辆车为例，其轨迹为如下的点序列：

T_i={<x₁,y₁,f₁>,…,<x_j,y_j,f_j>,…,<x_M,y_M,f_M>}

其中，各点均为三元组结构，x_j、y_j和f_j分别为第j个轨迹点的x轴、y轴图像坐标以及轨迹点所在的帧号，M为该轨迹中的点数。

对于各轨迹点<x_j,y_j,f_j>，首先根据坐标值判断其所在车道号k，对于本发明所述交通场景，计算公式如下：

$k=q,(x_q^l<x_j<x_q^r{|}_{y_q^l=y_j=y_q^r})$

其中，

为第q条车道左边界线上的点坐标，为第q条车道右边界线上的点坐标。

然后将该轨迹点绘入第k幅时空图的(t_j,d_j)坐标处，其中：

t_j＝f_j/F_s

d_j＝|Y_j-Y_O|

其中，F_s为摄像机帧频，Y_O为起始参考点的世界Y坐标，在此起始参考点取为第一级虚拟检测带中心线左端点；Y_j为当前轨迹点的世界Y坐标。各点的世界坐标值可由其图像坐标与摄像机标定参数计算得来。

（2）生成横向时空轨迹图。

对于第i辆车的轨迹点序列T_i={<x₁,y₁,f₁>,…,<x_j,y_j,f_j>,…,<x_M,y_M,f_M>}中的各轨迹点<x_j,y_j,f_j>,(1≤j≤M)，其在横向时空轨迹图中的坐标(t_j,d_j)满足：

t_j＝f_j/F_s

d_j＝|X_j-X_O|

其中，X_O为起始参考点的世界X坐标，X_j为当前轨迹点的世界X坐标，其余参数同上。

步骤5，如果处于车道修正阶段，则每新增γ条轨迹进行一次基于轨迹聚类法的车道位置自动修正，否则直接跳转至步骤2。

鉴于受摄像机拍摄角度及世界三维车辆到二维图像的投影变换影响，在步骤1的离线初始化阶段人工硬性设定的车道线往往与图像中车流实际情况存在一定偏差，如果不加以修正则会造成对某些车辆轨迹点所属车道号的判别错误。修正方法具体如下：考虑到实际交通场景中骑压车道线行驶的车辆轨迹对车道线修正效果的负面影响，本发明将粗糙集理论加以应用，提出了粗糙k-均值聚类法并结合轨迹的Hausdorff距离实现了车道位置的在线自动修正。包含如下步骤：

步骤5.1，在设置初始参数时确定聚类参数k及初始聚类中心。

（1）确定聚类参数k。

以车道中心线的位置确定为目标，聚类数目k即为车道数。对于确定的监控场景，车道数k为已知常数，因此可人工设定。

（2）确定初始聚类中心。

对于各聚类的初始聚类中心c_i，即每条车道的初始中心线，可根据步骤1中人工设定的感兴趣区域和各车道线初始位置数据计算得到，具体为：

$x_{c_i}=\frac{x_i^l+x_i^r}{2}{|}_{y_i^l=y_i^t=y_{c_i}\&Element;[y_{i,t},y_{i,b}]}$

其中，i为车道号且i∈{1,2,…,k}，c_i为第i条车道的中心线点集，

为c_i中的各点坐标；

为第i条车道左边界线上的点坐标，

为第i条车道右边界线上的点坐标，对于最两侧的车道其外边界即为感兴趣区域的外边界；y_i,t、y_i,d分别为当前车道的上、下顶点纵坐标。

步骤5.2，基于粗糙k-均值聚类实现车道分割线位置修正。

将每个车道聚类用粗糙集理论中的上、下近似集来描述，下近似集是上近似集的子集，用于保存肯定属于该车道类的轨迹样本，即在本车道内行驶的常规车辆，上下近似集之差称为边界区域，用于存放无法确定其明确归属车道的轨迹样本，即骑压车道线行驶的车辆轨迹。在每次更新聚类中心时，下近似集的轨迹被赋予较大的权重，边界区域的对象被赋予较小的权重。具体步骤如下：

（1）确定本次新增轨迹集中每条轨迹T_χ的上、下近似归属：

计算轨迹T_χ和各个聚类中心的Hausdorff距离h(T_χ,c_i),i＝1,2,...,k，公式如下：

$h(T_{\mathrm{\&chi;}},c_i)=\frac{1}{N_{\mathrm{\&chi;}}}\underset{a\&Element;T_{\mathrm{\&chi;}}}{\mathrm{\&Sigma;}}\underset{b\&Element;c_i}{\min }\left|\right|(x^a-x^b,y^a-y^b)\left|\right|$

其中，N_χ为轨迹T_χ的长度；

并找出与T_χ距离最短的类m，即

则：

其中，λ为判定阈值、

表示第m车道类的上近似集、

表示第m车道类的下近似集；

（2）更新聚类中心，公式如下：

$c_i=\left\{\begin{array}{cc}{w}_l\underset{T_{\mathrm{\&chi;}}\&Element;\underset{\&OverBar;}{C_i}}{\mathrm{\&Sigma;}}{T}_{\mathrm{\&chi;}}/\left|\underset{\&OverBar;}{C_i}\right|+(1-w_l)\underset{T_{\mathrm{\&chi;}}\&Element;(\stackrel{\&OverBar;}{C_i}-\underset{\&OverBar;}{C_i})}{\mathrm{\&Sigma;}}{T}_{\mathrm{\&chi;}}/|\stackrel{\&OverBar;}{C_i}-\underset{\&OverBar;}{C_i}|,& \stackrel{\&OverBar;}{C_i}\&NotEqual;\underset{\&OverBar;}{C_i}\\ \underset{T_{\mathrm{\&chi;}}\&Element;\underset{\&OverBar;}{C_i}}{\mathrm{\&Sigma;}}{T}_{\mathrm{\&chi;}}/\left|\underset{\&OverBar;}{C_i}\right|,& \stackrel{\&OverBar;}{C_i}=\underset{\&OverBar;}{C_i}\end{array}\right.$

其中，i∈{1,2,…,k}，权值w_l∈[0.5,1]、|·|表示基数。

（3）根据聚类后的车道中心线修正车道分割线位置，公式如下：

$x_{{\mathrm{Line}}_i}=\frac{x_{c_i}+x_{c_{i+1}}}{2}{|}_{y_{c_i}=y_{c_{i+1}}=y_{{\mathrm{Line}}_i}\&Element;[y_{{\mathrm{Line}}_i}^t,y_{{\mathrm{Line}}_i}^b]}$

其中，i为车道分割线标号且

及

为第i条车道线两侧的车道中心线上的点坐标；

和

分别为当前车道分割线的上、下顶点纵坐标。

（4）判断初始化阶段是否结束，判定条件为相邻两次修正的各聚类中心均无明显变化。

与现有技术相比，本发明具有以下明显优势：

（1）本发明采用三级虚拟检测带和局部背景变更新率动态建模的方法进行车辆特征点的提取，降低了特征点检测的计算量和活动阴影对精度的影响；

（2）本发明采用Adaboost智能分类器和运动一致性准则实现特征点的按车辆分组，提高了车辆检测的精度；

（3）在帧序列的图像金字塔中计算Kanade-Lucas-Tomasi稀疏光流并融合Uniform LocalBinary Pattern纹理特征进行特征点跟踪的方法，有效增强了对目标亮度变化的鲁棒性，充分提高了轨迹跟踪的精度；

（4）本发明采用结合轨迹的Hausdorff距离的粗糙k-均值聚类算法，能够很好地实现车道线的在线自动修正，避免车辆轨迹点所属车道号的错误判断，提高多车道时空轨迹图的生产精度。

附图说明

图1为交通监控摄像机安装实施例示意图；

图2为本发明所述方法的总体流程图；

图3为实施例中感兴趣区域及虚拟检测带示意图；

图4为车辆特征点的提取与分组方法的流程图；

图5为车辆特征点跟踪与分组修正方法的流程图；

图6为车辆轨迹提取与多个车道的时空轨迹图自动生成方法的流程图；

图7为基于粗糙k-均值聚类的车道线在线自动修正方法的流程图；

图8为实施例得到的各轨迹时空图的结果。

具体实施方式

下面结合附图和实例对本发明具体实施方式阐述如下。

本发明中，交通监控视频的获取通过采用固定安装方式的交通监控摄像机实现，摄像机的水平朝向与路面车流方向允许存在一定角度θ,(0≤θ＜90)。以北京市京哈高速北京段的视频监控场景为例，摄像机的架设如图1所示。不失一般性，下面以其中的B摄像机为例加以说明，该摄像机架设在厚俸桥上，距地面高度为8米，向东沿路面方向顺车流拍摄京哈高速北京段的南主路交通场景，场景中的车道数为4条、车道（分割）线为3条。摄像机拍摄帧频为25帧/秒，每帧图像大小为640×480像素，本发明所述跟踪方法采用隔帧执行方式。

本发明所述方法的总体流程图如图2所示，具体包括以下步骤：

步骤1，设置初始参数：首先，设置感兴趣区域为主路路面部分并根据路面标线手工指定3条车道线的初始位置；然后在感兴趣区域内靠近摄像机端设置3条垂直于车道的虚拟检测带，各条检测带宽度（垂直车道线方向的尺寸）均贯穿路面的感兴趣区域、高度（沿车道线方向的尺寸）为9个像素。如图3所示。另外，采用传统四点标定法对摄像机进行标定以获得单应矩阵：

$M=\left[\begin{array}{ccc}33.7333& 10.8956& 178.0000\\ -0.0000& -1.5419& 355.0000\\ -0.0000& 0.0360& 1.0000\end{array}\right]$

步骤2，逐帧读入道路交通监控摄像机的视频序列图像，并在帧图像中基于三级虚拟检测带进行车辆特征点的自动提取和按车辆分组。局部前景提取与活动阴影消除步骤中，对各虚拟检测带内的每个像素均采用5个高斯分布的叠加来进行建模、各像素颜色通道数取3。在自适应的动态背景更新速率α的计算过程中，交通畅通时的更新速率α₀取为0.003；参数β取为3；在计算场景中车辆的平均长度统计值ξ时，帧频F_s为25帧/秒，计入背景的最小阈值T取为0.8。

步骤3，特征点跟踪与分组修正。在稀疏光流跟踪过程中，如果发现第i帧的某个特征点在后继帧中的位置发生了变动，则本发明分别以变动前后的这两个点为中心，在各自的帧图像中分别计算其9×9邻域窗口的Uniform Local Binary Pattern局部二元图纹理特征值并进行匹配。若匹配成功则表明两点处的局部纹理特征一致，位置变动后的特征点有效；否则，表明该特征点在跟踪过程中发生了漂移，此时将其剔除。另外，在跟踪过程中，对少数的误分组特征点进行必要修正并对一些可能出现的干扰点进行筛选剔除过程中，特征点对系统中的系统噪声向量

和测量噪声向量η分别满足

～N(0,Q),P(η)～N(0,R)，Q和R分别为系统噪声和测量噪声的协方差矩阵。转移矩阵A、观测矩阵H为

$\begin{array}{cc}A=\left[\begin{array}{ccc}1& \mathrm{\&tau;}& {\mathrm{\&tau;}}^2/2\\ 0& 0& \mathrm{\&tau;}\\ 0& 0& 1\end{array}\right],& H=\left[\begin{array}{ccc}[1& 0& 0]\end{array}\right]\end{array}$

式中，采样周期τ取为1/25秒。

步骤4，车辆轨迹提取与多个车道的时空轨迹图自动生成。

步骤5，如果处于车道修正阶段，则每新增γ＝50条轨迹进行一次基于轨迹聚类法的车道位置自动修正，否则直接跳转至步骤2。聚类数目k即为车道数取为4；确定本次新增轨迹集中每条轨迹T_χ的上、下近似归属过程中，判定阈值λ取1～1.5，本实施例中基于实验取为1.2.；在第二步更新聚类中心过程中，权值w_l取为0.7。

本实施例中，为4条车道自动生成的驶过车辆的纵向时空轨迹图（TSD₁～TSD₄），以及一幅含所有车道的横向时空轨迹图（TLD）如图8所示。

Claims (1)

1.一种多车道交通时空轨迹图的自动生成方法，其特征在于包括以下步骤：

步骤1，设置初始参数；包括设置感兴趣区域、手工指定各车道线的初始位置、在车流上游位置设置三条垂直于车道的虚拟检测带1、2、3和对摄像机进行参数标定；

步骤2，在帧图像中基于三级虚拟检测带进行车辆特征点的自动提取和按车辆分组；

步骤2.1，基于局部变更新率背景建模提取局部前景并消除活动阴影；

（1）对三级虚拟检测带进行变更新率局部动态背景建模；

对各虚拟检测带内的每个像素，根据其在最近历史帧中的采样值{X¹,X²,…,X^t}采用K个高斯分布的叠加进行建模；然后，对新像素值X^t+1按ω_i,t/σ_i,t降序与该点各高斯分布进行匹配；对匹配成功的高斯分布更新如下：

$\left\{\begin{array}{c}{\mathrm{\&omega;}}_{i,t+1}=(1-\mathrm{\&alpha;}){\mathrm{\&omega;}}_{i,t}+\mathrm{\&alpha;}\\ {\mathrm{\&mu;}}_{i,t+1}=(1-\mathrm{\&rho;}){\mathrm{\&mu;}}_{i,t}+\mathrm{\&rho;}{X}^{t+1}\\ {\mathrm{\&sigma;}}_{i,t+1}^2=(1-\mathrm{\&rho;}){\mathrm{\&sigma;}}_{i,t}^2+\mathrm{\&rho;}{(X^{t+1}-{\mathrm{\&mu;}}_{i,t})}^2\\ \mathrm{\&rho;}=\mathrm{\&alpha;\&eta;}\left(X^{t+1}\right|{\mathrm{\&mu;}}_{i,t},{\mathrm{\&sigma;}}_{i,t})\end{array}\right.$

对未匹配的高斯分布，均值和方差保持不变，权值更新为ω_i,t+1＝(1-α)ω_i,t，α为Log-Sigmoid型动态自适应更新速率，计算公式如下：

$\mathrm{\&alpha;}=\frac{{\mathrm{\&alpha;}}_0}{1+e^{-\mathrm{\&beta;}(l-\mathrm{\&xi;})}}$

式中，α₀为交通畅通时的更新速率；β决定Log-Sigmoid型函数的陡度；l为由步骤3跟踪程序反馈而来的车体区域特征点组位置到虚拟检测带3的最小像素距离；ξ为车辆间距估计值，取为场景中车辆平均长度统计值，基于前2级虚拟检测带中心线生成的当前时刻之前的PVIs计算得到，公式如下：

$\mathrm{\&xi;}=\frac{1}{m}\underset{i=1}{\overset{m}{\mathrm{\&Sigma;}}}\frac{\mathrm{\&Delta;d}\&CenterDot;l_i}{F_s\&CenterDot;(t_i^2-t_i^1)}$

式中，m为两级PVIs中成功匹配的车辆Blob（前景块）个数，l_i为第i个Blob长度，

和

分别为第i辆车在两级PVIs中的Blob中心t坐标，Δd为两级虚拟带的距离，F_s为帧频；

对于与任何分布匹配均不成功的新像素值X^t+1，以X^t+1为均值、以一个较大值为方差构造一个新的高斯分布并代替ω_i,t/σ_i,t最小的分布来完成更新；对各高斯分布按ω_i,t/σ_i,t降序选取前B个分布作为当前点的背景模型，其中：

$B=\arg \underset{b}{\min }(\underset{k=1}{\overset{b}{\mathrm{\&Sigma;}}}{\mathrm{\&omega;}}_{i,t}>T)$

式中，T为计入背景的最小阈值；

（2）提取局部前景并消除活动阴影；

基于背景差分法分割各虚拟检测带上的活动前景，并采用前景像素值与对应位置处的背景像素值的局部纹理特征相似性原理有效去除前景中包含的车辆活动阴影；

步骤2.2，提取特征点并按车辆进行分组；

（1）提取特征点；

利用文献Tomasi方法，在三级虚拟检测带的局部前景图像序列中对驶过的车辆进行易跟踪特征点的逐级递增式提取，充分保障车辆前景特征点的提取效果；

（2）对特征点按车辆进行分组；

基于三级虚拟检测带生成三级前景PVIs并结合Adaboost智能分类器进行特征点的按车分组；具体方法如下：

首先，按下式方法分别将视频序列在三级虚拟检测带中心线上的前景像素信息以时间顺序依次排列生成三级前景PVIs：

$\begin{array}{c}{P}_t^k=X_k^t,k\&Element;\left\{\mathrm{1,2,3}\right\}\\ P_t^k={[P_{t,0}^k,P_{t,1}^k,\&CenterDot;\&CenterDot;\&CenterDot;,P_{t,N}^k]}^T\\ X_k^t=[X_{0,y_k}^t,X_{1,y_k}^t,\&CenterDot;\&CenterDot;\&CenterDot;,X_{N,y_k}^t]\end{array}$

其中，

为第k级PVI的t时刻列像素向量，

为第k级虚拟检测带中心线位置的t时刻像素向量，N为虚拟检测带中心线的像素级长度，y_k为第k级虚拟检测带中心线y坐标；

然后，在目标刚刚驶离第三级检测带时，利用Adaboost智能分类器在原始图像中进行车辆识别并根据识别结果对特征点进行分组；Adaboost智能分类器在原始图像中的施用区域由该目标在三级PVIs中的Blob（下面记为B^k,k∈{1,2,3}）及相互匹配情况来确定；

施用区域左下顶点图像坐标及区域宽度w^Adaboost计算如下：

$\begin{array}{c}{x}_{\mathrm{Lb}}^{\mathrm{Adaboost}}=x_{\mathrm{Lb}}^{B^3}\\ y_{\mathrm{Lb}}^{\mathrm{Adaboost}}=y^3\\ w^{\mathrm{Adaboost}}=w^{B^3}\end{array}$

其中，

为B³外接矩形的左下角顶点图像坐标，y³为第三级虚拟检测带中心线y坐标，

为B³外接矩形的宽度；

施用区域高度h^Adaboost的计算，分为如下三种情况：

第1种情况：目标在三级PVIs中形成的前景块B³、B²和B¹均相互匹配成功，则：

$h^{\mathrm{Adaboost}}=\frac{l^{B^3}\&CenterDot;|y^3-y^2|}{(t^3-t^2)\&CenterDot;F_s}$

其中，

为B³外接矩形的长度，y²为第二级虚拟检测带中心线y坐标，t³和t²分别是目标到达第三、二级虚拟检测带中心线的时刻值，F_s为帧频；

第2种情况：B³只与B²、B¹中的一个匹配成功，记为B'，则：

$h^{\mathrm{Adaboost}}=\frac{l^{B^3}\&CenterDot;|y^3-y^{\&prime;}|}{(t^3-t^{\&prime;})\&CenterDot;F_s}$

其中，y'为B'所在虚拟检测带中心线y坐标，t'为目标到达B'所在虚拟检测带中心线的时刻值；

第3种情况：B²与B¹匹配成功，但均不能与B³匹配，则：

$h^{\mathrm{Adaboost}}=\frac{|y^2-y^1|}{(t^2-t^1)\&CenterDot;F_s}\&CenterDot;(\max (l^{B^3},l^{B^2},l^{B^1})+\mathrm{\&Delta;h})$

其中，Δh为分离距离修正量，构造为

；

（3）同车特征点组的表示；

步骤（2）能将绝大多数的粘连车辆分离并确定各车辆的车体区域；对落在同一车体区域的特征点采用边加权动态图G＝{V,E,W}表述；其中，V为同车特征点集；E为V内各特征点间的优化边集，采用Delaunay三角剖分法构造得到；W为边权集，各边权值取为该边所连特征点对的运动特征相关度，初始化为车体区域长度l_vehicle，并在后继跟踪中动态修改；

步骤3，特征点跟踪与分组修正；

步骤3.1，特征点跟踪；

在后继帧的图像金字塔中计算Kanade-Lucas-Tomasi稀疏光流进行特征点跟踪，过程中若发现某个特征点在后继帧中的位置发生了变动，则分别以变动前后的这两个点为中心，在各自的帧图像中分别计算其n×n邻域窗口的Uniform Local Binary Pattern局部二元图纹理特征值并进行匹配；若匹配成功则表明两点处的局部纹理特征一致，位置变动后的特征点有效；否则，表明该特征点在跟踪过程中发生了漂移，此时将其剔除；

步骤3.2，基于运动特征相关度修正误分组特征点并剔除干扰点；

在跟踪过程中，采用状态空间转移模型来表述特征点对的运动关系并结合Kalman滤波对步骤2结果中存在的误分组特征点进行必要的分组修正，同时对一些可能出现的干扰点进行剔除；仅以x轴方向为例（y轴方向同理），对具有边连接关系的特征点对（i和j），以

表示其在t时刻的运动关系状态向量

，其中

和分别为两点间的x方向距离、距离变化速率和变化加速度，则此特征点对可表述为：

其中，

为t-1时刻特征点对的状态向量，

为t时刻的状态测量向量，

和η分别为系统噪声向量和测量噪声向量，A为转移矩阵，H为观测矩阵；t时刻状态预测值

及状态修正值

可由点对在上一时刻t-1的历史信息递推得到，具体如下：

$\begin{array}{c}{\hat{X}}_{i,j}^{\left(t\right)}=A{X}_{i,j}^{(t-1)}\\ X_{i,j}^{\left(t\right)}={\hat{X}}_{i,j}^{\left(t\right)}+K(Y_{i,j}^{\left(t\right)}-H{\hat{X}}_{i,j}^{\left(t\right)})\end{array}$

式中，K为Kalman修正矩阵；

为误差向量，其第一分量表征t时刻的点对距离测量值与预测值的差异，若该值为正，则表明点对可能出现分离趋势，此时将该边权值减1；若该边权值已减为0，则表明该点对不属于同一车辆，此时断开该点对的连接，实现特征点的分组修正；进而，如果某点已无任何边与其相连，则将其作为干扰点剔除；

步骤4，车辆轨迹提取与多个车道的时空轨迹图自动生成；

步骤4.1，辆车驶出感兴趣区域时提取其轨迹；

车辆驶出感兴趣区域时，基于下述加权平均式计算其V中各特征点的分布中心

$\begin{array}{c}\stackrel{\&OverBar;}{x}=\frac{1}{2n}\underset{1\&le;i\&le;n}{\mathrm{\&Sigma;}}{m}_i{x}_i\\ \stackrel{\&OverBar;}{y}=\frac{1}{2n}\underset{1\&le;i\&le;n}{\mathrm{\&Sigma;}}{m}_i{y}_i\end{array}$

其中，n为V的尺度，x_i和y_i为第i个特征点的图像坐标，m_i为第i个特征点的权重；

进而，计算V中各个特征点i与分布中心

的欧氏距离d_i,p，按照d_i,p升序选取前K个点构成分布中心的K-最近邻点集C并按下式选取第j个特征点的轨迹作为当前车辆轨迹：

$j=\underset{p\&Element;C}{\arg \max }{l}_p$

其中，l_p为K-最近邻点集C中特征点p的轨迹长度；

步骤4.2，根据各辆车的轨迹数据自动生成多个车道的时空轨迹图；

（1）生成纵向时空轨迹图；

为每条车道自动生成驶过车辆的纵向时空轨迹图；以第i辆车为例，其轨迹表述为点序列T_i={<x₁,y₁,f₁>,…,<x_j,y_j,f_j>,…,<x_M,y_M,f_M>}；其中，三元组结构分量x_j、y_j和f_j分别为第j个轨迹点的x轴、y轴图像坐标及轨迹点所在帧号，M为该轨迹中的点数；对各轨迹点<x_j,y_j,f_j>，首先根据坐标值判断其所在车道号k，对于本发明所述交通场景，计算公式为：

$k=q,(x_q^l<x_j<x_q^r{|}_{y_q^l=y_j=y_q^r})$

其中，

为第q条车道左边界线上的点坐标，

为第q条车道右边界线上的点坐标；

然后将该轨迹点绘入第k幅时空图的(t_j,d_j)坐标处：

t_j＝f_j/F_s

d_j＝|Y_j-Y_O|

其中，F_s为摄像机帧频，Y_O为起始参考点的世界Y坐标，在此起始参考点取为第一级虚拟检测带中心线左端点；Y_j为当前轨迹点的世界Y坐标；各点的世界坐标值可由其图像坐标与摄像机标定参数计算得来；

（2）生成横向时空轨迹图；

生成一幅含所有车道的横向时空轨迹图；对第i辆车的轨迹T_i中的各轨迹点<x_j,y_j,f_j>，其在横向时空轨迹图中的坐标(t_j,d_j)满足：

t_j＝f_j/F_s

d_j＝|X_j-X_O|

其中，X_O为起始参考点的世界X坐标，X_j为当前轨迹点的世界X坐标，其余参数同上；

步骤5，若处于车道修正阶段，则每新增γ条轨迹进行一次基于粗糙k-均值轨迹聚类的车道位置自动修正；否则，直接跳转至步骤2；

步骤5.1，在设置初始参数时确定聚类参数k及初始聚类中心；

（1）确定聚类参数k：

聚类数目k即为车道数，对于确定的监控场景k为已知常数，可人工设定；

（2）确定初始聚类中心：

对于各聚类的初始聚类中心c_i，即每条车道的初始中心线，根据步骤1中人工设定的感兴趣区域和各车道线初始位置数据计算得到，具体为：

$x_{c_i}=\frac{x_i^l+x_i^r}{2}{|}_{y_i^l=y_i^t=y_{c_i}\&Element;[y_{i,t},y_{i,b}]}$

其中，i为车道号且i∈{1,2,…,k}，c_i为第i条车道的中心线点集，

、

为c_i中的各点坐标；

为第i条车道左边界线上的点坐标，

为第i条车道右边界线上的点坐标，最两侧车道的外边界即为感兴趣区域外边界；y_i,t、y_i,d分别为当前车道上、下顶点纵坐标；

步骤5.2，基于粗糙k-均值聚类实现车道分割线位置修正；

（1）确定本次新增轨迹集中每条轨迹T_χ的上、下近似归属；

计算轨迹T_χ和各个聚类中心的Hausdorff距离h(T_χ,c_i),i＝1,2,...,k，公式如下：

$h(T_{\mathrm{\&chi;}},c_i)=\frac{1}{N_{\mathrm{\&chi;}}}\underset{a\&Element;T_{\mathrm{\&chi;}}}{\mathrm{\&Sigma;}}\underset{b\&Element;c_i}{\min }\left|\right|(x^a-x^b,y^a-y^b)\left|\right|$

其中，N_χ为轨迹T_χ的长度；

并找出与T_χ距离最短的类m，即

则：

其中，λ为判定阈值、表示第m车道类的上近似集、

表示第m车道类的下近似集；

（2）更新聚类中心，公式如下：

$c_i=\left\{\begin{array}{cc}{w}_l\underset{T_{\mathrm{\&chi;}}\&Element;\underset{\&OverBar;}{C_i}}{\mathrm{\&Sigma;}}{T}_{\mathrm{\&chi;}}/\left|\underset{\&OverBar;}{C_i}\right|+(1-w_l)\underset{T_{\mathrm{\&chi;}}\&Element;(\stackrel{\&OverBar;}{C_i}-\underset{\&OverBar;}{C_i})}{\mathrm{\&Sigma;}}{T}_{\mathrm{\&chi;}}/|\stackrel{\&OverBar;}{C_i}-\underset{\&OverBar;}{C_i}|,& \stackrel{\&OverBar;}{C_i}\&NotEqual;\underset{\&OverBar;}{C_i}\\ \underset{T_{\mathrm{\&chi;}}\&Element;\underset{\&OverBar;}{C_i}}{\mathrm{\&Sigma;}}{T}_{\mathrm{\&chi;}}/\left|\underset{\&OverBar;}{C_i}\right|,& \stackrel{\&OverBar;}{C_i}=\underset{\&OverBar;}{C_i}\end{array}\right.$

其中，i∈{1,2,…,k}，权值w_l∈[0.5,1]、|·|表示基数；

（3）根据聚类后的车道中心线修正车道分割线位置，公式如下：

$x_{{\mathrm{Line}}_i}=\frac{x_{c_i}+x_{c_{i+1}}}{2}{|}_{y_{c_i}=y_{c_{i+1}}=y_{{\mathrm{Line}}_i}\&Element;[y_{{\mathrm{Line}}_i}^t,y_{{\mathrm{Line}}_i}^b]}$

其中，i为车道分割线号且i∈{1,2,…,k-1}；、

、

及

为第i条车道线两侧的车道中心线上的点坐标；

和

分别为当前车道分割线的上、下顶点纵坐标。

Images