用matlab作地震波vsp图,利用VSP提高叠后地面地震资料分辨率

为提高地面地震的分辨率,特别是地面地震中深层的分辨率,人们通常采用反Q滤波、时变谱白化、振幅补偿、DMO(倾角时差)、小波高分辨处理等方法[1],而这些方法均依赖于精确的反褶积算子。子波反褶积是地震资料处理过程中最常用的一种提取反褶积算子的方法,该方法以某个子波作为地震波传播的已知子波,以子波的逆(即反子波)作为实际地震记录的反褶积算子[2]。但是,地下地质情况的差异往往导致地震波在不同地区传播特性的差异,于是用子波反褶积法所得到的地震子波往往与地下传播的真正地震子波存在较大的差异;其次,不同深度地层的地震子波经过上部地层的吸收和衰减作用后,已经不同于震源激发时产生的原始地震子波。因此仅利用受过严重干扰且经过衰减和吸收作用的地面地震资料无法得到高精度的反褶积算子[3]。VSP测井在潜水面以下激发、井下不同深度地层接收的探测手段为得到更加精确的各地层反褶积算子提供了一种切实可行的方法,又因VSP和地面地震所记录的都是地震反射波双程反射时间,故可利用VSP的子波来改善地面地震的分辨率[4]。对二维叠后地面地震资料进行反褶积处理并与其它反褶积方法所得结果进行的比较表明,该方法确能在一定程度上提高叠后处理剖面的频率和质量,特别对中深层更为有效。1方法原理设有一个地震子波w(t),存在一个滤波因子f(t),即反子波或反褶积算子,二者褶积存在以下关系:w(t)*f(t)=(t)(1)式(1)中:(t)为脉冲函数。对(1)式两边进行Z变换,可得:W(z)*F(z)=1(2)由(2)式得到子波的表达式为:W(z)=1/F(z)(3)经地震子波自相关的Z变换,有:Rw(z)=W(z)W(1/z)(4)式中:W(1/z)为W(z)的共扼。将(3)式代入(4)式,得Rw(z)*F(z)=W(1/z)(5)因为子波是一个实时函数,故有:rw(t)=rw(-t)Rw(z)=…+r2z-2+r1z-1+r0+r1z1+r2z2+…(6)考虑一个三点反滤波器(f0,f1,f2)的特殊情况:假设地震子波w(t)是最小相位的(因果的及可实现的),则它的反滤波器f(t)也是最小相位的,因此f(t)的Z变换为:F(z)=f0+f1z+f2z2(7)地震子波w(t)的自相关rw(t)的Z变换为Rw(z)=…+r2z-2+r1z-1+r0+r1z1+r2z2+…(8)地震子波w(t)的Z变换为:W(z)=w0+w1z+w2z2+…+wmzm(9)由此可得到:W(z)=w0+w1z+w2z2+…+wmzm(10)将(8),(9),(10)式代入(5)式,得:(r2z-2+r1z-1+r0+r1z1+r2z2)(f0+f1z+f2z2)=w0+w1z+w2z2(11)为求解(f0,f1,f2)需确定Z的各阶的系数:Z0的系数为:r0f0+r1f1+r2f2=w0(12)Z1的系数为:r1f0+r0f1+r1f2=0(13)Z2的系数为:r2f0+r1f1+r0f2=0(14)写成矩阵形式,系数方程为r0r1r2r1r0r2r2r1r0f0f1f2=w000(15)对f0作归一化处理,得:r0r1r2r1r0r2r2r1r01a1a2=L00(16)式中,a1=f1/f0;a2=f2/f0。对于矩阵中子波的自相关项ri,可根据反褶积的假设前提之一,即假设脉冲响应代表一个白化过程,再根据反射系数序列,得到:rx=r0rwrw=rx/r0(17)将(14)式矩阵中的子波自相关值用地震记录的自相关值来取代,然后通过解矩阵方程组,便可求得a1、a2的值。对于多个因子(f0,f1,f2,…,fn),(1