紀騰飛 楊鋮 趙毅

（1.天津大學，天津 300072；2.青島中科創(chuàng)電子技術有限公司，青島 266200）

引 言

5G和6G通信系統(tǒng)將毫米波頻段作為了主要頻段之一，使得毫米波通信機理研究、傳播特性分析等領域成為?近年來國內(nèi)外廣泛關注和研究的熱點[1-6].毫米波通信借助毫米波為傳輸信息的載體進行通信，其載波頻率覆蓋范圍為26.5～300 GHz，帶寬最高可達273.5 GHz[7]. 由于毫米波屬于極高頻段，干擾源很少，可以穩(wěn)定地在室內(nèi)或室外環(huán)境中傳播；同時在相同天線尺寸下毫米波的波束要比微波窄很多，可以分辨相距更近的小目標. 因此，毫米波通信是一種典型的具有高質量、恒定參數(shù)的無線傳輸信道的通信技術[8]. 室內(nèi)毫米波信號傳播環(huán)境復雜，多徑效應明顯，需要精確的場景模型預測信號覆蓋，以便提高信號服務質量和頻譜效率[9]. 室內(nèi)無線通信系統(tǒng)的研究與開發(fā)需要對信道特性[10-11]有廣泛的了解，射線追蹤及其改進模型是預測室內(nèi)信道特性的一種有效方法[11]，主要成果有：1) 文獻[12]提出了一種適用于辦公室會議環(huán)境的信道模型. 該模型考慮了傳輸損耗、信道聚類以及極化等影響，利用射線追蹤法驗證了信道的統(tǒng)計模型. 2) 文獻[13]提出了一種基于圖形處理器(graphics processing unit, GPU)的KDtree加速波束跟蹤方法. 該方法包含了一種基于插接器坐標的有效波束-三角形相交算法和一種擴展到波束跟蹤的KD-tree遍歷算法，克服了射線追蹤方法中的離散采樣效應，大大減少了計算時間. 3) 文獻[14]采用射線追蹤法檢測了室內(nèi)辦公環(huán)境中毫米波信道參數(shù)，并與空間交替廣義期望最大化(spatial-alternation generalized expectation-maximization，SAGE)估 計 算法處理的信道測量數(shù)據(jù)進行比較，發(fā)現(xiàn)在視距和非視距兩種情況下，優(yōu)勢路徑的模擬性能和實測性能都能達到很好的一致性. 4) 德國AWE Communications GMBH公司開發(fā)了WinProp室內(nèi)和室外傳播分析軟件. 該工具強大的傳播引擎包括經(jīng)驗和半經(jīng)驗模型以及三維射線跟蹤模型和獨有的優(yōu)勢路徑模型，其預測結果與實際測量結果能夠保持著很好的一致性，也使該工具成為了當今無線電波傳播與無線網(wǎng)絡規(guī)劃領域內(nèi)的標準工具[15]. 5)文獻[16]提出了一種快速射線追蹤算法，在對數(shù)字高程模型(digital elevation model，DEM)數(shù)據(jù)進行三角形剖分的基礎上，利用射線跟蹤方法對不規(guī)則地形進行了大量的覆蓋預測，分析了場分布和信道特性. 6)文獻[17]介紹了一種通用的射線追蹤模型，并根據(jù)我國的實際情況給出了相關參數(shù)的獲取方法. 仿真證明算法具有實用性強、計算精度高的特點. (7)文獻[18]提出了一種改進的空間分割射線追蹤算法跟蹤電磁波在海面上的傳播路徑，對無人機與船之間的移動信道進行了仿真和建模.

由于射線在復雜的室內(nèi)環(huán)境中傳播時間長，計算復雜度高，導致室內(nèi)毫米波研究受到了很大限制.因此，為了滿足工程需要，急需改善傳統(tǒng)的射線追蹤技術，提高算法效率. 國內(nèi)學者利用捷徑法[19]、鄰域輔助徑跡法[20]、優(yōu)化射線追蹤場[21]及改善接收球方法[22]等在一定程度上提升了射線追蹤效率. 為在保障預測精度的基礎上進一步提升算法效率，本文利用蟻群分布計算和啟發(fā)式搜索原理尋找電波傳播優(yōu)勢射線路徑，選取室內(nèi)傳播中能量較高的優(yōu)勢射線路徑代替所有射線路徑，用以提高射線追蹤的計算效率，降低算法復雜度.

1 基于蟻群算法的優(yōu)勢路徑確定

傳統(tǒng)的射線追蹤法是基于幾何光學和均勻衍射理論[7]通過射線尋跡過程，尋找從發(fā)射源發(fā)出并到達接收點的所有可能傳播路徑，計算反射、繞射和透射系數(shù)，將各條路徑的射線能量合并，從而得到接收點電平和傳播延時. 傳統(tǒng)的射線追蹤法在測試場景中發(fā)射大量射線，預測射線傳播的信道特性，計算效率低下，復雜度極高[23].

針對不同復雜的室內(nèi)環(huán)境，本文引入蟻群算法用于選取室內(nèi)傳播中能量較高的幾條優(yōu)勢射線路徑代替所有射線路徑，去除大量貢獻小的射線，大幅度降低了射線在室內(nèi)環(huán)境傳播的復雜度. 此外，在此基礎上，通過設置一定的誤差閾值，結合計算機多線程并行計算技術，在保證精度的前提下，可以大大提高射線追蹤的追蹤速度，減小計算復雜度. 快速射線追蹤與傳統(tǒng)射線追蹤對比如表1所示.

表1 快速射線追蹤與傳統(tǒng)射線追蹤對比Tab. 1 Comparison of fast ray tracing and traditional ray tracing processes

蟻群算法的研究模型源于對真實螞蟻覓食行為的模擬. 螞蟻在覓食過程中，會在所經(jīng)路徑上釋放出一種具有揮發(fā)性的物質——信息素，不同螞蟻個體通過感知信息素的存在及其強度指導自己的移動方向. 研究表明，螞蟻更傾向于選擇信息量較大的路徑，由此形成一種正反饋機制：最優(yōu)路徑上的信息量越來越大，其他路徑上的信息量則隨時間逐漸衰減.螞蟻個體間通過感知信息素交換路徑信息，最終整個蟻群在這種自組織作用下搜索出巢穴與食物源之間的最優(yōu)路徑.

根據(jù)蟻群算法中路徑信息量感知原理[24]，信息素濃度最高的路徑組合方案最優(yōu)，因此在蟻群算法中輸入發(fā)射點、反射點及接收點三維坐標，經(jīng)蟻群算法搜索后，在保證精度的前提下，可以輸出一部分優(yōu)勢 射線路徑替代所有傳播路徑.

1 .1 優(yōu)勢路徑搜索

1.1.2 射線方案禁忌表的建立

建立L個相互獨立的禁忌表：Chart1、Chart2、 ···、ChartL，分別存放射線路徑組合方案 φ1、φ2、 ···、 φL. 基于蟻群算法分布計算和啟發(fā)式搜索原理并行搜索L個禁忌表中的路徑組合方案，搜索結束后的每個禁忌表分別存放按接收射線傳輸衰減從低到高排列的射 線路徑組合方案φ′1、φ′2、 · ··、φ′L.

1.1.3 最優(yōu)射線路徑方案搜索

滿足 0＜ε≤εr且 0≤Δε≤0.01 dB的情況下，尋找優(yōu)勢射線數(shù)量最少的路徑組合方案，即為最終確定的優(yōu)勢射線路徑. 式中，εp和εp+1(0＜p＜L-1)分別為滿足約束條件下第p、p+1個禁忌表最優(yōu)路徑組合方案的均方誤差.

分析發(fā)現(xiàn)：傳統(tǒng)射線追蹤法發(fā)射角度間隔均勻的射線，通過一次次判交運算及是否到達接收點篩選符合條件的射線路徑，計算時間較長；而快速射線追蹤法根據(jù)蟻群算法確定到達接收點的射線路徑組合方案，避免大量無用或貢獻小的射線參與，進而省略判交運算時間，并結合計算機多線程并行計算技術，加速射線最優(yōu)路徑組合方案的篩選過程，從本質上提高了射線追蹤的追蹤速度，極大地降低了追蹤復 雜度.

1.2 射線傳播特性分析

射線在復雜的室內(nèi)場景中進行傳播時，將不可避免地與建筑面發(fā)生反射、繞射和透射[25-26]. 由于射線在不同電磁參數(shù)的建筑面衰減不同，需要考慮不同場景下射線的篩選情況. 當射線在室內(nèi)外之間的大場景傳播時，需要計算射線路徑上的自由空間傳播損耗，以及在建筑面產(chǎn)生的反射、繞射和透射損耗等. 然而，室內(nèi)小場景中穿透建筑的射線將會產(chǎn)生很大損耗，因此不考慮作為優(yōu)勢射線[27]. 選取優(yōu)勢射線，計算收發(fā)點間的自由空間傳播損耗、反射損耗、繞射損耗等，可以很好地預測室內(nèi)毫米波信道特性.

1.3 多線程并行計算

根據(jù)分布計算理念，在此采用四核CPU、Parpool并行池提供多核并行運算工具，將所有禁忌表中不同方案的路徑損耗計算任務分配給4個CPU同時運算，實現(xiàn)CPU的并行運算. 同時，每個CPU拿到分配的禁忌表方案后，逐個計算每種禁忌表中射線路徑組合方案的損耗值，實現(xiàn)禁忌表方案的多線程運算.采用多線程并行運算技術可以大幅度提高CPU的利用率[25]. 如快速射線追蹤法所示，蟻群算法確定射線路徑組合方案后，采用多線程并行計算加速不同組合方案的篩選情況，直到得到最優(yōu)射線組合方案. 利用多線程并行計算技術實現(xiàn)最優(yōu)路徑組合方案的加速篩選，其流程如圖1所示.

圖1 多線程并行計算流程Fig. 1 Multi-threading parallel computing process

2 仿真討論

基于上述理論，室內(nèi)毫米波信道特性預測時，首先根據(jù)室內(nèi)環(huán)境場景數(shù)據(jù)導入算法；基于此，利用蟻群算法確定優(yōu)勢射線路徑，設定誤差閾值，采用多線程并行加速；最后，在滿足誤差閾值的前提下，統(tǒng)計輸出優(yōu)勢射線的路徑損耗值及預測誤差.

2.1 場景設置

在此，利用AutoCAD軟件對室內(nèi)實測場景進行二維建模，統(tǒng)計室內(nèi)環(huán)境中的材質、高度及相對介電常數(shù)，并導入算法，放置收發(fā)天線(1個發(fā)射點，2個接收點)，進行建模仿真. 仿真室內(nèi)環(huán)境二、三維建模及收發(fā)機位置如圖2所示，室內(nèi)大小為13.8 m×12.1 m，并放有12張辦公桌子，6臺電腦，發(fā)射機位置(1，4.9)，2個接收機位置即兩個接收點分別為(6.225，4.9)、(11.289 5，4.9). 室內(nèi)建筑材質及電磁參數(shù)如表2所示，室內(nèi)收發(fā)天線參數(shù)設置如表3所示.

圖2 室內(nèi)二、三維環(huán)境模型及收發(fā)機位置Fig. 2 Indoor 2D, 3D environment model and transceiver location

表2 室內(nèi)建筑材質參數(shù)Tab. 2 Interior building material parameters

表3 室內(nèi)收發(fā)天線參數(shù)設置Tab. 3 Setting of indoor transmitting and receiving antenna parameters

2.2 優(yōu)勢射線路徑的確定

圖3 優(yōu)勢射線路徑的確定流程Fig. 3 Determination process of dominant ray path

圖4 室內(nèi)射線傳播路徑Fig. 4 Indoor ray propagation path

采用多線程并行計算技術分別計算兩個接收點禁忌表中射線路徑組合方案的路徑損耗值lPi與lPL，在滿足公式(6)、(7)的約束條件下，尋找最優(yōu)射線路徑組 合方案.

2.2.1 基于接收點1優(yōu)勢射線路徑規(guī)劃

Step1. 基于蟻群算法尋找接收點1最優(yōu)射線路徑組合方案：計算第一個接收點禁忌表Chart1、Chart2、 ···、Chart25的射線路徑組合方案的損耗值，并與WinProp傳統(tǒng)射線追蹤法預測結果進行對比，計算路徑損耗的均方誤差. 在本仿真案例中發(fā)現(xiàn)，在第一個接收點最優(yōu)射線路徑組合方案的篩選過程中，發(fā)現(xiàn)禁忌表Chart5中首個射線組合方案{L0，L11，L13，L12，L21}的均方誤差ε5=0.751 0 dB滿足不大于1 dB要求，與禁忌表Chart6中首個射線組合方案的均方誤差(ε6=0.750 1 dB)的偏差為0.000 9 dB，滿足收斂誤差不大于0.01 dB的要求.

Step2. 保存接收點1最優(yōu)射線路徑組合方案：由Step1可知，接收點1的最優(yōu)射線路徑組合方案為{L0，L11，L13，L12，L21}，優(yōu)勢射線路徑如圖5所示.

圖5 室內(nèi)第一個接收點優(yōu)勢射線路徑Fig. 5 The dominant ray path of the first receiving point in the room

設定誤差閾值εr=1 dB，對第一個接收點建立25個禁忌表：Chart1、Chart2、 ···、Chart25. 采用蟻群算法搜索第一個接收點的25個方案集合，并將搜索完成的方案集合放入對應禁忌表中：Chart1存放排列后的方案φ′1=[{L0}，{L11}，{L13}，{L14}，{L12}， ···，{L219}]；Chart2存放排列后的方案 φ′2=[{L0，L11}，{L0，L13}，{L0，L14}，{L0,L12}， ···， {L20，L19}]； ···；Chart25存放排列后的方案φ′25= [{L0，L11，L13， · ··，L219}].

同理，對第二個接收點建立23個禁忌表：Chart1、Chart2、 ···、Chart23. 分別存放蟻群算法搜索后的射線路徑組合方案：Chart1存放排列后的方案φ′1=[{L0}，{L11}，{L13}，{L14}，{L12}， ···，{L217}]；Chart2存放排列后的方案 φ′2=[{L0，L11}，{L0，L13}，{L0，L14}，{L0,L12}， ···， {L18，L17}]； ···；Chart23存放排列后的方Step3. 最優(yōu)射線路徑組合方案的確定：由最優(yōu)射線路徑組合方案獲知優(yōu)勢射線為1條直達路徑、3條一次反射路徑、1條二次反射路徑，然而從室內(nèi)電波傳播機理角度考慮，確定該5條優(yōu)勢射線路徑的篩選方法如圖6所示.

圖6 5條優(yōu)勢射線路徑的確定流程Fig. 6 Determination process of dominant 5 ray paths

導入場景數(shù)據(jù)后，先判斷室內(nèi)收發(fā)機之間的射線直達路徑是否存在，之后再對一次反射及二次反射路徑進行統(tǒng)計篩選. 如果直達路徑存在，最終的5條優(yōu)勢射線路徑應為1條直達路徑、一次反射中3條最短路徑、二次反射中1條最短路徑；若收發(fā)機之間的直達路徑不存在，則5條優(yōu)勢射線路徑為一次反射中4條最短路徑、二次反射中1條最短路徑.

Step4. 采用接收點1的最優(yōu)射線路徑組合方案驗證接收點2的場景：在未知接收點2最優(yōu)射線路徑組合方案的基礎上，利用上述接收點1確定的最優(yōu)方案{L0，L11，L13，L12，L21}直接計算第二個接收點，發(fā)現(xiàn)接收點2的均方誤差ε5=0.774 0 dB和收斂誤差Δε5=0.004 8 dB滿足要求，初步驗證了選取5條優(yōu)勢射線路徑替代全射線路徑分析在保證精度的前提下可以很好地達到預測效果. 由上述過程初步可以看出：由于僅有5條路徑參與計算，大大簡化了運算量 ，運算量降低了80%.

2.2.2 基于接收點2優(yōu)勢射線路徑規(guī)劃

Step1. 基于蟻群算法尋找接收點2最優(yōu)射線路徑組合方案：首先計算第二個接收點各禁忌表Chart1、Chart2、 ···、Chart22的射線路徑組合方案的損耗值，并與WinProp預測結果進行對比，發(fā)現(xiàn)禁忌表Chart4、Chart5和Chart6的首個路徑組合方案{L0，L11，L13，L12}、{L0，L11，L13，L12，L21}和{L0，L11，L13，L12，L21，L14}的 均 方 誤 差 分 別 為ε4=0.854 6 dB、ε5=0.774 1 dB和ε6=0.769 3 dB，收 斂 誤 差 分 別 為：Δε4=0.080 5 dB、Δε5=0.004 8 dB和Δε6=0.001 3 dB. 不難發(fā)現(xiàn)當取禁忌表Chart5首個射線路徑方案時，方可滿足誤差閾值不大于1 dB和收斂誤差不大于0.01 dB的要求.

Step2. 保存接收點2最優(yōu)射線路徑組合方案：由Step1可知，接收點2的最優(yōu)射線路徑組合方案為{L0，L11，L13，L12，L21}，優(yōu)勢射線路徑如圖7所示.

圖7 室內(nèi)第二個接收點優(yōu)勢射線路徑Fig. 7 The dominant ray path of the second receiving point in the room

Step3. 采用接收點2的最優(yōu)射線路徑組合方案驗證接收點1的場景：假定在未知接收點1最優(yōu)射線路徑組合方案的基礎上，利用接收點2的最優(yōu)射線路徑組合方案直接計算接收點1的均方誤差為0.751 0 dB，收斂誤差為0.000 9 dB，滿足誤差要求. 由上述過程同樣可以看出：由于僅有5條路徑參與計算 ，大大簡化了運算量，運算量降低了78%.

2.2.3 小結

在本文仿真案例中，收發(fā)機高度1.8 m高于室內(nèi)家具的建筑高度，因此兩個接收點的射線路徑均為收發(fā)機之間的直達射線路徑以及在墻面上反射的射線路徑. 在此基礎上，觀察2.2.1及2.2.2節(jié)兩接收點射線路徑組合方案的交叉驗證結果，發(fā)現(xiàn)禁忌表Chart5的最優(yōu)方案即為最優(yōu)射線路徑組合方案. 因此，在保證精度的前提下，選取5條優(yōu)勢射線替代全射線路徑分析可以在本質上改善射線追蹤算法，提高射線追蹤效率，降低復雜度，運算量降低了約80%.最后，統(tǒng)計兩個接收點各禁忌表最優(yōu)方案的路徑損耗誤差如圖8所示，發(fā)現(xiàn)誤差統(tǒng)計結果可以宏觀地驗 證上述結論.

圖8 室內(nèi)兩接收點各禁忌表最優(yōu)方案路徑損耗誤差Fig. 8 The path loss error of the optimal scheme of the tabu tables at 2 indoor receiving points

2.3 結果對比

在2.2節(jié)最優(yōu)射線路徑組合方案的基礎上，統(tǒng)計兩個接收點分別在28 GHz、32 GHz和38 GHz頻率下的路徑損耗如圖9所示：兩個接收點的路徑損耗值隨著頻率的增加而增加；由于第一個接收點距離發(fā)射機較近，所以射線在室內(nèi)傳播到達接收點的路徑損耗小于第二個接收點. 可以看出接收點1下本文方法較WinProp路徑損耗低，而接收點2兩種方法下路徑損耗相差較小，主要由于本文方法在挑選接收點1優(yōu)勢路徑時，除常規(guī)路徑外，還融入經(jīng)上方遠端一次反射路徑能量，體現(xiàn)出蟻群優(yōu)勢路徑搜索方法的優(yōu)勢.

圖9 兩個接收點在28 GHz、32 GHz和38 GHz下的路徑損耗對比Fig. 9 Path-loss comparison of 2 receiring points at 28 GHz, 32 GHz and 38 GHz frequency

近幾年國內(nèi)外研究學者提出了許多優(yōu)秀的射線追蹤改善方法：網(wǎng)格剖分加速法、基于KD-tree的快速剪枝法等. 基于本文接收點1仿真案例，網(wǎng)格剖分加速法[28-29]對某一條射線進行追蹤時，首先計算該條射線起點在二維平面的投影，然后根據(jù)該射線的傳播方向，不斷與其路徑沿線網(wǎng)格中的三角面(每個網(wǎng)格存放兩個三角面)進行判交運算，當遇見交點時停止判斷. 因此，網(wǎng)格剖分法得到表示向上取整)條射線路徑，雖然精度很高，但由于并未對優(yōu)勢路徑進行提取，即使采用GPU加速技術，仍需計算高達種射線路徑組合方案，計算復雜度將會很高. 基于KD-tree的快速剪枝法[24]將本文案例中的室內(nèi)場景進行網(wǎng)格剖分，并將三角面元信息存儲到KD-tree中，通過遍歷KD-tree樹，將不滿足要求的射線路徑進行快速剪枝，減少射線與三角面元判交計算的次數(shù). 因此，基于KD-tree的快速剪枝法相比于網(wǎng)格剖分法降低了射線的判交運算時間. 然而，最終得到的射線路徑仍然是本文中接收點1的25條射線路徑，需要計算的射線路徑組合方案為種，計算復雜度仍然很高. 相比于上述兩種方法，本文對兩個接收點眾多的射線路徑組合方案分類存放到禁忌表中，并利用蟻群算法搜索每個禁忌表中的射線路徑組合方案，依托多線程并行計算技術，最終得到滿足誤差閾值的最優(yōu)射線路徑組合方案，大幅降低了射線追蹤的計算時間及復雜度.

對比上述方法，統(tǒng)計快速射線追蹤法的計算效率和復雜度如表4所示(計算機仿真配置為Windows10系統(tǒng)、處理器i5-2430M、主頻2.4 G、運行內(nèi)存2 G).在追蹤時間方面：采用傳統(tǒng)射線追蹤法預測本案例室內(nèi)信道特性時耗時最長為1.929 0 s；雖然國內(nèi)外其他改善方法追蹤時間明顯降低，但仍然高于本文方法的追蹤時間0.402 8 s. 在提速效果方面，提速量為

表4 射線追蹤改善方法性能對比Tab. 4 Improved ray tracing method, calculation efficiency and complexity

式中：tΔ為 傳統(tǒng)射線追蹤法的追蹤時間；to為射線追蹤改善方法的追蹤時間. 本文方法提速效果378.9%，遠大于其他改善方法. 在計算復雜度方面：由于傳統(tǒng)射線追蹤法將室內(nèi)場景剖分為ξ2個單元格，需要發(fā)射條射線進行判交運算，因此復雜度最高為網(wǎng)格剖分法通過劃分三角面元尋找射線與面元交點，需要發(fā)射條射線進行判交運算，將復雜度降低至O(ξ|lg ξ|)；基于KD-tree的快速剪枝法將面元信息存儲到二叉樹結構中，遍歷滿足要求的射線路徑，只需發(fā)射ξ條射線路徑進行判交運算，復雜度為O(ξ)；而本文方法只需發(fā)射優(yōu)勢路徑進行判交運算，故復雜度最低為O(5). 對比不同射線追蹤改善方法，本文基于蟻群算法的優(yōu)勢射線路徑追蹤在保證精度的前提下，大幅提高了計算效率，降低了復雜度.

表5 不同測試頻率下的路徑損耗及均方根誤差Tab. 5 Path loss and prediction error at different test frequencies

3 結 論

針對毫米波在室內(nèi)傳播情況，基于蟻群分布計算和啟發(fā)式搜索原理，能夠尋找電波傳播優(yōu)勢射線路徑代替所有射線路徑. 文中依托多線程并行計算技術，預測了室內(nèi)28 GHz、32 GHz及38 GHz頻率下的路徑損耗. 經(jīng)仿真分析可以看出：1)在降低計算時間和復雜度方面，本文方法較網(wǎng)格剖分法及KD-tree快速剪枝法具有很大優(yōu)勢，在保證預測精度的前提下，可以大幅提高射線追蹤的計算效率，降低復雜度.2)文中方法基于蟻群算法選取5條優(yōu)勢射線路徑，可以在保證誤差的前提下很好地預測室內(nèi)毫米波傳播特性. 3)在確定性室內(nèi)場景下，利用蟻群算法搜索第一組收發(fā)機之間的優(yōu)勢射線路徑，作為相似位置在優(yōu)勢路徑選擇方案，能夠在保證路徑損耗預測結果準確性的基礎上，加速場景的整體預測過程.

上述成果提升了毫米波在室內(nèi)場景下的仿真效率，可為毫米波傳播分析、無線電業(yè)務規(guī)劃等方面提供參考依據(jù).