-
摘要: 高精度的地面LiDAR点云配准是空间目标三维表面拓扑重建的关键，针对待配准LiDAR点云和基准LiDAR点云存在位置、姿态和比例缩放差异的问题，提出了基于直线簇的地面LiDAR点云配准方法。首先，根据直线间相交、平行和异面的拓扑关系，分别对待配准和基准LiDAR点云的直线进行聚簇，构建直线簇；然后，分别将同名直线用Plücker坐标表示，通过待配准LiDAR点云的直线簇在空间中的螺旋缩放运动，使其与基准LiDAR点云的直线簇比例尺一致，且同名Plücker直线重合，构建基于直线簇的共线条件方程，实现了比例因子和相对位姿一体化解算。实验结果表明，直线簇的螺旋缩放增强了配准方程的几何约束性，提高了抗噪声能力，实现了高精度的地面LiDAR点云配准。Abstract: The high precision terrestrial LiDAR point clouds registration is the key to ensure 3D topological reconstruction of spatial objects. The reference LiDAR point clouds and the LiDAR point clouds to be registered have the issue of different position, attitude and scale. A registration method based on line cluster is proposed. Firstly, line cluster of the reference LiDAR point clouds and the LiDAR point clouds to be registered is respectively established according to the topological relation of lines. Then conjugate lines of point clouds are represented in Plücke coordinate. Through the spiral zooming motion of line cluster of the LiDAR point clouds to be registered, finally the scale of the two LiDAR point clouds is the same and the conjugate Plücke lines coincide with each other. The collinearity equations based on line cluster can solve the scale factor and the relative position and attitude simultaneously with high precision. The experimental results show that, the spiral and scaling of the line cluster enhances the geometric constraints of the registration equation, improves the ability of anti-noise, and realizes the high precision of terrestrial LiDAR point clouds registration.
-
Keywords:
- LiDAR /
- point clouds registration /
- line cluster /
- Plücke coordinate /
- spiral motion
-
守时实验室需要建立和维持一个准确、稳定、可靠的时间尺度作为时间基准[1-3]，其核心算法包括时间尺度算法[4]、钟差预测算法[5-6]、驾驭算法[7-9]等。原子钟是守时系统的核心部件。时间基准的性能不仅和算法设计有关，还和参与计算的各台原子钟的性能有关。每种算法需要针对各参与计算的原子钟的模型和频率稳定度来设计。因此分析原子钟的模型和稳定度具有重要意义。
目前的研究中，文献[10-13]只分析了指定平滑时间的稳定度，或只对时差观测量中的某一具体分量展开分析，没有具体分析时差观测量中各分量对Allan偏差的贡献，并通过估计值外推得到所有平滑时间的Allan偏差估计值；文献[12-14]对Kalman滤波器估计原子钟状态原理描述不清晰；当观测噪声过大、存在周期性波动时，无法使用斜率法准确直接估计原子钟噪声强度[15]；当周期性波动在时差中不明显时，目前估计方法较为复杂[12, 16-17]。
本文针对上述问题展开研究，探索原子钟模型和频率稳定度分析方法。大量实验表明，典型氢钟和铯钟的观测模型可以表示为[2, 18-20]：①确定性部分，用二次多项式（时差、频差和线性频漂）加周期性波动项表示；②随机性部分，即原子钟噪声，主要为频率白噪声（white frequency modulation noise，WFM）和频率随机游走噪声（walk random frequency modulation noise，RWFM）；③观测噪声，噪声类型为相位白噪声（white phase modulation noise，WPM）。本文详细分析了该模型各分量的Allan偏差表达式。
在此基础上，本文从最优估计和低通滤波器[7-9]两个角度描述Kalman滤波器（Kalman filter，KF）估计原子钟状态的原理；提出了综合KF状态估计的结果和Allan偏差图估计原子钟噪声和观测噪声强度的方法；提出了3种不同的估计线性频漂幅度的方法；结合原子钟随机微分方程模型，提出了综合Kalman滤波器状态估计的结果和对数Allan偏差图估计周期性波动周期和幅度的方法。
1 原子钟模型及频率稳定度
1.1 原子钟时差和时差量模型
典型的氢钟和铷钟，其确定性分量用二次多项式来表示；典型的铯钟，其确定性分量用一次多项式来表示。二次多项式模型[2]表示为：
$$ x\left( t \right) = {x_0} + {y_0}t + \frac{1}{2}d{t^2} + {\varepsilon _x}\left( t \right) $$ (1) 式中，等式右边前3项对应确定度分量，其中x0为初始时差，y0为初始频差，d为线性频漂；εx(t)对应时差的随机性部分，即噪声。一次多项式模型即式（1）中线性频漂d=0的情况。
某些时候，时差的确定性分量中还包含周期性波动分量。为简化分析，假设瞬时频差的周期性波动为标准正弦波或余弦波形式[19]，即：
$$ {y_s}\left( t \right) = A\cos \left( {2{\rm{ \mathit{ π} }}{f_0}t + \varphi } \right) $$ (2) 式中，A为幅度；f0为周期性波动的频率；φ为初始相位。
这时，频差周期性波动在时差上表现为：
$$ \begin{array}{*{20}{c}} {{x_s}\left( t \right) = \int_0^t {{y_s}\left( s \right){\rm{d}}s} = A\int_0^t {\cos \left( {2{\rm{ \mathit{ π} }}{f_0}s + \varphi } \right){\rm{d}}s} = }\\ {\frac{A}{{2{\rm{ \mathit{ π} }}{f_0}}}\sin \left( {2{\rm{ \mathit{ π} }}{f_0}s + \varphi } \right)\left| {_0^t} \right.} \end{array} $$ (3) 于是，二次多项式叠加周期性波动项的时差模型表示为：
$$ \begin{array}{*{20}{c}} {x\left( t \right) = {x_0} + {y_0}t + \frac{1}{2}d{t^2} + }\\ {\frac{A}{{2{\rm{ \mathit{ π} }}{f_0}}}\sin \left( {2{\rm{ \mathit{ π} }}{f_0}s + \varphi } \right)\left| {_0^t} \right. + {\varepsilon _x}\left( t \right)} \end{array} $$ (4) 通过大量观测发现，氢钟、铯钟时差中起主导作用的噪声是RWFM[2, 18, 20-21]和WFM。这两种噪声在时差上可以分别用维纳过程和积分维纳过程来建模[2, 18, 20-21]，即：
$$ {\varepsilon _x}\left( t \right) = {\sigma _1}{W_1}\left( t \right) + {\sigma _2}\int_0^t {{W_2}\left( s \right){\rm{d}}s} $$ (5) 式中，εx(t)为时差的随机性部分，即原子钟时差的噪声；W1(t)和W2(t)分别代表两个独立的维纳过程，并且有W(t)~N(0, t)，即每个维纳过程服从均值为0、方差为时间t的正态分布；σ1和σ2分别是这两个维纳过程的扩散系数，用于表明噪声的强度。
把式（5）代入式（4），时差可表示为：
$$ \begin{array}{*{20}{c}} {x\left( t \right) = {x_0} + {y_0}t + \frac{1}{2}d{t^2} + \frac{A}{{2{\rm{ \mathit{ π} }}{f_0}}}\sin \left( {2{\rm{ \mathit{ π} }}{f_0}s + \varphi } \right)\left| {_0^t} \right. + }\\ {{\sigma _1}{W_1}\left( t \right) + {\sigma _2}\int_0^t {{W_2}\left( s \right){\rm{d}}s} } \end{array} $$ (6) 时差的观测量表示为：
$$ \begin{array}{*{20}{c}} {z\left( t \right) = x\left( t \right) + \sigma \varepsilon \left( t \right) = {x_0} + {y_0}t + }\\ {\frac{1}{2}d{t^2} + \frac{A}{{2{\rm{ \mathit{ π} }}{f_0}}}\sin \left( {2{\rm{ \mathit{ π} }}{f_0}s + \varphi } \right)\left| {_0^t} \right. + }\\ {{\sigma _1}{W_1}\left( t \right) + {\sigma _2}\int_0^t {{W_2}\left( s \right){\rm{d}}s} + \sigma \varepsilon \left( t \right)} \end{array} $$ (7) 式中，ε(t)为WPM；σ用于表明观测噪声的强度。
1.2 各分量Allan方差的表达式
时域上通常用Allan方差来表征频率稳定度。Allan方差的定义如下[2, 19]：
$$ \sigma _y^2\left( \tau \right) = \frac{1}{2}E\left[ {{{\left( {\bar y\left( {{t_k} + \tau } \right) - \bar y\left( {{t_k}} \right)} \right)}^2}} \right] $$ (8) 式中，τ为平滑时间；y为平均频差，定义为：
$$ \bar y\left( {{t_k}} \right) = \frac{1}{\tau }\int_{{t_k} - \tau }^{{t_k}} {y\left( t \right){\rm{d}}t} = \frac{{x\left( {{t_k}} \right) - x\left( {{t_k} - \tau } \right)}}{\tau } $$ (9) 式中，y为瞬时频差；x即为式（6）中定义的时差。
实际上常用Allan方差的平方根，即Allan偏差来表征稳定度。文献[22]详细推导了扩散系数与Allan方差的关系，即：
$$ \sigma _y^2\left( \tau \right) = \sigma _1^2/\tau + \frac{1}{3}\sigma _2^2\tau $$ (10) 式中，σy2(τ)代表平滑时间为τ时的Allan方差。式（10）中等式右边第一项斜率为-1，第二项斜率为1，分别对应WFM和RWFM。这说明在对数Allan方差图中，WFM的斜率为-1，RWFM的斜率为1。在实际应用中，很容易通过对Allan方差拟合得到扩散系数的值。
文献[20]推导了WPM与Allan方差的关系：
$$ \sigma _y^2\left( \tau \right) = 3{\sigma ^2}/{\tau ^2} $$ (11) 式（11）表明，WPM在对数Allan方差图中斜率为-2。
把式（1）和式（9）代入式（8），得到线性频漂d与Allan方差的关系：
$$ \sigma _y^2\left( \tau \right) = \frac{1}{2}{d^2}{\tau ^2} $$ (12) 式（12）表明，d在对数Allan方差图中斜率为2。式（2）所示的周期性波动项和Allan方差的关系为[19]：
$$ \sigma _y^2\left( \tau \right) = {A^2}\frac{{{{\sin }^4}\left( {{\rm{ \mathit{ π} }}{f_0}\tau } \right)}}{{{{\left( {{\rm{ \mathit{ π} }}{f_0}\tau } \right)}^2}}} $$ (13) 综上，式（6）所示的时差的Allan方差表示为：
$$ \sigma _y^2\left( \tau \right) = \sigma _1^2/\tau + \frac{1}{3}\sigma _2^2\tau + \frac{1}{2}{d^2}{\tau ^2} + {A^2}\frac{{{{\sin }^4}\left( {{\rm{ \mathit{ π} }}{f_0}\tau } \right)}}{{{{\left( {{\rm{ \mathit{ π} }}{f_0}\tau } \right)}^2}}} $$ (14) 式（7）所示的时差观测量的Allan方差表示为：
$$ \begin{array}{*{20}{c}} {\sigma _y^2\left( \tau \right) = 3{\sigma ^2}/{\tau ^2} + \sigma _1^2/\tau + \frac{1}{3}\sigma _2^2\tau + }\\ {\frac{1}{2}{d^2}{\tau ^2} + {A^2}\frac{{{{\sin }^4}\left( {{\rm{ \mathit{ π} }}{f_0}\tau } \right)}}{{{{\left( {{\rm{ \mathit{ π} }}{f_0}\tau } \right)}^2}}}} \end{array} $$ (15) 式中，等式右边第1项为观测噪声；第2、3项为原子钟噪声；第4项为频漂；第5项为周期性波动在Allan方差中的分量。
由式（15）看出，WPM、WFM、RWFM和线性频漂在对数Allan方差图中的斜率分别为-2、-1、1和2，在对数Allan偏差图中的斜率分别为-1、-1/2、1/2和1。假如时差中存在较大幅度的周期性波动，那么Allan方差也存在明显的周期性波动，即Allan方差在中间某一平滑时间段内会凸起来。
但是，并不是每台钟的模型中都含有上述所有分量。例如大部分铯钟几乎没有线性频漂，大部分型号的地面钟都没有周期性波动。
2 时差模型及频率稳定度分析方法
本节首先从最优估计和低通滤波器两个角度描述了KF估计原子钟状态的原理；然后展示分析WPM、WFM、RWFM强度、线性频漂幅度、周期性波动的幅度和频率的方法。
2.1 使用Kalman滤波器估计原子钟状态
式（1）所示的原子钟模型用随机微分方程（stochastic differential equations，SDEs）闭合解的离散形式表示为[2, 18, 20-21]：
$$ \left\{ \begin{array}{l} x\left( {{t_{k + 1}}} \right) = x\left( {{t_k}} \right) + {x_2}\left( {{t_k}} \right)T + \frac{1}{2}d{T^2} + \\ \;\;\;\;\;\;\;\;\;\;\;\;\;{\sigma _1}\left( {{W_1}\left( {{t_{k + 1}}} \right) - {W_1}\left( {{t_k}} \right)} \right) + \\ \;\;\;\;\;\;\;\;\;\;\;\;\;{\sigma _2}\int_{{t_k}}^{{t_{k + 1}}} {\left( {{W_2}\left( s \right) - {W_2}\left( {{t_k}} \right)} \right){\rm{d}}s} \\ {x_2}\left( {{t_{k + 1}}} \right) = {x_2}\left( {{t_k}} \right) + dT + {\sigma _2}\left( {{W_2}\left( {{t_{k + 1}}} \right) - {W_2}\left( {{t_k}} \right)} \right) \end{array} \right. $$ (16) 式中，x和x2分别代表了两个状态分量，x和式（1）完全相同，代表时差，x2代表频差的一部分，即不含WFM的频差；T为时间间隔，T=tk+1-tk；其他符号的含义和§1.1相同。其中，噪声分量为：
$$ {J_k} = \left[ {\begin{array}{*{20}{c}} {{\sigma _1}\left( {{W_1}\left( {{t_{k + 1}}} \right) - {W_1}\left( {{t_k}} \right)} \right) + {\sigma _2}\int_{{t_k}}^{{t_{k + 1}}} {\left( {{W_2}\left( s \right) - {W_2}\left( {{t_k}} \right)} \right){\rm{d}}s} }\\ {{\sigma _2}\left( {{W_2}\left( {{t_{k + 1}}} \right) - {W_2}\left( {{t_k}} \right)} \right)} \end{array}} \right] $$ (17) 使用KF对上述两个状态分量进行估计，得到估计值记为$\hat x\left({{t_k}} \right)$和${\hat x_2}\left({{t_k}} \right)$，其步骤用以下5个方程表示：
$$ {\mathit{\boldsymbol{\hat s}}_{k, k - 1}} = \phi {\mathit{\boldsymbol{\hat s}}_{k - 1, k - 1}} $$ (18) $$ {P_{k, k - 1}} = \phi {P_{k - 1, k - 1}}{\phi ^{\rm{T}}} + \mathit{\boldsymbol{Q}} $$ (19) $$ {\mathit{\boldsymbol{K}}_k} = {P_{k, k - 1}}{H^{\rm{T}}}{\left( {H{P_{k, k - 1}}{H^{\rm{T}}} + R} \right)^{ - 1}} $$ (20) $$ {\hat s_{k, k}} = {\hat s_{k, k - 1}} + {K_k}\left( {{z_k} - H{{\hat s}_{k, k - 1}}} \right) $$ (21) $$ {\mathit{\boldsymbol{P}}_{k, k}} = \left( {\mathit{\boldsymbol{I}} - {\mathit{\boldsymbol{K}}_k}\mathit{\boldsymbol{H}}} \right){\mathit{\boldsymbol{P}}_{k, k - 1}} $$ (22) 式中各符号意义见文献[12, 21]。其中，H=[1 0]；$\phi $=$\left[ {\begin{array}{*{20}{c}} 1&T\\ 0&1 \end{array}} \right]$；过程噪声方差矩阵Q和观测噪声方差矩阵R分别表示为：
$$ \begin{array}{*{20}{c}} {\mathit{\boldsymbol{Q}} = {\rm{E}}\left[ {\left( {{J_k} - 0} \right){{\left( {{J_k} - 0} \right)}^{\rm{T}}}} \right] = }\\ {\left[ {\begin{array}{*{20}{c}} {\sigma _1^2T + \frac{1}{3}\sigma _2^2{T^3}}&{\frac{1}{2}\sigma _2^2{T^2}}\\ {\frac{1}{3}\sigma _2^2{T^3}}&{\sigma _2^2T} \end{array}} \right]} \end{array} $$ (23) $$ \mathit{\boldsymbol{R}} = {\sigma ^2} $$ (24) 当原子钟时差符合模型（16），Q和R的值符合${\sigma ^2}$、${\sigma _1}^2$和${\sigma _2}^2$时，KF估计得到的$\hat x\left({{t_k}} \right)$和${\hat x_2}\left({{t_k}} \right)$是最小均方意义下的最优估计。从频域角度，确定Q和R的值相当于确定KF的带宽[8-9]；KF作为低通滤波器，其作用就是滤除时差观测量中高频的WPM，$\hat x\left({{t_k}} \right)$中保留了低频的WFM和RWFM，${\hat x_2}\left({{t_k}} \right)$中只保留了最低频的RWFM[8-9]。
文献[12-14]将x2理解为频差，但笔者认为，频差是由时差的差分得到的，其中依然包含了WFM和RWFM，而x2中只包含RWFM。
由于$\hat x\left({{t_k}} \right)$中包含WFM和RWFM，而${\hat x_2}\left({{t_k}} \right)$中只包含RWFM，噪声分量更少，所以，从${\hat x_2}\left({{t_k}} \right)$中可以更有效地估计频漂、周期性分量等。当含有周期性波动项时，尽管原子钟模型不完全符合式（16），但从KF作为低通滤波器的角度[8-9]出发，KF依然可以有效滤除观测噪声，估计出$\hat x\left({{t_k}} \right)$和${\hat x_2}\left({{t_k}} \right)$。
2.2 使用斜率法和Kalman滤波器分析WPM、WFM和RWFM的强度
WPM、WFM和RWFM在对数Allan偏差图中的斜率分别为-1、-1/2、1/2，理论上可以通过斜率拟合出σ2、σ12和σ22的估计值，记为${\hat \sigma ^2}$、${\hat \sigma _1}^2$和${\hat \sigma _2}^2$。但是，实际上很多情况，尤其是远程比对的情况下，观测噪声很大，造成对数Allan偏差图中WFM分量淹没在WPM分量中。此外，周期性波动幅度很大时，Allan偏差图中间部分会向上凸起。这些都可能造成无法观察到WFM分量。针对这种情况，这里提出一种在${\hat \sigma _1}^2$无法准确获知的情况下，使用KF结合时差观测量z(t)和时差估计值$\hat x\left(t \right)$的Allan偏差，获取${\hat \sigma ^2}$、${\hat \sigma _1}^2$和${\hat \sigma _2}^2$的方法。步骤如下：
1) 计算式（7）所示z(t)的Allan偏差，斜率拟合${\hat \sigma _2}^2$和${\hat \sigma ^2}$；
2) 设置一个σ12的大致估计值，代入式（23）确定Q的值，使用KF估计得到Q的估计值，并画出$\hat x\left(t \right)$的Allan偏差图；理论上讲，$\hat x\left(t \right)$的噪声分量应该只包含WFM和RWFM；
3) 根据$\hat x\left(t \right)$的Allan偏差图，反复调整σ12的估计值，重复步骤2），观察$\hat x\left(t \right)$的Allan偏差图，确保噪声分量斜率为-1/2和1/2，即只含有WFM和RWFM，最终获取${\hat \sigma _1}^2$。
2.3 分析线性频漂的幅度
多种方法可以估计线性频漂$\hat d$。如直接对时差观测量z(t)进行最小二乘拟合（方法1），或对Allan偏差的斜率进行拟合（方法2），或对KF估计得到的${\hat x_2}\left(t \right)$的斜率进行最小二乘拟合（方法3）。
参照文献[5-6]，可以从理论上分析这些方法的估计不确定度，本文不展开分析。按照§2.1的分析，由于KF得到的${\hat x_2}\left(t \right)$中滤除了WFM，只含有RWFM，此外典型氢钟或铯钟WFM强度远大于RWFM强度，即σ12远大于σ22，所以从直观上理解，方法3的估计不确定度最小。
2.4 分析周期性波动分量
1）分析f0。按照§2.1方法，尽管式（16）没有对周期性波动建模，KF依然可以估计得到$\hat x\left(t \right)$和${\hat x_2}\left(t \right)$，由于${\hat x_2}\left(t \right)$中只含有RWFM，周期性波动在z(t)和$\hat x\left(t \right)$中不明显，在${\hat x_2}\left(t \right)$中很明显。可以通过观测较长时间段内${\hat x_2}\left(t \right)$波峰波谷的位置，得到波动频率的估计值，记为${\hat f_0}$。
2）分析A。由于噪声的存在，观察${\hat x_2}\left(t \right)$无法准确获取A的估计值，记为$\hat A$。本文采用如下方法：在$\hat x\left(t \right)$的Allan偏差图中，反复调整$\hat A$的值，观察式（13）平方根所示的周期性波动的Allan偏差与$\hat x\left(t \right)$的Allan偏差的吻合程度，最终得到$\hat A$。
3 实验分析
本文采用两台国产氢钟（分别记为Hm1和Hm2）相对于参考时间基准（记为Ref）的实测钟差（Hm1-Ref，Hm2-Ref）进行分析，以验证提出的原子钟模型和频率稳定度分析方法。其中Ref和国际协调世界时同步，稳定度远高于单台国产氢钟。
3.1 使用斜率法和KF分析WPM、WFM和RWFM的强度和线性频漂的幅度
取某一段长度约75 d的Hm1-Ref实测数据，作为时差观测量z(t)。按照§1.2的方法估计${\hat \sigma ^2}$、${\hat \sigma _1}^2$和${\hat \sigma _2}^2$，得到${\hat \sigma ^2}$=1×10-22 s2，${\hat \sigma _1}^2$=3×10-26 s和${\hat \sigma _2}^2$=1.2×10-33 s-1。其中，${\hat \sigma ^2}$=1×10-22 s2意味着观测噪声的标准差为0.01 ns或10 ps，符合测量设备的指标。
图 1展示了KF前后z(t)和$\hat x\left(t \right)$的Allan偏差，同时也画出了根据上述估计值，按照式（10）、（11）的平方根计算得到的WPM、WFM、RWFM部分的Allan偏差。综上，本实验表明§2的方法可以在观测噪声和周期性波动幅度很大时，有效估计出${\hat \sigma ^2}$、${\hat \sigma _1}^2$和${\hat \sigma _2}^2$。
![]()
图 1 z(t)和$\hat x\left(t \right)$的Allan偏差Figure 1. Allan Deviations of z(t) and $\hat x\left(t \right)$3.2 采用多种方法分析和验证线性频漂的幅度
本节分别对§2.3的3种方法进行验证。
方法1：对$\hat x\left(t \right)$分别用一阶多项式和二阶多项式进行最小二乘拟合，得到线性频漂的估计值为$\hat d$=-3.891×10-20 s/s2。图 2(a)和2(b)分别画出了拟合残差。对比图 2(a)和2(b)，可以明显看出Hm1中存在线性频漂。
![]()
图 2 $\hat x\left(t \right)$的一阶残差和二阶残差Figure 2. First Order and Second Order Residuals of $\hat x\left(t \right)$方法2：验证方法1。将方法1得到的$\hat d$值代入式（12）再平方根，计算得到Allan偏差的频漂部分。图 1中画出了Allan偏差的频漂部分，可知两种方法得到的$\hat d$值基本吻合。
方法3：图 3画出了KF估计得到的${\hat x_2}\left(t \right)$，可知${\hat x_2}\left(t \right)$中存在明显的线性频漂。对${\hat x_2}\left(t \right)$的斜率进行最小二乘拟合，得到$\hat d$=-3.86×10-20 s/s2，和方法1的结果基本吻合。
综上，§2.3的3种方法估计结果基本吻合。实际上，直观上分析，方法3估计不确定度较小，原因在于KF估计得到的${\hat x_2}\left(t \right)$中滤除了WFM，只含有RWFM。假如直接对时差一次差分或者二次差分，得到的序列中同样含有WFM和RWFM差分得到的噪声，由于WFM强度远大于RWFM，将导致无法直接在这些差分序列中观测到图 3中所示的周期性波动，周期性波动已经淹没于噪声中。
3.3 分析周期性波动项
采用§2.4的方法，观察图 3波峰波谷的位置，得到${\hat f_0}$=1/86 400 Hz，$\hat A$大约在1×10-14~2×10-14 s/s左右。反复比较式（13）开平方后所示的周期性波动和$\hat x\left(t \right)$的Allan偏差，最终得到$\hat A$=1.6×10-14 s/s。图 4画出了周期性波动项和$\hat x\left(t \right)$的Allan偏差，可见两者比较吻合。本实验表明了§2的方法可以有效估计周期性波动的周期和幅度。
![]()
图 4 周期性波动项和$\hat x\left(t \right)$的Allan偏差Figure 4. Allan Deviations of the Sinusoidal Component and $\hat x\left(t \right)$3.4 实例2分析
采用和§2相同的方法分析第2台国产氢钟Hm2的性能，得到该氢钟的参数估计值分别为：${\hat \sigma ^2}$=1×10-22 s2，${\hat \sigma _1}^2$=3×10-26 s，${\hat \sigma _2}^2$=8×10-33 s-1，$\hat d$=-3.8×10-20 s/s2，${\hat f_0}$=1/86 400 Hz和$\hat A$=1.0×10-14 s/s。图 5画出了KF前后z(t)和$\hat x\left(t \right)$的Allan偏差，以及通过参数估计值计算得到的WPM、WFM、RWFM、线性频漂、周期性波动项的Allan偏差。
![]()
图 5 z(t)和$\hat x\left(t \right)$的Allan偏差Figure 5. Allan Deviations of z(t) and $\hat x\left(t \right)$把§3.1~§3.4得到的参数估计值代入式（14）和式（15），得出τ < 1 d的实验结果和国产氢钟的说明书相符；当τ > 10 000 s时，观测噪声对Allan偏差的影响很小。
综上，本节实验结果验证了本文方法可以有效分析WPM、WFM、RWFM、线性频漂、周期性波动项各自的Allan偏差，以及总的Allan偏差，并通过这些估计值，拟合出任意平滑时间的Allan偏差估计值。
4 结语
本文展示了原子钟模型和频率稳定度分析方法，详细分析了原子钟时差观测量中的各分量，包括确定性部分（时差、频差、线性频漂和周期性波动项）、随机性部分（WFM、RWFM）和观测噪声（WPM）；分析了WPM、WFM、RWFM、线性频漂、周期性波动项在Allan偏差中的表达式，描述了KF用于估计原子钟状态的原理；提出了当在对数Allan偏差图中，WFM完全淹没于WPM时，使用KF估计WPM、WFM、RWFM强度的方法；提出了3种估计线性频漂幅度的方法和估计周期性波动周期和幅度的方法。通过两台国产氢钟的实测数据验证了本文方法的实用性。实际上，可以通过这些估计值拟合出任意平滑时间的Allan偏差估计值。本文提出的方法物理原理清晰，操作简便易行。该研究对于时间尺度、钟差预测、原子钟驾驭等算法具有重要意义。
-
![]()
图 4 噪声标准差对平移与缩放系数的影响
Figure 4. Influence of Noise Standard Deviation on the Pan and Scale Coefficient
![]()
图 6 噪声标准差对检查直线精度影响
Figure 6. Influence of Noise Standard Deviation on Accuracy of Check Line
表 1 不同方案模拟LiDAR点云配准结果
Table 1 Registration Results of Simulated LiDAR Point Clouds for Different Schemes
实验方案 X/m Y/m Z/m φ/(°) ω/(°) κ/(°) λ l4与l1 迭代不收敛 l4与l3 -0.9 -0.1 1 0 0 30 2 l4与l2 1 1 1 0 0 30 2 
下载: 导出CSV
表 2 实际LiDAR点云配准结果和精度
Table 2 The Registration Results and Accuracy of Actual LiDAR Point Clouds
方法 四元数法 Plücke直线 直线簇 直线簇(抽稀) XS/m -22.999 -22.993 -23.176 -23.047 YS/m 29.292 29.397 29.322 29.355 ZS/m -2.230 -2.314 -2.182 -2.298 φ/(°) 12.601 12.551 12.624 12.588 ω/(°) 1.0725 1.081 1.054 1.076 κ/(°) 28.900 28.902 28.900 28.862 λ 0.998 1 0.995 0.999 d/mm 8.7 5.7 3.0 3.2 θ/ (°) 0.030 0.024 0.024 0.026
下载: 导出CSV
-
[1] 杨荣华, 吕美英, 花向红.一种多站标靶点云整体配准算法[J].武汉大学学报·信息科学版, 2014, 39(10):1189-1193 http://ch.whu.edu.cn/CN/abstract/abstract3093.shtml Yang Ronghua, Lv Meiying, Hua Xianghong.An Algorithm for the Multiview Target Point Cloud Global Registration[J]. Geomatics and Information Science of Wuhan University, 2014, 39(10):1189-1193 http://ch.whu.edu.cn/CN/abstract/abstract3093.shtml
[2] 王永波, 杨化超, 刘燕华, 等.线状特征约束下基于四元数描述的LiDAR点云配准方法[J].武汉大学学报·信息科学版, 2013, 38(9):1057-1063 http://ch.whu.edu.cn/CN/abstract/abstract2751.shtml Wang Yongbo, Yang Huachao, Liu Yanhua, et al. Linear -Feature-Constrained Registration of LiDAR Point Cloud via Quaternion[J]. Geomatics and Information Science of Wuhan University, 2013, 38(9):1057-1063 http://ch.whu.edu.cn/CN/abstract/abstract2751.shtml
[3] 姚吉利, 贾象阳, 马宁, 等.地面激光扫描多站点云整体定向平差模型[J].测绘学报, 2014, 43(8):835-841 http://kns.cnki.net/KCMS/detail/detail.aspx?filename=chxb201408010&dbname=CJFD&dbcode=CJFQ Yao Jili, Jia Xiangyang, Ma Ning, et al. Overall Orientation Adjustment Model of Terrestrial Laser Scanning Multi-station Point Clouds[J].Acta Geodaetica et Cartographica Sinica, 2014, 43(8):835-841 http://kns.cnki.net/KCMS/detail/detail.aspx?filename=chxb201408010&dbname=CJFD&dbcode=CJFQ
[4] 王晏民, 石宏斌.应用三维线特征自动配准建筑物点云数据[J].测绘通报, 2014(9):23-25 http://kns.cnki.net/KCMS/detail/detail.aspx?filename=chtb201409007&dbname=CJFD&dbcode=CJFQ Wang Yanmin, Shi Hongbin. An Automatic Registration Method for Building Point Clouds Based on 3D Line Features[J]. Bulletin of Surveying and Mapping, 2014(9):23-25 http://kns.cnki.net/KCMS/detail/detail.aspx?filename=chtb201409007&dbname=CJFD&dbcode=CJFQ
[5] Jen J J, Tzu Y C. Registration of Ground-Based LiDAR Point Clouds by Means of 3D Line Features[J].Journal of the Chinese Institute of Engineers, 2008, 31(6):1031-1045 doi: 10.1080/02533839.2008.9671456
[6] Wang Yongbo, Yang Huachao, Liu Yanhua, et al. A Dual Quaternion-Based, Closed-Form Pairwise Registration Algorithm for Point Clouds[J]. ISPRS Journal of Photogrammetry and Remote Sensing, 2014, 94(4):63-69 https://www.sciencedirect.com/science/article/pii/S0924271614001087
[7] 寇媛, 徐景中.基于结构特征的遥感影像匹配[J].中国图象图形学报, 2013, 18(5):565-573 doi: 10.11834/jig.20130511 Kou Yuan, Xu Jingzhong. Algorithm forRemote Sensing Images Matching Based on the Structure Characteristics[J].Journal of Image and Graphics, 2013, 18(5):565-573 doi: 10.11834/jig.20130511
[8] Guan Y, Zhang H. Initial Registration for Point Clouds Based on Linear Features[C]. The 4th International Symposium on Knowledge Acquisition and Modeling, Sanya, China, 2011
[9] 徐景中, 寇媛, 袁芳, 等.基于结构特征的机载LiDAR数据与航空影像自动配准[J].红外与激光工程, 2013(12):42-52 http://kns.cnki.net/KCMS/detail/detail.aspx?filename=hwyj201312061&dbname=CJFD&dbcode=CJFQ Xu Jingzhong, Kou Yuan, Yuan Fang, et al. Auto-Registration of Aerial Imagery and Airborne LiDAR Data Based on Structure Feature[J]. Infrared and Laser Engineering, 2013(12):42-52 http://kns.cnki.net/KCMS/detail/detail.aspx?filename=hwyj201312061&dbname=CJFD&dbcode=CJFQ
[10] 郑德华, 岳东杰, 岳建平.基于几何特征约束的建筑物点云配准算法[J].测绘学报, 2008, 37(4):464-468 http://kns.cnki.net/KCMS/detail/detail.aspx?filename=chxb200804013&dbname=CJFD&dbcode=CJFQ Zhen Dehua, Yue Dongjie, Yue Jianping.Geometric Feature Constraint Based Algorithm for Building Scanning Point Cloud Registration[J]. Acta Geodaetica et Cartographica Sinica, 2008, 37(4):464-468 http://kns.cnki.net/KCMS/detail/detail.aspx?filename=chxb200804013&dbname=CJFD&dbcode=CJFQ
[11] Funda J, Taylor R H, Paul R P. On Homogeneous Transforms, Quaternions, and Computational Efficiency[J]. IEEE Transactions on Robotics and Automation, 1990, 1(6):382-388 http://citeseerx.ist.psu.edu/showciting?cid=972985
[12] Ronda J I, Gallego G, Valdes A. Camera Autocalibration Using Plücke Coordinates[C]. IEEE International Conference on Image Processing, Genova, 2005
[13] Sariyildiz E, Temeltas H. Solution of Inverse Kinematic Problem for Serial Robot Using Dual Quaterninons and PlückeCoordinates[C]. ASME International Conference on Advanced Intelligent Mechatronics, Singapore, 2009
[14] 盛庆红, 陈姝文, 费利佳, 等.基于Plücke直线的机载LiDAR点云与航空影像的配准[J].测绘学报, 2015, 44(7):761-767 doi: 10.11947/j.AGCS.2015.20140123 Sheng Qinghong, Chen Shuwen, Fei Lijia, et al. Registration of Aerial Image with Airborne LiDAR Data Based on Plücke Line[J]. Acta Geodaetica et Cartographica Sinica, 2015, 44(7):761-767 doi: 10.11947/j.AGCS.2015.20140123
[15] Lasenby J, Fitzgerald W J, Lasenby A N, et al. New Geometric Methods for Computer Vision:An Application to Structure and Motion Estimation[J]. International Journal of Computer Vision, 1998, 26(3):191-213 doi: 10.1023/A:1007901028047
-
期刊类型引用(21)
1. 李国俊，史丰丰，徐金锋，赵润，付桂涛，杨玉婷，李兆南，孙文龙，程梦飞，林勇昕. 基于钙离子光钟的160 d连续守时性能分析. 光学学报. 2025(04): 126-133 . 
百度学术
2. 崔海波，王伟，武建锋，李帅辰，李冲，王康. 机动状态下国产铯原子钟稳定性分析. 时间频率学报. 2024(01): 26-33 . 百度学术
3. 宋鑫，张秀荣，李卫东. 关联噪声对维纳过程和拉比振荡的影响. 山西大学学报(自然科学版). 2024(03): 578-582 . 百度学术
4. 李方能，梁益丰，许江宁，吴苗，朱兵. 基于鲁棒双参数指数平滑法的BDS卫星钟差预报. 中国惯性技术学报. 2024(07): 645-653 . 百度学术
5. 王芳敏，李汶林，戴鸿飞，李春怡，周建华，薛申辉，王波. A real-time performance improvement method for composite time scale. Chinese Physics B. 2024(09): 354-361 . 百度学术
6. 李宗源，袁海波，张虹，张继海，刘素芳，王一桁. 基于低通滤波器的原子钟时差测量数据降噪处理及其性能分析. 时间频率学报. 2024(03): 171-179 . 百度学术
7. 李昂，周铁中，薛潇博，易航，陈德好. 一种智能整定氢原子频标PID控制系统参数的方法. 宇航计测技术. 2024(05): 39-44 . 百度学术
8. 曹建峰，孔静，满海钧，鞠冰，张宇，刘荟萃. 天问一号星载USO的长期漂移标定. 武汉大学学报(信息科学版). 2024(12): 2181-2186 . 百度学术
9. 赵军华，张玉军，刘安斐，展昕，陈亮. 附加钟特性的原子钟预测不确定度分析. 测绘科学. 2024(11): 30-39 . 百度学术
10. 阚昊宇，胡志刚，吕逸飞，谢新，周仁宇，赵齐乐. 利用不同时间同步体制钟差评估北斗三号星载原子钟性能. 武汉大学学报(信息科学版). 2023(04): 604-610 . 百度学术
11. 李方能，梁益丰，许江宁，吴苗. BDS/GPS新型铷原子钟长期特性分析. 中国惯性技术学报. 2023(05): 452-461 . 百度学术
12. 李帅辰，武建锋，崔海波，方婧. 国产氢原子钟移动守时性能测试与分析. 导航定位学报. 2023(03): 38-44 . 百度学术
13. 郭文飞，朱萌萌，辜声峰，左鸿铭，陈金鑫. GNSS精密时频接收机时钟调控模型与参数设计方法. 武汉大学学报(信息科学版). 2023(07): 1126-1133 . 百度学术
14. 安景浩，王怀智，梁益丰，许江宁. BDS快速钟差预报方法对比分析. 舰船电子工程. 2023(03): 55-60 . 百度学术
15. 卢鋆，武建峰，袁海波，申建华，孟轶男，宿晨庚，陈颖. 北斗三号系统时频体系设计与实现. 武汉大学学报(信息科学版). 2023(08): 1340-1348 . 百度学术
16. 秦晓伟，杜二旺，赵春勃，王国永，杨涛. 基于频率源稳定度的钟差序列生成方法研究. 空间电子技术. 2023(04): 66-71 . 百度学术
17. 陈汗龙，董哲，周真帆，郑晓雪. 国产原子钟频率稳定度评估分析. 导航定位学报. 2023(06): 28-33 . 百度学术
18. 梁益丰，许江宁，吴苗，何泓洋. 北斗卫星钟差的CEEMDAN分解与周期项提取方法. 中国惯性技术学报. 2022(04): 476-484 . 百度学术
19. 张杰梁，陈鲤文，陈静，张煌辉，姜立斌. 守时系统钟组规模技术验证. 电子质量. 2022(10): 151-157 . 百度学术
20. 伍贻威，王世超. 不通过自由纸面时建立时间基准的方法与性能分析. 测绘学报. 2021(03): 343-355 . 百度学术
21. 惠恬，赵高长，苏军，龚莹莹. 原子钟数据的改进经验模态分解降噪. 河南科技大学学报(自然科学版). 2021(05): 45-50+56+7 . 百度学术
其他类型引用(11)
计量
- 文章访问数: 2522
- HTML全文浏览量: 284
- PDF下载量: 431
- 被引次数: 32
