全相位书中附录3遇到出现幅值NaN的问题
论坛中给的校正的判断部分是:%%%%%%%%%%%%%%%%%%%%%
dph=angle(F2)-angle(F1);
dph=mod(dph,2*pi);
if dph>pi
dph=dph-2*pi;
elseif dph<-pi
dph=dph+2*pi;
end
dk=dph/pi/2;
g=dk
h=2*pi*g.*(1-g.*g)./sin(pi*g);
AA=abs((h.^2).*F1)/2;
%%%%%%%%%%%%%%%%%%%%
这里如果g=0会出现NaN的结果
在丁康老师的书里好象叫什么死点的,这个我不太懂，希望各位高手能给点建议！！
下边是程序：
function[] = apFFT_Correct(N,Fs,f,A,Ph)
CNum = length(f);
t=-N+1:N*2-1; %apfft时移相位差法需要3N-1点
s = zeros(1,length(t));
for i=1:CNum
myt=A(i)*cos(2*pi*t*f(i)/Fs+Ph(i)*pi/180);
s=s+myt;
end
plot(s);
s=s(1:3*N-1);
win=hanning(N)';
win1=hann(N)';
win2=conv(win,win1);
win2=win2/sum(win2); %归一化卷积窗
w=pi*2;
s1=s(1:2*N-1); %取前2N-1点做FFT
s1=s1-ones(size(s1)).*mean(s1);%消直流
y1=s1.*win2;
y1a=y1(N:end)+;
F1=fft(y1a,N);
s2=s(1+N:3*N-1); %后2N-1点
s2=s2-ones(size(s2)).*mean(s2);
y2=s2.*win2;
y2a=y2(N:end)+;%全相位预处理
F2=fft(y2a,N);
a1=abs(F1);
L=fix(N/2)+1; %只查看一半的幅值
f=(0:L)*Fs/N;
for i=1:CNum
= max(a1(1:L));
k=Ind(i);
peak_F2(i)=F2(k);
peak_F1(i)=F1(k);
StartP=Ind(i)-2;
StopP=Ind(i)+2;
if Ind(i)-2<1
StartP=1;
end
if Ind(i)+2>L
StopP=L;
end
Indx=StartP:StopP;
a1(Indx)=zeros(1,length(Indx));
end
num=length(Ind);%%%%%设置所取的“谐波簇”个数
for m=1:num;
disp('振幅校正前')
peak_F2(m)
peak_F1(m)
= apFFT_Corret1(Ind(m),peak_F2(m),peak_F1(m));
end
disp('校正的频率')
fm=*Fs/N
disp('校正的相位')
PH=mod(PH*180/pi,360)
disp('校正的振幅')
AA
function = apFFT_Corret1(k0,F2,F1)
dph=angle(F2)-angle(F1);
dph=mod(dph,2*pi);
if dph>pi
dph=dph-2*pi;
elseif dph<-pi
dph=dph+2*pi;
end
dk=dph/pi/2;
g=dk
h=2*pi*g.*(1-g.*g)./sin(pi*g);
AA=abs((h.^2).*F1)/2;
fff=(k0+dk-1);
PH=angle(F1);
下边是测试程序：
%%%%%%%%%%%%%%%%%%校正程序测试
close all;
clear all;
clc;
%%%%%%%%%%%%%%%%%%%%%校正方法选择：0，时移位相位差法；1，fft/apfft法 %%%%%%%%%%%%%
method = 0;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
Fs = 200;%采样频率Hz
N = 256; %FFT计算的点数,实际采样点数比这个多2N-1/3N-1
f1=; %频率
a1=; %幅度
ph1 =; %相位
if method
fft_apfft(N,Fs,f1,a1,ph1)
else
apFFT_Correct(N,Fs,f1,a1,ph1);
end
%这里我模拟采样长度不是整树周期，Fs=200, N 固定255点。幅度校正出了个NaN 那个Fs = 200 N=256,上边打错了
回复 楼主 无限剑制 的帖子
f1=100, 在fs=200, N=256时, 折合频率=f1*N/fs=128, 即振幅谱峰值在N=256频谱的中点,不只振幅校正不对,相位也不对.在其它点, g=0(折合频率是整数)没有关系, 仍能校正振幅
[ 本帖最后由 zhwang554 于 2010-2-26 17:49 编辑 ] 哦，明白了，我这没满足采样定理，信号没学好。
感谢王老师的指点:lol 能把 apFFT_Correct 的程序分享下么？感谢
页:
[1]