CN105607092B - 基于tdoa和功率测量值的gnss欺骗干扰定位方法 - Google Patents
基于tdoa和功率测量值的gnss欺骗干扰定位方法 Download PDFInfo
- Publication number
- CN105607092B CN105607092B CN201610055511.2A CN201610055511A CN105607092B CN 105607092 B CN105607092 B CN 105607092B CN 201610055511 A CN201610055511 A CN 201610055511A CN 105607092 B CN105607092 B CN 105607092B
- Authority
- CN
- China
- Prior art keywords
- mrow
- msub
- msubsup
- mfrac
- mtd
- 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
Classifications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
- G01S19/00—Satellite radio beacon positioning systems; Determining position, velocity or attitude using signals transmitted by such systems
- G01S19/01—Satellite radio beacon positioning systems transmitting time-stamped messages, e.g. GPS [Global Positioning System], GLONASS [Global Orbiting Navigation Satellite System] or GALILEO
- G01S19/13—Receivers
- G01S19/35—Constructional details or hardware or software details of the signal processing chain
- G01S19/37—Hardware or software details of the signal processing chain
Landscapes
- Engineering & Computer Science (AREA)
- Signal Processing (AREA)
- Radar, Positioning & Navigation (AREA)
- Remote Sensing (AREA)
- Computer Networks & Wireless Communication (AREA)
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Other Investigation Or Analysis Of Materials By Electrical Means (AREA)
- Position Fixing By Use Of Radio Waves (AREA)
Abstract
本发明公开了一种基于TDOA和功率测量值的GNSS欺骗干扰定位方法,包括以下步骤:首先利用多个位置已知且静止的普通商用GNSS接收机组成一个欺骗干扰源定位系统,然后采用加权最小二乘法估计欺骗干扰源信号到达GNSS接收机的时间TOA;最后利用TOA估计值,采用最小二乘算法完成欺骗干扰源位置定位。该方法计算简单,且能避免伪距测量精度低以及接收信号功率弱导致的结果发散异常等问题,其定位精度逼近克拉美罗下限CRLB。
Description
技术领域
本发明涉及卫星导航技术领域,特别涉及一种基于差分到达时间(TimeDifference of Arrival,TDOA)和信号功率测量值的全球导航卫星系统(GlobalPositioning Satellite System,GNSS)欺骗干扰源定位方法。
背景技术
随着GNSS的发展,它在人们的生活中起着越来越重要的作用。然而由于到达接收机的GNSS信号很微弱,导致GNSS信号对带内干扰很脆弱,其中欺骗干扰是危害最大的一种。欺骗干扰的信号结构和功率等参数与真实信号相似,其目的是使接收机输出虚假的位置、时间结果,而不引起使用者注意。这将导致严重的后果,特别是对重要的基础设施,比如输电网络或数字通信网络使用的GNSS授时接收机被欺骗,时间被拉偏,将导致输电故障或通信中断。
因此,抗欺骗技术成为近来GNSS研究领域的一个研究热点。学者们提出了许多抗欺骗技术。带内功率监测、信号质量监测、天线阵技术、加密技术以及多接收机抗欺骗方法等是当前研究中较为流行的抗欺骗技术。这些方法主要集中于欺骗干扰的检测和抑制,对有欺骗干扰源定位的研究却很少。
距离测量值和功率测量值是GNSS接收机的两个基本测量值。因此,一些基于TDOA或能量的无源定位算法可以扩展到欺骗干扰源定位应用中。其中最大似然(MaximumLikelihood,ML)算法是一种较为吸引人的技术,因为它得出解的精度可以达到克拉美罗下限(Cramer-Rao lower bound,CRLB)。然而该算法需要迭代计算,并且需要一个好的初始解,否则计算结果会出现发散或者收敛到局部最优解。
发明内容
本发明的目的是针对上述已有技术的不足,提出一种基于差分到达时间TDOA和信号功率测量值的基于TDOA和功率测量值的GNSS欺骗干扰定位方法。
本发明的技术方案是:
一种基于TDOA和功率测量值的GNSS欺骗干扰定位方法,包括以下步骤:
S1:利用多个位置已知且静止的GNSS接收机组成一个欺骗干扰源定位系统;将伪距单差TDOA和与伪距单差TDOA对应的功率比的开方作为GNSS欺骗干扰源定位求解过程中的变量;
此欺骗干扰源定位系统由N个GNSS接收机组成,在柯西参考坐标系下,N个GNSS接收机分别位于已知的位置ri=[xi,yi,zi]T,此位置可以根据需求自由设定,欺骗干扰源位于so=[xo,yo,zo]T;第i个GNSS接收机接收到欺骗干扰信号xi(t)的模型为:
其中,t表示GNSS时间;pT为欺骗信号的发射功率;gi为第i个GNSS接收机的增益,包含天线增益和接收机前端处理损耗;为欺骗干扰源到第i个GNSS接收机之间的欧拉距离;F(t)为接收到的信号波形,由伪随机码和导航电文调制到载波上生成;τi为信号传播到第i个GNSS接收机的时延;ξi为高斯白噪声,均值为0,方差为N0B,N0为噪声功率谱密度,B为信号带宽;
欺骗干扰信号直视入射且符合空间自由传播模型,第i个GNSS接收机伪距测量值li的模型为:
式中c为信号传播速度;τf为干扰源模拟的虚假时延,dtr和dts分别表示GNSS接收机和欺骗干扰源的钟差,τa=τf+dtr-dts融合不同伪距测量值中相同的分量,εd,i为伪距测量噪声,其服从0均值高斯分布,方差为且各接收机伪距测量值li之间的噪声不相关;
对欺骗干扰信号进行解扩、相干积累可得:
式中Ii为相干累积的结果,T为相干累积间隔,F*(t)为接收到的信号波形F(t)的共轭;ξi′为式(1)噪声分量ξi积分后的结果,易知它是0均值高斯随机变量,方差为N0/T;
对欺骗干扰信号进行解扩和相干积累后结果Ii平方可得信号的功率测量值pi,其模型为:
其中,为信号功率测量噪声,忽略二次项后,它是一个0均值高斯噪声,方差为测量模型(2)、(4)中,作为中间变量,将测量值pi和li与待求解的干扰源位置so联系起来;gi能够通过对接收机校准获得,这里作为已知量;
由于伪距测量值和功率测量值中分别含有未知量τf,pT,因此本方法没有直接使用这些测量值进行定位求解,而采用了伪距单差TDOA和功率比的开方作为变量;
伪距单差TDOA模型,以第1个GNSS接收机为参考接收机可以表示为:
其中,ki1表示第i个GNSS接收机与第1个GNSS接收机的伪距单差TDOA,li为第i个GNSS接收机伪距测量值,l1为第1个GNSS接收机伪距测量值,为欺骗干扰源到第i个GNSS接收机之间的欧拉距离,为欺骗干扰源到第1个GNSS接收机之间的欧拉距离,εd,i为第i个GNSS接收机伪距测量噪声,εd,1为第1个GNSS接收机伪距测量噪声;
记k=[k21,...,kN1]T,其中k21,...,kN1由(5)式定义,则k的协方差Qk[i-1,j-1]为:
其中i,j=2,3,...,N,为第i个GNSS接收机伪距测量噪声方差,为第1个GNSS接收机伪距测量噪声方差;
与伪距单差TDOA对应的功率比的开方为:
式中,εp,i为第i个GNSS接收机信号功率测量噪声,εp,1为第1个GNSS接收机信号功率测量噪声,当信号的SNR足够大时,满足则对式(7)进行泰勒展开,并忽略2次及以上的项,可得:
记q=[q21,...,qN1]T,其中q21,...,qN1由(8)式定义,则q的协方差为:
S2:采用加权最小二乘法估计欺骗干扰源信号到达GNSS接收机的时间TOA,包括以下步骤:
S2.1:将式(5)、(8)分别移项可得相应的测量误差方程分别为:
式中Δki1,Δqi1分别表示TDOA测量误差和距离比测量误差;
S2.2:将式(11)等号两边同时乘以得:
方程(10)和(12)均为未知量的线性函数,他们的矩阵形式表示为:
e1=h1-G1do (13)
其中:
h1=[kT,01×N-1]T (15)
式中,01×N-1表示N-1维全0列向量,1N-1×1表示N-1维全1行向量,IN-1表示N维单位矩阵;
S2.3:加权矩阵W1定义为:
其中,E[·]-1表示将均值矩阵取逆,diag{}表示取矩阵主对角线上的元素,表示伪距单差向量k的协方差矩阵取逆,表示将功率比的开方q的协方差矩阵取逆,则可得使最小的加权最小二乘解为:
式中d表示TOA估计值;
S3:利用S2得到的TOA估计值,采用最小二乘算法完成欺骗干扰源位置定位,设欺骗干扰源初始位置为sg,迭代求解过程为式(19):
式中,
b=[d1-||s-r1||,...,dN-||s-rN||]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表示循环过程中获得的欺骗干扰源位置矩阵,r1表示第1个GNSS接收机的位置,rN表示第N个GNSS接收机的位置,d1表示第1个GNSS接收机的TOA估计值,dN表示第N个GNSS接收机的TOA估计值,Δs=[Δx,Δy,Δz]T表示循环当中的最小二乘解;当第m次循环的最小二乘解Δs(m)的模大于判决门限η时,循环结束,此时所获得s(m+1)的即为欺骗干扰源位置。
在步骤S1中,所有GNSS接收机共用同一个时钟源,不同接收机的钟差相同。
在步骤S1中,所述参考接收机可以为任意编号i的接收机。
与现有技术相比,本发明所具有的有益效果为:
本发明针对GNSS欺骗干扰抑制消除问题,提出了一种欺骗干扰源定位算法。利用商用级GNSS接收机的差分到达时间TDOA测量值和信号功率测量值,该算法分三步完成对干扰源的定位。首先搭建欺骗干扰源定位系统,步骤S2为闭式解析解,步骤S3的迭代次数少,因此算法的计算量小,可实现实时计算。仿真结果表明本文算法趋近克拉美罗下限CRLB。
附图说明
图1是本发明的流程图。
图2是本发明方法的欺骗干扰源定位场景空间示意图。
图3是本发明方法的欺骗干扰源定位近场仿真中的定位误差MSE;在近场干扰源定位性能仿真验证中,干扰源与参考接收机的距离为Rs=100m,俯仰角固定为βs=45°,方位角取值范围为[-30°,30°]。
图4是本发明方法的欺骗干扰源定位远场仿真中的定位误差MSE;在远场干扰源定位性能仿真验证中,干扰源与参考接收机的距离为Rs=300m,俯仰角固定为βs=45°,方位角取值范围为[-30°,30°]。
具体实施方式
以下结合附图对本发明的具体实施例进行详细描述,但不构成对本发明的限制。
一种基于TDOA和功率测量值的GNSS欺骗干扰定位方法,其流程如图1所示,包括以下步骤:
S1:利用多个位置已知且静止的普通商用GNSS接收机组成一个欺骗干扰源定位系统,如图2所示;此欺骗干扰源定位系统由N个GNSS接收机组成,在柯西参考坐标系下,N个GNSS接收机分别位于已知的位置ri=[xi,yi,zi]T,此位置可以根据需求自由设定,而且N个GNSS接收机共用同一个时钟源,不同接收机的钟差相同,欺骗干扰源位于so=[xo,yo,zo]T;第i个GNSS接收机接收到欺骗干扰信号xi(t)的模型为:
其中,t表示GNSS时间;pT为欺骗信号的发射功率;gi为第i个GNSS接收机的增益,包含天线增益和接收机前端处理损耗;为欺骗干扰源到第i个GNSS接收机之间的欧拉距离;F(t)为接收到的信号波形,由伪随机码和导航电文调制到载波上生成;τi为信号传播到第i个GNSS接收机的时延;ξi为高斯白噪声,均值为0,方差为N0B,N0为噪声功率谱密度,B为信号带宽;
欺骗干扰信号直视入射且符合空间自由传播模型,第i个GNSS接收机伪距测量值li的模型为:
式中c为信号传播速度;τf为干扰源模拟的虚假时延,dtr和dts分别表示GNSS接收机和欺骗干扰源的钟差,τa=τf+dtr-dts融合不同伪距测量值中相同的分量,εd,i为伪距测量噪声,其服从0均值高斯分布,方差为且各接收机伪距测量值li之间的噪声不相关;
对欺骗干扰信号的进行解扩、相干积累可得:
式中Ii为相干累积的结果,T为相干累积间隔,F*(t)为接收到的信号波形F(t)的共轭;ξi′为式(1)噪声分量ξi积分后的结果,易知它是0均值高斯随机变量,方差为N0/T;
对欺骗干扰信号进行解扩和相干积累后结果Ii平方可得信号的功率测量值pi,其模型为:
其中,为信号功率测量噪声,忽略二次项后,它是一个0均值高斯噪声,方差为测量模型(24)、(26)中,作为中间变量,将测量值pi和li与待求解的干扰源位置so联系起来;gi能够通过对接收机校准获得,这里作为已知量;
由于伪距测量值和功率测量值中分别含有未知量τf,pT,因此本方法没有直接使用这些测量值进行定位求解,而采用了伪距单差TDOA和功率比的开方作为变量;
伪距单差TDOA模型,以第1个GNSS接收机为参考接收机可以表示为:
其中,ki1表示第i个GNSS接收机与第1个GNSS接收机的伪距单差TDOA,li为第i个GNSS接收机伪距测量值,l1为第1个GNSS接收机伪距测量值,为欺骗干扰源到第i个GNSS接收机之间的欧拉距离,为欺骗干扰源到第1个GNSS接收机之间的欧拉距离,εd,i为第i个GNSS接收机伪距测量噪声,εd,1为第1个GNSS接收机伪距测量噪声;这里的参考接收机可以为任意编号i的接收机,不限于第一个GNSS接收机;
参考接收机可以为任意编号i的接收机,记k=[k21,...,kN1]T,其中k21,...,kN1由(27)式定义,则k的协方差为:
其中i,j=2,3,...,N,为第i个GNSS接收机伪距测量噪声方差,为第1个GNSS接收机伪距测量噪声方差;
与伪距单差TDOA对应的功率比的开方为:
式中,εp,i为第i个GNSS接收机信号功率测量噪声,εp,1为第1个GNSS接收机信号功率测量噪声,当信号的SNR足够大时,满足则对式(29)进行泰勒展开,并忽略2次及以上的项,可得:
记q=[q21,...,qN1]T,其中q21,...,qN1由(29)式定义,则q的协方差为:
S2:采用加权最小二乘法估计欺骗干扰源信号到达GNSS接收机的时间TOA,包括以下步骤:
A1:将式(27)、(30)分别移项可得相应的测量误差方程分别为:
式中Δki1,Δqi1分别表示TDOA测量误差和距离比测量误差;
A2:将式(33)等号两边同时乘以得:
方程(32)和(34)均为未知量的线性函数,他们的矩阵形式表示为:
e1=h1-G1do (35)
其中:
h1=[kT,01×N-1]T (37)
式中,01×N-1表示N-1维全0列向量,1N-1×1表示N-1维全1行向量,IN-1表示N维单位矩阵;
A3:加权矩阵W1定义为:
其中,E[·]-1表示将均值矩阵取逆,diag{}表示取矩阵主对角线上的元素,表示伪距单差向量k的协方差矩阵取逆,表示将功率比的开方q的协方差矩阵取逆,则可得使最小的加权最小二乘解为:
式中d表示TOA估计值;
S3:利用S2得到的TOA估计值,采用最小二乘算法完成欺骗干扰源位置定位,设欺骗干扰源初始位置为sg,迭代求解过程为式(41):
式中,
b=[d1-||s-r1||,...,dN-||s-rN||]T (43)
Δs=[Δx,Δy,Δz]T (44)
m表示迭代次数在while循环中的变量,η为判决门限,其值应大于克拉美罗下限CRLB,Δs(m)表示第m次循环的最小二乘解,Δs(m+1)表示第m+1次循环的最小二乘解,s(m)表示第m次循环所求得的欺骗干扰源位置,s(m+1)表示第m+1次循环所求得的欺骗干扰源位置,s表示循环过程中获得的欺骗干扰源位置矩阵,r1表示第1个GNSS接收机的位置,rN表示第N个GNSS接收机的位置,d1表示第1个GNSS接收机的TOA估计值,dN表示第N个GNSS接收机的TOA估计值,Δs=[Δx,Δy,Δz]T表示循环当中的最小二乘解;当第m次循环的最小二乘解Δs(m)的模大于判决门限η时,循环结束,此时所获得s(m+1)的即为欺骗干扰源位置。
本发明上述方法可以总结为,算法的计算由式(40)和(41)组成;加权矩阵W1如式(39)所示;因为Qq与距离有关,加权矩阵W1的计算需要知道干扰源位置;为了执行计算,可以先将W1置为单位矩阵,利用式(40)求解一个初始的粗略距离估计值;然后利用得到的估计值计算W1,并利用式(40)求解精确的距离值。
近场干扰源定位性能仿真验证结果如图3所示,解算误差的MSE大约比CRLB高1dB;远场干扰源定位性能仿真验证结果如图4所示,解算误差的MSE大约比CRLB高0.5dB。
Claims (6)
1.一种基于TDOA和功率测量值的GNSS欺骗干扰定位方法,其特征在于,包括以下步骤:
S1:利用多个位置已知且静止的GNSS接收机组成一个欺骗干扰源定位系统,将伪距单差TDOA和与伪距单差TDOA对应的功率比的开方作为GNSS欺骗干扰源定位求解过程中的变量;
此欺骗干扰源定位系统由N个GNSS接收机组成,在柯西参考坐标系下,N个GNSS接收机分别位于已知的位置ri=[xi,yi,zi]T,欺骗干扰源位于so=[xo,yo,zo]T;
伪距单差TDOA模型,以第1个GNSS接收机为参考接收机可以表示为:
<mrow>
<msub>
<mi>k</mi>
<mrow>
<mi>i</mi>
<mn>1</mn>
</mrow>
</msub>
<mo>=</mo>
<msub>
<mi>l</mi>
<mi>i</mi>
</msub>
<mo>-</mo>
<msub>
<mi>l</mi>
<mn>1</mn>
</msub>
<mo>=</mo>
<msubsup>
<mi>d</mi>
<mi>i</mi>
<mi>o</mi>
</msubsup>
<mo>-</mo>
<msubsup>
<mi>d</mi>
<mn>1</mn>
<mi>o</mi>
</msubsup>
<mo>+</mo>
<msub>
<mi>&epsiv;</mi>
<mrow>
<mi>d</mi>
<mo>,</mo>
<mi>i</mi>
</mrow>
</msub>
<mo>-</mo>
<msub>
<mi>&epsiv;</mi>
<mrow>
<mi>d</mi>
<mo>,</mo>
<mn>1</mn>
</mrow>
</msub>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>5</mn>
<mo>)</mo>
</mrow>
</mrow>
其中,ki1表示第i个GNSS接收机与第1个GNSS接收机的伪距单差TDOA,li为第i个GNSS接收机伪距测量值,l1为第1个GNSS接收机伪距测量值,为欺骗干扰源到第i个GNSS接收机之间的欧拉距离,为欺骗干扰源到第1个GNSS接收机之间的欧拉距离,εd,i为第i个GNSS接收机伪距测量噪声,εd,1为第1个GNSS接收机伪距测量噪声;
记k=[k21,...,kN1]T,其中k21,...,kN1由(5)式定义,则k的协方差Qk[i-1,j-1]为:
<mrow>
<msub>
<mi>Q</mi>
<mi>k</mi>
</msub>
<mo>&lsqb;</mo>
<mi>i</mi>
<mo>-</mo>
<mn>1</mn>
<mo>,</mo>
<mi>j</mi>
<mo>-</mo>
<mn>1</mn>
<mo>&rsqb;</mo>
<mo>=</mo>
<mfenced open = "{" close = "">
<mtable>
<mtr>
<mtd>
<mrow>
<msubsup>
<mi>&sigma;</mi>
<mrow>
<mi>d</mi>
<mo>,</mo>
<mn>1</mn>
</mrow>
<mn>2</mn>
</msubsup>
<mo>,</mo>
</mrow>
</mtd>
<mtd>
<mrow>
<mi>i</mi>
<mo>&NotEqual;</mo>
<mi>j</mi>
</mrow>
</mtd>
</mtr>
<mtr>
<mtd>
<mrow>
<msubsup>
<mi>&sigma;</mi>
<mrow>
<mi>d</mi>
<mo>,</mo>
<mn>1</mn>
</mrow>
<mn>2</mn>
</msubsup>
<mo>+</mo>
<msubsup>
<mi>&sigma;</mi>
<mrow>
<mi>d</mi>
<mo>,</mo>
<mi>i</mi>
</mrow>
<mn>2</mn>
</msubsup>
<mo>,</mo>
</mrow>
</mtd>
<mtd>
<mrow>
<mi>i</mi>
<mo>=</mo>
<mi>j</mi>
</mrow>
</mtd>
</mtr>
</mtable>
</mfenced>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>6</mn>
<mo>)</mo>
</mrow>
</mrow>
其中i,j=2,3,...,N,为第i个GNSS接收机伪距测量噪声方差,为第1个GNSS接收机伪距测量噪声方差;
与伪距单差TDOA对应的功率比的开方为:
<mrow>
<msub>
<mi>q</mi>
<mrow>
<mi>i</mi>
<mn>1</mn>
</mrow>
</msub>
<mo>=</mo>
<msqrt>
<mfrac>
<mrow>
<msub>
<mi>p</mi>
<mn>1</mn>
</msub>
<mo>/</mo>
<msub>
<mi>g</mi>
<mn>1</mn>
</msub>
</mrow>
<mrow>
<msub>
<mi>p</mi>
<mi>i</mi>
</msub>
<mo>/</mo>
<msub>
<mi>g</mi>
<mi>i</mi>
</msub>
</mrow>
</mfrac>
</msqrt>
<mo>=</mo>
<mfrac>
<msubsup>
<mi>d</mi>
<mi>i</mi>
<mi>o</mi>
</msubsup>
<msubsup>
<mi>d</mi>
<mn>1</mn>
<mi>o</mi>
</msubsup>
</mfrac>
<msqrt>
<mrow>
<mo>(</mo>
<mn>1</mn>
<mo>+</mo>
<mfrac>
<mrow>
<msubsup>
<mi>d</mi>
<mn>1</mn>
<mrow>
<mi>o</mi>
<mn>2</mn>
</mrow>
</msubsup>
<msub>
<mi>&epsiv;</mi>
<mrow>
<mi>p</mi>
<mo>,</mo>
<mn>1</mn>
</mrow>
</msub>
</mrow>
<mrow>
<msub>
<mi>p</mi>
<mi>T</mi>
</msub>
<msub>
<mi>g</mi>
<mn>1</mn>
</msub>
</mrow>
</mfrac>
<mo>)</mo>
<msup>
<mrow>
<mo>(</mo>
<mn>1</mn>
<mo>+</mo>
<mfrac>
<mrow>
<msubsup>
<mi>d</mi>
<mi>i</mi>
<mrow>
<mi>o</mi>
<mn>2</mn>
</mrow>
</msubsup>
<msub>
<mi>&epsiv;</mi>
<mrow>
<mi>p</mi>
<mo>,</mo>
<mi>i</mi>
</mrow>
</msub>
</mrow>
<mrow>
<msub>
<mi>p</mi>
<mi>T</mi>
</msub>
<msub>
<mi>g</mi>
<mi>i</mi>
</msub>
</mrow>
</mfrac>
<mo>)</mo>
</mrow>
<mrow>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msup>
</mrow>
</msqrt>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>7</mn>
<mo>)</mo>
</mrow>
</mrow>
式中,pi为第i个GNSS接收机的功率测量值,gi为第i个GNSS接收机的增益,pT为欺骗信号的发射功率,εp,i为第i个GNSS接收机信号功率测量噪声,εp,1为第1个GNSS接收机信号功率测量噪声,当信号的SNR足够大时,满足则对式(7)进行泰勒展开,并忽略2次及以上的项,可得:
<mrow>
<msub>
<mi>q</mi>
<mrow>
<mi>i</mi>
<mn>1</mn>
</mrow>
</msub>
<mo>&ap;</mo>
<mfrac>
<msubsup>
<mi>d</mi>
<mi>i</mi>
<mi>o</mi>
</msubsup>
<msubsup>
<mi>d</mi>
<mn>1</mn>
<mi>o</mi>
</msubsup>
</mfrac>
<msqrt>
<mrow>
<mn>1</mn>
<mo>+</mo>
<mfrac>
<mrow>
<msubsup>
<mi>d</mi>
<mn>1</mn>
<mrow>
<mi>o</mi>
<mn>2</mn>
</mrow>
</msubsup>
<msub>
<mi>&epsiv;</mi>
<mrow>
<mi>p</mi>
<mo>,</mo>
<mn>1</mn>
</mrow>
</msub>
</mrow>
<mrow>
<msub>
<mi>p</mi>
<mi>T</mi>
</msub>
<msub>
<mi>g</mi>
<mn>1</mn>
</msub>
</mrow>
</mfrac>
<mo>-</mo>
<mfrac>
<mrow>
<msubsup>
<mi>d</mi>
<mi>i</mi>
<mrow>
<mi>o</mi>
<mn>2</mn>
</mrow>
</msubsup>
<msub>
<mi>&epsiv;</mi>
<mrow>
<mi>p</mi>
<mo>,</mo>
<mi>i</mi>
</mrow>
</msub>
</mrow>
<mrow>
<msub>
<mi>p</mi>
<mi>T</mi>
</msub>
<msub>
<mi>g</mi>
<mi>i</mi>
</msub>
</mrow>
</mfrac>
</mrow>
</msqrt>
<mo>&ap;</mo>
<mfrac>
<msubsup>
<mi>d</mi>
<mi>i</mi>
<mi>o</mi>
</msubsup>
<msubsup>
<mi>d</mi>
<mn>1</mn>
<mi>o</mi>
</msubsup>
</mfrac>
<mrow>
<mo>(</mo>
<mn>1</mn>
<mo>-</mo>
<mfrac>
<mn>1</mn>
<mn>2</mn>
</mfrac>
<mo>(</mo>
<mrow>
<mfrac>
<mrow>
<msubsup>
<mi>d</mi>
<mn>1</mn>
<mrow>
<mi>o</mi>
<mn>2</mn>
</mrow>
</msubsup>
<msub>
<mi>&epsiv;</mi>
<mrow>
<mi>p</mi>
<mo>,</mo>
<mn>1</mn>
</mrow>
</msub>
</mrow>
<mrow>
<msub>
<mi>p</mi>
<mi>T</mi>
</msub>
<msub>
<mi>g</mi>
<mn>1</mn>
</msub>
</mrow>
</mfrac>
<mo>-</mo>
<mfrac>
<mrow>
<msubsup>
<mi>d</mi>
<mi>i</mi>
<mrow>
<mi>o</mi>
<mn>2</mn>
</mrow>
</msubsup>
<msub>
<mi>&epsiv;</mi>
<mrow>
<mi>p</mi>
<mo>,</mo>
<mi>i</mi>
</mrow>
</msub>
</mrow>
<mrow>
<msub>
<mi>p</mi>
<mi>T</mi>
</msub>
<msub>
<mi>g</mi>
<mi>i</mi>
</msub>
</mrow>
</mfrac>
</mrow>
<mo>)</mo>
<mo>)</mo>
</mrow>
<mo>=</mo>
<mfrac>
<msubsup>
<mi>d</mi>
<mi>i</mi>
<mi>o</mi>
</msubsup>
<msubsup>
<mi>d</mi>
<mn>1</mn>
<mi>o</mi>
</msubsup>
</mfrac>
<mo>-</mo>
<mfrac>
<mn>1</mn>
<mn>2</mn>
</mfrac>
<mrow>
<mo>(</mo>
<mfrac>
<mrow>
<msubsup>
<mi>d</mi>
<mn>1</mn>
<mi>o</mi>
</msubsup>
<msubsup>
<mi>d</mi>
<mi>i</mi>
<mi>o</mi>
</msubsup>
<msub>
<mi>&epsiv;</mi>
<mrow>
<mi>p</mi>
<mo>,</mo>
<mi>i</mi>
</mrow>
</msub>
</mrow>
<mrow>
<msub>
<mi>p</mi>
<mi>T</mi>
</msub>
<msub>
<mi>g</mi>
<mn>1</mn>
</msub>
</mrow>
</mfrac>
<mo>-</mo>
<mfrac>
<mrow>
<msubsup>
<mi>d</mi>
<mi>i</mi>
<mrow>
<mi>o</mi>
<mn>3</mn>
</mrow>
</msubsup>
<msub>
<mi>&epsiv;</mi>
<mrow>
<mi>p</mi>
<mo>,</mo>
<mi>i</mi>
</mrow>
</msub>
</mrow>
<mrow>
<msub>
<mi>p</mi>
<mi>T</mi>
</msub>
<msub>
<mi>g</mi>
<mi>i</mi>
</msub>
<msubsup>
<mi>d</mi>
<mn>1</mn>
<mi>o</mi>
</msubsup>
</mrow>
</mfrac>
<mo>)</mo>
</mrow>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>8</mn>
<mo>)</mo>
</mrow>
</mrow>
记q=[q21,...,qN1]T,其中q21,...,qN1由(8)式定义,则q的协方差为:
<mrow>
<msub>
<mi>Q</mi>
<mi>q</mi>
</msub>
<mo>&lsqb;</mo>
<mi>i</mi>
<mo>-</mo>
<mn>1</mn>
<mo>,</mo>
<mi>j</mi>
<mo>-</mo>
<mn>1</mn>
<mo>&rsqb;</mo>
<mo>=</mo>
<mfenced open = "{" close = "">
<mtable>
<mtr>
<mtd>
<mrow>
<mfrac>
<mrow>
<msubsup>
<mi>d</mi>
<mn>1</mn>
<mrow>
<mi>o</mi>
<mn>2</mn>
</mrow>
</msubsup>
<msubsup>
<mi>d</mi>
<mi>i</mi>
<mi>o</mi>
</msubsup>
<msubsup>
<mi>d</mi>
<mi>j</mi>
<mi>o</mi>
</msubsup>
<msubsup>
<mi>&sigma;</mi>
<mrow>
<mi>p</mi>
<mo>,</mo>
<mn>1</mn>
</mrow>
<mn>2</mn>
</msubsup>
</mrow>
<mrow>
<mn>4</mn>
<msubsup>
<mi>p</mi>
<mi>T</mi>
<mn>2</mn>
</msubsup>
<msubsup>
<mi>g</mi>
<mn>1</mn>
<mn>2</mn>
</msubsup>
</mrow>
</mfrac>
<mo>,</mo>
</mrow>
</mtd>
<mtd>
<mrow>
<mi>i</mi>
<mo>&NotEqual;</mo>
<mi>j</mi>
</mrow>
</mtd>
</mtr>
<mtr>
<mtd>
<mrow>
<mfrac>
<mrow>
<msubsup>
<mi>d</mi>
<mn>1</mn>
<mrow>
<mi>o</mi>
<mn>2</mn>
</mrow>
</msubsup>
<msubsup>
<mi>d</mi>
<mi>i</mi>
<mrow>
<mi>o</mi>
<mn>2</mn>
</mrow>
</msubsup>
<msubsup>
<mi>&sigma;</mi>
<mrow>
<mi>p</mi>
<mo>,</mo>
<mn>1</mn>
</mrow>
<mn>2</mn>
</msubsup>
</mrow>
<mrow>
<mn>4</mn>
<msubsup>
<mi>p</mi>
<mi>T</mi>
<mn>2</mn>
</msubsup>
<msubsup>
<mi>g</mi>
<mn>1</mn>
<mn>2</mn>
</msubsup>
</mrow>
</mfrac>
<mo>+</mo>
<mfrac>
<mrow>
<msubsup>
<mi>d</mi>
<mi>i</mi>
<mrow>
<mi>o</mi>
<mn>6</mn>
</mrow>
</msubsup>
<msubsup>
<mi>&sigma;</mi>
<mrow>
<mi>p</mi>
<mo>,</mo>
<mi>i</mi>
</mrow>
<mn>2</mn>
</msubsup>
</mrow>
<mrow>
<mn>4</mn>
<msubsup>
<mi>p</mi>
<mi>T</mi>
<mn>2</mn>
</msubsup>
<msubsup>
<mi>g</mi>
<mi>i</mi>
<mn>2</mn>
</msubsup>
<msubsup>
<mi>d</mi>
<mn>1</mn>
<mrow>
<mi>o</mi>
<mn>2</mn>
</mrow>
</msubsup>
</mrow>
</mfrac>
<mo>,</mo>
</mrow>
</mtd>
<mtd>
<mrow>
<mi>i</mi>
<mo>=</mo>
<mi>j</mi>
</mrow>
</mtd>
</mtr>
</mtable>
</mfenced>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>9</mn>
<mo>)</mo>
</mrow>
</mrow>
S2:采用加权最小二乘法估计欺骗干扰源信号到达GNSS接收机的时间TOA;
S2.1:将式(5)、(8)分别移项可得相应的测量误差方程分别为:
<mrow>
<msub>
<mi>&Delta;k</mi>
<mrow>
<mi>i</mi>
<mn>1</mn>
</mrow>
</msub>
<mo>=</mo>
<msub>
<mi>k</mi>
<mrow>
<mi>i</mi>
<mn>1</mn>
</mrow>
</msub>
<mo>-</mo>
<mrow>
<mo>(</mo>
<msubsup>
<mi>d</mi>
<mi>i</mi>
<mi>o</mi>
</msubsup>
<mo>-</mo>
<msubsup>
<mi>d</mi>
<mn>1</mn>
<mi>o</mi>
</msubsup>
<mo>)</mo>
</mrow>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>10</mn>
<mo>)</mo>
</mrow>
</mrow>
1
<mrow>
<msub>
<mi>&Delta;q</mi>
<mrow>
<mi>i</mi>
<mn>1</mn>
</mrow>
</msub>
<mo>=</mo>
<msub>
<mi>q</mi>
<mrow>
<mi>i</mi>
<mn>1</mn>
</mrow>
</msub>
<mo>-</mo>
<mfrac>
<msubsup>
<mi>d</mi>
<mi>i</mi>
<mi>o</mi>
</msubsup>
<msubsup>
<mi>d</mi>
<mn>1</mn>
<mi>o</mi>
</msubsup>
</mfrac>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>11</mn>
<mo>)</mo>
</mrow>
</mrow>
式中Δki1,Δqi1分别表示TDOA测量误差和距离比测量误差;
S2.2:将式(11)等号两边同时乘以得:
<mrow>
<msubsup>
<mi>d</mi>
<mn>1</mn>
<mi>o</mi>
</msubsup>
<msub>
<mi>&Delta;q</mi>
<mrow>
<mi>i</mi>
<mn>1</mn>
</mrow>
</msub>
<mo>=</mo>
<msubsup>
<mi>d</mi>
<mn>1</mn>
<mi>o</mi>
</msubsup>
<msub>
<mi>q</mi>
<mrow>
<mi>i</mi>
<mn>1</mn>
</mrow>
</msub>
<mo>-</mo>
<msubsup>
<mi>d</mi>
<mi>i</mi>
<mi>o</mi>
</msubsup>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>12</mn>
<mo>)</mo>
</mrow>
</mrow>
式(10)和(12)均为未知量的线性函数,它们的矩阵形式表示为:
e1=h1-G1do (13)
其中:
<mrow>
<msub>
<mi>e</mi>
<mn>1</mn>
</msub>
<mo>=</mo>
<msup>
<mrow>
<mo>&lsqb;</mo>
<msub>
<mi>&Delta;k</mi>
<mn>21</mn>
</msub>
<mo>,</mo>
<mo>...</mo>
<mo>,</mo>
<msub>
<mi>&Delta;k</mi>
<mrow>
<mi>N</mi>
<mn>1</mn>
</mrow>
</msub>
<mo>,</mo>
<msubsup>
<mi>d</mi>
<mn>1</mn>
<mi>o</mi>
</msubsup>
<msub>
<mi>&Delta;q</mi>
<mn>21</mn>
</msub>
<mo>,</mo>
<mo>...</mo>
<mo>,</mo>
<msubsup>
<mi>d</mi>
<mn>1</mn>
<mi>o</mi>
</msubsup>
<msub>
<mi>&Delta;q</mi>
<mrow>
<mi>N</mi>
<mn>1</mn>
</mrow>
</msub>
<mo>&rsqb;</mo>
</mrow>
<mi>T</mi>
</msup>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>14</mn>
<mo>)</mo>
</mrow>
</mrow>
h1=[kT,01×N-1]T (15)
<mrow>
<msub>
<mi>G</mi>
<mn>1</mn>
</msub>
<mo>=</mo>
<mfenced open = "[" close = "]">
<mtable>
<mtr>
<mtd>
<mrow>
<mo>-</mo>
<msub>
<mn>1</mn>
<mrow>
<mi>N</mi>
<mo>-</mo>
<mn>1</mn>
<mo>&times;</mo>
<mn>1</mn>
</mrow>
</msub>
</mrow>
</mtd>
<mtd>
<msub>
<mi>I</mi>
<mrow>
<mi>N</mi>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msub>
</mtd>
</mtr>
<mtr>
<mtd>
<mrow>
<mo>-</mo>
<mi>q</mi>
</mrow>
</mtd>
<mtd>
<msub>
<mi>I</mi>
<mrow>
<mi>N</mi>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msub>
</mtd>
</mtr>
</mtable>
</mfenced>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>16</mn>
<mo>)</mo>
</mrow>
</mrow>
式中,01×N-1表示N-1维全0列向量,1N-1×1表示N-1维全1行向量,IN-1表示N维单位矩阵;
S2.3:加权矩阵W1定义为:
<mrow>
<msub>
<mi>W</mi>
<mn>1</mn>
</msub>
<mo>=</mo>
<mi>E</mi>
<msup>
<mrow>
<mo>&lsqb;</mo>
<msub>
<mi>e</mi>
<mn>1</mn>
</msub>
<msubsup>
<mi>e</mi>
<mn>1</mn>
<mi>T</mi>
</msubsup>
<mo>&rsqb;</mo>
</mrow>
<mrow>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msup>
<mo>=</mo>
<mi>d</mi>
<mi>i</mi>
<mi>a</mi>
<mi>g</mi>
<mo>{</mo>
<msubsup>
<mi>Q</mi>
<mi>&gamma;</mi>
<mrow>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msubsup>
<mo>,</mo>
<mfrac>
<mn>1</mn>
<msubsup>
<mi>d</mi>
<mn>1</mn>
<mrow>
<mi>o</mi>
<mn>2</mn>
</mrow>
</msubsup>
</mfrac>
<msubsup>
<mi>Q</mi>
<mi>q</mi>
<mrow>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msubsup>
<mo>}</mo>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>17</mn>
<mo>)</mo>
</mrow>
</mrow>
其中,E[·]-1表示将均值矩阵取逆,diag{}表示取矩阵主对角线上的元素,表示伪距单差向量k的协方差矩阵取逆,表示将功率比的开方q的协方差矩阵取逆,则可得使最小的加权最小二乘解为:
<mrow>
<mi>d</mi>
<mo>=</mo>
<msup>
<mrow>
<mo>(</mo>
<msubsup>
<mi>G</mi>
<mn>1</mn>
<mi>T</mi>
</msubsup>
<msub>
<mi>W</mi>
<mn>1</mn>
</msub>
<msub>
<mi>G</mi>
<mn>1</mn>
</msub>
<mo>)</mo>
</mrow>
<mrow>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msup>
<msubsup>
<mi>G</mi>
<mn>1</mn>
<mi>T</mi>
</msubsup>
<msub>
<mi>W</mi>
<mn>1</mn>
</msub>
<msub>
<mi>h</mi>
<mn>1</mn>
</msub>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>18</mn>
<mo>)</mo>
</mrow>
</mrow>
式中d表示TOA估计值;
S3:利用S2得到的TOA估计值,采用最小二乘算法完成欺骗干扰源位置定位;
设欺骗干扰源初始位置为sg,迭代求解过程为式(19):
<mrow>
<mtable>
<mtr>
<mtd>
<mrow>
<mi>m</mi>
<mo>=</mo>
<mn>0</mn>
<mo>,</mo>
<msup>
<mi>s</mi>
<mrow>
<mo>(</mo>
<mn>0</mn>
<mo>)</mo>
</mrow>
</msup>
<mo>=</mo>
<msub>
<mi>s</mi>
<mi>g</mi>
</msub>
</mrow>
</mtd>
</mtr>
<mtr>
<mtd>
<mfenced open = "[" close = "">
<mtable>
<mtr>
<mtd>
<mtable>
<mtr>
<mtd>
<mrow>
<mi>w</mi>
<mi>h</mi>
<mi>i</mi>
<mi>l</mi>
<mi>e</mi>
</mrow>
</mtd>
<mtd>
<mrow>
<mo>|</mo>
<mo>|</mo>
<msup>
<mi>&Delta;s</mi>
<mrow>
<mo>(</mo>
<mi>m</mi>
<mo>)</mo>
</mrow>
</msup>
<mo>|</mo>
<mo>|</mo>
<mo>></mo>
<mi>&eta;</mi>
</mrow>
</mtd>
</mtr>
</mtable>
</mtd>
</mtr>
<mtr>
<mtd>
<mrow>
<msup>
<mi>&Delta;s</mi>
<mrow>
<mo>(</mo>
<mi>m</mi>
<mo>+</mo>
<mn>1</mn>
<mo>)</mo>
</mrow>
</msup>
<mo>=</mo>
<msup>
<mrow>
<mo>(</mo>
<msubsup>
<mi>G</mi>
<mn>2</mn>
<mi>T</mi>
</msubsup>
<msub>
<mi>G</mi>
<mn>2</mn>
</msub>
<mo>)</mo>
</mrow>
<mrow>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msup>
<msubsup>
<mi>G</mi>
<mn>2</mn>
<mi>T</mi>
</msubsup>
<mi>b</mi>
</mrow>
</mtd>
</mtr>
<mtr>
<mtd>
<mrow>
<msup>
<mi>s</mi>
<mrow>
<mo>(</mo>
<mi>m</mi>
<mo>+</mo>
<mn>1</mn>
<mo>)</mo>
</mrow>
</msup>
<mo>=</mo>
<msup>
<mi>s</mi>
<mrow>
<mo>(</mo>
<mi>m</mi>
<mo>)</mo>
</mrow>
</msup>
<mo>+</mo>
<msup>
<mi>&Delta;s</mi>
<mrow>
<mo>(</mo>
<mi>m</mi>
<mo>+</mo>
<mn>1</mn>
<mo>)</mo>
</mrow>
</msup>
</mrow>
</mtd>
</mtr>
<mtr>
<mtd>
<mrow>
<mi>m</mi>
<mo>=</mo>
<mi>m</mi>
<mo>+</mo>
<mn>1</mn>
</mrow>
</mtd>
</mtr>
</mtable>
</mfenced>
</mtd>
</mtr>
</mtable>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>19</mn>
<mo>)</mo>
</mrow>
</mrow>
式中,
<mrow>
<msub>
<mi>G</mi>
<mn>2</mn>
</msub>
<mo>=</mo>
<msup>
<mrow>
<mo>&lsqb;</mo>
<mfrac>
<mrow>
<mo>(</mo>
<mi>s</mi>
<mo>-</mo>
<msub>
<mi>r</mi>
<mn>1</mn>
</msub>
<mo>)</mo>
</mrow>
<mrow>
<mo>|</mo>
<mo>|</mo>
<mi>s</mi>
<mo>-</mo>
<msub>
<mi>r</mi>
<mn>1</mn>
</msub>
<mo>|</mo>
<mo>|</mo>
</mrow>
</mfrac>
<mo>,</mo>
<mo>...</mo>
<mo>,</mo>
<mfrac>
<mrow>
<mo>(</mo>
<mi>s</mi>
<mo>-</mo>
<msub>
<mi>r</mi>
<mi>N</mi>
</msub>
<mo>)</mo>
</mrow>
<mrow>
<mo>|</mo>
<mo>|</mo>
<mi>s</mi>
<mo>-</mo>
<msub>
<mi>r</mi>
<mi>N</mi>
</msub>
<mo>|</mo>
<mo>|</mo>
</mrow>
</mfrac>
<mo>&rsqb;</mo>
</mrow>
<mi>T</mi>
</msup>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>20</mn>
<mo>)</mo>
</mrow>
</mrow>
b=[d1-||s-r1||,...,dN-||s-rN||]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表示循环过程中获得的欺骗干扰源位置矩阵,r1表示第1个GNSS接收机的位置,rN表示第N个GNSS接收机的位置,d1表示第1个GNSS接收机的TOA估计值,dN表示第N个GNSS接收机的TOA估计值,Δs=[Δx,Δy,Δz]T表示循环当中的最小二乘解;当第m次循环的最小二乘解Δs(m)的模大于判决门限η时,循环结束,此时所获得s(m+1)的即为欺骗干扰源位置。
2.根据权利要求1所述的基于TDOA和功率测量值的GNSS欺骗干扰定位方法,其特征在于,在步骤S1中,第i个GNSS接收机接收到欺骗干扰信号xi(t)的模型为:
<mrow>
<msub>
<mi>x</mi>
<mi>i</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<mfrac>
<msqrt>
<mrow>
<msub>
<mi>p</mi>
<mi>T</mi>
</msub>
<msub>
<mi>g</mi>
<mi>i</mi>
</msub>
</mrow>
</msqrt>
<msubsup>
<mi>d</mi>
<mi>i</mi>
<mi>o</mi>
</msubsup>
</mfrac>
<mi>F</mi>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>-</mo>
<msub>
<mi>&tau;</mi>
<mi>i</mi>
</msub>
<mo>)</mo>
</mrow>
<mo>+</mo>
<msub>
<mi>&xi;</mi>
<mi>i</mi>
</msub>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>1</mn>
<mo>)</mo>
</mrow>
</mrow>
其中,t表示GNSS时间;pT为欺骗信号的发射功率;gi为第i个GNSS接收机的增益,包含天线增益和接收机前端处理损耗;为欺骗干扰源到第i个GNSS接收机之间的欧拉距离;F(t)为接收到的信号波形,由伪随机码和导航电文调制到载波上生成;τi为信号传播到第i个GNSS接收机的时延;ξi为高斯白噪声,均值为0,方差为N0B,N0为噪声功率谱密度,B为信号带宽。
3.根据权利要求2所述的基于TDOA和功率测量值的GNSS欺骗干扰定位方法,其特征在于,在步骤S1中,所有GNSS接收机共用同一个时钟源,不同接收机的钟差相同。
4.根据权利要求3所述的基于TDOA和功率测量值的GNSS欺骗干扰定位方法,其特征在于,在步骤S1中,欺骗干扰信号直视入射且符合空间自由传播模型。
5.根据权利要求4所述的基于TDOA和功率测量值的GNSS欺骗干扰定位方法,其特征在于,在步骤S1中,第i个GNSS接收机伪距测量值li的模型为:
<mrow>
<mtable>
<mtr>
<mtd>
<mrow>
<msub>
<mi>l</mi>
<mi>i</mi>
</msub>
<mo>=</mo>
<msubsup>
<mi>d</mi>
<mi>i</mi>
<mi>o</mi>
</msubsup>
<mo>+</mo>
<msub>
<mi>c&tau;</mi>
<mi>f</mi>
</msub>
<mo>+</mo>
<mi>c</mi>
<mrow>
<mo>(</mo>
<msub>
<mi>dt</mi>
<mi>r</mi>
</msub>
<mo>-</mo>
<msub>
<mi>dt</mi>
<mi>s</mi>
</msub>
<mo>)</mo>
</mrow>
<mo>+</mo>
<msub>
<mi>&epsiv;</mi>
<mrow>
<mi>d</mi>
<mo>,</mo>
<mi>i</mi>
</mrow>
</msub>
</mrow>
</mtd>
</mtr>
<mtr>
<mtd>
<mrow>
<mo>=</mo>
<msubsup>
<mi>d</mi>
<mi>i</mi>
<mi>o</mi>
</msubsup>
<mo>+</mo>
<msub>
<mi>c&tau;</mi>
<mi>a</mi>
</msub>
<mo>+</mo>
<msub>
<mi>&epsiv;</mi>
<mrow>
<mi>d</mi>
<mo>,</mo>
<mi>i</mi>
</mrow>
</msub>
</mrow>
</mtd>
</mtr>
</mtable>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>2</mn>
<mo>)</mo>
</mrow>
</mrow>
式中c为信号传播速度;τf为干扰源模拟的虚假时延,dtr和dts分别表示GNSS接收机和欺骗干扰源的钟差,τa=τf+dtr-dts;εd,i为伪距测量噪声,其服从0均值高斯分布,方差为且各GNSS接收机伪距测量值li之间的噪声不相关。
6.根据权利要求5所述的基于TDOA和功率测量值的GNSS欺骗干扰定位方法,其特征在于,在步骤S1中,对欺骗干扰信号进行解扩、相干积累可得:
<mrow>
<mtable>
<mtr>
<mtd>
<mrow>
<msub>
<mi>I</mi>
<mi>i</mi>
</msub>
<mo>=</mo>
<mfrac>
<mn>1</mn>
<mi>T</mi>
</mfrac>
<msubsup>
<mo>&Integral;</mo>
<mn>0</mn>
<mi>T</mi>
</msubsup>
<msub>
<mi>x</mi>
<mi>i</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
<msup>
<mi>F</mi>
<mo>*</mo>
</msup>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
<mi>d</mi>
<mi>t</mi>
</mrow>
</mtd>
</mtr>
<mtr>
<mtd>
<mrow>
<mo>=</mo>
<mfrac>
<msqrt>
<mrow>
<msub>
<mi>p</mi>
<mi>T</mi>
</msub>
<msub>
<mi>g</mi>
<mi>i</mi>
</msub>
</mrow>
</msqrt>
<msubsup>
<mi>d</mi>
<mi>i</mi>
<mi>o</mi>
</msubsup>
</mfrac>
<mo>+</mo>
<msubsup>
<mi>&xi;</mi>
<mi>i</mi>
<mo>&prime;</mo>
</msubsup>
</mrow>
</mtd>
</mtr>
</mtable>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>3</mn>
<mo>)</mo>
</mrow>
</mrow>
式中Ii为相干累积的结果,T为相干累积间隔,F*(t)为接收到的信号波形F(t)的共轭;ξ′i为式(1)噪声分量ξi积分后的结果,其是0均值高斯随机变量,方差为N0/T;
对欺骗干扰信号进行解扩和相干积累后结果Ii平方可得信号的功率测量值pi,其模型为:
<mrow>
<msub>
<mi>p</mi>
<mi>i</mi>
</msub>
<mo>=</mo>
<mfrac>
<msub>
<mi>g</mi>
<mi>i</mi>
</msub>
<msubsup>
<mi>d</mi>
<mi>i</mi>
<mrow>
<mi>o</mi>
<mn>2</mn>
</mrow>
</msubsup>
</mfrac>
<msub>
<mi>p</mi>
<mi>T</mi>
</msub>
<mo>+</mo>
<msub>
<mi>&epsiv;</mi>
<mrow>
<mi>p</mi>
<mo>,</mo>
<mi>i</mi>
</mrow>
</msub>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>4</mn>
<mo>)</mo>
</mrow>
</mrow>
3
其中,为信号功率测量噪声,忽略二次项后,它是一个0均值高斯噪声,方差为
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201610055511.2A CN105607092B (zh) | 2016-01-27 | 2016-01-27 | 基于tdoa和功率测量值的gnss欺骗干扰定位方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN201610055511.2A CN105607092B (zh) | 2016-01-27 | 2016-01-27 | 基于tdoa和功率测量值的gnss欺骗干扰定位方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN105607092A CN105607092A (zh) | 2016-05-25 |
CN105607092B true CN105607092B (zh) | 2017-10-31 |
Family
ID=55987168
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN201610055511.2A Active CN105607092B (zh) | 2016-01-27 | 2016-01-27 | 基于tdoa和功率测量值的gnss欺骗干扰定位方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN105607092B (zh) |
Families Citing this family (13)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN106353720B (zh) * | 2016-09-04 | 2019-09-24 | 中国人民解放军海军航空大学 | 基于tdoa/groa的多站连续定位模型 |
CN106772456B (zh) * | 2017-01-12 | 2019-05-10 | 清华大学 | 一种基于多用户协作的转发式欺骗源的定位方法 |
CN107193023B (zh) * | 2017-04-18 | 2020-10-13 | 中国铁建电气化局集团第二工程有限公司 | 一种具有闭式解的高精度北斗卫星系统单点定位方法 |
CN107024616A (zh) * | 2017-05-05 | 2017-08-08 | 郑州云海信息技术有限公司 | 一种提高电磁骚扰源诊断效率的方法及系统 |
CN108196270B (zh) * | 2017-12-31 | 2021-09-21 | 南京理工大学 | 一种双基线欺骗干扰检测算法的性能分析方法 |
CN110879380B (zh) * | 2019-08-19 | 2022-03-22 | 中国电波传播研究所(中国电子科技集团公司第二十二研究所) | 一种基于手机的gnss干扰源定位方法 |
CN110361761A (zh) * | 2019-08-28 | 2019-10-22 | 上海无线电设备研究所 | 一种生成式gnss欺骗干扰方法 |
CN110673168B (zh) * | 2019-09-05 | 2021-09-03 | 清华大学 | 异步多用户联合欺骗信号检测方法及装置 |
CN111522031B (zh) * | 2020-04-28 | 2022-06-07 | 中国南方电网有限责任公司超高压输电公司 | 针对gnss授时应用的多接收机欺骗检测方法 |
CN111751846B (zh) * | 2020-05-22 | 2021-02-26 | 中南民族大学 | 一种无人机载的卫星导航干扰信号检测方法 |
CN112152651B (zh) * | 2020-08-17 | 2021-08-13 | 西安交通大学 | 5g系统面向gnss接收机的干扰源定位方法、存储介质及设备 |
CN113109843B (zh) * | 2021-04-15 | 2022-02-18 | 中国人民解放军63812部队 | 基于双接收机伪距双差的欺骗信号检测抑制方法及装置 |
CN113640839B (zh) * | 2021-06-30 | 2024-05-28 | 湖南天熠电子科技有限公司 | 基于aoa/tdoa的gnss欺骗干扰辐射源定位方法 |
Family Cites Families (5)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
WO2014039100A1 (en) * | 2012-09-10 | 2014-03-13 | Research In Motion Limited | Resource block indication and allocation for in-device coexistence interference avoidance |
CN103727960B (zh) * | 2013-01-30 | 2017-10-03 | 中国人民解放军海军航空工程学院 | 一种基于drfm的无线电高度表干扰信号产生方法 |
CN103760532B (zh) * | 2014-01-13 | 2016-02-03 | 西安电子科技大学 | 干扰机与接收机联合组网的sar欺骗干扰方法 |
CN103941267A (zh) * | 2014-03-17 | 2014-07-23 | 中国民航大学 | 一种结合去噪和doa估计的卫星导航欺骗式干扰抑制方法 |
CN104391305A (zh) * | 2014-11-01 | 2015-03-04 | 中国民航大学 | 基于欺骗式干扰doa估计的卫星导航欺骗式干扰抑制方法 |
-
2016
- 2016-01-27 CN CN201610055511.2A patent/CN105607092B/zh active Active
Also Published As
Publication number | Publication date |
---|---|
CN105607092A (zh) | 2016-05-25 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN105607092B (zh) | 基于tdoa和功率测量值的gnss欺骗干扰定位方法 | |
CN104714244B (zh) | 一种基于抗差自适应Kalman滤波的多系统动态PPP解算方法 | |
Lam et al. | LoRa-based localization systems for noisy outdoor environment | |
CN104316903B (zh) | 一种三站时差定位性能试验评估方法 | |
Ma et al. | TOA localization in the presence of random sensor position errors | |
CN105334495B (zh) | 一种无线网络中基于信号到达时间的非视距稳健定位方法 | |
CN103954931B (zh) | 一种远场和近场混合信号源的定位方法 | |
CN103270801B (zh) | 用于定位用户设备的方位的方法和系统 | |
CN107193023B (zh) | 一种具有闭式解的高精度北斗卫星系统单点定位方法 | |
KR20110112829A (ko) | 일반화된 에러 분포를 사용하는 위치 추정 방법 | |
CN110174643A (zh) | 一种无需噪声功率信息的基于到达时间差的定位方法 | |
US9213100B1 (en) | Bearing-only tracking for horizontal linear arrays with rapid, accurate initiation and a robust track accuracy threshold | |
CN104749554A (zh) | 一种基于递归秩损的幅相误差校准和波达方向估计方法 | |
CN102196559A (zh) | 基于tdoa定位的通道时延误差消除方法 | |
CN105898865A (zh) | 非线性非高斯条件下基于ekf和pf的协同定位方法 | |
CN103616685B (zh) | 基于图像特征的isar图像几何定标方法 | |
CN107870314A (zh) | 基于极化敏感阵列的完备电磁分量加权融合测向优化方法 | |
CN110389326A (zh) | 一种接收站误差下多站多外辐射源雷达运动目标定位方法 | |
CN108089147A (zh) | 改进的短波单位定位方法 | |
CN110286396A (zh) | 一种基于电离层延迟先验信息和时空变化信息双约束的非差非组合ppp方法 | |
Ren et al. | A Novel Indoor Positioning Algorithm for Wireless Sensor Network Based on Received Signal Strength Indicator Filtering and Improved Taylor Series Expansion. | |
CN115524661A (zh) | 电离层高度与目标位置联合优化的短波时差定位方法 | |
CN111198387A (zh) | 一种抗欺骗干扰的空时采样导航定位方法 | |
CN104076324A (zh) | 一种未知信源数高精度波达方向估计方法 | |
CN104105049A (zh) | 一种减少传声器使用数量的房间冲激响应函数测量方法 |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
C06 | Publication | ||
PB01 | Publication | ||
C10 | Entry into substantive examination | ||
SE01 | Entry into force of request for substantive examination | ||
GR01 | Patent grant | ||
GR01 | Patent grant |