[发明专利]一种采用频变的复数相关系数估算瞬时慢度的方法

权利要求书

1.一种采用频变的复数相关系数估算瞬时慢度的方法，特点是具体实施步骤如下：

1）采集得到二维地震剖面；

2）对地震剖面进行频谱分析并确定地震剖面的主频f_c；对地震剖面进行F-K谱分析并确定地震剖面的最小和最大视速度v_min和v_max；

3）在时间方向，把地震道投影到以地震剖面的主频f_c为频率的富里叶变换的基函数上，得到投影后的地震剖面A(x,t,f_c)：

A(x,t,fc)=a(x,t)e-2πifct]]>

式中，a（x,t）是二维地震剖面，x和t分别为地震道的空间和时间坐标，f_c为地震剖面的主频；A(x,t,f_c)为投影后的地震剖面；

4）对投影后的地震剖面A(x,t,f_c)先按空间道序再按时间采样点顺序逐点按步骤5）-18）估算瞬时慢度；

5）根据步骤2)确定的地震剖面的主频f_c，用下式计算以采样点为中心的沿时间方向的处理时窗的长度：

wt=1fcΔt]]>

式中，f_c为地震剖面的主频，Δt为地震剖面的时间采样率，w_t为处理时窗的长度；

6）根据步骤2)确定的地震剖面的最小和最大视速度v_min和v_max，计算出慢度扫描的最少值和最大值；

p_min=1/v_max

p_max=1/v_min

式中，p_min和p_max表示慢度扫描的最少值和最大值，v_min和v_max是由步骤2）确定的地震剖面的最小和最大视速度；

7）根据慢度扫描个数，计算慢度扫描的增量：

dp=（p_max-p_min）/（n_p-1）

式中，n_p为用户给出的慢度扫描个数；dp为慢度扫描的增量；p_min和p_max表示慢度扫描的最少值和最大值；

8）以当前处理地震道为中心道，确定参与瞬时慢度估算的相邻地震道数w_x；

9）从慢度扫描的最少值开始到慢度扫描的最大值为止，按慢度扫描增量计算扫描慢度p_i：

p_i=（i-1）*dp（i=1,2,…,n_p）

其中，i表示扫描慢度的编号，p_i是第i个扫描慢度，dp表示由步骤7）计算的慢度增量；

10）以当前处理地震道和处理时窗的中心点为基点，对每一扫描慢度按照下式计算其对应的相邻地震道上的处理时窗的中心点时间τ_i,j和时窗对应的信号的延迟时Δτ_i,j；

τ_i,j=τ₀+（x_j-x₀）*p_i  （-w_x/2≤j≤w_x/2,i=1,2…..n_p）

Δτ_i,j=（x_j-x₀）*p_i（-w_x/2≤j≤w_x/2，i=1,2…..n_p）

式中，i表示扫描慢度的编号；j表示相邻地震道的编号，0编号表示当前处理道；p_i表示第i个扫描慢度。x₀表示当前处理地震道的空间作标，x_j表示第j个相邻地震道的空间作标；τ₀表示当前处理地震道的处理时窗的中心点时间；τ_i,j表示当用慢度p_i进行扫描时第j个相邻地震道的处理时窗的中心点时间；Δτ_i,j表示当用慢度p_i进行扫描时第j个相邻地震道与当前处理道的处理时窗中心点之间的时间延迟时；w_x表示由步骤8）定义的相邻地震道数；n_p表示由步骤7）定义的慢度扫描个数；

11）对每一扫描慢度根据步骤10）计算的当前道及其相邻道的处理时窗的中心点时间以及步骤5）定义的时窗长度，提取当前道及其相邻道的处理时窗内的采样点的数值

B~j,n(pi)=A(xj,mi,j+n,fc)]]>

m_i,j=τ_i,j/Δt

1≤i≤n_p，-w_x/2≤j≤w_x/2，-w_t/2≤n≤w_t/2

式中，i表示扫描慢度的编号；j表示相邻地震道的编号，0序号表示当前处理道；n表示处理时窗内的采样点的编号，0编号表示处理时窗的中心采样点；p_i表示第i个扫描慢度；τ_i,j表示当用慢度pi进行扫描时第j个相邻地震道的处理时窗的中心点时间；Δt为地震剖面的时间采样率；m_i,j表示当用慢度p_i进行扫描时第j个相邻地震道的处理时窗的中心点对应的采样点号；A(x_j,m_i,j+n,fc)表示当用慢度p_i进行扫描时第j个相邻地震道的处理时窗内的编号为n的采样点的数值；表示当用慢度p_i进行扫描时从第j个相邻地震道的处理时窗内提取的编号为n的采样点的数值；w_x表示步骤8）定义的相邻地震道数；n_p表示步骤7）定义的慢度扫描个数。w_t表示步骤5）定义的时窗长度；f_c表示由步骤2）确定的地震剖面的主频；

12）对每一扫描慢度，根据10）计算的相邻道处理时窗对应的信号延迟时，对提取的处理时窗内的采样点的数值进行相位校正；

Bj,n(pi)=B~j,n(pi)e-2πifcΔτi,j]]>

式中，i表示扫描慢度的编号，j表示相邻地震道的编号，0序号表示当前处理道；n表示处理时窗内的采样点的编号，0编号表示处理时窗的中心采样点；p_i表示第i个扫描慢度；表示当用慢度p_i进行扫描时从第j个相邻地震道的处理时窗内提取的编号为n的采样点的数值；B_j,n(p_i)表示经过相位校正后的当用慢度p_i进行扫描时从第j个相邻地震道的处理时窗内提取的编号为n的采样点的数值；Δτ_i,j表示当用慢度p_i进行扫描时第j个相邻地震道与当前处理道的处理时窗中心点之间的时间延迟时；f_c表示由2）确定的地震剖面的主频；

13）对每一扫描慢度，计算相位校正后的每个相邻地震道的处理时窗内的信号的能量E_j(p_i)：

Ej(pi)=Σn=-wt/2wt/2Bj,n(pi)B‾j,n(pi)]]>

式中，i表示扫描慢度的编号，j表示相邻地震道的编号，0序号表示当前处理道；n表示处理时窗内的采样点的编号，0编号表示处理时窗的中心采样点；p_i表示第i个扫描慢度；B_j,n(p_i)表示经过相位校正后的当用慢度p_i进行扫描时从第j个相邻地震道的处理时窗内提取的编号为n的采样点的数值；为B_j,n(p_i)的共轭；E_j(p_i)表示当用慢度p_i进行扫描时第j个相邻地震道处理时窗内信号的能量；w_t表示由步骤5）定义的时窗长度；

14）对每一扫描慢度,计算参与慢度估算的各地震道处理时窗内的信号的复数互相关系数C_j,k(p_i)；

Cj,k(pi)=Σn=-wt/2wt/2Bj,n(pi)B‾k,n(pi)/Ej(pi)Ek(pi)]]>

-w_x/2≤j≤w_x/2，k>j

式中，i表示扫描慢度的编号；j和k表示相邻地震道的编号，0序号表示当前处理道；n表示处理时窗内的采样点的编号，0编号表示处理时窗的中心采样点；p_i表示第i个扫描慢度；B_j,n(p_i)表示经过相位校正后的当用慢度p_i进行扫描时从第j个相邻地震道的处理时窗内提取的编号为n的采样点的数值；B_k,n(p_i)表示经过相位校正后的当用慢度p_i进行扫描时从第k个相邻地震道的处理时窗内提取的编号为n的采样点的数值；表示B_k,n(p_i)的共轭；E_j(p_i)表示当用慢度p_i进行扫描时第j个相邻地震道处理时窗内的信号的能量；E_k(p_i)表示当用慢度p_i进行扫描时第k个相邻地震道处理时窗内信号的能量；C_j,k(p_i)表示当用慢度p_i进行扫描时第j个和第k个相邻地震道的处理时窗内信号的复相关系数；w_t表示由步骤5）定义的时窗长度；w_x表示由步骤8）定义的相邻地震道数；

15）对每一扫描慢度,把其对应的各地震道处理时窗内信号的复数互相关系数进行求和，得到该慢度的总复数相关系数C(p_i)：

C(pi)=Σj=-wx/2wx/2Σk>jwx/2Cj,k(pi)]]>

式中，i表示扫描慢度的编号；j和k表示相邻地震道的编号，0序号表示当前处理道；p_i表示第i个扫描慢度；C_j,k(p_i)表示当用慢度p_i进行扫描时第j个和第k个相邻地震道的处理时窗内信号的复相关系数；C(p_i)表示当用慢度p_i进行扫描时所有相邻道的处理时窗内的信号的复数相关系数的和；w_x表示由步骤8）定义的相邻地震道数；

16）对每一扫描慢度,计算该慢度对应的总复数相关系数的实部和虚部；

RC(p_i)=real(C(p_i))

IC(p_i)=image(C(p_i))

式中，i表示扫描慢度的编号；p_i表示第i个扫描慢度；C(p_i)表示当用慢度p_i进行扫描时所有相邻道的处理时窗内的信号的复数相关系数的和；RC(p_i)表示C(p_i)的实部；IC(p_i)表示C(p_i)的虚部；

17）根据步骤16）的计算结果，在步骤6）慢度扫描的最小值和最大值之间寻找使复数相关系数虚部为零的绝对值最小的扫描慢度q₁：

IC(q₁)=0

式中，IC(p_i)表示C(p_i)的虚部。C(p_i)表示当用慢度p_i进行扫描时所有相邻道的处理时窗内的信号的复数相关系数的和；q₁表示使复数相关系数虚部为零的绝对值最小的扫描慢度；

18）在慢度0和使复数相关系数虚部为零的绝对值最小的零点之间，寻找使复数相关系数实部为最大的扫描慢度q₂，得到当前处理点的瞬时慢度；

RC(q₂)=maxRC(p_i)

式中，RC(p_i)表示C(p_i)的虚部；C(p_i)表示当用慢度p_i进行扫描时所有相邻道的处理时窗内的信号的复数相关系数的和，q₂表示使复数相关系数实部为最大的扫描慢度。

2.根据权利要求1所述的方法，特点是步骤8）所述的相邻地震道数w_x少于9大于3。