CN107748734A - One kind parsing empirical mode decomposition method - Google Patents
One kind parsing empirical mode decomposition method Download PDFInfo
- Publication number
- CN107748734A CN107748734A CN201711049787.0A CN201711049787A CN107748734A CN 107748734 A CN107748734 A CN 107748734A CN 201711049787 A CN201711049787 A CN 201711049787A CN 107748734 A CN107748734 A CN 107748734A
- Authority
- CN
- China
- Prior art keywords
- mrow
- imf
- msub
- mode
- decomposition
- 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
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
- G06F17/10—Complex mathematical operations
- G06F17/14—Fourier, Walsh or analogous domain transformations, e.g. Laplace, Hilbert, Karhunen-Loeve, transforms
- G06F17/141—Discrete Fourier transforms
- G06F17/142—Fast Fourier transforms, e.g. using a Cooley-Tukey type algorithm
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F30/00—Computer-aided design [CAD]
- G06F30/20—Design optimisation, verification or simulation
Abstract
The invention discloses one kind to parse empirical mode decomposition method, comprises the following steps：S1, by empirical mode decomposition method primary signal is sieved, obtain n IMF component and 1 residual volume；S2, by resolution modalities decomposition algorithm IMF components are handled, then according to certain rule reconstruct mode.The present invention is for the existing screening stop criterion problem of empirical mode decomposition, propose the screening stop criterion based on valid data section, for mode aliasing, propose a kind of parsing Empirical mode decomposition, the mode obscured mode occurs carries out interpretive model decomposition again, finally reconstructs these mode.The parsing empirical mode decomposition method of the present invention can effectively suppress mode aliasing, and the time decomposed is shorter, and method execution efficiency is higher.
Description
Technical field
The present invention relates to a kind of parsing-empirical mode decomposition method.
Background technology
Empirical mode decomposition (EMD) is a kind of widely used adaptive signal processing method, for empirical mode
Decompose and stop criterion problem is sieved present in implementation process, it is proposed that the screening stop criterion based on valid data section, and with
Traditional screening stop criterion compares, and new screening criterion avoids the influence of end effect, improves empirical modal
The accuracy of decomposition.For mode aliasing, it is proposed that parsing-Empirical mode decomposition, mould occurs to empirical mode decomposition
The mode that state is obscured carries out interpretive model decomposition again, finally reconstructs these mode, and the mode of obtained each rank mode is obscured bright
It is aobvious to be less than empirical mode decomposition.Using set forth herein method and existing empirical mode decomposition main improved method-totality
Average Empirical mode decomposition compares, and comparative result shows, this method and population mean empirical mode decomposition method one
Sample, it can effectively suppress mode aliasing.But the effect that this method suppression mode is obscured is more obvious, execution efficiency
It is higher.
When empirical mode decomposition (Empirical Mode Decomposition, EMD) method is based on signal local feature
Between yardstick, can the signal decomposition of complexity into limited individual intrinsic mode function (Intrinsic Mode Function, referred to as
IMF) sum, each intrinsic mode function represent a natural mode of vibration of original signal, then smoothly intrinsic to these
Mode function carries out analyzing and processing can be to extract the characteristic information of signal.
Empirical mode decomposition since the advent of the world just causes the very big concern of numerous scholars.Be widely used in economics,
The fields such as biomedical and engineering science, but there is also some shortcomings, mainly include, mode mixing (Mode Mixing), end points
The stop condition (Stop Criteria) and mistake envelope of effect (End Effects) and Intrinsic mode function (IMF) criterion,
Owe Inclusion etc..Empirical mode decomposition, which has these defects, will restrict the use of its further genralrlization.
Condition due to judging Intrinsic mode function in empirical mode decomposition method is that the average value of lower envelope line thereon is
Zero, but in true decomposable process, influence due to the interference of spline interpolation, end effect and using frequency etc.
Factor, the envelope average value for the Intrinsic mode function component that empirical mode decomposition obtains can not possibly be zero.Therefore, it is necessary to define
The Intrinsic mode function criterion problem of the stopping criterion of one feasible decomposable process in engineering, i.e. empirical mode decomposition is also referred to as
Sieve stop criterion problem.
Mode obscures mode mixing and refers to include difference in an IMF (Intrinsic Mode Function)
Great characteristic time scale, or similar characteristic time scale are distributed in different IMF, cause 2 adjacent IMF ripples
Shape aliasing, influences each other, it is difficult to recognizes.
Mode obscure be EMD a major defect, once mode obscures generation, finally give Intrinsic mode function component
Physical significance just less clearly, influence the accuracy of signal decomposition, seriously hinder its application in engineering.
For the mode confounding issues of empirical mode decomposition, scholars propose some improved methods.Foremost is Huang
The population mean empirical mode decomposition method (Ensemble Empirical Mode Decomposition, EEMD) of proposition.Should
Method has the statistical property of frequency-flat distribution using white Gaussian noise, makes the signal after addition white Gaussian noise in different chis
There is continuity on degree, therefore solve the problems, such as mode mixing on some scale.But it be present, for example, addition
The selection of the amplitude of the number of white noise and added white noise has a very strong subjectivity, and method adaptivity is poor.In addition, work as
Number and the amplitude selection of the white noise of addition are reasonable, are suppressing artificially cause low frequency again while high frequency mode is obscured
Mode is obscured.Moreover, the algorithm of overall part mean decomposition method is complex, program runtime is longer.
The content of the invention
It is an object of the invention to overcome the deficiencies of the prior art and provide a kind of mode obscured mode occurs to be carried out again
Interpretive model decomposes, and finally reconstructs these mode.Parsing-empirical mode decomposition method of the present invention can effectively suppress mode
Aliasing, the time of decomposition is shorter, the higher parsing-empirical mode decomposition method of method execution efficiency.
The purpose of the present invention is achieved through the following technical solutions：A kind of parsing-empirical mode decomposition method, including
Following steps：
S1, by empirical mode decomposition method primary signal is sieved, obtain n IMF component and 1 remnants
Amount；
S2, by resolution modalities decomposition algorithm IMF components are handled, reconstruct mode.
Further, the step S1 includes following sub-step：
S11, all maximum points and minimum point for determining primary signal x (t), are fitted respectively using cubic spline curve
Go out the upper and lower envelope of extreme point, and ensure that signal x (t) between upper and lower envelope, then calculates the office of upper and lower envelope
Portion's average, is denoted as m1；
S12, calculate signal x (t) and the difference h of local mean value1(t)：
h1(t)=x (t)-m1(1)；
S13, judge h1(t) condition that IMF is set up can be met；If not satisfied, then by h1(t) as new primary signal,
Repeat step S11 and S12, until h1(t) condition that IMF is set up is met；
S14, make C1(t)=h1(t), C1(t) first obtained IMF component is as decomposed；
S15, by C1(t) separated from signal x (t), obtained residual signal is designated as：r1(t)=x (t)-C1(t)；
S16, by r1(t) as new primary signal, repeat step S11~S15, n limited circulation is carried out, obtains n
Individual IMF components and 1 residual volume rn(t)；I.e.：
Primary signal x (t) is as：Wherein, Ci(t) i-th of Intrinsic mode function point is represented
Amount, rn(t) it is referred to as residual error function.
Further,
Used in the step S16 based on the screening stop criterion of valid data section to calculate cycle-index n, concrete operations
Method is：As Intrinsic mode function Cn(t) when meeting following condition, circulation is stopped：
In formula,Represent Cn(t) effective section of upper lower envelope Mean curve, δ represent deviation threshold, and span is
0.01~0.1.
Further, the step S2 includes following sub-step：
S21, the n Intrinsic mode function component that step S1 is calculated are designated as IMF successively1(t)~IMFn(t)；
To IMF1(t) component does Fast Fourier Transform (FFT), obtains its amplitude spectrum, and boundary frequency f is determined by amplitude spectrumb1；
To IMF1(t) interpretive model decomposition is carried out, by boundary frequency fb1Decomposition obtains two signal c1(t) andc1(t) represent
IMF1(t) correction component,It is c1(t) residual components；I.e.：
S22, generalIt is added to IMF2(t) in, the IMF after being updated2(t), it is denoted as IMF2 *(t)；To IMF2 *(t) enter
Row Fast Fourier Transform (FFT), obtains IMF2 *(t) amplitude spectrum, and boundary frequency f is determined by amplitude spectrumb2；To IMF2 *(t) enter
Row interpretive model decomposes, by boundary frequency fb2Obtain c2(t) andWherein, c2(t) it is exactly IMF2 *(t) correction component,It is c2(t) residual components；I.e.：
S23, to except IMFn(t) S21's and S22's the operation of all Intrinsic mode function component repeat steps beyond；
S24, for IMFn(t) following operate is carried out：WillIt is added to IMFn(t) in, IMF is obtained* n(t), to IMF* n
(t) do Fast Fourier Transform (FFT) and draw its amplitude spectrum, determine boundary frequency fbk；To IMF* n(t) interpretive model decomposition is carried out,
By boundary frequency fbkObtain signal cn(t) andI.e.：
WillIt is added to residual error function rn(t) in, final residual error is obtained, is denoted as vn(t)。
The beneficial effects of the invention are as follows：
1st, the present invention is directed to mode aliasing existing for empirical mode decomposition, it is proposed that a kind of parsing-empirical modal point
Solution, the mode obscured mode occurs carry out interpretive model decomposition, finally reconstruct these mode again.Parsing-warp of the present invention
Mode aliasing can effectively be suppressed by testing mode decomposition method, and the time decomposed is shorter, and method execution efficiency is higher.
2nd, existing intrinsic mode function screening stop criterion problem of the present invention for empirical mode decomposition, it is proposed that base
Stop criterion is sieved in valid data section, is contrasted with traditional screening stop criterion, is terminated using the screening of the present invention
The result for the empirical mode decomposition that criterion obtains is more accurate.
Brief description of the drawings
Fig. 1 is the time-domain diagram of emulation signal x (t) and the two composition component used in checking test (one)；
Fig. 2 is using the decomposition result figure of the loop termination criterion based on valid data section in checking test (one)；
Fig. 3 is using the decomposition result figure of the Huang loop termination criterions proposed in checking test (one)；
Fig. 4 is the vibration signal time-domain diagram of the gear tooth breakage failure used in checking test (two)；
Fig. 5 is first Intrinsic mode function component of vibration signal parsing-empirical mode decomposition of gear tooth breakage failure
Figure；
Fig. 6 is the time domain beamformer of emulation signal x (t) used in checking test (three)；
Fig. 7 is the emulation signal empirical mode decomposition result figure of checking test (three)；
Fig. 8 is the emulation signal population mean empirical mode decomposition result figure of checking test (three)；
Fig. 9 is emulation signal resolution-empirical mode decomposition result figure of checking test (three)；
Figure 10 is the time-domain diagram of rotor fault signal used in checking test (four)；
Figure 11 is rotor fault signal local mean value decomposition result；
Figure 12 is the result of rotor fault signal resolution-empirical mode decomposition.
Embodiment
Technical scheme is further illustrated below in conjunction with the accompanying drawings.
The invention discloses a kind of parsing-empirical mode decomposition method, comprise the following steps：
S1, by empirical mode decomposition method primary signal is sieved, obtain n IMF component and 1 remnants
Amount；Empirical mode decomposition is the process of one " screening "；First by the minimum component " sieve " of extrema elimination characteristic dimension out, so
" sieve " afterwards and go out larger characteristic dimension component, the maximum component " sieve " of characteristic dimension out to the last.That is, experience
Mode decomposition obtains Intrinsic mode function component C1, C2... ... average frequency gradually reduce.This step specifically include with
Lower sub-step：
S11, all maximum points and minimum point for determining primary signal x (t), are fitted respectively using cubic spline curve
Go out the upper and lower envelope of extreme point, and ensure that signal x (t) between upper and lower envelope, then calculates the office of upper and lower envelope
Portion's average, is denoted as m1；
S12, calculate signal x (t) and the difference h of local mean value1(t)：
h1(t)=x (t)-m1(1)；
S13, judge h1(t) condition that IMF is set up can be met；If not satisfied, then by h1(t) as new primary signal,
Repeat step S11 and S12, until h1(t) condition that IMF is set up is met；
S14, make C1(t)=h1(t), C1(t) first obtained IMF component is as decomposed；
S15, by C1(t) separated from signal x (t), obtained residual signal is designated as：r1(t)=x (t)-C1(t)；
S16, by r1(t) as new primary signal, repeat step S11~S15, n limited circulation is carried out, obtains n
Individual IMF components and 1 residual volume rn(t)；I.e.：
Primary signal x (t) is as：Wherein, Ci(t) i-th of Intrinsic mode function point is represented
Amount, rn(t) it is referred to as residual error function.
R is worked as in above-mentioned steps S16 cyclic processn(t) loop termination when being a monotonic function, but actual decomposable process
In, work as rn(t) when meeting monotonic function condition, the number of circulation is often very big, decomposes the Intrinsic mode function finally given
Component number will be a lot.Found from the experiment of substantial amounts of empirical mode decomposition, some last that empirical mode decomposition obtains
Component is all false Intrinsic mode function component mostly, does not have tangible physical significance.
The stopping criterion for the decomposable process that Huang is proposed is by limiting the big of the standard deviation between 2 adjacent IMF components
It is small to realize, pass through iteration threshold SdTo control screening number.SdIt is defined as：
In formula：T is the time span of signal；hk(t), hk-1(t) it is two adjacent IMF components in decomposable process；Work as Sd
Meet to stop circulation during 0.2~0.3 scope.
It is determined that rational iteration threshold is extremely important, if threshold value is too small, it will amount of calculation is greatly increased, finally in acquisition
It is also meaningless to report mode function component.If threshold value is too big, it will be difficult to meets the decomposition result of Intrinsic mode function condition.
In most cases effect is obvious to the screening stop criterion that Huang is defined.But it there is also one
A little problems.Work as h it can be seen from formula (3)k(t) be 0 when, the iteration threshold of gained will be a uncertain value, if again with
It is obviously improper that it makees screening criterion.In addition, the screening stop criterion does not take into full account the influence of end effect, mistake has been used
After genuine end-point data, final decomposition result error can be caused larger.Under many circumstances, the end points effect of empirical mode decomposition
Should be obvious, even if taking some way to suppress end effect, end effect is still existing.
For SdThe problem of criterion is present, terminated using the screening based on valid data section in step S16 of the present invention
Criterion calculates cycle-index n, and so-called effectively section refers to the remaining data section for removing two-end-point distortion data section.The criterion is filled
Point consider end effect to sieve criterion influence, only with valid data section data as calculating iteration threshold basis,
Avoid the influence of empirical mode decomposition end effect.The concrete operation method of this step is：As Intrinsic mode function Cn(t) it is full
During foot row condition, stop circulation：
In formula,Represent Cn(t) effective section of upper lower envelope Mean curve, δ represent deviation threshold, and span is
0.01~0.1.
Checking test (one)；Illustrate the loop termination criterion of the present invention in empirical mode decomposition by an emulation signal
Application in method, and compared with the criterion of N.E.Huang proposition.Emulating signal x (t) expression formula is：
X (t)=x1(t)+x2(t),t∈[0,1] (5)
x1(t)=(1+0.5sin (2 π * 4t)) * 2cos (2 π * 80t), x2(t)=sin (2 π * 24t) (6)
It is by x to emulate signal x (t)1And x (t)2(t) the two compositions are formed.Emulate signal and two composition components such as
Shown in Fig. 1.Result using the empirical mode decomposition of loop termination criterion proposed by the present invention is as shown in Figure 2.Carried using Huang
The empirical mode decomposition result of the stop criterion gone out is as shown in Figure 3.Comparison diagram 1, Fig. 2 and Fig. 3 are understood, utilize the circulation of the present invention
The component C that stop criterion obtains1The component imf obtained with the Huang stop criterions proposed1All relatively desired result x1
(t), C2And imf2All relatively desired result x2(t), no matter this explanation uses any screening stop criterion, empirical modal
Decomposition can be by the component of the two in primary signal x1And x (t)2(t) separate exactly.But C1Two-end-point near
Curve ratio imf1Be more nearly ideal curve x1(t), C2The neighbouring curve ratio imf of two-end-point2Be more nearly ideal curve x2
(t), this explanation uses loop termination criterion proposed by the present invention, more accurate than the result of empirical mode decomposition.
Checking test (two)：Fig. 4 show the Gearbox vibration signal of a generation broken teeth failure.The number of gear teeth is z=75,
Modulus m=2, turn frequency f=n2/ 60 ≈ 13.6Hz, now using Empirical mode decomposition decomposition vibration signal, empirical mode decomposition side
The screening end condition of method uses criterion proposed by the present invention.The final result of decomposition is 5 Intrinsic mode function component imf1,
imf2,…imf5。imf1Bulk information comprising gear distress, first to imf1Analyzed, imf1Time domain it is as shown in Figure 5.From
1st Intrinsic mode function component imf1Waveform in it can be seen that obvious modulation signature, the cycle T of modulating wave=
0.074s, frequency are approximately 13.6Hz, just equal with turn frequency of failure gear.Illustrate using proposed by the present invention based on effective
The empirical modal method of the loop termination criterion of data segment can decompose the letter that is out of order exactly to the Gearbox vibration signal of reality
Breath.
S2, by resolution modalities decomposition algorithm IMF components are handled, reconstruct mode；For empirical mode decomposition
Mode aliasing, this paper presents parsing-Empirical mode decomposition, using resolution modalities decomposition method to empirical mode decomposition method
The mode that mode is obscured occurs to be decomposed again, then reconstructs these mode, largely reduces mode aliasing.
Resolution modalities decomposition method (Analytical Mode Decomposition, AMD) can be isolated from signal
The harmonic components of each frequency band, this method principle are as follows：
Signal x0(t) by boundary frequency fbIt is decomposed into two parts：One fast changed signal xfAnd a tempolabile signal x (t)s
(t)：
x0(t)=xf(t)+xs(t) (7)
Corresponding Fourier transformation is designated as X respectivelyf(ω)、Xs(ω), the two is non-overlapping on frequency band.Respectively to cos (2 π
fbt)·x0(t) with sin (2 π fbt)·x0(t) make Hilbert transform, obtain：
H[cos(2πfbt)·x0(t)]=xs(t)·H[cos(2πfbt)]+cos(2πfbt)·H[xf(t)] (8)
H[sin(2πfbt)·x0(t)]=xs(t)·H[sin(2πfbt)]+sin(2πfbt)·H[xf(t)] (9)
X is derived by above formulafAnd x (t)s(t) computational methods are：
xs(t)=sin (2 π fbt)·H[cos(2πfbt)·x0(t)]-cos(2πfbt)·H[sin(2πfbt)·x0
(t)] (10)
xf(t)=x0(t)-xs(t) (11)。
The step S2 specifically includes following sub-step：
S21, the n Intrinsic mode function component that step S1 is calculated are designated as IMF successively1(t)~IMFn(t)；
To IMF1(t) component does Fast Fourier Transform (FFT), obtains its amplitude spectrum, and boundary frequency f is determined by amplitude spectrumb1；
To IMF1(t) interpretive model decomposition is carried out, by boundary frequency fb1Decomposition obtains two signal c1(t) andc1(t) represent
IMF1(t) correction component,It is c1(t) residual components；I.e.：
S22, generalIt is added to IMF2(t) in, the IMF after being updated2(t), it is denoted as IMF2 *(t)；To IMF2 *(t) enter
Row Fast Fourier Transform (FFT), obtains IMF2 *(t) amplitude spectrum, and boundary frequency f is determined by amplitude spectrumb2；To IMF2 *(t) enter
Row interpretive model decomposes, by boundary frequency fb2Obtain c2(t) andWherein, c2(t) it is exactly IMF2 *(t) correction component,It is c2(t) residual components；I.e.：
S23, to except IMFn(t) S21's and S22's the operation of all Intrinsic mode function component repeat steps beyond；
S24, for IMFn(t) following operate is carried out：WillIt is added to IMFn(t) in, IMF is obtained* n(t), to IMF* n
(t) do Fast Fourier Transform (FFT) and draw its amplitude spectrum, determine boundary frequency fbk；To IMF* n(t) interpretive model decomposition is carried out,
By boundary frequency fbkObtain signal cn(t) andI.e.：
WillIt is added to residual error function rn(t) in, final residual error is obtained, is denoted as vn(t)。
Checking test (three)：For the validity and advantage of the parsing-empirical mode decomposition method that prove out, using experience
Parsing-Empirical mode decomposition of Mode Decomposition, population mean Empirical mode decomposition and the present invention is believed emulation respectively
Number decomposed, and contrast their decomposition result.
Emulating signal x (t) is：
X (t)=x1(t)+x2(t)+x3(t)；t∈[0,1] (15)
x2(t)=2cos (60 π t)；x3(t)=2.5cos (28 π t)；t∈[0,1] (17)
It is as shown in Figure 6 to emulate signal.The result of empirical mode decomposition is as shown in Figure 7.From figure 7 it can be seen that the knot decomposed
Fruit mode aliasing is than more serious, first mode C1(t), high frequency intermittent signal and 30Hz sinusoidal signal are obscured one
Rise；Second mode C2(t), together with 30Hz sinusoidal signals are obscured with 14Hz sinusoidal signal；3rd mode with the 4th
Mode is false Intrinsic mode function component.Population mean times N=30 of population mean Empirical mode decomposition, noise amplitude
The standard deviation of 0.01 times of primary signal is set to, the result of population mean empirical mode decomposition is imf1、imf2、imf3、imf4, such as
Shown in Fig. 8.The result of parsing-empirical mode decomposition is as shown in Figure 9.In program operation, the operation of these three methods is recorded
Time (6 run times), as shown in table 1.In table 1, EMD represents Empirical mode decomposition, and EEMD represents population mean experience
Mode Decomposition, AEMD represent parsing-Empirical mode decomposition.
1 three kinds of method run time (unit s) of table
| Method | For the first time | Second | For the third time | 4th time | 5th time | 6th time |
| EMD | 1.3988 | 1.3860 | 1.3930 | 1.3846 | 1.3820 | 1.3940 |
| EEMD | 41.4875 | 44.8383 | 44.1800 | 42.8480 | 43.2321 | 44.024 |
| AEMD | 1.4027 | 1.3832 | 1.3843 | 1.3934 | 1.3880 | 1.3848 |
By the decomposition result for contrasting these three methods, it is found that population mean empirical mode decomposition and parsing experience
Mode Decomposition, which can suppress mode, to be obscured, but the effect for parsing empirical mode decomposition is better than population mean empirical mode decomposition
Method.It can be found by contrasting these three method resolving times, the run time of parsing-Empirical mode decomposition will be significantly lower than always
The average Empirical mode decomposition of body, it is approximate with the run time of empirical modal method.Therefore, parsing-empirical mode decomposition has
Certain advantage.
Checking test (four)：Figure 10 is the vibration signal of Rotor Rubbing Fault.The rotating speed of rotor by input motor control, if
The rotating speed of rotor is 300r/min.Displacement transducer is horizontally mounted, and the vibration of partial rub failure occurs for measuring rotor
Signal.Setting sensor sample frequency 8kHz.This is decomposed using improved method-local mean value method of Empirical mode decomposition to shake
Dynamic signal, as a result as shown in figure 11.It is as shown in figure 12 using parsing-empirical mode decomposition of the present invention, decomposition result.Comparison diagram
11 understand with Figure 12, and the result mode of local mean value method is obscured more seriously, and the mode of parsing empirical mode decomposition result is mixed
Confuse weaker, illustrate that parsing empirical modal can suppress mode and obscure, have certain practical value to analysis rotor oscillation signal.
One of ordinary skill in the art will be appreciated that embodiment described here is to aid in reader and understands this hair
Bright principle, it should be understood that protection scope of the present invention is not limited to such especially statement and embodiment.This area
Those of ordinary skill can make according to these technical inspirations disclosed by the invention various does not depart from the other each of essence of the invention
The specific deformation of kind and combination, these deform and combined still within the scope of the present invention.
Claims (4)
1. a kind of parsing-empirical mode decomposition method, it is characterised in that comprise the following steps：
S1, by empirical mode decomposition method primary signal is sieved, obtain n IMF component and 1 residual volume；
S2, by resolution modalities decomposition algorithm IMF components are handled, reconstruct mode.
2. a kind of parsing-empirical mode decomposition method according to claim 1, it is characterised in that the step S1 includes
Following sub-step：
S11, all maximum points and minimum point for determining primary signal x (t), pole is fitted using cubic spline curve respectively
It is worth the upper and lower envelope of point, and ensures signal x (t) between upper and lower envelope, then calculates the local equal of upper and lower envelope
Value, is denoted as m1；
S12, calculate signal x (t) and the difference h of local mean value1(t)：
h1(t)=x (t)-m1(1)；
S13, judge h1(t) condition that IMF is set up can be met；If not satisfied, then by h1(t) as new primary signal, repeat
Step S11 and S12, until h1(t) condition that IMF is set up is met；
S14, make C1(t)=h1(t), C1(t) first obtained IMF component is as decomposed；
S15, by C1(t) separated from signal x (t), obtained residual signal is designated as：r1(t)=x (t)-C1(t)；
S16, by r1(t) as new primary signal, repeat step S11~S15, n limited circulation is carried out, obtains n IMF
Component and 1 residual volume rn(t)；I.e.：
<mrow>
<mfenced open = "{" close = "">
<mtable>
<mtr>
<mtd>
<mrow>
<msub>
<mi>r</mi>
<mn>2</mn>
</msub>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<msub>
<mi>r</mi>
<mn>1</mn>
</msub>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
<mo>-</mo>
<msub>
<mi>C</mi>
<mn>2</mn>
</msub>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
</mrow>
</mtd>
</mtr>
<mtr>
<mtd>
<mrow>
<msub>
<mi>r</mi>
<mn>3</mn>
</msub>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<msub>
<mi>r</mi>
<mn>2</mn>
</msub>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
<mo>-</mo>
<msub>
<mi>C</mi>
<mn>3</mn>
</msub>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
</mrow>
</mtd>
</mtr>
<mtr>
<mtd>
<mtable>
<mtr>
<mtd>
<mo>.</mo>
</mtd>
</mtr>
<mtr>
<mtd>
<mo>.</mo>
</mtd>
</mtr>
<mtr>
<mtd>
<mo>.</mo>
</mtd>
</mtr>
</mtable>
</mtd>
</mtr>
<mtr>
<mtd>
<mrow>
<msub>
<mi>r</mi>
<mi>n</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<msub>
<mi>r</mi>
<mrow>
<mi>n</mi>
<mo>-</mo>
<mn>1</mn>
</mrow>
</msub>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
<mo>-</mo>
<msub>
<mi>C</mi>
<mi>n</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>2</mn>
<mo>)</mo>
</mrow>
</mrow>
Primary signal x (t) is as：Wherein, Ci(t) i-th of Intrinsic mode function component, r are representedn
(t) it is referred to as residual error function.
3. a kind of parsing-empirical mode decomposition method according to claim 2, it is characterised in that adopted in the step S16
With calculating cycle-index n based on the screening stop criterion of valid data section, concrete operation method is：As Intrinsic mode function Cn
(t) when meeting following condition, circulation is stopped：
<mrow>
<mo>-</mo>
<mi>&delta;</mi>
<mo><</mo>
<msub>
<mover>
<mi>m</mi>
<mo>~</mo>
</mover>
<mi>k</mi>
</msub>
<mo><</mo>
<mi>&delta;</mi>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>3</mn>
<mo>)</mo>
</mrow>
</mrow>
In formula,Represent Cn(t) effective section of upper lower envelope Mean curve, δ represent deviation threshold, span be 0.01~
0.1。
4. a kind of parsing-empirical mode decomposition method according to claim 2, it is characterised in that the step S2 includes
Following sub-step：
S21, the n Intrinsic mode function component that step S1 is calculated are designated as IMF successively1(t)~IMFn(t)；
To IMF1(t) component does Fast Fourier Transform (FFT), obtains its amplitude spectrum, and boundary frequency f is determined by amplitude spectrumb1；To IMF1
(t) interpretive model decomposition is carried out, by boundary frequency fb1Decomposition obtains two signal c1(t) andc1(t) IMF is represented1(t)
Correction component,It is c1(t) residual components；I.e.：
<mrow>
<msub>
<mi>IMF</mi>
<mn>1</mn>
</msub>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<msub>
<mi>c</mi>
<mn>1</mn>
</msub>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
<mo>+</mo>
<msubsup>
<mi>c</mi>
<mn>1</mn>
<mo>^</mo>
</msubsup>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>3</mn>
<mo>)</mo>
</mrow>
<mo>;</mo>
</mrow>
S22, generalIt is added to IMF2(t) in, the IMF after being updated2(t), it is denoted as IMF2 *(t)；To IMF2 *(t) carry out fast
Fast Fourier transformation, obtains IMF2 *(t) amplitude spectrum, and boundary frequency f is determined by amplitude spectrumb2；To IMF2 *(t) solved
Mode Decomposition is analysed, by boundary frequency fb2Obtain c2(t) andWherein, c2(t) it is exactly IMF2 *(t) correction component,It is
c2(t) residual components；I.e.：
<mrow>
<msubsup>
<mi>IMF</mi>
<mn>2</mn>
<mo>*</mo>
</msubsup>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<msub>
<mi>c</mi>
<mn>2</mn>
</msub>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
<mo>+</mo>
<msubsup>
<mi>c</mi>
<mn>2</mn>
<mo>^</mo>
</msubsup>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>5</mn>
<mo>)</mo>
</mrow>
<mo>;</mo>
</mrow>
S23, to except IMFn(t) S21's and S22's the operation of all Intrinsic mode function component repeat steps beyond；
S24, for IMFn(t) following operate is carried out：WillIt is added to IMFn(t) in, IMF is obtained* n(t), to IMF* n(t) do
Fast Fourier Transform (FFT) simultaneously draws its amplitude spectrum, determines boundary frequency fbk；To IMF* n(t) interpretive model decomposition is carried out, by demarcating
Frequency fbkObtain signal cn(t) andI.e.：
<mrow>
<msubsup>
<mi>IMF</mi>
<mi>n</mi>
<mo>*</mo>
</msubsup>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
<mo>=</mo>
<msub>
<mi>c</mi>
<mi>n</mi>
</msub>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
<mo>+</mo>
<msubsup>
<mi>c</mi>
<mi>n</mi>
<mo>^</mo>
</msubsup>
<mrow>
<mo>(</mo>
<mi>t</mi>
<mo>)</mo>
</mrow>
<mo>-</mo>
<mo>-</mo>
<mo>-</mo>
<mrow>
<mo>(</mo>
<mn>6</mn>
<mo>)</mo>
</mrow>
</mrow>
WillIt is added to residual error function rn(t) in, final residual error is obtained, is denoted as vn(t)。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201711049787.0A CN107748734B (en) | 2017-10-31 | 2017-10-31 | Analytic-empirical mode decomposition method |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201711049787.0A CN107748734B (en) | 2017-10-31 | 2017-10-31 | Analytic-empirical mode decomposition method |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN107748734A true CN107748734A (en) | 2018-03-02 |
| CN107748734B CN107748734B (en) | 2021-08-27 |
Family
ID=61252821
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201711049787.0A Active CN107748734B (en) | 2017-10-31 | 2017-10-31 | Analytic-empirical mode decomposition method |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN107748734B (en) |
Cited By (11)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN108508268A (en) * | 2018-04-02 | 2018-09-07 | 山西大学 | A kind of novel digital type flickermeter decomposed based on time-frequency domain |
| CN108596215A (en) * | 2018-04-04 | 2018-09-28 | 中国航发湖南动力机械研究所 | Multi-modal signal resolution separation method, device, equipment and storage medium |
| CN109567743A (en) * | 2018-10-24 | 2019-04-05 | 加康康健有限公司 | A kind of signal reconfiguring method based on EMD, device, terminal device and storage medium |
| CN109582913A (en) * | 2018-11-20 | 2019-04-05 | 国家海洋局第海洋研究所 | A kind of method of empirical mode decomposition screening iterative process stop criterion |
| CN109799535A (en) * | 2019-03-14 | 2019-05-24 | 中船海洋探测技术研究院有限公司 | A kind of filtering method of full tensor magnetic gradient detection and localization data |
| CN109859174A (en) * | 2019-01-09 | 2019-06-07 | 东莞理工学院 | A kind of OLED defect inspection method based on empirical mode decomposition and regression model |
| CN110096673A (en) * | 2019-04-29 | 2019-08-06 | 河北工业大学 | A kind of EMD improved method suitable for signal decomposition |
| CN110441063A (en) * | 2019-06-12 | 2019-11-12 | 祝思宁 | A kind of method of monitoring, diagnosing large high-speed armature spindle crackle |
| CN111537989A (en) * | 2020-03-25 | 2020-08-14 | 中国电子科技集团公司第二十九研究所 | Method for extracting signal micro Doppler modulation component based on empirical mode decomposition |
| CN111553308A (en) * | 2020-05-11 | 2020-08-18 | 成都亿科康德电气有限公司 | Reconstruction method of partial discharge signal of power transformer |
| CN114936375A (en) * | 2022-06-02 | 2022-08-23 | 北京理工大学 | Asymmetric integrated optical encryption system based on two-dimensional empirical mode decomposition |
Citations (15)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20100074496A1 (en) * | 2008-09-23 | 2010-03-25 | Industrial Technology Research Institute | Multi-dimensional empirical mode decomposition (emd) method for image texture analysis |
| US20120101781A1 (en) * | 2010-10-20 | 2012-04-26 | Huafan University | Operation circuit and method thereof |
| CN102778357A (en) * | 2012-08-15 | 2012-11-14 | 重庆大学 | Mechanical failure feature extracting method based on optimal parameter ensemble empirical mode decomposition (EEMD) |
| CN103617684A (en) * | 2013-12-12 | 2014-03-05 | 威海北洋电气集团股份有限公司 | Interference type optical fiber perimeter vibration intrusion recognition algorithm |
| CN103870694A (en) * | 2014-03-18 | 2014-06-18 | 江苏大学 | Empirical mode decomposition denoising method based on revised wavelet threshold value |
| CN104006961A (en) * | 2014-04-29 | 2014-08-27 | 北京工业大学 | Cycloid bevel gear fault diagnosis method based on empirical mode decomposition and cepstrum |
| CN104007315A (en) * | 2014-04-29 | 2014-08-27 | 昆明理工大学 | Improved empirical mode decomposition processing method |
| CN105893773A (en) * | 2016-04-20 | 2016-08-24 | 辽宁工业大学 | Vibration signal frequency characteristic extraction method based on EEMD technology, CMF technology and WPT technology |
| CN105938542A (en) * | 2016-03-16 | 2016-09-14 | 南京大学 | Empirical-mode-decomposition-based noise reduction method for bridge strain signal |
| CN106019102A (en) * | 2016-06-27 | 2016-10-12 | 国网北京市电力公司 | Signal de-noising method and apparatus |
| US20170116155A1 (en) * | 2015-10-22 | 2017-04-27 | National Central University | System and method of conjugate adaptive conjugate masking empirical mode decomposition |
| CN106844293A (en) * | 2017-02-03 | 2017-06-13 | 中国铁道科学研究院 | A kind of adaptive decoupling method of modal overlap problem in empirical mode decomposition |
| CN106874549A (en) * | 2017-01-10 | 2017-06-20 | 西安理工大学 | A kind of discrete distribution parabolic equation method in the arrowband of high-precision forecast ASF |
| CN107064752A (en) * | 2017-03-22 | 2017-08-18 | 北京航空航天大学 | A kind of distinguished number of aviation fault electric arc detection |
| CN107255563A (en) * | 2017-06-28 | 2017-10-17 | 石家庄铁道大学 | Realize gear-box mixed fault signal blind source separation method |
-
2017
- 2017-10-31 CN CN201711049787.0A patent/CN107748734B/en active Active
Patent Citations (16)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US20100074496A1 (en) * | 2008-09-23 | 2010-03-25 | Industrial Technology Research Institute | Multi-dimensional empirical mode decomposition (emd) method for image texture analysis |
| US20120101781A1 (en) * | 2010-10-20 | 2012-04-26 | Huafan University | Operation circuit and method thereof |
| CN102778357A (en) * | 2012-08-15 | 2012-11-14 | 重庆大学 | Mechanical failure feature extracting method based on optimal parameter ensemble empirical mode decomposition (EEMD) |
| CN103617684A (en) * | 2013-12-12 | 2014-03-05 | 威海北洋电气集团股份有限公司 | Interference type optical fiber perimeter vibration intrusion recognition algorithm |
| CN103870694A (en) * | 2014-03-18 | 2014-06-18 | 江苏大学 | Empirical mode decomposition denoising method based on revised wavelet threshold value |
| CN104006961A (en) * | 2014-04-29 | 2014-08-27 | 北京工业大学 | Cycloid bevel gear fault diagnosis method based on empirical mode decomposition and cepstrum |
| CN104007315A (en) * | 2014-04-29 | 2014-08-27 | 昆明理工大学 | Improved empirical mode decomposition processing method |
| US20170116155A1 (en) * | 2015-10-22 | 2017-04-27 | National Central University | System and method of conjugate adaptive conjugate masking empirical mode decomposition |
| CN106610918A (en) * | 2015-10-22 | 2017-05-03 | 中央大学 | Empirical mode decomposition method and system for adaptive binary and conjugate shielding network |
| CN105938542A (en) * | 2016-03-16 | 2016-09-14 | 南京大学 | Empirical-mode-decomposition-based noise reduction method for bridge strain signal |
| CN105893773A (en) * | 2016-04-20 | 2016-08-24 | 辽宁工业大学 | Vibration signal frequency characteristic extraction method based on EEMD technology, CMF technology and WPT technology |
| CN106019102A (en) * | 2016-06-27 | 2016-10-12 | 国网北京市电力公司 | Signal de-noising method and apparatus |
| CN106874549A (en) * | 2017-01-10 | 2017-06-20 | 西安理工大学 | A kind of discrete distribution parabolic equation method in the arrowband of high-precision forecast ASF |
| CN106844293A (en) * | 2017-02-03 | 2017-06-13 | 中国铁道科学研究院 | A kind of adaptive decoupling method of modal overlap problem in empirical mode decomposition |
| CN107064752A (en) * | 2017-03-22 | 2017-08-18 | 北京航空航天大学 | A kind of distinguished number of aviation fault electric arc detection |
| CN107255563A (en) * | 2017-06-28 | 2017-10-17 | 石家庄铁道大学 | Realize gear-box mixed fault signal blind source separation method |
Non-Patent Citations (9)
| Title |
|---|
| HAN, J. 等: "Gear Fault Diagnosis Based on Improved EMD Method and the Energy Operator Demodulation Approach", 《JOURNAL OF CHANGSHA UNIVERSITY OF SCIENCE & TECHNOLOGY(NATURAL SCIENCE)》 * |
| PEIMING SHI 等: "Vibration Analysis as a Diagnosis Tool for Health Monitoring of Industrial Machines", 《SHOCK AND VIBRATION》 * |
| UDDIN MD BURHAN 等: "A new machine learning approach to select adaptive IMFs of EMD", 《2016 2ND INTERNATIONAL CONFERENCE ON ELECTRICAL, COMPUTER & TELECOMMUNICATION ENGINEERING (ICECTE)》 * |
| VAN TUAN DO 等: "Adaptive Empirical Mode Decomposition for Bearing Fault Detection", 《JOURNAL OF MECHANICAL ENGINEERING》 * |
| 林近山: "基于EEMD和Hilbert变换的齿轮箱故障诊断", 《机械传动》 * |
| 汤宝平 等: "基于形态奇异值分解和经验模态分解的滚动轴承故障特征提取方法", 《机械工程学报》 * |
| 蒋婷 等: "基于EMD与分批估计的动态称量快速融合方法", 《仪器仪表学报》 * |
| 郑近德 等: "非平稳信号分析的广义解析模态分解方法", 《电子学报》 * |
| 雷亚国 等: "自适应总体平均经验模式分解及其在行星齿轮箱故障检测中的应用", 《机械工程学报》 * |
Cited By (15)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN108508268A (en) * | 2018-04-02 | 2018-09-07 | 山西大学 | A kind of novel digital type flickermeter decomposed based on time-frequency domain |
| CN108596215A (en) * | 2018-04-04 | 2018-09-28 | 中国航发湖南动力机械研究所 | Multi-modal signal resolution separation method, device, equipment and storage medium |
| CN108596215B (en) * | 2018-04-04 | 2020-11-17 | 中国航发湖南动力机械研究所 | Multi-modal signal analysis and separation method, device, equipment and storage medium |
| CN109567743A (en) * | 2018-10-24 | 2019-04-05 | 加康康健有限公司 | A kind of signal reconfiguring method based on EMD, device, terminal device and storage medium |
| CN109582913A (en) * | 2018-11-20 | 2019-04-05 | 国家海洋局第海洋研究所 | A kind of method of empirical mode decomposition screening iterative process stop criterion |
| CN109582913B (en) * | 2018-11-20 | 2023-08-22 | 国家海洋局第一海洋研究所 | Method for empirical mode decomposition screening iteration process termination criterion |
| CN109859174A (en) * | 2019-01-09 | 2019-06-07 | 东莞理工学院 | A kind of OLED defect inspection method based on empirical mode decomposition and regression model |
| CN109799535B (en) * | 2019-03-14 | 2020-10-30 | 中船海洋探测技术研究院有限公司 | Filtering method for full-tensor magnetic gradient positioning detection data |
| CN109799535A (en) * | 2019-03-14 | 2019-05-24 | 中船海洋探测技术研究院有限公司 | A kind of filtering method of full tensor magnetic gradient detection and localization data |
| CN110096673A (en) * | 2019-04-29 | 2019-08-06 | 河北工业大学 | A kind of EMD improved method suitable for signal decomposition |
| CN110096673B (en) * | 2019-04-29 | 2023-03-14 | 河北工业大学 | EMD improvement method suitable for signal decomposition |
| CN110441063A (en) * | 2019-06-12 | 2019-11-12 | 祝思宁 | A kind of method of monitoring, diagnosing large high-speed armature spindle crackle |
| CN111537989A (en) * | 2020-03-25 | 2020-08-14 | 中国电子科技集团公司第二十九研究所 | Method for extracting signal micro Doppler modulation component based on empirical mode decomposition |
| CN111553308A (en) * | 2020-05-11 | 2020-08-18 | 成都亿科康德电气有限公司 | Reconstruction method of partial discharge signal of power transformer |
| CN114936375A (en) * | 2022-06-02 | 2022-08-23 | 北京理工大学 | Asymmetric integrated optical encryption system based on two-dimensional empirical mode decomposition |
Also Published As
| Publication number | Publication date |
|---|---|
| CN107748734B (en) | 2021-08-27 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN107748734A (en) | One kind parsing empirical mode decomposition method | |
| Feng et al. | Recent advances in time–frequency analysis methods for machinery fault diagnosis: A review with application examples | |
| CN103226649B (en) | The integrated noise reconstruct ensemble empirical mode decomposition method of machinery early stage and combined failure | |
| CN102360502B (en) | Automatic baseline correction method | |
| CN105928702B (en) | Variable working condition box bearing method for diagnosing faults based on form PCA | |
| CN103870694B (en) | Empirical mode decomposition denoising method based on revised wavelet threshold value | |
| CN104236905A (en) | Bearing fault diagnosis method | |
| CN109241915B (en) | Intelligent power plant pump fault diagnosis method based on vibration signal stability and non-stationarity judgment and feature discrimination | |
| CN104133404B (en) | A kind of signal processing method and device | |
| CN105354377A (en) | Method for determining fluctuation wind induced vibration load of power transmission tower | |
| CN103441762A (en) | ADC dynamic parameter testing method based on Blackman window three-spectrum-line interpolation | |
| CN103983849B (en) | A kind of Electric Power Harmonic Analysis method of real-time high-precision | |
| CN104006961A (en) | Cycloid bevel gear fault diagnosis method based on empirical mode decomposition and cepstrum | |
| Xie et al. | Fast-varying AM–FM components extraction based on an adaptive STFT | |
| CN104215459A (en) | Bearing fault diagnosis method | |
| CN102636693A (en) | Harmonic analysis algorithm combining fast Fourier transform (FFT) and nonlinear least square | |
| CN102201240A (en) | Harmonic noise excitation model vocoder based on inverse filtering | |
| CN103795411A (en) | SFDR testing method based on five-maximum-sidelobe-damping-window three-spectral-line interpolation | |
| CN105266784A (en) | Oscillography-based blood pressure measuring method and device | |
| CN105510066B (en) | A kind of rotatory mechanical system method for diagnosing faults based on adaptive noise reduction algorithm | |
| Chen et al. | Adaptive scale decomposition and weighted multikernel correntropy for wheelset axle box bearing diagnosis under impact interference | |
| CN103440226A (en) | EMD (Empirical Mode Decomposition) endpoint effect suppression method based on HMM (Hidden Markov Model) correction and neural network extension | |
| CN107944405A (en) | A kind of cubic spline part mean decomposition method based on extreme value point calibration | |
| CN108959739A (en) | A kind of analysis method and device of the pressure fluctuation of hydroenergy storage station transient process | |
| CN107561420A (en) | A kind of cable local discharge signal characteristic vector extracting method based on empirical mode decomposition |
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 |