FFT频谱相频特性仿真结果的问题
本帖最后由 guodongliang 于 2011-5-23 16:53 编辑用FFT进行相频特性研究，如果横坐标都是模拟频率（Hz），matlab程序如下：
clc
clear
%---------------------- non-uniform sampling
N=1000;
n = 1:N;
f0=100; f1=200; f2=400;% the frequecy of signal
fs=1000; % the sampling frequency
t=n/fs;
x=sin(2*pi*f0*t)+sin(2*pi*f1*t+pi/2)+sin(2*pi*f2*t+pi/4); % signal
Xk = fft(x);
F=Xk(1:N/2+1); f=fs*(0:N/2)/N;
subplot(311);
plot(t,x);
title('sampled signal')
xlabel('Time / s')
subplot(312);
plot(f,abs(F))
xlabel('Frequency / Hz')
ylabel('DFT magnitude')
subplot(313);
ph = angle(F)/pi*180;%degree
plot(f,ph);
xlabel('Frequency / Hz')
ylabel('DFT Phase / Degree')
% find the phase
=sort(abs(F),'descend');
k=4;
num(1:k)
val(1:k)
ph(val(1:k))
结果如图所示：
问题是：相频的结果（如同噪声）对吗？
最后一段代码是查出幅频图上最大4个幅频对应的频点，然后查对应的相位：
ans =
500.0000500.0000500.0000 0.0000
ans =
201 401 101 400
ans =
72.0000 99.0000-54.0000 90.3364
结果与所设定的初相（0、pi/2 、pi/4）不吻合。
请高手指点。
怎么没人回复呢 本帖最后由 songzy41 于 2011-5-25 05:48 编辑
guodongliang 发表于 2011-5-23 15:32 http://www.chinavib.com/static/image/common/back.gif
用FFT进行相频特性研究，如果横坐标都是模拟频率（Hz），matlab程序如下：
clc
clear
把程序做了一些修改，程序下面列着。主要修改内容为：
1，时间变量t的设置不正确，改为t=(n-1)/fs;，即初始时刻t=0，这样才有可能求出初始相角（t=0时刻的相角）；
2，从angle得到的相角是针对余弦信号的，与正弦信号差90度；
3，因为3个主要频率对应的谱线分别为101，201，401条谱线相接从这几条谱线中读入计算，其答案为（其中直接计算了信号的幅值）
k=101Am=1.00000 Freq=100.00000 Angle=0.00000
k=201Am=1.00000 Freq=200.00000 Angle=90.00000
k=401Am=1.00000 Freq=400.00000 Angle=45.00000
程序为：
clc
clear
%---------------------- non-uniform sampling
N=1000;
n = 1:N;
f0=100; f1=200; f2=400;% the frequecy of signal
fs=1000; % the sampling frequency
t=(n-1)/fs;
x=sin(2*pi*f0*t)+sin(2*pi*f1*t+pi/2)+sin(2*pi*f2*t+pi/4); % signal
Xk = fft(x);
F=Xk(1:N/2+1); f=fs*(0:N/2)/N;
subplot(311);
plot(t,x);
title('sampled signal')
xlabel('Time / s')
subplot(312);
plot(f,abs(F))
xlabel('Frequency / Hz')
ylabel('DFT magnitude')
subplot(313);
ph = angle(F)/pi*180+90;%degree
plot(f,ph);
xlabel('Frequency / Hz')
ylabel('DFT Phase / Degree')
fprintf('k=101Am=%6.5f Freq=%6.5f Angle=%6.5f\n',abs(F(101))*2/N,f(101),ph(101));
fprintf('k=201Am=%6.5f Freq=%6.5f Angle=%6.5f\n',abs(F(201))*2/N,f(201),ph(201));
fprintf('k=401Am=%6.5f Freq=%6.5f Angle=%6.5f\n',abs(F(401))*2/N,f(401),ph(401)); 感谢大师！
还有一个疑问：
这个模拟信号只有3个频率分量，并无噪声，幅度频谱是只有3个点有频率的，也几乎无频谱泄露，为何相频曲线在所有频点处都有相位呢？ 是计算误差？ guodongliang 发表于 2011-5-25 11:41 static/image/common/back.gif
感谢大师！
还有一个疑问：
这个模拟信号只有3个频率分量，并无噪声，幅度频谱是只有3个点有频率的，也几 ...
这是由于计算误差造成的。首先计算机数字的表示都是用有限字长表示的，例如pi，我们知着它是一个无限长的无理数，但在计算机中只能用有限长度表示。pi=3.1415926...，而在计算机中只能表示为3.14159，后面的0.0000026...就是有限字长的截断误差了。由pi计算出的三角函数，例如sqrt(2)/2，sqrt(3)/2等等原本也应为无理数，现在也按有限字长表示了，就是数的表示中就有误差，用这些带有误差的数参与运算，使误差累计。当然这些误差很小，信号的幅值计算中都看不出来。但在相位计算出就不一样了，例如FFT后X(i)=1.2e-9+j2.4e-10，它的幅值(模值)就很小，但相角angle(1.2e-9+j*2.4e-10)=0.1974就不小了，所以相频曲线在-pi至pi之间振荡是正常的。 感谢大师！
以后多多指教！ {:{39}:}学习下！
页:
[1]