Seismic characterization of vertical fractures described as general linear-slip interfaces

Fluid flow in many hydrocarbon reservoirs is controlled by aligned fractures which make the medium anisotropic on the scale of seismic wavelength. Applying the linear‐slip theory, we investigate seismic signatures of the effective medium produced by a single set of ‘general’ vertical fractures embedded in a purely isotropic host rock. The generality of our fracture model means the allowance for coupling between the normal (to the fracture plane) stress and the tangential jump in displacement (and vice versa). Despite its low (triclinic) symmetry, the medium is described by just nine independent effective parameters and possesses several distinct features which help to identify the physical model and estimate the fracture compliances and background velocities. For example, the polarization vector of the vertically propagating fast shear wave S₁ and the semi‐major axis of the S₁‐wave normal‐moveout (NMO) ellipse from a horizontal reflector always point in the direction of the fracture strike. Moreover, for the S₁‐wave both the vertical velocity and the NMO velocity along the fractures are equal to the shear‐wave velocity in the host rock. Analysis of seismic signatures in the limit of small fracture weaknesses allows us to select the input data needed for unambiguous fracture characterization. The fracture and background parameters can be estimated using the NMO ellipses from horizontal reflectors and vertical velocities of P‐waves and two split S‐waves, combined with a portion of the P‐wave slowness surface reconstructed from multi‐azimuth walkaway vertical seismic profiling (VSP) data. The stability of the parameter‐estimation procedure is verified by performing non‐linear inversion based on the exact equations.     

Evaluation of combination rules for peak response calculation in three‐component seismic analysis

The responses, r_e, given by several multicomponent combination rules used in seismic codes for determining peak responses to three ground motion components are evaluated for elastic systems and compared with the critical response r_cr; this is defined as the largest response for all possible incident angles of the seismic components and obtained by means of the CQC3‐rule when a principal seismic component is vertical, or the GCQC3‐rule when it departs from the vertical direction. The combination rules examined are the SRSS‐, 30%‐, 40%‐ and IBC‐rules, considering different alternatives for the design horizontal spectrum. Assuming that a principal seismic component is along the vertical direction, the upper and lower bounds of the ratio r_e/r_cr for each combination rule are determined as a function of the spectral intensity ratio of the horizontal seismic components and of the responses to one seismic component acting alternately along each structural axis. Underestimations and overestimations of the critical response are identified for each combination rule and each design spectrum. When a component departs from the vertical direction, the envelopes of the bounds of the ratio r_e/r_cr for each combination rule are calculated, considering all possible values of the spectral intensity ratios. It is shown that the inclination of a principal component with respect to the vertical axis can significantly reduce the values of r_e/r_cr with respect to the case when the component is vertical. Copyright © 2003 John Wiley & Sons, Ltd.     

THE MAGNETIC ANOMALY OF 3D SOURCES HAVING ARBITRARY GEOMETRY AND MAGNETIZATION

Magnetic anomalies of complicated 3D sources can be calculated by using a combination of analytical and numerical integration. Two surfaces and the magnetization parameters (the amplitudes of the induced and remanent components and the direction cosines) of the source can be defined by arbitrary functions or by discrete data points in a plane. When combined with a polynomial magnetization function in the direction of the third axis, 3D magnetization distribution can also be modelled. The method gives very general equations for anomaly calculation. It can be used for direct modelling of sources interpreted by seismic or other methods and also for interactive interpretation with fast computers. It is possible to calculate anomalies of, for example, intrusives or folded sedimentary beds whose surfaces are functions of horizontal coordinates and which have polynomial magnetization variations in the vertical direction due to gravitational differentiation and arbitrarily varying magnetization in the horizontal direction due to regional metamorphosis. If the distribution of magnetization parameters in the vertical direction cannot be described satisfactorily by polynomials, models can be used whose surfaces are functions of the vertical coordinate and which can then have any arbitrary magnetization distribution in the vertical direction.     

Migration velocity analysis for tilted transversely isotropic media

Tilted transversely isotropic formations cause serious imaging distortions in active tectonic areas (e.g., fold‐and‐thrust belts) and in subsalt exploration. Here, we introduce a methodology for P‐wave prestack depth imaging in tilted transversely isotropic media that properly accounts for the tilt of the symmetry axis as well as for spatial velocity variations. For purposes of migration velocity analysis, the model is divided into blocks with constant values of the anisotropy parameters ε and δ and linearly varying symmetry‐direction velocity V_P0 controlled by the vertical (k_z) and lateral (k_x) gradients. Since determination of tilt from P‐wave data is generally unstable, the symmetry axis is kept orthogonal to the reflectors in all trial velocity models. It is also assumed that the velocity V_P0 is either known at the top of each block or remains continuous in the vertical direction. The velocity analysis algorithm estimates the velocity gradients k_z and k_x and the anisotropy parameters ε and δ in the layer‐stripping mode using a generalized version of the method introduced by Sarkar and Tsvankin for factorized transverse isotropy with a vertical symmetry axis. Synthetic tests for several models typical in exploration (a syncline, uptilted shale layers near a salt dome and a bending shale layer) confirm that if the symmetry‐axis direction is fixed and V_P0 is known, the parameters k_z, k_x, ε and δ can be resolved from reflection data. It should be emphasized that estimation of ε in tilted transversely isotropic media requires using nonhyperbolic moveout for long offsets reaching at least twice the reflector depth. We also demonstrate that application of processing algorithms designed for a vertical symmetry axis to data from tilted transversely isotropic media may lead to significant misfocusing of reflectors and errors in parameter estimation, even when the tilt is moderate (30°). The ability of our velocity analysis algorithm to separate the anisotropy parameters from the velocity gradients can be also used in lithology discrimination and geologic interpretation of seismic data in complex areas.     

Diffraction imaging by prestack reverse‐time migration in the dip‐angle domain

Prestack image volumes may be decomposed into specular and non‐specular parts by filters defined in the dip‐angle domain. For space‐shift extended image volumes, the dip‐angle decomposition is derived via local Radon transform in depth and midpoint coordinates, followed by an averaging over space‐shifts. We propose to employ prestack space‐shift extended reverse‐time migration and dip‐angle decomposition for imaging small‐scale structural elements, considered as seismic diffractors, in models with arbitrary complexity. A suitable design of a specularity filter in the dip‐angle domain rejects the dominant reflectors and enhances diffractors and other non‐specular image content. The filter exploits a clear discrimination in dip between specular reflections and diffractions. The former are stationary at the specular dip, whereas the latter are non‐stationary without a preferred dip direction. While the filtered image volume features other than the diffractor images (for example, noise and truncation artefacts are also present), synthetic and field data examples suggest that diffractors tend to dominate and are readily recognisable. Averaging over space‐shifts in the filter construction makes the reflectors? rejection robust against migration velocity errors. Another consequence of the space‐shift extension and its angle‐domain transforms is the possibility of exploring the image in a multiple set of common‐image gathers. The filtered diffractions may be analysed simultaneously in space‐shift, scattering‐angle, and dip‐angle image gathers by means of a single migration job. The deliverables of our method obviously enrich the processed material on the interpreter's desk. We expect them to further supplement our understanding of the Earth's interior.     

Logarithm of short‐time Fourier transform for extending the seismic bandwidth

Improving seismic resolution is essential for obtaining more detailed structural and stratigraphic information. We present a new algorithm to increase seismic resolution with a minimum of user‐defined parameters. The algorithm inherits useful properties of both the short‐time Fourier transform and the cepstrum to smooth and broaden the frequency spectrum at each translation of the spectral decomposing window. The key idea is to replace the amplitude spectrum with its logarithm in each window of the short‐time Fourier transform. We describe the mathematical formulation of the algorithm and its testing on synthetic and real seismic data to obtain broader frequency spectra and thus enhance the seismic resolution.     

SEISMIC HORIZON DETECTION USING IMAGE PROCESSING ALGORITHMS1

Seismic sections used in interpretation are actually images. We often superimpose colour-coded pictures of seismic attributes on seismic sections. Thus, it seems straightforward to use image processing algorithms to enhance the quality of the seismic images. From an image processing point of view, seismic horizons can be thought of as edges on the seismic image. We present a novel approach to detecting seismic horizons, which includes the 2D median filtering of the instantaneous phase attribute and applying an edge detection algorithm. The resulting edge magnitude picture provides a skeletonized image of the seismic section, on which the structural and stratigraphic patterns can be better recognized.     

场地对地震动放大作用的方向相关性研究

现有场地对地震动的放大作用分析中,常采用观测得到正交方向的HVSR和SBSR的矢量合成结果表征场地放大作用,而不考虑放大作用的方向相关性。为研究不同地震作用下场地放大作用的方向相关性,采用多向HVSR和多向SBSR及谱比峰值方向相关分布概率,以11次地震事件获得的地表、井下基岩记录作为数据,开展响嘡场地效应台阵场地对地震动放大作用的方向相关性分析。研究结果表明,常用的矢量合成法分析结果与多向谱比上包络曲线法分析结果存在差异,前者存在忽略部分峰值的可能;多向谱比的极坐标等高线图可较直观地体现场地对地震动放大作用的方向相关性特征,但不同多向谱比表现存在差异;多向谱比峰值方向分布概率可较好地定量分析放大作用峰值在不同方向上的差异,但不同地震事件存在差异,不同的多向谱比存在差异。可采用多向谱比上包络曲线法分析场地放大作用,除可避免矢量合成法引起的误估外,也可更好地分析实际场地频率-方向相关的放大作用。     

A short-spread configuration for mapping dipping horizons in reflection seismic surveys

Summary In areas of steeply dipping seismic horizons, the conventional split-spread seismic reflection method may yield poor seismograms. In order to develop a suitable field technique to improve the record quality in such areas, an analysis of the effect of spherical divergence, reflection coefficient, and absorption on the amplitude of the signal reflected from a dipping interface has been carried out. It has been found that for a dipping interface the reflected signal amplitude in many cases may be higher for the larger angles of incidence at the reflecting bed than for the normal incidence. Thus, a suitable inline-offset spread in the updip direction is likely to increase the recorded reflected signal amplidute and improve the signal-to-noisre ratio in areas of dipping reflectors. The design of such a shot-spread configuration has been developed in this paper, and a field method using this configuration is suggested for continuous mapping of subsurface horizons.     

谱分解技术在储层预测中的应用

目的：在传统的地震资料解释中存在着复杂地质体的反射层位准确拾取难、沿单一解释层位提取属性不客观等不足。文中介绍一种适用于三维地震数据体的短时窗谱分解和调谐三维体技术。方法：该技术针对薄层单元地质体计算属性，它利用有限带宽地震子波内各种频率分量特性。结果：精细刻划复杂地质体内部地层反射特征，从而客观揭示复杂储层的纵横向变化趋势。结论：以缝洞型碳酸盐岩溶洞实例结果说明方法的有效性。     

球坐标系下谱元法三维地震波场模拟

谱元法已成为区域性乃至大陆性尺度地震波场模拟的重要工具.对于区域或大陆尺度层析成像而言,地球曲率不可忽略,此时模拟地震波传播采用球坐标系更为合适.本文从球坐标系下弹性波动方程弱形式出发,基于球坐标系变分原理给出了球坐标系下求解三维地震波方程的谱元法.另一方面,计算Fréchet敏感核是进行全波形反演的关键,本文借助伴随原理,推导了全波走时层析成像三维Fréchet敏感核表达式.为了验证球坐标系下谱元法的精度,我们将数值模拟结果与normal mode方法得到的解析解在1-D PREM模型下进行了对比.同时,我们将此方法应用到华北克拉通区域,以期获得地球内部结构精确成像.基于3-D全球径向各向异性地幔模型S362ANI和3-D地壳模型Crust1.0,我们建立了华北克拉通初始3-D背景模型,并将数值模拟结果与实际观测台站记录波形资料进行对比分析,利用互相关方法提取走时残差,最后给出了Fréchet敏感核在3-D空间中的分布,这些工作为下一步开展球坐标系下三维大尺度全波形反演奠定了基础.     

The rotated Cartesian coordinate method to remove the axial singularity of cylindrical coordinates in finite‐difference schemes for elastic and viscoelastic waves

When modelling the propagation of 3D non‐axisymmetric elastic and viscoelastic waves in cylindrical coordinates using the finite‐difference time‐domain method, a mathematical singularity occurs due to the presence of terms in the elastic and viscoelastic wave equations. For many years, this issue has been impeding the accurate numerical solution near the axis. In this work, we propose a simple but effective method for the treatment of this numerical singularity problem. By rotating the Cartesian coordinate system around the z‐axis in cylindrical coordinates, the numerical singularity problems in both 2D and 3D cylindrical coordinates can be removed. This algorithm has three advantages over the conventional treatment techniques: (i) the excitation source can be directly loaded at , (ii) the central difference scheme with second‐order accuracy is maintained, and (iii) the stability condition at the axis is consistent with the finite‐difference time‐domain in Cartesian coordinates. This method is verified by several 3D numerical examples. Results show that the rotating the Cartesian coordinate method is accurate and stable at the singularity axis. The improved finite‐difference time‐domain algorithm is also applied to sonic logging simulations in non‐axisymmetric formations and sources.     

Characterization of dipping fractures in a transverselyisotropic background

Although it is believed that natural fracture sets predominantly have near‐vertical orientation, oblique stresses and some other mechanisms may tilt fractures away from the vertical. Here, we examine an effective medium produced by a single system of obliquely dipping rotationally invariant fractures embedded in a transversely isotropic with a vertical symmetry axis (VTI) background rock. This model is monoclinic with a vertical symmetry plane that coincides with the dip plane of the fractures. Multicomponent seismic data acquired over such a medium possess several distinct features that make it possible to estimate the fracture orientation. For example, the vertically propagating fast shear wave (and the fast converted PS‐wave) is typically polarized in the direction of the fracture strike. The normal‐moveout (NMO) ellipses of horizontal reflection events are co‐orientated with the dip and strike directions of the fractures, which provides an independent estimate of the fracture azimuth. However, the polarization vector of the slow shear wave at vertical incidence does not lie in the horizontal plane – an unusual phenomenon that can be used to evaluate fracture dip. Also, for oblique fractures the shear‐wave splitting coefficient at vertical incidence becomes dependent on fracture infill (saturation). A complete medium‐characterization procedure includes estimating the fracture compliances and orientation (dip and azimuth), as well as the Thomsen parameters of the VTI background. We demonstrate that both the fracture and background parameters can be obtained from multicomponent wide‐azimuth data using the vertical velocities and NMO ellipses of PP‐waves and two split SS‐waves (or the traveltimes of PS‐waves) reflected from horizontal interfaces. Numerical tests corroborate the accuracy and stability of the inversion algorithm based on the exact expressions for the vertical and NMO velocities.     

Dreamlet source‐receiver survey sinking prestack depth migration

Survey sinking migration downward continues the entire surface observed multi‐shot data to the subsurface step by step recursively. Reflected energy from reflectors at current depth appear at zero time and zero offset in the extrapolated wavefield. The data (seismic records) of t > 0 at this depth are equivalent to the data acquired by a survey system deployed at this depth. This is the reason to name the process ‘survey sinking’. The records of negative time need not to be further propagated since they carry no information to image structures beneath the new survey system. In this paper, we combine survey sinking with dreamlet migration. The dreamlet migration method decomposes the seismic wavefield and one‐way wave propagator by complete time‐space localized bases. The localization on time gives flexibility on time‐varying operations during depth extrapolation. In dreamlet survey sinking migration, it only keeps the data for imaging the structures beneath the sunk survey system and gets rid of the data already used to image structures above it. The deeper the depth is, the shorter is the valid time records of the remaining data and less computation is needed for one depth step continuation. For data decomposition, in addition to time axis, dreamlet survey sinking also decomposes the data for source and receiver gathers, which is a fully localized decomposition of prestack seismic data. A three‐scatter model is first used to demonstrate the computational feature and principle of this method. Tests on the two‐dimensional SEG/EAGE salt model show that with reduced data sets the proposed method can still obtain good imaging quality on complex geology structures and a strong velocity contrast environment.     

Determinations of directions of the mean stress field in Sichuan-Yunnan region from a number of focal mechanism solutions

Introduction With the model of double couple point source, both orientations of the T, B, and P axes of a single focal mechanism solution and parameters of the two nodal planes may be determined. The focal mechanism solution under ″the model of the best double couple″ can also be given by seis- mic moment tensor, such as those provided by Harvard CMT solutions. The orientations of the T, B, and P axes from one earthquake are associated with the released stress, but they cannot be re- gar…     

由多个地震震源机制解求川滇地区平均应力场方向

基于震源断层面解的空间取向和断层滑动方向,写出相应力轴张量在地理坐标系中的表达式,进而给出计算平均力轴张量及主值的方法,即通过求解相应的本征方程得到. 对使用多个震源机制解的T,B,P轴参数计算平均应力场的方法,以用滑动方向拟合法反演富蕴、唐山地区平均应力场数据进行验证,二者结果一致. 选择具有地震构造意义的地块或地震带内大量地震的震源机制解研究区域平均应力场. 根据川滇13个地震带（区）的256次中强地震的震源机制解,给出了各带（区）应力张量的定量分析结果. 该方法算法简便,使由大量地震震源机制解资料分析构造应力场的方法走向定量化.      

多尺度快速匹配追踪多域联合地震反演方法

多尺度快速匹配追踪多域联合地震反演是一种通过地震数据多尺度分解的迭代反演方法.与此同时,在快速匹配追踪算法中引入低频模型约束,有效提高了收敛精度,使反演结果具有丰富的高低频信息.首先通过对大尺度地震资料进行反演得到低频背景.在此基础上,采用中尺度与小尺度地震数据进行逐级迭代用以获得高频数据,因而有效缓解了常规反演方法对于初始模型精度的依赖.最后利用理论模型与实际地震数据进行测试,通过与常规时间域反演方法的反演结果进行对比可以看出,本文方法在地层连续变化处依然可以对变化地层进行精确刻画,且在纵向分辨率提升的同时保持了较好的横向连续性.     

2D inversion of refraction traveltime curves using homogeneous functions

A method using simple inversion of refraction traveltimes for the determination of 2D velocity and interface structure is presented. The method is applicable to data obtained from engineering seismics and from deep seismic investigations. The advantage of simple inversion, as opposed to ray‐tracing methods, is that it enables direct calculation of a 2D velocity distribution, including information about interfaces, thus eliminating the calculation of seismic rays at every step of the iteration process. The inversion method is based on a local approximation of the real velocity cross‐section by homogeneous functions of two coordinates. Homogeneous functions are very useful for the approximation of real geological media. Homogeneous velocity functions can include straight‐line seismic boundaries. The contour lines of homogeneous functions are arbitrary curves that are similar to one another. The traveltime curves recorded at the surface of media with homogeneous velocity functions are also similar to one another. This is true for both refraction and reflection traveltime curves. For two reverse traveltime curves, non‐linear transformations exist which continuously convert the direct traveltime curve to the reverse one and vice versa. This fact has enabled us to develop an automatic procedure for the identification of waves refracted at different seismic boundaries using reverse traveltime curves. Homogeneous functions of two coordinates can describe media where the velocity depends significantly on two coordinates. However, the rays and the traveltime fields corresponding to these velocity functions can be transformed to those for media where the velocity depends on one coordinate. The 2D inverse kinematic problem, i.e. the computation of an approximate homogeneous velocity function using the data from two reverse traveltime curves of the refracted first arrival, is thus resolved. Since the solution algorithm is stable, in the case of complex shooting geometry, the common‐velocity cross‐section can be constructed by applying a local approximation. This method enables the reconstruction of practically any arbitrary velocity function of two coordinates. The computer program, known as godograf , which is based on this theory, is a universal program for the interpretation of any system of refraction traveltime curves for any refraction method for both shallow and deep seismic studies of crust and mantle. Examples using synthetic data demonstrate the accuracy of the algorithm and its sensitivity to realistic noise levels. Inversions of the refraction traveltimes from the Salair ore deposit, the Moscow region and the Kamchatka volcano seismic profiles illustrate the methodology, practical considerations and capability of seismic imaging with the inversion method.     

利用岩性解释方法圈定岩浆岩侵入煤层范围

岩性地震反演技术是根据钻孔测井数据纵向分辨率很高的有利条件,对井旁地震资料进行约束;并在此基础上对孔间地震资料进行反演,推断煤系地层岩性在平面上的变化情况.谱分解技术是利用薄层调谐体离散频率特性,通过分析复杂岩层内陷频谱变化和局部相位的不稳定性,识别薄地层横向分布特征.本文根据崔庄煤矿三采区的地震资料,综合利用地震反演技术和谱分解技术等岩性解释方法圈定岩浆岩侵入煤层的范围.     

OBTAINING INTERVAL VELOCITIES FROM STACKING VELOCITIES WHEN DIPPING HORIZONS ARE INCLUDED *

Common-depth-point stacking velocities may differ from root-mean-square velocities because of large offset and because of dipping reflectors. This paper shows that the two effects may be treated separately, and proceeds to examine the effect of dip. If stacking velocities are assumed equal to rms velocities for the purpose of time to depth conversion, then errors are introduced comparable to the difference between migrated and unmigrated depths. Consequently, if the effect of dip on stacking velocity is ignored, there is no point in migrating the resulting depth data. For a multi-layered model having parallel dip, a formula is developed to compute interval velocities and depths from the stacking velocities, time picks, and time slope of the seismic section. It is shown that cross-dip need not be considered, if all the reflectors have the same dip azimuth. The problem becomes intractable if the dips are not parallel. But the inverse problem is soluble: to obtain, stacking velocities; time picks, and time slopes from a given depth and interval velocity model. Finally, the inverse solution is combined with an approximate forward solution. This provides an iterative method to obtain depths and interval velocities from stacking velocities, time picks and time slopes. It is assumed that the dip azimuth is the same for all reflectors, but not necessarily in the plane of the section, and that the curvature of the reflecting horizons is negligible. The effect of onset delay is examined. It is shown that onset corrections may be unnecessary when converting from time to depth.