盧 舟， 劉鐘中， 孫秀斌， 賈岸斌

（1.常德市氣象局，湖南 常德 415000；2.成都信息工程大學(xué) 電子工程學(xué)院，成都 610225；3.中國氣象局大氣探測重點(diǎn)實驗室，成都 610225）

0 引言

極端天氣中的閃電活動是一種災(zāi)害性天氣現(xiàn)象，以高電壓、大電流、強(qiáng)電磁輻射等為主要物理特征，常常造成嚴(yán)重的人員傷亡和經(jīng)濟(jì)損失［1］。隨著微電子、計算機(jī)、電力以及航空航天等技術(shù)的不斷發(fā)展，云閃引起的災(zāi)害越來越多，造成的危害也越來越大［2］。云閃過程中輻射的甚高頻段信號比地閃輻射的低頻段信號所包含的信息豐富很多，可利用云閃對地閃進(jìn)行預(yù)警。因此，研究精度高實時性好的云閃探測系統(tǒng)對云閃放電通道進(jìn)行定位就顯得非常緊要。

由于云閃多發(fā)生在云內(nèi)且具有瞬時性、劇烈性和頻繁性，要對其進(jìn)行準(zhǔn)確實時的監(jiān)測定位仍然有較大難度。近年來，隨著計算機(jī)技術(shù)、電子技術(shù)和相關(guān)理論的不斷發(fā)展，云閃定位隨之得到快速發(fā)展，許多理論方法在實際工程應(yīng)用中也逐步得到實現(xiàn)。根據(jù)張衛(wèi)杰等［3］所著文章中對數(shù)據(jù)采集系統(tǒng)的硬件設(shè)計可知，用于云閃定位的硬件系統(tǒng)可由云閃VHF信號接收天線、前端信號調(diào)理模塊、A/D采樣模塊和FPGA＋DSP數(shù)據(jù)處理模塊等組成。本文采用MUSIC測向算法作為DSP中的解算算法，對其在云閃定位上的應(yīng)用進(jìn)行了分析，提出了優(yōu)化MUSIC算法的方法，并用Matlab仿真了云閃空間譜圖。

1 閃電定位方法

閃電是發(fā)生于云間、云地間或者云體內(nèi)的瞬態(tài)強(qiáng)放電過程，大致可將其分為云閃和地閃兩類。閃電定位主要是通過測量閃電輻射的聲、光、電磁場信息，然后計算其相關(guān)參數(shù)來確定閃電的空間位置，進(jìn)而推測其活動趨勢［4］。云閃定位和地閃定位原理相近，但前者是利用云閃電磁脈沖輻射源發(fā)出的甚高頻信號，后者則是利用閃電輻射的低頻和甚低頻電磁波。目前，閃電定位方法主要包括定向法、時差定位法以及時差測向聯(lián)合法，總體趨于向時差測向聯(lián)合法（TDOA）發(fā)展，該方法是通過多個子站既記錄閃電輻射源信號到達(dá)時間差又解算其來波信號方向信息來確定閃電空間位置。采用時差測向聯(lián)合法的多站云閃定位系統(tǒng)空間分布圖如圖1所示：

圖1 多站云閃定位系統(tǒng)空間分布圖

在多站云閃定位系統(tǒng)中，各探測子站的定位誤差主要來源于觀測時間誤差、布站位置誤差以及算法誤差。本文采用MUSIC算法作為時差測向聯(lián)合法中的測向算法來測算云閃輻射源信號的方位角和仰角，具有很高的分辨率和實時性，可以提高測向精度、減少定位誤差。

2 基于均勻圓陣的二維MUSIC算法

2.1 均勻圓陣模型

均勻圓陣（UCA）是一種具有全向探測能力的陣列，它能同時獲取多個云閃輻射源的方位角和俯仰角，從而對云閃進(jìn)行定位。它由N個相同的全向陣元組成，按等間隔分布在半徑為R的圓周上，圓周中心O為球坐標(biāo)系的原點(diǎn)和云閃入射信號的相位參考點(diǎn)，如圖 2 所示［5］［6］。 圖中，θi和 φi分別為云閃輻射源入射信號的仰角和方位角，i＝1，2，3…。若有P個遠(yuǎn)場窄帶信號源S（t）入射到此均勻圓陣中，則第 i個信源的入射角度為（θi，φi）。

圖2 探測系統(tǒng)中的均勻圓陣模型

方向矢量 a（θi，φi）是入射信號（θi，φi）對應(yīng)的陣列響應(yīng)，其表達(dá)式為［7］：

上式中：γm＝2πn/N（n＝0，1，…，N-1）；R 為圓陣半徑。

事實上，均勻圓陣所接收到的信號和均勻線陣所接收到的信號是相同的，但是二者對比起來均勻圓陣的陣因子 A＝［a （θ1，φ1） a （θ2，φ2.） … a（θP，φP）］中有入射信號的方位角和仰角兩個因子，比線陣的陣因子包含的方位信息多。采用MUSIC算法實現(xiàn)測向?qū)﹃嚵薪Y(jié)構(gòu)沒有嚴(yán)格要求，因此可以直接使用均勻圓陣對來波的方位角和仰角進(jìn)行測量。

2.2 二維MUSIC算法原理

MUSIC算法是通過對均勻圓陣天線所接收到的云閃輻射源信號的協(xié)方差矩陣進(jìn)行特征分解，然后由信號特征值、噪聲特征值分別得到信號子空間VS和噪聲子空間VN，再利用它們的正交特性對譜函數(shù)進(jìn)行譜峰搜索，最終完成對云閃信號的方位估計［8］。

假設(shè)P個窄帶信號入射到均勻圓陣天線中，陣元數(shù)為M，陣元間夾角為固定值，而且陣元數(shù)大于此時的信號數(shù)。若這些信號的波達(dá)方向分別為（θ1，φ1），（θ2，φ2），...，（θP，φP），陣因子 A＝［a（θ1，φ1），a（θ2，φ2），...，a（θP，φP）］M×P，n（t）為噪聲矢量，則陣列輸出信號可用下式表示［9］：

相應(yīng)的輸入陣列協(xié)方差矩陣為：

其中，RSS＝E［S（t）SH（t）］是信號相關(guān)矩陣，σn2是噪聲方差。

對2-3式中的RX進(jìn)行特征值分解得到M個特征值，對應(yīng) M 個特征向量 VK（K＝1，2，…，M）。 由較大的P個特征值對應(yīng)信號子空間，較小的M-P個特征值對應(yīng)噪聲子空間。把噪聲特征向量作為列構(gòu)造成噪聲子空間VN：

對應(yīng)于信號分量的導(dǎo)向矢量空間與VN正交，可得到空間譜搜索函數(shù)Pmusic（θ）如下式：

通過上式可知，分母中的 VN和 a（θ，φ）正交，所以Pmusic（θ，φ）的值在云閃信源處會有一個很尖的譜峰，此處即為云閃信源的位置。當(dāng)云閃輻射出多個信號時，可以通過搜索空間譜中的波峰來估計各信號的來波角度信息。二維MUSIC算法處理流程圖如圖3所示：

圖3 MUSIC算法簡要流程

3 算法優(yōu)化

在云閃定位系統(tǒng)中采用MUSIC算法來解算云閃輻射源的方位角和仰角等閃電參數(shù)具有很大的優(yōu)勢，該方法具有實時性好、測向精度高、能多信號同時測向等優(yōu)點(diǎn)。采用MUSIC算法實現(xiàn)云閃定位主要包括空間譜函數(shù)構(gòu)造、信源數(shù)估計以及譜峰搜索三個步驟，下面詳細(xì)介紹對這三個部分的優(yōu)化。

3.1 空間譜函數(shù)優(yōu)化

對均勻圓陣輸出信號的協(xié)方差矩陣進(jìn)行優(yōu)化處理，充分利用求得的特征值對應(yīng)的信號子空間和噪聲子空間包含的信息。

優(yōu)化空間譜函數(shù)構(gòu)造方法的實現(xiàn)步驟如下：

（1）采集云閃輸入信號樣本 X（i），i＝1，2，…N，計算輸入?yún)f(xié)方差矩陣Rx，如公式（2-3）所示，構(gòu)造新的協(xié)方差矩陣RX如式（3-1）：

上式中I為M×M的反向單位矩陣，R^*X為R^X的共軛矩陣。

（2）對3-1式中新構(gòu)造的協(xié)方差矩陣RX進(jìn)行特征分解重構(gòu)矩陣得：RX＝UΛH，可得對應(yīng)的信號子空間RXS1和噪聲子空間RXN1如下式：

（3）用低秩矩陣代替滿秩矩陣RX，再次構(gòu)造新的協(xié)方差矩陣RXX并對其進(jìn)行特征值分解，重復(fù)步驟（2），再一次構(gòu)造新的信號子空間和噪聲子空間［10］：

（4）對兩次得到的信號子空間和噪聲子空間進(jìn)行平均，如式（9）：

（5）根據(jù)US和UN構(gòu)造新的空間譜函數(shù)，如下式：

從上式可知，該空間譜函數(shù)充分利用了信號子空間和噪聲子空間所包含的信息，提高了MUSIC算法的穩(wěn)定性和在非理想情況下的分辨性能。根據(jù)上式進(jìn)行譜峰搜索，可以得到準(zhǔn)確性高、精度好的入射云閃輻射源信號方位角θ和仰角φ。

3.2 信號源數(shù)估計優(yōu)化

在MUSIC算法中估計信號源的來波方向之前需要先估計出信源數(shù)，否則會出現(xiàn)虛警或漏警。一般選用基于信息論準(zhǔn)則的信源數(shù)估計方法來完成該步驟，該方法在有高斯白噪聲的情況下能獲得較好的估計性能。但是在實際的云閃輻射源探測環(huán)境中，天線陣接收到的噪聲往往是空間色噪聲，這樣便不能準(zhǔn)確獲取陣列協(xié)方差矩陣的噪聲特征值，從而造成該方法的估計性能較差。為了解決這個問題，采用基于總體最小二乘擬合法來估計云閃信號源數(shù)。

基于總體最小二乘法（TLS）的信源數(shù)估計方法步驟如下：

（1） 算出均勻圓陣接收信號的輸出協(xié)方差矩陣：

（2）根據(jù)MUSIC算法原理可知，對協(xié)方差矩陣進(jìn)行特征值分解可以得到M個特征值，將其按從大到小的順序排列如下［11］：

（3）對式3-7中后面的k個特征值進(jìn)行擬合，擬合點(diǎn)數(shù)i（2≤i≤M）從所得特征值個數(shù)M到2逐漸減小，依次可以獲得對應(yīng)的擬合偏差 ε（1），ε（2），…，ε（M-1）。

（4）假設(shè)云閃入射信號源數(shù)為 n（n＝1，2，…，M），因為ε（k）為參與擬合數(shù)目的倒數(shù)擬合，所以待估計的信號源數(shù)目為。按式＝argmin{ε（k）＋p（k）}構(gòu)造目標(biāo)估計函數(shù)，式中罰函數(shù) p（k）＝1/2k（2M-k）logL，L為快拍數(shù)，能夠使該式取得最小值的k值就是所要估計的云閃信號源數(shù)。

3.3 譜峰搜索優(yōu)化

由于譜峰搜索在MUSIC算法中耗時較多，影響了云閃定位的實時性，所以在保證一定測向精度的情況下，提出了一種精度高實時性好的搜索方法。（1）采用一個較大的角度搜索步長（如 2°），在［0°，180°］的范圍里盡可能尋找出所有極大值出現(xiàn)的大概位置。假設(shè)粗略搜索到M個極大值，在沒有噪聲和干擾信號的情況下，M值等于信源數(shù)P。但在實際探測環(huán)境中存在無法避免的噪聲以及干擾信號，所以會有多于信源數(shù)P的極大值出現(xiàn)。此時需要對M和信號源估計數(shù)P的大小進(jìn)行比較，若M＝P則直接進(jìn)入細(xì)搜，否則，將所有峰值依次從大到小排序，選出其中P個最大的極大值，這些值就對應(yīng)著信號到達(dá)方向。 （2）在每個經(jīng)粗選的極大值（i，j）附近采用較小的角度搜索步長（如0.1°）進(jìn)行搜索，并分別與周圍八個點(diǎn)（i-1，j）、（i＋1，j）、（i，j-1）、（i，j＋1）、（i＋1，j＋1）、（i＋1，j-1）/（i-1，j-1）、（i-1，j＋1）比較；如果該極大值大于周圍八個點(diǎn)，則其為有效譜峰，（i，j）即為來波信號的方位角和仰角。（3）如果譜峰搜索結(jié)果還未達(dá)到精度要求，那么返回（2）粗略搜索出的峰值附近繼續(xù)減小角度搜索步長，重復(fù)進(jìn)行搜索，直到結(jié)果滿足云閃定位的精度要求［12］。

采用這種一次粗搜全方位多次細(xì)搜的優(yōu)化方法，可以避免在要求較高精度的爬山搜索算法中不必要的計算量，并且可以快速確定云閃信號來波方向的大概角度，然后再通過多次細(xì)搜找到滿足精度要求的方位角、仰角。優(yōu)化后的譜峰搜索方法可以極大的減少搜索運(yùn)算時間，同時提高了云閃定位結(jié)果的估計精度。搜索算法流程圖如圖4所示：

圖4 譜峰搜索流程圖

4 算法仿真及性能分析

假設(shè)均勻圓陣半徑滿足 r≤1/［4sin（π/M）］，陣元個數(shù) M＝6，快拍數(shù)：1024，信號源數(shù) P＝3；信號的方位角分別為（80°、100°、120°），俯仰角分別為（40°、50°、60°）， 分別在信噪比取 20dB、30dB 的情況下進(jìn)行二維MUSIC譜估計，如圖5（見彩頁）、圖6（見彩頁）所示：

假設(shè)均勻圓陣半徑滿足 r≤1/［4sin（π/M）］，陣元個數(shù) M＝8，快拍數(shù)：1024，信號源數(shù) P＝3；信號的方位角分別為（80°、100°、120°），俯仰角分別為（35°、40°、45°），信噪比均為 30dB，二維 MUSIC 算法優(yōu)化前和優(yōu)化后的二維MUSIC算法譜估計如圖7、圖8所示：

通過多次二維分級譜峰搜索獲得云閃輻射源的估計方位角、仰角平均值，如表一所示。

表1 兩種情況下二維估計仿真結(jié)果

通過對比圖六和圖七的空間譜圖可以看出，如果其他條件相同，增大信噪比可以使基于均勻圓陣的二維MUSIC算法DOA估計性能得到提高，對三個信號的分辨能力也得到增強(qiáng)；圖7（見彩頁）和圖8（見彩頁）的仿真條件考慮到云閃發(fā)生的連續(xù)性，針對信源之間的角度間隔可能會比較小進(jìn)行了設(shè)置；將優(yōu)化前和優(yōu)化后的二維MUSIC算法進(jìn)行仿真對比發(fā)現(xiàn)，優(yōu)化后的二維MUSIC算法譜峰更加陡峭，表明對來波方向估計的分辨力得到進(jìn)一步提升。從表一中的方位角和仰角估計數(shù)據(jù)可知，優(yōu)化后的二維MUSIC算法的DOA估計精度得到了很大的提高，并且對多個相隔較近的云閃輻射源信號能夠同時測向。

5 結(jié)論與展望

本文簡要介紹了閃電定位基本方法和均勻圓陣建模模型，重點(diǎn)概述了MUSIC算法中各個重要步驟的優(yōu)化方法；對算法中構(gòu)造空間譜函數(shù)的方法進(jìn)行了優(yōu)化，充分利用了信號和噪聲所包含的信息，提高了算法的分辨率；針對實際的空間噪聲環(huán)境，提出采用基于總體最小二乘擬合法來估計云閃信號源數(shù)；針對MUSIC算法中譜峰搜索耗時較多影響定位實時性的問題，提出了一種精度高實時性好的譜峰搜索方法。最后給出了不同信噪比情況下基于二維MUSIC算法的空間譜估計結(jié)果，對優(yōu)化前和優(yōu)化后的算法性能進(jìn)行了分析。通過對比分析，可以看到基于均勻圓陣二維優(yōu)化MUSIC算法作為云閃定位系統(tǒng)測向算法的優(yōu)越性，使用這種優(yōu)化的測向方法可以有效提高測向結(jié)果的實時性和精度。將這種方法應(yīng)用于估計云閃輻射源的來波方向，再結(jié)合多片并行高速數(shù)字信號處理器，可以實時解算多個云閃信號的方位角和仰角，準(zhǔn)確獲取云閃的空間位置，為云閃定位系統(tǒng)的實現(xiàn)提供了一種更優(yōu)的方法。

參考文獻(xiàn)：

［1］李云敏，孫秀斌，楊崧令，等.云閃定位算法及誤差分析［J］.電子測量與儀器學(xué)報，2012，26（10）：917-922.

［2］張義軍，周秀驥.雷電研究的回顧和進(jìn)展［J］.應(yīng)用氣象學(xué)報，2006，17（6）：829-834.

［3］張衛(wèi)杰，候孝民.高速大容量數(shù)據(jù)采集系統(tǒng)設(shè)計與實現(xiàn)［J］.電子測量與儀器學(xué)報，2005，19（4）：51-55.

［4］張義軍，孟青，馬明，等.閃電探測技術(shù)發(fā)展和資料應(yīng)用［J］.應(yīng)用氣象學(xué)報，2006，17（5）：611-620.

［5］ Mathews C P， ZoltowskiM D.Eigenstructure Techniques for 2-D Angle Estimation with Uniform Circular Array［J］.IEEE trans on Sp r，1994，42（9）：2395-2404.

［6］ 宋士瓊，陳旗，胡生亮.基于九元均勻圓陣的MUSIC測向方法研究［J］.計算機(jī)仿真，2008，25（3）：123-125.

［7］ 張小飛，汪飛，陳偉華.陣列信號處理的理論與應(yīng)用（第2 版）［M］.北京，國防工業(yè)出版社，2013.

［8］王永良，陳 輝，彭應(yīng)寧，等.空間譜估計理論與算法［M］.北京：清華大學(xué)出版社，2005：82-139.

［9］ SchmidtR.O.Multiple emitterlocation and signal parameterestimation ［J］.IEEE Trans.OnAP， 1986，34（3）：276-280.

［10］趙謙，董民，梁文娟.DOA估計算法的一種修正MUSIC算法的研究［J］.計算機(jī)工程與應(yīng)用，2012，48（10）：102-105.

［11］藍(lán)曉宇.提高空間譜估計分辨率的超分辨測向算法研究［D］.哈爾濱，哈爾濱工程大學(xué)，2012.

［12］郭躍，王宏遠(yuǎn)，陳思捷，等.陣列測向中的精確高速并行譜峰搜索算法 ［J］. 微電子學(xué)與計算機(jī)，2007，24（12）：50-54.

［13 潘軍，鐘一帆，蒙劍.雷電技術(shù)在現(xiàn)代建筑中的應(yīng)用［J］.氣象研究與應(yīng)用，2008，29（3）：49-54.

［14］鄒興奮，劉志輝，李可娟.建筑物防雷接地裝置結(jié)構(gòu)探析［J］.氣象研究與應(yīng)用，2010，31（S2）.

［15］黃劍釗.新型國家自動站現(xiàn)場總線CAN的研究及維護(hù)方法［J］.氣象研究與應(yīng)用，2016，37（2）：87-89.

［16］張娟，林卓宏，陳巧淑，等.自動氣象站觀測場防雷接地制式的技術(shù)分析［J］.氣象研究與應(yīng)用，2012，33（4）69-71.

［17］姚家釗，陳華宣.基層氣象臺站防雷問題淺析及解決方法［J］.氣象研究與應(yīng)用，2009，30（1）：78-81.

［18］黃仁立，羅曉軍，周開春.防城港天氣雷達(dá)樓雷擊防護(hù)等級及防護(hù)設(shè)計［J］.氣象研究與應(yīng)用，2014，35（4）：111-114.

［19］黎錦雷，韋菊，楊玉靜.新型自動氣象站故障分析與排除［J］.氣象研究與應(yīng)用，2015，36（4）：102-104.

［20］傅俊霖，黃君健，何肖珍.防雷裝置接地電阻的認(rèn)識和探討［J］.氣象研究與應(yīng)用，2008，29（2）：68-69.

［21］韓建海，吳松.一次風(fēng)電場雷災(zāi)的特點(diǎn)及原因分析［J］.氣象研究與應(yīng)用，2016，37（1）：121-123.

［22］羅天龍.電涌保護(hù)器的簡易檢測方法［J］.氣象研究與應(yīng)用，2013，34（2）：92-94.

［23］楊仲江，盧燕，曹書華.用于防雷工程的電涌保護(hù)器的測試研究［J］.氣象研究與應(yīng)用，2007，28（4）：54-57.

［24］ 譚惠冰，杜建德，梁偉漢，等.電涌保護(hù)器（SPD）在低壓電氣系統(tǒng)中的設(shè)計［J］.氣象研究與應(yīng)用，2013，34（4）：92-97.

［25］林政，黎梓華，唐雷.淺談如何利用法拉第籠原理防護(hù)雷電電磁脈沖［J］.氣象研究與應(yīng)用，2009，30（1）：83-84.