2. 北京川速微波科技有限公司, 北京 100080
2. Beijing TransMicrowave Science and Technology Company Limited, Beijing 100080, China
在智能交通系统中, 车道检测是一个长期的研究热点。车道检测包括车道线的检测、道路边界的检测以及车辆可通行区域的检测等。目前, 基于视觉的检测技术[1-2]由于摄像机获取信息量大、成本低等优势应用最为广泛。但摄像机拍摄的图片极易受到光照和天气等外部环境的影响, 对环境条件要求较为苛刻。近年来, 随着雷达探测技术[3-5]的发展, 研究人员开始采用毫米波雷达和激光雷达来代替或者辅助摄像机。雷达不受光照和恶劣天气等环境因素影响, 并且具有探测范围广、测距精度高等优点。
史鹏波[6]利用雷达数据并采用了一种双阈值的方法来提取道路边界点, 但需要预先确定两个阈值, 缺乏自适应性; Xu[7]将雷达获取的数据点分成若干个区域, 计算该区域的随机密度来检测路边, 方法简单, 但计算量大要计算每个区域的协方差矩阵; Han等[8]采用阈值分割和综合概率数据关联滤波器(Integrated Probabilistic Data Association Filter, IPDAF)算法来检测和跟踪道路边沿; 吴维一等[9]采用改进的迭代自组织数据分析算法(Iterative Self Organizing Data Analysis Techniques Algorithm, ISODATA)对雷达数据进行聚类, 虽然算法具有一定的自组织性和启发性, 但还是需要给出先验的最小样本数目和长度约束。
本文提出了一种基于统计和密度特征的核聚类算法(Kernel Clustering algorithm based on Statistical and Density Features, K-CSDF), 从毫米波雷达获取的车辆数据中提取道路信息。无论雷达采用正装还是侧装的方式(见图 1), 该算法都可以对车道进行智能划分, 不需要人为地测量雷达的摆角以及雷达安装位置到车道中心的距离等信息, 提高了人工操作的简便性以及车道划分的准确率。K-CSDF的流程见图 2。
本文的实验载体是北京川速微波科技有限公司的多目标交通测速雷达, 主要用于在测速卡口对车辆进行超速抓拍。该雷达系统主要包括三部分:相机、雷达和补光灯。其中雷达是系统的核心设备, 它能捕获到车辆并触发相机对车辆进行抓拍。
多目标交通雷达采用频移键控(Frequency Shift Keying, FSK)体制, 利用多普勒频移对目标进行测速, 利用不同发射频率的相位差对目标进行测距, 并通过一发两收的天线设计来测量目标的角度。雷达对目标的测距和测角公式如下:
$R = \frac{{c \cdot \Delta \varphi }}{{4{\rm{\pi }}\left( {{f_1} - {f_2}} \right)}}$ | (1) |
$\theta {\kern 1pt} {\kern 1pt} {\rm{ = }}{\kern 1pt} {\kern 1pt} \arcsin \left( {\frac{{\lambda \cdot \Delta {\varphi ^\prime }}}{{2{\rm{\pi }}d}}} \right)$ | (2) |
其中:R为雷达到目标的距离; c为光速; Δφ为同一接收天线两个不同发射频率f1和f2的相位差; θ为雷达天线法向与目标的夹角; λw为雷达发射电磁波的波长; Δφ′为两个不同接收天线的同一发射频率的相位差; d为两个接收天线之间的距离。
在一段时间内, 雷达识别出的车辆目标的行驶轨迹分布如图 3所示, 为了直观起见, 已将雷达获取的车辆目标的极坐标距离信息(R, θ)转换成直角坐标信息(x, y), 即图 3中每个点坐标(x, y)表示车辆的距离信息, y的正负代表车辆行驶的方向。将每个样本点对应的幅度记为z, 所有点的幅度按由小到大排序得到车辆信号的能量分布, 记为q(z)。
对雷达获取的车辆数据, 首先利用车辆目标幅度信息的统计特征对数据进行阈值处理, 剔除掉部分异常数据。以监测来向车为例, 取y<0的数据进行分析, 如图 4所示。
通常, 雷达照射范围内的车辆反射信号很强, 但同时也存在邻近车道的车辆产生的干扰信号, 图 4中“鬼影区”就是受天线的测角范围所限, 由非监测区域的干扰目标所产生的干扰信号, 因此阈值处理的目的就是去掉“鬼影区”的异常数据。
为了提取出有效的车辆信息, 将车辆信号的能量分布q(z)的上分位数α定义如下:
$\alpha = \frac{1}{{{N_q}}}\sum\limits_{z > {z_\alpha }} {q\left( z \right)} $ | (3) |
其中:α表示能量高于zα的样本点的百分比, 0<α<1;Nq表示样本总数。q(z)的统计分布如图 5所示, 呈现出“双峰”特性, 低峰值处表示“鬼影区”的样本分布, 高峰值处表示监测车道区域的样本分布。
对于任一α值, 可以将q(z)分成两组, 取两组数据之间的方差达到最大时的α值作为保留数据的百分比, 即保留能量较高的α·Nq个点, 剔除能量较低的(1-α)Nq个点。假设两组数据的均值分别为λ1和λ2, 分别对应于α和1-α, 则样本总体均值λ为:
$\lambda = \alpha {\lambda _1} + \left( {1 - \alpha } \right){\lambda _2}$ | (4) |
两组数据之间的方差δ2定义如下:
$\begin{array}{l} {\delta ^2}\left( \alpha \right) = \alpha {\left( {{\lambda _1} - \lambda } \right)^2} + \left( {1 - \alpha } \right){\left( {{\lambda _2} - \lambda } \right)^2} = \\ \;\;\;\;\;\;\;\;\;\;\alpha \left( {1 - \alpha } \right){\left( {{\lambda _1} - {\lambda _2}} \right)^2} \end{array}$ | (5) |
在0~1之间改变α, 便能求得式(5) 取得最大值时的α, 此时阈值的取值为能量分布q(z)中的第(1-α)Nq」个点对应的幅值。经实验验证, 即使q(z)的统计分布无明显的“双峰”特性, 这种方法也能很好地剔除“鬼影区”的异常数据。
2.2 基于样本密度特征的动态半径提取对样本数据的特征提取对后续算法以及最终结果具有直接的影响, 更好的特征能够降低模型的复杂度并提高车道划分的准确性。
假设经过上述处理后的数据样本集合为X={x(1) , x(2) , …, x(m); x(i)∈Rn}。其中:x(i)是一个n维的向量, 代表第i个样本的n维信息, m表示样本的数量。
将第i个样本点的局部密度[10]定义如下:
${\rho _i} = \sum\limits_j {\chi \left( {{d_{ij}} - {d_c}} \right)} $ | (6) |
其中:
$\chi \left( x \right) = \left\{ {\begin{array}{*{20}{c}} {1,\;\;\;x < 0}\\ {0,\;\;\;x \ge 0} \end{array}} \right.$ | (7) |
dij表示样本x(i)与x(j)之间的距离, dc是截断距离(Cut-off distance)。
由定义可知, ρi表示与样本x(i)距离小于dc的样本点的个数, 当ρi大于N时, 该样本被视为有效样本。其中dc和N是超参数, 需要人为指定, 参数设置的不同可能会导致结果的较大差异。为了降低算法的参数敏感性, 把dc看作一个变量, 将第i个样本点的动态半径定义如下:
${\tau _i} = \mathop {\min }\limits_{{\rho _i} > N} \left( {{d_c}} \right)$ | (8) |
其中:τi表示样本x(i)达到密度N所需要的最小半径(如图 6所示), τ的值越小, 表明该样本点越可能为有效点; τ的值越大, 表明该样本点越可能为噪声点。
本章主要介绍K-CSDF的第二步:对提取出的有效点进行聚类。首先分析了传统K均值算法的不足, 然后通过引入核函数改进了原有算法, 最后给出了K-CSDF的聚类实现流程。
3.1 K均值算法的不足K均值算法[11]是一种简单、高效的动态聚类算法, 其时间复杂度接近线性, 因此在工业中有广泛的应用。
K均值算法采用迭代的思想, 利用最小误差平方和准则来判断失真函数(Distortion function)是否收敛, 定义失真函数如下:
${J_{c,\mathit{\boldsymbol{\mu }}}} = \frac{1}{m}{\sum\limits_{i = 1}^m {\left\| {{\mathit{\boldsymbol{x}}^{\left( i \right)}} - {\mathit{\boldsymbol{\mu }}_{{c^{\left( i \right)}}}}} \right\|} ^2}$ | (9) |
其中:
${\mathit{\boldsymbol{\mu }}_j} = \frac{1}{{{l_j}}}\sum\limits_{\mathit{\boldsymbol{x}} \in class\;j} \mathit{\boldsymbol{x}} ;j = 1,2,...,k$ | (10) |
μj是第j类样本的均值向量, 代表第j类的聚类中心; lj表示第j类样本的数量; 下标c(i)表示第i个样本的类别标签; k表示类别数。当Jc, μ取得最小值时的聚类就是误差平方和准则下的最优结果。
在本文的应用场景中, 直接采用K均值算法来进行聚类是不合适的。由于该算法采用欧氏距离来定义样本间的相似性, 并用均值来更新聚类中心, 只有当类内样本分布为超球状或接近超球状时, 才能取得较好的效果。另一种距离度量方法是采用闵可夫斯基距离(Minkowski distance), 设两个样本为(a1, a2, …, an)和(b1, b2, …, bn), 则它们之间的闵可夫斯基距离定义为:
$d = {\left( {\sum\limits_{i = 1}^n {{{\left| {{a_i} - {b_i}} \right|}^p}} } \right)^{1/p}}$ | (11) |
图 7显示了当p取不同值时, 样本逼近聚类中心的趋势。但采用闵可夫斯基距离也不能解决所有问题, 更一般的距离或相似性度量方式可以通过引入核函数[12]的方法来实现。
通过上述分析可知, 设计一个有效的核函数是K-CSDF第二步的关键。在本文应用场景中, 样本数据是由车辆在车道内行驶而产生的, 因此样本数据的特点是集中在相应的主轴方向, 即车道中心线的方向。因此, 定义主轴核函数如下:
${K_j}\left( {{\mathit{\boldsymbol{x}}^{\left( i \right)}},{\mathit{\boldsymbol{U}}_j}} \right) = {\mathit{\boldsymbol{U}}_j}^{\rm{T}}{\mathit{\boldsymbol{x}}^{\left( i \right)}}$ | (12) |
其中:Uj是样本类内离散度矩阵Sj的最大特征值所对应的特征向量。
${\mathit{\boldsymbol{S}}_j} = \sum\limits_{\mathit{\boldsymbol{x}} \in class\;j} {\left( {\mathit{\boldsymbol{x}} - {\mathit{\boldsymbol{\mu }}_j}} \right){{\left( {\mathit{\boldsymbol{x}} - {\mathit{\boldsymbol{\mu }}_j}} \right)}^{\rm{T}}}} $ | (13) |
相应地, 可以将样本与核函数之间的距离ζ定义如下:
$\xi \left( {{\mathit{\boldsymbol{x}}^{\left( i \right)}},{K_j}} \right) = {\mathit{\boldsymbol{\eta }}^{\rm{T}}}\mathit{\boldsymbol{\eta }}$ | (14) |
$\mathit{\boldsymbol{\eta }} = \left( {{\mathit{\boldsymbol{x}}^{\left( i \right)}} - {\mathit{\boldsymbol{\mu }}_j}} \right) - {\mathit{\boldsymbol{U}}_j}{\mathit{\boldsymbol{U}}_j}^{\rm{T}}\left( {{\mathit{\boldsymbol{x}}^{\left( i \right)}} - {\mathit{\boldsymbol{\mu }}_j}} \right)$ | (15) |
在定义了代表不同类的核函数以及样本与核函数之间的距离之后, 就可以参照K均值算法来构造K-CSDF的聚类部分, 具体流程如下。
输入 样本集合{x(1) , x(2) , …, x(m); x(i)∈Rn}
初始化 将样本x(i)初始化成k类; 随机初始化每类的核Kj。
Repeat
对每个样本x(i)利用式(14) 计算出它到初始核的距离ζ, 取最小距离并将样本归为c(i)类;
对所有归为c(i)类的样本利用式(13) 求出其Sj和Uj, 并更新核Kj
Until
首先, 给出K-CSDF特征提取的实验仿真流程, 见图 9。经过K-CSDF第一步的特征提取后通常可以得到2000个左右的有效样本点, 取100个有效点进行聚类, 最终的聚类结果见图 10。从图 10中可以看出算法识别出的三条车道中心线是由每一类的样本在核方向上的投影产生的。
其次, 分别取样本点数量为100、500、1000、2000对聚类效果进行分析。此时, 聚类结果的实验仿真如图 11所示, 通过对100组实际采集的数据进行分析, 结果表明即使只取100个样本点, 算法仍能很好地识别出车道中心线。表 1给出了图 11的聚类结果所对应的聚类中心μj, 样本类内离散度矩阵Sj和其最大特征值所对应的最大特征向量Uj的数值。
此外, 在实际采集数据时很容易遇到的一个问题就是, 各个车道的过车数量可能并不均匀, 这会导致采集到的样本分布不均且有明显的“间断”。为了衡量聚类效果, 定义评价指标eva如下:
$eva={\left( {{\beta }_{i}}-{{\alpha }_{i}} \right)}/{\max \left( {{\alpha }_{i}},{{\beta }_{i}} \right)}\;$ | (16) |
其中:αi表示第i个样本到此类的其他样本的平均距离; βi表示第i个样本分别到其他各类样本的平均距离中的最小值; eva值在-1~1范围内, 越接近1表明聚类效果越好。如图 12所示, 仿真结果能识别出车道中心线并且大部分样本的eva值都大于0.6, 表明算法能很好地适应这种情况。
最后, 考虑到国内该类产品仍处于研发阶段并涉及商业机密, 关于多目标交通雷达的车道划分理论研究还未见公开的文献以及实际的工程结果, 因此将本文提出的K-CSDF和另外两种具有代表性的聚类算法在算法用时以及车道划分的准确率上进行对比:
1) 高斯混合模型(Gaussian Mixture Model, GMM):相对于K均值算法强制地将每个样本分给某个类, GMM算法给出的是样本分到每个类的概率, 因而又称作软聚类(soft clustering)。图 13是GMM的聚类结果以及样本分类的后验概率。
2) 自组织映射神经网络(Self-Organizing Maps, SOM):采用2×3的拓扑网络将样本分成6类, 训练次数取200次。SOM网络的聚类结果给出了6个类的中心, 如图 14所示, 分别取左、中、右拓扑结构的两个聚类中心的连线作为识别出的车道中心线。
实验仿真显示三种聚类算法在大多数情况下都能达到90%以上的分车道正确率。其中, 本文算法和SOM算法可以给出车道中心线, 而GMM算法不能; 在样本分布不均衡的情况下, 本文算法也具有很好的鲁棒性, 仍能达到95%以上的正确率, GMM算法可以达到90%左右的正确率, 而SOM算法无法正确分类; 在算法用时上, 以取1000个样本点为例, GMM算法用时最快, 在0.2 s左右, 本文算法需要0.8 s左右, SOM算法需要2.5 s左右。具体对比见表 2和图 15。
测试设备:PC、雷达、相机和三脚架等。
测试地点:天桥。
测试步骤:
1) 将雷达用三脚架正装, 并连接相机和PC, 然后用雷达采集车辆数据5 min(约20辆车)。
2) 通过上位机发送命令, 执行车道划分算法, 通过聚类结果的聚类中心和特征向量对车道进行拟合, 上位机界面见图 16。
3) 将雷达设置成工作状态, 对车辆进行正常抓拍并保存原始数据、抓拍照片和视频用于统计分析。
4) 对雷达进行侧装, 重复上述三个步骤。
对三组测试结果分别统计分车道正确率, 并对每辆过车取10帧数据进行单帧分析, 如图 17所示。路测统计结果见表 3。
表 3中, 第1组数据为正装时采集, 第2、3组为侧装时采集, 可以看出正装的车道划分正确率要稍高于侧装时的正确率, 但从总体上, 该方法在两种安装方式下都可以达到95%以上的车道划分正确率, 可以满足实际应用的需求。
5 结语本文提出了K-CSDF用于车道划分, 该算法主要包括对原始数据的特征提取以及基于核的相似性的动态聚类两步。从实验仿真和路测结果可以看出, 该方法在保证实时性的同时, 可以达到95%以上的分车道正确率;当只取100个样本点进行拟合时, 算法也具有很好的鲁棒性。但本文只对监测3个车道的情况进行了分析, 对于监测更多车道以及更复杂环境下的道路情况, 还需要后续的研究。