%《數(shù)字信號(hào)處理教程——MATLAB釋義與實(shí)現(xiàn)》第四章例5.4.3程序hc543
% 非單位圓上螺旋線Chirp的頻譜分析
% 電子工業(yè)出版社出版  陳懷琛編著 2004年9月
%
p1=0.9*exp(j*pi/3);p2=0.9*exp(-j*pi/3);
%a=conv([1,-p1],[1,-p2]);b=1;
%[r,p,k]=residuez(b,a)
M=64,n=0:M-1;dw=pi/M;                         % 規(guī)定點(diǎn)數(shù)，求頻率分辨率
x=0.9.^n.*(cos(pi*n/3)+2*0.2887*sin(pi*n/3)); % 給定x序列
subplot(3,2,1),stem(n,x,'.')                  % 畫序列曲線
H=freqz(x,1,M);                               % 求序列頻率響應(yīng)
subplot(3,2,2),plot([0:M-1]*dw,abs(H))        % 畫H的幅頻特性
x1=x.*(0.9).^-n; H1=freqz(x1,1,64);           % 求r=0.9的中心圓上的H(z)，先修改x
subplot(3,2,3),plot([0:M-1]*dw,abs(H1))       % 畫H1的幅頻特性
A=1;W0=0.9^(-3/M),W=W0*exp(-j*dw);            % 確定螺旋線軌跡參數(shù)A和W
z = A * W.^(-(0:M-1));                        % 寫出螺旋線軌跡方程
subplot(3,2,4),plot(z),grid on                % 畫出螺旋線軌跡
axis([-1.2,1.2,-1.2,1.2]),axis equal
hold on,plot(p1,'x')                        % 把序列極點(diǎn)畫在螺旋線軌跡圖上
G=czt(x,M,W,A);                               % 調(diào)用czt函數(shù)
subplot(3,2,5),plot([0:M-1]*dw,abs(G(1:M)))   % 畫出czt函數(shù)的幅特性
subplot(3,2,6),plot([0:M/2-1]*dw,abs(G(1:M/2)))  % 畫前半段czt的幅特性益3顯示峰值
set(gcf,'color','w');                           % 置圖形背景色為白