日志
全相位比值法和Hann窗生成公式
热度 1 ||
在许多文献和Matlab中的hann(N)窗的生成公式为
hann(N)=sin(pi*(0:N-1)/(N-1)).^2 (1)
这样的生成公式产生的窗是对称窗。它的生成公式中分毋为N-1,对N阶fft不是最好,虽在泄漏等性质上区别不大,但在用於校正时误差大。
在许多文献和Matlab中的hanning(N)窗的生成公式为
hanning(N)=sin(pi*(1:N)/(N+1)).^2 (2)
这样的生成公式产生的窗是对称窗。它的生成公式中分毋为N+1,对N阶fft不是最好,虽在泄漏等性质上区别不大,但在校正时误差大。
若改为另一种Hann窗生成公式(本文写为n-hann)
n-hann(N)=sin(2*pi(1/2:N-1/2)/N).^2 (3)
这样的生成公式产生的窗是对称窗, 它的生成公式中分毋为N,对N阶fft合适,虽在泄漏等性质上区别不大,但在用於校正时改善明显。
也不是hann窗的生成公式都除以N-1,也有除N的. 在振动论坛上zhaoyixu的”学写程序－比值校正法”帖子中(本文写为振-hann)http://www.chinavib.com/thread-65243-1-1.html
振-hann=0.5-0.5*cos(2*pi*(0:N-1)/N)
=sin(pi*(0:N-1)/N).^2 (4)
这样的生成公式生成的是不对称窗, 它的生成公式中分毋为N,,对N阶fft合适,在校正时改善明显。由於是不对称窗,当振-hann和振-hann卷积窗用於apfft时,apfft的水平相位特性消失,相位特性是一条斜线.但振-hann和反振-hann窗(sin(pi*(N-1:-1:0)/N).^2)的卷积窗用於apfft时,apfft的水平相位特性保持。
最近汪意到在matlab中, 还有
hann(N,'periodic '),它的生成公式中分毋为N,它和上面的振-hann相同,它用於DFT/FFT。
hann(N,’symmetric'),它的生成公式中分毋为N-1,它和上面的hann相同,它用於filter design。
过去apfft使用matlab中的hann窗时没有注意到还有指令hann(N,'periodic '),但它的不对称性,用於apfft的卷积窗时还不如n-hann对称窗。
文献介绍说hann窗的频谱只有3根,主线为1,左右为1/2。
式(1)的8个取样值为(Matlab中的 hann(8)' )
hann(8)= [0 0.1883 0.6113 0.9505 0.9505 0.6113 0.1883 0]
其DFT的谱为
Ha=[3.5000 2.0800 0.2135 0.0541 0 0.0541 0.2135 2.0800];
式(2)的8个取样值为(Matlab中的 hanning(8)' )
hanning(8) =[0.1170 0.4132 0.7500 0.9698 0.9698 0.7500 0.4132 0.1170];
其DFT的谱为
Hb=[4.5000 1.8337 0.1080 0.0304 0 0.0304 0.1080 1.8337
式(3)的8个取样值为
n-hann(8) =[0.0381 0.3087 0.6913 0.9619 0.9619 0.6913 0.3087 0.0381];
其DFT的谱为
Hc=[4.00000 2.0000 0.0000 0.0000 0 0.0000 0.0000 2.0000]
Hc只有3根谱线,2:1的比值也对,而Ha和Hb不符合,所以从频谱看,hann窗的生成公式应为式(3)。
式(4)的8个取样值为
振-hann=[0 0.1464 0.5000 0.8536 1.0000 0.8536 0.5000 0.1464];
振-hann是不对称窗,其DFT的谱为
Hd =[4.0000 2.0000 0 0.0000 0 0.0000 0 2.0000];
由於是N等分点取样,Hd频谱和Hc一样.
任意N阶Hc和Hd都只有3条谱线,这是n-hann和振-hann校正性能好的原因。
在apfft/apfft校正法中,无噪时采用hanning和hanning卷积窗相位和频率精度都很高,可达10(=10),但振幅校正精度只有10(-3),采用hann和hanning卷积窗可达校正精度10(-6),还比相位和频率精度低,现采n-hann(N)=sin(pi*(1/2:N-1/2)/N).^2的卷积窗振幅校正精度也可达10(-9)(N=1024)。
如<数字信号全相位谱分析与滤波技术>一书附录三的apfft/a[pfft校正法程序中原调用matlab中的hanning窗(第5句), 改为win=sin(pi*(1/2:N-1/2)/N).^2的卷积,且第2句N=256改N=1024
close all;clc;clear all;
N=1024;
t=-N+1:N*2-1; f1=39.3;A1=1;ph1=60;
s=1*cos(1*(2*pi*t*f1/N+ph1*pi/180));
win=sin(pi*(1/2:N-1/2)/N).^2;
win2=conv(win,win);
win2=win2/sum(win2);
w=pi*2;
s1=s(1:2*N-1);
y1=s1.*win2;
y1a=y1(N:end)+[0 y1(1:N-1)];
Out1=fft(y1a,N);
a1=abs(Out1);
p1=mod(phase(Out1),2*pi);
s2=s(1+N:3*N-1);
y2=s2.*win2;
y2a=y2(N:end)+[0 y2(1:N-1)];
Out2=fft(y2a,N);
a2=abs(Out2);
p2=mod(phase(Out2),2*pi);
g=mod((p2-p1)/pi/2,1);
h=2*(1-g.*g)./sinc(g); aaa=abs((h.^2).*a1)/2;
rr=round(f1);
fff=floor(f1)+g(rr+1)
aaa=aaa(floor(f1)+1)
ppp=p1(rr+1)*180/pi
运行结果:
fff=39.30000000000006
aaa =0.9999999999990419
ppp =59.99999999998399
振幅校正精度10(-12),
可见用n-hann(N)的卷积窗振幅校正精度明显改善。
在apfft/apfft校正法振幅校正中,先用hanning和hanning卷积窗,校正精度仅10(-3),后来偶然发现用hanning和hann卷积窗,校正精度可达10(-6),提高了3-4个数量级。当时不知原因. 现在从hanning和hann窗的生成公式知,一个分毋为N-1,, 一个分毋为N+1,合用后二者平均在N等分上取点,在代入hanning窗插值公式时更合适。
不同N时测试结果如下:
表1 加n-hann和n-hann卷积窗apfft/apfft校正法
2048 39.3000000000001 1.0000000000003 59.999999999984
1024 39.3000000000001 0.999999999999505 59.999999999984
512 39.3000000000001 0.999999999986809 59.9999999999839
256 39.3000000000001 0.99999999978369 59.9999999999817
128 39.3000000000015 0.999999996540251 59.9999999996312
64 19.3000000000846 0.999999944943427 59.9999999792749
32 9.30000000438993 0.999999133831681 59.9999989246667
表2 加hann和hanning卷积窗apfft/apfft 校正法
2048 39.3000000000001 0.999999971337887 59.9999999999844
1024 39.3000000000001 0.999999885351212 59.9999999999854
512 39.3 0.999999541415475 59.9999999999893
256 39.3 0.999998165847839 60.0000000000004
128 39.2999999999998 0.999992666381374 60.0000000000485
64 19.2999999999898 0.999970713387453 60.0000000024893
32 9.29999999963711 0.999883621953122 60.0000000888927
表3 加振-hann和反振-hann卷积窗apfft/apfft 校正法
2048 39.3000000000001 1.00000000000041 59.999999999984
1024 39.3000000000001 1.00000000000132 59.999999999984
512 39.3000000000001 1.00000000001583 59.9999999999841
256 39.3 1.00000000024797 59.9999999999864
128 39.3000000000006 1.00000000396609 59.9999999998605
64 19.3000000000263 1.00000006356964 59.9999999935615
32 9.30000000080167 1.00000102050807 59.999999803628
表1加n-hann和n-hann卷积窗apfft/apff校正法 和 表2加hann和hanning卷积窗apfft/apff校正法比,频率校正精度相同,相位校正精度基本相同,振幅校正精度表1明显改善。说明在apfft/apfft中加n-hann卷积窗改善振幅校正精度有效。
表3加振-hann和反振-hann卷积窗apfft/apfft校正法 和 表1加n-hann和n-hann卷积窗apfft/apfft校正法比,频率校正精度,相位校正精度和振幅校正精度都相同。说明不对称的振-hann窗也可用於apfft。
全相位比值校正法
加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窗全相位比值校正法
2048 39.2999999999999 0.999999999999727 59.999999999984
1024 39.2999999999993 0.999999999998464 59.999999999984
512 39.2999999999897 0.999999999978253 59.9999999999839
256 39.2999999998362 0.999999999654862 59.9999999999817
128 39.2999999973784 0.999999994475913 59.9999999996312
64 19.2999999579652 0.999999911374186 59.9999999792749
32 9.29999932175938 0.999998568710431 59.9999989246668
表5 加hann窗全相位比值校正法
2048 39.2997989691631 1.0007044986664 59.9999999999814
1024 39.2995978306405 1.001408982995 59.9999999999792
512 39.2991952305092 1.00281790824099 59.9999999999765
256 39.2983887377612 1.00563558247314 59.9999999999818
128 39.2967705820544 1.01127020743653 59.9999999999848
64 19.2935136004499 1.02253645403005 59.9999999991149
32 9.28691721629751 1.04505664338407 59.9999999614661
表4 加n-hann窗fft比值校正法
2048 39.2999993667011 1.00000013265508 60.0000897758125
1024 39.2999993668401 1.0000001326292 60.0000897559346
512 39.299999369199 1.0000001321903 60.0000894184815
256 39.2999994175024 1.00000012324823 60.0000825036076
128 39.3000020151823 0.999999668955788 59.9997086611727
64 19.3000145685806 0.999998113251098 59.9978626792359
32 9.30009447285494 0.999995405304193 59.9856736074141
表7 加hann窗fft比值校正法
2048 39.3002002767745 0.999647817013724 59.9375635653786
1024 39.3004011133139 0.999295363095037 59.8750505713501
512 39.3008025605412 0.998590041246001 59.7500640917667
256 39.3016043984791 0.997177769998095 59.5002706862587
128 39.3032003918317 0.994347334695689 59.00214393244
64 19.306364923688 0.988661620229383 58.0111564495078
32 9.31256111932966 0.97718611895961 56.0545743160083
表8为N=1024时不同频偏下加n-hann窗全相位比值校正法的测试结果,频率从39.0到39.9,频率校正精度,相位校正精度和振幅校正精度都好於10(-10)。不出现apfft/apfft校正法当频偏量绝对值为0.5时频率估计是会出现相差一个频率分辨率的问题。
表8 不同频偏下加n-hann窗全相位比值校正法 (N=1024)
39 1 60.0000000000001
39.0999999999995 0.999999999999739 59.9999999999979
39.1999999999993 0.999999999999136 59.9999999999921
39.2999999999993 0.999999999998464 59.999999999984
39.3999999999996 0.999999999997957 59.9999999999759
39.5 0.999999999997772 59.9999999999723
39.6000000000004 0.999999999997963 59.999999999978
39.7000000000007 0.999999999998472 59.9999999999856
39.8000000000008 0.999999999999142 59.999999999993
39.9000000000005 0.999999999999741 59.9999999999982
apfft比值法程序
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
全部作者的其他最新日志
发表评论 评论 (3 个评论)
- 回复 zhwang554
- 最近发现加n-hann窗全相位比值校正法,其频率校正精度,相位校正精度和振幅校正精度都很高,可以和apfft/apfft校正法相比。但apfft/apfft校正法须3N-1个取样,而全相位比值校正法只须2N-1个取样,克服了apfft/apfft法的一个缺点。现在看来单频或非密集频谱信号全相位比值校正法也是一种不错的选择,过去窗没有正确使用。
“全相位比值法和Hann窗生成公式”
http://www.chinavib.com/home-space-uid-62061-do-blog-id-18820.html
将apfft/apfft程序中的hanning(N)窗改为生成函数hann(N)=sin(2*pi(1/2:N-1/2)/N).^2,可进一步改善振幅校正精度
"apfft谱分析中的振幅校正精度"
http://www.chinavib.com/home-space-uid-62061-do-blog-id-17986.html
