0 引言

铁路既有线的实际线形参数是工务部门进行轨道维护的基准。既有铁路在长期运营和维护过程中，其轨道线形会不断发生变化，因此需要对既有铁路的轨道线形进行定期复测^[1]。基于复测数据，通过合理的优化方法对轨道平面线形进行分段并估计最佳线形参数，有利于在最小施工量的前提下指导轨道维修达到平面最优线形。

目前轨道实际线形参数的检测，主要采用全站仪配合轨检小车获取轨道实际线形三维坐标数据^[2]。针对这种检测方法，文献^[2]提出依据相邻30 m弦之间夹角变化进行平面线形概略分段，该方法要求轨道自身形变小且测量精度高，因此实用性差。文献^[3]提出基于圆曲线最小二乘拟合方法优化平面线形，该方法需要人工根据正矢变化判断分界点，而且在大半径曲线优化时容易产生病态结果。文献^[4]提出利用三次样条曲线拟合测点曲率求得初始值，再迭代优化平面线形参数的方法，该方法需要人工判断初始分段点，并且在优化结果不满意时需人工多次修改分段点，自动化程度低。文献^[5]提出基于正交距离最短的平面线形拟合方法，该方法只拟合直线和圆曲线段，未考虑缓和曲线对施工量的影响。文献^[6]提出以夹直线的最小二乘拟合为切入点优化平面曲线，该方法需要依据超高变化人工识别分段点。文献^[7]提出采用三次样条函数与稳健估计相结合的方法识别分段点，该方法对采样精度要求高。文献^[8]提出平差法，文献^[9]提出方向加速法，但均是在分段点已知情况下优化平面线形参数，并未提出如何对轨道平面线形进行分段。另一方面，基于全站仪的轨道检测速度慢、需要建立铁路工程基桩控制网(Control network for railway Piles,CPIII)并长期维护，增加了轨道维护工作的成本和复杂度。

文献^[10-11]分别提出利用陀螺仪测量角度的高精确性配合里程计检测轨道实际线形的方法，说明了惯性器件检测轨道线形的可行性，但并没有提出相应的线形参数优化算法。文献^[12]提出了一种基于惯性器件的轨道检测新方法，该方法通过量测惯导角度数据检测轨道平面线形参数，具有检测效率高、精度高、结果重复性好等优点，为轨道检测领域开辟了新的发展方向；该文献提出根据方向角变量曲线人工分段，再用最小二乘法拟合线形参数的方法，该方法因存在人工介入数据处理过程，有人工误判进而影响参数估计结果的可能，并且其得出的并不是优化的轨道线形，不能保证施工量最小。

因此，基于新检测方法得到的包含轨道不平顺及量测噪声的轨道方向角、里程数据，通过合理的优化方法对轨道平面线形进行分段并估计最佳线形参数是一项十分重要的工作。为了解决上述问题，本文将轨道平面最佳线形确定建模成优化问题，提出首先根据定长曲率曲线最小二乘拟合斜率变化对轨道进行概略分段，然后使用组合迭代算法进行精确分段并拟合最佳线形参数，最后通过仿真算例及应用实例验证了本算法的有效性。

1 轨道平面最佳线形估计问题建模

轨道线形设计中最常见的类型为对称基本型，该线形由直线-缓和曲线-圆曲线-缓和曲线-直线组合而成，且前后缓和曲线长度相同。分析轨道线形沿里程的方向角特征有：直线段方向角保持不变，圆曲线段方向角线性变化，缓和曲线段方向角呈二次抛物线变化，如图 1所示。

本文将轨道最佳线形估计问题建模为以各量测点方向角与拟合方向角偏差平方和最小为目标函数，通过优化轨道线形参数确定轨道最佳线形。考虑到沿里程的轨道测点曲率变化规律明显，且直线拟合相比高次曲线拟合误差小、算法简单，因此本文从曲率-里程曲线入手对轨道进行分段，先拟合轨道测点曲率，然后积分得到拟合方向角，再优化目标函数以确定轨道最佳线形参数。

目标函数：min(f(ds,X)) = ∑Δ_i²，其中ds是自变量，其意义为轨道每一量测点里程及该里程处拟合曲率值的有序集合{(m₁,k₁),(m₂,k₂),…,(m_n,k_n)}，对ds沿里程进行积分可得S,其意义为轨道每一量测点里程及该里程处拟合方向角的有序集合{(m₁,α₁),(m₂,α₂),…,(m_n,α_n)}，由S及ds可判别轨道线形，计算轨道线形参数(线形长度、圆曲线半径等)，线形起、终点处的里程，因此可以得到全线轨道设计参数；X为量测数据里程、方向角的有序集合，X={(m₁,α₁),(m₂,α₂),…,(m_n,α_n)}；用Δ_i = α_i-α_i′表示测点方向角与拟合方向角的偏差。目标函数以最小化各量测点方向角与拟合方向角的偏差平方和为目标，其解为轨道各量测点最佳方向角序列，由其可计算轨道最佳线形参数。

约束条件应满足轨道线形曲率变化规律要求，以及铁路轨道设计规范所规定的各类参数最小值要求，且对称基本型线形应满足前后缓和曲线长度相等。因此有约束条件：

$\left\{ \begin{align} & {{k}_{Z}}=0 \\ & {{k}_{ci}}={{\varepsilon }_{i}} \\ & {{k}_{h(i-1)}}=a\times {{\varepsilon }_{i}}+{{b}_{(i-1)}} \\ & {{k}_{h(i-1)}}=-a\times {{\varepsilon }_{i}}+{{b}_{(i+1)}} \\ & R\ge {{R}_{\min }} \\ & {{l}_{R}}\ge {{l}_{{{R}_{\min }}}} \\ & l\ge {{l}_{\min }} \\ & {{l}_{former}}={{l}_{latter}} \\ \end{align} \right.$

式中：k_Z=0表示直线段的曲率为0; k_ci=ε_i表示第i段圆曲线的曲率为一定值ε_i，k_h(i－1)=a×ε_i+b_(i－1)表示第i段圆曲线的前缓和曲线上各测点曲率在0到ε_i间线性变化，k_h(i+1)=－a×ε_i+b_(i+1)表示第i段圆曲线的后缓和曲线上各测点曲率在ε_i到0间线性变化；R_min、l_Rmin、l_min分别为设计规范规定的最小曲线半径、最小圆曲线长度和最小缓和曲线长度；l_former、l_latter分别表示同一圆曲线前后两段缓和曲线的长度。

2 模型复杂性分析

轨道最佳线形参数确定模型的复杂性在于不同线形参数之间具有耦合作用，某一线形参数的微小变化都会使得线路整体发生改变，所以需要不断改变自变量ds中的有序集合{(m₁,k₁),(m₂,k₂),…,(m_n,k_n)}，每次改变后重新计算目标函数值min (f(ds,X)) = ∑Δ_i²直到最优。因此线路长度越长，其解空间规模越大，处理起来越复杂。

穷举法是处理此类搜索最优解问题的一种常用方法。定义一段完整的直线-缓和曲线-圆曲线-缓和曲线-直线线形组合为一基本线元组。一基本线元组中的任一测点在满足设计规范前提下都有可能成为分段点，因此在测点数量庞大的情况下，穷举法需要遍历的分段点组合规模相当大。假设在M个测点中只存在一组基本线元，随着M的不断增大，穷举法所要遍历的分段点组合种类呈指数型增长，如图 2所示。

实际上轨道检测线路长、测点多，因此使用穷举法解决该问题计算量庞大，时效性差。本文提出先将轨道分为多个基本线元组，然后在之前分段结果的基础上进行组合迭代得到每一基本线元组的最佳轨道参数，最终合并为整条轨道的最佳线形参数。在识别分段点时以每一基本线元组内必定存在的4个分段点(直缓(ZH)点、缓圆(HY)点、圆缓(YH)点、缓直(HZ)点)为一组进行分段，然后每次优化包含该组分段点在内的一段数据{(m_i,α_i),…,(m_ZH,α_ZH),…,(m_HY,k_HY),…,(m_YH,k_YH),…,(m_HZ,k_HZ),…,(m_i+n,k_i+n)}在当前基本线元组最优的基础上再对下一个基本线元组进行优化，相邻两个基本线元组之间设有一定交集以保证每一个基本线元组的最优性，且避免穷举法指数增长的遍历次数。

3 轨道最佳线形参数确定过程

针对曲率-里程曲线，先基于最小二乘拟合斜率变化规律提取出概略分段点，将轨道分为多个基本线元组，然后根据轨道线形曲率变化规律拟合出基于该分段点的轨道线形，再进行组合迭代，最终取使目标函数最优的一组线形参数为轨道最佳线形参数。

3.1 基于轨道曲率曲线最小二乘拟合斜率变化的概略分段

量测数据得到的曲率-里程曲线中包含测量误差和轨道不平顺，通过摒弃粗差和中值滤波能够得到较为平顺的曲率-里程曲线。对图 1一例中曲率-里程曲线中的测点i前后两段数据X={(m_i－n1,k_i－n1),…,(m_i－1,k_i－1),(m_i,k_i)}和X={(m_i,k_i),(m_i+1,k_i+1),…,(m_i+n1,k_i+n1)}分别应用最小二乘直线拟合算法进行拟合，计算前一段拟合曲线与后一段拟合曲线的斜率差如图 3所示。

由图 3可知，轨道曲率曲线基于最小二乘拟合的斜率变化存在如下规律：其斜率差在各分段点处有极值，其中ZH点和HZ点处有极小值，HY点和YH点处有极大值。对于圆曲线弯曲方向相反的情况，拟合斜率差变化规律相反，其中ZH点和HZ点处有极大值，HY点和YH点处有极小值，如图 4所示。

因此，根据轨道线形曲率特征以及曲率曲线基于最小二乘拟合斜率变化规律可以对线形进行概略分段得到多组分段点坐标，将轨道分为多段基本线元组。

3.2 基于组合迭代的轨道线形优化算法

针对每一基本线元组的分段点，先根据轨道线形曲率特征拟合线形，然后组合迭代至该基本线元组最优，再处理下一基本线元组。

轨道线形拟合应如下进行：确定起点、一组分段点以及终点横坐标(M₀,M₁,M₂,M₃,M₄,M_e)后，结合图 1~2中轨道平面线形曲率的特征，应按照以下原则确定其曲率值：

$\begin{align} & k({{M}_{i}})= \\ & \left\{ \begin{matrix} 0\text{ }i=0,1,4,e \\ \omega =E(k({{M}_{2}}),k({{M}_{2}}+1),\cdots ,k({{M}_{3}}))\text{ }i=2,3 \\ \end{matrix} \right. \\ \end{align}$

其中:直线段曲率为0;圆曲线段曲率为1/R。

将坐标确定的各点逐次相连后得到拟合的曲率-里程曲线，沿里程积分后得到拟合的方向角-里程曲线，由此可以计算每一分段点组合下的目标函数值。

组合迭代算法是将概略分段点一定阈值范围内的分段点组合作为搜索域，遍历搜索域内所有组合确定轨道最优线形的算法，步骤如下：

步骤1 确定迭代进行的搜索域范围，代入概略分段点轨道线形为初始目标解(ds,X)，计算初始目标函数值f(ds,X)。

步骤2 在搜索域内改变分段点组合，重新确定轨道线形(ds_n,X_n)，计算目标函数值f(ds_n,X_n)，若f(ds_n,X_n) ＜f(ds,X)，则接受新解; 否则接受原解。

步骤3 判断是否遍历搜索域内的所有分段点组合，是则算法终止;否则重复步骤2。

完整的算法流程如图 5所示，通过本算法可以得到整条轨道的最佳线形参数。

4 实验算例

引入均方根误差(Root Mean Square Error,RMSE)评价轨道线形拟合优劣程度，RMSE越小表示拟合程度越好：

$RMSE=\sqrt{\left( \sum{{{\Delta }_{i}}^{2}} \right)/N}$

其中:Δ_i = α_i-α_i′表示测点方向角与拟合方向角的偏差；N表示参与拟合的测点个数。

4.1 仿真算例

设计仿真实验，以0.6 m采样间隔仿真一条包含两段基本线元组的轨道线形，长900 m，有测点1501个。取最小缓和曲线长度、最小圆曲线长度均为30 m，惯导测角误差取不超过±0.05°的随机值，轨道不平顺取不超过±10 mm的随机值，分别应用穷举法、文献^[12]算法和本文算法对两段基本线元组进行处理。

本文算法首先根据定长曲率曲线最小二乘拟合斜率变化对轨道进行概略分段，取定长为30 m，最小二乘拟合斜率变化曲线经滤波后如图 6所示。

图 6中曲线的极值点即为本算法得到的概略分段点。基于概略分段结果使用组合迭代算法进行精确分段并拟合最佳线形参数，得到的线形优化结果如图 7所示。

文献^[12]根据带噪声的方向角变量曲线(如图 8所示)人工判断出两组分段点，基于这两组分段点得到线形参数拟合结果，与穷举法及本算法的结果比较如图 9。对图 9中方向角拟合结果曲线中矩形框一段放大后如图 10所示。各算法线形分段及参数估计结果见表 1，线形参数估计评价指标列表见表 2。

结果分析：由于量测值中包含噪声和轨道不平顺，因此得到的是有噪声干扰的最佳曲线而非真实的设计曲线，穷举法拟合出的轨道线形参数能够使得目标函数最小，由表 1中的线形分段及参数估计结果可知，本文算法与穷举法结果更相近，文献^[12]分段点2结果次之，文献^[12]分段点1结果与穷举法偏差较大。

从穷举法、文献^[12]算法、本文算法线形参数估计指标列表表 2中可以看出，目标函数值有∑Δ_{i 穷举法}²＜∑Δ_{i 本文算法}²＜∑Δ_{i 文献^[5]分段点1}²＜∑Δ_{i 文献^[5]分段点2}²，RMSE指标有RMSE_穷举法 ＜RMSE_本文算法 ＜RMSE_{文献^[5]分段点1} ＜RMSE_{文献^[5]分段点2},因此在本算例中，穷举法结果最优，然后依次是本算法、文献^[12]分段点2、文献^[12]分段点1。图 9所示结果放大图也表明本算法与穷举法结果相近，优于文献^[12]方法的结果。

与文献^[12]相比，由于文献^[12]方法不存在优化过程，其结果优劣完全取决于人工判断分段点的准确性，分段点判断越精确线形参数拟合结果越好。但由于方向角变量曲线存在噪声(如图 8所示)，人工往往很难从其中准确分辨出分段点，因此线形参数拟合结果具有不确定性。由于第二组分段点精确性优于第一组，因此基于第二组分段点得到的线形参数拟合结果更优。本算法经迭代优化得到的分段点精确性一般高于人工判断，因此本算法较文献^[12]更优。

与穷举法相比，本算法的RMSE指标仅比穷举法大4.93%，说明本算法计算结果的准确性与可行性。另一方面，穷举法的遍历次数高于4×10⁷次，本算法在仿真中设置搜索域为±20个测点因此遍历分段点组合8000个，计算量仅为穷举法的0.02%，在计算量上本算法远远优于穷举法。

4.2 应用算例

应用本算法对西安地铁3号线右线K18+554~K21+000段轨道按轨枕间距(平均为0.6 m)采样，该段轨道共包含圆曲线5段、缓和曲线10段、直线段6段，测点4030个，依本文算法进行线形优化，得到最终线形优化结果如图 11，线形分段结果及线形参数见表 3,表中缓和曲线长度为该圆曲线前后两段缓和曲线长度的平均值。

图 11中优化曲线与量测曲线基本重合，说明两者贴合度好；表 3中列出各直线段及圆曲线段起、终点处的里程值、各圆曲线的缓和曲线长度及圆曲线半径。结合数据算得目标函数值为∑Δ_i ²= 21.3015，RMSE为0.073，反映了量测噪声以及轨道不平顺。该应用算例表明本算法可以很好地对实际轨道中的基本线元组进行分段，并能够有效地从量测数据中提取出轨道最佳线形参数。

5 结语

本文提出先基于曲率曲线最小二乘拟合斜率变化对轨道线形概略分段，再使用组合迭代算法精确分段并优化线形的轨道平面最佳线形估计算法，可以获得与实际量测曲线最为贴合的符合设计规范的轨道线形参数。应用本算法可以处理轨道既有线线形优化问题，得到调拨量最小前提下的轨道平面线形最佳参数，也为后续处理轨道不平顺垫定了基础。本算法由量测值拟合轨道曲率再积分成方向角，积分起点选择不佳会造成误差随积分累积而导致结果不理想，应进一步改进轨道线形拟合方法以减少噪声对轨道线形确定的影响。