CN106053059A - Envelope analysis method based on intrinsic time scale decomposition filtering - Google Patents

Envelope analysis method based on intrinsic time scale decomposition filtering Download PDF

Info

Publication number
CN106053059A
CN106053059A CN201610492071.7A CN201610492071A CN106053059A CN 106053059 A CN106053059 A CN 106053059A CN 201610492071 A CN201610492071 A CN 201610492071A CN 106053059 A CN106053059 A CN 106053059A
Authority
CN
China
Prior art keywords
signal
envelope
time scale
component
data
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
Application number
CN201610492071.7A
Other languages
Chinese (zh)
Other versions
CN106053059B (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 CN201610492071.7A priority Critical patent/CN106053059B/en
Publication of CN106053059A publication Critical patent/CN106053059A/en
Application granted granted Critical
Publication of CN106053059B publication Critical patent/CN106053059B/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/02Gearings; Transmission mechanisms
    • G01M13/028Acoustic or vibration analysis
    • 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 an envelope analysis method based on intrinsic time scale decomposition filtering. The method includes that, firstly, the original signal is decomposed by the intrinsic time scale decomposition method, and then the noise component and the trend item in the decomposition result are eliminated by the data rearrangement and substitution operation, the signal after the primary filtering is analyzed by the spectral kurtosis method to obtain the center frequency and bandwidth of an optimal filter, the filter filters the signal after the primary filtering secondarily, the signal after the secondary filtering is subjected to envelope analysis by means of the rational spline iteration smoothing envelope analysis method, and at the end, the fault type of a rotary machine can be determined based on the spectrum envelope. The method is suitable for dealing with complex rotary machine fault signals, can accurately determine the fault type of the rotary machine, is good in noise resistance and robustness, and is convenient for engineering application.

Description

A kind of based on the interior envelope Analysis Method grasping time scale decomposition filtering
Technical field
The present invention relates to condition monitoring for rotating machinery and fault diagnosis field, be specifically related to one and grasp time scale based on interior Decompose 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: the most existing Envelope Analysis technology or directly primary signal is analyzed, or only to original Signal is analyzed after simply filtering again, and the most existing method is easily subject to the dry of noise, trend and other composition Disturb, thus cause the analysis precision of prior art relatively low;The most existing Envelope Analysis technology is to be transformed to basis with Hilbert, And Hilbert conversion requires that analyzed signal must be the narrow band signal of simple component, otherwise the frequency modulating section of signal will The amplitude envelope analysis result of signal to be polluted, but signal the most to be analyzed the most strictly meets the bar of simple component and arrowband Part, so may result in prior art and easily occurs erroneous judgement problem because precision is the highest;3. the envelope spectrum obtained by traditional method There is end effect.
Summary of the invention
The problem to be solved in the present invention is for above not enough, proposes a kind of based on the interior bag grasping time scale decomposition filtering Method analyzed by network, after using the envelope Analysis Method of the present invention, has analysis result accuracy and degree of accuracy is high, and can be exactly The advantage detecting rotating machinery fault type.
For solving above technical problem, the technical scheme that the present invention takes is as follows: one is decomposed based on interior time scale of grasping The envelope Analysis Method of filtering, it is characterised in that comprise the following steps:
Step 1: utilize acceleration transducer to measure the vibration signal x(k of rotating machinery with sample frequency fs), (k=1,2, ..., N), N is the length of sampled signal;
Step 2: grasp time scale decomposition algorithm in employing by signal x(k) resolve into n component and a trend term sum, i.e., wherein, ciK () represents and is grasped the i-th component that time scale decomposition algorithm obtains, r by interiorn(k) Represent and grasped, by interior, the trend term that time scale decomposition algorithm obtains;
Step 3: to ciK () performs reordering operations and substitutes operation, data c that rearranged operation obtainsi shuffleK () represents, Data c are obtained after substituting operationi FTranK () represents;
Step 4: to ci(k), ci shuffle(k) and ci FTranK () performs multi-fractal respectively and removes trend fluction analysis (Multifractal Detrended Fluctuation Analysis, MFDFA), obtains generalized Hurst index curve, ci The generalized Hurst index curve H of (k)iQ () represents;ci shuffleThe generalized Hurst index curve H of (k)i shuffle(q) table Show;ci FTranThe generalized Hurst index curve H of (k)i FTranQ () represents;
Step 5: if Hi(q) and Hi shuffle(q) or Hi(q) and Hi FTranQ the relative error between () is less than 5%, or Hi (q), Hi shuffle(q) and Hi FTranQ () three do not change with q, then abandon the c of correspondencei(k) component;
Step 6: to remaining ci(k) component sue for peace, by this and be designated as signal rearranged and substitute filtered result xf1(k);
Step 7: to xf1K () performs spectrum kurtosis analysis, obtain the mid frequency f corresponding to signal kurtosis maximum0And bandwidth B;
Step 8: according to mid frequency f0With bandwidth B to xf1K () carries out bandpass filtering, obtain xf2(k);
Step 9: to signal xf2K () performs rational spline iteration smoothed envelope and analyzes, obtain signal envelope eov(k);
Step 10: the signal envelope eov(k to obtaining) perform discrete Fourier transform obtain envelope spectrum, according to envelope spectrum feature Frequency judges the fault type of machine.
A kind of prioritization scheme, in described step 2 in grasp time scale decomposition algorithm and comprise the following steps:
1) for signal xt, (t=1,2 ..., N), define an operatorFor extracting low frequency background signal, it may be assumed that
WhereinIt is background signal,It is an intrinsic rotational component, it is assumed thatIt is One real-valued signal,Represent xtThe moment corresponding to local extremum, define for convenience
If xtCertain interval has steady state value, it is contemplated that neighbouring signal also exists fluctuation, and we still believe that xt? Extreme value is comprised, at this moment on this intervalIt it is the right endpoint in this interval;
For convenience, definition,;AssumeWith?On be defined, xk?On be defined, in intervalOn continuous threshold point between define the extraction of piecewise linear background signal Operator, it may be assumed that
Wherein
Here parameterIt is a linear gain,, in this example
2) one intrinsic rotational component extraction operator of definition, it may be assumed that
Wherein,It it is the component that the 1st time Breaking Recurrently obtains;It it is the background signal that the 1st time Breaking Recurrently obtains;The 1st In secondary iteration, xtRepresent x (k) in step 2 described in claim 1;
3) againAs new data, repeat the above steps, the intrinsic rotational component that frequency reduces successively can be isolated, until background signal becomes dullness;
So xtWhole catabolic process can be written as:
WhereinRepresent ith iteration and decompose the component obtained,Represent ith iteration and decompose the background signal obtained.
Further, in described step 3, data rearrangement operation comprises the following steps:
Upset component c at randomiPutting in order of (k).
Further, in described step 3, the operation of data replacement comprises the following steps:
1) to component ciK () performs discrete Fourier transform, it is thus achieved that component ciThe phase place of (k);
2) component c is replaced with one group of pseudo-independent same distribution number being positioned in (-π, π) intervaliThe original phase of (k);
3) frequency domain data after phase place substitutes is performed inverse discrete Fourier transform and obtain data ci IFFTK (), asks for number According to ci IFFTThe real part of (k).
Further, in described step 4, MFDFA method comprises the following steps:
1) structurex (k )(k =1,2,…,N) profileY (i):
x (k) represent the 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 nonoverlappingN s Segment length issData, due to data lengthNIt is generally not capable of Divide exactlys, can not utilize so one piece of data can be remained;
In order to make full use of the length of data, then from the opposite direction of data with identical length segmentation, obtain the most altogether 2N s Segment data;
3) utilize the polynomial trend of the every segment data of least square fitting, then calculate the variance of every segment data:
y v (i) it is the of matchingvThe trend of segment data, if the polynomial trend of matching ismRank, then remember that this goes trend process For (MF-) DFAm;In this example, m=1;
4) the is calculatedqThe meansigma methods of rank wave function:
5) ifx (k) there is self-similarity characteristics, thenqThe meansigma methods of rank wave functionF q (s) and time scalesIt Between there is power law relation:
F q (s )~s H(q)
WhenqWhen=0, the formula in step 4) dissipates, at this momentH(0) determined by logarithmic mean process defined in following formula:
6) are taken the logarithm in the formula both sides in step 5) can obtain ln [F q (s )]=H (q )ln(s )+c (cFor constant), by This can obtain the slope of straight lineH (q )。
Further, the spectrum kurtosis method in described step 7 comprises the following steps:
1) one cut-off frequency of structure isf c =0.125+εLow pass filterh (n);ε> 0, in this examplef c =0.3;
2) based onh (n) structure 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) signalWarph 0(n )、h 1 (n) filtering and down-sampled after resolve into low frequency partAnd HFS, the down-sampled factor is 2, then shaping filter tree after successive ignition filters, and kth layer has 2 k Individual frequency band, whereinRepresent in wave filter tree thekOn layeriThe output signal of individual 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 mid frequency of individual wave filterf ki And bandwidthB k It is respectively
f ki =(i +2-1)2-k -1
B k =2-k-1
5) each filter results is calculated(i=0,…,2k-1) kurtosis
6) all of spectrum kurtosis is collected, obtain the spectrum kurtosis that signal is total.
Further, the analysis of the rational spline iteration smoothed envelope in described step 9 method comprises the following steps:
1) signal calculatedz (k) absolute valuez (k) local extremum;In the 1st iteration,z (k) representing right will Seek x in step 9 described in 1f2(k);
2) rational spline curve matching Local Extremum is used to obtain envelope eov1(k);
3) rightz (k) be normalized and obtain
4) the 2nd iteration:z 1(k) again as new data, repeat above-mentioned steps 1) ~ 3), obtain
5) ith iteration:z i-1(k) again as new data, repeat above-mentioned steps 1) ~ 3), obtain
6) ifnSecondary iteration obtainsz n (k) amplitude less than or equal to 1, then iterative process stops, and finally obtains letter Numberz (k) envelope be
The present invention uses above technical scheme, compared with prior art, the invention have the advantages that
1) grasp time scale decomposition in utilizing primary signal is decomposed, then utilize the rearrangement of data and substitute operation eliminating Noise therein and trend component, the only useful component in stick signal component, thus avoid noise and trend component pair The impact of Envelope Analysis result, analysis result accuracy and degree of accuracy are high.
2) utilize rational spline iteration smoothed envelope to analyze method to be kept completely separate with frequency modulating section by signal envelope, energy Enough avoid the frequency modulating section impact on signal envelope analysis result, thus improve the precision of Envelope Analysis.
3) fault type of rotating machinery can be detected exactly.
4) there is end effect in the envelope spectrum obtained by traditional method, and the envelope spectrum obtained by the present invention it can be avoided that End effect.
The present invention will be further described with embodiment below in conjunction with the accompanying drawings.
Accompanying drawing explanation
Accompanying drawing 1 is the flow chart of the inventive method in the embodiment of the present invention.
Accompanying drawing 2 carries out preliminary exposition to signal show for using low pass filter and high pass filter in the embodiment of the present invention It is intended to.
Accompanying drawing 3 is the schematic diagram using tree-shaped filter construction quickly to calculate spectrum kurtosis in the embodiment of the present invention.
Accompanying drawing 4 is the bearing vibration signal in the embodiment of the present invention with inner ring fault.
Accompanying drawing 5 is to use tradition envelope Analysis Method to inner ring faulty bearing vibration signal in the embodiment of the present invention Analysis result.
Accompanying drawing 6 is the present invention analysis result to inner ring faulty bearing vibration signal in the embodiment of the present invention.
Accompanying drawing 7 is the bearing vibration signal in the embodiment of the present invention with outer ring fault.
Accompanying drawing 8 is to use tradition envelope Analysis Method to outer ring faulty bearing vibration signal in the embodiment of the present invention Analysis result.
Accompanying drawing 9 is the present invention analysis result to outer ring faulty bearing vibration signal in the embodiment of the present invention.
Detailed description of the invention
Embodiment, as shown in Figure 1, Figure 2, Figure 3 shows, a kind of based on the interior envelope Analysis Method grasping time scale decomposition filtering, It is characterized in that, comprise the following steps:
Step 1: utilize acceleration transducer to measure the vibration signal x(k of rotating machinery with sample frequency fs), (k=1,2, ..., N), N is the length of sampled signal;
Step 2: grasp time scale decomposition algorithm in employing by signal x(k) resolve into n component and a trend term sum, i.e., wherein, ciK () represents and is grasped the i-th component that time scale decomposition algorithm obtains, r by interiorn(k) Represent and grasped, by interior, the trend term that time scale decomposition algorithm obtains;
Step 3: to ciK () performs reordering operations and substitutes operation, data c that rearranged operation obtainsi shuffleK () represents, Data c are obtained after substituting operationi FTranK () represents;
Step 4: to ci(k), ci shuffle(k) and ci FTranK () performs multi-fractal respectively and removes trend fluction analysis (Multifractal Detrended Fluctuation Analysis, MFDFA), obtains generalized Hurst index curve, ci The generalized Hurst index curve H of (k)iQ () represents;ci shuffleThe generalized Hurst index curve H of (k)i shuffle(q) table Show;ci FTranThe generalized Hurst index curve H of (k)i FTranQ () represents;
Step 5: if Hi(q) and Hi shuffle(q) or Hi(q) and Hi FTranQ the relative error between () is less than 5%, or Hi (q), Hi shuffle(q) and Hi FTranQ () three do not change with q, then abandon the c of correspondencei(k) component;
Step 6: to remaining ci(k) component sue for peace, by this and be designated as signal rearranged and substitute filtered result xf1(k);
Step 7: to xf1K () performs spectrum kurtosis analysis, obtain the mid frequency f corresponding to signal kurtosis maximum0And bandwidth B;
Step 8: according to mid frequency f0With bandwidth B to xf1K () carries out bandpass filtering, obtain xf2(k);
Step 9: to signal xf2K () performs rational spline iteration smoothed envelope and analyzes, obtain signal envelope eov(k);
Step 10: the signal envelope eov(k to obtaining) perform discrete Fourier transform obtain envelope spectrum, according to envelope spectrum feature Frequency judges the fault type of machine.
Grasp time scale decomposition algorithm in step 2 to comprise the following steps:
1) for signal xt, (t=1,2 ..., N), define an operatorFor extracting low frequency background signal, it may be assumed that
WhereinIt is background signal,It is an intrinsic rotational component, it is assumed thatIt is One real-valued signal,Represent xtThe moment corresponding to local extremum, define for convenience
If xtCertain interval has steady state value, it is contemplated that neighbouring signal also exists fluctuation, and we still believe that xt? Extreme value is comprised, at this moment on this intervalIt it is the right endpoint in this interval;
For convenience, definition,;AssumeWith?On be defined, xk?On be defined, in intervalOn continuous threshold point between define a piecewise linear background signal take out Take operator, it may be assumed that
Wherein
Here parameterIt is a linear gain,, in this example
2) one intrinsic rotational component extraction operator of definition, it may be assumed that
Wherein,It it is the component that the 1st time Breaking Recurrently obtains;It it is the background signal that the 1st time Breaking Recurrently obtains;The 1st In secondary iteration, xtRepresent x (k) in step 2 described in claim 1;
3) againAs new data, repeat the above steps, the intrinsic rotational component that frequency reduces successively can be isolated, until background signal becomes dullness;
So xtWhole catabolic process can be written as:
WhereinRepresent ith iteration and decompose the component obtained,Represent ith iteration and decompose the background signal obtained.
In step 3, data rearrangement operation comprises the following steps:
Upset component c at randomiPutting in order of (k).
In step 3, the operation of data replacement comprises the following steps:
1) to component ciK () performs discrete Fourier transform, it is thus achieved that component ciThe phase place of (k);
2) component c is replaced with one group of pseudo-independent same distribution number being positioned in (-π, π) intervaliThe original phase of (k);
3) frequency domain data after phase place substitutes is performed inverse discrete Fourier transform and obtain data ci IFFTK (), asks for number According to ci IFFTThe real part of (k).
In step 4, MFDFA method comprises the following steps:
1) structurex (k )(k =1,2,…,N) profileY (i):
x (k) represent the 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 nonoverlappingN s Segment length issData, due to data lengthNIt is generally not capable of Divide exactlys, can not utilize so one piece of data can be remained;
In order to make full use of the length of data, then from the opposite direction of data with identical length segmentation, obtain the most altogether 2N s Segment data;
3) utilize the polynomial trend of the every segment data of least square fitting, then calculate the variance of every segment data:
y v (i) it is the of matchingvThe trend of segment data, if the polynomial trend of matching ismRank, then remember that this goes trend process For (MF-) DFAm;In this example, m=1;
4) the is calculatedqThe meansigma methods of rank wave function:
5) ifx (k) there is self-similarity characteristics, thenqThe meansigma methods of rank wave functionF q (s) and time scalesIt Between there is power law relation:
F q (s )~s H( q)
WhenqWhen=0, the formula in step 4) dissipates, at this momentH(0) determined by logarithmic mean process defined in following formula:
6) are taken the logarithm in the formula both sides in step 5) can obtain ln [F q (s )]=H (q )ln(s )+c (cFor constant), by This can obtain the slope of straight lineH (q)。
Spectrum kurtosis method in step 7 comprises the following steps:
1) one cut-off frequency of structure isf c =0.125+εLow pass filterh (n);ε> 0, in this examplef c =0.3;
2) based onh (n) structure 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) signalWarph 0(n )、h 1 (n) filtering and down-sampled after resolve into low frequency partAnd HFS, the down-sampled factor is 2, then shaping filter tree after successive ignition filters, and kth layer has 2 k Individual frequency band, whereinRepresent in wave filter tree thekOn layeriThe output signal of individual 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 mid frequency of individual wave filterf ki And bandwidthB k It is respectively
f ki =(i +2-1)2- k-1
B k =2- k-1
5) each filter results is calculated(i=0,…,2k-1) kurtosis
6) all of spectrum kurtosis is collected, obtain the spectrum kurtosis that signal is total.
Rational spline iteration smoothed envelope in step 9 is analyzed method and is comprised the following steps:
1) signal calculatedz (k) absolute valuez (k) local extremum;In the 1st iteration,z (k) representing right will Seek x in step 9 described in 1f2(k);
2) rational spline curve matching Local Extremum is used to obtain envelope eov1(k);
3) rightz (k) be normalized and obtain
4) the 2nd iteration:z 1(k) again as new data, repeat above-mentioned steps 1) ~ 3), obtain
5) ith iteration:z i-1(k) again as new data, repeat above-mentioned steps 1) ~ 3), obtain
6) ifnSecondary iteration obtainsz n (k) amplitude less than or equal to 1, then iterative process stops, and finally obtains letter Numberz (k) envelope be
Test 1, utilizes the bearing vibration data with inner ring fault to test the performance of algorithm of the present invention Card.
Experiment bearing used is 6205-2RS JEM SKF, utilizes electric discharge machining method working depth on bearing inner race For 0.2794mm, width be the groove of 0.3556mm to simulate bearing inner race fault, this experiment load is about 0.7457kW, drives Motor turns frequency and is about 29.5Hz, and bearing inner race fault characteristic frequency is about 160Hz, and sample frequency is 4.8KHz, during signal sampling A length of 1s.
The inner ring fault-signal collected is as shown in Figure 4.
Initially with traditional envelope Analysis Method, the signal shown in Fig. 4 is analyzed, the analysis result obtained such as Fig. 5 Shown in.From fig. 5, it can be seen that the fault signature of bearing is blanked completely, the most traditional envelope Analysis Method can not be effectively Extract the fault signature of bearing;Additionally, from fig. 5, it can be seen that the left end point of envelope spectrum also exists abnormal high level, this explanation is by passing The envelope spectrum that system method obtains also exists end effect.
Use method proposed by the invention that signal shown in Fig. 4 is analyzed, the analysis result obtained such as Fig. 6 institute Show.From fig. 6, it can be seen that the spectral line corresponding to 160Hz and 320Hz is apparently higher than other spectral line, the two frequency correspondence respectively 1 frequency multiplication of bearing inner race fault characteristic frequency and 2 frequencys multiplication, may determine that bearing has inner ring fault accordingly;Can from Fig. 6 Go out, the present invention envelope spectrum obtained does not has end effect.
Showing through many experiments, in the case of load and fault dimensional depth are constant, the present invention can reliable recognition Minimum inner ring fault dimension width is about 0.22 mm, and traditional method can the minimum inner ring fault dimension width of reliable recognition Being about 0.53mm, precision improves 58.5%.
Test 2, utilizes the bearing vibration data with outer ring fault to test the performance of algorithm of the present invention Card.
Experiment bearing used is 6205-2RS JEM SKF, utilizes electric discharge machining method working depth on bearing outer ring For 0.2794mm, width be the groove of 0.5334mm to simulate bearing outer ring fault, this experiment load is about 2.237 kW, drives Motor turns frequency and is about 28.7Hz, and bearing outer ring fault characteristic frequency is about 103Hz, and sample frequency is 4.8KHz, during signal sampling A length of 1s.
The outer ring fault-signal collected is as shown in Figure 7.
Initially with traditional envelope Analysis Method, the signal shown in Fig. 7 is analyzed, the analysis result obtained such as Fig. 8 Shown in.From figure 8, it is seen that the fault signature of bearing is blanked completely, the most traditional envelope Analysis Method can not be effectively Extract the fault signature of bearing;Additionally, from figure 8, it is seen that the left end point of envelope spectrum also exists abnormal high level, this explanation is by passing The envelope spectrum that system method obtains also exists end effect.
Use method proposed by the invention that signal shown in Fig. 7 is analyzed, the analysis result obtained such as Fig. 9 institute Show.From fig. 9, it can be seen that the spectral line corresponding to 103Hz and 206Hz is apparently higher than other spectral line, the two frequency correspondence respectively 1 frequency multiplication of bearing outer ring fault characteristic frequency and 2 frequencys multiplication, may determine that bearing has outer ring fault accordingly;Can from Fig. 9 Go out, the present invention envelope spectrum obtained does not has end effect.
Showing through many experiments, in the case of load and fault dimensional depth are constant, the present invention can reliable recognition Minimum outer ring fault dimension width is about 0.30mm, and traditional method can reliable recognition minimum outer ring fault dimension width about For 0.68mm, precision improves 55.9%.
According to result of the test, think after analysis:
1) traditional envelope Analysis Method directly carries out Envelope Analysis to primary signal, or to after merely through simple process Primary signal carries out Envelope Analysis, different from traditional envelope Analysis Method, and the present invention decomposes first with interior time scale of grasping Primary signal is decomposed, then utilizes the rearrangement of data and substitute operation eliminating noise therein and trend component, only Useful component in stick signal component, thus avoid the impact on Envelope Analysis result of noise and trend component, improve Accuracy and precision.
2) traditional envelope Analysis Method is transformed to basis with Hilbert, and Hilbert conversion requires analyzed letter Number must be the narrow band signal of simple component, otherwise the frequency modulating section of signal will pollute the Envelope Analysis result of signal, but It is the signal the most to be analyzed condition that the most strictly meets simple component and arrowband, so may result in prior art because of precision not High and erroneous judgement problem easily occur, different from tradition envelope Analysis Method, the present invention utilizes rational spline iteration smoothed envelope to divide Signal envelope is kept completely separate by analysis method with frequency modulating section, it is possible to avoid frequency modulating section to signal envelope analysis result Impact, thus improve the precision of Envelope Analysis.
3) fault type of rotating machinery can be detected exactly.
4) there is end effect in the envelope spectrum obtained by traditional method, and the envelope spectrum obtained by the present invention it can be avoided that End effect.
5) each step effect:
1st) step: gather vibration signal;
2nd) step: primary signal is resolved into the form of different component sum, some of which component correspondence noise and trend term, some Component correspondence useful signal;
3rd) ~ 5) step: the signal that obtains above-mentioned decomposition performs reordering operations and substitutes operation, reject noise component(s) therein and Trend term, only retains useful signal;
6th) step: remaining useful signal is sued for peace, should and as signal rearranged and substitute filtered result xf1(k);
7th) step: to filtered signal xf1K () performs spectrum kurtosis analysis, obtain center frequency corresponding at signal maximum kurtosis Rate f0And bandwidth B;
8th) step: according to mid frequency f0With bandwidth B to xf1K () carries out bandpass filtering, obtain signal xf2(k);
9th) step: signal calculated xf2Envelope eov (k) of (k);
10th) step: eov (k) is performed discrete Fourier transform and obtains envelope spectrum, judge the failure classes of bearing according to envelope spectrum Type.
One skilled in the art would recognize that above-mentioned detailed description of the invention is exemplary, be to make ability Field technique personnel can be better understood from 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, each fall within protection scope of the present invention.

Claims (7)

1. one kind is decomposed the envelope Analysis Method of filtering based on interior time scale of grasping, it is characterised in that comprise the following steps:
Step 1: utilize acceleration transducer to measure the vibration signal x(k of rotating machinery with sample frequency fs), (k=1,2, ..., N), N is the length of sampled signal;
Step 2: grasp time scale decomposition algorithm in employing by signal x(k) resolve into n component and a trend term sum, i.e., wherein, ciK () represents and is grasped the i-th component that time scale decomposition algorithm obtains, r by interiorn(k) Represent and grasped, by interior, the trend term that time scale decomposition algorithm obtains;
Step 3: to ciK () performs reordering operations and substitutes operation, data c that rearranged operation obtainsi shuffleK () represents, Data c are obtained after substituting operationi FTranK () represents;
Step 4: to ci(k), ci shuffle(k) and ci FTranK () performs multi-fractal respectively and removes trend fluction analysis (Multifractal Detrended Fluctuation Analysis, MFDFA), obtains generalized Hurst index curve, ci The generalized Hurst index curve H of (k)iQ () represents;ci shuffleThe generalized Hurst index curve H of (k)i shuffle(q) table Show;ci FTranThe generalized Hurst index curve H of (k)i FTranQ () represents;
Step 5: if Hi(q) and Hi shuffle(q) or Hi(q) and Hi FTranQ the relative error between () is less than 5%, or Hi (q), Hi shuffle(q) and Hi FTranQ () three do not change with q, then abandon the c of correspondencei(k) component;
Step 6: to remaining ci(k) component sue for peace, by this and be designated as signal rearranged and substitute filtered result xf1(k);
Step 7: to xf1K () performs spectrum kurtosis analysis, obtain the mid frequency f corresponding to signal kurtosis maximum0And bandwidth B;
Step 8: according to mid frequency f0With bandwidth B to xf1K () carries out bandpass filtering, obtain xf2(k);
Step 9: to signal xf2K () performs rational spline iteration smoothed envelope and analyzes, obtain signal envelope eov(k);
Step 10: the signal envelope eov(k to obtaining) perform discrete Fourier transform obtain envelope spectrum, according to envelope spectrum feature Frequency judges the fault type of machine.
The most according to claim 1 a kind of based on the interior envelope Analysis Method grasping time scale decomposition filtering, its feature exists In, in described step 2 in grasp time scale decomposition algorithm and comprise the following steps:
1) for signal xt, (t=1,2 ..., N), define an operatorFor extracting low frequency background signal, it may be assumed that
WhereinIt is background signal,It is an intrinsic rotational component, it is assumed thatIt is one Individual real-valued signal,Represent xtThe moment corresponding to local extremum, define for convenience
If xtCertain interval has steady state value, it is contemplated that neighbouring signal also exists fluctuation, and we still believe that xtAt this Extreme value is comprised, at this moment on individual intervalIt it is the right endpoint in this interval;
For convenience, definition,;AssumeWith?On be defined, xk?On be defined, in intervalOn continuous threshold point between define a piecewise linear background signal take out Take operator, it may be assumed that
Wherein
Here parameterIt is a linear gain,
2) one intrinsic rotational component extraction operator of definition, it may be assumed that
Wherein,It it is the component that the 1st time Breaking Recurrently obtains;It it is the background signal that the 1st time Breaking Recurrently obtains;
3) againAs new data, repeat the above steps, the intrinsic rotational component that frequency reduces successively can be isolated, until background signal becomes dullness;
So xtWhole catabolic process can be written as:
WhereinRepresent ith iteration and decompose the component obtained,Represent ith iteration and decompose the background signal obtained.
The most according to claim 1 a kind of based on the interior envelope Analysis Method grasping time scale decomposition filtering, its feature exists In, in described step 3, data rearrangement operation comprises the following steps:
Upset component c at randomiPutting in order of (k).
The most according to claim 1 a kind of based on the interior envelope Analysis Method grasping time scale decomposition filtering, its feature exists In: in described step 3, the operation of data replacement comprises the following steps:
1) to component ciK () performs discrete Fourier transform, it is thus achieved that component ciThe phase place of (k);
2) component c is replaced with one group of pseudo-independent same distribution number being positioned in (-π, π) intervaliThe original phase of (k);
3) frequency domain data after phase place substitutes is performed inverse discrete Fourier transform and obtain data ci IFFTK (), asks for data ci IFFTThe real part of (k).
The most according to claim 1 a kind of based on the interior envelope Analysis Method grasping time scale decomposition filtering, its feature exists In: in described step 4, MFDFA method comprises the following steps:
1) structurex (k )(k =1,2,…,N) profileY (i):
x (k) represent the 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 nonoverlappingN s Segment length issData, from the opposite direction of data with identical Length segmentation, obtains 2N s Segment data;
3) utilize the polynomial trend of the every segment data of least square fitting, then calculate the variance of every segment data:
y v (i) it is the of matchingvThe trend of segment data, if the polynomial trend of matching ismRank, then remember that this goes trend process For (MF-) DFAm
4) the is calculatedqThe meansigma methods of rank wave function:
5) ifx (k) there is self-similarity characteristics, thenqThe meansigma methods of rank wave functionF q (s) and time scalesIt Between there is power law relation:
F q (s )~s H(q )
WhenqWhen=0, the formula in step 4) dissipates, at this momentH(0) determined by logarithmic mean process defined in following formula:
6) are taken the logarithm in the formula both sides in step 5) can obtain ln [F q (s )]=H (q )ln(s )+c (cFor constant), thus The slope of straight line can be obtainedH (q )。
The most according to claim 1 a kind of based on the interior envelope Analysis Method grasping time scale decomposition filtering, its feature exists In: the spectrum kurtosis method in described step 7 comprises the following steps:
1) one cut-off frequency of structure isf c =0.125+εLow pass filterh (n);
2) based onh (n) structure 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) signalWarph 0(n )、h 1 (n) filtering and down-sampled after resolve into low frequency partAnd HFS, the down-sampled factor is 2, then shaping filter tree after successive ignition filters, and kth layer has 2 k Individual frequency band, whereinRepresent in wave filter tree thekOn layeriThe output signal of individual wave filter,i =0,…,2k-1,0≤k≤K-1;
4) in decomposition treekOn layeriThe mid frequency of individual wave filterf ki And bandwidthB k It is respectively
f ki =(i +2-1)2-k -1
B k =2- k-1
5) each filter results is calculated(i=0,…,2k-1) kurtosis
6) all of spectrum kurtosis is collected, obtain the spectrum kurtosis that signal is total.
The most according to claim 1 a kind of based on the interior envelope Analysis Method grasping time scale decomposition filtering, its feature exists In, the rational spline iteration smoothed envelope in described step 9 is analyzed method and is comprised the following steps:
1) signal calculatedz (k) absolute valuez (k) local extremum;
2) rational spline curve matching Local Extremum is used to obtain envelope eov1(k);
3) rightz (k) be normalized and obtain
4) the 2nd iteration:z 1(k) again as new data, repeat above-mentioned steps 1) ~ 3), obtain
5) ith iteration:z i-1(k) again as new data, repeat above-mentioned steps 1) ~ 3), obtain
6) ifnSecondary iteration obtainsz n (k) amplitude less than or equal to 1, then iterative process stops, and finally obtains letter Numberz (k) envelope be
CN201610492071.7A 2016-06-29 2016-06-29 It is a kind of based on it is interior grasp time scale decompose filtering envelope Analysis Method Expired - Fee Related CN106053059B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201610492071.7A CN106053059B (en) 2016-06-29 2016-06-29 It is a kind of based on it is interior grasp time scale decompose filtering envelope Analysis Method

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201610492071.7A CN106053059B (en) 2016-06-29 2016-06-29 It is a kind of based on it is interior grasp time scale decompose filtering envelope Analysis Method

Publications (2)

Publication Number Publication Date
CN106053059A true CN106053059A (en) 2016-10-26
CN106053059B CN106053059B (en) 2018-03-27

Family

ID=57165987

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201610492071.7A Expired - Fee Related CN106053059B (en) 2016-06-29 2016-06-29 It is a kind of based on it is interior grasp time scale decompose filtering envelope Analysis Method

Country Status (1)

Country Link
CN (1) CN106053059B (en)

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN108152363A (en) * 2017-12-21 2018-06-12 北京工业大学 A kind of defect of pipeline recognition methods for pressing down end intrinsic time Scale Decomposition
CN109977627A (en) * 2019-05-10 2019-07-05 江南大学 A kind of networking Multi-sensor Fusion fault detection method

Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2013088431A (en) * 2011-10-13 2013-05-13 General Electric Co <Ge> Method and system for automatically detecting rolling bearing fault
CN104677632A (en) * 2015-01-21 2015-06-03 大连理工大学 Rolling bearing fault diagnosis method using particle filtering and spectral kurtosis
CN104792528A (en) * 2014-01-22 2015-07-22 中国人民解放军海军工程大学 Adaptive optimal envelope demodulation method

Patent Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2013088431A (en) * 2011-10-13 2013-05-13 General Electric Co <Ge> Method and system for automatically detecting rolling bearing fault
CN104792528A (en) * 2014-01-22 2015-07-22 中国人民解放军海军工程大学 Adaptive optimal envelope demodulation method
CN104677632A (en) * 2015-01-21 2015-06-03 大连理工大学 Rolling bearing fault diagnosis method using particle filtering and spectral kurtosis

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
林近山等: "多重分形去趋势波动分析在滚动轴承损伤程度识别中的应用", 《中国机械工程》 *

Cited By (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN108152363A (en) * 2017-12-21 2018-06-12 北京工业大学 A kind of defect of pipeline recognition methods for pressing down end intrinsic time Scale Decomposition
CN108152363B (en) * 2017-12-21 2021-06-25 北京工业大学 Pipeline defect identification method based on restrained end intrinsic time scale decomposition
CN109977627A (en) * 2019-05-10 2019-07-05 江南大学 A kind of networking Multi-sensor Fusion fault detection method
CN109977627B (en) * 2019-05-10 2023-05-16 江南大学 Networked multi-sensor fusion fault detection method

Also Published As

Publication number Publication date
CN106053059B (en) 2018-03-27

Similar Documents

Publication Publication Date Title
CN106198015A (en) The VMD of a kind of rolling bearing, spectrum kurtosis and smooth iteration envelope Analysis Method
CN106096313A (en) A kind of envelope Analysis Method decomposed based on singular spectrum and compose kurtosis
CN106168538A (en) The ITD of a kind of rolling bearing, spectrum kurtosis and smooth iteration envelope Analysis Method
CN106198013A (en) A kind of envelope Analysis Method based on empirical mode decomposition filtering
CN106096200A (en) A kind of envelope Analysis Method based on wavelet decomposition with spectrum kurtosis
CN106153339B (en) A kind of envelope Analysis Method based on the filtering of variation Mode Decomposition
CN106096198A (en) A kind of envelope Analysis Method based on variation Mode Decomposition with spectrum kurtosis
CN105954031A (en) Envelopment analysis method based on singular spectrum decomposition filtering
CN106096199A (en) The WT of a kind of rolling bearing, spectrum kurtosis and smooth iteration envelope Analysis Method
CN106053069A (en) SSD, spectral kurtosis and smooth iteration envelope analysis method of antifriction bearing
CN106198009A (en) The EMD of a kind of rolling bearing, spectrum kurtosis and smooth iteration envelope Analysis Method
CN106053060A (en) Envelope analysis method based on non-linear mode decomposition filtering
CN106153333A (en) A kind of envelope Analysis Method based on wavelet decomposition filtering
CN106198012A (en) A kind of envelope Analysis Method decomposed based on local mean value and compose kurtosis
CN106198010A (en) A kind of envelope Analysis Method decomposing filtering based on local mean value
CN106053059A (en) Envelope analysis method based on intrinsic time scale decomposition filtering
CN105954030B (en) It is a kind of based on it is interior grasp time scale decompose and spectrum kurtosis envelope Analysis Method
CN106053061A (en) Envelope analysis method based on non-linear mode decomposition and spectrum kurtosis
CN106198014A (en) A kind of envelope Analysis Method based on empirical mode decomposition with spectrum kurtosis
CN106198017A (en) The LMD of a kind of rolling bearing, spectrum kurtosis and smooth iteration envelope Analysis Method
CN106096201A (en) The EEMD of a kind of rotating machinery and smoothed cubic spline envelope Analysis Method
CN106198018A (en) The EEMD of a kind of rotating machinery and smooth iteration envelope Analysis Method
CN105973603A (en) EEMD and rational spline smooth envelope analysis method for rotating machine
CN106198016A (en) The NMD of a kind of rolling bearing, spectrum kurtosis and smooth iteration envelope Analysis Method
CN106124200A (en) The ELMD of a kind of rotating machinery 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: 20180327

Termination date: 20210629

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