本文展示了双树复小波变换 (Dual-Tree Complex Wavelet Transforms，DTCWT) 在信号、图像和轴承故障诊断方面的应用。 DTCWT 为两个独立的两通道滤波器组，在实际应用用，不能随意选择两棵树中使用的尺度小波滤波器。第一棵树 {h0 ,h1 } 的低通（尺度）和高通（小波）滤波器生成一个尺度函数和小波，另一棵树是由第一棵树的尺度函数的近似希尔伯特变换以及相应的高通滤波器生成的小波函数组成，记作{g0，g1 }。因此由两棵树形成的复值尺度函数和小波函数是近似解析的，因此，DTCWT比DWT具有更小的移位方差和更大的方向选择性。DTCWT的冗余度明显小于未抽取的DWT 的冗余度。本文验证了 DTCWT 的近似移位不变性，DTCWT在2-D和3-D数据中的方向选择性，以及DTCWT变换在图像去噪等方面的应用。
双树复小波变换DTCWT基本理论
双树复小波变换DTCWT采用二叉树结构的两路滤波器组进行信号的分解和重构，第一棵树生成实部，第二棵生成虚部，合理设计实、虚部树低通滤波器，满足半采样延迟条件，具有近似平移不变性。两树的滤波器采样频率相同，但是它们之间的延迟恰好是一个采样间隔，这样虚部树中第１层的二抽取恰好采到实部树中二抽取所丢掉的采样值，在获得了复小波变换的平移不变性的同时避免了大量的计算并且具有容易实现的优势。下图为３层双树复小波的分解和重构过程。
３层双树复小波的分解和重构过程
虚线上方实部树变换的小波系数和尺度系数可由式（2）和（3）计算
相关的参考文献如下：
[1]Huang Tongyuan,Xu Jia,Yang Yuling,Han Baoru. Robust Zero-Watermarking Algorithm for Medical Images Using Double-Tree Complex Wavelet Transform and Hessenberg Decomposition[J]. Mathematics,2022,10(7).
[2]Zhou Yilu,Fu Xiaojin. Image Denoising Based on Dual-tree Complex Wavelet Transform and Convolutional Neural Network[J]. Journal of Physics: Conference Series,2021,1995(1).
[3]Lei Wang,Zhiwen Liu,Hongrui Cao,Xin Zhang. Subband averaging kurtogram with dual-tree complex wavelet packet transform for rotating machinery fault diagnosis[J]. Mechanical Systems and Signal Processing,2020,142(C).DTCWT的近似移位不变性
DWT受到移位变化的影响，这意味着输入信号或图像的微小偏移将导致DWT系数中跨尺度的信号或图像能量分布发生显着变化。 DTCWT是近似移位不变的。为了在测试信号中证明这一点，构建两个移位的离散时间脉冲，长度为128 个样本。一个信号在样本60处具有单位脉冲，而另一个信号在样本64处具有单位脉冲，两个信号显然都具有单位能量（l2)范数。
kronDelta1 = zeros(128,1);
kronDelta1(60) = 1;
kronDelta2 = zeros(128,1);
kronDelta2(64) = 1;
将DWT的扩展模式设置为周期性。使用长度为14的小波和尺度滤波器获得两个信号的DWT 和DTCWT，直到第3 级，提取第3级的细节系数进行比较。
origmode = dwtmode('status','nodisplay');
dwtmode('per','nodisp')
J = 3;
[dwt1C,dwt1L] = wavedec(kronDelta1,J,'sym7');
[dwt2C,dwt2L] = wavedec(kronDelta2,J,'sym7');
dwt1Cfs = detcoef(dwt1C,dwt1L,3);
dwt2Cfs = detcoef(dwt2C,dwt2L,3);
[dt1A,dt1D] = DTCWT(kronDelta1,'Level',J,'FilterLength',14);
[dt2A,dt2D] = DTCWT(kronDelta2,'Level',J,'FilterLength',14);
dt1Cfs = dt1D{3};
dt2Cfs = dt2D{3};
绘制第3级两个信号的DWT 和DTCWT系数的绝对值，并计算系数的能量（平方l2范数），以相同的比例绘制系数。信号中的 4 个样本移位导致第3级 DWT 系数的能量发生了显着变化，第 3 级DTCWT 系数中的能量仅变化了大约 3%。
figure
subplot(1,2,1)
stem(abs(dwt1Cfs),'markerfacecolor',[0 0 1])
title({'DWT';['Squared 2-norm = ' num2str(norm(dwt1Cfs,2)^2,3)]},...
'fontsize',10)
ylim([0 0.4])
subplot(1,2,2)
stem(abs(dwt2Cfs),'markerfacecolor',[0 0 1])
title({'DWT';['Squared 2-norm = ' num2str(norm(dwt2Cfs,2)^2,3)]},...
'fontsize',10)
ylim([0 0.4])
figure
subplot(1,2,1)
stem(abs(dt1Cfs),'markerfacecolor',[0 0 1])
title({'Dual-tree CWT';['Squared 2-norm = ' num2str(norm(dt1Cfs,2)^2,3)]},...
'fontsize',10)
ylim([0 0.4])
subplot(1,2,2)
stem(abs(dwt2Cfs),'markerfacecolor',[0 0 1])
title({'Dual-tree CWT';['Squared 2-norm = ' num2str(norm(dt2Cfs,2)^2,3)]},...
'fontsize',10)
ylim([0 0.4])
为了证明近似移位不变性在实际数据中的作用，我们分析心电图 (ECG) 信号。 ECG 信号的采样间隔为 1/180 秒，数据取自 Percival & Walden[4]，第 125 页（数据最初由华盛顿大学的 William Constantine 和 Per Reinhall 提供）。为方便起见，我们将数据从 t=0 开始。
[4]Percival, Donald B., and Andrew T. Walden. Wavelet Methods for Time Series Analysis. Cambridge Series in Statistical and Probabilistic Mathematics. Cambridge ; New York: Cambridge University Press, 2000.
load wecg
dt = 1/180;
t = 0:dt:(length(wecg)*dt)-dt;
figure
plot(t,wecg)
xlabel('Seconds')
ylabel('Millivolts')
相隔约0.7 秒的较大的正值峰是心律R波。首先，使用带有 Farras近似对称滤波器的临界采样DWT分解信号。绘制原始信号以及第2级和第3级小波系数。之所以选择第2级和第3级系数，是因为对于给定的采样频率，R波在这些尺度中被分离的较为显著。
figure
J = 6;
[df,rf] = dtfilters('farras');
[dtDWT1,L1] = wavedec(wecg,J,df(:,1),df(:,2));
details = zeros(2048,3);
details(2:4:end,2) = detcoef(dtDWT1,L1,2);
details(4:8:end,3) = detcoef(dtDWT1,L1,3);
subplot(3,1,1)
stem(t,details(:,2),'Marker','none','ShowBaseline','off')
title('Level 2')
ylabel('mV')
subplot(3,1,2)
stem(t,details(:,3),'Marker','none','ShowBaseline','off')
title('Level 3')
ylabel('mV')
subplot(3,1,3)
plot(t,wecg)
title('Original Signal')
xlabel('Seconds')
ylabel('mV')
进行DTCWT分析，只取第2级和第3级双树系数的实部，如下图所示
临界采样DWT和双树小波变换DTCWT都将 ECG 波形的一个重要特征定位到相似的尺度。小波在一维信号中的一个重要应用是获得按尺度的方差分析。按理说，这种方差分析不应该对输入信号的循环移位敏感。不幸的是，严格采样的DWT并非如此。为了证明这一点，我们将 ECG 信号循环移位4个样本，采用严格采样的DWT分析未移动和移动的信号，并计算跨尺度的能量分布。
wecgShift = circshift(wecg,4);
[dtDWT2,L2] = wavedec(wecgShift,J,df(:,1),df(:,2));
detCfs1 = detcoef(dtDWT1,L1,1:J,'cells');
apxCfs1 = appcoef(dtDWT1,L1,rf(:,1),rf(:,2),J);
cfs1 = horzcat(detCfs1,{apxCfs1});
detCfs2 = detcoef(dtDWT2,L2,1:J,'cells');
apxCfs2 = appcoef(dtDWT2,L2,rf(:,1),rf(:,2),J);
cfs2 = horzcat(detCfs2,{apxCfs2});
sigenrgy = norm(wecg,2)^2;
enr1 = cell2mat(cellfun(@(x)(norm(x,2)^2/sigenrgy)*100,cfs1,'uni',0));
enr2 = cell2mat(cellfun(@(x)(norm(x,2)^2/sigenrgy)*100,cfs2,'uni',0));
levels = {'D1';'D2';'D3';'D4';'D5';'D6';'A6'};
enr1 = enr1(:);
enr2 = enr2(:);
table(levels,enr1,enr2,'VariableNames',{'Level','enr1','enr2'})
注意到第3级和第4级的小波系数显示原始信号和移位信号之间的能量变化约为 3%。接下来，我们使用DTCWT进行分析，结果如下所示。
由结果可得，DTCWT变换为原始信号及其循环移位版本生成按比例的一致性方差。
图像处理中的方向选择性
我们绘制一下具有 4阶消失矩 sym4 的 Daubechies 最小不对称相位小波的 LH、HL 和 HH滤波器
由图可以看出LH和HL小波滤波器分别具有清晰的水平和竖直方向，最右边的HH小波滤波器混合了 45 度和-45 度方向，产生了伪影。这种方向的混合是由于使用了实值可分离滤波器。HH实值可分离滤波器在2维频率平面的所有四个高频角都有通带。DTCWT通过使用近似解析的小波实现方向选择性，在DTCWT中，实部和虚部有6个子带。6个子带的实部是通过将列滤波后的输出相加，然后对两棵树中的输入图像进行行滤波而形成的。6个子带的虚部是通过减去列滤波和行滤波的输出而形成的。应用于列和行的滤波器可能来自相同的滤波器 {h0，h1 } 或 {g0，g1}，或者来自不同的滤波器 {h0，g1},{g0，h1}。下面绘制对应于DTCWT实部和虚部的12个小波滤波器的方向，图的第1行是6个小波滤波器的实部，第2行是虚部。
二维边缘表示
DTCWT的近似解析性和方向选择性在图像边缘的表示方面提供了优于标准2维DWT的性能。为了说明这一点，我们使用临界采样的2-D DWT和 2-D DTCWT来分析测试图像。首先，导入图像
将图像分解到第 4 级，并根据第 4 级细节系数重构图像。
接下来，分析一个具有双曲边的八边形。
load woctagonHyperbolic
figure
imagesc(woctagonHyperbolic)
colormap gray
title('Octagon with Hyperbolic Edges')
axis equal
axis off
同样将图像分解到第 4 级，并根据第 4 级细节系数重构图像
注意到两幅图像的二维 DTCWT 中均不存在二维严格采样 DWT 中明显的伪影。 DTCWT 再现了直线和曲线奇点。
图像降噪
由于能够在独立的子带中分离不同的方向，DTCWT在图像降噪等应用中通常优于标准的可分离DWT。使用软阈值处理对严格采样的 DWT 和 DTCWT 进行降噪比较。对于每个阈值，均显示均方根 (RMS) 误差和峰值信噪比 (PSNR)。
DTCWT 在 RMS 误差和 PSNR 方面优于标准 DWT。接下来，简单地看下图像的降噪效果
DTCWT 提供的 PSNR 比标准 2-D DWT 高出近 4 dB
3-D空间中的方向选择性
当将小波分析扩展到更高维时，使用可分离DWT在二维空间中观察到的伪影会加剧。DTCWT 能以最小的冗余保持3-D空间中的方向选择性。在3-D空间中，DTCWT中有 28 个小波子带。为了演示3-DDTCWT 的方向选择性，可视化3-D DTCWT 和可分离 DWT小波的3-D等值面图。首先，分别可视化两个双树子带的实部和虚部。
等值面图的红色部分表示小波从零开始的正偏移，而蓝色表示负偏移。可以清楚地看到DTCWT的实部和虚部在空间中的方向选择性。现在可视化双树子带之一，将实部和虚部图一起绘制为一个等值面。
上图展示了实部和虚部在空间中相互转换，这反映了双树复小波的虚部是实部的近似希尔伯特变换。接下来，在 3-D 中可视化正交小波的等值面以进行比较
看下DTCWT和DWT的方向混叠现象
与 2-D情况类似，3-D DTCWT的方向选择特性通常降噪效果更好。为了证明这一点，考虑一个由 16 个切片组成的 MRI 数据集。原始数据集中加入了标准差为10的高斯噪声。
load MRI3D
montage(reshape(noisyMRI,[128 128 1 16]),'DisplayRange',[])
计算一下SNR
20*log10(norm(origMRI(:),2)/norm(origMRI(:)-noisyMRI(:),2))
使用 DTCWT 和 DWT 将 MRI 数据集降噪，共分解到4级，这两种情况下都使用相似的小波滤波器长度。将生成的 SNR 绘制为阈值的函数。
figure
montage(reshape(imrecDTCWT,[128 128 1 16]),'DisplayRange',[])
title('DTCWT Denoised Volume')
figure
montage(reshape(imrecDWT,[128 128 1 16]),'DisplayRange',[])
title('DWT Denoised Volume')
