气枪震源激发条件对走时变化观测结果的影响
周青云 陈俊磊
摘要：基于4支单枪容量2000 in3的气枪组成的阵列，2017年2月在宾川地震信号发射台大银甸水库内开展了不同激发气枪压力、沉放深度和水平位移情况下的对比试验。利用发射台周边9个流动宽频带数字地震仪与震源附近3个参考台的数据，使用互相关时延检测的方法计算了不同激发条件下的气枪信号的走时变化，并结合激发条件进行了分析。研究表明：（1）不同气枪压力或沉放深度下，激发的气枪信号的波形和互相关系数存在明顯差异：（2）气枪压力从9MPa升为15MPa，走时变化可达0.06s：（3）走时变化的可能原因是主频段（3～5Hz）中的相对高频成分发生变化，相对低频成分未变化，二者耦合在一起形成了伪走时变化：（4）由于气枪信号不同时窗走时变化不尽相同，技术手段难以消除压力的影响，所以建议未来气枪激发时固定激发条件；利用已有资料做走时变化相关研究时选用震幅相近的数据。
关键词：气枪震源；走时变化；激发条件
中图分类号：P315.3 文献标识码：A 文章编号：1000-0666（2018）02-0264-09
0 引言
板块运动会导致地下介质应力状态发生改变，地下介质应力状态的改变与地震的孕育和发生密切相关，研究地壳介质应力随时间的变化是研究地震规律的有效途径。地震波是获取地球内部结构和介质状态变化等重要而常见的手段，人们对地球内部的知识绝大部分来源于对地震波的研究（陈颙，朱日祥，2005）。随着数字地震台数量的不断增加、资料的逐渐积累和新方法的引入，利用重复地震和背景噪声等方法研究地球介质变化取得了重要进展。
重复地震是指发生在断层同一位置上具有近乎相似的震级和震源机制解，并且在同一台站上波形记录高度相似的一组天然地震（Nadeau et al，1995； Schaff，Beroza，2004）。利用重复地震，王鹏等（2016）研究了山东乳山震群前后的地壳介质波速变化。随着地震背景噪声成像技术的发展（Yao et al，2006：Bensen et al，2010：Zheng et al，2017；郑定昌，王俊，2017），背景噪声也被用于地下介质变化的监测研究中。王俊等（2016）利用背景噪声研究了芦山M7.0地震前后震源区波速变化时空特征。重复地震和背景噪声的方法是被动的等待和利用地震信号，前者的空间重复性和空间分布局限，后者时间分辨率较低并且受噪声源影响较大。
近年来，全国建立了多个气枪震源基地，利用气枪震源主动向地下发射重复性很高的地震波（王伟涛等，2017）。气枪震源具有激发间隔短（最小激发间隔低于10min）、波形之间相关系数高（近场气枪信号相关系数大于0.99）等特点（陈蒙，2014），有效弥补了被动源时间分辨率不足、空间和波形重复性不够等不足之处。魏芸芸等（2016）利用新疆气枪信号分析了新疆2次5级左右地震前后波速变化。
利用主动源信号开展走时变化的研究中，检测到传播路径长为5～50km时的走时变化大多是10-3～10-2s数量级（杨微等，2010；刘自凤等，2015；魏芸芸等，2016），大气压力变化、潮汐应力等因素的影响比走时变化小1个数量级左右（Wang et al，2008； Niu et al，2008；杨微等，2010）。这些研究考虑的都是传播路径上影响因素变化导致的走时变化，没有分析震源变化对走时变化的影响。本文利用2017年2月宾川地震信号发射台集中激发期间发射台周边流动台的数据，计算了不同气枪激发条件下的走时变化并分析了可能的成因，为以后的气枪激发工作提出了建议。
1 实验概况
1.1 宾川地震信号发射台
宾川地震信号发射台位于云南省宾川县城以西的大银甸水库内，建成于2011年，2012年正式投入使用。现拥有4支容量为2000in3的BoltLonglife型气枪及相关配套设施，4支气枪沿水平方向排列成正方形阵列沉放于水面以下10m。气枪激发时，4支气枪同时向水库喷射高压气体产生地震波。配合气枪发射台建立的流动地震台观测台网主要使用Reftec130B和GuralpCMG-40T地震仪，采样率为100Hz（王彬等，2015）。在气枪震源南西架设了3个参考台，将参考台记录的波形当作是气枪子波，与接收台站信号反卷积可以去除震源影响。3个参考台中CKT2更远一些，CKT4与CKTO位置相同但使用了低灵敏度的摆（图1）。已有研究表明，水库气枪阵列单枪信号传播距离超过35km，叠加后波形可见距离超过200km；并中激发的信号单枪传播距离达7km，叠加后波形可见距离超过27km（王彬等，2016）。
1.2 实验概况
2016年2月14-20日宾川地震信号发射台进行的气枪主动源实验中，气枪阵列总共激发161次。笔者就气枪压力、沉放深度和水平位置进行了不同情况下的测试，其中：气枪压力分别为9MPa、12MPa和15MPa（图2a）：沉放深度（气枪阵列到水面的距离）分别为8m、10m、12m（图2b）：水平位置分别为0m、-7m、7m（图2c）。实验共耗时约7d，为保证气枪信号信噪比，尽量选在夜间开展实验。
2 资料的选取及处理
2.1 资料的选取
当研究时间范围内气枪激发次数足够多时，可以利用叠加的方法来提高信号的信噪比，增大气枪信号的接收范围，提高走时变化计算的精度等（刘自凤等，2014；魏芸芸等，2016），但叠加会降低时间分辨率。考虑到本文研究对象的特征及气枪激发次数的限制，结合各流动台记录到单枪波形的信噪比，最终选取气枪周边10km范围内的8个流动台站及较远的1个流动台（图1）。
2.2 资料处理
2.2.1 利用参考台数据提取激发时刻
气枪控制台在气枪激发的过程中会生成激发日志记录了激发时刻，但是由于枪控设备时钟校准、机械延迟、时间记录精度等因素的影响，枪控日志中的激发时刻与真实激发时刻存在较大差异。因此，笔者利用气枪旁边参考台（CKTO）记录到的气枪信号初至时刻当作气枪激发时刻。首先，去除参考台连续记录中的非零均值和长周期线性趋势；再对信号进行主频段3～5Hz（栾奕等，2016）的Butterworth带通滤波；最后，选取参考波形窗口，利用互相关检测方法提取激发时刻及截取实验波形记录。
震中距较远的台站信噪比较低，个别波形可能无法识别。为提高后续计算结果的可信度，我们利用互相关检测的方法剔除了一些信噪比较低的波形：将3～5Hz带通滤波后各台站的161枪信号线性叠加后与每一枪计算相关系数，剔除掉相关系数小于0.6的波形。
2.2.2 利用互相关时延检测技术提取走时变化
利用互相关时延检测技术可以得到亚采样级别（10-3s）的走时变化（Wang et al，2008；劉自凤等，2015）。首先对参考台和流动台的气枪信号进行3～5Hz的带通滤波，降低噪声干扰；随后将同一枪的流动台数据与参考台数据反卷积，对齐到时并去除震源影响，得到台站到参考台的格林函数；将各流动台的格林函数进行叠加，叠加结果作为模板信号；将单枪格林函数与模板干涉，对干涉结果余弦插值后提取得到走时变化。
3 不同激发条件对波形及重复性的影响
3.1 对波形的影响
笔者在实验中改变了气枪压力、沉放深度和水平位移这3个影响因素。理论上气枪压力越大，释放的能量越大，波形的振幅越大；沉放深度和波形的振幅存在一定的关联；平行于气枪水平移动方向的的台站气枪到时会有更明显变化。图3为53265台在不同气枪压力下互相关对齐后的波形，可以看出气枪压力越大，激发地震波的振幅越大（图3a）。图3b～d是波形局部放大图，图中黑色倒三角为波形最低点，到时差异肉眼可见，压力越小，到时越早。激发时间上，3枪仅相差约4h，可以认为外部环境因素及地壳介质状态基本不变，“波速变化”的原因应仅来自于震源处。图4展示的是不同沉放深度的波形图，可以看出与不同压力的波形图具有相似的结果，沉放深度越大，到时越早。
3.2 对互相关系数的影响
相较于炸药、夯车等主动震源，气枪震源优势之一是气枪震源的重复性很好，10km以内台站记录到的信号可以达到0.9以上，震源附近可达到0.99以上（Wang et al.，2012；陈蒙，2014；王彬等，2016）。宾川地震信号发射台气枪阵列激发地震波的主频段大约为3-5Hz（刘自凤等，2012；陈佳等，2016h秕奕等，2016），因此我们将53265台的气枪进行3～5Hz的Butterworth带通滤波后，两两求互相关系数，结果如图5所示。从图5可以看出：不同激发条件相关系数之间形成不同的条带，压力和沉放深度对相关系数均有影响；相关系数条带受压力影响大，受沉放深度影响较小；对角线上形成一系列相关系数大于0.99的正方形（图5中A1、A2、A3等），表明同一种激发条件下激发的波形相关系数很高；高压和低压下激发的波形相关系数可能低于0.80，最低可达到0.73（图5中B1、 B2等）：相关系数超过0.9的占88.7%。
4 不同激发条件对走时变化计算的影响
4.1 走时变化检测原理
如果x（t），y（t）是能量有限信号，则它们的互相关函数可以定义为：式中：“*”表示复共扼。显然，互相关函数是2个信号之间时延：的函数。对于地震仪记录到的2个波形x（t）和y（t），如果二者波形高度相似，只存在时间上的延迟，则可通过计算2个波形的互相关系数函数来得到时延τ。x（t）和y（t）的时延相关系数定义为：式中：T为地震波持续时间。当相关系数取得最大值时，对应的τ就是2个信号的时延。
本文使用的流动台站的采样频率是100Hz，求得的相关函数的采样间隔为0.01s，相关函数的峰值并不会恰好位于采样点上，而计算相关函数时又是一次移动一个采样点，因而真实的走时变化可能和测得的走时变化存在一定差异。因此，对相关函数进行余弦插值，使新得到的相关函数满足精度要求。
4.2 数据处理流程
根据气枪源附近流动台数据的完整性和单枪信号的信噪比，笔者选择实验期间9个流动台的数据和气枪震源附近的3个参考台进行分析（图1）。数据处理流程如下（陈佳等，2017）：（1）去除相关系数较小的气枪信号记录，阈值取0.6；（2）数据预处理：3～5Hz带通滤波，去趋势，去均值；（3）将流动台气枪数据与参考台数据反卷积得到格林函数，以达到对齐到时和压缩子波的目的；（4）格林函数线性叠加作为模板，将每一条格林函数与模板进行互相关时延检测，得到相关函数；（5）对相关函数作余弦插值；（6）提取目标窗口的走时变化。
4.3 走时变化计算结果
利用前述方法处理了9个台站的数据，得到2017年2月14-17日最大振幅附近走时变化图。图6为53265台的走时变化图，其它8个台站的走时变化与53265台一致。从图6可以看出：（1）走时变化和气枪压力高度呈正相关：（2）CKT2对气枪压力的变化更加敏感：（3）难以分辨出沉放深度对走时变化的影响：（4）水平位移对走时变化有影响：（5）气枪气压改变可能导致的走时变化最大可达0.06s（53278台）。
波在水中的传播速度约为1.5km/s，当气枪信号人射角为60°时，水平视速度约为1.3km/s，气枪水平移动14m，走时变化应约为0.01s，与图6显示结果相近。与水平移动相比，沉放深度的改变除影响传播路径的几何响应外，还会影响到震荡气泡与水体的耦合（陈惠芳等，2016）。本次实验设计较不合理，沉放深度改变的同时，其它因素也发生了改变，未能剥离出沉放深度的影响。
5 分析与讨论
5.1 全时窗的走时变化
从以上分析可以看出，选取振幅峰值附近的时间窗，所有台站的走时变化在特定时窗内均呈现与气枪压力呈正相关的特征，表明空间域内这种相关性均存在。那么时间域内，气枪信号的各震相的走时变化是否均存在与气枪压力正相关的特征？为此笔者将53265台0～10.2s的气枪信号以0.3s为窗长，分割为34个窗口，按照前述方法计算了走时变化。为减少影响因素和增加可信度，我们着重关注沉放深度和水平位移均不变、仅气枪压力改变的第70～154枪之间的走时变化。气枪气压剧烈改变前后（第114枪）走时变化的趋势可以简单的分为上升、不变和下降3类，因此将时间上相邻、趋势上一致的走时变化曲线划分为一类，结果见图7。34个窗口共分为10类，其中：走时增加的2类（图7e，i，窗长共3.3 s）：走时减小的3类（图7b，d，g，窗长共2.1s）：走时基本不变的5类（窗长共4.8s）。分析表明，峰值附近震相的走时变化与气枪压力正相关，其它时间窗走时变化与气枪压力变化并非呈现绝对的正相关，也有可能不相关和负相关。
5.2 干涉模板的选择
前面的研究是将格林函数叠加后作为干涉模板计算得到的走时变化结果。为研究不同干涉模板对走时变化计算的影响，笔者使用了前一枪的格林函数及不同枪压下某一枪的格林函数作为干涉模板，计算了53265台时间窗为3.3～3.6s的走时变化（图8）。从图8a可以看出，将前一枪格林函数作为干涉模板，其结果等于将格林函数叠加作为干涉模板的结果的微分，即：
Δt（PreOne）=d[Δt（Stack）]（3）
从图8b可以看出，使用不同激发压力下的格林函数作为模板，其走时变化的结果除平台高度外基本相同。以上分析表明，选择不同的干涉模板计算出来的走时变化基本一致，即干涉模板的选择不是走时变化的成因。
5.3 反卷积的影响
将流动台记录到的气枪信号与参考台的信号反卷积得到的格林函数具有对齐到时、压制噪音的作用，但是会弱化震源信号不重复带来的影响。笔者采用2.2节中的数据处理流程但不进行反卷积，利用波形记录提取走时变化，结果见图9。从图9中可以看出，反卷积后格林函数提取的走时变化幅度约0.02s，远小于波形记录提取的走时变化幅度0.15s，表明反卷积会弱化震源不重复的影响。波形记录提取的走时变化与激发压力的相关系数较低，3个周期分别为0.30、0.29和0.23，而格林函数的相关系数为0.92。综合分析表明，使用格林函数计算走时变化会弱化震源不重复性的影响，但能提高计算结果和激发条件的相关程度。
5.4 成因推测
气枪产生的地震近场子波可分为压力脉冲和气泡脉冲，前者源于高压气体喷出瞬间的高频（>20Hz）脉冲，后者源于气泡在水中的低频（<10Hz）震荡（夏季等，2016），本文使用的是气泡震荡产生的低频段信号。有学者认为气枪气压越大，气泡震荡周期越大，主频越低（何汉漪，2001；陈蒙，2014），也有学者发现气枪压力对主脉冲与第1个气泡脉冲的时间间隔影响显著，对主频影响不显著（夏季等，2016）。本文研究结果与后者一致。
笔者计算了9MPa和15MPa气枪信号的归一化傅氏谱及时频图，将二者的时频谱相减，发现：（1）傅氏谱的峰值均为3.45Hz：（2）3～3.7Hz频段、4.1～4.3频段二者基本一致，3.7～4.1Hz频段、4.3～5Hz频段差异较大（图10a）；（3）二者的归一化时频图4.3Hz以下形态基本一致，4.3Hz以上形态存在差异（图10b～c）；（4）图10d中存在竖向条纹。这说明气枪压力的改变影响的并不是信号的主频，而是更高的一些频段，因此我們认为走时变化的成因是：气枪压力改变导致气枪子波主频段中相对高频成分发生变化，变化后的高频成分与未变化的低频成分耦合在一起，形成了伪走时变化。
5.5 建议
利用50km范围内台站接收到的主动源信号监测到的走时变化一般在10-2～10-3s数量级，例如：魏芸芸等（2016）观测到沙湾M5.0地震前后存在0.05～0.06s的走时变化；杨微等（2010）观测到绵竹M5.6地震前后存在5～9ms的走时变化；刘自凤等（2015）观测到宾川主动源无地震期间走时变化波动幅度±0.02s；王宝善等（2016）利用宾川主动源计算的走时日变化数量级为10-3s。这些走时变化与气枪气压导致的走时变化大小相当或更小，因此走时变化的研究中需去除气枪压力变化的影响。
通过前面的分析可以看出，气枪震源的主频段（3～5Hz）受气枪压力影响明显且无明显规律，接收台站的位置、用于反卷积的参考台的位置、选择的时间窗口的位置和长度对走时变化计算均有影响，因此难以给出一个可靠的去除气枪压力变化影响的方法。建议后续激发实验中，选用恒定的激发压力；对于已有的数据，建议使用震幅基本一致的气枪信号。
6 结论
本文利用2017年2月宾川地震信号发射台集中激发的数据，研究了不同激发压力、沉放深度和水平位移条件下气枪信号的走时变化特征，取得了以下研究结果：
（1）不同压力和沉放深度下的波形存在肉眼可见的差别；相同激发条件下激发的地震波的主频段相关系数可达到0.99，不同激发条件下的相关系数最低为0.73。
（2）当激发条件发生变化时，即使波形相关系数大于0.9，也有可能产生一些伪走时变化，这种走时变化反应了震源的变化而非介质或外界因素的变化，反卷积并不能完全去掉震源的影响。气枪气压变化对走时变化影响较大，气枪信号震幅峰值附近震相的走时变化与压力变化成正比，对于1～15km范围内的台站，最大可达0.06s；无法分辨出沉放深度对走时变化的影响，但根据其对波形的影响可以推测出对走时的影响可能达到10-2s级别；水平位移对走时变化的影响主要体现在传播几何路径改变上。
（3）气枪信号不同震相的走时变化与压力变化的相关性存在差异：峰值附近的走时变化与气压变化正相关，其它震相有可能正相关，也可能不相关或负相关。据傅氏谱及时频谱推测，其成因可能是气枪压力改变后，信号主频段中的相对高频成分发生了变化，变化后的高频成分与未变化的低频成分耦合在一起，形成了伪走时变化。
（4）气枪压力变化导致的走时变化最大可达0.06s，与已有研究中地震导致的走时变化相当或更大，因此相关研究中应考虑激发压力的影响。建议以后气枪的激发压力恒定；已有数据研究走时变化时，建议使用震幅基本一致的气枪信号。
两位审稿老师给出了宝贵意见，云南省地震局主动源创新团队各成员为开展实验所付出的艰辛努力，本文使用了中国地震局地球物理研究所王宝善老师提供的走时变化计算程序，在此一并表示感谢。
参考文献：
陈惠芳，林彬华，金星，等.2016，水库大容量气枪震源激发条件优化实验研究[J].中国地震，32（2）：241-248.
陈佳，李孝宾，杨军，等.2016.云南宾川大容量气枪震源波形频谱特征分析[J].中国地震，32（2）：216-221.
陳佳，叶泵，高琼，等.2017.利用气枪震源信号研究2016年云龙MS5.0地震前后波速变化特征[J].地震研究，40（4）：550-556.
陈蒙.2014.利用水库大容量非调制气枪阵列进行区域尺度地下结构探测和监测[D].北京：中国地震局地球物理研究所.
陈颙，朱日祥.2005.设立“地下明灯研究计划”的建议[J].地球科学进展，20（5）：485-489.
何汉漪.2001.海上高分辨率地震技术及其应用[M].北京：地质出版社，52-77
刘自凤，苏有锦，王宝善，等.2015.宾川主动源地震波走时变化分析方法研究[J].地震研究，38（4）：591-597
栾奕，杨宏峰，王宝善，等.2016.大容量气枪主动源波形资料处理（一）：云南宾川[J].中国地震，32（2）：305-318.
王宝善，葛洪魁，王彬，等.2016.利用人工重复震源进行地下介质结构及其变化研究的探索和进展[J].中国地震，32（2）：168-179.
王彬，李孝宾，刘自凤，等.2016.宾川地震信号发射台的震源系统、观测系统和观测结果[J].中国地震，32（2）：193-201.
王彬，吴国华，苏有锦，等.2015.宾川地震信号发射台的选址、建设及初步观测结果[J].地震研究，38（1）：1-6.
王俊，郑定昌，郑江蓉，等.2016：利用背景噪声自相关研究芦山M7.0地震震源区地壳相对波速的时空变化特征[J].地震地质，38（1）：15-168
王鹏，郑建常，谭毅培.2016.利用重复地震研究山东乳山地区地壳介质波速变化[J].地震学报，38（5）：728-738.
王伟涛，王宝善，蒋生森，等.2017.利用气枪震源探测大陆浅部的地震学研究回顾与展望[J].地震研究，40（4）：514-524.
魏芸芸，王海涛，苏金波，等.2016.新疆2次中强地震前气枪震源反射波震相走时异常变化初步研究[J].中国地震，32（2）：270-281.
夏季，金星，蔡辉腾，等.2016.大容量气枪震源子波时频特性及其影响因素[J].中国地震，32（2）：249-260.
杨微，葛洪魁，王宝善，等.2010.由精密控制人工震源观测到的绵竹5.6级地震前后波速变化[J].地球物理学报，53（5）：1149-1157.
郑定昌，王俊.2017.基于背景噪声的川滇地区勒夫波层析成像[J].地震学报，39（5）：633-647.
Bensen G D，Ritzwoller M H，Barmin M P，et at.2010.Processing seismicambient noise data to obtain rdiphle broad-band surface wave dis-persion measurements[J].Geophysical Journal of the Royal Astrv-nomibal Society，169（3）：1239-1260.
Nadeau R M，Foxall W，McEVIlly T V.1995.clustering and periodic。-currence of microearthquakes on the San Andreas faults at Parkfield，California[J].Science，267（5179）：503-507.
Niu F，Silver P G，Daley T M，et al.2008.Preseismic velocity changes ob-served from active source monitoring at the Parkfield SAFOD drillsite[J].Nature，454：204-208，doi：10.1038/nature07111.
Schaff D P，Beroza G C.2004.Coseismic and postseismic velocity changesmeasured by repeating earthquakes[J] .journal of Geophysical Re-search，109，B10302，doi：10.1029/2004JB003011.
Wang Rs，Ge H K，Yang W，et al.2012.Transmitting seismic station mo-nitors fault zone at depth[J].Eos Transactions American Geophysi-cal Union，93（5）：49-50.
Wang Bs，Zhu P，Chen Y，et al.2008.Continuous subsurface velocitymeasurement with coda wave interferometry[J].Geophys Res，113，B12313，doi：10.1029/2007JB005023.
Yao H，Vander Hilst RD，De Hoop M V.2006.Surface -wave arraytomography in SE Tibet from ambient seismic noise and two -stationanalysis -I，Phase velocity maps[J].Geophysical Journal Interna-
Zheng D C，Saygin E，Cummins P，et al.2017.Tranadimensional Bayesianseismic ambient noise tomography across SE Tibet[J].Journal of A-sian Earth Sciences，134：86 -93.
摘要：基于4支单枪容量2000 in3的气枪组成的阵列，2017年2月在宾川地震信号发射台大银甸水库内开展了不同激发气枪压力、沉放深度和水平位移情况下的对比试验。利用发射台周边9个流动宽频带数字地震仪与震源附近3个参考台的数据，使用互相关时延检测的方法计算了不同激发条件下的气枪信号的走时变化，并结合激发条件进行了分析。研究表明：（1）不同气枪压力或沉放深度下，激发的气枪信号的波形和互相关系数存在明顯差异：（2）气枪压力从9MPa升为15MPa，走时变化可达0.06s：（3）走时变化的可能原因是主频段（3～5Hz）中的相对高频成分发生变化，相对低频成分未变化，二者耦合在一起形成了伪走时变化：（4）由于气枪信号不同时窗走时变化不尽相同，技术手段难以消除压力的影响，所以建议未来气枪激发时固定激发条件；利用已有资料做走时变化相关研究时选用震幅相近的数据。
关键词：气枪震源；走时变化；激发条件
中图分类号：P315.3 文献标识码：A 文章编号：1000-0666（2018）02-0264-09
0 引言
板块运动会导致地下介质应力状态发生改变，地下介质应力状态的改变与地震的孕育和发生密切相关，研究地壳介质应力随时间的变化是研究地震规律的有效途径。地震波是获取地球内部结构和介质状态变化等重要而常见的手段，人们对地球内部的知识绝大部分来源于对地震波的研究（陈颙，朱日祥，2005）。随着数字地震台数量的不断增加、资料的逐渐积累和新方法的引入，利用重复地震和背景噪声等方法研究地球介质变化取得了重要进展。
重复地震是指发生在断层同一位置上具有近乎相似的震级和震源机制解，并且在同一台站上波形记录高度相似的一组天然地震（Nadeau et al，1995； Schaff，Beroza，2004）。利用重复地震，王鹏等（2016）研究了山东乳山震群前后的地壳介质波速变化。随着地震背景噪声成像技术的发展（Yao et al，2006：Bensen et al，2010：Zheng et al，2017；郑定昌，王俊，2017），背景噪声也被用于地下介质变化的监测研究中。王俊等（2016）利用背景噪声研究了芦山M7.0地震前后震源区波速变化时空特征。重复地震和背景噪声的方法是被动的等待和利用地震信号，前者的空间重复性和空间分布局限，后者时间分辨率较低并且受噪声源影响较大。
近年来，全国建立了多个气枪震源基地，利用气枪震源主动向地下发射重复性很高的地震波（王伟涛等，2017）。气枪震源具有激发间隔短（最小激发间隔低于10min）、波形之间相关系数高（近场气枪信号相关系数大于0.99）等特点（陈蒙，2014），有效弥补了被动源时间分辨率不足、空间和波形重复性不够等不足之处。魏芸芸等（2016）利用新疆气枪信号分析了新疆2次5级左右地震前后波速变化。
利用主动源信号开展走时变化的研究中，检测到传播路径长为5～50km时的走时变化大多是10-3～10-2s数量级（杨微等，2010；刘自凤等，2015；魏芸芸等，2016），大气压力变化、潮汐应力等因素的影响比走时变化小1个数量级左右（Wang et al，2008； Niu et al，2008；杨微等，2010）。这些研究考虑的都是传播路径上影响因素变化导致的走时变化，没有分析震源变化对走时变化的影响。本文利用2017年2月宾川地震信号发射台集中激发期间发射台周边流动台的数据，计算了不同气枪激发条件下的走时变化并分析了可能的成因，为以后的气枪激发工作提出了建议。
1 实验概况
1.1 宾川地震信号发射台
宾川地震信号发射台位于云南省宾川县城以西的大银甸水库内，建成于2011年，2012年正式投入使用。现拥有4支容量为2000in3的BoltLonglife型气枪及相关配套设施，4支气枪沿水平方向排列成正方形阵列沉放于水面以下10m。气枪激发时，4支气枪同时向水库喷射高压气体产生地震波。配合气枪发射台建立的流动地震台观测台网主要使用Reftec130B和GuralpCMG-40T地震仪，采样率为100Hz（王彬等，2015）。在气枪震源南西架设了3个参考台，将参考台记录的波形当作是气枪子波，与接收台站信号反卷积可以去除震源影响。3个参考台中CKT2更远一些，CKT4与CKTO位置相同但使用了低灵敏度的摆（图1）。已有研究表明，水库气枪阵列单枪信号传播距离超过35km，叠加后波形可见距离超过200km；并中激发的信号单枪传播距离达7km，叠加后波形可见距离超过27km（王彬等，2016）。
1.2 实验概况
2016年2月14-20日宾川地震信号发射台进行的气枪主动源实验中，气枪阵列总共激发161次。笔者就气枪压力、沉放深度和水平位置进行了不同情况下的测试，其中：气枪压力分别为9MPa、12MPa和15MPa（图2a）：沉放深度（气枪阵列到水面的距离）分别为8m、10m、12m（图2b）：水平位置分别为0m、-7m、7m（图2c）。实验共耗时约7d，为保证气枪信号信噪比，尽量选在夜间开展实验。
2 资料的选取及处理
2.1 资料的选取
当研究时间范围内气枪激发次数足够多时，可以利用叠加的方法来提高信号的信噪比，增大气枪信号的接收范围，提高走时变化计算的精度等（刘自凤等，2014；魏芸芸等，2016），但叠加会降低时间分辨率。考虑到本文研究对象的特征及气枪激发次数的限制，结合各流动台记录到单枪波形的信噪比，最终选取气枪周边10km范围内的8个流动台站及较远的1个流动台（图1）。
2.2 资料处理
2.2.1 利用参考台数据提取激发时刻
气枪控制台在气枪激发的过程中会生成激发日志记录了激发时刻，但是由于枪控设备时钟校准、机械延迟、时间记录精度等因素的影响，枪控日志中的激发时刻与真实激发时刻存在较大差异。因此，笔者利用气枪旁边参考台（CKTO）记录到的气枪信号初至时刻当作气枪激发时刻。首先，去除参考台连续记录中的非零均值和长周期线性趋势；再对信号进行主频段3～5Hz（栾奕等，2016）的Butterworth带通滤波；最后，选取参考波形窗口，利用互相关检测方法提取激发时刻及截取实验波形记录。
震中距较远的台站信噪比较低，个别波形可能无法识别。为提高后续计算结果的可信度，我们利用互相关检测的方法剔除了一些信噪比较低的波形：将3～5Hz带通滤波后各台站的161枪信号线性叠加后与每一枪计算相关系数，剔除掉相关系数小于0.6的波形。
2.2.2 利用互相关时延检测技术提取走时变化
利用互相关时延检测技术可以得到亚采样级别（10-3s）的走时变化（Wang et al，2008；劉自凤等，2015）。首先对参考台和流动台的气枪信号进行3～5Hz的带通滤波，降低噪声干扰；随后将同一枪的流动台数据与参考台数据反卷积，对齐到时并去除震源影响，得到台站到参考台的格林函数；将各流动台的格林函数进行叠加，叠加结果作为模板信号；将单枪格林函数与模板干涉，对干涉结果余弦插值后提取得到走时变化。
3 不同激发条件对波形及重复性的影响
3.1 对波形的影响
笔者在实验中改变了气枪压力、沉放深度和水平位移这3个影响因素。理论上气枪压力越大，释放的能量越大，波形的振幅越大；沉放深度和波形的振幅存在一定的关联；平行于气枪水平移动方向的的台站气枪到时会有更明显变化。图3为53265台在不同气枪压力下互相关对齐后的波形，可以看出气枪压力越大，激发地震波的振幅越大（图3a）。图3b～d是波形局部放大图，图中黑色倒三角为波形最低点，到时差异肉眼可见，压力越小，到时越早。激发时间上，3枪仅相差约4h，可以认为外部环境因素及地壳介质状态基本不变，“波速变化”的原因应仅来自于震源处。图4展示的是不同沉放深度的波形图，可以看出与不同压力的波形图具有相似的结果，沉放深度越大，到时越早。
3.2 对互相关系数的影响
相较于炸药、夯车等主动震源，气枪震源优势之一是气枪震源的重复性很好，10km以内台站记录到的信号可以达到0.9以上，震源附近可达到0.99以上（Wang et al.，2012；陈蒙，2014；王彬等，2016）。宾川地震信号发射台气枪阵列激发地震波的主频段大约为3-5Hz（刘自凤等，2012；陈佳等，2016h秕奕等，2016），因此我们将53265台的气枪进行3～5Hz的Butterworth带通滤波后，两两求互相关系数，结果如图5所示。从图5可以看出：不同激发条件相关系数之间形成不同的条带，压力和沉放深度对相关系数均有影响；相关系数条带受压力影响大，受沉放深度影响较小；对角线上形成一系列相关系数大于0.99的正方形（图5中A1、A2、A3等），表明同一种激发条件下激发的波形相关系数很高；高压和低压下激发的波形相关系数可能低于0.80，最低可达到0.73（图5中B1、 B2等）：相关系数超过0.9的占88.7%。
4 不同激发条件对走时变化计算的影响
4.1 走时变化检测原理
如果x（t），y（t）是能量有限信号，则它们的互相关函数可以定义为：式中：“*”表示复共扼。显然，互相关函数是2个信号之间时延：的函数。对于地震仪记录到的2个波形x（t）和y（t），如果二者波形高度相似，只存在时间上的延迟，则可通过计算2个波形的互相关系数函数来得到时延τ。x（t）和y（t）的时延相关系数定义为：式中：T为地震波持续时间。当相关系数取得最大值时，对应的τ就是2个信号的时延。
本文使用的流动台站的采样频率是100Hz，求得的相关函数的采样间隔为0.01s，相关函数的峰值并不会恰好位于采样点上，而计算相关函数时又是一次移动一个采样点，因而真实的走时变化可能和测得的走时变化存在一定差异。因此，对相关函数进行余弦插值，使新得到的相关函数满足精度要求。
4.2 数据处理流程
根据气枪源附近流动台数据的完整性和单枪信号的信噪比，笔者选择实验期间9个流动台的数据和气枪震源附近的3个参考台进行分析（图1）。数据处理流程如下（陈佳等，2017）：（1）去除相关系数较小的气枪信号记录，阈值取0.6；（2）数据预处理：3～5Hz带通滤波，去趋势，去均值；（3）将流动台气枪数据与参考台数据反卷积得到格林函数，以达到对齐到时和压缩子波的目的；（4）格林函数线性叠加作为模板，将每一条格林函数与模板进行互相关时延检测，得到相关函数；（5）对相关函数作余弦插值；（6）提取目标窗口的走时变化。
4.3 走时变化计算结果
利用前述方法处理了9个台站的数据，得到2017年2月14-17日最大振幅附近走时变化图。图6为53265台的走时变化图，其它8个台站的走时变化与53265台一致。从图6可以看出：（1）走时变化和气枪压力高度呈正相关：（2）CKT2对气枪压力的变化更加敏感：（3）难以分辨出沉放深度对走时变化的影响：（4）水平位移对走时变化有影响：（5）气枪气压改变可能导致的走时变化最大可达0.06s（53278台）。
波在水中的传播速度约为1.5km/s，当气枪信号人射角为60°时，水平视速度约为1.3km/s，气枪水平移动14m，走时变化应约为0.01s，与图6显示结果相近。与水平移动相比，沉放深度的改变除影响传播路径的几何响应外，还会影响到震荡气泡与水体的耦合（陈惠芳等，2016）。本次实验设计较不合理，沉放深度改变的同时，其它因素也发生了改变，未能剥离出沉放深度的影响。
5 分析与讨论
5.1 全时窗的走时变化
从以上分析可以看出，选取振幅峰值附近的时间窗，所有台站的走时变化在特定时窗内均呈现与气枪压力呈正相关的特征，表明空间域内这种相关性均存在。那么时间域内，气枪信号的各震相的走时变化是否均存在与气枪压力正相关的特征？为此笔者将53265台0～10.2s的气枪信号以0.3s为窗长，分割为34个窗口，按照前述方法计算了走时变化。为减少影响因素和增加可信度，我们着重关注沉放深度和水平位移均不变、仅气枪压力改变的第70～154枪之间的走时变化。气枪气压剧烈改变前后（第114枪）走时变化的趋势可以简单的分为上升、不变和下降3类，因此将时间上相邻、趋势上一致的走时变化曲线划分为一类，结果见图7。34个窗口共分为10类，其中：走时增加的2类（图7e，i，窗长共3.3 s）：走时减小的3类（图7b，d，g，窗长共2.1s）：走时基本不变的5类（窗长共4.8s）。分析表明，峰值附近震相的走时变化与气枪压力正相关，其它时间窗走时变化与气枪压力变化并非呈现绝对的正相关，也有可能不相关和负相关。
5.2 干涉模板的选择
前面的研究是将格林函数叠加后作为干涉模板计算得到的走时变化结果。为研究不同干涉模板对走时变化计算的影响，笔者使用了前一枪的格林函数及不同枪压下某一枪的格林函数作为干涉模板，计算了53265台时间窗为3.3～3.6s的走时变化（图8）。从图8a可以看出，将前一枪格林函数作为干涉模板，其结果等于将格林函数叠加作为干涉模板的结果的微分，即：
Δt（PreOne）=d[Δt（Stack）]（3）
从图8b可以看出，使用不同激发压力下的格林函数作为模板，其走时变化的结果除平台高度外基本相同。以上分析表明，选择不同的干涉模板计算出来的走时变化基本一致，即干涉模板的选择不是走时变化的成因。
5.3 反卷积的影响
将流动台记录到的气枪信号与参考台的信号反卷积得到的格林函数具有对齐到时、压制噪音的作用，但是会弱化震源信号不重复带来的影响。笔者采用2.2节中的数据处理流程但不进行反卷积，利用波形记录提取走时变化，结果见图9。从图9中可以看出，反卷积后格林函数提取的走时变化幅度约0.02s，远小于波形记录提取的走时变化幅度0.15s，表明反卷积会弱化震源不重复的影响。波形记录提取的走时变化与激发压力的相关系数较低，3个周期分别为0.30、0.29和0.23，而格林函数的相关系数为0.92。综合分析表明，使用格林函数计算走时变化会弱化震源不重复性的影响，但能提高计算结果和激发条件的相关程度。
5.4 成因推测
气枪产生的地震近场子波可分为压力脉冲和气泡脉冲，前者源于高压气体喷出瞬间的高频（>20Hz）脉冲，后者源于气泡在水中的低频（<10Hz）震荡（夏季等，2016），本文使用的是气泡震荡产生的低频段信号。有学者认为气枪气压越大，气泡震荡周期越大，主频越低（何汉漪，2001；陈蒙，2014），也有学者发现气枪压力对主脉冲与第1个气泡脉冲的时间间隔影响显著，对主频影响不显著（夏季等，2016）。本文研究结果与后者一致。
笔者计算了9MPa和15MPa气枪信号的归一化傅氏谱及时频图，将二者的时频谱相减，发现：（1）傅氏谱的峰值均为3.45Hz：（2）3～3.7Hz频段、4.1～4.3频段二者基本一致，3.7～4.1Hz频段、4.3～5Hz频段差异较大（图10a）；（3）二者的归一化时频图4.3Hz以下形态基本一致，4.3Hz以上形态存在差异（图10b～c）；（4）图10d中存在竖向条纹。这说明气枪压力的改变影响的并不是信号的主频，而是更高的一些频段，因此我們认为走时变化的成因是：气枪压力改变导致气枪子波主频段中相对高频成分发生变化，变化后的高频成分与未变化的低频成分耦合在一起，形成了伪走时变化。
5.5 建议
利用50km范围内台站接收到的主动源信号监测到的走时变化一般在10-2～10-3s数量级，例如：魏芸芸等（2016）观测到沙湾M5.0地震前后存在0.05～0.06s的走时变化；杨微等（2010）观测到绵竹M5.6地震前后存在5～9ms的走时变化；刘自凤等（2015）观测到宾川主动源无地震期间走时变化波动幅度±0.02s；王宝善等（2016）利用宾川主动源计算的走时日变化数量级为10-3s。这些走时变化与气枪气压导致的走时变化大小相当或更小，因此走时变化的研究中需去除气枪压力变化的影响。
通过前面的分析可以看出，气枪震源的主频段（3～5Hz）受气枪压力影响明显且无明显规律，接收台站的位置、用于反卷积的参考台的位置、选择的时间窗口的位置和长度对走时变化计算均有影响，因此难以给出一个可靠的去除气枪压力变化影响的方法。建议后续激发实验中，选用恒定的激发压力；对于已有的数据，建议使用震幅基本一致的气枪信号。
6 结论
本文利用2017年2月宾川地震信号发射台集中激发的数据，研究了不同激发压力、沉放深度和水平位移条件下气枪信号的走时变化特征，取得了以下研究结果：
（1）不同压力和沉放深度下的波形存在肉眼可见的差别；相同激发条件下激发的地震波的主频段相关系数可达到0.99，不同激发条件下的相关系数最低为0.73。
（2）当激发条件发生变化时，即使波形相关系数大于0.9，也有可能产生一些伪走时变化，这种走时变化反应了震源的变化而非介质或外界因素的变化，反卷积并不能完全去掉震源的影响。气枪气压变化对走时变化影响较大，气枪信号震幅峰值附近震相的走时变化与压力变化成正比，对于1～15km范围内的台站，最大可达0.06s；无法分辨出沉放深度对走时变化的影响，但根据其对波形的影响可以推测出对走时的影响可能达到10-2s级别；水平位移对走时变化的影响主要体现在传播几何路径改变上。
（3）气枪信号不同震相的走时变化与压力变化的相关性存在差异：峰值附近的走时变化与气压变化正相关，其它震相有可能正相关，也可能不相关或负相关。据傅氏谱及时频谱推测，其成因可能是气枪压力改变后，信号主频段中的相对高频成分发生了变化，变化后的高频成分与未变化的低频成分耦合在一起，形成了伪走时变化。
（4）气枪压力变化导致的走时变化最大可达0.06s，与已有研究中地震导致的走时变化相当或更大，因此相关研究中应考虑激发压力的影响。建议以后气枪的激发压力恒定；已有数据研究走时变化时，建议使用震幅基本一致的气枪信号。
两位审稿老师给出了宝贵意见，云南省地震局主动源创新团队各成员为开展实验所付出的艰辛努力，本文使用了中国地震局地球物理研究所王宝善老师提供的走时变化计算程序，在此一并表示感谢。
参考文献：
陈惠芳，林彬华，金星，等.2016，水库大容量气枪震源激发条件优化实验研究[J].中国地震，32（2）：241-248.
陈佳，李孝宾，杨军，等.2016.云南宾川大容量气枪震源波形频谱特征分析[J].中国地震，32（2）：216-221.
陳佳，叶泵，高琼，等.2017.利用气枪震源信号研究2016年云龙MS5.0地震前后波速变化特征[J].地震研究，40（4）：550-556.
陈蒙.2014.利用水库大容量非调制气枪阵列进行区域尺度地下结构探测和监测[D].北京：中国地震局地球物理研究所.
陈颙，朱日祥.2005.设立“地下明灯研究计划”的建议[J].地球科学进展，20（5）：485-489.
何汉漪.2001.海上高分辨率地震技术及其应用[M].北京：地质出版社，52-77
刘自凤，苏有锦，王宝善，等.2015.宾川主动源地震波走时变化分析方法研究[J].地震研究，38（4）：591-597
栾奕，杨宏峰，王宝善，等.2016.大容量气枪主动源波形资料处理（一）：云南宾川[J].中国地震，32（2）：305-318.
王宝善，葛洪魁，王彬，等.2016.利用人工重复震源进行地下介质结构及其变化研究的探索和进展[J].中国地震，32（2）：168-179.
王彬，李孝宾，刘自凤，等.2016.宾川地震信号发射台的震源系统、观测系统和观测结果[J].中国地震，32（2）：193-201.
王彬，吴国华，苏有锦，等.2015.宾川地震信号发射台的选址、建设及初步观测结果[J].地震研究，38（1）：1-6.
王俊，郑定昌，郑江蓉，等.2016：利用背景噪声自相关研究芦山M7.0地震震源区地壳相对波速的时空变化特征[J].地震地质，38（1）：15-168
王鹏，郑建常，谭毅培.2016.利用重复地震研究山东乳山地区地壳介质波速变化[J].地震学报，38（5）：728-738.
王伟涛，王宝善，蒋生森，等.2017.利用气枪震源探测大陆浅部的地震学研究回顾与展望[J].地震研究，40（4）：514-524.
魏芸芸，王海涛，苏金波，等.2016.新疆2次中强地震前气枪震源反射波震相走时异常变化初步研究[J].中国地震，32（2）：270-281.
夏季，金星，蔡辉腾，等.2016.大容量气枪震源子波时频特性及其影响因素[J].中国地震，32（2）：249-260.
杨微，葛洪魁，王宝善，等.2010.由精密控制人工震源观测到的绵竹5.6级地震前后波速变化[J].地球物理学报，53（5）：1149-1157.
郑定昌，王俊.2017.基于背景噪声的川滇地区勒夫波层析成像[J].地震学报，39（5）：633-647.
Bensen G D，Ritzwoller M H，Barmin M P，et at.2010.Processing seismicambient noise data to obtain rdiphle broad-band surface wave dis-persion measurements[J].Geophysical Journal of the Royal Astrv-nomibal Society，169（3）：1239-1260.
Nadeau R M，Foxall W，McEVIlly T V.1995.clustering and periodic。-currence of microearthquakes on the San Andreas faults at Parkfield，California[J].Science，267（5179）：503-507.
Niu F，Silver P G，Daley T M，et al.2008.Preseismic velocity changes ob-served from active source monitoring at the Parkfield SAFOD drillsite[J].Nature，454：204-208，doi：10.1038/nature07111.
Schaff D P，Beroza G C.2004.Coseismic and postseismic velocity changesmeasured by repeating earthquakes[J] .journal of Geophysical Re-search，109，B10302，doi：10.1029/2004JB003011.
Wang Rs，Ge H K，Yang W，et al.2012.Transmitting seismic station mo-nitors fault zone at depth[J].Eos Transactions American Geophysi-cal Union，93（5）：49-50.
Wang Bs，Zhu P，Chen Y，et al.2008.Continuous subsurface velocitymeasurement with coda wave interferometry[J].Geophys Res，113，B12313，doi：10.1029/2007JB005023.
Yao H，Vander Hilst RD，De Hoop M V.2006.Surface -wave arraytomography in SE Tibet from ambient seismic noise and two -stationanalysis -I，Phase velocity maps[J].Geophysical Journal Interna-
Zheng D C，Saygin E，Cummins P，et al.2017.Tranadimensional Bayesianseismic ambient noise tomography across SE Tibet[J].Journal of A-sian Earth Sciences，134：86 -93.
