学写程序之四(源代码)---比值校正法之时域平移法

zhaoyixu

第四贴，欢迎大家拍砖头啊，还请各位多多指教。参考自丁康老师《离散频谱校正技术》一书和版主yangzj的相关贴子。
~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
% Fs:采样频率  N：做谱点数  L：平移点数
clear all;clc  
Fs=1024;N=1024;L=100;
t =(0:N+L-1)/Fs;
tt=0:N-1;
hanning=0.5-0.5*cos(2*pi*tt/N);
windowtype=input('请选择加窗类型1.矩形窗2.汉宁窗');
x=10.343*cos(2*pi*298.30453*t+240*pi/180);%.*hanning;%L+N点时间序列
%取出两段时间序列
for i=1:N
    x1(i)=x(i);
    x2(i)=x(i+L);
end
if windowtype==1
    y1=fft(x1);%对第一段信号加矩形窗做FFT变换
    y2=fft(x2);%对第二段信号加矩形窗做FFT变换
elseif windowtype==2
    y1=fft(x1.*hanning);%对第一段信号加汉宁窗做FFT变换
    y2=fft(x2.*hanning);%对第二段信号加汉宁窗做FFT变换
else
    error('选择有误,请重新选择');
end
Y1=abs(y1(1:N/2)/N*2);%第一段信号幅值
Y2=abs(y2(1:N/2)/N*2);%第二段信号幅值
f=(1:N/2)*Fs/N;
subplot(211);
   if windowtype==1
        plot(f,Y1);
        xlabel('f');ylabel('A');title('加矩形窗校正前');grid on
    elseif windowtype==2
        plot(f,2*Y1);
        xlabel('f');ylabel('A');title('加汉宁窗校正前');grid on
   end
[Y1Amax,k1]=max(Y1);
[Y2Amax,k2]=max(Y2);
phase1=angle(y1(k1));%第一段峰值所对应的相位
phase2=angle(y2(k2));%第二段峰值所对应的相位
   if windowtype==1
        A_uncorrect=Y1Amax
    elseif windowtype==2
        A_uncorrect=Y1Amax*2
   end
f_uncorrect=(k1-1)*Fs/N
phase_uncorrect=phase1*180/pi
delt=mod(phase2-phase1-2*pi*(k1-1)*L/N,2*pi);
%将delt调整到（-pi,pi）之间
if delt<-pi
    delt1=delt+2*pi;
elseif delt>pi
    delt1=delt-2*pi;
else delt1=delt;
end
deltf=delt1/(2*pi*L/N);
f_correct=(k1-1+deltf)*Fs/N
phase_correct=(phase1-deltf*pi)*180/pi        
Y_correct=zeros(1,N/2);
if windowtype==1
    A_correct=Y1Amax/sinc(deltf)
    Y_correct(k1)=A_correct;
    f(k1)=f_correct;
    subplot(212);stem(f,Y_correct);grid on
    xlabel('f');ylabel('A');title('加矩形窗校正后');
elseif windowtype==2
    A_correct=2/sinc(deltf)*Y1Amax*(1-deltf^2)
    Y_correct(k2)=A_correct;
    f(k2)=f_correct;
    subplot(212);stem(f,Y_correct);grid on
    xlabel('f');ylabel('A');title('加汉宁窗校正后');
end
~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
运行结果：
理论值：
幅值：10.343    频率：298.30453        相位：240
加矩形窗
校正前                                                                校正后
A_uncorrect =8.830978191997250                   A_correct =10.362781259858098
f_uncorrect =298                                               f_correct =2.983068229594135e+002
phase_uncorrect =-65.283647874590685        phase_correct =-1.205117805690294e+002
加汉宁窗
校正前                                                               校正后
A_uncorrect =9.739025234720254                   A_correct =10.342999982203294   
f_uncorrect =298                                               f_correct =2.983045299993504e+002
phase_uncorrect =-65.184599901407665        phase_correct =-1.199999997844760e+002

[ 本帖最后由 zhaoyixu 于 2008-6-10 00:28 编辑 ]

yangzj

题目写错了,是相位差法不是比值法

zhaoyixu

原帖由 yangzj 于 2008-6-10 23:42 发表
题目写错了,是相位差法不是比值法

太不小心了貌似改不了了

cdccdq

看了很有收获，谢谢