小波變換第五次作業(yè)

1、小波變換第五次作業(yè)專業(yè)：信息與通信工程 學號：406130714098 姓名：徐標1. 設計一 CQMFB，低通濾波器 來自一半帶濾波器。該半帶濾波器的長度為47，通帶截止頻率，試給出，的幅頻響應，單位抽樣響應。2. 產(chǎn)生一信號，它由兩個正弦加白噪聲組成，一個在低頻，一個在高頻，正弦的頻率及與白噪聲的信噪比自己給定。試用所設計的濾波器組對該信號進行分解和重建。比較重建后的效果。設計思路：參照課本181頁（1） 首先設計一個半帶濾波器，N=47，=0.42。根據(jù)第六章半帶濾波器的設計思路，先要用Chebyshew最佳一致逼近法設計一個單帶濾波器G（z），令其通帶截止頻率為2=0.84，長度為2J

2、=24。由此單帶濾波器，可得半帶濾波器=，可以通過時域對g（n）作二倍的插值，并令插值后的序列的中心點位0.5。結果如下：（2）對半帶濾波器進行處理，得到幅頻響應非負的半帶濾波器。 方法：令中間過度的濾波器：，假定為零相位，實現(xiàn)上式的簡單辦法是令：，再令，則是一個半帶濾波器，是非負的。產(chǎn)生半帶濾波器的幅頻響應： 零極分析如圖：由上面左圖圖可以看出，共有46對極零點，其中11個零點在單位圓內，11個零點在單位圓外，其余24個零點在單位圓上，對 做譜分解，因為對于CQFMB而言，譜分解的最佳選擇是使成為最小相位系統(tǒng)，取單位圓內的11個零點以及12個單位圓上的零點賦予 ，從而構造出符合要求的，并得到

3、其幅頻響應、單位抽樣響應，同時，根據(jù) 之間的對應關系，很容易得到其他3個濾波器的幅頻響應、單位抽樣響應。結果如下：和的對數(shù)幅頻響應如下圖 ：下面產(chǎn)生一個信號來檢驗剛才所設計的濾波器組的實際性能，產(chǎn)生一個信號 ，包含兩個正弦信號，頻率分別為 ，還包含一個最大幅值 0.2 的高斯白噪。計算得 的信噪比 SNR=19.6040dB。抽樣頻率 200Hz下，用剛才所設計濾波器組進行處理，得到信號如下程序：clear all; close all; N0=512; %First step: To design a one-band filter G(z) by Chebyshev approximati

4、on % set the cutoff frequency of the G(w); F=0 0.84 1 1; A=1 1 0 0; N=23; B=firpm(N,F,A) a=1; wf=0:pi/N0:pi*(N0-1)/N0; wff=0:1/N0:(N0-1)/N0; Gw=freqz(B,a,wf); Ew=exp(i*N*wf/2); Gr=real(Gw.*Ew); gn=impz(B,a,N); %plot the response curve figure(1); subplot(2,1,1) plot(wff,real(Gr);grid on; title(single

5、-band filt G(jw); subplot(2,1,2) stem(1:N),gn,filled) title(single-band filt g(n) hold on plot(0:N),zeros(N+1),b);grid on; % To obtain half-band filter P(z) from G(z), % F(z)=G(z2)+(N-1)/2/2 N1=47; s=length(B); B2=zeros(1,2*s-1); for k=1:s, B2(k*2-1)=B(k); end, N2=length(B2); B2(N+1)=1+B2(N+1); B2=B

6、2/2; Pw=freqz(B2,a,wf); Ew=exp(i*N*wf); Pr=real(Pw.*Ew); pn=impz(B2,a,N1); figure(2); subplot(2,1,1) plot(wff,Pr);grid on; title(half-band filt H(jw); subplot(2,1,2) stem(pn,filled);grid on; title(half-band filt h(n) ) ;hold on plot(0:N1),zeros(N1+1),b);grid on; %Modify P(z) and prepare to do spectr

7、um decomposition minP=min(Pr); B2(N+1)=B2(N+1)+abs(minP); B3=B2*0.5/(0.5+abs(minP); Pw0=freqz(B3,a,wf); Ew=exp(i*N*wf); Pr0=real(Pw0.*Ew); figure; plot(wff,Pr0);grid; title(Impulse Response of P(z); %Spectrum Analysis A0=zeros(1,N1); A0(1)=1; figure(3); subplot(121); zplane(B3,A0); title(the zeros a

8、nd poles before resolve, P(z) Z,P,K=tf2zp(B3,A0); Z1=sort(Z); ll=length(Z1)/2; ZZ=Z1(1:ll); temp=1; for m=1:ll ; temp=(-ZZ(m)*temp; end;if imag(temp)0.0001 temp=real(temp); else disp(Spectrum Decomposition Failed!); end KK=sqrt(K/temp); l2=length(ZZ); PP=zeros(l2,1); H0Num,H0Den=zp2tf(ZZ,PP,KK); sub

9、plot(122) zplane(H0Num,H0Den); title(the zeros and poles after resolve, H_0(z) H1Num=qmf(H0Num,1); G0Num=-wrev(H0Num); G1Num=qmf(G0Num); %H0 H0w=freqz(H0Num,H0Den,wf); h0n=impz(H0Num,H0Den); figure(4); subplot(2,1,1) plot(wff,abs(H0w);grid on; title(Amplitude-Frequency Response of H_0(z); subplot(2,

10、1,2) stem(h0n,filled);grid on; title(Impulse Response of H_0(z) ;hold on; plot(0:N),zeros(N+1),b);grid on; %H1 H1w=freqz(H1Num,H0Den,wf); %Ew=exp(i*N*wf); %H1r=real(H1w.*Ew); h1n=impz(H1Num,H0Den); figure(5); subplot(2,1,1) plot(wff,abs(H1w);grid on; title(Amplitude-Frequency Response of H_1(z); sub

11、plot(2,1,2) stem(h1n,filled);grid on; title(Impulse Response of H_1(z) hold on; plot(0:N),zeros(N+1),b);grid on; %G0 G0w=freqz(G0Num,H0Den,wf); g0n=impz(G0Num,H0Den); figure(6); subplot(2,1,1) plot(wff,abs(G0w);grid on; title(Amplitude-Frequency Response of G_0(z); subplot(2,1,2) stem(g0n,filled);gr

12、id on; title(Impulse Response of G_0(z) hold on; plot(0:N),zeros(N+1),b);grid on; %G1 G1w=freqz(G1Num,H0Den,wf); g1n=impz(G1Num,H0Den); figure(7); subplot(2,1,1) plot(wff,abs(G1w);grid on; title(Amplitude-Frequency Response of G_1(z); subplot(2,1,2) stem(g1n,filled);grid on; title(Impulse Response o

13、f G_1(z) hold on; plot(0:N),zeros(N+1),b);grid on; %features of filter group H0,w= freqz(h0n,a); H1,w= freqz(h1n,a); absH0=abs(H0); absH1=abs(H1); ah0=20*log10(absH0); ah1=20*log10(absH1); figure(8) plot(wff,ah0,k-,wff,ah1,b-); grid on; xlabel(w/pi);ylabel(|H_0(ejw)|/dB,|H_1(ejw)|/dB); sumh=absH0.*a

14、bsH0+absH1.*absH1; sum=10*log10(sumh); save H_filters h0n disp(Press any key to process a signal using the filters.); pause; clear all; close all; load H_filters h0=h0n; h1=qmf(h0,1); g0=-wrev(h0); g1=qmf(g0); fs=200;f1=8;f2=30; step=-15; N=1000;n=1:N; Sig1=2*sin(2*pi*f1*n/fs); Sig2=2*sin(2*pi*f2*n/

15、fs); P1=var(Sig1+Sig2) Nos=0.2*randn(N,1); P2=var(Nos); snr=10*log10(P1/P2) Sig=Sig1+Sig2+Nos; %signal v0(n) and v1(n) after decompositon and decimation x0=filter(h0,1,Sig); v0=decimate(x0,2); %low pass x1=filter(h1,1,Sig); v1=decimate(x1,2); % high pass %reconstructed signal by filter g0(n) g1(n) u0=interp(v0,2); u1=interp(v1,2); x_recon=filter(g0,1,u0)+filter(g1,1,u1