基于LiDAR點云的地理電磁環(huán)境建模及仿真研究

王亞偉 鐘選明 周亮 廖成

(西南交通大學電磁場與微波技術(shù)研究所,成都 610031)

引 言

隨著信息化社會的發(fā)展,人們對通信質(zhì)量的要求越來越高,無線電波作為信息傳遞的一種重要載體,其傳播特性越來越受到人們的重視,準確地預(yù)測電波傳播特性能夠為基站的選址、網(wǎng)絡(luò)的優(yōu)化、電臺的配置等提供有效的理論依據(jù)和工程指導[1]. 從20世紀70年代開始,國內(nèi)外相關(guān)學者對電波傳播特性預(yù)測展開了大量的研究工作,建立了很多的電波預(yù)測模型,如射線跟蹤法、時域有限差分法等,但這些模型在處理復(fù)雜地理環(huán)境中的電波傳播問題時難以在精度和效率之間達到一個很好的平衡. 由波動方程推導而來的拋物方程方法(Parabolic Equation Method, PEM)能夠快速準確地處理復(fù)雜地表邊界和復(fù)雜大氣結(jié)構(gòu)對電波傳播的影響[2-3],在求解大區(qū)域復(fù)雜環(huán)境中的電波傳播問題時得到了廣泛的應(yīng)用.

不規(guī)則的地形,多樣化的地表覆蓋物等是影響電波傳播特性的主要因素. 目前,已有許多國內(nèi)外學者基于拋物方程對地理環(huán)境中的電波傳播問題展開了研究,Holm[4]、郭建炎[5]等人將森林視為有耗的介質(zhì)層,采用二維拋物方程研究了不規(guī)則地形中森林覆蓋環(huán)境下的電波傳播問題,并與實測結(jié)果對比,驗證了模型的正確性;張青洪[6]等人將森林看成空氣和水的混合物,建立了森林等效介電常數(shù)模型,對拋物方程的森林模型進行了改進,這些研究一般考慮的地表環(huán)境比較單一,電磁建模也比較粗糙. 白瑞杰[7]等人基于數(shù)字高程模型對傳播區(qū)域的電波特性做了預(yù)測分析,由于數(shù)字高程模型數(shù)據(jù)描述的是地形起伏的狀況,缺少森林、植被、建筑等地表覆蓋物的信息,難以滿足城市小區(qū)域電波傳播預(yù)測的電磁建模要求. 因此,研究更加接近真實場景下的電磁建模很有必要.

LiDAR點云是由機載激光雷達系統(tǒng)通過發(fā)射和接收激光脈沖獲取的地表密集的高精度的三維點坐標,是對地形、地物特征的一種高精度的描述[8]. 通過對其分析處理,不僅可以獲得地形起伏的信息,還可以獲得植被、建筑等地物信息,將其應(yīng)用于復(fù)雜地理環(huán)境的電磁建模,將會使電波傳播的預(yù)測更加準確,因此本文研究了基于LiDAR點云的地理電磁環(huán)境建模. 采用數(shù)學形態(tài)學方法對格網(wǎng)化后的LiDAR點云做了濾波等處理,獲取地形、地物等信息,結(jié)合相應(yīng)物質(zhì)的電磁參數(shù),進行了地理環(huán)境的電磁建模,在此基礎(chǔ)上,利用三維拋物方程模型對地理環(huán)境中的電波傳播特性進行了仿真分析.

1 三維拋物方程模型

在笛卡爾直角坐標系中,假定電磁場的時諧因子為e-iωt,電場或磁場分量ψ滿足如下形式的標量Helmholtz方程:

(1)

式中:n為媒質(zhì)的折射指數(shù);k0為自由空間中的波數(shù).定義波函數(shù)u(x,y,z)=ψe-ik0x(x,y,z),將其代入式(1),并進行因式分解,可以得到

(2)

式(2)中,Q為偽微分算子,

(3)

僅考慮前向傳播,對Q作Feit-Fleck型近似,可得到如下形式的三維拋物方程:

(4)

式(4)的數(shù)值解法主要為有限差分 (Finite Difference, FD[9]) 法和分步傅里葉變換 (Split-step Fourier Transform, SSFT[10]) 法,它們都是步進迭代算法. 其中FD算法的計算步長受波長嚴格限制,主要用于計算目標的雷達散射截面(Radar Cross Section, RCS)以及城鎮(zhèn)小尺度電波傳播特性的預(yù)測. SSFT算法引入了快速傅里葉變換(Fast Fourier Transform, FFT)技術(shù),在傳播方向上可以取較大的步長,能夠快速地預(yù)測復(fù)雜地理環(huán)境中的電波傳播特性,其解的形式如下[10]:

u(x+Δx,y,z)=eik0Δx(n-2)F-1{eikxΔxF[u(x,y,z)]}.

(5)

式中:kx表示k0在x方向上的分量;F和F-1分別表示二維傅里葉正逆變換.知道初始場,利用式(5)結(jié)合邊界條件可迭代遞推求得整個空間中的場強分布.

電磁波在空間遠距離傳輸時,森林的影響是一個不容忽視的因素. 目前在研究森林對電波傳播的影響時,通常把森林等效為均勻有耗的介質(zhì)層[4-5],其等效的折射系數(shù)可表示成

(6)

式中:εr,σ分別表示森林的等效相對介電常數(shù)和電導率. 在迭代式(5)的過程中,通過式(6)不斷地修正森林中的折射系數(shù),則可處理森林中的電波傳播問題. 通常森林中的等效折射系數(shù)趨近1,滿足拋物方程的近似條件.

2 地理電磁環(huán)境建模

機載LiDAR是一種新型的主動式遙感技術(shù),它通過測量激光脈沖的傳播時間,結(jié)合POS系統(tǒng)提供的定位定姿數(shù)據(jù),能夠直接、快速獲取地面的三維坐標[11]. LiDAR獲取的點云文件記錄了一系列被掃描物體表面的三維離散點坐標,其不規(guī)則地分布在三維空間中. 不同的地形、地物對應(yīng)不同的激光腳點分布,反應(yīng)到點云數(shù)據(jù)就會呈現(xiàn)出不同的幾何拓撲形態(tài). 通常情況下,地面點云的高程變化比較平緩,點的分布比較均勻;樹木點云分布比較散亂,一般在邊緣處會有高程突變;建筑物點云整體上幾何形態(tài)較為規(guī)整,高程變化也較為一致. 水體對激光信號(波長為1 064 nm)的吸收,則會造成大面積的水體空洞區(qū). 依據(jù)這些地物呈現(xiàn)出的不同特征,采用一定的技術(shù)手段可以對LiDAR點云數(shù)據(jù)進行分類.

2.1 點云濾波

點云分類中重要的一步是從點云數(shù)據(jù)分離出地面點和地物點,也就是點云濾波. 目前國內(nèi)外學者提出許多種的濾波方法,其中具有代表性的算法有移動曲面擬合法、迭代三角網(wǎng)加密法、數(shù)學形態(tài)學方法等[10]. 為了適應(yīng)拋物方程計算傳播迭代步進的需要,本文采用了一種基于網(wǎng)格操作的數(shù)學形態(tài)學濾波算法[13].

數(shù)學形態(tài)學是一門建立在集合論的學科,目前已被廣泛應(yīng)用于各類圖像分析和處理中,其基本操作有四種,分別為腐蝕、膨脹、開和閉運算. 對于LiDAR觀測點p(x,y,z),膨脹(Dilation)運算和腐蝕(Erosion)運算分別定義為:

(7)

(8)

點(xp,yp,zp)表示在窗口(也稱結(jié)構(gòu)元)大小為w時點p的鄰域點,理論上,窗口形狀可以是任意的,實際處理中,我們通常選擇簡單的形狀,如矩形、圓形、菱形等. 由式(7)、(8)可以推斷出,膨脹運算和腐蝕運算分別能夠獲得鄰域窗口內(nèi)高程最大值和最小值. 基于膨脹和腐蝕運算,可以得到用于處理LiDAR數(shù)據(jù)的開運算和閉運算,開運算是先對給定的數(shù)據(jù)執(zhí)行腐蝕運算,然后執(zhí)行膨脹運算,閉運算則剛好相反.

數(shù)學形態(tài)學濾波的主要思想是將地物(樹木、建筑物等)通過腐蝕運算腐蝕至地表,然后通過膨脹運算進行恢復(fù)[14]. 圖1是采用線形結(jié)構(gòu)元對點云數(shù)據(jù)的某一縱剖面執(zhí)行開運算濾波的示意圖. 當有多個點落入窗口中時,執(zhí)行開運算后,將保留高程值最小的數(shù)據(jù)點. 由于樹木的尺寸小于結(jié)構(gòu)元尺寸,經(jīng)過腐蝕后首先被去除,對于尺寸大于結(jié)構(gòu)元的建筑物,則會在膨脹運算中得到重建.如果某點濾波前后的高程差在閾值范圍內(nèi),則判定該點為地面點,否則判定為地物點,并標記. 可以看出,開運算能夠很好分離出地物點以及地面點,并且能夠獲得較為平滑的地形.

圖1 LiDAR點云開運算濾波示意圖

結(jié)構(gòu)元尺寸的選擇決定了最終的濾波效果,當采用較小的濾波窗口時,尺寸大于窗口的地物(如建筑物)得不到很好的分離,采用較大的濾波窗口在分離地物的同時,也會過度地平整更多的地形特征,造成錯誤的分類. 由此可見,采用固定大小的濾波窗口很難對點云進行很好的分類,因此本文采用了線性漸進增大的濾波窗口.

2.2 點云分類

圖2為某一地區(qū)的地理影像,圖3為對應(yīng)的原始LiDAR點云高程值圖,區(qū)域中存在著大量的樹木,右半部分有一個明顯的河流空洞區(qū).區(qū)域大小約為694 m×2100 m,點云高程差約為172 m,點云間距約1.2 m,共有1 717 762個離散點.

圖2 地理影像圖

圖3 原始LiDAR點云高程值圖

針對該區(qū)域的特點,采用以下流程將區(qū)域中的物質(zhì)分為樹木、地面、水三種類型.

1) 剔除噪點. 采用k鄰近查詢球[11]方法,剔除極低點、孤立點等高程異常值.

2) 點云格網(wǎng)化. 使用一個二維格網(wǎng)覆蓋在LiDAR數(shù)據(jù)點上,格網(wǎng)間距通常略小于LiDAR點云密度,如果一個格網(wǎng)單元里有多個數(shù)據(jù)點,取格網(wǎng)單元中高程的最小數(shù)據(jù)點代表該格網(wǎng)單元的值,若網(wǎng)格單元中沒有數(shù)據(jù)點的則不做處理[14].

3) 用二值圖像表示格網(wǎng)后的點云數(shù)據(jù),有值的網(wǎng)格單元高程值用“1”表示,否則用“0”表示. 用合適的“結(jié)構(gòu)元”對二值圖像執(zhí)行閉運算. 經(jīng)過閉運算后,小面積的空白區(qū)域?qū)惶畛?剩下的空白就是大面積的水體造成空洞區(qū),并標記.

4) 對非水體區(qū)域執(zhí)行數(shù)學形態(tài)學濾波,判斷濾波前后高程差閾值是否滿足條件,分離出樹木點和地面點.

5) 由地面點通過線性插值獲得地形.

圖4和圖5分別為最終獲得的地形以及分類結(jié)果.從圖4中可以看出,在整個區(qū)域中,樹木被很好地分離,而地形特征也保持相對圓滿. 將圖5對照原始地理影像圖,可以看出,各種地物得到了很好的分類. 將二者結(jié)合起來,對各種地物賦上相匹配的電磁參數(shù),就可以對上述地理環(huán)境進行電磁建模.

圖4 地形圖

圖5 點云分類結(jié)果

3 電磁環(huán)境仿真

利用三維拋物方程模型預(yù)測上節(jié)所述地理環(huán)境下的電波傳播特性.仿真參數(shù)設(shè)置如下:天線頻率為900 MHz;高斯方向圖3 dB帶寬為10°;水平極化波;天線架設(shè)在橫向距離為0 m,縱向距離為347 m,距離地表高度為45 m的地理環(huán)境處;大氣條件為理想大氣條件. 受限于數(shù)據(jù)源,最大傳播距離為2 100 m,dy=dz=0.5 m,傳播方向上的步長取dx=1 m.

圖6顯示了源所在的地形剖面上傳播因子的分布情況,分別選取了(a)單一地表為黏性干土(介電常數(shù)為5,電導率為0.021 S/m)、(b) 真實的地理環(huán)境兩種情況,其中地理環(huán)境中水的相對介電常數(shù)和電導率為81和0.22 S/m,樹木的等效介電常數(shù)和電導率分別為1.004和180 μS/m[5]. 圖6(a)中由于地形的存在,在地形后方出現(xiàn)了電磁波的盲區(qū)以及繞射現(xiàn)象,同時也可以看到由于地形的反射出現(xiàn)了明顯的干涉條紋. 圖6(b)中由于樹木引起的遮擋,散射和吸收會產(chǎn)生較大的路徑損耗,導致空間電磁波的整體損耗要大于沒有樹木時的情形,特別是在森林的內(nèi)部,電磁波的損耗明顯增大.

(a) 地表為干土

(b) 真實地理環(huán)境下圖6 傳播因子分布偽彩圖

圖8 兩種算法下傳播因子隨高度的變化

圖7為采用FD算法[15]計算得到的上述地理環(huán)境中傳播因子的分布情況.與圖6(b)對比,可以看出兩者總體趨勢大致相同,采用三維拋物方程能夠考慮到電波的橫向繞射,更加真實地反應(yīng)出地理環(huán)境中電波傳播特性. 圖8顯示了主輻射方向上,傳播距離為2 100 m處,采用FD算法和SSFT算法得到的傳播因子隨高度的變化,兩者吻合較好,說明了本文算法的正確性及有效性,SSFT算法相對于FD算法擁有更高的效率,能夠快速地預(yù)測環(huán)境中的電波傳播特性.

圖9中顯示的是主輻射方向上,傳播距離一定時,有樹木和無樹木、有水無水的場景下傳播因子隨高度的變化情況.圖9(a)是傳播距離為2 000 m處,對比有樹木和無樹木情況下的傳播因子.可以看出:在高度小于50 m的樹木內(nèi)部,電波的衰減較大;在遠離樹木的區(qū)域,電波的衰減相對較小. 圖9(b)是傳播距離為1 920 m處,有無水時傳播因子隨高度的變化曲線.從圖中可以看出,兩者相差很小,其主要原因是,水體的范圍只有幾十米,遠離輻射源,電磁波成掠入射接觸水體,并且水體前面存在著樹木的遮擋作用,總體上水體對電磁波的傳播影響較小.

(a) 傳播距離為2 000 m

(b) 傳播距離為1 920 m圖9 不同傳播距離處傳播因子隨高度的變化

圖10是在上述地理環(huán)境中,發(fā)射天線分別架設(shè)在距離地面高度35 m、45 m和55 m的情況下,距離發(fā)射天線2 100 m處的傳播因子隨高度的變化情況. 通過比較可以看出,在一定接收高度下,電波的衰減隨著發(fā)射天線架設(shè)高度的增加而減小,因此適當把天線架設(shè)在高處,可以減小到達接收點的信號衰減,擴大無線信號的覆蓋范圍.

圖10 不同天線高度下傳播因子隨高度的變化

4 結(jié) 論

本文研究了基于LiDAR點云的PEM電波傳播模型,針對傳統(tǒng)地理電磁環(huán)境建模的不足,研究了LiDAR點云數(shù)據(jù)中地物分類,對傳播區(qū)域的地理電磁環(huán)境進行了更為準確的建模,并采用拋物方程模型對有地表覆蓋物的環(huán)境下電波傳播特性進行了仿真分析. 相對于傳統(tǒng)的電磁建模方法,本方法能夠考慮到地表覆蓋物對電波的影響,從而更真實地反映出電波的傳播特性. 下一步的工作重點是在現(xiàn)有的基礎(chǔ)上,研究包含建筑物、樹木的復(fù)雜城區(qū)的LiDAR點云分類,并結(jié)合大氣環(huán)境對區(qū)域中的電波傳播特性進行預(yù)測.

[1] DUMONT N, WATSON R J, PENNOCK S R. Use of the parabolic equation propagation model to predict TV white space availability[C]//Antennas and Propagation Conference. Loughborough, November 8-9, 2010: 353-356.

[2] LEONTOVICH M, FOCK V. Solution of the problem of propagation of electromagnetic waves along the earth’s surface by the method of parabolic equation[J]. Journal of physics USSR, 1946(7): 557-573.

[3] DONOHUE D J, KUTTLER J R. Propagation modeling over terrain using the parabolic wave equation[J]. IEEE transactions on antennas & propagation, 2000, 48(2): 260-277.

[4] HOLM P, WAERN A. Wave propagation over a forest edge parabolic equation modelling vs. GTD modelling[C]//IEEE International Symposium on Electromagnetic Compatibility. Istanbul, May 11-16, 2003: 764-767.

[5] 郭建炎, 王劍瑩, 龍云亮, 等. 基于拋物方程法的部分森林覆蓋山區(qū)電波傳播分析[J]. 電波科學學報, 2008, 23(6): 1045-1050.

GUO J Y, WANG J Y, LONG Y L, et al. Analysis of radio propagation in partly forested terrain environment using parabolic equation approach[J]. Chinese journal of radio science, 2008, 23(6): 1045-1050. (in Chinese)

[6] 張青洪, 廖成, 盛楠, 等. 森林環(huán)境電波傳播拋物方程模型的改進研究[J]. 物理學報, 2013, 62(20): 204101.

ZHANG Q H, LIAO C, SHENG N, et al. Improved study on parabolic equation model for radio wave propagation in forest[J]. Acta physica sinica, 2013, 62(20):204101. (in Chinese)

[7] 白瑞杰, 廖成, 張青洪,等. 復(fù)雜地理環(huán)境的圖像分割及電波傳播特性[J]. 強激光與粒子束, 2015, 27(10): 63-67.

BAI R J, LIAO C, ZHANG Q H, et al. Image segmentation of complex geographical environment and wave propagation characteristics[J]. High power laser and particle beams, 2015, 27(10): 63-67. (in Chinese)

[8] 陳松堯, 程新文. 機載LIDAR系統(tǒng)原理及應(yīng)用綜述[J]. 測繪工程, 2007, 16(1): 27-31.

CHEN S Y, CHENG X W. The principle and application of airborne LIDAR[J]. Engineering of surveying and mapping, 2007, 16(1): 27-31. (in Chinese)

[9] LEVY M. Parabolic equation methods for electromagnetic wave propagation[M]. London: IEE Press, 2000: 35-40.

[10] HARDIN R H, TAPPERT F D. Application of the split-step Fourier method to the numerical solution of nonlinear and variable coefficient wave equations[J]. SIAM review, 1973, 15(423): 423.

[11] 成曉倩, 馬洪超, 趙紅強,等. 一種基于區(qū)域增長的機載LIDAR濾波[J]. 測繪科學, 2010, 35(1): 61-63.

CHEN X Q, MA H C, ZHAO H Q, et al. A filtering of LIDAR based on region growing[J]. Science of surveying and mapping, 2010, 35(1): 61-63. (in Chinese)

[12] ZHANG W, QI J, WAN P, et al. An easy-to-use airborne LiDAR data filtering method based on cloth simulation[J]. Remote sensing, 2016, 8(6): 501.

[13] PINGEL T J, CLARKE K C, MCBRIDE W A. An improved simple morphological filter for the terrain classification of airborne LIDAR data[J]. ISPRS journal of photogrammetry & remote sensing, 2013, 77(1): 21-30.

[14] ZHANG K, CHEN S C, WHITMAN D, et al. A progressive morphological filter for removing nonground measurements from airborne LIDAR data[J]. IEEE transactions on geoscience & remote sensing, 2003, 41(4): 872-882.

[15] DING W, WANG K, LONG Y. Study of the electromagnetic waves propagation over the improved fractal sea surface based on parabolic equation method[J]. International journal of antennas and propagation, 2016(3): 1-7.