陸 鑫，陳志敏,2

(1.上海電機(jī)學(xué)院 電子信息學(xué)院，上海 201306；2.東南大學(xué) 毫米波國(guó)家重點(diǎn)實(shí)驗(yàn)室，南京 210096)

0 引言

生命探測(cè)雷達(dá)是現(xiàn)代雷達(dá)技術(shù)和生物醫(yī)學(xué)工程技術(shù)結(jié)合的產(chǎn)物，生命探測(cè)雷達(dá)在軍事反恐、災(zāi)難營(yíng)救、醫(yī)療檢測(cè)等領(lǐng)域有廣泛應(yīng)用[1-3]，在航空航天領(lǐng)域也扮演重要的角色[4-5]。探測(cè)雷達(dá)通過(guò)電磁波穿透非金屬遮擋介質(zhì)等障礙物，探測(cè)到人體的生命特征，如呼吸、心跳、肢體動(dòng)作等,通過(guò)處理雷達(dá)回波信號(hào)，可以實(shí)現(xiàn)對(duì)飛行員的呼吸心跳的監(jiān)測(cè)。

雷達(dá)接收機(jī)接收到的回波信號(hào)除了含有生命體的特征信號(hào)，還包括非金屬介質(zhì)反射的強(qiáng)回波，以及周?chē)h(huán)境帶來(lái)的干擾和噪聲。另一方面雷達(dá)接收信號(hào)具有非平穩(wěn)、隨機(jī)性強(qiáng)等特點(diǎn)，從中提取微弱的生命信號(hào)非常困難。其次，雷達(dá)接收信號(hào)經(jīng)預(yù)處理，再經(jīng)A/D轉(zhuǎn)換后數(shù)據(jù)量較大，難以滿足現(xiàn)場(chǎng)實(shí)時(shí)、快速探測(cè)的需求。因此，現(xiàn)有的生命探測(cè)雷達(dá)仍存在許多技術(shù)難點(diǎn)需要解決。

首先，針對(duì)回波信號(hào)具有非平穩(wěn)性、非周期性和隨機(jī)性強(qiáng)等特點(diǎn)，可以采用卡爾曼濾波法、小波變換、經(jīng)驗(yàn)?zāi)B(tài)分解(empirical mode decomposition，EMD)和變分模態(tài)分解(variationalmode decomposition，VMD)等方法來(lái)處理。其中，卡爾曼濾波以最小均方誤差為最佳估計(jì)準(zhǔn)則來(lái)尋求遞推估計(jì)的算法，能夠?qū)r(shí)變系統(tǒng)、非平穩(wěn)信號(hào)、多維信號(hào)進(jìn)行處理，但當(dāng)噪聲為統(tǒng)計(jì)特性未知的有色噪聲或系統(tǒng)具有不確定性時(shí)，該方法會(huì)造成極大的誤差[6]；小波變換通過(guò)多尺度分析把信號(hào)分解為不同的頻率分量，但存在小波基函數(shù)的選擇問(wèn)題，適應(yīng)性較差[7]；而EMD算法在分解時(shí)會(huì)因循環(huán)遞歸迭代產(chǎn)生混疊，且分解的層數(shù)是隨機(jī)的，進(jìn)行分解時(shí)無(wú)法手動(dòng)確定不同模態(tài)的頻率范圍，缺乏理論上的證明，存在端點(diǎn)效應(yīng)問(wèn)題[8-10]。

VMD作為一種新的自適應(yīng)信號(hào)處理方法，運(yùn)算效率高，可克服EMD中的模態(tài)混疊問(wèn)題，實(shí)現(xiàn)信號(hào)的準(zhǔn)確分離，利用其自身具有的維納濾波特性可獲得更優(yōu)的噪聲濾除效果。該算法重點(diǎn)在于如何確定固有模態(tài)分量(intrinsic mode function，IMF)的分解層數(shù)k和懲罰因子α的個(gè)數(shù)。若k值過(guò)小會(huì)產(chǎn)生模態(tài)混疊，過(guò)大會(huì)導(dǎo)致過(guò)度分解；若懲罰因子α取值過(guò)大，則會(huì)造成分解的模態(tài)函數(shù)頻帶過(guò)窄，丟失有用信息，反之，會(huì)造成頻帶過(guò)寬而攜帶干擾信息。針對(duì)該問(wèn)題，本文提出一種基于包絡(luò)熵的粒子群(particle swarm optimization，PSO)算法對(duì)VMD算法進(jìn)行改進(jìn)以確定上述兩個(gè)參數(shù)，實(shí)現(xiàn)提取雷達(dá)接收機(jī)回波信號(hào)中生命特征的目的。

本文核心內(nèi)容主要包括4部分，第1部分給出了生命探測(cè)雷達(dá)的系統(tǒng)模型；第2部分提出一種基于包絡(luò)熵的PSO-VMD生命探測(cè)算法，并給出算法的具體實(shí)現(xiàn)；第3部分對(duì)本文所提出的算法進(jìn)行實(shí)驗(yàn)驗(yàn)證和仿真對(duì)比分析；第4部分主要是對(duì)本文工作的總結(jié)。

1 生命探測(cè)雷達(dá)回波信號(hào)建模

由于超寬帶微功率脈沖雷達(dá)穿透能力強(qiáng)，抗干擾能力強(qiáng)，探測(cè)靈敏度高，能探測(cè)呼吸和體動(dòng)等微弱的生命特征，因此被廣泛應(yīng)用于生命探測(cè)領(lǐng)域研究[11]。本文采用超寬帶雷達(dá)體制對(duì)生命信號(hào)進(jìn)行研究，重點(diǎn)關(guān)注人體生命信號(hào)的探測(cè)和辨識(shí)，不考慮人員的定位和數(shù)量?？紤]到人體呼吸和心跳引起的胸腔運(yùn)動(dòng)具有一定的周期性，可以用單頻正弦信號(hào)來(lái)近似表示。為簡(jiǎn)化分析，假設(shè)呼吸心跳引起的胸腔振動(dòng)是等距離振動(dòng)，即只需要考慮呼吸心跳的頻率特征，則呼吸心跳引起的胸腔振動(dòng)模型Z(t)可表示為：

Z(t)=Arsin(2πf1t+φ1)+Ahsin(2πf2t+φ2)

(1)

其中，Ar、Ah分別為人體呼吸和心跳引起的胸腔和脈搏運(yùn)動(dòng)幅度，f1、f2分別為呼吸和心跳的頻率，φ1、φ2分別為呼吸和心跳的初始相位。

對(duì)于頻率為fc的電磁波，將生命探測(cè)雷達(dá)與飛行器內(nèi)各種障礙目標(biāo)之間的距離記為Rn(n=1,2,…,N)，雷達(dá)與飛行員之間的距離記為R0，那么來(lái)自障礙物以及飛行員的回波信號(hào)可以表示為：

(2)

其中，An表示第n個(gè)目標(biāo)回波信號(hào)的幅度，fc為載波頻率，B表示飛行員回波幅度，Z(t)表示飛行員引起的胸腔距離變化。

那么，疊加噪聲后得到雷達(dá)回波信號(hào)R(t)為：

(3)

最終，根據(jù)式(3)，我們需要從得到的回波信號(hào)中，提取呼吸和心跳的頻率。接下來(lái)提出一種基于PSO-VMD的生命探測(cè)算法對(duì)接收信號(hào)進(jìn)行處理。

2 基于PSO-VMD的生命探測(cè)算法

2.1 信號(hào)預(yù)處理

生命探測(cè)雷達(dá)回波中包含了強(qiáng)背景噪聲與雜波等干擾，為消除干擾影響，有效提取生命體征信號(hào)，需要對(duì)接收信號(hào)進(jìn)行預(yù)處理[13-15]，具體如下：

(4)

2.2 PSO-VMD重構(gòu)雷達(dá)回波信號(hào)

在VMD的分解過(guò)程中，若分解層數(shù)k取值過(guò)大，則會(huì)造成過(guò)度分解，出現(xiàn)虛假模態(tài)；反之，會(huì)導(dǎo)致分解不充分，產(chǎn)生模態(tài)混疊。同時(shí)，若懲罰因子α取值過(guò)大，則會(huì)造成分解的模態(tài)函數(shù)頻帶過(guò)窄，丟失有用信息；反之，會(huì)造成頻帶過(guò)寬，攜帶干擾信息。對(duì)此，本文采用粒子群算法選擇包絡(luò)熵作為適應(yīng)度函數(shù)，來(lái)確定以上兩個(gè)參數(shù)的最優(yōu)組合，達(dá)到改進(jìn)VMD算法的目的。

2.2.1 VMD算法

VMD是一種新的自適應(yīng)信號(hào)處理方法，對(duì)非線性、非平穩(wěn)信號(hào)的處理具有明顯的優(yōu)勢(shì)。該方法運(yùn)算效率高，利用迭代搜索變分模型最優(yōu)解確定每個(gè)分解的分量中心頻率及帶寬。VMD算法可以抑制大部分的噪聲，有效降低信號(hào)的非平穩(wěn)性，并將生命體征信號(hào)分解為一系列的IMF分量。

VMD算法的實(shí)現(xiàn)核心主要分為構(gòu)造變分問(wèn)題和求解變分模型2個(gè)步驟[16-17]，具體如下。

1)構(gòu)造變分問(wèn)題

VMD算法是將原始信號(hào)f(t)分解為k個(gè)具有特定稀疏性的相互獨(dú)立的固有模態(tài)函數(shù)。首先，對(duì)uk(t)進(jìn)行希爾伯特變換，得到分解信號(hào)的單邊頻譜[(δ(t)+j/πt)uk(t)]。隨后引入e-jωkt對(duì)每一個(gè)模態(tài)中心頻率進(jìn)行調(diào)整，實(shí)現(xiàn)將頻譜調(diào)制到基帶。最后，計(jì)算每個(gè)解調(diào)信號(hào)梯度的平方L2范數(shù)，估算uk(t)寬度。受約束的變分問(wèn)題可描述為：

(5)

式中：δ(t)為狄利克來(lái)函數(shù)；uk，ωk分別為IMF具有有限帶寬的模態(tài)分量和中心頻率；f(t)為原始信號(hào)。

2)求解變分模型問(wèn)題

通過(guò)引入懲罰因子α和拉格朗日乘法算子λ，將約束性變分問(wèn)題變?yōu)榉羌s束性變分問(wèn)題。增廣拉格朗日表達(dá)式為：

(6)

最優(yōu)解采用二次懲罰和拉格朗日乘數(shù)將上訴約束問(wèn)題轉(zhuǎn)換為非約束問(wèn)題，并用交替方向乘子法(alternating direction method of multipliers，ADMM)求解這個(gè)非約束問(wèn)題，通過(guò)迭代更新最終得到信號(hào)分解的所有模態(tài)。分解的所有模態(tài)中有包含主要信號(hào)的模態(tài)和包含噪聲的模態(tài)。

2.2.2 PSO算法

本文采用粒子群優(yōu)化算法來(lái)確定VMD的懲罰因子α和分解層數(shù)k。粒子群算法是一種智能優(yōu)化算法，相比于其他算法，搜索速度快、算法簡(jiǎn)單、效率高。通過(guò)粒子狀態(tài)初始化，不斷尋優(yōu)和迭代比較個(gè)體最優(yōu)和全局最優(yōu)來(lái)調(diào)整粒子對(duì)環(huán)境的適應(yīng)度，搜索到最優(yōu)的區(qū)域的最優(yōu)解，PSO算法具體描述如下：

對(duì)于一個(gè)D維空間存在n個(gè)粒子組成的種群X，其中第i個(gè)粒子以速度為Vi=(vi1,vi2,…,viD)飛行，粒子在空間的位置為Xi=(xi1,xi2,…,xiD)，粒子的最佳位置記為pi=(pi1,pi2,…,piD)。其種群的最佳位置為pg=(pg1,pg2,…,pgD)。作為它的第d+1維的速度Vid+1與位置Xid+1的更新表達(dá)式為：

(7)

其中，vid為粒子速度：i∈[1,n]；d+1∈[1,D]，c1與c2為加速常數(shù)，ω為慣性權(quán)重，r1，r2是0到1的任意隨機(jī)數(shù)。為了確保搜索準(zhǔn)確，跳出局部獲取全局的最優(yōu)解，要將粒子的速度與位置控制在一定的范圍內(nèi)。

2.2.3 包絡(luò)熵

VMD算法雖然避免了EMD算法存在的模態(tài)混疊和端點(diǎn)效應(yīng)問(wèn)題，但由于VMD算法中的懲罰因子α和分解層數(shù)k是未給定的參數(shù)，而兩個(gè)參數(shù)的選擇又對(duì)VMD分解的影響很大。由于兩個(gè)參數(shù)一般是依靠人們的經(jīng)驗(yàn)確定，缺乏嚴(yán)格的證明，因此制約了VMD算法在各種工程中的應(yīng)用。信號(hào)的熵值反映了它的隨機(jī)程度和復(fù)雜程度，信號(hào)的周期性越明顯，復(fù)雜程度越低，熵值越小，反之非周期性噪聲干擾越多，信號(hào)越復(fù)雜，熵值越大[18]。因此本文采用包絡(luò)熵作為粒子群算法的適應(yīng)度函數(shù)。生命信號(hào)經(jīng)過(guò)VMD算法分解后的IMF分量的包絡(luò)熵可表示為：

(8)

式中，i為原信號(hào)分解得到的IMF層數(shù)；pi,j為ai(j)的歸一化形式；ai(j)為信號(hào)IMF經(jīng)希爾伯特變換后得到的包絡(luò)信號(hào)。當(dāng)信號(hào)經(jīng)VMD分解得到IMF分量包含的周期性生命體征信號(hào)越多，包絡(luò)熵較小；反之，周期性信號(hào)越少越稀疏，包絡(luò)熵值則越大。因此，為了得到周期性較好的生命體征信號(hào)，本文采用包絡(luò)熵的最小值作為PSO算法適應(yīng)度函數(shù)的判別方式。

2.2.4 基于PSO-VMD的參數(shù)尋優(yōu)算法

PSO-VMD算法引入包絡(luò)熵作為粒子群算法的適應(yīng)度函數(shù)，VMD分解計(jì)算適應(yīng)度函數(shù)值，通過(guò)對(duì)粒子速度與位置的更新，找到全局最佳的粒子速度與位置，確定VMD分解的分解層數(shù)和懲罰因子兩個(gè)參數(shù)。PSO-VMD算法確定兩個(gè)參數(shù)的搜索過(guò)程如下：

1)初始化PSO的參數(shù)并將包絡(luò)熵Ei的最小值作為粒子群算法的適應(yīng)度函數(shù)；

2)隨機(jī)產(chǎn)生粒子種群中的粒子位置(k,α)和粒子速度vid；

3)VMD分解計(jì)算信號(hào)在不同粒子位置(k,α)下每個(gè)粒子位置對(duì)應(yīng)的Ei值；

4)通過(guò)不斷對(duì)比Ei的大小更新個(gè)體極值和種群全局極值；

5)利用式(7)更新粒子的移動(dòng)速度和位置；

6)循環(huán)步驟3～5，當(dāng)?shù)螖?shù)達(dá)到預(yù)設(shè)最大值結(jié)束循環(huán)并輸出最佳粒子位置(k,α)。

2.3 生命探測(cè)雷達(dá)處理方法

生命探測(cè)雷達(dá)接收到回波信號(hào)，首先濾除環(huán)境噪聲，再經(jīng)低噪放、混頻解調(diào)、濾波等信號(hào)預(yù)處理后，采用PSO-VMD算法進(jìn)行信號(hào)處理，以確定最佳的VMD分解懲罰因子和分解層數(shù)參數(shù)組合。其次，將所得到的參數(shù)組合作為VMD分解的特定參數(shù)，對(duì)回波信號(hào)進(jìn)行分解，得到一系列IMF分量。接著對(duì)各個(gè)IMF分量進(jìn)行頻譜分析，選擇能夠反映生命體征信號(hào)的IMF分量重構(gòu)信號(hào)。最后，對(duì)重構(gòu)的呼吸心跳信號(hào)進(jìn)行傅里葉變換，完成對(duì)呼吸心跳等生命特征信號(hào)頻率的提取，最終實(shí)現(xiàn)對(duì)人體呼吸心跳的監(jiān)測(cè)。具體生命探測(cè)流程如圖1所示：

圖1 生命探測(cè)雷達(dá)信號(hào)處理流程Fig.1 Block diagram of life detection radar

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

正常人的呼吸頻率為每分鐘12～44次，頻率為0.2～0.8 Hz，呼吸引起的胸腔運(yùn)動(dòng)幅度約為1～30 mm；心跳頻率每分鐘60～120次，頻率為1～2 Hz，心跳引起的胸腔運(yùn)動(dòng)幅度約為0.1～10 mm。本文假定f1=0.3 Hz，Ar=8 mm，φ1=0°；f2=1.3 Hz，Ah=2 mm，φ2=60°，雷達(dá)距離障礙物的距離為R1=3 m，距離飛行員R0=5.1 m，呼吸心跳的引起的胸腔振動(dòng)模型Z(t)可表示為：

Z(t)=0.008sin(2π·0.3t)+
0.002sin(2π·1.3t+π/3)

(9)

雷達(dá)生命體征信號(hào)為：

(10)

其中，A1表示雷達(dá)穿透第一個(gè)障礙物回波信號(hào)的幅度，B表示飛行員回波幅度。

當(dāng)雷達(dá)回波經(jīng)過(guò)正交解調(diào)，得到同相分量I路基帶信號(hào)和正交分量Q路基帶信號(hào)：

(11)

其中，a是接收信號(hào)的幅度，fc為發(fā)射信號(hào)載頻，Δφ(t)為總的剩余相位噪聲。本文取a=1，fc=1 GHz，不考慮總的剩余相位噪聲，對(duì)Q路基帶信號(hào)進(jìn)行分析。生命體征信號(hào)和雷達(dá)回波信號(hào)的模擬仿真圖如圖2所示：

圖2 模擬仿真圖Fig.2 Simulation diagram

將回波信號(hào)用于優(yōu)化VMD參數(shù)的粒子群算法參數(shù)設(shè)置為kmim=3，kmax=10，迭代速度為1；αmin=500，αmax=3 000，迭代速度為10；初始種群規(guī)模sizepop=30，迭代次數(shù)設(shè)置maxgen=20。

用本文提出的PSO-VMD算法，取實(shí)驗(yàn)5組數(shù)據(jù)得到的兩個(gè)最優(yōu)參數(shù)如表1所列：

表1 PSO-VMD方法得到的兩個(gè)參數(shù)組合

根據(jù)表1可以發(fā)現(xiàn)最優(yōu)分解層數(shù)主要集中在8層，因此采用在分解層數(shù)8層下最接近5組懲罰因子平均數(shù)的2 374為最優(yōu)的懲罰因子數(shù)。最后，選擇參數(shù)組合(2 374，8)作為此類(lèi)生命信號(hào)進(jìn)行VMD分解的最優(yōu)參數(shù)設(shè)置。生命信號(hào)采用PSO-VMD算法得到參數(shù)VMD分解下的各個(gè)IMF分量時(shí)域和頻域圖如圖3所示：

圖3 各IMF分量時(shí)域及對(duì)應(yīng)的頻域圖Fig.3 The time and frequency domain of each IMF component

由圖3可以看出，在8個(gè)IMF分量中，IMF1包含了比較豐富的低頻的生命體征信息，因此，將IMF1作為重構(gòu)的生命信號(hào)，并對(duì)其進(jìn)行頻譜分析，結(jié)果如圖4所示：

圖4 重構(gòu)信號(hào)頻譜圖Fig.4 Spectrum diagram of reconstructed signal

從圖4中可以看出，重構(gòu)的生命信號(hào)中包含了仿真假定的心跳和呼吸頻率，說(shuō)明采用本文所提出的PSO-VMD算法可以有效檢測(cè)出微弱的生命體信號(hào)。

在相同噪聲條件下，進(jìn)一步仿真對(duì)比VMD算法和EMD算法性能。EMD算法重構(gòu)生命體信號(hào)與PSO-VMD算法類(lèi)似，不同在于PSO-VMD只需要選擇低頻段的固有模態(tài)分量(IMF1)即可重構(gòu)信號(hào)，而EMD要根據(jù)每個(gè)固有模態(tài)分量頻譜低頻分量選擇總能量占比大于70%來(lái)重構(gòu)生命體信號(hào)[19-20]，顯然要比PSO-VMD繁瑣。圖5給出了相同噪聲環(huán)境下，VMD和EMD算法提取呼吸心跳頻率的效果。可以看出，VMD算法和EMD算法都可以準(zhǔn)確提取呼吸頻率，然而EMD方法提取的心跳頻率淹沒(méi)在噪聲諧波中，無(wú)法準(zhǔn)確得出心跳頻率特征。這是因?yàn)镋MD算法因復(fù)雜噪聲影響包絡(luò)線的計(jì)算，導(dǎo)致計(jì)算不斷擴(kuò)大誤差，出現(xiàn)了模態(tài)混疊的現(xiàn)象，而VMD算法具有較好的噪聲魯棒性，復(fù)雜噪聲對(duì)其影響不明顯。所以在相同的噪聲條件下，PSO-VMD算法提取呼吸心跳頻率比EMD算法更準(zhǔn)確有效。

圖5 兩種方法的提取心跳呼吸頻率對(duì)比圖Fig.5 Comparison of extraction heartbeat and respiratory rate of the two methods

4 結(jié)論

本文介紹了一種基于PSO-VMD參數(shù)優(yōu)化的生命探測(cè)方法。首先，采用包絡(luò)熵作為PSO算法的適應(yīng)度函數(shù)，有效提取雷達(dá)接收機(jī)回波中呼吸心跳頻率的生命體征信號(hào)。其次，采用PSO算法調(diào)整參數(shù)重構(gòu)信號(hào)，尋找VMD最優(yōu)的分解層數(shù)k和懲罰因子α，以避免因k和α不準(zhǔn)確帶來(lái)的模態(tài)模糊等問(wèn)題。最后將本文所提出的算法與EMD算法進(jìn)行了仿真對(duì)比，結(jié)果顯示，PSO-VMD算法重構(gòu)信號(hào)提取呼吸心跳頻率更加準(zhǔn)確、快速和高效。今后將進(jìn)一步深入研究效果更佳的優(yōu)化算法，同時(shí)結(jié)合壓縮感知、稀疏重構(gòu)等先進(jìn)算法進(jìn)一步提高探測(cè)性能。