BDS偽距單站定位計算及精度分析

干堅定，周 適，段太生，郭 平，晏 勇，劉立正

（1.中鐵二局第五工程有限公司，四川 成都 610091；2.中鐵二局集團有限公司，四川 成都 610031）

北斗衛(wèi)星導航系統(tǒng)（BDS）最后一顆全球組網(wǎng)衛(wèi)星（GEO衛(wèi)星）已于2020年6月23日9：43在西昌衛(wèi)星發(fā)射中心成功發(fā)射，這是由我國自主研發(fā)的衛(wèi)星導航系統(tǒng)。對于全球導航衛(wèi)星系統(tǒng)（GNSS）而言，GPS雖然仍是當前技術最先進和最成熟的導航定位系統(tǒng)，但BDS通過近3年來高密度的衛(wèi)星發(fā)射，目前已具備和GPS共同作用、兼容并存的實力。從衛(wèi)星數(shù)量和衛(wèi)星類型上看，北斗的衛(wèi)星數(shù)量和亞太地區(qū)的衛(wèi)星分布比GPS更具優(yōu)勢，國內大量專業(yè)人士對BDS進行了相關深入研究。BDS即將組網(wǎng)完成，極大地增強了國人的自信心和民族自豪感，也激發(fā)了眾多相關從業(yè)人員學習和研究BDS的興趣。

1 相鄰歷元相位平滑偽距的數(shù)學模型

楊元喜等[1]提出，對北斗所有衛(wèi)星進行統(tǒng)計，北斗B1I和B2I偽距野外測量的精度約33cm，B1和B2載波相位野外測量的精度約2mm。由于載波相位觀測的精度遠高于偽距觀測[2]，因此采用偽距值進行單站定位計算時，有必要對偽距值進行平滑處理。偽距平滑的成熟算法包括載波相位平滑偽距算法和基于卡爾曼濾波的偽距平滑算法[3]。文章采用簡易數(shù)學公式計算，BDS共設為3個頻率，由于外業(yè)使用的接收機只采集B1、B2兩個頻率的觀測數(shù)據(jù)，因此文章中的數(shù)學模型是針對北斗雙頻（B1、B2）展開的。對BDS不同頻率分別列出偽距平滑公式，由于相鄰歷元時間較短，因此偽距中包含的噪聲（電離層延遲、對流層延遲、多路徑效應）可以忽略不計，整周模糊度N也可消去，具體公式如下[4]：

同理，可得B2頻率下的各參數(shù)含義。利用相位平滑后的偽距進行后續(xù)單站定位計算，若監(jiān)測到載波相位發(fā)生整周跳變，則相關歷元不采用相位平滑的偽距值，而采用原始偽距觀測值。

2 BDS單站單歷元定位的數(shù)學模型

在t歷元時刻，接收機共接收到n顆北斗衛(wèi)星，通過北斗廣播星歷文件可計算每顆北斗衛(wèi)星的空間坐標，選取高度角最高的衛(wèi)星作為參考衛(wèi)星r，其余衛(wèi)星記作衛(wèi)星s，將偽距觀測方程線性化后，通過星間單差（其他衛(wèi)星和參考衛(wèi)星作差），可列出B1頻率的偽距觀測方程，此時接收機鐘差可消去，并忽略多路徑效應誤差，具體公式如下[5]：

式中：x0、y0、z0為測站點初始坐標，初始坐標可設為(0,0,0)；(xr,yr,zr)為參考衛(wèi)星r坐標；(xs,ys,zs)為其他衛(wèi)星s坐標。由于t歷元時刻共接收n顆衛(wèi)星，由式（2）可列出n-1個方程。

同理可列出B2頻率的星間單差偽距觀測方程，具體公式如下：

采用不同頻率觀測值的線性無電離層組合，由于北斗B1頻率值為1561.098MHz，B2頻率值為1207.14MHz，可寫出北斗B1和B2頻率的波長比值，具體公式如下：

因此式（2）乘以7632，式（4）乘以-5902，相加后移項可得到如下公式：

式（6）已消去大部分電離層延遲，方程左邊為常數(shù)項，其中對流層延遲，可采用Saastamoinen模型進行改正計算，對流層改正模型如下[6]：

式中：z為衛(wèi)星天頂距，rad；p為大氣壓，hPa；T為開氏溫度，K；e為水氣分壓，hPa。

式中：a0、a1、a2分別為衛(wèi)星鐘的鐘差、鐘速、鐘漂。

式（6）可簡寫為如下公式[8]：

1個歷元接收到n顆北斗衛(wèi)星，可列出n-1個誤差方程，式（9）中，常數(shù)項矩陣為L(n-1)×1，系數(shù)陣為B(n-1)×3，改正數(shù)矩陣為V(n-1)×1，未知數(shù)矩陣為x3×1。構建誤差方程式后，根據(jù)最小二乘原理，可計算未知數(shù)x，具體公式如下[9]：

由于測站初值設為(0,0,0)，通過式（10）計算得出的x還不是最終結果，需進行迭代運算，因不知道迭代次數(shù)，C#程序中可采用while循環(huán)語句，循環(huán)結束的判別標志可令未知數(shù)x＜1cm，根據(jù)實測數(shù)據(jù)驗證，一般情況下迭代次數(shù)不超過5次即可滿足條件。

式（10）的權陣P可通過衛(wèi)星高度角不同進行定權，高度角越大權值越高。式（6）的觀測值為B1和B2頻率下的偽距觀測值，B1、B2兩個頻率的偽距值是獨立的。觀測值組合中誤差可設為σ，其取值與衛(wèi)星高度角有關。觀測值方差陣D可用如下公式表示：

權陣P可通過方差陣D求逆運算得到，代入式（10）可計算未知參數(shù)x。

3 BDS和GPS組合單站單歷元定位的數(shù)學模型

在BDS單衛(wèi)星系統(tǒng)下單歷元單站定位的數(shù)學模型中采用星間單差、雙頻觀測值無電離層組合的模式進行計算。若采用BDS和GPS組合衛(wèi)星系統(tǒng)，則星間單差計算方式有兩種，一種是采用GPS和BDS各系統(tǒng)選擇自己的參考衛(wèi)星分開組合星間單差；另一種是僅選擇一個參考衛(wèi)星組混合單差[10]。考慮到組混合單差需要引入系統(tǒng)間偏差，算法較為復雜，故采用第一種計算方式，即BDS和GPS分別各選擇系統(tǒng)內的一顆衛(wèi)星作為參考衛(wèi)星，按照上述公式進行計算。由于GPS雙頻觀測值（不采用GPS第3個民用頻率L5的數(shù)據(jù)進行計算，因L5觀測值目前部分GPS衛(wèi)星還沒有，僅采用L1、L2兩個頻率）L1頻率值為1575.42MHz，L2頻率值為1227.60MHz，易得GPS兩個頻率的波長關系，77λ1=60λ2。將BDS單站定位計算公式中式（6）、式（11）、式（12）中的系數(shù)763替換成77，590替換成60，即可得到GPS單站定位的公式。

現(xiàn)將兩者組合，設1個歷元共接收n顆北斗衛(wèi)星、m顆GPS衛(wèi)星的觀測數(shù)據(jù)，可按照上述公式寫出組合衛(wèi)星系統(tǒng)的誤差方程式，具體公式如下：

除未知參數(shù)x不變外，其余矩陣的行數(shù)有所擴充，北斗衛(wèi)星觀測值可列出n-1個方程，GPS衛(wèi)星觀測值可列出m-1個方程，組合形成n+m-2個方程。

此時方差陣D可用分塊矩陣表示，具體公式如下：

式中：DBDS為BDS觀測值的方差陣，采用式（12）計算；DGPS為GPS觀測值的方差陣，將式（12）的北斗雙頻組合系數(shù)換成GPS的組合系數(shù)即可；x為未知參數(shù)，按式（10）計算。

4 BDS單站多歷元定位的數(shù)學模型

文章采用單站定位的數(shù)學模型未知參數(shù)為3個，即測站點的未知數(shù)(dx,dy,dz)；而傳統(tǒng)單站定位的數(shù)學模型未知參數(shù)一般有4個，包括接收機鐘差。文章通過星間單差（其他衛(wèi)星和參考衛(wèi)星觀測方程作差）消去了接收機鐘差，不將其作為未知參數(shù)，由于是靜態(tài)觀測，測站點坐標不會改變，而接收機鐘差可能會發(fā)生變化，故這里消去了接收機鐘差，再利用多歷元觀測數(shù)據(jù)計算單站定位的坐標時，便于采用法方程疊加的原理進行計算，具體公式如下：

式中：Xk-1為前k-1個歷元計算的測站點坐標值；Xk為k個歷元計算的測站點坐標值；Bi、Pi、Li分別為第i個歷元所形成的系數(shù)陣、權陣和常數(shù)項矩陣，不同歷元觀測的北斗衛(wèi)星數(shù)量可能不同，因此各歷元的系數(shù)陣、權陣、常數(shù)項矩陣的行數(shù)可能不同。

式（15）是單獨BDS多歷元單站定位的數(shù)學模型，通過單歷元GPS和BDS組合衛(wèi)星系統(tǒng)的單站定位計算模型，易知GPS和BDS組合衛(wèi)星系統(tǒng)下多歷元單站定位的計算公式，亦采用法方程累加的原理，構建的各矩陣由于觀測方程的增加而需相應擴充。

當計算出測站點坐標時，精度指標也可同時計算得出。即可通過相應公式計算單位權中誤差σ0，求出測站在N、E、U方向的點位精度MN、ME、MU，平面精度因子HDOP值，高程精度因子VDOP值及空間精度因子PDOP值等精度指標。點位精度MN、ME、MU表征計算得到的測站點各方向上的精度，DOP值表征衛(wèi)星分布的合理性，DOP值與測站和觀測衛(wèi)星構成的單位多面體的體積成反比，DOP值越小，則相應的解算精度越高[8]。

5 實測數(shù)據(jù)的坐標及精度分析

利用2019年5月在山東菏澤地區(qū)某控制點觀測約70min的靜態(tài)數(shù)據(jù)，采樣間隔10s，高度截止角設為15°，共采集357個歷元數(shù)據(jù)。對每個歷元進行單站定位計算，只采用BDS觀測值，計算得到的測站空間坐標X、Y、Z坐標如圖1～圖3所示。

由圖1～圖3可知，單獨利用BDS偽距觀測值進行單站定位計算，X坐標變化幅度為7.01m，（X最大值為-2291357.8834m，最小值為-2291364.9011m），Y坐標變化幅度為10.21m，（Y最大值為4670024.9080m，最小值為4670014.6962m），Z坐標變化幅度為4.52m（Z最大值為3678370.6212m，最小值為3678366.1014m）。其中Y坐標變化幅度最大，這是因為中國大部分地區(qū)Y坐標（空間直角坐標）主要影響高程方向，衛(wèi)星定位高程方向的精度一般比平面精度低。各歷元測站X坐標值大部分落在-2291360～-2291362m，各歷元測站Y坐標值大部分落在4670017～4670021m，各歷元測站Z坐標值大部分落在3678367～3678369m。結合圖1～圖3分析，測站定位坐標在第150～180個歷元時間段內突變較大，其余時刻測站坐標各分量變化不大。參考RINEX數(shù)據(jù)O文件給出的測站近似坐標，與BDS各歷元計算的單站定位坐標平均值進行較差比較，如表1所示。

由表1可知，O文件提供的近似坐標與利用BDS觀測值數(shù)據(jù)計算各歷元單站定位坐標進行比較，坐標差別不大，各分量較差均在1m以下。

現(xiàn)將單獨利用BDS和單獨利用GPS計算測站各歷元單站定位坐標進行較差比較，如圖4～圖6所示。

圖1 BDS偽距觀測值計算每個歷元測站點X坐標（單位：m）

圖2 BDS偽距觀測值計算每個歷元測站點Y坐標（單位：m）

圖3 BDS偽距觀測值計算每個歷元測站點Z坐標（單位：m）

表1 RINEX數(shù)據(jù)的O文件測站近似坐標 單位：m

由圖4～圖6可知，單獨采用BDS系統(tǒng)和單獨采用GPS系統(tǒng)進行單站定位，各歷元計算測站點的坐標較差差別不大，大部分歷元X坐標較差在2m以下，Y坐標較差在4m以下，Z坐標較差在3m以下。將BDS和GPS各歷元單站定位坐標取平均值進行比較，如表2所示。

由表2可知，在分別利用BDS和GPS系統(tǒng)計算各歷元單站定位坐標平均值時，分量較差均在2m以下。在BDS系統(tǒng)下單站定位在N、E、U方向的點位精度MN、ME、MU如圖7～圖9所示。

由圖7～圖9可知，BDS系統(tǒng)下單站定位精度為米級，大部分歷元的E方向的ME和N方向的MN精度值在1.5m以下。該實例數(shù)據(jù)的E方向點位精度最高，U方向點位精度最低，在第170～180個歷元計算U方向的精度最低，MU出現(xiàn)峰值。從圖1～圖3各歷元單站定位計算的坐標各分量X、Y、Z值也可以看出，在第170～180個歷元區(qū)段內，計算坐標值有明顯差異，這與圖9歷元U方向上的點位精度MU變化特征是相符的。

由于單站定位的精度在米級，文章也只采用了廣播星歷計算衛(wèi)星坐標，再按照數(shù)學模型進行單站定位計算。從BDS和GPS計算的精度指標結果分析，單獨采用BDS和GPS系統(tǒng)進行單站定位精度相當，均能達到米級[11-13]，若采用BDS和GPS混合衛(wèi)星系統(tǒng)進行單站定位計算，由于廣播星歷計算的衛(wèi)星坐標精度約為2m[13]，單站定位精度也只能達到米級，但在個別時段衛(wèi)星數(shù)量較少時進行單站定位解算，此時混合衛(wèi)星系統(tǒng)下單站定位的解算精度能得到更好的保證。若需要提高定位精度，達到厘米級的定位精度[14]，衛(wèi)星星歷需采用精密星歷，觀測值需采用載波相位觀測值[15]，并準確計算整周模糊度，各種誤差分別采用各自的誤差改正模型進行改正計算。

圖4 BDS和GPS計算各歷元測站點X坐標較差（單位：m）

圖5 BDS和GPS計算各歷元測站點Y坐標較差（單位：m）

圖6 BDS和GPS計算各歷元測站點Z坐標較差（單位：m）

表2 BDS和GPS各歷元單站定位坐標的平均值較差比較表 單位：m

圖7 BDS計算各歷元單站定位的N方向的點位精度MN（單位：m）

圖8 BDS計算各歷元單站定位的E方向的點位精度ME（單位：m）

圖9 BDS計算各歷元單站定位的U方向的點位精度MU（單位：m）

6 結論

文章對BDS利用偽距觀測值進行單站定位的數(shù)學模型進行了介紹，并通過C#編程實現(xiàn)計算過程。通過實例數(shù)據(jù)分析，得出以下結論：（1）BDS和GPS單獨采用一種衛(wèi)星系統(tǒng)進行單站定位的精度相當，精度都能達到米級。若采用BDS、GPS混合衛(wèi)星系統(tǒng)進行單站定位計算，由于廣播星歷計算衛(wèi)星坐標的精度約為2m，因此BDS、GPS混合衛(wèi)星系統(tǒng)下單站定位的坐標精度仍為米級。（2）采用文章的數(shù)學模型僅使用BDS進行單站定位計算，將計算結果與RINEX的O文件中提供的近似坐標進行比較，文章實例數(shù)據(jù)計算的坐標較差在1m以下。BDS單站定位計算的點位精度，平面的點位精度MN、ME高于高程方向的點位精度MU。（3）采用多歷元靜態(tài)觀測數(shù)據(jù)進行單站定位計算時，由于數(shù)學模型消去了接收機鐘差未知參數(shù)，僅保留測站的3個未知參數(shù)，在各歷元測站點進行靜態(tài)觀測并保持不動，便于采用法方程累加的原理計算多歷元單站定位的測站坐標。