2010, Vol. 53
2. 中国地震局地震预测研究所, 北京 100036
, MA Yu-Qian1, WANG Huan-Yu1, LU Hong1, ZHANG Xue-Min2, WANG Ping1, SHI Feng1, XU Yan-Bing1, HUANG Jian-Ping2, MENG Xiang-Cheng1, WANG Hui1, ZHAO Xiao-Yun1
2. Institute of Earthquake Science, China Earthquake Administration, Beijing 100036, China
随着人们对地震电磁现象研究的逐渐深入，认识到地震前兆可能以电磁波发射等方式耦合于空间电离层、辐射带等区域.为此，已经建立了多个理论模型对其进行描述.所涉及的领域包括化学、声学以及电磁学等等，尤其是地震对空间的电磁影响，地震学家发现震前数小时至几十小时前有低频电磁信号异常现象，对地震短期预报工作具有非常重要的意义，因此，地震电磁现象越来越多地被人们关注.为了研究地震前兆在空间的电磁表现形式，电磁卫星应运而生，这些卫星上的探测器通常以测量电离层空间的电场、磁场和离子密度为主，并装有高能粒子探测器.至今已有一些关于发现了地震电磁异常现象的报道，1981年8月24日，Inter-cosmos-Bulgaria 1300卫星发现了对应于（经度181°，纬度51.5°）的5.2级地震的极低频（约8Hz）波段出现电磁扰动现象[1]；Parrot等[2]对GEOS-2卫星的统计分析表明，在距离卫星经度20°范围内的极低频ELF电磁信号增强与4.7级以上地震相关系数为0.54；Larkina[3]对INTERKOSMOS-19卫星的统计所得到的电磁异常与地震相关性更高，达到了0.8以上.也有关于地震前数小时有小时时标上的带电粒子涨高事例的报道[4].我国是一个地震多发国家，分析和研究地震先兆与空间电磁现象的关联，亦是提高地震预报准确度的重要手段之一.
法国的DEMETER卫星（Detection of Electro-Magnetic Emissions Transmitted from Earthquake Regions）是一个典型的电磁卫星[5~7]，主要测量空间多个频率上的电场、甚低频磁场、电子和等离子体的密度和温度，以及高能电子的通量和能谱（IDP）等.其科学目标是研究地震、火山及人类活动所引起的电离层现象，研究全球电离层电磁环境，为短期预测地震提供科学依据.该卫星自2004年6月29日发射以来，运行至今，这给我们分析2008年汶川地震的可能关联提供了机会.一年多来，国内已相继有一些研究成果发表[8~10]，但是至今还没有关于高能粒子方面的结果报道，我们从DEMETER卫星数据中心下载了高能粒子探测器IDP（Instrument for Detection of Particles）的相关数据，从了解IDP测量结果入手，分析了高能电子通量与能谱的日变化和短时标变化，发现在地震前后确实存在变化，本文将首次给出这些结果并给以讨论.
2 DEMETER卫星与高能粒子探测器DEMETER卫星总重量130kg，轨道高度发射初期为710km，2006年调整为665km，轨道倾角98°，卫星轨道为准太阳同步轨道，自北向南时（称作降轨），对应于当地时间的白天.自南向北（升轨）时对应于当地时间的夜晚，轨道周期102.86min[5].
DEMETER卫星上载有5类科学仪器：电场探测器（ICE），可测量从DC/ULF（0~15 Hz）到HF（10kHz~3.175 MHz）4个不同频率范围内电场的3个分量；磁场探测器（IMSC），由3个正交磁探测线圈组成，频率测量范围从DC-VLF；两个Langmuir探针（ISL），用于测量电子浓度、温度和卫星电压；离子分析仪（IAP），测量等离子体成分、离子浓度、离子温度及等离子体漂移速度；以及高能粒子探测器（IDP）用于测量高能电子的通量和能谱[5].
高能粒子探测器IDP的结构如图 1a所示[6]，其主体是一个具有一定准直范围的硅探测器，仪器的半张角16°，硅探测器厚度1mm，直径25mm，有效面积490mm2，探测能量在70~2340keV范围内的高能电子，该探测器可有两种工作模式：巡查模式和暴模式，分别具有4s和1s的采样时间精度，并且分别给出128道和256道的能谱测量，能道在测量范围内均匀划分，巡查模式对应能道宽度17.8keV，暴模式对应能道宽度8.9keV[5, 6].
IDP的最大几何因子为1.16cm2·sr.IDP的有效几何因子曲线如图 1b所示，可以看出，从800keV开始IDP的探测效率随着能量升高迅速下降.IDP被安装在DMEMTER卫星的侧面，在卫星飞行过程中，由于地球的自转，对于同一星下点区域而言，会造成升轨时入射窗口朝西，而降轨时朝东的效果.由于升轨和降轨朝向相反，两个轨道段的电子通量的差别能够在一定程度上反映电子在东西经度方向漂移的影响.
DEMETER地面数据中心可提供的数据有4个级别[7, 11, 12].其中，0级为观测原始数据；0'级为快速检查数据，用于测试原始数据，快速判断仪器工作状态和数据质量；1级数据包括物理数据以及由轨道参数、姿态参数和地震事例等构成的工程数据.
汶川地震发生在2008年5月12日14:27:58.7（UT），震级8.0，震中位置在（30.97°N，103.57°E），对应磁壳层参数L=1.225.我们下载了2008年全年的1级（Level-1）IDP观测数据供分析使用.IDP的1级数据文件中提供了两种工作模式下每能道的入射电子通量（counts·cm-2·s-1·sr-1·keV-1），为了分析方便我们把所有数据都归一到探测器接收的每能道计数率（counts·s-1）.两者之间只相差一个常数，即几何因子（1.16cm2·sr）与能道宽度（暴模式8.9keV，普查模式17.8keV）的积.
3 分析与结果 3.1 高能电子通量的东西向漂移我们研究了2008年5月份汶川附近区域的电子通量分布，按升轨和降轨分别作出北半球L：1~1.2（纬度：0°N~28°N），经度70°E~130°E的电子计数率的统计分布，能量划分为90~600 keV，600~1000keV和1000~2342keV三个范围，分别求出两个时标4s时标（对90~600keV能段）和28s时标（对于600~1000keV和1000~2342keV能段）上的平均计数率并进行统计.90~600keV的统计结果如图 2所示.
|
图 2 2008年5月，90~600keV的电子计数分布（L：1~1.2，经度：70°N~130°N） （a）为降轨情况，（b）为升轨情况.图中直方图为电子计数率频数分布，平滑曲线为高斯拟合结果. Fig. 2 The distribution of the 90~600 keV electron in May 2008 (L:1~1.2, Lat.:70°N~130°N). (a) Downward, (b) Upward.The histogram is the frequency distribution of the electron count rate, and the smooth line is the Gassian fit result. |
各个能段的高斯拟合参数详见表 1.
由表 1可知，在L：1~1.2（纬度：0°N~28°N），经度70°E~130°E范围内，在我国汶川及以下偏低纬度地区，因为地处正磁区，电子通量很低.在高于600keV的两个能区的电子通量都不存在升、降轨的差别，在大于1000keV相对偏高的计数率主要来自能量高端进入探测器的质子，很高的标准偏差值说明事例的稀少因而涨落很大.90~600keV能量范围内的电子升轨比降轨的平均计数率高出（1.18±0.04）/s，相对于降轨时的15.6%，但是与拟合的标准偏差值相比，这个差别是很小的，至少在短时标测量的情况下这个差别会被埋没在统计涨落中.
|
表 1 不同能量范围的电子计数率分布高斯拟合参数 Table 1 The parameters of Gaussian for the electron count rate in different energy range |
为了避免地域不同造成的差别，我们研究了同纬度范围内更小经度区间的计数率分布，经度在100°E~110°E的计数率中心值与表 1在误差范围内一致，而拟合的标准偏差减小为1.78±0.13，与平均值的东西差异相比仍然是高的.
了解升轨和降轨情况下电子通量的变化，目的是考察电子径向漂移的情况.原则上，昼夜情况下辐射带的高度因受太阳影响不同，白天会下降，在相同卫星轨道高度下，电子平均计数率低于夜晚，这个因素与地震无关；而由于地球的自转，电子有自西向东方向的经向漂移，此漂移量有可能是地震造成的甚低频（VLF）电磁波的变化与耦合而引起的变化.但是由于汶川地区位于正磁区，在中低纬度上，总的电子通量过低，东西差别很小，这一参量在本工作区间内用来做汶川地震的前兆分析会有困难，至少在短时标上是不可能的，需要比IDP更灵敏的探测器.
3.2 电子能谱与汶川地震的可能关联 3.2.1 按轨道寻找本文研究了2008年5月12日汶川地震前后4~6月三个月内过汶川附近的轨道在不同L值范围内累积的能谱的变化.以过汶川的轨道为中心，在东西10°范围共有54条升轨轨道，全部选出以供分析，如图 3所示，图中五角星标记位置为汶川，5个黑三角为5个VLF发射站的位置，五角星附近的平滑斜线为卫星轨道.
|
图 3 汶川区域轨道选择 五角星为汶川，五角星附近的黑曲线为卫星轨道 Fig. 3 Selected orbits nearby Wenchuan area The pentagram is Wenchuan, the black smooth lines nearby the pentagram are orbits |
如图 4所示，为2008年5月1日、12日和20日所选轨道在L：1.8~1.9（47°N~49°N）的能谱，可以看出5月1日和5月20日的能谱与5月12日的相比在低能端有明显的长高.
|
图 4 地震当天（5月12日（b））与地震前后（5月1日（a），5月20日（c））所选汶川区域轨道在北半球L：1.8~1.9范围内的能谱 Fig. 4 The spectrum of selected orbit range in the northern hemisphere L:1.8~1.9 in May 12 (b) and before and after ((a):May 1, (c):May 20) |
为了与5月12日当天过汶川的轨道比较，以5月12日过汶川的轨道能谱为准，分别求出其他轨道相同区域的能谱与其之差.发现在L：1.8~1.9北半球区域，地震期间与地震前后的能谱差有明显不同.图 5为汶川东西10°轨道，L：1.8~1.9北半球区域能谱差从2008年4月初到6月底的演化，纵向每个颜色条代表一条轨道与5月12日当天过汶川轨道的能谱差.由于能谱差主要位于低能段，所以图中能量范围只给到300keV以下.
|
图 5 汶川附近轨道在北半球L：1.8~1.9范围内的能谱减去5月12日同区域能谱所得到的能谱差 颜色表示能谱差的数值. Fig. 5 The spectrum difference of selected orbit range in the northern hemisphere L:1.8~1.9 The color represents the value of spectrum difference. |
从图 5中可以看出，在L：1.8~1.9范围内，能谱差在100~200keV能量范围内普遍存在一个峰，并且，该峰从5月5日开始基本消失，直到5月20日之后才逐渐恢复，并且，能谱差在这个范围内也基本都为正值，南半球L：1.8~1.9范围内的演化与此类似.此外，本文还研究了在L：1.6~1.7以及1.7~1.8范围内的能谱演化，L：1.6~1.7范围的能谱差在地震前后没有明显的变化；L：1.7~1.8范围内，北半球的能谱差在5月12日前后存在一个下降，而且，绝大多数时候的能谱差都为正值，这表明在这个范围内5月12日前后地震时期的通量最低.
这里的能谱差当中位于100~200keV能量范围的峰来自于VLF发射站的影响[11].
为了考察电子能谱变化的相关量，搜索了震中区域上空2000km范围内的升轨轨道，并提取了超低频400Hz的磁场频谱密度，为反映汶川震中上空的变化，计算了每日25°N~35°N区间范围内的均值，为比较演化趋势，把对应图 5当中能谱差的计数统计与之画在一起（注：能谱差为0或负值的点在对数坐标上未标出）.时间序列的比较如图 6所示.
|
图 6 对应图 5中高能电子计数率的时间演化（a）与2008年4月9日至6月23日400 Hz磁场在25°N~35°N的日均值序列（b） Fig. 6 The count rate of the spectrum difference of Fig. 5a from April 9 to June 23, 2008 and the 400 Hz magnetic field (b) in 25°N~35°N |
从图 6当中可以看出，400 Hz磁场在纬度25°N~35°N间的平均频谱值的时间演化与电子能谱的时间总体演化趋势较为一致，尤其是在5月12日地震前后，5月1日至5月20日，计数率和400Hz平均频谱值以及Kp指数都呈先下降后上升的演化趋势.据此本文认为：地震期间汶川附近的地磁扰动在一定程度上影响到汶川以北L：1.8~1.9区域的高能电子的通量和能谱.
3.2.2 按电子沉降带寻找首先作出了2008年全年（经度在60°~180°）范围内，每1°×1°范围内的电子平均计数率分布（图 7），图中白框所圈的带是由VLF发射站引起的电子沉降带[11].在卫星的轨道高度上，这条沉降带的L范围大约以1.7为中心分布.把白框所圈区域从左至右按1~4顺序编号分为4个区域（每个区域L值从1.6~1.8，经度范围占20°），我们进一步研究了这4个区域按轨道的电子能谱的时间演化.
|
图 7 中国及附近区域的2008年600keV以下电子计数率分布 图中黑色五角星为汶川位置，从汶川向右上方连出的虚线为主震破裂带方向，5个黑三角分别为文献[8]中的5个VLF信号发射站位置.空白地区无数据. Fig. 7 The distribution of electron ( < 600keV) count rate in 2008 (L < 2.5) The pentagram is Wenchuan, the dashed from Wenchuan to top right direction is the direction of the zone of fracture of the principal earthquake, the five triangles are the position of the five VLF station in [8].There is no data in the blank. |
为了更为直观地从数值上反映汶川地震前后的电子分布的变化情况，本文逐天统计了2008年5月份第2、3两个区域按轨道累积的电子能谱，如图 8所示.从8a可以看出，5月12日地震前第2区域（L：1.7±0.1，经度：120°±10°）能谱有连续涨高，在5月7日和5月8日连续两天达到最高.按照90~600keV的总计数率的值如图 9a所示，可明显看出5月7~8日的计数率涨高，并且达到这段时间的最高值（约22000个/s）.从图 8b和9b可以看出，3号区（L：1.7±0.1，经度：140°±10°）从5月6~16日约10天内，电子能谱一直保持在高水平状态，从该能段累积的总计数率上看，这些天电子计数率多数都达到10000个/s以上.1号区域（L：1.7±0.1，经度：100°±10°）与图 5中的能谱演化所处的位置南北紧邻，演化趋势基本一致；4号区（L：1.7±0.1，经度：160°±10°）计数率在震发当天有所提高，但与前后对比并不显著.
|
图 8 VLF北半球电子沉降带上4个区域当中两个区域的能谱在2008年5月份的演化 （a）区域为L：1.7±0.1，经度：120°±10°；（b）区域经度：140°±10°.图中红线标记的是5月12日. Fig. 8 The evolvement of electron spectrum in the two of the four area in the VLF electron precipitation belt in May 2008 (a) Second sub-area:L:1.7±0.1, lognitude:120°±10°; (b) Third sub-area:L:140°±10°.The day marked by red line is May 12. |
|
图 9 分别对应图 8两个区域的2008年5月份的相应能量范围内的电子计数率的逐天演化 （a）L：1.7±0.1，经度：120°±10°，（b）经度：140°±10°.图中红线标记的是5月12日. Fig. 9 The evolvement of electron count rate in the same area with the two sub-areas of Fig. 8 in May 2008 The day marked by red line is May 12. |
2008年太阳处于平静时期，5月份没有大的太阳活动，如M级以上耀斑，大的日冕物质抛射事件，磁暴等[13].表明图 8和9中显示的地震前后的电子通量涨高不是太阳活动的影响.
第2、3两区域电子能谱和计数率异常出现的开始和持续时间的特征都有所不同，第2区的异常出现于5月7日和8日两天，而第三区的异常则从5月6日开始一直持续到5月18日.这一现象是否与地震发生的特点存在着某种关联呢?我们得知地震的破裂面“从震中汶川县开始破裂，并且破裂以3.1km/s的平均速度向北偏东49°方向传播，破裂长度约300km，破裂过程总持续时间近120s，地震的主要能量于前80s内释放”[9]，在图 7中标记了地震破裂面的传播方向，其延长线恰好位于分析的第2区.这一符合完全是偶然的还是预示了某种关联的存在，需要进一步的研究，了解这个方向上的磁场的变化可能有助于得出结论.
3.3 电子计数率的短时标变化 3.3.1 与离子密度短时标变化的关联对应于曾中超等[10]提到的4次离子密度与地震关联的观测异常，我们按照文中指出的时间和轨道，对IDP数据作了考察，发现只有1次存在电子计数率的明显长高，发生的时间完全一致，但是该涨高主要出现在高能段（1000keV以上），而且仅出现在4s的时标上，从图 10中即可看出，这个涨高在28s时标上就变得不明显了，考证一下该事例的能谱，高能段的计数很低，且出现在能段的最高端，造成峰值的统计误差很大，因而不能作为一次与电子通量的相关事例.
|
图 10 IDP电子计数率异常（图（c）为1000~2342keV）与IAP异常图（d）的时间比较，图（d）图来自文献[10]. Fig. 10 Compare of the abnormity of electron count rate in IDP (arrowhead in the left third plot, 1000~2342keV) with the IAP abnormity (arrowhead in the right (c) plot).Right plot comes from the literature [10]. |
本工作未在与文献[10]中的其他3次异常的比较中发现IDP电子计数率的异常变化.
3.3.2 与汶川余震的可能关联本文研究了从2008年5月12日到5月16日期间汶川地区上空电子通量随时间的变化在28s时标上与4级以上的余震发生的时间之间可能存在的关联.这段时间范围内，共有12次大于4级的余震落在了IDP观测时间范围内.
在这几天卫星进入以过汶川轨道为中心，东西经度20°（相比图 3的轨道区域东西各扩展10°）区的每个夜晚的半轨中，我们选出了有余震发生的轨道进行分析，作出每28s时标上的3个能段的电子平均计数率随时间的变化.
对这12次余震期的升轨数据进行统计，共发现5次与地震存在关联的电子计数率超出本底的变化，其中有3次电子计数率升高，两次电子计数率降低.从统计结果看，计数率变化幅度都偏小，多数在3~4σ，并主要发生在高于1 MeV的低统计量区.如果用高于5σ作为标准，都无法判定其存在关联.此外，由于高纬度区的高电子计数及造成的短暂的可用时间段，用半轨数据寻找关联存在较多局限，对于频繁余震的情况更是困难.
4 结论用DEMETER卫星高能电子探测器的观测数据，对在2008年5月份我国地区上空28°N以下地区升轨和降轨时的电子通量做出统计，发现在这个纬度上电子通量很低，在90~600keV能区升轨时的计数率比降轨时高出15.6%，其差别小于各自测量时的1倍标准偏差.说明电子通量的由西向东的经向漂移及昼夜差别相对于其统计涨落是个小量.
分析研究了2008年5月份过汶川地区上空的54条轨道上测量的电子能谱，在磁壳参数L < 2范围内按L值分段累计能谱，并比较于5月12日能谱的不同，发现在L：1.8~1.9区域内低能（90~600keV）区，地震前后的电子能谱，有明显时间演化，对应能谱差的计数率与400 Hz磁场在25°N~35°N间的平均频谱值的时间演化有明显的正相关.
VLF基站所引起的北半球电子沉降带与地震破裂带方向相交处的计数率在震前以及地震期间有显著长高，而且在不同经度区域涨高特征有所不同.
在短时标上，与文献[10]所提到的4次与地震相关异常对应的IDP未发现显著异常.在28s时标上研究了汶川地震后余震与临近汶川轨道上电子计数率变化的关联，未发现有统计显著性的关联.
致谢感谢法国DEMETER卫星数据中心提供数据下载服务.
| [1] | Chmyrev V M, Isaev N V, et al. Observation by space-borne detectors of electric fields and hydromagnetic waves in the ionosphere over an earthquake centre. Phys. Earth Planet. Inter. , 1989, 57: 110-114. DOI:10.1016/0031-9201(89)90220-3 |
| [2] | Parrot M, Lefeuvre F, Corcuff F Y, Godefory P. Observations of emissions at the time of earthquakes in the Kerguelen 1985 Modelling the ionospheric disturbance caused by an explosion on the ground. Ann. Geophys. , 1985, 3: 695-704. |
| [3] | Larkina V I, Migulin W, Molchanov O A, et al. Some statistical results on very low frequency radio wave emissions in the upper ionosphere over earthquake zones. Phys. Earth Planet. Inter. , 1989, 57: 100-109. DOI:10.1016/0031-9201(89)90219-7 |
| [4] | Aleksandrin S Yu, Galper A M, Grishantzeva L A, et al. High-energy charged particle bursts in the near-Earth space as earthquake precursors. Ann. Geophys. , 2003, 21: 597-602. DOI:10.5194/angeo-21-597-2003 |
| [5] | Cussac T, Clair M A, Ultré-Guerard P, et al. The DEMETER microsatellite and ground segment. Planetary and Space Science , 2006, 54: 413-427. DOI:10.1016/j.pss.2005.10.013 |
| [6] | Lasoutte D, Brochot J Y, Carvalho D De, et al. The DEMETER Science Mission Center. Planetary and Space Science , 2006, 54: 428-440. DOI:10.1016/j.pss.2005.10.014 |
| [7] | Sauvaud J A, Moreau T, Maggiolo R, et al. High-energy electron detection onboard DEMETER:The IDP spectrometer, descripsion and first results on the inner belt. Planetary and Space Science , 2006, 54: 502-511. DOI:10.1016/j.pss.2005.10.019 |
| [8] | 何宇飞, 杨冬梅, 陈化然, 等. DEMETER卫星探测到可能与汶川地震有关的地面VLF发射站信号的信噪比变化. 中国科学(D辑) , 2009, 39(4): 403–412. He Y F, Yang D M, Chen H R, et al. DEMETER satellite detected the transformation of the Signal-to-Noise of VLF signal from ground station, which may related with Wenchuan Earthquake. Science in China (Ser.D) (in Chinese) , 2009, 39(4): 403-412. |
| [9] | 2008年四川汶川8级地震的成因分析, 中国地震局地质研究所, http://www.csi.ac.cn/sichuan/beijing04.htm The analysis of the cause of Wenchuan magnitude 8 earthequake in 2008, Institute of Geology, China Earthquake Administration, http://www.csi.ac.cn/sichuan/beijing04.htm |
| [10] | 曾中超, 张蓓, 方广有, 等. 利用DEMETER卫星数据分析汶川地震前的电离层异常.. 地球物理学报 , 2009, 52(1): 11–19. Zeng Z C, Zhang B, Fang G Y, et al. The analysis of ionospheric variations before Wenchuan earthquake with Demeter data. Chinese J. Geophys. (in Chinese) , 2009, 52(1): 11-19. |
| [11] | Sauvaud J A, Maggiolo R, Jacquey C, et al. Radiation belt electron precipitation due to VLF transmitters:Satellite observations. Geophys. Res. Lett , 2008, 35: L09101. |
| [12] | DEMETER data server, http://demeter.cnrs-orleans.fr |
| [13] | News and information about the Sun-Earth environment, http://spaceweather.com |