在北京王府井书店,西单图书大厦,中关村图书大厦都有买,各地大书店也一样,或网上买 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编辑过] | 
