%《數(shù)字信號(hào)處理教程——MATLAB釋義與實(shí)現(xiàn)》第七章例7.4.2a程序hc742a
% 過(guò)渡帶設(shè)一個(gè)樣本時(shí)的試湊程序
% 電子工業(yè)出版社出版  陳懷琛編著 2004年9月
%
N=input('N= (N必須為奇數(shù)）');wc=pi/3;           % 給出原始數(shù)據(jù)
N=N+mod(N+1,2);                                 % 使N成為奇數(shù) 
N1=fix(wc/(2*pi/N));N2=N-2*N1-1;        % 劃分通帶和阻帶的樣本點(diǎn)數(shù)
T1=input('過(guò)渡帶樣本T1的值');              % 
A=[ones(1,N1+1),T1,zeros(1,N2-2),T1,ones(1,N1)];	% 符幅特性樣本序列
theta = - pi*[0:N-1]*(N-1)/N;			% 相位特性樣本序列
H=A.*exp(j*theta);				% 頻率特性樣本序列
h=real(ifft(H));				% 反變換求出脈沖序列，去掉運(yùn)算誤差造成的虛部
wp1=2*pi/N*fix(wc/(2*pi/N));ws1=wp1+2*2*pi/N;
[db,mag,pha,grd,w] = myfreqz(h,[1]);    % 檢驗(yàn)設(shè)計(jì)出的濾波器的頻率響應(yīng)
[Ar,ww,type,L0]=amplres(h);           % 檢驗(yàn)設(shè)計(jì)出的濾波器的符幅特性
dw = 2*pi/1000;                     % myfreqz函數(shù)中的頻率分辨率
Rp = -min(db(1:fix(wp1/dw)+1))      % 實(shí)際的通帶波動(dòng)
As = -round(max(db(fix(ws1/dw)+1:501)))  % 最小阻帶衰減
l = 0:N-1; wl = (2*pi/N)*l;
wdl=[0,wc,wc,2*pi-wc,2*pi-wc,2*pi]/pi;Adl=[1,1,0,0,1,1];
subplot(2,2,1);plot(wl(1:N)/pi,A(1:N),'.',wdl,Adl); 
axis([0,1,-0.1,1.1]); title('頻率樣本')
xlabel(''); ylabel('A(k)')
set(gca,'XTickMode','manual','XTick',chop([0,2/N,4/N,1],2)) % 繪制刻度線
set(gca,'YTickMode','manual','YTick',[0,1]); grid
subplot(2,2,2); stem(l,h,'.');
title('脈沖響應(yīng)');ylabel('h(n)');text(N+1,-0.1,'n')
subplot(2,2,3); plot(ww/pi,Ar,wl(1:N)/pi,A(1:N),'.');
axis([0,1,-0.2,1.2]); title('符幅響應(yīng)')
xlabel('頻率（單位：pi）'); ylabel('Ar(w)')
set(gca,'XTickMode','manual','XTick',chop([0,2/N,4/N,1],2)) % 繪制刻度線
set(gca,'YTickMode','manual','YTick',[0,1]); grid
subplot(2,2,4);plot(w/pi,db); axis([0,1,-50,10]); grid
title('幅度響應(yīng)'); xlabel('頻率（單位：pi）');
ylabel('分貝數(shù)');
set(gca,'XTickMode','manual','XTick',chop([0,2/N,4/N,1],2)) % 繪制刻度線
set(gca,'YTickMode','Manual','YTick',[-As;0]);
set(gca,'YTickLabelMode','manual','YTickLabels',['As';' 0'])
set(gcf,'color','w');                % 置圖形背景色為白