|
楼主 |
发表于 2015-3-15 23:56:26
|
显示全部楼层
%-- 2015/3/15 23:55 --%
start = 200000;
length =10000000;
ender =start+length;
dt(1:ender,1)=data(2:ender+1,4)-data(1:ender,4);%deltaT计算
offsetX = 0;offsetY=0;offsetZ=0;
sample=500*10; %500HZ的采样 计算前10S的均值作为零偏
for n=start:start+sample %估计零偏
offsetX=offsetX+data(n,1);
offsetY=offsetY+data(n,2);
offsetZ=offsetZ+data(n,3);
end
offsetX=offsetX/sample;offsetY=offsetY/sample;offsetZ=offsetZ/sample;
BiasX=0;BiasY=0;BiasZ=0;allanvarX=0;allanvarY=0;allanvarZ=0;
BX(1:ender,1)=0;BY(1:ender,1)=0;BZ(1:ender,1)=0;
for n=start:ender-1 %对数据进行积分,获取零偏的值
BiasX=BiasX+((data(n,1)-offsetX)*dt(n,1));
BX(n,1)=BiasX;
BiasY=BiasY+((data(n,2)-offsetY)*dt(n,1));
BY(n,1)=BiasY;
BiasZ=BiasZ+((data(n,3)-offsetZ)*dt(n,1));
BZ(n,1)=BiasZ;
allanvarX = allanvarX+( (data(n+1,1) - data(n,1))^2 ) ;%计算艾伦方差
allanvarY = allanvarY+( (data(n+1,2) - data(n,2))^2 ) ;
allanvarZ = allanvarZ+( (data(n+1,3) - data(n,3))^2 ) ;
end
allanvarX=allanvarX/(2*(ender-1));%单位为°/s, °/hr直接*3600
allanvarY=allanvarY/(2*(ender-1));
allanvarZ=allanvarZ/(2*(ender-1));
BiasX=BiasX/10000;%deltaT单位为100us,因此这里需要/10000将零偏变为°/h
BiasY=BiasY/10000;
BiasZ=BiasZ/10000;
time = (data(ender,4)-data(start,4))/36000000;%数据总时长
HrBiasX=BiasX/time;%平均每小时零偏
HrBiasY=BiasY/time;
HrBiasZ=BiasZ/time;
附上带注释的matlab数据处理程序,添加了艾伦方差的计算 |
|