基于波形相似特征的微地震事件初至拾取及全局校正

微地震事件初至拾取是井下微地震监测数据处理的关键步骤之一.初至误差的存在会使微地震震源定位结果产生较大偏差,进而影响后续的压裂裂缝解释.通常初至拾取过程对所有的微地震事件选择相同的特征函数并采用一致的拾取参数进行统一处理,然而当事件的能量、震源机制、传播路径以及背景噪声等存在明显差异时,所得初至拾取结果差别显著.为了提高微地震事件初至拾取标准一致性,本文提出基于波形相似特征的初至拾取及全局校正方法.该方法首先利用互相关函数对每个事件内的各道记录进行时差校正,得到初始初至信息并形成叠加道,再对所有事件的叠加道进行全局互相关得到事件间初至相对校正量,最终初至结果可以通过各个事件的初始初至信息与其相对校正量相加得到.方法将所有微地震事件初至结果作为一个整体处理,从而能够克服常规方法初至拾取标准一致性差的缺陷.实际资料处理结果表明,相比于常规方法,该方法可以有效提高事件初至拾取和定位结果的一致性.     

Improvements in mining induced microseismic source locations at the Lucky Friday mine using an automated whole-waveform analysis system

For years, severe rockburst problems at the Lucky Friday mine in northern Idaho have been a persistent safety hazard and an impediment to production. An MP250 based microseismic monitoring system, which uses simple voltage threshold picking of first arrivals, has been used in this mine since 1973 to provide source locations and energy estimates of seismic events. Recently, interest has been expressed in developing a whole waveform microseismic monitoring system for the mine to provide more accurate source locations and information about source characteristics. For this study, we have developed a prototype whole-waveform microseismic monitoring system based on a 80386 computer equipped with a 50 kHz analog-digital convertor board. The software developed includes a data collection program, a data analysis program, and an event detection program. Whole-waveform data collected and analyzed using this system during a three-day test have been employed to investigate sources of error in the hypocenter location process and to develop an automatic phase picker appropriate for microseismic events.Comparison of hypocenter estimates produced by the MP250 system to those produced by the whole-waveform system shows that significant timing errors are common in the MP250 system and that these errors caused a large part of the scatter evident in the daily activity plots produced at the mine. Simulations and analysis of blast data show that analytical control over the solutions is strongly influenced by the array geometry. Within the geophone array, large errors in the velocity model or moderate timing errors may result in small changes in the solution, but outside the array, the solution is very sensitive to small changes in the data.Our whole-waveform detection program picks event onset times and determines event durations by analysis of a segmented envelope function (SEF) derived from the microseismic signal. The detection program has been tested by comparing its arrival time picks to those generated by human analysis of the data set. The program picked 87% of the channels that were picked by hand with a standard error of 0.75 milliseconds. Source locations calculated using times provided by our entire waveform detection program were similar to those calculated using hand-picked arrival times. In particular, they show far less scatter than source locations calculated using arrival times based on simple voltage threshold picking of first arrivals.     

Hypocenter location of mixed events in real-time processing systems

Automatic pickings in earthquake real-time monitoring systems often contain noise bursts and/or phases of different event(s) occurring almost simultaneously. Typically, a locator uses these picks as P and S waves arrival times coming from a single event and, therefore, should be complemented by a distinctive phase association logic. The method we propose manages to automatically associate data related to different events and eliminates the influence of spoiled data from single events. The method is based on “network beamforming”, a robust and stable algorithm, which utilizes a hypocenter grid search for the stack maximum of a set of complex exponents applied to the P phase readings. The algorithm separates the residual outliers and then uses them for location. If successful, a hypocenter is established for the interfering event. The solutions obtained are overall robust and independent from the estimate of origin times. The preliminary epicenter for the grid search is provided by the intersection of perpendicular bisectors in the modified “arrival order algorithm” or by the modified “T^now” algorithm, which uses non-arrival information. We applied this method to automatic first arrival phase readings of 915 events registered by the Hi-net Japan seismic network and our results are statistically promising. Here, we present two interesting and complicated examples.     

Seismic phase picking in China Seismic Array using a deep convolutional neuron network*

Seismic phase picking is the preliminary work of earthquake location and body-wave travel time tomography. Manual picking is considered as the most accurate way to access the arrival times but time consuming. Many automatic picking methods were proposed in the past decades, but their precisions are not as high as human experts especially for events with low ratio of signal to noise and later arrivals. As the increasing deployment of large seismic array, the existing methods can not meet the requirements of quick and accurate phase picking. In this study, we applied a phase picking algorithm developed on the base of deep convolutional neuron network (PickNet) to pick seismic phase arrivals in ChinArray-Phase III. The comparison of picking error of PickNet and the traditional method shows that PickNet is capable of picking more precise phases and can be applied in a large dense array. The raw picked travel-time data shows a large variation deviated from the traveltime curves. The absolute location residual is a key criteria for travel-time data selection. Besides, we proposed a flowchart to determine the accurate location of the single-station earthquake via dense seismic array and phase arrival picked by PickNet. This research expands the phase arrival dataset and improves the location accuracy of single-station earthquake.     

The Influence of Path Corrections and a Three-dimensional Global <Emphasis Type="Italic">P</Emphasis>-wave Velocity Model on Seismic Event Location in Kazakhstan

We relocate 81 large nuclear explosions that were detonated at the Balapan and Degelen Mountain subregions of the Semipalatinsk test site in Kazakhstan during the years 1978 to 1989. The absolute locations of these explosions are available, as well as very accurate estimates of their origin times. This ground truth information allows us to perform a detailed analysis of location capability. We use a sparse network of stations with highly accurate first arrival picks measured using a waveform cross-correlation method. These high quality data facilitate very accurate location estimates with only a few phases per event. We contrast two different approaches: 1) a calibration-based approach, where we achieve improved locations by using path corrections, and 2) a model-based approach, where we achieve improved locations by relocating in a recently published global 3-D P-velocity model. Both methods result in large improvements in accuracy of the obtained absolute locations, compared to locations obtained in a 1-D reference earth model (ak135). The calibration-based approach gives superior results for this test site, in particular when arrival times from regional stations are included. Estimated locations remain well within a 1000 km² region surrounding the ground truth locations when the path corrections for the Balapan and Degelen Mountain subregions are interchanged, but even for the short separation between these two regions, we find variations in the path corrections that cause systematic mislocations. The model-based approach also results in substantially reduced mislocation distances and has the distinct advantage that it is, in principle, transportable to other source regions around the world.     

Locating microseismic events using borehole data

Constraining microseismic hypocentres in and around hydrocarbon reservoirs and their overburdens is essential for the monitoring of deformation related to hydraulic fracturing, production and injection and the assessment of reservoir security for CO₂ and wastewater storage. Microseismic monitoring in hydrocarbon reservoirs can be achieved via a variety of surface and subsurface acquisition geometries. In this study we use data from a single, subsurface, vertical array of sensors. We test an existing technique that uses a 1D velocity model to constrain locations by minimizing differential S‐to‐P arrival times for individual sensors. We show that small errors in either arrival time picks or the velocity model can lead to large errors in depth, especially near velocity model discontinuities where events tend to cluster. To address this issue we develop two methods that use all available arrival times simultaneously in the inversion, thus maximizing the number of potential constraints from to N, where N is the number of phase picks. The first approach minimizes all available arrival time pairs whilst the second approach, the equal distance time (EDT) method defines the hypocentre as the point where the maximum number of arrival time surfaces intersect. We test and compare the new location procedures with locations using differential S‐to‐P times at each individual sensor on a microseismic data set recorded by a vertical array of sensors at the Ekofisk reservoir in the North Sea. Specifically, we test each procedure's sensitivity to perturbations in measured arrival times and the velocity model using Monte Carlo analysis. In general, location uncertainties increase with increasing raypath length. We show that errors in velocity model estimates are the most significant source of uncertainty in source location with these experiments. Our tests show that hypocentres determined by the new procedures are less sensitive to erroneous measurements and velocity model uncertainties thus reducing the potential for misinterpretation of the results.     

Improved microseismic event location by inclusion of a priori dip particle motion: a case study from Ekofisk

Microseismic monitoring in petroleum settings provides insights into induced and naturally occurring stress changes. Such data are commonly acquired using an array of sensors in a borehole, providing measures of arrival times and polarizations. Events are located using 1D velocity models, P‐ and S‐wave arrival times and the azimuths of P‐wave particle motions. However in the case of all the sensors being deployed in a vertical or near‐vertical borehole, such analysis leads to an inherent 180° ambiguity in the source location. Here we present a location procedure that removes this ambiguity by using the dip of the particle motion as an a priori information to constrain the initial source location. The new procedure is demonstrated with a dataset acquired during hydraulic fracture stimulation, where we know which side of the monitoring well the events are located. Using a 5‐step location procedure, we then reinvestigate a microseismic data set acquired in April 1997 at the Ekofisk oilfield in the North Sea. Traveltimes for 2683 candidate events are manually picked. A noise‐weighted analytic‐signal polarization analysis is used to estimate the dip and azimuth of P‐wave particle motions. A modified t‐test is used to statistically assess the reliability of event location. As a result, 1462 events are located but 627 are deemed to be statistically reliable. The application of a hierarchal cluster analysis highlights coherent structures that cluster around wells and inferred faults. Most events cluster at a depth of roughly 3km in the Ekofisk chalk formation but very little seismicity is observed from the underlying Tor chalk formation, which is separated from the Ekofisk formation by an impermeable layer. We see no evidence for seismicity in the overburden but such events may be too distant to detect. The resulting picture of microseismicity at Ekofisk is very different from those presented in previous studies.     

Reflection-constrained 2D and 3D non-hyperbolic moveout analysis using particle swarm optimization

The objective of moveout parameter inversion is to derive sets of parameter models that can be used for moveout correction and stacking at each common midpoint location to increase the signal-to-noise ratio of the data and to provide insights into the kinematic characteristics of the data amongst other things. In this paper, we introduce a data-driven user-constrained optimization scheme that utilizes manual picks at a point on each reflector within a common midpoint gather to constrain the search space in which an optimization procedure can search for the optimal parameter sets at each reflection. The picks are used to create boundary curves which can be derived approximately via an optimization technique or analytically via the derivation of an analytical bounds function. In this paper, we derive analytical forms of bounds functions for four different moveout cases. These are normal moveout, non-hyperbolic moveout, azimuthally dependent normal moveout and azimuthally dependent non-hyperbolic moveout. The optimization procedure utilized here to search for the optimal moveout parameters is the particle swarm optimization technique. However, any metaheuristic optimization procedure could be modified to account for the constraints introduced in this paper. The technique is tested on two-layer synthetic models based on three of the four moveout cases discussed in this paper. It is also applied to an elastic forward modelled synthetic model called the HESS model, and finally to real 2D land data from Alaska. The resultant stacks show a marked improvement in the signal-to-noise ratio compared to the raw stacks. The results for the normal moveout, non-hyperbolic moveout and azimuthally dependent normal moveout tests suggest that the method is viable for said models. Results demonstrate that our method offers potential as an alternative to conventional parameter picking and inversion schemes, particularly for some cases where the number of parameters in the moveout approximation is 2 or greater.     

微地震事件初至拾取SLPEA算法

微地震事件初至拾取是微地震数据处理的关键步骤之一.实际微地震监测资料中存在大量低信噪比事件,而传统方法对这些事件的应用效果并不理想.为了克服传统方法抗噪性弱的缺点,本文通过综合地震信号与环境噪声在振幅、偏振以及统计特征等方面的存在的差异,设计了一种针对低信噪比微地震事件的初至拾取方法——SLPEA算法.为了检验本文方法的可行性和有效性,分别对模型数据和实际资料进行了处理,并将处理结果与传统方法及手工拾取的结果进行了对比.分析表明,利用本文方法得到的初至到时与手工拾取结果的绝对误差平均值仅为1.33×10~(-3)s,小于3个采样点;方差为3.21×10~(-6)s~2;初至到时在手工拾取结果±0.005s误差范围内的个数占总数的95.8%.这些参数值均优于传统方法的同类参数,证明了本文方法的可靠性.     

STA/LTA—AIC算法对地震P波震相拾取稳定性影响

选取区域地震台网记录的地震波形数据,使用STA/LTA算法与STA/LTA—AIC算法,进行地震P波震相初至到时自动拾取,对地方震及震中距较大的震相进行P波震相拾取效果分析,发现:STA/LTA算法对于地方震P波震相识别精度较高,与STA/LTA—AIC算法拾取的P波震相初至到时相差不大;震中距变大后,STA/LTA算法对P波拾取位置相对于最佳位置向后延迟,STA/LTA—AIC算法有效矫正了STA/LTA算法拾取位置的延迟问题,与人工拾取位置差别可忽略不计。     

基于三分向台站波形的重复地下爆炸相关检测

在将相关检测方法应用于三分向台站记录数据时不能采用台阵数据检测时所使用的基于相关系数束慢度估计的虚假触发筛查方法来控制误检测.为此,本文根据重复事件的震中位置本身固定,各台站记录到的重复事件信号之间的到时差与相应参考事件的信号到时差基本相同的特性,根据两个以上三分向台站的相关检测触发到时差筛查相关检测虚假触发,从而解决了相关检测方法在应用于三分向台站数据时虚假触发过多的问题.利用新疆的三个三分向台站一个月的连续数据对该方法进行测试的结果表明,该方法能在接近零误检率的情况下对重复地下爆炸进行检测.     

IDENTIFICATION OF SEISMIC REFLECTIONS USING SINGULAR VALUE DECOMPOSITION*

Singular value decomposition (SVD) is applied to the identification of seismic reflections by using two different models: the impulse response model, where a seismic trace is assumed to consist of a known signal pulse convolved with a reflection coefficient series plus noise, and the delayed pulse model, where the seismic signal is assumed to consist of a small number of delayed pulses of known shape and with unknown amplitudes and arrival times. SVD clearly shows how least-squares estimation of the reflection coefficients may become unstable, since a division by the singular values is required. Two methods for stabilizing this procedure are investigated. The inverse of the singular values may be replaced by zeros when they are less than a given threshold. This is called the SVD cut-off method. Alternatively, we may use ridge regression which in filter design corresponds to assuming white noise. Statistical methods are used to compute an optimal SVD cut-off level and also to compute an optimal weighting parameter in ridge regression. Numerical studies indicate that the use of SVD cut-off or ridge regression stabilizes the least-squares procedure, but that the results are inferior to maximum-likelihood estimation where the noise is assumed to be filtered white noise. For the delayed pulse model, we use a linearization procedure to iteratively update the estimates of both the reflection amplitudes and the arrival times. In each step, the optimal SVD cut-off method is used. Confidence regions for the estimated reflection amplitudes and arrival times are also computed. Synthetic data examples demonstrate the effectiveness of this method. In a real data example, the maximum-likelihood method assuming an impulse response model is first used to obtain initial estimates of the number of reflections and their amplitudes and traveltimes. Then the iterative procedure is used to obtain improved estimates of the reflection amplitudes and traveltimes.     

Experiments on Exploration of Shallow Fine Structures and the Construction of the 1-D Velocity Model in the Pingtan Island, Fujian

112 short-period seismographs were set up in the 400km² area of Pingtan Island and its surrounding areas in Fujian. The combined observations of the airgun source and ambient noise source were carried out using a dense array to receive the 387 airgun signals excited around the island and one month of continuous ambient noise recording. The 1-D P-wave and S-wave shallow velocity model of Pingtan Island is obtained by the inversion of the airgun body wave''s first arrival time data, and the reliability of the velocity model is verified by using the surface wave phase velocity dispersion curve, which can provide initial model for subsequent 3-D imaging. The experimental results show that this experiment is a successful demonstration of local scale green non-destructive detection, which can provide basic data for shallow surface structure research and strong vibration simulation of the Pingtan Island.     

Detection and Analysis of Near-Surface Explosions on the Kola Peninsula

Seismic and infrasonic observations of signals from a sequence of near-surface explosions at a site on the Kola Peninsula have been analyzed. NORSAR’s automatic network processing of these events shows a significant scatter in the location estimates and, to improve the automatic classification of the events, we have performed full waveform cross-correlation on the data set. Although the signals from the different events share many characteristics, the waveforms do not exhibit a ripple-for-ripple correspondence and cross-correlation does not result in the classic delta-function indicative of repeating signals. Using recordings from the ARCES seismic array (250 km W of the events), we find that a correlation detector on a single channel or three-component station would not be able to detect subsequent events from this source without an unacceptable false alarm rate. However, performing the correlation on each channel of the full ARCES array, and stacking the resulting traces, generates a correlation detection statistic with a suppressed background level which is exceeded by many times its standard deviation on only very few occasions. Performing f-k analysis on the individual correlation coefficient traces, and rejecting detections indicating a non-zero slowness vector, results in a detection list with essentially no false alarms. Applying the algorithm to 8 years of continuous ARCES data identified over 350 events which we confidently assign to this sequence. The large event population provides additional confidence in relative travel-time estimates and this, together with the occurrence of many events between 2002 and 2004 when a temporary network was deployed in the region, reduces the variability in location estimates. The best seismic location estimate, incorporating phase information for many hundreds of events, is consistent with backazimuth measurements for infrasound arrivals at several stations at regional distances. At Lycksele, 800 km SW of the events, as well as at ARCES, infrasound is detected for most of the events in the summer and for few in the winter. At Apatity, some 230 km S of the estimated source location, infrasound is detected for most events. As a first step to providing a Ground Truth database for this useful source of infrasound, we provide the times of explosions for over 50 events spanning 1 year.     

东祁连山北缘断裂带基于深度学习的密集台阵地震事件快速检测与定位研究

为监测东祁连山北缘断裂带附近的地震活动性,布设包含240台短周期地震仪的面状密集台阵,进行约30 d的连续观测。首先使用基于深度学习的多台站地震事件检测算法（CNNDetector）进行地震事件检测,然后使用震相拾取网络（PhaseNet）对地震事件进行P波和S波到时拾取,其次使用震相关联算法（REAL）进行震相关联及初定位,最后使用双差定位（hypoDD）进行地震重定位,最终的精定位地震目录中共有517个地震。在密集台阵观测期间,中国地震台网正式地震目录中共有39个位于台阵内的地震事件,相比而言,密集台阵检测到大量小于0级的地震。因此通过布设密集台阵,可提高活动断裂微地震活动性的监测能力。与历史地震空间分布相比,密集台阵地震精定位分布具有较好的一致性,表现出更明显的线性分布特征。基于地震分布,发现研究区域存在与地表断层迹线走向不同的隐伏活跃断裂。     

PRACTICAL ASPECTS OF THE DETERMINATION OF 3-D STACKING VELOCITIES*

Proper stacking of three-dimensional seismic CDP-data generally requires the knowledge of normal moveout velocities in all source-receiver directions contributing to a CDP-gather. The azimuthal variation of the stacking velocities mainly depends on the dip of the seismic interfaces. For a single dipping plane a simple relation exists between the dip and the azimuthal variation of NMO-velocity. Varying strike and dip of subsequent reflectors, however, result in a complex dependency of the seismic parameters. Reliable information on the spatial distribution of the normal moveout (NMO)-velocity can be derived from a wavefront curvature estimation using a 3-D ray-tracing technique. These procedures require additional information, e.g. reflection time gradients or depth maps to show interval velocities between leading interfaces. Moreover, their application to an extended 3-D data volume is restricted by high costs. The need for a routine 3-D procedure resulted in a special data selection to create pseudo 2-D profiles and to apply existing velocity estimation routines to these profiles. At least three estimates in different directions are necessary to derive the full azimuthal velocity variation, characterized by the large and the small main axis and the orientation of the velocity ellipse. Errors are estimated by means of computer models. Stacking velocities obtained by mathematical routines (least-squares fit) and by seismic standard routines (NMO-correction and correlation) are compared. Finally, a general 3-D velocity procedure using cross-correlation of preliminarily NMO-corrected traces is proposed.     

3种远震P波到时拾取方法的比较及其参数优化

分析了自回归赤池信息准则（AR-AIC）、高阶统计量（HOS）和累积和（CUSUM）等3种到时拾取方法中参数对远震P波到时估计的影响,以450个远震P波信号为样本集,参考人工拾取到时,以网格搜索方式确定了每种方法的最优参数。之后重新选取信噪比处于[2,20]区间的100个远震P波信号,用确定最优参数后的3种到时拾取方法估计其P波初至时间,并比较了3种方法对低信噪比远震P波的拾取准确度。结果表明,AR-AIC方法和CUSUM方法对低信噪比远震P波的拾取准确度要优于HOS方法,CUSUM方法的计算速度最快,HOS方法由于其原理的限制更适用于信噪比较大、初动较尖锐的信号。      

区域震相初至估计

本在地震数据自动化处理中，给出一种基于自回归模型的Akaike information criteria(AIC)算法和信号平均幅值比的混合方法来估计地震信号的初至．用信号的AIC曲线和平均幅值比曲线构造一种叠加曲线，再进行类似于坐标旋转的校正，可以准确估计低信噪比记录中信号的初至，尤其对于震相类型比较复杂的后续震相(如S波、Lg波)的初至估计结果很好．通过对中国数字地震台网乌鲁木齐台记录到的23次天然地震中P波、S波和Lg波的初至估计，与人工分析结果相比，P波初至估计的均方误差为0．71s，后续震相(S波、Lg波)的均方误差为1．64s，优于传统AIC算法的估计结果．     

MAXIMUM LIKELIHOOD ESTIMATION OF SEISMIC REFLECTION COEFFICIENTS*

A seismic trace is modeled as a moving average (MA) process both in signal and noise: a signal wavelet convolved with a reflection coefficient series plus colored random noise. Seismic reflection coefficients can be estimated from seismic traces using suitable estimation algorithms if the input wavelet is known and vice versa. The maximum likelihood (ML) algorithm is used to estimate the system order and the reflection coefficients. The system order is related to the arrival time of the latest signal in a complex seismic reflection event. The least-squares (LS) method does not provide such information. The ML algorithm makes assumptions only about the Gaussian nature of the noise. It is better suited for seismic applications since the LS method inherits the white noise assumption. The Gauss-Newton (G-N) and Newton-Raphson (N-R) optimization algorithms are used to obtain the ML and the LS estimates. Reflection coefficient estimations are affected by the choice of sampling rate of seismic data. Theoretically, the optimum choice in system identification is the Nyquist rate. Experience with synthetic data confirms the theory. In practice, good estimates of reflection coefficients are possible only up to certain pulse separations (or, equivalently, orders). This is mostly due to numerical problems with the optimization algorithms used and partly due to the limited bandwidth of seismic signals. Good estimates from data simulated using three airgun array pulses recorded with 6–128 Hz filter setting are possible up to about 40.0 ms pulse separations. Successful estimations from pinchout and thin layer simulations and well controlled offshore “bright-spots” are given.     

Ground truth seismic events and location capability at Degelen mountain, Kazakhstan

We utilized nuclear explosions from the Degelen Mountain sub-region of the Semipalatinsk Test Site (STS), Kazakhstan, to assess seismic location capability directly. Excellent ground truth information for these events was either known or was estimated from maps of the Degelen Mountain adit complex. Origin times were refined for events for which absolute origin time information was unknown using catalog arrival times, our ground truth location estimates, and a time baseline provided by fixing known origin times during a joint hypocenter determination (JHD). Precise arrival time picks were determined using a waveform cross-correlation process applied to the available digital data. These data were used in a JHD analysis. We found that very accurate locations were possible when high precision, waveform cross-correlation arrival times were combined with JHD. Relocation with our full digital data set resulted in a mean mislocation of 2 km and a mean 95% confidence ellipse (CE) area of 6.6 km² (90% CE: 5.1 km²), however, only 5 of the 18 computed error ellipses actually covered the associated ground truth location estimate. To test a more realistic nuclear test monitoring scenario, we applied our JHD analysis to a set of seven events (one fixed) using data only from seismic stations within 40° epicentral distance. Relocation with these data resulted in a mean mislocation of 7.4 km, with four of the 95% error ellipses covering less than 570 km² (90% CE: 438 km²), and the other two covering 1730 and 8869 km² (90% CE: 1331 and 6822 km²). Location uncertainties calculated using JHD often underestimated the true error, but a circular region with a radius equal to the mislocation covered less than 1000 km² for all events having more than three observations.