[发明专利]基于TDOA和功率测量值的GNSS欺骗干扰定位方法

权利要求书

1.一种基于TDOA和功率测量值的GNSS欺骗干扰定位方法，其特征在于，包括以下步骤：

S1：利用多个位置已知且静止的GNSS接收机组成一个欺骗干扰源定位系统，将伪距单差TDOA和与伪距单差TDOA对应的功率比的开方作为GNSS欺骗干扰源定位求解过程中的变量；

此欺骗干扰源定位系统由N个GNSS接收机组成，在柯西参考坐标系下，N个GNSS接收机分别位于已知的位置r_i＝[x_i,y_i,z_i]^T，欺骗干扰源位于s^o＝[x^o,y^o,z^o]^T；

伪距单差TDOA模型，以第1个GNSS接收机为参考接收机可以表示为：

ki1=li-l1=dio-d1o+ϵd,i-ϵd,1---(5)]]>

其中，k_i1表示第i个GNSS接收机与第1个GNSS接收机的伪距单差TDOA，l_i为第i个GNSS接收机伪距测量值，l₁为第1个GNSS接收机伪距测量值，为欺骗干扰源到第i个GNSS接收机之间的欧拉距离，为欺骗干扰源到第1个GNSS接收机之间的欧拉距离，ε_d,i为第i个GNSS接收机伪距测量噪声，ε_d,1为第1个GNSS接收机伪距测量噪声；

记k＝[k₂₁,...,k_N1]^T，其中k₂₁,...,k_N1由(5)式定义，则k的协方差Q_k[i-1,j-1]为：

Qk[i-1,j-1]=σd,12,i≠jσd,12+σd,i2,i=j---(6)]]>

其中i,j＝2,3,...,N，为第i个GNSS接收机伪距测量噪声方差，为第1个GNSS接收机伪距测量噪声方差；

与伪距单差TDOA对应的功率比的开方为：

qi1=p1/g1pi/gi=diod1o(1+d1o2ϵp,1pTg1)(1+dio2ϵp,ipTgi)-1---(7)]]>

式中，p_i为第i个GNSS接收机的功率测量值，g_i为第i个GNSS接收机的增益，p_T为欺骗信号的发射功率，ε_p,i为第i个GNSS接收机信号功率测量噪声，ε_p,1为第1个GNSS接收机信号功率测量噪声，当信号的SNR足够大时，满足则对式(7)进行泰勒展开，并忽略2次及以上的项，可得：

qi1≈diod1o1+d1o2ϵp,1pTg1-dio2ϵp,ipTgi≈diod1o(1-12(d1o2ϵp,1pTg1-dio2ϵp,ipTgi))=diod1o-12(d1odioϵp,ipTg1-dio3ϵp,ipTgid1o)---(8)]]>

记q＝[q₂₁,...,q_N1]^T，其中q₂₁,...,q_N1由(8)式定义，则q的协方差为：

Qq[i-1,j-1]=d1o2diodjoσp,124pT2g12,i≠jd1o2dio2σp,124pT2g12+dio6σp,i24pT2gi2d1o2,i=j---(9)]]>

S2：采用加权最小二乘法估计欺骗干扰源信号到达GNSS接收机的时间TOA；

S2.1：将式(5)、(8)分别移项可得相应的测量误差方程分别为：

Δki1=ki1-(dio-d1o)---(10)]]> 1

Δqi1=qi1-diod1o---(11)]]>

式中Δk_i1,Δq_i1分别表示TDOA测量误差和距离比测量误差；

S2.2：将式(11)等号两边同时乘以得：

d1oΔqi1=d1oqi1-dio---(12)]]>

式(10)和(12)均为未知量的线性函数，它们的矩阵形式表示为：

e₁＝h₁-G₁d^o (13)

其中：

e1=[Δk21,...,ΔkN1,d1oΔq21,...,d1oΔqN1]T---(14)]]>

h₁＝[k^T,0_1×N-1]^T (15)

G1=-1N-1×1IN-1-qIN-1---(16)]]>

式中，0_1×N-1表示N-1维全0列向量，1_N-1×1表示N-1维全1行向量，I_N-1表示N维单位矩阵；

S2.3：加权矩阵W₁定义为：

W1=E[e1e1T]-1=diag{Qγ-1,1d1o2Qq-1}---(17)]]>

其中，E[·]^-1表示将均值矩阵取逆，diag{}表示取矩阵主对角线上的元素，表示伪距单差向量k的协方差矩阵取逆，表示将功率比的开方q的协方差矩阵取逆，则可得使最小的加权最小二乘解为：

d=(G1TW1G1)-1G1TW1h1---(18)]]>

式中d表示TOA估计值；

S3：利用S2得到的TOA估计值，采用最小二乘算法完成欺骗干扰源位置定位；

设欺骗干扰源初始位置为s_g，迭代求解过程为式(19)：

m=0,s(0)=sgwhile||Δs(m)||>ηΔs(m+1)=(G2TG2)-1G2Tbs(m+1)=s(m)+Δs(m+1)m=m+1---(19)]]>

式中，

G2=[(s-r1)||s-r1||,...,(s-rN)||s-rN||]T---(20)]]>

b＝[d₁-||s-r₁||,...,d_N-||s-r_N||]^T(21)

Δs＝[Δx,Δy,Δz]^T (22)

m表示迭代次数在while循环中的变量，η为判决门限，其值应大于克拉美罗下限CRLB，Δs^(m)表示第m次循环的最小二乘解，Δs^(m+1)表示第m+1次循环的最小二乘解，s^(m)表示第m次循环所求得的欺骗干扰源位置，s^(m+1)表示第m+1次循环所求得的欺骗干扰源位置，s表示循环过程中获得的欺骗干扰源位置矩阵，r₁表示第1个GNSS接收机的位置，r_N表示第N个GNSS接收机的位置，d₁表示第1个GNSS接收机的TOA估计值，d_N表示第N个GNSS接收机的TOA估计值，Δs＝[Δx,Δy,Δz]^T表示循环当中的最小二乘解；当第m次循环的最小二乘解Δs_(m)的模大于判决门限η时，循环结束，此时所获得s^(m+1)的即为欺骗干扰源位置。