[发明专利]基于时频域目标函数的盐下构造最小二乘逆时偏移方法

权利要求书

1.一种基于时频域目标函数的盐下构造最小二乘逆时偏移方法，其步骤是：

步骤1.使用MATLAB软件及其地震数据处理软件Crews工具包对地震数据进行预处理，并将预处理后的地震数据作为时频域最小二乘逆时偏移的数据输入；

步骤2.构建初始速度模型，作为波动方程正演模拟的模型输入；

步骤3.定义观测系统，提取震源子波；利用初始速度模型进行波动方程正演模拟，存储正传波场，并与观测数据的反传波场互相关，获得逆时偏移成像结果；在速度模型上正演模拟对应的声波方程为：

其中x,z为坐标轴，s表示慢度场，u为地震波场，f为震源子波，t为时间；

步骤4.计算扰动介质产生的扰动地震波场，扰动介质表示为：

将慢度场的平方s²表示为背景慢度的平方与扰动慢度平方Δs²的叠加；对应的总地震波场同样分解为背景地震波场u₀和扰动地震波场u_s的叠加，即：

u＝u₀+u_s

其中背景地震波场u₀与总地震波场u都满足波动方程：

当扰动地震波场u_s＜＜u₀时，用背景地震波场代替总地震波场u₀≈u₀+u_s，即：

步骤5.利用Gabor变换获得观测数据和模拟数据的时频域振幅相位信息；其中观测数据和模拟数据Gabor变换为：

其中h(τ-t)表示高斯窗函数；t表示时间；τ为中间变量；ω为角频率，u(t)和d(t)分别表示时间域模拟数据和观测数据；和表示时频域模拟数据和观测数据，F_h[·]表示对地震数据作用的Gabor变换算子；

步骤6.构建基于时频域振幅相位的最小二乘逆时偏移目标函数：

其中i表示虚数；ns为震源数目；nr为检波器数目；模型扰动和表示模拟数据和观测数据在时频域的相位信息；和表示观测数据和模拟数据在时频域振幅信息；ε∈[0,1]表示权重因子，用来控制目标函数中振幅和相位的比例；当ε＝0时，目标函数为纯相位目标函数，完全忽略了波形振幅的影响；当ε＝1时，目标函数将退化为常规最小二乘目标函数；求目标函数关于模型参数的偏导数，获得目标函数对模型扰动的更新量：

步骤7.利用Gabor反变换将时频域的梯度表达式转换到时间域中，其中为了对公式做进一步的简化，地震观测数据和模拟数据Gabor逆变换表示为：

其中表示Gabor逆变换算子；因此，目标函数的梯度进一步在时空域简化为：

步骤8.定义伴随震源，并将其反传至模型空间；对应的伴随震源为：

步骤9.利用零延迟互相关方法获得最小二乘逆时偏移的梯度；根据波动方程的表达式，

Born正演模拟用矩阵的形式来表示：

A_su_s＝A₀u₀

其中计算其关于模型扰动的偏导数：

同样将梯度写成矩阵的形式为：

步骤10.利用L-BFGS优化算法计算模型扰动的更新方向；迭代公式如下：

m_k+1＝m_k-α_kH_kg_k

其中m_k为第k步模型扰动的更新量，α_k为步长，H_k为近似海森矩阵的逆，为模型扰动的梯度；在L-BFGS优化更新中，只需要保存少数的向量对，用于更新Hessian矩阵，其更新公式如下：

H_k+1＝V_k^TH_kV_k+ρ_ks_ks_k^T

V_k＝I-ρ_ky_ks_k^T

s_k＝m_k+1-m_k，y_k＝g_k+1-g_k

其中H_k+1是根据向量对{s_k,y_k}和H_k计算得到；H_kg_k的乘积通过梯度g_k与向量对{s_k,y_k}之间一系列向量的内积与向量的和来获得；其中的近似Hessian矩阵的逆矩阵H_k需满足以下更新公式：

步骤11.判断是否满足终止条件，若满足则输出基于时频域目标函数的最小二乘逆时偏移结果；若不满足终止条件，将当前的成像结果继续作为下一次循环的初始扰动模型，直到满足终止条件。

2.根据权利要求1所述的基于时频域目标函数的盐下构造最小二乘逆时偏移方法，其特征是：采用有限差分法进行波动方程正演模拟。