應(yīng)用獨(dú)立分量分析的激光測(cè)風(fēng)雷達(dá)湍流頻譜分解

李一越，胡姝玲

(北京航空航天大學(xué) 儀器科學(xué)與光電工程學(xué)院，北京 100191)

1 引 言

大氣湍流是大氣運(yùn)動(dòng)的一種形式，不同于由規(guī)則的大氣分子運(yùn)動(dòng)形成的層流，大氣湍流屬于不規(guī)則運(yùn)動(dòng)，且運(yùn)動(dòng)的劇烈程度通常大于層流。在地表附近，由于受地形起伏、熱量交換等的影響，近地表的大氣運(yùn)動(dòng)往往存在強(qiáng)度不等的湍流。這些湍流會(huì)導(dǎo)致光信號(hào)在傳播過(guò)程中質(zhì)量下降，具體表現(xiàn)為光束漂移、光束擴(kuò)展、光強(qiáng)閃爍及到達(dá)角起伏等。

在風(fēng)力發(fā)電領(lǐng)域中，湍流作為一種常見(jiàn)的大氣效應(yīng)，對(duì)風(fēng)電系統(tǒng)的危害主要包括風(fēng)能利用率下降、發(fā)電量降低、風(fēng)機(jī)軸疲勞及載荷增大等[1]。如新疆大型陸上風(fēng)電廠中[2]，由于湍流特別是風(fēng)力發(fā)電機(jī)尾流的影響范圍達(dá)到3倍葉片直徑以上，即300 m以上，由此導(dǎo)致的最大風(fēng)速損失率為50%。由于傳統(tǒng)激光測(cè)風(fēng)雷達(dá)無(wú)法將風(fēng)機(jī)尾流同當(dāng)?shù)刂行娘L(fēng)速區(qū)分開(kāi)，使后方風(fēng)機(jī)在受到干擾時(shí)無(wú)法做出正確調(diào)整，因此年平均風(fēng)能損失率約為10%。若不能將尾流及時(shí)分離出來(lái)，則需要根據(jù)當(dāng)?shù)仫L(fēng)場(chǎng)情況適當(dāng)增加風(fēng)力發(fā)電機(jī)間距離，這同樣會(huì)造成額外的風(fēng)能浪費(fèi)。同時(shí)，如果能夠及時(shí)觀測(cè)尾流等湍流對(duì)風(fēng)場(chǎng)的影響，則可以及時(shí)調(diào)節(jié)不同風(fēng)機(jī)的參數(shù)，保證系統(tǒng)整體的工作效率最高。

在風(fēng)力發(fā)電機(jī)機(jī)組中，由于尾流問(wèn)題導(dǎo)致不同位置風(fēng)機(jī)的實(shí)際可利用風(fēng)速的大小和方向均有所區(qū)別。若在風(fēng)速分析中將湍流與中心風(fēng)速進(jìn)行混合分析，則所獲得的風(fēng)速穩(wěn)定性低，同時(shí)需要部署多臺(tái)測(cè)風(fēng)雷達(dá)用于整個(gè)風(fēng)機(jī)組陣列的風(fēng)場(chǎng)測(cè)量，以保證綜合利用率最大。但是這種做法會(huì)明顯提高系統(tǒng)的成本，同時(shí)降低系統(tǒng)的穩(wěn)定性(由于需要多臺(tái)測(cè)風(fēng)雷達(dá)同時(shí)保持正常工作狀態(tài))。因此，需要對(duì)湍流進(jìn)行研究并設(shè)計(jì)一種能在風(fēng)機(jī)組陣列的任意或部分位置，在獲得先驗(yàn)的湍流信息后，能夠根據(jù)即時(shí)風(fēng)速信息，對(duì)整個(gè)陣列風(fēng)場(chǎng)進(jìn)行估計(jì)的方法。

在20世紀(jì)70年代至90年代，激光雷達(dá)成為大氣湍流探測(cè)領(lǐng)域的一種重要手段，以美國(guó)海洋大氣中心和俄羅斯大氣環(huán)境中心等一部分組織和機(jī)構(gòu)對(duì)大氣湍流的性質(zhì)和特點(diǎn)開(kāi)展了一系列的相關(guān)研究[3]。在湍流測(cè)量分析領(lǐng)域的相關(guān)研究中，2008年Smalikho等[4]利用徑向速度譜來(lái)確定湍流的耗散率。2016年Matsuo等[5]利用神經(jīng)網(wǎng)絡(luò)對(duì)湍流信號(hào)進(jìn)行分析，克服了傳統(tǒng)VPP、VAD算法在湍流分析領(lǐng)域的缺陷。徐春鳳等[6]針對(duì)湍流大氣中激光光強(qiáng)的閃爍特性進(jìn)行了研究。2017年，Zhai 等[7]利用激光雷達(dá)對(duì)由于地面高度起伏和摩擦造成的近地表湍流進(jìn)行了研究。

本文基于風(fēng)力發(fā)電機(jī)組中常見(jiàn)的近地表湍流和渦流進(jìn)行研究，在結(jié)合現(xiàn)有non-Kolmogorov湍流模型對(duì)湍流的脈動(dòng)風(fēng)速特點(diǎn)進(jìn)行分析后，提出一種利用獨(dú)立分量分析(Independent Component Analysis，ICA)對(duì)常見(jiàn)的湍流大氣中回波信號(hào)的相位噪聲譜進(jìn)行分離的方法。

2 激光測(cè)風(fēng)雷達(dá)的設(shè)計(jì)參數(shù)與湍流

頻譜特點(diǎn)

2.1 激光雷達(dá)的設(shè)計(jì)參數(shù)

在風(fēng)力發(fā)電機(jī)組中，激光測(cè)風(fēng)雷達(dá)往往針對(duì)風(fēng)機(jī)前方50～80 m處的近場(chǎng)風(fēng)速信息進(jìn)行探測(cè)。在這一探測(cè)范圍內(nèi)，采用連續(xù)激光雷達(dá)可以獲得更高的探測(cè)精度和處理速度。本系統(tǒng)中采用的是由窄線寬半導(dǎo)體激光器作為激光光源，經(jīng)放大器后由收發(fā)一體式透鏡組雷達(dá)天線發(fā)射并在目標(biāo)內(nèi)獲得粒子的后向散射信號(hào)，再經(jīng)過(guò)Mach-Zehnder干涉系統(tǒng)將其多普勒頻移信息提取并進(jìn)行處理，最終獲得相關(guān)的風(fēng)速解算信號(hào)，其原理圖及設(shè)計(jì)參數(shù)如圖1所示。

圖1 激光測(cè)風(fēng)雷達(dá)原理Fig.1 Block diagram of wind lidar system

圖中，LD為半導(dǎo)體激光器，ISO為光隔離器，EDFA為激光放大器，AOM為聲光移頻器，circulator為環(huán)形器，antenna為雷達(dá)天線，PD為光電探測(cè)器，ESA為電頻譜儀。在測(cè)量精度方面，雷達(dá)系統(tǒng)采用定時(shí)轉(zhuǎn)盤速度校準(zhǔn)方法及回光信噪比監(jiān)測(cè)方法，確保測(cè)量精度在0.5 m/s，系統(tǒng)穩(wěn)定時(shí)最優(yōu)精度為0.1 m/s。

表1 激光測(cè)風(fēng)雷達(dá)的設(shè)計(jì)參數(shù)Tab.1 Parameter designed for windlidar

2.2 風(fēng)機(jī)附近風(fēng)場(chǎng)湍流模型

通常情況下在風(fēng)力發(fā)電機(jī)組之間，以尾流為主導(dǎo)的混合湍流中往往包含一些較為復(fù)雜的氣流擾動(dòng)因素，如地形湍流等。大多數(shù)情況下，隨著氣象條件的穩(wěn)定，尾流在湍流中的影響會(huì)愈發(fā)明顯，因此在以尾流為主的大氣湍流的研究中，常用基于大氣凍結(jié)湍流假設(shè)條件下的non-Kolmogorov模型對(duì)湍流特性進(jìn)行建模。其中，大氣凍結(jié)湍流假設(shè)在時(shí)間尺度不超過(guò)100 ms的前提條件下可以認(rèn)為是有效的[8]，其過(guò)程表示為：

φ(z，t)=φ(z-ωtd，t-td)，

(1)

其中：td為系統(tǒng)延遲，ω為橫向風(fēng)。這一過(guò)程表明，在特定的系統(tǒng)時(shí)延下，湍流風(fēng)場(chǎng)中的某一點(diǎn)空間特性可以認(rèn)為保持不變。在這一假設(shè)的基礎(chǔ)上，可以對(duì)湍流風(fēng)場(chǎng)采用non-Kolmogorov模型來(lái)描述它的空間分布特性。這一模型具有可變的大氣指數(shù)α值，且其功率譜函數(shù)包括non-Kolmogorov、廣義von Karman譜等[9]。

其功率譜密度函數(shù)如下：

(2)

由式(2)建立湍流信號(hào)仿真模型，并得到如圖 所示的仿真結(jié)果，通過(guò)調(diào)節(jié)折射率結(jié)構(gòu)常數(shù)及湍流尺寸可以獲得不同的湍流信號(hào)分布情況?？梢钥闯?，隨著湍流頻率分布范圍的擴(kuò)大以及湍流整體功率的提高，信號(hào)的分布范圍逐漸擴(kuò)大，同時(shí)信噪比下降，頻譜峰值穩(wěn)定程度下降，最終被噪聲淹沒(méi)。

(a)低湍流強(qiáng)度頻譜 (a)Moderate turbulence spectrum

(b)中湍流強(qiáng)度頻譜 (b)Medium turbulence spectrum

(c)高湍流強(qiáng)度頻譜 (c)Severe turbulence spectrum圖2 不同湍流強(qiáng)度的仿真頻譜Fig.2 Simulated power spectra in different turbulence levels

圖2表明當(dāng)湍流增強(qiáng)時(shí)，系統(tǒng)的回光信號(hào)由于被湍流信道調(diào)制而出現(xiàn)整體質(zhì)量下降。但同時(shí)，目標(biāo)區(qū)域附近的湍流本身也是粒子運(yùn)動(dòng)所產(chǎn)生的真實(shí)后向散射信號(hào)，湍流的實(shí)時(shí)強(qiáng)度及其變化規(guī)律也是風(fēng)場(chǎng)所需要的重要數(shù)據(jù)。因此，針對(duì)湍流信號(hào)的功率譜提取方法也顯得更加重要。

3 獨(dú)立分量分析模型的構(gòu)建

獨(dú)立分量分析的概念可以追溯到20世紀(jì)90年代，由法國(guó)學(xué)者Jutten等首次提出了獨(dú)立分量分析(Independent Component Analysis，ICA)概念，該方法在90年代中期獲得國(guó)際信號(hào)處理界的廣泛關(guān)注。

假設(shè)利用3個(gè)探測(cè)器記錄信號(hào)xi(t)(i=1，2，3)為3個(gè)信號(hào)來(lái)源si(t)(i=1，2，3)的一個(gè)線性加權(quán)[11]，則它可用矩陣形式表示為：

(3)

其中hij(i，j=1，2，3)為信道權(quán)系數(shù)，其矩陣組合為信道模型。當(dāng)信號(hào)來(lái)源增多時(shí)，可以利用更多的探測(cè)器對(duì)信號(hào)進(jìn)行探測(cè)與辨識(shí)，同樣用矩陣形式表示，則有：

(4)

其展開(kāi)式為：

(5)

這種多源信號(hào)混合分析問(wèn)題被稱為“雞尾酒會(huì)”問(wèn)題，如圖3所示。在忽略時(shí)延、非線性效應(yīng)等因素的影響后，系統(tǒng)可視為最簡(jiǎn)單的混合系統(tǒng)-線性瞬時(shí)混合系統(tǒng)。而構(gòu)建獨(dú)立分量分析模型的目的就是僅利用觀測(cè)數(shù)據(jù)估計(jì)出或恢復(fù)出未知源信號(hào)的大部分或全部特征。

圖3 多源信號(hào)的“雞尾酒會(huì)”模型Fig.3 Cocktail party effect of multisource signal

為了獲得信號(hào)源s的一個(gè)良好估計(jì)y=Wx，引入目標(biāo)函數(shù)-損失函數(shù)ρ(y，W)，并設(shè)其數(shù)學(xué)期望為：

R(W)=E{ρ(y，W)}，

(7)

其中R(W)為風(fēng)險(xiǎn)函數(shù)。由于相互獨(dú)立性的測(cè)度可以有多個(gè)衡量標(biāo)注，在信號(hào)中最常用的如KL散度。其定義為：

(8)

在ICA問(wèn)題中，q(y)是y的邊緣概率密度的乘積，可以表示為：

(9)

其中pi(yi)即為所求的邊緣概率密度，即：

(10)

(11)

其中ρ(y，W)=log(py(y)/q(y))。因此，可以獲得一個(gè)便于求解的獨(dú)立性測(cè)度R(W)。

為使代價(jià)函數(shù)(如R(W))最小化，可以利用對(duì)隨機(jī)梯度搜尋最優(yōu)解的思想獲得獨(dú)立分量分析的自然梯度模型：

(12)

(13)

其中激活函數(shù)向量f(y)=[f1(y1)，f2(y2)，f3(y3)，…，fn(yn)]T，為列向量，其單個(gè)分量為：

(14)

其中qi(yi)為源信號(hào){si}的近似概率密度函數(shù)。

(15)

為了配合算法的實(shí)時(shí)性，結(jié)合式(13)后，式(15)可改寫為實(shí)時(shí)自然梯度算法如式(16)所示：

(16)

利用式(16)便可以從混合信號(hào)中通過(guò)迭代逐步分離出不同類型的信號(hào)，進(jìn)而還原出整個(gè)信號(hào)的特點(diǎn)。

4 湍流頻譜估計(jì)與中心風(fēng)速分離

在利用獨(dú)立分量分析對(duì)湍流進(jìn)行頻譜分解前，首先要對(duì)獲得的信號(hào)進(jìn)行空域解處理。在這一步預(yù)處理中，分子的熱運(yùn)動(dòng)和粒子的布朗運(yùn)動(dòng)被去相關(guān)化[13]。由于在最終的頻譜分析中，這兩者并不是本文關(guān)注的對(duì)象，因此這一處理會(huì)大幅降低二者對(duì)后續(xù)頻譜分離與風(fēng)速估計(jì)的影響，從而提高分離與估計(jì)的準(zhǔn)確度。

(17)

由輸入信號(hào)變?yōu)榘谆幚硇盘?hào)的過(guò)程則可由式(18)表示：

(18)

其中Q矩陣稱為白化矩陣。當(dāng)信號(hào)在信道傳輸過(guò)程中出現(xiàn)混疊時(shí)，白化信號(hào)矩陣Q由n×n變?yōu)閚×m(m

圖4 白化處理矩陣結(jié)構(gòu)Fig.4 Framework of whiten matrix

因此，經(jīng)由白化矩陣Q處理之后的信號(hào)矩陣如下：

(19)

顯然，對(duì)這一白化矩陣Q左乘任何一個(gè)正交矩陣，式(19)結(jié)果不變，因此矩陣Q對(duì)信號(hào)處理而言并非唯一解。

由于零均值信號(hào)的相關(guān)矩陣Rxx=E{xxT}通常為對(duì)稱正定矩陣，因此，根據(jù)矩陣?yán)碚摽煞纸鉃椋?/p>

(20)

其中Vx為正交矩陣，Λx=diag{λ1，λ2，λ3，…，λn}為對(duì)角矩陣，對(duì)應(yīng)特征值由大到小排列，且由于其為正定矩陣，故所有特征值均大于等于0。

因此，白化矩陣Q的一種求法可以表示為：

(21)

也可以表示為：

(22)

假設(shè)兩個(gè)信號(hào)s1和s2是相互獨(dú)立的，如圖5(a)所示，將兩信號(hào)合成后，其結(jié)果如圖5(b)所示。此時(shí)兩信號(hào)x1和x2含有相關(guān)成分，若直接帶入后續(xù)處理中，則會(huì)產(chǎn)生混疊影響估計(jì)質(zhì)量，而通過(guò)預(yù)白化處理后可以達(dá)到如圖5(c)所示的效果，圖中y1與y2分量與s1和s2分量在能量分布上均可以視為獨(dú)立的。與s1和s2的區(qū)別僅為分布位置的變化，而與信號(hào)自身的獨(dú)立性無(wú)關(guān)。

綜上，白化預(yù)處理可以使具有相關(guān)性的混合信號(hào)解相關(guān)，并削弱非相關(guān)信號(hào)對(duì)混合信號(hào)的影響，在此基礎(chǔ)上進(jìn)行的頻譜分離則具有更高的穩(wěn)定性和準(zhǔn)確度。

(a)獨(dú)立信號(hào)成分(a)Independent signals(b)獨(dú)立成分混合(b)Mixed signals(c)混合成分白化(c)Whiten signals圖5 混合獨(dú)立分量與白化處理示意圖Fig.5 Schematic diagram of mixed and whiten signals

5 測(cè)量實(shí)驗(yàn)與結(jié)果

基于之前的理論分析，對(duì)ICA在風(fēng)速處理中的實(shí)際表現(xiàn)進(jìn)行了測(cè)試，在26 ℃恒溫地下實(shí)驗(yàn)室搭建了激光測(cè)風(fēng)雷達(dá)光路并結(jié)合電路部分，如圖6所示。圖7所示為激光測(cè)風(fēng)雷達(dá)系統(tǒng)樣機(jī)，該樣機(jī)主要包括基于伽利略望遠(yuǎn)鏡結(jié)構(gòu)設(shè)計(jì)的雙目激光雷達(dá)天線，基于Mach-Zehnder干涉系統(tǒng)設(shè)計(jì)的延時(shí)自外差光纖測(cè)量光路，經(jīng)光電轉(zhuǎn)換后轉(zhuǎn)變?yōu)镸Hz級(jí)電信號(hào)，再經(jīng)數(shù)字電路簡(jiǎn)單處理后可得原始信號(hào)(Raw Signal)，將原始信號(hào)傳遞至上位機(jī)，即可利用ICA方法對(duì)原始信號(hào)進(jìn)行湍流分離與風(fēng)場(chǎng)估計(jì)。

圖6 測(cè)風(fēng)雷達(dá)光路部分Fig.6 Optical system of wind lidar

圖7 測(cè)風(fēng)雷達(dá)系統(tǒng)樣機(jī)Fig.7 Prototype of lidar

圖8 混合風(fēng)速測(cè)量結(jié)果Fig.8 Raw signal in wind measurement

首先利用轉(zhuǎn)盤進(jìn)行速度校準(zhǔn)，通過(guò)控制轉(zhuǎn)盤速度及測(cè)量位置在0～3 m/s，每隔0.5 m/s采樣10組對(duì)回光信號(hào)進(jìn)行校準(zhǔn)，其結(jié)果表2所示。

表2 轉(zhuǎn)盤轉(zhuǎn)速與回光頻率關(guān)系Tab.2 Relationship between disk speed andfrequency of backscatter light

在校準(zhǔn)中發(fā)現(xiàn)，1 MHz以下部分波動(dòng)劇烈，頻率校準(zhǔn)方差較大且與實(shí)際速度關(guān)聯(lián)性弱，因此考慮它為光路內(nèi)反饋或電磁干擾等形成的噪聲。因此利用高通濾波器對(duì)噪聲進(jìn)行濾除，隨后利用擋板障礙物與風(fēng)扇組合構(gòu)建復(fù)雜風(fēng)場(chǎng)環(huán)境。每獲得一組1 000×1 024個(gè)有效風(fēng)速數(shù)據(jù)點(diǎn)(約10 min)，之后對(duì)測(cè)風(fēng)雷達(dá)進(jìn)行一次定速轉(zhuǎn)盤速度校準(zhǔn)，在校準(zhǔn)結(jié)果誤差小于5%時(shí)繼續(xù)測(cè)量，否則停機(jī)進(jìn)行進(jìn)一步校驗(yàn)并放棄前一組數(shù)據(jù)?；旌巷L(fēng)速測(cè)量結(jié)果如圖8所示。在獲得有效數(shù)據(jù)后對(duì)復(fù)雜風(fēng)場(chǎng)回光信號(hào)進(jìn)行預(yù)白化處理，得到的白化信號(hào)如圖9所示。在白化過(guò)程中，雙目信號(hào)中具有相關(guān)性的部分被逐步分離，通過(guò)解相關(guān)運(yùn)算，強(qiáng)相關(guān)部分(可認(rèn)為是湍流或風(fēng)速的“中心”)被保留，而弱相關(guān)部分(可認(rèn)為是湍流或風(fēng)速的“邊緣”)則轉(zhuǎn)化為本底噪聲。這樣的信號(hào)在ICA處理中速度更快且收斂性更好。

圖9 原始信號(hào)與白化及濾波處理后的信號(hào)Fig.9 Raw and whiten signals

利用梯度下降的ICA方法可以將雙目信號(hào)分解為相互獨(dú)立的信號(hào)，在頻譜中表現(xiàn)為兩個(gè)信號(hào)峰，如圖9所示，其中ICA Processed曲線即為ICA處理結(jié)果，其頻譜中的獨(dú)立峰參數(shù)如表3中所示。

表3 雙目信號(hào)ICA方法在風(fēng)速分離處理中頻譜獨(dú)立峰的平均參數(shù)

Tab.3 Main parameters of two peaks measured by ICA processing

如表3所示，雖然在實(shí)際的風(fēng)速信號(hào)處理中，ICA方法的表現(xiàn)不如模擬信號(hào)，但是同樣可以獲得較為穩(wěn)定且信噪比滿足要求的頻譜峰值。結(jié)合風(fēng)速的實(shí)際情況與包含湍流的風(fēng)場(chǎng)特性，以高斯函數(shù)和正弦函數(shù)作為先驗(yàn)信息中的參考進(jìn)行獨(dú)立分量的分離，得到ICA Processed 曲線對(duì)應(yīng)的分離結(jié)果，如圖10所示。將一組中所有數(shù)據(jù)分為50個(gè)小組進(jìn)行50次ICA分離，其結(jié)果穩(wěn)定性如表4所示。由此可以看出，在可用數(shù)據(jù)量減少時(shí)，ICA方法處理中雖然會(huì)出現(xiàn)少量錯(cuò)誤結(jié)果，但正確辨識(shí)率依舊能夠保證在90%以上。

圖10 原始信號(hào)、白化信號(hào)與ICA方法處理后的信號(hào)Fig.10 Raw， whiten and ICA processed signals

表4 50次ICA分離結(jié)果統(tǒng)計(jì)
Tab.4 ICA separation results for 50 tests

分離結(jié)果次 數(shù)正確識(shí)別率/%完全分離46中心風(fēng)速與渦流混合192.0未解混2無(wú)法辨識(shí)的混合結(jié)果1

圖11 ICA方法所獲得風(fēng)速頻譜峰值的1 s穩(wěn)定性曲線Fig.11 Two peaks values of wind frequency spectra measured by ICA processing in 1 s

對(duì)一組中所有數(shù)據(jù)進(jìn)行統(tǒng)計(jì)可以得到頻域穩(wěn)定性情況，如圖11所示。低頻區(qū)域曲線穩(wěn)定性較差，而上方曲線穩(wěn)定性較好，因此可以初步判斷上方曲線為中心風(fēng)速，而下方曲線為湍流中心的可能性較大。若考慮到低頻區(qū)域會(huì)受到更大的噪聲干擾時(shí)，可以結(jié)合主成份分析或自適應(yīng)濾波等方法作為輔助信號(hào)處理方法進(jìn)一步判斷中心風(fēng)速的頻率。

6 結(jié) 論

本文應(yīng)用ICA方法，對(duì)雙目激光測(cè)風(fēng)雷達(dá)回光信號(hào)的頻譜進(jìn)行分解，可以從中分離出兩個(gè)獨(dú)立信號(hào)，得到信噪比分別為25.93和31.01 dB的湍流信號(hào)和中心風(fēng)速信號(hào)估計(jì)，通過(guò)觀測(cè)其頻域穩(wěn)定性，認(rèn)為其中(2.59±0.05)MHz為中心風(fēng)速，能夠?qū)崿F(xiàn)在湍流大氣中捕捉穩(wěn)定中心風(fēng)速信號(hào)的功能。