多径干扰是全球导航卫星系统(global navigation satellite system, GNSS)的主要误差来源之一。在诸如测量型接收机、参考站、监测站等高精度应用中,由于其他误差被更加有效地抑制,多径干扰逐渐成为最主要的误差来源。因此,对于此类应用,针对站点的多径干扰的评估对于其他利用该站点差分信息的接收机意义重大。近些年来,随着北斗导航卫星系统(BeiDou navigation satellite system, BDS)的不断发展,地球同步轨道(geosynchronous earth orbit, GEO)卫星在北斗导航卫星系统中扮演着重要角色,同时也给多径效应带来了新特征。针对GEO卫星的多径误差特征的观测最早出现在星基增强系统中,Wanninger等在对欧洲地球静止导航重叠服务(European geostationary navigation overlay service, EGNOS)系统的GEO卫星的载波相位观测中,发现一个变化缓慢的由地面多径导致的误差[1]。在对广域增强系统(wide area augmentation system, WAAS)GEO卫星的伪距观测中,研究者们也发现了由多径导致的固定偏差,该偏差无法消除,被称为Standing Multipath现象[2],这是由于GEO卫星相对于地面观测点静止不动,因此多径导致的偏差几乎不变。北斗系统除了具有上述GEO卫星的多径特征以外,其多径变化周期还呈现时快时慢的特征[3]。高扬等通过分析卫星轨道参数和多径衰落频率之间的关系,提出了轨道多径衰落因子,成功地从理论上解释了GEO卫星的多径误差变化特征的机理[4]。
虽然GEO卫星产生的慢变多径的机理得到了解释,但现有的多径检测技术在运用到GEO卫星时性能都有所下降。常用的多径检测技术包含两大类:1)以多径信号参数估计为基础的估计方法,其代表为多径估计延迟锁定环(multipath estimating delay lock loop, MEDLL)[5]和多径消除技术(multipath mitigation technique, MMT)[6]。这两种方法都基于最大似然估计对接收信号中的直达信号和多径信号的参数进行估计,从而检测多径信号的存在。该类方法虽然对多径的慢变不敏感,但其性能受到多径模型、信号的极限分辨率等估计理论性能上限的制约,对短时延多径的检测能力十分有限。2)基于测量域的多径检测技术最为常用。该技术利用受多径影响的观测量的特征来检测多径。基于载噪比(carrier to noise ratio, CNR)起伏的多径检测方法利用接收信号的幅度受多径干扰后会有波动的现象来检测多径[7]。然而,对于GEO卫星,由于其多径相对于直达信号的变化十分缓慢,起伏周期较长,因此不易检测慢变多径。基于码减载波(code minus carrier, CMC)观测量的技术利用多径对伪距和载波观测量的影响不同,通过二者之间差值的波动来检测多径[8]。然而,若多径衰落频率过低,二者的差值会引入电离层时延的影响,从而导致针对GEO卫星的多径检测性能下降。基于双频伪距、载波相位观测量的方法[9-11]虽然可以消除电离层的影响,但其检测原理仍依赖于多径偏差在一段时间内的平均偏差为零,而GEO卫星的慢变多径使其检测性能受限,并且无法实时反映多径干扰的情况。
综上所述,传统的测量域多径检测技术都是根据中地球轨道(medium earth orbit, MEO)卫星的多径变化特征设计的,而对于GEO卫星多径的慢变特征,尚无有效的、基于测量域的多径检测方法。
上述GEO卫星多径检测困难问题将直接影响地基增强系统中参考站的选址;同时,对于已经受到GEO卫星多径干扰的参考站,其对多径干扰的检测能力有限。此外,现有的测量域多径检测方法受限于GEO卫星的多径衰落频率,其检测周期较长,无法及时反映当前时段多径干扰的强度。另外,提升对GEO卫星多径的检测能力有助于提升对抗多径天线、算法的性能评估。
针对GEO卫星多径检测问题,本文考虑卫星和接收机天线的运动,推导了多径衰落频率和接收机锁相环路相关器输出值之间的关系,并提出了基于由接收机天线圆周运动带来的接收信号幅度振荡的多径检测方法。由于天线运动所带来的多径衰落频率远大于卫星运动所导致的,因此接收信号的多径衰落频率将完全由天线的运动决定。通过控制天线运动的参数可以产生一个特殊的多径衰落频率,进而利用在接收机环路相关器输出幅度的频谱中出现的明显特征,快速、灵敏地检测GEO卫星的多径干扰。
1 多径衰落频率与多径环境多径衰落频率由多径信号相对于直达信号的时延变化或者载波相对相位的变化来确定,可以表示为[12]
$ {f_{{\rm{fade}}}}\left( t \right) = \frac{1}{{2\pi }}\frac{{{\rm{d}}\left[{\Delta \varphi \left( \mathit{t} \right)} \right]}}{{{\rm{d}}\mathit{t}}} = \frac{1}{\lambda }\frac{{{\rm{d}}\left[{\Delta L\left( \mathit{t} \right)} \right]}}{{{\rm{d}}\mathit{t}}}. $ | (1) |
其中:ΔL(t)为多径相对于直达信号的额外传播距离,Δφ(t)为由额外传播距离带来的相位偏差,λ为信号的载波波长。ΔL(t)与接收机天线、反射体与卫星的相对位置有关,本文利用文[12]中的任意反射体位置模型对多径衰落频率进行分析。
图 1显示了卫星、反射体(面)、接收机天线三者的位置p(t)以及速度v(t)。在反射体是散射点的情况下,反射体是静止的;而当反射体是平面时,反射点将会在反射面上移动,因此本文引入反射体的运动速度vd(t)来统一静态散射点模型和反射面模型。
$ \Delta L\left( t \right) = \left| {\mathit{\boldsymbol{m}}\left( t \right)} \right|\cdot\left[{1-\frac{{\mathit{\boldsymbol{m}}{{\left( t \right)}^{\rm{T}}}\mathit{\boldsymbol{k}}\left( t \right)}}{{\left| {\mathit{\boldsymbol{m}}\left( t \right)} \right|\left| {\mathit{\boldsymbol{k}}\left( t \right)} \right|}}} \right]. $ | (2) |
其中:m(t)=pu(t)-pd(t)为由反射体位置pd(t)指向用户位置pu(t)的多径向量,k(t)=ps(t)-pu(t)为由卫星位置ps(t)指向用户位置的直达信号向量。将式(2)代入式(1)中,可以得到多径衰落频率和模型中反射场景的关系,
$ \begin{array}{l} {{\mathit{f}}_{{\rm{fade}}}} = \frac{1}{\lambda }\left[{\frac{{\mathit{\boldsymbol{m}}{{\left( t \right)}^{\rm{T}}}}}{{\left| {\mathit{\boldsymbol{m}}\left( t \right)} \right|}}-\frac{{\mathit{\boldsymbol{k}}{{\left( t \right)}^{\rm{T}}}}}{{\left| {\mathit{\boldsymbol{k}}\left( t \right)} \right|}}} \right]\cdot\left[{\mathit{\boldsymbol{v}}_{\rm{u}}}\left( t \right) - {\mathit{\boldsymbol{v}}_{\rm{d}}}\left( t \right)\right] - \\ \;\;\;\;\;\;\;\;\;\frac{1}{\lambda }\mathit{\boldsymbol{m}}{\left( t \right)^{\rm{T}}}\cdot{\mathit{\boldsymbol{P}}_{ \bot k(t)}}\cdot\frac{{\left[{{\mathit{\boldsymbol{v}}_{\rm{s}}}\left( t \right)-{\mathit{\boldsymbol{v}}_{\rm{u}}}\left( t \right)} \right]}}{{\left| {\mathit{\boldsymbol{k}}\left( t \right)} \right|}}. \end{array} $ | (3) |
其中:
$ {\mathit{\boldsymbol{P}}_{ \bot k(t)}} = \mathit{\boldsymbol{I}} - \frac{{\mathit{\boldsymbol{k}}\left( t \right)\mathit{\boldsymbol{k}}{{\left( t \right)}^{\rm{T}}}}}{{\mathit{\boldsymbol{k}}{{\left( t \right)}^{\rm{T}}}\mathit{\boldsymbol{k}}\left( t \right)}} $ | (4) |
为投影至k(t)的正交空间的投影矩阵,vu(t)和vs(t)分别为用户和卫星的运动速度。式(3)等式右侧的第1项为由用户运动引起的多径衰落,第2项为由卫星运动引起的多径衰落。对于GEO卫星,vs(t)通常很小而距离矢量k(t)很大,因此衰落频率的第1项远大于第2项,第2项可以忽略。式(4)可简化为式(5),以表示由用户运动所带来的多径衰落,
$ {{\mathit{f}}_{{\rm{fade}}}} = \frac{1}{\lambda }\left[{\frac{{\mathit{\boldsymbol{m}}{{\left( t \right)}^{\rm{T}}}}}{{\left| {\mathit{\boldsymbol{m}}\left( t \right)} \right|}}-\frac{{\mathit{\boldsymbol{k}}{{\left( t \right)}^{\rm{T}}}}}{{\left| {\mathit{\boldsymbol{k}}\left( t \right)} \right|}}} \right]\cdot\left[{\mathit{\boldsymbol{v}}_{\rm{u}}}\left( t \right) - {\mathit{\boldsymbol{v}}_{\rm{d}}}\left( t \right)\right]. $ | (5) |
当接收机天线按图 1中的虚线所示的方式以转动角速度为ω、半径为r做圆周运动时,遍历任意的直达和多径方向矢量k(t)及m(t),可以得到给定运动参数下的多径衰落频率的上界,
$ {f_{{\rm{up}}}} = \frac{{2r\omega }}{\lambda }. $ | (6) |
该值仅与圆周运动参数相关,是天线圆周运动下多径衰落频率的重要参数。
2 天线运动与相关器幅度振荡由式(3)可知,当天线和反射体均静止时,多径衰落频率完全由卫星运动决定, 对于GEO卫星,其多径衰落频率极小,该低衰落频率的多径干扰对码相位测量值带来的偏差如图 2中的虚线所示,表现为幅度高、变化周期长。天线运动将大大提升多径衰落频率,此时多径带来的码相位偏差将受到码环滤波器的影响,如图 2中的实线所示,其幅度、变化周期均减小。该原理也被用于多径抑制[13]。
不同于利用转动天线来消除多径偏差,其设定运动参数使多径衰落频率高于跟踪环路带宽[13],本文所提出的多径检测方法必须保证接收机锁相环路(phase lock loop, PLL)带宽大于式(6)所示的最大衰落频率,从而保证环路能完全跟上由多径衰落带来的接收机信号Doppler变化。在此,本文使用常用的Costas锁相环跟踪信号,用一个对准支路相关器的输出(IP(t))作为PLL的输入[14]。基于含有一条直达路径和一条非直达路径的双径模型,接收信号经下变频和相关(匹配接收)之后,相关器的输出可以表示为
$ \begin{array}{l} {I_{\rm{P}}}\left( t \right) + {\rm{j}}{{Q}_{\rm{P}}}\left( t \right) = {A_0}{{\rm{e}}^{{\rm{j}}{\varepsilon _{{\varphi _0}}}}}R(t - {\varepsilon _{{\tau _0}}}) + \\ \;\;\;{A_1}{{\rm{e}}^{{\rm{j}}({\varepsilon _{{\varphi _0}}} + \Delta \varphi \left( t \right))}}R\left( {t - {\varepsilon _{{\tau _0}}} - \Delta \tau \left( t \right)} \right). \end{array} $ | (7) |
其中:A0和A1分别为直达信号和多径信号的幅度;R(t)为归一化互相关函数;εφ0和ετ0分别为接收信号和本地生成信号的载波相位和码相位之间的跟踪误差;Δφ(t)和Δτ(t)分别为多径信号相对于直达信号的载波相位和码相位偏差。
当PLL环路锁定时,跟踪误差εφ0和ετ0均较小且可被忽略,因此对准支路同相相关器输出所表示的接收信号幅度为
$ \begin{array}{l} \;\;\;\;\;\;\;\;\;{I_{\rm{P}}}\left( t \right) = {A_0}{\rm{cos}}({\varepsilon _{{\varphi _0}}})R({\varepsilon _{{\tau _0}}}) + \\ \;\;\;\;\;{A_1}{\rm{cos}}({\varepsilon _{{\varphi _0}}} + \Delta \varphi \left( t \right))R({\varepsilon _{{\tau _0}}} - \Delta \tau \left( t \right)) \approx \\ {A_0} + {A_1}{\rm{cos}}\left[{2\pi {f_{{\rm{fade}}}}t + \Delta \varphi \left( 0 \right)} \right]R(\Delta L\left( t \right)/{\lambda _{\rm{c}}}). \end{array} $ | (8) |
其中λc表示伪随机码码片长度。因为λc≫λ,由ffade产生的幅度振荡远高于由R(ΔL(t)/λc)产生的幅度振荡,并且初相位Δφ(0)与振荡频率没有直接关系。因此,相关器输出IP(t)与衰落频率ffade之间的关系可以简化为
$ {I_{\rm{P}}}\left( t \right) = {A_0} + {{A'}_1}{\rm{cos}}\left( {2\pi {f_{{\rm{fade}}}}t} \right). $ | (9) |
其中A′1 =A1·R(ΔL(t)/λc)。多径相对于直达信号的时延越短,R(ΔL(t)/λc)越接近于1,则A′1越大。因此, 多径时延越短,振荡幅度越大。当有多径干扰存在时,IP(t)会出现小于fup的频谱分量;反之,当没有多径干扰存在时,IP(t)表现为A0均值加白噪声的分布。因此,通过控制天线运动参数以及转动平面的角度,可以控制fup使幅度振荡易于在IP(t)的频谱中显示出来,从而检测多径。
3 多径检测器由式(5)可知,随着天线的运动,ffade并非固定不变,其所有可能值均落在0~fmax的“衰落频带”内,其中fmax<fup且由具体的反射场景几何关系决定。天线运动除了在上述衰落频段内会产生频谱峰值以外,还可能在2个频点产生峰值:圆周运动的频率ω/(2π)和2倍的圆周运动频率ω/π。这两个频率都描述了多径衰落频率变化的周期,为圆周运动带来的固有频率。图 3为一个典型情况下相关器输出幅度振荡的时域图和频谱的示例。参数fup、fmax、ω/(2π)、ω/π以及衰落频带之间的关系在图 3中已标出;信号多径功率比(signal to multipath ratio, SMR)为20 dB,ω=2π/10 rad/s,r=1 m。
利用相关器输出幅度振荡来检测多径干扰比基于CNR和CMC的方法更有优势。图 3b中的虚线为多径不存且未有其他干扰时相关器输出幅度的频谱,为白噪声;而相同场景下基于CNR和CMC的检测量均为非白噪声,会为检测带来额外干扰,即使接收机是静止的。
从图 3b中可以看出,拥有较高能量的频谱分量集中在fmax附近,因此可以通过减小转动半径r来使fmax和固有频率重叠,形成谐振以获得更高的谱峰。但是,对于那些fmax已经小于ω/(2π)的多径,增大转动半径才能使fmax落入衰落频带从而提升可检测的多径的范围。事实上,增大转动半径有助于提升对衰落频率低的多径的检测性解,但会减弱对衰落频率高的多径的检测性解,反之亦然。因此,本文建议最小的转动半径必须满足fup≥2ω/π也即r≥λ/π=0.06 m。转动角频率的选择要保证由运动带来的多径衰落频率能明显区别于其他机制带来的,同时fup小于PLL环路带宽。以下为本文所提出的多径检测器的实现步骤:
步骤1 对以下相关器的输出的变换应用快速Fourier变换(fast Fourier transformation, FFT):
步骤2 在频率从ω/(2π)至fup=2rω/λ的衰落频带内找到频谱的峰值,其幅度记为AS。
步骤3 在步骤2中所述的衰落频带以外的噪声频谱中选定噪声门限并记为AN。AN的选择应满足大于99%的噪声频带的频谱功率同时小于余下的1%。
步骤4 如果检测量Γ=AS/AN≥γ,γ=1.5,可判定存在多径干扰。
4 仿真本文使用Spirent GSS9000 GNSS信号发生器来设定运动多径场景。如图 1所示,假定有一个无穷大的与xoz平面平行的反射面,其到原点的距离为d。仿真中的两个场景的参数如表 1所示。场景A为短时延多径场景,天线的运动使其多径衰弱频率较低;场景B为极短时延多径场景,天线运动使其多径衰落频率较高。在两个场景中,ω=2π/10 rad/s,r=1 m,SMR为30 dB,CNR为50 dB·Hz, θAzm=90°。PLL环路带宽为10 Hz,环路积分时间为10 ms。用于多径检测的数据长度为20 s。
场景 | θElv/(°) | d/m | fmax/Hz | |
A | 85 | 160 | 93.0 | 0.576 |
B | 45 | 3 | 14.2 | 4.670 |
图 4a和4b分别展示了两个场景下的检测结果。可以看到,本文提出的检测方法能很好地检测多径干扰。在衰落频率更低从而使fmax更接近于转动频率的场景A中,检测器的性能更优,这也验证了第3节中设计检测器时的分析。
基于表 1所示的两个多径场景,本文也测试了基于CNR和CMC的多径检测方法。基于CNR的检测方法虽然和本文所提的基于幅度振荡的方法有相类似的地方,但在这两个场景中,由于多径信号强度过低,CNR方法都无法有效检测到多径。基于CMC的方法无法在场景B中检测到极短时延多径。在场景A中,如图 5a所示,CMC检测量虽然有一定周期起伏,但和无多径情况因噪声带来的起伏相似,二者差异较小;从图 5b频域上的结果比较也能看出,当多径存在时,低频部分的能量有所增强,但依然和无多径时区分较小。因此,基于CMC的检测方法不论在时域上还是频域上,其短时延多径检测效果都不如本文方法。在由运动造成的多径衰落频率一致时,随着多径时延的增加,如式(9)所示,基于幅度振荡的检测方法性能逐渐下降,而基于CMC的方法通过增大相关器间隔,性能逐渐提升。对于中长时延的多径(相对多径时延大于500 ns或150 m),基于CMC的检测方法开始优于幅度振荡方法,因此本文所提出的基于幅度振荡的检测量可以代替CNR检测量,与CMC多径检测方法互为补充。此外,通过天线运动,本文解决了GEO卫星慢变多径检测困难问题并将多径的检测时间缩短至20 s。本文方法同样适用于其他轨道类型的卫星的多径检测。
5 结论
GEO卫星慢变多径的特征给多径检测带来了诸多挑战。传统的基于测量域的多径检测技术大多基于MEO卫星的多径衰落特征设计,在应用于GEO卫星多径检测时性能下降。针对上述问题,本文分析了多径衰落与天线运动之间的关系,进而得到天线圆周运动带来的多径衰落导致的接收信号幅度振荡与衰落频率之间的关系。通过控制天线转动参数,可以大幅提升原来GEO卫星慢变多径的多径衰落频率,从而使接收到的受多径干扰的信号的幅度振荡在频谱上易于分辨。基于上述分析,本文提出了一种基于接收机天线圆周运动的利用接收信号幅度振荡的多径检测方法。该方法能提升接收机对GEO卫星慢变多径干扰的检测性能,同时还能缩短检测时间。最后,基于GNSS信号源的实验验证了本文所提方法的有效性。此方法能有效检测出低于直达信号30 dB的短时延多径干扰,且其性能优于基于载噪比(CNR)和码减载波(CMC)的多径检测方法。本文方法大大地缩短了多径检测时间,解决了GEO卫星多径检测的问题。
[1] | WANNINGER L, WALLSTAB-FREITAG S. Combined processing of GPS, GLONASS, and SBAS code phase and carrier phase measurements[C]//Proceedings of the 20th International Technical Meeting of the Satellite Division of the Institute of Navigation (ION GNSS 2007). Fort Worth, USA, 2007: 866-875. |
[2] | SCHEMPP T, BURKE J, RUBIN A. WAAS benefits of GEO ranging[C]//Proceedings of the 21st International Technical Meeting of the Satellite Division of the Institute of Navigation (ION GNSS 2008). Savannah, USA, 2008: 1903-1910. |
[3] |
冯晓超, 金国平, 范建军, 等.
GNSS接收机伪距测量中的多径效应试验分析[J]. 现代电子技术, 2013, 36(5): 77–81.
FENG X C, JIN G P, FAN J J, et al. Experimentation and analysis of multipath effect in pseudo-range measurement of GNSS receiver[J]. Modern Electronics Technique, 2013, 36(5): 77–81. (in Chinese) |
[4] | GAO Y, YAO Z, CUI X, et al. Analysing the orbit influence on multipath fading in global navigation satellite systems[J]. IET Radar, Sonar and Navigation, 2014, 8(1): 65–70. |
[5] | SOKHANDAN N, CURRAN J T, BROUMANDAN A, et al. An advanced GNSS code multipath detection and estimation algorithm[J]. GPS Solutions, 2016, 20(4): 627–640. DOI:10.1007/s10291-015-0475-z |
[6] | WEILL L R. Achieving theoretical bounds for receiver-based multipath mitigation using Galileo OS signals[C]//ION GNSS 19th International Technical Meeting of the Satellite Division. Fort Worth, USA, 2006: 1035-1047. |
[7] | STRODE P R, GROVES P D. GNSS multipath detection using three-frequency signal-to-noise measurements[J]. GPS Solutions, 2016, 20(3): 399–412. DOI:10.1007/s10291-015-0449-1 |
[8] | BEITLER A, TOLLKUEHN A, PLATTNER D B. CMCD: Multipath detection for mobile GNSS receivers[C]//Proceedings of the 2015 International Technical Meeting of the Institute of Navigation. Dana Point, USA, 2015: 455-464. |
[9] | SHALLBERG K, SHLOSS P, ALTSHULER E L. WAAS measurement processing, reducing the effects of multipath[C]//Proceedings of the 14th International Technical Meeting of the Satellite Division of the Institute of Navigation (ION GPS 2001). Salt Lake City, USA, 2001: 2334-2340. |
[10] | SHALLBERG K, SHENG F. WAAS measurement processing: Current design and potential improvements[C]//Proceedings of 2008 IEEE/ION Position, Location and Navigation Symposium. Monterey, USA, 2008: 253-262. |
[11] | WU X L, ZHOU J H, WANG G, et al. Multipath error detection and correction for GEO/IGSO satellites[J]. Science China:Physics, Mechanics and Astronomy, 2012, 55(7): 1297–1306. |
[12] | IRSIGLER M. Characterization of multipath phase rates in different multipath environments[J]. GPS Solutions, 2010, 14(4): 305–317. DOI:10.1007/s10291-009-0155-y |
[13] | VAN DEN BREKEL B J H, VAN NEE D J R. GPS multipath mitigation by antenna movement[J]. Electronic Letters, 1992, 28(25): 2286–2288. |
[14] | KAPLAN E D, HEGARTY C J. Understanding GPS:Principles and applications[M]. 2nd ed, Boston, USA: Artech House Publishers, 2006: 164-173. |