CN112611310B - 一种磁偶极子目标测距测向方法 - Google Patents
一种磁偶极子目标测距测向方法 Download PDFInfo
- Publication number
- CN112611310B CN112611310B CN202011461502.6A CN202011461502A CN112611310B CN 112611310 B CN112611310 B CN 112611310B CN 202011461502 A CN202011461502 A CN 202011461502A CN 112611310 B CN112611310 B CN 112611310B
- Authority
- CN
- China
- Prior art keywords
- magnetic
- point
- sign
- target
- gradient tensor
- 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
Images
Classifications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01B—MEASURING LENGTH, THICKNESS OR SIMILAR LINEAR DIMENSIONS; MEASURING ANGLES; MEASURING AREAS; MEASURING IRREGULARITIES OF SURFACES OR CONTOURS
- G01B7/00—Measuring arrangements characterised by the use of electric or magnetic techniques
- G01B7/004—Measuring arrangements characterised by the use of electric or magnetic techniques for measuring coordinates of points
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01B—MEASURING LENGTH, THICKNESS OR SIMILAR LINEAR DIMENSIONS; MEASURING ANGLES; MEASURING AREAS; MEASURING IRREGULARITIES OF SURFACES OR CONTOURS
- G01B7/00—Measuring arrangements characterised by the use of electric or magnetic techniques
- G01B7/02—Measuring arrangements characterised by the use of electric or magnetic techniques for measuring length, width or thickness
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01C—MEASURING DISTANCES, LEVELS OR BEARINGS; SURVEYING; NAVIGATION; GYROSCOPIC INSTRUMENTS; PHOTOGRAMMETRY OR VIDEOGRAMMETRY
- G01C1/00—Measuring angles
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01R—MEASURING ELECTRIC VARIABLES; MEASURING MAGNETIC VARIABLES
- G01R33/00—Arrangements or instruments for measuring magnetic variables
- G01R33/02—Measuring direction or magnitude of magnetic fields or magnetic flux
- G01R33/0206—Three-component magnetometers
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01R—MEASURING ELECTRIC VARIABLES; MEASURING MAGNETIC VARIABLES
- G01R33/00—Arrangements or instruments for measuring magnetic variables
- G01R33/02—Measuring direction or magnitude of magnetic fields or magnetic flux
- G01R33/022—Measuring gradient
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
- G06F17/10—Complex mathematical operations
- G06F17/15—Correlation function computation including computation of convolution operations
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
- G06F17/10—Complex mathematical operations
- G06F17/16—Matrix or vector computation, e.g. matrix-matrix or matrix-vector multiplication, matrix factorization
-
- Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
- Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
- Y02A—TECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE
- Y02A90/00—Technologies having an indirect contribution to adaptation to climate change
- Y02A90/30—Assessment of water resources
Landscapes
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Engineering & Computer Science (AREA)
- Mathematical Physics (AREA)
- Computational Mathematics (AREA)
- Data Mining & Analysis (AREA)
- Mathematical Analysis (AREA)
- Mathematical Optimization (AREA)
- Theoretical Computer Science (AREA)
- Pure & Applied Mathematics (AREA)
- Condensed Matter Physics & Semiconductors (AREA)
- Algebra (AREA)
- Databases & Information Systems (AREA)
- Software Systems (AREA)
- General Engineering & Computer Science (AREA)
- Computing Systems (AREA)
- Radar, Positioning & Navigation (AREA)
- Remote Sensing (AREA)
- Measuring Magnetic Variables (AREA)
Abstract
本发明公开了一种磁偶极子目标测距测向方法,使用测量轴对称配置的三轴磁强计十字阵列测量阵列中心处的磁偶极子磁梯度张量,依据磁矩大小与磁性目标姿态无关的特点可以得到被探测目标的磁矩大小范围,由该点磁梯度张量的比值以及目标磁矩大小范围反演出被探测目标与十字阵列中心的距离和十字阵列测量坐标系下的磁偶极子目标的方位角。本发明所提的方法无需磁梯度张量测量***对磁性目标的磁场矢量的测量,也无需磁梯度张量测量***移动以测量空间多个点的磁梯度张量值;本发明所提的方法能使用三轴磁强计十字阵列作为磁梯度张量测量***,所使用的三轴磁强计相对较少,降低了测量***的复杂性和重量,便于轻型平台搭载使用。
Description
技术领域
本发明涉及一种磁偶极子目标测距测向方法,特别是一种基于磁梯度张量分量比值的磁偶极子目标测距测向方法,属于磁性目标探测及其参数反演技术领域
背景技术
基于磁异常场的被动式目标探测技术具有隐蔽性好、能连续侦测、效率高、使用简单可靠、反应迅速等优点,还不受天文、气象和水文等的影响,可作为声学探测手段的有益补充,在水下目标探测与追踪、矿产勘探、体内微型诊疗装置的无创定位等方面有着重要的应用价值。
作为共轭抑制背景磁场的差分测量方式,磁梯度张量的测量结果受地磁场的倾角和偏角的影响小,所构成的张量不变量无需特殊处理就可很好地描述磁场源。因此,基于磁梯度张量的磁性目标参数反演技术得到了关注与研究。
Wiegert等提出了张量约缩量STAR的定位方法,张量约缩量为球等值面的不变量(Roy Wiegert,J.Oeschger.Generalized magnetic gradient contraction basedmethod for detection,localization and discrimination of underwater mines andunexploded ordnance[J].Mine Warfare&Ship Self-Defence,2010,2:1325-1332)。Clark利用磁梯度张量特征值计算归一化磁源强度,结合梯度张量确定磁源位置的唯一解(DavidA Clark.New methods for interpretation of magnetic vector and gradient tensordata I:Eigenvector analysis and the normalized source strength[J].Explorationgeophysics,2012,43(4):267-282)。Teixeira等利用张量矩阵最小特征值对应的特征向量与相对位置矢量之间的正交关系,提出一种张量欧拉反褶积和磁场梯度张量特征分析的联合反演方法,增强了反演算法的收敛性(Francisco Curado Teixeira,AntónioPascoal.Magnetic navigation and tracking of underwater vehicles[J].IFACProceedings Volumes,2013,46(33):239-244)。Lee等用磁梯度张量导出的两个标量参数来表征目标的磁矩和位置,再由梯度法求解位置并实验验证这种方法(Kok-Meng Lee,MinLi.Magnetic tensor sensor for gradient-based localization of ferrous objectin geomagnetic field[J].IEEE Transactions on Magnetics,2016,52(8):1-10)。
海军工程大学提出了一种基于等效衰减磁矩的定位方法,分析了地磁场起伏、三轴磁强计噪声及测量基线长等对定位结果的影响(贾文抖,林春生,孙玉绘,翟国君.地磁场中磁目标的等效衰减磁矩定位方法[J].兵工学报,2018,39(2):283-289)。海军航空工程学院提出了一种基于正六面体磁梯度张量***的单点定位改进方法,消除了张量约缩量STAR法的***误差(吕俊伟,迟铖,于振涛,毕波,宋庆善.磁梯度张量不变量的椭圆误差消除方法研究[J].物理学报,2015,64(9):52-59)。陆军工程大学研究了一种两点磁梯度张量定位算法并进行了铁磁体定位实验(Gang Yin,Yingtang Zhang,Zhining Li,Hongbo Fan andGuoquan Ren.Detection of ferromagnetic target based on mobile magneticgradient tensor system[J].Journal of Magnetism and Magnetic Material,2016,402(6):1-7)。
通常的磁梯度张量定位方法需要测量磁性目标的磁场矢量或者测量磁性目标在空间多个点的磁梯度张量。一方面,在地球背景场存在的条件下不易测量磁性目标的磁场矢量;另一方面,要想测量空间多个点的磁梯度张量需要移动磁梯度张量测量***或者用更多的三轴磁强计构建合适的测量阵列,这增加了测量***的复杂性。
发明内容
针对上述现有技术,本发明要解决的技术问题是提供一种基于磁梯度张量分量比值的磁偶极子目标测距测向方法,所使用的三轴磁强计相对较少,降低测量***的复杂性和重量,便于轻型平台搭载使用。
为解决上述技术问题,本发明的一种磁偶极子目标测距测向方法,包括以下步骤:
步骤1、将四个同规格的三轴磁强计按照测量轴对称配置方式安装在由非磁性材料制作而成的十字基座上,构成三轴磁强计十字阵列;四个三轴磁强计位置点分别为A、B、C和D点,其中A点和C点关于y轴对称,B点和D点关于x轴对称;
步骤2、读取四个三轴磁强计的测量输出,计算十字阵列中心处的磁梯度张量的独立分量值a、b、d、e和f,具体为:
式中,δx和δy分别为三轴磁强计的xs轴和ys轴偏离zs轴的位置误差,Lx为十字阵列的点A和点C之间的长度,Ly为点B和点D之间的长度,BAx、BAy和BAz为点A处的磁性体磁场的三分量,BBx、BBy和BBz为点B处的磁性体磁场的三分量,BCx、BCy和BCz为点C处的磁性体磁场的三分量,BDx、BDy和BDz为点D处的磁性体磁场的三分量;
步骤3、设位置坐标向量X=[x,y,z]T∈R3,构造关于X的方程式,建立优化目标函数F(X),具体为:
F(X)=|f1(X)|+|f2(X)|+|f3(X)|+|f4(X)|
kax=(3r2-5x2)x,kay=(r2-5x2)y,kaz=(r2-5x2)z,kbx=(r2-5y2)x,kby=(3r2-5y2)y,kbz=(r2-5y2)z,kdx=kay,kdy=kbx,kdz=-5xyz,kex=kdz,key=kbz,kez=y(r2-5z2),kfx=kaz,kfy=kdz,kfz=(r2-5z2)x,ηab=a/b,ηad=a/d,ηae=a/e和ηaf=a/f,x、y和z为磁梯度张量的测量点相对于磁偶极子的空间坐标,r为磁梯度张量的测量点与磁偶极子之间的距离;
式中,h1(X)=xv31+yv32+zv33,v31、v32和v33是磁梯度张量矩阵的第三个特征值对应的特征向量元素;
步骤5、计算得到磁偶极子磁矩的三分量mx、my和mz,具体为:
式中,矩阵[Kij]的元素分别为K11=(3r2-5x2)x、K12=(r2-5x2)y、K13=(r2-5x2)z、K21=(r2-5y2)x、K22=(3r2-5y2)y、K23=(r2-5y2)z、K31=K12、K32=K21、K33=-5xyz、K41=K33、K42=K23、K43=(r2-5z2)y、K51=K13、K52=K33和K53=(r2-5z2)x;
步骤6、求解关于k的多项式,得到k的六个根,具体为:
A6k6+A5k5+A4k4+A3k3+A2k2+A1k+A0=0
式中,A6=d2(a+2b)-e2(a-b)+2def,A5=-2d[(a-b)(a+2b)+(d2+e2+f2)],A4=(a-b)2(a+2b)+d2(4a-7b)+(f2-2e2)(a-b)+6def,A3=-4d[(a-b)2+(-d2+e2+f2)],A2=(a-b)2(2a+b)+d2(4b-7a)+(2f2-e2)(a-b)+6def,A1=2d[(a-b)(2a+b)-(d2+e2+f2)],A0=d2(2a+b)+f2(a-b)+2def。
步骤7、舍弃k的复数根,根据k的剩余实数根计算q值,具体为:
步骤8、根据十字对称阵列的中心点位置处于偶极子磁性目标的上方或下方以及建立的测量坐标系确定z的正负,并根据k、q和t的值计算z值,具体为:
其中t是磁偶极子磁场大小,当十字对称阵列的中心点位置处于偶极子磁性目标的上方,z为负,当十字对称阵列的中心点位置处于偶极子磁性目标的下方,z为正。
分别计算出x和y,具体为:
x=kqz
y=qz
得到位置坐标的实数集合{xl,yl,zl|l=1,2,…,Nr},Nr为k的实数根的个数;
步骤9、利用{xl,yl,zl|l=1,2,…,Nr}计算出误差集合{Δl},具体为:
步骤10、从{Δl}寻找最小值对应的l,即
则k的值选取为:
k=kL
然后得到q和z的值;
步骤11、确定x和y的符号,具体为:
式中,sign(·)表示求符号值,磁偶极子目标相当于阵列中心点P的位置坐标xQ和yQ的符号分别为sign(xQ)=-sign(x)和sign(yQ)=-sign(y);
计算磁偶极子目标的方位角ψ∈[0,2π),具体为:
当sign(xQ)为“+”且sign(yQ)为“+”时,ψ=arctan(1/k);
当sign(xQ)为“+”且sign(yQ)为“-”时,ψ=2π+arctan(1/k);
当sign(xQ)为“-”且sign(yQ)为“+”或“-”时,ψ=π+arctan(1/k);
其中,当|x|-0≤A时且sign(yQ)为“+”时,ψ=π/2,;当|x|-0≤A时且sign(yQ)为“-”时,ψ=3π/2,A为给定值,A≤10-4;
步骤12、利用距离值r计算得到磁偶极子磁场和磁矩的大小,具体为:
本发明的有益效果:本发明提出一种基于磁梯度张量分量比值的磁偶极子目标测距测向方法,使用测量轴对称配置的三轴磁强计十字阵列测量阵列中心处的磁偶极子磁梯度张量,依据磁矩大小与磁性目标姿态无关的特点可以得到被探测目标的磁矩大小范围,由该点磁梯度张量的比值以及目标磁矩大小范围反演出被探测目标与十字阵列中心的距离和十字阵列测量坐标系下的磁偶极子目标的方位角。本发明所提的方法无需磁梯度张量测量***对磁性目标的磁场矢量的测量,也无需磁梯度张量测量***移动以测量空间多个点的磁梯度张量值;本发明所提的方法能使用三轴磁强计十字阵列作为磁梯度张量测量***,所使用的三轴磁强计相对较少,降低了测量***的复杂性和重量,便于轻型平台搭载使用。
相比于现有的基于磁梯度张量的磁偶极子参数反演方法,本发明提出一种基于磁梯度张量分量比值的磁偶极子目标测距测向方法,无需磁梯度张量测量***测量磁性目标的磁场矢量,也无需磁梯度张量测量***移动以测量空间多个点的磁梯度张量值,拓展了方法的应用对象;本发明所提的方法能使用三轴磁强计十字阵列作为磁梯度张量测量***,所使用的三轴磁强计相对较少,降低了测量***的复杂性和重量,便于轻型平台搭载使用,具有很好的实际应用价值;所提的方法通过磁梯度张量分量之间的比值运算消减了介质磁导率的影响。
附图说明
图1是三轴磁强计敏感轴对称配置的十字阵列;
图2(a)是测距误差随三轴磁强计噪声标准差的变化曲线;
图2(b)是测向误差随三轴磁强计噪声标准差的变化曲线;
图3(a)是测距误差随磁矩大小的初始估算误差的变化曲线;
图3(b)是测向误差随磁矩大小的初始估算误差的变化曲线。
具体实施方式
下面结合附图对本发明具体实施方式做进一步说明。
本发明的一种基于磁梯度张量分量比值的磁偶极子目标测距测向方法的实现原理具体为:
三轴磁强计的三个正交测量轴通常不相交于一点,也就是说三轴磁强计的xs轴与其zs轴之间存在δx的位置偏差,三轴磁强计的ys轴与其zs轴之间存在δy的位置偏差。为了减少测量误差,采用如图1所示的三轴磁强计测量轴对称配置的十字阵列测量阵列中心处的磁梯度张量的独立分量值,四个三轴磁强计的zs轴与十字阵列坐标系的z轴一致,且其分别放置在点A、B、C和D处。
点A和点C之间的长度为Lx,点B和点D之间的长度为Ly,则磁梯度张量矩阵的5个独立分量a、b、d、e和f的测量值分别为
式中,BAx、BAy和BAz为点A处的磁性体磁场的三分量,BBx、BBy和BBz为点B处的磁性体磁场的三分量,BCx、BCy和BCz为点C处的磁性体磁场的三分量,BDx、BDy和BDz为点D处的磁性体磁场的三分量。
对于磁偶极子磁性体而言,其磁梯度张量矩阵的5个独立分量a、b、d、e和f的表达式分别为
式中,μ为磁导率,x、y和z为磁梯度张量的测量点相对于磁偶极子的空间坐标,r为磁梯度张量的测量点与磁偶极子之间的距离,mx、my和mz为磁偶极子磁矩的三分量。
根据式(2),求得独立分量之间的比值分别为
式中,kax=(3r2-5x2)x,kay=(r2-5x2)y,kaz=(r2-5x2)z,kbx=(r2-5y2)x,kby=(3r2-5y2)y,kbz=(r2-5y2)z,kdx=kay,kdy=kbx,kdz=-5xyz,kex=kdz,key=kbz,kez=y(r2-5z2),kfx=kaz,kfy=kdz,kfz=(r2-5z2)x。
由于mx、my和mz不可能全部为零,因此由式(3)可得
式(4)可改写为
式中,X=[x,y,z]T∈R3为阵列中心点相当于磁偶极子的位置坐标向量。
式(5)可改写为
式(6)可改写为
式(7)可改写为
磁偶极子的磁梯度张量场还存在式(12)所示的关系。
h1(X)=xv31+yv32+zv33=0 (12)
式中,v31、v32和v33是磁梯度张量矩阵的第三个特征值对应的特征向量元素。
根据探测目标的类型和归一化强度u确定探测距离的范围,即有
通过求解式(14)所示的约束最优化问题计算磁偶极子目标与三轴磁强计十字阵列中心之间的距离。
式中,F(X)=|f1|+|f2|+|f3|+|f4|。
再由式(15)计算出mx、my和mz。
式中,矩阵[Kij]的元素分别为K11=(3r2-5x2)x、K12=(r2-5x2)y、K13=(r2-5x2)z、K21=(r2-5y2)x、K22=(3r2-5y2)y、K23=(r2-5y2)z、K31=K12、K32=K21、K33=-5xyz、K41=K33、K42=K23、K43=(r2-5z2)y、K51=K13、K52=K33和K53=(r2-5z2)x。
求解式(16)所示的关于k的多项式,得到k的六个根。
A6k6+A5k5+A4k4+A3k3+A2k2+A1k+A0=0 (16)
式中,A6=d2(a+2b)-e2(a-b)+2def,A5=-2d[(a-b)(a+2b)+(d2+e2+f2)],A4=(a-b)2(a+2b)+d2(4a-7b)+(f2-2e2)(a-b)+6def,A3=-4d[(a-b)2+(-d2+e2+f2)],A2=(a-b)2(2a+b)+d2(4b-7a)+(2f2-e2)(a-b)+6def,A1=2d[(a-b)(2a+b)-(d2+e2+f2)],A0=d2(2a+b)+f2(a-b)+2def。
舍弃k的复数根,将k的剩余实数根代入式(17)计算q值。
根据十字对称阵列的中心点位置处于偶极子磁性目标的上方或下方以及建立的测量坐标系,判断z的正负,以此确定式(18)中的正负号,并将k、q和t的值代入式(18)计算z值。
再根据式(19)和式(20)分别计算出
x=kqz (19)
y=qz (20)
得到位置坐标的实数集合{xl,yl,zl|l=1,2,…,Nr},Nr为k的实数根的个数,将这些空间坐标的实数解代入式(21)计算出误差集合{Δl}。
从{Δl}寻找最小值对应的l,即
则k的值选取为
k=kL (23)
由此也确定出q和z的值。
根据式(24)确定x和y的符号,即
式中,sign(·)表示求符号值。磁偶极子目标相当于阵列中心点P的位置坐标xQ和yQ的符号分别为sign(xQ)=-sign(x)和sign(yQ)=-sign(y)。
由表1所描述的计算方法计算出磁偶极子目标的方位角ψ∈[0,2π)。
表1方位角ψ的计算公式
将距离值r代入式(25)和式(26)还可计算出磁偶极子磁场大小t和磁矩大小m。
结合图1,本发明具体实施方式包括以下步骤:
步骤1、将四个同规格的三轴磁强计按照测量轴对称配置方式安装在由非磁性材料制作而成的十字基座上,构成如图1所示的三轴磁强计十字阵列。
步骤2、读取四个三轴磁强计的测量输出,根据式(1)计算十字阵列中心处的磁梯度张量的独立分量值a、b、d、e和f。
式中,δx和δy分别为三轴磁强计的xs轴和ys轴偏离zs轴的位置误差,Lx为十字阵列的点A和点C之间的长度,Ly为点B和点D之间的长度,BAx、BAy和BAz为点A处的磁性体磁场的三分量,BBx、BBy和BBz为点B处的磁性体磁场的三分量,BCx、BCy和BCz为点C处的磁性体磁场的三分量,BDx、BDy和BDz为点D处的磁性体磁场的三分量。
步骤3、设位置坐标向量X=[x,y,z]T∈R3,构造关于X的方程式,按式(2)建立优化目标函数F(X)。
F(X)=|f1(X)|+|f2(X)|+|f3(X)|+|f4(X)| (2)
kax=(3r2-5x2)x,kay=(r2-5x2)y,kaz=(r2-5x2)z,kbx=(r2-5y2)x,kby=(3r2-5y2)y,kbz=(r2-5y2)z,kdx=kay,kdy=kbx,kdz=-5xyz,kex=kdz,key=kbz,kez=y(r2-5z2),kfx=kaz,kfy=kdz,kfz=(r2-5z2)x,ηab=a/b,ηad=a/d,ηae=a/e和ηaf=a/f。
式中,h1(X)=xv31+yv32+zv33,v31、v32和v33是磁梯度张量矩阵的第三个特征值对应的特征向量元素。
步骤5、按式(5)计算出mx、my和mz。
式中,矩阵[Kij]的元素分别为K11=(3r2-5x2)x、K12=(r2-5x2)y、K13=(r2-5x2)z、K21=(r2-5y2)x、K22=(3r2-5y2)y、K23=(r2-5y2)z、K31=K12、K32=K21、K33=-5xyz、K41=K33、K42=K23、K43=(r2-5z2)y、K51=K13、K52=K33和K53=(r2-5z2)x。
步骤6、求解式(6)所示的关于k的多项式,得到k的六个根。
A6k6+A5k5+A4k4+A3k3+A2k2+A1k+A0=0 (6)
式中,A6=d2(a+2b)-e2(a-b)+2def,A5=-2d[(a-b)(a+2b)+(d2+e2+f2)],A4=(a-b)2(a+2b)+d2(4a-7b)+(f2-2e2)(a-b)+6def,A3=-4d[(a-b)2+(-d2+e2+f2)],A2=(a-b)2(2a+b)+d2(4b-7a)+(2f2-e2)(a-b)+6def,A1=2d[(a-b)(2a+b)-(d2+e2+f2)],A0=d2(2a+b)+f2(a-b)+2def。
步骤7、舍弃k的复数根,将k的剩余实数根代入式(7)计算q值。
步骤8、根据十字对称阵列的中心点位置处于偶极子磁性目标的上方或下方以及建立的测量坐标系,判断z的正负,以此确定式(8)中的正负号,并将k、q和t的值代入式(8)计算z值。
再根据式(9)和式(10)分别计算出
x=kqz (9)
y=qz (10)
得到位置坐标的实数集合{xl,yl,zl|l=1,2,…,Nr},Nr为k的实数根的个数。
步骤9、将{xl,yl,zl|l=1,2,…,Nr}代入式(11)计算出误差集合{Δl}。
步骤10、从{Δl}寻找最小值对应的l,即
则k的值选取为
k=kL (13)
由此也确定出q和z的值。
步骤11、根据式(14)确定x和y的符号,即
式中,sign(·)表示求符号值。磁偶极子目标相当于阵列中心点P的位置坐标xQ和yQ的符号分别为sign(xQ)=-sign(x)和sign(yQ)=-sign(y)。
由表1所描述的计算方法计算出磁偶极子目标的方位角ψ∈[0,2π)。
表1方位角ψ的计算公式
步骤12、将距离值r代入式(15)和式(16)得到磁偶极子磁场和磁矩的大小。
为检验本发明提出的测距测向方法的可行性,进行如下数值仿真实验。磁偶极子的磁矩三分量分别为mx=8×108A·m、my=6×108A·m和mz=3×108A·m,其位置坐标真值为xQ=180m、yQ=220m和zQ=-150m,三轴磁强计十字阵列在x和y方向上的基线长分别为Lx=0.8m和Ly=0.8m,测量轴的位置偏差为δx=2mm和δy=2mm。三轴磁强计在每个轴的噪声是相互独立的高斯过程,其均值为0,标准差为σm。磁偶极子磁矩大小的初始估算误差为ρ。
仿真实例一:设磁偶极子磁矩大小的初始估算误差为10%,即和在不同的三轴磁强计噪声标准差的情况下,进行磁偶极子目标测距与测向的50次蒙特卡洛仿真实验。得到三轴磁强计十字阵列对磁偶极子目标的测距与测向误差随三轴磁强计噪声标准差的变化曲线如图2(a)和2(b)所示,其中图2(a)为σm从0nT增加到0.8nT时的距离测量误差曲线,图2(b)为σm从0nT增加到0.8nT时的方向角测量误差曲线。由图2(a)和2(b)可知,距离测量误差随σm的增加而呈增大趋势,但增大的幅度不大;方向角测量误差随σm的增加而呈线性增大趋势。
仿真实例二:设三轴磁强计噪声标准差σm为0.1nT,在不同的磁偶极子磁矩大小的初始估算误差的情况下,进行磁偶极子目标测距与测向的50次蒙特卡洛仿真实验。得到三轴磁强计十字阵列对磁偶极子目标的测距与测向误差随磁偶极子磁矩大小的初始估算误差的变化曲线如图3(a)和3(b)所示,其中图3(a)为ρ从2%增加到20%时的距离测量误差曲线,图3(b)为ρ从2%增加到20%时的方向角测量误差曲线。由图3(a)和3(b)可知,在ρ较大时距离测量误差随ρ的增加而呈线性增大趋势,在ρ小时距离测量误差随ρ的增加而几乎不变;方向角测量误差随ρ的增加而变化不大。
以测量轴对称配置的三轴磁强计十字阵列为磁梯度张量测量装置,进行了磁偶极子磁性目标的测距测向的蒙特卡洛仿真实验。仿真实验给出了在不同的三轴磁强计噪声标准差的情况下三轴磁强计十字阵列对磁偶极子目标的测距与测向误差随三轴磁强计噪声标准差的变化曲线,结果表明距离测量误差随三轴磁强计测量噪声标准差的增加而呈增大趋势,但增大的幅度不大;方向角测量误差随轴磁强计测量噪声标准差的增加而呈线性增大趋势。仿真实验给出了在不同的磁偶极子磁矩大小的初始估算误差的情况下,三轴磁强计十字阵列对磁偶极子目标的测距与测向误差随磁偶极子磁矩大小的初始估算误差的变化曲线,结果表明在磁矩大小的初始估算误差较大时距离测量误差随初始估算误差的增加而呈线性增大趋势,在磁矩大小的初始估算误差小时距离测量误差随初始估算误差的增加而几乎不变;方向角测量误差随初始估算误差的增加而变化不大。仿真实验证明了本发明方法的可行性。
相比于一般的基于磁梯度张量的磁偶极子参数反演方法,本发明提出的方法无需磁梯度张量测量***测量磁性目标的磁场矢量,也无需磁梯度张量测量***移动以测量空间多个点的磁梯度张量值,拓展了方法的应用对象;本发明所提的磁梯度张量测量方案所需的三轴磁强计相对较少,构造简单,降低了测量***的复杂性和重量,便于轻型平台搭载使用,三轴磁强计测量轴对称配置也降低了磁梯度张量差分式测量的误差,具有很好的实际应用价值;磁梯度张量分量之间的比值运算也消减了介质磁导率的影响。
Claims (1)
1.一种磁偶极子目标测距测向方法,其特征在于,包括以下步骤:
步骤1、将四个同规格的三轴磁强计按照测量轴对称配置方式安装在由非磁性材料制作而成的十字基座上,构成三轴磁强计十字阵列;四个三轴磁强计位置点分别为A、B、C和D点,其中A点和C点关于y轴对称,B点和D点关于x轴对称;
步骤2、读取四个三轴磁强计的测量输出,计算十字阵列中心处的磁梯度张量的独立分量值a、b、d、e和f,具体为:
式中,δx和δy分别为三轴磁强计的xs轴和ys轴偏离zs轴的位置误差,Lx为十字阵列的点A和点C之间的长度,Ly为点B和点D之间的长度,BAx、BAy和BAz为点A处的磁性体磁场的三分量,BBx、BBy和BBz为点B处的磁性体磁场的三分量,BCx、BCy和BCz为点C处的磁性体磁场的三分量,BDx、BDy和BDz为点D处的磁性体磁场的三分量;
步骤3、设位置坐标向量X=[x,y,z]T∈R3,构造关于X的方程式,建立优化目标函数F(X),具体为:
F(X)=|f1(X)|+|f2(X)|+|f3(X)|+|f4(X)|
kax=(3r2-5x2)x,kay=(r2-5x2)y,kaz=(r2-5x2)z,kbx=(r2-5y2)x,kby=(3r2-5y2)y,kbz=(r2-5y2)z,kdx=kay,kdy=kbx,kdz=-5xyz,kex=kdz,key=kbz,kez=y(r2-5z2),kfx=kaz,kfy=kdz,kfz=(r2-5z2)x,ηab=a/b,ηad=a/d,ηae=a/e和ηaf=a/f,x、y和z为磁梯度张量的测量点相对于磁偶极子的空间坐标,r为磁梯度张量的测量点与磁偶极子之间的距离;
式中,h1(X)=xv31+yv32+zv33,v31、v32和v33是磁梯度张量矩阵的第三个特征值对应的特征向量元素;
步骤5、计算得到磁偶极子磁矩的三分量mx、my和mz,具体为:
式中,矩阵[Kij]的元素分别为K11=(3r2-5x2)x、K12=(r2-5x2)y、K13=(r2-5x2)z、K21=(r2-5y2)x、K22=(3r2-5y2)y、K23=(r2-5y2)z、K31=K12、K32=K21、K33=-5xyz、K41=K33、K42=K23、K43=(r2-5z2)y、K51=K13、K52=K33和K53=(r2-5z2)x;
步骤6、求解关于k的多项式,得到k的六个根,具体为:
A6k6+A5k5+A4k4+A3k3+A2k2+A1k+A0=0
式中,A6=d2(a+2b)-e2(a-b)+2def,A5=-2d[(a-b)(a+2b)+(d2+e2+f2)],A4=(a-b)2(a+2b)+d2(4a-7b)+(f2-2e2)(a-b)+6def,A3=-4d[(a-b)2+(-d2+e2+f2)],A2=(a-b)2(2a+b)+d2(4b-7a)+(2f2-e2)(a-b)+6def,A1=2d[(a-b)(2a+b)-(d2+e2+f2)],A0=d2(2a+b)+f2(a-b)+2def;
步骤7、舍弃k的复数根,根据k的剩余实数根计算q值,具体为:
步骤8、根据十字对称阵列的中心点位置处于偶极子磁性目标的上方或下方以及建立的测量坐标系确定z的正负,并根据k、q和t的值计算z值,具体为:
其中t是磁偶极子磁场大小,当十字对称阵列的中心点位置处于偶极子磁性目标的上方,z为负,当十字对称阵列的中心点位置处于偶极子磁性目标的下方,z为正;
分别计算出x和y,具体为:
x=kqz
y=qz
得到位置坐标的实数集合{xl,yl,zl|l=1,2,…,Nr},Nr为k的实数根的个数;
步骤9、利用{xl,yl,zl|l=1,2,…,Nr}计算出误差集合{Δl},具体为:
步骤10、从{Δl}寻找最小值对应的l,即
则k的值选取为:
k=kL
然后得到q和z的值;
步骤11、确定x和y的符号,具体为:
式中,sign()表示求符号值,磁偶极子目标相当于阵列中心点P的位置坐标xQ和yQ的符号分别为sign(xQ)=-sign(x)和sign(yQ)=-sign(y);
计算磁偶极子目标的方位角ψ∈[0,2π),具体为:
当sign(xQ)为“+”且sign(yQ)为“+”时,ψ=arctan(1/k);
当sign(xQ)为“+”且sign(yQ)为“-”时,ψ=2π+arctan(1/k);
当sign(xQ)为“-”且sign(yQ)为“+”或“-”时,ψ=π+arctan(1/k);
其中,当|x|-0≤A时且sign(yQ)为“+”时,ψ=π/2,;当|x|-0≤A时且sign(yQ)为“-”时,ψ=3π/2,A为给定值,A≤10-4;
步骤12、利用距离值r计算得到磁偶极子磁场和磁矩的大小,具体为:
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202011461502.6A CN112611310B (zh) | 2020-12-11 | 2020-12-11 | 一种磁偶极子目标测距测向方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202011461502.6A CN112611310B (zh) | 2020-12-11 | 2020-12-11 | 一种磁偶极子目标测距测向方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN112611310A CN112611310A (zh) | 2021-04-06 |
CN112611310B true CN112611310B (zh) | 2022-09-27 |
Family
ID=75233601
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN202011461502.6A Active CN112611310B (zh) | 2020-12-11 | 2020-12-11 | 一种磁偶极子目标测距测向方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN112611310B (zh) |
Families Citing this family (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN115542225B (zh) * | 2022-10-12 | 2023-05-19 | 中国科学院国家空间科学中心 | 一种提升悬丝扭秤装置磁性测量精度的校正方法 |
CN116804773B (zh) * | 2023-08-21 | 2023-11-07 | 崂山国家实验室 | 基于高阶偏导磁梯度张量的水下磁目标定位方法及装置 |
Citations (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN103630139A (zh) * | 2013-12-17 | 2014-03-12 | 哈尔滨工程大学 | 一种基于地磁梯度张量测量的水下载体全姿态确定方法 |
CN104567871A (zh) * | 2015-01-12 | 2015-04-29 | 哈尔滨工程大学 | 一种基于地磁梯度张量的四元数卡尔曼滤波姿态估计方法 |
Family Cites Families (9)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
FR2560388B1 (fr) * | 1980-03-17 | 1986-08-01 | Commissariat Energie Atomique | Gradientmetre a couche magnetique mince |
US6841994B1 (en) * | 2004-03-01 | 2005-01-11 | The United States Of America As Represented By The Secretary Of The Navy | Magnetic anomaly sensing system for detection, localization and classification of magnetic objects |
US7038458B1 (en) * | 2004-10-12 | 2006-05-02 | The United States Of America As Represented By The Secretary Of The Navy | Magnetic anomaly homing system and method using rotationally invariant scalar contractions of magnetic gradient tensors |
US10338261B2 (en) * | 2015-09-16 | 2019-07-02 | Raytheon Company | Measurement of magnetic field gradients |
RU164969U1 (ru) * | 2016-03-27 | 2016-09-27 | Общество С Ограниченной Ответственностью "Феникс" | Узел датчиков для диагностики технического состояния подземных трубопроводов |
CN109814163B (zh) * | 2019-02-28 | 2020-09-01 | 中国科学院遥感与数字地球研究所 | 一种基于扩展补偿模型的航磁张量数据抑噪方法及*** |
CN110333536B (zh) * | 2019-07-22 | 2021-04-02 | 哈尔滨工程大学 | 一种测距线性定位算法 |
CN112050799B (zh) * | 2020-08-19 | 2022-11-18 | 哈尔滨工程大学 | 一种基于磁梯度张量缩并量之比的测距定位方法 |
CN112050800B (zh) * | 2020-08-19 | 2022-11-18 | 哈尔滨工程大学 | 一种基于日字型三轴磁强计对称配置平面阵列的磁梯度张量定位方法 |
-
2020
- 2020-12-11 CN CN202011461502.6A patent/CN112611310B/zh active Active
Patent Citations (2)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN103630139A (zh) * | 2013-12-17 | 2014-03-12 | 哈尔滨工程大学 | 一种基于地磁梯度张量测量的水下载体全姿态确定方法 |
CN104567871A (zh) * | 2015-01-12 | 2015-04-29 | 哈尔滨工程大学 | 一种基于地磁梯度张量的四元数卡尔曼滤波姿态估计方法 |
Also Published As
Publication number | Publication date |
---|---|
CN112611310A (zh) | 2021-04-06 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
US6841994B1 (en) | Magnetic anomaly sensing system for detection, localization and classification of magnetic objects | |
Storms et al. | Magnetic field navigation in an indoor environment | |
CN112611310B (zh) | 一种磁偶极子目标测距测向方法 | |
US20020171427A1 (en) | Magnetic anomaly sensing system and methods for maneuverable sensing platforms | |
US7603251B1 (en) | Magnetic anomaly sensing system for detection, localization and classification of a magnetic object in a cluttered field of magnetic anomalies | |
US7688072B1 (en) | Portable magnetic sensing system for real-time, point-by-point detection, localization and classification of magnetic objects | |
CN107272069A (zh) | 基于磁异常梯度的磁性目标追踪方法 | |
Springmann et al. | Magnetic sensor calibration and residual dipole characterization for application to nanosatellites | |
US7038458B1 (en) | Magnetic anomaly homing system and method using rotationally invariant scalar contractions of magnetic gradient tensors | |
Lee et al. | Magnetic tensor sensor for gradient-based localization of ferrous object in geomagnetic field | |
Riwanto et al. | Particle swarm optimization with rotation axis fitting for magnetometer calibration | |
Wiegert et al. | Generalized magnetic gradient contraction based method for detection, localization and discrimination of underwater mines and unexploded ordnance | |
Li et al. | Magnetic object positioning method based on tensor spacial invariant relations | |
CN113624240A (zh) | 一种基于磁感应强度与特征矢量的位姿识别方法和装置 | |
Wiegert et al. | Improved magnetic STAR methods for real-time, point-by-point localization of unexploded ordnance and buried mines | |
Lin et al. | Improvement and omnidirectional analysis of magnetic gradient tensor invariants method | |
Wang et al. | Magnetoresistive sensor error compensation method using geometry-constraint contour scaling | |
Liu et al. | Two-step calibration method for three-axis magnetic sensor error based on particle swarm optimization | |
Li et al. | An efficient method for tri-axis magnetometer calibration | |
Wang et al. | Joint use of scalar and vector magnetometers for underwater magnetic target localization | |
Menezes Filho et al. | Accurate and computational-efficient analytical solutions for the extended two-step magnetometer calibration | |
CN114234958A (zh) | 一种基于磁场特征值的磁信标定向方法、存储介质及设备 | |
Včelák et al. | Long-range magnetic tracking system | |
Liu et al. | A hybrid calibration method for a three-axis magnetometer in limited-range attitudes | |
Xu et al. | Simulation Analysis of Magnetic Gradient Full‐Tensor Measurement System |
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 |