本发明是一种基于指数衰减的时频域相位反演方法，具体是对地震数据进行短时傅里叶变换处理，提取时频域相位信息，并加入指数衰减项，进行地震波速度建模的一种方法。本方法可以缓解全波形反演中存在的主要问题：因缺失低频信息导致的周期跳跃，并为后续的常规全波形反演提供一个高精度速度建模结果。
背景技术：
1、20世纪80年代，lailly和tarantola将clearbout提出的偏移成像原理重塑为一个最小二乘局部优化问题，率先提出了全波形反演(fwi)理论。常规fwi主要通过l2范数目标函数拟合模拟数据和观测数据，使其残差最小并更新相应的地下介质参数，以获得高精度的反演结果。同时，由于直接求取fréchet导数计算量庞大，tarantola指出可以利用伴随状态法直接求取速度模型的梯度，即正演模拟波场与残差反传波场的零延迟互相关，有效降低了时间域fwi的耗时。然而，由于常规基于l2范数目标函数的fwi的强非线性特征，导致其反演精度在很大程度上取决于地震数据的低频信息、初始模型的选取和震源估计准确与否等。但是实际地震勘探中难以获得4hz以下的低频数据，且低频数据易受到噪声污染。这就会导致常规fwi梯度朝着错误的方向进行更新，出现严重的周期跳跃问题，无法获得高精度的反演结果。
2、为了克服常规fwi目标函数的强非线性特征，许多方法被相继提出。bunks于1995年提出使用地震数据低频到高频的多尺度反演策略，采用分频的反演策略可以有效缓解全波形反演中的周期跳跃问题。与之类似，在pratt率先在频率域中应用fwi方法，频率域全波形反演有分频多尺度反演的天然优势，仅需要离散的频率值就可以反演至高分辨率速度模型。luo等提出了基于波动方程的走时反演方法；使用互相关函数计算观测数据和模拟数据的走时差，避免手动拾取地震事件的走时，相比于走时层析成像可以有效提高走时反演的精度。然而，fwi对于低频数据的高度依赖和野外实际采集数据质量较低等问题，严重限制了fwi在实际勘探中的应用。地震数据低频重建是缓解fwi周期跳跃的有效途径，等提出地震数据的包络方法，指出包络信号可以有效降低fwi目标函数的非线性，获得地下速度模型的大尺度构造信息。wu等在数字信号理论中详细分析了利用包络可以解调出地震数据中缺失的低频信号，进而反演平滑的地下大尺度构造，并提出基于l2范数目标函数的包络平方、包络对数的目标函数。xiong等提出了改进的包络反演方法，详细分析了传统的包络反演所存在的问题，即包络伴随震源对反射信息的压制，导致伴随震源中丢失了无法匹配的反射波信息，改进的包络反演很好地解决了这一问题。oh和alkhalifah等结合信号包络和全局互相关目标函数，在降低振幅误差和常规l2范数目标函数强非线性特征的同时，利用互相关包络反演模型的宏观构造，并成功应用于北海obc数据。
3、相比于振幅信息，相位信息包含着地震波传播的运动学特征，能够揭示地震波的传播路径和地下介质等信息。sun和schuster提出了时间域相位反演目标函数，忽略地震数据振幅的变化，侧重于利用地震波相位信息并进行优化反演。相位反演可以在初始模型存在较大误差的情况下，通过迭代和优化逐步修正和完善得到较为可靠的反演结果。kim将时间域相位反演方法推广至频率域，构建频率域纯相位目标函数，基于频率域反演框架的相位反演方法可以分频反演并且更高效的得到较好的初始模型。luo等利用hilbert变换提取地震信号的瞬时相位信息，并用于时间域波形反演。hu等发展了时频域加权的相位互相关波形反演方法，利用短时傅里叶变换可以充分考虑地震数据的局部特性，并得到较好的反演结果，验证了相位反演的巨大优势。
4、由于信号相位信息的固有特性，纯相位目标函数容易受到噪声影响，难以获得地下速度模型的长波长结构。为此，本发明在时频域相位反演的研究基础上，提出基于指数衰减的时频域相位反演方法。该方法可以很好的压制时频域相位反演带来的噪声影响，提高反演精度。在反演过程中，通过不断调整衰减因子的值以达到分层多尺度的反演效果，压制噪声干扰，提高反演精度，并为后续常规fwi提供了较好的初始速度建模。
技术实现思路
1、本发明的主要技术创新点为：通过对地震数据进行短时傅里叶变换提取时频域相位信息，并加入衰减因子构建基于指数衰减的时频域相位反演目标函数，获得地下速度模型高精度反演结果。该方法可以压制时频域相位反演目标函数带来的噪声干扰，在反演过程中调整衰减因子，从模型浅部反演至模型深部，形成分层多尺度策略，以恢复速度模型的低波数分量。
2、本发明是基于指数衰减的时频域相位反演方法，主要步骤为：
3、步骤1.使用crews工具包对采集的地震数据进行预处理，将处理后的地震数据作为观测数据输入至基于指数衰减的时频域相位反演方法中。
4、步骤2.根据工区地质学信息，测井信息和走时信息等，设定反演的初始速度模型，作为基于指数衰减的时频域相位反演方法的初始模型。
5、步骤3.读取观测系统，输入震源子波。利用初始速度模型进行正演模拟，存储每一个采样时间节点的波场值用于后续互相关计算速度模型的更新梯度。
6、为了简化整个物理过程，加快反演效率，本次研究基于常密度声波方程，
7、
8、其中u表示波场值，z,x为纵向和横向坐标，t为时间，f(t)为震源子波。
9、步骤4.对地震数据的时频变换。对地震数据进行短时傅里叶变换处理，本次研究采用gabor变换。对时间域数据的gabor变换定义为：
10、
11、u(τ)表示时间域地震数据，h表示高斯窗函数，对应的t和ω表示时频谱中时间和频率的坐标，表示gabor后地震数据的时频域表达，fh[·]表示对数据的gabor变换算子。gabor逆变换定义为：
12、
13、其中表示对数据的gabor逆变换算子。
14、步骤5.计算地震信号时频域相位信息，地震数据时频域信号可以表示为：
15、
16、则时频域地震数据的指数相位表示为：
17、
18、其中表示地震数据的时频域包络。
19、步骤6.构建基于指数衰减的时频域相位反演目标函数。相位信息包含着地震波传播的运动学特征，能够揭示地震波的传播路径和地下介质等信息。时频域相位信息更关注地震信号的局部相位变化。此外，在时频域纯相位目标函数中添加衰减因子，可以有效缓解目标函数的非线性特征，为常规全波形反演提供一个较好的初始速度模型，缓减反演过程中的周期跳跃问题。
20、常规fwi的目标函数可以表示为：
21、
22、其中u和d分别表示模拟地震数据和实际观测数据。ns和nr分别表示震源炮点数和检波器个数，v表示速度模型参数。
23、时频域相位反演目标函数为：
24、
25、和分别表示模拟数据和观测数据的时频域相位。
26、基于指数衰减的时频域相位反演目标函数可以表示为：
27、
28、其中a为引入的衰减因子
29、步骤7.计算目标函数关于模型速度参数的偏导数和伴随震源。首先是常规fwi对于速度参数的偏导数可以表示为：
30、
31、则对应的伴随震源为rs＝(u-d)，即为残差记录并用于反传至模型空间。时频域相位反演目标函数对速度参数的偏导数可以表示为：
32、
33、其中re[·]表示取地震数据的实部，*表示复共轭。根据链式法则进一步化简可得：
34、
35、其中可进一步化简得到：
36、
37、其中表示gabor逆变换算子，则常规时频域相位反演目标函数的伴随震源可以表示为：
38、
39、依据上述推导，基于指数衰减的时频域相位反演目标函数对速度参数的偏导数可以表示为：
40、
41、其中re[·]表示取地震数据的实部，*表示复共轭。根据链式法则进一步化简可得：
42、
43、其中可进一步化简为：
44、
45、其中表示gabor逆变换算子，基于指数衰减的时频域相位反演目标函数的伴随震源可以表示为：
46、
47、步骤8.计算速度模型更新的方向。常规fwi目标函数的梯度可以用正演模拟波场与残差反传波场的零延迟互相关直接求得，表示为：
48、
49、其中nt表示时间域采样数，表示正演模拟波场对时间的二阶偏导，ub表示残差反传波场。
50、步骤9.利用优化算法对速度参数进行更新迭代。本发明选用最速下降方法，更新迭代公式可以表示为：
51、
52、其中αk表示更新步长，vk+1表示通过上一步速度模型vk更新得到的速度参数。
53、表示速度模型的更新梯度。
54、步骤10.输出反演结果。判断每一步的反演结果是否满足对目标函数设置的终止条件，若满足，输出反演结果；若不满足，将当前的反演结果作为下一次反演的初始速度模型，继续迭代直至满足终止条件；反演过程以目标函数作为约束条件，迭代过程需要满足目标函数值连续下降。
55、本发明提出了一种基于指数衰减的时频域相位反演方法，用于获得地下高精度速度建模。该方法可以压制常规时频域相位反演方法带来的相位噪声干扰，提高速度建模的精度。首先将地震数据变换至时频域，构建时频域相位反演目标函数，并引入衰减因子，以达到分层多尺度反演的效果。然后，推导基于指数衰减的时频域相位反演目标函数的梯度和伴随震源。最后，使用局部优化算法更新速度模型，为常规全波形反演提供一个较好的初始模型。本发明使用marmousi速度模型进行数值测试并验证本方法可以为常规fwi提供一个较好的初始速度建模，获得高精度反演结果。
1.本发明是基于指数衰减的时频域相位反演方法，主要步骤为：
