楊 敢，陳佳雯，戴炳哲，姚年鵬，蔣 堯，張其林

(南京信息工程大學氣象災害教育部重點實驗室/氣候與環(huán)境變化國際合作聯(lián)合實驗室/氣象災害預報預警與評估協(xié)同創(chuàng)新中心/中國氣象局氣溶膠與云降水重點開放實驗室，南京 210044)

0 引言

閃電作為大氣對流活動中劇烈的放電現(xiàn)象，所產(chǎn)生的強烈電磁輻射能量可能會對電子電器設備造成破壞作用。同時，我們也可以利用雷電產(chǎn)生的電磁場對雷電進行精準定位、雷電流參數(shù)反演以及防雷工程的設計等。而雷電電磁波在近地面?zhèn)鞑ミ^程中必然會受到土壤電導率、土壤色散以及高低起伏的不規(guī)則地形等因素的影響，導致電磁波波形峰值、上升沿時間和半峰值時間等都被改變，從而降低了雷擊點的定位精度。因此，研究雷電電磁波傳播特性有重要意義。

目前對雷電電磁波的影響因素，比如土壤電導率[1-3]、土壤各向異性[4-5]和土壤色散效應[6-7]等已有較多研究。另外，我國地形類型眾多，在海面上雷電電磁場的傳播特征也有一些研究[8-9]。而對于山地地形，張源源[10]研究了不規(guī)則地表粗糙程度跟電磁場衰減的關(guān)系。Li等[11]發(fā)現(xiàn)雷電電磁波在傳播路徑上遇到錐型山體時，電磁波的波形峰值大小以及波形形狀會發(fā)生改變。而且路徑上山越高，影響越大。接著Li等[12]測量了雷擊山頂高塔的雷電流實際數(shù)據(jù)，把它和雷電探測網(wǎng)觀測數(shù)據(jù)反演出來的雷電流進行對比，發(fā)現(xiàn)反演出來的雷電流幅值達到實測的1.8倍，即山體地形對電場有明顯的增強作用。Khosravi等[13]建立了雷擊山頂高塔模型，模擬發(fā)現(xiàn)雷擊山體會使得電磁場幅值增大。Hou等[14]也建立二維FDTD山體模型，研究了雷擊山頂時近距離和遠距離電磁場的變化規(guī)律，發(fā)現(xiàn)垂直電場在近距離受山體的屏蔽作用下會被削弱，在遠距離因山腳與平地的反射會增強。對于真實地形對雷電電磁波的影響，Schulz和Diendorfer等[15]用真實地形的傳播路徑對模式進行修正，進行修正后的閃電定位點會更接近真實雷擊點。Li等[12]比較了terrain-envelope method和tight-terrain-fit method兩種方法對閃電定位的修正作用。

總的來說，以上研究涉及實際山體地形的雷電電磁波傳播特性較少，且大部分都局限在孤立山體。雷電電磁場沿復雜起伏地形傳播時，會以繞射和反射的方式經(jīng)過山體，會引起電磁場到達測站時間的增加，另外電磁場也會出現(xiàn)強度的變化，從而使根據(jù)雷電電磁場波形進行的閃電定位產(chǎn)生誤差。因此如何將地形地貌帶來的誤差進行修正，提高閃電的探測效率和定位精度，是個亟需解決的問題。

因此筆者在二維柱坐標系下FDTD模型加入了真實地形，研究了真實地形對切向磁場波形以及閃電定位的影響。在考慮到現(xiàn)有的TOA時差定位算法的缺陷，即忽略了地形對電磁波的到達時間延遲作用，提出了不同步長的地形包絡法，對定位精度進行優(yōu)化。

1 研究內(nèi)容和方法

1.1 測站和雷擊點位置

烏魯木齊市位于新疆中部，地處天山北麓、準噶爾盆地南緣，地勢起伏較大，山地面積廣闊。另外烏魯木齊地屬中溫帶大陸性干旱氣候，閃電主要發(fā)生在每年的4月-9月，7月達到峰值。雷暴空間分布總體表現(xiàn)為山區(qū)多于平原、南部多北部少的特征[16]。本研究把實際安裝在烏魯木齊周邊區(qū)域的5個測站位置應用于FDTD模式中。測站共5個，分別為柴窩鋪站、阿什里站、呼圖壁站、雷達站和阜康站。最大基線長度接近130 km。

表1 測站坐標Table 1 Station coordinates

表2 雷擊點坐標Table 2 Lightning coordinates

圖1為烏魯木齊周邊地區(qū)測站和雷擊點位置圖，圖中5個測站位置用三角形表示。在雷電探測網(wǎng)區(qū)域選取了4個模擬雷擊點L1-L4，用十字表示。其中3個點設置在探測網(wǎng)西南方向的天山山脈上，一個點設在北部盆地區(qū)域作為對比。選取這些雷擊點原因在于：1)點在方圓幾公里區(qū)域內(nèi)為最高處，也就是假設閃電更趨向于擊中山頂，其中L1和L2海拔高度達到4 000 m；2)點處于測站探測范圍內(nèi)不同的位置，點之間距離在20 km以上，減少模擬結(jié)果的偶然性。

圖1 測站和雷擊點位置圖Fig.1 Location map of stations and lightning

在圖1中，L3相比L1和L2更接近探測網(wǎng)的中心，且海拔高度更低。L1、L2和L3這3個雷擊點到測站的傳播路徑上海拔高度基本都呈現(xiàn)下降的趨勢。L1和L2雖然處在探測網(wǎng)不同位置，但路徑中都會遇到高于雷擊點的山峰，L3到測站路徑上海拔高度基本一直降低。L4處于盆地區(qū)域，周圍地勢平坦，從L4到每個測站的路徑高度起伏不大，可以作為其他點的參考。

1.2 FDTD模式簡介

時域有限差分(FDTD)是把麥克斯韋方程組的微分方程轉(zhuǎn)化為差分形式，將電場和磁場在空間格點化和在時間離散化的一種電磁場計算方法。本研究采用二維柱坐標系下的FDTD模型來模擬雷電通道產(chǎn)生的電磁場沿真實地形路徑的傳播特征。另外因為采用磁場波形進行閃電定位，下面給出切向磁場的表達式[17]：

(1)

圖2為FDTD模式示意圖，模式中網(wǎng)格寬度Δr×Δz取10 m×10 m，滿足數(shù)值穩(wěn)定性條件，同時也與插值后的地形網(wǎng)格數(shù)據(jù)相匹配。時間步長 Δt=16.7 ns，符合Courant 穩(wěn)定性條件。土壤電導率σg=0.001 S/m，相對介電常數(shù)εr=10。空氣電導率σair=0 S/m，相對介電常數(shù)εair=10。邊界條件采用CPML吸收邊界。閃電回擊通道放置于計算區(qū)域的最左側(cè)山頂?shù)纳厦?，通道高度H=7 500 m。

圖2 模型示意圖Fig.2 Model diagram

高程數(shù)據(jù)采用美國太空總署(NASA)和美國國防部國家測繪局(NIMA)聯(lián)合測量的，公開的SRTM(Shuttle Radar Topography Mission)V4.1地形高程數(shù)據(jù)，數(shù)據(jù)的分辨率為0.000 83°。利用SRTM數(shù)據(jù)先計算得到雷擊點到測站的地形剖面，再將其線性插值成分辨率為10 m的數(shù)據(jù)后放入FDTD模式中。

1.3 地閃回擊模型

本研究使用Heidler雙指數(shù)函數(shù)[18]來模擬繼后回擊下雷電流通道底部的基電流，用式(2)表示。典型的繼后回擊的電流參數(shù)見表3。

表3 繼后回擊電流參數(shù)典型值Table 3 Typical subsequent return stroke current parameter

(2)

其中，

(3)

(4)

式(2)中，i01和i02分別為擊穿電流和電暈電流峰值，式(3)、(4)中η1、η2分別為擊穿電流和電暈電流峰值的修正因子。τ11和τ12分別為擊穿電流和電暈電流的波形上升時間,τ12和τ22分別為擊穿電流和電暈電流的波形下降時間。

回擊工程模型采用MTLE模型[19]，即通道中回擊電流沿通道高度以指數(shù)形式衰減。通道中t時刻z′高度處的電流分布如下：

i(z′,t)=e-z′/λi(0,t-z′/v)

(5)

式中,e-z′/λ為指數(shù)衰減系數(shù)，λ為衰減因子，取2 000 m，v為回擊速度，取1.5×108m/s。

1.4 波形到達時間計算方法

時差定位算法(Time of Arrival, TOA)是目前常用的一種閃電定位算法，在TOA定位算法中，因為測站間的距離是固定的，而閃電到達不同測站會有一個時間差，利用該時間差可以在每兩個測站中構(gòu)造出一條雙曲線，多個測站能夠構(gòu)造出的多條雙曲線。曲線相交的點即認為閃電信號輻射點的位置。如果測站有多于四站，可以得到三維空間下的雙曲面。如果多于五站即出現(xiàn)多個非線性方程，可以通過優(yōu)化算法求解輻射信號點的三維空間坐標位置以及信號發(fā)生時間。但TOA算法里是假定電磁波以光速沿直線路徑傳播到測站，在實際傳播中，起伏地形會導致電磁波尤其是高頻分量的衰減以及信號到達時間的延長，電磁場波形形狀以及波形峰值大小必定會改變，進而產(chǎn)生定位誤差。

而對于TOA定位算法本身，不同的波形到達時間計算方法同樣可能產(chǎn)生誤差。只有先根據(jù)波形確定每個站的電磁波到達時間，從而獲得每個站之間的時間差，才能進行定位。而如何定義到達時間在波形的位置，筆者簡要介紹兩種使用到的方法即10%峰值法以及逐峰法，并加入互相關(guān)法一起討論對比。

如圖3所示，10%峰值法時間tt可以用以下公式得到：

圖3 到達時間示意圖Fig3 Diagram of arrival time

(6)

其中，tpeak是波形峰值時間，即逐峰法定義的時間。t10是波形上升沿達到峰值強度的10%的時間，tt即10%峰值法定義的時間。

互相關(guān)法是對不同測站接收到的同一個信號波形做互相關(guān)計算，求互相關(guān)函數(shù)在時間窗口內(nèi)取得最大值時兩個信號波形的時間差。所以，相比逐峰法和10%峰值法確定的是某一個時間點，互相關(guān)法計算的是幾段波形間的平均時間差。本研究選擇雷達站作為中心站，式(7)為互相關(guān)函數(shù)表達式。

(7)

其中，T為時間窗口長度，f(t),g(t)為該時間窗口中兩個觀測子站接收到的閃電脈沖輻射源磁信號。

1.5 地形包絡法

雖然對閃電定位普遍采用的是TOA時差法，但顯然這種算法并沒有考慮實際地形對閃電電磁波產(chǎn)生的路徑延長問題，即信號到達每個測站的時間都會比理想情況下略長，因此本研究用不同步長的地形包絡法來近似估計地形引起的延長時間。

如圖4所示，地形包絡即用連續(xù)線段將地形包圍起來。考慮到所有測站和雷擊點的距離大小，每隔5 km、10 km、20 km和40 km步長距離取每段地形的最高點并連線形成包絡線。把這些地形包絡線長度作為電磁波在傳播中的近似傳播距離，它們大致反映出了信號傳播路徑隨山體高度的變化。然后用包絡線長度減去水平距離后除以電磁波的傳播速度作為電磁波到達時間的誤差訂正值。假設真實路徑下用包絡法得到測站至雷擊點路徑長度為l1，水平路徑為l2，可以得到誤差訂正時間:

圖4 地形包絡法示意圖Fig.4 Diagram of terrain envelope method

Δt=(l1-l2)/c

(8)

把原各個測站的信號到達時間分別減去對應的誤差訂正時間，重新進行雷擊點定位。

1.6 雷擊點到測站真實地形

圖5和圖6為雷擊點到測站的真實地形圖。

圖5 雷擊點L1和L2到測站地形圖Fig.5 Topographic map of lightning L1 and L2 to the station

圖6 雷擊點L3和L4到測站地形圖Fig.6 Topographic map of lightning L3 and L4 to the station

可以看到，雷擊點L1、L2和L3傳播路徑在海拔高度上總體都呈現(xiàn)下降的趨勢，垂直落差達到3 000 m左右，主要在路徑的前半段有比較復雜的地形起伏。而且L1和L2路徑上會遇到高于雷擊點的幾個山峰，起伏變化程度也比L3大。L4由于取點在盆地區(qū)域，到每個測站的路徑都較為平坦，只有到柴窩鋪和阿什里站的后一部分會遇到山地地形。路徑的高程數(shù)據(jù)都線性插值到10 m的分辨率，與模式中網(wǎng)格大小相匹配。

2 模擬結(jié)果與討論

2.1 平地與真實的波形對比

如圖7和圖8為不同雷擊點到每個測站平地與真實地形的波形對比。實線和虛線分別平地和真實地形的模擬結(jié)果。山地地形明顯改變了波形的上升沿時間和波形峰值。上升沿時間增加范圍在1～5 μs，這主要是因為電磁波在傳播中受土壤有限電導率和色散效應影響，高頻電磁場存在更大的衰減。同樣峰值時間也都是滯后于平地的，但峰值的大小變化卻不盡相同。

圖7 雷擊點L1和L2到測站磁場波形Fig.7 Magnetic field waveform from lightning L1 and L2 to the station

圖8 雷擊點L3和L4到測站磁場波形Fig.8 Magnetic field waveform from lightning L3 and L4 to the station

對于雷擊點L1和雷擊點L2，柴窩鋪站的峰值增大了，而其他4個站的峰值都有略微減小，這也說明了地形對電磁波影響的復雜性。峰值減小可能是因為在長距離傳播過程中，電磁波能量受有限電導率影響不斷衰減。而峰值增強的一個可能原因是雷擊山頂時，山腳的反射增強了測站接收到的電磁場；另一個可能原因是在路徑上的特殊地形結(jié)構(gòu)導致電磁波在山峰間多次反射。而對于雷擊點L3，除了阿什里站，其他4個站的峰值都變大了，變化最大的呼圖壁站峰值增加了57%。波形上升沿時間變化也比較明顯，比如雷達站延后了2 μs左右。而雷擊點L4，因為傳播路徑比較平坦，和平地相比起伏變化差距不大，只有柴窩鋪站和阿什里站能看出波形的一些區(qū)別。由此可以得到，雷電電磁波在山體傳播時，無論用哪種波形到達時間計算方法，都可能產(chǎn)生幾個微秒的時間誤差，進而產(chǎn)生上千米的定位誤差。

2.2 定位點的偏差及修正

圖9-圖12為在互相關(guān)法、10%峰值法和逐峰法3種不同波形到達時間算法下，4個雷擊點在平地、真實地形以及不同步長地形包絡法修正后的定位誤差對比圖。對于平坦地表，雷擊點在10%峰值法下的定位效果最好，定位誤差都在100多米，其他兩種到達時間算法誤差較大。這說明在平坦開闊的地形環(huán)境下，10%峰值法定位效果最好，其次是互相關(guān)法，再是逐峰法。但在真實地形下，10%峰值法定位效果卻未必最好。比如真實地形下L1和L3用逐峰法定位效果就更好,而L2用互相關(guān)法更好。

圖9 雷擊點L1定位誤差Fig.9 Location error of lightning L1

圖10 雷擊點L2定位誤差Fig.10 Location error of lightning L2

圖11 雷擊點L3定位誤差Fig.11 Location error of lightning L3

圖12 雷擊點L4定位誤差Fig.12 Location error of lightning L4

對比L1、L2和L3在真實地形和平坦地形條件下的定位誤差可以看出，真實地形較大地增加了定位誤差。L2在10%峰值法下增加了3 800 m的誤差，L1也增加了一千多米的誤差，這說明了復雜地形對閃電定位的影響不可忽略。而L4變化不大，相差一百多米，這是因為L4到每個測站路徑都較為平坦。

4種不同步長的地形包絡法對雷擊點的修正效果不同。10 km步長的地形包絡法的修正效果最好，它對每個點都成功修正，即降低了雷擊點的定位誤差。其中L1在逐峰法下修正的最多，達到510 m，而10%峰值法下修正467 m。L2在每種時間到達法下都修正了將近400 m。L3也修正了幾十到上百米。

但包絡法的修正效果并不是絕對的。在互相關(guān)法和10%峰值法算法下，對L1應用20 km和40 km步長包絡法修正，以及對L3應用5 km步長包絡法修正時反而增大了定位誤差。即在一定范圍內(nèi)，步長取的越短，包絡線更貼近地表，包絡線更長，對地形引起的時間延遲量修正的越多。但如果步長取的過短，又會過高地估計時間延遲量，反而增加了誤差。

L4作為對比，因為只有到柴窩鋪站和阿什里站的路徑上有山，山體較少，所以用地形包絡法修正后的定位誤差基本沒有變化。

3 結(jié) 論

1)相比平地，真實的起伏山地地形會增加閃電切向磁場波形的上升沿時間，也會改變波形的峰值大小。本研究中切向磁場波形峰值變化最大的增加了57%，上升沿時間增加范圍在1～5 μs，進而導致閃電定位誤差變大。相比平地，雷擊點L2在真實地形路徑下定位誤差增加3 800 m。

2)本研究提供了一種可以修正閃電定位誤差的方法，即不同步長的地形包絡法，它在一定程度上可以減少由地形帶來的閃電定位誤差。在考慮了烏魯木齊周邊實際地形下，雷擊點L1最大修正了510 m。而且不同步長下的地形包絡法修正效果不同，10 km步長的包絡法修正效果最好。