跑道表面存在水、冰或雪覆盖会严重降低道面摩擦系数，显著影响飞机刹车性能^[1]. 根据ICAO 于2019年发布的《Circular 355》（https://news.ncac.mn/uploads/bookSubject/2022-05/6279e8f7346ae.pdf），当水或干雪厚度超过3 mm时，飞机的制动和方向控制能力将急剧下降^[2]. 因此，机场须对道面覆盖物进行毫米级精度的深度检测，这对保障飞行安全至关重要^[3].

超声波技术在道面状态检测方面已有广泛应用. Gudra等^[4]使用40 kHz超声波研究冰雪参数，通过反射系数获取深度信息. Bystrov等^[5]基于超声回波信号分析不同道面状态（干、湿、冰、雪）的特性. 然而实际应用中，环境因素会显著影响测量精度. 例如Campbell的SR50A传感器在强风条件下误差可达3 cm^[6]，Herman等^[7]的研究也表明渡越时间（time-of-flight，ToF）精度会受到环境因素如温度、湿度的影响.

ToF估计的准确性直接影响深度测量结果^[8]. Liu等^[9]提出的双阈值法虽然改进了起振点检测，但在低信噪比条件下仍存在局限. Sunol等^[10]提出的互相关法通过动态参考信号提高了精度，但未充分考虑信噪比影响. Xiong等^[11]采用的参数估计法在SNR<10 dB时表现良好，但仅适用于特定信号模型. 此外，Zhang等^[12]开发了基于小波变换的阈值检测方法，但弱信号下的ToF分辨能力依然不足.

对于信号滤波技术，卡尔曼滤波因其最优递推估计能力展现出独特优势. Akhlaghi等^[13]提出的自适应方法改进了过程噪声协方差矩阵(${\boldsymbol{Q}}$)和观测噪声协方差矩阵(R)估计, 然而${\boldsymbol{Q}}$和${\boldsymbol{R}}$的比值选择不当依然可能导致估计误差增大. Wang等^[14]提出ARMA模型方法自适应调整增益矩阵，减少了相位滞后的问题，但该方法的滤波方法会使信号失去一定带宽，导致信号失真. Wu等^[15]采用的变分贝叶斯方法通过噪声协方差矩阵估计增强了滤波稳定性，然而不符合指数衰减的信号模型滤波效果不足. Sui等^[16]提出自适应移动平均技术的混合卡尔曼滤波器，对于阶跃式信号有显著的信号提取能力，对周期性信号平滑程度过高导致信号失真. Panda等^[17]提出的deep Koopman Kalman filter采用三阶段训练策略解决训练速度和梯度问题，然而该方法仍存在相位滞后且对于不规则信号滤波能力不足的问题.

综上，卡尔曼滤波方法依然存在以下问题: 1) 在低信噪比条件下，滤波对信号前沿的噪声抑制较差; 2)滤波后信号平滑度较差，影响过零点位置估计; 3)滤波后信号失真严重，影响特征点判断; 4)可能存在明显的相位滞后问题. 这些问题将影响起振点位置和ToF估计准确性，进而降低覆盖物深度计算精度. 因此，本研究提出结合总变差和贝叶斯估计的自适应卡尔曼滤波方法(TV-BAKF)，动态调整$Q$和$R$，提升滤波效果. 通过模拟混合指数模型信号并添加不同功率噪声，分析信号前沿波形提取准确性以及卡尔曼滤波参数的动态变化，并在低温环境舱中测试5种覆盖物 (水、冰、雪浆、湿雪、干雪) 的超声回波信号，并通过对比不同滤波方法，验证本研究提出的TV-BAKF方法在低信噪比条件下的ToF估计准确性.

在超声波检测冰雪水深度时，外界环境因素（包括风、温度波动和机械振动）以及覆盖物自身的介质特性都会影响回波信号的信噪比. 其中，风和温度不仅会改变波速，还会导致信号振幅变化和波形畸变；机械振动则会产生高频噪声，干扰原始信号并降低信噪比. 此外，覆盖物本身若具有松散多孔或表面粗糙的特性，也会减弱回波信号的能量，进一步降低信噪比，这种低信噪比环境会对信号起振点的准确识别和ToF的精确估计带来挑战. 现有方法（如小波去噪）虽能提升整体信噪比和平滑波形，但对起振点附近的噪声抑制效果有限，难以准确提取信号前沿波形特征，导致起振点估计误差较大，ToF精度较低.

卡尔曼滤波通过动态融合预测和测量信息抑制噪声，其性能由过程噪声${{Q}}$和观测噪声${{R}}$决定. 如图1所示为不同${{Q}}$/${{R}}$组合的滤波结果影响. 其中，${A_0}$、${A_{\max }}$、$\Delta {\varphi _i}$和SNR分别表示信号的幅值、最大幅值、相位滞后、信噪比. 以图1信号为例，当${{Q}}$较大${{R}}$较小时，滤波器侧重观测数据，能保留波形细节但去噪效果不足；反之${{Q}}$较小R较大时，信噪比提升，但会导致显著相位滞后（如${{Q}}$=0.05、${{R}}$=2时滞后5个采样点）和64%的幅值衰减. 图1中的热力图揭示了信噪比、相位滞后和幅值间的互斥关系，可以看出，相位滞后影响ToF计算精度，而过度平滑会削弱信号前沿特征.

卡尔曼滤波通过动态融合状态预测值与实测值，实现系统状态的最优估计. 其核心表达式如下：

(1)$ \boldsymbol{x}_{k-1}=\boldsymbol{A} \hat{\boldsymbol{x}}_{k-1}^-+\boldsymbol{B} \hat{\boldsymbol{u}}_{k-1}^- \;,$

(2)$ \boldsymbol{P}_k^-=\boldsymbol{A} \boldsymbol{P}_{k-1} \boldsymbol{A}^{\mathrm{T}}+\boldsymbol{Q} \;.$

式中：${\boldsymbol{A}}$为状态矩阵，${\boldsymbol{B}}$为控制矩阵，$ {\boldsymbol{Q}} $为过程噪声协方差，$ {{\boldsymbol{x}}_{k - 1}} $为$ k - 1 $时刻真实值，$ \hat {\boldsymbol{x}}_{k - 1}^ - $为$ k - 1 $时刻的先验状态估计值，$ \hat {\boldsymbol{u}}_{k - 1}^ - $为$ k - 1 $时刻控制输入量，$ {\boldsymbol{P}}_k^ - $为$ k $时刻先验误差协方差，$ {{\boldsymbol{P}}_{k - 1}} $为$ k - 1 $时刻后验误差协方差. 则卡尔曼滤波状态更新方程为

(3)$ \boldsymbol{K}_k=\frac{\boldsymbol{P}_k^- \boldsymbol{H}_k^{\mathrm{T}}}{\boldsymbol{H}_k \boldsymbol{P}_k^- \boldsymbol{H}_k^{\mathrm{T}}+\boldsymbol{R}}\;, $

(4)$ \hat{\boldsymbol{x}}_k=\hat{\boldsymbol{x}}_k^-+\boldsymbol{K}_k\left(\boldsymbol{y}_k-\boldsymbol{H}_k \hat{\boldsymbol{x}}_k^-\right)\;, $

(5)$ \boldsymbol{P}_k=\left(\boldsymbol{I}-\boldsymbol{K}_k \boldsymbol{H}_k\right) \boldsymbol{P}_k^-\;. $

式中：$ {{\boldsymbol{K}}_k} $为$ k $时刻卡尔曼增益系数，${\boldsymbol{H}}_k$为增益矩阵，$ {\boldsymbol{R}} $为观测噪声协方差，$ {\hat {\boldsymbol{x}}_k} $为$ k $时刻后验状态估计值，$ \hat {\boldsymbol{x}}_k^ - $为$ k $时刻先验状态估计值，$ {{\boldsymbol{y}}_k} $为$ k $时刻的观测值，${{\boldsymbol{P}}_k}$、$ {\boldsymbol{P}}_k^ - $分别为后验和先验估计方差，${\boldsymbol{ I }}$为单位矩阵.

总变差法（total variation denoising，TV）常用于去除信号中的噪声并平滑信号. 设有一维连续函数$f(x)$，其在区间$[a,b]$上的全变分可定义为$ x $在$ [a, b] $上的弧长，表达式如下：

(6)$ V_a^b\left( f \right) = \int_a^b {\left| {f'\left( x \right)} \right|{\mathrm{d}}x}\;. $

对于一维离散信号的全变分，设某信号序列为$\left\{ {{y_i}} \right\},i = 1,2,3, \cdots ,n$，其全变分的定义为

(7)$ V\left( y \right) = \sum\limits_{i = 1}^n {\left| {{y_{i+1}} - {y_i}} \right|}\;. $

在一维信号降噪中，总变差法的目标是使得去噪后的信号${\gamma _i}$尽可能与观测信号$ {y_i} $差异最小化，并保持${\gamma _i}$的平滑度. $ \varGamma $表示最终通过总表差法处理获得信号的期望状态.

(8)$ \varGamma = \mathop {\min }\limits_y\; E\left( {y,\gamma } \right)+\lambda E\left( {y,\gamma } \right)\;. $

式中：$ E\left( {y,\gamma } \right) $表征观测信号与处理后信号的相似程度（数据误差项）；参数$\lambda \in \left[ {0,1.0} \right]$，用于调节$ E\left( {y,\gamma } \right) $和$ V\left( \gamma \right) $的权重.

(9)$ E\left( {y,\gamma } \right) = \dfrac{1}{2}\sum\limits_i^{} {{{\left( {{y_i} - {\gamma _i}} \right)}^2}}\;. $

对于一维信号的去噪，本研究采用改进的自适应卡尔曼滤波. 由于是对一维信号处理，因此卡尔曼滤波参数，如${\boldsymbol{Q}}$、${\boldsymbol{R}}$、${\boldsymbol{P}}$、${\boldsymbol{Q}}$，将变为1×1的标量形式. 自适应卡尔曼滤波器通过融合残差与新息矩阵，适应信号的动态变化. 其$ k $时刻先验残差$ \varepsilon _k^ - $表达式如下：

(10)$ \varepsilon _k^ - = \hat x_k^ - - {y_k}\;. $

式中：$ {y_k} $为$ k $时刻信号的观测值. $ k $时刻先验新息 $ \hat d_k^ - $表达式为

(11)$ \hat d_k^ - = {y_k} - {H_k}\hat x_k^ - = {H_k}{x_k}+{v_k} - {H_k}\left( {{x_k}+\varepsilon _k^ - } \right). $

式中：$ {v_k} $为$ k $时刻观测噪声. 则观测噪声期望值为

(12)$ E\left[ {\hat v_k^ - \hat v_k^{ - }} \right] = {R_k}+{H_k}P_k^ - {H_k}. $

由于$ E\left[ {\hat v_k^ - \hat v_k^{ - }} \right] = E\left( {{y_k}+H_k\varepsilon _k^ - } \right){\left( {{y_k} - H_k\varepsilon _k^ - } \right)} $，则观测噪声协方差$ {R_k} $为

(13)$ {R_k} = {\alpha _k}{R_{k - 1}}+\left( {1 - {\alpha _k}} \right)\left( {\varepsilon _k^ - \varepsilon _k^{ -} - H_kP_k^ - {H_k}} \right). $

式中：$ {\alpha _k} $为修正$ {R_k} $的遗忘因子.

$ {Q_k} $为过程噪声协方差，可以表示为

(14)$ \begin{split} {Q_k} =& E\left( {{x_k} - \left( {{A_{k - 1}}{x_{k - 1}}+{B_{k - 1}}{u_{k - 1}}} \right)} \right) \times \\ &{\left( {{x_k} - \left( {{A_{k - 1}}{x_{k - 1}}+{B_{k - 1}}{u_{k - 1}}} \right)} \right)} .\end{split} $

式中：$ {u_{k - 1}} $、$ {A_{k - 1}} $、$ {B_{k - 1}} $分别为$ k - 1 $时刻的控制输入、$ k - 1 $到$ k $时刻状态转移和$ k - 1 $时刻控制输入. 对于$ \hat x_k^ - $到$ {\hat x_k} $的状态转移，有

(15)$ {\hat x_k} - \hat x_k^ - = {\varepsilon _k}+{\omega _{k - 1}} - {A_{k - 1}}{\varepsilon _{k - 1}}. $

式中：$ {\varepsilon _k} $、$ {\varepsilon _{k - 1}} $分别为$ k $和$ k - 1 $时刻的残差后验值，$ {\omega _{k - 1}} $为$ k - 1 $时刻过程噪声. 对式（15）进行协方差运算，则有

(16)$ \begin{split} &E\left( {\left( {{A_{k - 1}}{\varepsilon _{k - 1}}} \right){{\left( {{K_k}H_k\left( {{A_{k - 1}}\left( {{x_{k - 1}} - {{\hat x}_k}} \right)+{\omega _{k - 1}}} \right)} \right)}}} \right) = \\ &\qquad - {A_{k - 1}}{P_{k - 1}}A_{k - 1}{H_k}K_k .\end{split}$

因此，$ E\left( {{\varepsilon _k}w_k} \right) = - {Q_k}+{K_k}H_k{Q_k} $，有

(17)$ E\left( {{{\hat x}_k} - \hat x_k^ - } \right){\left( {{{\hat x}_k} - \hat x_k^ - } \right)} = {P_k}+{Q_k}+{A_{k - 1}}{P_{k - 1}}A_{k - 1}- 2{P_k}. $

又由于$ {\hat x_k} = \hat x_k^ - +{K_k}\left( {{y_k} - {H_k} \hat x_k^ - } \right) $，结合新息的协方差，则$ {Q_k} $可表示为

(18)$ {Q_k} = {\beta _k}{Q_{k - 1}}+\left( {1 - {\beta _k}} \right)\left( {{K_k}{d_k}d_kK_k+{P_k} - {P_{k - 1}}} \right). $

式中：$ {\beta _{{k}}} $为修正$ {Q_k} $的遗忘因子.

采用总变差（TV）与贝叶斯估计相结合的方法，是为了应对现有常见滤波算法在超声波回波信号处理中的局限性. 小波变换的时频局部化特导致其在信号突变处难以准确区分噪声和真实信号，从而产生周期跳跃误差；巴特沃斯滤波器陡峭的截止特性容易导致信号过度平滑，高斯滤波器则存在边缘模糊的问题，两者都会带来明显的相位偏移，且参数设置不当将导致波形过度平滑，影响ToF的估计. 相比之下，TV方法通过分析信号的局部梯度特征，能够动态获得噪声分布情况，为后续处理提供更准确的噪声评估依据.

在参数适应性方面，传统自适应卡尔曼滤波算法对初始参数敏感，不同的参数选择可能导致信号特征丢失或噪声残留的问题. 贝叶斯估计通过概率化的处理方式，可以根据实际观测数据动态调整参数更新的权重，使得滤波过程对初始设置的依赖性有所降低. TV方法与贝叶斯估计的结合，一方面通过梯度分析提供了更可靠的噪声评估，另一方面通过概率框架增强了参数调整的适应性.

本研究提出的TV-BAKF算法包含3个步骤.

1）采用高斯-牛顿法优化窗口内信号的$ \varGamma $，求解经全变分处理后的信号$ {\gamma _j} $(j表示滑动窗口序列号)，以获取差异因子. 设窗口长度为$ w $，信号中某序列点为$ {y_k} $，预处理段为$ {y_k} $前$ w $个信号点，即$ {y_j} \in \left[ {{y_k},{y_{k+w}}} \right] $，则$ \varGamma $的优化区段为

(19)$ \varGamma = \mathop {\min }\limits_{{y_j} \in \left[ {{y_{k - w}},{y_k}} \right]} E\left( {{y_k},{\gamma _k}} \right)+\lambda V\left( {{\gamma _k}} \right), $

(20)$ E\left( {{y_k},{\gamma _k}} \right) = \frac{1}{2}\sum\limits_j^{} {{{\left( {{y_k} - {\gamma _k}} \right)}^2}} ,$

(21)$ V\left( {{\gamma _k}} \right) = \sum\limits_{l=k - w}^w {\left| {{\gamma _{l+1}} - {\gamma _l}} \right|}. $

式中：$ \gamma_l $为滑动窗口内全变分处理后信号.

采用高斯-牛顿法求解预处理后的最终信号$ {\gamma _k} $（迭代过程中的信号记为$ {\hat \gamma _k} $），通过对目标函数Γ进行一阶泰勒展开实现求解：

(22)$ \varGamma \left( {{{\hat \gamma }_k}+\Delta {{\hat \gamma }_k}} \right) \approx f\left( {{{\hat \gamma }_k}} \right)+J{\left( {{{\hat \gamma }_k}} \right)}\Delta {\hat \gamma _k}+o\left( {\Delta {{\hat \gamma }_k}} \right). $

解$ {\hat \gamma _k} $可以通过海森矩阵$ {\boldsymbol{D}} $和雅各布矩阵$ {\boldsymbol{J}} $迭代求解，并输出处理后信号序列号：$[ {\gamma _{k - w}},{\gamma _{k - w}}+ $$ 1, \cdots ,{\gamma _k} $]

(23)$ \begin{split} & {\boldsymbol{D}}\left(\hat{\gamma}_k\right) \Delta \hat{\gamma}_k={\boldsymbol{J}}; \\& \text {s.t. } {\boldsymbol{D}}\left(\hat{\gamma}_k\right) = {\boldsymbol{J}}\left(\hat{\gamma}_k\right) {\boldsymbol{J}}\left(\hat{\gamma}_k\right)^{\mathrm{T}} , \;{\boldsymbol{G}}\left(\hat{\gamma}_k\right) = -{\boldsymbol{J}}\left(\hat{\gamma}_k\right) \varGamma\left(\hat{\gamma}_k\right).\end{split} $

式中：G为梯度矩阵.

2)求解差分因子$ {\theta _k} $. 为了表征窗口内信号平滑度，并基于分段信噪比自适应调整卡尔曼滤波参数，引入差分因子$ {\theta _k} $修正$ {Q_k} $和$ {R_k} $，$ {\theta _k} $可以表示为

(24)$ {\theta _k} = \frac{{V\left( {{\gamma _k}} \right)}}{{V\left( {{y_k}} \right)}} = \frac{{\displaystyle\sum\nolimits_{k - w}^w {\left| {{\gamma _{k+1}} - {\gamma _k}} \right|} }}{{\displaystyle\sum\nolimits_{k - w}^w {\left| {{y_{k+1}} - {y_k}} \right|} }}. $

式中：$ V\left( {{\gamma _k}} \right) < V\left( {{y_k}} \right), 0 < {\theta _k} < 1 .0$. 当$ {\theta _k} $接近0时，噪声能量较大，须适当增大$ {R_k} $，减小$ {Q_k} $；当$ {\theta _k} $接近1.0时，$ {Q_k} $和$ {R_k} $的调节相反. 则$ {\theta _k} $修正$ {Q_k} $和$ {R_k} $的方法如下：

(25)$ {R_k} = {\alpha _k}{R_k}+\left( {1 - {\alpha _k}} \right)\left[ {\frac{1}{{{\theta _k}}} \left( {\varepsilon _k^ - \varepsilon _k^{ - } - H_kP_k^ - {H_k}} \right)} \right] ，$

(26)$ {Q_k} = {\beta _k}{Q_{k - 1}}+\left( {1 - {\beta _k}} \right)\left[ {{\theta _k} \left( {{K_k}{d_k}d_kK_k+{P_k} - {P_{k - 1}}} \right)} \right]. $

3）基于贝叶斯估计计算遗忘因子$ {\alpha _k} $和$ {\beta _k} $. 对于$ {\alpha _k} $和$ {\beta _k} $，假设其先验分布为

(27)$ p\left( {{\alpha _k}} \right) = {\text{Beta}}\left( {{\alpha _k};a,b} \right) = \frac{{\alpha _k^{a - 1}{{\left( {1 - {\alpha _k}} \right)}^{b - 1}}}}{{B\left( {a,b} \right)}}.$

式中：$ a、b $为控制先验分布形状的超参数. 其取值对自适应卡尔曼滤波器的动态特性具有一定的影响. 若$ a $和$ b $取值均较小(如$ a = 1、b = 1 $)，初始不确定性过大，易致滤波器对噪声过度敏感而产生振荡；当$ a $和$ b $的取值相差较大时 (如$ a = 5,b = 1 $)，初始确定性过强，过度依赖历史信息而抑制新数据，导致相位滞后；当$ a $和$ b $均取适中值时，允许滤波器根据实际观测证据动态调整历史与新数据权重，既保持稳定性又能快速适应变化，因此，超参数选用$ a = 2,b = 2 $. 则基于残差构建的似然函数$ p\left( {{\gamma _k}\left| {{\alpha _k}} \right.} \right) $可以表示为

(28)$ p\left( {{\gamma _k}\left| {{\alpha _k}} \right.} \right) \propto \exp \;\left( { - \frac{1}{2}\varepsilon _k\sum\nolimits_\varepsilon ^{ - 1} {{\varepsilon _k}} } \right). $

对于后验分布，$ p\left( {{\alpha _k}\left| {{\gamma _k}} \right.} \right) $表示为

(29)$ p\left( {{\alpha _k}\left| {{\gamma _k}} \right.} \right) \propto \exp\; \left( { - \dfrac{1}{2}d_k\sum\nolimits_e^{ - 1} {{d_k}} } \right) \cdot \alpha _k^{a - 1}{\left( {1 - {\alpha _k}} \right)^{b - 1}}. $

若似然函数与Beta分布共轭，则后验分布仍持Beta分布形式，此时，有

(30)$ p\left( {{\gamma _k}\left| {{\alpha _k}} \right.} \right) \propto \alpha _k^{{n_1}}{\left( {1 - {\alpha _k}} \right)^{{n_2}}}. $

式中：$ {n_1} $表示支持$ {R_{k - 1}} $的证据量，$ {n_2} $表示支持$ {R_k} $的证据量. 此时，后验分布为

(31)$ p\left( {{\gamma _k}\left| {{\alpha _k}} \right.} \right) = {\text{Beta}}\left( {{\alpha _k};a+{n_1},b+{n_2}} \right). $

对于每个事件步$ k $，$ {n_1} $和$ {n_2} $可表示为

(32)$ {n}_{1}=\mathrm{exp}\left(-\dfrac{1}{2}{\epsilon }_{k}{\displaystyle {\sum }_{d}^{-1}{\epsilon }_{k}}\right)\text{，}{n}_{2}=1-{n}_{1}. $

以后验分布均值作为权重$ \alpha $的估计值，则有

(33)$ {\alpha _k} = \frac{{a+{n_1}}}{{a+b+{n_1}+{n_2}}}. $

对于$ {Q_k} $，与$ {R_k} $类似，区别在于似然函数$ p\left( {{\gamma _k}\left| {{\beta _k}} \right.} \right) $可基于新息$ {d_k} $设计，其定义为

(34)$ p\left( {{\gamma _k}\left| {{\beta _k}} \right.} \right) \propto \exp \;\left( { - \frac{1}{2}d_k\sum\nolimits_d^{ - 1} {{d_k}} } \right). $

最终，将遗忘因子$ {\alpha _k} $和$ {\beta _k} $代入式(25)和(26)，对$ {Q_k} $和$ {R_k} $进行更新. 综上，本研究提出的TV-BAKF算法框架如图2所示.

为了验证TV-BAKF方法在噪声回波信号估计中的准确性，基于MATLAB构建符合低频空气耦合探头窄带特性的混合指数模型噪声信号：

(35)$ S\left( t \right) = {A_0}{\left( {\frac{{t - \tau }}{T}} \right)^m}\exp \;\left( { - \frac{{t - \tau }}{T}} \right) \sin \;\left( {2{\text{π}} {f_0}\left( {t - \tau } \right)+\phi } \right) .$

式中：$ t $为采样点时刻，m、$ T $为传感器特性参数(决定波形)，$ \tau $为时延，$ {f_0} $为中心频率，$ \phi $为初始相位(一般为0).

通过向基础信号逐级添加噪声，模拟实际冰雪水回波信号，结果如图3所示.

为了验证所提出的信号滤波方法提取冰雪覆盖物回波信号起振点的能力，在课题组自主研发的低温环境实验舱中开展了冰雪水超声波回波信号测试. 该实验舱配备了造雪设备、水冷却系统、风扇及温/湿度控制系统，可调节舱内温度和湿度，并通过控制造雪机进水与进气比例来改变雪的含水率. 结合温湿度控制与流动水雾化系统，实验舱成功制造出干雪、湿雪和雪浆等多种雪样. 课题组还通过中国北方实地雪样的光谱测试验证了模拟雪样的真实性，基于光谱分析，课题组研发的低温环境舱制造的雪和真实世界雪特性一致. 测试系统采用40 kHz空气探头、示波器、放大器和信号发生器，其中探头距离覆盖物为20 cm. 低温环境舱、各类覆盖物样本及试验装置如图4所示.

由于用于计算污染物深度的ToF参数基于超声波在空气中的传播速度，而温湿度变化会影响波速，为了获得精确的计算结果，采用以下公式对不同环境条件下的波速进行修正^[17]：

(36)$ v\left( {\theta ,h_{\mathrm{m}}} \right) = 331+0.6 \theta +{k_{\mathrm{m}}} {h_{\mathrm{m}}}. $

式中：$v\left( {\theta,h_{\mathrm{m}}} \right)$为超声波在空气温度 $\theta $和相对湿度 ${h_{\mathrm{m}}}$下的传播速度；${k_{\mathrm{m}}}$为湿度相关系数, 取值为0.012.

针对冰雪水介质表面的复杂反射特性（雪表面粗糙多孔、冰表面裂缝）、低温环境及高精度ToF测量对信噪比的要求，首先进行超声探头比选，如图5所示. 其中，t为时间. 通过低温环境箱试验，测试了不同频率、指向性和阻抗层材料的探头. 发现较高频率探头(200 kHz)信号因波长小于表面不平整度，导致强散射，回波幅值低(<0.1V)且信噪比较低；40 kHz波长($\lambda $≈8.5 mm)满足瑞利准则，在雪、冰、水表面均能产生信噪比≥12 dB的清晰且更为稳定的回波. 同时，高指向性探头对平整冰/水表面虽好，但在疏松多孔的粗糙雪面易导致回波能量过低甚至丢失；中等指向性(40°~50°)探头则能兼顾不同介质，避免回波信号能量过低导致低信噪比. 对于换能材料，压电陶瓷凭借高灵敏度、可控余振时间(<50 μs)避免干扰以及优良的温度稳定性，能有效捕捉更为微弱的干/湿雪面反射信号. 综上，最终选定40 kHz, 指向性50°的空气耦合探头.

对R(观测噪声协方差)、Q(过程噪声协方差) 、P(误差协方差) 、K(卡尔曼增益系数)进行分析. TV-BAKF对模拟指数信号施加不同噪声功率后，滤波效果如图6所示. 结果表明，该方法在不同信噪比条件下均能有效平滑信号，显著减少毛刺干扰. 特别是信号前沿的前3个周期波形清晰可辨，即使在14 dB的低信噪比情况下，仍能准确提取信号前沿特征.

TV-BAKF通过总差分预处理技术能够获取信号的噪声分布特征并实现参数动态调整. 在信号处理过程中，Q和R会根据信号特征进行协同调节：在信噪比较低的零线及信号前沿段，Q自适应性降低使滤波器更信任预测结果，同时R相应增大以抑制测量噪声影响；在波包段则呈现相反趋势，Q增大而R自适应性减小以提高观测权重，这种自适应调节在起振点前后表现尤为明显（如黑色虚线标示），特别是在低信噪比区域（信噪比<20 dB）调节幅度明显. 实验结果表明，该机制使TV-BAKF在14~40 dB的信噪比范围内均能保持良好的滤波性能，有效抑制噪声干扰.

通过P的变化趋势可以看出，在信号零线段P明显增大，而在波包前沿段则显著下降，特别是在第3个周期位置P降至1.0以下，这表明系统状态估计精度明显提升. 与之相对应的是， K呈现出与P相反的变化规律：在高信噪比区域K维持在较低水平，说明此时滤波器更倾向于采用预测值；而在信号前3个周期之后，K明显增大，表明测量值的权重显著提高. 值得注意的是，在低信噪比条件下，K整体降低，这种自适应调整使滤波器更依赖系统状态估计值，从而确保在强噪声干扰下仍能保持滤波结果的稳定性.

如图7所示展示了贝叶斯估计动态调整遗忘因子α 和 β的结果. 在信号零线位置，由于信号变化平缓且观测噪声占主导，α 和 β趋近于1.0，滤波器更信任预测值；而在波包段信号变化显著时，α 和 β明显减小，使滤波器更依赖当前观测数据. 信噪比的影响表现为：低信噪比时α 和 β趋近1.0（依赖历史信息），高信噪比时趋近0（信任当前观测）.

值得注意的是，α 的变化幅度明显大于β的. 这是因为α用于调节观测噪声协方差 R，由于观测噪声易受外部环境干扰，须扩展其调节区间以更好地权衡历史信息与实时观测数据. 因此，通过贝叶斯估计动态调整遗忘因子，卡尔曼滤波能够自适应地平衡历史信息和当前观测数据，从而更好地适应系统状态和噪声特性的变化.

通过设置4种差异较大的初始参数组合(即Q=10, R=0.1; Q=0.1, R=10; Q=10, R=10; Q=0.01, R=0.01）对滤波结果进行比较分析，结果如图8所示. 可以看出，无论初始Q、R设定如何，TV-BAKF方法都能通过捕捉当前信号波形特征及噪声分布，动态调整Q、R至最优范围. 也可以看出，在4种Q、R参数组合下，卡尔曼增益在波包段的峰值均接近0.45，而P在零线噪声段的最高值均接近80，在波包段则均下降至接近0. 因此，分析其参数动态变化可见，在低信号幅值区段，Q均自适应地减小而R相应增大，以更好抑制测量噪声并保持稳定；在高信号能量区段，Q增大而R减小，以增强对信号动态变化的适应性. 由此可见，TV-BAKF方法对初始参数敏感性低，能够依据信号本身的动态特性和噪声水平进行实时调整，因而对初始参数的选择具有较高鲁棒性.

为了直观地比较传统卡尔曼滤波和本研究提出的优化卡尔曼滤波方法的滤波效果，图9展示了采用不同卡尔曼滤波方法所得的Lissajous图，该图能够直观评价滤波后信号保真能力，如相位滞后程度与波形平滑度. 在Lissajous图中，相位滞后现象严重，表现为椭圆曲线更宽，而椭圆的圆滑则表征信号平滑程度.

如图9所示为不同卡尔曼滤波方法(及Q、R组合)的Lissajous图，其中${A_\gamma }$和${A_x}$分别表示去噪后信号和原信号的幅值. 可以看出，传统EKF方法的性能受Q、R参数显著影响：增大Q使滤波器更信任观测值，虽减小相位滞后但导致噪声残留和平滑度下降；增大R则增强对预测值的依赖，虽改善平滑度却加剧相位滞后，揭示出相位保真度与平滑度之间的矛盾. 传统AKF方法通过自适应调整Q、R在一定程度上缓解了这一矛盾，但在低信噪比条件下因参数修正幅度过大，反而导致相位滞后增大且平滑度降低. 相比之下，TV-BAKF方法展现出显著优势：在不同信噪比条件下均能维持较小的椭圆曲线半径，即便在14 dB时相位滞后仅轻微增加，同时在低信噪比环境下椭圆曲线依旧保持圆滑，表明滤波后信号有较高的平滑度.

进一步地，本研究对目前先进的滤波方法和TV-BAKF方法的滤波效果进行评估，并分别从均方误差(MSE)、平滑度($ {J_{{\text{smooth}}}} $)、相位保真度（$ {J_{{\text{phase}}}} $）、幅值保真度($ {J_{{\text{amp}}}} $) 4个维度评价滤波效果和波形保持能力. 4个指标表达式如下：

(37)$ {\mathrm{MSE}} = \frac{1}{N}{\sum\nolimits_{i = 1}^N {\left( {{y_{{i}}} - {\gamma _{{i}}}} \right)} ^2}, $

(38)$ {J_{{\text{smooth}}}} = \frac{1}{N}{\sum\nolimits_{i = 1}^N {\left( {\frac{{{{\mathrm{d}}^2}{\gamma _{{i}}}}}{{{\mathrm{d}}{\gamma _{{i}}}^2}}} \right)} ^2}, $

(39)$ {J_{{\text{phase}}}} = \sqrt {\frac{1}{N}{{\sum\nolimits_{i = 1}^N {\left( {{\phi _{{y_{{i}}}}} - {\phi _{{\gamma _{{i}}}}}} \right)} }^2}}, $

(40)$ J_{\mathrm{amp}}=\dfrac{\displaystyle{\sum}_{i=1}^N\left(y_{{i}}-\gamma_{{i}}\right)^2}{\displaystyle{\sum}_{i=1}^N y_{{i}}^2}.$

式中：$ {\gamma _{{i}}} $、$ {y_{{i}}} $分别表示滤波后及滤波前信号，$ \phi $表示相位.

如表1所示列举了5种滤波方法（改进小波变换^[18]、改进的S-滤波^[19]、UAKF^[20]、VBAKF^[15]以及TV-BAKF）的滤波效果. 5种滤波方法的对比分析显示，在信噪比较高（25 dB）时，改进小波变换的均方误差最小，表现最佳；然而在低信噪比（14 dB）下，其均方误差超过0.5，滤波性能下降. 改进S-G滤波在低信噪比（14 dB）下误差最大，达到约0.82，2种改进卡尔曼滤波方法（UAKF和VBAKF）的均方误差均超过0.4，而TV-BAKF方法的误差始终稳定在0.35以下，展现出其良好的低信噪比下抗噪能力. 从平滑度指标来看，UAKF在低信噪比下平滑效果最显著，VBAKF平滑度最高，但其相位保真度超过0.20 rad，存在明显相位滞后；相比之下，改进小波变换和TV-BAKF方法的相位保真度均低于0.15 rad. 在幅值保真度方面，UAKF和VBAKF在20 dB时低于0.90，在14 dB时进一步降至0.85以下，表明其过度平滑导致波形失真. 虽然改进小波变换在保真度上略优，但其低信噪比下的误差和平滑度表现较差. 综合来看，TV-BAKF在噪声抑制能力与信号特征保持之间实现了良好平衡.

根据ICAO《Circular 355》GRF要求，重点研究了水、冰、雪浆、湿雪和干雪5种污染物状态，同时对冰水、冰雪、水雪混合覆盖物进行了测试. 通过分析不同覆盖物的回波信号特征发现：水和冰的回波信号具有较高的幅值和信噪比，波形与标准指数模型吻合度好；湿雪和雪浆则表现出幅值降低，波形畸变、波形前沿畸变严重，及信号拖尾延长等特征；干雪由于多孔结构和粗糙表面，导致信号能量大幅衰减，信噪比显著下降. 混合覆盖物的信号特征介于其组成成分之间，其信噪比和幅值分布包含在水和干雪的区间范围内. 不同覆盖物(包括水/冰/雪浆/湿雪/干雪和混合覆盖物)的信噪比和波形特征研究流程如图10所示. 进一步分析表明，冰、水及冰水混合物的信号特征集中分布在高信噪比区域，而各类雪覆盖物则分散在低信噪比区域，其中干雪的信号最为离散. 基于研究低信噪比条件下ToF估计精度的核心目标，并响应GRF规范需求，本研究选取水、冰、雪浆、湿雪和干雪5种典型GRF覆盖物作为研究对象.

选取某次检测中的低信噪比回波信号进行分析. 首先通过深度测量和波形前沿分析确定信号前若干周期的起振点位置作为基准. 随后采用双阈值法，对比3种去噪方法（小波、AKF和TV-BAKF）处理后的起振点估计结果. 双阈值法通过设置2个电压阈值，获取信号与阈值交叉后的首个零交叉点作为特征点，进而估计起振点位置.

如图11所示展示了3种滤波方法的处理效果. 小波去噪在波包段效果较好，但信号前沿和零线位置的噪声滤除不足，可能导致双阈值法误判. AKF方法虽然使信号前沿更平滑，但过度依赖预测值导致幅值畸变，影响特征点判断. 相比之下，TV-BAKF通过总变差预处理和贝叶斯动态调整Q、R值，显著提升了零线与波包段的分割准确性，在保持波形平滑度的同时有效减少了相位滞后，展现出良好的信号前沿波形特征提取能力.

在4.1.4节已经通过仿真测试对比分析了改进小波变换、改进S-G滤波、AUKF、VBAKF和TV-BAKF方法的滤波效果. 接下来，利用上述5种滤波方法，对5种道面覆盖物的超声波回波信号进行滤波，采用双阈值方法估计起振点位置，并计算覆盖物深度. 如图12所示展示了5种道面覆盖物下，3种信号处理方法在不同信噪比时的深度误差(${D_{{\text{error}}}}$)的分布，图13以小提琴图形式对比了5种方法的误差状况. 可以看出，水面表面平滑，回波信号信噪比高，5种去噪方法的深度误差均不超过0.55 mm，其中改进小波去噪精度最优. 冰面因表面裂纹，冰颗粒导致信噪比略低但仍高于22 dB，各方法误差均小于0.65 mm，但AUKF和VBAKF因相位滞后问题误差稍高. 雪浆在信噪比为18~24 dB时误差随信噪比下降而增大，改进小波和VBAKF的误差接近0.80 mm，改进S-G和AUKF的误差超过0.88 mm，TV-BAKF的误差仍小于0.55 mm. 湿雪信噪比更低（14~18 dB），改进小波和S-G滤波因噪声抑制不足误差超过1.25 mm，AUKF因过度平滑误差达1.32 mm，VBAKF因波形畸变误差超1.20 mm，而TV-BAKF误差小于0.70 mm. 干雪因表面松散多孔导致信噪比极低（＜16 dB），改进小波和S-G滤波误差超1.25 mm，AUKF误差接近2.00 mm，VBAKF超1.15 mm，TV-BAKF虽因轻微相位滞后误差增至0.96 mm，但仍优于其他方法. 总体而言，本研究提出的TV-BAKF方法，对比的其他方法，对于低信噪比的回波信号的信号前沿波形提取更为准确，同时相比AUKF和VBAKF，该方法有效改善了相位滞后和波形特征提取问题，同时减小了信号畸变，更有利于后续双阈值法、互相关法对深度的准确计算.

提出结合总变差和贝叶斯估计的自适应卡尔曼滤波方法(TV-BAKF)，旨在提升对道面覆盖物(水/冰/雪)回波信号的ToF估计精度.

(1) TV-BAKF方法利用总变差法检测滑动窗口信号内的噪声分布，结合贝叶斯估计计算差分因子以反映信号变化率与噪声水平的相对关系，并借助残差矩阵和创新矩阵构建的自适应调整机制，结合差分因子加权动态调整Q和R.

(2) TV-BAKF方法在处理低信噪比信号时表现出色，尤其是在信号的前沿位置和零线位置. 信号模拟实验表明，该方法能够准确提取低信噪比信号的波形，保持信号的平滑度，同时减少相位滞后，提高信号的起振点估计精度.

(3) TV-BAKF方法在处理较高信噪比(水、冰)覆盖物回波信号时，计算误差和小波及VBAKF方法差异较小；而对于松散且表面粗糙(干/湿雪)覆盖物回波信号，TV-BAKF相比其他方法仍能够对ToF进行准确估计，并能够将深度精度提升至1 mm以内.

本研究还须进一步优化TV-BAKF方法中滑动窗口尺寸的选择，以及其在非平稳噪声（如突发干扰）场景下的适应性表现. 针对现有不足，后续研究将重点优化滑动窗口的动态调整策略，在保证滤波精度的同时降低计算复杂度；同时，通过改进噪声协方差矩阵的更新机制，增强算法对极端噪声条件的鲁棒性，进一步提升滤波稳定性.