宇宙塵埃軌跡探測(cè)器信號(hào)提取方法及反演精度影響因素分析

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

（1.北京衛(wèi)星環(huán)境工程研究所； 2.可靠性與環(huán)境工程技術(shù)重點(diǎn)實(shí)驗(yàn)室；3.航天東方紅衛(wèi)星有限公司：北京 100094； 4.中國(guó)航天科技集團(tuán)有限公司，北京 100048）

0 引言

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

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

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

宇宙塵埃大都分布在亞微米到厘米量級(jí)，速度在每秒幾千米到幾十千米不等，DTS 的探測(cè)精度也因塵埃性質(zhì)不同而顯示差異。本文基于典型的DTS結(jié)構(gòu)，針對(duì)塵埃軌跡和電荷量的測(cè)量建立DTS數(shù)值仿真模型，研究從DTS 金屬絲信號(hào)提取塵埃粒子速度矢量和電荷量的方法，重點(diǎn)分析塵埃速度、帶電量和入射方向?qū)λ惴ň鹊挠绊?，有助于探索適用于DTS 探測(cè)的塵埃范圍。

1 DTS 結(jié)構(gòu)及工作原理

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

圖1 DTS 結(jié)構(gòu)示意[5]Fig.1 Structure schematic of DTS[5]

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

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

圖3 DTS 數(shù)值模型Fig.3 Numerical model of DTS

探測(cè)器金屬絲的誘導(dǎo)電荷量與塵埃粒子的入射位置、角度、探測(cè)器尺寸和構(gòu)型有關(guān)。Auer 等[1]給出了當(dāng)帶電量為Q的球形塵埃顆粒穿過DTS時(shí)，金屬絲i的誘導(dǎo)電荷qi的近似表達(dá)式

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

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

圖4 Auer 在塵埃實(shí)驗(yàn)室記錄到的DTS 信號(hào)[1]Fig.4 DTS signals recorded in the dust laboratory by Auer[1]

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

粒子速度的大小為[5]

粒子入射角為[5]

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

2 塵埃參數(shù)提取方法

2.1 計(jì)算流程

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

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

參數(shù)提取的具體步驟為：

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

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

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

式中：χ2為仿真信號(hào)誤差；N為數(shù)據(jù)總個(gè)數(shù)；Dm為第m個(gè)實(shí)測(cè)數(shù)據(jù)；Sm為第m個(gè)仿真數(shù)據(jù)。

重復(fù)上述隨機(jī)選取參數(shù)P即χ2計(jì)算103次，χ2值隨P變化，χ2最小值對(duì)應(yīng)的參數(shù)集為P1,min。

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

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

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

2.2 數(shù)值仿真與結(jié)果分析

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

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

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

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

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

參數(shù)初值P0與最優(yōu)估計(jì)Popt對(duì)應(yīng)的金屬絲信號(hào)分別如圖8(a)和(b)所示。

圖8 不同參數(shù)估計(jì)對(duì)應(yīng)的金屬絲信號(hào)Fig.8 Metal wire signals corresponding to different parameter estimates

圖中藍(lán)色表示模擬的實(shí)際觀測(cè)數(shù)據(jù)，紅色表示根據(jù)P0和Popt反演出的金屬絲信號(hào)。可以看到，經(jīng)過6 次迭代后，Popt能復(fù)現(xiàn)實(shí)際的觀測(cè)信號(hào)，說明時(shí)間估計(jì)的誤差只影響塵埃速度大小的反演結(jié)果，對(duì)塵埃速度方向及DTS 金屬絲信號(hào)的反演結(jié)果沒有影響。

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

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

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

表3 不同θy 對(duì)應(yīng)的塵埃軌跡Table 3 Dust trajectories for different values of θy

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

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

其中符號(hào)a代表參數(shù)Q或v。

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

圖9 反演誤差均值及其絕對(duì)值均值隨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 反演誤差均值及其絕對(duì)值均值隨θy 的變化Fig.10 Variations of mean inversion error and mean absolute inversion error with θy

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

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

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

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

3.4 控制變量為QNR 的情況

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

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

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

4 結(jié)束語

塵埃穿過DTS 會(huì)誘導(dǎo)金屬絲產(chǎn)生電信號(hào)，通過對(duì)觀測(cè)信號(hào)和仿真信號(hào)進(jìn)行擬合可提取塵埃的電荷量與軌跡信息。本文對(duì)DTS 探測(cè)塵埃速度矢量進(jìn)行了數(shù)值仿真，結(jié)果表明：

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

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

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