CN114509805B - 矢量凸集投影多分量三维地震数据重建方法及装置 - Google Patents

矢量凸集投影多分量三维地震数据重建方法及装置 Download PDF

Info

Publication number
CN114509805B
CN114509805B CN202110529379.5A CN202110529379A CN114509805B CN 114509805 B CN114509805 B CN 114509805B CN 202110529379 A CN202110529379 A CN 202110529379A CN 114509805 B CN114509805 B CN 114509805B
Authority
CN
China
Prior art keywords
component
complex
reconstruction
quaternion
frequency
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.)
Active
Application number
CN202110529379.5A
Other languages
English (en)
Other versions
CN114509805A (zh
Inventor
高建军
李帆
王赟
芦俊
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
China University of Geosciences Beijing
Original Assignee
China University of Geosciences Beijing
Priority date (The priority date 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 date listed.)
Filing date
Publication date
Application filed by China University of Geosciences Beijing filed Critical China University of Geosciences Beijing
Priority to CN202110529379.5A priority Critical patent/CN114509805B/zh
Publication of CN114509805A publication Critical patent/CN114509805A/zh
Application granted granted Critical
Publication of CN114509805B publication Critical patent/CN114509805B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

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
    • 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
    • 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
    • G01V1/00—Seismology; Seismic or acoustic prospecting or detecting
    • G01V1/28—Processing seismic data, e.g. for interpretation or for event detection
    • G01V1/36—Effecting static or dynamic corrections on records, e.g. correcting spread; Correlating seismic signals; Eliminating effects of unwanted energy
    • 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/624—Reservoir parameters

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

本发明提供了一种矢量凸集投影多分量三维地震数据重建方法及装置，其中，该方法包括：将三或四分量三维地震数据的各分量数据变换至频率域；在频率域提取各分量正频率部分的单频率切片数据；将提取的正频率部分的单频率切片数据矢量联合表示为复四元数；基于复四元数傅里叶变换，利用矢量凸集投影迭代方法对单频率切片对应的复四元数矢量联合重建，得到正频率部分单频率切片的重建结果；对正频率部分单频率切片的重建结果中的各分量数据取共轭，得到负频率部分单频率切片的重建结果；根据正频率部分和相应负频率部分的重建结果，得到频率域所有频率切片的重建结果，并变换至时间域。通过上述方案能实现三或四分量数据联合重建而且能提高重建效率。

Description

矢量凸集投影多分量三维地震数据重建方法及装置
技术领域
本发明涉及地质勘探技术领域，尤其涉及一种矢量凸集投影多分量三维地震数据重建方法及装置。
背景技术
在陆上地震数据的采集阶段，受障碍物、采集禁区及采集成本因素的影响，采集到的地震数据通常存在稀疏或者不规则缺失问题。在海上地震数据的采集阶段，采用 海洋拖缆(Tower Streamer，简称TC)采集地震数据时会受到羽状漂移影响，采用海 底电缆(Ocean Bottom Cable，简称OBC)采集数据时会受海流和潮汐的影响，使用 海底节点地震仪(Ocean Bottom Node，简称OBN)采集数据时会受到海底崎岖复杂 地形和仪器成本因素的影响，这些影响使得采集到的地震数据也存在稀疏采样或者不 规则采样问题。另外，在室内地震资料处理阶段，由于剔除废炮和废道，以及数据处 理域转换也会引起地震道的不规则缺失。不规则缺失地震数据会对常规多道处理技术 的正确运行产生不良影响，降低地震资料品质，因此，需要对不规则地震数据进行规 则化重建处理。
此外，随着油气勘探难度的日益加大，单分量纵波勘探已不能满足现今勘探需求，多波多分量勘探因加入了横波信息，能有效降低储层预测的多解性，提高储层流体识 别和预测精度，因而，受到了石油工业界的普遍重视，多分量地震数据采集技术也逐 渐进入常规工业生产流程。在多分量地震数据采集时，陆地上使用三分量检波器进行 记录，即水平方向的X分量、Y分量以及垂向的Z分量，而海上使用四分量检波器 进行记录，即X分量、Y分量、Z分量和压力分量。
像单分量采集一样，多分量数据在采集过程中也会面临不规则缺失道问题，也同样需要进行规则化插值重建。然而，目前已有的重建方法大都是针对单分量数据的规 则化重建问题研发的，针对多分量数据重建的方法非常稀缺。当前对多分量缺失数据 的重建普遍采用单分量重建的方法，即逐个分量用单分量的方法各自进行重建。这种 做法破坏了各个分量之间的相互关系以及地震波波场结构的完整性，因此迫切需要开 发能适应于多分量数据联合重建的新方法。
在多分量地震数据重建方面，有人提出矢量凸集投影重建算法实现了对X分量 和Y分量的二分量数据联合重建。该方法在频率域将X分量的实部和虚部、Y分量 的实部和虚部联合在一起组建实数四元数，然后利用实数傅里叶变换和凸集投影方法 实现二分量不规则缺失道数据的规则化重建，但是受实数四元数理论的限制，该方法 无法实现对X、Y和Z三分量数据的联合重建，并且无法实现对X和Y分量具有不 同缺失模式数据的重建。有人对基于K-SVD(Singular Value Decomposition，SVD， 奇异值分解)稀疏字典学习算法的单分量地震数据重建方法进行改进，提出了基于 K-SVD算法的三分量地震数据重建方法。该方法将三个分量X、Y和Z数据串联成 一个长的列向量来实现三分量数据的联合，并定义了新的能适用于三分量数据运算的 内积表达式来实现基于K-SVD算法的三分量二维地震数据重建。由于每一次对字典 矩阵更新都需要做K次SVD分解，导致计算量大，使得该方法很难扩展到三维以及 更高维三分量数据的重建。有人提出基于四元数矩阵降秩的三分量三维地震数据重建 方法，该方法借助实四元数傅里叶变换将用于单分量重建的奇异谱分析算法(Singular Spectrum Analysis，SSA，奇异谱分析)改造为实数四元数奇异谱分析算法(Quaternionic Singular Spectrum Analysis，QSSA，四元数奇异谱分析)，利用实数四 元数SVD实现矩阵降秩，进而实现三分量三维不规则缺失地震数据重建。受奇异谱 分析算法适用条件的限制，该方法只适应于线性同相轴数据重建，对于弯曲同相轴数 据需分窗口进行重建。有人提出基于四元数稀疏紧致集方法的三分量二维地震数据重 建方法。该紧致集方法通过四元数奇异值分解(Quaternionic Singular Value Decomposition，QSVD，四元数奇异值分解)和稀疏正则化来实现重建。但由于QSVD 计算量大，成本高，使得该方法目前还停留在二维数据重建阶段。此外，上述方法都 存在一个共同的缺陷，即对多分量数据中每个分量具有不同缺失模式的数据无法实现 重建，只能对各分量保持相同缺失模式的多分量数据进行重建。
目前随着陆上地震勘探三分量数据采集技术和海上四分量数据采集技术的普遍推广，三分量和四分量数据的规则化重建问题已显得十分必要和迫切。但是，目前矢 量凸集投影重建算法存在只能重建二分量数据，无法重建三分量数据，以及只能对各 分量保持相同缺失模式数据进行重建，无法对各分量具有不同缺失模式数据进行重建 的缺陷。
发明内容
鉴于此，本发明提供了一种矢量凸集投影多分量三维地震数据重建方法及装置，以解决现有技术存在的一个或多个缺陷。
为了达到上述目的，本发明采用以下方案实现：
根据本发明实施例的一个方面，提供了一种矢量凸集投影多分量三维地震数据重建方法，包括：
对时间域的多分量三维地震数据的各分量数据做关于时间变量的普通傅里叶变换，得到频率域的多分量三维地震数据的各分量数据；其中，所述多分量为三分量或 四分量；
从频率域的多分量三维地震数据的各分量数据中分别提取正频率部分的频率切片数据；
利用复四元数框架将从各分量数据提取的正频率部分的单频率切片数据进行矢量联合表示，得到相应正频率部分的单频率切片对应的复四元数；
基于复四元数傅里叶变换，利用矢量凸集投影迭代方法对单频率切片对应的复四元数进行矢量联合重建，得到相应正频率部分的单频率切片的多分量三维地震数据的 重建结果；
对正频率部分的单频率切片的多分量三维地震数据的重建结果中的各分量数据取共轭，得到相应负频率部分的单频率切片的多分量三维地震数据重建结果；
根据正频率部分的单频率切片的多分量三维地震数据的重建结果和相应负频率部分的单频率切片的多分量三维地震数据的重建结果，得到整个频率域的多分量三维 地震数据的重建结果；
对整个频率域的多分量三维地震数据的重建结果中的各分量数据分别做普通傅里叶反变换，得到时间域的多分量三维地震数据的重建结果。
根据本发明实施例的另一个方面，还提供了一种计算机设备，包括存储器、处理器及存储在存储器上并可在处理器上运行的计算机程序，所述处理器执行所述程序时 实现上述任一实施例所述方法的步骤。
根据本发明实施例的另一个方面，还提供了一种计算机可读存储介质，其上存储有计算机程序，该程序被处理器执行时实现上述任一实施例所述方法的步骤。
本发明实施例的矢量凸集投影多分量三维地震数据重建方法、计算机设备及计算机可读存储介质，采用复数四元数傅里叶变换方法，在频率域能将三分量数据或四分 量数据联合成复数四元数形式的矢量数据，然后，借助复数四元数傅里叶变换和凸集 投影算法实现三分量地震数据的联合重建。而且采用复数四元数傅里叶变换能够保持 地震数据在频率域正负频率切片的共轭对称性，在重建时只需要对正频率切片进行迭 代重建，负频率部分通过对正频率部分的重建结果取共轭对称即可得到，减少了数据 重建时的计算量，提高了重建计算效率。
附图说明
为了更清楚地说明本发明实施例或现有技术中的技术方案，下面将对实施例或现有技术描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅 是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提 下，还可以根据这些附图获得其他的附图。在附图中：
图1是本发明一实施例的矢量凸集投影多分量三维地震数据重建方法的流程示意图；
图2是本发明一具体实施例中的X分量的原始完整数据和随机缺失数据的示意图；
图3是本发明一具体实施例中的Y分量的原始完整数据和随机缺失数据的示意图；
图4是本发明一具体实施例中的Z分量的原始完整数据和随机缺失数据的示意图；
图5是利用单分量凸集投影方法和本发明实施例的矢量凸集投影重建X分量的结果示意图；
图6是利用单分量凸集投影方法和本发明实施例的矢量凸集投影重建X分量的重建误差示意图；
图7是利用单分量凸集投影方法和本发明实施例的矢量凸集投影重建Y分量的结果示意图；
图8是利用单分量凸集投影方法和本发明实施例的矢量凸集投影重建Y分量的重建误差示意图；
图9是利用单分量凸集投影方法和本发明实施例的矢量凸集投影重建Z分量的结果示意图；
图10是利用单分量凸集投影方法和本发明实施例的矢量凸集投影重建Z分量的重建误差示意图；
图11是利用单分量重建方法得到的X、Y和Z三个分量重建结果与原始真实的X、 Y和Z三个分量地震波波形对比的矢端图；
图12是利用本发明一实施例的重建方法得到的X、Y和Z三个分量重建结果与原 始真实的X、Y和Z三个分量地震波波形对比矢端图；
图13是不同缺失度下单分量方法、二分量方法及本发明一实施例的三分量方法的X分量重建质量因子曲线；
图14是不同缺失度下单分量方法、二分量方法及本发明一实施例的三分量方法的Y分量重建质量因子曲线；
图15是不同缺失度下单分量方法、二分量方法及本发明一实施例的三分量方法的Z分量重建质量因子曲线；
图16是缺失炮对应的共炮点道集利用单分量凸集投影方法重建的结果；
图17是缺失炮对应的共炮点道集利用本发明一实施例的多分量矢量凸集投影方法重建的结果。
具体实施方式
为使本发明实施例的目的、技术方案和优点更加清楚明白，下面结合附图对本发明实施例做进一步详细说明。在此，本发明的示意性实施例及其说明用于解释本发明， 但并不作为对本发明的限定。
需要预先说明的是，下述实施例或示例的描述或其中所提及的特征可以用相同或类似的方式，与其他实施例或示例中的特征组合，或替换其他实施例或示例中的特征， 以形成可能的实施方式。另外，本文所使用的术语“包括/包含”是指特征、要素、 步骤或组件的存在，但并不排除还存在一个或多个其他特征、要素、步骤或组件。
基于傅里叶变换的凸集投影重建一般包括以下步骤：
(1)将t-x-y域三维观测数据Dobs(t,x,y)做关于时间t的普通傅里叶变换到频 率域，得到f-x-y域三维数据
Figure BDA0003067454800000061
简记为
Figure BDA0003067454800000062
(2)设定重建的频带宽度，在该频带内对每一个频率切片fn对应的数据
Figure BDA00030674548000000622
使用凸集投影迭代公式
Figure BDA0003067454800000063
l＝1,2,...,Niter进行迭代，实现缺失道的恢 复重建；
(3)对每一个频率切片fn重复步骤(2)进行重建；
(4)对频率f做傅里叶反变换到时间域完成重建。
基于凸集投影方法进行二分量地震数据重建，可以利用实数四元数工具将频率域两个分量数据的实数实部和实数虚部进行矢量联合组建为实四元数，然后利用凸集投 影算法对实数四元数进行重建，进而实现对两分量数据的矢量联合重建。
以X分量数据DX(t,x1,x2)和Y分量数据DY(t,x1,x2)的重建为例，基于凸集投影 理论进行二分量地震数据重建的方法可包括以下步骤：
(1)首先将DX(t,x1,x2)和DY(t,x1,x2)做关于时间t的傅里叶变换到频率域得
Figure BDA0003067454800000064
和
Figure BDA0003067454800000065
简记为
Figure BDA0003067454800000066
和
Figure BDA0003067454800000067
(2)因在频率域，X分量数据
Figure BDA0003067454800000068
为复数，设
Figure BDA0003067454800000069
的实部为
Figure BDA00030674548000000610
虚部为
Figure BDA00030674548000000611
Figure BDA00030674548000000612
Y分量数据
Figure BDA00030674548000000613
也为复数，设
Figure BDA00030674548000000614
的实部为
Figure BDA00030674548000000615
虚部为
Figure BDA00030674548000000624
Figure BDA00030674548000000616
构建实数四元数q，
Figure BDA00030674548000000618
实现对X和Y两分量数据的矢量联合， 其中，
Figure BDA00030674548000000619
和
Figure BDA00030674548000000623
均为实数矩阵；
(3)利用实数四元数傅里叶变换和凸集投影算法实现二分量数据重建，迭代重建的表达式为：
Figure BDA00030674548000000620
l＝1,2,...,Niter，其中，ql表示第 l次迭代的重建结果，α表示加权因子，I表示全1矩阵，S表示采样算子，
Figure BDA00030674548000000621
表示 实数四元数傅里叶反变换，Tl表示第l次重建的阈值算子，Fq表示实数四元数傅里叶 变换，ql-1表示第l-1次迭代的重建结果，Niter表示最大迭代次数；
(4)迭代停止时，得到
Figure BDA0003067454800000071
从重建结果ql中提取出
Figure BDA0003067454800000072
和
Figure BDA0003067454800000073
构建X分量重建数据
Figure BDA0003067454800000074
Figure BDA0003067454800000077
和Y分量数据
Figure BDA0003067454800000075
Figure BDA0003067454800000076
(5)对每个频率切片fn重复步骤(2)-(4)完成频率域重建；
(6)对频率f做傅里叶反变换到时间域得到时间域重建数据两个分量的重建结果DX(t,x1,x2)和DY(t,x1,x2)。
但是，该方法存在如下缺陷：(1)采用实数四元数工具，利用二分量数据的实部 和虚部构建实数四元数，实现二分量数据的联合表示，无法对三分量或者更多分量数 据的实部和虚部进行四元数表示，也即无法实现对三分量或者更多分量数据的联合重 建；(2)无法实现不同分量数据在不同缺失模式下的重建，当两个分量各自缺失度不 同或者在相同缺失度下缺失道位置不同时，该方法无法实现重建。
另外，对于三分量数据X、Y和Z，在使用凸集投影方法进行重建时，需要对每 个分量X、Y和Z分别单独进行重建。该方法存在如下缺陷：(1)属于单分量重建， 每次只能重建一个分量，缺少多分量数据的矢量联合机制，无法同时实现三个分量数 据的矢量联合重建。对三个分量数据进行重建时只能逐个分量逐一进行单独重建。
针对目前矢量凸集投影重建技术存在的只能重建二分量数据无法同时重建三分量和四分量数据且二分量数据需要保持相同缺失模式的问题，以及采用实四元数傅里 叶变换重建导致的正负频率切片丧失共轭对称性问题，本发明提出一种矢量凸集投影 多分量三维地震数据重建方法。本发明实施例的方法基于复数四元数傅里叶变换进行 矢量凸集投影重建，可以有效解决目前矢量凸集投影重建存在的不足。
图1是本发明一实施例的矢量凸集投影多分量三维地震数据重建方法的流程示意图，如图1所示，该实施例的矢量凸集投影多分量三维地震数据重建方法可包括以 下步骤S110～步骤S170。
下面将对步骤S110至步骤S170的具体实施方式进行详细说明。
步骤S110：对时间域的多分量三维地震数据的各分量数据做关于时间变量的普通傅里叶变换，得到频率域的多分量三维地震数据的各分量数据；其中，所述多分量 为三分量或四分量。
该步骤S110中，可以利用布置在地面不同位置的检波器分别采集一定时间内的三维地震数据。三分量三维地震数据可以是三分量勘探中三分量检波器所采集的共炮 点域X、Y及Z分量数据。四分量三维地震数据可以是四分量勘探中三分量检波器所 采集的共炮点域X、Y、Z分量数据和单分量检波器采集的P分量数据。
步骤S120：从频率域的多分量三维地震数据的各分量数据中分别提取正频率部分的频率切片数据。
该步骤S120中，可以按一定频率间隔从频率连续的频率域的多分量三维地震数据的各分量数据中提取得到单频率切片的数据。正频率部分的频率切片数据和负频率 部分的频率切片数据是互为共轭关系(复数共轭)。
步骤S130：利用复四元数框架将从各分量数据中提取的正频率部分的单频率切片数据进行矢量联合表示，得到相应正频率部分单频率切片对应的复四元数。
该步骤S130中，利用复四元数框架构建的复四元数(复数四元数)也称超复数， 包含一个实部与三个虚部，其形式可以表示为q＝a+bi+cj+dk，其中，当a,b,c,d为实数，则称为实四元数，当a,b,c,d为普通复数则称为复四元数，即 q＝areal+aimagI+(breal+bimagI)i+(creal+cimagI)j+(dreal+dimagI)k。当a＝0时， 四元数q也被称为纯四元数。i,j,k为虚数单位，i,j,k的关系为 i2＝j2＝k2＝ijk＝-1以及j＝k＝-ji,jk＝i＝-kj,ki＝j＝-ik。
从几何意义上讲，虚数单位i,j,k对应于三维空间中直角坐标系下的X、Y和Z 轴，但是，纯四元数不限于利用正交基i,j,k来表示，还可以利用任意方向三个正交 基μ,ν,ξ来表示，这类似于将三维空间中X、Y和Z轴在保持正交结构不变的同时又 以一定角度旋转后形成的坐标轴。式 q＝areal+aimagI+(breal+bimagI)i+(creal+cimagI)j+(dreal+dimagI)k在广义基(1,μ,ν,ξ)下可以表示为q′＝a′+b′μ+c′ν+d′ξ，其中，μ，ν和ξ是正交的单位 纯四元数，满足μ⊥ν⊥ξ，ξ＝μν，且μνξ＝-1，μ，ν和ξ与i,j,k之间的关系可以 用一个3×3的复正交矩阵表示：
Figure BDA0003067454800000081
其中，μl,νl,ξl,l＝1,2,3为普通复数。
在以上复正交变换矩阵式的作用下，式q＝a+bi+cj+dk和式 q′＝a′+b′μ+c′ν+d′ξ满足如下对应关系：
a′＝a，b′＝bμ1+cμ2+dμ3，c′＝bν1+cν2+dν3，d′＝bξ1+cξ2+dξ3。
在一些实施例中，可以根据不同分量的单频率切片数据构建得到复四元数矩阵。例如，图1所示的方法，还可包括步骤：S170，将得到的各分量单频率切片对应的复 四元数按照检波器在地面不同位置处的横坐标和纵坐标两个空间变量维度表示成复 四元数矩阵，以用于矢量联合重建。该步骤S170中，复四元数矩阵的行和列中一个 可以对应检波器在地面的横坐标x1，另一个可以对应检波器在地面的纵坐标x2。
在一些实施例中，对于三分量和四分量数据均可以构建复四元数。
示例性地，在多分量为三分量的情况下，上述步骤S130，即，利用复四元数框 架将从各分量数据提取的正频率部分的单频率切片数据进行矢量联合表示，得到相应 正频率部分的单频率切片对应的复四元数，具体可包括步骤：S1311，在多分量为三 分量的情况下，将从各分量数据提取的正频率部分的单频率切片数据中的每个分量作 为复四元数框架中的一个虚部，得到相应正频率部分的单频率切片对应的复四元数。 例如，构建的复四元数表示为q＝a+bi+cj+dk，其中，b,c,d分别对应分别为各 分量数据。在其他实施例中三个分量也可以与实部和两个虚部对应。
示例性地，在多分量为四分量的情况下，上述步骤S130，即，利用复四元数框 架将从各分量数据提取的正频率部分的单频率切片数据进行矢量联合表示，得到相应 正频率部分的单频率切片对应的复四元数，具体可包括步骤：S1321，在多分量为四 分量的情况下，将从各分量数据提取的正频率部分的单频率切片数据中任选其中一个 分量作为复四元数框架中的实部，其余三个分量中的每个分量作为复四元数框架中的 一个虚部，得到相应正频率部分的单频率切片对应的复四元数。例如，构建的复四元 数表示为q＝a+bi+cj+dk，其中，a,b,c,d分别对应分别为各分量数据。
步骤S140：基于复四元数傅里叶变换，利用矢量凸集投影迭代方法对单频率切 片对应的复四元数进行矢量联合重建，得到相应正频率部分的单频率切片的多分量三 维地震数据的重建结果。
该步骤S140中，矢量凸集投影迭代方法可以基于现有重建表达式做所需改进得到。重建表达式中使用到复四元数傅里叶变换，还可以用到复四元数傅里叶反变换。 在前述步骤构建得到复四元数矩阵的情况下，重建表达式中相应量可以表示为矩阵形 式。
在进行矢量联合重建过程中，可以根据不同分量的缺失道模式即采样算子的情况，可以采用不同的重建表达式进行重建。
在一些实施例中，在多分量三维地震数据的各分量数据的采样算子相同的情况下，上述步骤S140，即，基于复四元数傅里叶变换，利用矢量凸集投影迭代方法对 单频率切片对应的复四元数进行矢量联合重建，得到相应正频率部分的单频率切片的 多分量三维地震数据的重建结果，具体可包括步骤：S141，在多分量三维地震数据的 各分量数据的采样算子相同的情况下，基于包含复四元数傅里叶正变换和复四元数傅 里叶反变换的第一重建表达式，利用矢量凸集投影迭代方法对单频率切片对应的复四 元数进行矢量联合重建，得到相应正频率部分的单频率切片的多分量三维地震数据的 重建结果的矢量联合表示，从重建结果的矢量联合表示中提取得到重建结果的各分量 数据。
上述步骤S141中，第一重建表达式可以为：
Figure BDA0003067454800000101
其中，Ql(fn,x1,x2)表示单频率切片对应的复四元数的第l次重建结果， l＝1,2,...,Niter，Niter表示最大迭代次数，fn表示第n个频率切片，x1和x2表 示检波器在地面上的两个空间坐标变量，α表示加权因子，Qobs(fn,x1,x2)表示 单频率切片对应的初始复四元数，I表示全1矩阵，S表示采样算子，
Figure BDA0003067454800000102
表示复 四元数傅里叶反变换，Tl表示第l次重建的阈值算子，FQ表示复四元数傅里叶变 换，Ql-1(fn,x1,x2)表示单频率切片对应的复四元数的第l-1次重建结果， Q0(fn,x1,x2)＝Qobs(fn,x1,x2)。
其中，Qobs(fn,x1,x2)、
Figure BDA0003067454800000103
等均体现了利用各分量进 行联合重建，而不是一个分量一个分量单独进行重建。对于构建得到复四元数矩阵的 情况，I为全1矩阵，
Figure BDA0003067454800000104
为相应四元矩阵。
在另一些实施例中，在多分量三维地震数据的各分量数据的采样算子存在不同的情况下，上述步骤S140，即，基于复四元数傅里叶变换，利用矢量凸集投影迭代方 法对单频率切片对应的复四元数进行矢量联合重建，得到相应正频率部分的单频率切 片的多分量三维地震数据的重建结果，具体可包括步骤：S142，在多分量三维地震数 据的各分量数据的采样算子存在不同的情况下，基于包含复四元数傅里叶正变换和复 四元数傅里叶反变换的第二重建表达式，利用矢量凸集投影迭代方法对单频率切片对 应的复四元数进行矢量联合重建，得到相应正频率部分的单频率切片的多分量三维地 震数据的重建结果的各分量数据。
上述步骤S142中，第二重建表达式可以为：
Figure BDA0003067454800000111
其中，
Figure BDA0003067454800000112
表示单频率切片的三维地震数据的分量υ的第l次重建结果，υ取三分量中的分量x、y、z或四分量中的分量p、x、y、z，l＝1,2,...,Niter， Niter表示最大迭代次数，fn表示第n个频率切片，x1和x2表示检波器在地面上的 两个空间坐标变量，α表示加权因子，
Figure BDA0003067454800000113
表示单频率切片的三维地 震数据的分量υ的初始数据，I表示全1矩阵，Sυ表示分量υ的采样算子，
Figure BDA0003067454800000114
表 示复四元数傅里叶反变换，Tl表示第l次重建的阈值算子，FQ表示复四元数傅里 叶变换，Ql-1(fn,x1,x2)表示单频率切片对应的复四元数的第l-1次重建结果， Q0(fn,x1,x2)＝Qobs(fn,x1,x2)，
Figure BDA0003067454800000115
表示从
Figure BDA0003067454800000116
取出分量υ。
其中，
Figure BDA0003067454800000117
体现了利用各分量进行联合重建，而不 是一个分量一个分量单独进行重建。
具体实施例中，上述各实施例的阈值算子可以为线性、指数或数据驱动的阈值模型。可以根据实际情况进行选择。
进一步的实施例中，上述步骤S140，即，基于复四元数傅里叶变换，利用矢量 凸集投影迭代方法对单频率切片对应的复四元数进行矢量联合重建，得到相应正频率 部分的单频率切片的多分量三维地震数据的重建结果，包括：S143，对单频率切片对 应的复四元数的迭代重建结果做复四元数傅里叶变换。具体地，可以是上述步骤S141 或步骤S142中包含该步骤S143。如第一重建表达式和第二重建表达式中的 FQQl-1(fn,x1,x2)。
具体地，该步骤S143，即，对单频率切片对应的复四元数的迭代重建结果做复 四元数傅里叶变换，具体可包括步骤：
S1431，将单频率切片对应的复四元数由单位正交基的表示形式变换至由包含一个复数型实部和三个复数型虚部的广义基表示形式，得到以广义基表示的第一复四元 数；
S1432，取广义基表示的单频率切片对应的复四元数的复数型实部的实数实部和实数虚部，得到第一实数实部和第一实数虚部；取广义基表示的单频率切片对应的复 四元数的第一复数型虚部的实数实部和实数虚部，得到第二实数实部和第二实数虚 部；取广义基表示的单频率切片对应的复四元数的第二复数型虚部的实数实部和实数 虚部，得到第三实数实部和第三实数虚部；取广义基表示的单频率切片对应的复四元 数的第三复数型虚部的实数实部和实数虚部，得到第四实数实部和第四实数虚部；
S1433，以第一实数实部作为实部并以第二实数实部作为虚部构建第一普通复数；以第一实数虚部作为实部并以第二实数虚部作为虚部构建第二普通复数；以第三实数 实部作为实部并以第四实数实部作为虚部构建第三普通复数；以第三实数虚部作为实 部并以第四实数虚部作为虚部构建第四普通复数；
S1434，对第一普通复数、第二普通复数、第三普通复数及第四普通复数分别做 普通傅里叶变换，得到第一普通复数的普通傅里叶变换结果、第二普通复数的普通傅 里叶变换结果、第三普通复数的普通傅里叶变换结果及第四普通复数的普通傅里叶变 换结果；
S1435，以第一普通复数的普通傅里叶变换结果的实数实部作为实部并以第二普通复数的普通傅里叶变换结果的实数实部作为虚部，构建第五普通复数；以第一普通 复数的普通傅里叶变换结果的实数虚部作为实部并以第二普通复数的普通傅里叶变 换结果的实数虚部作为虚部，构建第六普通复数；以第三普通复数的普通傅里叶变换 结果的实数实部作为实部并以第四普通复数的普通傅里叶变换结果的实数实部作为 虚部，构建第七普通复数；以第三普通复数的普通傅里叶变换结果的实数虚部作为实 部并以第四普通复数的普通傅里叶变换结果的实数虚部作为虚部，构建第八普通复 数；
S1436，以第五普通复数作为实部，以第六普通复数、第七普通复数及第八普通 复数作为三个虚部，构建以广义基表示的第二复四元数，并将以广义基表示的第二复 四元数变换至以单位正交基表示的复数四元数，作为单频率切片对应的复四元数的迭 代重建结果的复四元数傅里叶变换的结果。
举例而言，复四元数傅里叶变换可包括以下步骤：
(1)根据前述复正交矩阵以及式q＝a+bi+cj+dk和式q′＝a′+b′μ+c′ν+d′ξ满 足的对应关系，可以将复四元数q转化成基为(1,μ,ν,ξ)的q′形式， q′＝a′+b′μ+c′ν+d′ξ；
(2)对q′＝a′+b′μ+c′ν+d′ξ中的a′，b′，c′，d′取实部
Figure RE-GDA0003112777830000131
和取虚部
Figure RE-GDA0003112777830000132
并构建复数，得到
Figure RE-GDA0003112777830000133
Figure RE-GDA0003112777830000134
并对q′1,q′2,q′3,q′4做普通傅里叶变换，得到变换后的数据，记作 Q′1,Q′2,Q′3,Q′4；
(3)对Q′1,Q′2,Q′3,Q′4分别取实部
Figure RE-GDA0003112777830000135
和取虚部
Figure RE-GDA0003112777830000136
并构建复数，得到
Figure RE-GDA0003112777830000137
Figure RE-GDA0003112777830000138
构建复四元数Q′＝Q1+Q2μ+Q3ν+Q4ξ,利用复正交矩阵以及 式q＝a+bi+cj+dk和式q′＝a′+b′μ+c′ν+d′ξ满足的对应关系将Q′转化成 基为(1,i,j,k)的形式，得到q＝a+bi+cj+dk的复四元数傅里叶变换的结果： Q＝A+Bi+Cj+Dk，其中，A,B,C,D与Q1,Q2,Q3,Q4的对应关系如下，
A＝Q1
B＝Q2μ1+Q3ν1+Q4ξ1
C＝Q2μ2+Q3ν2+Q4ξ2
D＝Q2μ3+Q3ν3+Q4ξ3
通过以上步骤(1)～(3)即可完成复四元数傅里叶变换过程。
在另一些实施例中，复四元数傅里叶反变换可以是类似于上述步骤S1431～S1436的逆过程。
步骤S150：对正频率部分的单频率切片的多分量三维地震数据的重建结果中的各分量数据取共轭，得到相应负频率部分的单频率切片的多分量三维地震数据重建结 果。
该步骤S150，对于得到的重建结果为各分量联合结果的情况，可以先提取分量 再进行后续步骤，如上述采样算子相同的情况。对于得到的重建结果为分量的情况下， 可直接进行后续步骤。得到的正频率部分的单频切片为复数，取复数共轭可以得到负 频率部分的单频切片的重建结果。
步骤S160：根据正频率部分的单频率切片的多分量三维地震数据的重建结果和相应负频率部分的单频率切片的多分量三维地震数据的重建结果，得到整个频率域的 多分量三维地震数据的重建结果。
该步骤S160中，可以得到整个频率域所有频率切片的多分量三维地震数据的重建结果。
步骤S170：对整个频率域的多分量三维地震数据的重建结果中的各分量数据分别做普通傅里叶反变换，得到时间域的多分量三维地震数据的重建结果。
该步骤S170，普通傅里叶变换/普通傅里叶反变换主要是相对于本发明实施例中的复四元数傅里叶变换/复四元数傅里叶反变换而言的。普通傅里叶反变换是指常见 的傅里叶反变换，可直接称为傅里叶反变换。相应的，普通傅里叶变换是指常见的傅 里叶变换，可直接称为傅里叶变换。
另外，本发明实施例还提供了一种计算机设备，包括存储器、处理器及存储在存储器上并可在处理器上运行的计算机程序，所述处理器执行所述程序时实现上述任一 实施例所述矢量凸集投影多分量三维地震数据重建方法的步骤。
本发明实施例还提供了一种计算机可读存储介质，其上存储有计算机程序，该程序被处理器执行时实现上述任一实施例所述矢量凸集投影多分量三维地震数据重建 方法的步骤。
本发明实施例的矢量凸集投影多分量三维地震数据重建方法、计算机设备及计算机可读存储介质，采用复数四元数傅里叶变换方法，在频率域能将三分量数据或四分 量数据联合成复数四元数形式的矢量数据，然后，借助复数四元数傅里叶变换和凸集 投影算法实现三分量地震数据的联合重建。而且采用复数四元数傅里叶变换能够保持 地震数据在频率域正负频率切片的共轭对称性，在重建时只需要对正频率切片进行迭 代重建，负频率部分通过对正频率部分的重建结果取共轭对称即可得到，减少了数据 重建时的计算量，提高了重建计算效率。进一步，本发明实施例的方法能够实现不同 分量各自具有不同缺失模式的重建，即当每个分量的缺失度不同或者相同缺失度下每 个分量的缺失道位置不同时多个分量数据的联合重建。
下面结合一个具体实施例对上述方法进行说明，然而，值得注意的是，该具体实施例仅是为了更好地说明本申请，并不构成对本申请的不当限定。
在一具体实施例中，以三分量数据重建为例，基于复数四元数傅里叶变换的矢量凸集投影重建技术实现多分量三维不规则缺失地震数据的规则化重建的具体实施方 式。实现三分量三维不规则缺失地震数据的规则化重建具体可包括以下步骤:(1)对 三分量三维地震数据做关于时间t的傅里叶变换，将地震数据由时间域变换到频率f 域；(2)借助复四元数工具，对正频率部分的任一单个频率切片数据进行联合，构建 复四元数矩阵，实现三分量数据的矢量联合；(3)运用复四元数傅里叶变换和凸集投 影迭代技术实现三分量数据的矢量联合重建，得到三分量数据重建结果；(4)对正频 率部分的每一个频率切片按照步骤(2)和(3)进行重建，并对正频率切片取共轭得 到负频率切片的重建结果；(5)对频率f做傅里叶反变换到时间t域，得到时间域的 重建数据，完成重建过程。
举例而言，可以在共炮点域、共接收点域或者共中心点-偏移距域实施重建过程。以三分量共炮点域数据为例，实施方式可包括以下过程：
(1)设Dx(t,x1,x2)，Dy(t,x1,x2)和Dz(t,x1,x2)分别表示三分量勘探所采集的共炮点域X、Y和Z分量数据，将三分量数据沿时间t维做普通傅里叶变换到频率f域， 得到频率域三分量数据
Figure BDA0003067454800000151
(2)提取三个分量正频率部分的任意单频切片数据
Figure BDA0003067454800000152
Figure BDA0003067454800000153
利用复四元数框架，将三个分量的频率切片数据进行矢量联合表示， 记为Q(fn,x1,x2)，Q(fn,x1,x2)可以表示为：
Figure BDA0003067454800000154
其中，对于四分量数据的重建而言，可将第四个分量的单频率切片数据替换上式右边为0处的部分。
为了书写简洁，我们忽略频率变量fn，空间变量x1和x2(检波器的在地面的横 纵坐标，一个检波器对应多个频率切片/采样点数据)，则Q(fn,x1,x2)的表达式可以简 写为：
Figure BDA0003067454800000161
需要指出的是，
Figure BDA0003067454800000162
和
Figure BDA0003067454800000163
均为复数，因此，构建的四元数矩阵Q为一个复数 型四元数矩阵；
(3)设
Figure BDA0003067454800000164
和
Figure BDA0003067454800000165
为频率域原始观测的不规则采样数据，S为三个分量的 共同采样算子，其元素为0和1,0表示缺失道，1表示已知道。记
Figure BDA0003067454800000166
借助复数四元数傅里叶变换和凸集投影迭代算法进行 不规则数据的规则化重建，第l次重建表达式为：
Figure BDA0003067454800000167
其中，Ql表示第l次迭代重建结果，
Figure BDA0003067454800000168
Q0＝Qobs，I表示 全1矩阵，FQ和
Figure BDA0003067454800000169
分别表示空间二维复四元数傅里叶正反变换，Tl表示阈值算子， Tl有线性、指数和数据驱动三种阈值模型可以选择。α表示加权因子，0＜α≤1，用 于控制原始已知道加入重建数据的比重，进而实现去噪重建，Niter表示最大迭代次数。 若三个分量数据
Figure BDA00030674548000001610
和
Figure BDA00030674548000001611
各自的缺失模式不同，即采样算子S取不同值，分别为 Sx，Sy和Sz，此时三分量数据的矢量凸集投影重建表达式可以写为：
Figure BDA00030674548000001612
其中，υ分别取p，x，y和z，
Figure BDA00030674548000001613
表示取出四 元数数组中υ所对应分量的数据。
(4)对正频率部分的每一个频率切片fn按照步骤(2)和(3)进行重建，再对 正频率部分取共轭得到负频率部分的重建结果，进而得到整个频带内的三分量数据重 建结果
Figure BDA00030674548000001614
和
Figure BDA00030674548000001615
(5)对频率f做普通傅里叶反变换到时间t域，得到时间域重建数据Dx(t,x1,x2)，Dy(t,x1,x2)和Dz(t,x1,x2)，完成重建。
举例而言，设计一个三层水平层状地质模型，各层的厚度h，密度ρ，纵波速度Vp，和横波速度VS等参数设置如表1所示。设置199条测线，每条测线199道，每道含时间 采样点2001个，时间采样间隔1ms。然后，以模型中心点为炮点进行激发，合成三分 量共炮点地震数据Dx(t,x1,x2)，Dy(t,x1,x2)和Dz(t,x1,x2)，如图2(a)，图3(a)和图 4(a)所示。x1和x2分别表示检波器在空间Inline和Crossline方向的位置。
表1水平地层模型物理参数表
参数 h(m) ρ(kg/m<sup>3</sup>) V<sub>p</sub>(m/s) V<sub>s</sub>(m/s)
第一层 600 2000 2000 1000
第二层 600 2100 2500 1250
第三层 600 2500 2800 1500
随机剔除地震道，使得X、Y和Z分量的缺失度达到60％，但保持缺失道的位置在三个分量中互不相同，形成具有不同缺失模式的三分量数据，记为
Figure BDA0003067454800000171
Figure BDA0003067454800000172
和
Figure BDA0003067454800000174
三个分量各自的缺失度记为Sx，Sy和Sz。图2中的(a) 图和(b)图分别显示了X分量原始完整数据Dx(t,x1,x2)和60％随机缺失道数据
Figure BDA0003067454800000175
图3中的(a)图和(b)图分别显示了Y分量原始完整单炮记录Dy(t,x1,x2) 和随机缺失道数据
Figure BDA0003067454800000176
图4中的(a)图和(b)图分别显示Z分量原始完整 单炮记录Dz(t,x1,x2)和缺失道数据
Figure BDA0003067454800000177
对三分量数据做关于时间t的普通傅 里叶变换到频率域，然后用三分量数据矢量凸集投影重建方法进行迭代重建，图5中 的(a)图和(b)图分别显示了X分量用单分量凸集投影方法的重建结果和采用本实 施例的矢量凸集投影重建结果。图6中的(a)图和(b)图分别显示了X分量单分量凸 集投影重建的误差和本实施例的矢量凸集投影重建的误差。图7中的(a)图和(b) 图分别显示了Y分量采用单分量凸集投影方法的重建结果和本实施例的矢量凸集投影 重建结果。图8中的(a)图和(b)图分别显示了Y分量单分量凸集投影重建的误差和 本实施例矢量凸集投影重建的误差。图9中的(a)图和(b)图分别显示了Z分量采用 单分量凸集投影方法的重建结果和本实施例的矢量凸集投影重建结果。图10中的(a)图和(b)图分别显示了Z分量单分量凸集投影重建的误差和本实施例的矢量凸集投影 重建的误差。从图5中的(a)和(b)，图7中的(a)和(b)以及图9中的(a)和(b) 可以看出无论单分量凸集投影重建还是矢量凸集投影重建，不规则缺失道数据都被有 效恢复，同相轴恢复了连续性。但是，对比残差剖面图6中的(a)和(b)，图8中的 (a)和(b)以及图10中的(a)和(b)可以看出单分量重建的差剖面中有效信号的残留 量大，能量痕迹明显，矢量凸集投影重建结果的有效信号残留量小，差剖面干净。这 表明本实施例的矢量凸集投影可以对具有各自不同缺失模式的三分量数据进行有效 重建，而且重建效果优于单分量重建方法。
为了进一步显示本实施例的三分量重建技术具有良好的保持地震波波场结构完整性的优势，在图2的(a)图，图3的(a)图，图4的(a)图，图5的(a)图和(b) 图，图7的(a)图和(b)图以及图9的(a)图和(b)图中取出InLine＝50和CrossLine＝45 交叉点处的单道，截取时间窗在[557,706]ms处的地震道波形，绘制矢端图。图11 显示了使用单分量重建技术得到的X、Y和Z三个分量重建结果与原始真实的X、Y和Z 三个分量地震波波形对比的矢端图，图12显示了使用本实施例重建技术得到的X、Y 和Z三个分量重建结果与原始真实的X、Y和Z三个分量地震波波形对比矢端图，从图11 上可以看到用单分量重建技术得出的波形矢端图(虚线)与原始真实数据的矢端图轮 廓(实线)存在明显偏差，而从图12上可以看到采用本发明的矢量凸集投影重建出的 地震波波形矢端图(虚线)与原始真实数据矢端图重合较好，这表明本实施例的重建 技术对地震波波场结构保持较好，经该技术重建后的数据能有效保持各分量之间的相 互关系，获得比单分量重建技术更为精细和准确的重建结果。
以表1合成的三分量模型数据为基础，令X、Y和Z分量分别缺失10％，20％,...,80％， 然后利用单分量凸集投影技术，二分量凸集投影技术和本实施例的矢量凸集投影技术 对不同缺失度下的三个分量数据进行重建，同时定义重建质量因子Q，
Figure BDA0003067454800000181
Drecon表示重建数据，Dtrue表示原始真实 数据，三种方法的重建质量因子Q曲线如图13～15所示。从图上可以看出在缺失度为 10％-60％时，本实施例的三分量矢量凸集投影重建技术的重建效果优于单分量凸集投 影重建技术和二分量凸集投影重建技术。在70％-80％的严重缺失状态下本发明的重建 效果不低于二分量凸集投影重建技术，但均优于单分量凸集投影重建技术。
最后将本实施例应用于某区块三维三分量VSP(垂直地震剖面)实际资料插值重建处理。该VSP资料用取walk-around方式采集，井筒位于圆心处，在井下布置26道检 波器接收。选取地表水平方向26炮，垂直方向32炮的矩形数据区域进行插值重建，该 区域内含有7个缺失炮。分别采用单分量凸集投影方法和本发明的矢量凸集投影技术 对共接收点道集进行重建，共重建26次。最后，从重建的26个共接收点道集中，提取 出缺失炮对应的共炮点道集。单分量和三分量矢量凸集投影重建结果如图16和图17 所示。其中，图16中(a)图为R分量，(b)图为T分量，(c)图为Z分量。图17中(a) 图为R分量，(b)图为T分量，(c)图为Z分量。对比单分量和多分量重建结果可以看 出，在本实施例的矢量凸集投影重建结果中地震波的波组特征明显，同相轴连续清晰， 并与附近已知炮的共炮点道集中地震波的波组结构特征相似，而在单分量重建结果中 同向轴的连续性较差，波组特征不明显。这充分说明了本实施方法的先进性。
本实施例中，对三分量不规则缺失地震数据采用基于复数四元数傅里叶变换的矢量凸集投影重建技术进行规则化重建，该方法不仅能够实现三个或者四分量数据的联 合重建而且能够实现三个或者四个分量在各自不同缺失模式下进行联合重建。另外， 该方法采用的是复数四元数傅里叶变换不同于实数四元数傅里叶变换，复数四元数傅 里叶变换能够保持地震数据在频率域正负频率切片的共轭对称性，在重建时只需要对 正频率切片进行迭代重建，负频率部分通过对正频率部分的重建结果取共轭对称即可 得到，不再需要对负频率切片进行重建，因此，减少了计算量，提高了重建计算效率。
在本说明书的描述中，参考术语“一个实施例”、“一个具体实施例”、“一些 实施例”、“例如”、“示例”、“具体示例”、或“一些示例”等的描述意指结合 该实施例或示例描述的具体特征、结构、材料或者特点包含于本发明的至少一个实施 例或示例中。在本说明书中，对上述术语的示意性表述不一定指的是相同的实施例或 示例。而且，描述的具体特征、结构、材料或者特点可以在任何的一个或多个实施例 或示例中以合适的方式结合。各实施例中涉及的步骤顺序用于示意性说明本发明的实 施，其中的步骤顺序不作限定，可根据需要作适当调整。
本领域内的技术人员应明白，本发明的实施例可提供为方法、系统、或计算机程序产品。因此，本发明可采用完全硬件实施例、完全软件实施例、或结合软件和硬件 方面的实施例的形式。而且，本发明可采用在一个或多个其中包含有计算机可用程序 代码的计算机可用存储介质(包括但不限于磁盘存储器、CD-ROM、光学存储器等) 上实施的计算机程序产品的形式。
本发明是参照根据本发明实施例的方法、设备(系统)、和计算机程序产品的流 程图和/或方框图来描述的。应理解可由计算机程序指令实现流程图和/或方框图中的 每一流程和/或方框、以及流程图和/或方框图中的流程和/或方框的结合。可提供这些 计算机程序指令到通用计算机、专用计算机、嵌入式处理机或其他可编程数据处理设 备的处理器以产生一个机器，使得通过计算机或其他可编程数据处理设备的处理器执 行的指令产生用于实现在流程图一个流程或多个流程和/或方框图一个方框或多个方 框中指定的功能的装置。
这些计算机程序指令也可存储在能引导计算机或其他可编程数据处理设备以特定方式工作的计算机可读存储器中，使得存储在该计算机可读存储器中的指令产生包 括指令装置的制造品，该指令装置实现在流程图一个流程或多个流程和/或方框图一 个方框或多个方框中指定的功能。
这些计算机程序指令也可装载到计算机或其他可编程数据处理设备上，使得在计算机或其他可编程设备上执行一系列操作步骤以产生计算机实现的处理，从而在计算 机或其他可编程设备上执行的指令提供用于实现在流程图一个流程或多个流程和/或 方框图一个方框或多个方框中指定的功能的步骤。
以上所述的具体实施例，对本发明的目的、技术方案和有益效果进行了进一步详细说明，所应理解的是，以上所述仅为本发明的具体实施例而已，并不用于限定本发 明的保护范围，凡在本发明的精神和原则之内，所做的任何修改、等同替换、改进等， 均应包含在本发明的保护范围之内。

Claims (8)

1.一种矢量凸集投影多分量三维地震数据重建方法，其特征在于，包括：
对时间域的多分量三维地震数据的各分量数据做关于时间变量的普通傅里叶变换，得到频率域的多分量三维地震数据的各分量数据；其中，所述多分量为三分量或四分量；
从频率域的多分量三维地震数据的各分量数据中分别提取正频率部分的频率切片数据；
利用复四元数框架将从各分量数据中提取的正频率部分的单频率切片数据进行矢量联合表示，得到相应正频率部分的单频率切片对应的复四元数；
其中，在多分量为三分量的情况下，将从各分量数据提取的正频率部分的单频率切片数据中的每个分量作为复四元数框架中的一个虚部，得到相应正频率部分的单频率切片对应的复四元数；在多分量为四分量的情况下，将从各分量数据提取的正频率部分的单频率切片数据中任选其中一个分量作为复四元数框架中的实部，其余三个分量中的每个分量作为复四元数框架中的一个虚部，得到相应正频率部分的单频率切片对应的复四元数；基于复四元数傅里叶变换，利用矢量凸集投影迭代方法对单频率切片对应的复四元数进行矢量联合重建，得到相应正频率部分的单频率切片的多分量三维地震数据的重建结果；
对正频率部分的单频率切片的多分量三维地震数据的重建结果中的各分量数据取共轭，得到相应负频率部分的单频率切片的多分量三维地震数据重建结果；
根据正频率部分的单频率切片的多分量三维地震数据的重建结果和相应负频率部分的单频率切片的多分量三维地震数据的重建结果，得到整个频率域的多分量三维地震数据的重建结果；
对整个频率域的多分量三维地震数据的重建结果中的各分量数据分别做普通傅里叶反变换，得到时间域的多分量三维地震数据的重建结果。
2.如权利要求1所述的矢量凸集投影多分量三维地震数据重建方法，其特征在于，基于复四元数傅里叶变换，利用矢量凸集投影迭代方法对单频率切片对应的复四元数进行矢量联合重建，得到相应正频率部分的单频率切片的多分量三维地震数据的重建结果，包括：
在多分量三维地震数据的各分量数据的采样算子相同的情况下，基于包含复四元数傅里叶正变换和复四元数傅里叶反变换的第一重建表达式，利用矢量凸集投影迭代方法对单频率切片对应的复四元数进行矢量联合重建，得到相应正频率部分的单频率切片的多分量三维地震数据的重建结果的矢量联合表示，从重建结果的矢量联合表示中提取得到重建结果的各分量数据；
第一重建表达式为：
Figure FDA0003962608210000021
其中，Ql(fn,x1,x2)表示单频率切片对应的复四元数的第l次重建结果，l＝1,2,...,Niter，Niter表示最大迭代次数，fn表示第n个频率切片，x1和x2表示检波器在地面坐标的两个空间变量，α表示加权因子，Qobs(fn,x1,x2)表示单频率切片对应的初始复四元数，I表示全1矩阵，S表示采样算子，
Figure FDA0003962608210000022
表示复四元数傅里叶反变换，Tl表示第l次重建的阈值算子，FQ表示复四元数傅里叶变换，Ql-1(fn,x1,x2)表示单频率切片对应的复四元数的第l-1次重建结果，Q0(fn,x1,x2)＝Qobs(fn,x1,x2)。
3.如权利要求1所述的矢量凸集投影多分量三维地震数据重建方法，其特征在于，基于复四元数傅里叶变换，利用矢量凸集投影迭代方法对单频率切片对应的复四元数进行矢量联合重建，得到相应正频率部分的单频率切片的多分量三维地震数据的重建结果，包括：
在多分量三维地震数据的各分量数据的采样算子存在不同的情况下，基于包含复四元数傅里叶正变换和复四元数傅里叶反变换的第二重建表达式，利用矢量凸集投影迭代方法对单频率切片对应的复四元数进行矢量联合重建，得到相应正频率部分的单频率切片的多分量三维地震数据的重建结果的各分量数据；
第二重建表达式为：
Figure FDA0003962608210000023
其中，
Figure FDA0003962608210000024
表示单频率切片的三维地震数据的分量υ的第l次重建结果，υ取三分量中的分量x、y、z或四分量中的分量p、x、y、z，l＝1,2,...,Niter，Niter表示最大迭代次数，fn表示第n个频率切片，x1和x2表示检波器在地面坐标的两个空间变量，α表示加权因子，
Figure FDA0003962608210000031
表示单频率切片的三维地震数据的分量υ的初始数据，I表示全1矩阵，Sυ表示分量υ的采样算子，
Figure FDA0003962608210000032
表示复四元数傅里叶反变换，Tl表示第l次重建的阈值算子，FQ表示复四元数傅里叶变换，Ql-1(fn,x1,x2)表示单频率切片对应的复四元数的第l-1次重建结果，Q0(fn,x1,x2)＝Qobs(fn,x1,x2)，
Figure FDA0003962608210000033
表示从
Figure FDA0003962608210000034
中取出分量υ。
4.如权利要求1所述的矢量凸集投影多分量三维地震数据重建方法，其特征在于，基于复四元数傅里叶变换，利用矢量凸集投影迭代方法对单频率切片对应的复四元数进行矢量联合重建，得到相应正频率部分的单频率切片的多分量三维地震数据的重建结果，包括：对单频率切片对应的复四元数的迭代重建结果做复四元数傅里叶变换；
对单频率切片对应的复四元数的迭代重建结果做复四元数傅里叶变换，包括：
将单频率切片对应的复四元数由单位正交基的表示形式变换至由包含一个复数型实部和三个复数型虚部的广义基表示形式，得到以广义基表示的第一复四元数；
取广义基表示的单频率切片对应的复四元数的复数型实部的实数实部和实数虚部，得到第一实数实部和第一实数虚部；取广义基表示的单频率切片对应的复四元数的第一复数型虚部的实数实部和实数虚部，得到第二实数实部和第二实数虚部；取广义基表示的单频率切片对应的复四元数的第二复数型虚部的实数实部和实数虚部，得到第三实数实部和第三实数虚部；取广义基表示的单频率切片对应的复四元数的第三复数型虚部的实数实部和实数虚部，得到第四实数实部和第四实数虚部；
以第一实数实部作为实部并以第二实数实部作为虚部构建第一普通复数；以第一实数虚部作为实部并以第二实数虚部作为虚部构建第二普通复数；以第三实数实部作为实部并以第四实数实部作为虚部构建第三普通复数；以第三实数虚部作为实部并以第四实数虚部作为虚部构建第四普通复数；
对第一普通复数、第二普通复数、第三普通复数及第四普通复数分别做普通傅里叶变换，得到第一普通复数的普通傅里叶变换结果、第二普通复数的普通傅里叶变换结果、第三普通复数的普通傅里叶变换结果及第四普通复数的普通傅里叶变换结果；
以第一普通复数的普通傅里叶变换结果的实数实部作为实部并以第二普通复数的普通傅里叶变换结果的实数实部作为虚部，构建第五普通复数；以第一普通复数的普通傅里叶变换结果的实数虚部作为实部并以第二普通复数的普通傅里叶变换结果的实数虚部作为虚部，构建第六普通复数；以第三普通复数的普通傅里叶变换结果的实数实部作为实部并以第四普通复数的普通傅里叶变换结果的实数实部作为虚部，构建第七普通复数；以第三普通复数的普通傅里叶变换结果的实数虚部作为实部并以第四普通复数的普通傅里叶变换结果的实数虚部作为虚部，构建第八普通复数；
以第五普通复数作为实部，以第六普通复数、第七普通复数及第八普通复数作为三个虚部，构建以广义基表示的第二复四元数，并将以广义基表示的第二复四元数变换至以单位正交基表示的复数四元数，作为单频率切片对应的复四元数的迭代重建结果的复四元数傅里叶变换的结果。
5.如权利要求2或3所述的矢量凸集投影多分量三维地震数据重建方法，其特征在于，阈值算子为线性、指数或数据驱动的阈值模型。
6.如权利要求1所述的矢量凸集投影多分量三维地震数据重建方法，其特征在于，还包括：
将得到的各分量单频率切片对应的复四元数按照检波器在地面不同位置处的横坐标和纵坐标两个空间变量维度表示成复四元数矩阵，以用于矢量联合重建。
7.一种计算机设备，包括存储器、处理器及存储在存储器上并可在处理器上运行的计算机程序，其特征在于，所述处理器执行所述程序时实现如权利要求1至6任一项所述方法的步骤。
8.一种计算机可读存储介质，其上存储有计算机程序，其特征在于，该程序被处理器执行时实现如权利要求1至6任一项所述方法的步骤。
CN202110529379.5A 2021-05-14 2021-05-14 矢量凸集投影多分量三维地震数据重建方法及装置 Active CN114509805B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202110529379.5A CN114509805B (zh) 2021-05-14 2021-05-14 矢量凸集投影多分量三维地震数据重建方法及装置

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202110529379.5A CN114509805B (zh) 2021-05-14 2021-05-14 矢量凸集投影多分量三维地震数据重建方法及装置

Publications (2)

Publication Number Publication Date
CN114509805A CN114509805A (zh) 2022-05-17
CN114509805B true CN114509805B (zh) 2023-02-28

Family

ID=81547636

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202110529379.5A Active CN114509805B (zh) 2021-05-14 2021-05-14 矢量凸集投影多分量三维地震数据重建方法及装置

Country Status (1)

Country Link
CN (1) CN114509805B (zh)

Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN105425301A (zh) * 2016-01-08 2016-03-23 东华理工大学 一种频率域三维不规则地震数据重建方法
CN106842321A (zh) * 2015-12-04 2017-06-13 中国石油天然气股份有限公司 地震数据重建方法和装置
CN108983287A (zh) * 2018-10-23 2018-12-11 东华理工大学 一种基于凸集投影算法的曲波变换反假频地震数据重建方法
CN109164483A (zh) * 2018-08-29 2019-01-08 中国科学院地球化学研究所 多分量地震数据矢量去噪方法及多分量地震数据矢量去噪装置
CN112445649A (zh) * 2019-08-30 2021-03-05 中国石油化工股份有限公司 地震缺失数据恢复方法、计算机存储介质和计算机系统

Family Cites Families (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US8103453B2 (en) * 2008-05-21 2012-01-24 Bp Corporation North America Inc. Method of seismic data interpolation by projection on convex sets
US8352190B2 (en) * 2009-02-20 2013-01-08 Exxonmobil Upstream Research Company Method for analyzing multiple geophysical data sets
CN106646612B (zh) * 2016-12-20 2018-11-30 中国地质大学（北京） 基于矩阵降秩的地震数据重建方法
CN110764135A (zh) * 2018-07-26 2020-02-07 中国石油化工股份有限公司 不规则地震数据全频带重建方法

Patent Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN106842321A (zh) * 2015-12-04 2017-06-13 中国石油天然气股份有限公司 地震数据重建方法和装置
CN105425301A (zh) * 2016-01-08 2016-03-23 东华理工大学 一种频率域三维不规则地震数据重建方法
CN109164483A (zh) * 2018-08-29 2019-01-08 中国科学院地球化学研究所 多分量地震数据矢量去噪方法及多分量地震数据矢量去噪装置
CN108983287A (zh) * 2018-10-23 2018-12-11 东华理工大学 一种基于凸集投影算法的曲波变换反假频地震数据重建方法
CN112445649A (zh) * 2019-08-30 2021-03-05 中国石油化工股份有限公司 地震缺失数据恢复方法、计算机存储介质和计算机系统

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
"Irregular seismic data reconstruction based on exponential threshold model of POCS method";Gao jianjun 等;《Applied Geophysics》;20100930;第7卷(第3期);第229-238页 *
"Vector reconstruction of multicomponent seismic data";Aaron Stanton 等;《Geophysics》;20130831;第78卷(第4期);第V131-V145页 *

Also Published As

Publication number Publication date
CN114509805A (zh) 2022-05-17

Similar Documents

Publication Publication Date Title
CN105319581B (zh) 一种高效的时间域全波形反演方法
Jiang et al. A convolutional autoencoder method for simultaneous seismic data reconstruction and denoising
US20100054082A1 (en) Reverse-time depth migration with reduced memory requirements
CN109738952B (zh) 基于全波形反演驱动的被动源直接偏移成像方法
Sun et al. Deep learning for low-frequency extrapolation of multicomponent data in elastic FWI
Ovcharenko et al. Multi-task learning for low-frequency extrapolation and elastic model building from seismic data
CN113885077B (zh) 一种基于深度学习的多震源地震数据分离方法
CN105277985A (zh) 一种基于图像处理的ovt域地震数据规则化方法
CN114839673B (zh) 多震源高效采集波场分离方法、分离系统及计算机设备
CA3024240A1 (en) System and method for real-time passive seismic event localization
CN113269818A (zh) 一种基于深度学习的地震数据纹理特征重建方法
Kong et al. Deep prior-based seismic data interpolation via multi-res u-net
Zhang et al. VelocityGAN: Subsurface velocity image estimation using conditional adversarial networks
CN114325821A (zh) 基于3d-snacnn网络的叠前地震资料中强散射噪声压制方法及系统
López et al. Spectral gap-based seismic survey design
CN114509805B (zh) 矢量凸集投影多分量三维地震数据重建方法及装置
CN116719086B (zh) 基于点扩散函数的稀疏海底四分量数据高分辨率成像方法
Greiner et al. Cross-streamer wavefield interpolation using deep convolutional neural network
CN106950597A (zh) 基于三边滤波的混合震源数据分离方法
Trad Five dimensional seismic data interpolation
Rusmanugroho et al. 3D velocity model building based upon hybrid neural network
US20230066911A1 (en) Method for reconstructing at least one trace in a seismic image
CN115166823A (zh) 一种地震纵横波数据分离方法、模型构建方法及系统
CN111239828B (zh) 基于最优双曲线积分路径叠加的多次波压制方法
CN111830560B (zh) 一种基于降秩算法的地震数据重建方法

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