CN106198010B - A kind of envelope Analysis Method that filtering is decomposed based on local mean value - Google Patents

A kind of envelope Analysis Method that filtering is decomposed based on local mean value Download PDF

Info

Publication number
CN106198010B
CN106198010B CN201610492072.1A CN201610492072A CN106198010B CN 106198010 B CN106198010 B CN 106198010B CN 201610492072 A CN201610492072 A CN 201610492072A CN 106198010 B CN106198010 B CN 106198010B
Authority
CN
China
Prior art keywords
envelope
signal
data
mean value
local mean
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.)
Expired - Fee Related
Application number
CN201610492072.1A
Other languages
Chinese (zh)
Other versions
CN106198010A (en
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.)
Weifang University
Original Assignee
Weifang University
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 Weifang University filed Critical Weifang University
Priority to CN201610492072.1A priority Critical patent/CN106198010B/en
Publication of CN106198010A publication Critical patent/CN106198010A/en
Application granted granted Critical
Publication of CN106198010B publication Critical patent/CN106198010B/en
Expired - Fee Related legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • GPHYSICS
    • G01MEASURING; TESTING
    • G01MTESTING STATIC OR DYNAMIC BALANCE OF MACHINES OR STRUCTURES; TESTING OF STRUCTURES OR APPARATUS, NOT OTHERWISE PROVIDED FOR
    • G01M13/00Testing of machine parts
    • G01M13/04Bearings
    • G01M13/045Acoustic or vibration analysis

Landscapes

  • Physics & Mathematics (AREA)
  • Acoustics & Sound (AREA)
  • General Physics & Mathematics (AREA)
  • Measurement Of Mechanical Vibrations Or Ultrasonic Waves (AREA)
  • Complex Calculations (AREA)

Abstract

The invention discloses a kind of envelope Analysis Methods that filtering is decomposed based on local mean value, this method decomposes original signal first with part mean decomposition method, then the noise component(s) and trend term in decomposition result are excluded using the rearrangement of data and replacement operation, then filtered signal for the first time is analyzed using spectrum kurtosis method again, obtain the centre frequency and bandwidth of optimal filter, then second of filtering is carried out again to filtered signal for the first time using the wave filter, then Envelope Analysis is carried out to second of filtered signal using rational spline iteration smoothed envelope analysis method, the fault type of rotating machinery is finally determined according to envelope spectrum.The present invention is suitable for the complicated rotating machinery fault signal of processing, can accurately determine the fault type of rotating machinery, have good noise immunity and robustness, convenient for engineer application.

Description

A kind of envelope Analysis Method that filtering is decomposed based on local mean value
Technical field
The present invention relates to condition monitoring for rotating machinery and fault diagnosis field, and in particular to one kind is decomposed based on local mean value The envelope Analysis Method of filtering.
Background technology
Envelope Analysis technology is widely used in the fault diagnosis of gear and rolling bearing.Existing Envelope Analysis technology has Three defects below:1. existing Envelope Analysis technology directly analyzes original signal or only to original Signal analyzed again, therefore existing method is easily done by noise, trend and other ingredients after simply filtering It disturbs, it is relatively low so as to cause the analysis precision of the prior art;2. existing Envelope Analysis technology is based on Hilbert is converted, And the signal that Hilbert transformation requirements are analyzed must be the narrow band signal of simple component, otherwise the frequency modulating section of signal will The amplitude envelope analysis result of signal is polluted, but signal to be analyzed at present does not meet the item of simple component and narrowband strictly Part may result in the prior art and erroneous judgement problem be susceptible to due to precision is not high in this way;3. the envelope spectrum obtained by conventional method There are end effects.
Invention content
The problem to be solved in the present invention is for a kind of above insufficient, envelope point that filtering is decomposed based on local mean value of proposition Analysis method after envelope Analysis Method using the present invention, has that analysis result accuracy and accuracy are high, and can accurately detect The advantages of going out rotating machinery fault type.
For solution more than technical problem, the technical solution that the present invention takes is as follows:One kind decomposes filtering based on local mean value Envelope Analysis Method, which is characterized in that include the following steps:
Step 1:The vibration signal x of rotating machinery is measured with sample frequency fs using acceleration transducer(k), (k=1, 2, …,N), N is the length of sampled signal;
Step 2:Using local mean value decomposition algorithm by signal x(k)The sum of n component and trend term are resolved into, i.e.,, wherein, ci(k)Represent i-th of the component obtained by local mean value decomposition algorithm, rn(k)It represents The trend term obtained by local mean value decomposition algorithm;
Step 3:To ci(k)It performs reordering operations and substitutes and operate, it is rearranged to operate obtained data ci shuffle(k)Table Show, data c is obtained after substituting operationi FTran(k)It represents;
Step 4:To ci(k)、ci shuffle(k)And ci FTran(k)Multi-fractal is performed respectively removes trend fluction analysis (Multifractal Detrended Fluctuation Analysis, MFDFA), obtain generalized Hurst index curve, ci (k)Generalized Hurst index curve Hi(q)It represents;ci shuffle(k)Generalized Hurst index curve Hi shuffle(q)Table Show;ci FTran(k)Generalized Hurst index curve Hi FTran(q)It represents;
Step 5:If Hi(q)With Hi shuffle(q)Or Hi(q)With Hi FTran(q)Between relative error be less than 5% or Person Hi(q) 、Hi shuffle(q)And Hi FTran(q)Three does not change with q, then abandons corresponding ci(k)Component;
Step 6:To remaining ci(k)Component is summed, and by this and to be denoted as signal rearranged and substitute filtered result xf1 (k);
Step 7:To xf1(k)Spectrum kurtosis analysis is performed, the centre frequency f corresponding to signal kurtosis maximum is obtained0And band Wide B;
Step 8:According to centre frequency f0With bandwidth B to xf1(k)Bandpass filtering is carried out, obtains xf2(k);
Step 9:To signal xf2(k)The analysis of rational spline iteration smoothed envelope is performed, obtains signal envelope eov(k);
Step 10:To obtained signal envelope eov(k)It performs discrete Fourier transform and obtains envelope spectrum, according to envelope spectrum Characteristic frequency judges the fault type of machine.
A kind of prioritization scheme, local mean value decomposition algorithm includes the following steps in the step 2:
1)Calculate local mean value function:Determine original signalx(k) all Local Extremumn i , calculate two neighboring pole Value pointn i Withn i+1Average valuem i , i.e.,
By the average value of all two neighboring extreme pointsm i It is connected, is then carried out using rolling average method smooth with broken line Processing, obtains local mean value functionm 11(k) ;In this example, the smooth step-length in rolling average method is set as 5;
2)Estimate the envelope value of signal:Using Local Extremumn i Calculate envelope estimated valuea i
Equally, by all two neighboring envelope estimated valuesa i It is connected with broken line, is then put down using rolling average method Sliding processing, obtains envelope estimation functiona 11(k);
3)By local mean value functionm 11(k) from original signalx(k) in separate, obtain
h 11(k)= x(k)- m 11(k);
4)Withh 11(k) divided by envelope estimation functiona 11(k) so as to righth 11(k) be demodulated, it obtains
It is desirable thats 11(k) be a pure FM signal, i.e. its envelope estimation functiona 12(k) meeta 12(k)=1;Ifs 11(k) condition is not satisfied, then wills 11(k) as initial data repeatedly more than iterative process m times, until obtaining a pure tune Frequency signals 1m(k), i.e.,s 1m(k) satisfaction -1≤s 1m(k)≤1, its envelope estimation functiona 1(m+1)(k) meeta 1(m+1)(k)=1; Therefore have
In formula
The condition of iteration ends is
In practical applications, a variation Δ can be set, when meeting 1- Δs≤a1m(k)≤1+ Δs when, iteration is whole Only;Variation Δ=0.01 in this example;
5)All envelope estimation functions multiplication generated in iterative process can be obtained envelope signal
6)By envelope signala 1 (k) and pure FM signals 1m (k) can be obtained by being multiplied by 1st point of original signal Amount
c 1 (k)=a 1(k)s 1m(k);
7)By the 1st componentc 1 (k) from original signalx(k) in separate, obtain a new signalr 1 (k), it willr 1 (k) as new data repetition above step, cyclenIt is secondary, untilr n (k) for until monotonic function
So far, by original signalx(k) be decomposed intonA component and a monotonic functionr n (kThe sum of), i.e.,
Further, data rearrangement operation includes the following steps in the step 3:
Upset component c at randomi(k)Put in order.
Further, the operation of data replacement includes the following steps in the step 3:
1)To component ci(k)Discrete Fourier transform is performed, obtains component ci(k)Phase;
2)It is located at the pseudo- independent same distribution number in (- π, π) section with one group to replace component ci(k)Original phase;
3)Inverse discrete Fourier transform is performed to the frequency domain data after phase substitutes and obtains data ci IFFT(k), ask for Data ci IFFT(k)Real part.
Further, MFDFA methods include the following steps in the step 4:
1)Constructionx(k)(k=1,2,…,N) profileY(i):
x(k) represent c in step 4 described in claim 1i(k)Or ci shuffle(k)Or ci FTran(k);
2)By signal profileY(i) be divided into it is nonoverlappingN s Segment length issData, due to data lengthNIt generally can not be whole It removess, so the remaining one piece of data of meeting cannot utilize;
In order to make full use of the length of data, then it is obtained with identical length segmentation, such one from the negative direction of data 2N s Segment data;
3)Using polynomial trend of the least square fitting per segment data, the variance per segment data is then calculated:
y v (i) for fitting thevThe trend of segment data, if the polynomial trend of fitting ismRank then remembers that this goes trend process For(MF-)DFAm;In this example, m=1;
4)Calculate theqThe average value of rank wave function:
5)Ifx(k) there are self-similarity characteristics, thenqThe average value of rank wave functionF q (s) and time scalesBetween There are power law relations:
F q (s)~s H(q)
WhenqWhen=0, step 4)In formula diverging, at this momentH(0) by logarithmic mean process defined in following formula come really It is fixed:
6)To step 5)In formula both sides take the logarithm can obtain ln [F q (s)]=H(q)ln(s)+c (cFor constant), thus The slope of straight line can be obtainedH(q)。
Further, the spectrum kurtosis method in the step 7 includes the following steps:
1)Constructing a cutoff frequency isf c =0.125+εLow-pass filterh(n);ε>0, in this examplef c =0.3;
2)It is based onh(n) construction passband be [0;0.25] quasi- low-pass filterh 0(n) and passband be [0.25;0.5] Quasi- high-pass filterh 1 (n),
3)SignalThroughh 0(n)、h 1 (n) filter and resolve into low frequency part after down-sampledAnd high frequency section, the down-sampled factor is 2, then the shaping filter tree after successive ignition filters, kth layer have 2 k A frequency band, whereinIt represents the in wave filter treekOn layeriThe output signal of a wave filter,i=0,…,2k- 1,0≤k≤K-1, this example Middle K=8;c 0(n) represent x in step 7 described in claim 1f1(k);
4)In decomposition treekOn layeriThe centre frequency of a wave filterf ki And bandwidthB k Respectively
f ki =(i+2-1)2-k-1
B k =2-k-1
5)Calculate each filter results(i=0,…,2k- 1) kurtosis
6)All spectrum kurtosis are summarized, obtain the total spectrum kurtosis of signal.
Further, the rational spline iteration smoothed envelope analysis method in the step 9 includes the following steps:
1)Calculate signalz(k) Jue Dui Zhi ∣z(k) the local extremum of ∣;In the 1st iteration,z(k) represent claim X in 1 step 9f2(k);
2)Envelope eov is obtained using rational spline curve fitting Local Extremum1(k);
3)It is rightz(k) be normalized to obtain
4)2nd iteration:z 1(k) new data is re-used as, repeat above-mentioned steps 1)~3), obtain
5)Ith iteration:z i-1(k) new data is re-used as, repeat above-mentioned steps 1)~ 3) it, obtains
If 6)nWhat secondary iteration obtainedz n (k) amplitude be less than or equal to 1, then iterative process stop, finally obtaining Signalz(k) envelope be
The present invention is using above technical scheme, and compared with prior art, the present invention has the following advantages:
1) original signal is decomposed using local mean value decomposition, is then excluded using the rearrangement of data and replacement operation Noise and trend component therein, the only useful component in stick signal component, so as to avoid noise and trend component pair The influence of Envelope Analysis result, analysis result accuracy and accuracy are high.
2) signal envelope and frequency modulating section are kept completely separate using rational spline iteration smoothed envelope analysis method, energy Influence of the frequency modulating section to signal envelope analysis result is enough avoided, so as to improve the precision of Envelope Analysis.
3) fault type of rotating machinery can be accurately detected.
4) envelope spectrum obtained by conventional method is there are end effect, and can be avoided by the envelope spectrum that the present invention obtains End effect.
The present invention will be further described with reference to the accompanying drawings and examples.
Description of the drawings
Attached drawing 1 is the flow chart of the method for the present invention in the embodiment of the present invention.
Attached drawing 2 is to carry out preliminary exposition to signal using low-pass filter and high-pass filter in the embodiment of the present invention to show It is intended to.
Attached drawing 3 is the schematic diagram for quickly calculating spectrum kurtosis in the embodiment of the present invention using tree-shaped filter construction.
Attached drawing 4 is the bearing vibration signal for having in the embodiment of the present invention inner ring failure.
Attached drawing 5 is to inner ring faulty bearing vibration signal in the embodiment of the present invention using traditional envelope Analysis Method Analysis result.
Attached drawing 6 is the present invention in the embodiment of the present invention to the analysis result of inner ring faulty bearing vibration signal.
Attached drawing 7 is the bearing vibration signal for having in the embodiment of the present invention outer ring failure.
Attached drawing 8 is to outer ring faulty bearing vibration signal in the embodiment of the present invention using traditional envelope Analysis Method Analysis result.
Attached drawing 9 is the present invention in the embodiment of the present invention to the analysis result of outer ring faulty bearing vibration signal.
Specific embodiment
Embodiment, as shown in Figure 1, Figure 2, Figure 3 shows, a kind of envelope Analysis Method that filtering is decomposed based on local mean value, including Following steps:
Step 1:The vibration signal x of rotating machinery is measured with sample frequency fs using acceleration transducer(k), (k=1, 2, …,N), N is the length of sampled signal;
Step 2:Using local mean value decomposition algorithm by signal x(k)The sum of n component and trend term are resolved into, i.e.,, wherein, ci(k)Represent i-th of the component obtained by local mean value decomposition algorithm, rn(k)It represents The trend term obtained by local mean value decomposition algorithm;
Step 3:To ci(k)It performs reordering operations and substitutes and operate, it is rearranged to operate obtained data ci shuffle(k)Table Show, data c is obtained after substituting operationi FTran(k)It represents;
Step 4:To ci(k)、ci shuffle(k)And ci FTran(k)Multi-fractal is performed respectively removes trend fluction analysis (Multifractal Detrended Fluctuation Analysis, MFDFA), obtain generalized Hurst index curve, ci (k)Generalized Hurst index curve Hi(q)It represents;ci shuffle(k)Generalized Hurst index curve Hi shuffle(q)Table Show;ci FTran(k)Generalized Hurst index curve Hi FTran(q)It represents;
Step 5:If Hi(q)With Hi shuffle(q)Or Hi(q)With Hi FTran(q)Between relative error be less than 5% or Person Hi(q) 、Hi shuffle(q)And Hi FTran(q)Three does not change with q, then abandons corresponding ci(k)Component;
Step 6:To remaining ci(k)Component is summed, and by this and to be denoted as signal rearranged and substitute filtered result xf1 (k);
Step 7:To xf1(k)Spectrum kurtosis analysis is performed, the centre frequency f corresponding to signal kurtosis maximum is obtained0And band Wide B;
Step 8:According to centre frequency f0With bandwidth B to xf1(k)Bandpass filtering is carried out, obtains xf2(k);
Step 9:To signal xf2(k)The analysis of rational spline iteration smoothed envelope is performed, obtains signal envelope eov(k);
Step 10:To obtained signal envelope eov(k)It performs discrete Fourier transform and obtains envelope spectrum, according to envelope spectrum Characteristic frequency judges the fault type of machine.
Local mean value decomposition algorithm includes the following steps in step 2:
1)Calculate local mean value function:Determine original signalx(k) all Local Extremumn i , calculate two neighboring pole Value pointn i Withn i+1Average valuem i , i.e.,
By the average value of all two neighboring extreme pointsm i It is connected, is then carried out using rolling average method smooth with broken line Processing, obtains local mean value functionm 11(k) ;In this example, the smooth step-length in rolling average method is set as 5;
2)Estimate the envelope value of signal:Using Local Extremumn i Calculate envelope estimated valuea i
Equally, by all two neighboring envelope estimated valuesa i It is connected with broken line, is then put down using rolling average method Sliding processing, obtains envelope estimation functiona 11(k);
3)By local mean value functionm 11(k) from original signalx(k) in separate, obtain
h 11(k)= x(k)- m 11(k);
4)Withh 11(k) divided by envelope estimation functiona 11(k) so as to righth 11(k) be demodulated, it obtains
It is desirable thats 11(k) be a pure FM signal, i.e. its envelope estimation functiona 12(k) meeta 12(k)=1;Ifs 11(k) condition is not satisfied, then wills 11(k) as initial data repeatedly more than iterative process m times, until obtaining a pure tune Frequency signals 1m(k), i.e.,s 1m(k) satisfaction -1≤s 1m(k)≤1, its envelope estimation functiona 1(m+1)(k) meeta 1(m+1)(k)=1; Therefore have
In formula
The condition of iteration ends is
In practical applications, a variation Δ can be set, when meeting 1- Δs≤a1m(k)≤1+ Δs when, iteration is whole Only;Variation Δ=0.01 in this example;
5)All envelope estimation functions multiplication generated in iterative process can be obtained envelope signal
6)By envelope signala 1 (k) and pure FM signals 1m (k) can be obtained by being multiplied by 1st point of original signal Amount
c 1(k)=a 1(k)s 1m(k);
7)By the 1st componentc 1 (k) from original signalx(k) in separate, obtain a new signalr 1 (k), it willr 1 (k) as new data repetition above step, cyclenIt is secondary, untilr n (k) for until monotonic function
So far, by original signalx(k) be decomposed intonA component and a monotonic functionr n (kThe sum of), i.e.,
Data rearrangement operation includes the following steps in step 3:
Upset component c at randomi(k)Put in order.
Data substitute operation and include the following steps in step 3:
1)To component ci(k)Discrete Fourier transform is performed, obtains component ci(k)Phase;
2)It is located at the pseudo- independent same distribution number in (- π, π) section with one group to replace component ci(k)Original phase;
3)Inverse discrete Fourier transform is performed to the frequency domain data after phase substitutes and obtains data ci IFFT(k), ask for Data ci IFFT(k)Real part.
MFDFA methods include the following steps in step 4:
1)Constructionx(k)(k=1,2,…,N) profileY(i):
x(k) represent c in step 4 described in claim 1i(k)Or ci shuffle(k)Or ci FTran(k);
2)By signal profileY(i) be divided into it is nonoverlappingN s Segment length issData, due to data lengthNIt generally can not be whole It removess, so the remaining one piece of data of meeting cannot utilize;
In order to make full use of the length of data, then it is obtained with identical length segmentation, such one from the negative direction of data 2N s Segment data;
3)Using polynomial trend of the least square fitting per segment data, the variance per segment data is then calculated:
y v (i) for fitting thevThe trend of segment data, if the polynomial trend of fitting ismRank then remembers that this goes trend process For(MF-)DFAm;In this example, m=1;
4)Calculate theqThe average value of rank wave function:
5)Ifx(k) there are self-similarity characteristics, thenqThe average value of rank wave functionF q (s) and time scalesBetween There are power law relations:
F q (s)~s H(q)
WhenqWhen=0, step 4)In formula diverging, at this momentH(0) by logarithmic mean process defined in following formula come really It is fixed:
6)To step 5)In formula both sides take the logarithm can obtain ln [F q (s)]=H(q)ln(s)+c (cFor constant), thus The slope of straight line can be obtainedH(q)。
Spectrum kurtosis method in step 7 includes the following steps:
1)Constructing a cutoff frequency isf c =0.125+εLow-pass filterh(n);ε>0, in this examplef c =0.3;
2)It is based onh(n) construction passband be [0;0.25] quasi- low-pass filterh 0(n) and passband be [0.25;0.5] Quasi- high-pass filterh 1 (n),
3)SignalThroughh 0(n)、h 1 (n) filter and resolve into low frequency part after down-sampledAnd high frequency section, the down-sampled factor is 2, then the shaping filter tree after successive ignition filters, kth layer have 2 k A frequency band, whereinIt represents the in wave filter treekOn layeriThe output signal of a wave filter,i=0,…,2k- 1,0≤k≤K-1, this example Middle K=8;c 0(n) represent x in step 7 described in claim 1f1(k);
4)In decomposition treekOn layeriThe centre frequency of a wave filterf ki And bandwidthB k Respectively
f ki =(i+2-1)2-k-1
B k =2-k-1
5)Calculate each filter results(i=0,…,2k- 1) kurtosis
6)All spectrum kurtosis are summarized, obtain the total spectrum kurtosis of signal.
Rational spline iteration smoothed envelope analysis method in step 9 includes the following steps:
1)Calculate signalz(k) Jue Dui Zhi ∣z(k) the local extremum of ∣;In the 1st iteration,z(k) represent claim X in 1 step 9f2(k);
2)Envelope eov is obtained using rational spline curve fitting Local Extremum1(k);
3)It is rightz(k) be normalized to obtain
4)2nd iteration:z 1(k) new data is re-used as, repeat above-mentioned steps 1)~3), obtain
5)Ith iteration:z i-1(k) new data is re-used as, repeat above-mentioned steps 1)~ 3) it, obtains
If 6)nWhat secondary iteration obtainedz n (k) amplitude be less than or equal to 1, then iterative process stop, finally obtaining Signalz(k) envelope be
Experiment 1, tests the performance of algorithm of the present invention using the bearing vibration data with inner ring failure Card.
Bearing used in experiment be 6205-2RS JEM SKF, using electric discharge machining method on bearing inner race working depth The groove for being 0.3556mm for 0.2794mm, width simulates bearing inner race failure, this experiment load is about 0.7457kW, driving It is about 29.5Hz that motor, which turns frequency, and bearing inner race fault characteristic frequency is about 160Hz, sample frequency 4.8KHz, during signal sampling A length of 1s.
Collected inner ring fault-signal is as shown in Figure 4.
Signal shown in Fig. 4 is analyzed using traditional envelope Analysis Method first, obtained analysis result such as Fig. 5 It is shown.From fig. 5, it can be seen that the fault signature of bearing is blanked completely, therefore traditional envelope Analysis Method cannot be effectively Extract the fault signature of bearing;In addition, from fig. 5, it can be seen that the left end point of envelope spectrum has abnormal high level, this explanation is by tradition There is end effects for the envelope spectrum that method obtains.
Signal shown in Fig. 4 is analyzed using method proposed by the invention, obtained analysis result such as Fig. 6 institutes Show.From fig. 6, it can be seen that the spectral line corresponding to 160Hz and 320Hz is apparently higher than other spectral lines, the two frequencies correspond to respectively 1 frequency multiplication and 2 frequencys multiplication of bearing inner race fault characteristic frequency may determine that bearing has inner ring failure accordingly;It can from Fig. 6 Go out, the envelope spectrum obtained by the present invention does not have end effect.
Show that the present invention is capable of reliable recognition in the case where loading and failure dimensional depth being constant through many experiments Minimum inner ring failure dimension width is about 0.19mm, and conventional method is capable of the minimum inner ring failure dimension width of reliable recognition about For 0.53mm, precision improves 64.2%.
Experiment 2, tests the performance of algorithm of the present invention using the bearing vibration data with outer ring failure Card.
Bearing used in experiment be 6205-2RS JEM SKF, using electric discharge machining method on bearing outer ring working depth The groove for being 0.5334mm for 0.2794mm, width simulates bearing outer ring failure, this experiment load is about 2.237 kW, driving It is about 28.7Hz that motor, which turns frequency, and bearing outer ring fault characteristic frequency is about 103Hz, sample frequency 4.8KHz, during signal sampling A length of 1s.
Collected outer ring fault-signal is as shown in Figure 7.
Signal shown in Fig. 7 is analyzed using traditional envelope Analysis Method first, obtained analysis result such as Fig. 8 It is shown.From figure 8, it is seen that the fault signature of bearing is blanked completely, therefore traditional envelope Analysis Method cannot be effectively Extract the fault signature of bearing;In addition, from figure 8, it is seen that the left end point of envelope spectrum there is abnormal high level, this explanation is by passing There is end effects for the envelope spectrum that system method obtains.
Signal shown in Fig. 7 is analyzed using method proposed by the invention, obtained analysis result such as Fig. 9 institutes Show.From fig. 9, it can be seen that the spectral line corresponding to 103Hz and 206Hz is apparently higher than other spectral lines, the two frequencies correspond to respectively 1 frequency multiplication and 2 frequencys multiplication of bearing outer ring fault characteristic frequency may determine that bearing has outer ring failure accordingly;It can from Fig. 9 Go out, the envelope spectrum obtained by the present invention does not have end effect.
Show that the present invention is capable of reliable recognition in the case where loading and failure dimensional depth being constant through many experiments Minimum outer ring failure dimension width is about 0.29mm, and conventional method is capable of the minimum outer ring failure dimension width of reliable recognition about For 0.68mm, precision improves 57.4%.
According to result of the test, think after analysis:
1) traditional envelope Analysis Method directly carries out original signal Envelope Analysis or to merely through simple process Original signal afterwards carries out Envelope Analysis, different from traditional envelope Analysis Method, and the present invention is decomposed first with local mean value Original signal is decomposed, then excludes noise and trend component therein using the rearrangement of data and replacement operation, only Useful component in stick signal component so as to avoid the influence of noise and trend component to Envelope Analysis result, improves Accuracy and precision.
2) traditional envelope Analysis Method is based on Hilbert is converted, and the letter that Hilbert transformation requirements are analyzed Number must be the narrow band signal of simple component, otherwise the frequency modulating section of signal will pollute signal Envelope Analysis as a result, but It is the condition that signal to be analyzed does not meet simple component and narrowband strictly at present, may result in the prior art in this way because of precision not High and be susceptible to erroneous judgement problem, different from traditional envelope Analysis Method, the present invention is divided using rational spline iteration smoothed envelope Signal envelope and frequency modulating section are kept completely separate by analysis method, can avoid frequency modulating section to signal envelope analysis result Influence, so as to improve the precision of Envelope Analysis.
3) fault type of rotating machinery can be accurately detected.
4) envelope spectrum obtained by conventional method is there are end effect, and can be avoided by the envelope spectrum that the present invention obtains End effect.
5)Each step effect:
1) step:Acquire vibration signal;
2) step:Original signal is resolved into the form of different component sums, some of which component corresponds to noise and trend term, Some components correspond to useful signal;
3) ~ 5) step:Is performed by reordering operations and is substituted for the signal that above-mentioned decomposition obtains and is operated, rejects noise therein point Amount and trend term only retain useful signal;
6) step:Remaining useful signal is summed, using this and as signal it is rearranged and substitute filtered result xf1 (k);
7) step:To filtered signal xf1(k) spectrum kurtosis analysis is performed, corresponding center at signal maximum kurtosis is obtained Frequency f0And bandwidth B;
8) step:According to centre frequency f0With bandwidth B to xf1(k) bandpass filtering is carried out, obtains signal xf2(k);
9) step:Calculate signal xf2(k) envelope eov (k);
10) step:Discrete Fourier transform is performed to eov (k) and obtains envelope spectrum, the failure of bearing is judged according to envelope spectrum Type.
One skilled in the art would recognize that above-mentioned specific embodiment is only exemplary, it is to make ability Field technique personnel can be better understood from the content of present invention, should not be understood as limiting the scope of the invention, as long as According to technical solution of the present invention improvements introduced, protection scope of the present invention is each fallen within.

Claims (7)

1. a kind of envelope Analysis Method that filtering is decomposed based on local mean value, which is characterized in that include the following steps:
Step 1:The vibration signal x of rotating machinery is measured with sample frequency fs using acceleration transducer(k), k=1, 2, …, N, N are the length of sampled signal;
Step 2:Using local mean value decomposition algorithm by vibration signal x(k)The sum of n component and trend term are resolved into, i.e.,, wherein, ci(k)Represent i-th point obtained by local mean value decomposition algorithm Amount, rn(k)Represent the trend term obtained by local mean value decomposition algorithm;
Step 3:To ci(k)It performs reordering operations and substitutes and operate, it is rearranged to operate obtained data ci shuffle(k)It represents, Data c is obtained after substituting operationi FTran(k)It represents;
Step 4:To ci(k)、ci shuffle(k)And ci FTran(k)Multi-fractal is performed respectively removes trend fluction analysis (Multifractal Detrended Fluctuation Analysis, MFDFA), obtain generalized Hurst index curve, ci (k)Generalized Hurst index curve Hi(q)It represents;ci shuffle(k)Generalized Hurst index curve Hi shuffle(q)Table Show;ci FTran(k)Generalized Hurst index curve Hi FTran(q)It represents;
Step 5:If Hi(q)With Hi shuffle(q)Or Hi(q)With Hi FTran(q)Between relative error be less than 5% or Hi (q) 、Hi shuffle(q)And Hi FTran(q)Three does not change with q, then abandons corresponding ci(k)Component;
Step 6:To remaining ci(k)Component is summed, and by this and to be denoted as signal rearranged and substitute filtered result xf1(k);
Step 7:To xf1(k)Spectrum kurtosis analysis is performed, the centre frequency f corresponding to signal kurtosis maximum is obtained0And bandwidth B;
Step 8:According to centre frequency f0With bandwidth B to xf1(k)Bandpass filtering is carried out, obtains xf2(k);
Step 9:To signal xf2(k)The analysis of rational spline iteration smoothed envelope is performed, obtains signal envelope eov(k);
Step 10:To obtained signal envelope eov(k)It performs discrete Fourier transform and obtains envelope spectrum, according to envelope spectrum signature Frequency judges the fault type of machine.
A kind of 2. envelope Analysis Method that filtering is decomposed based on local mean value according to claim 1, which is characterized in that institute Local mean value decomposition algorithm in step 2 is stated to include the following steps:
1)Calculate local mean value function:Determine vibration signalx(k) all Local Extremumn i , calculate two neighboring extreme pointn i Withn i+1Average valuem i , i.e.,
By the average value of all two neighboring extreme pointsm i It is connected with broken line, is then smoothed using rolling average method, Obtain local mean value functionm 11(k) ;In this example, the smooth step-length in rolling average method is set as 5;
2)Estimate the envelope value of signal:Using Local Extremumn i Calculate envelope estimated valuea i
Equally, by all two neighboring envelope estimated valuesa i It is connected with broken line, is then smoothly located using rolling average method Reason, obtains envelope estimation functiona 11(k);
3)By local mean value functionm 11(k) from vibration signalx(k) in separate, obtain
h 11(k)= x(k)- m 11(k);
4)Withh 11(k) divided by envelope estimation functiona 11(k) so as to righth 11(k) be demodulated, it obtains
It is desirable thats 11(k) be a pure FM signal, i.e. its envelope estimation functiona 12(k) meeta 12(k)=1;Ifs 11 (k) condition is not satisfied, then wills 11(k) as initial data repeatedly more than iterative process m times, until obtaining a pure frequency modulation letter Numbers 1m(k), i.e.,s 1m(k) satisfaction -1≤s 1m(k)≤1, its envelope estimation functiona 1(m+1)(k) meeta 1(m+1)(k)=1;Therefore Have
In formula
The condition of iteration ends is
In practical applications, a variation Δ is set, when meeting 1- Δs≤a1m(k)≤1+ Δs when, iteration ends;In this example Variation Δ=0.01;
5)All envelope estimation functions generated in iterative process are mutually obtained envelope signal at convenience
6)By envelope signala 1 (k) and pure FM signals 1m (k) the 1st component of original signal is mutually obtained at convenience
c 1 (k)= a 1 (k)s 1m(k);
7)By the 1st componentc 1 (k) from vibration signalx(k) in separate, obtain a new signalr 1 (k), it willr 1 (k) as new data repetition above step, cyclenIt is secondary, untilr n (k) for until monotonic function
So far, by vibration signalx(k) be decomposed intonA component and a monotonic functionr n (kThe sum of), i.e.,
A kind of 3. envelope Analysis Method that filtering is decomposed based on local mean value according to claim 1, which is characterized in that institute Data rearrangement operation in step 3 is stated to include the following steps:
Upset component c at randomi(k)Put in order.
4. a kind of envelope Analysis Method that filtering is decomposed based on local mean value according to claim 1, it is characterised in that:Institute It states data in step 3 and substitutes to operate and include the following steps:
1)To component ci(k)Discrete Fourier transform is performed, obtains component ci(k)Phase;
2)It is located at the pseudo- independent same distribution number in (- π, π) section with one group to replace component ci(k)Original phase;
3)Inverse discrete Fourier transform is performed to the frequency domain data after phase substitutes and obtains data ci IFFT(k), ask for data ci IFFT(k)Real part.
5. a kind of envelope Analysis Method that filtering is decomposed based on local mean value according to claim 1, it is characterised in that:Institute MFDFA methods in step 4 are stated to include the following steps:
1)Construct xy(k) profileY(i), k=1,2 ..., N:
,
,
xy(k) c in step 4 described in claim 1 is representedi(k)Or ci shuffle(k)Or ci FTran(k);
2)By signal profileY(i) be divided into it is nonoverlappingN s Segment length issData, due to data lengthNIt generally can not divide exactlys, So the remaining one piece of data of meeting cannot utilize;
In order to make full use of the length of data, then from the negative direction of data 2 are obtained with identical length segmentation, such oneN s Section Data;
3)Using polynomial trend of the least square fitting per segment data, the variance per segment data is then calculated:
y v (i) for fitting thevThe trend of segment data, if the polynomial trend of fitting ismRank then remembers that this goes the trend process to be (MF-)DFAm;In this example, m=1;
4)Calculate theqThe average value of rank wave function:
5)If xy(k) there are self-similarity characteristics, thenqThe average value of rank wave functionF q (s) and time scalesBetween exist Power law relation:
F q (s)~s H(q)
WhenqWhen=0, step 4)In formula diverging, at this momentH(0) it is determined by logarithmic mean process defined in following formula:
6)To step 5)In formula both sides take the logarithm can obtain ln [F q (s)]=H(q)ln(s)+c,cFor constant, thus to obtain straight line SlopeH(q)。
6. a kind of envelope Analysis Method that filtering is decomposed based on local mean value according to claim 1, it is characterised in that:Institute The spectrum kurtosis method stated in step 7 includes the following steps:
1)Constructing a cutoff frequency isf c =0.125+εLow-pass filterh(n);ε>0, in this examplef c =0.3;
2)It is based onh(n) construction passband be [0;0.25] quasi- low-pass filterh 0(n) and passband be [0.25;0.5] standard High-pass filterh 1 (n),
3)SignalThroughh 0(n)、h 1 (n) filter and resolve into low frequency part after down-sampledAnd high frequency section, the down-sampled factor is 2, then the shaping filter tree after successive ignition filters, kth layer have 2 k A frequency band, whereinIt represents the in wave filter treekOn layeriThe output signal of a wave filter,i=0,…,2k- 1,0≤k≤K-1, this K=8 in example;c 0(n) represent x in step 7 described in claim 1f1(k);
4)In decomposition treekOn layeriThe centre frequency of a wave filterf ki And bandwidthB k Respectively
f ki =(i+2-1)2-k-1
B k =2-k-1
5)Calculate each filter results i=0,…,2k- 1 kurtosis
6)All spectrum kurtosis are summarized, obtain the total spectrum kurtosis of signal.
A kind of 7. envelope Analysis Method that filtering is decomposed based on local mean value according to claim 1, which is characterized in that institute The rational spline iteration smoothed envelope analysis method stated in step 9 includes the following steps:
1)Calculate signalz(k) Jue Dui Zhi ∣z(k) the local extremum of ∣;In the 1st iteration,z(k) represent claim 1 institute State x in step 9f2(k);
2)Envelope eov is obtained using rational spline curve fitting Local Extremum1(k);
3)It is rightz(k) be normalized to obtain
4)2nd iteration:z 1(k) new data is re-used as, repeat above-mentioned steps 1)~3), obtain
5)Ith iteration:z i-1(k) new data is re-used as, repeat above-mentioned steps 1)~ 3) it, obtains
If 6)nWhat secondary iteration obtainedz n (k) amplitude be less than or equal to 1, then iterative process stop, finally obtaining signalz (k) envelope be
CN201610492072.1A 2016-06-29 2016-06-29 A kind of envelope Analysis Method that filtering is decomposed based on local mean value Expired - Fee Related CN106198010B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201610492072.1A CN106198010B (en) 2016-06-29 2016-06-29 A kind of envelope Analysis Method that filtering is decomposed based on local mean value

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201610492072.1A CN106198010B (en) 2016-06-29 2016-06-29 A kind of envelope Analysis Method that filtering is decomposed based on local mean value

Publications (2)

Publication Number Publication Date
CN106198010A CN106198010A (en) 2016-12-07
CN106198010B true CN106198010B (en) 2018-06-22

Family

ID=57462462

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201610492072.1A Expired - Fee Related CN106198010B (en) 2016-06-29 2016-06-29 A kind of envelope Analysis Method that filtering is decomposed based on local mean value

Country Status (1)

Country Link
CN (1) CN106198010B (en)

Families Citing this family (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN107389343B (en) * 2017-07-24 2019-07-30 潍坊学院 The instantaneous Frequency Estimation method examined based on Roberts operator and t
CN108489719B (en) * 2018-04-09 2019-11-05 常州湖南大学机械装备研究院 A kind of combined failure of rotating machinery diagnostic method based on the unusual spectral factorization of G-P
CN113505703A (en) * 2021-07-13 2021-10-15 天津工业大学 Spectral signal denoising method for uneven noise distribution

Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
EP2581725A2 (en) * 2011-10-13 2013-04-17 General Electric Company Methods and systems for automatic rolling-element bearing fault detection
CN105004523A (en) * 2015-08-04 2015-10-28 潍坊学院 Rolling bearing state monitoring method based on weighted similarity measure
CN105067262A (en) * 2015-08-04 2015-11-18 潍坊学院 Rolling bearing state monitoring method
CN105092239A (en) * 2014-05-09 2015-11-25 潍坊学院 Method for detecting early stage fault of gear
CN105588717A (en) * 2015-12-10 2016-05-18 潍坊学院 Gearbox fault diagnosis method

Patent Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
EP2581725A2 (en) * 2011-10-13 2013-04-17 General Electric Company Methods and systems for automatic rolling-element bearing fault detection
CN105092239A (en) * 2014-05-09 2015-11-25 潍坊学院 Method for detecting early stage fault of gear
CN105004523A (en) * 2015-08-04 2015-10-28 潍坊学院 Rolling bearing state monitoring method based on weighted similarity measure
CN105067262A (en) * 2015-08-04 2015-11-18 潍坊学院 Rolling bearing state monitoring method
CN105588717A (en) * 2015-12-10 2016-05-18 潍坊学院 Gearbox fault diagnosis method

Also Published As

Publication number Publication date
CN106198010A (en) 2016-12-07

Similar Documents

Publication Publication Date Title
CN106198015B (en) A kind of VMD of rolling bearing, spectrum kurtosis and smooth iteration envelope Analysis Method
CN106096313B (en) A kind of envelope Analysis Method based on unusual spectral factorization and spectrum kurtosis
CN106153339B (en) A kind of envelope Analysis Method based on the filtering of variation Mode Decomposition
CN106198013B (en) A kind of envelope Analysis Method based on empirical mode decomposition filtering
CN106096200B (en) A kind of envelope Analysis Method based on wavelet decomposition and spectrum kurtosis
CN106168538B (en) A kind of ITD of rolling bearing, spectrum kurtosis and smooth iteration envelope Analysis Method
CN106096198B (en) A kind of envelope Analysis Method based on variation Mode Decomposition and spectrum kurtosis
CN106198009B (en) A kind of EMD of rolling bearing, spectrum kurtosis and smooth iteration envelope Analysis Method
CN105954031B (en) A kind of envelope Analysis Method based on unusual spectral factorization filtering
CN106096199B (en) A kind of WT of rolling bearing, spectrum kurtosis and smooth iteration envelope Analysis Method
CN106053069B (en) A kind of SSD of rolling bearing, spectrum kurtosis and smooth iteration envelope Analysis Method
CN106198010B (en) A kind of envelope Analysis Method that filtering is decomposed based on local mean value
CN106198012B (en) A kind of envelope Analysis Method for decomposing and composing kurtosis based on local mean value
CN106153333B (en) A kind of envelope Analysis Method based on wavelet decomposition filtering
CN106053060B (en) A kind of envelope Analysis Method that filtering is decomposed based on nonlinear model
CN106198017B (en) A kind of LMD of rolling bearing, spectrum kurtosis and smooth iteration envelope Analysis Method
CN106198014B (en) A kind of envelope Analysis Method based on empirical mode decomposition and spectrum kurtosis
CN106198018B (en) A kind of EEMD of rotating machinery and smooth iteration envelope Analysis Method
CN106053061B (en) A kind of envelope Analysis Method for decomposing and composing kurtosis based on nonlinear model
CN106096201B (en) A kind of EEMD and smoothed cubic spline envelope Analysis Method of rotating machinery
CN106053059B (en) It is a kind of based on it is interior grasp time scale decompose filtering envelope Analysis Method
CN105954030B (en) It is a kind of based on it is interior grasp time scale decompose and spectrum kurtosis envelope Analysis Method
CN105973603B (en) The EEMD and rational spline smoothed envelope analysis method of a kind of rotating machinery
CN106124200B (en) A kind of ELMD of rotating machinery and smooth iteration envelope Analysis Method
CN106198016B (en) A kind of NMD of rolling bearing, spectrum kurtosis and smooth iteration envelope Analysis Method

Legal Events

Date Code Title Description
C06 Publication
PB01 Publication
C10 Entry into substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant
CF01 Termination of patent right due to non-payment of annual fee

Granted publication date: 20180622

Termination date: 20210629

CF01 Termination of patent right due to non-payment of annual fee