到2009年为止，在国际陨石坑资料库中地球上已发现的陨石撞击坑数量为176个，而中国境内的撞击坑在其中为零^[1].我国从20世纪70年代起陆续在地质、地球化学领域开展了撞击构造的研究工作，至今已提出了太湖、沈阳东南郊区、上河湾、多伦、香港、海南白沙、河北矾山和罗圈里等8个撞击构造^[2~13]，但均未被国际同行确认.国际学术界将陨石撞击坑判别的诊断性标准明确为撞击坑岩石造岩矿物石英和长石中发育击变面状页理（planar deformation features，PDFs）.中国陨石坑研究未被国际同行确认的主要原因是未能获得证实陨石撞击成因的关键证据，值得注意的是，陈鸣2007年在岫岩陨石坑（也称罗圈里撞击坑）发现了具陨石撞击起源的诊断性标志的击变面状页理（PDFs），这在我国陨石撞击坑研究中具有重要意义^[14].此外，国际上常采用钻探和地球物理方法揭示撞击构造深部及整体特征，比如美国最大的Chesapeake湾撞击构造就是通过钻探和地震勘探手段而得以发现和深入研究的.我国撞击构造的研究除岫岩坑外尚未见到钻探和地球物理方法的具体应用，致使研究手段比较单一，研究深度不够，因而对撞击构造不能取得整体的认识.借鉴国外的研究方法，在进一步加强地质、地球化学研究的基础上，运用钻探和地球物理方法，特别是地震勘探方法，对我国的撞击构造进行研究，将会使我国在该领域的研究水平提升到一个新高度^[15].

岫岩陨石坑位于我国辽宁省岫岩县，坐落在辽东半岛北部低山丘陵地区，为一简单碗形坑.赵树明^[12]和覃功炯^[11]分别在1992和1998年提出该坑可能为陨石撞击起源.陈鸣2000年对该坑开展了岩石和矿物冲击变质的调查，发现了含“击变面状页理”（PDFs）的石英颗粒，该发现为该坑的成因研究提供了重要线索^[14].鉴于对撞击构造的地球物理研究在我国还是空白，2008年在国家自然科学基金项目“罗圈里撞击坑高精度地震探测研究”的资助下，中国地震局地球物理勘探中心首次在中国对陨星撞击构造进行了二维地震反射和三维地震折射探测研究.此外，陈鸣于2009年在该坑中心位置实施了地质钻探，钻孔岩芯提供了坑内中心位置处地质和地层结构资料，为岩石矿物冲击变质研究提供了关键样品^[16].至此对岫岩陨石坑先后进行了地质、地球化学、钻探和地震勘探方法的综合性的研究，是中国进行综合研究的第一个陨石坑，这将为中国陨石坑后续的研究奠定基础.

本文主要从地震波衰减的角度来研究该坑的三维Q值结构.Q值与岩石的物理状态（如压力、温度、流体饱和度、流体的黏滞系数等）密切相关，是反映地震波衰减的一个重要参数.衰减的微观机制相当复杂，但在宏观上，可主要归结为“内摩擦”引起的固有衰减、介质的非均匀性所引起的散射衰减和层间耦合机制所引起的力学结构效应三大类.这三方面的因素所引起的衰减总效应在介质的品质因数Q值上得以体现.

岫岩陨石坑（图 1）形成在下元古界变质岩系中，该岩系由浅粒岩、角闪岩、大理岩、片麻岩和变质玄武岩等组成，坑唇顶部和两侧山坡普遍被厚度达数米的风化土壤、残积和坡积物覆盖，局部有零星基岩（主要为浅粒岩）出露，坑底被第四系湖泊相沉积覆盖.陨石坑形态保存良好，平面形状近似五边形到似圆形，形态呈碗状，直径约为1800 m.现在的坑底到坑唇山峰最大高程差约为200 m，平均150 m，坑底中央区海拔高度约140 m.

图 1

Fig. 1

图 1 罗圈里撞击坑地质图(据覃功炯等, 2001)[11] 1.第四纪冲击物 Quaternary Period alluvium; .河流水塘 River and pond; 放射状裂隙 Fractures; .坑外裂隙 Fractures beyond crater; 5.角砾岩块 Breccia block; 6.地层产状 Stratum occurrence; 7.撞击溶岩 Impact magma; 8.铁陨石 Iron meteorite; 9.坑边及震裂片理 Schistosity; 10.震裂维 Shatter cones; 11.下元古界辽河群 Lower Proterozoic; 12 片理产状 Schistosity occurrence. Fig. 1 Luoquanli crater geological map (according to Qin Gong-Jiong, 2001)[11]

根据本项目研究目的和所要达到的科学目标，针对撞击坑的客观情况并经实地踏勘后拟定设计的方阵为5 km²，边长长度为2.24 km，布设观测点200个，撞击坑内观测点距130 m左右，撞击坑外观测点距250~300 m左右，沿方阵布置9个炮点，采用单井或多井组合爆破激发方式，单炮药量为10~20 kg，观测距离最大约为3 km，大多在2 km范围，见图 2.

图 2

Fig. 2

图 2 台站和炮点位置分布图 圆圈表示台站, 五星表示炮点. Fig. 2 Station and shot point locations Circles stand for station, five-stars stand for shot.

地震波波形提供了地壳和上地幔内部物质的成分和物态等信息，利用它不仅可以确定速度和界面结构，而且可以确定Q值结构，该结构对破碎带反应敏感，在破碎带内，对地震波能量吸收和散射作用较强，大多表现为低Q值^{[17, 18]}.

目前，利用地震体波获得Q值结构主要有三种方法：频谱比法，Atkinson方法和Q值层析成像方法^[19~21].前两种方法侧重研究Q与频率的关系，而把与空间位置的关系置于次要地位.第三种方法主要研究Q与空间位置的关系，而把与频率的关系置于次要地位，这就避免了同时研究Q值与频率和空间位置关系的难度.在第三种方法中，从利用地震波的信息角度而言，有主要利用质心频率（centroid frequency）的^[22]，有主要利用振幅的^{[23, 24]}，有利用波形的^{[17, 25, 26]}，从利用信息量的角度而言，波形所含信息量最多.

本文利用人工地震P波波形和基于吸收特征时间t^*的三维Q值层析成像方法，首先求取吸收特征时间t^*，然后反演三维Q值结构.

首先根据地震折射P波走时选择波形，使其尽量少包含其他震相的波形和干扰，使之得到主要反映该震相的振幅谱，称其为实测振幅谱.然后，根据爆破地震的P波远场位移谱与t^*的关系式（1），计算理论振幅谱并与实测振幅谱进行最小二乘拟合，得到观测吸收特征时间t_obs^*.最后利用t_obs^*反演三维Q值结构.

对人工地震，利用人工爆炸源的震源谱表达式^[27]得到的P波远场位移谱为：

(1)

式中，f为频率，Ω₀为零频振幅，f_c为拐角频率，t^*由式（2）定义.

t^*通常是频率和空间位置的函数.实验室以及实际地震测量都表明^[28~30]：对于固体，Q值在地震观测的短周期频段内基本上与频率无关.在忽略地震波频散并采用定常Q模式（即Q与频率无关）之后，三维介质中t^*的定义为：

(2)

式中，Q为品质因数，V为波速，r为射线路径.

Q的定义^[31]为1/Q=ΔE/2πE，其中，ΔE为地震波在一个周期内或一个波长内的运动过程中波动能量的衰减，E为波动能量的最大值.Q值反映了地震波能量衰减的大小，其值越大，能量衰减越小，反之能量衰减越大.

三维Q值模型离散网格节点上Q值的扰动经射线路径而引起地震波t^*的扰动为^{[17, 32]}：

(3)

Q值的层析成像问题归解为迭代求解下式：

(4)

式中，A是由δt^*/δQ组成的Frechet矩阵，Δm是Q值残差构成的向量，Δd是观测的t_obs^*与计算的t_cal^*构成的残差向量.

层析成像过程为：给定初始Q值模型和初始阻尼值λ进行迭代，在迭代过程中自动选择合适的阻尼因子，为了使Q值的修正量Δm更合理，在每次迭代时重新计算新的Frechet偏导矩阵，如此不断迭代，直到t^*残差小于预先给定的可接收的值（本文为0.003）为止.

初始Q值模型的确定基于以下几方面：①参考相关地区浅层工程物探和石油物探方面得到的地表层的Q值；②参考实验室确定的相关岩性的Q值；③参考该陨石坑的三维速度结构；④根据上述三方面的资料设计几个初始模型，通过试验从中选择一个比较好的.

初始阻尼因子的确定主要基于经验和试验，迭代过程中根据目标函数的值自动修改阻尼因子，若值变大就增大阻尼，否则就减小阻尼.

图 2为岫岩陨石坑三维地震折射观测系统，圆圈为接收台站的位置，约200台三分量便携式地震仪^[33]、五角形为炮点位置，共9个炮点.台站间距100~200 m.共获得约1800个地震记录，每个地震记录包含三个分量，两个水平分量和一个垂直分量，水平分量记录的是S波，垂直分量记录的是P波，从中我们精选了约1700个P波震相可靠的地震记录，用于反演该陨石坑的三维Q值结构.图 3给出了SP2的垂直分量记录.

图 3

Fig. 3

图 3 SP2的垂直分量记录 Fig. 3 Vertical component records of SP2

资料处理分四步：首先确定P波震相及做频谱分析用的时窗段，为了提高频谱的分辨率并减少漏能，根据前人的工作和我们的试验，频谱的窗长取1.5~2 s；之后对其作频谱分析，求出振幅谱，称其为实测振幅谱，并做自由表面和仪器频率特性校正；然后根据式（1）计算理论振幅谱，与校正后的实测振幅谱用非线性最小二乘法进行拟合得到观测吸收特征时间t_obs^*；最后反演三维Q值结构.

我们用SeisProc软件包^[34]读取了1700多个P波走时并按1.5~2s的窗函数长度，确定作频谱分析用的时窗段.然后用SpectraFit频谱分析软件对这些地震记录作频谱分析，求出实测振幅谱.地震记录的最高频率为100 Hz，最低频率为0.5 Hz.为减少谱的“泄漏”，采用了谱边辩比较小的汉宁窗，其窗函数长度尽量是信号周期的整数倍，以避免人为产生的附加频率成分.对求出的实测振幅谱进行了自由表面和仪器频率特性的校正.利用式（1）计算基于点源爆破模型的远场位移谱，该谱包括了固有衰减、散射衰减、几何扩散、源幅射图形等因素.利用SpectraFit软件将计算的理论振幅谱与校正后的实测振幅谱利用非线性最小二乘方法进行迭代拟合.因一次爆破只有一个拐角频率，故对同一次爆破的所有台站取统一的f_c，这样将进一步减少t^*的不确定性，提高其精度.我们计算的t^*，其误差范围为0.0004~0.0018 s，大部分t^*的误差为0.001s左右.图 4给出了部分实测振幅谱和理论振幅谱的拟合图.

图 4

Fig. 4

图 4 振幅谱拟合图 粗线表示理论振幅谱, 细线表示观测振幅谱.图上面的数0.24 km等为炮检距. Fig. 4 Amplitude spectrum titting Thick lines stand for theoretical amplitude spectra, thin lines represent observed amplitude spectra.

我们对Zelt92正演程序进行了改进，使其既可以计算地震波的理论走时t，又可以计算地震波的吸收特征时间t^*.在利用SpectraFit软件获得t_obs^*值后，利用改进的Zelt92正演和反演程序进行三维Q值反演.其原理跟走时反演类似，走时反演是给定初始速度模型通过正演程序计算理论走时，再通过反演程序用理论走时拟合观测走时.Q值反演是在已知速度模型（在Q值反演过程中不变）的条件下，给定初始Q值模型（图 6a给出了三维Q值初始模型沿不同深度的水平切面图），通过改进的Zelt92正演程序由式（2）经射线追踪计算理论吸收特征时间t_cal^*，再通过反演程序用t_cal^*与t_obs^*进行拟合（图 5给出了t^*曲线拟合图），最后得到三维Q值模型.图 6b给出了三维Q值结果模型沿不同深度的水平切面. 图 7给出了9炮的射线分布.

图 6

Fig. 6

图 6 三维Q值初始模型(a)及最终模型(b)沿不同深度的水平切面 Fig. 6 Horizontal sections of the initial (a) and final (b) Q-value models at different depths

图 5

Fig. 5

图 5 9炮的t*曲线拟合图 十字表示理论的tcal*, 圆圈表示观测的tobs*. Fig. 5 The curve fittings of t* shots The crosses stand for theoretical tcal*, and the circles indicate the observed tobs*.

图 7

Fig. 7

图 7 9炮人工爆破沿不同视角的射线分布图 Fig. 7 Ray distribution at the different perspectives

从图 6a可知，初始Q值模型中，深度从0.3~0.7 km为均匀层，但从图 6b的最终Q值模型可知，深度从0.3~0.7 km仍有不同程度的低Q值.这表明在陨石坑内的破碎角砾岩层（深度从0.1~0.3 km）下面仍然有不同程度的破碎层，只是其破碎程度较破碎角砾岩层的程度低了许多.

岫岩陨石坑地震数据处理结果表明：陨石坑的三维Q值结构中，不同深度的水平切片上的Q值与该坑地质构造和岩性结构有一定的相关性.陈鸣2009年在该坑中心位置实施了地质钻探，钻孔深度为307 m，钻孔揭露了坑内上部107 m厚的湖泊相沉积物以及下部厚度为188 m的角砾岩堆积透镜体，坑内填充的角砾岩是不同程度冲击变质岩石的混合堆积物^[16].这一钻探结果与图 6b所示的三维Q值结构有一定的相关性：坑内从地表到100 m深处的Q值层对应坑内上部107 m厚的湖泊相沉积物，该沉积物在地下100 m深处的平面分布面积大约是以坑中心为圆心500 m为半径所构成的圆.深度100~300 m的Q值层对应坑内下部厚度为188 m的角砾岩堆积透镜体.从图 6b中的深度为100~300 m的3个小图可知，深度为200 m处的低Q值面积最大，半径约为700 m，而深度300 m处的低Q值面积半径减小为300 m左右.深度300~700 m的Q值层对应角砾岩堆积透镜体下面的下元古界变质岩系的、程度较低的撞击破碎层.从深度400 m处的Q值分布看，中心附近有一较高的Q值突起，可能是撞击坑内的中央峰.

从图 7的沿两个不同方向的射线分布来看，岫岩坑地表下面约250 m、500 m和750 m深处为3个波阻抗分界面.250 m深处的波阻抗界面与2009年陈鸣地质钻孔的260 m深处的陨石冲击分界面相对应：钻孔204~260 m深度区间由弱冲击变质浅粒岩、片麻岩和少量角闪岩混合碎屑组成；而钻孔260~295 m深度区间含大量的石质角砾岩和含熔体多相角砾岩，该层位岩石显示强烈的冲击变质特征并含有大量岩石熔体玻璃碎片^[16].

由于该陨石坑内沉积物的速度较其他一些地方的地表低速带的速度偏低，因而对地震波的衰减更大，故所得Q值较其他一些地方的地表低速带的Q值偏低.

本研究在700 m深度范围内能够得到Q值随深度的变化，与三方面因素有关.其一与高分辨观测系统有关，其二与高质量的观测数据有关，其三与所用的反演方法有关.

就目前三维Q值反演方法而言，基于吸收特征时间t^*的反演方法得到的结果从介质非均匀角度看应该说更接近实际一些，该方法能得到与三维非均匀速度结构有一定对应关系的结果，这是其他方法不容易得到的，尽管该方法所得的Q值仍有一定的波动范围.

从吸收特征时间t^*的定义以及其正反演过程可以看出，其原理跟走时反演类似.走时反演通过迭代过程不断修改速度模型用理论走时拟合观测走时，最终得到非均匀速度结构.Q值反演是在已知速度模型的条件下（在Q值反演过程中保持不变），通过迭代过程不断修改Q值模型用理论t_cal^*拟合观测t_obs^*，最终得到非均匀Q值结构.

若三维走时反演方法得到的速度结构是三维非均匀的，那么，从物理角度而言，反演得到的Q值结构也应该是三维非均匀的，而与Q值的大小以及所处范围的大小关系不大.

基于人工爆炸源的P波远场位移谱（式（1））的非线性最小二乘拟合方法和人机交互软件，可方便快捷地计算地震体波（P波和S波）资料的振幅谱，通过拟合理论振幅谱和观测振幅谱可得到表征地震波衰减的吸收特征时间t^*.

基于吸收特征时间t^*的三维Q值层析成像方法，对确定陨石坑的衰减结构是有效且可行的.

衷心感谢审稿专家和编辑部各位老师的辛勤劳动！