宇宙塵埃軌跡探測器信號提取方法及反演精度影響因素分析

馮桃君，焦子龍,2，姜利祥,2，鄭慧奇，姜海富,2，劉宇明,2，王 鵬，劉學超，李 濤

（1.北京衛星環境工程研究所； 2.可靠性與環境工程技術重點實驗室；3.航天東方紅衛星有限公司：北京 100094； 4.中國航天科技集團有限公司，北京 100048）

0 引言

微米和亞微米級的宇宙塵埃粒子通常來自彗星、小行星、月球、星體間的碰撞碎片和行星際物質[1-2]。宇宙塵埃暴露在等離子體、太陽紫外輻射、高能粒子、電磁場等組成的復雜空間環境中，運動過程中會攜帶電荷，因此塵埃軌跡是行星和星際磁場、太陽和行星引力、等離子體阻力及輻射壓力等因素共同作用的結果[3-4]。Xie 等[5]綜述并研究宇宙塵埃的運動方式發現：小行星塵埃由于Poynting-Robertson 效應會以螺旋的方式緩慢接近太陽；β-流星體主要受到太陽輻射壓的作用，以雙曲線軌跡逃逸出太陽系；Ulysses 號航天器上的塵埃探測器發現木星塵埃粒子受到磁層電場的加速并以大于100 km/s 的速度逃逸出木星。并且指出，宇宙塵埃的軌跡和電荷信息有助于確定粒子起源，辨別不同起源的粒子，揭示帶電塵埃粒子與行星際空間環境的相互作用。

對宇宙塵埃粒子進行探測的方法很多是利用塵埃探測器的原位測量[6-8]，高速粒子在探測器內部發生碰撞電離，可通過電離產生的總電荷量和電脈沖信號的上升時間分別推算出粒子質量和碰撞速度；并通過對粒子飛行時間的分析得出粒子的化學成分。該項技術已在Vega、Wind、Cassini、STEREO等航天器上用于探測塵埃的質量、成分和速度大小[9-11]，但速度方向只能通過統計數據分析得到，要求探測器進行全方位掃描，限制了碰撞探測器的應用。Grün 等[12]提出了塵埃天文學概念，涉及宇宙塵埃粒子的起源以及粒子化學成分、電荷、速度、軌跡的測量。塵埃軌跡探測器（dust trajectory sensor,DTS）是塵埃天文學的關鍵探測組件，當帶電粒子穿過DTS 內部的金屬絲陣列時，會產生誘導電荷，可通過分析金屬絲輸出的電信號提取粒子的電荷量和軌跡信息[1,13]。DTS 與不同的儀器結合可構成塵埃望遠鏡[14]、離子質量分析器[15]和靜電質量探測器[1]等，適用于不同速度、不同尺寸的塵埃顆粒探測。

DTS 在幾十年前就已被提出[5]，其幾何構型已經過不斷優化。Auer 等[1-2]利用COULOMB 軟件對不同幾何構型的DTS 進行數值仿真，分析金屬絲的幾何參數、間距及數量等對電荷分布、信號靈敏度的影響，獲得了適用于速度在2～5 km/s 的微米級塵埃探測的DTS 結構參數；并得出DTS 的電荷噪聲比（charge to noise ratio, QNR）≥6.25 才能從背景噪聲中探測到塵埃信號。Xie 等[5]同樣利用COULOMB 軟件，發現探測器的QNR 越小電荷量探測誤差分布范圍越大，當QNR=10 時，塵埃電荷量探測誤差為-1.5%～5%。為提高儀器對塵埃電荷量的探測靈敏度，Li 等[16]提出了分段式低電容金屬絲陣列結構，使塵埃電荷量探測閾值降至0.2 fC。Voronov 等[17]介紹了一種基于網格面和金屬絲陣列面交替的DTS 設計，不僅可以測量帶電微流星體的電荷量、軌跡，還能計算微流星體的質量。

宇宙塵埃大都分布在亞微米到厘米量級，速度在每秒幾千米到幾十千米不等，DTS 的探測精度也因塵埃性質不同而顯示差異。本文基于典型的DTS結構，針對塵埃軌跡和電荷量的測量建立DTS數值仿真模型，研究從DTS 金屬絲信號提取塵埃粒子速度矢量和電荷量的方法，重點分析塵埃速度、帶電量和入射方向對算法精度的影響，有助于探索適用于DTS 探測的塵埃范圍。

1 DTS 結構及工作原理

DTS 一般由中間4 個金屬絲陣列面和分別位于頂部和底部的2 個屏蔽柵網組成，每個金屬絲陣列面包含7～16 條相互平行的金屬絲，兩相鄰平面的金屬絲方向正交。如圖1[5]所示，兩屏蔽柵網接地；每條金屬絲都獨自與一個電荷靈敏放大器（charge sensitive amplifier, CSA）連接，CSA 通道的輸出信號被瞬態記錄器以一定頻率采集。當帶電塵埃粒子穿過DTS 時，會產生誘導電荷。根據同一平面金屬絲信號特性可確定粒子穿過當前平面的一維坐標（x1、y2、x3、y4），如圖2[5]所示，進而可推導出塵埃粒子的速度矢量和電荷量[1,5]。

圖1 DTS 結構示意[5]Fig.1 Structure schematic of DTS[5]

圖2 帶電塵埃粒子在DTS 中的穿越軌跡示意[5]Fig.2 Trajectory schematic of charged particles travelling in DTS[5]

本文研究所用DTS 模型含4 個平行于xy平面的金屬絲陣列面，每個平面包含7 條金屬絲，相鄰平面間隔40 mm，金屬絲直徑0.4 mm，長140 mm，間隔20 mm； DTS 的整體尺寸為160 mm×160 mm×200 mm（長×寬×高），如圖3 所示。平面編號從上到下依次增加，平面1、3 的金屬絲平行于y軸，平面2、4 的金屬絲平行于x軸。模型坐標原點位于結構的幾何中心，圖中以紅點標記。

圖3 DTS 數值模型Fig.3 Numerical model of DTS

探測器金屬絲的誘導電荷量與塵埃粒子的入射位置、角度、探測器尺寸和構型有關。Auer 等[1]給出了當帶電量為Q的球形塵埃顆粒穿過DTS時，金屬絲i的誘導電荷qi的近似表達式

式中：ri為塵埃到金屬絲i的垂直距離；j(1/rj)代表求和遍歷了DTS 所有的金屬絲。塵埃的位置(xp,yp,zp)隨時間變化，因此塵埃到金屬絲的距離和金屬絲的誘導電荷均為時間的函數。塵埃電荷量可由所有金屬絲誘導電荷相加近似得到。

Auer 等[1]在地面試驗中采集到的隨時間變化的DTS 輸出信號如圖4 所示，圖中信號來自兩金屬絲陣列平面，每陣列包含16 條金屬絲。輸出信號的形狀反映金屬絲與塵埃粒子的距離，信號強度反映金屬絲帶電量，信號峰值表明塵埃粒子近距離地從某金屬絲旁穿過該金屬絲所在平面。

圖4 Auer 在塵埃實驗室記錄到的DTS 信號[1]Fig.4 DTS signals recorded in the dust laboratory by Auer[1]

假設高速塵埃粒子在DTS 內的軌跡為一條直線，由塵埃粒子進、出DTS 的時間（tin、tout）和位置坐標（(xin, yin, zin)、(xout, yout, zout)）可以確定塵埃粒子的速度和方向。通過式(2)可由4 個一維坐標確定塵埃粒子進、出DTS 的坐標(xin, yin)、(xout, yout) [5]：

粒子速度的大小為[5]

粒子入射角為[5]

其中，θx、θy分別為塵埃軌跡在xz、yz平面的投影與z軸的夾角。

2 塵埃參數提取方法

2.1 計算流程

為提取塵埃粒子的速度矢量和電荷量，在已知探測器幾何構型的前提下，表征塵埃軌跡和電荷量需要7 個獨立參量：塵埃電荷量（Q）；塵埃進、出DTS 的時刻（tin、tout）；塵埃穿過4 個金屬絲陣列平面的一維坐標（x1、y2、x3、y4），故設參數矢量為P=(Q,tin,tout,x1,y2,x3,y4)。文獻[5]提出了一種基于探測數據與仿真數據擬合的參數提取方法，其反演流程如圖5 所示。

圖5 參數提取反演流程Fig.5 Inversion process for parameter extraction

參數提取的具體步驟為：

1）初始化參數矢量P0。塵埃初始電荷量Q可設置為DTS 所有金屬絲電荷量之和；x1、y2、x3、y4根據式(1)估算得到。文獻[5]中并未介紹tin和tout的取值方法，本文設DTS 的采樣頻率為r（單位為MS/s，即106sample/s），對應的信號采樣個數為n，并將第1 個采樣信號的時刻t1設置為0 作為參考，則tin∈[-1/r, 0]、tout∈[(n-1)/r,n/r]，單位為μs，tin和tout的初值在各自范圍內均勻隨機取值。

2）初始化參數的誤差范圍ΔP0。ΔQ為±10%，Δt為±1/r，Δx或Δy為±10 mm。

3）數據擬合。從P0±ΔP0的范圍中隨機選取一組P=(Q,tin,tout,x1,y2,x3,y4)，根據式(1)以及DTS采樣頻率仿真計算DTS 各金屬絲的采樣信號，然后計算仿真信號相對實測信號的誤差，

式中：χ2為仿真信號誤差；N為數據總個數；Dm為第m個實測數據；Sm為第m個仿真數據。

重復上述隨機選取參數P即χ2計算103次，χ2值隨P變化，χ2最小值對應的參數集為P1,min。

4）優化參數矢量P1。由于χ2也是單一參數的函數，依次改變P1,min中第l（l=1, 2, ···, 7）個參數P(l)的值，使χ2最小的參數值成為該參數的新估計值P1(l)。

5）縮小參數的誤差范圍ΔP1。ΔP1取2ΔP0/3和3 |P1-P0|中的大者。

6）重復步驟3）～5），直到χ2的變化小于0.1%。

2.2 數值仿真與結果分析

為驗證2.1 節塵埃參數提取方法的有效性，針對典型算例進行計算分析。塵埃預設軌跡為入射點坐標(13, -22, 100)、出射點坐標(-7, 38, -100)，入射角θx=-5.7°、θy=16.7°，z方向速度分量vz=5 km/s，塵埃穿過各金屬絲平面的一維坐標為x1=9 mm，y2=2 mm，x3=1 mm，y4=26 mm，涵蓋了塵埃與金屬絲之間的近距離、中等距離和遠距離。另外，金屬絲輸出信號用誘導電荷與塵埃電荷之比表示，設塵埃電荷為1 C。設置QNR=10，采樣頻率為r=10 MS/s，塵埃穿過DTS 的時間為t=40 μs，因此每個金屬絲有n=400 個采樣信號。設塵埃進入DTS 的時間tin=0，則tout=40 μs，參數真值Ptrue=(1, 0, 40, 9, 2, 1,26)，本文用下標“opt”和“true”分別表示參數的最優估計值和真值。離塵埃近的金屬絲將聚集大部分誘導電荷，因此在信號分析中，為節省計算時間只對每個平面中離塵埃最近的2 條金屬絲的信號進行擬合，這樣用于擬合的信號一共有3200 個。實際的觀測信號包含了電荷靈敏放大器產生的隨機噪聲，為模擬真實信號，在式(1)計算出的電荷信號中添加正態分布的隨機噪聲，其均值為μ=0，標準差σ=1/QNR。根據上述仿真條件，隨塵埃位置zp變化的無噪聲金屬絲模擬信號如圖6(a)所示，圖中虛線表示各平面對應的z軸坐標，頂端數字為平面序號。從圖中可看出，當塵埃穿過平面時與金屬絲距離很近（如y2=2 mm，x3=1 mm），該金屬絲信號幅值驟升，隨后下降；與此同時，塵埃另一側相距較遠的金屬絲信號出現下凹。當塵埃穿過平面時與兩側金屬絲距離相差不大（如x1=9 mm），兩金屬絲信號形狀相似。圖6(b)所示是疊加噪聲后的金屬絲信號，用于模擬DTS 的實測信號。

圖6 塵埃粒子從(13, -22, 100)到(-7, 38, -100)穿越DTS的模擬信號Fig.6 Simulated signals generated by dust particles travelling from (13, -22, 100) to (-7, 38, -100) through DTS

按上述塵埃軌跡與電荷信號的分析方法，經過6 次迭代得到各參數的最優估計Popt=(0.99, -0.05,41.27, 9.04, 2.06, 0.92, 26.07)，χ2隨迭代次數的變化如圖7 所示。參數估計偏差ΔP=Popt-Ptrue=(-0.01,-0.05, 1.27, 0.04, 0.06, -0.08, 0.07)，塵埃電荷和軌跡的反演結果如表1 所示。其中，速度的反演結果誤差最大，被低估了7.63%，其原因可能是在獨立參數的估計中算法能獲得精度較高的塵埃位置，但塵埃穿過DTS 的時間估值偏差較大（1.32 μs），導致反演的塵埃速度減小。

表1 塵埃電荷與軌跡的反演結果Table 1 Inversion results of the dust charge and trajectory

圖7 χ2 隨迭代次數的變化Fig.7 Variation of χ2 with iteration steps

參數初值P0與最優估計Popt對應的金屬絲信號分別如圖8(a)和(b)所示。

圖8 不同參數估計對應的金屬絲信號Fig.8 Metal wire signals corresponding to different parameter estimates

圖中藍色表示模擬的實際觀測數據，紅色表示根據P0和Popt反演出的金屬絲信號?？梢钥吹?，經過6 次迭代后，Popt能復現實際的觀測信號，說明時間估計的誤差只影響塵埃速度大小的反演結果，對塵埃速度方向及DTS 金屬絲信號的反演結果沒有影響。

3 反演精度的影響因素分析

分別改變用于仿真的塵埃速度分量vz、入射角度θy、塵埃電荷量Q和探測器電荷噪聲比QNR，分析它們對反演精度的影響。各參數取值如表2 所示，當改變其中1 個參量時，其他3 個參量按表2中粗體數字取值。根據采用的DTS 結構，塵埃穿過DTS 的最大入射角為38.6°，仿真設置θx=0°，θy在0°～30°范圍內變化，表3 列出了與表2 中θy對應的預設塵埃軌跡（進、出DTS 的位置坐標）。

表2 各控制變量的取值Table 2 Values for control variables

表3 不同θy 對應的塵埃軌跡Table 3 Dust trajectories for different values of θy

按表2 依次改變各控制變量數值，針對每一組控制變量（vz、θy、Q、QNR），重復圖5 的反演流程100 次，每次的模擬信號疊加不同的噪聲，會反演得到100 組不同的塵埃電荷量（Q）及軌跡參數（vz、θy、θx）。塵埃入射角度的反演誤差Δθ為反演最優值與真值之差，即

電荷量和速度大小的反演誤差Δa以相對誤差表征，即

其中符號a代表參數Q或v。

計算得到100 組反演結果的反演誤差后，分析這100 組反演誤差的均值及誤差絕對值的均值隨表2 中各控制變量的變化，結果如圖9～圖12 所示。其中，誤差絕對值的均值（圖中藍線）反映反演結果相對 “真”值的整體偏離幅度，誤差均值（圖中粉線）反映反演結果相對“真”值的整體偏離方向（被高估或被低估）。

圖9 反演誤差均值及其絕對值均值隨vz 的變化Fig.9 Variations of mean inversion error and mean absoluteinversion error with vz

3.1 控制變量為vz 的情況

反演誤差隨vz的變化如圖9 所示。

3.2 控制變量為軌跡傾角θy 的情況

反演誤差隨θy的變化如圖10 所示。

圖10 反演誤差均值及其絕對值均值隨θy 的變化Fig.10 Variations of mean inversion error and mean absolute inversion error with θy

3.3 控制變量為塵埃電荷量Q 的情況

反演誤差隨Q的變化如圖11 所示。

圖11 反演誤差均值及其絕對值均值隨Q 的變化Fig.11 Variations of mean inversion error and mean absolute inversion error with Q

從圖11 可以看出，塵埃電荷量的變化對各參數反演誤差沒有顯著影響。保持在0.4%左右，在0.2%左右，和接近0.15°。在反演結果偏離方向上，塵埃電荷量和速度大小整體被低估0.1%左右，軌跡傾角相對于“真”值同樣沒有明顯的偏移方向。

3.4 控制變量為QNR 的情況

反演誤差隨QNR 的變化如圖12 所示。

圖12 反演誤差均值及其絕對值均值隨QNR 的變化Fig.12 Variations of mean inversion error and mean absolute inversion error with QNR

通過上述分析可知，本文所提出的方法在其考慮的參數范圍內，對塵埃電荷量的平均反演精度優于1.5%，速度大小的平均反演精度優于0.8%，軌跡傾角的平均反演精度優于0.6°。

4 結束語

塵埃穿過DTS 會誘導金屬絲產生電信號，通過對觀測信號和仿真信號進行擬合可提取塵埃的電荷量與軌跡信息。本文對DTS 探測塵埃速度矢量進行了數值仿真，結果表明：

1）塵埃電荷量和軌跡的反演誤差隨塵埃速度vz的增大而增大，隨探測器QNR 的增大而減小，基本不受軌跡傾角變化和塵埃帶電量的影響。

2）塵埃電荷量反演精度優于1.5%，速度大小反演精度優于0.8%，速度方向反演精度優于0.6°。

本研究通過DTS 探測獲得了高精度的塵埃電荷量和軌跡信息，并對探測精度影響因素進行了仿真分析，可為DTS 研制與科學探測提供參考。需要說明的是：由于反演過程需要7 個獨立參數，本文算法中通過隨機計算103次來確定χ2的最小值很大可能是局部最小而不是全局最??；本文用于仿真的DTS 結構為較優結構，并未針對其他結構類型（如不同的柵網平面間距、金屬絲間距等）進行研究；本文對塵埃顆粒形狀進行了簡化（以球形代替），且用于計算探測器信號的公式與實際情況間存在一定偏差。后續將重點開展金屬絲誘導電荷信號的理論和試驗研究，依據更符合實際情況的信號修正本文的計算方法，進一步開展反演精度分析。