【作者】 巩向博；

【导师】 韩立国；

【作者基本信息】 吉林大学 ， 地球探测与信息技术， 2011， 博士

【摘要】 采用地震方法进行深部金属矿的探测,可以弥补重磁电传统方法探测深度不足的缺点,为第二找矿空间提供了物探技术可能性。但金属矿区地质条件往往复杂多变；地下矿体形态复杂、倾角大、尺度小；矿体与围岩之间的波阻抗差异较小；且矿体内部常呈现非均匀体性质,使得后续地震资料处理面临很多难题。本文针对金属矿地震数据的这些特点,研究了一系列提高地震分辨率和信噪比的处理技术,有些技术是首次提出,有些技术是应用研究,旨在综合运用以提高金属矿地震数据的分辨率和信噪比,从而为我国深部矿产资源的探测与开发提供技术保障。金属矿地震勘探中的分辨率分为垂直分辨率和水平分辨率,一般认为垂直分辨率为调谐厚度——“四分之一波长”,水平分辨率为第一菲涅尔带半径。现代叠前深度偏移技术从波动方程出发,或为积分成像或为延拓成像,可以将菲涅尔带收敛为点,从而大幅提高地震记录中的分辨率,理论上深度偏移成像结果的分辨率仅与空间采样率有关。我们首先对二维复杂Marmousi模型进行了三种叠前深度偏移方法的对比研究,其中有Kirchhoff积分偏移、单程波延拓偏移、逆时偏移。通过对比发现,逆时偏移的成像精度最高,对标准Marmousi模型中的大断裂、小断块、陡倾角断层、高速楔状体、背斜构造、目标储层成像准确清晰,因而适用于大倾角、小尺度、速度陡变的金属矿地质模型。我们对逆时偏移成像技术进行了深入研究,总结了目前逆时偏移面临的问题与难点：1、由于时间一致性成像条件多解性,互相关成像条件将折射波、回折波、逆散射波错误互相关而造成低频假象,恶化了浅层的成像效果；2、高计算量与高内存开销；3、数值计算中面临的所有问题。同时也总结了解决这些问题与难点的方案：1、采用多种成像条件综合应用及偏移后剖面多种滤波方法,以压制低频假象；2、平面波偏移技术,多步偏移策略以及采用高性能科学计算技术,包括多节点并行计算技术、GPU通用计算技术等,以应对高计算量、高内存开销问题；3、采用高阶有限差分求解波动方程、吸收边界等技术提高数值计算精度。我们针对金属矿地震勘探的特点,选用了上述消除假象、提高计算效率的一些技术与方法,主要有高阶差分延拓波场、归一化成像条件、成像结果滤波、高性能并行计算等方法,并进行了金属矿模型的逆时偏移试算。以湖北某铜矿和铜铁矿共生的隐伏矿床模型为例,进行了逆时偏移成像,成像后的结果可以清晰的分辨细长隐伏在矿床中的共生铜矿体,铜铁矿与铜矿体的分界面已经可以明显分辨,许多细节特征都可以被刻画出来。叠前深度偏移技术相对于叠后和时间偏移,对速度的依赖更为严重,速度误差引起的成像误差已经超过偏移算法本身,研究一种效率高、精度高的速度建模方法十分重要。我们研究了上文介绍的三种叠前深度方法对速度的敏感性问题,给出了速度误差引起深度误差的解析式,分析了在深度偏移剖面和偏移后共成像点道集(CIG)出现“微笑”和“皱眉”现象。通过水平层状模型、高速体模型以及复杂Marmousi模型在相同速度误差下的深度误差,CIG同相轴误差及弯曲程度的分析,得出三种叠前深度偏移方法在简单地质模型下对速度敏感性一致,在复杂地质模型下速度敏感性相似；在简单地质构造情况下,偏移速度分析(MVA)中偏移方法的选择对速度分析准确性的影响不大。根据这个结论,我们提出了联合速度建模方法,这种方法集成了Kirchhoff积分偏移计算效率高、逆时偏移计算精度高的优点。建模的流程如下：第一步,常规速度分析方法估计最佳初始速度模型；第二步,初始速度模型下使用Kirchhoff偏移进行成像；第三步,抽取CIG或者ADCIG(角度域CIG),根据CIG的剩余曲率(深度误差)计算速度扰动,进而更新速度模型；第四步,新速度模型下重复二、三步,将主要构造速度界面估计准确,一般需迭代运算2-3次；第五步,使用逆时偏移进行成像；第六步,抽取CIG或者是ADCIG(角度域CIG),根据CIG的剩余曲率(深度误差)计算速度扰动,进而更新速度模型；第七步,新速度模型下重复五、六步得到最终偏移速度模型,需迭代运算1-2次,将小尺度复杂构造速度界面估计准确,输出最终偏移速度模型；第八步,使用逆时偏移在最终偏移速度模型下进行成像,输出偏移结果。以某蘑菇状高速金属矿模型为例,进行了联合速度建模方法试算,采用逐层分析策略,每层选择适量的速度分析点,其余点通过内插和平滑得到,更新速度时采用剩余曲率分析(RCA)方法,根据CIG同相轴的弯曲曲率反算速度误差。经联合速度建模后,有效地得到了高质量金属矿成像剖面及金属矿体速度模型。金属矿地震数据另一特点是信噪比较低,噪声较多且强,有效信号相对较弱,地震数据信噪比低会直接降低地震剖面的分辨率,我们给出了地震数据信噪比与分辨率之间的关系,分析得出金属矿地震数据中保持信噪比大于2是必要的。Radon变换是将地震数据通过一定的路径积分,变换到另一个域中,使得t-x域内耦合的信号和噪声在Radon域内得到分离,经噪声切除后反变换,得到去噪后的地震数据。根据积分路径不同,可将Radon变换分为线性、抛物、多项式、双曲、相移双曲、椭圆Radon变换；其中线性、抛物、多项式Radon变换算子具有时不变特性,可以在频率域求解,这样做的好处是避免了求解t-x域所有慢度采样构成的大矩阵求逆问题,从而提高了计算速度,算子矩阵还具有Toeplitz结构,可以用莱文森递推算法快速求解。基于贝叶斯原理的高分辨率Radon变换技术,运用了频域空间稀疏约束算法进行迭代计算,利用加权矩阵将前一次计算结果联系起来,更新加权矩阵以得到Radon域的高分辨率结果,这种算子矩阵可以使用双共轭梯度算法求解。地下介质往往是各向异性的,在Radon变换的积分路径中考虑各向异性参数,符合真实地下介质模型,对于实际数据处理去噪可以收到更好的效果,我们从非双曲时差公式出发,定义了各向异性Radon正变换公式,相比于常规Radon变换公式,增加了各向异性参数(非椭圆率)η,这个参数可以用Thomson各向异性参数ε、δ表示,重新写成τ关于t的形式,得到各向异性反变换公式。Radon域信噪分离时,直接切除噪声方法容易损失有效信号,采用二维蒙版滤波技术可以自动拾取聚焦能量剪切线外侧的能量残余,利用Radon域截断效应的轨迹拾取剪切线外侧的能量残余,这样可以将干扰噪声更好的分离,使有效波的能量更好的保留。使用加入高斯随机噪声的模拟记录,对常规Radon变换、高分辨率Radon变换和二维蒙版滤波去噪进行了对比分析,可以得出后者有更好的压制随机噪声效果；对某中间放炮记录进行了多种方法对比研究,分别使用了FK域滤波、高分辨率线性Radon变换、高分辨率抛物Radon变换以及二维蒙版滤波去噪,比较得出采用Radon变换结合蒙版滤波去噪的结果最佳；对庐枞实际金属矿地震数据进行了去噪处理,通过去噪结果和差值剖面分析,Radon变换及蒙版滤波技术压制了强面波及强线性干扰波,大幅提高了地震数据的信噪比。对金昌镍铜矿区地震试验数据进行了特色处理研究。我们根据此地区地质调绘资料制作了镍铜矿的地质模型,其中镍铜矿体呈50-60。倾斜层状分布在围岩之中。在此地质模型下使用高阶差分波动方程进行了数值模拟及波场特征分析,倾斜层引起的倾斜同相轴由于层距微小产生相互干涉现象,其中混杂耦合了来自尖灭点、突变点引起的绕射波能量,使得波场特征十分复杂。模拟数据分别进行了常规地震处理和叠前深度偏移处理,前者得到的叠加剖面矿脉被严重拉伸,细长矿体整体向左上方移动,绕射波和散射波仍严重影响叠后的剖面质量,造成解释上的误差；后者得到的偏移图像有很高的分辨率,已将矿脉、矿体位置正确还原,可以清晰分辨矿脉边界,结合地质资料容易进行地质解译。针对金昌镍铜矿区观测系统方式及地震数据特点,我们使用了自主知识产权的处理软件系统及多种自主开发的特色处理模块对实际数据进行了处理,其中包括：高分辨率速度分析技术,以提高速度谱拾取的精度；Gabor变换地表一致性反褶积技术,以消除由于近地表风化严重,地震检波器与基岩耦合不佳的影响,并提高了地震资料的频带；Curvelet域组合变换法压制随机噪声技术,以提高地震记录信噪比；Radon变换及蒙版自适应切除滤波去噪技术,压制强线性干扰、提取弱反射信号。对比常规处理后的叠加剖面,综合使用这些特色处理技术,可以大幅提高金属矿地震数据的分辨率和信噪比,获得高质量的地震叠加剖面,为指示深部金属矿体的延伸和分布形态提供了技术支持。更多还原

【Abstract】 Using seismic method for deep mineral exploration, can make up for the shortcoming of probing depth of the traditional exploration method, such as gravity, magnetic and electrical prospecting, it provides the necessary techniques for the second prospecting space. However, metallic mines geological conditions are often complex and undulate; ore bodies underground present complex structure, steep dip, small scale; the impedance has little difference between the ore body and the surrounding rock; and in the ore body often present the heterogeneous nature. The subsequent seismic data processing faces many problems. In this paper, with these characteristics of metallic mines, we study a series of seismic processing techniques to improve the resolution and signal-to-noise ratio of the data, some of these technologies are first proposed, some of these technologies are applied research. With comprehensive application, it provides the technical support for deep exploration and development of mineral resources.The resolutions of seismic exploration are the vertical resolution and horizontal resolution. It is generally recognized that the vertical resolution is the tuning thickness "quarter wavelength", the horizontal resolution is the first Fresnel zone radius. Pre-stack depth migration focuses the Fresnel zone on a point, and then it greatly improves the resolution of seismic data. The resolutions of pre-stack depth migration only have relation to the spatial sampling rate in theory. We compare three pre-stack migration methods in two-dimensional complex Marmousi model, which are Kirchhoff integral migration, one-way wave-equation migration and reverse-time migration. By comparison, the lines in the image of Kirchhoff integral migration are thicker, target reservoir is not clear, there are some blind areas and caustic areas, and the imaging is not very well; the image of one-way wave equation migration has higher resolution, but the imaging of high dip angle faults and small-scale fault blocks are indistinct; the reverse time migration has the best imaging result and the resolution is high, many features of details are retained, the imaging of three faults, fault blocks, anticline, high-velocity anomaly and target reservoir are clearer, and the locations are also very accurate. If we neglect the imaging accuracy, the locations of target reservoir obtained by three methods are same correspondent with the speed model. We made a further study on imaging technique of reverse time migration and summed up the present problem and difficulty of reverse time migration we stated:1. Due to the multiple solutions of time consistency imaging conditions, cross-correlation imaging condition cause low-frequency false by mistakenly correlating the refraction wave, diving wave, inverse scattering wave with received wavefield,which deteriorates the imaging quality of shallow;2. The cost of computation and memory is high;3; All problems which numerical calculation faces. Also I summarized to solve these problems and difficulties of the program:1. I comprehensively apply a variety of imaging conditions and the various filtering methods of migrated section, in order to suppress low-frequency false;2.Plane wave migration, multi-step migration strategy and the use of high-performance scientific computing, including multi-node parallel computing, GPU general computing technologies are applied in order to cope with the high computation and high memory overhead;3. I use higher order finite difference solving wave equation, absorbing boundary and other technology to improve the calculation precision. Aiming at the characters of metal mine seismic exploration, certain technologies mentioned above are selected to eliminate the pseudomorph and improve calculation efficiency, mainly including high order wavefield continuation、normalized imaging conditions、the filter of imaging results、the methods such as high-performance parallel computing, and reverse-time migration trial calculation of metal mine model was carried out. Taking the concealed mineral deposit model accreted of copper mine and delafossite in Hubei province for example, reverse-time migration imaging formation was done. Slender symbiotic copper ore body concealed in the mineral deposit and the interface of delafossite and copper mine could be clearly figured out by the results of imaging, and many detail characters could be depicted.Pre-stack depth migration respect to post-stack and time migration is more serious dependence on the velocity. The error of migration caused by velocity error has been more than caused by migration algorithm. The study on high efficiency and high precision velocity mode-building method is important. We studied the three pre-stack depth migration methods on velocity sensitive issues, and given analytic equation between the velocity error and the depth error. We analysis the "Smile" and "frown" phenomenon in the depth migrated section and common image gathers (CIG). Through the analysis of horizontal layered model, high-velocity model and Marmousi model, we give below conclusion:the migration depths of three methods with the same velocity perturbation are the same, in other words, they have the same depth error, namely the similar velocity error sensitivity. Consider the imaging effect and computational efficiency, using the advantages of above two imaging methods, we proposed the combined velocity model-building method which based on the Kirchhoff integral migration and reverse time migration, it combines the efficient computing of Kirchhoff integral migration and the high accuracy imaging of reverse time migration. The following eight steps describe the process of combined velocity model-building method:First, Establish the "initial velocity model" by general velocity analysis; Second, Use Kirchhoff integral migration to image section with "initial velocity model"; Third, Form the CIG or ADCIG of the initial velocity model, use the relationship between the imaging depth error and formation parameters error to update the velocity model (RCA strategy layer by layer); Fourth, Repeat steps 2 to 3, the process of updating velocity must be iterative calculation 2-3 times, it will be form "middle velocity model" with major structures, we call it is "part MVA"; Fifth, Use reverse time migration to image section with "middle velocity model"; Sixth, Same to step 3; Seventh, Repeat steps 5 to 6, the process of updating velocity must be iterative calculation 1-2 times, it will be form "final velocity model" with the complex structures of steep dip angle and small scale; Eighth, Use reverse time migration to image section with "final velocity model", as final imaging results. Taking a mushroom-shaped high-speed metal mining model as an example, the velocity analysis using RCA at each layer, we should choose the appropriate amount of control point; the residual points are available through interpolation and smoothing. We use the smooth function between adjacent layers, increase the control point near the complex structure and the acute velocity regions to improve the accuracy of the velocity analysis. The combined velocity modeling can effectively get high quality imaging section and velocity model of metal mine.Another characteristic of mineral mining seismic data is low signal-to-noise rate, the noise is strong, and effective signal is relatively weak. Low signal-to-noise ratio of seismic data will directly reduce the resolution of seismic section, we give the relationship between the resolution and the signal-to-noise ratio, analysis of seismic data obtained in the metallic mine to maintain signal-to-noise ratio greater than 2 is necessary. According path can be integral to Radon transform respectively linear, parabolic, polynomial, hyperbolic, elliptic transform. Phase shift hyperbolic Dix hyperbolic equation is far offset error for the larger substitute. Compared with conventional transformation, far offset hyperbolic reflection energy is better focused. Linear, parabolic, polynomial Radon transform operator has invariance when, in the frequency domain can be achieved. Time-the spatial domain into the frequency of the problem-to solve the spatial domain, the benefits of doing so is the speed of Fourier transform is very fast and can be carried out in each frequency component Radon transform, and avoid solving all the time--offset the speed of the large sample matrix problems, and enhance the computing speed. Toeplitz operator matrix structure of (uniform sampling parameters), can be used Levinsion fast recursive algorithm to solve. High resolution Radon transform is a frequency domain space sparse bound algorithm, in the inversion iteration process, in accordance with the results of the previous iteration, through Bayesian principles will be weighted matrix and the results of the previous iteration link be new weighted matrix; weighted matrix and then solve the equation, the frequency domain by the sparse solution. As operator matrix weighted matrix does not have the presence of Toeplitz structure, but still Hermite matrix, we can use the Cholesky decomposition method, LU decomposition method, the faster algorithm is conjugate gradient algorithm. In foreign, some scholars have begun studying time-space domain to direct calculation precision hyperbolic Radon transform algorithm, the algorithm used in least squares method to the most rapid decline in law. Conjugate gradient method fast convergence, the process for stability and the high accuracy solution, which has greatly improved its calculation speed. Subsurface media is often anisotropic, in order to consistent with the real subsurface media, anisotropy parameters should be considered in the integral path of Radon transform. For actual data processing, de-noising can receive better results, based on non-hyperbolic travel time offset function, anisotropic Radon transform formula has been defined, Compared to general Radon transform formula, anisotropy parametersη(anisotropy parameters) has been added, this parameters can be indicated by anisotropy parameters epsilon and delta, re-writtenτin the form of t, we obtain anisotropic inverse transform formula. We use simulation seismic records which added Gaussian random noise, and compared conventional Radon transform, high-resolution Radon transform and 2D mask filtering de-noising, we concluded that 2D mask filtering de-noising has a better effect of random noise suppression. By a variety of methods of a comparative study, which contains FK filtering,high-resolution linear Radon transform, high-resolution parabolic Radon transform and 2D mask filtering, on a middle shooting seismic record, we concluded that Radon transform combined with mask filtering has a best result of filtering and de-noising. De-noising a actual metal mining seismic data,and comparing de-noising results and difference profile, we concluded that Radon transform and mask filtering can suppress strong surface waves and strong linear interference, and then substantially raise signal to noise ratio of seismic data.We did characteristic processing and research on seismic test data of Jinchang Nickel-Copper Mining Area. According to geological mapping data of this area, we made nickel-copper geological model, where the layered nickel-copper orebody with dip was distributed in surrounding rock at 50～60 degrees. Then we did numerical simulating and wave-field characteristic analysis with high-order difference wave-equation in this model. Because layer spacing is tiny, the inclined events caused by inclined layers lead to mutual interferences, which are confounded and coupled with diffraction energy from pinch outs and catastrophe points, that’s why the wave-field is very complicated. We did conventional seismic processing and pre-stack depth migration processing on the test data, respectively. Mineral vein from the post-stack section of the former is stretched seriously, slim lined orebody entirely moves upward, diffracted wave and scattered wave seriously affect the quality of the post-stack section, this will cause error of interpretation; the migration image from the latter has higher resolution and restores the locations of mineral vein and orebody, the mineral vein boundary can be distinguished clearly, it’s easier to process geologic interpretation combined with geologic data. Based on the observation system and seismic data of Jinchang Nickel-Copper Mining Area, we processed the real data using processing software system with self-owned intellectual property rights and several self-developped characteristic processing modules, including Velocity Analysis Technique of High-resolution to improve accuracy of velocity spectrum picking; Surface-consistent Deconvolution Technique of Gabor Transformation to remove uncoupled effect between seismometer and bedrock caused by near surface weathering and to enhance the frequency band of seismic data; Combined Transformation Method of Suppression of Random Noise in Curvelet field to enhance seismic signal/noise ratio; Radon Transformation and Adaptive Mute Mask Filtering and Denoising Technique to suppress strong linear interference and extract weak reflection signals. Compared with post-stack section of conventional processing, synthetically using these characteristic processing techniques can greatly improve resolution and signal/noise ratio of seismic data of metallic orebody,obtain seismic stacked section of high quality and provide technical support for extension and distribution morphology of deep metallic orebody.更多还原