CN107607065A - 一种基于变分模态分解的冲击回波信号分析方法 - Google Patents
一种基于变分模态分解的冲击回波信号分析方法 Download PDFInfo
- Publication number
- CN107607065A CN107607065A CN201710864339.XA CN201710864339A CN107607065A CN 107607065 A CN107607065 A CN 107607065A CN 201710864339 A CN201710864339 A CN 201710864339A CN 107607065 A CN107607065 A CN 107607065A
- Authority
- CN
- China
- Prior art keywords
- hilbert
- echo signal
- time
- spectrum
- 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.)
- Pending
Links
Abstract
本发明公开了一种基于变分模态分解的冲击回波信号分析方法，包括1）根据现场测量环境要求，选择测量参数和设备工作参数，采集冲击回波信号；2）根据冲击回波信号，设定变分模态分解中的参数；3）利用变分模态分解方法将采集的冲击回波信号分解为若干个本征模量函数，求解变分问题；4）利用Hilbert变换，得出本征模量函数的Hilbert时频谱；5）将不同频率的Hilbert时频谱在时间域进行积分得到最终的边际谱。本发明通过冲击回波信号的边际谱识别混凝土构件的缺陷，在强噪声影响下能很好的压制噪声干扰，相对于传统Fourier变换具有更高的分辨率。
Description
技术领域
本发明涉及一种基于变分模态分解的冲击回波信号分析方法，属于土木结构工程检测技术领域。
背景技术
目前，采用冲击回波法对混凝土构件厚度与缺陷检测有了一定的研究，一般采用Fourier变换求解出振动信号的频谱特征，再根据主频和表观速度计算出混凝土构件的厚度。但实际的冲击回波信号并非平稳信号，却是非平稳的信号，在实际工程中冲击回波不但从底板反射而且从缺陷处发生反射，使得Fourier 变换频谱常出现多重峰值，因而采用传统的Fourier变换很难正确识别。Hilber-Huang变换将非平稳信号分解为有限个本征模态函数，然后对本征模态函数进行Hilbert变换，再构造出信号的Hilbert振幅谱，最后通过对Hilbert谱积分得出边际谱，边际谱相对于Fourier频谱更加准确反映信号的实际频谱，但是 Hilber-Huang变换经验模态分解法常出现模态混叠和端点效应等问题。
发明内容
本发明所要解决的技术问题是克服现有技术的缺陷，提供一种基于变分模态分解的冲击回波信号分析方法，根据变分模态分解和Hibert变换方法原理，提出通过冲击回波信号的边际谱识别混凝土构件的缺陷的检测方法，该方法在强噪声影响下能很好的压制噪声干扰，相对于传统Fourier变换具有更高的分辨率。
为解决上述技术问题，本发明提供一种基于变分模态分解的冲击回波信号分析方法，包括以下步骤：
1)根据现场测量环境要求，选择测量参数和设备工作参数，采集冲击回波信号；
2)根据冲击回波信号，设定变分模态分解中的参数；
3)利用变分模态分解方法将采集的冲击回波信号分解为n个不同尺度的本征模量函数；
4)利用Hilbert变换，得出本征模量函数的Hilbert时频谱；
5)将不同频率的Hilbert时频谱在时间域进行积分得到最终的边际谱。
前述的步骤1)中，测量参数和设备工作参数包括采样频率，采样点及震源。
前述的步骤2)中，变分模态分解中的参数包括拉格朗日乘子及中心频率。
前述的步骤3)中，变分模态分解包括小波变换、短时傅立叶变换和经验模态分解。
前述的步骤3)的分解过程具体如下：
3-1)假设将采集的冲击回波信号分解为n个本征模量函数：
其中，f(t)为冲击回波信号，uk(t)表示第k个本征模量函数，t表示t时刻；
3-2)对每个本征模量函数进行Hilbert变换得到解析信号：
其中，δ(t)是狄利克来函数，j是虚数符号；
3-3)通过移频方式，将各解析信号的频谱变换到基带上：
其中，ωk表示中心频率；
3-4)计算解析信号的欧式距离，计算各本征模量函数的带宽，变分问题如下：
3-5)利用二次罚函数项和Lagrange乘子将变分问题转化为无约束问题：
其中，α为惩罚因子，λ(t)为Lagrange乘子，表示范数，<·>表示内积运算；
3-6)利用乘子交替方向算法求取方程(5)的无约束变分问题，最终得到多个不同尺度的本征模量函数。
前述的步骤4)中，本征模量函数的Hilbert时频谱表示为：
其中，P为柯西主值，τ为时间变量。
前述的步骤5)，求解边际谱，具体如下：
5-1)由本征模量函数的Hilbert时频谱和本征模量函数构造出解析信号，形式为：
A[uk(t)]＝uk(t)+iH[uk(t)] (7)
其中，A[uk(t)]表示解析信号；
5-2)求解解析信号的振幅和相位：
其中，ak(t)表示振幅，θk(t)表示相位；
5-3)对相位进行求导，计算得出解析信号的瞬时频率：
其中，ωk(t)表示瞬时频率。
5-4)将振幅ak(t)表示在联合的时频平面上，得出本征模量函数uk(t)的Hilbert时频谱：
其中，Hk(t,ω)表示Hilbert时频谱，ω为频率；
5-5)对冲击回波信号f(t)进行整体Hilbert谱分析，Hilbert时频谱表示为：
5-6)根据式(12)的Hilbert时频谱进一步定义出边际谱：
其中，h(ω)表示边际谱。
本发明所达到的有益效果：
1)本发明针对非稳定频率的冲击回波信号，采用边际谱识别混凝土构件的缺陷，在强噪声影响下能很好的压制噪声干扰，相对于传统Fourier变换具有更高的分辨率。
2)本发明首次提出用变分模态分解法对冲击回波信号进行分析，该方法相对现有的多尺度分析方法有更完善的数学基础及实际效果，该方法采用温纳滤波理论不仅能分解出不同尺度多模态，而且能压制噪声。
3)本发明并不仅限于对冲击回波信号的特征分析，也可以用于其它振动检测信号分析，例如桩基低应变检测信号的分析。
附图说明
图1为本发明方法的流程图；
图2为冲击回波检测混凝土底板几何模型；
图3为底板模拟信号，(a)为原数值模拟信号，(b)为加入白噪声的强干扰信号；
图4为对强干扰信号处理后的频谱图，(a)为Fourier变换频谱，(b)为VMD边际频谱。
具体实施方式
下面对本发明作进一步描述。以下实施例仅用于更加清楚地说明本发明的技术方案，而不能以此来限制本发明的保护范围。
本发明的基于变分模态分解的冲击回波信号分析方法，如图1所示，包括以下步骤：
1、根据现场测量环境要求，选择合适的测量参数和设备工作参数，如设置采样频率，采样点及震源，采集冲击回波信号。
2、根据冲击回波信号，设定变分模态分解中的参数；参数包括拉格朗日乘子及中心频率。
3、利用变分模态分解方法将采集的冲击回波信号分解为若干个不同尺度的本征模量函数，求解变分问题，具体如下：
3-1)假设将冲击回波信号分解为n个本征模量函数：
其中，f(t)为冲击回波信号，uk(t)表示第k个本征模量函数，t表示t时刻；
变分模态分解，可以是传统的信号分析方法，包括小波变换、短时傅立叶变换和经验模态分解等方法。
3-2)对每个本征模量函数进行Hilbert变换得到其解析信号：
其中，δ(t)是狄利克来函数，j是虚数符号；
3-3)通过移频方式，将各解析信号的频谱变换到基带上：
其中，ωk表示中心频率；
3-4)计算上述解析信号的欧式距离，估计各模态带宽，其变分问题如下：
3-5)利用二次罚函数项和Lagrange乘子将其转化为无约束问题：
式中，α为惩罚因子，λ(t)为Lagrange乘子，表示范数，<·>表示内积运算；
3-6)再利用乘子交替方向算法求取方程(5)无约束变分问题，最终得到多个不同尺度的本征模量函数。
4、利用Hilbert变换，得出本征模量函数的Hilbert时频谱；变换形式表示为：
其中，P为柯西主值，τ为时间变量。
5、将不同频率的Hilbert时频谱在时间域进行积分得到最终的边际谱，具体过程如下：
5-1)由本征模量函数的Hilbert时频谱式(6)和本征模量函数构造出解析信号，其形式为：
A[uk(t)]＝uk(t)+iH[uk(t)] (7)
A[uk(t)]表示解析信号。
5-2)对公式(7)的解析信号，求出其振幅和相位：
其中，ak(t)表示振幅，θk(t)表示相位。
5-3)对相位进行求导，计算得出解析信号的瞬时频率：
其中，ωk(t)表示瞬时频率。
5-4)将振幅ak(t)表示在联合的时频平面上，得出本征模量函数uk(t)的Hilbert时频谱：
其中，Hk(t,ω)表示Hilbert时频谱，ω为频率。
5-5)对f(t)进行整体Hilbert谱分析，Hilbert时频谱可表示为：
5-6)根据式(12)的Hilbert时频谱进一步定义出边际谱：
其中，h(ω)表示边际谱。
实施例
冲击回波是基于瞬态P波进行的结构无损检测，构件的厚度按下列公式计算：
其中，T为混凝土结构构件的厚度，Vp为P波的波速，ωc为振幅谱主频。
P波波速Vp与结构的杨氏模量E，泊松比ν及结构密度ρ存在下列关系：
以一个混凝土底板作为本次试验对象，混凝土底板尺寸1.50m×1.50m，厚度为0.20m。杨氏模量E，泊松比ν及结构密度ρ分别为37.8GPa，2500kg/m3，0.2。根据公式(15)可得出纵波波速4099m/s。冲击回波厚度根据公式(14)计算，Fourier主频理论值应为10247Hz。
对于冲击回波表面响应加速度值，根据模型参数通过有限元软件进行仿真。试验中将混凝土底板中心作为冲击球捶击点，冲击点与加速度接收点如图2所示。图中，带箭头为冲击点，另一点为加速度接收点。冲击球的冲击点与加速度检波器之间距离为0.06m，加速度采样频率取1MHz，采样长度取2ms。冲击球荷载的由下列函数确定：
P(t)＝(sin(πt/40))3 0≤t≤40μs (16)
通过上述参数设置，模拟出底板表面接收点的加速度随时间的变化情况，如图3(a)所示，取初始信号的信噪比SNR＝-8，将原数值模拟信号图3(a)加入白噪声，得到加入强干扰后的结果(图3(b))。
将含有强干扰的信号分别利用Fourier和VMD变换，求解出Fourier频谱和边际谱，变换结果如图3 所示，图4(a)为Fourier变换频谱，图4(b)为VMD边际频谱。
根据图4(a)分析结果，在强白噪声背景下采用Fourier变换得到的频谱图中显示多个振动主峰，呈锯齿状，不包含单一的主频率峰值，主频率特征不明显。而在图4(b)采用VMD方法得出的边际谱，背景噪声对主频峰值的识别影响较小，即使在强白噪声背景下主频峰值也很明显。因此，VMD方法得出的边际谱能很好的压制干扰的影响，相对传统的Fourier频谱具有显著优势。
以上所述仅是本发明的优选实施方式，应当指出，对于本技术领域的普通技术人员来说，在不脱离本发明技术原理的前提下，还可以做出若干改进和变形，这些改进和变形也应视为本发明的保护范围。
Claims (7)
1.一种基于变分模态分解的冲击回波信号分析方法，其特征在于，包括以下步骤：
1)根据现场测量环境要求，选择测量参数和设备工作参数，采集冲击回波信号；
2)根据冲击回波信号，设定变分模态分解中的参数；
3)利用变分模态分解方法将采集的冲击回波信号分解为n个不同尺度的本征模量函数；
4)利用Hilbert变换，得出本征模量函数的Hilbert时频谱；
5)将不同频率的Hilbert时频谱在时间域进行积分得到最终的边际谱。
2.根据权利要求1所述的一种基于变分模态分解的冲击回波信号分析方法，其特征在于，所述步骤1)中，测量参数和设备工作参数包括采样频率，采样点及震源。
3.根据权利要求1所述的一种基于变分模态分解的冲击回波信号分析方法，其特征在于，所述步骤2)中，变分模态分解中的参数包括拉格朗日乘子及中心频率。
4.根据权利要求1所述的一种基于变分模态分解的冲击回波信号分析方法，其特征在于，所述步骤3)中，变分模态分解包括小波变换、短时傅立叶变换和经验模态分解。
5.根据权利要求1所述的一种基于变分模态分解的冲击回波信号分析方法，其特征在于，所述步骤3)的分解过程具体如下：
3-1)假设将采集的冲击回波信号分解为n个本征模量函数：
<mrow>
<mi>f</mi>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<munderover>
<mo>&Sigma;</mo>
<mrow>
<mi>k</mi>
<mo>=</mo>
<mn>1</mn>
</mrow>
<mi>n</mi>
</munderover>
<msub>
<mi>u</mi>
<mi>k</mi>
</msub>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>1</mn>
<mo>)</mo>
</mrow>
</mrow>
其中，f(t)为冲击回波信号，uk(t)表示第k个本征模量函数，t表示t时刻；
3-2)对每个本征模量函数进行Hilbert变换得到解析信号：
<mrow>
<mo>(</mo>
<mi>&delta;</mi>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
<mo>+</mo>
<mfrac>
<mi>j</mi>
<mrow>
<mi>&pi;</mi>
<mi>t</mi>
</mrow>
</mfrac>
<mo>)</mo>
<mo>*</mo>
<msub>
<mi>u</mi>
<mi>k</mi>
</msub>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>2</mn>
<mo>)</mo>
</mrow>
</mrow>
其中，δ(t)是狄利克来函数，j是虚数符号；
3-3)通过移频方式，将各解析信号的频谱变换到基带上：
<mrow>
<mo>&lsqb;</mo>
<mrow>
<mo>(</mo>
<mi>&delta;</mi>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
<mo>+</mo>
<mfrac>
<mi>j</mi>
<mrow>
<mi>&pi;</mi>
<mi>t</mi>
</mrow>
</mfrac>
<mo>)</mo>
</mrow>
<mo>*</mo>
<msub>
<mi>u</mi>
<mi>k</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
<mo>&rsqb;</mo>
<msup>
<mi>e</mi>
<mrow>
<mo>-</mo>
<msub>
<mi>j&omega;</mi>
<mi>k</mi>
</msub>
<mi>t</mi>
</mrow>
</msup>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>3</mn>
<mo>)</mo>
</mrow>
</mrow>
其中，ωk表示中心频率；
3-4)计算解析信号的欧式距离，计算各本征模量函数的带宽，变分问题如下：
<mrow>
<mfenced open = "{" close = "">
<mtable>
<mtr>
<mtd>
<mrow>
<munder>
<mi>min</mi>
<mrow>
<mo>{</mo>
<msub>
<mi>u</mi>
<mi>k</mi>
</msub>
<mo>}</mo>
<mo>,</mo>
<mo>{</mo>
<msub>
<mi>&omega;</mi>
<mi>k</mi>
</msub>
<mo>}</mo>
</mrow>
</munder>
<mo>{</mo>
<munder>
<mo>&Sigma;</mo>
<mi>k</mi>
</munder>
<mo>|</mo>
<mo>|</mo>
<msub>
<mo>&part;</mo>
<mi>t</mi>
</msub>
<mo>&lsqb;</mo>
<mrow>
<mo>(</mo>
<mi>&delta;</mi>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
<mo>+</mo>
<mfrac>
<mi>j</mi>
<mrow>
<mi>&pi;</mi>
<mi>t</mi>
</mrow>
</mfrac>
<mo>)</mo>
</mrow>
<mo>*</mo>
<msub>
<mi>u</mi>
<mi>k</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
<mo>&rsqb;</mo>
<msup>
<mi>e</mi>
<mrow>
<mo>-</mo>
<msub>
<mi>j&omega;</mi>
<mi>k</mi>
</msub>
<mi>t</mi>
</mrow>
</msup>
<mo>|</mo>
<msup>
<mo>|</mo>
<mn>2</mn>
</msup>
<mo>}</mo>
</mrow>
</mtd>
</mtr>
<mtr>
<mtd>
<mtable>
<mtr>
<mtd>
<mrow>
<mi>s</mi>
<mo>.</mo>
<mi>t</mi>
<mo>.</mo>
</mrow>
</mtd>
<mtd>
<mrow>
<munder>
<mo>&Sigma;</mo>
<mi>k</mi>
</munder>
<msub>
<mi>u</mi>
<mi>k</mi>
</msub>
<mo>=</mo>
<mi>f</mi>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
</mrow>
</mtd>
</mtr>
</mtable>
</mtd>
</mtr>
</mtable>
</mfenced>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>4</mn>
<mo>)</mo>
</mrow>
</mrow>
3-5)利用二次罚函数项和Lagrange乘子将变分问题转化为无约束问题：
<mrow>
<mtable>
<mtr>
<mtd>
<mrow>
<mi>L</mi>
<mrow>
<mo>(</mo>
<mo>{</mo>
<msub>
<mi>u</mi>
<mi>k</mi>
</msub>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
<mo>}</mo>
<mo>,</mo>
<mo>{</mo>
<msub>
<mi>&omega;</mi>
<mi>k</mi>
</msub>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
<mo>}</mo>
<mo>,</mo>
<mi>&lambda;</mi>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
<mo>)</mo>
</mrow>
<mo>=</mo>
<mi>&alpha;</mi>
<munder>
<mo>&Sigma;</mo>
<mi>k</mi>
</munder>
<mo>|</mo>
<mo>|</mo>
<msub>
<mo>&part;</mo>
<mi>t</mi>
</msub>
<mo>&lsqb;</mo>
<mrow>
<mo>(</mo>
<mi>&delta;</mi>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
<mo>+</mo>
<mfrac>
<mi>j</mi>
<mrow>
<mi>&pi;</mi>
<mi>t</mi>
</mrow>
</mfrac>
<mo>)</mo>
</mrow>
<mo>*</mo>
<msub>
<mi>u</mi>
<mi>k</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
<mo>&rsqb;</mo>
<msup>
<mi>e</mi>
<mrow>
<mo>-</mo>
<msub>
<mi>j&omega;</mi>
<mi>k</mi>
</msub>
<mi>t</mi>
</mrow>
</msup>
<mo>|</mo>
<msubsup>
<mo>|</mo>
<mn>2</mn>
<mn>2</mn>
</msubsup>
<mo>+</mo>
</mrow>
</mtd>
</mtr>
<mtr>
<mtd>
<mrow>
<mo>|</mo>
<mo>|</mo>
<mi>f</mi>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
<mo>-</mo>
<munder>
<mo>&Sigma;</mo>
<mi>k</mi>
</munder>
<msub>
<mi>u</mi>
<mi>k</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
<mo>|</mo>
<msubsup>
<mo>|</mo>
<mn>2</mn>
<mn>2</mn>
</msubsup>
<mo>+</mo>
<mo><</mo>
<mi>&lambda;</mi>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
<mo>,</mo>
<mi>f</mi>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
<mo>-</mo>
<munder>
<mo>&Sigma;</mo>
<mi>k</mi>
</munder>
<msub>
<mi>u</mi>
<mi>k</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
<mo>></mo>
</mrow>
</mtd>
</mtr>
</mtable>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>5</mn>
<mo>)</mo>
</mrow>
</mrow>
其中，α为惩罚因子，λ(t)为Lagrange乘子，表示范数，<·>表示内积运算；
3-6)利用乘子交替方向算法求取方程(5)的无约束变分问题，最终得到多个不同尺度的本征模量函数。
6.根据权利要求5所述的一种基于变分模态分解的冲击回波信号分析方法，其特征在于，所述步骤4)中，本征模量函数的Hilbert时频谱表示为：
<mrow>
<mi>H</mi>
<mo>&lsqb;</mo>
<msub>
<mi>u</mi>
<mi>k</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
<mo>&rsqb;</mo>
<mo>=</mo>
<mfrac>
<mn>1</mn>
<mi>&pi;</mi>
</mfrac>
<mi>P</mi>
<msubsup>
<mo>&Integral;</mo>
<mrow>
<mo>-</mo>
<mi>&infin;</mi>
</mrow>
<mrow>
<mo>+</mo>
<mi>&infin;</mi>
</mrow>
</msubsup>
<mfrac>
<mrow>
<msub>
<mi>u</mi>
<mi>k</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>&tau;</mi>
<mo>)</mo>
</mrow>
</mrow>
<mrow>
<mi>t</mi>
<mo>-</mo>
<mi>&tau;</mi>
</mrow>
</mfrac>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>6</mn>
<mo>)</mo>
</mrow>
</mrow>
其中，P为柯西主值，τ为时间变量。
7.根据权利要求6所述的一种基于变分模态分解的冲击回波信号分析方法，其特征在于，所述步骤5)，求解边际谱，具体如下：
5-1)由本征模量函数的Hilbert时频谱和本征模量函数构造出解析信号，形式为：
A[uk(t)]＝uk(t)+iH[uk(t)] (7)
其中，A[uk(t)]表示解析信号；
5-2)求解解析信号的振幅和相位：
<mrow>
<msub>
<mi>a</mi>
<mi>k</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<msqrt>
<mrow>
<msubsup>
<mi>u</mi>
<mi>k</mi>
<mn>2</mn>
</msubsup>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
<mo>+</mo>
<msup>
<mi>H</mi>
<mn>2</mn>
</msup>
<mo>&lsqb;</mo>
<msub>
<mi>u</mi>
<mi>k</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
<mo>&rsqb;</mo>
</mrow>
</msqrt>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>8</mn>
<mo>)</mo>
</mrow>
</mrow>
<mrow>
<msub>
<mi>&theta;</mi>
<mi>k</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<mi>a</mi>
<mi>r</mi>
<mi>c</mi>
<mi>t</mi>
<mi>a</mi>
<mi>n</mi>
<mfrac>
<mrow>
<mi>H</mi>
<mo>&lsqb;</mo>
<msub>
<mi>u</mi>
<mi>k</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
<mo>&rsqb;</mo>
</mrow>
<mrow>
<msub>
<mi>u</mi>
<mi>k</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
</mrow>
</mfrac>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>9</mn>
<mo>)</mo>
</mrow>
</mrow>
其中，ak(t)表示振幅，θk(t)表示相位；
5-3)对相位进行求导，计算得出解析信号的瞬时频率：
<mrow>
<msub>
<mi>&omega;</mi>
<mi>k</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<mfrac>
<mn>1</mn>
<mn>2</mn>
</mfrac>
<mfrac>
<mrow>
<msub>
<mi>d&theta;</mi>
<mi>k</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
</mrow>
<mrow>
<mi>d</mi>
<mi>t</mi>
</mrow>
</mfrac>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>10</mn>
<mo>)</mo>
</mrow>
</mrow>
其中，ωk(t)表示瞬时频率。
5-4)将振幅ak(t)表示在联合的时频平面上，得出本征模量函数uk(t)的Hilbert时频谱：
<mrow>
<msub>
<mi>H</mi>
<mi>k</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>,</mo>
<mi>&omega;</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<mfenced open = "{" close = "">
<mtable>
<mtr>
<mtd>
<mrow>
<msub>
<mi>a</mi>
<mi>k</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
<mo>,</mo>
<mi>&omega;</mi>
<mo>=</mo>
<msub>
<mi>&omega;</mi>
<mi>k</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
</mrow>
</mtd>
</mtr>
<mtr>
<mtd>
<mrow>
<mn>0</mn>
<mo>,</mo>
<mi>&omega;</mi>
<mo>&NotEqual;</mo>
<msub>
<mi>&omega;</mi>
<mi>k</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
</mrow>
</mtd>
</mtr>
</mtable>
</mfenced>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>11</mn>
<mo>)</mo>
</mrow>
</mrow>
其中，Hk(t,ω)表示Hilbert时频谱，ω为频率；
5-5)对冲击回波信号f(t)进行整体Hilbert谱分析，Hilbert时频谱表示为：
<mrow>
<mi>H</mi>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>,</mo>
<mi>&omega;</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<munderover>
<mo>&Sigma;</mo>
<mrow>
<mi>k</mi>
<mo>=</mo>
<mn>1</mn>
</mrow>
<mi>n</mi>
</munderover>
<msub>
<mi>H</mi>
<mi>k</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>,</mo>
<mi>&omega;</mi>
<mo>)</mo>
</mrow>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>12</mn>
<mo>)</mo>
</mrow>
</mrow>
5-6)根据式(12)的Hilbert时频谱进一步定义出边际谱：
<mrow>
<mi>h</mi>
<mrow>
<mo>(</mo>
<mi>&omega;</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<msubsup>
<mo>&Integral;</mo>
<mrow>
<mo>-</mo>
<mi>&infin;</mi>
</mrow>
<mrow>
<mo>+</mo>
<mi>&infin;</mi>
</mrow>
</msubsup>
<mi>H</mi>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>,</mo>
<mi>&omega;</mi>
<mo>)</mo>
</mrow>
<mi>d</mi>
<mi>t</mi>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>13</mn>
<mo>)</mo>
</mrow>
</mrow>
其中，h(ω)表示边际谱。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201710864339.XA CN107607065A (zh) | 2017-09-22 | 2017-09-22 | 一种基于变分模态分解的冲击回波信号分析方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201710864339.XA CN107607065A (zh) | 2017-09-22 | 2017-09-22 | 一种基于变分模态分解的冲击回波信号分析方法 |
Publications (1)
| Publication Number | Publication Date |
|---|---|
| CN107607065A true CN107607065A (zh) | 2018-01-19 |
Family
ID=61061682
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201710864339.XA Pending CN107607065A (zh) | 2017-09-22 | 2017-09-22 | 一种基于变分模态分解的冲击回波信号分析方法 |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN107607065A (zh) |
Cited By (12)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN108596215A (zh) * | 2018-04-04 | 2018-09-28 | 中国航发湖南动力机械研究所 | 多模态信号解析分离方法、装置、设备及存储介质 |
| CN108875170A (zh) * | 2018-06-05 | 2018-11-23 | 天津大学 | 一种基于改进变分模态分解的噪声源识别方法 |
| CN108919203A (zh) * | 2018-05-21 | 2018-11-30 | 西安电子科技大学 | 一种雷达有源干扰的识别方法及系统 |
| CN109323832A (zh) * | 2018-09-12 | 2019-02-12 | 温州大学 | 一种冷镦成型机模具冲击状态的监测方法 |
| CN109359633A (zh) * | 2018-12-10 | 2019-02-19 | 西北工业大学 | 基于希尔伯特-黄变换和小波脊线的信号联合分类方法 |
| CN109948249A (zh) * | 2019-03-19 | 2019-06-28 | 东南大学 | 基于变分模态分解的斜拉索时变索力识别方法 |
| CN110082436A (zh) * | 2019-04-25 | 2019-08-02 | 电子科技大学 | 一种基于变分模态的高提离电磁超声信号消噪方法 |
| CN110319995A (zh) * | 2019-08-14 | 2019-10-11 | 清华大学 | 火工冲击响应数据时频谱图分析方法 |
| CN111401168A (zh) * | 2020-03-06 | 2020-07-10 | 上海神添实业有限公司 | 一种无人机的多层雷达特征提取与选择方法 |
| CN112924655A (zh) * | 2021-01-28 | 2021-06-08 | 北京中水科海利工程技术有限公司 | 一种渡槽高排架结构混凝土质量的检测方法 |
| CN115436486A (zh) * | 2021-11-12 | 2022-12-06 | 江西理工大学 | 一种基于Hilbert变换的岩体破裂声发射信号伪维格纳分布二次时频分析方法 |
| CN116013240A (zh) * | 2023-01-07 | 2023-04-25 | 广西大学 | 基于变分模态分解与数字滤波的钢管混凝土信号降噪方法 |
Citations (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN103034757A (zh) * | 2012-12-02 | 2013-04-10 | 中国科学院电工研究所 | 基于经验模态分解的风电场时频域建模方法 |
| CN105160077A (zh) * | 2015-08-07 | 2015-12-16 | 中国电力科学研究院 | 一种高压直流线路电晕电流变化趋势提取方法 |
| CN105758644A (zh) * | 2016-05-16 | 2016-07-13 | 上海电力学院 | 基于变分模态分解和排列熵的滚动轴承故障诊断方法 |
| CN106053597A (zh) * | 2016-05-20 | 2016-10-26 | 湖南科技大学 | 一种基于hht特征提取的钢管混凝土脱空缺陷的检测方法 |
-
2017
- 2017-09-22 CN CN201710864339.XA patent/CN107607065A/zh active Pending
Patent Citations (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN103034757A (zh) * | 2012-12-02 | 2013-04-10 | 中国科学院电工研究所 | 基于经验模态分解的风电场时频域建模方法 |
| CN105160077A (zh) * | 2015-08-07 | 2015-12-16 | 中国电力科学研究院 | 一种高压直流线路电晕电流变化趋势提取方法 |
| CN105758644A (zh) * | 2016-05-16 | 2016-07-13 | 上海电力学院 | 基于变分模态分解和排列熵的滚动轴承故障诊断方法 |
| CN106053597A (zh) * | 2016-05-20 | 2016-10-26 | 湖南科技大学 | 一种基于hht特征提取的钢管混凝土脱空缺陷的检测方法 |
Non-Patent Citations (3)
| Title |
|---|
| （俄罗斯）勃罗斯库列亚柯夫: "《数与多项式》", 31 January 2016 * |
| 曾向阳: "《智能水中目标识别》", 31 March 2016 * |
| 郑近德 等: "广义变分模态分解方法及其在变工况齿轮故障诊断中的应用", 《振动工程学报》 * |
Cited By (19)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN108596215B (zh) * | 2018-04-04 | 2020-11-17 | 中国航发湖南动力机械研究所 | 多模态信号解析分离方法、装置、设备及存储介质 |
| CN108596215A (zh) * | 2018-04-04 | 2018-09-28 | 中国航发湖南动力机械研究所 | 多模态信号解析分离方法、装置、设备及存储介质 |
| CN108919203A (zh) * | 2018-05-21 | 2018-11-30 | 西安电子科技大学 | 一种雷达有源干扰的识别方法及系统 |
| CN108875170A (zh) * | 2018-06-05 | 2018-11-23 | 天津大学 | 一种基于改进变分模态分解的噪声源识别方法 |
| CN109323832A (zh) * | 2018-09-12 | 2019-02-12 | 温州大学 | 一种冷镦成型机模具冲击状态的监测方法 |
| CN109359633B (zh) * | 2018-12-10 | 2021-12-07 | 西北工业大学 | 基于希尔伯特-黄变换和小波脊线的信号联合分类方法 |
| CN109359633A (zh) * | 2018-12-10 | 2019-02-19 | 西北工业大学 | 基于希尔伯特-黄变换和小波脊线的信号联合分类方法 |
| CN109948249A (zh) * | 2019-03-19 | 2019-06-28 | 东南大学 | 基于变分模态分解的斜拉索时变索力识别方法 |
| CN110082436A (zh) * | 2019-04-25 | 2019-08-02 | 电子科技大学 | 一种基于变分模态的高提离电磁超声信号消噪方法 |
| CN110082436B (zh) * | 2019-04-25 | 2022-01-11 | 电子科技大学 | 一种基于变分模态的高提离电磁超声信号消噪方法 |
| CN110319995B (zh) * | 2019-08-14 | 2020-12-15 | 清华大学 | 火工冲击响应数据时频谱图分析方法 |
| CN110319995A (zh) * | 2019-08-14 | 2019-10-11 | 清华大学 | 火工冲击响应数据时频谱图分析方法 |
| CN111401168A (zh) * | 2020-03-06 | 2020-07-10 | 上海神添实业有限公司 | 一种无人机的多层雷达特征提取与选择方法 |
| CN111401168B (zh) * | 2020-03-06 | 2023-11-17 | 上海神添实业有限公司 | 一种无人机的多层雷达特征提取与选择方法 |
| CN112924655A (zh) * | 2021-01-28 | 2021-06-08 | 北京中水科海利工程技术有限公司 | 一种渡槽高排架结构混凝土质量的检测方法 |
| CN115436486A (zh) * | 2021-11-12 | 2022-12-06 | 江西理工大学 | 一种基于Hilbert变换的岩体破裂声发射信号伪维格纳分布二次时频分析方法 |
| CN115436486B (zh) * | 2021-11-12 | 2024-04-16 | 江西理工大学 | 一种基于Hilbert变换的岩体破裂声发射信号伪维格纳分布二次时频分析方法 |
| CN116013240A (zh) * | 2023-01-07 | 2023-04-25 | 广西大学 | 基于变分模态分解与数字滤波的钢管混凝土信号降噪方法 |
| CN116013240B (zh) * | 2023-01-07 | 2023-10-31 | 广西大学 | 基于变分模态分解与数字滤波的钢管混凝土信号降噪方法 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN107607065A (zh) | 一种基于变分模态分解的冲击回波信号分析方法 | |
| CN101887119B (zh) | 基于子带anmf海杂波中动目标检测方法 | |
| CN107329932B (zh) | 基于非线性调频分量分解的时频域模态参数辨识方法 | |
| CN112101245B (zh) | 基于频域窗函数的短时傅里叶变换机械冲击特征提取方法 | |
| CN111665489B (zh) | 一种基于目标特性的线谱提取方法 | |
| CN101627322A (zh) | 可控震源采集方法 | |
| CN104793253A (zh) | 基于数学形态学的航空电磁数据去噪方法 | |
| CN101893698B (zh) | 噪声源测试分析方法及其装置 | |
| CN101674520A (zh) | 振动式拾音器参数分析测试方法和装置 | |
| CN109340586B (zh) | 一种供水管道泄露的检测方法及系统 | |
| CN104112072A (zh) | 基于小波阈值去噪的主成分分析的工作模态参数识别方法 | |
| CN104680541A (zh) | 基于相位一致性的遥感图像质量评价方法 | |
| CN106534014A (zh) | 一种多分量lfm信号的精确检测与分离方法 | |
| CN105637331A (zh) | 异常检测装置、异常检测方法及计算机可读取的存储介质 | |
| CN104215456A (zh) | 一种基于平面聚类和频域压缩感知重构的机械故障诊断方法 | |
| Wei et al. | Incipient fault diagnosis of rolling element bearing based on wavelet packet transform and energy operator | |
| CN106199185A (zh) | 一种基于连续对数扫频的线性脉冲响应测量方法及系统 | |
| Avanesians et al. | Wave separation: Application for arrival time detection in ultrasonic signals | |
| CN109581516A (zh) | 曲波域统计量自适应阈值探地雷达数据去噪方法及系统 | |
| CN106548031A (zh) | 一种结构模态参数识别方法 | |
| CN110057918B (zh) | 强噪声背景下的复合材料损伤定量识别方法及系统 | |
| CN113189641B (zh) | 一种两道多模式瑞利波地下探测系统及方法 | |
| CN102590348A (zh) | 一种检测物体局部阻抗变化的方法及系统 | |
| Li et al. | Magnetotelluric signal-noise separation method based on SVM–CEEMDWT | |
| CN113188797A (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 | ||
| RJ01 | Rejection of invention patent application after publication | ||
| RJ01 | Rejection of invention patent application after publication |
Application publication date: 20180119 |