TDOA-GDOP计算

TDOA-GDOP计算

背景

TDOA定位是一种利用时间差进行定位的方法。通过测量信号到达监测站的时间，可以确定信号源的距离。利用信号源到各个监测站的距离(以监测站为中心，距离为半径作圆)，就能确定信号的位置。

一个目标点：$E(x,y,z)$

四个测向站：$S_0(x_0,y_0,z_0),S_1(x_1,y_1,z_1),S_2(x_2,y_2,z_2),S_3(x_3,y_3,z_3),$

（设定$S_0$站点为主站）

基本参数计算

各测向站与目标点的距离：
$$r_i=\sqrt{(x-x_i{)}^2+(y-y_i{)}^2+(z-z_i{)}^2}i=0,1,2,3$$
各测向站与主站的距离：
$$△r_i=r_i-r_0=c\ast △t_{i}i=1,2,3$$
($△t_i$表示在接收信号过程中第$i$个测向站与主站的接收的时间差，$c$为光速)

误差分析

假设目标定位误差为：$(dx,dy,dz)$

各站点的测量误差为：$(d{x}_i,d{y}_i,d{z}_i)i=0,1,2,3$

(上述各分量之间不相关)

根据误差传递原理：$△r_i=r_i-r_{0}i=1,2,3$

等式两边同时微分：
$$d△r_i=(F_{ix}-F_{0x})dx+(F_{iy}-F_{0y})dy+(F_{iz}-F_{0z})dz+(k_0-k_i)i=1,2,3$$

$$\{\begin{array}{l}{F}_{ix}=\frac{\partial r_i}{\partial x}=-\frac{\partial r_i}{\partial x_i}=\frac{x-x_i}{r_i}\\ F_{iy}=\frac{\partial r_i}{\partial y}=-\frac{\partial r_i}{\partial y_i}=\frac{y-y_i}{r_i}\\ F_{iz}=\frac{\partial r_i}{\partial z}=-\frac{\partial r_i}{\partial z_i}=\frac{z-z_i}{r_i}\\ k_i=F_{ix}d{x}_i+F_{iy}d{y}_i+F_{iz}d{z}_i\end{array}i=1,2,3$$

故$\mathbf{F}$矩阵如下：
$$\mathbf{F}=\left[\begin{array}{ccc}\frac{x-x_1}{r_1}-\frac{x-x_0}{r_0}& \frac{y-y_1}{r_1}-\frac{y-y_0}{r_0}& \frac{z-z_1}{r_1}-\frac{z-z_0}{r_0}\\ \frac{x-x_2}{r_2}-\frac{x-x_0}{r_0}& \frac{y-y_2}{r_2}-\frac{y-y_0}{r_0}& \frac{z-z_2}{r_2}-\frac{z-z_0}{r_0}\\ \frac{x-x_3}{r_3}-\frac{x-x_0}{r_0}& \frac{y-y_3}{r_3}-\frac{y-y_0}{r_0}& \frac{z-z_3}{r_3}-\frac{z-z_0}{r_0}\end{array}\right]$$
将同时微分的结果表示为矩阵形式：
$$\mathbf{d}△\mathbf{r}=\mathbf{F}\cdot \mathbf{dr}+\mathbf{dS}$$
$\mathbf{d}△\mathbf{r}$表示TDOA估计引入的误差，同时联立（2）式：
$$\mathbf{d}△\mathbf{r}=\left[\begin{array}{c}d△r_1\\ d△r_2\\ d△r_3\end{array}\right]=\left[\begin{array}{c}c\ast d(t_1-t_0)\\ c\ast d(t_2-t_0)\\ c\ast d(t_3-t_0)\end{array}\right]$$
$\mathbf{dr}$表示待求目标点的定位误差：
$$\mathbf{dr}=\left[\begin{array}{c}dx\\ dy\\ dz\end{array}\right]$$
$$\mathbf{dS}=\left[\begin{array}{c}{k}_0-k_1\\ k_0-k_2\\ k_0-k_3\end{array}\right]$$

进一步求解

求解目标：目标点定位误差$\mathbf{dr}$

利用伪逆法：
$$\mathbf{dr}=(\mathbf{F}{\mathbf{F}}^T{)}^{-1}{\mathbf{F}}^T(\mathbf{d}△\mathbf{r}-\mathbf{dS})$$
​ 由此可知，目标辐射源的定位误差与站址布局方式、接收站站址误差以及时差估计误差有关。由TDOA的公式可知：各TDOA估计都和主站相关，也就是说TDOA估计中含有相同的误差信息，因此可得出结论：测量误差在$△r_i$处的值是彼此相关的。现假设修正后的测量误差为零均值，且站址误差为恒定值，而其各误差矢量之间以及矢量各元素之间均不相关。在该情况下，定位误差的协方差矩阵可写成如下形式：
$${\mathbf{P}}_{\mathbf{dr}}=E[\mathbf{dr}\cdot {\mathbf{dr}}^T]$$
令$\mathbf{C}=(\mathbf{F}{\mathbf{F}}^T{)}^{-1}{\mathbf{F}}^T=(c_{ij}{)}_{3\ast 3}$，因此公式可化为：
P d r = E [ d r ⋅ d r T ] = C { E [ d △ r ⋅ d △ r T ] + E [ d S ⋅ d S T ] } C T \boldsymbol{P_{dr}} = E[\boldsymbol{dr}·{\boldsymbol{dr}}^T] = \boldsymbol{C}\{ E[\boldsymbol{d\triangle r} · \boldsymbol{{d\triangle r}}^T]+E[\boldsymbol{dS}·{\boldsymbol{dS}}^T]\}\boldsymbol{C}^T Pdr​=E[dr⋅drT]=C{E[d△r⋅d△rT]+E[dS⋅dST]}CT
其中：
$$E[\mathbf{d}△\mathbf{r}\cdot {\mathbf{d}△\mathbf{r}}^T]=\left[\begin{array}{ccc}{\sigma }_{△r_1}^2& {\eta }_{12}{\sigma }_{△r_1△r_2}& {\eta }_{13}{\sigma }_{△r_1△r_3}\\ {\eta }_{12}{\sigma }_{△r_1△r_2}& {\sigma }_{△r_2}^2& {\eta }_{23}{\sigma }_{△r_2△r_3}\\ {\eta }_{13}{\sigma }_{△r_1△r_3}& {\eta }_{23}{\sigma }_{△r_2△r_3}& {\sigma }_{△r_3}^2\end{array}\right]$$
本质上而言，${\sigma }_{△r_i}$是测量误差$△r_i$的标准差，参考（2）实质上就是对时间测量的误差，即：
$$E[\mathbf{d}△\mathbf{r}\cdot {\mathbf{d}△\mathbf{r}}^T]=c^2\ast \left[\begin{array}{ccc}{\sigma }_{△t_1}^2& {\eta }_{12}{\sigma }_{△t_1△t_2}& {\eta }_{13}{\sigma }_{△t_1△t_3}\\ {\eta }_{12}{\sigma }_{△t_1△t_2}& {\sigma }_{△t_2}^2& {\eta }_{23}{\sigma }_{△t_2△t_3}\\ {\eta }_{13}{\sigma }_{△t_1△t_3}& {\eta }_{23}{\sigma }_{△t_2△t_3}& {\sigma }_{△t_3}^2\end{array}\right]$$

令：
$$E[\mathbf{dS}\cdot {\mathbf{dS}}^T]=\left[\begin{array}{c}{k}_0-k_1\\ k_0-k_2\\ k_0-k_3\end{array}\right]\ast \left[\begin{array}{c}{k}_0-k_1{k}_0-k_2{k}_0-k_3\end{array}\right]$$

对上述公式（13）-（15）的部分变量进行注解：

${\sigma }_{△r_i}$：是测量误差$△r_i$的标准差

${\eta }_{ij}$：是$△r_i$和$△r_j$的相关系数（或者描述$△t_i$和$△t_j$），计算公式如下：
$${\eta }_{ij}=\frac{cov(△r_i,△r_j)}{{\sigma }_{△r_i}\cdot {\sigma }_{△r_j}}=\frac{cov(△t_i,△t_j)}{{\sigma }_{△t_i}\cdot {\sigma }_{△t_j}}$$

得出答案(各个站各个分量标准差是一致的情况)

假设各测向站在各个分量上的标准差大小一样，故设：
$${\sigma }_{x_i}^2={\sigma }_{y_i}^2={\sigma }_{z_i}^2={\sigma }_s^{2}i=0,1,2,3$$
同时由于之前的推导，可以得知：
$$F_{ix}^2+F_{iy}^2+F_{iz}^2=1i=0,1,2,3$$

即：
$$E[\mathbf{dS}\cdot {\mathbf{dS}}^T]=\left[\begin{array}{ccc}2{\sigma }_s^2& {\sigma }_s^2& {\sigma }_s^2\\ {\sigma }_s^2& 2{\sigma }_s^2& {\sigma }_s^2\\ {\sigma }_s^2& {\sigma }_s^2& 2{\sigma }_s^2\end{array}\right]$$
设
$$E[\mathbf{d}△\mathbf{r}\cdot {\mathbf{d}△\mathbf{r}}^T]+E[\mathbf{dS}\cdot {\mathbf{dS}}^T]=({\alpha }_{ij}{)}_{3\ast 3}$$
即：
$${\alpha }_{ij}=\begin{array}{c}\{\begin{array}{l}{c}^2\ast {\sigma }_{△t_i}^2+2\ast {\sigma }_s^{2}i=j\\ c^2\ast {\eta }_{ij}{\sigma }_{△r_i△r_j}+{\sigma }_s^{2}i\ne j\end{array}\end{array}$$
联立上述式子：
P d r = C { E [ d △ r ⋅ d △ r T ] + E [ d S ⋅ d S T ] } C T = C [ σ △ r 1 2 + 2 σ s 2 η 12 σ △ r 1 △ r 2 + σ s 2 η 13 σ △ r 1 △ r 3 + σ s 2 η 12 σ △ r 1 △ r 2 + σ s 2 σ △ r 2 2 + 2 σ s 2 η 23 σ △ r 2 △ r 3 + σ s 2 η 13 σ △ r 1 △ r 3 + σ s 2 η 23 σ △ r 2 △ r 3 + σ s 2 σ △ r 3 2 + 2 σ s 2 ] C T \boldsymbol{P_{dr}} = \boldsymbol{C}\{ E[\boldsymbol{d\triangle r} · \boldsymbol{{d\triangle r}}^T]+E[\boldsymbol{dS}·{\boldsymbol{dS}}^T]\}\boldsymbol{C}^T\\ =\boldsymbol{C} \begin{bmatrix} {\sigma^2_{\triangle r_1}}+2\sigma^2_s & \eta_{12} \sigma_{\triangle r_1 \triangle r_2}+\sigma^2_s & \eta_{13} \sigma_{\triangle r_1 \triangle r_3}+\sigma^2_s \\ \eta_{12} \sigma_{\triangle r_1 \triangle r_2}+\sigma^2_s &{\sigma^2_{\triangle r_2}}+2\sigma^2_s & \eta_{23} \sigma_{\triangle r_2 \triangle r_3}+\sigma^2_s \\ \eta_{13} \sigma_{\triangle r_1 \triangle r_3}+\sigma^2_s & \eta_{23} \sigma_{\triangle r_2 \triangle r_3}+\sigma^2_s&{\sigma^2_{\triangle r_3}}+2\sigma^2_s \\ \end{bmatrix}\boldsymbol{C}^T Pdr​=C{E[d△r⋅d△rT]+E[dS⋅dST]}CT=C ​σ△r1​2​+2σs2​η12​σ△r1​△r2​​+σs2​η13​σ△r1​△r3​​+σs2​​η12​σ△r1​△r2​​+σs2​σ△r2​2​+2σs2​η23​σ△r2​△r3​​+σs2​​η13​σ△r1​△r3​​+σs2​η23​σ△r2​△r3​​+σs2​σ△r3​2​+2σs2​​ ​CT
$$=\mathbf{C}\left[\begin{array}{ccc}{c}^2\ast {\sigma }_{△t_1}^2+2{\sigma }_s^2& c^2\ast {\eta }_{12}{\sigma }_{△t_1△t_2}+{\sigma }_s^2& c^2\ast {\eta }_{13}{\sigma }_{△t_1△t_3}+{\sigma }_s^2\\ c^2\ast {\eta }_{12}{\sigma }_{△t_1△t_2}+{\sigma }_s^2& c^2\ast {\sigma }_{△t_2}^2+2{\sigma }_s^2& c^2\ast {\eta }_{23}{\sigma }_{△t_2△t_3}+{\sigma }_s^2\\ c^2\ast {\eta }_{13}{\sigma }_{△t_1△t_3}+{\sigma }_s^2& c^2\ast {\eta }_{23}{\sigma }_{△t_2△t_3}+{\sigma }_s^2& c^2\ast {\sigma }_{△t_3}^2+2{\sigma }_s^2\end{array}\right]{\mathbf{C}}^T$$
最终的求解结果：
$$\mathbf{GDOP}=\sqrt{trace({\mathbf{P}}_{\mathbf{dr}})}$$

得出答案(各个站各个分量标准差不一致的情况)

$$E[\mathbf{dS}\cdot {\mathbf{dS}}^T]=\left[\begin{array}{c}{k}_0-k_1\\ k_0-k_2\\ k_0-k_3\end{array}\right]\ast \left[\begin{array}{c}{k}_0-k_1{k}_0-k_2{k}_0-k_3\end{array}\right]$$
$$=\left[\begin{array}{ccc}(k_0-k_1{)}^2& (k_0-k_1)(k_0-k_2)& (k_0-k_1)(k_0-k_3)\\ (k_0-k_1)(k_0-k_2)& (k_0-k_2{)}^2& (k_0-k_2)(k_0-k_3)\\ (k_0-k_1)(k_0-k_3)& (k_0-k_2)(k_0-k_3)& (k_0-k_3{)}^2\end{array}\right]$$

但是那个情况是各站点本身的位置的误差，一般不进行讨论。