CN106052694A - 基于重力矢量及其梯度张量对单个运动物体进行定位跟踪的方法 - Google Patents

Abstract

本发明提供了一种基于重力矢量及其梯度张量对单个运动物体进行定位跟踪的方法，包括以下步骤：(1)在一个重力基站，按一定时间间隔，观测由质量为M的单个运动物体引起的剩余重力矢量及其梯度张量，得到一系列时刻的剩余重力矢量g及其梯度张量T。(2)根据观测数据，利用解析公式，精确计算单个运动目标体的坐标位置、质量和速度矢量。本发明可应用于民用或军事领域，对单个物体运动轨迹和质量的精确实时监测。

Description

基于重力矢量及其梯度张量对单个运动物体进行定位跟踪的 方法

技术领域

本发明涉及一种地球重力学方法，特别涉及一种基于重力矢量及其梯度张量对单个运动物体进行定位跟踪的方法。

背景技术

在地球重力学中，重力梯度张量除用于勘探地球地质结构构造，寻找油气、矿产资源外，还有其他多种应用。例如，利用卫星重力梯度，高精度恢复中短波的地球重力场；利用高精度重力梯度张量，为舰艇或巡航导弹进行导航；2010年，Majid Beiki and LaustB.Pedersen(Majid Beiki and Laust B.Pedersen，Eigenvector analysis of gravitygradient tensor to locate geologic bodies，Geophysics,2010)应用重力梯度张量特征值与特征向量、垂直方向的重力场计算了地下单一异常体的位置与质量，并被用于重力梯度张量剖面数据的定性解释。2014年，Lockerbie(N A Lockerbie，The location ofsubterranean voids using tensor gravity gradiometry，Class.Quantum Grav.2014)基于2个不同测点处的重力梯度张量特征向量唯一确定了单一目标体的位置，并应用于地下掩体范围的大致圈定。

现有技术中，必须基于2个不同测点处的重力梯度张量确定单一目标体位置测量，且只能应用于地下掩体范围的大致圈定。因此，有必要设计一种利用一个基站的剩余重力异常及其梯度张量进行单个运动物体定位和跟踪的解析方法。

发明内容

本发明所解决的技术问题是，针对现有技术的不足，提供了一种基于重力矢量及其梯度张量对单个运动物体进行定位跟踪的方法，利用一个重力观测基站观测重力矢量和重力梯度张量异常，以确定单一运动物体的位置和质量，并跟踪其运动轨迹。

本发明的技术方案为：

一种基于重力矢量及其梯度张量对单个运动物体进行定位跟踪的方法，包括以下步骤：

步骤1、布置重力观测基站：

根据观测目标确定重力观测基站的位置，并在重力观测基站安装重力仪和重力梯度仪；

步骤2、计算质量为M的单个运动物体引起的剩余重力矢量g(t)和剩余重力梯度张量T(t)：

定义空间直角坐标系(x,y,z)；利用重力仪和重力梯度仪按一定时间间隔记录空间直角坐标系(x,y,z)下重力观测基站的重力矢量g_total(t)和重力梯度张量T_total(t)，根据公式(1)计算质量为M的单个运动物体引起的剩余重力矢量g(t)和剩余重力梯度张量T(t)：

$\begin{array}{c}T\left(t\right)=T_{\mathrm{total}}\left(t\right)-T_{\mathrm{tital}}\left(t_0\right)\\ g\left(t\right)=g_{\mathrm{total}}\left(t\right)-g_{\mathrm{total}}\left(t_0\right)\end{array}---\left(1\right)$

上式，T_total(t₀)，g_total(t₀)为单个运动物体对重力观测基站不产生影响时观测到的重力矢量及重力梯度张量，t₀称为时间基点；

步骤3、由剩余重力梯度张量和剩余重力矢量g(t)的垂直分量g_z(t)，根据公式(2)计算单个运动物体在t时刻的质心坐标位置(运动轨迹)和质量：

$\left\{\begin{array}{c}{D}_{xz}\left(t\right)=\frac{T_{xy}\left(t\right)}{T_{yz}\left(t\right)}\\ D_{yz}\left(t\right)=\frac{T_{xy}\left(t\right)}{T_{xz}\left(t\right)}\\ A\left(t\right)={\[D_{xz}^2\left(t\right)+D_{yz}^2\left(t\right)+1\]}^{\frac{1}{2}}\\ d_z\left(t\right)=\frac{A{\left(t\right)}^2-3}{A{\left(t\right)}^2}\frac{g_z\left(t\right)}{T_{zz}\left(t\right)}\\ x^{\′}\left(t\right)=x\left(t\right)-D_{xz}\left(t\right)d_z\left(t\right)\\ y^{\′}\left(t\right)=y\left(t\right)-D_{yz}\left(t\right)d_z\left(t\right)\\ z^{\′}\left(t\right)=z\left(t\right)-d_z\left(t\right)\\ M\left(t\right)=d_z^2\left(t\right)\·g_z\left(t\right)\·A{\left(t\right)}^3\left(\frac{1}{-G}\right)\end{array}\right.,---\left(2\right)$

上式中，T_yz(t)≠0,T_xz(t)≠0，x(t),y(t),z(t)为t时刻重力观测基站在x、y、z轴上的坐标位置，G为引力常数，x′(t),y′(t),z′(t)为t时刻单个运动物体质心在x、y、z轴上的坐标位置，M(t)为t时刻单个运动物体质量。

当T_xy(t)≠0,T_yz(t)≠0时，所述步骤3中，采用公式(3)和公式(4)所示的第一和第二变种公式计算单个运动物体在t时刻的质心坐标位置(运动轨迹)和质量：

$\left\{\begin{array}{c}{D}_{zy}\left(t\right)=\frac{T_{xz}\left(t\right)}{T_{xy}\left(t\right)},D_{xy}\left(t\right)=\frac{T_{xz}\left(t\right)}{T_{yz}\left(t\right)},\\ A\left(t\right)={\[D_{zy}^2\left(t\right)+D_{xy}^2\left(t\right)+1\]}^{\frac{1}{2}},\\ d_y\left(t\right)=\frac{A{\left(t\right)}^2-3D_{zy}^2\left(t\right)}{A{\left(t\right)}^2{D}_{zy}\left(t\right)}\frac{g_z\left(t\right)}{T_{zz}\left(t\right)},\\ x^{\′}\left(t\right)=x\left(t\right)-D_{xy}\left(t\right)d_y\left(t\right),\\ y^{\′}\left(t\right)=y\left(t\right)-d_y\left(t\right),\\ z^{\′}\left(t\right)=z\left(t\right)-D_{zy}\left(t\right)d_y\left(t\right),\\ M\left(t\right)=d_y^2\left(t\right)\·g_z\left(t\right)\·A{\left(t\right)}^3(-\frac{1}{{\mathrm{GD}}_{zy}\left(t\right)})\end{array}\right..---\left(3\right);$

$\left\{\begin{array}{c}{D}_{zy}\left(t\right)=\frac{T_{xz}\left(t\right)}{T_{xy}\left(t\right)},D_{xy}\left(t\right)=\frac{T_{xz}\left(t\right)}{T_{yz}\left(t\right)},\\ A\left(t\right)={\[D_{zy}^2\left(t\right)+D_{xy}^2\left(t\right)+1\]}^{\frac{1}{2}},\\ d_y\left(t\right)=\frac{A{\left(t\right)}^2-3}{A{\left(t\right)}^2}\frac{g_y\left(t\right)}{T_{yy}\left(t\right)},\\ x^{\′}\left(t\right)=x\left(t\right)-D_{xy}\left(t\right)d_y\left(t\right),\\ y^{\′}\left(t\right)=y\left(t\right)-d_y\left(t\right),\\ z^{\′}\left(t\right)=z\left(t\right)-D_{zy}\left(t\right)d_y\left(t\right),\\ M\left(t\right)=|d_y^2\left(t\right)\·g_y\left(t\right)\·A{\left(t\right)}^3\left(\frac{1}{-G}\right)|\end{array}\right.---\left(4\right).$

当T_xy(t)≠0,T_xz(t)≠0时，所述步骤3中，采用公式(5)和公式(6)所示的第三和第四变种公式计算单个运动物体在t时刻的质心坐标位置(运动轨迹)和质量：

$\left\{\begin{array}{c}{D}_{zx}\left(t\right)=\frac{T_{yz}\left(t\right)}{T_{xy}\left(t\right)},D_{yx}\left(t\right)=\frac{T_{yz}\left(t\right)}{T_{xz}\left(t\right)},\\ A\left(t\right)={\[D_{zx}^2\left(t\right)+D_{yx}^2\left(t\right)+1\]}^{\frac{1}{2}},\\ d_x\left(t\right)=\frac{A{\left(t\right)}^2-3D_{zx}^2\left(t\right)}{A{\left(t\right)}^2{D}_{zx}\left(t\right)}\frac{g_z\left(t\right)}{T_{zz}\left(t\right)},\\ x^{\′}\left(t\right)=x\left(t\right)-d_x\left(t\right),\\ y^{\′}\left(t\right)=y\left(t\right)-D_{yx}\left(t\right)d_x\left(t\right),\\ z^{\′}\left(t\right)=z\left(t\right)-D_{zx}\left(t\right)d_x\left(t\right),\\ M\left(t\right)=d_x^2\left(t\right)\·g_z\left(t\right)\·A{\left(t\right)}^3\left(\frac{1}{-{\mathrm{GD}}_{zx}\left(t\right)}\right)\end{array}\right.---\left(5\right);$

$\left\{\begin{array}{c}{D}_{zx}\left(t\right)=\frac{T_{yz}\left(t\right)}{T_{xy}\left(t\right)},D_{yx}\left(t\right)=\frac{T_{yz}\left(t\right)}{T_{xz}\left(t\right)},\\ A\left(t\right)={\[D_{zx}^2\left(t\right)+D_{yx}^2\left(t\right)+1\]}^{\frac{1}{2}},\\ d_x\left(t\right)=\frac{A{\left(t\right)}^2-3}{A{\left(t\right)}^2}\frac{g_x\left(t\right)}{T_{xx}\left(t\right)},\\ x^{\′}\left(t\right)=x\left(t\right)-d_x\left(t\right),\\ y^{\′}\left(t\right)=y\left(t\right)-D_{yx}\left(t\right)d_x\left(t\right),\\ z^{\′}\left(t\right)=z\left(t\right)-D_{zx}\left(t\right)d_x\left(t\right),\\ M\left(t\right)=|d_x^2\left(t\right)\·g_x\left(t\right)\·A{\left(t\right)}^3\left(\frac{1}{-G}\right)|\end{array}\right.---\left(6\right).$

当g_z≠0时，所述步骤3中，采用公式(7)所示的第五变种公式计算单个运动物体在t时刻的质心坐标位置(运动轨迹)和质量：

$\left\{\begin{array}{c}{D}_{xz}\left(t\right)=\frac{g_x\left(t\right)}{g_z\left(t\right)},D_{yz}\left(t\right)=\frac{g_y\left(t\right)}{g_z\left(t\right)},\\ A\left(t\right)={\[D_{xz}^2\left(t\right)+D_{yz}^2\left(t\right)+1\]}^{\frac{1}{2}},\\ d_z\left(t\right)=\frac{A{\left(t\right)}^2-3}{A{\left(t\right)}^2}\frac{g_z\left(t\right)}{T_{zz}\left(t\right)},\\ x^{\′}\left(t\right)=x\left(t\right)-D_{xz}\left(t\right)d_z\left(t\right),\\ y^{\′}\left(t\right)=y\left(t\right)-D_{yz}\left(t\right)d_z\left(t\right),\\ z^{\′}\left(t\right)=z\left(t\right)-d_z\left(t\right),\\ M\left(t\right)=d_z^2\left(t\right)\·g_z\left(t\right)\·A{\left(t\right)}^3\left(\frac{1}{-G}\right)\end{array}\right.---\left(7\right).$

当g_y≠0时，所述步骤3中，当采用公式(8)和公式(9)所示的第六和第七变种公式计算单个运动物体在t时刻的质心坐标位置(运动轨迹)和质量：

$\left\{\begin{array}{c}{D}_{zy}\left(t\right)=\frac{g_z\left(t\right)}{g_y\left(t\right)},D_{xy}\left(t\right)=\frac{g_x\left(t\right)}{g_y\left(t\right)},\\ A\left(t\right)={\[D_{zy}^2\left(t\right)+D_{xy}^2\left(t\right)+1\]}^{\frac{1}{2}},\\ d_y\left(t\right)=\frac{A{\left(t\right)}^2-3D_{zy}^2\left(t\right)}{A{\left(t\right)}^2{D}_{zy}\left(t\right)}\frac{g_z\left(t\right)}{T_{zz}\left(t\right)},\\ x^{\′}\left(t\right)=x\left(t\right)-D_{xy}\left(t\right)d_y\left(t\right),\\ y^{\′}\left(t\right)=y\left(t\right)-d_y\left(t\right),\\ z^{\′}\left(t\right)=z\left(t\right)-D_{zy}\left(t\right)d_y\left(t\right),\\ M\left(t\right)=d_y^2\left(t\right)\·g_z\left(t\right)\·A{\left(t\right)}^3\left(\frac{1}{-{\mathrm{GD}}_{zy}\left(t\right)}\right)\end{array}\right.---\left(8\right);$

$\left\{\begin{array}{c}{D}_{zy}\left(t\right)=\frac{g_z\left(t\right)}{g_y\left(t\right)},D_{xy}\left(t\right)=\frac{g_x\left(t\right)}{g_y\left(t\right)},\\ A\left(t\right)={\[D_{zy}^2\left(t\right)+D_{xy}^2\left(t\right)+1\]}^{\frac{1}{2}},\\ d_y\left(t\right)=\frac{A{\left(t\right)}^2-3}{A{\left(t\right)}^2}\frac{g_y\left(t\right)}{T_{yy}\left(t\right)},\\ x^{\′}\left(t\right)=x\left(t\right)-D_{xy}\left(t\right)d_y\left(t\right),\\ y^{\′}\left(t\right)=y\left(t\right)-d_y\left(t\right),\\ z^{\′}\left(t\right)=z\left(t\right)-D_{zy}\left(t\right)d_y\left(t\right),\\ M\left(t\right)=|d_y^2\left(t\right)\·g_y\left(t\right)\·A{\left(t\right)}^3\left(\frac{1}{-G}\right)|\end{array}\right.---\left(9\right).$

当g_x≠0时，采用公式(10)和公式(11)所示的第八和第九变种公式计算单个运动物体在t时刻的质心坐标位置(运动轨迹)和质量：

$\left\{\begin{array}{c}{D}_{zx}\left(t\right)=\frac{g_z\left(t\right)}{g_x\left(t\right)},D_{yx}\left(t\right)=\frac{g_y\left(t\right)}{g_x\left(t\right)},\\ A\left(t\right)={\[D_{zx}^2\left(t\right)+D_{yx}^2\left(t\right)+1\]}^{\frac{1}{2}},\\ d_x\left(t\right)=\frac{A{\left(t\right)}^2-3D_{zx}^2\left(t\right)}{A{\left(t\right)}^2{D}_{zx}\left(t\right)}\frac{g_z\left(t\right)}{T_{zz}\left(t\right)},\\ x^{\′}\left(t\right)=x\left(t\right)-d_x\left(t\right),\\ y^{\′}\left(t\right)=y\left(t\right)-D_{yx}\left(t\right)d_x\left(t\right),\\ z^{\′}\left(t\right)=z\left(t\right)-D_{zx}\left(t\right)d_x\left(t\right),\\ M\left(t\right)=d_x^2\left(t\right)\·g_z\left(t\right)\·A{\left(t\right)}^3\left(\frac{1}{-{\mathrm{GD}}_{zx}\left(t\right)}\right)\end{array}\right.---\left(10\right);$

$\left\{\begin{array}{c}{D}_{zx}\left(t\right)=\frac{g_z\left(t\right)}{g_x\left(t\right)},D_{yx}\left(t\right)=\frac{g_y\left(t\right)}{g_x\left(t\right)},\\ A\left(t\right)={\[D_{zx}^2\left(t\right)+D_{yx}^2\left(t\right)+1\]}^{\frac{1}{2}},\\ d_x\left(t\right)=\frac{A{\left(t\right)}^2-3}{A{\left(t\right)}^2}\frac{g_x\left(t\right)}{T_{xx}\left(t\right)},\\ x^{\′}\left(t\right)=x\left(t\right)-d_x\left(t\right),\\ y^{\′}\left(t\right)=y\left(t\right)-D_{yx}\left(t\right)d_x\left(t\right),\\ z^{\′}\left(t\right)=z\left(t\right)-D_{zx}\left(t\right)d_x\left(t\right),\\ M\left(t\right)=|d_x^2\left(t\right)\·g_x\left(t\right)\·A{\left(t\right)}^3\left(\frac{1}{-G}\right)|\end{array}\right.---\left(11\right).$

由单个运动物体的质心坐标位置(运动轨迹)，根据公式(12)计算单个运动物体相对于重力观测基站的运动速度V＝(V_x,V_y,V_z)：

$\begin{array}{c}{V}_x=\frac{x^{\′}\left(t_{i+1}\right)-x^{\′}\left(t_i\right)}{t_{i+1}-t_i}\\ V_y=\frac{y^{\′}\left(t_{i+1}\right)-y^{\′}\left(t_i\right)}{t_{i+1}-t_i}\\ V_z=\frac{z^{\′}\left(t_{i+1}\right)-z^{\′}\left(t_i\right)}{t_{i+1}-t_i}\end{array},---\left(12\right)$

其中，t_i+1第i+1个时刻，t_i为第i个时刻。

所述单个运动物体位于太空、地表或水面以下、地表或水面以上，或位于地表面或水面。

所述重力观测基站固定在地表面或地表以上某个位置，或位于船只、飞机或卫星等移动平台上。

设重力观测基站的移动速度为V_s，则单个运动物体的绝对速度为：

V_绝对速度＝V_s+V， (13)。

本发明原理为：

假设质量异常体为V，剩余密度为ρ，观察点(x,y,z)的剩余重力梯度张量满足下列积分方程：

$\begin{array}{c}{T}_{xx}=-G\underset{V}{\∫\∫\∫}\mathrm{\ρ}\frac{r^2-3{(x-x^{\′})}^2}{r^5}dv\\ T_{xy}=-G\underset{V}{\∫\∫\∫}\mathrm{\ρ}\frac{-3(x-x^{\′})(y-y^{\′})}{r^5}dv\\ T_{xz}=-G\underset{V}{\∫\∫\∫}\mathrm{\ρ}\frac{-3(x-x^{\′})(z-z^{\′})}{r^5}dv\\ T_{yy}=-G\underset{V}{\∫\∫\∫}\mathrm{\ρ}\frac{r^2-3{(y-y^{\′})}^2}{r^5}dv\\ T_{yz}=-G\underset{V}{\∫\∫\∫}\mathrm{\ρ}\frac{-3(y-y^{\′})(z-z^{\′})}{r^5}dv\\ T_{zz}=-G\underset{V}{\∫\∫\∫}\mathrm{\ρ}\frac{r^2-3{(z-z^{\′})}^2}{r^5}dv\end{array}---\left(14\right)$

其中T_xx+T_yy+T_zz＝0，(x′,y′,z′)∈V，r＝|(x,y,z)-(x′,y′,z′)|，上式中只存在5个独立的张量分量。当物体与观测点有一定距离时，复杂的质量异常体可以用规则的球体来代替。当观测点位于球体的外部时，式(14)中的重力张量异常表示式可进一步简化为：

$\begin{array}{c}{T}_{xx}=-GM\frac{r^2-3{(x-x^{\′})}^2}{r^2}\\ T_{xy}=-GM\frac{-3(x-x^{\′})(y-y^{\′})}{r^5}\\ T_{xz}=-GM\frac{-3(x-x^{\′})(z-z^{\′})}{r^5}\\ T_{yz}=-GM\frac{-3(y-y^{\′})(z-z^{\′})}{r^5}\\ T_{zz}=-GM\frac{r^2-3{(z-z^{\′})}^2}{r^5}\end{array}---\left(15\right)$

同样地，在观测点，可以利用重力仪测得剩余重力矢量g的水平分量g_x、g_y，剩余重力矢量g的垂直分量g_z：

$\begin{array}{c}{g}_x=-GM\frac{x-x^{\′}}{r^3}\\ g_y=-GM\frac{y-y^{\′}}{r^3}\\ g_z=-GM\frac{z-z^{\′}}{r^3}\end{array}---\left(16\right)$

以T_xy为参考分量，可得：

$\begin{array}{c}{D}_{xz}=\frac{T_{xy}}{T_{yz}}=-\frac{x-x^{\′}}{z-z^{\′}}=\frac{d_x}{d_z}\\ D_{yz}=\frac{T_{xy}}{T_{xz}}=-\frac{y-y^{\′}}{z-z^{\′}}=\frac{d_y}{d_z}\end{array}---\left(17\right)$

然后，剩余重力矢量g的垂直分量g_z可为：

$g_z=-{\mathrm{GMA}}_z=-GM\frac{d_z}{{\[D_{xz}^2+D_{yz}^2+1\]}^{\frac{3}{2}}{d}_z^3}=-GM\frac{1}{A^3{d}_z^2},---\left(18\right)$

其中，定义对于T_zz来说，可改写为：

$T_{zz}=-{\mathrm{GMA}}_{zz}=-GM\frac{d_x^2+d_y^2+d_z^2-3d_z^2}{{(d_x^2+d_y^2+d_z^2)}^{\frac{5}{2}}}=-GM\frac{(A^2-3)}{A^5{d}_z^3},---\left(19\right)$

利用T_zz与g_z的比值，可以精确获得d_z＝z-z′，

$d_z=\frac{A^2-3}{A^2}\frac{g_z}{T_{zz}}---\left(20\right)$

进一步可以获得d_x,d_y，

$\begin{array}{c}{d}_x=D_{xz}{d}_z\\ d_y=D_{yz}{d}_z\end{array},---\left(21\right)$

和物体质量M，

$M=d_z^2\·g_z\·A^3\left(\frac{1}{-G}\right).---\left(22\right)$

物体质心坐标与质量计算公式归类为：

$\begin{array}{c}{D}_{xz}=\frac{T_{xy}}{T_{yz}}\\ D_{yz}=\frac{T_{xy}}{T_{xz}}\\ A={\[D_{xz}^2+D_{yz}^2+1\]}^{\frac{1}{2}}\\ d_z=\frac{A^2-3}{A^2}\frac{g_z}{T_{zz}}\\ x^{\′}=x-D_{xz}{d}_z\\ y^{\′}=y-D_{yz}{d}_z\\ z^{\′}=z-d_z\\ M=d_z^2\·g_z\·A^3\left(\frac{1}{-G}\right)\end{array}---\left(23\right)$

其中，x,y,z为观测点坐标位置，G为引力常数，x′,y′,z′为物体质心的坐标位置，M为物体质量。公式(23)即为上述定位公式(2)。

重新定义公式(17)如下：

$\begin{array}{c}{D}_{zy}=\frac{T_{xz}}{T_{xy}}=\frac{z-z^{\′}}{y-y^{\′}}=\frac{d_z}{d_y}\\ D_{xy}=\frac{T_{xz}}{T_{yz}}=\frac{x-x^{\′}}{y-y^{\′}}=\frac{d_x}{d_y}\end{array}---\left(24\right)$

并结合公式(18)和公式(19)，即可获得上述定位公式(3)和公式(4)。

重新定义公式(17)如下：

$\begin{array}{c}{D}_{zx}=\frac{T_{yz}}{T_{xy}}=\frac{z-z^{\′}}{x-x^{\′}}=\frac{d_z}{d_x}\\ D_{yx}=\frac{T_{yz}}{T_{xz}}=\frac{y-y^{\′}}{x-x^{\′}}=\frac{d_y}{d_x}\end{array}---\left(25\right)$

并结合公式(18)和公式(19)，即可获得上述定位公式(5)和公式(6)。

重新定义公式(17)如下：

$\begin{array}{c}{D}_{xz}=\frac{g_x}{g_z}=\frac{x-x^{\′}}{z-z^{\′}}=\frac{d_x}{d_z}\\ D_{yz}=\frac{g_y}{g_z}=\frac{y-y^{\′}}{z-z^{\′}}=\frac{d_y}{d_z}\end{array}---\left(26\right)$

并结合公式(18)和公式(19)，即可获得上述定位公式(7)。

重新定义公式(17)如下：

$\begin{array}{c}{D}_{zy}=\frac{g_z}{g_y}=\frac{z-z^{\′}}{y-y^{\′}}=\frac{d_z}{d_y}\\ D_{xy}=\frac{g_x}{g_y}=\frac{x-x^{\′}}{y-y^{\′}}=\frac{d_x}{d_y}\end{array}---\left(27\right)$

并结合公式(18)和公式(19)，即可获得上述定位公式(8)和公式(9)。

重新定义公式(17)如下：

$\begin{array}{c}{D}_{zx}=\frac{g_z}{g_x}=\frac{z-z^{\′}}{x-x^{\′}}=\frac{d_z}{d_x}\\ D_{yx}=\frac{g_y}{g_x}=\frac{y-y^{\′}}{x-x^{\′}}=\frac{d_y}{d_x}\end{array}---\left(28\right)$

并结合公式(18)和公式(19)，即可获得上述定位公式(10)和公式(11)。

有益效果：

本方法给出了利用一个重力基站的重力矢量及其梯度张量监测一个运动物体的坐标位置和质量的核心算法，可以实时地监测运动物体的位置和质量，精度高。通过观测剩余重力异常及其梯度张量随时间的变化，即可求出物体的运动轨迹。运动物体可以位于外太空、地表面及以下，如地下洞穴及其中运动的汽车、人员，常规潜艇及核潜艇，航空母舰及各种水面舰艇等；也可以位于空气中，如飞机、导弹等。本发明可应用于民用或军事领域，对单个运动物体的轨迹和质量进行精确实时监测。

附图说明

图1为本发明利用一个重力观测基站进行单体运动轨迹识别的示意图。

图2为利用本发明进行外太空未知运动质量体(如未知彗星)的轨迹识别装置图。在卫星上布置重力梯度张量仪与重力仪，实时监测重力矢量及其梯度张量随时间的变化，通过本发明的公式(1)-(13)计算出未知运动质量体的坐标和质量，然后通过通信系统返回跟踪结果给地表监测站。

图3为利用本发明进行地下未知运动质量体的轨迹识别装置图。在地表上布置重力梯度张量仪与重力仪，实时监测重力矢量及其梯度张量随时间的变化，通过本发明的公式(1)-(13)计算出未知运动质量体的坐标和质量，然后通过通信系统返回给地表监测站。

图4为利用本发明进行水下未知运动质量体的轨迹识别装置图。在海面船舰上布置重力梯度张量仪与重力仪，实时监测重力矢量及其梯度张量随时间的变化，通过本发明的公式(1)-(13)计算出未知运动质量体的坐标和质量，然后通过通信系统返回给船舰监测站。

图5为实施例1中运动质量体的真实运动轨迹与基站的位置。

图6为实施例1中基站记录的50个时刻的剩余重力梯度张量。

图7是实施例1中运动质量体轨迹的监测结果。

图8是实施例1中运动质量体轨迹的监测结果的相对误差图。

具体实施方式

以下结合附图和具体实施方式对本发明作进一步的说明。

本发明公开了一种基于重力矢量及其梯度张量对单个运动物体进行定位跟踪的方法，包括以下步骤：

(1)重力观测基站布置设计：确定观测目标及仪器灵敏度范围，确定重力观测基站的合理布置；

(2)剩余重力矢量及其梯度张量的计算：实时记录重力观测基站的重力矢量与重力梯度张量，以便实时计算出重力矢量与重力梯度张量的异常变化，即剩余重力矢量g及其梯度张量T。

(3)实时坐标与质量的计算与监测：整理基站采集数据，利用公式(1)-(13)所示的计算方法，实时确定单一运动目标体的坐标、质量和运动速度等。

实施例1：

本实施例为进行水下运动质量体轨迹定位的实例。运动质量体的质量为4×10⁶kg，运动速度为40节(74080米每小时)，运动轨迹为的抛物型方程，y＝0.1875x²，z＝100，50个时间间隔，时间间隔为5分钟，Δt＝1,2,...,50。图5展示了运动质量体的运动轨迹与重力观测基站的位置。运动质量体起始点x坐标为-4000米，重力观察基站位于(0，-200米，0)处(用符号“+”标识)。

利用图4所示的装置进行水下未知运动质量体的轨迹识别在海面船舰上布置重力梯度张量仪与重力仪，实时监测重力矢量及其梯度张量随时间的变化，共记录50个时刻的剩余重力梯度张量及剩余重力矢量的垂直分量。图6为基站记录的50个时刻的剩余重力梯度张量，显然地，从剩余重力梯度张量随时间的变化曲线上无法鉴别出运动质量体的位置和轨迹。

通过本发明的公式(1)-(13)计算出未知运动质量体的坐标和质量，然后通过通信系统返回给船舰监测站。得到如图7所示的运动质量体的轨迹的监测结果，结果表明识别出的轨迹、质量与真实运动轨迹、质量完全重合。图8是运动质量体运动轨迹的监测结果的相对误差图，结果表明识别出的轨迹、质量的相对误差小于0.0001％,均可以忽略不计，证明了本发明方法的有效性。

Claims (10)

1.一种基于重力矢量及其梯度张量对单个运动物体进行定位跟踪的方法，其特征在于，包括以下步骤：

步骤1、布置重力观测基站：

根据观测目标确定重力观测基站的位置，并在重力观测基站安装重力仪和重力梯度仪；

步骤2、计算质量为M的单个运动物体引起的剩余重力矢量g(t)和剩余重力梯度张量T(t)：

定义空间直角坐标系(x,y,z)；利用重力仪和重力梯度仪按一定时间间隔记录空间直角坐标系(x,y,z)下重力观测基站的重力矢量g_total(t)和重力梯度张量T_total(t)，根据公式(1)计算质量为M的单个运动物体引起的剩余重力矢量g(t)和剩余重力梯度张量T(t)：

$\begin{array}{c}T\left(t\right)=T_{total}\left(t\right)-T_{total}\left(t_0\right)\\ g\left(t\right)=g_{total}\left(t\right)-g_{total}\left(t_0\right)\end{array}---\left(1\right)$

上式，T_total(t₀)，g_total(t₀)为单个运动物体对重力观测基站不产生影响时观测到的重力矢量及重力梯度张量，t₀称为时间基点；

步骤3、由剩余重力梯度张量和剩余重力矢量g(t)的垂直分量g_z(t)，根据公式(2)计算单个运动物体在t时刻的质心坐标位置(运动轨迹)和质量：

$\{\begin{array}{c}{D}_{xz}\left(t\right)=\frac{T_{xy}\left(t\right)}{T_{yz}\left(t\right)}\\ D_{yz}\left(t\right)=\frac{T_{xy}\left(t\right)}{T_{xz}\left(t\right)}\\ A\left(t\right)={\[D_{xz}^2\left(t\right)+D_{yz}^2\left(t\right)+1\]}^{\frac{1}{2}}\\ d_z\left(t\right)=\frac{A{\left(t\right)}^2-3}{A{\left(t\right)}^2}\frac{g_z\left(t\right)}{T_{zz}\left(t\right)}\\ x^{\′}\left(t\right)=x\left(t\right)-D_{xz}\left(t\right)d_z\left(t\right)\\ y^{\′}\left(t\right)=y\left(t\right)-D_{yz}\left(t\right)d_z\left(t\right)\\ z^{\′}\left(t\right)=z\left(t\right)-d_z\left(t\right)\\ M\left(t\right)=d_z^2\left(t\right)\·g_z\left(t\right)\·A{\left(t\right)}^3\left(\frac{1}{-G}\right)\end{array},---\left(2\right)$

上式中，T_yz(t)≠0,T_xz(t)≠0，x(t),y(t),z(t)为t时刻重力观测基站在x、y、z轴上的坐标位置，G为引力常数，x′(t),y′(t),z′(t)为t时刻单个运动物体质心在x、y、z轴上的坐标位置，M(t)为t时刻单个运动物体质量。

2.权利要求1所述的基于重力矢量及其梯度张量对单个运动物体进行定位跟踪的方法，其特征在于：当T_xy(t)≠0,T_yz(t)≠0时，所述步骤3中，采用公式(3)和公式(4)所示的第一和第二变种公式计算单个运动物体在t时刻的质心坐标位置和质量：

$\left\{\begin{array}{c}{D}_{zy}\left(t\right)=\frac{T_{xz}\left(t\right)}{T_{xy}\left(t\right)},D_{xy}\left(t\right)=\frac{T_{xz}\left(t\right)}{T_{yz}\left(t\right)},\\ A\left(t\right)={\[D_{zy}^2\left(t\right)+D_{xy}^2\left(t\right)+1\]}^{\frac{1}{2}},\\ d_y\left(t\right)=\frac{A{\left(t\right)}^2-3D_{zy}^2\left(t\right)}{A{\left(t\right)}^2{D}_{zy}\left(t\right)}\frac{g_z\left(t\right)}{T_{zz}\left(t\right)},\\ x^{\′}\left(t\right)=x\left(t\right)-D_{xy}\left(t\right)d_y\left(t\right),\\ y^{\′}\left(t\right)=y\left(t\right)-d_y\left(t\right),\\ z^{\′}\left(t\right)=z\left(t\right)-D_{zy}\left(t\right)d_y\left(t\right),\\ M\left(t\right)=d_y^2\left(t\right)\·g_z\left(t\right)\·A{\left(t\right)}^3\left(\frac{1}{-{\mathrm{GD}}_{zy}\left(t\right)}\right)\end{array}\right..---\left(3\right);$

$\left\{\begin{array}{c}{D}_{zy}\left(t\right)=\frac{T_{xz}\left(t\right)}{T_{xy}\left(t\right)},D_{xy}\left(t\right)=\frac{T_{xz}\left(t\right)}{T_{yz}\left(t\right)},\\ A\left(t\right)={\[D_{zy}^2\left(t\right)+D_{xy}^2\left(t\right)+1\]}^{\frac{1}{2}},\\ d_y\left(t\right)=\frac{A{\left(t\right)}^2-3}{A{\left(t\right)}^2}\frac{g_y\left(t\right)}{T_{yy}\left(t\right)},\\ x^{\′}\left(t\right)=x\left(t\right)-D_{xy}\left(t\right)d_y\left(t\right),\\ y^{\′}\left(t\right)=y\left(t\right)-d_y\left(t\right),\\ z^{\′}\left(t\right)=z\left(t\right)-D_{zy}\left(t\right)d_y\left(t\right),\\ M\left(t\right)=|d_y^2\left(t\right)\·g_y\left(t\right)\·A{\left(t\right)}^3\left(\frac{1}{-G}\right)|\end{array}\right.---\left(4\right).$

3.根据权利要求1所述的基于重力矢量及其梯度张量对单个运动物体进行定位跟踪的方法，其特征在于：当T_xy(t)≠0,T_xz(t)≠0时，所述步骤3中，采用公式(5)和公式(6)所示的第三和第四变种公式计算单个运动物体在t时刻的质心坐标位置和质量：

$\left\{\begin{array}{c}{D}_{zx}\left(t\right)=\frac{T_{yz}\left(t\right)}{T_{xy}\left(t\right)},D_{yx}\left(t\right)=\frac{T_{yz}\left(t\right)}{T_{xz}\left(t\right)},\\ A\left(t\right)={\[D_{zx}^2\left(t\right)+D_{yx}^2\left(t\right)+1\]}^{\frac{1}{2}},\\ d_x\left(t\right)=\frac{A{\left(t\right)}^2-3D_{zx}^2\left(t\right)}{A{\left(t\right)}^2{D}_{zx}\left(t\right)}\frac{g_z\left(t\right)}{T_{zz}\left(t\right)},\\ x^{\′}\left(t\right)=x\left(t\right)-d_x\left(t\right),\\ y^{\′}\left(t\right)=y\left(t\right)-D_{yx}\left(t\right)d_x\left(t\right),\\ z^{\′}\left(t\right)=z\left(t\right)-D_{zx}\left(t\right)d_x\left(t\right),\\ M\left(t\right)=d_x^2\left(t\right)\·g_z\left(t\right)\·A{\left(t\right)}^3\left(\frac{1}{-{\mathrm{GD}}_{zx}\left(t\right)}\right)\end{array}\right..---\left(5\right);$

$\left\{\begin{array}{c}{D}_{zx}\left(t\right)=\frac{T_{yz}\left(t\right)}{T_{xy}\left(t\right)},D_{yx}\left(t\right)=\frac{T_{yz}\left(t\right)}{T_{xz}\left(t\right)},\\ A\left(t\right)={\[D_{zx}^2\left(t\right)+D_{yx}^2\left(t\right)+1\]}^{\frac{1}{2}},\\ d_x\left(t\right)=\frac{A{\left(t\right)}^2-3}{A{\left(t\right)}^2}\frac{g_x\left(t\right)}{T_{xx}\left(t\right)},\\ x^{\′}\left(t\right)=x\left(t\right)-d_x\left(t\right),\\ y^{\′}\left(t\right)=y\left(t\right)-D_{yx}\left(t\right)d_x\left(t\right),\\ z^{\′}\left(t\right)=z\left(t\right)-D_{zx}\left(t\right)d_x\left(t\right),\\ M\left(t\right)=|d_x^2\left(t\right)\·g_x\left(t\right)\·A{\left(t\right)}^3\left(\frac{1}{-G}\right)|\end{array}\right.---\left(6\right).$

4.根据权利要求1所述的基于重力矢量及其梯度张量对单个运动物体进行定位跟踪的方法，其特征在于：当g_z≠0时，所述步骤3中，采用公式(7)所示的第五变种公式计算单个运动物体在t时刻的质心坐标位置和质量：

$\left\{\begin{array}{c}{D}_{xz}\left(t\right)=\frac{g_x\left(t\right)}{g_z\left(t\right)},D_{yz}\left(t\right)=\frac{g_y\left(t\right)}{g_z\left(t\right)},\\ A\left(t\right)={\[D_{xz}^2\left(t\right)+D_{yz}^2\left(t\right)+1\]}^{\frac{1}{2}},\\ d_z\left(t\right)=\frac{A{\left(t\right)}^2-3}{A{\left(t\right)}^2}\frac{g_z\left(t\right)}{T_{zz}\left(t\right)},\\ x^{\′}\left(t\right)=x\left(t\right)-D_{xz}\left(t\right)d_z\left(t\right),\\ y^{\′}\left(t\right)=y\left(t\right)-D_{yz}\left(t\right)d_z\left(t\right),\\ z^{\′}\left(t\right)=z\left(t\right)-d_z\left(t\right),\\ M\left(t\right)=d_z^2\left(t\right)\·g_z\left(t\right)\·A{\left(t\right)}^3\left(\frac{1}{-G}\right)\end{array}\right.---\left(7\right).$

5.根据权利要求1所述的基于重力矢量及其梯度张量对单个运动物体进行定位跟踪的方法，其特征在于：当g_y≠0时，所述步骤3中，当采用公式(8)和公式(9)所示的第六和第七变种公式计算单个运动物体在t时刻的质心坐标位置和质量：

$\left\{\begin{array}{c}{D}_{zy}\left(t\right)=\frac{g_z\left(t\right)}{g_y\left(t\right)},D_{xy}\left(t\right)=\frac{g_x\left(t\right)}{g_y\left(t\right)},\\ A\left(t\right)={\[D_{zy}^2\left(t\right)+D_{xy}^2\left(t\right)+1\]}^{\frac{1}{2}},\\ d_y\left(t\right)=\frac{A{\left(t\right)}^2-3D_{zy}^2\left(t\right)}{A{\left(t\right)}^2{D}_{zy}\left(t\right)}\frac{g_z\left(t\right)}{T_{zz}\left(t\right)},\\ x^{\′}\left(t\right)=x\left(t\right)-D_{xy}\left(t\right)d_y\left(t\right),\\ y^{\′}\left(t\right)=y\left(t\right)-d_y\left(t\right),\\ z^{\′}\left(t\right)=z\left(t\right)-D_{zy}\left(t\right)d_y\left(t\right),\\ M\left(t\right)=d_y^2\left(t\right)\·g_z\left(t\right)\·A{\left(t\right)}^3\left(\frac{1}{-{\mathrm{GD}}_{zy}\left(t\right)}\right)\end{array}\right.---\left(8\right);$

$\left\{\begin{array}{c}{D}_{zy}\left(t\right)=\frac{g_z\left(t\right)}{g_y\left(t\right)},D_{xy}\left(t\right)=\frac{g_x\left(t\right)}{g_y\left(t\right)},\\ A\left(t\right)={\[D_{zy}^2\left(t\right)+D_{xy}^2\left(t\right)+1\]}^{\frac{1}{2}},\\ d_y\left(t\right)=\frac{A{\left(t\right)}^2-3}{A{\left(t\right)}^2}\frac{g_y\left(t\right)}{T_{yy}\left(t\right)},\\ x^{\′}\left(t\right)=x\left(t\right)-D_{xy}\left(t\right)d_y\left(t\right),\\ y^{\′}\left(t\right)=y\left(t\right)-d_y\left(t\right),\\ z^{\′}\left(t\right)=z\left(t\right)-D_{zy}\left(t\right)d_y\left(t\right),\\ M\left(t\right)=|d_y^2\left(t\right)\·g_y\left(t\right)\·A{\left(t\right)}^3\left(\frac{1}{-G}\right)|\end{array}\right.---\left(9\right).$

6.权利要求1所述的基于重力矢量及其梯度张量对单个运动物体进行定位跟踪的方法，其特征在于：当g_x≠0时，采用公式(10)和公式(11)所示的第八和第九变种公式计算单个运动物体在t时刻的质心坐标位置和质量：

$\left\{\begin{array}{c}{D}_{zx}\left(t\right)=\frac{g_z\left(t\right)}{g_x\left(t\right)},D_{yx}\left(t\right)=\frac{g_y\left(t\right)}{g_x\left(t\right)},\\ A\left(t\right)={\[D_{zx}^2\left(t\right)+D_{yx}^2\left(t\right)+1\]}^{\frac{1}{2}},\\ d_x\left(t\right)=\frac{A{\left(t\right)}^2-3D_{zx}^2\left(t\right)}{A{\left(t\right)}^2{D}_{zx}\left(t\right)}\frac{g_z\left(t\right)}{T_{zz}\left(t\right)},\\ x^{\′}\left(t\right)=x\left(t\right)-d_x\left(t\right),\\ y^{\′}\left(t\right)=y\left(t\right)-D_{yx}\left(t\right)d_x\left(t\right),\\ z^{\′}\left(t\right)=z\left(t\right)-D_{zx}\left(t\right)d_x\left(t\right),\\ M\left(t\right)=d_x^2\left(t\right)\·g_z\left(t\right)\·A{\left(t\right)}^3\left(\frac{1}{-{\mathrm{GD}}_{zx}\left(t\right)}\right)\end{array}\right.---\left(10\right);$

$\left\{\begin{array}{c}{D}_{zx}\left(t\right)=\frac{g_z\left(t\right)}{g_x\left(t\right)},D_{yx}\left(t\right)=\frac{g_y\left(t\right)}{g_x\left(t\right)},\\ A\left(t\right)={\[D_{zx}^2\left(t\right)+D_{yx}^2\left(t\right)+1\]}^{\frac{1}{2}},\\ d_x\left(t\right)=\frac{A{\left(t\right)}^2-3}{A{\left(t\right)}^2}\frac{g_x\left(t\right)}{T_{xx}\left(t\right)},\\ x^{\′}\left(t\right)=x\left(t\right)-d_x\left(t\right),\\ y^{\′}\left(t\right)=y\left(t\right)-D_{yx}\left(t\right)d_x\left(t\right),\\ z^{\′}\left(t\right)=z\left(t\right)-D_{zx}\left(t\right)d_x\left(t\right),\\ M\left(t\right)=|d_x^2\left(t\right)\·g_x\left(t\right)\·A{\left(t\right)}^3\left(\frac{1}{-G}\right)|\end{array}\right.---\left(11\right).$

7.权利要求1～6所述的基于重力矢量及其梯度张量对单个运动物体进行定位跟踪的方法，其特征在于：由单个运动物体的质心坐标位置，根据公式(12)计算单个运动物体相对于重力观测基站的运动速度V＝(V_x,V_y,V_z)：

$\begin{array}{c}{V}_x=\frac{x^{\′}\left(t_{i+1}\right)-x^{\′}\left(t_i\right)}{t_{i+1}-t_i}\\ V_y=\frac{y^{\′}\left(t_{i+1}\right)-y^{\′}\left(t_i\right)}{t_{i+1}-t_i}\\ V_z=\frac{z^{\′}\left(t_{i+1}\right)-z^{\′}\left(t_i\right)}{t_{i+1}-t_i}\end{array}---\left(12\right)$

其中，t_i+1第i+1个时刻，t_i为第i个时刻。

8.权利要求1～6所述的基于重力矢量及其梯度张量对单个运动物体进行定位跟踪的方法，其特征在于：单个运动物体位于太空、地表或水面以下、地表或水面以上，或位于地表面或水面。

9.权利要求1～6所述的基于重力矢量及其梯度张量对单个运动物体进行定位跟踪的方法，其特征在于：重力观测基站固定在地表面或地表以上某个位置，或位于移动平台上。

10.权利要求7所述的基于重力矢量及其梯度张量对单个运动物体进行定位跟踪的方法，其特征在于：设重力观测基站的移动速度为V_s，则单个运动物体的绝对速度为：

V_绝对速度＝V_s+V， (13)。