CN102620729B - 一种机抖激光陀螺惯性测量单元数字滤波器设计方法 - Google Patents
一种机抖激光陀螺惯性测量单元数字滤波器设计方法 Download PDFInfo
- Publication number
- CN102620729B CN102620729B CN201210117027.XA CN201210117027A CN102620729B CN 102620729 B CN102620729 B CN 102620729B CN 201210117027 A CN201210117027 A CN 201210117027A CN 102620729 B CN102620729 B CN 102620729B
- Authority
- CN
- China
- Prior art keywords
- digital filter
- linear phase
- low
- centerdot
- laser gyroscope
- 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.)
- Expired - Fee Related
Links
Landscapes
- Gyroscopes (AREA)
Abstract
一种机抖激光陀螺惯性测量单元数字滤波器设计方法,其特征在于滤波器包括带阻部分和低通部分;首先根据机抖激光陀螺机抖噪声特点设计线性相位IIR数字滤波器带阻部分;然后在时域内采用带遗忘因子的最小二乘参数估计使线性相位IIR数字滤波器低通部分逼近线性相位FIR低通滤波器,获得近似线性相位;最后建立频域内线性相位IIR数字滤波器低通部分目标函数,利用DFP法优化目标函数,修正时域设计结果,使线性相位IIR数字滤波器低通部分具有所期望的频率特性。其优点是相位近似线性,运算量低,群延时小;降低了滤波器对IMU硬件的要求,提高了IMU测量实时性,保证了时延补偿精度,为IMU高精度实时测量奠定了基础。
Description
技术领域
本发明属于惯性技术领域,涉及一种机抖激光陀螺惯性测量单元(IMU)数字滤波器设计方法,可以应用于位置姿态测量系统(Position and OrientationSystem,POS)和捷联惯性导航系统(Strap down Navigation System,SINS),以及SINS/GPS(Global Position System,全球定位系统)组合导航系统的数据预处理。
背景技术
惯性测量单元(Inertial Measurement Unit,IMU)是POS、SINS、SINS/GPS等惯性测量系统的核心部件,惯性测量系统是以牛顿力学定律为基础,利用IMU测量载体的线运动和角运动参数,在给定初始条件下,由计算机推算得到载体的位置速度姿态参数。它是一种自主式的位置速度姿态测量单元,完全依靠机载设备自主地完成测量任务,和外界不发生任何光、电联系,具有隐蔽性好,工作不受气象条件限制的优点,因此在航空、航天和航海领域得到广泛的使用。IMU主要包含陀螺仪模块、加速度计模块、数据采集模块和数据处理模块,陀螺仪的精度很大程度上决定了IMU的测量精度。由于激光陀螺具有高动态、低成本、高可靠、性能价格比高的特点,在军用、民用方面被广泛应用。目前全世界应用的激光陀螺大部分是美国霍尼尔公司和利顿公司的机抖激光陀螺,而在国内激光陀螺主要是国防科技大学和航空618所的机抖激光陀螺,由于机抖激光陀螺的突出优点与广泛应用,机抖激光陀螺IMU已经成主要发展方向之一。
激光陀螺由于“闭锁效应”的存在,通常采用机械抖动偏频,这种激光陀螺即为机抖激光陀螺,机抖抖动偏频是在激光陀螺敏感轴方向输入偏置角速率,使其工作点全部或大部分从锁区偏置出来,因此机抖激光陀螺输出的信号包含抖动角速率分量,需要进行抖动解调,剥离由机械抖动引起的抖动角速率分量,保留IMU载体的运动角速率分量。抖动解调的方法有光学补偿法,整周期计数法,误差信号修正法和高速采样低通滤波方法。其中高速采样低通滤波方法解调精度高、设计简单,是目前主要的解调方式之一,该方法的关键是低通数字滤波器的设计。从数字滤波器的冲激响应来看,滤波器可分有限长单位冲激响应(FIR)滤波器和无限长单位冲激响应(IIR)滤波器。机抖激光陀螺IMU信号滤波采用FIR滤波器能够获得严格的线性相位,但由于机抖激光陀螺机械抖动噪声强度远远超出有用信号强度,要想达到理想的幅频特性,FIR滤波器阶数将会很高,运算量大,对IMU硬件系统提出了苛刻的要求,而且IMU输出的群延时也与FIR滤波器阶数成比例,因此FIR滤波器使得IMU输出数据延时大。机抖激光陀螺IMU采用IIR滤波器可以获得良好的噪声抑制效果,而且IIR滤波器阶数低,IIR滤波器包括巴特沃斯滤波器、切比雪夫I型滤波器、切比雪夫II型滤波器和椭圆滤波器,按照这些准则设计的IIR滤波器,不仅相位非线性严重,使得不同频率的IMU测量信号输出延时差别大,导致IMU测量精度下降,而且其群延时也很大,降低了IMU测量的实时性。
发明内容
本发明的技术解决问题是:克服现有技术的不足,提供一种机抖激光陀螺IMU数字滤波器设计方法,该方法设计的数字滤波器具有相位近似线性,运算量小,群延时小,群延时波动小的优点,降低了滤波器对IMU硬件系统的要求,提高了IMU测量实时性和时延补偿精度,为IMU高精度实时测量打下基础。
本发明技术解决方案为:一种机抖激光陀螺IMU数字滤波器设计方法,其特点在于包括下列步骤:
(1)根据机抖激光陀螺IMU性能要求,确定线性相位IIR数字滤波器带阻部分和低通部分的指标参数:带阻部分阻带衰减参数Astop1,中心频率f0;低通部分阻带衰减参数Astop2,低通部分通带纹波Apass2,低通部分群延时波动Gdw;线性相位IIR数字滤波器整体阻带衰减参数Astop。并确定机抖激光陀螺抖动噪声峰值频率fn和强度An以及机抖激光陀螺抖频波动范围frange;
(2)采用通带最平坦准则,以中心频率f0=fn,设计线性相位IIR数字滤波器带阻部分,获得其传递函数分子系数as1,as2,…,asn和分母系数bs0,bs1,bs2,…,bsn,其传递函数
(3)在机抖激光陀螺的抖频调节范围为3frange~5frange的情况下测试Hs(z),采集经Hs(z)滤波后的信号进行功率谱分析,如果抖动噪声峰值低于0.7An,调节Hs(z)的中心频率f0重新进行步骤(2)、步骤(3);
(4)在时域内,以线性相位FIR低通滤波器与线性相位IIR数字滤波器低通部分Hl(z)的误差为量测误差,利用最小二乘估计使Hl(z)相位线性化,获得Hl(z)初始分子系数a1(0),a2(0),…,an(0)和分母系数b0(0),b1(0)…bn(0),即初始
(5)在频域内,建立优化Hl(z)的目标函数J,J的相频误差评价函数权系数为Kp,幅频误差评价函数过渡带和阻带权系数为KM2,KM3;以a1(0),a2(0),…,an(0),b0(0),b1(0),b2(0)…bn(0)为初值,采用多变量无约束非线性极小值问题求解方法DFP(Daridon-Fletcher-Powell)法对J进行优化,获得频域优化后的Hl(z)的分子和分母系数al1(n),al2(n),…,aln(n),bl0(n),bl1(n)…bln(n);
(6)通过单位冲激响应测试法,测试频域优化后的Hl(z),并采集经Hl(z)滤波后机抖激光陀螺信号,进行功率谱分析,如果Hl(z)的通带纹波大于Apass2,同时减小Kp,KM2,KM3,如果Hl(z)的群延时波动大于Gdw,增大权系数Kp,如果机抖激光陀螺抖动噪声衰减低于Astop1,增大权系数KM3,重新进行步骤(6);
(7)利用公式H(z)=Hs(z)Hl(z)求线性相位IIR数字滤波器传递函数H(z),并采集滤波后信号进行功率谱分析,如果H(z)阻带内噪声衰减低于Astop,增大Astop1、Astop2,重新进行,求得H(z)的分子系数b0(n),b1(n),b2(n)…bn(n)和分母系数a1(n),a2(n),…,an(n),最终获得线性相位IIR数字滤波器传递函数 完成线性相位IIR数字滤波器设计。
本发明的原理是:机抖激光陀螺信号包含被测角速率信息,机械抖动噪声和其它噪声。其中机械抖动噪声具有频率集中,噪声强度大。一般单独采用FIR低通滤波器阶数会很高,其主要用于抑制抖动噪声,但抖动噪声频带很窄,高阶的FIR滤波器却将整个阻带信号全部衰减,FIR滤波器实际利用率低。因此本发明采用二阶或三阶线性相位IIR数字滤波器带阻部分,抑制机抖激光陀螺抖动噪声,该线性相位IIR数字滤波器带阻部分在低频段幅值波动和相位非线性比线性相位IIR数字滤波器低通部分小一个数量级以上,可以忽略,对整个数字滤波器通带内的幅相特性没有影响。为进一步保证整个数字滤波器的滤波器效果,设计一个线性相位FIR低通滤波器,由于IMU不要求数字滤波器具有严格的线性相位,滤波器非线性只要对IMU测量精度影响可以忽略即可,因此考虑到如果能够设计一个线性相位的IIR滤波器低通部分,整个滤波器运算量会进一步降低。但是求解IIR滤波器最优问题是一个多变量无约束优化问题,这种问题在参数较多时很难求得优化解。本发明采用时域设计与频域设计的优化方法,首先在时域内,利用最小二乘法IIR滤波器低通部分优化跟踪所期望频率特性的线性相位FIR低通滤波器的单位冲激响应,从而获得这个多变量无约束优化问题的初值,优化范围已经缩小,然后在频域中采用DFP法求解多变量无约束优化问题。通过这种时域与频域结合设计的线性相位IIR数字滤波器低通部分相位接近线性,而且相位线性度可以通过频域设计时权重系数调节,同时该滤波器阶数相对于具有同等幅频特性效果的FIR滤波器阶数降低超过一半,这种近似线性相位滤波器进一步降低了整个数字滤波器的运算量,而且保证了整个频段的幅频特性。
本发明与现有技术相比的优点在于:本发明兼具了常规IIR滤波和FIR滤波器不能同时具有的低阶小运算量和近似线性相位的优点,而且群延时小;在机抖激光陀螺IMU中应用,降低了数字滤波器对IMU硬件的要求,提高了IMU测量实时性,保证了IMU时延补偿精度,为IMU高精度实时测量奠定了基础。
附图说明
图1为本发明的具体实施的POS系统组成框图;
图2为本发明的算法应用的硬件组成框图;
图3为本发明的线性IIR滤波器设计的流程框图;
图4为本发明的线性相位IIR数字滤波器低通部分时域设计原理图;
图5为本发明的递推最小二乘估计的解算流程图。
具体实施方式
如图1所示,为本发明方法实施对象IMU所在的POS系统,IMU包括了激光陀螺仪模块、石英加速度计模块、数据采集电路模块和数据处理模块。本发明设计的滤波器是数据预处理模块的主要组成部分。如图2所示为滤波器所在硬件系统平台,图中标明了信号流向,滤波器程序将在图2中DSP最小子系统中实现。如图3所示,具体实施方法如下:
(1)通过4KHz采集机抖激光陀螺输出信号,采用功率谱密度分析方法,确定机抖激光陀螺抖动噪声峰值频率fn和强度An以及机抖激光陀螺抖频波动范围frange。线性相位IIR数字滤波器带阻部分和低通部分的指标参数:带阻部分阻带衰减参数Astop1=0.7An,中心频率f0=fn;低通部分阻带衰减参数Astop2=60dB,低通部分通带纹波Apass2=0.01dB,以IMU时延要求的波动范围的0.9倍为滤器低通部分群延时波动范围Gdw;线性相位IIR数字滤波器整体阻带衰减参数Astop应满足通带内信号强度超过阻带内噪声强度一个数量级以上;
(2)采用通带最平坦准则设计IIR滤波器带阻部分Hs(z),通过选择阶数为二阶或三阶的带阻部分,并选用对称于f0的带阻参数,使带阻部分的相位非线性相对于低通部分低一个数量级以上;
(3)在机抖激光陀螺的抖频调节范围为3frange~5frange的情况下测试Hs(z)。以4KHz采样频率采集Hs(z)滤波后的机抖激光陀螺2s信号,采集三组分析其中功率谱密度,求其均值,如果Hs(z)阻带衰减不到0.7An,将f0上移1Hz重新进行(2)(3);再将Hs(z)中心频率移动3frange~5frange,如果Hs(z)的阻带衰减不到0.7As,将f0上移0.2Hz重新进行步骤(2)、步骤(3);
(4)线性相位IIR数字滤波器低通部分时域设计。如图4所示线性相位IIR数字滤波器低通部分时域设计原理图,根据等纹波设计准则设计线性相位FIR低通滤波器Hf(z);以线性相位IIR数字滤波器低通部分H(z)输入输出组成状态方程量测矩阵,建立最小二乘参数估计器;然后在高斯白噪声激励下,利用最小二乘估计器使线性相位IIR数字滤波器低通部分优化跟踪上述线性相位FIR低通滤波器,获得时域设计的线性相位IIR数字滤波器低通部分,具体步骤如下:
(a)以待求系数a1,a2,…,an,b0,b1…bn建立Hl(z)的模型:
采用等纹波设计准则设计线性相位FIR低通滤波器其中h0,h1,…,hN-1为线性相位FIR低通滤波器传递函数系数,FIR低通滤波器阶数约为Hl(z)阶数的一倍,通带纹波应比Hl(z)大一个数量级以上,阻带衰减应比Hl(z)少一个数量级;
(b)建立最小二乘估计方程。以a1,a2,…,an,b0,b1,b2…bn组成最小二乘估计的估计状态X(k)=[a1,a2,…,an,b0,b1,…,bn]T,以Hl(z)输入x(k)和输出y(k)组成最小二乘估计量测矩阵H(k)=[-y(k-1),…,-y(k-n),x(k),…,x(k-n)],建立最小二乘估计方程:
Z(k)=H(k)X(k)+V(k) (2)
其中Z(k)为k时刻量测量,H(k)为k时刻状态转移矩阵,X(k)为k时刻系统状态量,V(k)为k时刻量测噪声矩阵;
采用带遗忘因子的最小二乘递推方法,算法如下:
K(k)=P(k-1)HT(k)[H(k)P(k-1)HT(k)+r(k)]-1 (3)
其中为k时刻系统状态估计量,为k-1时刻系统状态估计量,K(k)为k时刻增益矩阵,Z(k)为k时刻量测量,H(k)为k时刻状态转移矩阵,P(k)为k时刻系统估计量的协方差矩阵,P(k-1)为k-1时刻系统估计量的协方差矩阵,r(k)为k时刻遗忘因子,I为单位矩阵,最小二乘估计遗忘因子采用常值遗忘因子r(k)=λ,为兼顾最小二乘收敛速度和灵敏度取0.980<λ<0.998。最小二乘解算流程如图5所示;
(c)利用(b)中的最小二乘算法在高斯白噪声激励下,取最小二乘估计初值P(0)=rI,r>0,使Hl(z)逼近(a)中的线性相位FIR低通滤波器Hf(z),以获得具有线性相位的Hl(z),其系数为a1(0),a2(0),…,an(0),b0(0),b1(0),b2(0)…bn(0);
(5)线性相位IIR数字滤波器低通部分频域设计。建立线性相位IIR数字滤波器低通部分与两个待逼近FIR低通滤波器单位冲激响应的幅值和相位误差的评价函数,以时域设计结果为初值,采用DFP法对J进行优化,以获得具有所期望频率特性的线性相位IIR数字滤波器;然后将优化后的IIR滤波器用于机抖激光陀螺抖频抑制,测试线性相位IIR数字滤波器低通部分滤波效果。具体步骤如下:
(a)在频域内建立两个线性相位FIR低通滤波器与Hl(z)的幅值误差和相位误差的评价函数:
其中为线性相位IIR数字滤波器低通部分的传递函数,Hdm(ejω)和Hdp(ejω)为两个线性相位FIR低通滤波器传递函数;JM1,JM2,JM3分别为通带,过渡带和阻带的幅值误差的评价函数,JP为相位误差评价函数;lm0,lm1,lm2,lm3为幅频误差评价区间的分界点,lp0,lp1为通带相频特性的评价区间点;
(b)建立线性相位IIR数字滤波器低通部分的最优准则:
J=aJP1+β1JM1+β2JM2+β3JM3 (8)
评价函数JM1,JM2,JM3和JP有相当大差距,为使权系数直观反映各部分的相对权重,取
(9)
其中,X(n)=[a1(0),a2(0),…,b0(0),b1(0),b2(0)…bn(0)]T为时域设计得到的Hl(z)的时域初值,Kp,KM2,KM3是直观反映IIR滤波器幅相特性各部分权重的权系数;
(c)设计待逼近的两个FIR低通滤波器Hdm(ejω)和Hdp(ejω)。FIR低通滤波Hdm(ejω)具有待求线性相位IIR数字滤波器低通部分期望的幅频特性,其阶数为Hl(z)阶数的3~5倍,FIR低通滤波Hdp(ejω)具有所期望的相频特性,其阶数不超过Hl(z)阶数的2倍;
(d)根据Hl(z)幅值和相位权重,确定J权系数Kp,KM2,KM3,以a1(0),a2(0),…,an(0),b0(0),b1(0),b2(0)…bn(0)为初值利用DFP法求目标函数J极小值,以获得具有所期望频率特性的线性相位IIR数字滤波器低通部分分子系数al1(n),al2(n),…,aln(n)和分母系数bl0(n),bl1(n)…bln(n),传递函数
(e)通过单位冲激响应测试法,测试优化后的Hl(z),并将优化后的Hl(z)用于机抖激光陀螺抖频抑制,以4KHz采集2s滤波后信号,进行功率谱密度,如果Hl(z)的通带纹波大于Apass2,同时减小Kp,KM2,KM3至原来的0.9倍,如果Hl(z)的群延时波动大于Gdw,增大权系数Kp至原来的1.05倍,如果机抖激光陀螺抖动噪声衰减低于Astop1,增大权系数KM3至原来的1.05倍,重新进行(6);
(7)利用公式(10)求线性相位IIR数字滤波器H(z):
H(z)=Hs(z)Hl(z) (10)
将H(z)用于机抖激光陀螺噪声抑制,并以4KHz采集2s滤波后信号进行功率谱分析,如果H(z)阻带内噪声衰减低于Astop,将Astop1和Astop2增大5dB,重新进行步骤(2)-步骤(7),求得H(z)的分子系数b0(n),b1(n),b2(n)…bn(n)和分母系数a1(n),a2(n),…,an(n),最终获得线性相位IIR数字滤波器传递函数 完成线性相位IIR数字滤波器设计;在具有相同机抖激光陀螺噪声抑制效果的情况下本发明设计的线性相位IIR数字滤波器比线性相位FIR低通滤波器阶数低,运算量减小,群延时可减少50%以上,为IMU高精度实时测量奠定了良好的基础。
本发明未详细阐述部分属于本领域公知技术。
显然,对于本领域的普通技术人员来说,参照上文所述的实施例还可能做出其它的实施方式。本发明中的实施例都只是示例性的、而不是局限性的。所有的在本发明的权利要求技术方案的本质之内的修改都属于其所要求保护的范围。
Claims (7)
1.一种机抖激光陀螺惯性测量单元数字滤波器设计方法,其特征在于包括下列步骤:
(1)根据机抖激光陀螺惯性测量单元性能要求,确定线性相位IIR数字滤波器带阻部分和低通部分的指标参数,所述指标参数包括:带阻部分阻带衰减参数Astop1,中心频率f0;低通部分阻带衰减参数Astop2,低通部分通带纹波Apass2,低通部分群延时波动Gdw;线性相位IIR数字滤波器整体阻带衰减参数Astop;并确定机抖激光陀螺抖动噪声峰值频率fn和强度An以及机抖激光陀螺抖频波动范围frange;
(2)采用通带最平坦准则,以中心频率f0=fn,设计线性相位IIR数字滤波器带阻部分,获得其传递函数分母系数as1,as2,…,asn和分子系数bs0,bs1,bs2,…,bsn,其传递函数
(3)在机抖激光陀螺的抖频调节范围为3frange~5frange的情况下测试Hs(z),采集经Hs(z)滤波后的信号进行功率谱分析,如果抖动噪声峰值低于0.7An,调节Hs(z)的中心频率f0重新进行步骤(2)、步骤(3);
(4)在时域内,以线性相位FIR低通滤波器与线性相位IIR数字滤波器低通部分Hl(z)的误差为量测误差,利用最小二乘估计使Hl(z)相位线性化,获得Hl(z)初始分母系数a1(0),a2(0),…,an(0)和分子系数b0(0),b1(0)…bn(0),即初始
(5)在频域内,建立优化Hl(z)的目标函数J,J的相频误差评价函数权系数为Kp,幅频误差评价函数过渡带和阻带权系数为KM2,KM3;以a1(0),a2(0),…,an(0),b0(0),b1(0)…bn(0)为初值,采用DFP法对J进行优化,获得频域优化后的Hl(z)的分母和分子系数al1(n),al2(n),…,aln(n),bl0(n),bl1(n)…bln(n);
(6)通过单位冲激响应测试法,测试频域优化后的Hl(z),并采集经Hl(z)滤波后机抖激光陀螺信号,进行功率谱分析,如果Hl(z)的通带纹波大于Apass2,同时减小Kp,KM2,KM3,如果Hl(z)的群延时波动大于Gdw,增大权系数Kp,如果机抖激光陀螺抖动噪声衰减低于Astop1,增大权系数KM3,重新进行步骤(6);
(7)利用公式H(z)=Hs(z)Hl(z)求线性相位IIR数字滤波器传递函数H(z),并采集滤波后信号进行功率谱分析,如果H(z)阻带内噪声衰减低于Astop,增大Astop1、Astop2,重新进行步骤(2)-步骤(7),求得H(z)的分母系数a1(n),a2(n),…,an(n)和分子系数b0(n),b1(n)…bn(n),最终获得线性相位IIR数字滤波器传递函数 完成了线性相位IIR数字滤波器设计。
2.根据权利要求1所述的机抖激光陀螺惯性测量单元数字滤波器设计方法,其特征在于:所述步骤(1)中,通过4KHz~10KHz的采样频率采集机抖激光陀螺输出信号,利用功率谱分析方法,确定机抖激光陀螺抖动噪声峰值频率fn和强度An。
3.根据权利要求1所述的机抖激光陀螺惯性测量单元数字滤波器设计方法,其特征在于:所述步骤(1)中,以惯性测量单元时延要求波动范围的0.9倍为线性相位IIR数字滤波器低通部分群延时波动Gdw。
4.根据权利要求1所述的机抖激光陀螺惯性测量单元数字滤波器设计方法,其特征在于:所述步骤(1)中,线性相位IIR数字滤波器整体阻带衰减参数Astop应满足通带内信号强度超过阻带内噪声强度一个数量级以上。
5.根据权利要求1所述的机抖激光陀螺惯性测量单元数字滤波器设计方法,其特征在于:所述步骤(2)中,采用通带最平坦准则设计IIR数字滤波器带阻部分时,通过选择阶数为二阶或三阶的带阻部分,并选用对称于f0的带阻参数,使带阻部分的相位非线性相对于低通部分低一个数量级以上。
6.根据权利要求1所述的机抖激光陀螺惯性测量单元数字滤波器设计方法,其特征在于:所述步骤(4)中利用最小二乘估计使Hl(z)相位线性化,获得Hl(z)初始系数a1(0),a2(0),…,an(0),b0(0),b1(0)…bn(0)步骤如下:
(a)以待求系数a1,a2,…,an,b0,b1…bn建立Hl(z)模型,设计线性相位FIR低通滤波器其中h0,h1,…,hN-1为线性相位FIR低通滤波器传递函数系数,FIR低通滤波器阶数约为Hl(z)阶数的一倍,通带纹波应比Hl(z)大一个数量级以上,阻带衰减应比待求Hl(z)少一个数量级;
(b)建立最小二乘估计方程
以a1,a2,…,an,b0,b1…bn组成最小二乘估计的估计状态X(k)=[a1,a2,…,an,b0,b1,…,bn]T,以Hl(z)输入x(k)和输出y(k)组成最小二乘估计量测矩阵H(k)=[-y(k-1),…,-y(k-n),x(k),…,x(k-n)],建立最小二乘估计方程:
Z(k)=H(k)X(k)+V(k)
其中Z(k)为k时刻量测量,H(k)为k时刻状态转移矩阵,X(k)为k时刻系统状态量,V(k)为k时刻量测噪声矩阵;
采用带遗忘因子的最小二乘递推方法,算法如下:
K(k)=P(k-1)HT(k)[H(k)P(k-1)HT(k)+r(k)]-1
其中为k时刻系统状态估计量,为k-1时刻系统状态估计量,K(k)为k时刻增益矩阵,P(k)为k时刻系统估计量的协方差矩阵,P(k-1)为k-1时刻系统估计量的协方差矩阵,r(k)为k时刻遗忘因子,I为单位矩阵,最小二乘估计遗忘因子采用常值遗忘因子r(k)=λ,为兼顾最小二乘收敛速度和灵敏度取0.980<λ<0.998;
(c)利用步骤(b)中的最小二乘算法在高斯白噪声激励下,取最小二乘估计状态量初值估计量协方差矩阵初值P(0)=rI,r>0为常数,I为单位矩阵,使Hl(z)逼近(a)中的线性相位FIR低通滤波器Hf(z),以获得具有线性相位的Hl(z),其系数为a1(0),a2(0),…,an(0),b0(0),b1(0)…bn(0)。
7.根据权利要求1所述的机抖激光陀螺惯性测量单元数字滤波器设计方法,其特征在于:所述步骤(5)中建立优化Hl(z)的目标函数J步骤如下,
(a)在频域内建立两个线性相位FIR低通滤波器与Hl(z)的幅值误差和相位误差的评价函数:
其中为线性相位IIR数字滤波器低通部分的传递函数,Hdm(ejω)和Hdp(ejω)为两个线性相位FIR低通滤波器传递函数;JM1,JM2,JM3分别为通带,过渡带和阻带的幅值误差的评价函数,JP为相位误差评价函数;lm0,lm1,lm2,lm3为幅频误差评价区间的分界点,lp0,lp1为通带相频特性的评价区间点;
(b)建立线性相位IIR数字滤波器低通部分的最优准则:
J=aJP1+β1JM1+β2JM2+β3JM3
评价函数JM1,JM2,JM3和JP有相当大差距,为使权系数直观反映各部分的相对权重,取:
其中,X(n)=[a1(0),a2(0),…,b0(0),b1(0),b2(0)…bn(0)]T为时域设计得到的Hl(z)的初值,Kp,KM2,KM3是直观反映IIR数字滤波器幅相特性各部分权重的权系数;
(c)设计待逼近的两个FIR低通滤波器Hdm(ejω)和Hdp(ejω),FIR低通滤波Hdm(ejω)具有待求线性相位IIR数字滤波器期望的幅频特性,其阶数为Hl(z)阶数的3~5倍;FIR低通滤波Hdp(ejω)具有所期望的相频特性,其阶数不超过Hl(z)阶数的2倍;
(d)根据Hl(z)幅值和相位权重,确定J权系数Kp,KM2,KM3,以a1(0),a2(0),…,an(0),b0(0),b1(0),b2(0)…bn(0)为初值利用DFP法求目标函数J极小值,以获得具有所期望频率特性的线性相位IIR数字滤波器低通部分分母系数al1(n),al2(n),…,aln(n)和分子系数bl0(n),bl1(n)…bln(n),传递函数
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201210117027.XA CN102620729B (zh) | 2012-04-19 | 2012-04-19 | 一种机抖激光陀螺惯性测量单元数字滤波器设计方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201210117027.XA CN102620729B (zh) | 2012-04-19 | 2012-04-19 | 一种机抖激光陀螺惯性测量单元数字滤波器设计方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN102620729A CN102620729A (zh) | 2012-08-01 |
CN102620729B true CN102620729B (zh) | 2014-12-31 |
Family
ID=46560796
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201210117027.XA Expired - Fee Related CN102620729B (zh) | 2012-04-19 | 2012-04-19 | 一种机抖激光陀螺惯性测量单元数字滤波器设计方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN102620729B (zh) |
Families Citing this family (16)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN103152013B (zh) * | 2013-03-27 | 2016-03-30 | 北京众谱达科技有限公司 | 车载在用无线电设备检测系统滤波器单元 |
CN104154916B (zh) * | 2013-08-30 | 2018-11-30 | 北京航天发射技术研究所 | 一种基于激光陀螺捷联惯组的车载定位设备 |
CN104655136B (zh) * | 2015-02-17 | 2017-08-29 | 西安航天精密机电研究所 | 一种适用于激光捷联惯性导航系统的多凹点fir滤波方法 |
CN106840157A (zh) * | 2015-12-07 | 2017-06-13 | 上海新跃仪表厂 | 一种光纤陀螺惯性测量装置角频率特性实现方法 |
CN106840193A (zh) * | 2015-12-07 | 2017-06-13 | 上海新跃仪表厂 | 一种惯性测量线角耦合抑制方法 |
CN105547294B (zh) * | 2016-01-14 | 2018-02-23 | 中国人民解放军国防科学技术大学 | 二频机抖激光陀螺惯性测量单元最优安装配置的评估方法 |
CN109073381B (zh) * | 2016-05-11 | 2022-03-22 | 株式会社村田制作所 | 具有力反馈能力的副感测回路 |
CN107389980B (zh) * | 2017-07-25 | 2019-08-02 | Oppo广东移动通信有限公司 | 一种加速度计控制方法、加速度计控制装置及智能终端 |
CN110896303B (zh) * | 2018-09-12 | 2024-04-05 | 浙江菜鸟供应链管理有限公司 | 滤波方法和滤波装置以及存储介质 |
CN112378418A (zh) * | 2020-10-30 | 2021-02-19 | 哈尔滨理工大学 | 一种陀螺信号高阶低通滤波及滞后补偿方法 |
CN114838720A (zh) * | 2021-02-02 | 2022-08-02 | 湖南二零八先进科技有限公司 | 一种用于激光陀螺的自适应抖动剥除方法和装置 |
CN112965365B (zh) * | 2021-02-23 | 2023-03-31 | 浙江中智达科技有限公司 | Pid控制回路的模型辨识方法、装置、系统及存储介质 |
CN113595527B (zh) * | 2021-07-30 | 2023-10-20 | 国光电器股份有限公司 | 一种滤波参数确定方法、滤波方法及相关装置 |
CN113459733B (zh) * | 2021-08-06 | 2023-07-21 | 常州禄安汽车科技有限公司 | 一种判断车辆四轮同时缺气的方法 |
CN115133906A (zh) * | 2022-07-04 | 2022-09-30 | Oppo广东移动通信有限公司 | 一种滤波器设计方法及装置、存储介质 |
CN117879540B (zh) * | 2024-03-12 | 2024-06-18 | 西南应用磁学研究所(中国电子科技集团公司第九研究所) | 基于改进卡尔曼滤波的磁罗盘传感器自适应信号滤波方法 |
Citations (4)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN101281036A (zh) * | 2008-05-15 | 2008-10-08 | 哈尔滨工程大学 | 一种基于fpga的机抖激光陀螺抖动解调装置及解调方法 |
JP2010141780A (ja) * | 2008-12-15 | 2010-06-24 | Audio Technica Corp | Iirフィルタ設計方法 |
US20110153995A1 (en) * | 2009-12-18 | 2011-06-23 | Electronics And Telecommunications Research Institute | Arithmetic apparatus including multiplication and accumulation, and dsp structure and filtering method using the same |
CN102169127A (zh) * | 2010-12-22 | 2011-08-31 | 北京航空航天大学 | 一种具有自适应能力的机抖激光陀螺实时去抖方法 |
-
2012
- 2012-04-19 CN CN201210117027.XA patent/CN102620729B/zh not_active Expired - Fee Related
Patent Citations (4)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN101281036A (zh) * | 2008-05-15 | 2008-10-08 | 哈尔滨工程大学 | 一种基于fpga的机抖激光陀螺抖动解调装置及解调方法 |
JP2010141780A (ja) * | 2008-12-15 | 2010-06-24 | Audio Technica Corp | Iirフィルタ設計方法 |
US20110153995A1 (en) * | 2009-12-18 | 2011-06-23 | Electronics And Telecommunications Research Institute | Arithmetic apparatus including multiplication and accumulation, and dsp structure and filtering method using the same |
CN102169127A (zh) * | 2010-12-22 | 2011-08-31 | 北京航空航天大学 | 一种具有自适应能力的机抖激光陀螺实时去抖方法 |
Non-Patent Citations (3)
Title |
---|
IIR Digital Filter Design With New Stability;Aimin Jiang;《IEEE TRANSACTIONS ON SIGNAL PROCESSING》;20090331;全文 * |
JP特开2010141780A 2010.06.24 * |
激光陀螺POS惯性数据滤波及时延补偿;钟麦英;《中国惯性技术学报》;20111231;第19卷(第6期);全文 * |
Also Published As
Publication number | Publication date |
---|---|
CN102620729A (zh) | 2012-08-01 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN102620729B (zh) | 一种机抖激光陀螺惯性测量单元数字滤波器设计方法 | |
Nassar et al. | Modeling inertial sensor errors using autoregressive (AR) models | |
Arató | Linear stochastic systems with constant coefficients: a statistical approach | |
CN100368774C (zh) | 光纤陀螺惯测装置快速启动和精度保证的工程实现方法 | |
CN107036589B (zh) | 一种用于mems陀螺仪的角度测量系统及其方法 | |
CN102243080A (zh) | 高精度光纤陀螺带温度补偿的信号检测方法及装置 | |
CN108955727B (zh) | 一种光纤线圈性能评价方法 | |
CN105486225A (zh) | 一种抑制光强波动噪声的相位解调装置及解调方法 | |
CN114218750A (zh) | 基于数字滤波器的星载微推力器推力响应时间测量方法 | |
CN102055434A (zh) | 一种应用于惯性器件中数字滤波器的设计方法 | |
CN108196242A (zh) | 基于边沿检测的激光雷达计时方法及数据处理单元 | |
CN111879348A (zh) | 一种惯性仪表性能地面测试系统效能分析方法 | |
CN110726852A (zh) | 一种mems加速度计温度补偿方法 | |
CN109521488A (zh) | 用于卫星重力梯度数据的arma最优滤波模型构建方法 | |
CN106153029A (zh) | 二频机抖激光陀螺抖动信号抵消装置 | |
CN110888142B (zh) | 基于mems激光雷达测量技术的航天器隐藏目标点测量方法 | |
CN102589551A (zh) | 一种基于小波变换的船用光纤陀螺信号实时滤波方法 | |
CN104501832A (zh) | 一种改进型实用惯性传感器降噪装置 | |
CN105180946A (zh) | 基于宽频测量的卫星高精度姿态确定方法及系统 | |
CN107576323B (zh) | 基于fir和lms自适应滤波组合型光纤陀螺滤波方法 | |
Weiquan et al. | Research of the inertial navigation system with variable damping coefficients horizontal damping networks | |
Wang et al. | Design and Implementation of Strong Tracking Combined Filtering Algorithm for MEMS Gyroscope | |
CN210198392U (zh) | 一种新型mems谐振式陀螺仪测控装置 | |
Xia et al. | Adaptive Kalman filtering based on higher-order statistical analysis for digitalized silicon microgyroscope | |
Bolduc et al. | Determination of the seismograph phase response from the amplitude response |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
C06 | Publication | ||
PB01 | Publication | ||
C10 | Entry into substantive examination | ||
SE01 | Entry into force of request for substantive examination | ||
C14 | Grant of patent or utility model | ||
GR01 | Patent grant | ||
CF01 | Termination of patent right due to non-payment of annual fee |
Granted publication date: 20141231 Termination date: 20190419 |
|
CF01 | Termination of patent right due to non-payment of annual fee |