CN106908853B - 基于相关分析与经验模分解的捷联式重力仪误差矫正方法 - Google Patents
基于相关分析与经验模分解的捷联式重力仪误差矫正方法 Download PDFInfo
- Publication number
- CN106908853B CN106908853B CN201710152576.3A CN201710152576A CN106908853B CN 106908853 B CN106908853 B CN 106908853B CN 201710152576 A CN201710152576 A CN 201710152576A CN 106908853 B CN106908853 B CN 106908853B
- Authority
- CN
- China
- Prior art keywords
- result
- modal components
- strapdown
- empirical mode
- mode decomposition
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Active
Links
Classifications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V7/00—Measuring gravitational fields or waves; Gravimetric prospecting or detecting
- G01V7/02—Details
- G01V7/06—Analysis or interpretation of gravimetric records
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V7/00—Measuring gravitational fields or waves; Gravimetric prospecting or detecting
- G01V7/16—Measuring gravitational fields or waves; Gravimetric prospecting or detecting specially adapted for use on moving platforms, e.g. ship, aircraft
Abstract
本发明属于重力测量领域,公开了基于相关分析与经验模分解的捷联式重力仪误差矫正方法。通过等效带通滤波器研究相应频段内的重力仪误差以提高方法的精度,利用经验模分解得到测量结果的本征模态分量,利用本征模态分量与伪比力的相关性来进行误差分离,并且利用迭代的方式避免经验模分解中的模态混叠,本发明在不影响捷联式航空重力测量分辨率的情况下提高了捷联式航空重力测量的精度,对于捷联式重力仪的精度提升具有较大意义。
Description
技术领域:
本发明涉及一种捷联式航空重力仪误差矫正方法,特别涉及了基于相关分析与经验模分解的捷联式重力仪误差矫正方法。
背景技术:
地球重力场的信息在地球物理、大地测量以及大地动力学等领域具有非常重要的作用。在所有测量地球重力场的手段当中,航空重力测量是一种获取高精度大范围地球重力场的高效率手段。航空重力仪分为平台式与捷联式重力仪,捷联式重力仪同平台式重力仪最大的不同在于捷联式重力仪使用数学方式而不是实际的稳定平台来获取重力场信息。捷联式重力仪具有结构简单,体积较小,成本低以及功耗低的优势,并且前期研究工作表明捷联式航空重力仪在获取高精度高分辨率重力场方面具有同平台式重力仪相当甚至更高的表现,所以捷联式重力仪已经成为航空重力测量的趋势所在。捷联式航空重力仪面对的一个问题是,它容易受到气流以及飞行器运动特性的影响。这是由于捷联式重力仪直接固定在舱室内部没有物理平台用于隔离气流以及载体运动对捷联式重力仪的影响,这类由于气流以及载体动态特性导致的误差被称为动态性误差。
经验模分解是一种自适应的信号分解手段,通过经验模分解可以得到信号的本征模态分量。经验模分解可以用于海洋与地震数据分析、机械振动信号分析、声音处理、纹理分析、SAR影像滤波以及气象等领域,并取得了较好的结果。从目前情况来看,直接针对捷联式重力仪动态性误差的研究工作以及利用经验模分解对重力测量信号进行处理的文献较少。文献(Zhao L.et al.Airborne Gravity Data Denoising Based on Empirical ModeDecomposition:A Case Study for SGA-WZ Greenland Test Data.ISPRS Int.J.Geo-Inf.2015,4,2205-2218)应用了经验模分解对实际测量得到的捷联式重力仪数据进行了处理,得到了同航空重力测量中通常使用的低通滤波器相接近的效果,但是对测量结果精度的提高并没有显著帮助。从重力测量理论来看,动态性误差的源头不能通过卡尔曼估计的方式得到,但是动态性误差与比力是强相关的。因此,实现一种基于相关性分析与经验模分解的捷联式航空重力仪动态性误差矫正方法具有很高的实际应用价值。
发明内容:
本发明针对捷联式航空重力仪受外界气流以及载体动态特性影响,进而产生相关的动态性误差的问题,提出了一种基于相关分析与经验模分解的捷联式重力仪误差矫正方法,进而减小动态性对航空重力测量的影响,提升航空重力测量的精度。
主要的技术方案如下:
基于相关分析与经验模分解的捷联式重力仪误差矫正方法,包括以下步骤:
步骤一,对实际测量得到的重力数据进行重力解算并且经过低通滤波得到特定分辨率下的测量结果:
捷联式重力测量的模型为
其中,和为由GPS测量得到的相对地球的加速度与速度,fb是捷联式惯导系统测量得到的在载体系下的比力信息,是由载体系到导航系的姿态变化矩阵,是地球相对于惯性系的转动角速度,是导航系相对于地球的转动角速度,γn是在导航系下的正常重力值,dgn是重力异常结果在导航系下的表示;
步骤二,利用等效带通滤波器对测量结果进行滤波并得到对应频段测量结果,其中等效带通滤波器的结果相当于两个低通滤波器的差值;
步骤三,利用经验模分解对重力测量的带通滤波结果进行分解,得到其本征模态分量;
步骤四,利用GPS信息计算伪比力,之后同样利用步骤二中的等效带通滤波器对其进行带通滤波;
伪比力的计算公式为
其中,IF是所求的不包含惯性测量误差的伪比力;
步骤五,计算步骤三中得到的带通滤波结果的本征模态分量与伪比力带通滤波结果的相关性系数,其中相关性系数最大的本征模态分量为由动态性误差构成的本征模态分量,将其从测量结果中减去以得到中间结果;
步骤六,由于经验模分解中可能存在模态混淆问题,对中间结果再次进行带通滤波,并对结果进行经验模分解,得到中间结果经带通滤波后的本征模态分量;
步骤七,计算新的本征模态分量同伪比力带通滤波结果的相关性系数,依旧认为相关性系数最大的本征模态分量为由动态性误差构成的本征模态分量应从测量结果中减去,得到最终结果;
在本发明中,通过以上七个步骤,便可以实现基于相关分析与经验模分解的捷联式重力仪误差矫正方法,在不改变分辨率的情况下分离了动态性误差。
本发明具有以下优点:
1)直接针对捷联式重力仪中的动态性误差进行矫正,实现误差分离;
2)不改变分辨率的情况下提高测量精度。
附图说明:
图1为某次重力测量实验中的高度变化情况;
图2为飞行实验中测线1、5、6的结果图;
图3为带通滤波器过后测线1、5、6的结果图;
图4为经验模分解的流程图;
图5为测线1的各个本征模态分量;
图6为矫正后飞行实验中测线1、5、6的结果图;
图7为整个发明的流程图。
具体实施方式:
下面结合附图和某次实际航空重力测量实验对本发明中的方法做进一步详细阐述,实验的组织实施情况如表1,在该次飞行实验中包含6条重复测线(以下统称飞行实验),其中第5、6条飞行测线为起伏飞行以验证动态性误差,飞行实验的高度变化如图1所示。
表1
步骤一,对实际测量得到的重力数据进行重力解算并且经过低通滤波得到3km分辨率下的测量结果;
捷联式重力测量的模型为
其中,和为由GPS测量得到的相对地球的加速度与速度,fb是捷联式惯导系统测量得到的在载体系下的比力信息,是由载体系到导航系的姿态变化矩阵,是地球相对于惯性系的转动角速度,是导航系相对于地球的转动角速度,γn是在导航系下的正常重力值,dgn是重力异常结果在导航系下的表示;
对飞行实验的数据进行重力解算,得到6条重复测线的测量精度为1.73mGal/3km(1mGal=1*10-5m/s2),其中测线1、5、6的精度对比如图2。
步骤二,利用等效带通滤波器对测量结果进行滤波并得到对应频段测量结果,其中等效带通滤波器的结果相当于两个低通滤波器的差值;
将数据经过带通滤波器,图3给出了经过带通滤波器后测线1、5、6的结果。
步骤三,利用经验模分解对重力测量的带通滤波结果进行分解,得到其本征模态分量:
经验模态分解的步骤如下所示(图4为其流程图):
(1):x(t)为将要进行经验模分解的信号,确定x(t)的所有极大值和极小值,分别对极大值点和极小值点进行插值(例如采用三次样条进行插值计算),构造的上下包络线xup(t)和xlow(t),计算上下包络线的均值m1(t)=(xup(t)+xlow(t))/2;
(2):把原信号序列x(t)减去该均值m1(t)可得到一个去掉低频的新数据序列g1(t),即g1(t)=x(t)-m1(t);
(3):判断g1(t)是否是一个本征模态分量(intrinsic mode function,IMF):
如果g1(t)符合IMF的定义条件,是一个IMF,则提取g1(t)作为第一个IMF分量,令imf1(t)=g1(t),并求原信号与imf1之间的差值r1(t),r1(t)=x(t)-imf1(t)。
如果g1(t)不是一个IMF,则将g1(t)视为一个新的信号序列,重复第一步和第二步,求其包络均值及g1(t)与m11(t)间的差值g1(t);对g1(t)重复上述过程n次,直到符合IMF的定义条件,则令为x(t)的第一个IMF分量,并求原信号与imf1之间的差值r1(t),即r1(t)=x(t)-imf1(t);
(4):将r1(t)作为一个新的“原始”信号,重复步骤1至3,提取第二个本征模态函数分量imf2,令r2(t)=r1(t)-imf2(t);将r2(t)作为一个新的“原始”信号,抽取第三个本征模态函数分量imf3;以此类推,直到第K次的余项rK=r(K-1)-imfk满足终止条件,则停止迭代,x(t)的经验模态分解完成。
经以上步骤,此时原始信号x(t)可以被表示为以上分解得到的所有本征模态函数和一个残余趋势项之和:
式中imfk(t)表示第k个本征模态函数,每个本征模态函数imf1,imf2,……,imfk表示原始信号从小尺度到大尺度的时间特征。
对各条测线进行经验模分解,得到其本征模态分量,图5给出了测线1的各个本征模态分量。
步骤四,利用式(2)计算伪比力,之后利用步骤二中相同的带通滤波器对其进行带通滤波;
步骤五,计算步骤三中得到的本征模态分量与伪比力带通滤波结果的相关性系数,其中相关性系数最大的本征模态分量为由动态性误差构成的本征模态分量,将其从测量结果中减去以得到中间结果。
测线1中本征模态分量与伪比力带通滤波结果的相关性系数如表2所示,其中第2个本征模态分量的相关性系数最大,所以认为第2个本征模态分量为由动态性误差构成的本征模态分量,在测量结果中去除,得到中间结果;
表2
本征模态分量编号 | 1 | 2 | 3 | 4 | 5 | 6 |
相关性系数 | -0.2 | -0.56 | -0.13 | -0.09 | 0.04 | 0.01 |
步骤六,由于经验模分解中可能存在模态混淆问题,对中间结果再次利用步骤二中带通滤波器进行带通滤波,并对结果进行经验模分解,得到新的本征模态分量。
步骤七,计算新的本征模态分量与伪比力带通滤波结果的相关性系数,依旧认为相关性系数最大的本征模态分量为由动态性误差构成的本征模态分量应从测量结果中减去,得到最终结果,图7为整个算法流程图。
经过矫正后,飞行实验的精度从1.73mGal/3km提升到1.27mGal/3km,图6给出了提升后测线1、5、6的精度。以上结果充分说明了本发明中算法的突出效果。
以上仅是本发明的优选实施方式,本发明的保护范围并不仅局限于上述实施例,凡属于本发明思路下的技术方案均属于本发明的保护范围。应当指出,对于本技术领域的普通技术人员来说,在不脱离本发明原理前提下的若干改进和润饰,应视为本发明的保护范围。
Claims (1)
1.基于相关分析与经验模分解的捷联式重力仪误差矫正方法,其特征在于包含以下步骤:
步骤一,对实际测量得到的重力数据进行重力解算并且经过低通滤波得到特定分辨率下的测量结果:
捷联式重力测量的模型为
其中,和为由GPS测量得到的相对地球的加速度与速度,fb是由捷联式惯导系统测量得到的在载体系下的比力信息,是由载体系到导航系的姿态变化矩阵,是地球相对于惯性系的转动角速度,是导航系相对于地球的转动角速度,γn是在导航系下的正常重力值,dgn是重力异常结果在导航系下的表示;
步骤二,利用等效带通滤波器对测量结果进行滤波并得到对应频段测量结果,其中等效带通滤波器的结果相当于两个低通滤波器的差值;
步骤三,利用经验模分解对重力测量的带通滤波结果进行分解,得到其本征模态分量;
步骤四,利用GPS信息计算伪比力,之后同样利用步骤二中的等效带通滤波器对其进行带通滤波;
伪比力的计算公式为
其中,IF是所求的不包含惯性测量误差的伪比力;
步骤五,计算步骤三中得到的带通滤波结果的本征模态分量同伪比力带通滤波结果的相关性系数,其中相关性系数最大的本征模态分量为由动态性误差构成的本征模态分量,将其从测量结果中减去以得到中间结果;
步骤六,由于经验模分解中可能存在模态混淆问题,对中间结果再次进行带通滤波,并对结果进行经验模分解,得到中间结果带通滤波后的本征模态分量;
步骤七,计算新的本征模态分量同伪比力带通滤波结果的相关性系数,依旧认为相关性系数最大的本征模态分量为由动态性误差构成的本征模态分量应从测量结果中减去,得到最终结果。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201710152576.3A CN106908853B (zh) | 2017-03-15 | 2017-03-15 | 基于相关分析与经验模分解的捷联式重力仪误差矫正方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201710152576.3A CN106908853B (zh) | 2017-03-15 | 2017-03-15 | 基于相关分析与经验模分解的捷联式重力仪误差矫正方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN106908853A CN106908853A (zh) | 2017-06-30 |
CN106908853B true CN106908853B (zh) | 2018-09-18 |
Family
ID=59186451
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201710152576.3A Active CN106908853B (zh) | 2017-03-15 | 2017-03-15 | 基于相关分析与经验模分解的捷联式重力仪误差矫正方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN106908853B (zh) |
Families Citing this family (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN109001829B (zh) * | 2018-07-12 | 2019-11-05 | 中国人民解放军国防科技大学 | 一种捷联式水下动态重力测量仪 |
CN110632674B (zh) * | 2019-10-28 | 2021-04-23 | 中国自然资源航空物探遥感中心 | 一种航空重力测量数据的弱信息提取方法 |
Family Cites Families (4)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN103917850B (zh) * | 2011-10-25 | 2016-01-20 | 中国人民解放军国防科学技术大学 | 一种惯性导航系统的运动对准方法 |
CN102538781B (zh) * | 2011-12-14 | 2014-12-17 | 浙江大学 | 基于机器视觉和惯导融合的移动机器人运动姿态估计方法 |
CN103529482B (zh) * | 2013-10-25 | 2016-05-11 | 中国人民解放军国防科学技术大学 | 一种高精度确定载体动态加速度的方法 |
CN106052686B (zh) * | 2016-07-10 | 2019-07-26 | 北京工业大学 | 基于dsptms320f28335的全自主捷联惯性导航系统 |
-
2017
- 2017-03-15 CN CN201710152576.3A patent/CN106908853B/zh active Active
Also Published As
Publication number | Publication date |
---|---|
CN106908853A (zh) | 2017-06-30 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN110632674B (zh) | 一种航空重力测量数据的弱信息提取方法 | |
US20140081595A1 (en) | Gravity gradiometer survey techniques | |
CN109814163B (zh) | 一种基于扩展补偿模型的航磁张量数据抑噪方法及系统 | |
Cai et al. | The first airborne scalar gravimetry system based on SINS/DGPS in China | |
CN106908853B (zh) | 基于相关分析与经验模分解的捷联式重力仪误差矫正方法 | |
US20220075091A1 (en) | Post-compensation method for motion errors of rotating accelerometer gravity gradiometer | |
CN109443342A (zh) | 新型自适应卡尔曼无人机姿态解算方法 | |
CN105628189B (zh) | 一种直升机机体振动信号的自适应处理方法 | |
CN110907208A (zh) | 一种运载火箭的模态试验方法 | |
CN112257186B (zh) | 针对小型四旋翼飞行器气动参数的时域辨识方法 | |
WO2012001388A2 (en) | Gravity survey data processing | |
CN108983795A (zh) | 一种三轴姿态校正方法及设备 | |
CN109282786A (zh) | 用于固定翼无人机的测量高度方法及装置 | |
Wu et al. | A Kalman filter approach based on random drift data of fiber optic gyro | |
RU2382988C1 (ru) | Бесплатформенная инерциальная система ориентации на "грубых" чувствительных элементах | |
CN107290794A (zh) | 一种时间域航空电磁探测系统接收线圈运动噪声的数值仿真方法 | |
CN106596014A (zh) | 一种直升机舱内飞行振动环境模拟试验方法 | |
CN112051595A (zh) | 利用dgps位置信息求解载体运动加速度的后向差分滤波方法 | |
CN103217715A (zh) | 多尺度规则网格层析反演静校正方法 | |
CN106323280A (zh) | 用于bds和sins导航定位系统的滤波器和滤波方法 | |
CN104748747B (zh) | 一种适用于装甲车的定位定向装置及方法 | |
Burdette et al. | Using an array of in-canopy sensors for guided airdrop system state estimation | |
CN101706283B (zh) | 高度通道的跨音速激波干扰自适应抑制方法 | |
Wu et al. | Real-time correction for sensor's dynamic error based on DSP | |
CN206479651U (zh) | 一种具有姿态记录功能的地空tem接收系统 |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
PB01 | Publication | ||
PB01 | Publication | ||
SE01 | Entry into force of request for substantive examination | ||
SE01 | Entry into force of request for substantive examination | ||
GR01 | Patent grant | ||
GR01 | Patent grant |