一种基于Huber范数的非均一采样地震数据鲁棒重构方法

专利名称：一种基于Huber范数的非均一采样地震数据鲁棒重构方法

专利类型：发明专利

专利申请号：CN202410313200.6

专利申请（专利权）人：哈尔滨工业大学
权利人地址：黑龙江省哈尔滨市南岗区西大直街92号

专利发明（设计）人：马坚伟,于四伟,徐英杰,陈尧

专利摘要：本发明提出一种基于Huber范数的非均一采样地震数据鲁棒重构方法，属于地震数据处理技术领域。包括：步骤一、对重构问题进行建模，根据地震数据采样点的坐标信息计算二维重心拉格朗日插值算子；具体为：获取地震数据的观测系统信息，通过实际非均一网格采样点和规则网格采样点的坐标信息计算插值算子，利用插值算子，重构目标函数；步骤二、设置最大迭代次数参数和阈值相关参数，计算阈值函数；步骤三、计算迭代误差ek和伪地震数据pk，更新当前数据yk；通过迭代阈值算法进行循环求解xk+1，直到最大迭代次数。解决现有技术中含有随机和异常值噪声情况下重构非均一采样的数据处理效果不好的技术问题。本发明能够起到在重构的同时压制随机噪声和异常值。

主权利要求：
1.一种基于Huber范数的非均一采样地震数据鲁棒重构方法，其特征在于，包括以下步骤：步骤一、对重构问题进行建模，根据地震数据采样点的坐标信息计算二维重心拉格朗日插值算子；
具体为：获取地震数据的观测系统信息，通过实际非均一网格采样点和规则网格采样点的坐标信息计算插值算子，利用插值算子，重构目标函数；
计算插值算子为：
其中，t＝x‑xi,u＝y‑yi,(x,y)是规则网格采样点，(xi,yi)是实际的非均一网格采样点，ωi,j为重心拉格朗日的权重系数：利用插值算子W，将重构目标函数写作：其中，d为观测数据，C为曲波变换，x为原始完整地震数据，ρc表示Huber范数：参数c＝1.345σ，其中，σ为绝对偏差中值，计算公式为：median表示中值函数，及选取该数据的中间值，ρc(x)的导数定义为：步骤二、设置最大迭代次数参数和阈值相关参数，计算阈值函数；
最大迭代次数参数为：30‑80次，最小阈值为：0.1‑0.4，最大阈值为：0.5‑0.9；
计算阈值函数为：
c(k‑1)/(N‑1)
λk＝τmaxe ,c＝ln(τmin/τmax),k＝1,2,...N(7)其中，τmax和τmin为人为设定的最大和最小阈值；
步骤三、计算迭代误差ek和伪地震数据pk，更新当前数据yk；通过迭代阈值算法进行循环求解xk+1，直到最大迭代次数；
迭代误差ek和伪地震数据pk的计算为：ek＝d‑Wxk‑1(8)
其中，xk‑1为第k‑1迭代出的重构数据，设x1＝x0＝d；
xk+1通过阈值函数更新：
T
yk＝xk‑A(Axk‑pk)(11)对公式(8)‑(12)的更新迭代，直到最大迭代次数为止。
2.一种电子设备，其特征在于，包括存储器和处理器，存储器存储有计算机程序，所述的处理器执行所述计算机程序时实现权利要求1所述的一种基于Huber范数的非均一采样地震数据鲁棒重构方法的步骤。
3.一种计算机可读存储介质，其特征在于，其上存储有计算机程序，所述计算机程序被处理器执行时实现权利要求1所述的一种基于Huber范数的非均一采样地震数据鲁棒重构方法。 说明书 : 一种基于Huber范数的非均一采样地震数据鲁棒重构方法技术领域[0001] 本发明涉及一种基于Huber范数的非均一采样地震数据鲁棒重构方法，属于地震数据处理技术领域。背景技术[0002] 目前，压缩感知(CS)采集在一个或两个空间方向上实现不规则采样，并在处理阶段重建缺失的样本，因此在工业界引起了广泛关注，并且也逐步在实际采集中运用这种方式。CS采集的主要挑战在于缺失数据的重构，特别是重构的准确性和鲁棒性。此外，现场施工条件可能会破坏采样几何，导致从设计的网格中产生缺失样本和偏离网格的样本，这会影响重构的准确性。压缩感知技术的关键思想是，在信号具有稀疏性的前提下，通过较少的测量就能恢复出完整的信号。这种方法在信号处理、图像处理和数据传输等领域都有重要应用。在理想情况下，CS通过对信号进行非均匀采样，降低了对数据采集的要求，从而减少了数据的存储和传输成本。然而，在实际应用中，由于各种因素如设备限制、环境干扰等，实际采集的样本可能会偏离理想的采样位置，产生所谓的“偏离网格”样本。这些偏离网格的样本，我们将其称之为非均一采样数据，它们会增加信号重构的难度，并且采样的数据中通常都含有随机噪声和异常值噪声，这使得构建一种针对非均一采样地震数据的鲁棒重构算法尤为重要。[0003] 传统的地震数据重构算法中主要有两类针对非均一采样数据设计的方法。一类是基于稀疏变换的方法，常用的有非均一傅里叶变换和非均一曲波变换。该类方法通过调整对稀疏变换的计算过程，适应非均匀采样的特点，从而能够更准确地分析和重构地震数据。还有一类方法是基于插值算子的方法，该类方法通过将地震数据的采样算子改进为插值算子，来实现对于非均一采样地震数据的重构。这两类方法在处理具有高噪声水平或含有异常值的数据时，在实际应用时重构效果都会受到限制。[0004] 近年来，基于Huber范数的鲁棒算法研究越来越多，该类算法通常设计得更加灵活和适应性强，能够应对各种不同的信号处理场景和挑战。而在地震数据重构中，现有的算法不能够很好的在含有随机和异常值噪声情况下重构非均一采样的数据。因此，研究一种能同时处理非均一采样数据并消除噪声和异常值干扰的地震数据鲁棒重构算法是很有必要的。发明内容[0005] 在下文中给出了关于本发明的简要概述，以便提供关于本发明的某些方面的基本理解。应当理解，这个概述并不是关于本发明的穷举性概述。它并不是意图确定本发明的关键或重要部分，也不是意图限定本发明的范围。其目的仅仅是以简化的形式给出某些概念，以此作为稍后论述的更详细描述的前序。[0006] 鉴于此，为解决现有技术中存在的含有随机和异常值噪声情况下重构非均一采样的数据处理效果不好的技术问题，本发明提供一种基于Huber范数的非均一采样地震数据鲁棒重构方法。[0007] 方案一、一种基于Huber范数的非均一采样地震数据鲁棒重构方法，包括以下步骤：[0008] 步骤一、对重构问题进行建模，根据地震数据采样点的坐标信息计算二维重心拉格朗日插值算子；[0009] 具体为：获取地震数据的观测系统信息，通过实际非均一网格采样点和规则网格采样点的坐标信息计算插值算子，利用插值算子，重构目标函数；[0010] 步骤二、设置最大迭代次数参数和阈值相关参数，计算阈值函数；[0011] 步骤三、计算迭代误差ek和伪地震数据pk，更新当前数据yk；通过迭代阈值算法进行循环求解xk+1，直到最大迭代次数。[0012] 优选的，步骤一具体是，计算插值算子为：[0013][0014] 其中，t＝x‑xi,u＝y‑yi,(x,y)是规则网格采样点，(xi,yi)是实际的非均一网格采样点，ωi,j为重心拉格朗日的权重系数：[0015][0016] 利用插值算子W，将重构目标函数写作：[0017][0018] 其中，d为观测数据，C为曲波变换，x为原始完整地震数据，ρc表示Huber[0019] 范数：[0020][0021] 参数c＝1.345σ，其中，σ为绝对偏差中值，计算公式为：[0022][0023] median表示中值函数，及选取该数据的中间值，ρc(x)的导数定义为：[0024][0025] 优选的，计算阈值函数为：[0026] λk＝τmaxec(k‑1)/(N‑1),c＝ln(τmin/τmax),k＝1,2,...N(7)[0027] 其中，τmax和τmin为人为设定的最大和最小阈值。[0028] 优选的，迭代误差ek和伪地震数据pk的计算为：[0029] ek＝d‑Wxk‑1(8)[0030][0031] 其中，xk‑1为第k‑1迭代出的重构数据，设x1＝x0＝d；[0032] xk+1通过阈值函数更新：[0033][0034] yk＝xk‑AT(Axk‑pk)(11)[0035][0036] 对公式(8)‑(12)的更新迭代，直到最大迭代次数为止。[0037] 方案二、一种电子设备，包括存储器和处理器，存储器存储有计算机程序，所述的处理器执行所述计算机程序时实现方案一所述的一种基于Huber范数的非均一采样地震数据鲁棒重构方法的步骤。[0038] 方案三、一种计算机可读存储介质，其上存储有计算机程序，所述计算机程序被处理器执行时实现方案一所述的一种基于Huber范数的非均一采样地震数据鲁棒重构方法。[0039] 本发明的有益效果如下：[0040] (1)本发明考虑地震数据采样的实际坐标信息，通过构建二维重心拉格朗日插值算子对非均一采样地震数据进行重构，在含有随机和异常值噪声情况下重构非均一采样的数据处理上效果很好；[0041] (2)本发明将Huber范数加入到重构模型中，使得该方法具有鲁棒性能够起到在重构的同时压制随机噪声和异常值的作用。[0042] (3)本发明利用迭代阈值算法进行模型的求解，整个过程较为简单，参数容易设定，且方法的适应范围广。附图说明[0043] 此处所说明的附图用来提供对本发明的进一步理解，构成本发明的一部分，本发明的示意性实施例及其说明用于解释本发明，并不构成对本发明的不当限定。在附图中：[0044] 图1为一种基于Huber范数的非均一采样地震数据鲁棒重构方法流程图；[0045] 图2为大小为260*30*30的模拟地震数据示意图，其中，a为模拟地震数据图，b为含有随机噪声和异常值的地震数据图；[0046] 图3为大小为260*30*30的模拟地震数据使用鲁棒凸集投影方法进行重构的结果和误差图，其中，a为重构结果图，b为误差图；[0047] 图4为大小为260*30*30的模拟地震数据使用本发明方法进行重构的结果和误差图，其中，a为重构结果图，b为误差图；[0048] 图5为大小为512*128*128的模拟地震数据示意图，其中，a为模拟地震数据图，b为含有随机噪声和异常值的地震数据图；；[0049] 图6为大小为512*128*128的模拟地震数据使用本发明方法进行重构的结果和误差图，其中，a为重构结果图，b为误差图；[0050] 图7为大小为512*128*128的模拟地震数据使用鲁棒凸集投影方法进行重构的结果和误差图，其中，a为重构结果图，b为误差图；[0051] 图8为一个单道的重构数据和频谱对比图，其中，a为单道重构数据对比图，b为单道重构数据的频谱对比图。具体实施方式[0052] 为了使本发明实施例中的技术方案及优点更加清楚明白，以下结合附图对本发明的示例性实施例进行进一步详细的说明，显然，所描述的实施例仅是本发明的一部分实施例，而不是所有实施例的穷举。需要说明的是，在不冲突的情况下，本发明中的实施例及实施例中的特征可以相互组合。[0053] 实施例1：参照图1‑图8说明本实施方式，一种基于Huber范数的非均一采样地震数据鲁棒重构方法，包括以下步骤：[0054] 步骤一、对重构问题进行建模，根据地震数据采样点的坐标信息计算二维重心拉格朗日插值算子；[0055] 具体为：获取地震数据的观测系统信息，通过实际非均一网格采样点和规则网格采样点的坐标信息计算插值算子，利用插值算子，重构目标函数；[0056] 步骤一具体是，计算插值算子为：[0057][0058] 其中，t＝x‑xi,u＝y‑yi,(x,y)是规则网格采样点，(xi,yi)是实际的非均一网格采样点，ωi,j为重心拉格朗日的权重系数：[0059][0060] 利用插值算子W，将重构目标函数写作：[0061][0062] 其中，d为观测数据，C为曲波变换，x为原始完整地震数据，ρc表示Huber[0063] 范数：[0064][0065] 参数c＝1.345σ，其中，σ为绝对偏差中值，计算公式为：[0066][0067] median表示中值函数，及选取该数据的中间值，ρc(x)的导数定义为：[0068][0069] 步骤二、设置最大迭代次数参数N和阈值相关参数，计算阈值函数；[0070] 最大迭代次数参数N的范围为：通常设置为30‑80次，在本实施例的实验中设置为50次。[0071] 最大和最小阈值的范围为：通常设置最小阈值为0.1‑0.4之间，最大阈值为0.5‑0.9之间。在本实施例的实验中设置最小和最大阈值分别为0.2和0.8。[0072] 计算阈值函数为：[0073] λk＝τmaxec(k‑1)/(N‑1),c＝ln(τmin/τmax),k＝1,2,...N(7)[0074] 其中，τmax和τmin为人为设定的最大和最小阈值。[0075] 步骤三、计算迭代误差ek和伪地震数据pk，更新当前数据yk；通过迭代阈值算法进行循环求解xk+1，直到最大迭代次数。[0076] 迭代误差ek和伪地震数据pk的计算为：[0077] ek＝d‑Wxk‑1(8)[0078][0079] 其中，xk‑1为第k‑1迭代出的重构数据，设x1＝x0＝d；[0080] xk+1通过阈值函数更新：[0081][0082] yk＝xk‑AT(Axk‑pk)(11)[0083][0084] 对公式(8)‑(12)的更新迭代，直到最大迭代次数为止。[0085] 以下使用鲁棒凸集投影方法与本发明对50％非均一采样的3维模拟地震数据进行插值操作，分别得到重构地震数据。根据信噪比对重构地震数据进行比较；[0086] 信噪比计算公式为：[0087][0088] 其中，X为完整地震数据，Xrec为重构地震数据。[0089] 图2是大小为260*30*30的模拟地震数据；图5是大小为512*128*128的模拟地震数据，其中，X轴方向道间距＝20m,Y轴方向道间距＝20m，时间采样点t＝0.004s；[0090] 图3‑图4是大小为260*30*30的模拟地震数据使用鲁棒凸集投影方法和本发明方法进行重构的结果和误差图；[0091] 通过图3‑图4可以明显的看到本发明方法的重构误差更小，几乎没有信号泄露，重构的同相轴连续性也更强，重构后的数据信噪比为15.97dB高于鲁棒凸集投影方法的11.83dB；[0092] 图6‑图7是大小为512*128*128的模拟地震数据使用鲁棒凸集投影方法和本发明方法进行重构的结果和误差图；[0093] 通过图6‑图7可以看出可以明显的看到本发明方法的重构误差主要由随机噪声构成，信号泄露极少。而于鲁棒凸集投影方法表现出有限的有效性，其重构误差明显。[0094] 为了进一步的比较，图8为两种方法对其中一个单道的重构数据和频谱对比图，图中表明本发明方法重构的数据和频谱与原始数据更为接近，误差更小。[0095] 实施例2：本发明的计算机装置可以是包括有处理器以及存储器等装置，例如包含中央处理器的单片机等。并且，处理器用于执行存储器中存储的计算机程序时实现上述的一种基于Huber范数的非均一采样地震数据鲁棒重构方法的步骤。[0096] 所称处理器可以是中央处理单元(CentralProcessingUnit，CPU)，还可以是其他通用处理器、数字信号处理器(DigitalSignalProcessor，DSP)、专用集成电路(ApplicationSpecificIntegratedCircuit，ASIC)、现成可编程门阵列(Field‑ProgrammableGateArray，FPGA)或者其他可编程逻辑器件、分立门或者晶体管逻辑器件、分立硬件组件等。通用处理器可以是微处理器或者该处理器也可以是任何常规的处理器等。[0097] 所述存储器可主要包括存储程序区和存储数据区，其中，存储程序区可存储操作系统、至少一个功能所需的应用程序(比如声音播放功能、图像播放功能等)等；存储数据区可存储根据手机的使用所创建的数据(比如音频数据、电话本等)等。此外，存储器可以包括高速随机存取存储器，还可以包括非易失性存储器，例如硬盘、内存、插接式硬盘，智能存储卡(SmartMediaCard,SMC)，安全数字(SecureDigital,SD)卡，闪存卡(FlashCard)、至少一个磁盘存储器件、闪存器件、或其他易失性固态存储器件。[0098] 实施例3：计算机可读存储介质实施例。[0099] 本发明的计算机可读存储介质可以是被计算机装置的处理器所读取的任何形式的存储介质，包括但不限于非易失性存储器、易失性存储器、铁电存储器等，计算机可读存储介质上存储有计算机程序，当计算机装置的处理器读取并执行存储器中所存储的计算机程序时，可以实现上述的一种基于Huber范数的非均一采样地震数据鲁棒重构方法的步骤。[0100] 所述计算机程序包括计算机程序代码，所述计算机程序代码可以为源代码形式、对象代码形式、可执行文件或某些中间形式等。所述计算机可读介质可以包括：能够携带所述计算机程序代码的任何实体或装置、记录介质、U盘、移动硬盘、磁碟、光盘、计算机存储器、只读存储器(ROM，Read‑OnlyMemory)、随机存取存储器(RAM，RandomAccessMemory)、电载波信号、电信信号以及软件分发介质等。需要说明的是，所述计算机可读介质包含的内容可以根据司法管辖区内立法和专利实践的要求进行适当的增减，例如在某些司法管辖区，根据立法和专利实践，计算机可读介质不包括电载波信号和电信信号。[0101] 尽管根据有限数量的实施例描述了本发明，但是受益于上面的描述，本技术领域内的技术人员明白，在由此描述的本发明的范围内，可以设想其它实施例。此外，应当注意，本说明书中使用的语言主要是为了可读性和教导的目的而选择的，而不是为了解释或者限定本发明的主题而选择的。因此，在不偏离所附权利要求书的范围和精神的情况下，对于本技术领域的普通技术人员来说许多修改和变更都是显而易见的。对于本发明的范围，对本发明所做的公开是说明性的，而非限制性的，本发明的范围由所附权利要求书限定。

专利地区：黑龙江

专利申请日期：2024-03-19

专利公开日期：2024-11-29

专利公告号：CN118169758B