2. 中国科学院国家天文台, 北京 100101;
3. 中国科学院月球与深空探测重点实验室, 北京 100101
2. National Astronomical Observatories, Chinese Academy of Sciences, Beijing 100101, China;
3. Key Laboratory of Lunar and Deep Space Exploration, Chinese Academy of Sciences, Beijing 100101, China
美国著名工程师卡尔·央斯基在1930年用他自己建造的“旋转木马”首次观测到来自银河系的射电辐射[1],自此展开了对射电天文的广泛研究。天文学家的研究发现,不同天体因为辐射机制的不同产生不同波段的电磁辐射。因此,对宇宙进行全波段观测是一件非常有意义的事。但是,地球大气的存在使很多射电波段的信号特别是频率低于10 MHz的甚低频波段难以到达地面。然而,低频波段的探测对于研究天体的演化与起源有至关重要的作用。为了弥补甚低频波段的观测空白,许多国家或组织发射卫星到太空进行观测,比如太阳长期探测卫星(WIND)和尤利西斯号卫星(Ulysses)。
由于地球以及众多在轨卫星的电磁干扰,空间低频射电观测也面临挑战。因为月球自身可以遮挡来自地球方向的射电干扰,所以月球背面的射电宁静区是十分理想的观测地点。欧洲航天局于2015年1月提出在2025年向月球背面发射一颗载有低频射电频谱仪着陆器的月球背面探测任务(Farside Explorer)计划[2]。我国提出嫦娥四号探月计划,嫦娥四号是世界首颗在月球背面软着陆和巡视探测的航天器。2018年5月21日,嫦娥四号中继星“鹊桥”成功发射,着陆器和巡视器将于年底发射,着陆器上搭载了一台科学载荷低频射电频谱仪,是人类在月球背面首次实现对宇宙和太阳的低频段(0.1~40 MHz)观测。
根据着陆器整器在微波暗室的电磁兼容性试验结果,虽然来自外界的信号干扰可以很好地被屏蔽,但是着陆器本身搭载的电子设备以及其它的载荷对低频观测存在电磁干扰,由于无法对平台进行工程上的调整,所以只能通过数据处理的方式消除平台的干扰。目前降低信号的噪声处理方法主要有3种:谱减法、维纳滤波及自适应滤波方法。本文结合低频射电频谱仪的科学目标,对比、分析信号增强的效果,得到最优化的方法,为低频射电频谱仪在轨探测任务的数据处理提供依据。
1 低频射电频谱仪的工作原理以及噪声来源低频射电频谱仪系统由4根接收天线(3根5 m的长天线A,B,C和1根20 cm的短天线D)、前置放大器以及电子学单元组成。3根长天线A,B,C相互垂直[3-4] (如图 1),接收空间电磁信号在三个方向上的分量,根据电磁波的传播理论计算接收信号的强度以及方向。短天线D与长天线B,C处于同一个平面,根据电磁波的传播理论以及不同长度天线的增益不同,A,B,C长天线可以同时接收来自宇宙的远场信号和来自着陆器本身的近场信号,而天线D因天线短且离着陆器近,主要接收来自近场的着陆器噪声信号[3]。
嫦娥四号搭载的低频射电频谱仪的工作频率包括低频(100 KHz~2 MHz)和高频(1~40 MHz)两个频段,频谱仪的主要技术参数如表 1。
项目 | 设计参数 | 实测结果 |
工作频率 | 100 KHz~40 MHz | 100 KHz~40 MHz |
接收机灵敏度(nV/ |
6 | 8.1 |
频率分辨率 | 低频:5 KHz | 4.88 KHz |
高频:100 KHz | 97.66 KHz | |
最大数据率/bps | 4.5 M | 4.34 M |
功耗/W | 24 | 24 |
在低频射电频谱仪工作的同时, 着陆器上其他部分载荷以及一些电子模块也处于运行状态,比如电源控制器、数传调制器、发射段RE102、数传控制电源测量单元(Source Measure Unit, SMU)等。图 2给出了月球表面的电磁波通量密度,天线附近区域的电场强度和接收机输出的电压值根据(1)式计算:
$ {U_{{\rm{out}}}} = \frac{{{G_{{\rm{preamp}}}}{G_{{\rm{Receiver}}}}I{L_{{\rm{coax}}}}}}{{1 + {Z_\alpha }/{Z_{{\rm{preamp}}}}}}\overrightarrow {{h_{\rm{e}}}} \vec E, $ | (1) |
其中,Uout为接收机的输出电压值;Gpreamp和GReceiver分别为前置放大器增益和接收机增益;ILcoax为高频电缆插损;Zα和Zpreamp分别为天线阻抗和前置放大器阻抗;
(1) 式中各项参数根据低频射电频谱仪生产方提供的定标报告获取,结合各项参数以及(1)式将接收机测得的电压数据反演出空间电场以及通量密度。图 2给出了各天体到达月球表面的通量密度,图中的着陆器噪声是根据中国空间技术研究院微波暗室中对嫦娥四号着陆器整器进行的电磁兼容性试验测得工况三(着陆器电源控制器主份、数管管理单元主份、数传调制器主份、载荷电控箱主份、低频射电频谱仪等依次加电连续工作)天线A接收的高频段数据反演得到的。从图 2可以发现,着陆器的噪声通量密度大于普通的太阳爆发产生的通量密度,与太阳爆发峰值通量密度相当。所以不对接收到的信号进行降噪处理,难以反演出太阳爆发的真实强度。
2 低频射电频谱仪有效信号提取方法由于低频射电频谱仪的3根主要天线A,B,C收到的目标信号都包含来自着陆器自身的噪声信号,所以在对信号进行提取之前需要通过降噪算法抑制噪声并最大限度地分离出有用的信号。用x(n)表示观测到的信号,s(n)表示有用信号,v(n)表示噪声信号。它们之间的关系表示为
$ x\left( n \right) = s\left( n \right) + v\left( n \right). $ | (2) |
在本文的模拟仿真中,将一组强度与太阳爆发强度接近的模拟信号作为有用信号s(n),用电磁兼容性试验测得的工况三的着陆器噪声作为噪声信号v(n),通过谱减法、维纳滤波及自适应滤波的方法进行仿真处理。
2.1 谱减法谱减法是信号增强的一种常见方法,在频域上用信号加噪声的幅值减去噪声的幅值得到有用信号的幅值。文[5]于1979年改进了该方法,称为改进的谱减法。改进的谱减法在谱减的过程中使用噪声的平均功率,有效地减小剩余的“音乐噪声”。改进的谱减法用公式表示为
$ {\left| {S\left( n \right)} \right|^\gamma } = \left\{ {\begin{array}{*{20}{c}} \begin{array}{l} {\left| {X\left( n \right)} \right|^\gamma } - aD\left( n \right);\\ bD\left( n \right); \end{array}&\begin{array}{l} {\left| {X\left( n \right)} \right|^\gamma } \ge aD\left( n \right)\\ {\left| {X\left( n \right)} \right|^\gamma } < aD\left( n \right) \end{array} \end{array}} \right., $ | (3) |
其中,D(n)为噪声的平均幅值,
$ D\left( n \right) = \frac{1}{N}\sum\limits_{i = 1}^N {{{\left| {{D_i}\left( n \right)} \right|}^\gamma }} , $ | (4) |
(3) 式、(4)式中,γ为1时,幅值相减;γ为2时,功率相减。a(a ≥ 1)为过减因子,b(b > 0)为补偿因子,通过调节两个参数减少周期性的误差。结合低频射电频谱仪的结构特点以及对科学数据的要求,用长天线A,B,C收到信号的平均幅值减去短天线D收到信号的平均幅值获取信号的幅值。
由于天线A,B,C与天线D的长度及材质不一样导致其增益不同,所以在进行谱减时需要补充不同增益带来的信号强度差异。具体来说,先取一小段时间内的信号计算它们的对消系数,然后在后续的处理过程中将天线D收到信号的平均幅值乘以对消系数获取噪声幅值。结合低频射电频谱仪的结构特点,运用改进的谱减法进行降噪处理,以天线A和天线D为例,用公式表示为
$ \left| {S\left( f \right)} \right| = \left\{ {\begin{array}{*{20}{c}} \begin{array}{l} \left| {{V_{\rm{A}}}\left( f \right)} \right| - \left| {P\left( f \right){V_{\rm{D}}}\left( f \right)} \right|;\\ b\left| {P\left( f \right){V_{\rm{D}}}\left( f \right)} \right|; \end{array}&\begin{array}{l} \left| {{V_{\rm{A}}}\left( f \right)} \right| \ge \left| {P\left( f \right){V_{\rm{D}}}\left( f \right)} \right|\\ \left| {{V_{\rm{A}}}\left( f \right)} \right| < \left| {P\left( f \right){V_{\rm{D}}}\left( f \right)} \right| \end{array} \end{array}} \right., $ | (5) |
其中,VA (f)为长天线A收到的信号进行傅里叶变换后在相应频点的平均值;VD (f)为短天线D收到的信号在相应频点的平均值;γ取1;P (f)取太阳宁静时一小段时间内4个接收机接收到的信号分别算出:
$ {P_{\rm{A}}}\left( f \right) = \frac{{{V_{\rm{A}}}\left( f \right)}}{{{V_{\rm{D}}}\left( f \right)}}. $ | (6) |
处理流程如图 3。
下面3种数据处理方法以天线A与天线D接收到的数据进行仿真为例,谱减法仿真结果如图 4。
图 4展示了快速傅里叶变换后的频域信号,横轴表示频率(单位为MHz),纵轴为幅度(单位为dBuV)。图 4中黑色曲线是带有噪声的信号,蓝色曲线是经过改进谱减法处理后的信号,红色曲线是有用信号。从图 4可以看出,该方法在整个波段对噪声进行了有效抑制,在10 MHz以下频段内噪声幅度下降了20 dB左右,在其它频段噪声降低10 dB左右,且处理后的信号波形基本能还原原本的有用信号,在背景噪声较低的频段(20~40 MHz)处理效果更好。
2.2 维纳滤波从(2)式可知,进行信号处理的目标是得到不含噪声的有用信号s(n)。但是在实际信号处理过程中,求得的有用信号并不完全等于s(n),只是s(n)的逼近值或估计值,因此对信号的处理可以看作是对s(n)的估计,处理信号就是找到一个最佳的估计器。维纳滤波是一个寻找最佳估计器来估计信号的方法[6-7]。
图 5为维纳滤波的原理示意图,用公式表示为
$ y\left( n \right) = x\left( n \right) * h\left( n \right), $ | (7) |
其中,h(n)为维纳滤波器的滤波参数;y(n)为对有用信号s(n)的估计值。
根据最小误差准则使ε=E[{s(n)-y(n)}2]达到最小可以算出最佳滤波参数,
$ h\left( n \right) = \frac{{E\left[ {s\left( n \right) * x\left( n \right)} \right]}}{{E\left[ {{x^2}\left( n \right)} \right]}}, $ | (8) |
将上式进行傅里叶变换可以导出,
$ H\left( k \right) = \frac{{{P_{{\rm{sx}}}}\left( k \right)}}{{{P_{\rm{x}}}\left( k \right)}}, $ | (9) |
其中,Psx(k)为s(n)与x(n)的互功率谱密度;Px(k)为x(n)的功率谱密度。由于有用信号s(n)与噪声信号v(n)互不相关,即有Psv(k)=0,则可以得到:
$ {P_{{\rm{sx}}}}\left( k \right) = {P_{\rm{s}}}\left( k \right), $ | (10) |
$ {P_{\rm{x}}}\left( k \right) = {P_{\rm{s}}}\left( k \right) + {P_{\rm{d}}}\left( k \right). $ | (11) |
将(10)式和(11)式代入(9)式可以得到:
$ H\left( k \right) = 1 - \frac{1}{{\gamma \left( k \right)}}. $ | (12) |
(12) 式中,γ(k)为傅里叶变换后对应频点处带噪信号的功率谱与噪声功率谱的比值[8]。根据天线D和天线A接收到的数据计算出γ(k),从而计算出最佳的滤波参数。图 6为仿真结果。
从图 6可以看出,维纳滤波能够降低噪声,但是相对于谱减法的降噪效果要差一些,基本能还原有用信号的波形,但是处理后的信号在背景噪声比较大的波段(0~20 MHz)仍难以分辨有用信号。
2.3 自适应滤波不同于维纳滤波,自适应滤波器根据收到的信号不断地自动调节自身的参数,最大限度地降低噪声信号[9]。自适应滤波已经广泛应用在系统识别、信道均衡、信号增强以及信号预测方面。本文主要运用信号增强的功能,结构示意图如图 7。自适应滤波算法至少需要两套接收设备,其中一套用来接收带有噪声的信号,一套用来接收噪声,然后通过自适应算法求得有用信号,而低频射电频谱仪正好有这样两套接收设备。
图 7为一个自适应降噪滤波器的结构示意图。其中,x(n)指接收机收到的信号(包括有用信号s(n)以及噪声信号v(n));v1(n)为另一个接收机收到的噪声信号;e(n)为误差信号。它们之间的关系表示为
$ e\left( n \right) = x\left( n \right) - \sum\limits_{i = 1}^N {\left[ {w\left( i \right){v_1}\left( {n - i} \right)} \right]} , $ | (13) |
其中,w(n)为滤波器的加权系数,该系数通过最小均方误差准则或者最小二乘法评估E[e2(n)]来调整,使滤波器工作在最佳状态[10-12]。下面利用最小二乘法进行仿真。
从图 8可以发现,利用自适应滤波方法仿真的结果比较差,降噪效果不明显,对信号的还原度不高。在噪声信号中一些幅度比较大的频点处(如1.73 MHz,2.32 MHz,2.9 MHz等),经过自适应滤波处理后噪声信号幅度没有下降。
3 低频射电频谱仪背景噪声消除结果分析前文通过3种不同的数据处理方法对着陆器噪声进行了仿真处理,其中自适应滤波方法是在时域上的处理方法,改进的谱减法和维纳滤波则是在频域上的处理方法。表 2给出这3种信号增强方法对信号信噪比提升效果的对比。
数据处理方法 | 改进谱减法 | 维纳滤波 | 自适应滤波 |
处理前信噪比/dB | -27.64 | -27.64 | -27.64 |
处理后信噪比/dB | 1.99 | -4.64 | -25.35 |
信噪比提升/dB | 29.63 | 23.00 | 2.29 |
从表 2可以发现,改进谱减法和维纳滤波的处理结果明显好于自适应滤波,且其中改进谱减法的信噪比提升略大于维纳滤波。自适应滤波降噪能力受输入端的噪声信号v1(n)与带噪信号中的噪声信号相关度影响很大。在低频射电频谱仪设计中用于接收信号的长天线A,B,C和用于接收着陆器噪声信号的短天线D的参数不同,仿真信号中有用信号的信噪比非常低等因素导致自适应滤波降噪的效果变差。长天线A,B,C与短天线D收到噪声信号的不同可以通过评估短天线与长天线在不同频点处的增益差异来抵消,所以在频域上处理的两种算法在降噪方面取得比较好的效果。
为了验证低频射电频谱仪的性能以及数据处理方法,在中国科学院国家天文台密云观测站食堂楼顶安装低频射电频谱仪及接收天线进行外场试验,并且在距离接收天线100 m左右的宿舍楼里用任意波形发生器连接天线发射模拟信号。着陆器噪声信号采用电磁兼容性工况三的试验结果,由一个放在接收天线下方1 m长的小天线连接波形发生器发出。图 9为试验天线A接收到的一组数据利用改进的谱减法处理后的结果,其中发射的有用信号为1.8 MHz点频信号。
图 9截取了1.7~1.9 MHz的信号,从图中可以看出,经过改进谱减法处理后,着陆器的噪声被有效地抑制,幅度下降10 dB,而有用信号的幅度并没有降低,信噪比从-25.23 dB提升到-2.60 dB。通过改进谱减法处理后的信号信噪比提升了-22.63 dB,说明改进的谱减法在外场实验中对着陆器背景噪声信号的降低也是有效的。
4 总结在月球背面进行极低频(0.1~40 MHz)探测非常有科学意义[13],同时由于该频段信号微弱,信噪比低,所以提取有用信号十分困难。本文利用电磁兼容性测试的着陆器噪声以及嫦娥四号搭载的低频射电频谱仪对月球背面接收到的目标信号强度进行了仿真实验。利用改进的谱减法、维纳滤波以及自适应滤波这3种常见的信号增强算法对仿真信号进行提取,其中改进的谱减法对信号的提取效果显著,经处理后信噪比有明显提升,可以直观地从处理后的波形中看到目标信号波形。改进的谱减法对在密云进行外场试验的数据处理也有很好的效果。因此,将改进的谱减法作为低频射电频谱仪主要的预处理方法,实现对太阳爆发较为准确的观测。低频射电频谱仪在月球背面进行探测的过程中,温度等环境的变化可以引起对消系数的改变,所以在实际探测过程中需要结合环境情况,每隔一段时间计算一次对消系数。
[1] | JANSKY K G. Directional studies of atmospherics at high frequencies[C]//Proceedings of the Institute of Radio Engineers. 1932, 20: 1920-1932. |
[2] | MIMOUN D, WIECZOREK M A, ALKALAI L, et al. Farside explorer:unique science from a mission to the farside of the moon[J]. Experimental Astronomy, 2012, 33(2/3): 529–585. |
[3] | 纪奕才, 赵博, 方广有, 等. 在月球背面进行低频射电天文观测的关键技术研究[J]. 深空探测学报, 2017, 4(2): 150–157 |
[4] | 施硕彪, 董亮, 高冠男, 等. 米波太阳射电频谱仪的科学目标和技术方案[J]. 天文研究与技术——国家天文台台刊, 2011, 8(3): 229–235 |
[5] | BOLL S. Suppression of acoustic noise in speech using spectral subtraction[J]. IEEE Transactions on Acoustics Speech & Signal Processing, 1979, 27(2): 113–120. |
[6] | BROWN R G, HWANG P Y C. Introduction to random signals and applied Kalman filtering:with MATLAB exercises and solutions[M]. Hoboken: John Wiley & Sons, 1997. |
[7] | WIENER N. Extrapolation, interpolation, and smoothing of stationary time series:with engineering applications[M]. Hoboken: John Wiley & Sons, 1949. |
[8] | 宋知用. MATLAB语音信号分析与合成[M]. 北京: 北京航天航空大学出版社, 2017. |
[9] | 耿妍, 张端金. 自适应滤波算法综述[J]. 信息与电子工程, 2008, 6(4): 315–320 DOI: 10.3969/j.issn.1672-2892.2008.04.018 |
[10] | DINIZ P S R. Adaptive filtering: algorithms and practical implementation[M]. 4th ed.[S.l.]: Kluwer Academic Publishers, 2008. |
[11] | 沈福民. 自适应信号处理[M]. 西安: 西安电子科技大学出版社, 2001. |
[12] | 王树棠, 王蕴红. 自适应滤波器阶数与稳态均方误差[J]. 电路与系统学报, 1997(1): 65–67 |
[13] | 梅丽, 苏彦, 周建锋. 极低频射电天文观测现状与未来发展[J]. 天文研究与技术, 2018, 15(2): 127–149 |