空間目標監視雷達大氣折射修正技術研究

唐新杰 劉默然 喬建東 周晨

（武漢大學電子信息學院, 武漢 430072）

0 引 言

無線電波在大氣層中傳播時，由于大氣介質不均勻分布的影響產生折射效應，使得電波射線發生彎曲、傳播速度小于光速，導致雷達等測控系統在定位、導航時產生一定的偏差[1]。隨著技術的發展和各種應用場景的需要，大氣電波折射誤差逐漸成為制約航天測定軌、空間目標監視等遠程探測系統精度的重要因素，因此需要對測控系統探測目標測量值的折射誤差予以修正。

現有的雷達折射誤差修正方法可以分為射線描跡法、近似修正法和新型修正法三種，其中射線描跡法因為其精度高、實用性強等優點在研究和工程實踐中都有廣泛的應用[2]。折射誤差修正精度不僅和選用的算法有關，還取決于選取的大氣折射率剖面的準確度，學者們對此進行了大量的研究，通過模型或一些探測方法對大氣折射率剖面進行建模。王寧等[3]引入徑向基函數(radial basis function, RBF)神經網絡算法反演大氣折射率，驗證了利用微波輻射計對流層電波折射進行修正的可行性；徐艷等[4]對中近程雷達應用場景中高精度目標定位誤差進行了研究分析；張瑜等[5]針對下墊面復雜地區的三維大氣結構，采用直接探測法和全國大氣剖面模型組合的方法，有效提高了電波射線折射率分布的準確度；喬江等[6]對比研究了折射率剖面模型和映射函數在對流層折射誤差修正中的精度表現；李川川等[7]結合微波輻射計和GNSS接收機研制的電波折射誤差修正設備，實現了S波段雷達實時、高精度折射誤差修正；劉琨等[8]提出了適用于高軌道目標折射誤差修正的空間投影法和自適應網格法。

對于航天目標的雷達探測，電波傳播主要受對流層和電離層兩方面的影響。而對于P波段雷達而言，由于其頻率相對較低，受電離層的影響更為顯著，也明顯比對流層影響更為嚴重。同時，低仰角相對于高仰角的誤差更大，而低仰角意味著水平距離更遠，其環境參數與探測點的環境參數差別也更大，因此對低仰角的誤差修正就更為困難，這也一直是相關領域的難點之一。本文結合NRLMSISE-00大氣模型對海拔60 km以下的對流層折射率剖面進行建模，并利用TEC數據同化的方式建立電離層現報模型，得到更精確的電離層實時信息，完成電離層折射率剖面的建模。在此基礎上，利用射線描跡法對P波段雷達探測數據進行修正。與基于國際參考電離層（International Reference Ionosphere, IRI）模型的電離層折射率剖面相比，利用TEC數據同化現報模型得到的電離層折射率剖面可以極大提高修正精度。

1 建模方法

1.1 折射誤差計算方法

射線描跡法從斯奈爾定律和費馬原理出發，基于數學幾何關系嚴格推導而來，具有較高的精度，因此成為工程中廣泛應用的修正方法[9]。由于精度高，該方法也常被用來與其他的修正模型進行對比驗證。圖1給出了球面分層電波傳播軌跡示意圖，C為地心，O為測控站位置，T為目標位置。

圖1 球面分層電波射線軌跡幾何圖形Fig.1 Geometry of spherical layered radio ray

高軌道遠程目標一般沉浸于電離層之中的空間環境，由斯奈爾定理可以推導出雷達探測獲得的目標視在距離Re計算式為

式中：rI為電離層底的高度；r0為測控站處的地心半徑，r0=a+h0，a為地球半徑，h0為測控站海拔高度；n為電波射線軌跡上變化的折射指數；r為目標至地心之間的距離；n0為地表處的折射指數；θ0為目標視在仰角。

測控站與目標之間地心張角φ的計算式為

目標真實仰角α0的計算式為

目標真實距離R0的計算式為

由上述式子計算得到距離折射誤差?R和仰角折射誤差ε分別為

1.2 基于電離層現報同化模型的電離層折射率剖面

數據同化按照理論基礎分類主要分為兩類：一類是基于統計學的估計理論，如最優插值、卡爾曼濾波、集合卡爾曼等；一類是基于變分理論的方法[10]。卡爾曼濾波算法是常見的線性濾波和預測方法中的一種，是由匈牙利科學家Rudolph E.Kalman在1960年提出的。它是最早出現的順序同化算法，一般認為它是順序同化算法的理論基礎。在近年的研究中，卡爾曼濾波算法越來越多地被應用于電離層數據同化中。

根據卡爾曼濾波理論，對觀測資料的基本同化過程可表示為：

式中：為t+1時刻的背景場，由背景電離層模型給出；M為狀態轉換矩陣；為t時刻的分析場，即同化后得到的電子密度；為t+1時刻的背景場誤差協方差矩陣；為t時刻的分析場誤差協方差矩陣；Q為模型誤差方差矩陣；為t+1時刻的分析場；K為增益矩陣，起到觀測數據對背景場的調整作用；Y為觀測場；H為觀測算子，代表觀測量與模型參量之間的關系；R為觀測場誤差協方差矩陣；為t+1時刻的分析場誤差協方差矩陣。

本文采用目前被廣泛應用于電離層研究的經驗模型-IRI模型作為背景模型，版本為IRI2016。IRI由國際無線電科學聯盟和空間研究委員會共同研發和維護，其數據來源于全球180多個電離層探測站的資料以及衛星的資料，主要模擬了電子密度、離子成分、離子溫度、電子溫度等一系列電離層參量的全球分布[11]。背景場的經緯度范圍選取以雷達測站為中心的探測覆蓋范圍及周邊區域，緯度為15°N~55°N，經度為70°E~140°E，水平格點分辨率設定為1°×1°，高度范圍為60~1 500 km， 步長為10 km，因此背景場誤差協方差矩陣Pb的大小為419 184 ×419 184。參與三維數據同化的觀測數據主要為中國地殼運動觀測網絡(Crustal Movement Observation Network of China, CMONOC)的261個GNSS監測臺站提供的電離層TEC數據，GNSS接收機的采樣間隔為30 s。TEC是GNSS信號傳播路徑對路徑上電子密度的線積分，所以觀測算子H即為GNSS信號傳播路徑在電離層柵格中穿越的長度。

基于統計的協方差模型在估計背景誤差協方差矩陣時會引入虛假相關，即在物理上不相關或空間上距離較遠的狀態變量之間存在虛假相關。為了降低虛假相關對同化過程的影響，提升同化效率，本文應用了協方差局地化(covariance localization，CL)方法[12]。CL方法采用一個基于距離的相關系數ρ與背景誤差協方差矩陣作舒爾積，來替代原有的背景協方差矩陣。相關函數采用常用的GC (Gaspari-Cogn)五階分段多項式函數，表達式如下：

式中：z為物理空間中兩個格點之間的距離；c為長度尺度。

1.3 基于空間格點化獲取的射線路徑上折射率剖面

傳統的射線描跡法忽略水平方向的不均勻性，往往采用垂直方向上的空間環境剖面來代表不同仰角下觀測條件，即傳統的射線描跡法適用于下墊面平坦、大氣水平相對均勻的地區，而對于高目標低仰角的目標觀測而言，電波實際傳輸路徑較長，跨越的區域范圍更大，經歷的空間環境也更加多變復雜，因而不再滿足水平均勻的假設條件[13]。

考慮到大氣不均勻性產生的影響以及三維空間格點化大氣高度分布形式，根據射線傳播路徑獲取射線路徑上位于不同空間柵格對應的空間環境參數值，進而計算得到更加真實的射線傳播路徑上的折射率剖面[14]。

三維空間中的電波射線的路徑可以通過大地坐標系（L，B，H）和站心地平直角坐標系（E，N，U）描述，在站心地平直角坐標系中以雷達測站位置作為站心P0建立坐標系，電波射線軌跡以站心為起點向外延展到達目標，電波射線路徑上的每個位置P可以用（EP，NP，UP）來表示，如圖2所示。

圖2 地心空間直角坐標系和站心地平直角坐標系Fig.2 Geocentric space and station center ground horizontal Cartesian coordinate system

NRLMSISE-00大氣模型是由美國海軍研究實驗室設計研發的全球大氣經驗模型，描述了從地面到熱層高度范圍內（0~1 000 km）的中性大氣密度、溫度等大氣物理性質，是目前使用最多的大氣模型之一[15]。該模型在長時間的觀測數據基礎上建立并不斷更新，主要數據源為火箭探測數據、衛星遙感數據和非相干散射雷達數據等，通過采用低階球諧函數擬合大氣性質隨經緯度、年周期、半年周期、地方時的變化而建立的[16]。

本文使用NRLMSISE-00模型來描述和構建雷達電波射線經過的對流層大氣環境，通過公式(12)可以計算得到對流層的折射率剖面：

式中：P為壓強，hPa；T為大氣溫度，K；e為水汽壓，hPa。

NRLMSISE-00模型和IRI模型是基于大地坐標系來計算獲取全球空間環境參數的，為了保持空間匹配，NRLMSISE-00模型也按照經緯度1°×1°的柵格空間分辨率把雷達測站視野范圍內的三維空間格點化處理，高度范圍為0~60 km，步長為1 km。通過站心地平直角坐標系和大地坐標系的轉換，可以將電波射線路徑經過的每個位置投影到通過模型構建的空間柵格中，進而獲取射線路徑上相關的空間環境參數P、T、e、ne關于高度的分布（如示意圖3所示），在此基礎上進一步計算得到射線傳輸路徑上的折射率剖面。電離層折射指數n與空間環境中電子密度的關系可表示為

圖3 射線路徑上不同空間柵格中環境參數選取Fig.3 Selection of environmental parameters in different spatial grids on the ray path

式中：ne為電子密度，el/m3；f為電波頻率，Hz。

站心地平直角坐標系轉換為空間直角坐標系的公式為[17]

空間直角坐標系轉換為大地坐標系的公式為

式中：N為卯酉圈的半徑，為第一偏心率，為橢球長半軸和短半軸。

2 結果分析

2.1 電子密度剖面同化模型效果比較

圖4給出了2015-03-17T9:00LT的數據同化實例，其中圖4(a)為背景模式的區域TEC分布結果，圖4(b)為進行數據同化后的TEC分布結果。可以看出，數據同化方法實現了觀測數據和背景模型有效融合，數據同化后的區域TEC分布相比IRI模型直接獲取的結果展現了更為豐富的細節，更加符合實際。

圖4 2015-03-17T9:00LT數據同化前后區域TEC分布的對比Fig.4 Comparison of regional ionospheric TEC maps using data assimilation on 2015-03-17T9:00LT

電離層的折射率高度分布依賴于電子密度這一關鍵的環境參量，同時電子密度剖面可以有效地描述電離層垂直結構的細節。為了驗證三維同化模型數據同化的效果，選取處于同化區域內的北京昌平站的電離層測高儀SAO格式觀測數據對同化前后的電子密度剖面進行對比，SAO觀測數據來自空間環境地基綜合監測網。Digisonde系列數字測高儀內嵌的軟件，可自動度量標定電離圖得到電離層特征參數和電子密度剖面，生成SAO數據[18]。圖5給出了北京昌平站上空的電子密度剖面對比結果，可以直觀地看出，在大多數時刻同化后的電子密度剖面與真實的觀測值更為接近，電子密度峰值與測高儀探測結果相差更小，說明同化算法可以有效地把TEC觀測數據融合到背景場，改善提升背景場的數據精度。

圖5 2015-12-01北京昌平站上空的電子密度剖面同化結果對比Fig.5 Comparison of electron density profile assimilation results above Changping station, Beijing on 2015-12-01

2.2 電波折射誤差修正效果比較

為了驗證本文提出的利用空間格點化柵格獲取電波射線路徑上折射率剖面以及結合IRI模型數據同化進行電波折射誤差修正的效果，選取某臺工作于P波段的雷達在年積日（day of year, DOY）分別為71、72、73對高軌道衛星目標的三次跟蹤實測記錄作為實驗數據，目標在雷達跟蹤視野范圍內仰角歷經了從大到小的變化，俯仰范圍為4°~80°，南北跨越緯度約23°，東西跨越經度約42°。使用GPS精密軌道數據作為測量基準來確定和比較不同方法對各種距離仰角條件下目標的折射誤差修正精度。

本文利用射線描跡法對比三種電離層折射率剖面構建方法對修正結果的影響：1）垂直剖面法，即通過IRI模型構建的水平均勻的電離層折射率剖面；2）射線路徑剖面法，即通過IRI模型構建的隨傳播路徑變化的水平不均勻的電離層折射率剖面；3）數據同化法，即利用電離層同化現報模型構建的電離層折射率剖面。基于以上的折射率剖面，利用射線描跡法可以計算得到折射誤差修正量和剩余殘差，雷達測距測角誤差分別為雷達探測視在距離和視在仰角與對應目標精密軌道值之差，剩余殘差為經修正后的測距測角誤差。

圖6～8分別給出了三次觀測跟蹤中三種不同方法應用于折射誤差修正中對距離誤差修正的表現。可以看出，隨著觀測仰角降低，三種方法計算的誤差修正量均隨之增大，在低仰角處基于模型數據同化獲取射線路徑剖面計算得到的距離折射誤差修正量與雷達測距折射誤差有更好的一致性，相應地，修正殘差結果（子圖b）顯示利用數據同化法進行距離折射誤差修正的殘差明顯更小，并且修正殘差伴隨觀測仰角變化的趨勢更加平緩，很大程度上消除了觀測仰角的影響，表明了數據同化方法的有效性。

圖6 DOY=71時三種方法距離折射誤差修正結果Fig.6 Correction results of distance refraction error by 3 methods when DOY=71

圖7 DOY=72時三種方法距離折射誤差修正結果Fig.7 Correction results of distance refraction error by 3 methods when DOY=72

圖8 DOY=73時三種方法距離折射誤差修正結果Fig.8 Correction results of distance refraction error by 3 methods when DOY=73

圖9～11分別給出了三次觀測跟蹤中三種不同方法應用于折射誤差修正中對仰角誤差修正的表現。其低仰角情況下誤差更為顯著的結論對于仰角誤差而言仍然成立，最大仰角誤差可達0.4°。對比三種方法在仰角誤差修正中的差別可以發現，三種方法的修正結果依然是數據同化法更優秀，射線路徑剖面法次之，垂直剖面法最差，但是不同于距離誤差，三種方法對于仰角誤差的修正效果差距并沒有那么明顯，或者說數據同化方法對仰角誤差的修正相對于其他兩種方法效果并不顯著，這個問題還有待更深入的研究。

圖9 DO Y=71時三種方法仰角折射誤差修正結果Fig.9 Correction results of elevation angle refraction error by 3 methods when DOY=71

圖10 DOY=72時三種方法仰角折射誤差修正結果Fig.10 Correction results of elevation angle refraction error by 3 methods when DOY=72

圖11 DOY=73時三種方法仰角折射誤差修正結果Fig.11 Correction results of elevation angle refraction error by 3 methods when DOY=73

為了量化對比三種方法的修正結果，對利用上述三種方法進行雷達電波折射誤差修正的殘差進行統計分析，結果見表1、2。其中，數據同化方法距離折射誤差修正平均殘差分別為15.43、16.97、13.64 m，相比于未進行同化的射線路徑剖面方法的修正精度提升了41.2%、35.5%、44.6%，這表明，在三次雷達觀測修正中，通過數據同化方法進行電波折射誤差修正能極大地提高修正精度。同時，數據同化法的修正結果均方根誤差(root mean square error,RMSE)也最小，這表明修正后的殘差更為“平坦”，也正如前文所說，對于低仰角的誤差修正得到了很大改善。而對于仰角誤差而言，如表2，其提升結果有限。

表1 不同修正方法距離折射誤差修正殘差統計分析Tab.1 Statistical analysis of the correction residuals of distance refraction error

表2 不同修正方法仰角折射誤差修正殘差統計分析Tab.2 Statistical analysis of the correction residuals of elevation refraction error

3 結 論

本文針對傳統折射誤差修正方法忽略大氣水平不均勻性而采用測站單點位置垂直方向的折射率剖面進行電波折射誤差修正的局限性，提出了空間格點化大氣和電離層模型同時結合數據同化獲取電波射線路徑上電子密度剖面的方法。基于某臺P波段雷達的外場實測數據和精密軌道數據對三種方法應用于折射誤差修正的精度進行了比較分析，實驗結果表明，通過空間格點化獲取射線路徑上的折射率剖面來進行電波折射誤差修正相比于利用大氣和電離層模型直接獲取測站單點的折射率剖面來進行修正有一定的改善，但仍有較大的修正剩余誤差，難以滿足工程應用中的精度需求。基于空間格點化大氣和電離層模型結合數據同化的方法可以有效地增強電離層電子密度分布的準確性，并顯著降低雷達測距測角的折射誤差修正殘差。研究結果驗證了利用大氣和電離層模型數據同化方法提升電波折射誤差修正精度是一個可靠有效的途徑。