第4版程序集_数字信号处理第四版所有例程_

%《数字信号处理(第4版）》第7章图7.3.1和图7.3.2绘图程序：Fig731-2.m % 西安电子科技大学出版社出版 高西全 丁玉美 合著 2016年 % 用频率采样法设计FIR低通滤波器 clear;close all %============================================================================================ %采样点数N=13 N=15;wc=pi/3; Np=fix(wc/(2*pi/N));Ns=N-2*Np-1; %Np+1为通带[0,wc]上采样点数,Ns为阻带[wc,2*pi-wc]上采样点数 Ak=[ones(1,Np+1),zeros(1,Ns),ones(1,Np)]; %幅度采样向量A(k) thetak=-pi*(N-1)*(0:N-1)/N; %相位采样向量θ(k) Hk=Ak.*exp(j*thetak); %构造频域采样向量H(k) hn=real(ifft(Hk)); %h(n)=IDFT[H(k)] %============================================================================================== %绘图7.3.1 Aw=[1,1,0,0,1,1];wk1=[0,wc,wc,2*pi-wc,2*pi-wc,2*pi]/pi; subplot(2,2,1);plot(wk1,Aw);axis([0,2,-0.2,1.2]); xlabel('ω/π');ylabel('H_g(k)');title('图7.3.1(a) 频域幅度采样') hold on wk2=[0:N-1]*2/N; plot(wk2,Ak,'.'); n=0:N-1; subplot(2,2,2);stem(n,hn,'.');axis([0,15,-0.2,0.5]) xlabel('n');ylabel('h(n)');title('图7.3.1(b)') Hw=fft(hn,1024); wk=2*pi*[0:1023]/1024; Hgw=Hw.*exp(j*wk*(N-1)/2); %===================================================== %计算N=15时的通带最大衰减rp15和阻带最小衰减rs15 rp15=max(20*log10(abs(Hgw))) hgmin=min(real(Hgw)) rs15=20*log10(abs(hgmin)) %======================================================== subplot(2,2,3);plot(wk/pi,real(Hgw));axis([0,2,-0.2,1.2]); hold on; plot(wk1,Aw,':');plot(wk2,Ak,'.'); xlabel('ω/π');ylabel('Hg(ω)');title('图7.3.1(c)频域幅度内插波形') subplot(2,2,4);plot(wk/pi,20*log10(abs(Hgw)));axis([0,1,-30,3]);grid on xlabel('ω/π');ylabel('20lg|Hg(ω)|');title('图7.3.1(d)');set(gca,'Xtick',0:0.2:1,'Ytick',-40:10:10) %============================================================================================ figure(2) %图7.3.2 plot(wk/pi,Hgw);axis([0,1,-0.2,1.2]);title('图7.3.2 N=15和N=75时的幅度内插波形') hold on; plot(wk1,Aw,':');plot(wk2,Ak,'o'); %============================================================================================= %采样点数N=65 N=75; Np=fix(wc/(2*pi/N));Ns=N-2*Np-1; %Np+1为通带[0,wc]上采样点数,Ns为阻带[wc,2*pi-wc]上采样点数 Ak=[ones(1,Np+1),zeros(1,Ns),ones(1,Np)]; %幅度采样向量A(k) wk2=[0:N-1]*2/N; plot(wk2,Ak,'.'); thetak=-pi*(N-1)*(0:N-1)/N; %相位采样向量θ(k) Hk=Ak.*exp(j*thetak); %构造频域采样向量H(k) hn=real(ifft(Hk)); %h(n)=IDFT[H(k)] Hw=fft(hn,1024); wk=2*pi*[0:1023]/1024; Hgw=Hw.*exp(j*wk*(N-1)/2); %===================================================== %计算N=75时的通带最大衰减rp15和阻带最小衰减rs15 rp75=max(20*log10(abs(Hgw))) hgmin=min(real(Hgw)) rs75=20*log10(abs(hgmin)) %======================================================== plot(wk/pi,Hgw); xlabel('ω/π');ylabel('Hg(ω)')