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

（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，韓建國(guó)3

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

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

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

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

1 利用天空偏振光測(cè)量載體方位

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

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

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

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

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

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

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

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

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

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

①使用

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

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

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

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

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

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

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

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

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

1）個(gè)體初始化。

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

2）解的期望值估計(jì)。

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

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

計(jì)算實(shí)時(shí)噪聲／信號(hào)比

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

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

4）群初始化。

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

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

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

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

6）步驟重復(fù)。

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

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

4.1 實(shí)驗(yàn)描述

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

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

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

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

2）近似值計(jì)算。

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

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

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

α'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的上、下限曲線(xiàn)以及它們所包含的φ＾'0j曲線(xiàn)）。

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

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

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

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

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

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

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

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

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

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

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

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

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

5 實(shí)際測(cè)量實(shí)驗(yàn)

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

圖4 過(guò)渡值曲線(xiàn)Fig.4 Transition value curves

5.2 測(cè)量值序列實(shí)驗(yàn)

圖5 偏振測(cè)量?jī)x輸出Fig.5 The polarization instrument output

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

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

6 結(jié)束語(yǔ)

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

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

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

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