CN108345031A - 一种弹性介质主动源和被动源混采地震数据全波形反演方法 - Google Patents
一种弹性介质主动源和被动源混采地震数据全波形反演方法 Download PDFInfo
- Publication number
- CN108345031A CN108345031A CN201810024765.7A CN201810024765A CN108345031A CN 108345031 A CN108345031 A CN 108345031A CN 201810024765 A CN201810024765 A CN 201810024765A CN 108345031 A CN108345031 A CN 108345031A
- Authority
- CN
- China
- Prior art keywords
- source
- component
- seismic
- passive
- virtual
- 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 42
- 239000012530 fluid Substances 0.000 title abstract 2
- 238000012545 processing Methods 0.000 claims abstract description 6
- 238000007781 pre-processing Methods 0.000 claims description 12
- 238000005065 mining Methods 0.000 claims description 8
- 230000005284 excitation Effects 0.000 claims description 5
- 238000012937 correction Methods 0.000 claims description 3
- 238000010606 normalization Methods 0.000 claims description 3
- 230000003068 static effect Effects 0.000 claims description 3
- 238000012217 deletion Methods 0.000 abstract 1
- 230000037430 deletion Effects 0.000 abstract 1
- 238000003384 imaging method Methods 0.000 description 11
- 238000010276 construction Methods 0.000 description 7
- 238000011160 research Methods 0.000 description 7
- 230000008901 benefit Effects 0.000 description 5
- 238000004364 calculation method Methods 0.000 description 5
- 230000000694 effects Effects 0.000 description 5
- 230000008569 process Effects 0.000 description 3
- 238000011084 recovery Methods 0.000 description 3
- 239000004065 semiconductor Substances 0.000 description 3
- NAWXUBYGYWOOIX-SFHVURJKSA-N (2s)-2-[[4-[2-(2,4-diaminoquinazolin-6-yl)ethyl]benzoyl]amino]-4-methylidenepentanedioic acid Chemical compound C1=CC2=NC(N)=NC(N)=C2C=C1CCC1=CC=C(C(=O)N[C@@H](CC(=C)C(O)=O)C(O)=O)C=C1 NAWXUBYGYWOOIX-SFHVURJKSA-N 0.000 description 2
- 238000004458 analytical method Methods 0.000 description 2
- 238000004422 calculation algorithm Methods 0.000 description 2
- 238000002939 conjugate gradient method Methods 0.000 description 2
- 230000007547 defect Effects 0.000 description 2
- 238000000605 extraction Methods 0.000 description 2
- 238000005286 illumination Methods 0.000 description 2
- 238000004519 manufacturing process Methods 0.000 description 2
- 238000012067 mathematical method Methods 0.000 description 2
- 239000002184 metal Substances 0.000 description 2
- 239000003921 oil Substances 0.000 description 2
- 238000003325 tomography Methods 0.000 description 2
- FGUUSXIOTUKUDN-IBGZPJMESA-N C1(=CC=CC=C1)N1C2=C(NC([C@H](C1)NC=1OC(=NN=1)C1=CC=CC=C1)=O)C=CC=C2 Chemical compound C1(=CC=CC=C1)N1C2=C(NC([C@H](C1)NC=1OC(=NN=1)C1=CC=CC=C1)=O)C=CC=C2 FGUUSXIOTUKUDN-IBGZPJMESA-N 0.000 description 1
- 238000010793 Steam injection (oil industry) Methods 0.000 description 1
- 230000005540 biological transmission Effects 0.000 description 1
- 230000015572 biosynthetic process Effects 0.000 description 1
- FFGPTBGBLSHEPO-UHFFFAOYSA-N carbamazepine Chemical compound C1=CC2=CC=CC=C2N(C(=O)N)C2=CC=CC=C21 FFGPTBGBLSHEPO-UHFFFAOYSA-N 0.000 description 1
- 230000000295 complement effect Effects 0.000 description 1
- 239000000295 fuel oil Substances 0.000 description 1
- 230000005484 gravity Effects 0.000 description 1
- 230000006872 improvement Effects 0.000 description 1
- 238000002347 injection Methods 0.000 description 1
- 239000007924 injection Substances 0.000 description 1
- 238000005259 measurement Methods 0.000 description 1
- 238000013508 migration Methods 0.000 description 1
- 230000005012 migration Effects 0.000 description 1
- 238000012544 monitoring process Methods 0.000 description 1
- 230000001902 propagating effect Effects 0.000 description 1
- 238000005070 sampling Methods 0.000 description 1
- 238000004088 simulation Methods 0.000 description 1
- 238000001228 spectrum Methods 0.000 description 1
- 238000010183 spectrum analysis Methods 0.000 description 1
- 238000012360 testing method Methods 0.000 description 1
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/282—Application of seismic models, synthetic seismograms
-
- 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/30—Analysis
- G01V1/307—Analysis for determining seismic attributes, e.g. amplitude, instantaneous phase or frequency, reflection strength or polarity
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01V—GEOPHYSICS; GRAVITATIONAL MEASUREMENTS; DETECTING MASSES OR OBJECTS; TAGS
- G01V2210/00—Details of seismic processing or analysis
- G01V2210/60—Analysis
- G01V2210/62—Physical property of subsurface
- G01V2210/622—Velocity, density or impedance
- G01V2210/6222—Velocity; travel time
Landscapes
- Engineering & Computer Science (AREA)
- Remote Sensing (AREA)
- Physics & Mathematics (AREA)
- Life Sciences & Earth Sciences (AREA)
- Acoustics & Sound (AREA)
- Environmental & Geological Engineering (AREA)
- Geology (AREA)
- General Life Sciences & Earth Sciences (AREA)
- General Physics & Mathematics (AREA)
- Geophysics (AREA)
- Geophysics And Detection Of Objects (AREA)
Abstract
本发明涉及一种弹性介质主动源和被动源混采地震数据全波形反演方法。首先对采集到的被动源弹性多分量地震数据进行地震干涉处理，得到四种不同模式的虚拟炮多分量记录；然后对不同模式的虚拟炮记录两两分组，分别进行弹性波多分量不依赖子波的反演；最后采用交叉选取的方法，分别用高质量Vp和Vs宏观模型作为缺失低频主动源弹性波多分量全波形反演的初始模型，进行反演迭代。被动源弹性多分量地震记录与主动源地震记录串联处理既能克服反演的低频缺失系列问题，也能得到准确的细节和深部构造信息，能更好地刻画纵、横波速度模型的宏观和细节构造信息，并且反演的跳周问题被大大缓解，最终结果既有准确的宏观构造，也有较精确的细节构造信息。
Description
技术领域：
本发明涉及一种地震数据全波形反演方法，尤其是弹性介质主动源和被动源混采地震数据全波形反演方法。
背景技术：
在常规地震采集中，被动源地震信号往往由于其能量弱、不易提取等原因而被忽视。近年来，对被动源地震勘探方法的深入研究发现其具有一系列常规主动源地震勘探不具有的优点。为了更精确地刻画地下介质的细节和深部构造信息，必须充分利用地表接收到的多震源多分量波场信息。本发明主要利用主动源和被动源弹性波多分量地震数据的联合应用提高速度模型构建的精度。
主动源地震采集是指采用人工激发的震源产生地震波，在地表或地下一定深度布设检波器接收记录。被动源地震采集是指在地表或地下一定深度布设检波器，接收由非人工震源激发产生的地震波。近年来，对被动源地震勘探方法的研究使得被动源地震数据的重要性越来越受到重视。被动源地震勘探至少有以下优点：(1)受地表环境影响较小，可在城市、自然保护区等不适合进行主动源勘探的地区展开；(2)携带真实的地下介质信息，用于缺失信息补偿时，比纯数学方法更具物理依据；(3)可进行灵活重构，且重构后虚拟震源位置确定已知；(4) 对深部照明更充分，来自目标层附近的震源信号携带有丰富的构造参数信息；(5) 频带宽，可用于补偿常规主动源缺失的频带信息。尽管被动源地震勘探存在如上所述的优点，但是其能量弱、信噪比低等缺点也限制了其应用效果。针对此，本发明提出了主动源和被动源弹性波地震数据的多分量混合采集和联合应用。利用被动源地震数据辅助和补偿主动源地震数据，从而提高弹性波全波形反演的效果。
对主动源和被动源地震数据的重构以及联合成像已有一些研究。Berkhout 和Verschuur(2009.Integrated imaging with active and passive seismic data.79thSEG Annual Meeting,Expanded Abstracts,1592-1596.)提出了一种主动源和被动源联合成像的理论框架。Vesnaver等(2010.Joint 3D processing of active and passiveseismic data.Geophysical Prospecting, 58:831-844.)模拟了CO2注入区的三维主动源和被动源地震数据混合采集，并进行了层析成像研究，发现被动源可以作为新的照明震源，并且能提高薄气储层的成像效果。Forgues等(2011.Simultaneous active/passiveseismic monitoring of steam assisted heavy oil production.81th SEG AnnualMeeting,Expanded Abstracts,4082-4086.)在一个注蒸汽重油开采区同时采集了主动源和被动源地震数据，并成功监测到蒸汽注入期间储层的时移变换。张盼等提出了主动源和被动源地震数据的插值重构方法，并将联合地震数据应用于叠前偏移成像，验证了主被动源地震数据在成像中的互补作用。Zhang等 (2015.Passive-source multitaper-spectralmethod based low-frequency data reconstruction for active seismicsources.Applied Geophysics,12(4):585-597.)利用被动源的低频信息成功重构了主动源地震数据的缺失低频信息，并发现利用低频重构后的地震数据能提高深部成像效果。Alali等(2016.Merging active and passive seismic reflection data withinterferometry by multidimensional deconvolution.86th SEG Annual Meeting,Expanded Abstracts,5182-5186.)采用多维反褶积方法对主动源和被动源地震数据进行了融合，成功达到了数据拓频的目的。另外，一些针对浅层面波勘探的主被动源地震数据联合应用研究也已展开。Baradello等(2015. Elastic tomographic inversion of thenear surface by active and passive seismic.Near-Surface Asia PacificConference,Expanded Abstracts, 280-283.)利用被动源面波地震数据补偿主动源面波地震数据，成功得到了浅部介质的纵横波速度比。李欣欣(2017.主动源与被动源瑞利波联合成像方法研究[博士论文].西安:长安大学.)研究了主动源与被动源瑞利波联合成像方法，使瑞利波成像的精度和深度均得到了提高。Hayashi(2017.Application of activeand passive surface wave methods to shallow to deep S-wave velocityestimation.87th SEG Annual Meeting,Expanded Abstracts, 5202-5206.)将主动源和被动源联合面波勘探方法应用于实际数据，成功反演了地下横波速度结构。
地震波在地下介质中传播时，不仅有纵波，还有横波和各种转换波。因此，弹性介质比声波介质更接近实际复杂的地下介质情况。弹性波全波形反演理论和方法最初由Tarantola(1986.A strategy for nonlinear elastic inversion of seismicreflection data.Geophysics,51(10):1893-1903.)和Mora(1987. Nonlinear two-dimensional elastic inversion of multioffset seismic data. Geophysics,52(9):1211-1228.)提出和发展。Crase等(1990.Robust elastic nonlinear waveforminversion:Application to real data. Geophysics,55(5):527-538.)最早将弹性波全波形反演方法应用于实际数据。近年来，随着计算能力的提升，弹性波全波形反演的方法研究和实际应用得到了快速的发展。Prieux等(2013.Building starting models for fullwaveform inversion from wide-aperture data by stereotomography.GeophysicalProspecting,61(s1):109-137.)通过立体层析成像构建初始模型，对一套宽角数据进行了频域弹性波全波形反演计算。Vigh等(2014.Elastic full-waveform inversionapplication using multicomponent measurements of seismic datacollection.Geophysics,79(2):R63-R77.) 对墨西哥湾的一套海底电缆四分量地震数据进行了三维弹性波全波形反演计算。孙宏宇等(2015.基于可视性分析与能量补偿的金属矿弹性波全波形反演.地球物理学报，58(12):4605-4616.)针对金属矿速度建模问题，提出了基于可视性分析与能量补偿的弹性波全波形反演方法。王毓伟等(2016.降低弹性波全波形反演强烈非线性的分步反演策略.石油地球物理勘探，51(2):288-294.) 提出了综合利用包络目标函数和多尺度等策略的弹性波全波形反演分步反演策略。Zhang等(2016.Robustsource-independent elastic full-waveform inversion in the timedomain.Geophysics,81(2):R29-R44.)通过对时间域褶积型不依赖子波的目标函数进行加窗，提出了一种稳健的不依赖子波的弹性波全波形反演方法。张盼等(2018.基于各向异性全变分约束的多震源弹性波全波形反演.地球物理学报.doi:10.6038/cjg2018L0122.)将各向异性全变分约束引入到弹性波全波形反演中，用以压制多震源同时反演时的串扰噪声。
在现有技术中没有人提出针对弹性波主动源和被动源多分量地震数据联合全波形反演的方法和策略。主动源和被动源地震数据都是真实的经过地下介质传播的波场记录，比利用数学方法得到的记录更有物理依据。现有技术中也没有涉及到如何挖掘弹性波被动源多分量数据的低频信息，并将其用于全波形反演。
发明内容：
本发明的目的是针对上述现有技术的不足，提供一种弹性介质主动源和被动源混采地震数据全波形反演方法。
本发明的目的是通过以下技术方案实现的：
本发明的思想是：首先，对弹性波被动源多分量地震记录进行干涉重构，得到四种不同模式的虚拟震源；其次，采用相当于Z分量为虚拟震源的两组虚拟炮记录作为观测记录，进行低频弹性波不依赖子波的反演，得到的结果中保留纵波速度模型；再次，采用相当于X分量为虚拟震源的两组虚拟炮记录作为观测记录，进行低频弹性波不依赖子波的反演，得到的结果中保留横波速度模型；最后，以得到的纵、横波速度模型作为初始速度，进行高频主动源弹性波全波形反演，得到的纵、横波速度模型为最终结果。
本发明所述的弹性介质主动源和被动源混采地震数据全波形反演方法是通过MATLAB平台实现的。
一种弹性介质主动源和被动源混采地震数据全波形反演方法，包括以下步骤：
a、安装MATLAB软件平台，要求采用MATLAB R2016a及以上版本。并且已配备并行工具包；
b、对被动源地震数据进行预处理，主要包括能量归一化和去噪处理，得到高质量的多分量被动源地震数据；
c、应用基于互相关的被动源地震干涉法，重构不同模式的虚拟炮记录。重构公式为：
其中，xi和xj代表两个检波器位置，上角标p代表检波器分量，上角标q代表震源分量，对于二维情况，均可以代表x分量或z分量。Gp,q(xj,xi,t)代表在xi位置的q分量震源激发，在xj位置接收的p分量格林函数，t代表时间。*表示褶积运算，S(t)表示虚拟震源子波，ρ代表密度，vp表示纵波速度，dp(xj,t)表示在xj接收的p分量地震记录，dq(xi,t)表示在xi接收的q分量地震记录，<g>表示空间整体平均。公式(1)可以简写为
其中，代表重构的虚拟炮记录，i和j分别表示虚拟震源和虚拟检波器的位置，q和p分别表示震源和检波器的分量；
重构后得到四种类型的虚拟炮记录，分别是和
d、构建初始速度模型，采用估计的最大速度和最小速度进行线性插值，得到由浅到深递增的梯度模型作为初始速度模型；
e、构建弹性波多分量不依赖子波的目标函数如下：
其中，u表示合成地震记录，上角标p代表检波器分量，上角标q代表震源分量， k表示参考道索引；
f、采用虚拟炮记录和作为多分量观测记录，选取低频子波作为模拟子波，进行不依赖子波的弹性波多分量全波形反演，更新迭代后，选取结果中的纵波速度模型Vp；
g、采用虚拟炮记录和作为多分量观测记录，选取低频子波作为模拟子波，进行不依赖子波的弹性波多分量全波形反演，更新迭代后，选取结果中的横波速度模型Vs；
h、对主动源地震记录进行预处理，主要包括静校正、去噪、去多次波，得到高质量的主动源地震数据；
i、根据主动源地震数据估计震源子波；
j、以第f步的Vp模型和第g步的Vs模型作为初始模型，进行主动源弹性波全波形反演的更新迭代。达到迭代停止条件后，输出最终反演结果，即纵波速度模型和横波速度模型。
有益效果：
本发明避免了对弹性波被动源多分量地震数据的震源子波估计和震源定位的问题，并能有效提取和利用被动源多分量地震数据的低频信息，构建出纵、横波速度模型的宏观结构。此宏观结构对于克服低频缺失主动源弹性波全波形反演的跳周问题具有重要的意义。此外，高频主动源弹性波全波形反演串联在低频被动源弹性波全波形反演之后，能够得到更多的细节构造信息。
本发明与现有技术相比具有的实质性特点：
1.本发明利用了被动源弹性波多分量地震数据，并且将被动源多分量地震数据重构得到的所有模式的虚拟炮记录应用到弹性波全波形反演中，充分挖掘了被动源波场携带的地下弹性参数信息。
2.提出了对主动源和被动源多分量弹性波混采地震数据的联合处理，即充分利用采集到的一切与地下介质参数相关的波场信息。不再将被动源地震信号作为噪声而去除，也不将被动源地震信号作为主体，而是将其作为主动源多分量弹性波地震记录的辅助和补充。
3.通过将弹性波地震干涉法和弹性波多分量不依赖子波算法进行联合，有效避免了弹性波被动源多分量地震数据的震源定位和子波估计难题。
4.本发明充分利用了被动源地震数据频带较宽的特点。通过选取低频模拟子波，有效地提取被动源多分量地震波场的低频信息。并将其用于弹性波全波形反演，得到地下纵、横波速度模型的宏观结构。使得整个反演过程更加稳定，结果更加可靠，不再依赖于常规采集中的主动源低频地震信息。
5.采用交叉反演方式。本发明利用Z分量为震源的虚源记录得到纵波速度模型，采用X分量为震源的虚源记录得到横波速度模型。这种交叉组合很好的汇集了反演结果中的高精度结果，能保证下一步反演更好的收敛。
本发明公开的全波形反演方法是目前精度最高的弹性介质主被动源地震数据联合速度建模方法。是在无需震源定位和子波估计的情况下，通过被动源数据反演得到弹性介质的宏观背景构造，之后通过主动源全波形反演可以得到细节构造信息。
附图说明：
图1一种弹性介质主动源和被动源混采地震数据全波形反演方法流程图。
图2真实速度模型。
(a)真实Vp(纵波速度)模型；
(b)真实Vs(横波速度)模型。
图3初始速度模型。
(a)初始Vp模型；
(b)初始Vs模型。
图4缺失低频主动源子波。
(a)子波波形；
(b)子波频谱。
图5被动震源分布图。
图6弹性波被动源干涉重构记录。
(a)重构记录Z分量；
(b)重构记录X分量；
(c)主动源参考记录Z分量；
(d)主动源参考记录X分量。
图7弹性介质主动源和被动源混采地震数据全波形反演结果。
(a)常规主动源Vp反演结果；
(b)常规主动源Vs反演结果；
(c)Z分量虚源被动源Vp反演结果；
(d)Z分量虚源被动源Vs反演结果；
(e)X分量虚源被动源Vp反演结果；
(f)X分量虚源被动源Vs反演结果；
(g)主被动源串联Vp反演结果；
(h)主被动源串联Vs反演结果。
图8反演结果抽道对比。
(a)Vp结果抽道对比；
(b)Vs结果抽道对比。
具体实施方式：
下面结合附图和实例对本发明进一步的详细说明。
本发明所述的一种弹性介质主动源和被动源混采地震数据全波形反演方法，包括以下步骤：
a、在win7或Linux***下安装MATLAB软件平台，要求采用MATLAB R2016a及以上版本，并且已配备并行工具包(Parallel Computing Toolbox)；
b、对被动源地震数据进行预处理，主要包括能量归一化和去噪处理，得到高质量的多分量被动源地震数据；
c、应用基于互相关的被动源地震干涉法，重构不同模式的虚拟炮记录，重构公式为：
其中，xi和xj代表两个检波器位置，上角标p代表检波器分量，上角标q代表震源分量，对于二维情况，均可以代表x分量或z分量。Gp,q(xj,xi,t)代表在xi位置的q分量震源激发，在xj位置接收的p分量格林函数，t代表时间，*表示褶积运算，S(t)表示虚拟震源子波，ρ代表密度，vp表示纵波速度，dp(xj,t)表示在xj接收的p分量地震记录，dq(xi,t)表示在xi接收的q分量地震记录，<g>表示空间整体平均；公式(1)可以简写为
其中，代表重构的虚拟炮记录，i和j分别表示虚拟震源和虚拟检波器的位置，q和p分别表示震源和检波器的分量；
重构后得到四种类型的虚拟炮记录，分别是和
d、构建初始速度模型，采用估计的最大速度和最小速度进行线性插值，得到由浅到深递增的梯度模型作为初始速度模型；
e、构建不依赖子波的弹性波多分量全波形反演的目标函数如下：
其中，u表示合成地震记录，上角标p和q分别表示检波器分量和震源分量，k表示参考道索引；
公式(3)关于模型参数m的导数为
其中，
利用伴随状态法求取梯度时，伴随源为：
其中，表示互相关运算，rj',i为除参考道之外，其他检波器位置处的伴随源，r”j,i为参考道位置的伴随源；利用链式法则，可以得到纵横波速度的计算式如下：
目标函数关于拉梅常数的导数为
其中，ux和uz分别表示正传波场的x和z分量，和分别表示伴随波场的x 和z分量；
f、采用虚拟炮记录和作为多分量观测记录，选取低频子波作为模拟子波，进行不依赖子波的弹性波全波形反演。用低频子波作为震源模拟得到正传波场，用公式(5)和公式(6)作为伴随源，反传得到伴随波场。利用公式(7) —公式(10)计算梯度，用共轭梯度法更新模型。更新迭代后，选取结果中的纵波速度模型Vp；
g、采用虚拟炮记录和作为多分量观测记录，选取低频子波作为模拟子波，进行不依赖子波的弹性波全波形反演。反演过程与第f步类似，只是所有的观测记录换为和且震源作用在x分量上。更新迭代后，选取结果中的横波速度模型Vs；
h、对主动源地震记录进行预处理，主要包括静校正、去噪、去多次波，得到高质量的主动源地震数据；
i、根据主动源地震数据估计震源子波；
j、以第f步的Vp模型和第g步的Vs模型作为初始模型，进行主动源弹性波全波形反演的更新迭代。梯度计算仍然采用公式(7)—公式(10)，但是计算伴随波场的伴随源为x分量的记录残差和z分量的计算残差，即和采用共轭梯度法更新模型。达到迭代停止条件后，输出最终反演结果，即纵波速度模型和横波速度模型。
实施例1：
本发明的整体流程如图1所示。
采用的真实Vp(纵波速度)模型如图2a所示，采用的真实Vs(横波速度) 模型如图2b所示。采用的初始Vp模型如图3a所示，采用的初始Vs模型如图 3b所示。主动源采集时，震源在地表放炮，震源数量为24，左侧第一个震源位于x＝50m位置处，炮间距为80m。检波器布设在地表0m-1910m范围内，总共192个检波器，检波器间距为10m。主动源子波采用缺失低频的雷克子波，如图4所示。图4a为20Hz主频的雷克子波去除9Hz以下低频成分后的子波。图4b为震源子波的频谱分析，可见，9Hz以下的频率成分几乎为零。被动源采集时，200个被动震源随机分布于地下，震源位置如图5所示。震源子波采用频带范围为0-15Hz的随机噪声序列。检波器布设与主动源采集时相同。
首先模拟弹性波被动源地震数据。200个噪声震源随机激发地震信号，采样间隔为1ms，记录总时长为1200s。得到被动源原始记录后，进行基于互相关的被动源地震干涉重构，得到重构的虚拟炮记录。其中，虚拟炮点位置在x＝850m 的X分量激发，Z分量接收的一个虚拟炮记录如图6a所示。虚拟炮点位置在 x＝850m的Z分量激发，X分量接收的一个虚拟炮记录如图6b所示。为对比方便，采用10Hz主频的雷克子波模拟了相同炮点位置和相同激发接收分量的弹性波主动源多分量记录，分别如图6c和图6d所示。由于被动源重构记录的直达波信息时不准确的，因此，主动源参考记录均为去除了直达波之后的结果。对比被动源重构记录和主动源参考记录可知，虽然被动源地震记录中含有较多的干扰，但是仍可见许多的有效反射信息。
下面进行弹性波全波形反演测试。采用图4a所示的缺失低频子波，直接进行常规的主动源弹性波全波形反演。采用如图3所示的模型作为初始模型，迭代 100次的Vp和Vs结果分别如图7a和图7b所示。可见，由于子波主频较高且缺失低频，结果的浅部有明显的局部极小值，且深部反演效果较差。采用不依赖子波的算法对被动源重构数据进行反演。这里采用的模拟子波为8Hz主频的雷克子波。首先，对Z分量激发的被动源多分量虚拟炮记录进行反演，迭代100 次的Vp和Vs结果分别如图7c和图7d所示。可见，Vp结果体现了较好的宏观背景速度构造信息。Vs结果中含有一些干扰信息，但是基本上也体现了Vs 模型的大尺度构造信息。然后，对X分量激发的被动源多分量虚拟炮记录进行反演，迭代100次的Vp和Vs结果分别如图7e和图7f所示。可见，Vp结果构造信息十分不明显，Vs结果体现了大量的大尺度构造信息。对比两组被动源反演结果可知：虚拟子波的影响在反演中被成功的消除；Z分量激发的被动源虚拟炮记录反演得到的Vp结果较好，X分量激发的被动源虚拟炮记录反演得到的Vs结果较好。因此，本发明提出了将两组结构进行交叉组合，采用图7c和图 7f分别作为下一步反演的Vp和Vs的初始模型。以图7c和图7f为初始模型，采用图4a所示的缺失低频子波进行常规主动源弹性波全波形反演。迭代200次的Vp和Vs结果分别如图7g和图7h所示。与图7a和图7b所示的常规反演结果对比可知，浅部的跳周现象被克服，模型的主体构造清晰，而且深部细节构造也得到了较好的恢复。对图7所示的反演结果进行抽道对比，Vp和Vs的抽道结果分别如图8a和图8b所示。可见，常规全波形反演对模型的浅部信息有一定的恢复，但是没能得到深部的速度结构信息。本发明的主动源和被动源串联全波形反演方法不仅能得到较为准确的浅部速度构造，还能得到较好的深部速度结构信息。
Claims (1)
1.一种弹性介质主动源和被动源混采地震数据全波形反演方法，其特征在于，包括以下步骤：
a、安装MATLAB软件平台，要求采用MATLAB R2016a及以上版本，并且已配备并行工具包；
b、对被动源地震数据进行预处理，主要包括能量归一化和去噪处理，得到高质量的多分量被动源地震数据；
c、应用基于互相关的被动源地震干涉法，重构不同模式的虚拟炮记录，重构公式为：
其中，xi和xj代表两个检波器位置，上角标p代表检波器分量，上角标q代表震源分量，对于二维情况，均代表x分量或z分量，Gp,q(xj,xi,t)代表在xi位置的q分量震源激发，在xj位置接收的p分量格林函数，t代表时间，*表示褶积运算，S(t)表示虚拟震源子波，ρ代表密度，vp表示纵波速度，dp(xj,t)表示在xj接收的p分量地震记录，dq(xi,t)表示在xi接收的q分量地震记录，<g>表示空间整体平均，公式(1)简写为：
其中，代表重构的虚拟炮记录，i和j分别表示虚拟震源和虚拟检波器的位置，q和p分别表示震源和检波器的分量；
重构后得到四种类型的虚拟炮记录，分别是和
d、构建初始速度模型，采用估计的最大速度和最小速度进行线性插值，得到由浅到深递增的梯度模型作为初始速度模型；
e、构建不依赖子波的弹性波多分量全波形反演的目标函数如下：
其中，u表示合成地震记录，上角标p代表检波器分量，上角标q代表震源分量，k表示参考道索引；
f、采用虚拟炮记录和作为多分量观测记录，选取低频子波作为模拟子波，进行不依赖子波的弹性波多分量全波形反演，更新迭代后，选取结果中的纵波速度模型Vp；
g、采用虚拟炮记录和作为多分量观测记录，选取低频子波作为模拟子波，进行不依赖子波的弹性波多分量全波形反演，更新迭代后，选取结果中的横波速度模型Vs；
h、对主动源地震记录进行预处理，主要包括静校正、去噪、去多次波，得到高质量的主动源地震数据；
i、根据主动源地震数据估计震源子波；
j、以第f步的Vp模型和第g步的Vs模型作为初始模型，进行主动源弹性波全波形反演的更新迭代，达到迭代停止条件后，输出最终反演结果，即纵波速度模型和横波速度模型。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201810024765.7A CN108345031B (zh) | 2018-01-11 | 2018-01-11 | 弹性介质主动源和被动源混采地震数据全波形反演方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201810024765.7A CN108345031B (zh) | 2018-01-11 | 2018-01-11 | 弹性介质主动源和被动源混采地震数据全波形反演方法 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN108345031A true CN108345031A (zh) | 2018-07-31 |
| CN108345031B CN108345031B (zh) | 2020-01-17 |
Family
ID=62960452
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201810024765.7A Active CN108345031B (zh) | 2018-01-11 | 2018-01-11 | 弹性介质主动源和被动源混采地震数据全波形反演方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN108345031B (zh) |
Cited By (18)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN109407151A (zh) * | 2018-12-18 | 2019-03-01 | 吉林大学 | 基于波场局部相关时移的时间域全波形反演方法 |
| CN109407152A (zh) * | 2018-12-18 | 2019-03-01 | 吉林大学 | 基于零均值归一化互相关目标函数的时间域全波形反演方法 |
| CN109507726A (zh) * | 2017-09-15 | 2019-03-22 | 中国石油化工股份有限公司 | 时间域弹性波多参数全波形的反演方法及*** |
| CN109738952A (zh) * | 2019-01-24 | 2019-05-10 | 吉林大学 | 基于全波形反演驱动的被动源直接偏移成像方法 |
| CN111103621A (zh) * | 2019-12-09 | 2020-05-05 | 东华理工大学 | 一种主动源共成像点叠加多道面波分析方法 |
| CN111505714A (zh) * | 2020-04-16 | 2020-08-07 | 吉林大学 | 基于岩石物理约束的弹性波直接包络反演方法 |
| CN111665555A (zh) * | 2019-03-07 | 2020-09-15 | 中普宝信(北京)科技有限公司 | 拉梅参数反演方法 |
| CN111781639A (zh) * | 2020-06-04 | 2020-10-16 | 同济大学 | 针对obs多分量数据的炮检互易弹性波全波形反演方法 |
| CN112987099A (zh) * | 2021-04-19 | 2021-06-18 | 吉林大学 | 基于多震源卷积神经网络的低频地震数据重构方法 |
| CN113391351A (zh) * | 2021-06-10 | 2021-09-14 | 中国地质科学院 | 一种基于被动源地震波场分析提取矿集区结构的方法 |
| CN113614578A (zh) * | 2019-03-07 | 2021-11-05 | 沙特***石油公司 | 使用双向虚拟源重建的地震勘测 |
| CN113866821A (zh) * | 2021-09-26 | 2021-12-31 | 吉林大学 | 一种基于照明方向约束的被动源干涉偏移成像方法和*** |
| CN114002740A (zh) * | 2020-07-28 | 2022-02-01 | 北京大学 | 一种低信噪比地震初至信号自动增强方法及*** |
| CN114779335A (zh) * | 2022-03-31 | 2022-07-22 | 吉林大学 | 基于各向异性全变分约束的弹性波直接包络反演方法 |
| CN115184986A (zh) * | 2022-06-28 | 2022-10-14 | 吉林大学 | 不依赖震源的全局包络互相关全波形反演方法 |
| CN116338779A (zh) * | 2023-05-11 | 2023-06-27 | 中国铁路设计集团有限公司 | 用于密集线性台阵的交火式主被动源联合探测方法 |
| CN116660979A (zh) * | 2023-06-02 | 2023-08-29 | 哈尔滨工程大学 | 一种基于Kaiser时窗积分的OBN资料全波形反演方法 |
| WO2024016574A1 (zh) * | 2022-07-19 | 2024-01-25 | 中国石油天然气集团有限公司 | 地下介质速度的确定方法、装置、设备 |
Citations (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20100219980A1 (en) * | 1998-08-28 | 2010-09-02 | Marathon Oil Company | Method and system for performing operations and for improving production in wells |
| CN103105624A (zh) * | 2011-11-11 | 2013-05-15 | 中国石油集团川庆钻探工程有限公司地球物理勘探公司 | 基于数据库技术的纵横波时差定位方法 |
| CN106594170A (zh) * | 2016-11-08 | 2017-04-26 | 浙江大学 | 一种主被动混合减震文物减震浮放平台的控制方法 |
-
2018
- 2018-01-11 CN CN201810024765.7A patent/CN108345031B/zh active Active
Patent Citations (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20100219980A1 (en) * | 1998-08-28 | 2010-09-02 | Marathon Oil Company | Method and system for performing operations and for improving production in wells |
| CN103105624A (zh) * | 2011-11-11 | 2013-05-15 | 中国石油集团川庆钻探工程有限公司地球物理勘探公司 | 基于数据库技术的纵横波时差定位方法 |
| CN106594170A (zh) * | 2016-11-08 | 2017-04-26 | 浙江大学 | 一种主被动混合减震文物减震浮放平台的控制方法 |
Non-Patent Citations (1)
| Title |
|---|
| PAN ZHANG EL: ""Passive seismic full waveform inversion based on seismic interferometry and source independent method"", 《2017 SEG INTERNATIONAL EXPOSITION AND 87TH ANNUAL MEETING》 * |
Cited By (25)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN109507726A (zh) * | 2017-09-15 | 2019-03-22 | 中国石油化工股份有限公司 | 时间域弹性波多参数全波形的反演方法及*** |
| CN109407152A (zh) * | 2018-12-18 | 2019-03-01 | 吉林大学 | 基于零均值归一化互相关目标函数的时间域全波形反演方法 |
| CN109407151A (zh) * | 2018-12-18 | 2019-03-01 | 吉林大学 | 基于波场局部相关时移的时间域全波形反演方法 |
| CN109738952A (zh) * | 2019-01-24 | 2019-05-10 | 吉林大学 | 基于全波形反演驱动的被动源直接偏移成像方法 |
| CN109738952B (zh) * | 2019-01-24 | 2020-05-08 | 吉林大学 | 基于全波形反演驱动的被动源直接偏移成像方法 |
| CN113614578A (zh) * | 2019-03-07 | 2021-11-05 | 沙特***石油公司 | 使用双向虚拟源重建的地震勘测 |
| CN111665555A (zh) * | 2019-03-07 | 2020-09-15 | 中普宝信(北京)科技有限公司 | 拉梅参数反演方法 |
| CN111103621A (zh) * | 2019-12-09 | 2020-05-05 | 东华理工大学 | 一种主动源共成像点叠加多道面波分析方法 |
| CN111505714A (zh) * | 2020-04-16 | 2020-08-07 | 吉林大学 | 基于岩石物理约束的弹性波直接包络反演方法 |
| CN111505714B (zh) * | 2020-04-16 | 2021-05-25 | 吉林大学 | 基于岩石物理约束的弹性波直接包络反演方法 |
| CN111781639A (zh) * | 2020-06-04 | 2020-10-16 | 同济大学 | 针对obs多分量数据的炮检互易弹性波全波形反演方法 |
| CN114002740A (zh) * | 2020-07-28 | 2022-02-01 | 北京大学 | 一种低信噪比地震初至信号自动增强方法及*** |
| CN112987099A (zh) * | 2021-04-19 | 2021-06-18 | 吉林大学 | 基于多震源卷积神经网络的低频地震数据重构方法 |
| CN113391351B (zh) * | 2021-06-10 | 2024-02-13 | 中国地质科学院 | 一种基于被动源地震波场分析提取矿集区结构的方法 |
| CN113391351A (zh) * | 2021-06-10 | 2021-09-14 | 中国地质科学院 | 一种基于被动源地震波场分析提取矿集区结构的方法 |
| CN113866821A (zh) * | 2021-09-26 | 2021-12-31 | 吉林大学 | 一种基于照明方向约束的被动源干涉偏移成像方法和*** |
| CN113866821B (zh) * | 2021-09-26 | 2022-08-02 | 吉林大学 | 一种基于照明方向约束的被动源干涉偏移成像方法和*** |
| CN114779335A (zh) * | 2022-03-31 | 2022-07-22 | 吉林大学 | 基于各向异性全变分约束的弹性波直接包络反演方法 |
| CN115184986A (zh) * | 2022-06-28 | 2022-10-14 | 吉林大学 | 不依赖震源的全局包络互相关全波形反演方法 |
| CN115184986B (zh) * | 2022-06-28 | 2023-09-29 | 吉林大学 | 不依赖震源的全局包络互相关全波形反演方法 |
| WO2024016574A1 (zh) * | 2022-07-19 | 2024-01-25 | 中国石油天然气集团有限公司 | 地下介质速度的确定方法、装置、设备 |
| CN116338779B (zh) * | 2023-05-11 | 2023-08-08 | 中国铁路设计集团有限公司 | 用于密集线性台阵的交火式主被动源联合探测方法 |
| CN116338779A (zh) * | 2023-05-11 | 2023-06-27 | 中国铁路设计集团有限公司 | 用于密集线性台阵的交火式主被动源联合探测方法 |
| CN116660979A (zh) * | 2023-06-02 | 2023-08-29 | 哈尔滨工程大学 | 一种基于Kaiser时窗积分的OBN资料全波形反演方法 |
| CN116660979B (zh) * | 2023-06-02 | 2024-01-30 | 哈尔滨工程大学 | 一种基于Kaiser时窗积分的OBN资料全波形反演方法 |
Also Published As
| Publication number | Publication date |
|---|---|
| CN108345031B (zh) | 2020-01-17 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN108345031B (zh) | 弹性介质主动源和被动源混采地震数据全波形反演方法 | |
| Pan et al. | High-resolution characterization of near-surface structures by surface-wave inversions: from dispersion curve to full waveform | |
| CN103713315B (zh) | 一种地震各向异性参数全波形反演方法及装置 | |
| CN109669212B (zh) | 地震数据处理方法、地层品质因子估算方法与装置 | |
| CN104570082B (zh) | 一种基于格林函数表征的全波形反演梯度算子的提取方法 | |
| CN107765308B (zh) | 基于褶积思想与精确震源的重构低频数据频域全波形反演方法 | |
| CN101334483A (zh) | 一种在地震数据处理中衰减瑞雷波散射噪声的方法 | |
| CN108845351A (zh) | 一种vsp地震资料转换波全波形反演方法 | |
| CN111045077B (zh) | 一种陆地地震数据的全波形反演方法 | |
| CN104237945A (zh) | 一种地震资料自适应高分辨处理方法 | |
| Fei et al. | Full waveform inversion without low frequencies: A synthetic study | |
| Yin et al. | Improving horizontal resolution of high-frequency surface-wave methods using travel-time tomography | |
| CN102385066B (zh) | 一种叠前地震定量成像方法 | |
| CN110737018B (zh) | Vsp地震资料各向异性建模方法 | |
| CN110850469A (zh) | 一种基于克希霍夫积分解的地震槽波深度偏移的成像方法 | |
| CN110018516A (zh) | 一种解耦波场微地震逆时干涉定位方法 | |
| CN111505714B (zh) | 基于岩石物理约束的弹性波直接包络反演方法 | |
| Ernst | Modal elastic inversion | |
| CN111175822B (zh) | 改进直接包络反演与扰动分解的强散射介质反演方法 | |
| Trinh et al. | Data-oriented strategy and Vp/Vs model-constraint for simultaneous Vp and Vs reconstruction in 3D visco-elastic FWI: Application to the SEAM II Foothill dataset | |
| Hanafy et al. | Early arrival waveform inversion of shallow seismic land data | |
| Kurzmann et al. | Real data applications of seismic full waveform inversion | |
| CN116699686B (zh) | 基于余弦角差恒等式的浅地表瑞雷波瞬时能量反演方法 | |
| CN109471167A (zh) | 一种针对多震源缺失数据的波场重构反演方法 | |
| Li et al. | An explicit stable Q-compensated reverse time migration scheme for complex heterogeneous attenuation media |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant |