一種水平不均勻反常大氣折射環(huán)境反演方法

劉成國(guó) 周鴻錕 胡文韜 段開源 程潤(rùn)生 黃立峰郭相明 王紅光 蔡紅濤

（1. 湖北省射頻微波應(yīng)用工程技術(shù)研究中心, 武漢 430070；2. 武漢理工大學(xué)理學(xué)院, 武漢 430070；3. 武漢船舶通信研究所, 武漢 430205；4. 中國(guó)電波傳播研究所, 青島 266107；5. 武漢大學(xué)電子信息學(xué)院, 武漢 430072）

引 言

微波超視距反常傳播應(yīng)用必須獲得相應(yīng)的大氣折射率環(huán)境，為此產(chǎn)生了多種獲取方法[1]，包括蒸發(fā)波導(dǎo)剖面的氣海作用理論計(jì)算方法[2]等. 滿足空間連續(xù)覆蓋、電波傳播精度要求的環(huán)境參數(shù)必須能很好地描述其中存在的多種波導(dǎo)層結(jié)和分層水平不均勻等足夠精細(xì)的特征，然而上述方法并不能滿足這些要求. 近年來又有了很多新的方法[3]，其中反演方法受到越來越多的關(guān)注.

為了得到有效的新反演方法，需要分析應(yīng)用場(chǎng)景，充分融合反演理論、電波傳播計(jì)算、多源電磁環(huán)境與大氣環(huán)境數(shù)據(jù)來選擇適合的算法. 21 世紀(jì)后期開始，各種場(chǎng)合的反演技術(shù)不斷出現(xiàn). 如利用海雜波的反演技術(shù)，初步結(jié)果在50 km 內(nèi)，與直升機(jī)直接測(cè)量的結(jié)果相近[4]；利用全球衛(wèi)星導(dǎo)航系統(tǒng)反演大氣和蒸發(fā)波導(dǎo)[5-7]. 除此之外，還有用多個(gè)發(fā)射-接收電路上的點(diǎn)對(duì)點(diǎn)微波傳播測(cè)量進(jìn)行反演的方法[8]. 這些方法需要解決和場(chǎng)景密切相關(guān)的兩個(gè)核心問題.

第一個(gè)是電波傳播確定模式構(gòu)成的正問題. 正問題建立適用于復(fù)雜環(huán)境的電波傳播模型，預(yù)測(cè)給定系統(tǒng)在給定環(huán)境中的電波覆蓋情況. 目前拋物型方程[9]是廣泛使用的模型，形成了分步傅里葉算法[10]、有限差分法(finite differential method, FDM)[11]和有限元法(finite element method, FEM)[12]等具有代表性、使用廣泛的方法.

第二個(gè)是適用算法問題. 微波超視距反常傳播信息反演大氣折射環(huán)境時(shí)，涉及學(xué)科多、可選信息源多. 需要根據(jù)不同情景和大氣折射率環(huán)境的不同層次要求，選擇反演算法及其解決方案，特別是一直探索的水平不均勻復(fù)雜大氣折射率環(huán)境反演. 于是各種新發(fā)展智能算法的可用性不斷得到研究，如模擬退火算法和遺傳算法[4]、粒子群優(yōu)化算法[13]、蟻群優(yōu)化算法[14]、深度學(xué)習(xí)算法[15]及其各種混合算法[16]等等，以形成有效的大氣折射率環(huán)境反演方法.

本文研究微波超視距遠(yuǎn)距離傳播場(chǎng)景多信息的水平不均勻分層復(fù)雜大氣折射率環(huán)境反演問題，形成了一種基于雙向觀測(cè)數(shù)據(jù)的新反演方法，可為電波傳播研究等應(yīng)用提供滿足要求的海陸交界地區(qū)水平不均勻復(fù)雜大氣折射率環(huán)境信息. 該方法技術(shù)上只需要收發(fā)電路實(shí)現(xiàn)雙工功能，而天線架設(shè)不用變動(dòng)，大大降低了觀測(cè)電路工程難度，同時(shí)反演的大氣折射率環(huán)境信息能夠真實(shí)反應(yīng)海陸交界地區(qū)大氣折射率水平不均勻的情況.

1 原理要求

1.1 反演原理

反演技術(shù)用于難以測(cè)得的分布參數(shù)獲取問題.這些問題往往是不適定問題，可用一定的方法在所需空間由這些參數(shù)預(yù)測(cè)出可觀測(cè)量的值. 但是這些觀測(cè)量客觀上只能測(cè)得部分區(qū)域的值，造成反演的多解、收斂性等等不適定問題，所以需要研究盡可能用有限信息有效解決這些問題的方法[17]. 反演微波超視距反常傳播的大氣折射環(huán)境，空間分布的大氣折射指數(shù)n(r)是需要反演的參量. 可觀測(cè)量是空間中的電波參量，通常是功率P(f,r). 它是電波頻率f和空間位置r的4 維函數(shù)，實(shí)際中也只能獲得部分給定電波頻率和某些空間位置的測(cè)量數(shù)據(jù).

上述過程中折射指數(shù)的時(shí)空分布多解性、反演收斂性等不適定問題，需要通過合適精度正問題解決方法確定折射指數(shù)剖面的參數(shù)及其變化范圍，再選擇合適的優(yōu)化算法、目標(biāo)函數(shù)優(yōu)化準(zhǔn)則和搜索策略解決. 在選定的正問題求解方法中，給出所需的修正折射率M(M=(n-1)×106+103×h/6370，h是海拔高度，單位m)剖面模型搜索參數(shù)，利用歷史數(shù)據(jù)統(tǒng)計(jì)等方法獲得求解空間中折射率分布的先驗(yàn)規(guī)律(折射率剖面變化范圍的時(shí)空規(guī)律)，確定搜索參數(shù)的空間范圍，界定選定算法和目標(biāo)函數(shù)優(yōu)化準(zhǔn)則下搜索參數(shù)的搜索范圍，克服求解的不適定性，形成大氣折射率環(huán)境信息反演結(jié)果.

1.2 正問題

1.3 反演算法

目前的電磁環(huán)境和建立的電波電路中，都會(huì)設(shè)法獲得規(guī)定區(qū)域中某些位置確定頻率電波的功率，包括獲得大量數(shù)據(jù)的雷達(dá)海雜波問題、各種固定業(yè)務(wù)和移動(dòng)業(yè)務(wù)數(shù)據(jù)鏈路中確定位置的電波功率等等.具體情況中測(cè)量信號(hào)的維度和時(shí)空數(shù)據(jù)量決定了相應(yīng)的反演難度.

不同問題中環(huán)境的難度和復(fù)雜度制約了反演算法的選擇. 多個(gè)頻率、復(fù)雜下墊面、觀測(cè)點(diǎn)少的情況是高難度的代表. 復(fù)雜度大的大氣折射率環(huán)境有含波導(dǎo)層結(jié)的水平不均勻大氣環(huán)境(它支撐著微波超視距遠(yuǎn)距離傳播)等. 已有反演方法的研究表明，受限于實(shí)際使用場(chǎng)景，需要根據(jù)各種算法的特點(diǎn)，深入研究大氣的水平不均勻性反演中的參數(shù)敏感性、反演結(jié)果的準(zhǔn)確性與可靠性等問題[19]，有選擇地使用.

除了確定算法外，還需要確定解的搜索參數(shù)和搜索空間，以保證反演過程的收斂性和搜索速度等.根據(jù)反演空間的大氣折射率環(huán)境和微波超視距遠(yuǎn)距離反常傳播的先驗(yàn)知識(shí)，建立合適的先驗(yàn)?zāi)Ｐ停拍艽_定有效的搜索參數(shù)和搜索空間[13].

2 場(chǎng)景選擇和反演方法構(gòu)建

2.1 研究場(chǎng)景

本文研究波導(dǎo)層結(jié)類型多、剖面復(fù)雜、水平不均勻和出現(xiàn)概率高的跨海面固定電路上的水平不均勻大氣折射率環(huán)境的反演問題. 選擇文獻(xiàn)[20]中的微波超視距遠(yuǎn)距離傳播試驗(yàn)的頻段和線路構(gòu)成的場(chǎng)景，并使用這個(gè)場(chǎng)景在2019 年9 月測(cè)試的C 波段和X 波段電路測(cè)量數(shù)據(jù)驗(yàn)證反演效果.

試驗(yàn)中這個(gè)場(chǎng)景建立的4 條觀測(cè)線路如圖1 所示. 各線路的長(zhǎng)度為：線路1 約50 km，線路2 約70 km，線路3 約85 km，線路4 約150 km，天線的架設(shè)高度滿足各線路都出現(xiàn)波導(dǎo)傳播的超視距傳播線路. 本文針對(duì)線路1 開展反演研究，圖1 中標(biāo)出的①～④ 4 個(gè)點(diǎn)是不均勻大氣折射率反演的4 個(gè)位置. 試驗(yàn)期間4 條線路最南端的共同端設(shè)置C 波段和X 波段的收發(fā)設(shè)備開展觀測(cè)試驗(yàn)；線路1 的北端也設(shè)置C 波段和X 波段的收發(fā)設(shè)備，其中X 波段同時(shí)收發(fā)，C 波段單向測(cè)試，且其他線路測(cè)試時(shí)，該端不測(cè)試C 波段. 線路1 在試驗(yàn)期間每天開展測(cè)試，選擇超視距傳播有利的日期，其他線路只開展C 波段測(cè)試.

圖1 試驗(yàn)觀測(cè)場(chǎng)景和4 條傳播線路Fig. 1 Scenario of test and its 4 transmission lines

2.2 反演算法選擇

2.3 正問題的下墊面和天線初始場(chǎng)

正問題使用文獻(xiàn)[18]的方法，用拋物型方程的寬角有限差分進(jìn)行前向求解. 本文研究的海岸-海面-海岸下墊面是海陸混合下墊面，使大氣波導(dǎo)傳播模擬下墊面問題變得復(fù)雜. 本文研究中根據(jù)算法的數(shù)值計(jì)算網(wǎng)格和參數(shù)要求，提取了地理信息的海岸-海面-海岸剖面高程數(shù)據(jù)，按照試驗(yàn)場(chǎng)景對(duì)應(yīng)的海岸地面和海面情況，解決這一復(fù)合下墊面問題.

利用公開的地理信息數(shù)據(jù)得到下墊面需要的地面和海底高程數(shù)據(jù)，再根據(jù)Leontovich 阻抗邊界條件構(gòu)建跨海電路下墊面的海岸、水面邊界阻抗數(shù)據(jù).海岸陸地通過線性地形平移法處理描述成一定斜率地形上的邊界條件；海面部分，取用的海水相對(duì)介電常數(shù)為溫度T、鹽度S和頻率f的函數(shù)，表達(dá)式為

式中：ε∞取4.9；ε1，τ，σ 由溫度和鹽度確定；ω=2πf為電磁波的角頻率；ε0為真空中的介電常數(shù).

初始場(chǎng)設(shè)置通過對(duì)天線方向圖傅里葉逆變換得到，采用P. Valtr[21]給出的方法計(jì)算初始場(chǎng).

計(jì)算域的上邊界采用窗函數(shù)進(jìn)行計(jì)算. 所采用的窗函數(shù)形式為

式中：Z為整個(gè)計(jì)算域垂直方向的總離散格點(diǎn)數(shù)；hab為窗函數(shù)所起作用的吸收層格點(diǎn)數(shù)；ζ 為經(jīng)驗(yàn)系數(shù). ζ 與hab的具體值需要根據(jù)每次所計(jì)算頻段和計(jì)算域大小進(jìn)行調(diào)整.

2.4 大氣折射率環(huán)境先驗(yàn)?zāi)Ｐ?/h3>

研究中使用復(fù)雜大氣修正折射率M剖面及其沿水平方向的變化構(gòu)建先驗(yàn)?zāi)Ｐ? 根據(jù)波導(dǎo)環(huán)境大氣修正折射率剖面特征的研究結(jié)果，考慮場(chǎng)景構(gòu)建出如圖2 所示的復(fù)合修正大氣折射率剖面模型. 它由從地海面開始高度為hedh的蒸發(fā)波導(dǎo)層結(jié)、斜率為a1厚度為h1的大氣層結(jié)、斜率為a2厚度為h2的波導(dǎo)層結(jié)和斜率為a3的上層大氣4 層構(gòu)成. 這幾個(gè)參數(shù)和地海面的修正折射率M0一起構(gòu)建粒子群.

圖2 復(fù)合修正大氣折射率M 剖面模型Fig. 2 Hybrid M profile of the complex atmospheric environment

2.5 反演算法的實(shí)現(xiàn)

根據(jù)反演的原理和粒子群算法要求，本文構(gòu)建了如圖3 所示的基本代碼流程.

圖3 反演的粒子群算法基本流程Fig. 3 Basic flowchart of the particle swarm inversion algorithm

首先在搜索空間產(chǎn)生粒子群中的每個(gè)粒子，即先驗(yàn)?zāi)Ｐ椭械? 個(gè)參數(shù)；然后用它們產(chǎn)生的折射指數(shù)按照均勻性假設(shè)用于拋物型方程中的正問題計(jì)算，得到傳播損耗的預(yù)測(cè)值；再提取和實(shí)測(cè)值位置對(duì)應(yīng)的數(shù)據(jù)，得到適應(yīng)性函數(shù)的值，判斷是否滿足優(yōu)化準(zhǔn)則要求. 若不滿足，按照搜索規(guī)則產(chǎn)生新的參數(shù)作為粒子重新進(jìn)行計(jì)算和判斷，直至獲得滿足要求的粒子. 最后使用這些粒子參數(shù)用先驗(yàn)?zāi)Ｐ偷钠拭婧瘮?shù)構(gòu)建出折射指數(shù)的空間分布.

3 實(shí)現(xiàn)情況分析討論

針對(duì)選定場(chǎng)景及其在2019 年9 月進(jìn)行的C 波段和X 波段電路的試驗(yàn)測(cè)試，本文解決反演技術(shù)中的正問題計(jì)算和下墊面邊界問題，給出先驗(yàn)搜索空間的范圍和水平均勻與不均勻要求的問題，實(shí)現(xiàn)反演.

3.1 反演方法具體實(shí)現(xiàn)情況

1) 下墊面和天線初始場(chǎng)數(shù)據(jù)生成方法. 對(duì)于選擇的試驗(yàn)場(chǎng)景和電路，下墊面為海岸陸地和海面構(gòu)成的混合下墊面. 使用地理數(shù)據(jù)和海洋數(shù)據(jù)按2.3 節(jié)方法構(gòu)建下墊面的阻抗數(shù)據(jù). 圖4 分別是2019 年9 月14 日上午試驗(yàn)電路海域海面的溫度和鹽度分布. 將其用于式(7)計(jì)算海面海水復(fù)介電常數(shù)，進(jìn)而生成海面阻抗數(shù)據(jù).

圖4 試驗(yàn)附近海域溫度和鹽度分布Fig. 4 Sea surface temperature and salinity in the sea area near the test

式中，β 為波束寬度. 本文中，試驗(yàn)天線波束寬度為3°，天線架設(shè)高度為20 m.

2) 反演參數(shù)和適應(yīng)性函數(shù). 根據(jù)試驗(yàn)場(chǎng)景條件，形成基于物理原理分析反演參數(shù)物理特征確定粒子參數(shù)搜索空間的方法，用接收信號(hào)的來源情況構(gòu)建合適的適應(yīng)性函數(shù).

通過物理原理分析反演涉及參數(shù)的特征，確定出表1 所列的反演參數(shù)粒子和搜索區(qū)間. 比如c7對(duì)應(yīng)的參數(shù)a3，它是指數(shù)大氣的線性近似，通過物理分析，得到它的搜索范圍為

表1 反演參數(shù)粒子及其搜索區(qū)間Tab. 1 Particles of inversed parameter and its searching range

適用性函數(shù)根據(jù)試驗(yàn)電路1 的雙向收發(fā)測(cè)量做改進(jìn)后，雙向觀測(cè)數(shù)據(jù)的約束條件更好，特別是試驗(yàn)中觀察到的雙向傳播損耗不一致結(jié)果，體現(xiàn)出來能夠反演波導(dǎo)環(huán)境水平不均勻性的影響. 試驗(yàn)中的X 波段觀測(cè)中用連續(xù)波模式觀測(cè)，頻率功率測(cè)量數(shù)據(jù)可以直接使用，對(duì)其專門構(gòu)造了適應(yīng)性函數(shù)：

式中：下標(biāo)f 表示該鏈路傳播的正向；下標(biāo)r 表示該鏈路傳播的反向. 實(shí)際測(cè)量時(shí)對(duì)于正向傳播損耗發(fā)射天線約距海面20 m 高的位置，接收天線位于距海面2.5 m 高的位置，對(duì)于反向傳播損耗則反之.

3) 水平均勻大氣環(huán)境反演結(jié)果. 使用表1 的搜索條件和式(6)的適用性函數(shù)，實(shí)現(xiàn)水平均勻大氣的反演. 利用9 月14 日早上8:30 電路1 的C 波段測(cè)量數(shù)據(jù)反演，反演的M剖面與適應(yīng)度變化如圖5 所示，為一個(gè)由蒸發(fā)波導(dǎo)層結(jié)、一個(gè)懸空波導(dǎo)層結(jié)和上層大氣構(gòu)成的復(fù)合剖面.

圖5 反演的水平均勻大氣M 剖面及其適應(yīng)度變化Fig. 5 Inversed M profile and variation of fitness function with iteration for horizontal homogeneous atmosphere

4) 水平不均勻大氣環(huán)境反演結(jié)果. 電路設(shè)置4 個(gè)剖面，兩個(gè)在電路兩端，兩個(gè)分別距兩端15 km，相鄰兩個(gè)剖面之間的大氣參數(shù)由兩側(cè)的剖面插值得到. 參數(shù)搜索仍采用表1 方案，適用性函數(shù)采用式(15). 利用9 月14 日早上8:30 電 路1 的X 波段雙向觀測(cè)數(shù)據(jù)，反演得到的M剖面與適應(yīng)度變化如圖6所示. 可以看到，海面蒸發(fā)波導(dǎo)層結(jié)存在和海岸上變性為線性層結(jié)的實(shí)際情形. 這個(gè)變化和海岸-海面-海岸位置、海陸氣團(tuán)流動(dòng)、海岸海陸風(fēng)均影響局地修正大氣折射率剖面有關(guān). 剖面1 和剖面4 因?yàn)槲挥诤０渡?，下午海陸風(fēng)效應(yīng)使得陸地地面的氣流向海面流動(dòng)，產(chǎn)生一個(gè)陸性大氣的薄層緊貼地面，這一層上面是變性的海面氣流層，再上面是逐漸接近正常大氣的空氣，形成了這兩個(gè)剖面. 而剖面2 是一個(gè)典型的蒸發(fā)波導(dǎo)剖面，是陸上氣流影響不到的海面空氣變化，海面附近水汽急劇下降形成了這個(gè)剖面. 剖面3 雖然也在海上，但是更靠近北面的中國(guó)大面積陸地，北方大陸氣團(tuán)的陸上氣流會(huì)沿其比剖面2 附近更平坦的陸地地形流向離岸更遠(yuǎn)的海面上，其上出現(xiàn)類似剖面1 和剖面3 的變性海面氣流層. 由剖面3 和剖面4 的比較也可以看出同一北方大陸氣團(tuán)的陸上氣流沿著它們附近更平坦的陸地地形流向離岸更遠(yuǎn)的海面上的效果.

圖6 反演的水平不均勻大氣M 剖面及其適應(yīng)度變化Fig. 6 Inversed M profile and variation of fitness function with iteration for horizontal inhomogeneous atmosphere

3.2 反演效果分析

由3.1 節(jié)的結(jié)果可以看出，本文方法可實(shí)現(xiàn)微波超視距水平均勻和水平不均勻大氣折射率環(huán)境的反演情況. 在水平均勻的情況下，獲得了包含蒸發(fā)波導(dǎo)和波導(dǎo)層結(jié)的復(fù)合大氣折射率剖面，反演迭代不到20 次即迅速收斂，適應(yīng)性函數(shù)達(dá)到10-4量級(jí)；而水平不均勻的反演中，得到的剖面是蒸發(fā)波導(dǎo)和波導(dǎo)層結(jié)在不同位置出現(xiàn)的水平演變，符合大氣海陸環(huán)流實(shí)際產(chǎn)生的影響，但是收斂效果較差，適應(yīng)性函數(shù)接近200.

再分析剖面中產(chǎn)生超視距遠(yuǎn)距離傳播的大氣波導(dǎo)層結(jié)的反演結(jié)果. 收集了觀測(cè)當(dāng)天8:00 附近香港探空數(shù)據(jù)，把反演結(jié)果中波導(dǎo)層結(jié)相關(guān)的參數(shù)列在表2 中進(jìn)行對(duì)比.

表2 中數(shù)據(jù)表明，探空測(cè)得的是波導(dǎo)結(jié)構(gòu)為無基礎(chǔ)層貼地波導(dǎo)，其波導(dǎo)頂高78 m，波導(dǎo)強(qiáng)度為50.21 M 單位. 對(duì)于水平均勻大氣環(huán)境的反演，得到了蒸發(fā)波導(dǎo)層結(jié)與大氣波導(dǎo)層結(jié)復(fù)合環(huán)境，反演所得波導(dǎo)頂高度為67.3 m、波導(dǎo)強(qiáng)度為63.91 M 單位；對(duì)于水平非均勻大氣環(huán)境而言，靠近陸地可獲得約65 m 以上的波導(dǎo)高度，波導(dǎo)強(qiáng)度則因剖面情況而異，蒸發(fā)波導(dǎo)層結(jié)的情況有50 M 單位以上、線性層結(jié)最強(qiáng)約達(dá)到30 M 單位. 可以看出，得到的水平均勻和水平不均勻反演各剖面的貼地波導(dǎo)高度和探空測(cè)量的貼地波導(dǎo)高度誤差分別是13.7%、17.7%、33.0%、2.2%和14%；反演的水平不均勻平均貼地波導(dǎo)高度的誤差為16.6%. 其中剖面3 與4 和實(shí)測(cè)剖面的位置更近，將它們和探空的剖面比較更合理，說明本文貼地波導(dǎo)高度的反演誤差平均為9.4%.

表2 反演剖面參數(shù)和實(shí)測(cè)數(shù)據(jù)對(duì)比Tab. 2 Comparison of inversed profile parameterswith radio sounding

4 結(jié) 論

針對(duì)微波超視距遠(yuǎn)距離傳播的復(fù)雜大氣環(huán)境反演問題，開展了原理和具體方法的研究，形成了一種適用于微波超視距遠(yuǎn)距離傳播線路觀測(cè)場(chǎng)景的水平不均勻分層復(fù)雜大氣折射率環(huán)境反演方法，并闡明了海陸交界地區(qū)水平不均勻分布大氣折射率反演的原理和方法，重點(diǎn)說明了和場(chǎng)景需求結(jié)合需要解決的關(guān)鍵問題.

和常見觀測(cè)線路取的多點(diǎn)傳播數(shù)據(jù)不同，本文用雙向傳輸?shù)奈⒉ǔ暰噙h(yuǎn)距離傳播測(cè)量數(shù)據(jù)，建立了以粒子群優(yōu)化算法為核心的反演方法，實(shí)現(xiàn)了基于收發(fā)位置不變電路上傳播觀測(cè)數(shù)據(jù)的水平均勻和不均勻大氣環(huán)境的反演，有利于形成解決觀測(cè)工程難度的技術(shù)方法.

雖然反演結(jié)果和收集的探空數(shù)據(jù)的比較分析表明本文方法是有效的，但是海岸-海面-海岸下墊面的水平不均勻問題還需進(jìn)行深入的研究，進(jìn)而獲得更高的反演速度和更穩(wěn)定的收斂結(jié)果.