在北京王府井书店,西单图书大厦,中关村图书大厦都有买,各地大书店也一样,或网上买 http://www.xhsd.com.cn/books/views.asp?PLUCODE=712108133 下面是书中附录一的FFT和全相位FFT比较程序,运行后如下图,图中兰色为FFT,红色为全相位FFT,从相位谱图可见红色全相位FFT在20赫处相位45度,在28赫处相位90度, 全相位FFT可在信号非整数倍取样时测相位,如下例中20.4赫时在20赫处测相位45度,28.2赫在28赫测相位90度,为什么能测相位只好看书了 close all;clc;clear all; N=256; t=-N+1:N-1; y= cos(2*pi*t*20.4/N+pi/4)+0.001*cos(2*pi*t*28.2/N+pi/2);%正弦波初相为pi/4和pi/2 y1 = y(N:2*N-1);%取后N个数据为FFT输入数据 win = hanning(N)'; win1 = win/sum(win);%窗归1 y11= y1.*win1; y11_fft = fft(y11,N); a1 = abs(y11_fft);%apFFT振幅谱 p1 = mod(phase(y11_fft)*180/pi,360);%FFT相位谱(度) y2 = y(1:2*N-1);%取2N-1个数据为apFFT输入数据 winn = conv(win,win);%apFFT须要卷积窗 win2 = winn/sum(winn);%窗归1 y22= y2.*win2; y222=y22(N:end)+[0 y22(1:N-1)];%构成长N的FFT输入数据 y2_fft = fft(y222,N); a2 = abs(y2_fft);%apFFT振幅谱 p2=mod( phase(y2_fft)*180/pi,360);%apFFT相位谱 tt=0:N-1; subplot(211);plot(tt,10*log10(a1),'b',tt,10*log10(a2*1),'r'); title('amplitude spectrum'); ylim([-150,30]);xlim([0,N/2]); legend('fft','apfft'); xlabel('f');ylabel('db'); grid subplot(212);plot(tt,p1,'b',tt,p2,'r'); title('phase spectrum'); ylim([0,400]);xlim([0,N/2]); legend('fft','apfft'); xlabel('f');ylabel('度'); grid
[此贴子已经被作者于2009-6-16 19:03:36编辑过] |