地震勘探资料采集中精确定位干扰源位置的方法

摘要

本发明公开了一种地震勘探资料采集中精确定位干扰源位置的方法，包括如下步骤：在地震勘探野外采集现场，记录外界干扰源产生的干扰波；提取干扰波同一波前面对应的时间；计算三维空间中一个点震源即干扰源到观测点的距离、方位，实现对外界干扰源位置的定位。本发明通过精确定位干扰源位置，从而实现快速控制采集时外界干扰源产生的噪音，提高资料品质和施工效率。

说明书

技术领域

本发明涉及一种地震勘探资料采集中精确定位干扰源位置的方法，属于地震 勘探资料采集技术领域。

背景技术

长期以来，地震勘探资料采集都受到来自外界干扰源产生的环境噪音的影响。 定位外界干扰源位置只能采取人工实地寻访的方法，一个20条测线(测线间距 500m)、测线总长度200Km的二维地震勘探项目，或勘探面积为100Km²三维地 震勘探项目，采取人工实地寻访干扰源，需要排出15人次、5车次在崎岖陡峭的 山路上寻访，耗时一周左右方能完成；若工区交通不便完成时间则会相应更长。 这种方法不仅效率低、成本高，且存在明显的疏漏，难于发现断续工作或临时性 作业新增干扰源，如为灌溉临时启动的抽水机、工矿里间歇运行的机器等干扰。 那些疏漏的干扰源势必会影响采集资料的质量，正常采集作业不得不面对临时中 断，临时协停干扰，这不仅延长了正常采集作业时间，工期延后，降低了采集施 工效率；而且需要增加作业人员的数量、加大了车辆使用频率，造成人力、物力 的浪费。

在采集过程中，为了准确定位外界干扰源位置，往往要花费大量人力、物力 和时间，极大地制约了采集施工效率，同时外界干扰也严重影响了资料品质，但 目前没有一种定位外界干扰源位置的技术能解决采集施工过程面临的难题。

发明内容

本发明的目的在于克服现有地震勘探资料采集时定位外界干扰源位置存在的 上述问题，提供一种地震勘探资料采集中精确定位干扰源位置的方法，本发明通 过精确定位干扰源位置，从而实现快速控制采集时外界干扰源产生的噪音，提高 资料品质和施工效率。

为实现上述目的，本发明采用的技术方案如下：

一种地震勘探资料采集中精确定位干扰源位置的方法，其特征在于，包括如 下步骤：

a、在地震勘探野外采集现场，记录外界干扰源产生的干扰波；

b、提取干扰波同一波前面对应的时间；

c、计算三维空间中一个点震源即干扰源到观测点的距离、方位，实现对外界 干扰源位置的定位。

所述a步骤中，采用三维观测系统，布设多条排列进行观测，利用地震仪器 记录外界干扰源产生的干扰波。

所述b步骤中，波前是指波速分界面上，各质点开始震动的点的面，或者某 一时刻介质中各点开始振动形成的曲面为波前；在弹性介质中，已知T时刻的同 一波前面上的各个点，看作各个点从该时刻产生子波的新的点震源，子波的包络 面就是波T+ΔT时刻到达的新的波前面。

所述b步骤中，同一波前面对应的时间是指干扰源多次震动之中的某一次或 同一次震动产生的球面波，经过任何一个ΔT时刻后，被观测和记录。

所述b步骤中，外界干扰源产生干扰波，以干扰源为圆心，在弹性介质中以 球面波形式向前传播，在不同的地震道拾取干扰波同一波前面相对应的相同时间。

所述c步骤中，基于相同时间对应的地震道坐标拟合平面直角坐标系圆曲线， 采用最小二乘法拟合圆曲线，求取圆心坐标，求出三维空间中一个点震源到观测 点的距离、方位，实现外界干扰源位置的精确定位。

所述c步骤中，计算方法为：

最小二乘法拟合圆曲线：R²＝(x-A)²+(y-B)²

R²＝x²-2Ax+A²+y²-2By+B²

令

a＝-2A

b＝-2B

c＝A²+B²-R²

得圆曲线方程的另一个形式：

x²+y²+ax+by+c＝0    (1)

求出参数a，b，c，求得圆心半径的参数：

$A=\frac{a}{-2}$

$B=\frac{b}{-2}$

$R=\frac{1}{2}\sqrt{a^2+b^2-4c}$

平面直角坐标系中，设圆心坐标为O(x1，y1)，圆上一点P(x2，y2)，已知两点 坐标，计算两点距离及两点方向的方位角；

计算公式为：

S12＝sqr((x2-x1)²+(y2-y1)²)    (7)

A12＝arcsin((y2-y1)/S12)       (8)

S12为圆心至圆上一点P的距离；

A12为圆心至圆上一点P的坐标方位角；

或者A12＝arccos(Δx21/S12)*sgn(Δy21)+360°    (9)

式中sgn()为取符号函数；

由式(1)拟合圆，据式(7)、(8)或(9)求出三维空间中一个点震源到观 测点的距离、方位，实现外界干扰源位置的精确定位。

采用本发明的优点在于：

一、本发明通过精确、高效地定位干扰源位置，从而实现快速控制采集时外 界干扰源产生的噪音，保障了资料品质，提高施工效率。

二、与现有的人工实地寻方的方法相比，本发明投入的人力、物力和时间较 少，因而极大地降低了成本。

三、本发明有利于发现断续工作或临时性作业新增干扰源，探测干扰源不存 在任何的疏漏和盲区。

四、正常采集作业不会因新增干扰源临时中断，避免了临时协停干扰造成的 正常采集作业时间延长，确保了工期。

五、避免了临时协停干扰需要增加作业人员，降低了车辆使用频率，提高了 人力、物力的使用能效。

具体实施方式

以下对本发明作展开说明：

一种地震勘探资料采集中精确定位干扰源位置的方法，包括如下步骤：

a、在地震勘探野外采集现场，采用三维观测系统，布设多条排列进行观测， 利用地震仪器记录外界干扰源产生的干扰波；

b、提取干扰波同一波前面对应的时间。波前是波速分界面上，各质点开始震 动的点的面；或者说，某一时刻介质中各点刚好开始振动，这一曲面叫波前，也 叫波阵面。在弹性介质中，已知T时刻的同一波前面上的各个点，可以看作各个 点从该时刻产生子波的新的点震源，这些子波的包络面就是波T+ΔT时刻到达的 新的波前面。

外界干扰源产生干扰波，以干扰源为圆心，在弹性介质中以球面波形式向前 传播，同一波前面对应的时间是指干扰源多次震动之中的某一次(或同一次)震 动产生的球面波，经过任何一个ΔT时刻后，被观测并记录下来。

依据干扰波的动力学特征，主要是能量特性，一般包括振幅、频率、相位和 衰减等特性进行定量分析，正确识别干扰源多次震动之中的同一次震动产生的球 面波波前面。

依据地震波以球面波形式进行传播的理论，同一波前面对应的时间值具有等 时性。在不同的地震道拾取干扰波同一波前面相对应的相同时间。

c、计算三维空间中一个点震源(干扰源)到观测点的距离、方位。

基于相同时间对应的地震道坐标拟合平面直角坐标系圆曲线。采用最小二乘 法拟合圆曲线，求取圆心坐标。

最小二乘法(least squares analysis)是一种数学优化技术，它通过最小化误差 的平方和找到一组数据的最佳函数匹配。最小二乘法是用最简的方法求得一些绝 对不可知的真值，而令误差平方之和为最小。最小二乘法通常用于曲线拟合(least  squares fitting)。

最小二乘法拟合圆曲线：R²＝(x-A)²+(y-B)²

R²＝x²-2Ax+A²+y²-2By+B²

令

a＝-2A

b＝-2B

c＝A²+B²-R²

可得圆曲线方程的另一个形式：

x²+y²+ax+by+c＝0    (1)

只要求出参数a，b，c就可以求得圆心半径的参数：

$A=\frac{a}{-2}$

$B=\frac{b}{-2}$

$R=\frac{1}{2}\sqrt{a^2+b^2-4c}$

最小二乘法(least squares analysis)拟合圆曲线的公式推导如下：

样本集(X_i，Y_i)i∈(1，2，3...N)中点到圆心的距离为d_i：

d_i²＝(X_i-A)²+(Y_i-B)²

点(X_i，Y_i)到圆边缘的距离的平方与和半径平方的差为：

δ_i＝d_i²-R²＝(X_i-A)²+(Y_i-B)²-R²＝X_i²+Y_i²+aX_i+bY_i+c

令Q(a，b，c)为δ_i的平方和：

Q(a，b，c)＝∑δ_i²＝∑[(X_i²+Y_i²+aX_i+bY_i+c)]²

求参数a，b，c使得Q(a，b，c)的值最小值；

解：

平方差Q(a，b，c)大于0，因此函数存在大于或等于0的极小值，极 大值为无穷大。

F(a，b，c)对a，b，c求偏导，令偏导等于0，得到极值点，比较所有 极值点的函数值即可得到最小值。

$\frac{\partial Q(a,b,c)}{\partial a}=\Sigma 2({X_i}^2+{Y_i}^2+a{X}_i+b{Y}_i+c)X_i=0---\left(2\right)$

$\frac{\partial Q(a,b,c)}{\partial b}=\Sigma 2({X_i}^2+{Y_i}^2+a{X}_i+b{Y}_i+c)Y_i=0---\left(3\right)$

$\frac{\partial Q(a,b,c)}{\partial c}=\Sigma 2({X_i}^2+{Y_i}^2+a{X}_i+b{Y}_i+c)=0---\left(4\right)$

解这个方程组。

先消去c

(2)*N-(4)*∑X_i得：

N∑(X²+Y²+aX+bY+c)X-∑(x²+Y²+aX+bY+c)×∑X＝0

N∑(X²+Y²+aX+bY)X-∑(X²+Y²+aX+bY)×∑X＝0

(N∑X_i²-∑X_i∑X_i)a+(N∑X_iY_i-∑X_i∑Y_i)b

                                      (5)

+N∑X_i³+N∑X_iY_i²-∑(X_i²+Y_i²)∑X_i＝0

(3)*N-(4)*∑Y_i得：

N∑(X²+Y²+aX+bY+c)Y-∑(X²+Y²+aX+bY+c)×∑Y＝0

N∑(X²+Y²+aX+bY)Y-∑(X²+Y²+aX+bY)×∑Y＝0

(N∑X_iY_i-∑X_i∑Y_i)a+(N∑Y_i²-∑Y_i∑Y_i)b

                                          (6)

+N∑X_i²Y_i+N∑Y_i³-∑(X_i²+Y_i²)∑Y_i＝0

令

C＝(N∑X²-∑X∑X)

D＝(N∑XY-∑X∑Y)

E＝N∑X²+N∑XY-∑(X²+Y²)∑X  

G＝(N∑Y²-∑Y∑Y)

H＝N∑X²Y+N∑Y²-∑(X²+Y²)∑Y 

可解得：|

Ca+Db+E＝0

Da+Gb+H＝0

$a=\frac{\mathrm{HD}-\mathrm{EG}}{\mathrm{CG}-D^2}$

$b=\frac{\mathrm{HC}-\mathrm{ED}}{D^2-\mathrm{GC}}$

$c=-\frac{\Sigma ({X_i}^2+{Y_i}^2)+\mathrm{a\Sigma }{X}_i+\mathrm{b\Sigma }{Y}_i}{N}$

得A、B、R的估计拟合值：

$A=\frac{a}{-2}$

$B=\frac{b}{-2}$

$R=\frac{1}{2}\sqrt{a^2+b^2-4c}$

平面直角坐标系中，假设圆心坐标为O(x1，y1)，圆上一点P(x2，y2)，已知两 点坐标，计算两点距离及两点方向的方位角。

计算公式为：

S12＝sqr((x2-x1)²+(y2-y1)²)    (7)

A12＝arcsin((y2-y1)/S12)       (8)

S12为圆心至圆上一点P的距离；

A12为圆心至圆上一点P的坐标方位角。

按公式A12＝arcsin((y2-y1)/S12)计算出的方位角都要进行象限判断后加常数 才是真正的方位角。

也可通过下述计算公式进行：

A12＝arccos(Δx21/S12)*sgn(Δy21)+360°    (9)

式中sgn()为取符号函数，该公式只需加上条件(A12＞360°，A12＝A12-360 °)就可以计算出坐标方位角，不需要进行象限判断。

由式(1)拟合圆，据式(7)、(8)或(9)可求出三维空间中一个点震源(干 扰源)到观测点的距离、方位，从而实现外界干扰源位置的精确定位。