现代导航实验报告光纤陀螺静态测试Allan方差分析

现代导航测试实验报告光纤陀螺静态性能测试Allan方差分析姓名学号学校南京航空航天大学学院自动化学院专业自动化专业班级2014年11月一、 实验目的：1. 了解光学陀螺静态测试的过程。

2. 通过实验测试得到的数据，利用Allan 方差法分析其随机误差特性其随机噪声特性。

二、 实验原理：1. 光纤陀螺仪静态测试1) 静态测试方法：测试转台工作于静止状态，启动陀螺仪稳定工作状态后，以一定的频率采集陀螺仪的输出。

伺服控制测试设备原理图如下图：环境温度图表 1 伺服控制测试设备原理图2) 考虑地球自转带来的静态角速率被陀螺仪敏感的情况，需在输出角速率中去除地球自转角速率在实验所在地（南京：北纬32°03′）的分量：s L n iez /0032.0sin =Ω=ω其中地球自转角速率s rad /10292.75-⨯=Ω。

2. Allan 方差定义与计算Allan 方差法是在时域上对频域特性进行分析的一种方法，为评价光纤陀螺仪的各类误差(包括角度随机游走、零偏不稳定性、角速率随机游走、量化噪声和速率斜坡)特性提供了一种简便的手段．采用该方法，通过对陀螺输出数据构成的一个样本空间进行处理，就可以辨识出陀螺各项误差的系数。

计算Allan 方差的步骤如下所示：1)获取数据。

以固定的采样周期Ts ，采集光纤陀螺的输出角速率，共采样N 个点，得到长度为N 的样本空间。

2)动态分组，分成的每组数据个数是动态变化的。

将样本空间中每m(m=1，2，…，M ，M<N ／2)个数据分成一组，得到k 个独立的数组, 令k=[N ／m]且K=[N ／M]。

3)平均数据。

针对每组数据个数为m 的情况，对每组数据取平均值，即求群平均。

得到元素为群平均的随机变量集合，每一组的平均值为k j mm miim j j ,,2,1其中,1)(1)1( ==∑=+-ωωNm N m N m m m m +-+-++1ωωωωωωωωω,,,,,,,212212)(1ω))m图表 2 Allan 方差计算中的数组平均过程示意图4)计算特定相关时间的Allan 方差。

定义每个数组的持续时间s m mT =τ为相关时间，对于每个特定的相关时间m τ，Allan 方差是通过对相邻群平均之差的平方求平均再乘以1／2得到的：∑-=+--=11212))()(()1(21)(k j j j m m m k ωωτσ 5)计算Allan 方差。

对于不同的m ，分别取s m mT =τ,（s s s m MT T T ,2,=τ），得到在双对数坐标系中的m m ττσ-)(曲线，称为Allan 方差曲线。

3. Allan 方差的意义由Allan 方差的定义可以看出，它是光学陀螺稳定性的一个度量，它和影响陀螺仪性能的固有的随机过程统计特性有关。

Allan 方差与原始测量数据中噪声项的双边功率谱密度（PSD ）)(f S ω存在关系⎰∞=0242)()(sin )(4)(df f f f S A ππππτσω 上式说明：当通过一个传递函数为24)/()(sin ππππf f 的滤波器时，Allan 方差与陀螺仪输出地噪声总能量成正比。

由此，Allan 方差提供了一种方法，能够识别并量化存在于数据中的不同噪声项。

4. 随机噪声的Allan 方差5. Allan 方差的估计估计Allan 方差时，与实验所用的陀螺仪的类型和数据获取的环境有关，实验数据中可能存在各种成分的随机噪声，若噪声源是统计独立的，则Allan 方差可以表示成各类型误差的平方和。

如下：∑-===22222)()(n nntotalA A ττστσ由拟合函数)(τσA 可以求出n A ，再通过下面的计算可以得到量化误差（Q ）、角度随机游走（N ）、零偏不稳定性（B ）、角速率随机游走（K ）和速率斜坡（R ）的估计值，如下：rad A Q μ32-=h A N /601︒=- h A B /6643.00︒=h h A K /)(3601︒= h h A R /)(236002︒=为了估计的准确，陀螺仪输出数据的样本必须足够长，样本长度较短时，Allan 方差的可信度低，导致误差系数估计的可信度不高。

当样本很长、积分时间也较长时，速率斜坡系数R 才能辨识出来。

6. Allan 标准差拟合函数算法解释Allan 方差可以表示成各个类型误差的平方和：n nn total A A ττστσ∑-===22222)()(由于方差一般较小，拟合标准差可以提高拟合精度，上式可近似为()∑-==222/nn n A ττσ需根据前面得到的Allan 标准差)(τσA 与τ的关系求出二阶拟合多项式的系数n A 。

设采样时间[]n ττττ 10=相应Allan 标准差为[])()()()(10n A A AA τστστστσ =设拟合多项式的系数[]21012A A A A AA --=由于⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎣⎡⋅⎥⎦⎤⎢⎣⎡==-----=∑2/22/12/02/12/2210122/22)(τττττττσA A A A A A n nn A 所以[][]⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎣⎡⋅=⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎣⎡=----210122/22/12/02/12/210)()()()()()()()()(A A A A A TT T T T n A A A T A ττττττστστστσ)2,1,0,1,2()2/(.^210122102101212/102/11122/12022/1212112/11012/111112/10002/1010--=⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎣⎡⋅⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎣⎡=⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎣⎡⋅⎥⎥⎥⎥⎥⎥⎦⎤⎢⎢⎢⎢⎢⎢⎣⎡=------------i A A A A A i A A A A A n n n n n n ττττττττττττττττττττττττA i T ⋅=)2/(.^τ（i=-2,-1,0,1,2）令[]),2/(.^,)(,i X F Y X TA T===τστ则Y F \A =，此三步即拟合函数编程算法步骤。

三、 实验数据分析与处理对光纤陀螺的实测数据进行Allan 方差分析。

实验光纤陀螺去除地球自转角速率后的输出测试数据如图所示，取最小分组数为5组，计算Allan 方差，并绘制Allan 标准差相对相关时间的双对数曲线图，最后用二阶曲线拟合求出量化噪声源。

可以看到Allan 方差曲线图稍微有波动，且速率斜坡曲线并不能很好很清晰地表现出来，这是因为陀螺仪输出数据的样本不够长，各个系数还没有稳定。

当样本长度比较短时，Allan 方差的可信度低，从而会导致误差系数估计的可信度不高。

当样本很长，一般大于4h 时，各系数才趋于稳定，速率斜坡系数R 才能辨识出来。

图表 4测试数据误差曲线图图表 5 Allan 方差分析曲线图及拟合曲线最终求出Q =8.6019510-⨯rad μN=5.8278610-⨯h /B= 6.4166410-⨯h /K = 0.0032h h /)/(R =0.0052h h //。

四、 源程序%陀螺静态数据Allan 方差分析 clear; clc[W]=textread('静态.dat','%*f%f%*f%*f%*f');% 读取角速率数据W=W/3600-0.0032; %考虑地球自转角速率ω t=0.02; %采样时间 NN=length(W); %一共有NN 个点M=NN/5; %最少分成5组，每组M 个采样点for m=1:M; %分成的每组数据个数是变化的，先每1个分成1组，然后递增%最后每M 个（M<=NN/2）分成1组k=NN/m; %现在是每m 个1组的特定分法，共能分成k 组 for i=1:kb_aver(i)=mean(W(m*(i-1)+1:m*i)); %此K 组中每组取平均值,即求群平均 ends=0;for j=1:k-1s=s+(b_aver(j+1)-b_aver(j))^2;endsigma(m)=s/(k-1)/2; %Allan方差与m有关，而m是变化的，% m*t为相关时间，可以画出其双对数图形 endtau=(1:M)*t; %相关时间H=2; %拟合阶数X=tau';Y=sigma'; %共轭转置F=zeros(length(X),2*H+1);for i=1:2*H+1d=i-H-1;F(:,i)=X.^(d/2); %F矩阵的每一列为tau的n/2次方（n=-2，-1,0,1,2）endC=F\sqrt(Y); %Y左除F，C是二阶拟合多项式系数CQ=abs(C(1)/sqrt(3))*(1000000/180/3600);QN=abs(C(2)/1)/60;NB=abs(C(3))/(log(2)*2/pi);BK=abs(C(4)*sqrt(3))*60; KR=abs(C(5)*sqrt(2))*3600;Rsigma=sqrt(sigma);figure(1);loglog(tau,sigma); %拟合前的allan标准差曲线grid on; %打开分格线hold on;sigma_nihe=zeros(length(X),1);for i=-H:Hsigma_nihe=C(H+i+1).*X.^(i/2)+sigma_nihe; %sigma_nihe累加，%最终和为拟合的sigma方差endloglog(X,sigma_nihe,'r');xlabel('相关时间对数\tau');ylabel('Allan方差对数\sigma');title('Allan Analysis');figure(2);plot(1:NN,W); %原陀螺仪输出角速率xlabel('采样点');ylabel('输出信号');title('输出信号波形');。