CN1282424A - 叠加地震痕迹的方法 - Google Patents
叠加地震痕迹的方法 Download PDFInfo
- Publication number
- CN1282424A CN1282424A CN98812246A CN98812246A CN1282424A CN 1282424 A CN1282424 A CN 1282424A CN 98812246 A CN98812246 A CN 98812246A CN 98812246 A CN98812246 A CN 98812246A CN 1282424 A CN1282424 A CN 1282424A
- Authority
- CN
- China
- Prior art keywords
- centerdot
- noise
- stack
- vestige
- signal
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Granted
Links
- 238000000034 method Methods 0.000 title claims abstract description 70
- 230000002596 correlated effect Effects 0.000 claims description 9
- 238000009795 derivation Methods 0.000 claims description 4
- 239000000470 constituent Substances 0.000 abstract 1
- 238000005516 engineering process Methods 0.000 description 21
- 230000008569 process Effects 0.000 description 14
- 239000011159 matrix material Substances 0.000 description 12
- 238000004458 analytical method Methods 0.000 description 10
- 238000006073 displacement reaction Methods 0.000 description 10
- 238000005070 sampling Methods 0.000 description 10
- 230000008859 change Effects 0.000 description 7
- 238000012937 correction Methods 0.000 description 7
- 230000006872 improvement Effects 0.000 description 7
- 238000004364 calculation method Methods 0.000 description 6
- 238000001914 filtration Methods 0.000 description 5
- 238000012545 processing Methods 0.000 description 5
- 230000003068 static effect Effects 0.000 description 5
- 230000008901 benefit Effects 0.000 description 4
- 230000000875 corresponding effect Effects 0.000 description 4
- XLYOFNOQVPJJNP-UHFFFAOYSA-N water Substances O XLYOFNOQVPJJNP-UHFFFAOYSA-N 0.000 description 4
- 238000001514 detection method Methods 0.000 description 3
- 238000000926 separation method Methods 0.000 description 3
- 238000006467 substitution reaction Methods 0.000 description 3
- 239000000654 additive Substances 0.000 description 2
- 230000000996 additive effect Effects 0.000 description 2
- 238000003491 array Methods 0.000 description 2
- 238000006243 chemical reaction Methods 0.000 description 2
- 230000006835 compression Effects 0.000 description 2
- 238000007906 compression Methods 0.000 description 2
- 239000012141 concentrate Substances 0.000 description 2
- 230000000694 effects Effects 0.000 description 2
- 230000010355 oscillation Effects 0.000 description 2
- 239000004215 Carbon black (E152) Substances 0.000 description 1
- 208000027418 Wounds and injury Diseases 0.000 description 1
- 230000005540 biological transmission Effects 0.000 description 1
- 230000015556 catabolic process Effects 0.000 description 1
- 230000001427 coherent effect Effects 0.000 description 1
- 230000000052 comparative effect Effects 0.000 description 1
- 238000010276 construction Methods 0.000 description 1
- 230000008094 contradictory effect Effects 0.000 description 1
- 230000006378 damage Effects 0.000 description 1
- 230000007547 defect Effects 0.000 description 1
- 238000006731 degradation reaction Methods 0.000 description 1
- 238000012217 deletion Methods 0.000 description 1
- 230000037430 deletion Effects 0.000 description 1
- 238000010586 diagram Methods 0.000 description 1
- 238000005553 drilling Methods 0.000 description 1
- 239000002360 explosive Substances 0.000 description 1
- 229930195733 hydrocarbon Natural products 0.000 description 1
- 150000002430 hydrocarbons Chemical class 0.000 description 1
- 230000005764 inhibitory process Effects 0.000 description 1
- 208000014674 injury Diseases 0.000 description 1
- 229910052500 inorganic mineral Inorganic materials 0.000 description 1
- 238000003780 insertion Methods 0.000 description 1
- 230000037431 insertion Effects 0.000 description 1
- 238000004519 manufacturing process Methods 0.000 description 1
- 238000005259 measurement Methods 0.000 description 1
- 239000011707 mineral Substances 0.000 description 1
- 239000000203 mixture Substances 0.000 description 1
- 238000004321 preservation Methods 0.000 description 1
- 230000001915 proofreading effect Effects 0.000 description 1
- 230000001902 propagating effect Effects 0.000 description 1
- 230000009467 reduction Effects 0.000 description 1
- 230000011514 reflex Effects 0.000 description 1
- 230000003252 repetitive effect Effects 0.000 description 1
- 230000004044 response Effects 0.000 description 1
- 239000004576 sand Substances 0.000 description 1
- 230000011218 segmentation Effects 0.000 description 1
- 230000035945 sensitivity Effects 0.000 description 1
- 238000004088 simulation Methods 0.000 description 1
- 238000001228 spectrum Methods 0.000 description 1
- 230000002269 spontaneous effect Effects 0.000 description 1
- 230000005654 stationary process Effects 0.000 description 1
- 230000001360 synchronised effect Effects 0.000 description 1
- 238000012360 testing method Methods 0.000 description 1
- 230000007704 transition Effects 0.000 description 1
- 230000003313 weakening effect Effects 0.000 description 1
Images
Classifications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V1/00—Seismology; Seismic or acoustic prospecting or detecting
- G01V1/28—Processing seismic data, e.g. for interpretation or for event detection
- G01V1/36—Effecting static or dynamic corrections on records, e.g. correcting spread; Correlating seismic signals; Eliminating effects of unwanted energy
- G01V1/362—Effecting static or dynamic corrections; Stacking
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V2210/00—Details of seismic processing or analysis
- G01V2210/30—Noise handling
- G01V2210/32—Noise reduction
- G01V2210/322—Trace stacking
Landscapes
- Engineering & Computer Science (AREA)
- Remote Sensing (AREA)
- Physics & Mathematics (AREA)
- Life Sciences & Earth Sciences (AREA)
- General Physics & Mathematics (AREA)
- Environmental & Geological Engineering (AREA)
- Geology (AREA)
- General Life Sciences & Earth Sciences (AREA)
- Acoustics & Sound (AREA)
- Geophysics (AREA)
- Filters That Use Time-Delay Elements (AREA)
- Noise Elimination (AREA)
- Measurement Of Mechanical Vibrations Or Ultrasonic Waves (AREA)
- Radar Systems Or Details Thereof (AREA)
- Geophysics And Detection Of Objects (AREA)
- Indication And Recording Devices For Special Purposes And Tariff Metering Devices (AREA)
Abstract
一种叠加地震痕迹的方法,包括步骤:导出地震痕迹每一条的噪声分量的表示;由噪声分量的表示导出噪声分量之间的相关指示;及通过按从相关指示导出的比例把地震痕迹相结合产生一个叠加地震痕迹。该方法能用来产生与常规叠加地震痕迹相比具有改进信噪比的地震痕迹。该方法的最佳实施例有时称作“协方差”叠加,因为导出的加权因数取决于其构成痕迹的噪声分量的协方差。
Description
本发明涉及一种叠加地震痕迹的方法,特别是一种产生改进信噪比的叠加地震痕迹的方法。
地震勘探的基本原理是众所周知的。声波由地球地表下岩石的特征反射,并且反射的信号(痕迹)由灵敏的换能器拾取,并被记录供以后分析。在陆地上,地震源典型地由爆炸提供,并且地震痕迹由称作地震检波器的换能器检测。在水上,压缩空气的脉冲由所谓的“气枪”提供,并且地震痕迹由水中检波器检测。(然而,其他源是可能的)。记录痕迹然后由多种已知技术处理以提供地表下岩石的表示。地质学家然后能使用该表示确定碳氢化合物矿床的可能性。多种技术可用来分析地震痕迹,并且由这些痕迹产生地表下岩石的显示超出了本描述的范围。
本发明更确切地说,涉及在地震分析中称作“叠加”的最早阶段之一。将会理解地震痕迹的能级一般非常低,即使当使用非常高的振幅源时也是如此。另外,在采样地震痕迹时,多个一般不可避免的噪声源是存在的。道路上车辆、头顶通过的飞机、工厂、钻探及诸如雨和风之类的自然现象的噪声,有可能把大量噪声添加到采样痕迹上。因为地震痕迹的低电平,采样痕迹的信噪比可能易于是如此之差,以致于阻止或至少严重妨碍这些信号的进一步分析。的确,有可能在一些情况下,这些痕迹的噪声分量的振幅大于信号分量的振幅。叠加处理寻求通过把来自不同换能器的那些痕迹和/或在不同时间(或在整个时间上)从相同换能器得到的样本相结合,来增大地震痕迹的信噪比。一般说来,地震痕迹的希望分量一直认为具有高的相关度,而地震痕迹的噪声元素一直认为是不相关的。因而,相同地震痕迹的简单求和,如在多不同换能器处接收的那样,一般将提供信噪比的改进。
然而,这种原始技术没有考虑关于不同痕迹的相对噪声电平。在授予P.Embree标题为“相异地震记录叠加方法和系统”的美国专利No.3398396中,确定地震痕迹的每一条的噪声能量,并且由来自换能器每一个的地震痕迹的加权和导出叠加的痕迹。施加到每条地震痕迹上的权重与该特定地震痕迹的噪声能量成反比。因而,从地震痕迹每一条导出对加权和的影响,但来自具有较差信噪比的那些痕迹的影响,小于来自具有较好信噪比的痕迹的影响。
相异叠加技术具有如下主要缺点。该技术依赖于地震痕迹的每一条的噪声分量与其他痕迹任一条的噪声分量无关的假设。然而,在多种实际情况下,该假设实际上无效。作为结果,已经发现,相异叠加经常导致相当差的信噪比。噪声分量的相关可能归因于起源在相同的源处(例如一个工厂),或者由于相关是数学概念,所以偶然来自分离的源。
本发明的一个目的在于,提供一种能提供改进信噪比的叠加地震痕迹的方法。
根据本发明,提供一种叠加多条地震痕迹的方法,该方法包括:导出地震痕迹每一条的噪声分量的表示;由噪声分量的表示导出在噪声分量之间的相关指示;及通过按从相关指示导出的比例把多个地震痕迹相结合产生一个叠加地震痕迹。
因为该方法不依赖于在每条痕迹上的噪声是独立的假设,所以也能衰减相关的噪声列。已经发现这提供具有比先技术的技术好在2dB至7dB之间的的信噪比的叠加痕迹。
从多个地震传感器可以导出痕迹,但也能在两次或多扫描期间从单独的传感器(例如检波器)导出。在这种情况下,想要的信号在以后扫描上相同或非常接近,并且这使得本发明特别适用,如下面更详细讨论的那样。
将表明,用于地震痕迹的组合的最佳加权因数与各种噪声信号的交叉能量密度相关。在一个最佳实施例中,用来导出比例的地震痕迹的每一条的噪声分量的表示是噪声分量的能量。因而,该实施例不仅在进行相异叠加时使用零滞后自相关(等效于能量),而且也使用在噪声分量之间的零滞后交叉相关。
最好,地震痕迹的每一个的噪声分量的协方差用来为多条地震痕迹的组合提供一组权重。
本发明可以进一步包括把样本打散成多个时间窗口。一般地,随着减小窗口长度,将改进组合地震痕迹的信噪比。然而,这对于信息信号的一个到两个波长的下限成为增大计算强度。已经发现400毫秒的窗口持续时间在10至80 Hz的区域内提供与典型地震痕迹的良好平衡。
本发明可以应用在时域或频域中。在时域中进行计算的优点在于减小地震痕迹的预处理。然而,在频域中处理信号具有能优化噪声减小的优点。这对频域分析或在时域分析之前把信号分离成频带都是正确的。
即使在时域中,该技术也通常得益于把地震痕迹的每一条分离成多个不同的频带。以前提及的授予P.Embree的美国专利No.3398396表明把地震数据分离成频带,但仅与相异叠加有关。可是,最好使用正交反射滤波器(QMF)把信号分离成单独的频带。对于进一步的细节,见WO 9620451或GB9614148-A0。这样做的优点在于这样的事实:在地震痕迹中,仅在一定的频带内发现重复噪声信号,而另一个频带或诸带主要包含纯随机噪声。通过分开地震痕迹频谱,能把叠加技术选择和修整为噪声存在的类型。
例如,如果在痕迹之间的需要信号的变化很大,协方差叠加可能实际除去一些信号。当信号振幅很大时,这可能发生。相异叠加可能最好应用在这些具体情况下。
现在通过例子、参照附图将描述本发明,在附图中
图1表示布置在陆地上的一组地震传感器,
图2表示采样噪声两种组合痕迹的噪声能量相对于应用于其组合上的加权因数,
图3表示在计算用于两个痕迹的加权因数相对于信噪比的相对误差,
图4表示正交反射频率(QMF)分离和重新组合的示意图,
图5表示使用先有技术和本发明的那些技术处理来自一组检波器的输出的一些比较例子,及
图6是在对其可应用本发明的技术的地震传感器的一种布置的平面中的图。
图1表示多个地震传感器(在这种情况下为检波器)G和一个地震源,爆炸源E,表示地球表面下少数几个示范射线路径。在爆炸发生之后,在检波器处检测的信号被采样并且被存储以供进一步分析。在水上使用船和多个水中检波器也能进行该过程。连续震动法(vibroseis)是一种使用以后更详细讨论的振荡源的可选择技术。
首先,我们讨论构成本发明基础的一些数学理论。
我们引入包括T个样本的信号a的期望值E(a)为 和两个采样信号a和b的交叉能量密度E(ab)
其中星号“*”指示共轭复数值。
通常,a和b表示(实数)时域信号;然而,该方法也可应用于频域数据。当把(开窗口的)时域数据转换成到频域以优化噪声减小时,这是相关的,如将在以后详细描述的那样。包括复数信号是直接的,并且产生作为其实数相应部分的表示;因此,复数信号包括在随后的分析中。在实数分析的情况下,必须计算(ai×bi *)的实数部分。
开始,为了分析容易,我们考虑需要叠加以得到一个单输出痕迹的两条采样痕迹x1和x2。经两条痕迹的线性加权过程进行叠加;第T个输出样本是两条痕迹每一条的第T个样本的加权和。因而,与痕迹加权有关的滤波过程纯粹是一维的,换句话说,没有来自以前或以后样本值的影响。
痕迹可能是时域痕迹(来自地震痕迹的专用时间窗口),但是本发明也可以应用于频域数据。采样痕迹包含一个反射信号和一个噪声信号ni(i=1,2)。现在进行反射信号r对于所有痕迹都是相同的假设。因为传感器布置、在多次扫描期间来自单独传感器的痕迹的叠加或应用于痕迹的预处理,这将是实际中的情形。例如,间隔不大于20米远的传感器将返回满足该标准的痕迹。因此,然后能把采样信号x1和x2写为
x1=r+n1 S1
x2=r+n2 S2
在如下分析中要叠加的所有痕迹包含相同反射信号的这种基本假设适用于如下情况:
ⅰ.使用来自一起靠近布置、在成群组中的检波器的样本(称作“空间跳跃”)
ⅱ.除去振幅扰动和静态扰动(“扰动校正”)
ⅲ.根据其偏移差校正在痕迹之间的剩余正常时差(以后更详细地描述该过程)
ⅳ.分析在不同时刻来自单独检波器(或其他传感器)的痕迹
由这两个采样痕迹,我们需要建造单独输出痕迹。我们通过两条痕迹的加权求和或叠加形成一个输出痕迹y,其中加权因数用wi(i=1,2)指示:
y=w1x1+w2x2 S3
通过引入约束条件:
w1+w2=1 S4
我们保证信号能量对于所有组的加权因数是相同的。该约束条件在该分析的整个期间适用,从而由不同技术导出的叠加痕迹的信号能量值总是相同的。因而,信噪比的改进直接对应于噪声能量的减小。通过使噪声能量最小,信噪比是最佳的。
我们现在计算在输出信号y中存在的噪声能量。噪声能量N2由下式给出
在包含T个样本的整个窗口上计算能量值。以后讨论窗口长度T的影响和选择。由于在两个信号的噪声分量之间的交叉能量密度给出为:
E(n1n2)+E(n2n1)=2Re{E(n1n2)} S6
其中“Re”指示实数部分,所以有可能导出用于噪声能量的如下表达式:
使用加权因数之和是1的约束条件,我们能代入w2=1-w1,并且得到:
作为加权因数w1的函数的噪声能量的例子表示在图2中。这里测量数据的两条雨噪声痕迹用来表示噪声。
当w1=1,w2=0时,第一痕迹的噪声能量在其中w1=1的点处给出。因而第一痕迹的噪声能量是11×10-6。当w1=0,w2=1时,第二痕迹的噪声能量在其中w1=0的点处给出,约为3.6×10-6,与第一痕迹在相同的量级上。最佳加权因数与该曲线的最小值相对应。
噪声能量当其相对于w1的导数为零时为最小(能容易地证明曲线的极值点的确总是最小值,而绝不是最大值)。微分给出:
并且使用约束条件w2=1-w1,能导出用于两个加权因数的各自值:
对于任何种类的噪声,公式S10代表用于痕迹叠加的最佳加权因数。该技术称作“协方差叠加”,因为加权因数取决于其构成痕迹的噪声分量的协方差。诸如反比叠加或痕迹求和之类的其他技术将在输出痕迹中给出较大的(或者最好是相同量的)噪声能量。一个例子在图2中给出,其中表示用于协方差叠加、反比叠加和简单痕迹求和的总噪声能量。协方差叠加位于曲线上的最低点。该技术不要求加权因数是正的,只要其和等于一就行,他们根据地震数据能假定任何值。例如,图2表示当加权因数之一(w1)是-0.25的一个例子。因此另一个加权因数(w2)是1.25。
为了表明按照本发明的叠加技术对于反比叠加的优越性,进行如下比较。
通过把噪声协方差设置为0,从以上为协方差叠加计算的加权因数可以得到用于反比叠加的加权因数。这是从反比叠加基于的、在不同痕迹中的噪声不相关的假设得出的:
Re{E(n1n2)}=0 D1
把该假设插入在公式S10中,由此得出对于反比叠加,加权因数由下式给出:
能重写该公式,从而与反比叠加基于的原理的关系-借助于其噪声能量的倒数加权痕迹-是清楚的:
其中A是比例因数,引入以满足w1与w2之和等于1的条件(约束条件S4)。按如下计算A:
因而,当-并且仅当-要叠加的痕迹满足如下任一个时,由反比叠加提供的输出痕迹等效于由协方差叠加提供的:
Re{E(n1n2)}=0(即,噪声完全不相关) D5
或者
E(n2 1)=E(n2 2)(即,两条痕迹具有相同的噪声能量)D6
在实际情形中满足这两个条件的任一个的可能性是可忽略的。在所有其他情况下,协方差叠加技术的性能是优秀的。
为了完整起见,我们把协方差叠加与直接痕迹求和相比较。这里,我们简单地有
因而,如果-并且仅如果-两条痕迹具有相同的噪声能量,则直接痕迹求和等效于反比叠加。这等效于以上公式D6。因而,当噪声完全不相关、并且痕迹的噪声能量相等时,直接痕迹求和将给出与协方差叠加相同的性能。在所有其他情况下,协方差叠加技术的性能是优秀的。
叠加输出痕迹的噪声能量N2由以上公式S7给出
对于叠加技术的每一种下面是用于输出噪声能量的一般表达式。能证明对于两种技术的噪声能量由如下给出,
反比叠加
协方差叠加
现在讨论本发明提供的优于先有技术反比叠加技术的改进量。
我们现在计算在叠加之后对于信噪比发生的情况。叠加地震痕迹的最终目标是改进信噪(S/N)比。在叠加之后,信号能量S2由下式给出:
S2=E(r2) SN1
这对于先有技术的反比叠加和对于协方差叠加都是相同的,因为加权因数等于一。我们把在每条痕迹i上的原始S/N比,S/Ni定义为 并且把在叠加之后的S/N比定义为
在包含T个数据样本的有限窗口上计算S/N比。应该注意,为了方便起见,依据能量(而不是能量的平方根)定义S/N;这意味着必须按101og(S/N)而不是201og(S/N)计算S/N分贝(dB)。
我们现在比较对于反比叠加的S/N比与对于协方差叠加的S/N比的比值,并且我们证明该比值总小于1(条件是各个痕迹的噪声能量不相同,并且噪声不是完全不相的关)。在一些变换之后,它能表示为: 现在我们知道:
通过公式SN4和SN5之间的模拟,我们能导出如果(S/N)DIV小于(S/N)COV必须满足的简化表达式(注意,通过假定a>0,我们已经隐含地假定噪声不是完全不相关的。)
这等效于q<p,并且我们现在证明,对于由反比叠加导出的痕迹的S/N比总是小于对于由协方差叠加导出的的痕迹的S/N比。该不等式能写为
显然,只要E(N2 1)≠E(n2 2),二次表达式就总是正的,并且我们已经证明,对于不是完全不相关的且从痕迹到痕迹不同的噪声,协方差叠加总是导致相对于反比叠加的改进S/N比。
注意噪声不是完全不相关的条件不必归因于物理过程(例如,通过的汽车在相邻痕迹上产生相关的噪声)。条件E(n1n2)≠0是一个数学条件,意味着在T个样本的时间窗口上两条噪声痕迹的交叉相关不应该与零相同。仅由于纯粹的巧合,交叉相关才为零;在实际中,噪声在时间窗口上总是相关到一定的程度,即使来自每条痕迹的噪声纯粹是随机的也是如此。相同的理由适用于噪声痕迹不同能量值的条件。换句话说,噪声能量在不同的痕迹上绝不会准确地相同,以及绝不完全不相关。
至此,计算一直基于噪声能量能准确地从采样迹确定的假设。然而,在实际中,只可能得到畸变的噪声估计,因为信号能量也是存在的。结果,我们得到估计的加权因数
和
而不是真正的加权因数w1和w2,
和
对于协方差叠加由下式给出 而对于反比叠加由下式给出
如果我们把加权因数Δw的值中的误差定义为Δw=
w1-w1,对于协方差叠加我们得到:
注意信号能量不进入该公式;误差由信号与噪声与之间的相关确定。如果信号和噪声是不相关的,并且他们任一个具有零平均值,则误差是0。
对于反比叠加该情形相差很远。一般地说,误差由下式给出:
NE4
我们现在代入
xi=r+n1 NE5
但表达式变得相当长。现在讨论与反比叠加有关的误差。首先,反比叠加假定信号与噪声是不相关的,导致误差。当试图从数据导出加权因数时产生另外一个误差,并且该误差包括由在信号与噪声之间的互相关是非零的事实造成的第一分量。因为迹包含信号能量以及噪声,所以第二分量产生。这里为了简单起见忽略第一分量,并且这是一种有效的简化,因为如果对于协方差叠加这样做,则不会偏置加权因数。为了简化本计算中的问题,我们假定信号与噪声是完全不相关的,并且他们的任一个具有零平均:
E(mi)=0 NE6
注意,在该假设下,协方差叠加产生加权系数的完好估计。如果我们使用在痕迹i上的信噪比的定义,在公式SN1中给出S/Ni,则误差能写为S/N1和S/N2的函数 通过表达下式从以上表达式能计算相对误差Δw/w1:
因而,由于反射信号的存在在加权系数的计算中引入的误差现在取决于信噪比;即使信号与噪声是不相关的,对于反比叠加也将得到有偏估计。相对误差Δw/w1(×100%)表示在图3中,作为S/N1和S/N2的函数。
在理论上,当使用协方差叠加时,即使在每条痕迹上的噪声完全相关,那么噪声也将完全衰减噪声。只要噪声的振幅从痕迹到痕迹变化,同时信号所有痕迹上相同,则完全相关的噪声从迹到痕迹仍然不同,并且能完全除去。然而,在实际中,反射信号在所有痕迹上不是相同的。因此,即使噪声的振幅从痕迹到痕迹变化,完全相关的噪声也不可能与包含振幅扰动的希望反射信号区分开。为了明白这点,考虑如下理论例子。在假设无噪声的一天测量两个采样痕迹。噪声实际是零,但信号之一的反射元素在振幅中稍受扰动,其振幅是其他痕迹的振幅的(1+δ)倍。痕迹因而能写为
x1=r+0 AP1
x2=r+δr AP2
在这种情况下,协方差叠加将给出零的输出。这显然是不可接收的,所以现在论述振幅扰动的问题。为此目的,我们考虑如下两条痕迹:
x1=r+n1 AP3
x2=r+δr+n2 AP4
其次,我们进行2种假设:
ⅰ.信号与噪声是不相关的(为了数学上方便这样做,并且与扰动问题没有特别关系)
ⅱ.振幅扰动是小的;特别是,我们将利用近似δ2<<δ。
假设ⅰ.是合理的,因为信号扰动问题带有除去反射信号的严重危险。相反,相关信号与噪声的问题导致加权因数的偏差。在输出上的S/N比可以是非优化的,但有很小或没有严重损害反射的危险。
假设ⅱ.仅反映了我们需要应用认真扰动的事实。即使如此,如果我们不校正振幅扰动,假设AP5仍将适用于虚拟的所有实际情形。即使对于30%的振幅扰动,假设ⅱ.也仍是良好的近似。
组合之后的输出痕迹将是:
y=w1x1+w2x2=r(w1+w2(1+δ))+[w1n1+w2n2] AP5
应该简单等于r的信号输出现在等于
注意我们能识别项
作为残余信号能量与微分噪声能量的比值。在其中噪声实际为零的上个例子中,信号完全消失。该公式解释为什么这会发生。如果微分噪声能量是零,则用于r的加权因数也等于零(因为分母成为无限大)。
只有在各条痕迹上的S/N非常高时该问题才会出现,并且振幅扰动是显著的。例如,对于每个单独检波器上的20 dB S/N比:
10log(S2/N2)=20 or S2/N2=100 AP8
和10%的振幅扰动,信号输出等于0.5r(与1×r的希望输出相对)。注意这里的S/N比值指平均信号和噪声能量,而不是指波峰振幅。在一些情况下,波峰信号振幅可以非常高;但甚至在这种情况下,平均能量比值也处于大得多的中级处。尽管如果已经采用适当的预处理则其中出现问题的情形可能特别少见,但对问题的解法是存在的,从而协方差叠加仍然是一种通用的可适用技术。
保护信号的最简单方法是把阈值放在两条痕迹之间的校正系数ρ上。ρ由下式给出:
其中能量估计如以前那样从包含T个样本的数据窗口得到。如果ρ具有超过阈值的值,有可能对于较小的计算开销回复反比叠加;在这种情况下,痕迹的能量是如此强大,以致于我们不需要由协方差技术提供的改进噪声衰减。作为用于阈值的实际值,在0.8至0.9量级中的值是适当的。因为事实上,如果信号扰动与噪声相比较大,信号降低才发生,基于完整信号的该标准提供克服信号变坏的更适当保护。
在上面,表明用于协方差叠加的最佳加权因数能准确地由数据确定,条件是
E(mi)=0 DC1
换句话说:信号和噪声需要是不相关的。如果信号和噪声是独立的,并且另外,如果信号或噪声具有零平均值,则满足不相关信号和噪声的条件。在这些约束条件后面的原因是信号和噪声的直流分量总是完全相关的(注意:也能完全否定地相关)。这里假定最坏情况的情形,其中直流分量具有等同的符号,并且相关是肯定的。直流分量没有时间变化;用于在两条痕迹a与b之间完全相关的条件a=kb,因而总满足直流信号。
直流分量的问题不专用于协方差叠加;它同样好地适用于反比叠加。事实上,如已经表明的那样,即使信号和噪声是完全不相关的,则反比叠加仍不能得到加权因数的无偏估计,而协方差叠加能实现无偏估计。
通常,地震痕迹的确具有零平均值。这由接收检波器(探测仪器)的输出的电路的低灵敏度滤波器证明。然而,如果选择时间窗口,或者甚至时间窗口的频域表示,则一般直流分量将存在。
如果在加权因数和叠加痕迹的计算之前,除去或保持采样数据窗口的直流分量,则比较跟随在所发生的事件之间。考虑如下两个处理序列:
ⅰ.不要从痕迹除去平均值;计算加权因数,并且叠加
ⅱ.除去平均值,计算加权因数,叠加,及把一些种类的原始平均值的平均添加到叠加结果上。我们仍然需要直流分量的原因在于,一般地说,直流分量将在窗口之间变化;因而,局部它是直流分量,但对于整个地震痕迹的重新建造它不是。
一般地说,因为三种原因这些结果是不等同的:
a.在第一方法中,平均值影响加权因数的值;在第二方法中,它不会。因而,对于两种方法加权因数是不同的。
b.第二方法要求需要添加的原始平均值的定义。例如,它能是各条痕迹平均值的平均值。即使平均值不影响加权因数,两种方法仍产生不同的结果,除非在方法ⅱ中减去的平均值是两个平均值的加权和,其中加权因数等于在噪声减小算法中计算的加权因数。
c.如果信号和噪声是独立的,第二方法不再遭受漏入到加权系数计算中的信号能量之害(这对于协方差叠加是真实的;对于反比叠加,减小但不会完全消除信号能量泄漏)。因此,第二方法产生准确的加权因数,并因而产生更好的噪声减小,以防信号强大(然而,这准确地是其中增大噪声减小较不关键的情形)。
第一方法提供最好的噪声减小;加权因数也按比例缩小直流噪声能量。在第二方法中,不相关信号和噪声的假设更接近相似;然而,特别是如果窗口长度变小,并且平均值变成痕迹的有效部分,则平均值的非最佳加权对噪声减小是有害的。因而,第一方法(在考虑的数据窗口中没有直流除去)是最好的。当数据的平均值接近零时达到最好的性能。这对能采用的最小窗口尺寸施加一个下限。一般地说,整条痕迹的直流值实际上是由检测仪器的低灵敏度滤波器生成的零。信号样本的直流值当然等于样本值本身。随着样本数量增大,直流值向零衰减。几百毫秒量级的窗口长度提供良好的折衷。
因此建议不从数据除去直流分量,但使用具有这样一种长度的窗口,从而平均值与总痕迹能量相比较小。
我们现在考虑窗口长度的影响。必须记住,在包含T个时间样本的时间窗口上计算在这种分析中给出的加权因数(对应于时间t1,t2,…tT,其中ti=iΔt)
输入地震痕迹的开窗口,即细分,对协方差叠加的精度有重要的影响。在非常小的窗口上的叠加、计算负载及希望数据的波长之间有矛盾的要求。
有两个涉及输入痕迹开窗口的主要问题:
a.哪种窗口长度导致协方差叠加的最好性能?
换句话说,哪种窗口长度产生最好的噪声衰减?
b.对于给出的窗口长度,应该以什么样的间隔计算加权因数?
a.就第一问题而论:确定加权因数的准确度、和信号保存的程度肯定与窗口长度有关;但使用小窗口长度达到最高的噪声衰减。几百毫秒的窗口长度是真实的实际值。下面,我们概述该论述所基于的理论。
在第一情况下,在协方差叠加之前,由数据估计加权因数。在该过程中,我们必须假定信号和噪声是不相关的。如果以下条件成立则该假设是正确的:
Re{E(r(n2-n1))}=0 W1
为了估计W1是否合理,我们按如下进行。首先,合理的是假定反射信号和噪声是独立的:
Re{E(r(n2-n1))}=Re{E(r)}×Re{E(n2-n1)) W2
因而,由数据对协方差叠加加权因数的估计是准确的,条件是噪声或信号分量具有零平均值:
Re(E(r)}=01 or Re(E(n2-n1)}=0 W3
信号和噪声的平均值取决于所考虑的时间窗口。一般地说,窗口越长,平均值越小。例如,把整条迹作为窗口保证数据的平均值的确(实际上)是零,因为在检测仪器中通常存在的低灵敏度滤波器。仅使用1个样本作为窗口的另一极端情形,导致平均值等于一般不为零的单独样本值。注意,1个样本的窗口受另一个问题之害:在这种情况下没有区分信号和噪声的依据,并且协方差叠加将简单地产生零输出(如果样本具有相同的极性)。
在一个例子中,为了得到用于信号能量的近似零平均值,我们能把目的定为选择与地震小波持续时间有关的窗口长度。合理的是,假定具体的反射信号具有近似的零平均值。对于20-40Hz的典型反射主频率,我们得到25-50毫秒的小波持续时间。然后建议几倍该小波持续时间的最小窗口长度(因而几百毫秒),以解释延伸到考虑窗口外部的反射。没有上限,但一般地说,噪声减小随增大窗口长度而减弱。
总之,在第二种情况下,噪声衰减的程度随减小窗口长度而增大。原因在于,当在各条噪声痕迹之间的振幅差最大时,或者当痕迹具有高相关度时,协方差叠加工作得最好。选择大窗口趋于减小这些效果。因而,从该观点出发,推荐尽可能短的窗口长度。作为极端情况的一个例子,如果我们考虑纯噪声记录,1样本窗口对于协方差叠加导致完全噪声取消,而当然也导致信号取消。
b.就第二问题而论:首先,我们假定把对于给出窗口计算的加权因数分配给窗口中的中心样本。一般地说,我们计算每第N个样本的加权因数,由此对于每次以后计算向下移动窗口N个样本。如果我们指示窗口长度为L,并且窗口重叠量用O指示,其中O由下式给出
现在我们考虑要使用什么样的重叠O。最好具有99左右%的重叠(因而对于每个样本计算加权因数):这给出最高精度,并且加权因数的平稳状态保证信号扰动的适当处理,但以较高的计算负载成本。下面,我们概述该论述所基于的理论。
ⅰ.精度。一般地说,噪声序列是非静态的。这通过使用大窗口重叠(和短窗口长度)最好地描述。
ⅱ.计算负载。显然,计算负载随增大的窗口重叠增大。
ⅲ.信号扰动的存在。这些引起性能随长重叠变坏。
在理论上,我们假定反射信号在所有迹上都是相同的。因此,对于加权因数沿地震痕迹组中的痕迹的所有可能组合、和对于作为痕迹内时间函数的加权因数的所有变化,信号输出总是相同的(条件是加权因数之和是1)。然而,在实际中,从痕迹到痕迹的信号变化将出现,作为时间和痕迹数量的函数的强烈变化的加权因数因此产生地震波记录的不连续外观。选择100%重叠意味着用于每个新的加权因数计算的能量项几乎是相同的,每次以后的计算要求一个样本的删除和一个新样本的引入,同时保持所有的其他样本值,因而,变化不期望很大。尽管如此,如果新引入的样本代表一个尖峰,加权因数跳跃出现,一些合适的滤波应该用来防止这种情况,例如使用滑动平均、低通或类似的滤波器。注意,窗口长度也影响加权因数的变化:如果窗口长度较大,则只替换一个样本的影响相当小。
信号扰动问题不是协方差叠加特有的,它也出现在反比叠加中。加权系数的变化不仅取决于窗口重叠,它也取决于痕迹的数量。一般地说,在估计中使用的痕迹越多,加权系数的状态越平稳。
当噪声相关最大时,协方差叠加的性能对于短窗口长度是最好的。然而,为了完成理论上的假设(用于信号或噪声的零平均值),窗口长度最好应该很大。而且,重叠窗口和/或加权因数的平稳可以用来进攻反射信号扰动和非静态噪声特性。
通常,协方差叠加适用于时间窗口;在一个时间窗口内,噪声差和相关性能比在整个痕迹中更显著。我们能通过在几个频带中划分每个时间窗口能更进一步。这能进一步改进协方差叠加的噪声减小能力。考虑这样一个例子:其中数据具有在30与50Hz之间变化、相关的噪声分量,但也包含强大的但均匀的、在0与30Hz之间的不相关噪声分量。在这种情况下,反比叠加和协方差叠加都不能显著改进简单求和结果。原因是,两个叠加过程所需的(交叉)能量计算由该0-30Hz分量支配;关于30-50Hz噪声的所有信息丢失。如果我们能把协方差叠加用于不同的频带(例如,0-30Hz和30-60Hz),我们能够减弱变化、相关的30-60Hz噪声序列。一种实现这点的装置是所谓的正交反射滤波器(QMF),它能把时间序列分裂成N个序列,其中N个分裂时间序列的每一个对应于不同的频带。值得注意,如果噪声是非白的和/或噪声的相关性质随频率变化,则该方法是便利的。一般地说,在实际中一定程度地满足这些条件的两个。
图4表示在协方差叠加过程中如何集成QMF过程。在该例子中,我们假定把数据分裂成4个频带(F1、F2、F3、F4),在实际中能选择任何分裂,只要它是2的乘方就行(对于关于QMF的更多细节,读者参考以上专利文件)。开始点是两条输入痕迹x1和x2。其次,把每条痕迹划分成4部分,每部分代表使用QMF滤波器的不同频率范围。重新分类分裂痕迹,从而协方差叠加适用于每条痕迹的相对频带。最后,使用逆QMF技术建造输出数据,组合4个不同的频带。
为了简单起见,以前的描述仅集中在两条地震痕迹上,并且这种描述现在扩展到包括任意数量的地震痕迹。我们考虑包含M个检波器的检波器组。如果检波器靠近布置在一起,则已经除去静态和振幅扰动,而如果已经描述了在检波器之间的残余NMO(见后面),则这些检波器将接收相同的反射信号r,但他们一般将包括不同的噪声信号ni。然而,本发明也可应用于单独检波器,特别是对于连续震动法,并且可用于给出适当预处理的更大间隔的传感器。第i个检波器的检波器输出xi写为
xi=r+ni M1
其次,我们考虑痕迹xi的加权叠加y,其中加权因数由wi指示 在求和输出痕迹y中的平均信号能量S2由下式给出 在求和输出中的平均噪声能量N2由下式给出 信噪(S/N)比然后从下式得出
其次,我们确定加权因数wi,从而使S/N比最大,服从加权因数之和是一的约束条件:
因此,信号能量对于满足约束条件M6的加权因数的所有可能组合是相同的,并且问题是现在找到最小噪声能量。为此目的,约束条件M6写为
注意,我们选择根据其他加权系数(W1、W2、W3、W(M-1))表示。这种选择是任意的,并且有利于数学便利,仅把用于WM的约束条件M7代入公式M4产生
我们通过把噪声能量相对于加权系数的每一个的导数设置为零,使噪声能量最小
在一些变换之后,这产生 这必须满足k=1,2,…,M-1现在便利的是,按如下建造每条噪声痕迹与最后噪声痕迹的差
为了符号便利,我们也引入零痕迹n0(就是说,仅包含零的痕迹,分配后辍零),并且定义: 公式M10然后能写为 现在,为了更进一步简化,我们引入符号 并且公式M13最终成为 这能按如下写为矩阵公式 其中
是(M-1)×(M-1)矩阵,其元素由
给出,
w是包含加权系数的向量 并且
是其元素由下式给出的向量 加权因数然后简单地从下式得出 注意
是对称矩阵,并且使用Cholesky分解能高效地计算其逆矩阵。解矩阵公式产生加权系数M17的向量,然后作为最后步骤,经如下关系计算WM 能表明在求和输出中的生成噪声能量, 对于协方差叠加减小到
在公式M22中能使用j在区间[1,M]中的任意值以计算噪声;它将总是给出相同的结果。
至今,我们已经按噪声能量表示方差加权因数。在实际中,噪声是未知的,因为记录包含希望的反射和附加噪声。然而,矩阵
容易从数据建造。因为ni-nj = (r+ni)-(r+nj) = xi-xj NE1矩阵
能直接由矩阵
代替,从而我们得到 对于 并且从下式得到加权因数 剩下问题是由数据估计
。
的元素由下式给出 在表达式中的最后项容易从数据得到ni-nM=xi-xM NE6从而bi=Re{E(-nM(ni-nM))}=Re{E(-nM(xi-xM))} NE7
这给我们留下确定噪声nM的问题;这里我们简单地从原始数据估计噪声痕迹(注意:能进行其他近似)
如以前讨论的那样应该通过在痕迹之间限制最大允许相关系数,再加强信号保存。
另外,我们能使用不同的信号估计。例如,如果求和所有痕迹,则输出包含信号“r”外加所有痕迹的求和噪声。因为噪声相对于信号衰减,所以信号估计具有较小噪声。如果M个传感器之和称作SUM,则噪声估计是x(M)-SUM。为了改进噪声估计,我们一般需要信号的一些估计,从而能减去他。因此,代之以求和所有痕迹,我们能进行反比叠加或导出中间值等。的确,有可能使用方差叠加导出信号估计,然后重新进行协方差叠加,但这当然会增大需要的处理。
由协方差叠加施加的计算负载大于反比叠加要求的计算负载。然而,增大是按如下表示能提供的增大。
M个痕迹的反比叠加需要M个能量值的计算。如果我们取T个时间样本的情形,则需要进行总数2MT次运算(MT次乘和MT次加),以得到所有需要的比例因数。
相反,对于M条痕迹的协方差叠加,我们需要形成一个包含所有痕迹的交叉能量项的(M-1)×(M-1)矩阵。这要要K个能量值的计算,其中K由下式给出
运算的总次数因此等于2KT。由协方差叠加要求的过多运算次数因而是:
例如,如果输入是每条由2000个样本组成的5条痕迹,则反比叠加需要20,000次运算,而协方差叠加需要60,000次运算。另外,(M-1)×(M-1)矩阵需要求逆,并且逆矩阵乘以(M-1)向量。矩阵求逆和矩阵向量乘积的次数等于在单条痕迹中的时间窗口的次数。
得到的试验结果表示在如下表中:
噪声记录 | 单传感器配置 | 阵列配置 | 噪声类型 |
1 | 5星形,4米。在传感器之间 | 在10×10米栅格上的33个元件 | 汽车发动机 |
2 | 5星形,4米。在传感器之间 | 在10×10米栅格上的33个元件 | 汽车通过 |
3 | 5星形,4米。在传感器之间 | 在10×10米栅格上的33个元件 | 在路上的汽车 |
4 | 5星形,4米。在传感器之间 | 在10×10米栅格上的33个元件 | 在路上的卡车 |
5 | 5星形,4米。在传感器之间 | 在10×10米栅格上的33个元件 | 钻探操作 |
6 | 5星形,4米。在传感器之间 | 在10×10米栅格上的33个元件 | VVFB发动机 |
7 | 5星形,4米。在传感器之间 | 在10×10米栅格上的33个元件 | 飞机 |
8 | 5星形,4米。在传感器之间 | 在10×10米栅格上的33个元件 | 工厂 |
9 | 5星形,4米。在传感器之间 | 在10×10米栅格上的33个元件 | 靠近传播的脚步 |
10 | 5星形,1米。在传感器之间 | 在2×2米栅格上的15个元件 | 雨 |
11 | 阵列的前5个元素 | 18元件阵列 | 风 |
VVFB指由Geco-Prakla制造的地震振动器。为了比较目的,以七种不同的方式处理结果:
ⅰ.单独传感器输出;没有应用处理
ⅱ.5个传感器之和;直接求和
ⅲ.5个单传感器的反比叠加
ⅳ.5个单传感器的QMF反比叠加
ⅴ.5个单传感器的协方差叠加
ⅵ.5个单传感器的QMF协方差叠加
ⅶ.阵列响应;所有阵列元件的直接求和
在下表中给出应用于十一种噪声记录的七种技术的结果
ⅰ ⅱ ⅲ ⅳ ⅴ ⅵ ⅶ1 3.7dB 0.0dB -0.6dB -0.7dB -2.2dB -2.6 dB -1.9dB2 4.3dB 0.0dB -1.4dB -1.3dB -4.1dB -4.6 dB -1.3dB3 3.1dB 0.0dB -0.1dB -0.2dB -1.5dB -1.7 dB -1.3dB4 2.4dB 0.0dB -1.2dB -1.7dB -2.3dB -4.7 dB -7.6dB5 10.8dB 0.0dB -1.4dB -1.3dB -5.6dB -5.9 dB -4.3dB6 5.7dB 0.0dB -1.1dB -1.0dB -3.6dB -3.7 dB -1.3dB7 4.1dB 0.0dB -0.7dB -0.8dB -3.0dB -3.2 dB -1.7dB8 4.3dB 0.0dB -1.4dB -1.3dB -4.1dB -4.6 dB -1.3dB9 1.3dB 0.0dB -2.5dB -2.7dB -3.3dB -4.7 dB -4.1dB10 1.7dB 0.0dB -1.6dB -2.0dB -3.4dB -4.4 dB -1.4dB11 2.4dB 0.0dB -2.8dB -5.0dB -3.6dB -6.7 dB -8.3dB
QMF过程允许叠加之后的完整重新建造。在所示的例子中,数据窗口由QMF过程细分成8个频带。如以上讨论的那样,分裂成频带改进噪声衰减,因为清楚地说明作为频率函数的非均匀噪声状态。
用于每个噪声减小过程的总噪声能量表示在附图5中。把5个单传感器的直接痕迹求和选择为基准(0dB)减小过程。阵列结果不可能给出最佳衰减,因为在结果中,这些阵列的空间范围受到限制;尽管如此,阵列噪声减小值给出一些指示。在这些结果中的一个问题是来自附近工厂的几乎连续存在的(低频)工厂噪声,工厂噪声因此影响噪声减小的计算量。通过在f-t域中显示整个噪声记录(一种有助于噪声状态解释的技术),释放这种限制,然后能得到作为频率函数的衰减的良好抑制。
对于200毫秒或400毫秒窗口得到良好的结果。没有上限,但一般地说,窗口越长,衰减越低。
我们考虑由5个单独传感器组成的单独接收器组。单独传感器以星形图案成组,每个元件到组中心的距离等于Δx。配置表明在图6中的平面中。反射由单独传感器的每一个接收。为了简单起见,我们假定反射的行进时间t服从如下通常的时差(NMO)公式
其中t0指示零偏移行进时间,X指示源接收器距离,及Vrms指示叠加速度。在处理期间,行进时间校正(通常指NMO校正)用来校正偏移对行进时间的影响。为此目的,以上公式按其一阶Taylor近似展开。这产生
数量t-t0是NMO校正,并且指示为tNMO;而且,由t近似t0产生用于NMO校正的如下表达式
在该段中,我们对针对数据计算和应用NMO校正不感兴趣。这里,我们集中在这样的观察上:在组内的不同传感器具有不同的NMO位移,并且在数字组形成之前应该从数据中除去在NMO校正中的该变化(因而,不是NMO位移本身)。现在更详细地解释这点。
如果在一个源与跨过组的单独传感器之间的距离有变化,则每个单独传感器具有不同的反射到达时间。到达时间的不同假定仅由位移X引起(其他扰动,如近表面静态变化,这里不考虑,因为他们通过扰动校正处理,如已知的那样)。为了计算一个组中在单个传感器之间的反射到达时间差,我们需要计算tNMO相对于X的导数
如果我们假定把组中心选择为我们的基准位置,并且该组中心具有X的源接收器距离,我们因而必须借助于时间位移校正每个剩余的单独传感器,以描述源接收器偏移的差。相对组中心-单独传感器距离不是简单地在两个之间的距离,而是在线(即,在源与组中心之间的)分量。我们用ΔxNMO指示该相对距离。对于表示在图6中的各种传感器,下表概括在元件间隔Δx与至组中心的相对距离ΔxNMO之间的关系
传感器号 | 至组中心的距离 | 相对距离ΔxNMO |
S2 | Δx | -Δx |
S3 | Δx | 0 |
S4 | Δx | Δx |
S5 | Δx | 0 |
注意相对距离也能负的。通过设置下式得到我们然后需要应用的“残余”NMO校正ΔxNMO
dX=ΔxNMO NM5
这产生:
简单地求和所有信号是不允许的,因为作为传感器不同源-接收器距离的结果在反射信号之间存在的小时间位移。为了确定NMO影响的大小,我们需要估计在以上公式中的各项有多大。位移X一般在0与4千米之间变化。RMS速度VRMS一般在1500m/s(浅深度)至5000m/s(在几秒的深度处)变化。到达时间t在0与6秒之间变化。相对检波器间隔ΔxNMO在0与10米之间(对于单独传感器记录)、和在0与25米之间(对于常规检波器阵列)之间变化。
显然,最坏的情况(即,最大时间位移)在较大偏差、小到达时间、小叠加速度及至组中心的大距离期间出现。一些典型的例子如下。
检波器阵列,总长度50米,以1秒反射,叠加速度3000米/秒,偏移3000米。在阵列内的最大时间位移:8.3毫秒。
单独检波器图案,10米,在传感器之间的距离,以1秒反射,叠加速度3000米/秒,偏移3000米。在图案内的最大时间位移:3.3毫秒。
对于较深反射,我们呈现如下例子:
检波器阵列,总长度50米,以3秒反射,叠加速度5000米/秒,偏移3000米。在阵列内的最大时间位移:1毫秒。
单独检波器图案,10米,在传感器之间的距离,以3秒反射,叠加速度5000米/秒,偏移3000米。在图案内的最大时间位移:0.4毫秒。
这些时间位移是有效的,因为他们降低求和痕迹的质量。影响是频率相关的,并且其严重性随频率增大。当然,时间位移随减小的偏移线性地减小;而且,较小的单独传感器图案将减小时间位移。但是,作为用于这些时间位移的校正可能是完全可行的,对于单独传感器记录,我们应该利用该机会从数据中除去这些NMO扰动。对于其中叠加速度一般比对于P波(压缩波)小的剪力波记录,影响是最大的。
本发明特别适用于利用振荡源的连续震动法。这样一种源一般包括一块安装在卡车上与地球接触且振动的底板。长信号(一般为10秒的量级)由一个检波器接收,并且这些信号通过与输入信号的相关被处理,以把信号压缩到50-200毫秒的范围。因为在单独检波器处的接收信号在测量时间上几乎是相同的,连续震动法测量需要的竖直叠加理想地适用于协方差叠加。
总之,用于本发明的三种一般用途作为具体例子给出。
1.通过把检波器靠近放在一起(比如,高达约20米远)求和来自现场的检波器的痕迹,通过进行扰动校正和施加通常的时差校正除去在振幅和到达时间中的局部变化。
2.求和以后由同一振动器(阵列)发射和由同一检波器接收的几次扫描。这称作竖直叠加,并且原则上,不需要对接收信号的差进行校正。以1指示的处理能用来补偿在扫描之间、由例如振动器与地之间的较差接触造成的源强度差。由于以这种方式收集的数据具有非常差的信噪比,本发明的技术特别适用于该类型的数据,特别是在进行相关之前。
3.求和其反射信号源于地表下岩石上的相同反射点的发射-接收器对收集的痕迹(CMP收集)。诸痕迹由通常的时差已经校正。与情形2一样,可以选择应用情形1的处理。
在这些情形之间的共同点在于,痕迹在处理为适当的之后都具有非常类似的反射信号。
总之:
·对于估计的所有数据集,协方差叠加提供比反比叠加或直接痕迹求和更好的随机噪声衰减。
·QMF的应用改进协方差叠加的性能;改进量从非常小到约3dB变化。
·增大窗口尺寸会减小协方差叠加的性能。
为了减小加权因数计算时间,可以仅对于具体时间和/或频率计算加权因数,并且通过在这些计算的加权因数之间的内插,可以确定用于其他时间和/或频率的加权因数。例如,可以使用线性插入,以保证插入的加权因数总是添加到一。这种类型的方法也提供加权因数的平稳过渡,并且可以消除对上述加权因数平稳过程的需要。
上述实施例的改进在减弱地震痕迹中的尖峰方面特别有用。在地震痕迹中的尖峰例如可以归因于局部高振幅噪声(由车辆交通、脚步、电导墙网等引起的),或者归因于地震数据检测系统问题(由系统噪声、模数转换器误差、数据传输问题等引起)。对于多种类型的地震数据检测系统问题,“失去”或错误的地震数据值将分配零幅值,并且这些零幅值将产生与由上以讨论的尖峰噪声的常规源引入的问题非常类似的数据处理问题。
重要的是,在地震数据处理期间(并且最好在叠加处理期间)尽可能快和尽可能完全地从地震痕迹除去尖峰,因为当在频率和/或波数域中进行一定类型的运算时,在时域中短持续时间但高强度的尖峰成为污染到相邻样本或相邻痕迹中。一旦这种类型的噪声已经污染到相邻痕迹中的数据中,就变得使用常规噪声除去技术,如带通滤波、波数滤波、或F-F滤波,几乎不可能完全除去该噪声。
为了改进发明方法的除去尖峰的能力,可以把计算的加权因数分配到具有最大相对振幅的具体时间和/或频率。当在时间域中运算时,这可能涉及把计算的加权因数分配到与用于所有地震痕迹的窗口内的波峰信号振幅的样本有关的具体时间。这可能与以上讨论的把用于给定窗口的计算加权因数分配到窗口中心(即中点)相反。
如果地震痕迹窗口包含一个尖峰,为该具体地震痕迹窗口计算的加权因数与为其他地震痕迹窗口计算的加权因数相比较低。如果该尖峰位于窗口边缘之一附近,并且使用以上讨论的加权因数插入方案,则施加到具体尖峰样本或样本上的加权因数是用于尖峰位于其中的实际窗口的较小加权因数、和为最近相邻窗口导出的加权因数的平均值。然而,如果在高振幅尖峰的位置处精确地施加加权因数,则计算的较小的加权因数将直接施加到尖峰值本身上,因而更有效地衰减在输出痕迹中的尖峰。能明白,当最高振幅样本位于窗口的中点处时,该方法供给相同的结果,而当尖峰具有较大振幅并且更靠近窗口边缘之一定位时,该方法提供改进的结果。
在发明方法的通常用途中,白噪声将加到矩阵
(在公式NE2中描述的)对角线上,以增大加权因数确定过程的稳定性。添加的白噪声数值代表在加权因数确定过程期间在滤波器精度与滤波稳定性之间的折衷。
最好可能是也为计算的加权因数设置一个最大极限值,以防止弱或死地震痕迹(即具有非常低振幅的地震痕迹)不适当地影响叠加地震痕迹的组成。
Claims (16)
1.一种叠加多条地震痕迹的方法,包括:导出地震痕迹每一条的噪声分量的表示;由噪声分量的表示导出在噪声分量之间的相关指示;及通过按从相关指示导出的比例把多个地震痕迹相结合产生一个叠加地震痕迹。
2.根据权利要求1所述的方法,其中在噪声分量表示之间的相关指示是噪声分量的交叉能量密度。
3.根据权利要求2所述的方法,其中由使噪声分量交叉能量密度的微分最小导出该比例。
4.根据权利要求1至3任一项所述的方法,其中以迭代方式导出该比例。
5.根据权利要求1至4任一项所述的方法,其中把地震痕迹分离成多个时间窗口。
6.根据权利要求5所述的方法,其中时间窗口具有在50毫秒与1秒之间的长度。
7.根据权利要求1至6任一项所述的方法,其中在叠加之前把多条地震痕迹分离成在不同频率下的多条带。
8.根据权利要求7所述的方法,其中使用正交反射滤波器分离地震痕迹。
9.根据权利要求1至8任一项所述的方法,其中噪声分量的表示是噪声分量的噪声能量。
10.根据权利要求1至9任一项所述的方法,其中在地震痕迹中除去直流分量之后,导出该比例。
11.根据权利要求1至10任一项所述的方法,进一步包括导出在地震痕迹之间的相关系数,并且如果该相关系数大于一个阈值则把反比叠加用于痕迹。
12.根据权利要求11所述的方法,其中阈值在0.6与0.9之间。
13.根据权利要求1至12任一项所述的方法,其中由不同的传感器同时导出地震痕迹。
14.根据权利要求1至12任一项所述的方法,其中由连续震动法导出地震痕迹。
15.根据权利要求1至14任一项所述的方法,其中把由相关指示导出的比例分配到具体时间/或频率,并且通过在分配的导出比例之间内插来计算用于在这些具体时间和/或频率之间的中间点的比例。
16.根据权利要求15所述的方法,其中把所述导出的比例分配到其中所述地震痕迹具有最大相对振幅的具体时间和/或频率。
Applications Claiming Priority (2)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
GBGB9726928.6A GB9726928D0 (en) | 1997-12-19 | 1997-12-19 | Method of stacking seismic signals |
GB9726928.6 | 1997-12-19 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN1282424A true CN1282424A (zh) | 2001-01-31 |
CN1210579C CN1210579C (zh) | 2005-07-13 |
Family
ID=10823952
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CNB988122464A Expired - Fee Related CN1210579C (zh) | 1997-12-19 | 1998-12-18 | 叠加地震记录迹的方法 |
Country Status (9)
Country | Link |
---|---|
US (1) | US6535818B1 (zh) |
EP (1) | EP1040367A1 (zh) |
CN (1) | CN1210579C (zh) |
AU (1) | AU739793B2 (zh) |
CA (1) | CA2309140A1 (zh) |
EA (1) | EA002368B1 (zh) |
GB (2) | GB9726928D0 (zh) |
NO (1) | NO20003094L (zh) |
WO (1) | WO1999032903A1 (zh) |
Cited By (4)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN102478664A (zh) * | 2010-11-23 | 2012-05-30 | 中国石油天然气集团公司 | 一种有效信号无污染的空间采样间隔确定方法 |
CN104995534A (zh) * | 2013-05-02 | 2015-10-21 | 雪佛龙美国公司 | 增强三维地震解释中的平点的方法 |
CN106415322A (zh) * | 2014-02-13 | 2017-02-15 | 阿德洛克有限公司 | 识别反射信号的方法 |
CN107835962A (zh) * | 2015-05-13 | 2018-03-23 | 科诺科菲利浦公司 | 钻探数据的时间校正 |
Families Citing this family (18)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US6950371B2 (en) * | 2003-05-20 | 2005-09-27 | Chevron U.S.A. Inc. | Method for signal-to-noise ratio enhancement of seismic data using frequency dependent true relative amplitude noise attenuation |
US6912468B2 (en) * | 2003-08-14 | 2005-06-28 | Westerngeco, L.L.C. | Method and apparatus for contemporaneous utilization of a higher order probe in pre-stack and post-stack seismic domains |
US7230879B2 (en) | 2005-02-12 | 2007-06-12 | Chevron U.S.A. Inc. | Method and apparatus for true relative amplitude correction of seismic data for normal moveout stretch effects |
AU2006214688B2 (en) * | 2005-02-16 | 2010-06-17 | Exxonmobil Upstream Research Company | Estimating noise at one frequency by sampling noise at other frequencies |
JP4699113B2 (ja) * | 2005-07-05 | 2011-06-08 | パナソニック株式会社 | マルチキャリア通信装置、及びマルチキャリア通信方法 |
US7379386B2 (en) * | 2006-07-12 | 2008-05-27 | Westerngeco L.L.C. | Workflow for processing streamer seismic data |
US7869303B2 (en) * | 2007-08-14 | 2011-01-11 | Pgs Geophysical As | Method for noise suppression in seismic signals using spatial transforms |
US9348049B2 (en) * | 2012-01-05 | 2016-05-24 | Cgg Services Sa | Simultaneous joint estimation of the P-P and P-S residual statics |
US9341724B2 (en) | 2012-03-22 | 2016-05-17 | Westerngeco L.L.C. | Sweep sequence determination for overlapping sweeps |
US9726772B2 (en) * | 2013-10-17 | 2017-08-08 | Westerngeco L.L.C. | Iterative stacking of seismic image partitions |
WO2015083000A2 (en) * | 2013-12-05 | 2015-06-11 | Cgg Services Sa | Methods and systems of detecting a microseismic event using an iterative non-linear inversion algorithm |
US9857490B2 (en) * | 2013-12-30 | 2018-01-02 | Pgs Geophysical As | Methods and systems for optimizing generation of seismic images |
US10578759B2 (en) | 2014-02-10 | 2020-03-03 | Westerngeco L.L.C. | Quality control and preconditioning of seismic data |
WO2015128732A2 (en) | 2014-02-25 | 2015-09-03 | Cgg Services Sa | Subterranean formation monitoring using frequency domain weighted analysis |
CN106226818B (zh) * | 2016-04-01 | 2018-05-04 | 中国石油天然气股份有限公司 | 地震数据处理方法和装置 |
WO2020231287A1 (en) | 2019-05-13 | 2020-11-19 | Saudi Arabian Oil Company | Providing seismic images of the subsurface using enhancement of pre-stack seismic data |
CN113534259B (zh) * | 2021-07-09 | 2024-05-31 | 中石化石油工程技术服务有限公司 | 一种可控震源高效采集实时叠前时间偏移成像方法 |
CN113687421B (zh) * | 2021-08-23 | 2022-10-21 | 中国石油大学(北京) | 地震信号的数据处理方法、装置、电子设备及存储介质 |
Family Cites Families (9)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
OA02288A (fr) * | 1967-01-03 | 1970-05-05 | Texas Instruments Inc | Combinaison horizontale optimale de signaux d'exploration sismique. |
US3622967A (en) * | 1968-11-07 | 1971-11-23 | Mobil Oil Corp | Optimum stack |
US4561075A (en) * | 1982-12-29 | 1985-12-24 | Standard Oil Company | Method and apparatus for selecting inverse power weighting prior to vertically stacking acquired seismic data for suppressing noise |
US5392213A (en) * | 1992-10-23 | 1995-02-21 | Exxon Production Research Company | Filter for removal of coherent noise from seismic data |
US5572483A (en) * | 1995-07-25 | 1996-11-05 | Western Atlas International, Inc. | Method of reducing noise in seismic signals by adaptive filtering of a noise reference |
US5774417A (en) * | 1996-10-25 | 1998-06-30 | Atlantic Richfield Company | Amplitude and phase compensation in dual-sensor ocean bottom cable seismic data processing |
US5850622A (en) * | 1996-11-08 | 1998-12-15 | Amoco Corporation | Time-frequency processing and analysis of seismic data using very short-time fourier transforms |
GB9623945D0 (en) * | 1996-11-15 | 1997-01-08 | Geco Prakla Uk Ltd | Detection of ground roll cone |
US5978314A (en) * | 1997-03-21 | 1999-11-02 | Exxon Production Research Company | Method for determining seismic velocities |
-
1997
- 1997-12-19 GB GBGB9726928.6A patent/GB9726928D0/en not_active Ceased
-
1998
- 1998-12-18 EP EP98962557A patent/EP1040367A1/en not_active Withdrawn
- 1998-12-18 CN CNB988122464A patent/CN1210579C/zh not_active Expired - Fee Related
- 1998-12-18 GB GB9827995A patent/GB2332521B/en not_active Expired - Fee Related
- 1998-12-18 EA EA200000688A patent/EA002368B1/ru not_active IP Right Cessation
- 1998-12-18 CA CA002309140A patent/CA2309140A1/en not_active Abandoned
- 1998-12-18 WO PCT/GB1998/003819 patent/WO1999032903A1/en not_active Application Discontinuation
- 1998-12-18 US US09/581,417 patent/US6535818B1/en not_active Expired - Fee Related
- 1998-12-18 AU AU17693/99A patent/AU739793B2/en not_active Ceased
-
2000
- 2000-06-15 NO NO20003094A patent/NO20003094L/no not_active Application Discontinuation
Cited By (5)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN102478664A (zh) * | 2010-11-23 | 2012-05-30 | 中国石油天然气集团公司 | 一种有效信号无污染的空间采样间隔确定方法 |
CN102478664B (zh) * | 2010-11-23 | 2013-09-04 | 中国石油天然气集团公司 | 一种有效信号无污染的空间采样间隔确定方法 |
CN104995534A (zh) * | 2013-05-02 | 2015-10-21 | 雪佛龙美国公司 | 增强三维地震解释中的平点的方法 |
CN106415322A (zh) * | 2014-02-13 | 2017-02-15 | 阿德洛克有限公司 | 识别反射信号的方法 |
CN107835962A (zh) * | 2015-05-13 | 2018-03-23 | 科诺科菲利浦公司 | 钻探数据的时间校正 |
Also Published As
Publication number | Publication date |
---|---|
EA002368B1 (ru) | 2002-04-25 |
GB2332521A (en) | 1999-06-23 |
AU1769399A (en) | 1999-07-12 |
EP1040367A1 (en) | 2000-10-04 |
GB9827995D0 (en) | 1999-02-10 |
EA200000688A1 (ru) | 2000-12-25 |
GB2332521A9 (en) | |
NO20003094L (no) | 2000-08-11 |
GB2332521B (en) | 2000-03-08 |
WO1999032903A1 (en) | 1999-07-01 |
AU739793B2 (en) | 2001-10-18 |
GB9726928D0 (en) | 1998-02-18 |
CA2309140A1 (en) | 1999-07-01 |
US6535818B1 (en) | 2003-03-18 |
NO20003094D0 (no) | 2000-06-15 |
CN1210579C (zh) | 2005-07-13 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN1282424A (zh) | 叠加地震痕迹的方法 | |
CN1203324C (zh) | 自适应地震噪声和干扰衰减方法 | |
Wathelet et al. | Direct inversion of spatial autocorrelation curves with the neighborhood algorithm | |
Moya et al. | Inversion of source parameters and site effects from strong ground motion records using genetic algorithms | |
RU2116657C1 (ru) | Способ обработки сейсмических данных для подавления многократных переотражений | |
CN1503006A (zh) | 形成地下区域中代表一物理量分布的模型的方法 | |
CN1214255C (zh) | 齐发震动地震数据的处理 | |
Di Giulio et al. | Deriving wavefield characteristics and shear-velocity profiles from two-dimensional small-aperture arrays analysis of ambient vibrations in a small-size alluvial basin, Colfiorito, Italy | |
Langston et al. | The vertical component P-wave receiver function | |
RU2005108662A (ru) | Процедура поиска сигналов для системы определения местоположения | |
CN1101428A (zh) | 用于检测地震信号以便获得钻孔操作过程的垂直地震剖面图的方法和装置 | |
CN1832633A (zh) | 一种声源定位方法 | |
CN1331408A (zh) | 机动车检测装置和机动车检测方法 | |
US20030220744A1 (en) | Method for processing seismic data to attenuate multiples | |
WO1997003370A2 (en) | Calibrating vertical particle velocity and pressure detectors | |
CN1646941A (zh) | 海洋延时地震勘测 | |
EA001766B1 (ru) | Способ выделения сейсмических сигналов, а также определения и коррекции геометрических и статических ошибок в сейсмических данных | |
CN112068195B (zh) | 微地震p&s波匹配事件初至自动拾取方法和计算机存储介质 | |
CN104932015B (zh) | 构建地震数据的速度模型的方法 | |
JP3463677B2 (ja) | 震源位置の決定法 | |
EA005232B1 (ru) | Способ и устройство обработки сейсмических данных | |
CN1017094B (zh) | 带宽增大的地震采集系统 | |
CN1359473A (zh) | 组合三分量地震数据的方法和系统 | |
KAWASHIMA et al. | Attenuation of peak ground motions and absolute accleration response spectra of vertical earthquake ground motion | |
Kind et al. | Seismological evidence for a very sharp Sorgenfrei-Tornquist Zone in southern Sweden |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
C06 | Publication | ||
PB01 | Publication | ||
C10 | Entry into substantive examination | ||
SE01 | Entry into force of request for substantive examination | ||
C14 | Grant of patent or utility model | ||
GR01 | Patent grant | ||
C17 | Cessation of patent right | ||
CF01 | Termination of patent right due to non-payment of annual fee |
Granted publication date: 20050713 Termination date: 20111218 |