加hann窗全相位比值校正法和加hann窗fft比值校正法校正方法类同,只须将二个振幅比改为振幅开方比即可。这里加hann窗是关键,过去测试时,直接调用Matlab中的hann(N)窗,频率和振幅校正效果差,见表5加hann窗全相位比值校正法测试结果。表4为加n-hann窗全相位比值校正法,其频率校正精度,相位校正精度和振幅校正精度都很高,甚至可以和表1加n-hann卷积窗apfft/apfft校正法相比,相位和振幅校正精度都很高。但apfft/apfft校正法须3N-1个取样,而全相位比值校正法只须2N-1个取样,克服了apfft/apfft法的一个缺点。现在看来全相位比值校正法也是一种不错的选择,过去窗没有正确使用。对密集频谱,apfft比值法不如apfft/apffy时移法,因时移法只涉及峰线,比值法还涉及次峰线,所以apfft/apfft只须隔2条谱线,比值法须隔5条谱线。

表6加n-hann窗fft比值校正法和表7加hann窗fft比值校正法比,表6的频率校正精度,相位校正精度和振幅校正精度都明显好於表7,说明fft比值校正法加n-hann窗有效。看来zhaoyixu在”学写程序－比值校正法”帖子中己注意到不能调用Matlab中的hann(N)窗,而用公式产生 振-hann窗。

比较表4加n-hann窗全相位比值校正法和表6加n-hann窗fft比值校正法,因apfft的泄漏小,表4频率校正精度,相位校正精度和振幅校正精度都明显好於表6。虽apfft比值校正法须2N-1个取样,而fft比值校正法只须N个取样,但从表4和表6见,无噪时N阶apfft比值法的精度好於N+1阶fft比值法的精度。

表4   加n-hann窗全相位比值校正法204839.29999999999990.99999999999972759.999999999984

102439.29999999999930.99999999999846459.999999999984

51239.29999999998970.99999999997825359.9999999999839

25639.29999999983620.99999999965486259.9999999999817

12839.29999999737840.99999999447591359.9999999996312

6419.29999995796520.99999991137418659.9999999792749

329.299999321759380.99999856871043159.9999989246668

表5   加hann窗全相位比值校正法

204839.29979896916311.000704498666459.9999999999814

102439.29959783064051.00140898299559.9999999999792

51239.29919523050921.0028179082409959.9999999999765

25639.29838873776121.0056355824731459.9999999999818

12839.29677058205441.0112702074365359.9999999999848

6419.29351360044991.0225364540300559.9999999991149

329.286917216297511.0450566433840759.9999999614661表4   加n-hann窗fft比值校正法

204839.29999936670111.0000001326550860.0000897758125

102439.29999936684011.000000132629260.0000897559346

51239.2999993691991.000000132190360.0000894184815

25639.29999941750241.0000001232482360.0000825036076

12839.30000201518230.99999966895578859.9997086611727

6419.30001456858060.99999811325109859.9978626792359

329.300094472854940.99999540530419359.9856736074141

表7     加hann窗fft比值校正法

204839.3002002767745 0.99964781701372459.9375635653786

102439.30040111331390.99929536309503759.8750505713501

51239.30080256054120.99859004124600159.7500640917667

25639.30160439847910.99717776999809559.5002706862587

12839.30320039183170.99434733469568959.00214393244

6419.3063649236880.98866162022938358.0111564495078

329.312561119329660.9771861189596156.0545743160083

表8为N=1024时不同频偏下加n-hann窗全相位比值校正法的测试结果,频率从39.0到39.9,频率校正精度,相位校正精度和振幅校正精度都好於10(-10)。不出现apfft/apfft校正法当频偏量绝对值为0.5时频率估计是会出现相差一个频率分辨率的问题

表8 不同频偏下加n-hann窗全相位比值校正法 (N=1024)

39160.0000000000001

39.09999999999950.99999999999973959.9999999999979

39.19999999999930.99999999999913659.9999999999921

39.29999999999930.99999999999846459.999999999984

39.39999999999960.99999999999795759.9999999999759

39.50.99999999999777259.9999999999723

39.60000000000040.99999999999796359.999999999978

39.70000000000070.99999999999847259.9999999999856

39.8000000000008 0.99999999999914259.999999999993

39.90000000000050.99999999999974159.9999999999982apfft比值法程序

function XfCorrect=SpectrumCorrect(N,xf,CorrectNum,WindowType)

close all;clear all;clc

fs=1024;

N=1024;

t=(-N+1:N-1)/fs;

x=3*cos(2*pi*80*t+30*pi/180)+4*cos(2*pi*150.232*t+80*pi/180)+2*cos(2*pi*253.5453*t+240*pi/180);

w2=sin(pi*(1/2:N-1/2)/N).^2;

w22=conv(w2,w2);

xx2=x.*w22;

x2=xx2(N:end)+[0 xx2(1:N-1)];

xfw=fft(x2);

xfw=xfw(1:N/2);

XfCorrectW_apfft_1024=SpectrumCorrect(N,xfw,3,2);

XfCorrectW_apfft_1024(:,1)=XfCorrectW_apfft_1024(:,1)*fs/N;

XfCorrectW_apfft_1024

function XfCorrect=SpectrumCorrect(N,xf,CorrectNum,WindowType)

XfCorrect=zeros(CorrectNum,3);

for i=1:CorrectNum

A=abs(xf);

[Amax,index]=max(A);

phmax=angle(xf(index));

if (WindowType==2)

indsecL=sqrt(A(index-1))>sqrt(A(index+1));

df=indsecL.*(2*sqrt(A(index-1))-sqrt(Amax))./(sqrt(Amax)+sqrt(A(index-1)))-(1-indsecL).*(2*sqrt(A(index+1))-sqrt(Amax))./(sqrt(Amax)+sqrt(A(index+1)));

XfCorrect(i,1)=index-1-df;

XfCorrect(i,2)=(1-df.^2).^2/(sinc(df).^2)/N/N*Amax*8;

XfCorrect(i,3)=phmax*180/pi;

xf(index-4:index+4)=zeros(1,9);

end

XfCorrect(i,3)=mod(XfCorrect(i,3),360);

XfCorrect(i,3)=XfCorrect(i,3)-(XfCorrect(i,3)>180)*360;

end

转载本文请联系原作者获取授权，同时请注明本文来自王兆华科学网博客。

链接地址：http://blog.sciencenet.cn/blog-212920-567066.html

上一篇：4cfft插值公式证明1

下一篇：Hanning窗生成公式sin(pi*(1/2:N-1/2)/N).^2