基于天空偏振模式圖的抗噪聲方位角算法

關(guān)桂霞，韓愈章，吳敏華，韓建國

（1.首都師范大學(xué)信息工程學(xué)院，北京100048；2.Department of Computer Science，Duke University，Durham，North Carolina，U.S.A；3.北京化工大學(xué)信息科學(xué)與技術(shù)學(xué)院，北京100029）

基于天空偏振模式圖的抗噪聲方位角算法

關(guān)桂霞1，韓愈章2，吳敏華1，韓建國3

（1.首都師范大學(xué)信息工程學(xué)院，北京100048；2.Department of Computer Science，Duke University，Durham，North Carolina，U.S.A；3.北京化工大學(xué)信息科學(xué)與技術(shù)學(xué)院，北京100029）

針對運(yùn)動載體的相對方位角誤差問題，提出一種基于移動-分段執(zhí)行的卡爾曼濾波算法，以及具有可變初始群值域的專門化遺傳基因算法，提出的算法用來抑制信號源系統(tǒng)的有色噪聲和測量系統(tǒng)的白色噪聲對載體相對方位角的影響。該算法首先對利用判據(jù)解算出的相對方位角的近似解運(yùn)行移動-分段執(zhí)行的卡爾曼濾波算法，得出相對方位角的期望值，然后運(yùn)行遺傳算法得到相對方位角的最優(yōu)估值。計算結(jié)果的精確度和執(zhí)行過程的流暢性表明，該算法與天空偏振光傳感器相配合，能夠為當(dāng)前迅速發(fā)展的自主定位技術(shù)提供一個可靠的實(shí)現(xiàn)手段。

相對方位角；卡爾曼濾波算法；遺傳基因算法；最優(yōu)估值；偏振光傳感器；自主定位

大量的研究表明，基于天空偏振模式圖的仿生導(dǎo)航技術(shù)能夠為載體提供方位信息［1?9］。同時，相關(guān)文獻(xiàn)也提供了載體利用天空偏振模式圖獲取方位信息的手段［4?7］。然而，由測量噪聲和信息源本身的誤差對仿生偏振導(dǎo)航系統(tǒng)的影響不言而喻，國內(nèi)外學(xué)者在這方面都有一些研究［7?9］，相關(guān)文獻(xiàn)也從不同側(cè)面分析了誤差產(chǎn)生的原因并提出了解決方案，但將測量噪聲和信息源本身的誤差綜合考慮研究的較少。為此，本文研究了一種在線的動態(tài)方位角最優(yōu)估計算法，該算法與天空偏振光傳感儀表相配合，能夠為當(dāng)前迅速發(fā)展的自主定位技術(shù)提供一個可靠的實(shí)現(xiàn)手段。

1 利用天空偏振光測量載體方位

基于沙蟻的偏振導(dǎo)航原理，仿生學(xué)家嘗試著利用無機(jī)光電材料研制偏振敏感羅盤模型［7］。偏振敏感羅盤的基本組成是夾角為60°的3對偏振正交單元，3對偏振正交單元的輸出如下：

式中：φ為載體的相對方位角，d為載體觀測的偏振度，pi（φ）為偏振正交單元的輸出。利用式（1）即可解算載體運(yùn)行方向相對于太陽的方位角φ。

如果對式（1）進(jìn)行取指數(shù)變換，得到

當(dāng)i＝0，1，2時可以得到一個重要的數(shù)學(xué)特征：

理論上可以利用這一結(jié)論進(jìn)行精度校驗。

2 相對方位角的最優(yōu)估計算法

實(shí)際應(yīng)用中，以上的計算方法在有噪聲干擾的情況下結(jié)果不盡人意。以表1為例。

表1 理論值與含誤差的仿真測量值比較Table1 True values and measured values with errors

表1中，REP為p1（φ）、p2（φ）和p3（φ）的相對誤差；REφ為φ的相對誤差；“0%”、“15%”、“30%”分別對應(yīng)的是理論計算的p1（φ）、p2（φ）和p3（φ）的仿真值和人為添加15%、30%相對誤差之后的仿真值。表中數(shù)據(jù)顯示：p1（φ）、p2（φ）和p3（φ）相對誤差的增大，不僅使得φ的相對誤差增大，而且導(dǎo)致了判據(jù)f零值的不可實(shí)現(xiàn)。這就是說，實(shí)際上不可能建立起判據(jù)，也不可能獲得φ的精確結(jié)果?？紤]到仿生導(dǎo)航技術(shù)在實(shí)際應(yīng)用中所遇到的噪聲主要來源于現(xiàn)場測量系統(tǒng)（以白色為主）和被測量信號的信號源系統(tǒng)（以有色為主），而前者的信息主要靠對現(xiàn)場實(shí)際環(huán)境的噪聲狀況進(jìn)行實(shí)時分析而獲得；后者的信息具有規(guī)律性，可以通過理論分析和典型實(shí)驗分析來獲得。據(jù)此，本文提出了如下應(yīng)對手段：在采取不同辦法獲得這2個方面的信息之后，通過一個對上述二者都生效的手段來削弱噪聲影響，獲取所期待的最優(yōu)估計的近似值。此手段基于如下4個步驟：

1）在粗略進(jìn)行直接解方程的基礎(chǔ)上運(yùn)行一個專用的遺傳算法，該算法具有如下特點(diǎn)：

①使用

作為目標(biāo)函數(shù)，其中a1，a2和a3為權(quán)值（實(shí)數(shù)）；q為下面步驟中的各遺傳代的序列號。

②使用有控地隨時間變化的量值區(qū)域作為初始群的取值區(qū)（記作RIPB）。

2）運(yùn)行“基于在線建模的移動-分段卡爾曼算法”［10?11］獲取實(shí)時的相對方位角粗略解的期望值，用以作為白噪聲狀況分析的依據(jù)。

3）通過典型實(shí)驗分析得到相對方位角的粗略解的有色噪聲狀況（這一工作應(yīng)在事前完成）。

4）綜合步驟2、3的2個因素而確定信噪比并構(gòu)建RIPB，用以支持上述專用的遺傳算法［12?13］，從而實(shí)現(xiàn)最后生成所期待的相對方位角最優(yōu)近似值。

通過上述4個關(guān)鍵步驟，測量計算過程中遇到的純隨機(jī)干擾和系統(tǒng)性干擾將得到有效克服。

3 實(shí)現(xiàn)算法

按照上文所闡述的問題解決方案，針對經(jīng)歷時間k＝1，2，…，N作持續(xù)運(yùn)動的載體，建立在線估計載體自身相對于太陽方位的相對方位角和天空光偏振度｛φk，dk｝的算法（OAE算法）如下：

假設(shè)載體運(yùn)動過程中于特定時間k＝1，2，…，N由偏振傳感器的偏振正交單元輸出3個電壓值為。

1）個體初始化。

建立并粗略求解公式（1）的3個方程組｛（1），（2）｝，｛（2），（3）｝和｛（1），（3）｝，從中選擇出最合理的解（記作｛φ'0k，d'0k｝），選擇的標(biāo)準(zhǔn)是從式（3）演變出來的判據(jù)：“量值f 近似等于0”。最后所得結(jié)果將作為下文中執(zhí)行基于在線建模的移動-分段卡爾曼算法的依據(jù)。

2）解的期望值估計。

在隨時間移動的時間段上，即在時間段［j-T，j］（T：一個依賴于信號與噪聲狀況的偶數(shù)；j：取值于［k-T，k］）上運(yùn)行基于在線建模的移動-分段卡爾曼算法，獲得｛φ'0j，d'0j｝，j＝k-T，k-T+1，…，k的期望值，記作｛φ＾'0j，d＾'0j｝。

3）初始群取值區(qū)（RIPB）的構(gòu)成。

計算實(shí)時噪聲／信號比

然后構(gòu)建k鄰域的RIPB如下：

式中：實(shí)數(shù)α'1，α'2≥0.5Rφk和γ'1，γ'2≥0.5Rdk表述來自測量系統(tǒng)的白噪聲的噪聲／信號比；實(shí)數(shù)γ''1，γ''2，α''1和α''2表述來自信號源系統(tǒng)的有色噪聲的噪聲／信號比，可事先通過典型實(shí)驗理論分析確定。

4）群初始化。

以合適的密度（取決于實(shí)際需要）在步驟3所獲得的RIPB中選擇足夠數(shù)目的值作為φ和d的初始值，記作φ0和d0，從而構(gòu)成初始群p0。p0包括初始子群Gφ0（屬于量值φ0）和Gd0（屬于量值d0）。

5）最優(yōu)值確定。

使用式（4）所表述的適配函數(shù)，運(yùn)行專用的遺傳算法。每個時刻的運(yùn)算周期都起步于上一步所確定的初始群p0。

最后，經(jīng)過在k鄰域上運(yùn)行足夠次數(shù)的專用的遺傳算法周期，所期待的載體相對方位角和所在位置天空光偏振度｛φk，dk｝的最優(yōu)估計值（記作｛｝）將被產(chǎn)生。

6）步驟重復(fù)。

設(shè)置k＝k+1，重復(fù)步驟1～5，直到k＝N。

4 仿真實(shí)驗

4.1 實(shí)驗描述

本實(shí)驗的目的是在含白噪聲和有色噪聲混合干擾（噪聲／信號比達(dá)10%）的仿真測量值數(shù)組序列｝，k＝1，2，…，150的基礎(chǔ)上產(chǎn)生出｛φk，dk｝的最優(yōu)估計值，k＝1，2，…，150。實(shí)驗由如下4個步驟組成：

1）實(shí)際測量值序列仿真。

在（0，π）弧度區(qū)間內(nèi)選取平滑變化的數(shù)值作為載體運(yùn)動軌跡的相對方位角，并用數(shù)學(xué)方法產(chǎn)生相應(yīng)的天空光偏振度（弧度），構(gòu)成150對理論計算的真值｛φk，dk｝，k＝1，2，…，150（圖1中φk曲線）。在此基礎(chǔ)上用數(shù)學(xué)方法仿真出含白噪聲和有色噪聲混合干擾（噪聲／信號比約為10%）的150組偏振傳感單元輸出端測量數(shù)據(jù)｝，k＝1，2，…，150，簡稱“仿真測量值”。

圖1 不同運(yùn)算階段的數(shù)據(jù)曲線Fig.1 Data?curves in different stages

2）近似值計算。

使用數(shù)據(jù)序列｛p?1（k），p2?（k），p3?（k）｝建立并粗劣求解方程組｛（1），（2）｝，｛（2），（3）｝，和｛（1），（3）｝，使用判據(jù)“量值f 近似等于0”從解中選出最合理的近似解｛φ'0k，d'0k｝ （圖1中φ'0k曲線）。

｛JP3）期望值估計與構(gòu)成RIPB。

針對｛φ'0k，d'0k｝運(yùn)行基于在線建模的移動-分段卡爾曼算法，估計出｛φk，dk｝的盡可能精確的期望值｛φ＾'0j，d＾'0j｝（圖1中φ＾'0j曲線），并選取

α'1，α'2＝0.6Rφk， α''2，α''2＝0.3

γ'1，γ'2＝0.6Rdk， γ''2，γ''2＝0.3

以此構(gòu)成RIPB（圖2中的L1u，L1l-Gφ0的上、下限曲線以及它們所包含的φ＾'0j曲線）。

4）最優(yōu)值生成。

使用適配函數(shù)式（4）并遵從非線性約束

c＝0.008x1k+0.013sin（0.025x1k）-x2k-0.05，T＝10；在這些條件下運(yùn)行專用的遺傳算法，其中x1k、x2k為｛φk，dk｝的遺傳個體的過渡；算法的每個運(yùn)算周期從上述步驟3所建立的RIPB開始。最后，在k時刻（k＝11，12，…，150）通過足夠的運(yùn)算周期而獲得了所期待的相對方位角和偏振度｛φk，dk｝的最優(yōu)估計值｛｝（圖1中曲線）。

圖2 Gφ0的取值區(qū)Fig.2 Boundaries of Gφ0

4.2 實(shí)驗結(jié)果分析

比較仿真測量值與算法仿真結(jié)果（圖1中φk、，可計算出計算結(jié)果的噪聲／信號比，記作λφ：

符合弱噪聲標(biāo)準(zhǔn)。

另一方面，由圖3的曲線f可看出，從k＝11到k＝150的適配函數(shù)的值滿足f?1，這說明專用的遺傳算法發(fā)揮了有效功能，因此生成了較高精度的最優(yōu)值。

圖3 適配函數(shù)值曲線Fig.3 Fitness function values curve

通過基于在線建模的移動-分段卡爾曼算法而獲得的期望值估計結(jié)果｛φ＾'k，d＾'k｝的隨機(jī)擺動顯示出比粗略計算的值｛φ'k，d'k｝的隨機(jī)擺動有明顯縮小，但偏差的實(shí)時積分值：

上式結(jié)果明顯大于采用遺傳算法獲得的最終結(jié)果的實(shí)時積分值

這意味著盡管期望值估計曲線比較平滑，但較大的系統(tǒng)性誤差被基于在線建模的移動-分段卡爾曼算法保留了，這是因為卡爾曼濾波算法的功能是削弱白噪聲而保留有色噪聲。與此相反，專用的遺傳算法則將白噪聲和有色噪聲的干擾都削弱了，因為專用的遺傳算法優(yōu)化準(zhǔn)則是與信號性質(zhì)無關(guān)的適配函數(shù)。

值得注意的是，在建立式（5）、（6）時同時考慮了白噪聲和有色噪聲的影響，因而OAE算法能夠同時削弱白噪聲和有色噪聲的影響。

5 實(shí)際測量實(shí)驗

5.1 單組數(shù)據(jù)實(shí)驗

圖4 過渡值曲線Fig.4 Transition value curves

5.2 測量值序列實(shí)驗

圖5 偏振測量儀輸出Fig.5 The polarization instrument output

圖6 參數(shù)估計與實(shí)際值比較Fig.6 Estimated and practical parameters

針對上述50組測量值運(yùn)行OAE算法，結(jié)果生成了最優(yōu)的方位角參數(shù)值｛φ?（k），d?（k）｝，k＝1，2，…，50，（圖6中φ?（k）、d?（k））。

6 結(jié)束語

本文設(shè)計并實(shí)現(xiàn)了用于仿生偏振導(dǎo)航的OAE算法及其軟件，該算法利用移動-分段卡爾曼算法和專用的遺傳算法使得載體在有色噪聲和白噪聲同時干擾的環(huán)境下能夠在線地獲得盡可能精確的相對方位角最優(yōu)估計值。仿真實(shí)驗和實(shí)際測量實(shí)驗的結(jié)果表明，OAE算法不僅削弱了測量系統(tǒng)產(chǎn)生的白色噪聲對載體相對方位角的影響，而且削弱了信號源系統(tǒng)產(chǎn)生的有色噪聲對載體相對方位角的影響，這說明本文提出的算法對仿生偏振導(dǎo)航載體相對方位角的誤差補(bǔ)償是可行且有效的。

［1］MICHAEL L.Dynamic patterns of skylight polarization as clock and compass［J］.Theory Biol，1980，86：507?512.

［2］ABLE K P，ABLE M A.Calibration of the magnetic compass of a migratory bird by celestial rotation［J］.Nature，1990，347：378?380.

［3］LAMBRINOS L.Navigation in desert ants：the robotic solu?tion［J］.Robotic，2003：407?426.

［4］YAN Lei，GUAN Guixia，CHEN Jiabin，et al.The bionic o?rientation mechanism in the skylight polarization pattern［J］.Acta Scientiarum，Naturalium Universitatis，Pekinensis，2009（4）：616?620.

［5］CAREY N，STURZL W.An Insect?inspired omnidirectional vision systemincludingUV?Sensitivityandpolarization［C］／／IEEE International Conference on Computer Vision Workshops.Barcelona，Spain，2011：312?319.

［6］CHU Jinkui，ZHAO Kaichun.Study of angle measurement optoelectronic model on emulating polarization?sensitive com?pound eye of insect［J］.Nanoelectronic Device＆Technolo?ny，2005（12）：541?545.

［7］LAMBRINOS L，MOLLER R，LABHART T.A mobile robot employing insect strategies for navigation［J］.Robotics and Autonomous Systems，2000，30（1）：39?64.

［8］CHU Jinkui，ZHAO Kaichun，ZHANG Qiang.Construction and performance test of a novel polarization sensor for navi?gation［J］.Sensors and Actuators，2008，A148：75?82.

［9］CHU Jinkui，CHEN Wenjing，WANG Hongqing，at el.Mo?bile robot navigation tests with polarization sensors［J］.Op?tics and Precision Engineering，2011，19：2419?2426.

［10］MOHINDER S，GREWAL A，ANDREWS P.Kalman filte?ring theory and practice using MATLAB［M］.Hoboken N.L.：Wiley，2008：53?66，132?141.

［11］HAN Jianguo，SUN Cauhong，LI Yanqin.An on?line mod?eling based Kalman filtering process for time?interval?varia?ble sequences with application to astronomy surveying［J］.Journal Engineering and Electronics，2003，14（2）：52?56.

［12］HYUN C J，KIM Y K.A genetic algorithm for multiple ob?jective sequencing problems in mixed model assembly lines［J］.Computers＆Operations Research，1998（25）：7?8.［13］TSOU Mingcheng.Genetic algorithm for solving celestial navigation fix problems［J］.Polish Maritime Research，2012（3）：53?59.

The anti?noise algorithm of azimuth based on sky?lighted polarization pattern

GUAN Guixia1，HAN Yuzhang2，WU Minhua1，HAN Jianguo3
（1.Information Engineering College，Capital Normal University，Beijing 100048，China；2.Department of Computer Science，Duke University，Durham，North Carolina，U.S.A；3.School of Information Science and Technology，Beijing University of Chemical Tech?nology，Beijing 100029，China）

Considering the relative azimuth error of the moving object，a Kalman filter algorithm based on moving block execution and a specialized genetic algorithm with variable initial group of domain were proposed.It is to re?strain the influence brought by the colored noise of signal source system and the white noise of measurement system error.The expected value of the relative azimuth was calculated firstly by using Kalman filter algorithm based on moving block execution，and then the optimal estimation of the relative azimuth was got by running the genetic algo?rithm.The accuracy of the results and the smoothness of the execution show that the algorithm cooperating with the polarization instrument can provide a reliable realization means for the rapidly growing autonomous positioning tech?nology.

relative azimuth；Kalman filter algorithm；genetic algorithm；optimal estimation；polarization instru?ment；autonomous positioning

10.3969／j.issn.1006?7043.201311072

http：／／www.cnki.net／kcms／detail／23.1390.U.20150109.1706.024.html

TP721.1

A

1006?7043（2015）03?0404?05

2013?11?21.網(wǎng)絡(luò)出版時間：2015?01?09.

國家自然科學(xué)基金資助項目（61174220）；北京市教委基金資助項目（KM201210028002）.

關(guān)桂霞（1972?），女，副教授，博士.

關(guān)桂霞，E?mail：guanguixia＠sina.com.