查看: 3187|回复: 0|关注: 0

[我分享] 离散频谱校正相位差+单点FT法

[复制链接]

论坛优秀回答者

中级

858 麦片

财富积分


5001500


53

主题

1717

帖子

154

最佳答案
  • 关注者: 47
发表于 2017-6-30 19:54:48 | 显示全部楼层 |阅读模式
function resultCorrect=spectrumcorrectftphaseoffsetmethod(inputDate,N,L,fs)
%功能:离散频谱校正相位差+单点FT法,适用单频点信号校正,只采用了加汉宁窗的结果
%注意:信号的模型为Acos(2*pi*f*t+pha),注意t从0开始,注意信号长度和变换长度的关系,通常信号长度选为N+L
%输入:inputDate待分析数据,数据统一为行向量;N为变换点数,L平移点数;fs采样频率
%输出:resultCorrect校正后的频率,幅值,相位结果

resultCorrect=zeros(1,3);   
w=hann(N,'periodic');   %生成汉宁窗
x=inputDate(1:N).*w';       %加窗信号
y=inputDate(L+1:N+L).*w';   %时移加窗信号
X=fft(x);  
Y=fft(y);  
X=X(1:N/2)/N*4;
Y=Y(1:N/2)/N*4;             %单边复数谱,汉宁窗要乘以恢复系数2
magX=abs(X);
[~,maxIndex]=max(magX);   %最大值对应幅值和位置
maxAngleX=angle(X(maxIndex));   %最大幅值处对应的相位
maxAngleY=angle(Y(maxIndex));   %最大幅值处对应的相位
phaseOffset=maxAngleY-maxAngleX; %相位差
detla=phaseOffset-2*pi*L*(maxIndex-1)/N;
detla=mod(detla,2*pi);
if  detla<-pi
    detla=detla+2*pi;
end
if  detla>pi
    detla=detla-2*pi;
end
df=detla/(2*pi*L/N);
resultCorrect(1,1)=(maxIndex-1+df)*fs/N;       %频率校正结果,注意matlab下标是从1开始的
n=1:N;
% fftDateReal=sum(x.*cos(2*pi*(n-1)*df));   %细化实部
% fftDateImag =-sum(x.*sin(2*pi*(n-1)*df));  %细化虚部
% fftDateMag=sqrt(fftDateReal^2+fftDateImag^2)/N;
fftDateReal=sum(x.*cos(2*pi*(n-1)*resultCorrect(1,1)/fs));   %细化实部
fftDateImag =-sum(x.*sin(2*pi*(n-1)*resultCorrect(1,1)/fs));  %细化虚部
fftDateMag=4*sqrt(fftDateReal^2+fftDateImag^2)/N;
resultCorrect(1,2)=fftDateMag;        %细化幅值谱
resultCorrect(1,3)=maxAngleX-pi*df;           %校正相位结果
resultCorrect(1,3)=mod(resultCorrect(1,3),2*pi);
resultCorrect(1,3)=resultCorrect(1,3)-(resultCorrect(1,3)>pi)*2*pi;    %象限定在(-pi,pi]  
end

回复主题 已获打赏: 0 积分

举报

您需要登录后才可以回帖 登录 | 注册

本版积分规则

关闭

站长推荐上一条 /4 下一条

快速回复 返回顶部 返回列表