捷联惯导Matlab程序
惯性导航作业
惯性导航作业一、数据说明:1:惯导系统为指北方位的捷连系统。
初始经度为116.344695283度、纬度为39.975172度,高度h为30米。
初速度v0=[-9.993908270;0.000000000;0.348994967]。
2:jlfw中为600秒的数据,陀螺仪和加速度计采样周期分别为为1/100秒和1/100秒。
3:初始姿态角为[2 1 90](俯仰,横滚,航向,单位为度),jlfw.mat中保存的为比力信息f_INSc(单位m/s^2)、陀螺仪角速率信息wib_INSc(单位rad/s),排列顺序为一~三行分别为X、Y、Z向信息.4: 航向角以逆时针为正。
5:地球椭球长半径re=6378245;地球自转角速度wie=7.292115147e-5;重力加速度g=g0*(1+gk1*c33^2)*(1-2*h/re)/sqrt(1-gk2*c33^2);g0=9.7803267714;gk1=0.00193185138639;gk2=0.00669437999013;c33=sin(lat纬度);二、作业要求:1:可使用MATLAB语言编程,用MATLAB编程时可使用如下形式的语句读取数据:load D:\...文件路径...\jlfw,便可得到比力信息和陀螺仪角速率信息。
用角增量法。
2:(1) 以系统经度为横轴,纬度为纵轴(单位均要转换为:度)做出系统位置曲线图;(2) 做出系统东向速度和北向速度随时间变化曲线图(速度单位:m/s,时间单位:s);(3) 分别做出系统姿态角随时间变化曲线图(俯仰,横滚,航向,单位转换为:度,时间单位:s);以上结果均要附在作业报告中。
3:在作业报告中要写出“程序流程图、现阶段学习小结”,写明联系方式。
(注意程序流程图不是课本上的惯导解算流程,而是你程序分为哪几个模块、是怎样一步步执行的,什么位置循环等,让别人根据该流程图能够编出相应程序) (学习小结按条写,不用写套话) 4:作业以纸质报告形式提交,附源程序。
基于Matlab的SINS算法仿真
• 95•为了使实际SINS(捷联惯导系统)算法仿真程序编写更加方便,列写了一些与SINS 更新算法有关的函数计算公式。
对SINS 算法进行了仿真,仿真结果与理论分析结论吻合。
1.引言SINS 算法的核心部分是姿态更新解算,由于四元数算法的优良特性,它在工程上常被采用。
本文总结了与SINS 算法有关的函数计算公式,对SINS 算法进行了仿真,为SINS/GPS 组合导航的研究打下基础。
2.SINS算法本文选用地理坐标系为导航坐标系,记为n 系,载体坐标系记为b 系。
2.1 SINS更新算法SINS 更新的基本理念是,把前一时刻的姿态、速度和位置等导航参数作为初值,利用前一时刻至当前时刻的惯性器件采样输出,解算此时刻的导航参数,作为下一时刻SINS 解算初值,如此反复。
惯性器件采样经过误差补偿后获得等效旋转矢量和比力速度增量,再经过以下三步骤便可实现SINS 更新,计算公式为:(1)速度更新算法位置更新算法(3)姿态更新算法基于Matlab 的SINS算法仿真沈阳理工大学 王海生图1 陀螺漂移对导航精度的影响• 96•图2 捷联惯导算法仿真2.2 SINS速度与位置误差分析速度与位置是积分关系,二者紧密相连,在误差上也是递推的。
速度误差的递推关系式是:位置误差:3.SINS仿真与分析以ENU 坐标系为导航坐标系,在41.5纬度下仿真东、北、天三方向的初始姿态误差、初始位置误差和加速度计零偏误差对导航精度的影响,结果如图1所示。
由仿真图可知,东向速度精度受初始姿态偏差北向投影和初始北向零偏影响较大;经度精度受初始东向零偏影响较大;纬度精度受初始北向零偏和天向零偏影响较大。
而且由于地球自转和地球加速度的影响,速度和姿态都会发生略微偏差。
分别把初始平台的误差角设置为,对静态导航进行仿真,仿真1h 导航参数结果如图2所示。
图2中,速度误差和位置误差小图中,高度通道的误差为红色发散曲线,而水平平台误差、经纬度和水平速度误差呈现振荡趋势,方位平台误差变化比较小。
基于Matlab的低成本捷联惯导的仿真系统研究
基于Matlab的低成本捷联惯导的仿真系统研究作者:刘智平,邵小兵来源:《现代电子技术》2010年第12期摘要:低成本捷联惯导系统的研制和开发已成为导航领域的主要发展趋势之一,在选定了低成本的惯性器件后,系统预期的性能指标和最后样机的性能测试结果之间往往存在着一定的差距,如何降低开发成本并缩短开发时间,就成了低成本捷联惯导系统研制中的一个亟待解决的问题。
这里基于Matlab设计了一个仿真系统,将捷联惯导的软件和硬件模块化,该系统不仅可以完成对系统性能的预测和评估,也可完成对算法和器件选择的优化。
关键词:Matlab; 陀螺; 捷联惯导系统; 仿真系统中图分类号:TP274 文献标识码:A文章编号:1004-373X(2010)12-0078-03Simulation System of Low Cost Strapdown Inertial Navigation System Based on MatlabLIU Zhi-ping, SHAO Xiao-bing(Co llege of Computer, Xi’an Technological University, Xi’an 710032,China)Abstract:The research and development of low cost strapdown inertial navigation system has been the tendency in the field of navigation. The discrepancy between the expected performance and the test data was inevitable with the usage of inertial sensors. How to reduce the cost and shorten the development period is a hot issue in the development of low cost strapdown inertial navigation system.A simulation system based on MATLAB was designed and constructed, the software and hardware of strapdown inertial navigation were modularized, the optimization of various algorithms and selection of sensors should be more convenient and feasible.Keywords:Matlab; gyro; strapdown inertial navigation; simulation system0 引言捷联惯导系统由于体积小、成本低、可靠性好而大量应用在航空、航天、航海、兵器等领域[1-3]。
基于MATLAB_Simulink的捷联惯性导航系统仿真
设计与应用计算机测量与控制.2008.16(8) Computer Measurement &Control#1161#收稿日期:2007211206; 修回日期:2007212225。
基金项目:总装备部十一五惯性技术基金(51309060401)。
作者简介:齐 鑫(19822),男,陕西西安人,硕士研究生,主要从事惯性导航及检测技术方面的研究。
秦永元(19462),男,江苏省常熟市人,教授,博士生导师,主要从事惯性导航与容错组合导航系统,最优估计理论,GPS 信号处理的研究。
文章编号:167124598(2008)0821161203 中图分类号:T P27315;TH 873127 文献标识码:A基于MATLAB/Simulink 的捷联惯性导航系统仿真齐 鑫,秦永元,朱新颖,张芳芳(西北工业大学自动化学院,陕西西安 710072)摘要:介绍了仿真软件MAT LAB /Simulin k 及先进的仿真系统平台软件RT-Lab,设计了龙格-库塔积分模块,用于解决在使用定步长求解器及数字时钟情况下的积分精度不高的问题;建立了基于MAT LAB/Simu link 的捷联惯性导航仿真模型;并将系统分解成一系列功能相对独立的模块,如轨迹发生模块,捷联惯导解算模块;并进行了仿真验证,得到的姿态、速度、位置误差与理论值相符,证明了仿真模型建立正确、方法采用得当、有效。
关键词:RT-LAB;MAT LAB/Simu link;捷联惯导;仿真Simulation on Strap -down Inertial NavigationSystem Based on MATLAB/SimulinkQi Xin,Qin Yongyuan,Zhu Xinying,Zhang Fangfang(Automatization College of Northwest P olytechnic univer sity,Xi'an 710072,China)Abstract:Th e advanced simulation s ystem platform RT-LAB and software MATL AB /Simulin k are ex pou nded,th e Rung e-Kutta in 2tegrator models is designed ,and strap-d own inertial navigation simu lation models bas ed on M ATLAB/Sim ulink is built 1Th e system was divided in to a series of independen t functional units,such as track generation unit,strap-down inertial navigation unit 1And take the vali 2date,Comparing with the real data,the results sh ow that the m odel is sound and the m ethod is validity 1Key word s :RT-LAB;MATL AB /Simulin k;s trap-down inertial navigation;simulation0 引言MATLAB 是美国MathWorks 公司推出的功能强大的大型数学软件,它以优秀的数值计算和卓越的数据可视化能力被广泛应用于自动控制、数字信号处理、时间序列分析、动态系统仿真等各个领域。
惯性导航四元数法
基于Matlab 的捷联惯导算法设计及仿真1严恭敏西北工业大学航海学院,西安 (710072)E-mail :yangongmin@摘 要:根据圆锥误差补偿算法和划船误差补偿算法的研究成果,考虑到实际捷联惯导算法仿真程序编写的方便性,总结了一些与捷联惯导更新算法有关的函数的计算公式。
对圆锥误差补偿算法和捷联惯导算法进行了仿真,仿真结果和理论分析结论吻合。
在附录中给出了Matlab 的m 文件源程序代码,具有一定的参考价值。
关键词:捷联惯导;四元数;等效旋转矢量;Matlab ;算法;仿真 中图分类号:V249.31. 引言在捷联惯导系统中采用数学平台,姿态更新解算是捷联惯导系统算法的核心部分,由于四元数法算的优良特性,它在工程实际中经常被采用。
为了减小姿态计算的不可交换性误差,前人研究并建立了等效旋转矢量方程,高精度姿态更新解算的研究主要集中在等效旋转矢量方程的求解上,在圆锥运动环境下,许多研究者提出并完善了圆锥误差补偿算法。
基于圆锥误差补偿算法和划船误差补偿算法的等效原理,可将圆锥误差补偿算法移植到划船误差补偿算法中去,从而减少了划船误差推导的繁琐过程。
上述研究都已经比较成熟[1-6],本文根据这些研究结果,并考虑到实际仿真程序编写的方便性,总结了一些与捷联惯导算法有关的函数的计算公式或步骤,其中更详细的推导过程可见参考文献[7,8]。
最后,对圆锥误差补偿算法和捷联惯导算法进行了仿真。
附录中给出了Matlab 的m 文件源程序代码具有一定的参考价值。
2. 捷联惯导算法文中选取东-北-天(E-N-U )地理坐标系为导航坐标系,记为n 系;捷联惯组坐标系记为b 系。
2.1相关函数(1)四元数的共轭与乘积。
四元数q 可表示为k q j q i q q q v 32100+++=+=q q 。
用∗q表示q 的共轭四元数;21q q q =表示四元数q 是四元数1q 与2q 的乘积,四元数相乘是不可交换的。
imu惯导matlab解算
imu惯导matlab解算
在MATLAB中解算惯性导航(IMU)可以使用以下步骤:
1. 准备数据:你需要收集IMU传感器的数据,包括加速度计和陀螺仪的测量值。
这些数据通常以时间序列的形式提供。
2. 读取数据:使用MATLAB的文件读取函数(例如`readmatrix`)将数据读取到MATLAB的工作空间中。
3. 数据预处理:IMU数据通常需要进行预处理,例如去除噪声、校准传感器和对数据进行滤波。
可以使用MATLAB的信号处理工具箱进行这些预处理步骤。
4. 设计滤波器:在IMU数据预处理之后,可能需要设计一个滤波器来进一步减少噪声和提高数据质量。
MATLAB的信号处理工具箱提供了多种滤波器设计方法。
5. 解算姿态:使用IMU的加速度计和陀螺仪测量值,可以通过解算算法来估计姿态。
常见的解算算法包括加速度计积分法、四元数解算、互补滤波器等。
可以使用MATLAB中的数学函数和算法来实现这些解算算法。
6. 可视化结果:可以使用MATLAB的绘图函数(例如`plot`)将解算结果可视化,以便更好地理解IMU的姿态估计。
请注意,解算IMU的精确性和准确性可能需要对算法进行调优,并进行对比和评估。
IMU还可能受到环境因素和传感器误差的影响,因此在解算IMU数据时需要考虑到这些因素。
捷联惯导的解算程序
%====本程序为捷联惯导的解算程序(由惯性器件的输出解算出飞行器的位置、速度、姿态信息)======clear all;close all;clc;deg_rad=pi/180; %由度转化成弧度rad_deg=180/pi; %由弧度转化成度%-------------------------------从源文件中读入数据----------------------------------fid_read=fopen('IMUout.txt','r'); %path1_Den.dat 是由轨迹发生器产生的数据[AllDataNumofAllData]=fscanf(fid_read,'%g %g %g %g %g %g %g %g %g %g %g %g %g %g %g %g',[17 inf]);AllData=AllData';NumofEachData=round(NumofAllData/17);Time=AllData(:,1);longitude=AllData(:,2); %经度单位:弧度latitude=AllData(:,3); %纬度单位:弧度High=AllData(:,4); %高度单位:米Ve=-AllData(:,6); % 东向、北向、天向速度单位:米/妙Vn=AllData(:,5);Vu=AllData(:,7);fb_x=AllData(:,9); %比力(fx,fy,fz)fb_y=AllData(:,8); %指向右机翼方向为x正方向,指向机头方向为y正向,z轴与x轴和y轴构成右手坐标系单位:米/秒2fb_z=-AllData(:,10); %右前上pitch=AllData(:,11); %俯仰角(向上为正)单位:弧度head=-AllData(:,13); %偏航角(偏西为正)roll=AllData(:,12); %滚转角(向右为正)omigax=AllData(:,15); %陀螺输出(单位:弧度/秒,坐标轴的定义与比力的相同)omigay=AllData(:,14);omigaz=-AllData(:,16);%-------------------------------程序初始化--------------------------------------latitude0=latitude(1);longitude0=longitude(1); %初始位置High0=High(1);Ve0=Ve(1);Vn0=Vn(1); %初始速度Vu0=Vu(1);pitch0=pitch(1);head0=head(1); %初始姿态roll0=roll(1);TimeEach=0.005; %周期和仿真总时间TimeAll=(NumofEachData-1)*TimeEach;Omega_ie=0.7292115147E-4;%0.00007272205216643040; %地球自转角速度单位:弧度每妙g0=9.78;%------------------------------导航解算开始--------------------------------------%假设没有初始对准误差pitch_err0=pitch0+0*deg_rad;head_err0=head0+0*deg_rad;roll_err0=roll0+0*deg_rad;%初始捷联矩阵的计算《捷联惯导系统》P63 旋转顺序 head - pitch - roll %导航坐标系n为东北天方向载体坐标系b为右前上偏航角北偏西为正Tbn(1,1)=cos(roll_err0)*cos(head_err0)-sin(roll_err0)*sin(pitch_err0) *sin(head_err0);Tbn(1,2)=cos(roll_err0)*sin(head_err0)+sin(roll_err0)*sin(pitch_err0) *cos(head_err0);Tbn(1,3)=-sin(roll_err0)*cos(pitch_err0);Tbn(2,1)=-cos(pitch_err0)*sin(head_err0);Tbn(2,2)=cos(pitch_err0)*cos(head_err0);Tbn(2,3)=sin(pitch_err0);Tbn(3,1)=sin(roll_err0)*cos(head_err0)+cos(roll_err0)*sin(pitch_err0) *sin(head_err0);Tbn(3,2)=sin(roll_err0)*sin(head_err0)-cos(roll_err0)*sin(pitch_err0) *cos(head_err0);Tbn(3,3)=cos(roll_err0)*cos(pitch_err0);Tnb=Tbn';%位置矩阵的初始化《捷联惯导系统》P46 其中游动方位角 a=0 假使初始经纬度确知Cne(1,1) = - sin(longitude0);Cne(1,2) = cos(longitude0);Cne(1,3) = 0;Cne(2,1) = - sin(latitude0) * cos(longitude0);Cne(2,2) = - sin(latitude0) * sin(longitude0);Cne(2,3) = cos(latitude0);Cne(3,1) = cos(latitude0) * cos(longitude0);Cne(3,2) = cos(latitude0) * sin(longitude0);Cne(3,3) = sin(latitude0);Cen=Cne';%初始四元数的确定《捷联惯导系统》 P151-152 方法本身保证了q1^2+q2^2+q3^2+q4^2=1q(2,1) = sqrt(abs(1.0 + Tnb(1,1) - Tnb(2,2) - Tnb(3,3))) / 2.0; q(3,1) = sqrt(abs(1.0 - Tnb(1,1) + Tnb(2,2) - Tnb(3,3))) / 2.0; q(4,1) = sqrt(abs(1.0 - Tnb(1,1) - Tnb(2,2) + Tnb(3,3))) / 2.0;q(1,1) = sqrt(abs(1.0 - q(2,1) ^2 - q(3,1) ^2 - q(4,1) ^2));% 判断q(1,1)的符号flag_q11=cos(head_err0/2.0)*cos(pitch_err0/2.0)*cos(roll_err0/2.0)-sin(head_err0/2.0)*sin(pitch_err0/2.0)*sin(roll_err0/2.0);if (flag_q11 >0) %此时q(1,1)取正if (Tnb(3,2) < Tnb(2,3))q(2,1) = - q(2,1);endif (Tnb(1,3) < Tnb(3,1))q(3,1) = - q(3,1);endif (Tnb(2,1) < Tnb(1,2))q(4,1) = - q(4,1);endelse %此时q(1,1)取负或0q(1,1) = - q(1,1);if (Tnb(3,2) > Tnb(2,3))q(2,1) = - q(2,1);endif (Tnb(1,3) > Tnb(3,1))q(3,1) = - q(3,1);endif (Tnb(2,1) > Tnb(1,2))q(4,1) = - q(4,1);endend%-------------------------迭代推算用到的参数的初始化------------------------Wiee_e = 0;Wiee_n = 0;Wiee_u = Omega_ie;Wiee = [Wiee_e Wiee_n Wiee_u]'; %地球速率在地球系中的投影东-北-天Lat_err(1)=latitude0;Lon_err(1)=longitude0;High_err(1)=High0;Ve_err(1)=Ve0;Vn_err(1)=Vn0;Vu_err(1)=Vu0;pitch_err(1)=pitch_err0;head_err(1)=head_err0;roll_err(1)=roll_err0;Re=6378137.0;%6378245.0; %地球长轴《惯性导航系统》 P28e=0.0033528106647474807198455286185206; %地球扁率精确值ee=0.00669437999014131699614;%----------------------------迭代推算开始-----------------------------------for i=1:NumofEachData%----------------------------惯性仪表数据的获得------------------------Wibb(1,1)=omigax(i); %指向右机翼方向为x正方向,指向机头方向为y正向,z轴与x轴和y轴构成右手坐标系Wibb(2,1)=omigay(i); %单位:弧度/妙Wibb(3,1)=omigaz(i); %右前上fb(1,1)=fb_x(i); %指向右机翼方向为x正方向,指向机头方向为y正向,z轴与x轴和y轴构成右手坐标系fb(2,1)=fb_y(i); %单位:米/秒2fb(3,1)=fb_z(i); %右前上%--------计算在姿态矩阵和位置矩阵更新时用到的参数------------------ RM=Re*(1.0-2.0*e+3.0*e*Cne(3,3)^2)+High_err(i); %《捷联惯导系统》P233 P235RN=Re*(1.0+e*Cne(3,3)^2)+High_err(i);% RN=Re*(1-ee)/(sqrt(1-ee*sin(Lat_err(i))))^3+High_err(i);% RM=Re/sqrt(1-ee*sin(Lat_err(i)))+High_err(i);%实验当地重力加速度计算《捷联惯导系统》P150 《惯性导航系统》 P35g=g0*((1.0+0.0052884*Cne(3,3)^2)-0.0000059*(1-(1-2*Cne(3,3)^2)^2))*(1 .0-2.0*High_err(i)/Re);tmp_slat=sin(Lat_err(i))*sin(Lat_err(i));Wien = Cne * Wiee; %地球速率在导航系中的投影Wenn(1,1) = -Vn_err(i)/RM;Wenn(2,1) = Ve_err(i)/RN; % <<惯性导航系统>> P45 考虑了地球转动的影响.Wenn(3,1) = Ve_err(i)*tan(Lat_err(i))/RN; %计算Wenn(不太精确),更新速度和位置矩阵时用Winn=Wien+Wenn;Winb=Tbn*Winn;Wnbb=Wibb-Winb; %姿态速率在姿态更新时用到fn=Tnb*fb; % x-y-z 东-北-天% 速度的更新《捷联惯导系统》 P30 33 东-北-天difVe_err=fn(1,1)+(2*Wien(3,1)+Wenn(3,1))*Vn_err(i)-(2*Wien(2,1)+Wenn (2,1))*Vu_err(i);difVn_err=fn(2,1)-(2*Wien(3,1)+Wenn(3,1))*Ve_err(i)+(2*Wien(1,1)+Wenn (1,1))*Vu_err(i);difVu_err=fn(3,1)+(2*Wien(2,1)+Wenn(2,1))*Ve_err(i)-(2*Wien(1,1)+Wenn (1,1))*Vn_err(i)-g;Ve_err(i+1)=Ve_err(i)+difVe_err*TimeEach;Vn_err(i+1)=Vn_err(i)+difVn_err*TimeEach;Vu_err(i+1)=Vu_err(i)+difVu_err*TimeEach;High_err(i+1)=High_err(i)+Vu_err(i)*TimeEach;% 位置矩阵的实时更新《惯性导航系统》 P190Cne(1,1)=Cne(1,1)+TimeEach*(Wenn(3,1)*Cne(2,1)-Wenn(2,1)*Cne(3,1)); Cne(1,2)=Cne(1,2)+TimeEach*(Wenn(3,1)*Cne(2,2)-Wenn(2,1)*Cne(3,2)); Cne(1,3)=Cne(1,3)+TimeEach*(Wenn(3,1)*Cne(2,3)-Wenn(2,1)*Cne(3,3)); Cne(2,1)=Cne(2,1)+TimeEach*(-Wenn(3,1)*Cne(1,1)+Wenn(1,1)*Cne(3,1)); Cne(2,2)=Cne(2,2)+TimeEach*(-Wenn(3,1)*Cne(1,2)+Wenn(1,1)*Cne(3,2)); Cne(2,3)=Cne(2,3)+TimeEach*(-Wenn(3,1)*Cne(1,3)+Wenn(1,1)*Cne(3,3)); Cne(3,1)=Cne(3,1)+TimeEach*(Wenn(2,1)*Cne(1,1)-Wenn(1,1)*Cne(2,1)); Cne(3,2)=Cne(3,2)+TimeEach*(Wenn(2,1)*Cne(1,2)-Wenn(1,1)*Cne(2,2)); Cne(3,3)=Cne(3,3)+TimeEach*(Wenn(2,1)*Cne(1,3)-Wenn(1,1)*Cne(2,3));% Mat_Wenn(1,1)=0;% Mat_Wenn(1,2)=Wenn(3,1);% Mat_Wenn(1,3)=-Wenn(2,1); %Wenn的反对阵矩阵取负% Mat_Wenn(2,1)=-Wenn(3,1); %这里位置矩阵的及时修正为:dCne/dt=Mat_Wenn*Cne% Mat_Wenn(2,2)=0;% Mat_Wenn(2,3)=Wenn(1,1);% Mat_Wenn(3,1)=Wenn(2,1);% Mat_Wenn(3,2)=-Wenn(1,1);% Mat_Wenn(3,3)=0;%% Mat_Wenn=Mat_Wenn*Cne*TimeEach;% Cne=Cne+Mat_Wenn;Cen=Cne';% 计算经纬度Lat_err(i+1)=asin(Cne(3,3));Lon_err(i+1)=atan(Cne(3,2)/Cne(3,1)); %这是经度的主值if (Cne(3,1) < 0)if (Lon_err(i+1) > 0)Lon_err(i+1) = Lon_err(i+1) - pi;elseLon_err(i+1) = Lon_err(i+1) + pi;endend% 四元数的及时修正《惯性导航系统》 P194% Mat_Wnbb=[ 0, -Wnbb(1,1), -Wnbb(2,1), -Wnbb(3,1); % Wnbb(1,1), 0, Wnbb(3,1), -Wnbb(2,1); % Wnbb(2,1), -Wnbb(3,1), 0, Wnbb(1,1); % Wnbb(3,1), Wnbb(2,1), -Wnbb(1,1), 0];% q=q+Mat_Wnbb*q*TimeEach/2.0;q(1,1)=q(1,1)+TimeEach*(-Wnbb(1,1)*q(2,1)-Wnbb(2,1)*q(3,1)-Wnbb(3,1)* q(4,1))/2.0;q(2,1)=q(2,1)+TimeEach*(Wnbb(1,1)*q(1,1)+Wnbb(3,1)*q(3,1)-Wnbb(2,1)*q (4,1))/2.0;q(3,1)=q(3,1)+TimeEach*(Wnbb(2,1)*q(1,1)-Wnbb(3,1)*q(2,1)+Wnbb(1,1)*q (4,1))/2.0;q(4,1)=q(4,1)+TimeEach*(Wnbb(3,1)*q(1,1)+Wnbb(2,1)*q(2,1)-Wnbb(1,1)*q (3,1))/2.0;% 四元数归一化处理q_norm=sqrt(sum(q.*q));q=q/q_norm;% 计算姿态矩阵 TnbTnb(1,1) = q(1,1) ^2 + q(2,1) ^2 - q(3,1)^2 - q(4,1)^2;Tnb(1,2) = 2.0 * (q(2,1) * q(3,1) - q(1,1) * q(4,1));Tnb(1,3) = 2.0 * (q(2,1) * q(4,1) + q(1,1) * q(3,1));Tnb(2,1) = 2.0 * (q(2,1) * q(3,1) + q(1,1) * q(4,1));Tnb(2,2) = q(1,1)^2 - q(2,1)^2 + q(3,1)^2 - q(4,1)^2;Tnb(2,3) = 2.0 * (q(3,1) * q(4,1) - q(1,1) * q(2,1));Tnb(3,1) = 2.0 * (q(2,1) * q(4,1) - q(1,1) * q(3,1));Tnb(3,2) = 2.0 * (q(3,1) * q(4,1) + q(1,1) * q(2,1));Tnb(3,3) = q(1,1)^2 - q(2,1)^2 - q(3,1)^2 + q(4,1)^2;Tbn=Tnb';flag_pitch=asin(Tnb(3,2));flag_roll=atan(-Tnb(3,1)/Tnb(3,3));flag_head=atan(-Tnb(1,2)/Tnb(2,2));if(Tnb(3,3)<0)if(flag_roll<0)flag_roll=flag_roll+pi;endif(flag_roll>0)flag_roll=flag_roll-pi;endend% 偏航角范围 -180度——180度北偏西为正if(Tnb(2,2)<0)if(flag_head<0)flag_head=flag_head+pi;endif(flag_head>0)flag_head=flag_head-pi;endend% 姿态角更新pitch_err(i+1)=flag_pitch;head_err(i+1)=flag_head;roll_err(i+1)=flag_roll;% 解算完毕由对准结果、陀螺、加表的输出解算出载体的位置、速度、姿态%----------------------计算解算误差------------------ddLat(i)=(Lat_err(i)-latitude(i))*rad_deg; %纬度误差单位:度ddLog(i)=(Lon_err(i)-longitude(i))*rad_deg; %经度误差单位:度ddHigh(i)=High_err(i)-High(i); %高度误差单位:米ddVe(i)=Ve_err(i)-Ve(i);ddVn(i)=Vn_err(i)-Vn(i); % 速度误差单位:米/妙2ddVu(i)=Vu_err(i)-Vu(i);ddpitch(i)=(pitch_err(i)-pitch(i))*rad_deg*3600; %姿态误差单位:度ddhead(i)=(head_err(i)-head(i))*rad_deg*3600;ddroll(i)=(roll_err(i)-roll(i))*rad_deg*3600;endfclose(fid_read);%---------------------------绘图开始--------------------------------- figure(1)plot(Time,ddLog)ylabel('经度误差(度)'),xlabel('时间(秒)');figure(2)plot(Time,ddLat)ylabel('纬度误差(度)'),xlabel('时间(秒)');figure(3)plot(Time,ddHigh);ylabel('高度误差(米)'),xlabel('时间(秒)');figure(4)plot(Time,ddhead)ylabel('偏航角误差(角妙)'),xlabel('时间(秒)'); figure(5)plot(Time,ddpitch)ylabel('俯仰角误差(角妙)'),xlabel('时间(秒'); figure(6)plot(Time,ddroll);ylabel('滚转角误差(角妙)'),xlabel('时间(秒)'); figure(7)plot(Time,ddVe);ylabel('东向速度误差(米/秒)'),xlabel('时间(秒)'); figure(8)plot(Time,ddVn)ylabel('北向速度误差(米/秒)'),xlabel('时间(秒)'); figure(9)plot(Time,ddVu)ylabel('天向速度误差(米/秒)'),xlabel('时间(秒)');%------------------------------绘图结束-------------------------------。
基于MATLAB_Simulink的捷联惯性导航仿真研究
基于MATLAB_Simulink的捷联惯性导航仿真研究
蓝仁恩; 付战平
【期刊名称】《《军民两用技术与产品》》
【年(卷),期】2007(000)011
【摘要】利用MATLAB数值计算的强大功能和Simulink建模容易实现的特点,进行了捷联式惯性导航系统的仿真研究。
对四元数及姿态矩阵等的解算采用MATLAB编程实现,陀螺仪和加速计的建模用Simulink实现,将MATLAB程序与Simulink模型进行交互,得出仿真结果。
该方法降低了单纯应用MATLAB或Simulink进行仿真的复杂度,使仿真既可靠又易于实现。
【总页数】2页(P45-46)
【作者】蓝仁恩; 付战平
【作者单位】海军驻航天三院代表室北京100074; 海军航空工程学院青岛分院青岛266041
【正文语种】中文
【中图分类】V249.322
【相关文献】
1.基于Simulink的海底捷联惯性导航系统数学仿真 [J], 余洁;杨平
2.基于MATLAB/Simulink的捷联惯性导航系统仿真 [J], 齐鑫;秦永元;朱新颖;张芳芳
3.基于捷联惯性导航载体轨迹发生器的设计与仿真 [J], 罗宇锋;刘勇
4.基于MATLAB_Simulink的捷联惯性导航仿真研究 [J], 蓝仁恩; 付战平
5.基于MATLAB组件GUI的捷联惯性导航系统教学仿真平台设计 [J], 周雪梅;许德新;张晓宇;郭立东
因版权原因,仅展示原文概要,查看原文内容请购买。
一种基于Matlab的捷联惯导系统仿真轨迹发生器设计
一种基于Matlab的捷联惯导系统仿真轨迹发生器设计
文钢
【期刊名称】《舰船电子工程》
【年(卷),期】2016(036)005
【摘要】在捷联惯导系统的算法设计、精度仿真等前期验证工作中,通过软件编写的轨迹发生器程序发挥着降低试验成本,减少后期试验风险的重要作用.论文针对惯导系统的各类误差特点,提出了一种基于Matlab的捷联惯导系统仿真轨迹发生器.该轨迹发生器可以模拟生成惯导系统载体的各种姿态轨迹数据,另外还能在数据中加入各类惯导器件常见误差及噪声等,有效提高其仿真度.
【总页数】6页(P87-91,155)
【作者】文钢
【作者单位】海装重庆局重庆 400042
【正文语种】中文
【中图分类】TN219
【相关文献】
1.基于解析式轨迹的捷联惯导系统仿真研究 [J], 李建文;黄国荣;张宗麟
2.一种用于捷联惯导系统仿真的飞行轨迹发生方法 [J], 胡恒章;周荻;胡国辉
3.基于捷联惯导的舰船轨迹仿真发生器设计 [J], 张斌
4.基于捷联惯导的舰船轨迹仿真发生器设计 [J], 张斌
5.基于Matlab/Simulink的捷联惯导系统仿真运动轨迹生成器 [J], 谭红力;鲁小强;黄新生
因版权原因,仅展示原文概要,查看原文内容请购买。
捷联式惯性导航系统
1 绪论随着计算机和微电子技术的迅猛发展,利用计算机的强大解算和控制功能代替机电稳定系统成为可能。
于是,一种新型惯导系统--捷联惯导系统从20世纪60年代初开始发展起来,尤其在1969年,捷联惯导系统作为"阿波罗"-13号登月飞船的应急备份装置,在其服务舱发生爆炸时将飞船成功地引导到返回地球的轨道上时起到了决定性作用,成为捷联式惯导系统发展中的一个里程碑。
捷联式惯性导航(strap-down inertial navigation),捷联(strap-down)的英语原义是“捆绑”的意思。
因此捷联式惯性导航也就是将惯性测量元件(陀螺仪和加速度计)直接装在飞行器、舰艇、导弹等需要诸如姿态、速度、航向等导航信息的主体上,用计算机把测量信号变换为导航参数的一种导航技术。
现代电子计算机技术的迅速发展为捷联式惯性导航系统创造了条件。
惯性导航系统是利用惯性敏感器、基准方向及最初的位置信息来确定运载体的方位、位置和速度的自主式航位推算导航系统。
在工作时不依赖外界信息,也不向外界辐射能量,不易受到干扰破坏。
它完全是依靠载体自身设备独立自主地进行导航,它与外界不发生任何光、声、磁、电的联系,从而实现了与外界条件隔绝的假想的“封闭”空间内实现精确导航。
所以它具有隐蔽性好,工作不受气象条件和人为的外界干扰等一系列的优点,这些优点使得惯性导航在航天、航空、航海和测量上都得到了广泛的运用[1]1.1 捷联惯导系统工作原理及特点惯导系统主要分为平台式惯导系统和捷联式惯导系统两大类。
惯导系统(INS)是一种不依赖于任何外部信息、也不向外部辐射能量的自主式导航系统,具有隐蔽性好,可在空中、地面、水下等各种复杂环境下工作的特点。
捷联惯导系统(SINS)是在平台式惯导系统基础上发展而来的,它是一种无框架系统,由三个速率陀螺、三个线加速度计和微型计算机组成。
平台式惯导系统和捷联式惯导系统的主要区别是:前者有实体的物理平台,陀螺和加速度计置于陀螺稳定的平台上,该平台跟踪导航坐标系,以实现速度和位置解算,姿态数据直接取自于平台的环架;后者的陀螺和加速度计直接固连在载体上作为测量基准,它不再采用机电平台,惯性平台的功能由计算机完成,即在计算机内建立一个数学平台取代机电平台的功能,其飞行器姿态数据通过计算机计算得到,故有时也称其为"数学平台",这是捷联惯导系统区别于平台式惯导系统的根本点。
matlab惯导标定算法
matlab惯导标定算法
在MATLAB中进行惯导标定通常涉及以下几个步骤:
1. 数据采集:首先,你需要采集原始的惯导数据。
这通常涉及到使用IMU (惯性测量单元)设备来获取加速度计和陀螺仪的读数。
2. 数据预处理:在分析之前,你可能需要对数据进行一些预处理,例如滤波、去噪或平滑。
3. 算法实现:然后,你可以开始编写或应用一个标定算法。
以下是一个简化的标定流程示例,该流程可能包括:
姿态估计:使用加速度计和陀螺仪数据来估计设备的姿态(方向)。
误差计算:基于估计的姿态和已知的参考姿态,计算误差。
参数优化:使用优化算法(如梯度下降法)来调整IMU的内部参数(例如陀螺仪的偏置和标度因子),以最小化误差。
4. 验证与评估:完成算法后,需要在实际数据或模拟数据上验证其性能。
这可能涉及比较算法输出的姿态与已知的参考姿态,计算误差等。
5. 后处理:根据验证结果,可能需要对算法进行微调或改进。
请注意,具体的实现细节将取决于你的应用需求、IMU的特性和性能以及你的算法目标。
MATLAB是一个强大的工具,可以用来开发、测试和验证各种惯导标定算法。
如果你需要更具体的帮助,例如代码示例或特定算法的实现,请提供更多详细信息,我会很乐意帮助你。
基于MATLAB.Simulink的捷联惯性导航仿真
基于!"#$"%&’()*+(,-的捷联惯性导航仿真’./01234,5,6/.(0+708(90.(3,’()*+0.(3,:0;623,!"#$"%&’()*+(,-练军想杨壮志胡小平<国防科技大学长沙=>??@A BC 摘要D 阐述了EF G H F I &J K L M N K O P 的主要特点Q 设计了基于EF G H F I &J K L M N K O P 的捷联惯性导航仿真模型Q 并给出了仿真实例和结果Q 表明利用EF G H F I &J K L M N K O P 可以提高仿真效率R C关键词D EF G H F I &J K O L M N K O PQ 仿真Q 捷联惯性导航Q 传递对准"%’#S "T #U OV W K X Y Z Y [\Q V W [L Z K O]W Z \Z ]V [\K X V K ]X ^_EF G H F I &J K L M N K O PZ \[K N N M X V \Z V [‘a F X V \Z Y ‘^b OK O [\V K Z N O Z c K d Z V K ^OEF G e H F I &J K L M N K O PL ^‘[N K X‘[X K d O [‘Q Z O ‘V W [X K L M N Z V K ^O[f Z L Y N [Z O ‘\[X M N V Z \[d K c [OZ N X ^a F X K L M N Z V K ^O[__K ]K [O ]g_\^L EF G H F I&J K L M N K O P]Z Oh [^h V Z K O [‘a i j k lm S n ’EF G H F I &J K L M N K O P Q X K L M N Z V K ^O Q X V \Z Y ‘^b OK O [\V K Z N O Zc Kd Z V K ^O Q V \Z O X _[\Z N K d O L [O V仿真即建立系统的模型Q 并在模型上进行试验R 长期以来Q 仿真领域的研究重点是如何建立合理的仿真计算模型R 为了获得有用的仿真结果Q 研究者不得不将大量的精力投入到仿真算法的设计和仿真程序的编制o 调试上Q 而这一过程非常复杂和繁琐Q 且多为重复性劳动R 这在很大程度上阻碍了仿真技术的推广和应用R 近年来各种仿真软件的出现Q 上述不利情况得到了缓解R EF G H F I &J K L M N K O P 无疑是众多仿真软件中最成功的一个R 在EF G H F I &J K L M N K O P 中Q 建模过程变得非常简单Q 仿真过程可交互Q 还可以利用EF G H F I 各种工具对仿真结果进行分析和可视化处理等R 使用EF G H F I &J K L M N K O P 这一强大的仿真软件Q 可以轻松构建捷联惯性导航系统仿真模型p 使得对各种姿态算法Q 对准方案的分析研究大大简化q 也便于用户理解和掌握惯性导航的原理o 方法q 还可利用r G s<r [Z N G K L [s^\P X K W ^Y B 生成优化的高级语言代码Q 直接用于实时控制和快速原型设计Rt EF G H F I &J K L M N K O P 及其特点EF G H F I &J K L M N K O P 最大的优点是将用户从许多重复的代码编写工作中解脱出来Q让用户把精力从编程转向模型的构造o 试验以及结果的分析处理R EF G e H F I &J K L M N K O P 提供基本的功能模块Q 存放在J K L M N K O P 库中R 虽然每个模块对用户而言是完全透明的Q 但用户在使用时只需要知道模块的输入o 输出和模块的功能Q 而不必考虑模块内部具体如何实现的R 对于建模Q 用户所需完成的工作就是将这些基本模块根据自己的需要连接起来R 这就好比用户需要一台电脑Q 没有必要从单个的集成电路芯片开始设计Q 只需将厂商提供的主板o u v w o 显卡o 声卡o 电源等设备合理组装起来一样R 搭建好仿真模型后Q 用户就可以设定参数o进行仿真o 分析结果了R 对于具体的运行细节o 事件如何驱动等Q 用户不必关心R 总之Q EF G H F I &J K L M N K O P给用户提供一个非常友好的仿真环境REF G H F I &J K L M N K O P 的特点如下px 支持图形用户界面功能R EF G H F I &J K L M N K O P 的图形用户界面功能使建模的过程简化为鼠标的点击和拖拽Q 用户几乎不需要编写任何代码就能完成整个模型的构建R 这与编程语言的仿真软件包以及传统的EF G H F I 命令方式相比Q 既简单Q 又直观R y 模型层次性强R EF G H F I &J K L M N K O P 支持子系统功能Q 用户可以将联系比较紧密的几个模块<子系统B 组合成一个子系统<更高级子系统B Q 这样留在顶层的只有几个功能相对独立的高层子系统R 高层子系统与下一级子系统<模块B 之间的切换也很简单Q 双击鼠标左键或单击按钮即可R 此外Q 用户还可以通过点击工具栏中的按钮查看呈树型结构排列的各子系统的内嵌关系R 这使得用户对整个模型的结构一目了然Q 便于模型的设计o 修改Rz 模块封装性好R 对于功能比较独立的通用型模块Q EF G H F I &J K L M N K O P 可以对其进行封装R 这有利于将封装模块内部的各参数集中到一个对话框统一设置Q 还可避免其它使用者对模块不可预知的改动R{模块的可移植性强R功能相同或相近的模块Q 仅|}??}?~}?收到Q }??}?!}?改回||练军想Q 男Q >"@!年生Q在读硕士研究生Q 研究方向p 导航o 制导与控制R #}}#<总="}B 基于EF G H F I &J K L M N K O P 的捷联惯性导航仿真}??}年E=D F >G 取当地水平地理坐标系HH 北天东系作为导航坐标系"54比力"74速度">4角速度"C 4重力加速度"下标?94地球系相对惯性系"下标9:4导航系相对地球系"D 4姿态四元数">G 4载体相对于导航系的角速度"F4四元数相乘!$2$仿真模型的设计依据%&’(&)*+,-./,01作为仿真开发平台"利用其优秀的图形用户界面功能"按照自顶向下的顺序构建仿真模型"将整个捷联惯性导航模型分为四大模块4数据源#轨迹发生器#导航解算和结果比较"如图E 所示!数据源提供速度#加速度和姿态信息I 轨迹发生器依据理想的速度#加速度和姿态生成比力5和角速度>"作为加速度计和陀螺的采样值送给导航解算模块I 导航解算模块再根据比力方程和姿态微分方程解算出导航参数HH 速度#位置#姿态和高度I 最后"在结果比较模块中将理想的导航参数与导航解算得到的导航参数作比较"对仿真结果进行分析和评估!$2J 仿真应用实例K 引入陀螺和加速度计误差"观察导航解算精度!仿真结果表明"如果三个陀螺的常值漂移均为E L *M"仿真E M 后"位置偏差可达E N N 2E 1-!O 捷联惯导姿态解算模块分别采用单子样#双子样#三子样#四子样旋转矢量算法"评估子样数对姿态解算精度的影响!仿真结果表明"不考虑惯性器件误差和采样过程的量化误差时"遵循子样数越多精度越高的规律I 如果引入量化误差"过分增加子样数并不能提高精度!P 利用%&’(&)*+,-./,01模块可移植性强的特点"将图E 所示的惯性导航仿真模型稍作修改"扩充为主#子惯导两套导航解算分系统"再设计一个卡尔曼滤波器仿真模块"就可以进行导航传递对准的仿真研究!根据不同的传递对准匹配方案"分别选取不同的滤波状态和观测量"可以很容易比较各种方案的优劣!仿真结果表明"计算参数匹配的效果优于测量参数匹配I对于舰载武器传递对准"速度;角速率匹配不论是对准精度还是在对准速度上"明显优于单独的速度匹配或单独的角速率匹配!采用速度;角速率匹配传递对准"设置仿真的初始条件4初始位置4北纬Q R L "东经E E R L初始速度<-*S @47T 67U 67V 6R加速度<-*S =@4W T 6R 2=X Y S R 2E Z "W U 6R"W V 6R2E X Y S R 2E Z 主惯导欧拉角<[\]@4^_‘‘6R 2Q S ,0Z "a Wb 6R 2E S ,0Zc ?Zde 6R 2Q S ,0Z臂杆矢量<-@4f W W G 6g E R R R h 陀螺常值漂移4R 2E L *M 加速度计零偏<-*S =@4R 2R R Ri N 初始方位失准角4E L 初始水平失准角4R L上述Z 为时间"单位S !根据以上初始条件"j S 后"失准角可收敛在R 2k 角分以内!J 结束语与传统的编程语言用公式表达微分方程的仿真软件包相比"%&’(&)*+,-./,01可以大大简化建模过程"提高工作效率!而且"所建模型易于修改#移植和扩充"仿真结果便于分析#处理!利用%&’(&)*+,-./,01对捷联惯性导航系统进行仿真"得到了较好的仿真结果!参考文献E 薛定宇"陈阳泉2基于%&’(&)*+,-./,01的系统仿真技术与应用2北京4清华大学出版社"=R R =4E i =l==Q =范影乐"杨胜天"李轶2%&’(&)仿真应用详解2北京4人民邮电出版社"=R R E 4全书Q 陈桂明"张明照"戚红雨等2应用%&’(&)建模与仿真2北京4科学出版社"=R R E 4全书k 袁信"俞济祥"陈哲2导航系统2北京4航空工业出版社"E i i Q 4E m ln ioQ =o 第E m 卷第E R 期电脑开发与应用<总k i Q @。
基于MATLAB的捷联惯组测量数据分析平台
基于MATLAB的捷联惯组测量数据分析平台作者:陈春歌宋玉珍来源:《数字技术与应用》2020年第09期摘要:针对某型飞行器的运动轨迹参数的计算问题,建立了基于MATLAB的捷联惯组测量数据分析平台。
该平台采用MATLAB GUI图形化编程语言编制,功能设计兼具模块化和结构化,通过设定初始参数生成飞行器运动轨迹。
利用仿真试验数据结果表明,该平台实现了捷联惯组数据在发射坐标下的解算,可用于某飞行器试验数据事后处理和运动轨迹复现工作。
关键词:飞行器;捷联惯组;运动轨迹;MATLAB GUI中图分类号:TP273 文献标识码:A 文章编号:1007-9416(2020)09-0161-040 引言某飞行器惯性测量装置用于测量飞行器相对惯性空间转动的角速度和视加速度在载体坐标系三个轴上的投影[1]。
然而在实际应用中,由于测量所得角速率和加速度为载体坐标下数据,在事后数据处理中需将其转换到发射坐标系下进行分析,以复现出发射坐标下的弹道数据,并需要将数据进行存储再利用。
MATLAB GUI是基于MATLAB的图形化编程语言,基于MATLAB GUI的编程既能有效利用MATLAB的计算功能,又能通过GUI直观显示计算结果及曲线等,能更直观地分析飞行器飞行弹道特征及捷联惯组性能。
本文捷联惯组测量数据分析平台各部分功能相互独立,初始参数、原始测量数据读取部分根据飞行器的运动特征和初始环境,设定初始参数选取原则,并在弹道复现计算中进行了初始误差补偿;导航数据解算部分是分析平台的核心部分,对飞行器的姿态矩阵更新选取了两种算法:欧拉角法和四元数法,实现了两种算法下姿态解算的精度比较;该平台另一特色即是结合飞行器飞行弹道特征点,选取主要特征点进行解算分析。
通过设置不同初始参数,复现不同初始条件下各时刻载体的运动信息、姿态信息及弹道特征点信息,解算结果明确,便于工程技术人员进行事后数据分析。
1 捷联惯组测量数据分析建模1.1 选取坐标系捷联惯性测量装置固连于飞行器载体上,测量到的是载体坐标系轴向上的分量,飞行器运动轨迹的观测点在发射坐标系,初始运动参数的设置建立在艇体坐标系[1]。
捷联惯性导航总结
z
0 x
y
姿态矩阵计算方法
3 四元数微分方程及其解
1 qk 1 exp( )qk Mq k qk rk 2
rk ac as x i as y j as z k
(0.5 )2 (0.5 )4 ac cos( ) 1 2 2! 4!
tk b nb
t
v f dt
b tk
t
速度微分方程求解
2 速度积分算法
v f dt [2ω ωen ] v dt g1 dt
n n 0 0 n ie n n n 0
t
t
t
v
n k 1
v u gt
n k n
b
v
n l 1
[ I 2Ωie tl Ωen tl ]v
v f dt
b
t
1 tk 1 αk 1 α ωdt 2 tk
α ωdt
tk
tk
t
捷联导航算法
捷联惯导误差分析
主要误差源
1 初始对准误差
2 惯性器件误差
3 计算误差
误差传播方程:误差源与系统性能的相互关系的 方程,主要用来阐述各种类型误差随时间的传播 过程。
捷联惯导误差分析
sin( / 2) (0.5 ) (0.5 ) as 0.5(1 ) 3! 5!
2 4
姿态矩阵计算方法
4 圆锥误差
σ αk 1 αk 1
α ωdt
tk t
1 tk 1 αk 1 α ωdt 2 tk
姿态矩阵计算方法
5 姿态角计算
b nbz
捷联惯导系统的算法研究及其仿真实现(捷联惯导系统的发展趋势 初始对准技术的发展与研究现状)
捷联惯导系统的算法研究及其仿真实现Study and Simulation of Strapdown Inertial Navigation System1.1.3捷联惯导系统的发展趋势捷联式惯导系统是从20世纪60年代初开始发展起来的。
20世纪70年代以来,作为捷联系统的核心部件—惯性测量装置和计算机技术有了很大发展,而电子技术、计算机技术、现代控制理论的不断进步,为捷联惯性技术的发展创造了有利条件。
在硬件方面,新一代惯性器件如激光陀螺、光纤陀螺的成功研制,为捷联惯导的飞速发展打下了物质基础。
进入20世纪80-90年代,在航天飞机、宇宙飞船、卫星等民用领域及各种战略、战术导弹、军用飞机、反潜武器、作战舰艇等军事领域开始采用动力调谐式陀螺、激光陀螺和光纤式陀螺的捷联惯导系统。
其中激光陀螺和光纤式陀螺是捷联惯导系统的理想器件。
激光陀螺具有角速率动态范围宽、对加速度和震动不敏感、不需温控、启动时间特别短和可靠性高等优点。
激光陀螺惯导系统己在波音757/767、A310民机以及F-20战斗机上试用,精度达到 1.85km/h 的量级。
20世纪90年代,激光陀螺惯导系统估计占到全部惯导系统的一半以上,其价格与普通惯导系统差不多,但由于增加了平均故障间隔时间,其寿命期费用只有普通惯导系统的15%-20%。
光纤陀螺实际上是激光陀螺中的一种,其原理与环型激光陀螺相同,它克服了由激光陀螺闭锁带来的负效应,具有检测灵敏度和分辨率极高、启动时间极短、动态范围极宽、结构简单、零部件少体积小、造价低、可靠性高等优点。
采用光纤陀螺的捷联航姿系统已用于战斗机的机载武器系统及波音777飞机中。
波音777由于采用了光纤陀螺的捷联惯导系统,其平均故障间隔时间可高达20000h。
采用光纤陀螺的捷联惯导系统被认为是一种极有发展前途的导航系统。
而随着航空航天技术的发展及新型惯性器件关键技术的陆续突破,捷联惯导系统的可靠性、精度将会更高。
捷联惯性导航系统下运载体典型作战动作的航迹仿真
高摇-摇适用于当自身飞行方向和敌机飞行方 向夹角不太大时,攻击机几乎与目标机速度相同, 且没有剩余提前量进行后置滚转的情况。若离目标 很接近,敌机减速,意图甩至我机前方,此时迅速 拉起,避免被咬尾,并且将速度优势转化为高度优 势。
低摇-摇目的是增加接近速度和提前量,以便 切入对方内圈达成前置追踪。大角度弧线下滑,获 得速度优势,然后内切使攻击机能够获得有利的攻 击站位。
1.2 前置转弯
前置转弯的定义实质即为较早的转弯,是指攻 击机在前半球接近对手的情况下,在与对手交会之 前就开始进行转弯。前置转弯是一种非常有效的进 攻机动,转弯越早,潜在的好处越大。
在所有空战机动中,盘旋急转是最基本的动作, 是一个高过载机动。因为长时间的急转,会导致自 身飞行速度迅速降低,所以,为了避免敌方击中, 所以必须立刻进入其它动作。可以看到,上述的两 个动作,都可以看作盘旋急转的一种变形。同时, 拉起和俯冲也是构成基本动作的关键。所以,在后 续的仿真中,主要以盘旋、转弯、拉起作为主要分 析对象。
1 典型作战动作
攻击机作为将武器系统携带到一定高度并进 行发射的平台,其存在意义即为攻击消灭敌机。当 攻击机相对于敌机拥有高度优势时,意味着可攻可
第4期
刘俊妧等:捷联惯性导航系统下运载体典型作战动作的航迹仿真
·269·
退的速度优势;当攻击机直飞或者长时间维持相同 的机动动作,将大幅度增加敌机攻击成功概率[2]。 不同的机动动作对决定机载武器是否可以成功使 用的因素,如攻击机的射击距离、瞄准动作、攻击 机和目标机的相对位置、惯性系统精度等,将造成 不同的影响。
捷联惯性导航系统下运载体典型作战动作的航迹进行仿真。通过模型仿真,评估现役装备使用条
件下各种机动动作对载机惯导精度的影响,研究在不同机动条件下惯导系统精度是否满足实际作
捷联惯性导航系统算法程序汇总
捷联惯性导航系统算法
1.经典捷联惯性导航算法(毕卡逼近、旋转矢量、四阶龙格库塔算法),使用C语言编写,在
实际的系统中得到验证;
2.组合导航算法,包括:速度匹配、位置匹配、姿态角匹配等;
3.捷联惯性导航系统初始对准算法,粗对准方法:经典解析法、基于惯性系抗晃动基座解析
法,精对准方法:基于Kalman滤波的速度匹配、位置匹配精对准方法;
4.捷联惯导系统姿态算法研究,包括:四阶龙格库塔算法、旋转矢量算法,在典型圆锥运动
环境下对姿态解算算法系数进行优化;
5.利用Allan方差分析对光纤陀螺随机误差进行分析,为了抑制随机误差采用Kalman滤波
器对其进行滤波;
6.单轴旋转捷联惯导系统(SINS)多位置初始对准算法以及导航解算方法;
以上所有算法均采用C语言编写,且已经在实际的惯性导航系统中进行了充分的验证,如果需要交流,可以进一步进行联系!。
基于Matlab的轨迹发生器与惯导解算
由姿态转换矩阵便可计算航向角 ψ、 俯仰角 θ 和横滚角 φ。 写成矩阵形式,即:
o 0 q b q 1 = 1 ωnbx b q 2 2 ωnby b 3 ωnbz q
— 31 —
2016年信息与电脑9下-正文.indd 31
2016/12/19 8:59:06
计算机工程应用技术
更新方法为四元数姿态更新法。 2.1.1 四元数姿态更新法
信息与电脑 China Computer&Communication
2016 年第 18 期
2.1.3 位置更新方程 捷联惯导系统的位置更新是通过速度和角速度增量实现 的,位置更新微分方程为:
1 引言
水下航行器在海洋开发、科学研究和军事应用中都占有 重要地位,这些工作能否顺利完成、完成效果的好坏,很大 程度上依赖于水下航行器的导航与定位能力。 惯性是所有质量体的基本属性,因此,基于惯性原理工 作的惯导系统具有极好的自主性和隐蔽性,特别适合水下的 特殊环境。然而惯导系统的误差是随时间累积的,长时间的 水下导航误差必然是逐渐发散的。因此,无论是解算算法、 组合导航方式,最终目的都是抑制惯导系统的误差发散,提 高惯导系统的导航精度。然而这些研究若在真实的惯导系统 上进行,不仅耗费大量的人力、物力,研究结果也不易观察、 体现、总结。因此,众多研究者们将目光投向计算机仿真。 轨迹发生器的设计与使用应运而生,并在惯导系统的研究中 发挥非常重要的作用。轨迹发生器的设计可以模拟水下航行 器在水下运动的各类参数信息,如姿态、速度、加速度等, 极大提高了研究效率,又可在轨迹发生器的基础上进行惯导 解算的研究。 本文提出一种水下航行器轨迹生成的方法,通过这种方 法模拟航行器在水下运动的各个姿态过程,可以作为研究惯 导系统工作原理,乃至于组合导航系统的数据依据。并在此 基础上进行了惯导解算,分析了解算结果。
捷联惯导中微分方程的变步长解法
捷联惯导中微分方程的变步长解法
李斌;董慧颖;郝永平;王磊
【期刊名称】《科技创新导报》
【年(卷),期】2008(000)031
【摘要】Dormand-Prince方法是许多仿真软件(如matlab,Octave)默认的解微分方程的方法,该方法能在保证较高的精度前提下尽量减少运算时间.文中介绍了Dormand-Prince方法及一种变步长策略,给出了用C语言实现该方法的流程图,蕾方法用于捷联惯导微分方程的求解,取得了较高精度.
【总页数】1页(P236)
【作者】李斌;董慧颖;郝永平;王磊
【作者单位】沈阳理工大学,辽宁沈阳,100168;沈阳理工大学,辽宁沈阳,100168;沈阳理工大学,辽宁沈阳,100168;沈阳理工大学,辽宁沈阳,100168
【正文语种】中文
【中图分类】O182
【相关文献】
1.抛物型偏微分方程的变步长显式差分解法 [J], 秦学志;年四洪
2.一种新的捷联矩阵更新算法在无陀螺捷联惯导系统中的应用 [J], 赵龙;史震;马澍田
3.强跟踪扩展卡尔曼滤波及其在捷联惯导初始对准中的应用 [J], 吴苗;郭士荦;许江宁
4.EMD滤波在捷联惯导晃动基座初始对准中的应用 [J], 于晓雪; 张志鑫; 夏金桥
5.广义回归神经网络在冗余捷联惯导故障诊断中的应用研究 [J], 乔鹏超;孙湘钰;罗广地
因版权原因,仅展示原文概要,查看原文内容请购买。
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
捷联惯导程序,依据加表和陀螺仪的输出数据来求解飞行器的姿态clc;clear;format long; %设置数据精度为15位小数Data=importdata('temp.txt'); % 导入实验所采集的数据,以矩阵形式赋给Data变量,temp.txt 必须与该M文件在同一个文件夹中Px=Data(:,3); % Px,Py,Pz为陀螺仪的输出值Py=Data(:,4);Pz=Data(:,5);Nx=Data(:,6); % Nx,Ny,Nz为加速度计的输出值Ny=Data(:,7);Nz=Data(:,8);% 陀螺仪模型参数标定如下:Sx = -4.085903e-006 ; Sy = -4.085647e-006 ; Sz = -4.085170e-006 ;Mxy = 5.059527e-003 ; Mxz = -1.031103e-003 ; Myx = -3.355451e-003 ;Myz = 3.508468e-003 ; Mzx = -1.266671e-003 ; Mzy = -2.318244e-004 ;Dx = -2.009710e-006 ; Dy = 8.156346e-007 ; Dz = -5.749059e-007 ;GyroCali_A = [ 1 -Mxy -Mxz ; -Myx 1 -Myz ; -Mzx -Mzy 1 ];% 加速度计模型参数标定如下:Kx = 9.272930e-004 ; Ky = 9.065544e-004 ; Kz = 9.443748e-004 ;Ixy = 6.533872e-003 ; Ixz = 9.565992e-004 ; Iyx = -6.319376e-003 ;Iyz = -6.902339e-004 ; Izx = -1.144549e-003 ; Izy = -3.857963e-004 ;Bx = -3.400847e-002 ; By = -8.916341e-003 ; Bz = -9.947414e-003 ;AccCali_A = [1 -Ixy -Ixz ; -Iyx 1 -Iyz ; -Izx -Izy 1 ];Delta_t = 0.05; %采样时间为0.05秒Delta_Theta_x = 0;Delta_Theta_y = 0;Delta_Theta_z = 0; %定义陀螺仪输出的角度增量Delta_Vx = 0;Delta_Vy = 0;Delta_Vz = 0; %定义加速度计输出的速度增量L = zeros(1,12001);L(1)= 45.7328*pi/180 ; %纬度用L表示,纬度的初始值划为弧度形式,因为后面计算位置矩阵更新L(2)= 45.7328*pi/180 ; %时需要用到前两次的L值来计算当前L值,所以在此定义2个初始L值Lamda = 126.6287*pi/180 ; %经度用Lamda表示,经度的初始值划为弧度形式h = 136 ; %高度用h表示V = [ 0 ; 0 ; 0 ]; %导航坐标系中的东北天初始速度都为0Vx = 0; %方便后面的速度计算与速度更新Vy = 0;Vz = 0;Theta = 0;Gama = 0;Fai = 0; %初始姿态角(俯仰角/倾斜角/航向角)都为0,此处均为弧度Re = 6378254 ; Rp = 6356803 ;%定义地球的半长轴与半短轴e = (Re - Rp)/Re ; %定义旋转椭球扁率(椭球度)Wie = 15.04107/180*pi ; %定义地球自转角速度,地球坐标系相对于惯性坐标系的角速度Theta_Matrix = zeros(1,12000); %定义姿态角矩阵,供画图用Gama_Matrix = zeros(1,12000);Fai_Matrix = zeros(1,12000);L_Matrix = zeros(1,12001); %定义经纬度矩阵,供画图用,L的特殊性决定了其数据个数为12001L_Matrix(1) = 45.7328;Lamda_Matrix = zeros(1,12000);Ve_Matrix = zeros(1,12000); %定义速度矩阵,供画图用Vn_Matrix = zeros(1,12000);Vu_Matrix = zeros(1,12000);%以下计算捷联矩阵的初始值,捷联矩阵的初始值仅仅由Theta,Gama,Fai的初始值决定T = [ cos(Gama)*cos(Fai)-sin(Gama)*sin(Theta)*sin(Fai) -cos(Theta)*sin(Fai) sin(Gama)*cos(Fai)+cos(Gama)*sin(Theta)*sin(Fai) ;cos(Gama)*sin(Fai)+sin(Gama)*sin(Theta)*cos(Fai) cos(Theta)*cos(Fai) sin(Gama)*sin(Fai)-cos(Gama)*sin(Theta)*cos(Fai) ;-sin(Gama)*cos(Theta) sin(Theta) cos(Gama)*cos(Theta) ];%由捷联矩阵的初始值计算初始四元数值,为捷联矩阵的实时更新做准备if(T(3,2)-T(2,3)>0)Q1 = 0.5*sqrt(1+T(1,1)-T(2,2)-T(3,3));else if (T(3,2)-T(2,3)==0)Q1 = 0;else Q1 = -0.5*sqrt(1+T(1,1)-T(2,2)-T(3,3));end %求解Q1endif(T(1,3)-T(3,1)>0)Q2 = 0.5*sqrt(1-T(1,1)+T(2,2)-T(3,3));else if (T(1,3)-T(3,1)==0)Q2 = 0;else Q2 = -0.5*sqrt(1-T(1,1)+T(2,2)-T(3,3));end %求解Q2endif(T(2,1)-T(1,2)>0)Q3 = 0.5*sqrt(1-T(1,1)-T(2,2)+T(3,3));else if (T(2,1)-T(1,2)==0)Q3 = 0;else Q3 = -0.5*sqrt(1-T(1,1)-T(2,2)+T(3,3));end %求解Q3endQ0 = 0.5*sqrt(1-Q1*Q1-Q2*Q2-Q3*Q3); %求解Q0Q = [Q0 ; Q1 ; Q2 ; Q3]; %四元数初始值Q = Q / sqrt(Q0*Q0+Q1*Q1+Q2*Q2+Q3*Q3); %四元数的初始归一化,为得到最小漂移误差%以下求位置矩阵的初始值,通过位置矩阵更新后,反过来算运载体所在的经纬度%位置矩阵仅仅与经纬度有关系,Ce2n表示把地球坐标系转换为导航坐标系的转换矩阵Ce2n = [ -sin(Lamda) cos(Lamda) 0 ;-sin( L(1) )*cos(Lamda) -sin(L(1))*sin(Lamda) cos( L(1) );cos( L(1) )*cos(Lamda) cos( L(1) )*sin(Lamda) sin( L(1) ) ];%大循环,共执行12000次,实时更新捷联矩阵,速度矩阵,位置矩阵,保存作图所需数据for k = 1:12000;GyroCali_B = [Sx*Px(k)-Dx*Delta_t ; Sy*Py(k)-Dy*Delta_t ; Sz*Pz(k)-Dz*Delta_t ];Delta_Theta = GyroCali_A * GyroCali_B ; %计算陀螺仪输出的角度增量Delta_Theta_x = Delta_Theta(1);Delta_Theta_y = Delta_Theta(2);Delta_Theta_z = Delta_Theta(3);Delta_Theta_Module = sqrt( Delta_Theta_x * Delta_Theta_x + Delta_Theta_y * Delta_Theta_y + Delta_Theta_z * Delta_Theta_z );AccCali_B = [Kx*Nx(k)-Bx*Delta_t ; Ky*Ny(k)-By*Delta_t ; Kz*Nz(k)-Bz*Delta_t ];Delta_V = AccCali_A * AccCali_B ; %计算加速度计输出的速度增量Delta_Vx = Delta_V(1);Delta_Vy = Delta_V(2);Delta_Vz = Delta_V(3);Delta_V_Module = sqrt( Delta_Vx * Delta_Vx + Delta_Vy * Delta_Vy + Delta_Vz * Delta_Vz );%使用毕卡法求解四元数更新矩阵,即捷联矩阵Bika = zeros(4);Bika(1,1) = cos(0.5 * Delta_Theta_Module);Bika(1,2) = -Delta_Theta_x / Delta_Theta_Module * sin(0.5 * Delta_Theta_Module);Bika(1,3) = -Delta_Theta_y / Delta_Theta_Module * sin(0.5 * Delta_Theta_Module);Bika(1,4) = -Delta_Theta_z / Delta_Theta_Module * sin(0.5 * Delta_Theta_Module);Bika(2,1) = -Bika(1,2);Bika(2,2) = Bika(1,1);Bika(2,3) = -Bika(1,4);Bika(2,4) = Bika(1,3);Bika(3,1) = -Bika(1,3);Bika(3,2) = -Bika(2,3);Bika(3,3) = Bika(1,1);Bika(3,4) = -Bika(1,2);Bika(4,1) = -Bika(1,4);Bika(4,2) = -Bika(2,4);Bika(4,3) = -Bika(3,4);Bika(4,4) = Bika(1,1);Q = Bika * Q; % 每循环一次,更新一次四元素Q值,为求捷联矩阵Q = Q / sqrt(Q0*Q0+Q1*Q1+Q2*Q2+Q3*Q3); %四元数的归一化,为得到最小漂移误差Q0 = Q(1);Q1 = Q(2);Q2 = Q(3);Q3 = Q(4);%捷联矩阵的四元数表达式T = [ Q0*Q0+Q1*Q1-Q2*Q2-Q3*Q3 2*(Q1*Q2-Q0*Q3) 2*(Q1*Q3+Q0*Q2)2*(Q1*Q2+Q0*Q3) Q0*Q0-Q1*Q1+Q2*Q2-Q3*Q3 2*(Q2*Q3-Q0*Q1)2*(Q1*Q3-Q0*Q2) 2*(Q2*Q3+Q0*Q1) Q0*Q0-Q1*Q1-Q2*Q2+Q3*Q3 ];%*********************************************************************%********************求三个姿态角Theta,Gama和Fai ******************** %*********************************************************************Theta_Main = asin( T(3,2) );Gama_Main = atan( -T(3,1) / T(3,3));Fai_Main = atan( -T(1,2) / T(2,2));Theta = Theta_Main;if (T(3,3)>0)Gama = Gama_Main ;else if (T(3,3)<0 && Gama_Main > 0)Gama = Gama_Main + pi;else Gama = Gama_Main - pi; %此处用else实为不妥,不过为了程序的完善性,只能这样了endendif ( T(2,2)<0 )Fai = Fai_Main + pi ;else if (T(2,2)==0)Fai = pi/2;else if ( Fai_Main>0)Fai = Fai_Main;else Fai = Fai_Main + 2*pi ;endendend%以下存储姿态角到三个矩阵里面,为画图做准备Theta_Matrix(k) = Theta*180/pi; %作图用矩阵,以角度表示Gama_Matrix(k) = Gama*180/pi; %作图用矩阵,以角度表示if (Fai<2*pi)Fai_Matrix(k) = Fai*180/pi;else Fai_Matrix(k) = Fai*180/pi-360; %作图用矩阵,以角度表示end%到此为止,姿态角的求解完毕,以下先求速度%*********************************************************************%********************求飞行器相对于东北天的速度************************* %*********************************************************************Rm = Re*( 1-2*e+3*e*sin( L(k+1) )^2 );Rn = Re*( 1+e*sin(L(k+1))^2 );LL = 3/2*L(k+1) - 1/2*L(k) ;F = [ 0 -1/(Rm + h) 0 ; 1/(Rn + h) 0 0 ; tan( LL )/(Rn + h) 0 0];g = 9.7803+0.051799*sin(L(k+1))^2-0.94114e-006*h;G = [0;0;-g];Wen2n = F*V;Wie2n = [0 ; Wie*cos(L(k+1)); Wie*sin(L(k+1))];W = 2*Wie2n + Wen2n;W_X = [ 0 -W(3) W(2) ; W(3) 0 -W(1) ; -W(2) W(1) 0 ]; %此式中W_X为2*Wie+Wen2n 的反对称矩阵V = V + T*( Delta_V+ [ 0 -Delta_Theta_z Delta_Theta_y ; Delta_Theta_z 0 -Delta_Theta_x ; -Delta_Theta_y Delta_Theta_x 0 ]*Delta_V ) + Delta_t*(G-W_X*V);Vx = V(1);Vy = V(2);Vz = V(3);V e_Matrix(k) = V(1);Vn_Matrix(k) = V(2);Vu_Matrix(k) = V(3);%*********************************************************************%********************求飞行器所在的经纬度****************************** %*********************************************************************Epsilon = F*V*Delta_t;Ce2n =( eye(3) - [ 0 -Epsilon(3) Epsilon(2) ; Epsilon(3) 0 -Epsilon(1) ; -Epsilon(2) Epsilon(1) 0 ] )*Ce2n ; %位置矩阵实时更新%下面通过位置矩阵来实时更新经纬度L(k+2) = asin( Ce2n(3,3)); %由于L本来就是以矩阵形式定义的,下面定义一个把L用角度表示的矩阵L_Matrix(k+1) = L(k+2)*180/pi;Lamda_Main = atan( Ce2n(3,2)/Ce2n(3,1) ); %计算出来的L和Lamda都是弧度制的if (Ce2n(3,1)>0)Lamda = Lamda_Main;else if (Lamda_Main<0)Lamda = Lamda_Main + pi;else Lamda = Lamda_Main - pi;endendLamda_Matrix(k) = Lamda*180/pi; %作图用矩阵,以角度表示end%**************************************************************************%********************以下是画图程序*****************************************%**************************************************************************k=1:1:12000; %绘制三轴姿态变化图线-绿色figure(1);plot(k/20,Theta_Matrix(k),'g');xlabel('Time(second)');ylabel('Angle(degree)');title('Theta(俯仰角)');grid on;figure(2);plot(k/20,Gama_Matrix(k),'g');xlabel('Time(second)');ylabel('Angle(degree)');title('Gama(滚转角)');grid on;figure(3);plot(k/20,Fai_Matrix(k),'m');xlabel('Time(second)');ylabel('Angle(degree)');title('Fai(偏航角)');grid on;%绘制东北天个方向的速度变化曲线-红色figure(4);plot(k/20,Ve_Matrix(k),'r');xlabel('Time(second)');ylabel('Velocity(m/s)');title('Ve(东向速度)');grid on;figure(5);plot(k/20,Vn_Matrix(k),'r');xlabel('Time(second)');ylabel('Velocity(m/s)');title('Vn(北向速度)');grid on;figure(6);plot(k/20,Vu_Matrix(k),'r');xlabel('Time(second)');ylabel('Velocity(m/s)');title('Vu(天向速度)');grid on;%绘制飞行器所在经纬度曲线-蓝色figure(7);plot(k/20,L_Matrix(k),'b');xlabel('Time(second)');ylabel('Degree');title('Latitude L(纬度)');grid on;figure(8);plot(k/20,Lamda_Matrix(k),'b');xlabel('Time(second)');ylabel('Degree');title('Longitude Lamda(经度)');grid on;。