一种基于原子训练的大地电磁信号噪声压制方法及系统
技术领域
本发明属于大地电磁信号处理技术领域,具体涉及一种基于原子训练的大地电磁信号噪声压制方法及系统。
背景技术
大地电磁测深法(Magnetotelluric,MT)是上个世纪50年代提出来的一种天然大地电磁场源勘探方法,由于天然的大地电磁信号极易受到外界条件的干扰,且天然大地电磁信号极其微弱以及频带范围宽,导致野外采集的数据因地形结构,人文电磁环境等因素,不可避免地受到了各种电磁噪声的干扰,正是因为这些干扰的存在,导致大地电磁数据质量严重堪忧,从而制约了地质结构的可解释性,由此可见大地电磁信号去噪技术愈发重要。其中,受干扰的信号主要可以划分为两大类:非线性信号、非平稳信号,其干扰从时间域波形上来说又可大致划分为三大类:类充放电三角波干扰、类方波干扰、类尖脉冲干扰。
从大地电磁法诞生至今,获取高质量的大地电磁信号方式众多,不断涌现了很多大地电磁噪声压制方法。如远参考法、robust法、时间域去噪法、频率域去噪法、正交匹配追踪算法等,这些信号处理方法的出现在一定程度上都可以压制大地电磁信号中的强干扰,从而提升大地电磁有用信号质量,但这些信号处理方法都存在着一定的局限性,往往丢失大量的有用信号,不能根据信号及干扰的特征灵活处理,从而导致去噪后的信号与实际信号存在着比较大的偏差,对后续的电磁法反演带来更大的困扰。
其中,传统的正交匹配追踪算法使用单一类型的原子,原子与待处理的大地电磁信号之间不存在任何的关联,从而导致在使用该类型原子对大地电磁信号进行去噪时,不能够达到令人满意的效果,因此,利用传统的正交匹配追踪算法进行大地电磁信号噪声压制的去噪效果还有待进一步提高。
发明内容
本发明的目的是针对利用传统的正交匹配追踪算法进行噪声压制所存在去噪精度不高的技术问题,提供一种基于原子训练的大地电磁信号噪声压制方法及系统,其使用原子训练的方式,自适应的根据大地电磁信号训练原子,构建大地电磁信号与原子之间的关联关系,再使用训练出的原子对大地电磁信号进行去噪,解决因为原子与大地电磁信号无关而导致去噪效果不理想的问题。
一方面,本发明提供一种基于原子训练的大地电磁信号噪声压制方法,包括如下步骤:
S1:将待去噪的大地电磁信号进行分段,并对每段大地电磁信号分别进行步进分割构建出每段大地电磁信号对应的训练矩阵;
S2:利用每段大地电磁信号对应的训练矩阵对预设的初始字典分别进行原子训练得到每段大地电磁信号的原子字典;
S3:获取每段大地电磁信号对应的待去噪矩阵,将步骤S2中得到原子字典作为正交匹配追踪算法中信号稀疏表示的去噪字典并对同一段大地电磁信号对应的待去噪矩阵进行去噪得到编码矩阵;
S4:分别计算同一段大地电磁信号对应原子字典与编码矩阵的内积得到去噪矩阵,再对去噪矩阵进行复原操作得到去噪的大地电磁信号。
本发明针对大地电磁信号干扰的不确定性和种类的复杂性选用原子训练和正交匹配追踪算法相结合的方式,一方面,将大地电磁信号分段处理,二方面,将每段大地电磁信号转换成一个训练矩阵,然后将该训练矩阵作为原子训练的训练样本,训练出来一个与该段大地电磁信号相关联的原子字典,进而构建了原子字典与该段大地电磁信号的关联关系,再将该原子字典作为正交匹配追踪算法的去噪字典,再将该段大地电磁信号转换成为一个待去噪矩阵,并利用正交匹配追踪算法对该待去噪矩阵进行去噪得到去噪矩阵,重构大地电磁有用信号。其中,本发明通过构建每段大地电磁信号自适应的原子字典,使得正交匹配追踪去噪过程中的去噪字典与待去噪的信号是紧密相关的,能够得到更加匹配的原子,从而有效地去除了大地电磁信号中的大尺度强干扰,解决了原子种类单一的传统正交匹配追踪算法所存在无法完全剔除强干扰,同时丢失大量有用信息的问题。
进一步优选,步骤S2中原子训练的方法为:将所述预设的初始字典作为正交匹配追踪算法中信号稀疏表示的去噪字典,并对同一段大地电磁信号对应的训练矩阵进行去噪得到训练阶段的去噪矩阵;再将基于训练阶段的去噪矩阵与训练矩阵的误差得到待更新原子对应的残差,再对残差进行奇异值分解完成以原子更新。
进一步优选,步骤S2中获取任意一段大地电磁信号对应的原子字典的过程均执行如下步骤:
S21:将初始字典作为正交匹配追踪算法中信号稀疏表示的去噪字典,并对大地电磁信号段对应的训练矩阵进行去噪得到训练阶段的编码矩阵;
S22:基于初始字典、训练阶段的编码矩阵以及训练矩阵计算出各个原子对应的残差,并对残差进行奇异值分解,再基于分解结果对原子进行更新以及对训练阶段的编码矩阵进行更新;
其中,残差以及残差的奇异值分解如下:
G
k表示当前原子d
k对应的残差,d
k、d
j为初始字典中第k个原子、第j个原子,H为训练矩阵,
为训练阶段的编码矩阵中的第j个行向量;
G
k=Λ∑V
T,Λ、∑、V均为残差的奇异值分解后得到酉矩阵、半正定对角矩阵、特征向量矩阵,其中,矩阵Λ中第一个列向量表示当前原子更新后的值d
k,矩阵∑的第一个值与矩阵V的第一列相乘作为当前训练阶段的编码矩阵更新后的第k个行向量
S23:基于更新后的初始字典以及编码矩阵再返回步骤S21进行下一次迭代更新,直至达到预设迭代次数。
进一步优选,同一段大地电磁信号对应的训练矩阵步进分割时的步进长度小于待去噪矩阵的步进长度。
训练矩阵的步进长度小于待去噪矩阵的步进长度,可以理解为训练矩阵中重合信息更多,由于大地电磁信号包含着强干扰信息,重合信息更多可以更大程度的获取到大地电磁信号中的细节信息,细节信息越多得到的原子字典与大地电磁信号之间的关系更加紧密,进而处理待去噪矩阵时,去噪效果更佳。
进一步优选,训练矩阵对应的步进长度为1,待去噪矩阵对应的步进长度为32。
通过实验与论证,上述参数的选取可以带来较佳的去噪效果。训练矩阵的步进长度为1是最大程度可以获取到更多的细节信息,而待去噪矩阵的步进长度为32,是因为在构造待去噪矩阵的时候,不需要获取大量的细节信息,只需将信号进行稀疏表示得到去噪矩阵即可。使用细节信息过多的矩阵进行稀疏表示时,往往会导致引入系统误差,因为在进行稀疏表示时不能确保每次都选到最为匹配的原子,所以使用重叠32个点的方式即可以将信息进行稀疏表示也可以减少信号系统误差。
进一步优选,所述初始字典的线性代数矩阵A如下:
式中,W等于训练矩阵的行数,所述初始字典是将矩阵A转置,再对转置矩阵AT中的每一列进行归一化处理得到。
由于大地电磁信号中存在着大量的低频信号,而低频信号在大地电磁信号领域属于有用信号,因此我们需要将低频信号保留下来,而本发明使用的初始字典中存在着大量的缓变化型原子,缓变化型原子对应初始字典线性代数矩阵中全为1的一行和部分周期较小函数(频率小于5hz的函数),缓变化原子的存在可以使得大地电磁信号中的低频信息得以保存下来,而使用本发明提供的初始字典作为原子训练的字典可以使得训练出来的分为两大类:一类是可以保留低频信号的原子、一类是可以与强干扰相匹配的原子,与强干扰相匹配的原子可以剔除强干扰信号,与低频信号相匹配的原子可以保留大地电磁信号中的低频有用信号。
本发明中初始字典的长度与训练矩阵的行数相同,初始字典中原子个数为原子长度的2倍,进而构建一个过完备字典,可以得到更好的去噪效果。应当理解,上述矩阵A的形式是本发明的一种具体形式,可以不限制于该形式,本发明的初始字典只需要满足存在部分缓变化型原子均可以实现上述提到的效果。
进一步优选,步骤S3中获取任一段大地电磁信号对应的编码矩阵的过程均执行如下步骤:
S31:获取待去噪矩阵,并将同一段大地电磁信号对应的原子字典作为正交匹配追踪算法的去噪字典;
S32:依次对待去噪矩阵中每列信号进行处理得到每一列信号对应的稀疏逼近信号:
首先,从待去噪矩阵匹配的原子字典中选择出与当前处理的一列信号相匹配的最佳原子;
其次,按照如下公式计算出当前处理的一列信号对应的稀疏逼近信号x:
x=arg min||f-dif||
式中,f为待去噪矩阵中当前处理的一列信号,di表示与当前处理的一列信号相匹配的最佳原子,arg min|| ||表示为范数最小值的复数辐角值;
然后,循环迭代处理待去噪矩阵中下一列信号直至计算出每一列信号对应的稀疏逼近信号;
S33:利用每一列信号对应的稀疏逼近信号构建编码矩阵。
本发明不单只针对每段大地电磁信号分别训练出一个自适应的原子字典,同时也针对每段大地电磁信号的待去噪矩阵中每一列信号均从原子字典中匹配一个最佳原子进行去噪,进一步提高了去噪效果;正是基于上述多重因素,本发明所使用的正交匹配追踪算法不需要进行多次稀疏表示也可满足去噪精度需求,从而更不需要人为确定信号的稀疏度,极大程度降低了人为因素的干扰,同时提高了算法的运行效率。
进一步优选,按照如下公式从待去噪矩阵匹配的原子字典中选择出与当前处理的一列信号相匹配的最佳原子:
式中,dj表示原子字典中的第j个原子,arg max||表示为最大值的复数辐角值,N表示为原子字典中原子的个数,i表示为最佳原子在原子字典中的位置。
进一步优选,所述待去噪矩阵是对大地电磁信号进行步进分割得到,步骤S4中对去噪矩阵进行复原操作的过程如下:
将去噪矩阵中各列重合部分的信号进行均值处理,再将均值处理后的信号进行拼接得到去噪的大地电磁信号。
另一方面,本发明提供一种基于上述方法的系统,包括:
构建模块:用于将待去噪的大地电磁信号进行分段,并对每段大地电磁信号分别进行步进分割构建每段大地电磁信号对应的训练矩阵;以及还用于获取每段大地电磁信号对应的待去噪矩阵;
原子字典训练模块:用于利用每段大地电磁信号对应的训练矩阵对预设的初始字典分别进行原子训练得到每段大地电磁信号的原子字典;
去噪模块:用于将得到原子字典作为正交匹配追踪算法的去噪字典并对同一段大地电磁信号对应的待去噪矩阵进行去噪得到编码矩阵;以及分别计算同一段大地电磁信号对应的编码矩阵与原子字典的内积得到去噪矩阵;
复原模块:用于对去噪矩阵进行复原操作得到去噪的大地电磁信号。
有益效果
1、本发明提供的一种基于原子训练的大地电磁信号噪声压制方法,通过原子训练的方式可以得到一个与大地电磁信号相关联的原子字典,训练后得到的原子字典中的原子与用于训练的大地电磁信号相关联,从而提高了字典中原子与大地电磁信号的匹配度,最终提高了噪声压制的效果。尤其是相较于传统的冗余字典使用的是固定原子类型的方式,传统方式中由于冗余字典中的原子与待匹配的大地电磁信号之间不存在实质性的关联,所以对于后续的实验和研究之间还存在着不可控的因素存在,从而影响研究的结果,本发明完全可以克服传统方式所存在的问题。
2、本发明进一步的优选方案中,本发明利用步进分割的手段构建训练矩阵以及待去噪矩阵,一方面,利用大地电磁信号构建的训练矩阵中重合信息远多于待去噪矩阵中的重合信息,使得训练矩阵中包含更多的大地电磁信号的细节信息,得到的原子字典与大地电磁信号之间的紧密性更高,去噪效果更优。另一方面,基于步进分割的手段,使得待去噪矩阵中包含重合信息,因此,本发明针对重合部分进行均值处理,而均值处理后再拼接信号可以有效避免传统信号拼接中存在的阶梯现象,这是因为去噪矩阵中每段需要拼接的信号之间存在部分信号的重合,因此,通过均值处理可以使得信号段之间更加的连续,从而克服阶梯现象,使得还原的大地电磁信号更加可靠。
3、本发明进一步的优选方案中,利用原子训练后的原子字典对同一段大地电磁信号进行去噪时,对待去噪矩阵中每一列信号分别从原子字典中查找出最佳原子,使得信号与原子之间的匹配度更高,因此,本发明在正交匹配追踪过程中仅仅只需要一次稀疏分解即可,不需要人为确定信号的稀疏度,极大程度降低了人为因素的干扰,同时提高了算法的运行效率。
4、本发明所利用的初始字典中包含了大量缓变化型原子,使用该初始字典作为原子训练训练出来的原子分为两大类,一类是可以保留低频信号的原子,一类是可以与强干扰相匹配的原子,与强干扰相匹配的原子可以剔除强干扰信号,与低频信号相匹配的原子可以保留大地电磁信号中的低频有用信号。
附图说明
图1为本发明实施例提供的流程图;
图2中图(a)模拟充放电三角波信号的噪声压制效果图,图(b)为该段模拟充放电三角波信号训练出来的原子;
图3中图(a)模拟大尺度方波信号的噪声压制效果图,图(b)为该段模拟大尺度方波信号训练出来的原子;
图4中图(a)模拟脉冲信号的噪声压制效果图,图(b)为该段模拟脉冲信号训练出来的原子;
图5实测大地电磁信号噪声压制效果图,其中,图(a)、图(b)和图(c)为噪声压制图;
图6为实测点EL22203A的视电阻率曲线和相位曲线,其中,图(a)为原始的电阻率曲线图和相位曲线图。发明方法处理后的视电阻率曲线和相位曲线其中图(b)为电阻率曲线图和相位曲线图,其中Rxy表示电道Ex和磁道Hy的视电阻率,Ryx表示电道Ey和磁道Hx的视电阻率,Pxy表示电道Ex和磁道Hy的相位,Pyx表示电道Ey和磁道Hx的相位。
具体实施方式
下面将结合实施例对本发明做进一步的说明。
如图1所示,本发明提供的一种基于原子训练的大地电磁信号噪声压制方法,包括如下步骤:
步骤1:将采集到的大地电磁信号进行去均值处理,并对去均值处理后的大地电磁信号进行等长度分段;
其中,大地电磁信号的数据量很大,所以我们需要对大地电磁信号进行等长度分段,来减少信号处理时的计算量,从而达到提高效率的目的,同时,大地电磁信号干扰复杂,本发明分段处理后,每段信号自适应训练出一个原子字典,相较于整体处理,分段得到原子字典与信号段的细节信息更加相关以及匹配,从而可以提高去噪精度。本发明实施例中进行等长度分段,其他可行的实施例中,可以选择不等长度,本发明对此不进行具体的限定。同时,本实施例此步骤进行了均值处理,其他可行的实施例中,可以不进行均值处理,本发明对此不进行具体的限定。
步骤2:对每个等长度分段过后的大地电磁信号段ai(t)分别进行构造训练矩阵T;
大地电磁信号段表示为ai(t),其中t=1,2,3,…,n(n>1000),ai(t)中的i表示为等距离分段后的大地电磁信号的第i段(i=1,2,3,…,m),m表示为信号段数,n为每段大地电磁信号段中采样点个数。本发明实施例中,选择步进分割长度为64,步进长度为1,即每隔1个点划分一次的分段处理,其中每两小段之间信号中重合的数目为63;其他可行的实施例中,可以选择其他参数。本实施例中大地电磁信号段ai(t)构建的训练矩阵T如下所示:
步骤3:获取一个离散的字典,作为原子训练的初始字典;
本发明中优选初始字典的长度与训练矩阵的行数相同,初始字典中原子个数为原子长度的2倍,进而构建一个过完备字典,可以得到更好的去噪效果。
本发明实施例中,训练矩阵T的行数为64,获取的初始字典的线性代数矩阵A如下所示:
本实施例中,W等于64,因此矩阵A为128×64,将其进行转换得到初始字典的大小为64×128。,初始字典的内容如下所示:
步骤4:利用每段大地电磁信号对应的训练矩阵对预设的初始字典分别进行原子训练得到每段大地电磁信号的原子字典。
应当理解,不同段大地电磁信号均利用初始字典进行原子训练,训练后,每一段大地电磁信号得到一个原子字典。
本实施例中,初始字典中原子的更新方式如下:
步骤4.1:初始化操作:将大地电磁信号段的训练矩阵T用矩阵H表示,初始字典用D表示;
步骤4.2:使用初始字典D对大地电磁信号段对应的训练矩阵进行正交匹配追踪,得到训练阶段的编码矩阵B;
步骤4.3:将初始字典D中的原子(列)与稀疏矩阵B中的行相乘得到训练阶段的去噪矩阵DB;
步骤4.4:分别更新字典中原子,即基于步骤4.3中去噪矩阵与训练矩阵的误差计算出各个原子对应的残差;
误差计算的规则如下:
式中,
为Frobenius范数的平方,d
j为字典D中的第j个原子,残差
为稀疏矩阵B中的第j个行向量,d
k为字典D中的第k个原子,
稀疏矩阵B中的第k个行向量,K表示为字典中原子的个数,k表示为字典中要更新的原子位置;
步骤4.5:对得到的各个原子对应的残差进行奇异值分解,譬如,针对原子d
k,对其残差进行奇异值分解:G
k=Λ∑V
T,取出矩阵Λ的第一个列向量来表示d
k用于更新,取Σ矩阵的第一个值与V的第一列相乘来表示
用于更新,从而更新完了其中的一个原子,按照相同方式更新下一个原子,直至所有原子更新完毕;
其中,奇异值分解Gk=ΛΣVT的规则如下:
将矩阵G
k与其的转置相乘
得到一个方阵,求取方阵的特征值(线性代数中求解矩阵的特征值)
式中λ为
的特征值集合,V为特征向量,对特征值集合λ进行如下操作得到矩阵Λ:
特征值集合λ中的特征值按照从小到大的循序进行排列得到半正定对角矩阵∑;
式中,Λ为一个酉矩阵,酉矩阵中的向量时正交的,∑为半正定对角矩阵,除了对角线上的元素,其余全是0,VT为V的共轭转置是一个酉矩阵。
步骤4.6,字典完整的更新一次完成字点变成Dnew,跳转至步骤4.2将正交匹配追踪算法的字典进行更新D=Dnew开始下一次的循环,直至循环预设次数后完成原子训练,本实施例中,循环128次。
其中,图2、图3、图4中的(b)图分别表示为训练出来的三角波原子、方波原子、脉冲原子。
步骤5:利用每段大地电磁信号分别构建每段大地电磁信号对应的待去噪矩阵Q,将步骤4训练出来的字典作为正交匹配追踪算法的去噪字典,并对同一段大地电磁信号对应的待去噪矩阵Q进行正交匹配追踪处理得到编码矩阵U。
其中,用于构建待去噪矩阵Q的大地电磁信号与用于构建训练矩阵T的大地电磁信号必须是同一段。本实施例中,构建待去噪矩阵Q是步进分割长度为64,步进长度为32,即每隔32个点划分一次的分段处理,其中每两小段信号之间中重合的点数为32,将大地电磁信号段转变成为一个行数为64的待去噪矩阵Q。如图所示,大地电磁信号段ai(t)构建的待去噪矩阵Q如下所示:
针对任意一段大地电磁信号获取编码矩阵的过程如下:
步骤5.1:利用步骤4训练出来的原子字典作为正交匹配追踪算法的去噪字典;
步骤5.2:获取待去噪矩阵Q的列数,矩阵的列数即为循环迭代的次数;
步骤5.3:从原子字典中选择一个与待去噪矩阵Q中当前处理的一列信号相匹配的最佳原子;本实施例中,选取规则如下:
式中,dj表示原子字典中的第j个原子,arg max||表示为最大值的复数辐角值,N表示为原子字典中原子的个数,i表示为最佳原子在原子字典中的位置。
步骤5.4:求解大地电磁信号的稀疏逼近信号x,稀疏逼近信号的得到通过如下公式所得:
x=arg min||f-dif||
式中,f为待去噪矩阵中Q的当前处理的一列信号,di表示与当前处理的一列信号相匹配的最佳原子,arg min|| ||表示为范数最小值的复数辐角值;
步骤5.5:将得到的稀疏逼近信号储存到编码矩阵U中,并处理待去噪矩阵Q中的下一列信号,转跳到步骤5.3,直至计算完待去噪矩阵中的所有列得到编码矩阵U。
步骤6:计算步骤4训练出来的原子字典与步骤5得到的编码矩阵U的內积得到去噪矩阵S,并对去噪矩阵进行复原操作得到去噪信号。
本实施例中,去噪矩阵S中的内容可以表示为如下形式:
式中,M表示为去噪矩阵的列数。
由于去噪矩阵中存在重合部分,本实施例中对重合部分采取均值处理,譬如:矩阵中的第l列和第l+1列表示如下所示:
则,信号的重叠的部分进行如下操作:
均值处理后,再对信号进行拼接得到大地电磁信号。
其中,如图2、图3、图4中的(a)图分别表示模拟三角波干扰、模拟方波干扰,模拟脉冲干扰的电磁信号去噪的效果,如图2、图3、图4中的(b)图分别表示模拟三角波干扰、模拟方波干扰,模拟脉冲干扰训练出来的原子。如图5所示,(a)、(b)、(c)图为表示实测干扰的电磁信号去噪的效果图。
对比实测点(EL22203A)经本发明方法处理前后的视电阻率曲线来评价效果优劣。图6所示的(a)图为未经处理的大地电磁信号的视电阻率曲线和相位曲线,图6所示的(b)图为经过本发明方法处理后的视电阻率曲线和相位曲线,通过对比处理前和处理后的视电阻率曲线和相位曲线,本发明的方法有效地剔除了大地电磁信号中的强干扰信号,经本发明处理后的视电阻率曲线和相位更加平缓,有效地压制了近源干扰,对于大地电磁信号数据处理有着广阔的应用前景。
此外,基于上述方法,本发明还提供一种基于上述方法的系统,包括:
构建模块:用于将待去噪的大地电磁信号进行分段,并对每段大地电磁信号分别进行步进分割构建每段大地电磁信号对应的训练矩阵;以及还用于获取每段大地电磁信号对应的待去噪矩阵;
原子字典训练模块:用于利用每段大地电磁信号对应的训练矩阵对预设的初始字典分别进行原子训练得到每段大地电磁信号的原子字典;
去噪模块:用于将得到原子字典作为正交匹配追踪算法的去噪字典并对同一段大地电磁信号对应的待去噪矩阵进行去噪得到编码矩阵;以及分别计算同一段大地电磁信号对应的编码矩阵与原子字典的内积得到去噪矩阵;
复原模块:用于对去噪矩阵进行复原操作得到去噪的大地电磁信号。
应当理解:本发明各个实施例中的功能单元模块可以集中在一个处理单元中,也可以是各个单元模块单独物理存在,也可以是两个或两个以上的单元模块集成在一个单元模块中,可以采用硬件或软件的形式来实现。
需要强调的是,本发明所述的实例是说明性的,而不是限定性的,因此本发明不限于具体实施方式中所述的实例,凡是由本领域技术人员根据本发明的技术方案得出的其他实施方式,不脱离本发明宗旨和范围的,不论是修改还是替换,同样属于本发明的保护范围。