变形监测所获取的数据往往包含真值和误差(噪声)2个部分，对变形数据进行预处理可以有效地减小观测误差，提取变形特征，分析变形规律^[1]。在实际应用中，受仪器自身精度、周围环境和人为操作等多种因素的共同影响，噪声和真实信号相互叠加，使得变形数据呈现出一定的波动性。传统的离散小波变换(discrete wavelet transform，DWT)具有良好的时频特性，能较好地反映信号的局部特征，被广泛应用于变形分析中^[2]。众多学者从小波分解层数、阈值估计、阈值函数去噪等方面对传统小波去噪进行研究，取得了较好的效果^[3]。然而DWT仍存在几个方面的劣势：信号变换有平移敏感性、缺乏平移不变性、分解频率混叠、分解的方向有限等，主要原因是小波变换不能有效描述高维信息的线状奇异性或面状奇异性。为了克服离散小波变换存在的缺陷，Kingsbury等^[4]首次提出双树复小波变换(dual-tree complex wavelet transform，DTCWT)的概念；Selesnick等^[5]进一步提出双树复小波变换分解与重构算法。双树复小波变换具有近似平移不变性、良好的方向选择性、有限的数据冗余性、完全重构性和计算效率高等优良特性。

本文将双树复小波引入到变形监测数据处理中，探讨不同阈值函数下的去噪效果，通过仿真和实际算例，并与db10、sym6小波进行对比分析，验证其有效性和可行性。

1 去噪模型原理 1.1 双树复小波原理

在实际工程中，原始变形监测数据都是呈离散分布的。如果要对数据进行处理，需要将连续小波变换进行离散化，得到离散小波。双树复小波变换是采用二叉树结构(实部树a和虚部树b)的两路离散小波变换，树a生成实部，树b生成虚部。对于一维数据而言，其分解流程如图 1所示，其中h₀(n)、h₁(n)表示共轭正交滤波器对，g₀(n)、g₁(n)表示共轭积分滤波器对，↓2表示各点的采样。

根据DTCWT的构造方法，双树复小波可表示为：

$ \psi \left( t \right) = {\psi _h}\left( t \right) + {\rm{i}}{\psi _g}\left( t \right) $

(1)

式中，ψ_h(t)、ψ_g(t)为2个实小波，i为复数单位。

由于DTCWT由2个平行的小波变换组成，因此根据小波理论，图 1中树a上面的实部树小波变换的小波系数和尺度系数可由式(2)和(3)计算：

$ \begin{array}{*{20}{c}} {d_j^{{\mathop{\rm Re}\nolimits} }\left( n \right) = {2^{j/2}}\int_{ - \infty }^{ + \infty } {x\left( t \right){\psi _h}\left( {{2^j}t - n} \right){\rm{d}}t} }\\ {j = 1,2, \cdots ,J} \end{array} $

(2)

$ c_J^{{\mathop{\rm Re}\nolimits} }\left( n \right) = {2^{J/2}}\int_{ - \infty }^{ + \infty } {x\left( t \right){\psi _h}\left( {{2^J}t - n} \right){\rm{d}}t} $

(3)

同理，虚部树小波变换的小波系数和尺度系数可由式(4)和(5)计算：

$ d_j^{{\mathop{\rm Im}\nolimits} }\left( n \right) = {2^{j/2}}\int_{ - \infty }^{ + \infty } {x\left( t \right){\psi _g}\left( {{2^j}t - n} \right){\rm{d}}t} $

(4)

$ c_J^{{\mathop{\rm Im}\nolimits} }\left( n \right) = {2^{J/2}}\int_{ - \infty }^{ + \infty } {x\left( t \right){\psi _g}\left( {{2^J}t - n} \right){\rm{d}}t} $

(5)

由此可得DTCWT的小波系数和尺度系数为：

$ d_j^{\left( \psi \right)}\left( n \right) = d_j^{{\mathop{\rm Re}\nolimits} }\left( n \right) + {\rm{i}}d_j^{{\mathop{\rm Im}\nolimits} }\left( n \right) $

(6)

$ c_J^{\left( \psi \right)}\left( n \right) = c_J^{{\mathop{\rm Re}\nolimits} }\left( n \right) + {\rm{i}}d_J^{{\mathop{\rm Im}\nolimits} }\left( n \right) $

(7)

通过式(8)和(9)对DTCWT的小波系数和尺度系数进行重构：

$ \begin{array}{*{20}{c}} {{d_j}\left( t \right) = {2^{\left( {j - 1} \right)/2}} \cdot }\\ {\left[ {\sum\limits_{n = - \infty }^\infty {d_j^{{\mathop{\rm Re}\nolimits} }\left( n \right){\psi _h}\left( {{2^j}t - n} \right)} + \sum\limits_{k = - \infty }^\infty {d_j^{{\mathop{\rm Im}\nolimits} }\left( n \right){\psi _g}\left( {{2^j}t - k} \right)} } \right]} \end{array} $

(8)

$ \begin{array}{*{20}{c}} {{c_J}\left( t \right) = {2^{\left( {J - 1} \right)/2}} \cdot }\\ {\left[ {\sum\limits_{N = - \infty }^\infty {d_J^{{\mathop{\rm Re}\nolimits} }\left( n \right){\psi _g}\left( {{2^J}t - n} \right)} + \sum\limits_{k = - \infty }^\infty {c_j^{{\mathop{\rm Im}\nolimits} }\left( n \right){\psi _g}\left( {{2^J}t - k} \right)} } \right]} \end{array} $

(9)

因此，双树复小波变换可对分解后的重构信号表示为：

$ x\left( t \right) = {d_j}\left( t \right) + {c_J}\left( t \right) $

(10)

最后，原始信号x(t)的重构信号可表示为：

$ \hat x\left( t \right) = {T^{{\mathop{\rm Re}\nolimits} }} + {\rm{i}}{T^{{\mathop{\rm Im}\nolimits} }} $

(11)

式中，T^Re为小波分解过程中的实部树系数，T^Im为小波分解过程中的虚部树系数。

1.2 去噪流程和质量评估

去噪流程如下：1)对变形监测数据分别进行离散小波和双树复小波变换，得到一系列的小波系数ω_{j, k}，并对其进行对比分析；2)对小波系数分别进行多种阈值处理，得到新的小波系数$ {\mathit{\hat \omega }_{\mathit{j}, \mathit{k}}}$；3)对新的小波系数进行重构，得到去噪后的信号。

为了更好地评价变形去噪的效果，采用信噪比(SNR)、均方根误差(MSE)、平滑度(r)和相关系数(R)进行评估，具体计算公式见文献[6]。

相关系数(R)一般指原始信号(含噪信号)与去噪信号之间的相似度。由于在实际作业中变形的真值是未知的，对去噪效果的评估并不准确，一般只能作为辅助指标。然而在仿真数据中，由于真值已知，因此可直接用真值作为参考值与去噪后数据对比相似度，得出最准确的去噪评价。

2 算例分析 2.1 仿真实验

模拟一组真值已知的变形监测仿真数据进行分析。经MATLAB软件编程，生成一组长度为1 200的变形监测仿真数据，由3个振幅频率不相同的周期项和1个趋势项组成。其表达式为：

$ \begin{array}{*{20}{c}} {S = 3\sin \left( {2{\rm{ \mathsf{ π} }}t/500} \right){{\rm{e}}^{ - t/1\;000}} + 2\sin \left( {2{\rm{ \mathsf{ π} }}t/200} \right) + }\\ {2\sin \left( {2{\rm{ \mathsf{ π} }}t/50} \right) + t/1\;000} \end{array} $

(12)

由于变形监测数据中的噪声多以白噪声为主，因此对该仿真数据添加信噪比为db20的高斯白噪声。原始数据序列和含噪数据序列见图 2。

采用db10、sym6和双树复小波对含噪数据进行4层分解，利用SNR、MSE和R综合对比分析。鉴于该信号噪声较为简单，统一采用较为保守的极大极小原理计算阈值。在阈值函数方面，使用硬阈值、软阈值和文献[7]的改进阈值进行对比实验，各种小波分解效果见图 3。

由图 3可知，双树复小波中d₁和d₂层的高频信号比db10和sym6小波的图像更为密集，可见其对高频信号的分离能力强于传统离散小波。对比d₃和d₄层，db10和sym6小波都出现了较大程度的信号振荡与频率混叠现象，而双树复小波产生的频率混叠与振荡极少。对于a₄层，传统离散小波出现严重失真现象，这主要是由于a₄层信号失真部分已分离到d₄层，并在d₄层发生频率混叠现象。从图 3(a)与3(b)画圈部分发现，d₄层出现频率混叠的时间轴刚好与a₄层失真部位吻合。由于小波重构是基于a₄层信号进行的，出现失真会直接影响到信号重构的效果，可见双树复小波的分解效果较优。

各小波使用极大极小阈值估计原理改进阈值去噪的效果(图 4)。限于篇幅，选取去噪效果最好的方案进行分析。

由图 4可以看出，d₁~d₃层各小波去噪效果相同，采用较为保守的极大极小原理得出的阈值均可对各层小波系数完全置零，这主要是由于添加的高斯白噪声较为单一，容易被分离并去除。对于db10和sym6小波分解中标出的信号振荡、频率混叠处，均在d₄层去噪后保留了下来，去噪后图像依然与a₄层信号失真处吻合。各小波去噪重构后的去噪图像见图 5。

由图 5可知，db10小波的去噪图像在0~40处出现严重的信号变形，结合图 3(a)分析，这主要是由于信号分解时信噪分离不够彻底造成的。采用阈值估计和阈值函数均存在着一定的缺陷，对d₄层去噪后破坏了有用信号的结构，使得重构后出现失真。由于仿真数据的噪声较为单一，较为容易去除，使得sym6小波与双树复小波都与真实信号吻合得较好。

基于不同阈值函数3种小波去噪后的质量评价见表 1。对比硬阈值、软阈值和改进阈值发现，改进阈值的去噪效果要优于传统阈值函数，主要是由于改进阈值是在硬阈值、软阈值的基础上改进的。双树复小波针对各种阈值函数的去噪效果相同，说明其信噪分离到一种较为理想的状态，优于db10和sym6小波。可见，有效的信噪分离能在一定程度上避开或削弱阈值函数存在的固有缺陷。

2.2 实测数据

实际工作中，变形监测含有的噪声往往是随机且无固定规律的，将双树复小波引入到实际变形监测数据分析中，进一步验证其性能。选取文献[7]某市地铁施工建筑物变形监测中一组倾斜值的监测数据，见图 6。

由图 6看出，该建筑物的变形监测数据含有一定的噪声，带有一定的毛刺，总体呈逐渐向下倾斜的趋势。同算例1，采用3种小波对其进行4层分解，d₁~d₄分解结果见图 7。

由图 7可知，3种分解方法均能有效分离出变形数据中的不同频率信息。对比各分量发现，双树复小波分解得到的结果充分表现出周期性变化特征，较好地表现了细节部分的频率信息，可见双树复小波在分解信号过程中表现得更为彻底。同算例1，各层的去噪效果见图 8。

由图 8看出，极大极小原理得出的阈值在3种小波中都对信号d₁层的小波系数完全置零了，而d₃、d₄层的信号保留完好。db10小波由于其信噪分离得不彻底，导致d₂层含有噪声而对变形序列进行了去噪。重构后的变形序列见图 9。

由图 9可得，db10的去噪图像过于平滑，失去了很多细节成分，这主要是由于其将d₂层有用信号当作噪声去除了，导致图像失真。在极大极小原理计算出的阈值向量下，双树复小波和sym6小波都只对d₁层进行去噪。因此，去噪后的图像中差别不明显，需要借助评价指标来进行综合评价，见表 2。

由表 2可知，双树复小波的所有去噪指标均最优；db10小波的去噪结果过于平滑，失去了很多信号细节，导致信号失真，只在平滑度上较优。结合相关系数R来看，双树复小波对信号的还原度要优于其他2种传统离散小波，这主要是由于双树复小波拥有完全重构性等特点，使得重构后的去噪信号还原度更高，更为接近原始信号。综上，双树复小波不论在信号分离还是去噪效果方面，均优于db10小波和sym6小波。

3 结语

本文将双树复小波应用于变形监测数据去噪中，通过理论分析和实例表明：

1) 双树复小波的信噪分离能力优于传统离散小波，产生的频率混叠现象与信号震荡现象也较少；分解后的高频系数比传统离散小波更多地集中于低层次，低频近似系数能更为清晰地分离出来。

2) 双树复小波能够更好地提取变形中隐含的周期性波动和信号的低频近似细节。根据其特有的完全重构性，去噪后能更好地保留变形监测数据细节和周期波动，重构后得出的变形序列更为接近原始序列，且去噪后的曲线波动比传统离散小波更为平滑。

3) 当信噪分离效果较为理想时，各层阈值向量随着层数的增高而减小，直至接近或甚至为0。此时，不同阈值函数对信号去噪效果差别较小，在极端情况下得出的结果会完全相同，从而在一定程度上削弱了阈值函数的缺陷。该方法运用于变形监测数据降噪是可行、有效的。