多相分解法实现增采样

1、多相分解法实现增采% M120401% 演示用多相分解法实现增采样，说明使用窗函数对冲激响应序列进行截取的必要性。 % figure(1); set(gcf, 'color','w')for i = 1 : 13% 求理想冲激响应序列dd(i) = sin(pi * (i-1+0.000001) / 4) / (pi * (i-1+0.000001) / 4);end% 对冲激响应序列不加窗d = fliplr(dd(2:13), dd;for i = 1 : length(d)h(i) = d(i);end hidea = h;% 设置插值节点xup = 4,

2、0,0,0,6.6,0,0,0,5,0,0,0,6,0,0,0,6,0,0,0,7,0,0,0,8,0,0,0,7,0,0,0,7,0,0,0,7,0,0,0,6,0,0,0,. 6,0,0,0,5,0,0,0,5,0,0,0,4,0,0,0,4,0,0,0,6,0,0,0,7,0,0,0,8,0,0,0,6,0,0,0,5,0,0,0,5,0,0,0,6,0,0,0;subplot(311)grid; axis(0,70,0,10);H = stem(1:length(xup),xup); set(H,'markersize',2); xlabel('n',

3、'fontsize', 12); ylabel('x(n)', 'fontsize', 12); title(' (1)在给定序列的样点之间添加零值样点');% i 是行序号， % j 是列序号，代表子滤波器序号由和 j 求出冲激响应序列样点序号% 求变换矩阵的 4 个行向量 for i = 1 : 4 for j = 1 : 6：');for N = 4 : 15for i = 1 : 6hh(i, j) = h(-4 + i + 4 * j); end%-0.00000.0000-0.00001.0000-0.0000

4、0.0000%0.0818-0.12860.30010.9003-0.18010.1000%0.1273-0.21220.63660.6366-0.21220.1273%0.1000-0.18010.90030.3001-0.12860.0818% hh =取出 6 个节点，得到列向量enddis p('不加窗时的子滤波器冲激响应序列hhtemp2(i) = xup(4 * N + 9 - 4 * i); end % 每组插值节点确定 for j = 1 : 4for i = 1 : 64 个输出样点% j 代表子滤波器序号 ( 插值点序号 )% 取出 1 个行向量 ( 第 j 个子滤

5、波器的冲激响应序列 ) temp1(i) = hh(j, i);end% 求行向量与列向量的点积，得到第 j 个插值点a = 0;for k = 1 : 6a = a + temp1(k) * temp2(k);endxup(4 * N -8 + j) = a; endendsubplot(312)H = stem(1:length(xup),xup);set(H, 'markersize', 3);grid; axis(0,70,0,10);xlabel('n', 'fontsize', 12); ylabel('xup1(n)'

6、;, 'fontsize', 12);title(' ( 2)插值滤波后所得的增采样序列( 使用矩形窗 ) ');% % 使用 Hanning 窗w = hamming(25); for i = 1 : 25 h(i) = d(i) .* w(i);end hwindow = h % 设置插值节点% 对冲激响应序列加窗% i 是行序号，代表子滤波器序号% j 是列序号，由 i 和 j 求出冲激响应序列样点序号xup = 4,0,0,0,6.6,0,0,0,5,0,0,0,6,0,0,0,6,0,0,0,7,0,0,0,8,0,0,0,7,0,0,0,7,0,0,

7、0,7,0,0,0,6,0,0,0,. 6,0,0,0,5,0,0,0,5,0,0,0,4,0,0,0,4,0,0,0,6,0,0,0,7,0,0,0,8,0,0,0,6,0,0,0,5,0,0,0,5,0,0,0,6,0,0,0;4 个行向量% 求变换矩阵的 for i = 1 : 4 for j = 1 : 6 hw(i, j)= h(-4 + i + 4 * j);窗时的子滤波器冲激响应序列：');end end disp(' 加 Hanning hw% hw =%-0.00000.0000-0.00001.0000-0.00000.0000%0.0078-0.05410

8、.25970.8862-0.11870.0215%0.0180-0.11460.59740.5974-0.11460.0180%0.0215-0.11870.88620.2597-0.05410.0078for N =4 : 15for i = 1 : 6%取出 6 个节点，得到列向量temp2(i) = xup(4 * N + 9 - 4 * i); end % 每组插值节点确定 for j = 1 : 4for i = 1 : 64 个输出样点% j 代表子滤波器序号 ( 插值点序号% 取出 1 个行向量 ( 第 j 个子滤波器的冲激响应序列 ) temp1(i) = hw(j, i);e

9、nd% 求行向量与列向量的点积，得到第 j 个插值点a = 0;for k = 1 : 6a = a + temp1(k) * temp2(k); endxup(4 * N -8 + j) = a; endendsubplot(313)H = stem(1:length(xup),xup);set(H, 'markersize', 3);grid; axis(0,70,0,10);xlabel('n', 'fontsize', 12);ylabel('xup1(n)', 'fontsize', 12);title(

10、' ( 3)插值滤波后所得的增采样序列( 使用 Hanning 窗 ) ');% figure(2);set(gcf, 'color','w')H,w = freqz(hidea); subplot(211) plot(w,abs(H); grid;xlabel('w (rad. / sec)', 'fontsize', 12); title(' (4)实际低通滤波器的模频特性 H,w = freqz(hwindow);subplot(212) plot(w,abs(H); grid;xlabel('

11、;w (rad. / sec)', 'fontsize', 12);ylabel('|H(jw)|', 'fontsize', 12); ( 使用矩形窗 ) ');ylabel('|H(jw)|', 'fontsize', 12);卄*卄*卄*卄*卄*卄*卄*卄*卄*卄*卄*卄*卄*卄*卄*卄*卄*卄*卄*卄*卄*卄*卄*卄*读者若要看懂以上 MATLAB程序，请先参阅下面的文章，谢谢。米样率转换12.1 概述在DSP的许多实际应用中，系统往往由多个工作于不同采样率的子系统构成。例如，在数字电话、电

12、传打字、电报、 传真、语音、视频等电信系统中，需要用与带宽相适应的不同速率处理各种信号，使待处理信号既符合采样定理又可以减少 数据量。在声音光盘(CD )或数字音频磁带(DAT)播放器中，信号在进入 D/A变换器之前，为了避免使用技术指标很高 的后置模拟低通滤波器，需要提高采样率。种种客观需要使我们面临改变信号采样率的问题。DSP中，采用多种采样率的系统被称为多将信号从某个采样率转换到另一个采样率，这样的过程称为采样率转换。在 速率数字信号处理系统。采样率转换中，使采样率提高者，称为内插(interpolation)或增采样(up-sa mp ling)。使采样率降低者，称为抽取(decima

13、tion) 或减采样(down-sampling)。采样率转换可以用两种方式实现。第一种是通过 D/A转换器，把原有的数字信号转变为模拟信号，然后重新采样。第 二种方法是完全在数字域实现采样率转换。第一种方法能实现任意比率的采样率转换。但它的主要缺点是，在转变为模拟信号时，D/A转换器会引入信号失真，而在重新采样时，A/D转换器会引入由量化效应造成的信号失真。此外，在许多实际应用中，不可能采用这种方式。例如，在 CD播放器中，采样率为 44.1KHZ的数字音乐信号已存放在 CD中。在将这数字信号恢复为模拟信号以进行播放时，为了得 到高保真度，需要提高采样率。但曲终人散，人们不可能得到原声频信号

14、，以更高的采样率对其重新采样，从而不得不在数 字域实现采样率转换。第二种方法在数字域完成采样率转换的全部过程，避免了第一种方法的缺点。本章介绍在数字域中直接实现采样率转换的原理和方法。首先阐述增采样，然后介绍减采样。无论增采样或减采样，都 需要使用FIR数字滤波器对原序列进行低通滤波。为此，本章详细地介绍多相滤波技术。最后，介绍基于采样率转换的一种 被广泛应用的新技术-噪声整形技术。倖尽管可以直接在数字域实现采样率转换，但是，在研究采样率转换过程的频谱变化时，可以从采样定理出发，设想从原序列的频谱先还原出模拟信号的频谱，然后再按新采样率求出转换后的频谱。当然可以直接导出增采样或减采样后的信号频

15、谱，但从采样定理出发来理解频谱变化是非常容易的。为了叙述方便，这里重新表述采样定理如下。假设模拟信号xa(t)的频谱为Xa(f)，采样间隔为Ts (即采样频率为fs=1/Ts )，采样所得序列x(n)的频谱为X ( f)，则按照采样定理，有从幅值来看，X ( f ) = X a( f ) / TsX ( f )变成周期函数，周期为 fs =1 /Ts希望读者在学完本章后,熟悉增采样和减采样导致的频谱变化，理解低通数字滤波器的作用。 掌握多相分解技术。掌握采样率变换器的基本设计方法：直接法和多相分解法。12.2内插器的频谱与内插滤波器内插的作用是提高采样率，又称增采样，其倍数记为整数U (代表u

16、p)。图12.2.1以U =4为例说明内插过程。如图 所示，增采样器的输入序列x(n)工作于低采样率(简称低率)fs =1/Ts。其中，Ts为采样间隔。增采样器的作用是在x(n)的两个样点之间插入 U -1个零样点，使采样间隔变为 Ts'=Ts/U。所得序列记为xup(n')，它工作于高采样率(简称高率)fs'=Ufs。居于两个低率样点之间的 U -1个零值样点最终应取恰当值。最简单的办法是使样点之间的增量是常数,FIR内插滤波器，进行低率样点之间的插值。如图12.2.2所示。这种线性插值往往不能满足要求。所以增采样器一般要后接 后面将要说明，插值过程是数字低通滤波过程

17、。丿 upM)012043起'04 a 3如1零値内插軽知(那)FIR1Flip (幷)k-増采样ttz内插浦波器rJL1r p _/増采样因千U=A+1 -A1 -;空-7 >A*= A/4低率样点O插值样点图12.2.1增采样原理图图 12.2.2过于粗糙的插值方法其实，这里的内插滤波器与重建模拟信号所需的低通滤波器(见第有类似的功能。所不同者，前者是低通数字滤波器，在插值区间仅需要插入 个样点，所完成的是数字信号之间的转换；而后者则是低通模拟滤波器，理需要在插值区间插入“无限”多个值，以完成数-模转换。2.11max /rwc节)具 有 论图12.2.3示出模拟信号xa (

18、t)的频谱模值I X a( f ) I，其幅度设为 为基础，图12.2.4说明增米样前后的频谱变化。与频谱对应的图12.2.3模拟信号的频谱时域序列示于各子图的左边。首先，以fs =1/Ts对模拟信号xa(t)进行采样。根据采样定理，采样所得序列x(n)的频谱应如子图A所示。以因子U = 4进行增采样时，可以分两步走。第一步，设想在x(n)的每两个样点之间添加 U -1=3个零值样点，得序列x'up (n')，其频谱示于子图Box'up(n')每两个样点之间的时间间隔是Ts=Ts/4，即fs' = 1/Ts' = 4fs。容易证明x'up

19、(n')和x(n)的频谱完全相同，程序M120201验证了这一点【解说12.2.1】。第二步，在各个零值点插值。这两步合起来就相当于序列 x(n)先被转换成模拟信号 xa(t)，得频谱模值|Xa(f)|，然后以新采样率fs'=4fs对其进行采样，得到高率序列Xu p(n')，其频谱模值| Xu P ( f) 1示于子图C。读者运行程序 M120202，即可得到验证【解说 12.2.2】。比较子图B和C,可见为了在序列x'up (n')的零值样点处进行插值，得到序列xup (n')，需要滤去子图 B中用圆圈标出的镜像频谱块。为此，要使用低通滤波器对

20、序列XuP '(n')进行滤波。该滤波器的模频特性示于子图DoIJTCOIX河it'AAAAAAAA；,-A r4k4k+$4Jt J4+1 -Lt + SUZE±-A'-Q20 A'/S M!1旧0 -/LT-防L理撩删滤液器的戟止坳率D)菸厲£茂篠的频華鱼图12.2.4增采样前后的信号频谱及插值滤波器的作用3以上论述表明，增采样可以理解为插值过程，也可以理解为滤波过程。从时域看，增采样是在低率序列的样点之(anti-image filter)。增采样问题就变成怎样设计内插滤波器的问题了。间进行插值；从频域看，则是用低通滤波器滤掉镜

21、像频谱块。 镜像滤波器这样，所以这个滤波器称为内插滤波器(interpolation filter )，又称抗由子图A、D可见，滤波器的截止频率应为 fc = fs / 2。这里，fs是原来的采样率。从子图 D又可以看出，新采样率fs'相应的数字频率为2兀，故滤波器的数字截止频率为(12.3.1)在实践中，讲座10所介绍的级联积分器-梳状滤波器(cascade integrator -comb filters，CICF)非常适合于在进行插值 操作(提高采样率)之后，用作去像频滤波器。内插滤波器的实现方法有二，即直接法和多相分解法。下面在介绍直接法之后，将重点阐述多相分解法。直接法12.

22、3内插滤波器实现方法之一为了便于阅读，今将图12.2.1重新绘出。0 1 2 «(A) «人序列4 S H'"*（B）零H内插器鸯出序?1 高率a仁 増采样因子U=A儿/丹）Tfirrn曲十，04' B 刚插值浦£春出序列高率fZm图 12.3.1 U倍增采样系统的序列波形图12.3.1C中，实心样点表示来自输入序列的低率样点，空心样点表示零值内插操作产生的样点。零值内插器在相邻两 个低率样点等间隔地插入 U -1个零值样点（子图B）插值滤波器根据插值区间前、后各 M个低率样点（节点）求出各个 插值样点值。节点总数 2M的值越大，插值精度

23、就越高。图12.3.1中，n和n分另U作为低率和高率序列的序号。显然可见,Xup (nU)=x (n)(12.3.2)在序列 xup (n )中，xup (nU)与 Xup (nU+ U ）之间的U 1隔样点为零，即Xup (nU+i) =0,i =1,2,U 1(12.3.3)合并以上两式，得Xup(n') ix(n).n' = nU0,其它(12.3.4)和余数。由于Ts序号n'总可以写成n'= nU +i，其中，i的取值仅限于i =1,2,., U -1。换句话说，n和i分别是n'除以U所得的商 = UTs'，故若以秒来计量时间，贝U n所

24、相应的绝对时间为=n'Ts' = nUT s'+iTs' = nTs +iTs'i的最大可取值为U 1。因此，为了满足实时处理要求，两个低率样点之间的U -1个插值必须在（U -1）Ts'内完成。(12.3.5)由图12.2.4可以看出，插值滤波器是一个低通滤波器，其截止频率为fc = fs / 2 = fs'/(2U )若用数字频率©' = 2兀（f / fs'）来表示，则截止频率为国=2兀(fs/2) / fs'(12.3.6)子图C示出理想低通数字滤波器的幅频特性。图中,选倍增率U =4。滤波器的理

25、想冲激响应为这是一个用处j f.'i'md/m) “ H©')e心兀/U2兀U/UUe2兀sin(兀 m / U )d/m)=-3C < m(12.3.7)兀m / Usin Xsinc(X)=函数表出的无限长的非因果对称序列:X对称点在m =0 ,d/0) =1若选U =4，则在对称点 m =0的两侧，每经过U -1=3个样点，函数d处(m)有一个过零点:m = ±Uk (U、k为整数)(12.3.8)对d比(m)进行对称截取，使之在 m =0两侧各有M个过零点，得序列d (m):d/m),d (m) = «0,I m |WMU其它

26、(12.3.9)d (m)的长度为N =2MU +1(12.3.10)图12.3.2示出U =4,M =3时的序列d(m)。-U5151111II!1191,I1111II、1003>0 6366111 1111.¥ f.I"'0300'1 0.00】o'. 12730F f 申 0.03IS'i" A !厂i-0 10011 i k丘肪！cTo: : : ：-0.2122 :1 1 1 1 1O'010-1015ifi图12.3.2从理想冲激响应序列截取所得的对称序列d (m)对 Xup (n')式(12.3

27、.4)与 d(m)式(12.3.9)进行卷积，得UMyup(n') = Z d(nm)Xup(m)m =_UM(12.3.11)这样就得到图12.2.1的插值输出yup(n')。从式(12.3.11)来看，图12.2.1所示的整个系统来说，系统是线性、时变的【解说FIR插值滤波器运行于高率上，其滤波操作是线性、时不变的。120301】。但就图12.2.1上述直接算法的效率非常低，因为序列Xup(n')中，每两个低率样点之间有U 1个零值样点，而按式(12.3.11 )又不得不对大量零值样点进行无谓的乘法运算。为了克服这个缺点，需要采用下节所讲的多相分解法进行滤波。12.

28、4内插滤波器实现方法之二多相分解法从卷积概念和插值概念都可以导出多相分解法。本节将介绍第一种方法。12.4.1从卷积算法导出多相分解法下面用图12.4.1说明通过卷积实现插值的过程。(A)池囚：1 :1;約1: 1 :1:約 1: 1:1 : 1 1 1 : 1 : 1 : : 1：1 ；榕0勺：1 ： 1 ：思凹1:11 W9 S654 S2 10 1；34 ?6" S 9101112用冋憑帀：1:1 :沆 1 : 1:1 :诊悶：1 ： 1 ： 1 J 1 ：1 :1 :诊a环:丨：丨：緻勺：1:1 :帥习12 L110S -16413： Jf 0113t 41?6- S 910

29、 11： 12 m1*1V1V1*1t1 t1 11 11 11 11 11 1k-1kk+lk+2u程1111插值节克qff1 1 1 1 1 1 11待插ffiJ 111 j 11114k 4fc + -l4 + L2 卅4倉+ 1图 12.4.1通过低通滤波进行插值假设u=4, M 二=3，求出了图12.3.2所示的序列d(m) O为了简单起见,将该序列画成图12.4.1A的样子。A每个格子代表一个样点的对称点。在该处，d(m) =1。带斜线的格子表示 d(m)的过零点,符号“”标出序列d (m)即 d (m) =0。以":”标出的格子表示具有非零值的d( m)样点。将序列右移

30、一步，仍记为 d (m)，得子图B。子图C、D分别表示低率序列 x(n)和高率序列Xup (n')，带底纹的格子是有定义的样点，称为节点。已设U =4 , 故子图D中，每两个低率样点之间有 3个待插值的点第2章已经指出，卷积过程中，被卷积序列之一（X1 ）保持不动，而另一个序列（X2 ）在被反转之后，作前向或后 向移动。每次移动时，两序列对应点相乘并相加（做点积运算）。这样就得到卷积序列 X-12的一个样点。子图A表示序列d （m）的对称点已移到n =k （即n'=4k ）处。这里，为了简单起见，只对序列d （m）标出相对于对称点的下标。 而且，由于d （m）是对称序列，故负值

31、下标都改成正值（例如，将m = _2记为m = 2 ）。对子图A与D所示的两个序列求点积，得卷积序列的一个点：Xu p (4k +0) =Xup (4k +12) * d(12) +Xup (4 k +8)* d (8) +Xu p (4 k + 4) * d (4) + Xu p (4k)* d (0)+ xup (4k -4)* d(4) +xup (4k 8)* d(8)(12.4.1) 由于 d(0) =1，而 d(12)、d(8)和 d (4)均为 0,故Xup (4k) =x(4k)即在节点处，插值点就是节点本身，不受其它节点的影响。子图B表示序列d（m）右移一步。对子图 B与D所

32、示的两个序列求点积，得xup (4k+1) =Xup(4k +12)* d (1 1) + Xup (4k +8)* d(7) + Xup(4k +4)* d (3) + Xup(4k)* d (1)(12.4.2) + Xup(4k 4)* d(5) +Xup(4k -8)* d(9)同理可得Xu p(4k+ 2) =Xu p(4k +12)* d(10) +Xu p(4k +8)* d(6) +Xu p (4 k + 4) * d (2) + Xu p(4k)* d(2)+ Xup(4k 4)* d(6) +Xup(4k 8)* d(10)(12.4.3)Xup (4k+ 3) =Xup

33、(4k +12) *d(9) +xup(4k +8) * d(5) +xup(4k +4) * d(1) +xup(4k)* d(3)+ xup(4k 4)* d(7) +xup (4 k 8)* d (11)(12.3.4)由于 Xup (4 k +12) =x(k +3),，Xup (4k +8) =x(k +2)，Xup (4k + 4) = x(k +1)，x. p(4k) =x(k)， xu p(4k-4) =x(k 1)， xup (4k-8) =x(k 2)，故综合以上各式得x(k +3)Xup (4 k +0)d(12)d(8)d(4) d (0) d (4) d (8)

34、9;x(k +2)Xup (4 k +1)d(11) d (7)d(3) d (1) d (5) d (9)x(k +1)11t»Xup (4 k +2)d(10)d(6)d(2) d (2) d (6) d (10)x(k)Xup (4 k +3)d(9) d (5)d(1) d (3) d (7) d (11 )- *x(k -1)x(k 2)(12.4.5)由此可见，可以用 4个子滤波器通过实现插值，得到 Xup (4k +i), i= 0,1,2,3。这种算法称为多相分解法。式中,Xup (4k +0)、Xup(4k+1)、Xup(4k+2)和 Xup (4k +3)分别是子

35、滤波器 subfilter0、subfilter 1、subfilter 2和subfilter 3的输出。矩阵的4个行向量是相应子滤波器的冲激响应序列。每个子滤波器负责算出一个插值点。由于d (0) =1，d (12) =d (8) =d (4) =0，故在插值区间内，第 0个插值点就是原有的节点，其值与其它节点无关,即 Xup (4k) = x(k)。(12.4.5) 表明，输入序列是低率序列，输出插值序列则是高率序列。在多相分解法中，全部插值运算都在低率下完 成。而且可以采用并联运行的硬件同时计算各个插值，从而大大地提高效率。图12.4.2表示多相滤波器的工作情况。x(r + l)X菽)

36、低率M+1)+ -)即可根据这个规律求出各个子滤波器的冲激响应【解说12.4.1】。上面以U 子滤波器数目) 从式(12.4.5)12.4.5)。式中的矩阵的列数为 2M ，行数为U (即 )。12.4.2用多相滤波器实现插值=4和M =3为例，说明多相分解法的原理，导出式(，每一行是相应子滤波器的冲激响应序列(长度为2M可以看出矩阵元素下标的变化有图12.4.3所示的规律。图中，符号代表矩阵元素。如果给出U和M ,D 列 112CTxJW =4x3 = 41II d- 12丄 S丄0 行一Hir<Hr#2行！*；!*4丄斗 A/WlM 7图12.4.3式(12.4.5 )的矩阵元素下标

37、的变化对于增采样系统，必须注意两点：首先，第12.3节已经指出，增采样滤波器本身运行于高率，如式(12.3.11)所示，滤波操作是线性时不变的。但是,若采用分相技术，用 U个子滤波器实现插值，则系统输入、输出端运行于不同的采样率，从式(12.4.5)可见，整个系统的滤波操作虽是线性的， 但却是时变的，即：子滤波器冲激响应向量di (n); n= 0,1,，U-1的元素既与样点序号 n有关,又特别与子滤波器序号i (即插值点序号)有关。一般来说，输入序列x(n)与系统冲激响应序列卷积，求得第 i点的滤波值，如果冲激响应仅仅是样点序号n的函数，而与i无关，则系统是时不变的，否则是时变的【解说 12

38、0301】。其次，从式(12.4.5)可见，整个系统的滤波操作是非因果的，因为在某个低率样点x(k)后面插入高率样点时，除了使用先前的低率样点 x(k -2)、x(k -1)和x(k)夕卜，还需要使用“将来的”低率样点x(k +1)、x(k +2)和x(k +3)。12.2节已经指出，利用低通滤波器在频域中滤去镜像频谱块，就相当于在时域中实现插值(见图12.2.4)。是据此导出的。从时域插值概念也能导出同样的结果【解说12.4.2】。.这说明插值过程(在时域)相当于滤波过程式( 12.4.5)(在频域)。12.4.2多相分解滤波器演示程序M120401.M对多相分解法进行验证。部分运行结果示1

39、41r于图 12.4.4。程序首先生成一个高率序列Xup0 (n'):在原低率序列每两个样I .IIII I*1gji19»44» U 74«刘a低宰祥血之问砂辛恤轻ja点之间插入3个零值样点(子图A)其次，按式(12.3.9)求出序列d (m)，子滤波器的冲激响应序列，按式(12.4.8)求插值，得到子图 B1D|1Ml并计算各个La序列Xup1 (n')。但插值结果很不理想，包络线不平滑，有明显n' '<1ID»30ID50mro(St 低CTJdlW)输岀1D|_:RRRRRR1所示的插值的咼频分Bir%14a

40、 »44« UI门 心酒社器(加it) »山量。程序进一步作了改进：对序列d/m）加Hamming窗，得到冲激响应序列h（n'）。这时，式（12.4.8）中的矩阵行向量方法求出。于是，各子滤波器的冲激响应向量百i ;i =0,1,2,3可以表述如下:元素d （n'）应置换成h（n'）样点，并按此式进行插值，得到序图12.4.4用多相分解法实现增采样列 Xup2 (n')(子图 C)。由程序M110401.M的运行结果可以得出结论,简单地对理想低通滤波器的冲激响应序列d比（n'）进行截取，以此作为插值滤波器的冲激响应,不能得到

41、良好的插值效果。在实际应用中，应该采用适当的窗函数 （例如，Hamming窗）对d（n'）加窗。程序 M120401.M模频特性得到很大改善。还比较了简单截取与使用这是插值序列变得十分平滑的根本原因。Hamming窗所得滤波器的模频特性。结果示于图12.4.5。显然，加窗后,旧9)1图12.4.5低通滤波器的模频特性12.5内插滤波器的结构系统函数在阐述内插滤波器的结构之前，先说明多相滤波器的系统函数。为简单起见，取U =4 （即每两个低率样点之间插入3个样点，使采样率增为 4倍）和M =2 （在插值区间两侧各使用 2个低率样点）。在以下讨论中，假定对理想低通滤波器 的冲激响应序列加窗

42、。相应地，本节采用字符“ h ”以代替上节所用的字符“U =4，所以在多相分解的情况下，有4个子滤波器。又由于选12.5.1由于取d ”，表示冲激响应序列。M =2，故每个子滤波器的冲激响应序列长度为2M。各子滤波器的冲激响应用向量hi表示，下标i表示子滤波器序号,i =0,1,2,3。冲激响应样点用 hij表示,其中，下标i代表所属子滤波器的序号，j代表该样点在向量中的序号。系数组hjj ; i = 0 3, j = 0 3可按上节所讲=d 0 =h00 , h01 , h02 , h03 0, 0, 1, 0=d 1 =hio , hii, 02 , hi3 = -0.13, 0.30,

43、0.90, -0.18(12.5.1)=d 2 =h20,h21,h22,h23 =-0.21, 0.64, 0.64, 0.21=d3 = hgo ,山31,山32 ,山33 = 0.18, 0.90, 0.30,-0.13这样，在u = 4, M =2场合，可以用图12.5.1所示的4个子滤波器分别求出插值样点y 0 (n), yj n), y2( n)和y3 (n)。其中，y0( n)是原来的低率样点 x(n)，而y1( n), y2( n)和 氷n)是x(n)和x(n+1)之间的3个高率插值样点。图中,zin亠是低率延迟算子，延迟时间是低率采样间隔Ts =1/ fs。由于4个子滤波器是

44、并联的，故整个滤波器的系统函数为1231231H (z) =(hoo +ho1zi- +h02 zin + h03zin ) +00 +山1需 + h13zin ) zout1232_1+ (h20 +h21Zin +h22 zinh23 zin ) zout +(h3h31zin233h32 zinh33 zin )zou t(12.5.2)考虑到高率延迟算子zou1 = zj/4 (因为现在U =4，故相应的延迟时间等于Ts/4)，得H(Z)=(帕0 卄01 zout 卄02Zout 中 h03zout)+(h10 zout + h11 zout 十 02 zout 中 h13 zout

45、)261014371115+ (h20Zout+h21Zout+h22Zout+h23Zout )十(h30 Zout中 h31Zout + h32Zout+ h33Zout)15h(k)z 止k £(12.5.3)这是16抽头FIR滤波器的标准系统函数。式(12.5.2)被称为式(12.5.3)的多相分解(poly phase deco mp osition)12.5.2 结构由式(12.5.2.)可以得到U =4和M =2时的增采样多相结构如图12.5.1所示。低率输入工如%0I切换器T为0):訂卜.+>aut延旳g h TT丹30+十子謔液器延时27； /4子池破器3图1

46、2.5.1增采样的多相结构图12.5.1右端的相加器相当于一个切换器,由于切换，依次取出4个子滤波器的输出。延时H174-3二_L/ 4ZOUt = ZinAfi)，所以在输入序列的每个样点间隔内，该切换器进行测)12.5.2滤波器子采用独立延迟线12.5.3具有最少存储器的增采样分相结构增采样分相结构多相结构的形式多种多样。图12.5.1的各个子滤波器采用公共延迟线。 图12.5.2所示的子滤波器则采用独立的延迟线。将这种结构加以改变，得图12.5.3。它是一种具有最少存储器的结构。输入一个新的x(n)样点后，三个切换器同步地反时针旋转，切换4组滤波器系数，分别计算 4个y(n')样

47、点。12.5.3多相分解法的优点并行算法。其次，12313)计算插基于多相分解，整个滤波器被分解为U个并联的多相子滤波器，这样就能在软、硬件实现中，采用U分之一。多相分解法的运算效率比直接法高得多。因为作用于各个子滤波器的信号是低率样点，而用直接法按照式(值，所用的信号则是高率样点xuP (n')。容易证明，多相分解法的计算量仅为直接法的12.6抽取器的频谱与抽取滤波12.6.1抽取器的频谱与抽取滤波本节讨论按整数因子 D对信号x(n)进行减采样。减采样是增采样的逆过程。假定模拟信号xa(t)的频谱Xa( f)局限于频率区间一 fmax < f < fmax内(图12.6.

48、1)，其幅度为1。若采样率为fs，则所得序列xo (n)的频谱是Xa(f )的周期延拓,图 12.6.1模拟信号的频谱周期为fs ,而幅度为1/Ts = fs (图12.6.2A )。这里，fs >2fmax如果另以采样率fs = fs/ D对xa(t)进行采样,则序列x(n)的频谱的重复周期为fs (子图B),频谱幅度变为1/(DTs)。x(n)是从x。(n)每隔D个样点抽取一个而得到的。由于采样率降为fs =fs/D，故折叠频率变为 fs /(2D)。所以，如果模拟信号的最高频率仍为fmax，则有可能发生频谱混叠。子图 B示出频谱混叠情况。因此，在直接对序列x(n)进行减采样时，为了

49、避免混叠，必须采用低通数字滤波器，将x(n)的频谱截止频率减小到fmax兰f s /( 2d )，然后，才可以按D因子进行减采样而不致混叠。为此，这个滤波器又称为抗混叠滤波器(anti-aliasing filters )。子图C、D分别示出抗混叠滤波器的频率响应和无混叠的减采样器的频谱。假定模拟信号频谱幅值为1，如图12.6.1所示，那么，从子图 C、D，可以看出减采样前后的频谱幅度分别是1/Ts和s1 /(DTs)。子图C添加了数字频率标尺。由图可以看出，理想低通滤波器的截止频率应小于(12.6.1)图中，D =2。H；釆徉佰号童雀vyz<7v/vVKyMKLI、-II二Jl、亠I丄

50、二JJ0 g眞 切 / -HS査斯柯OP嶽菲ff僧号的样频率兀Z?=2有爲看i1：(Z迅加+177?° A 心 切 /*1他 械釆牌也号的頻腐z?=2),兀探*图12.6.2减采样前后的信号频谱及预滤波器的作用为了使减采样后的信号频谱不致混叠，一般减采样器应先进行低通滤波，削除可能混叠的那部分频谱，然后实施减采样操作。如图12.6.3所示。图12.6.3中，输入序列经过一个由冲激响应h(n)和频率响应H d )描述的低通滤波器。理想情况下，H d ()满足Hd®)珥0,I ©匡兀/D其它(1 2. 6. 2)讪)，1企圃匸波审呛)、蘇样1 7仇F叫2i班妙.高率L

51、:低率4寸:e 1 / = 7；图12.6.3减采样器原理图x(n)在范围因此，该滤波器削除了x(n)在区间H/D £©!；!的谱。这意味着在对信号的进一步处理中，仅仅对(12.6.3)(12.6.4)(1)减采样之前，D =1。根据采样定理，模频特性的幅值随着D值的增大而下降(D =1时，Ximax =1 ； D = 2I © |兀/ D的频率成分感兴趣。滤波器的输出v(n)按下式v(n)=艺 h(k)x(n -k)k=0算出。然后，用因子 D对v(n)减采样，生成y(m)，于是,y (m) = v(mD )3C=Z h(k)x(mD -k)k亠虽然，由式(12.6.3)描述的滤波操作是线性时不变的，但是，与滤波结合的减采样操作产生了一个时变系统。这一点可以很容易得到证明。假定x(n)产生y(m)。我们注意到，除非n0是D的整倍数，x(n -n0)并不意味着产生 y(n -no)。因此，在x(n)上的全部线性操作(先做预滤波，繼之以减采样)不是线性时不变的。程序M120601.演示减采样对信号频谱的影响。图12.6.4示出减采样导致的时域序列和模频特性的变化。图中，先后考 虑建采样因子D =1,2,3,4。对于此图，我们应该注意:时，X 2 max =0.5,)。(2)D =2时的序列长度比 D =1时减小了一半。(3)图