基于小波分析的虚拟消噪仪

一、引言

1807年，Fourier曾断言，任何一个周期为2的函数f(x)都是Fourier级数的和。绝大多数周期为2的连续函数满足Fourier的断言，但在1873年P.DuBois-Reymond构造了一个实变量x的2周期连续函数，它的Fourier级数在给定的点是发散的，这在当时的数学界引起了悍然大波，也引起了众多数学家的兴趣。在那时，数学界有三个未开发的领域：

（1） 修改函数的一般要求，并求得适合于Fourier级数的函数； （2） 修改Fourier级数收敛的定义；

（3） 求得另外的正交族，使得对于在三角族情形发生的P.DuBois-Reymond发散现象

在此时不至于发散。

这三个方向后来导致了重要的结果。其中，第三个方向导致了小波分析的产生。

自1986年以来，小波分析得到了迅速的发展，其应用也逐渐地变得越来越广泛，在某些领域中取得了一些超乎寻常的成果，特别是在噪声消除、特征信号的提取、图象处理等方面。 本文主要讨论小波分析在噪声消除方面的应用。 二、小波分析基础 1 Fourier变换

对于f(t)L空间的能量信号f(t)，它的Fourier变换定义为：

F()f(t)eitdt （8-1）

F()的逆Fourier变换定义为：

f(t)F()eitd （8-2）

实际应用中，为了计算Fourier变换及其逆Fourier变换，需要用数值积分，取f(t)在R上的离散点值来计算这个积分。设f(t)由采样得到，采样间隔为t，f(nt),n0,1,,N1为采样值，对应的离散Fourier变换和逆Fourier变换为：

F()N1kNN1k0f(nt)einti2nkN （8-3）

f(n)F(k)e （8-4）

2 短时Fourier变换

短时Fourier变换即时间信号加窗后的Fourier变换，其定义为：

wbF()eitf(t)w(tb)dt （8-5）

其中，w(t)为一个窗口函数，窗口函数w(t)的中心t*与半径w分别定义为：

t*1w2tw(t)2dt （8-6）

W1w2(tt*)w(t)dt21/2 （8-7）

这时，wbF()给出了时间信号在时间窗：

[t*bw,t*bw] （8-8）

的局部信息。

3 小波分析

f（t） 1 g（t） 0 g（t-b0） t 图8.1 窗口Fourier变换 如果把短时Fourier变换中的窗口函数W,b(t)替代为a,b(t)，其中：

a,b(t)a那么式（8-5）变为：

1/2tb （8-9） aWf(a,b)a该式即为小波变换定义式。

对应于（8-10）式，小波逆变换为：

1/2tbf(t)dt （8-10）

a1 f(t)C1tbWf(a,b)dbda （8-11） a2a比较式（8-5）与式（8-10），可以看到短时Fourier变换与小波变换之间的类似性，它

们都是函数f(t)与另一个具有两个指标函数族的内积。

对于(t)的一个典型的选择是：

(t)(1t2)exp(t2) （8-12）

2它是Gauss函数二阶导数，有时称这个函数为墨西哥帽函数，因为它象墨西哥帽的截面。墨西哥帽函数在时间与频率都有很好的局部化，函数图形如图8.2所示。

Ψ（t）

0 t

图8.2 墨西哥帽函数图形

短时Fourier变换与小波变换之间的不同，可由窗口函数的图形来说明，如图8.3所示。对于W,b，不管值的大小，具有同样的宽度。相比之下，a,b在高频（1/a相当于Fourier变换中的，a越大，频率越低）时很窄，低频时很宽的特性。因此，在很暂短的高频信号

上，小波变换能比窗口Fourier变换更好地进行“移近”观察。

a1a,b w(t) b0 0 x 0 t

Ψ（x） w,b 0 x

a,b 0 t

a1 b0 （a）

0 x

(b) 图8.3 (a)窗口Fourier变换函数W,b的形状 （b）小波a,b的形状

8.1.2 离散小波

m如果a，b都是离散值。这时，对于固定的伸缩步长a00，可选取aa0,mZ，

不失一般性，可假设a01或a01。在m0时，取固定的b0(b00)整数倍离散化b，选

mm取b0使(xnb0)覆盖整个实轴。选取aa0，其中m,n取遍整个整数域，而,bnb0a0a01，b00是固定的。于是，相应的离散小波函数族为：

m,n(t)am/20mxnb0a0m/2ma0a0xnb0 （8-13） ma0对应的离散小波变换系数为：

Cm,nf(t)0*m,n(t)dt

离散小波逆变换为：

f(t)CCm,nm,n(t)

其中，C为一常数。

8.1.3 小波级数

对应于Fourier级数的定义：

f(x)kF(kik0t )e0k0,1,, （8-14）

1其中，F(k0)TtTtf(x)eik0tdt。

同样可以定义小波级数：

f(x)j,kZZcj,kj,k(x)j,kZZdj,k~(x) （8-15） j,k~cj,kf,j,k其中，，称这两个无限级数为“小波级数”，并且是L2(R)收敛的。

dj,kf,j,k8.1.4 多分辨分析

1 多分辨分析的概念

如何由(x)L2(R)出发，使由k,n(x)张成L2(R)的闭子空间

VkclosL2(R)k,n(x):nZZ （8-16）

{(xn):nZZ}是V0的一个Riesz基，(x)称为尺度函数，这就是多分辨分析。

Vk，由于(x)V0V1，所以(x)可以用V1的基底设(x)生成一个多分辨分析1,n:nZZ表示。由于1,n:nZZ是V1的一个Riesz基，所以存在唯一l2序列pn，使

(x)npn(2xn) （8-17）

式（8-17）即为函数(x)的两尺度关系。系列pn称为两尺度序列。

对于模为1的复数z，引入如下记号

1 P(z)pnzn （8-18）

2n称为序列pn的符号。对式（8-17）两边作Fourier变换，则得到两尺度关系式

ˆ()P(z)ˆ()， zei/2 （8-19） 2同样地，由于(x)W0V1，所以存在唯一l2序列qn，使

(x)引入序列qn的符号

nqn(2xn) （8-20）

1 Q(z)qnzn （8-21）

2n对式（8-20）两边作Fourier变换，类似地得到

ˆ()，zei/2 （8-22） ˆ()Q(z) 22 分解算法与重构算法

由前所述可知，对于f(x)L(R)有唯一分解： f(x)2kgk(x)g1(x)g0(x)g1(x) （8-23）

其中，gk(x)Wk。令fk(x)Vk，则有：

fk(x)gk1(x)gk2(x) （8-24）

并且

fk(x)gk1(x)fk1(x) （8-25）

令

P(z)P(z)M(z) （8-26） Q(z)Q(z)在z1上，作函数

G(z)Q(z)P(z) ， H(z) （8-27）

detM(z)detM(z)则

MT(z)1 （8-28） H(z)H(z)对于符号G(z)，H(z)的序列{gn}，{hn}l1，存在如下的分解关系式 (2xl)若令an1{g2nl(xn)h2nl(xn)}，lZZ (8-29) 2nG(z)G(z)11gn，bnhn，则(8-29)式变成： 22n (2xl){al2n(xn)bl2n(xn)}，l0,1,2, (8-30)

为计算方便及以免产生混淆，有：

fk(x)jck,j(2kxj) （8-31） dk,j(2kj) （8-32）

gk(x)j在ck,j,dk,j中，k代表分解的“水平”。

2对于每个f(x)L(R)，固定NZZ，设fN是f在空间VN上的投影，有：

fNProjVNf （8-33） 可以把VN看作是“抽样空间”，而把fN看作f在VN上的“数据”（或者说测量采样值）。由于：

VN1WN1WN2WNMVNM （8-34） VNWN1所以，fN(x)有唯一分解

fN(x)gN1(x)gN2(x)gNMfNM （8-35）

对于固定的k，由{ck1,n}求{ck,n}，{dk,n}的算法称为分解算法。应用分解关系式(8-30)有

fk1(x)llck1,l(2k1xl)

nck1,l[{al2n(2kxn)bl2n(2kxn)}]{al2nck1,l}(2kn){bl2nck1,l}(2kxn)nlnl分解fk1(x)fk(x)gk(x)，得到：

{ck,nal2nck1,l}(2kn){dk,nbl2ndk1,l}(2kn)0

nlnl所以，由{k,n:nZZ}，{ k,n:nZZ}的线性无关性以及VkWk{0}，得到分解算法：

ck,nal2nck1,ll （8-36） dbck,nl2nk1,ll分解过程如图8.3所示。

cNdN1 dN2 dNM

cN1 cN2cNM

图8.3 分解过程

在实际计算中，假定取值点所对应的f(x)的水平为N，即

f(x)fN(x)

对于某个正数M(0MN)，信号由N水平分解到NM水平，即已知{cN,n}，求{dk,n}，

kN1,,NM及{cNM,n}。

同样地，固定k，由{ck,n}，{dk,n}求{ck1,n}的算法称为重构算法。应用两尺度关系有

fk(x)gk(x)ck,l(2kxl)dk,l(2kxl)llck,lpn(2k1x2ln)dk,lqn(2k1x2ln)lln(ck,lpn2ldk,lqn2l)(2nnllk1nxn)

{(pn2lck,lqn2ldk,l)}(2k1xn)因为fk(x)gk(x)fk1(x)，有fk1ck1(2k1xn)及{k1,n:nZZ}的线性无关性，得

n到重构算法：

ck1,n(pn2lck,lqn2ldk,l) （8-37）

l重构过程如图8.4所示。

dN1 dN2 dNM

cNcN1 cN2cNM

图8.4 重构过程

三、Matlab工具箱中小波分析函数

在一维信号的分析中，Matlab工具箱提供了许多小波分析功能函数。

四、基于小波分析的虚拟消噪仪

在信号的采集、传输与处理过程中，由于外界或电路内部因素干扰，使得信号被噪声污染，所处理的信号中夹杂着一些无用的噪声信号，通过各类电路或算法将噪声从所处理的信号中消除，称为信号的消噪，应用小波进行信号消噪处理是小波分析的一个重要应用之一。

4．1小波消噪原理

小波变换在信号消噪中的思想同傅立叶变换滤波思想相似，只不过傅立叶变换的数字滤波是等步长频谱滤波，而小波变换消噪则是二等分频谱滤波，只有进行小波包分解才能实现等步长频谱滤波。由于变换的基波不一样，经典的滤波效果和小波消噪的效果也不一样，在小波消噪处理中，选用的小波不同，消噪效果也不一样。

应用小波分析进行消噪主要涉及到小波的分解与重构，下面以一维信号为例来介绍小波消噪的原理。

含有噪声的一维信号可以表示成如下的形式：

s(i)f(i)e(i) ，i0,1,2,,n1

其中f(i)为真实信号；e(i)为高斯白噪声，噪声级为1；s(i)为含噪声的信号。

对信号s(i)进行消噪的目的就是要抑制信号中的噪声部分，从而在s(i)中恢复出真实信号f(i)。在实际工程中，有用信号通常表现为低频信号或是一些比较平稳的信号，而噪声信号则通常表现为高频信号。一般来说，一维信号的消噪算法可以分为三个步骤进行： 1) 对信号进行小波分解

选择一个小波并确定一个小波分解的层次N，然后对信号s进行N层小波分解。分解过程如图8.9a，图8.9b是对应于图8.9a的频谱图，分解算法见式8-36。 2)小波分解高频系数的阈值量化

对第一层到第N层的每一层高频系数选择一个阈值进行软阈值量化处理。噪声部分通常包含在CD1,CD2,CD3中，因而，可以以门限阈值等形式对小波系数进行处理，也就是令CD1,CD2,CD3部分的为零，或者适当地减小。 3）对信号进行重构

根据小波分解的第N层的低频系数和经过量化处理后的第一层到第N层的高频系数，小波系数进行处理完毕后，进行一维信号小波的重构，即可以达到消噪的目的，如图8.4所示。重构算法见公式8-37。

S CA1 CD1 CA2 CD2 CA3 CD3 （a） CA1 CA2 CA3 CD3 8CD2 4CD1 2 ω（b）

图8.9 小波消噪原理的分解示意图

在这三个步骤中，最关键的一步就是如何选取阈值和如何进行阈值的量化，从某种程度上来说，它直接关系到信号消噪的质量。

在这三个步骤中，最关键的一步就是如何选取阈值和如何进行阈值的量化，从某种程度上来说，它直接关系到信号消噪的质量。

8.3.2 虚拟小波消噪仪设计

虚拟小波消噪仪的设计包括仪器前面板设计、与Matlab接口设计及数据流程设计三部分的内容的设计。

1 虚拟小波消噪仪原理

虚拟小波消噪仪的设计思路为LabVIEW通过数据采集或仿真生成含有噪声的信号，通过仪器前面板设置消噪处理的参数，将参数通过LabVIEW与Matlab接口传递给Matlab相应的功能函数，完成信号分析与处理功能，最后将处理结果回传给LabVIEW进行显示，虚拟小波消噪仪的原理框图如图8.10所示。

传递参数 传递参数 LabVIEW LabVIEW Matlab系

接口 计算后参数 统 传递参数包括数组， 回传参数 与 matlab 相关参数 函数、功能

图8.10 虚拟小波消噪仪的原理框图

2 虚拟小波消噪仪设计要求与功能

设计如图8.11所示的仪器前面板，点击“信号频率”按钮，选择仿真正弦信号的频率；点击“采样频率”按钮，选择仿真正弦信号的采样频率；点击“消噪层数”按钮，选择小波消噪的层数；“含噪声波形”显示屏中显示含白噪声的正弦信号；点击“动态波形”按钮，动态显示波形。

程序运行后，调用Matlab工具箱中小波消噪函数WDEN（）对含白噪声的正弦信号进

行消噪处理，并将消噪后的正弦信号在“消噪后波形”显示屏中显示出来。

图8.11 虚拟小波消噪仪的仪器面板

图8.12 仪器面板及控件

5 运行虚拟仪器

运行虚拟仪器VI文件，运行结果如图8.14所示。

图8.14 运行结果