利用差分體像素CT由TWINS衛(wèi)星ENA成像反演暴時環(huán)電流離子三維分布

顏偉男, 馬淑英

武漢大學(xué)電子信息學(xué)院， 武漢 430072

?

利用差分體像素CT由TWINS衛(wèi)星ENA成像反演暴時環(huán)電流離子三維分布

顏偉男, 馬淑英*

武漢大學(xué)電子信息學(xué)院， 武漢 430072

本文將體像素CT方法應(yīng)用于TWINS衛(wèi)星高能中性原子立體成像數(shù)據(jù)，重建環(huán)電流離子通量強度三維分布.為了克服衛(wèi)星ENA成像可能存在的儀器偏差造成反演誤差增大甚至迭代發(fā)散無解的問題，發(fā)展了差分ENA體像素CT反演方法.采用環(huán)電流離子分布理論模型和實際衛(wèi)星觀測幾何構(gòu)形進行模擬檢驗，證明了該方法的可靠性和優(yōu)越性.利用此方法分析2012年7月15日磁暴(minDst=-133 nT)主相期間TWINS 兩顆衛(wèi)星的實測ENA數(shù)據(jù)(能量范圍4 keV～50 keV)，反演得到環(huán)電流離子微分通量隨L值/緯度/地方時的三維分布與能譜.所得離子通量呈現(xiàn)顯著不對稱的部分環(huán)電流特征，主要分布在低緯-赤道區(qū)磁午夜前后至黎明前數(shù)小時，L值約為3.5～6.5的區(qū)域；將反演得到的磁午夜后赤道區(qū)環(huán)電流離子通量能譜，與THEMIS衛(wèi)星當(dāng)?shù)販y量得到的此區(qū)域同時的能譜作比較，兩者符合得很好；證明本文發(fā)展的差分體像素CT是由多衛(wèi)星ENA二維圖像重建暴時環(huán)電流離子分布的有效方法.

ENA成像； 差分CT； TWINS衛(wèi)星； 環(huán)電流； 磁暴

1 引言

單電荷高能離子與較冷中性原子之間發(fā)生共振電荷交換而生成高能中性原子(ENA, Energetic Neutral Atoms)是太陽系空間等離子體中一種比較普遍的現(xiàn)象(Fahr et al., 2007).高能中性原子衛(wèi)星成像測量可廣泛應(yīng)用于地球和其他行星的磁層與電離層，以及太陽日冕與太陽風(fēng)等離子體的探測(Gruntman, 1997).對于地球內(nèi)磁層，環(huán)電流中的高能離子通過與中性冷地冕氫間的共振電荷交換生成高能中性原子，它們不再受磁場與電場的約束，以原有能量和動量幾乎沿直線軌跡從磁層逃逸；這些中性原子攜帶著發(fā)生電荷交換處的能量離子性態(tài)(成份、能量、方位等)的消息，是環(huán)電流能量離子的遠程信息使者；利用專門的探測儀器可在衛(wèi)星上對環(huán)電流ENA進行全球遙感成像.20世紀(jì)90年代后期，瑞典發(fā)射了首個ENA成像衛(wèi)星Astrid，用專門設(shè)計的成像儀，從1000 km的低軌道上對磁層進行了短期的掃描成像(Barabash et al., 1997; C:son Brandt et al.，1997, 2001).21世紀(jì)初，美國實施了首個對磁層進行多手段全面成像的IMAGE衛(wèi)星計劃(Burch, 2000，2003), 其中包括多能段ENA成像，取得有價值的數(shù)據(jù)，深化了人們對磁暴和亞暴全球過程的認(rèn)識(Perez et al., 2001; C:son Brandt et al., 2002; Mitchell et al., 2003; Pollock et al., 2003)，揭示了暴時非對稱環(huán)電流的最強離子通量會出現(xiàn)在午夜后，而不是傳統(tǒng)認(rèn)為的總是出現(xiàn)在午夜前等新現(xiàn)象(Fok et al., 2003).中國科學(xué)家在21世紀(jì)初即開展了衛(wèi)星ENA成像的模擬研究(Shen and Liu, 2002)，中國雙星Double Star(Liu et al., 2005)的極軌道星TC-2上攜載了高能段ENA成像儀，并由TC-2 ENA成像數(shù)據(jù)反演獲得赤道面環(huán)電流離子分布(Lu et al., 2008, 2010).

雖然ENA成像是空間等離子體的一種重要遙測診斷工具，但是從ENA二維圖像定量反演其生成離子源的三維空間(加上投擲角和能量共五維)分布是一個具有挑戰(zhàn)性的難題.IMAGE和TC-2衛(wèi)星在極區(qū)上空幾個地球半徑以外對內(nèi)磁層進行ENA全球成像，都是單顆星觀測，一次掃描所生成的二維圖像的觀察視角有限，這使得人們難以由此提取得到可靠的離子三維空間分布.自從1987年Roelof(Roelof, 1987)利用ISEE-1衛(wèi)星觀測到的ENA(它們是被原本用來測量帶電粒子的儀器觀測到的)，采用環(huán)電流赤道離子分布零階模型等簡化假設(shè)，首次反演得到粗略的環(huán)電流離子分布以來，已發(fā)展了多種反演方法以獲取ENA的離子源分布，包括使用參數(shù)化離子分布模型求解模型參數(shù)的方法(Roelof and Skinner, 2000; Perez et al., 2000, 2012)和仿照“有約束的線性反演”(Twomey，1977)直接解得赤道離子通量強度(C:son Brandt et al., 2002; DeMajistre et al., 2004; Lu et al., 2010)，以及利用環(huán)電流理論模型的數(shù)據(jù)融合(data assimulation)方法(Nakano et al., 2008)，等等.但至今所進行的反演都是對赤道面二維環(huán)電流離子分布的重建，為了使求解得以進行，大都采用絕熱不變假設(shè)，借此約束赤道面與其他緯度上離子分布之間的定量關(guān)系而將三維分布的求解降為二維.然而磁暴主相/亞暴期間的環(huán)電流離子分布往往會打破絕熱不變.美國TWINS(Two Wide-angle Imaging Neutral-atom Spectro-meters)是世界上首個對地球磁層進行立體成像的星座(McComas et al., 2009；Goldstein & McComas, 2013)，它利用相互交叉莫爾尼亞軌道上的兩顆衛(wèi)星，在北極上空有利位置，同時從不同角度對磁層ENA進行全球掃描成像；這為由多幅二維圖像利用體像素(Voxel)CT技術(shù)重建環(huán)電流離子能譜的三維空間分布提供了有利條件.

本文以下將介紹我們發(fā)展的磁層ENA發(fā)射CT反演方法，以及在此基礎(chǔ)上，為克服多衛(wèi)星ENA成像存在的儀器不一致性造成反演誤差增大甚至迭代發(fā)散無解的問題，而引入的差分ENA體像素CT反演方法；然后利用離子分布模型和TWINS衛(wèi)星實測成像幾何構(gòu)形，對差分反演方法進行數(shù)值模擬檢驗；并利用該方法由衛(wèi)星實測ENA數(shù)據(jù)重建2012年7月15日磁暴期間的環(huán)電流離子三維分布與能譜，將其與THEMIS衛(wèi)星當(dāng)?shù)販y量結(jié)果進行比較，驗證此方法的可靠性.

2 基本原理與方法

2.1 磁層ENA發(fā)射CT

CT探測與反演的基本思想可以簡要地概括為：從物理上，對于一個有限的空間區(qū)域，如果有某種波(包括輻射)從所有各種不同方向穿過該區(qū)域，波與該區(qū)域介質(zhì)相互作用的路徑積分效應(yīng)能夠被觀測到，波與介質(zhì)相互作用的物理規(guī)律已知，那么可以由所有這些路徑積分效應(yīng)，確定該區(qū)域介質(zhì)(或輻射源)的分布；從數(shù)學(xué)上來說，CT探測與反演分別對應(yīng)獲得Fredholm型積分方程組和求解該方程組，理論基礎(chǔ)是Radon變換與求逆(Radon，1917；Deans，1984)，及其三維推廣(Chiu et al.，1980；Tuy，1983；徐繼生和鄒玉華，2003).路徑積分值在CT領(lǐng)域稱為投影，投影與待重建的目標(biāo)分布之間，由Furiuer切片定理或稱投影定理(Kak and Slaney，1998；徐繼生和鄒玉華，2003)直接聯(lián)系起來.

磁層衛(wèi)星ENA成像測量到的高能中性原子強度(計數(shù))C，是到達儀器的ENA微分通量Jena被儀器響應(yīng)函數(shù)A加權(quán)的多重積分：

Ci=∫Δt∫ΔE∫β∫ηAi(η,β,E,t)JenasinηdηdβdEdt，

(1)

式中Ai是儀器因子，表示t時刻，第i個像素對能量為E，以仰角η和方位角β入射的ENA粒子的響應(yīng).

ENA微分通量Jena含有兩個部分，一是較高高度上能量離子(對于我們所研究的情況，主要是環(huán)電流離子)與稀薄地冕氫發(fā)生電荷交換生成的ENA，稱為高高度發(fā)射(High Altitude Emission)，記作JHAE；另一部分是到達外逸層底較濃密大氣中的高能離子與氧原子發(fā)生電荷交換生成的，稱為低高度發(fā)射(Low Altitude Emission)，記作JLAE，即

(2)

其中

(3)

式中σO表示高能離子與中性氧原子共振電荷交換截面，σstrip表示ENA被碰撞剝離又生成離子的交換截面.這部分涉及復(fù)雜的多重碰撞和多種類型的電荷交換，需要單獨專門處理.本文主要關(guān)心較高高度上環(huán)電流ENA發(fā)射，其微分通量可近似表示為

(4a)

其中

(4b)

上兩式中nH是冷的背景氫中性原子數(shù)密度，σH表示高能離子與中性氫原子共振電荷交換截面，其上標(biāo)表示離子的種類，Jion是沿視線方向的高能離子微分通量，s是從衛(wèi)星出發(fā)的視線，沿著由仰角η和方位角β決定的方向，積分的上限是視線首次和外逸層底相交的地方，如果沒有相交則理論上為+∞，但針對實際的探測源區(qū)，積分上限是有限的.

由(4a)式可見，到達儀器的ENA微分通量Jena是沿觀測視線上低能中性原子數(shù)密度與高能離子強度(通量)及其電荷交換截面的乘積的路徑積分.這樣，如果電荷交換截面和地冕氫分布已知，且假定高能離子源存在于有限的空間，背景地冕氫為光學(xué)薄的稀疏介質(zhì)，并進一步假定中性原子成像測量在待測區(qū)域之外的多點從不同角度同時進行，那末從原則上就可以采用CT反演技術(shù)重建高能離子源的分布.此情況可歸屬于文獻Bates(Bates et al.，1983)中所述遙感CT，在某種意義上類似醫(yī)學(xué)中的單光子發(fā)射CT(Knoll, 1983).

由于衛(wèi)星觀測的投影數(shù)據(jù)不完備，不適合采取直接依據(jù)投影定理導(dǎo)出的富里葉積分反演算法，需采用離散的有限級數(shù)展開法，在處理問題的最開始就進行離散化(Censor, 1983).將反演區(qū)域劃分為三維網(wǎng)格即體像素，式(4a)離散化，省略左端的上標(biāo)HAE，代之以觀測視線的序號i，得到

(i=1,2,…,M，j=1,2,…,N),

(5)

其中i、j分別是視線序號和網(wǎng)格序號，M為視線數(shù)目，N為網(wǎng)格數(shù)目，sij是第i條視線在第j個網(wǎng)格內(nèi)的截距.

上述離散化ENA微分通量方程組可表示為以下矩陣形式：

(6)

其中

(7)

稱為投影矩陣，它依賴于反演區(qū)域網(wǎng)格與觀測視線族之間的幾何關(guān)系.于是，反演問題轉(zhuǎn)化為對以Jion為未知數(shù)的代數(shù)方程組(6)的求解.本文利用代數(shù)重建法(Algebraic Reconstruction Technique)求解方程組，以得到反演區(qū)域離子通量的三維分布.ART是一種迭代算法，在每一次迭代過程中計算觀測值與上一輪迭代結(jié)果之差，然后根據(jù)這個差值對反演迭代結(jié)果進行反復(fù)修正，最終得到一個穩(wěn)定的收斂解.ART迭代公式如下：

(j=1,2,…,N)，

(8)

式中k表示第k次迭代，μk為松弛因子，取值在0到2之間，對含測量誤差的數(shù)據(jù)，適當(dāng)選取松弛因子的值可以改善重建圖像的質(zhì)量和迭代效率.一般情況下，迭代過程中松弛因子可以取一個定值1.由于ART算法并不能保證迭代結(jié)果為正，因而從實際的物理意義考慮，在ART算法中需要加上離子通量的非負(fù)約束.

2.2 差分體像素CT方法

隨著時間的推移，星載探測儀器性能會有所降低，從而探測結(jié)果可能存在因儀器而異的偏差，造成多衛(wèi)星測量的不一致性，這將導(dǎo)致相互矛盾的ENA測量代數(shù)方程組，使得反演誤差增大，甚至迭代反演不收斂.為了解決這個問題，本文引入差分體像素 CT方法，此方法的思想與電離層/等離子體層無線電射線CT中為消除多站點接收不同GPS信號所得TEC的不一致性而使用的差分TEC反演方法(Xiao et al., 2012)是相似的.我們假定這個偏差是加性誤差，并假定同一個探測儀器的偏差為一個常數(shù)，在這樣的簡化假設(shè)下，(5)式可表示為

(9)

式中εq為衛(wèi)星q的儀器偏差.對于某一衛(wèi)星，選擇一條合適的觀測視線m作為參照，用沿任意第l條視線測量到的ENA通量，減去沿參考視線m的ENA通量，得到第l條視線差分ENA通量，它不再含有這個偏差，從而消除了任意射線l上的儀器偏差，可用公式表示為

(10)

由(9)式得

(11)

此為差分體像素CT所依據(jù)的公式，其中用來進行反演的ENA通量變?yōu)椴罘滞?，投影矩陣亦作相?yīng)差分轉(zhuǎn)換.

2.3 反演區(qū)域及其體像素劃分

對于本文所用的磁層ENA成像衛(wèi)星觀測區(qū)域，ENA離子源區(qū)主要包括環(huán)電流離子源區(qū)和極光帶沉降離子源區(qū).我們假定，環(huán)電流離子源位于偶極地磁場中2

考慮到環(huán)電流粒子被束縛在磁殼上，并依賴于磁緯和磁地方時，我們采用地磁偶極子場模型，將上述反演區(qū)域依L值、偶極磁緯λ、以及磁地方時MLT劃分為三維網(wǎng)格(即體像素)，采用太陽磁坐標(biāo)系(Solar Magnetic Coordinates)，徑向地心距離r(以地球半徑RE為單位)與L值和偶極磁緯λ之間由下式聯(lián)系：

(12)

根據(jù)衛(wèi)星錐束掃描的視線分辨率，以及視線與反演區(qū)域的相對位置，合理確定體像素的大小.鑒于本文考慮的衛(wèi)星位于北半球極區(qū)上空，角度分辨率4°×4°，錐束掃描視線在北半球較密，在南半球較疏，因此環(huán)電流區(qū)的緯度網(wǎng)格采取不等間隔劃分.低緯赤道區(qū)域與較高緯度區(qū)不同，南北兩半球的較高緯度區(qū)又各不相同.在北緯30°以上較高緯度區(qū)，間隔為12.5°，而在南半球30°以上較高緯度區(qū)，間隔則為20°，南北緯30°以下赤道低緯區(qū)域間隔為15°.又鑒于緯度越高，相鄰L殼之間的空間距離越小，對于北緯30°以上環(huán)電流區(qū)，L值間隔取為1；對于其他緯度的環(huán)電流區(qū)，L值間隔為0.5.所有反演區(qū)域的磁地方時MLT的網(wǎng)格間隔均為1小時(經(jīng)度間隔15°).

對于磁緯55°～78°的北半球極區(qū)，不再對緯度進行劃分，磁地方時0—24 h劃分間隔與環(huán)電流區(qū)相同，為1小時，高度以5000 km為界劃分為高低兩個區(qū)域.

這樣劃分網(wǎng)格的結(jié)果，反演區(qū)域體像素總數(shù)為1608個，略少于有效視線(即觀測方程)數(shù)；對于所關(guān)注的環(huán)電流反演區(qū)，兩顆衛(wèi)星所處的位置比較有利時，每個體像素內(nèi)至少有數(shù)條，多至數(shù)十條不同方向的視線穿過.

2.4 其他簡化假設(shè)與模型/公式

2.4.1 離子通量投擲角分布

環(huán)電流離子通量嚴(yán)格來說是各向異性的，即隨投擲角而變，這使CT反演變得復(fù)雜；另一方面，衛(wèi)星當(dāng)?shù)販y量表明，對于本文所考慮的幾keV至幾十keV的中等能量離子，離子通量投擲角分布的各向異性相對較弱(Fok et al., 1996)，而且在擾動期間特別是大磁暴期間，離子分布各向異性有弱化的趨勢(Chen et al., 1998 )，因此本文采用離子通量投擲角分布為各向同性的假設(shè)，這也是以往由ENA反演環(huán)電流的經(jīng)典研究采用的假設(shè)(Roelof, 1987；C:son Brandt et al., 2001, 2002).各向異性投擲角分布情況下的處理將另文討論.

2.4.2 地冕氫原子密度模型

(13)

式中nH表示氫原子數(shù)密度，單位cm-3；r是地心距離，以地球半徑為單位；φ是地方時角，正午為零；θ是SM坐標(biāo)系下從z軸出發(fā)的極角.其中

(14)

H的單位為地球半徑，此量給出氫原子密度分布的晝夜不對稱性.

2.4.3 電荷交換截面

環(huán)電流區(qū)域能量離子主要有H+，O+，其他成分(如N+，He+，He++等)含量很低，均可以忽略.環(huán)電流區(qū)域高能H+，O+與地冕中性H發(fā)生電荷交換產(chǎn)生ENA，其電荷交換截面隨能量E變化，如式(15)和圖1所示(Lindsay and Stebbings,2005)：

圖1 H+和O+與H的電荷交換作用截面Fig.1 Cross sections for charge transfer of ions H+(black cross) and O+ (gray circle) with atom H

(15)

2.4.4 高能離子成分的簡化

如式(4b)所示，忽略次要離子成分后，電荷交換截面加權(quán)的能量離子通量包含氫離子和氧離子兩種成分，精確的反演應(yīng)先將氫ENA和氧ENA加以分離，然后分別反演得到氫離子和氧離子的分布.中等能量不同成分ENA的分離需要特殊的技術(shù)來實現(xiàn)，我們對這一技術(shù)的有關(guān)工作將另文介紹.本文通過使用一個氧和氫離子通量強度比率δ=JO/JH，將(4a)式簡化為

(16)

取δ=0,即假定環(huán)電流能量離子全部是氫離子，本文在數(shù)值模擬反演中取δ=0；在實測數(shù)據(jù)反演中，針對本文所考慮的大磁暴事件，采用δ=1/4.

3 差分CT的數(shù)值模擬

在本文這一節(jié)，采用環(huán)電流離子分布理論模型和實際衛(wèi)星觀測幾何構(gòu)形，對差分ENA體像素CT反演進行數(shù)值模擬，并就存在儀器偏差與不存在偏差情況，對差分與非差分CT反演效果進行比較，以檢驗差分CT方法的可靠性和優(yōu)勢.

圖2 數(shù)值模擬過程框圖Fig.2 Block diagram of numerical simulation

數(shù)值模擬的框圖如圖2所示.

3.1 離子微分通量模型

(17)

其中κ是kappa指數(shù)(Christon et al.,1991)，κ=5.5，e=(1+1/κ)κ+1≈2.962 ；L0=7.3是環(huán)電流離子注入?yún)^(qū)域外邊界，E0是注入離子特征能量，對于本文模擬所考慮的氫離子，此特征能量約為7 keV，J0是離子注入最大微分通量，取J0=2.0×105(s·sr·cm2·keV)-1(Shen and Liu，2002)；因子h(φ)反映隨地方時的分布：

(18)

考慮能量E=20 keV的環(huán)電流離子，L取值范圍2～7.5，磁緯 -55°<λ<55°，地方時角φ取值0°～360°.

3.2 衛(wèi)星位置與掃描視線

我們選取2012年7月15日磁暴主相期間16∶45—17∶00UT，美國TWINS衛(wèi)星(McComas et al., 2009)的位置和姿態(tài)及相應(yīng)掃描視線用于本文的模擬.表1給出SM坐標(biāo)系中兩顆衛(wèi)星在上述時間段大約中間時刻的位置，可見兩顆星都在北半球較高的中緯區(qū)，距地心5.3RE以遠，均處于本文考慮的反演區(qū)域以外；所處磁地方時相差約14個小時，比較接近同一個磁子午面；兩顆衛(wèi)星進行扇形錐束掃描的中軸線接近正交，有比較大的視角跨度.

圖3a和圖3b分別給出在所考慮的時間段的中點左右TWINS兩顆衛(wèi)星在磁子午面和赤道面內(nèi)的觀測視線覆蓋示意圖.圖3a中，綠色和藍色射線分別代表S1、S2的探測視線，紅色星號是探測視線與反演網(wǎng)格邊界的的交點.圖3b是兩顆衛(wèi)星探測視線在緯度-15°～15°網(wǎng)格內(nèi)的截距在赤道面上的投影，紅色為衛(wèi)星S1，藍色為衛(wèi)星S2.

表1 TWINS兩顆衛(wèi)星的位置(2012-07-15, 16∶45—17∶00 UT)Table 1 Positions of TWINS satellites (2012-07-15, 16∶45—17∶00 UT)

可以看出，由于兩顆衛(wèi)星所在位置相隔較遠，在幾乎相對的地方時扇區(qū)，因此覆蓋區(qū)域的射線交叉情況較好，有比較大的角度覆蓋，有利于CT反演.在昏側(cè)半球的環(huán)電流區(qū)域視線覆蓋情況好于晨側(cè)半球.由于兩顆衛(wèi)星均處于北半球，所以南半球高緯區(qū)域的視線覆蓋不理想.

3.3 模擬結(jié)果

3.3.1 ENA圖像的正演模擬

圖4給出在上文所述條件下，對TWINS兩顆衛(wèi)星觀測分別進行模擬計算得到的ENA成像.在模擬ENA成像中，能量取E=20 keV，因此由(15)式得到電荷交換截面σH=5×10-16cm2.圖像由魚眼視角給出，圖中的顏色值代表ENA微分通量大小，兩幅圖使用相同的色標(biāo).圖像中間的白色實線圓圈代表地球所在位置，紅色和黃色磁力線分別代表正午和黃昏.圖像外圈所標(biāo)示的角度代表儀器探測的方位角β，取值-90°～270°，徑向標(biāo)示的角度代表儀器仰角η，中心仰角為90°.

圖3 TWINS視線覆蓋示意圖.(a) 子午面；(b) 赤道面Fig.3 Illustration of projection of satellite line-of-sight on (a) Meridian plane and (b) Equatorial plane

圖4 魚眼視角下TWINS兩顆衛(wèi)星模擬得到的ENA圖像.紅色磁力線為正午，黃色為黃昏Fig.4 Simulated ENA images in view of ‘fish-eye’ from the two satellites of TWINS, in which the dipole magnetic field lines are shown, with red line for those at noon and yellow at dusk

3.3.2 由模擬的ENA圖像反演離子源分布

利用本文第2節(jié)所述的普通CT和差分CT方法，分別對3.3.1節(jié)正演模擬得到的ENA微分通量進行反演，ART迭代初始分布取為0，使用非負(fù)約束.

圖5給出普通ENA CT重建得到的反演區(qū)域離子微分通量分布(中間一列)和差分CT重建得到的離子分布(右邊一列)；為了比較，圖中左邊一列給出了正演模擬所采用的離子微分通量分布模型.普通CT和差分CT反演迭代次數(shù)相同.圖中上面一行是赤道區(qū)(-15°<λ<15°)離子微分通量隨L(2～7.5)和地方時的分布，左側(cè)為午夜，右側(cè)為正午；下面一行是子午面00∶00—12∶00MLT內(nèi)(左側(cè)為午夜)離子微分通量隨L(2～7.5)和緯度(-55°<λ<55°)的分布.所有圖像均使用相同的色標(biāo).

由圖5所示結(jié)果可見，在沒有加入儀器偏差的情況下，差分CT和普通CT都能很好地重建離子微分通量分布.其中赤道面的反演結(jié)果(圖5c和圖5e)在04∶00—08∶00 MLT扇區(qū)的外圈(L大于3)有明顯的誤差.結(jié)合圖3b可以發(fā)現(xiàn),這個區(qū)域只能被一顆衛(wèi)星探測到，從而沒有足夠的交叉視線，與CT反演需要的條件相差較遠.子午面的反演結(jié)果中(圖5d和圖5f)，南半球區(qū)域也存在觀測視線覆蓋不足情況，因此這些區(qū)域也存在較大誤差.迭代終止時兩種反演方法的總體相對均方根誤差都降到0.14左右(如圖6).可以說，沒有儀器偏差的情況下，兩種反演方法都可以有效地重建待反演區(qū)域的離子微分通量分布.

3.3.3 由存在儀器偏差的模擬ENA圖像反演離子源分布

隨著衛(wèi)星運行時間增加，探測儀器的性能會有所改變，測量參數(shù)出現(xiàn)偏移，比如儀器響應(yīng)因子的增大或減小，且難以及時進行在軌標(biāo)定，從而探測結(jié)果可能存在儀器偏差，特別是兩顆衛(wèi)星偏差不一致.這種偏差和不一致性造成反演涉及的代數(shù)方程組出現(xiàn)相互矛盾的方程，導(dǎo)致反演迭代誤差增大甚至迭代不收斂.我們對儀器偏差作一簡單的加性誤差模擬，在3.3.1節(jié)得到的兩顆衛(wèi)星的模擬ENA分布(圖4)中，分別加入以下偏差值：在TWINS-1的模擬ENA中減去其均值的十分之一；在TWINS-2中加上其均值的五分之一.使用相同的迭代條件，分別用普通CT和差分CT兩種方法進行重建，得到的反演結(jié)果如圖7所示.圖中左邊一列仍為離子分布模型，中間一列為利用普通CT方法重建得到的離子微分通量分布，右邊一列為利用差分CT方法得到的通量分布.

圖5 (a，b) 模擬的離子微分通量分布模型; (c，d) 利用普通CT反演得到的離子微分通量分布； (e，f) 利用差分CT反演結(jié)果.上面一行是赤道面內(nèi)分布，下面一行是00—12 MLT子午面內(nèi)分布Fig.5 Comparison of the ion flux distribution from (a，b) model, (c，d) classic CT retrieval and (e，f) differential CT retrieval. The top raw is for equatorial plane, the bottom raw is for meridian plane. Red parts of the earth indicate the dayside

由圖7可以發(fā)現(xiàn)，加入儀器偏差后，普通CT(圖7(c，d))反演得到的結(jié)果雜亂無章，無法重建反演區(qū)域離子微分通量分布，而差分CT(圖7(e，f))可以有效地消除儀器偏差的影響，得到接近原始模型(圖7(a，b))的結(jié)果，其結(jié)果也與圖5相近.通過比較兩者的相對均方根誤差隨著迭代次數(shù)的變化(圖8)，可以發(fā)現(xiàn)在加入儀器偏差后，普通CT迭代反演過程中，均方根誤差越來越大，迭代無法收斂；而差分CT反演結(jié)果在10次迭代內(nèi)快速收斂，其相對均方根誤差最后降到0.14左右.

圖6 無儀器偏差情況下反演的相對均方根誤差.(a) 普通CT反演； (b) 差分CT反演Fig.6 Relative root-mean-square error for reconstruction without instrument bias by (a) classic CT method and (b) differential CT method

圖7 對于具有儀器偏差的模擬ENA圖像，反演得到的離子微分通量分布. (a，b) 模擬的離子微分通量分布模型 ; (c，d) 利用普通CT反演得到的離子微分通量分布； (e，f) 利用差分CT反演結(jié)果.上面一行是赤道面內(nèi)分布，下面一行是00—12 MLT子午面內(nèi)分布Fig.7 The same as Fig.5 but the instrument biases are applied in ENA images which are further used to reconstruct the RC ion flux distribution

圖8 加入儀器偏差情況下反演的相對均方根誤差.(a) 普通CT反演；(b) 差分CT反演Fig.8 Relative root-mean-square error of inversion by (a) classic CT method and (b) differential CT method from images with instrument bias included

由此可見，在兩顆衛(wèi)星ENA圖像具有不一致的固定偏差時，普通CT反演難以重建反演區(qū)域離子微分通量分布，而差分CT方法能消除這個偏差帶來的影響，很好的重建得到反演區(qū)域離子微分通量分布，表明差分CT方法處理此類問題的優(yōu)越性.

4 實測數(shù)據(jù)反演

在這一節(jié)，我們由TWINS衛(wèi)星實測ENA數(shù)據(jù)，利用差分CT方法反演環(huán)電流離子微分通量強度的三維分布與能譜，并與THEMSIS衛(wèi)星的離子能譜當(dāng)?shù)販y量數(shù)據(jù)進行對比.

4.1 ENA數(shù)據(jù)

選取2012年7月15日磁暴主相期間，TWINS兩顆衛(wèi)星在大約15 min里(16∶45～17∶00 UT)約12輪掃描平均的ENA圖像，由其重建環(huán)電流離子分布.單次掃描可能會因干擾而不完整或者有較大噪音，所以選取多輪掃描的累計平均以提高圖像的質(zhì)量.如果采用每顆衛(wèi)星的12～13幅圖像來進行反演，雖然觀測視線的數(shù)目明顯增多，但每顆衛(wèi)星觀測視線的角度變化很小(總計約2°以內(nèi))，對于改進觀測視角覆蓋范圍的作用不大，同時又損失了用來進行反演的ENA圖像的質(zhì)量.二者權(quán)衡，采用多次掃描平均為宜.考慮能量為4, 8, 12, 16, 20, 25, 30, 50 keV的離子.磁暴期間Dst指數(shù)如圖9所示，其中灰色豎線表示進行環(huán)電流反演的時間段，處在臨近主相極大時.

圖9 磁暴期間Dst變化，灰色豎線是本文研究的時間段Fig.9 Dst index during the storm on July 14—16, 2012. The gray bar is for concerned time period

TWINS衛(wèi)星的位置如3.2節(jié)中的表1所示；這里取15 min時間段大約中點時刻衛(wèi)星所在的位置.在本文分析的時間段內(nèi)，TWINS兩顆衛(wèi)星的地心距離變化約0.2RE，磁緯變化約2°，磁地方時改變約12 min；此變化相對于本文反演區(qū)域及其網(wǎng)格大小來說很小，不足以影響反演結(jié)果，因此忽略衛(wèi)星的位置變化，以TWINS衛(wèi)星在中間時刻的位置來代表.此間衛(wèi)星掃描視線對反演區(qū)域的覆蓋情況如3.2節(jié)的圖3所示.所用原始數(shù)據(jù)是由美國西南研究所提供的經(jīng)過平滑的ENA微分?jǐn)?shù)通量二維圖像數(shù)據(jù)，取自cdaweb網(wǎng)站.參照McComas等(2012)給出的方法，對圖像進行消除背景，得到用來進行反演的來自兩顆衛(wèi)星的一對ENA圖像，圖10給出能量為20 keV時的一個示例，在儀器坐標(biāo)系中給出，徑向為仰角η，取值30°～90°；角向為方位角β，取值-90°～270°.在其他能量上做相同處理，得到分別在4, 8, 12, 16, 20, 25, 30, 50 keV上的TWINS兩顆衛(wèi)星的用來進行反演的圖像.

圖10 消除干擾背景后的TWINS ENA圖像Fig.10 TWINS ENA images with background interferences removed

圖11 差分CT重建的環(huán)電流區(qū)離子分布.(a) 赤道區(qū)H+微分通量強度隨L值和MLT的分布； (b) 01∶00 MLT子午面內(nèi)H+離子微分通量強度隨L值和緯度的分布.圖中左側(cè)為夜晚Fig.11 Reconstructed RC ion distribution at different energy levels by differential voxel CT method. (a) The distribution of H+ differential flux versus L-value and MLT in equatorial plane； (b) The distribution of H+ differential flux versus L-value and Mlat in meridian plane. Red parts of the earth indicate the dayside

4.2 反演結(jié)果

利用消除背景等預(yù)處理后的ENA數(shù)據(jù)，反演重建環(huán)電流離子微分通量的三維分布.如2.4節(jié)所述，在反演中假設(shè)氧離子和氫離子的通量強度比為常數(shù)0.25，離子通量分布為各向同性；依然采用ART迭代算法，初始分布取作零.

圖11a給出在幾個不同能量(4，8，16，30 keV)上，赤道區(qū)(-15°—15°N)H+離子微分通量隨L值和MLT的分布.圖11b給出在01—13 MLT子午面內(nèi)不同能量的離子通量隨L值和緯度的分布.圖中的空缺部分(白色網(wǎng)格)是由于缺少衛(wèi)星射線覆蓋以及數(shù)據(jù)不可用；所有各幅圖像都使用相同的色標(biāo).

從圖11可以看出，反演所得離子通量呈現(xiàn)出顯著不對稱的部分環(huán)電流特征，強的離子通量主要分布在磁午夜前后至黎明前數(shù)小時L值約3.5～6.5的區(qū)域，對于較低的能量(如4，8 keV)，在午夜前和黃昏扇區(qū)也有較強的通量；從子午面內(nèi)的分布可以看出，強的通量主要出現(xiàn)在低緯-赤道區(qū)，與環(huán)電流的一般特征相符.

4.3 與THEMIS衛(wèi)星當(dāng)?shù)販y量的比較

在本文研究的磁暴期間，THEMIS-D衛(wèi)星剛好處于環(huán)電流區(qū)域內(nèi)：r=5.31RE，MLT=01∶45，Mlat=15.07°S附近，該衛(wèi)星攜載的ESA (ElectroStatic Analyzer)和SST(Solid State Telescope)兩種儀器聯(lián)合，能夠探測6 eV到300 keV的離子能通量譜(其間20～30 keV有一小段能量空缺).圖12給出THEMIS-D在2012年7月15日16∶00到17∶00 UT期間ESA和SST測量到的離子能通量譜，橫坐標(biāo)是時間，縱坐標(biāo)是能量，顏色代表能通量強度.其中ESA數(shù)據(jù)的能量范圍是6 eV到20 keV，時間分辨率為97 s；SST數(shù)據(jù)的能量范圍是30 keV到300 keV，時間分辨率為3 s.

圖12 THEMIS-D衛(wèi)星實測離子能譜隨時間的變化.上半部為SST 測量(30 keV～300 keV)，下半部為ESA測量(6 eV～20 keV)Fig.12 Ion energy spectra measured in situ onboard THEMIS-D satellite during 16∶00—17∶00 UT on July 15, 2012. The top one is measured by SST (at 30 keV～300 keV) and the bottom is by ESA (6 eV～ 20 keV)

圖13 TWINS ENA反演得到的離子能通量譜(紅色粗點線) 與THEMIS-D實測結(jié)果的比較Fig.13 Comparison of the ion energy spectra between retrieved (red line) from TWINS ENA data and measured by THEMIS-D

我們將THEMIS-D實地探測結(jié)果與利用ENA圖像反演得到的離子能通量譜作比較，以驗證反演結(jié)果的可靠程度.首先將16∶45—17∶00 UT期間反演得到的在5.0

5 總結(jié)與討論

本文發(fā)展了由衛(wèi)星ENA二維圖像重建磁層環(huán)電流能量離子通量三維分布的CT反演方法，綜合考慮電荷交換物理機制與衛(wèi)星掃描視線幾何特征，比較合理地確定了反演區(qū)域及其體像素劃分方案；針對衛(wèi)星ENA成像可能存在的儀器偏差問題進一步發(fā)展了差分ENA 體像素CT方法，并通過正/反演數(shù)值模擬，證明差分CT可以比較有效地降低儀器整體偏差造成的大的反演誤差，避免迭代反演發(fā)散，比普通CT具有更寬的適用性和優(yōu)越性.將此差分ENA體像素CT方法應(yīng)用于2012年7月15日磁暴事件中TWINS兩顆衛(wèi)星的ENA數(shù)據(jù)，反演環(huán)電流離子微分通量的三維分布，得到了合理的結(jié)果，這是在不受絕熱不變條件限制下反演得到的環(huán)電流區(qū)離子三維分布，至今國內(nèi)外已有反演大都使用絕熱不變假設(shè)重建赤道面內(nèi)二維分布(C:son Brandt et al., 2002; DeMajistre et al., 2004; Lu et al., 2010; Perez et al., 2012).所得離子分布呈現(xiàn)顯著晝夜不對稱，夜側(cè)遠大于白天，這與磁暴主相期間的環(huán)電流為部分環(huán)電流的特征相呼應(yīng)；對于所分析的磁暴事件，較強的離子通量主要出現(xiàn)在低緯-赤道區(qū)磁午夜前后至黎明前數(shù)小時，這與傳統(tǒng)認(rèn)為的強的環(huán)電流離子通量只出現(xiàn)在黃昏至午夜前不同，但與IMAGE衛(wèi)星揭示的現(xiàn)象(Fok et al., 2003)相符；而且我們的反演結(jié)果還顯示，環(huán)電流強度的地方時分布依賴于離子的能量，在此次事件中較低能量(10 keV以下)的離子在午夜前和黃昏扇區(qū)也有較強的通量.將反演所得磁午夜后赤道區(qū)環(huán)電流離子通量能譜，與THEMIS衛(wèi)星當(dāng)?shù)販y量得到的此區(qū)域同時的能譜作比較，兩者符合得很好，證明本文發(fā)展的差分ENA體像素CT是由多衛(wèi)星ENA二維圖像重建暴時環(huán)電流離子分布的有效方法.

應(yīng)該指出，由于磁層成像衛(wèi)星探測到的ENA的離子源并非單一的被地磁場捕獲的環(huán)電流離子，還有近地等離子體片中的離子，以及到達較低高度大氣層的沉降離子(來自環(huán)電流和等離子體片)；與高能離子發(fā)生電荷交換的背景冷中性氣體既有光學(xué)薄的地冕氫，也有濃密的光學(xué)厚大氣層中的氧原子，后者涉及多重碰撞、往復(fù)的不同類型電荷交換與能量損失過程.因而下一步工作還需要更加合理細(xì)致地區(qū)分和確定反演區(qū)域以及相應(yīng)區(qū)域電荷交換過程的物理性質(zhì)與理論公式，以便對沉降離子源引起的低高度光學(xué)厚介質(zhì)中ENA發(fā)射做出合理修正，對近地等離子體片ENA發(fā)射源的影響作出評估.高空磁場與地冕氫分布模型對反演具有一定的影響，采用更接近實際的非偶極子磁場模型，以及利用TWINS衛(wèi)星的Ly-α成像反演得到的地冕氫分布，將會有助于改進反演結(jié)果.另外，盡管TWINS兩顆衛(wèi)星從不同角度掃描給出的兩幅ENA二維圖像(即層析技術(shù)的所謂“投影”)，提供了比單顆衛(wèi)星更為豐富的離子源分布信息，有可能進行CT 反演，但是投影數(shù)據(jù)還是很不完備的，還需要更為充分的合理先驗知識的約束以求得更為可靠的物理上合理的解，反演算法還有進一步優(yōu)化的空間；在條件允許情況下，發(fā)展ENA衛(wèi)星編隊CT遙測，應(yīng)是從根本上解決反演數(shù)據(jù)不完備的最好辦法.

致謝 感謝美國TWINS衛(wèi)星計劃團隊與PIs通過互聯(lián)網(wǎng)提供TWINS ENA數(shù)據(jù)；感謝J. Goldstein博士對本文工作的支持與有益建議.感謝美國國家宇航局通過SPDF提供THEMIS衛(wèi)星數(shù)據(jù).

Barabash S, C: son Brandt P, Norberg O, et al. 1997. Energetic neutral atom imaging by the Astrid microsatellite.AdvancesinSpaceResearch, 20(4-5): 1055-1060.

Bates R H T, Garden K L, Peters T M. 1983. Overview of computerized tomography with emphasis on future developments.Proc.oftheIEEE, 71(3): 356-372.Burch J L. 2000. Image mission overview.SpaceScienceReviews, 91(1-2): 1-14.

Burch J L. 2003. The first two years of image.SpaceScienceReviews, 109: 1-24.

C:son Brandt P, Barabash S, Norberg O, et al. 1997. ENA imaging from the Swedish micro satellite astrid during the magnetic storm of 8 February, 1995.AdvancesinSpaceResearch, 20(4-5): 1061-1066.

C:son Brandt P, Barabash S, Roelof E C, et al. 2001. Energetic neutral atom imaging at low altitudes from the Swedish microsatellite Astrid: Extraction of the equatorial ion distribution.JournalofGeophysicalResearch:SpacePhysics, 106(A11): 25731-25744.

C:son Brandt P, Demajistre R, Roelof E C, et al. 2002. IMAGE/high-energy energetic neutral atom: Global energetic neutral atom imaging of the plasma sheet and ring current during substorms.JournalofGeophysicalResearch:SpacePhysics, 107(A12): 1454, doi:10.1029/2002JA009307.

Censor Y. 1983. Finite series-expansion reconstruction methods.Proc.oftheIEEE, 71(3): 409-419.

Chen M W, Roeder J L, Fennell J F, et al. 1998. Simulations of ring current proton pitch angle distributions.JournalofGeophysicalResearch:SpacePhysics, 103(A1): 165-178.

Chiu M Y, Barrett H H, Simpson R G. 1980. Three-dimensional reconstruction from planar projections.JournaloftheOpticalSocietyofAmerica, 70(7): 755-762.

Christon S P, Williams D J, Mitchell D G, et al. 1991. Spectral characteristics of plasma sheet ion and electron populations during disturbed geomagnetic conditions.JournalofGeophysicalResearch:SpacePhysics, 96(A1): 1-22. Deans S R. 1984. Book-review—the radon transform and some of its applications.Science, 223(4635): 484.

DeMajistre R, Roelof E C, C:son Brandt P, et al. 2004. Retrieval of global magnetospheric ion distributions from high-energy neutral atom measurements made by the IMAGE/HENA instrument.JournalofGeophysicalResearch:SpacePhysics, 109(A4), doi: 10.1029/2003JA010322.Fahr H J, Fichtner H, Scherer K. 2007. Theoretical aspects of energetic neutral atoms as messengers from distant plasma sites with emphasis on the heliosphere.ReviewsofGeophysics, 45(4), doi: 10.1029/2006RG000214.

Fok M C, Moore T E, Greenspan M E. 1996. Ring current development during storm main phase.JournalofGeophysicalResearch:SpacePhysics, 101(A7): 15311-15322.

Fok M C, Moore T E, Wilson G R, et al. 2003. Global ena image simulations.SpaceScienceReviews, 109(1-4): 77-103. Goldstein J, McComas D J. 2013. Five years of stereo magnetospheric imaging by TWINS.SpaceScienceReviews, 180(1-4): 39-70. Gruntman M. 1997. Energetic neutral atom imaging of space plasmas.ReviewofScientificInstruments, 68(10): 3617-3656.

Kak A C, Slaney M. 1988. Principles of Computerized Tomographic Imaging. New York: IEEE press.

Knoll C F. 1983. Single-photon emission computed tomography.Proc.oftheIEEE, 71(3): 320-329, doi: 10.1109/PROC.1983.12590.

Lindsay B G, Stebbings R F. 2005. Charge transfer cross sections for energetic neutral atom data analysis.JournalofGeophysicalResearch:SpacePhysics, 110(A12), doi: 10.1029/2005JA011298. Liu Z X, Escoubet C P, Pu Z, et al. 2005. The Double Star mission.AnnalesGeophysicae, 23(8): 2707-2712.

Lu L, McKenna-Lawlor S, Barabash S, et al. 2008. Iterative inversion of global magnetospheric ion distributions using energetic neutral atom (ENA) images recorded by the NUADU/TC2 instrument.AnnalesGeophysicae, 26(6): 1641-1652. Lu L, McKenna-Lawlor S, Barabash S, et al. 2010. Comparisons between ion distributions retrieved from ENA images of the ring current and contemporaneous, multipoint ion measurements recorded in situ during the major magnetic storm of 15 May 2005.JournalofGeophysicalResearch, 115(A12), doi: 10.1029/2010JA015770.

McComas D J, Allegrini F, Baldonado J, et al. 2009. The two wide-angle imaging neutral-atom spectrometers (TWINS) NASA mission-of-opportunity.SpaceScienceReviews, 142(1-4): 157-231.

McComas D J, Buzulukova N, Connors M G, et al. 2012. Two wide-angle imaging neutral-atom spectrometers and interstellar boundary explorer energetic neutral atom imaging of the 5 April 2010 substorm.JournalofGeophysicalResearch, 117(A3), doi: 10.1029/2011JA017273.

Mitchell D G, C:son Brandt P, Roelof E C, et al. 2003. Global imaging of O+from IMAGE/HENA.SpaceScienceReviews, 109(1-4): 63-75.

Nakano S, Ueno G, Ebihara Y, et al. 2008. A method for estimating the ring current structure and the electric potential distribution using energetic neutral atom data assimilation.JournalofGeophysicalResearch:SpacePhysics, 113(A5), doi: 10.1029/2006JA011853.

Perez J D, Fok M C, Moore T E. 2000. Deconvolution of energetic neutral atom images of the earth′s magnetosphere.SpaceScience Reviews, 91(1-2): 421-436.

PerezJD,KozlowskiG,C:son-BrandtP,etal. 2001.InitialionequatorialpitchangledistributionsfrommediumandhighenergyneutralatomimagesobtainedbyIMAGE. Geophysical Research Letters, 28(6): 1155-1158.

PerezJD,GrimesEW,GoldsteinJ,etal. 2012.EvolutionofCIRstormon22July2009. Journal of Geophysical Research, 117(A9),doi: 10.1029/2012JA017572.

PollockCJ,C:son-BrandtP,BurchJL,etal. 2003.Theroleandcontributionsofenergeticneutralatom(ENA)imaginginmagnetosphericsubstormresearch. Space Science Reviews, 109(1-4): 155-182.

RadonJ. 1917. überdieBestimmungvonFunktionendurchirherIntegralwertel?ngsgewisserMannigfalligkeiten. Berichete S?chsische Akademie der Wissenschaften, 69: 262-277.RairdenRL,FrankLA,CravenJD. 1986.GeocoronalimagingwithDynamicsExplorer. Journal of Geophysical Research: Space Physics, 91(A12): 13613-13630.

RoelofEC. 1987.Energeticneutralatomimageofastorm-timeringcurrent. Geophysical Research Letters, 14(6): 652-655.

RoelofEC,SkinnerAJ. 2000.ExtractionofiondistributionsfrommagnetosphericENAandEUVimages. Space Science Reviews, 91(1-2): 437-459.

ShenC,LiuZX. 2002.PropertiesoftheneutralenergeticatomsemittedfromEarth′sringcurrentregion. Physics of Plasmas, 9(9) : 3984.TuyHK. 1983.Aninversionformulaforcone-beamreconstruction. SIAM Journal on Applied Mathematics, 43(3): 546-552.TwomeyS. 1977.Introductiontothemathematicsofinversioninremotesensingandindirectmeasurements. ∥Develop.Geomath.,vol. 3,NewYork:ElsevierScienceLtd.XiaoR,XuJS,MaSY,etal. 2012.AbnormaldistributionofionosphericelectrondensityduringNovember2004super-stormby3DCTreconstructionsfromIGSandLEO/GPSobservations. Sci. China Technol. Sci., 55(5): 1230-1239.XuJS,ZouYH. 2003.Reconstructionformulaoftime-dependent3-Dcomputerizedionospherictomography. Chinese J. Geophys. (inChinese), 46(4): 438-445,doi: 10.3321/j.issn:0001-5733.2003.04.002.

附中文參考文獻

徐繼生, 鄒玉華. 2003. 時變?nèi)S電離層層析成像重建公式. 地球物理學(xué)報, 46(4): 438-445, doi: 10.3321/j.issn:0001-5733.2003.04.002.

(本文編輯 胡素芳)

Reconstruction of 3-D distribution of storm-time RC ions from TWINS ENA images with differential voxel CT method

YAN Wei-Nan, MA Shu-Ying*

SchoolofElectronicInformation,WuhanUniversity,Wuhan430072，China

Onboard satellite imaging of ENA (Energetic Neutral Atom) produced through resonant charge exchange of energetic ions with cold neutral atoms is a useful tool for remote sensing of the space plasma. However, a quantitative retrieval of 3-D (amounting to 5-D if adding pitch angle and energy) distribution of parent ions from 2-D images of ENA is really a challenging task. ENA emission image data have been used to reconstruct global 2-D distributions of the earth magnetospheric ring current (RC) ions by various methods. These methods are mostly under assumption of retaining adiabatic invariant by which the 3-D problem is reduced to 2-D inversion in the equatorial plane. However, during great magnetic storms the adiabatic invariant law may probably be violated and some irregular constructions of ion distribution would come forth. It is thus motivated to develop a new method without the limitation of adiabatic invariant to retrieve the 3-D spatial distributions of storm-time RC energetic ions.A method of volume pixel (voxel) computerized tomography (CT) is initially explanted in the 3-D reconstruction of energetic RC ion distributions from ENA stereoscopic images measured on multi-satellites. To overcome the large error and divergence in inversion iteration caused by discordant instrument biases, a differential ENA voxel CT (DEVCT) method is developed. The performance of both the classical voxel CT method and the DEVCT method are demonstrated by numerical simulations using an analytic RC ion model along with real satellite observational configuration. The DEVCT method is then implemented to real ENA emission images measured by TWINS (Two Wide-angle Imaging Neutral-atom Spectrometers) during a major magnetic storm with minimumDst index about -133 nT. The TWINS Mission has a twain of spacecrafts flying on two widely-separated Molniya orbits, making ENA images simultaneously from largely different viewing angles out of the RC region. The concerned ions′ energies range 4～50 keV exampled at 8 different levels. The ion species of hydrogen is considered by taking a constant ion flux ratio of oxygen to hydrogen.The numerical simulations proved the reliability and effectiveness of the differential ENA voxel CT method, with total error (relative RMS) approaching as low as 0.14 for the modeled ion flux distribution. The retrieved distributions of the storm-time RC ion intensity (isotropic differential flux) from TWINS ENA observations seem corresponding to an asymmetric partial RC, being located mainly around the midnight favoring the post-midnight sectors ranging fromL=3.5 to 6.5 in the equatorial plane. This MLT feature contradicts with the conventional opinion that peak ion flux of storm-time enhanced RC would appear only in the pre-midnight sectors, while it supports once again that ring current enhancement could peaked at midnight-dawn during the main phase of the magnetic storm. The results also show that the RC ion distributions with MLT depend on the energy level. Furthermore, the resulted RC ion flux spectrum from remotely sensed ENA images is compared with the satellite in-situ measurements by THEMIS-D in the same time interval and within nearly the same spatial region, showing a very good consistent with each other. All the results cited above demonstrated that the DEVCT method developed in this study is an effective technique to reconstruct the storm-time RC ion distribution from multi-satellite ENA images.This study provides an effective new method for retrieving 3-D storm-time RC ion distributions from multi-satellite remotely sensed ENA images. However, it should be noted that the TWINS data is still poorly incomplete for making CT retrieval. The pitch angle anisotropy of the ion flux distribution should be taken into account in the next inversion algorithm. And the intensity ratio of oxygen to hydrogen should be accurately determined to obtain much more significant results on the RC dynamics.

ENA imaging; Differential CT; TWINS Mission; Ring current; Magnetic storm

10.6038/cjg20150701.

國家自然科學(xué)基金重點項目(41431073)資助.

顏偉男，1986年生，男，武漢大學(xué)博士研究生，研究方向空間探測與信息處理技術(shù).E-mail:wnyan@whu.edu.cn

*通訊作者 馬淑英，武漢大學(xué)教授，博士導(dǎo)師，主要研究領(lǐng)域空間物理與探測.E-mail:syma@whu.edu.cn

10.6038/cjg20150701

P352

2015-04-10，2015-06-03收修定稿

顏偉男, 馬淑英. 2015. 利用差分體像素CT由TWINS衛(wèi)星ENA成像反演暴時環(huán)電流離子三維分布.地球物理學(xué)報，58(7):2221-2234,

Yan W N, Ma S Y. 2015. Reconstruction of 3-D distribution of storm-time RC ions from TWINS ENA images with differential voxel CT method.ChineseJ.Geophys. (in Chinese),58(7):2221-2234,doi:10.6038/cjg20150701.