allan方差分析法.docx
《allan方差分析法.docx》由会员分享,可在线阅读,更多相关《allan方差分析法.docx(11页珍藏版)》请在冰豆网上搜索。
allan方差分析法
《现代导航测试技术》
实验报告
实验名称:
光纤陀螺随机误差的Allan方差分析法
班级:
0309103
学号:
030910309
姓名:
许莹
时间:
2012-12-17
一:
实验目的
由于光学陀螺的工作原理和环境干扰等原因,在光学陀螺输出信号中包含很多确定性和随机性的误差项。
光学陀螺的随机误差主要包括量化噪声、角速度随机游走、零偏不稳定性】角速率随机游走、速度斜坡和正弦分量,其中前三项误差被认为是其光学性能指标一部分。
对于这些随机误差,利用常规的分析方法,例如计算机样本均值和方差并不能揭示潜在的误差源,另一方面,虽然自相关函数和功率谱密度函数分别从时域和频域描述了随机误差的统计特性,但是在实际工作中通过这些函数加以分析将随机误差分离出来是很困难的。
Allan方差法是20世纪60年代由美国国家标准局DavidAllan提出的,它是一种基于时域的分析方法,不仅可以用来分析光学陀螺的误差特性,而且还可以应用于其他任何精密测量仪器.Allan方差法的主要特点是能非常容易对各种误差源及其对整个噪声统计特性的贡献进行细致的表征和辨识,而且便于计算,易于分离。
它提供了一种识别并量化存在于数据中的不同噪声项。
二:
实验原理与实验内容
1.Allan方差定义与计算
设以采样时间
0对陀螺仪输出角速率进行采样,共采样N个点,把所获得的N个数据分成K,每组包含M个采样点。
K=N/M,M
(N-1)/2
如图:
1,
2,…,
M
M+1,
M+2,…,
2M
N-M+1,
N-M+2,…,
N
K=1k=2k=K
每一组的持续时间
M=M
0,称之为相关时间,每一组的平均值为
k(M)=
k=1,2,3,…,K
Allan方差定义为:
2.Allan方差法最大的优点是可以简便细化分离、辨识光学陀螺的各项误差,同时确定各个误差项对总误差的贡献。
估计方差时,与实验所用的陀螺仪类型和数据类型获取的环境有关,实验数据中可能存在各种成分的随机噪声。
若各项噪声是统计独立的,则Allan方差可以表示成各类误差的平方和,如:
观察之后发现除去正弦噪声含有正弦分量外,其他都可以表示为
(n=-2,-1,0,1,2)。
因此上式可以简化为:
由于方差一般较小,拟合标准差则可以提高拟合精度,上式可以近似为:
其中
=m
0,
0为采样时间
在最小均方意义上,拟合函数
可以求出An,再通过以下计算可以得到量化误差Q,角度随机游走N,零偏不稳定性B,角速率随机游走K和速率斜坡R的估计值,为
R=A23600*sqrt
(2);
K=A160*sqrt(3);
B=A00.6643;
N=A-10;
Q=A-2qrt(3);
三:
实验过程(解算和数据拟合)
1.选取数据文件:
从计算机得到的数据文件中抽取一部分,数据特点是除陀螺仪绕Z轴输出变化以外都比较稳定,可以由Origin软件的图表显示。
如图:
选取其中一部分之后得到的z轴角速率输出曲线(50到-50之间变化)
2.源程序
A.将取得的数据分为5组,得到不同的总采样时间,再有Allan方差公式得到不同的方差,通过列出五个方程的方程组解出五个系数值。
程序为:
%dynamictestaboutz-rotate
%byallanvarianceanalysis
clear
clc
formatlonge
[t,w_x,w_y,w_z,a_x,a_y,a_z,b]=textread('E:
\allan\z-rotate.txt');
N=size(t,1);
sigma_squa=zeros(5,1);
w=zeros(5,5);
M=[floor(N/5),floor(N/4),floor(N/3),floor(N/2),N];
k=[5,4,3,2,1];
forj=1:
5
m=M(1,j)
g=N/m
fork=1:
g
sum=0
fori=1:
m
sum=sum+w_y((k-1)*m+i,1)
aver=sum/m
end
w(j,k)=aver
end
fori=1:
5
forj=1:
(g-1)
temp=0
temp=temp+(w(i,j+1)-w(i,j))^2
sigma_squa(i,1)=temp/(2*g-2)
end
end
end
sigma=sqrt(sigma_squa)
theta=0.02*M;
p=[theta(1,1)^-1,theta(1,2)^-1,theta(1,3)^-1,theta(1,4)^-1,theta(1,5)^-1
theta(1,1)^-0.5,theta(1,2)^-0.5,theta(1,3)^-0.5,theta(1,4)^-0.5,theta(1,5)^-0.5
11111
theta(1,1)^-0.5,theta(1,2)^-0.5,theta(1,3)^-0.5,theta(1,4)^-0.5,theta(1,5)^-0.5
theta(1,1)^1,theta(1,2)^1,theta(1,3)^1,theta(1,4)^1,theta(1,5)^1]
if(rank(p)<5)
h=pinv(p);
else
h=pinv(p)
A=abs(sigma'*h)
R=A(5)*3600*sqrt
(2)
K=A(4)*60*sqrt(3)
B=A(3)/0.6643
N=A
(2)/60
Q=A
(1)/sqrt(3)
xishu1=[QNBKR]'
x=0:
0.1:
6.40
y=A(1,1)*x.^-1+A(1,2)*x.^-0.5+A(1,3)+A(1,4)*x.^0.5+A(1,5)*x
plot(x,y,'r')
函数曲线为:
B:
将所取的数据文件分为11组,得到11个方差值,进行最小二乘法拟合,得到拟合函数。
程序为:
clear
clc
formatlonge
[t,w_x,w_y,w_z,a_x,a_y,a_z,b]=textread('E:
\allan\z-rotate.txt')
N=size(t,1)
sigma_squa1=zeros(11,1)
w=zeros(11,11)
M=[floor(N/11),floor(N/10),floor(N/9),floor(N/8),floor(N/7),floor(N/6),floor(N/5),floor(N/4),floor(N/3),floor(N/2),N]
k=[11,10,9,8,7,6,5,4,3,2,1]
forj=1:
11
m=M(1,j)
g=N/m
fork=1:
g
sum=0
fori=1:
m
sum=sum+w_y((k-1)*m+i,1)
aver=sum/m
end
w(j,k)=aver
end
fori=1:
11
forj=1:
(g-1)
temp=0
temp=temp+(w(i,j+1)-w(i,j))^2
sigma_squa1(i,1)=temp/(2*g-2)
end
end
end
y=sqrt(sigma_squa1)
x=0.02*M
plot(x,y,'o')
holdon
p=polyfit(x,y',3)
q=polyfit(x,y',2)
holdon
plot(x,polyval(p,x),'r')
holdon
plot(x,polyval(q,x),'g')
图像为:
四:
实验结果分析
1.采用A方法得到的各个系数值:
A=Columns1through3
1.785332140245362e+0011.599223756997718e+0011.546731144161181e+001
Columns4through5
1.599223756997741e+0018.767867579144998e-001
xishu1=[QNBKR]'
xishu1=1.030761991763551e+001
2.665372928329531e-001
2.328362402771611e+001
1.661962079874762e+003
4.463869407666557e+003
2.采用拟合函数
(1)分为11组,每组采样结果
M=2932354045536480106160320
(2)
结果:
sigma_squa1=6.958237942450183e-006
4.068707951998810e-006
2.148247678457314e-005
1.673889799999269e-005
2.790504838944067e-005
3.823311607431442e-006
2.455548450315828e-007
7.483690726323448e-009
5.078988766446280e-007
1.475849156097745e-001
0
(3)标准差y=
2.637847217419952e-003
2.017103852556633e-003
4.634919285658935e-003
4.091319835944470e-003
5.282522918969748e-003
1.955329027921245e-003
4.955349079848793e-004
8.650832749697250e-005
7.126702439730649e-004
3.841678221946426e-001
0
五:
实验感想
通过这次实验,我更加深入地理解了陀螺性能和物随机误差特性,并切身感受到了matlab的实用性。
在上实验课时,老师仔细地给我们演示了实验的步骤,并解释了数据的含义。
课后做时,在处理巨大的数据量,我感到了前所未有的压力。
之前的matlab有关练习中,我大部分都是借鉴别人的程序,自己并没完整的写过,所以我写这样一个数据处理程序还是花费很多的时间,先是读数据文件不知道选哪种系统,失败了好多次,这个开头就给我很大打击,我有时觉得自己真心不适合编程。
但是经过反思,我觉得自己太过浮躁了,现在这么多的空闲时间不多学一点必然落后别人越来越多。
然后慢慢开始写,遇到不会的XX,自己总结和记忆,最后把结果都得到了。
与此同时,还熟悉了mathtype工具的使用以及origin的使用,虽然都是一些基础的东西,但是在我进行数据分析和处理的过程都起到了很大的作用。
报告中所有的公式都是一个一个自己打出来,觉得好简单的,以前都没有用过。
在老师给的两个数据文件中,我用的是第二个数据文件,去了其中一部分,先是分成五组,用Allan公式得到的无法方差组成5个方程,方程组的解就是我们要找的5个系数。
然后是将数据分成11组,进行解算得到11个方差,再分别进行最小二乘法拟合和三次曲线拟合,得到不同的响应结果。
因为取得组数太少,所以拟合的效果不好,不够圆滑,有些脱离事实,我觉得可以多取几组,得到更多的方差进行拟合,效果会好许多,但是会增加计算量,计算速度大为减慢。
总而言之,我觉得这次经历对我具有很大的启发和激励作用。
我觉得学好编程还需要多加练习吧,很多东西都不懂,但是不去血永远都不会懂。
在解算A系数矩阵时,我弄了好久都不懂,最后才仔细研究才发现p是不可逆,在方程组中就是超定方程组,仿真结果一直都是参数错误什么的。
然后在进行函数拟合时,变量后面少了小数点,系统一直提示矩阵必须是方阵的错误。
六:
附件
1.采样数据文件
2.matlab-m文件
a.求解系数
b.拟合函数