谐波叠加法在道路谱重构中的改进与应用研究

钟岳， 张纬华， 徐峰， 方虎生

（1. 中国人民解放军32184部队，北京 100072；2. 陆军工程大学，南京 210007）

路面不平度是汽车行驶最主要的激励，对汽车舒适性、可靠耐久性等有直接影响，是汽车底盘部件磨损、疲劳失效的主要外因。采用虚拟试验来检验并预测样本的可靠性、安全性以及舒适性，可以减少汽车开发的周期从而节省开发成本。道路虚拟试验技术的关键是获取准确可靠的路面激励信号，即道路路面谱。道路谱是将实际路面的形状和特性参数用数字化的形式表示出来，以用于理论研究和试验，它表示的是路面表面偏离理想平面的幅度的一种时间域或空间域的模型。在工程实际中，重构路面的准确性是关键，结合实测路面数据进行路面重构能保证这一点。

目前，国内外对道路谱获取方法的研究主要有直接测量法、仿真模拟法和载荷谱迭代法。直接测量法测量精度相对较高。黄亮等［1］采用nCode 软件，对采集的路谱进行预处理得到准确可用的原始采集数据，将实测的较多工况路面载荷重构等效为部分工况进行后续的分析，取得了较好的成效，但是测量得到的只是单一道路轨迹。载荷谱迭代法测试仪器原理简单，操作方便，成本较低，但是误差源多，数据精度受迭代误差影响较大，周期长，迭代原理较复杂。仿真模拟法包括谐波叠加法、滤波白噪声法、Fourier逆变换生成法等。仿真模拟法不需要测试设备，成本很低，周期短，仿真结果可直接导入道路模拟系统。其中，采用滤波白噪声作为路面输入的算法具有简单易行的优点，然而其生成的时域路面随机模型仅能用于拟合型的路面谱结构，具有局限性。谐波叠加法是一种离散化数值模拟路面的方法，能满足路面不平度的随机各态特性，在道路谱获取方面有一定的优势。关于谐波叠加法在道路谱重构中的应用，相关学者进行了一些研究。聂彦鑫等［2］应用谐波叠加法对路面不平度随机序列进行仿真。结果表明通过谐波叠加法能很好地拟合路面不平度，且AR 模型谱估计具有拟合结果平滑且精度高的特点。刘大维等［3］基于谐波叠加法和三角网格法的基本原理，建立了三维随机路面数学模型。程国勇等［4］提出了基于谐波叠加法的道面不平度模拟方法，利用功率谱估计的方法对模拟路面谱和标准路面谱进行了比较，结果二者是一致的，表明谐波叠加法能很好地拟合路面不平度。此外，谐波叠加法在电力、风场谱重构等方面的应用也取得了较好的效果。例如，CALVO 等［5］提出了大正则系综的谐波叠加方法；何甜田等［6］提出了一种自适应谐波叠加的复调音乐主旋律提取算法，自适应叠加谐波构建显著函数；骆光杰等［7］推导了模拟脉动风速的谐波叠加法，利用Matlab数学平台编写了具备时空相关性的脉动风速时程模拟程序。由此可知，与直接测量法、载荷谱迭代法等相比，谐波叠加法具有显著的优势，为复杂势能图重构提供了理论框架。目前关于谐波叠加法的研究较多，但是其在道路谱重构中的应用研究以理论分析为主，缺乏对谐波叠加法不足的关注以及优化方案，对于谐波叠加法在道路谱重构中的应用效果也还有待检验。鉴于此，本文采用谐波叠加法进行道路谱重构，针对谐波叠加法的不足，提出谐波叠加法的改进方案，并应用道路谱重构，通过数值模拟仿真验证改进谐波叠加法在道路谱重构中的优势，试图做出新的探索。

1 基于谐波叠加法的道路谱重构

道路谱重构主要有两种思路，一是通过实际测试方法，记录路面不平度数据并根据测得数据进行重现；二是进行时域信号重构。在路面实测信号点预处理基础上，进行功率谱密度分析，得到基于谐波叠加法的路面功率谱密度，从而实现道路谱重构。

1.1 路面实测信号点预处理

本文以纵向距离为100 m 的普通砂石路面为对象进行高程数据采集，得到了相关原始数据，图1分别展示了砂石路面实车测试图和所得到的路面实测信号数据。

图1 砂石路面实车试验结果

在数据采集过程中，受到外部信号的干扰，由于随机干扰信号的频带较宽，导致采集的离散数据形成尖峰产生了较多的异常值点。为此，利用莱茵达准则对其进行异常值剔除。

将测量值中与平均值相差大于三倍标准差的点进行剔除，从而得到了相应的信号点，如图2所示。

图2 异常值剔除后的路面数据

与原始采集信号相比，经过莱茵达准则处理后，信号中的异常值明显减少，然而，可以发现，图中仍含有部分的趋势项数据，即信号中周期大于采样长度的频率部分，这是由于传感器等数据采集仪器在采集过程中出现的零点漂移、基础运动等引起了信号波形偏移。在试验中，利用最小二乘原理建立趋势项多项式的一般模型，从而进行趋势项去除，图3则展示了去趋势项后的结果。

图3 路面实测高程信号去趋势项后结果

在此基础上，选用小波阈值方式对其进行去噪，最终利用小波重构方式得到了所需的“纯净”信号。常用的阈值选择准则有4 种［8］，分别为无偏似然估计（rigrsure）原则、固定阈值（sqtwolog）原则、启发式阈值（heursure）原则以及极值阈值（minmax）原则，图4 则展示了4 类阈值原则下的去噪结果。

图4 4类小波阈值去噪方式结果对比

分别计算4种去噪结果后的信噪比（SNR），并比较大小，见表1。

表1 四种阈值原则去噪后的信噪比

根据信噪比结果对比，初步得出结论：基于rigrsure原则的小波去噪方式效果最佳。

1.2 基于谐波叠加法的路面不平度

谐波叠加法源于谱估计的周期延拓，是一种以离散谱逼近目标随机过程的模型，具有较严谨的数学推导过程。

首先，利用平稳随机过程将不同空间频率范围内的路面位移功率谱密度Gq(n)转化为路面不平度方差。

然后，通过中心频率代替整个频段方式得到了路面不平度方差在离散空间下的求和方法。

对于频率为nmid_k(k= 1，2，……，m)的各频带，通过正弦波函数叠加方式从而得到了空间路面的不平度，如式（1）所示。

式中：x为路面纵向距离；θk为［0，2π］区间内均匀分布的随机数；Gq(nmid_k)为各频段的中心频率处的功率谱密度；Δnk为频段数。

考虑到车辆在行驶时产生的左右轮横向路面不平度差异性，将谐波叠加法拓展至三维空间，从而得到了三维空间条件下的路面不平度函数，如式（2）所示。

式中：x为纵向距离；y为横向距离；θk(x，y)为［0，2π］内均匀分布的二维随机矩阵；F(x，y)表示由路面横纵向位置决定的位置参数，设为F(x，y) =(x2+y2)12。

2 基于改进谐波叠加法的道路谱重构

2.1 三次样条插值法

谐波叠加算法模型是一个平稳的各态历经、零均值的随机过程，广泛应用于路面谱［9］、风场谱重构［10］、复杂桥梁结构方程［11］求解等多个数值模型问题。

样条插值，是指在某一节点上，用一组已测量的点，在各节点间隔上构造三次多项式函数，并在一定的条件下，得到的曲线或曲面的一次导数和二次导数都是存在的，也就是说，拟合出的曲面或曲面为连续且光滑。因此，在实际应用中，样条插值是一种较普遍的拟合逼近方法。相较于传统的谐波叠加或特征参数描述等数值模拟方法，引入三次样条插值函数的最大优势在于在不影响模拟精度的前提下，通过压缩道路互功率谱矩阵的正定平方分解（Cholesky）次数来减少计算内存，大大地缩短模拟计算所花费的时间，从而提升道路谱不平度模型的拟合时间。因此，采用三次样条插值技术对谐波叠加法进行改进，以提高道路谱重构精度。

卵石路、砂石路等不规则路面的侧向特性具有很强的随机性，同时原始数据不够，分辨率低。因此，利用三次样条插值法，对采集得到的高程数据进行插值，如图5 所示。将多个测量点的高程资料以平均路径S为水平间距。通过三次样条插值，获得了整个道路的高程资料，图中K（1）为传感器采集的第1 组数据，根据采集的高程数据采用Matlab 编程软件进行5 mm 间隔横向插值，插值出第二组数据K（2），同理可以插值出第3 组数据K（3），进而得到整条路面的数据点信息。使建立的数字路面不平度与实际路面不平度一致。

图5 三次样条插值方法

2.2 中心频率点的互功率谱获取

根据中心频带方差表示法，得到了频率为nmid_k(k= 1，2，……，m)的各个频带，再根据式（3）计算方式可以得到各中心频带之间的互功率谱密度。

式中：Gii(n)为第i个中心点的功率谱密度；Gjj(n)为第j个中心点的功率谱密度；cohij为两点之间的相干函数。

陈水生等［12］曾根据汽车左右轮之间的相干性拟合出基于轮距参数的相干函数模型，如式（4）所示。

式中：(n1，n2)为有效频带范围；B为左右轮轮距；ρ为经验系数，通常取为1。

根据式（4），可以得到关于各频带中心频率点的相干函数表达式，如式（5）所示。

将式（5）代入至式（3）中，从而得到了关于频带中心点间的互功率谱模型。

2.3 谐波叠加函数获取

通过轮距之间的相干函数得到了各频带中心的互功率谱密度模型。将其转化为矩阵形式，从而得到了各频带中心频率处的互功率谱密度矩阵，如式（6）所示。

利用Cholesky 分解可以关于矩阵G(n)的下三角矩阵，Cholesky分解如式（7）所示。

式中：LT(n)为L(n)的转置矩阵，且为下三角矩阵，L(n)如式（8）所示。

相较于基于插值技术的风速场谐波叠加，此类方法最大的优势在于获得了更多具有“插值”属性的离散点，因此，下一步骤将直接对这m个中心频率点进行多项式插值并获得相对应的下三角矩阵L(n)，令L′(nmid_k)=Mk，s=(nmid_d-nmid_l)/k。从而可得到互谱分解矩阵L′(nmid_k)关于s的插值方程，如式（9）所示。

由式（9）可知，若确定了谱密度Cholesky 分解矩阵L′(nmid_k)的二阶导数Mk，则此p元方程可以进行求解，根据三次样条函数性质可知，内部的m- 2 个中心频率点处一阶导数连续，从而确立了m- 2个约束方程，而针对两边界处的中心频率点，试验中进行了两类边界条件下的插值结果对比。

边界条件1：两边界处的中心频率点nmid_l和nmid_m处的二阶导数连续，即：

边界条件2（拉格朗日边界条件）：令两边界处中心频率点的斜率与通过前后各4 个离散点的三次函数在端点处的斜率相同，可表示为：

为验证样条插值技术试验结果的可靠性，试验同时引入了多项式插值技术，多项式插值线性插值的推广，是利用更高阶的多项式对原始的离散中心频率点进行函数拟合，而对于m个中心频率点，则建立m- 1次多项式，如式（12）所示。

根据式（7）～（10）得到了拟合后的中心离散点互谱分解矩阵Lij(nmid_n)，再代入式（6）得到了互功率谱密度矩阵Gij(nmid_n)，最后将其代入式（2）得到了关于路面不平度的插值谐波叠加函数。

3 基于改进谐波叠加法的道路谱重构仿真

为说明以上基于改进谐波叠加法的道路谱重构算法的有效性，选择了普通砂石路面的采样路段进行虚拟重构。

3.1 参数设置

选取了长度为100 m，宽度为4 m 的普通砂石路面模型。图6 为采用韦尔奇法对图3 进行功率谱分析所得到的曲线，通过观察该路面的功率谱曲线可以得到砂石路面不平度信号的主频为0.117 2 Hz、功率谱密度峰值为（1 589×10-6）m3。

图6 原始路面功率谱密度曲线

通过消除趋势项、去噪处理等方法，获得一维的自相关系数。基于道路的等向性特征，结合坐标转动转换的基本原理，利用一维自相关函数进行插补，得到了二维自相关函数矩阵。重构路面一维时域序列如图7所示。

图7 重构路面一维时域序列

砂石路三维模型如图8所示。

图8 砂石路三维模型

3.2 路面谱重构流程

根据普通砂石路面得到的道路谱重构模型，尝试将其推广到各类强化试验道路。为制定更加具有泛化性的道路谱重构模型，试验中利用Matlab/GUI建立了基于第1 类边界条件下三次样条插值技术的谐波叠加重构模型，其流程如图9 所示，其中N表示横纵向距离上的采样频率。

图9 三次样条插值谐波叠加法的随机路面重构流程

3.3 路面谱重构仿真

结合三次样条插值谐波叠加法的随机路面重构流程，对基于改进谐波叠加法的道路谱重构效果进行仿真分析。

首先，选取了长度为100 m 的普通砂石路作为数值模拟对象，并根据不同的中心频率划分范围分别抽取了10，20，30，40，50 个中心频率点进行比较分析，表2 则比较了选取不同数量的中心频率点时，各类模型重构方法的时间。

表2 四类方法模拟路面不平度函数运行时间 单位：s

根据不同的路面不平度函数重构方法运行时间结果可以发现，传统的谐波叠加法运行时间最多，多项式插值技术次之，样条插值技术所耗时间最少，并且随着模拟点数的不断增长，插值技术的时间缩短更加明显，从而验证了功率谱密度经Cholesky 分解后在优化插值叠加函数效率中的可行性；而比较表中两类边界条件下的样条插值技术运行时间发现，其运行时间随模拟点数的增长较相似，所以仅仅根据重构时间难以得出可靠性结论，需要对四类方法所得到的路面不平重构函数的准确性进行分析。

然后，运用相关系数法在级分布图中，对实测路面原始数据功率谱和重构路面功率谱做验证分析。为提升各类重构方法功率谱与原始路面的功率谱密度的直观差异度，将A 级至H 级道路谱［13］作为参考从而得到了原始路面的位移功率谱密度-空间频率图，如图10所示。

图10 砂石路原始位移功率谱密度与谐波叠加重构数据对比图

对其余三类方法所得到的重构功率谱密度函数进行同理变换，得到了三类重构谱与原始谱在空间频率上的变化对比，如图11所示。

图11 三类插值重构法与原始数据功率谱密度对比图

利用式（15）分别计算谐波叠加法以及三类插值重构谱与原始谱的相关系数r，得到结果见表3。

表3 4类重构谱与原始谱的相关系数比较

根据图10～11 所示的功率谱密度-空间频率对比曲线可以发现：相较于原始的谐波叠加重构法，三类插值重构法在与原始功率谱数值的相关程度上均获得了较大的提升，且大大提升了主频、峰值与原始谱的距离，从而验证了插值重构方法在道路谱重构模拟工程中的可行性。

同时将三类插值重构法重构精度与直接测量法、载荷谱迭代法、谐波叠加法以及滤波白噪声法等进行对比。重构精度计算式为：

式中：xi为原始道面不平度数据样本；yi为重构路面不平度数据样本；m为正整数。

结果见表4。

表4 道路谱重构精度对比

由表3 可知，三类插值重构法的重构精度明显优于现有方法，具有显著的优势。

在此基础上，利用Matlab可以做出上一节中利用改进谐波叠加法重构获得的道路不平度模拟信号离散序列的功率谱密度曲线，将其与目标信号的功率谱密度曲线进行对比可以验证信号重构的效果。GB/T7031—2005 是描述道路谱特性的最新国家标准，给出了从位移功率谱出发的典型道路的近似分类方法并提供了常规道路的分类。进一步利用标准路面谱生成模拟路面后，对生成的路面进行评估，分析本文模拟路面的功率谱与标准路面谱的接近程度［16］。标准功率谱曲线与重构信号的功率谱对比，如图12所示。

图12 标准功率谱与重构信号功率谱

由图12 可知，改进谐波叠加法重构的功率谱能很好地逼近目标谱，较好地符合了国标中规定的标准功率谱，满足了重构要求。

最后，利用标准路面谱生成模拟路面后，对生成的路面进行评估，分析模拟路面的功率谱与标准路面谱的接近程度。具体是采用经典的周期图法和现代AR 参数模型进行功率谱估计［17-18］，验证改进谐波叠加法重构的有效性。

周期图法估计的改进谐波叠加法重构功率谱密度如图13所示。

图13 周期图法估计的功率谱密度

由图13 可知，采用周期图法估计的功率谱得到的仿真结果与理论基本上一致。谱线的分辨率较高，但是离散度较大，资料波动较大，方差较大，特别是在高频部分，曲线较粗糙，信息十分丰富。

AR 参数模型估计功率谱与理论谱对比如图14所示。

图14 AR参数模型估计的功率谱密度

由图14 可知，AR 参数模型估算的功率频谱与理论频谱具有很好的一致性，而且曲线的平滑程度明显高于传统频谱估算，其处理结果的方差比经典谱估计方法的小。

无论是采用经典周期图法谱估计，还是现代AR 参数模型谱估计，得到的结果与标准路面理论谱具有很好的一致性，从而验证了本文所提出的改进的谐波叠加算法的正确性和实用性，能较好地符合路面谱重建的需求。同时，可以根据不同的参数设定，使模拟结果更理想。

4 结论

本文对谐波叠加法在道路谱重构中的应用进行研究，结合谐波叠加法在道路谱重构中的优势与存在的不足，采用三次样条插值技术对谐波叠加法进行改进，提出改进的道路谱重构方法，并进行仿真检验，主要结论如下。

1）针对传统谐波叠加法重构运行时间较长以及精度较低这一现象，本文引入了样条插值函数中的正定平方分解（Cholesky）这一方法，通过压缩道路互功率谱矩阵的Cholesky 次数来减少计算内存，利用Cholesky 分解得到了互谱密度分解矩阵，提升了道路谱模型的重构拟合时间，提高了道路谱重构精度。

2）利用Matlab/GUI 建立基于第一类边界条件下三次样条插值技术的谐波叠加重构模型，比较插值重构模型与原始道路谱模型的相关系数，并对常见道路谱重构方法的精度进行对比。研究发现，相较于传统的谐波叠加重构方法，其平均重构运行模拟时间可减少37.22%，精度提升可达14.1%。与直接测量法、载荷谱迭代法、谐波叠加法以及滤波白噪声法、Fourier逆变换生成法等进行对比，改进的谐波叠加法道路谱重构精度在90%以上，具有显著优势。并且，采用经典周期图法和现代AR 参数模型谱估计，基于改进谐波叠加法的道路谱重构模拟仿真重构的路面功率谱与标准路面谱能很好地契合，得到更加满意的仿真效果。

总之，本文提出的基于谐波叠加法的道路谱重构方法为车辆虚拟道路试验提供了一种比较有效的谱重构模型方法，不仅在道路谱重构中有重要的应用价值，同时在随机风场模拟、牵引供电系统、车辆荷载模拟等方面也具有广阔的应用前景。