基于改進(jìn)蝴蝶搜索算法的DGPS 整周模糊度快速解算

尚俊娜，羅照旺

（杭州電子科技大學(xué) 通信工程學(xué)院，杭州 310018)

在差分全球定位系統(tǒng)（Differential Global Positioning System,DGPS）定位方式中，載波相位觀測(cè)值[1]是否能夠達(dá)到厘米級(jí)甚至毫米級(jí)的定位結(jié)果，關(guān)鍵在于能否正確快速地解算出整周模糊度。通常解算整周模糊度的關(guān)鍵步驟為：首先通過(guò)最小二乘估計(jì)得到模糊度的浮點(diǎn)解[2]和協(xié)方差矩陣；其次，利用去相關(guān)操作處理協(xié)方差矩陣，從而降低模糊度間的相關(guān)性；最后通過(guò)搜索算法得到模糊度固定解[3-5]。

近年來(lái)，大量學(xué)者在模糊度解算方法上不斷探索，在提高搜索效率和降低模糊度之間的相關(guān)性方面進(jìn)行了研究，尤其在搜索算法方面采用了現(xiàn)代啟發(fā)式算法。劉新華等[6]提出了一種改進(jìn)蟻群算法（Improved Ant Colony Optimization,IACO），在傳統(tǒng)蟻群算法的基礎(chǔ)上引入自反饋因子，改善了傳統(tǒng)蟻群的尋優(yōu)性能，將整周模糊度解算成功率提升了3%。然而，蟻群算法容易缺乏初始信息素，從而導(dǎo)致模糊度解算時(shí)間較長(zhǎng)。歐陽(yáng)利等[7]提出了一種基于短基線整周模糊度求解的改進(jìn)雞群優(yōu)化算法（Improved Chicken Swarm Optimization,ICSO），對(duì)模糊度解算過(guò)程進(jìn)行優(yōu)化，提高了搜索速度和求解成功率，但是在基線不固定時(shí)或者在長(zhǎng)基線環(huán)境下存在搜索空間變大、搜索效率下降和模糊度搜索速度變慢的問(wèn)題。針對(duì)模糊度解算時(shí)間長(zhǎng)、效率低的問(wèn)題，張波等[8]提出了改進(jìn)粒子群（Improved Particle Swarm Optimization,IPSO）的模糊度解算優(yōu)化算法，通過(guò)引入sin 函數(shù)對(duì)慣性權(quán)重進(jìn)行改進(jìn)，縮短了模糊度解算時(shí)間。上述算法在模糊度解算問(wèn)題上雖然取得了一定的性能提升，但現(xiàn)代啟發(fā)式算法自身所固有的尋優(yōu)精度不高、收斂速度慢、易于早熟收斂等問(wèn)題還有待進(jìn)一步解決，這使得在提高模糊度搜索效率和模糊度解算成功率等問(wèn)題上仍然有進(jìn)一步的優(yōu)化空間。

為了提高模糊度固定解的搜索成功率，本文提出了一種改進(jìn)的蝴蝶搜索算法（Improved Butterfly Optimization Algorithm,IBOA）來(lái)求解整周模糊度。首先，引入動(dòng)態(tài)切換概率來(lái)權(quán)衡蝴蝶搜索算法（Butterfly Optimization Algorithm,BOA）中全局搜索和局部搜索的比例，從而提高全局搜索能力；其次，在BOA 算法的香味系數(shù)中加入一個(gè)自適應(yīng)權(quán)重，在全局搜索和局部搜索階段引入新的位置更新策略，不僅提升了全局尋優(yōu)精度，還改善了跳出局部最優(yōu)的能力；最后，通過(guò)實(shí)驗(yàn)對(duì)比IBOA 算法與LAMBDA 算法解算的整周模糊度。

1 載波相位模型分析

載波相位觀測(cè)量是通過(guò)接收機(jī)捕獲到的相位與自身生成的相位之差得到的，載波相位觀測(cè)量的觀測(cè)精度比偽距觀測(cè)量更高，但其中包含著未知的整周模糊度，如何快速準(zhǔn)確地求解這個(gè)未知數(shù)是實(shí)現(xiàn)高精度載波相位相對(duì)定位的關(guān)鍵所在。

1.1 模糊度最小二乘估計(jì)模型

在實(shí)現(xiàn)高精度載波相位相對(duì)定位時(shí)，通常對(duì)站星間載波相位雙差觀測(cè)方程進(jìn)行處理[9]。假設(shè)接收機(jī)q和r同步觀測(cè)衛(wèi)星s和k，那么在接收機(jī)q和r之間求得單差觀測(cè)值后，再在衛(wèi)星s和k之間作差求得雙差觀測(cè)值，其方程為：

將式(1)進(jìn)行載波相位雙差模型線性化，其數(shù)學(xué)模型為：

其中，y為載波相位觀測(cè)值；N為整周模糊度向量；b為兩接收機(jī)之間的基線向量；A、B均為設(shè)計(jì)矩陣；ε為觀測(cè)噪聲；D(y)和Qy分別為觀測(cè)量y的方差陣和協(xié)方差陣，是其單位先驗(yàn)方差。

在GPS 領(lǐng)域，一般采用最小二乘方法或Kalman濾波方法[10]求解式(2)。由于式(2)中不僅包含未知的實(shí)數(shù)參數(shù)，還包含未知的整數(shù)信息，因此，式(2)被稱為混合整數(shù)模型[11]。通過(guò)最小二乘估計(jì)準(zhǔn)則求解式(2)，有：

如果先不考慮整周模糊度N的整數(shù)特性，利用標(biāo)準(zhǔn)最小二乘方法來(lái)解算，可以得到整周模糊度N和基線向量b各自的浮點(diǎn)解及其協(xié)方差矩陣和：

在解算出N和b的浮點(diǎn)解后，再考慮整周模糊度的整數(shù)特性，式(3)可以進(jìn)一步展開：

由于b為實(shí)數(shù)向量，所以當(dāng)，且取最小值時(shí)取得最小值。因此，式(3)的最小問(wèn)題就轉(zhuǎn)化為：

1.2 模糊度降相關(guān)

由式(4)得到的N和b的浮點(diǎn)解存在非常強(qiáng)的相關(guān)性，若不采用降相關(guān)方法處理就直接進(jìn)行整周模糊度的搜索，會(huì)造成模糊度搜索時(shí)間過(guò)長(zhǎng)。為了獲得更好的搜索效率，在搜索前需要先進(jìn)行模糊度降相關(guān)操作。本文采用Z 變化[12]的降相關(guān)處理方法：

其中，z和Q分別為N和經(jīng)整數(shù)變化后的向量，為經(jīng)整數(shù)變化后的浮點(diǎn)解。注意的是，為保證N和z之間的變化是整數(shù)變換，z? Zn×1,Z?Zn×n，Z和Z-1中的元素均為整數(shù)。

式(6)轉(zhuǎn)換成一個(gè)新的整數(shù)最小二乘問(wèn)題：

通過(guò)Z 變化降低協(xié)方差矩陣之間的相關(guān)性，這樣能夠減少一些搜索時(shí)間，由于Z 變化自身的運(yùn)算量增加了模糊度解算的整體計(jì)算時(shí)間，因此，想進(jìn)一步提升搜索效率，需要使用一個(gè)高效的搜索算法。

2 IBOA 算法求解整周模糊度

模糊度搜索的本質(zhì)是在搜索空間中找到J(z)的最小值。為了進(jìn)一步提升模糊度的搜索效率，本文提出了一種IBOA 算法進(jìn)行模糊度搜索。BOA 算法是一種模擬覓食行為的新型優(yōu)化技術(shù)[13-15]，蝶群中的每個(gè)蝴蝶通過(guò)感知其他蝴蝶散發(fā)出的香味進(jìn)行位置更新，并且每個(gè)蝴蝶的位置是解空間中的一個(gè)解。

2.1 BOA 算法的基本原理

2.1.1 初始化種群位置

假設(shè)在一個(gè)D維搜索空間中，存在n只蝴蝶，第i只蝴蝶在D維空間中的位置為：

其中，1≤i≤n，1≤d≤D，xi,d是第i只蝴蝶在d維空間的位置。

為了更好地解算整周模糊度，搜索空間的邊界條件可以設(shè)置為：

其中，Xmin和Xmax分別表示搜索空間的下限和上限。

蝴蝶Xi的第j維參數(shù)的初始化公式為：

其中，R是[0,1]之間的隨機(jī)數(shù)。

蝴蝶在進(jìn)行位置更新時(shí)，需要保證位置更新后的向量元素處于預(yù)先設(shè)定的搜索空間之內(nèi)。若超出范圍，則需要進(jìn)行修正。修正方法為：

2.1.2 更新蝴蝶位置

蝴蝶位置更新分為兩個(gè)階段：在全局搜索階段，蝴蝶個(gè)體會(huì)朝著產(chǎn)生香味最濃的蝴蝶進(jìn)行位置更新；在局部搜索階段，當(dāng)蝴蝶個(gè)體感知不到其它蝴蝶的香味時(shí)，就會(huì)隨機(jī)移動(dòng)。在迭代階段，為了控制全局搜索和局部搜索的占有比例，傳統(tǒng)的BOA 算法采用一個(gè)切換概率p=0.8的固定值與每次迭代產(chǎn)生的一個(gè)隨機(jī)數(shù)rand進(jìn)行判斷。如果rand＜p，則進(jìn)行全局搜索；反之則進(jìn)行局部搜索。蝴蝶的位置更新為：

其中，c表示感覺(jué)模態(tài)；α表示功率指數(shù)，反映蝴蝶香味吸收程度的變化。α和c的取值范圍都在[0,1]之間。I表示刺激強(qiáng)度，在整周模糊度搜索[14,15]中將I替換成目標(biāo)函數(shù)J(z)，其中z用蝴蝶位置Xi表示：

2.1.3 搜索空間的構(gòu)建

BOA 算法尋優(yōu)本質(zhì)是在備選空間中搜索能夠使目標(biāo)函數(shù)J(Xi)取得最小值的一組解向量，所以搜索空間中必須包含正確的模糊度整數(shù)解。因此，在DGPS基線定位中可以利用基線長(zhǎng)度構(gòu)造搜索空間[16,17]。搜索空間為：

基線長(zhǎng)度可以約束搜索空間，通過(guò)多次測(cè)量取平均值來(lái)減小基線長(zhǎng)度的測(cè)量誤差。例如1.0 m 長(zhǎng)的基線，L1 載波波長(zhǎng)為19.03 cm，則。其中，l為基線長(zhǎng)度，[l/λ]表示向上取整。因此，邊界條件確定為xmax,j=[zj]+6，xmin,j=[zj] -6。通過(guò)基線長(zhǎng)度來(lái)約束搜索空間，一般基線長(zhǎng)度測(cè)量誤差比較小，如果出現(xiàn)基線長(zhǎng)度測(cè)量誤差大的情況，會(huì)造成搜索空間變大和候選值冗余的問(wèn)題。因此，下文提出的IBOA算法搜索速度較快，即使搜索空間變大了，對(duì)解算效率影響很小。

2.2 IBOA 算法

在整周模糊度解算問(wèn)題上傳統(tǒng)的BOA 算法還存在很多問(wèn)題，比如收斂速度慢和容易陷入局部最優(yōu)值等。為了解決這些問(wèn)題，本文從搜索策略上對(duì)傳統(tǒng)BOA 算法進(jìn)行改進(jìn)，提出一種IBOA 算法。該算法引入非線性遞減函數(shù)P改進(jìn)切換概率，有效地平衡了全局搜索和局部搜索的比例，提升了算法的收斂性；引入自適應(yīng)權(quán)重S對(duì)香味公式進(jìn)行改進(jìn)，增強(qiáng)算法的尋優(yōu)能力。

2.2.1 動(dòng)態(tài)切換概率策略

在BOA 算法中，采用了一個(gè)固定值切換概率p與每次迭代產(chǎn)生的一個(gè)隨機(jī)數(shù)進(jìn)行比較，從而實(shí)現(xiàn)全局搜索和局部搜索之間的切換比例。如果p的取值越大，進(jìn)行全局搜索的概率就越大；反之，局部搜索的概率就越大。為了更有效地平衡全局搜索和局部搜索的比例，本文基于cos 函數(shù)的權(quán)重遞減方法[18]對(duì)BOA算法的開關(guān)概率進(jìn)行改進(jìn)，改進(jìn)后的開關(guān)概率P是一個(gè)遞減函數(shù)，前期函數(shù)下降速度慢，后期函數(shù)下降速度快，保證了前期全局搜索可以充分進(jìn)行，后期局部搜索也可以更好地二次精細(xì)尋優(yōu)。

其中，t和Tmax分別為當(dāng)前迭代次數(shù)和最大迭代次數(shù)。

2.2.2 自適應(yīng)權(quán)重和新的位置更新策略

為了減小BOA 中香味系數(shù)fi,d的隨機(jī)性和盲目性[19]，在香味系數(shù)fi,d中加入自適應(yīng)權(quán)重S，使香味分布范圍越來(lái)越小，算法能夠越來(lái)越仔細(xì)地搜索。當(dāng)搜索過(guò)程開始時(shí)，初始香味分布范圍很大。隨著迭代過(guò)程的進(jìn)行，香味分布范圍逐漸減小，保證了收斂速度和跳出局部極值點(diǎn)的能力。

新的香味系數(shù)fnew：

理想的尋優(yōu)過(guò)程包括前期的高搜索能力和后期的高精度能力。對(duì)于全局搜索階段，IBOA 算法通過(guò)調(diào)整全局最優(yōu)位置g*與當(dāng)前迭代的每個(gè)蝴蝶位置之間的差異來(lái)進(jìn)行位置更新。因此，IBOA 算法使每個(gè)個(gè)體圍繞當(dāng)前迭代的最優(yōu)位置有一個(gè)更加仔細(xì)和徹底的搜索機(jī)制。全局搜索階段新的位置更新公式為：

對(duì)于局部搜索階段，IBOA 算法使用每次迭代中最優(yōu)位置g*和最差位置w*的平均值來(lái)更新每只蝴蝶的下一個(gè)位置，減少隨機(jī)選擇兩只蝴蝶的盲目性。新的局部搜索階段的位置更新公式為：

由于整周模糊度具有整數(shù)特性，IBOA 算法的位置更新公式需要進(jìn)行四舍五入：

其中，round(?)是四舍五入函數(shù)。

2.2.3 IBOA 算法的流程圖

IBOA 算法的流程如圖1 所示，輸入值為經(jīng)過(guò)Z變化后的浮點(diǎn)解和協(xié)方差矩陣Q，通過(guò)IBOA 算法搜索后輸出整周模糊度固定解。

圖1 IBOA 算法搜索流程圖Fig.1 Flow chart of the search of the IBOA algorithm

3 實(shí)驗(yàn)與結(jié)果分析

3.1 仿真實(shí)驗(yàn)

為了驗(yàn)證IBOA 算法解算整周模糊度的性能，利用最小二乘估計(jì)得到模糊度浮點(diǎn)解和協(xié)方差矩陣，如式(25)所示，再根據(jù)經(jīng)典算例[20]進(jìn)行仿真處理。

為了有更好的搜索效率，在搜索前先進(jìn)行降相關(guān)處理和矩陣變化，即對(duì)進(jìn)行 Z 變化（z=ZTN），運(yùn)算步驟及結(jié)果為和經(jīng)過(guò)處理后的浮點(diǎn)解和協(xié)方差為：

為了驗(yàn)證IBOA 算法的有效性，與IACO 算法、ICSO 算法、IPSO 算法和傳統(tǒng)BOA 算法進(jìn)行對(duì)比。將降相關(guān)處理后的數(shù)據(jù)進(jìn)行整周模糊度解算，五種算法的種群數(shù)量都為80，最大迭代次數(shù)Tmax=30。表1 為五種算法的參數(shù)設(shè)置，圖2 是以J(Xi)為目標(biāo)函數(shù)的五種算法解算整周模糊度的收斂過(guò)程，表2 是五種算法的最優(yōu)解演變過(guò)程。

表1 五種算法參數(shù)設(shè)置Tab.1 Parameter setting of five algorithms

表2 最優(yōu)解的演變Tab.2 Evolution of the optimal solution

圖2 五種算法的收斂曲線Fig.2 Convergence curves of five algorithms

從圖2 可知，在解算整周模糊度的收斂過(guò)程中，BOA 算法在第9 次迭代時(shí)進(jìn)入局部最優(yōu)，在第15 次迭代時(shí)才跳出，最后仍然沒(méi)有達(dá)到最優(yōu)解。IPSO 算法在第18 次迭代時(shí)達(dá)到最優(yōu)解，目標(biāo)函數(shù)值J(Xi)達(dá)到最小，不再變化；ICSO 算法在第19 次迭代時(shí)達(dá)到最優(yōu)解，目標(biāo)函數(shù)值J(Xi)達(dá)到最小，不再變化；IACO 在第15 次迭代時(shí)達(dá)到最優(yōu)解，目標(biāo)函數(shù)值J(Xi)達(dá)到最小，不再變化；IBOA 算法在第5 次迭代時(shí)就達(dá)到最優(yōu)解，目標(biāo)函數(shù)值J(Xi)達(dá)到最小，且不再發(fā)生變化。IPSO、ICSO、IACO 和IBOA 算法都能夠?qū)さ阶顑?yōu)解，但是IBOA 收斂效果比其他三種算法要好。此外，IBOA 算法具有不錯(cuò)的全局尋優(yōu)能力，解決了BOA 算法無(wú)法跳出局部最優(yōu)的問(wèn)題，進(jìn)一步提升了整周模糊度解算的效率。表2 的最優(yōu)解演變過(guò)程證明了IBOA 算法求解整周模糊度固定解的有效性。

3.2 實(shí)驗(yàn)分析

本文實(shí)驗(yàn)分為兩次測(cè)試：第一次測(cè)試使用武漢大學(xué)IGS（International GPS service for Geodynamics）數(shù)據(jù)中心下載的觀測(cè)值文件進(jìn)行測(cè)試，第二次測(cè)試是在蘇州市科技城中科地理信息產(chǎn)業(yè)園進(jìn)行。由于兩次測(cè)試的數(shù)據(jù)量不同，為了方便進(jìn)行比較，統(tǒng)一從各自數(shù)據(jù)量中按順序選取其中1000 個(gè)歷元進(jìn)行解算。

武漢大學(xué)IGS 數(shù)據(jù)中心是由國(guó)際地球動(dòng)力學(xué)服務(wù)機(jī)構(gòu)組織授權(quán)建立，可免費(fèi)下載國(guó)際IGS 基準(zhǔn)站觀測(cè)數(shù)據(jù)、GNSS 衛(wèi)星精密軌道和精密誤差等數(shù)據(jù)服務(wù)。本次測(cè)試數(shù)據(jù)來(lái)自香港基站hkws 和hksl，其中hkws基站作為流動(dòng)站，hksl 基站作為基準(zhǔn)站，基線長(zhǎng)度為42.7 km。下載兩個(gè)基站同一天的數(shù)據(jù)文件，使用RTKLIB（Real-Time Kinematic Library）程序進(jìn)行數(shù)據(jù)解析。采用傳統(tǒng)LAMBDA、BOA、IPSO、ICSO、IACO和IBOA 算法進(jìn)行整周模糊度解算，一天的數(shù)據(jù)總計(jì)2347 歷元，選取其中1000 個(gè)歷元進(jìn)行解算，對(duì)六種算法的成功率和解算時(shí)間進(jìn)行對(duì)比，如表3 所示。

表3 六種算法解算結(jié)果Tab.3 Solution results of six algorithms

第二次測(cè)試的實(shí)驗(yàn)數(shù)據(jù)是使用兩個(gè)GPS/BDS 雙模天線進(jìn)行不同距離的短基線解算，基線長(zhǎng)度分別為0.43 m 和3.45 m，可以從接收模塊中得到載波相位和偽距等數(shù)據(jù)，然后利用串口把數(shù)據(jù)傳輸給計(jì)算機(jī)，最后在MATLAB 上運(yùn)行六種算法程序進(jìn)行整周模糊度解算。實(shí)驗(yàn)總共采集了兩次數(shù)據(jù)，第一次采集的是GPS數(shù)據(jù)，總共觀測(cè)了3875 個(gè)歷元；第二次采集的是BDS數(shù)據(jù)，總共觀測(cè)了3627 個(gè)歷元。兩次數(shù)據(jù)都取前1000個(gè)歷元進(jìn)行解算，表4 和表5 是六種算法分別使用GPS 和BDS 數(shù)據(jù)的解算結(jié)果。

表4 六種算法使用GPS數(shù)據(jù)的解算結(jié)果Tab.4 Solution results of six algorithms using GPS data

表5 六種算法使用BDS數(shù)據(jù)的解算結(jié)果Tab.5 Solution results of six algorithms using BDS data

由表3～表5 可知，長(zhǎng)基線的解算成功率比短基線的低很多，因?yàn)榛€長(zhǎng)度的增加使得模糊度范圍增大，導(dǎo)致解算成功率降低。本文提出的IBOA 算法在解算整周模糊度的性能上依然有一定的提升，解算時(shí)間比LAMBDA、IPSO、ICSO 和IACO 四種算法都要快，成功率也優(yōu)于其他五種算法。其中，IBOA 算法成功率比LAMBDA 算法平均高5.07%，比IPSO 算法平均高5.77%，比ICSO 算法平均高4.90%，比BOA 算法平均高24.13%，比IACO 算法平均高3.37%。雖然BOA算法解算平均時(shí)間在六種算法中最少，但其解算成功率在六種算法中最低。IBOA 算法在GPS 和BDS 雙模模糊度解算結(jié)果上比較接近，因此將IBOA 算法應(yīng)用于整周模糊度解算，在理論研究和實(shí)際應(yīng)用上都有很大的優(yōu)勢(shì)，適用于GPS 差分和BDS 差分。

4 結(jié)論

本文主要研究了IBOA 算法在DGPS 整周模糊度求解中的應(yīng)用。首先，采用動(dòng)態(tài)切換概率平衡全局搜索和局部搜索的比例，改善了傳統(tǒng)BOA 算法搜索能力弱的缺點(diǎn)；其次，引入自適應(yīng)權(quán)重使香味分布因子范圍越來(lái)越小，同時(shí)通過(guò)調(diào)整最優(yōu)解與每個(gè)解之間的差值來(lái)更新位置，改善了傳統(tǒng)BOA 算法易于陷入局部最優(yōu)的問(wèn)題；最后，采用IBOA 算法求解整周模糊度固定解。仿真結(jié)果表明，IBOA 算法比BOA、IPSO、ICSO 和IACO 算法的收斂速度更快、效果更好。實(shí)驗(yàn)結(jié)果表明IBOA 算法的解算成功率更高，也進(jìn)一步驗(yàn)證了IBOA 算法在GPS 和BDS 雙模模糊度固定中的有效性和應(yīng)用可行性。