CN103499443B  A kind of gear distress is without key phase angular domain average computation order analysis method  Google Patents
A kind of gear distress is without key phase angular domain average computation order analysis method Download PDFInfo
 Publication number
 CN103499443B CN103499443B CN201310416504.7A CN201310416504A CN103499443B CN 103499443 B CN103499443 B CN 103499443B CN 201310416504 A CN201310416504 A CN 201310416504A CN 103499443 B CN103499443 B CN 103499443B
 Authority
 CN
 China
 Prior art keywords
 signal
 frequency
 omega
 time
 domain average
 Prior art date
Links
 238000004458 analytical method Methods 0.000 title claims abstract description 29
 230000001133 acceleration Effects 0.000 claims abstract description 17
 238000010845 search algorithm Methods 0.000 claims abstract description 16
 238000001914 filtration Methods 0.000 claims abstract description 11
 238000000034 method Methods 0.000 claims description 21
 238000005070 sampling Methods 0.000 claims description 10
 230000021615 conjugation Effects 0.000 claims description 3
 235000009808 lpulo Nutrition 0.000 claims description 3
 238000003745 diagnosis Methods 0.000 abstract description 9
 238000005516 engineering processes Methods 0.000 abstract description 6
 238000004064 recycling Methods 0.000 abstract 1
 238000010586 diagram Methods 0.000 description 6
 238000001228 spectrum Methods 0.000 description 5
 230000000694 effects Effects 0.000 description 4
 230000001360 synchronised Effects 0.000 description 4
 230000005540 biological transmission Effects 0.000 description 3
 230000000875 corresponding Effects 0.000 description 3
 230000000737 periodic Effects 0.000 description 3
 238000006243 chemical reaction Methods 0.000 description 2
 238000005260 corrosion Methods 0.000 description 2
 238000001514 detection method Methods 0.000 description 2
 229910000831 Steel Inorganic materials 0.000 description 1
 238000004364 calculation method Methods 0.000 description 1
 239000000969 carrier Substances 0.000 description 1
 238000010276 construction Methods 0.000 description 1
 238000000605 extraction Methods 0.000 description 1
 238000009499 grossing Methods 0.000 description 1
 238000003780 insertion Methods 0.000 description 1
 238000009434 installation Methods 0.000 description 1
 238000005259 measurement Methods 0.000 description 1
 239000000203 mixture Substances 0.000 description 1
 239000003638 reducing agent Substances 0.000 description 1
 239000010959 steel Substances 0.000 description 1
 230000001629 suppression Effects 0.000 description 1
Abstract
The invention provides a kind of gear distress without key phase angular domain average computation order analysis method, first the vibration acceleration signal of collection is carried out calculating its timefrequency distributions by smoothed pseudo wigner ville disstribution after low pass protects phase filtering, then the instantaneous frequency of gear case rotating shaft is estimated by Viterbi optimum route search algorithm, recycling key signal estimation model carries out pointbypoint integration to instantaneous frequency and obtains estimating key signal, finally combine angularly resampling and angle domain average technology carries out calculating order analysis to vibration acceleration signal, the order obtained based on instantaneous Frequency Estimation is composed.This order spectrogram fully can reflect the characteristic information of gearbox fault.The present invention combines the advantage of smoothed pseudo wigner ville disstribution, Viterbi optimum route search algorithm, angle domain average technology, calculating order analysis, effectively can carry out fault diagnosis to the gear case under variable speed operating condition.
Description
Technical field
The invention belongs to gear distress analytical approach field, be specifically related to a kind of gear distress without key phase angular domain average computation order analysis method.
Background technology
Traditional gear case vibrationtesting and the analysis of test signal, test are carried out under constant rotational speed, constant duty, and gear case in real work runs under variable speed, variable working condition load often, the such as gear case of wind power generator group, Helicopter gearboxes, car deceleration device, crane etc.For crane, handling steel ladle is responsible for by crane, the process that during work, gear case experience raising speed promotes, reduction of speed is toppled over.At present, under gear case variable speed operating condition condition, in calculating order analysis method used, the key signal of rotating shaft be carry out gearbox fault calculate order analysis institute requisite.But, for some casing relative closures (rotating machinery is nondisconnectable) or can not, at axle head sensor installation gear case, be difficult to extract key signal.Particularly some gear reduction design time inconsiderate or assembling time compact conformation, cannot fitting key phase pulse transducer, to be difficult under variable rotate speed realization and to wait periodic sampling, to cause in practical work process to gear case vibrationtesting and to the analysis of test signal, test is more difficult completes.Therefore; research effective method; for the gear case under variable rotate speed and rotary speed information cannot gather when monitor the running status of gear case; arrange rational predetermined spare part, shut down replacing construction; before maintenance of equipment is arranged in the generation of equipment failure accident, reduces burst shutdown event to greatest extent, give full play to equipment effectiveness; extension device serviceable life, significant to the security of national product, reliability.
Summary of the invention
The object of the present invention is to provide a kind of gear distress without key phase angular domain average computation order analysis method, the method can solve the difficulties carrying out fault diagnosis under variable rotate speed in the not collectable situation of key signal.
To achieve these goals, the technical solution used in the present invention is:
A kind of gear distress, without key phase angular domain average computation order analysis method, comprises the following steps:
1) instantaneous Frequency Estimation, specifically comprises following three steps:
11) gather vibration acceleration signal, low pass is carried out to the vibration acceleration signal gathered and protects phase filtering, obtain low frequency signal x (t) reflecting rotor turns;
The timefrequency distributions of low frequency signal x (t) 12) calculated by smoothed pseudo wigner ville disstribution, obtains time frequency distribution map;
In the time frequency distribution map obtained, 13) utilize Viterbi optimum route search algorithm to estimate the instantaneous frequency (InstantaneousFrequency, IF) of gear case rotating shaft; .
2) order based on instantaneous Frequency Estimation is composed, and specifically comprises following four steps:
21) utilizing key signal estimation model to step 13) instantaneous frequency that obtains carries out pointbypoint integration, obtains estimating key signal;
22) utilize the estimation key signal that obtains to carry out angularly resampling to vibration acceleration signal, obtain the signal after angularly resampling;
23) angle domain average process is carried out to the signal after angularly resampling, obtain the signal after angle domain average process;
24) Fast Fourier Transform (FFT) (FFT) is carried out to the signal after angular domain average treatment, obtain reflecting that the order based on instantaneous Frequency Estimation of gearbox fault characteristic information is composed, gear distress is analyzed.
Described step 12) concrete operations be:
Step 11) the smoothed pseudo wigner ville disstribution SPW (t, ω) of low frequency signal x (t) that obtains is by shown in formula (1):
Namely time frequency distribution map is obtained by formula (1); In formula (1), * represents conjugation, and t represents the time, and ω represents frequency, and h (τ) and g (u) is window function, and τ, u represent the shift factor of time shaft and frequency axis respectively.
Described window function comprises Gauss window and Hamming window.
Described step 13) concrete operations be:
Get the time series [n on the smoothed pseudo wigner ville disstribution of x (t)
_{1}, n
_{2}], time point n ∈ [n
_{1}, n
_{2}], the frequency values on each time point n is SPW
_{n}, SPW
_{n}by the array after descending sort for being designated as SPW '
_{n};
From 0 s, press increment 1 successively by order from big to small and distinguish assignment again, the timefrequency distributions sequence after assignment be designated as f (n, SPW '
_{n});
Instantaneous frequency is the path making formula (2) minimum:
Wherein, g (ω
_{n}, ω
_{n+1}) be continuity in order to meet frequency change and the penalty factor that larger frequency change of adjusting the distance adds.
Described penalty factor g (ω
_{n}, ω
_{n+1}) such as formula shown in (3):
Wherein δ is the maximum expected value of the Instantaneous frequency variations of 2 continuity points.
Described step 21) concrete operations be:
Work under a rotating shaft is in variable rotate speed, rotation instantaneous frequency is f (t); Elapsed time t
_{n}, rotating shaft turns over n circle, and rotating shaft turns over angle such as formula shown in (4):
To f (t) pointbypoint integration, record integrated value is the t of integer
_{n}, at t
_{n}the key signal of a step is inserted in position, is namely the estimation key signal obtained by instantaneous frequency.
Described step 22) concrete operations be:
Timing node t when obtaining angularly θ resampling according to formula (5), carries out angularly resampling, obtains the signal after angularly resampling;
B in formula (5)
_{0}, b
_{1}, b
_{2}obtain according to formula (6);
Wherein
Described step 23) concrete operations be:
Burst after angularly resampling is designated as
the burst y after angle domain average process is obtained by the formula formula (7) of angle domain average process
_{n},
Wherein, P is the hop count signal averaging after angularly resampling be divided into, and N is each section of identical sampling number.
Relative to prior art, beneficial effect of the present invention is:
The gear distress that this patent proposes is without key phase angular domain average computation order analysis method method, the rotary speed information of rotating shaft can be extracted from vibration acceleration signal, effectively solve the problem being difficult to diagnosis without gear case fault under variable speed operating condition time key signal (key signal can not gather), namely can without the diagnosis completed during key signal and under variable rotate speed gearbox fault.The method combines smoothed pseudo wigner ville disstribution, Viterbi optimum route search algorithm, key signal estimation, angle domain average process and calculates the advantage of order analysis, has following differences in the significant advantage of classic method:
1) in instantaneous Frequency Estimation, present invention incorporates the timefrequency focusing that smoothed pseudo wigner ville disstribution is good, and Viterbi optimum route search algorithm follows the feature of freerunning frequency change continuity principle, can estimate the instantaneous frequency of gear case rotating shaft effectively accurately.
2) compose in calculating at the order based on instantaneous Frequency Estimation, utilize key signal estimation model to carry out pointbypoint integration to instantaneous frequency and estimated key signal accurately, then combine angularly resampling technique and angle domain average treatment technology and calculating order analysis is carried out to vibration acceleration signal, effectively resist the interference of noise and unrelated signal, obtained reflecting that the order based on instantaneous Frequency Estimation of gearbox fault characteristic information is composed.
3) on the whole, method provided by the invention combines smoothed pseudo wigner ville disstribution, Viterbi optimum route search algorithm, key signal estimation, angle domain average treatment technology and calculates the advantage of order analysis, utilize the present invention effectively can carry out fault diagnosis to the gear case that cannot gather key signal under variable rotate speed, improve the accuracy rate of fault diagnosis.
Accompanying drawing explanation
Fig. 1 (a) is the schematic diagram of Viterbi optimum route search algorithm sweep forward;
Fig. 1 (b) is the schematic diagram of Viterbi optimum route search algorithm sweep backward;
Fig. 1 (c) is the schematic diagram in the local optimum path of Viterbi optimum route search algorithm;
Fig. 1 (d) is the schematic diagram in the global optimum path of Viterbi optimum route search algorithm;
Fig. 2 is the schematic diagram of the rotational model of axle;
Fig. 3 (a) is emulation vibration signal for faster figure;
Fig. 3 (b) is the interpolation resampling figure to emulation vibration signal for faster;
Fig. 3 (c) is the estimation key signal figure of emulation vibration signal for faster;
Fig. 3 (d) carries out the signal graph after angularly resampling to estimation key signal;
Fig. 4 tests the vibration acceleration signal figure recorded;
Fig. 5 vibrates the signal graph after signal for faster lowpass filtering to gained;
Fig. 6 is the SPWVD time frequency distribution map of signal after lowpass filtering, and wherein a is that 1X turns frequently, and b is that 2X turns frequently;
Fig. 7 is the instantaneous frequency figure utilizing Viterbi optimum route search algorithm search to obtain;
Fig. 8 is the estimation key signal figure utilizing key signal estimation model to obtain;
Fig. 9 is the signal graph after angularly resampling;
Figure 10 is the order spectrogram based on instantaneous Frequency Estimation;
Figure 11 is that the epicyclic gearbox centre wheel flank of tooth evenly puts corrosion figure.
Embodiment
Below in conjunction with accompanying drawing, content of the present invention is described in further detail.
The invention discloses a kind of gear distress without key phase angular domain average computation order analysis method (AngledomainaverageAndNonkeyphasorCalculatedOrderTrack ing, ANCOT), specifically implement according to the following steps:
1) instantaneous Frequency Estimation.Low pass is carried out to the vibration acceleration signal gathered and protects phase filtering, obtain the low frequency signal reflecting rotor turns.Then the timefrequency distributions of low frequency signal is calculated by smoothed pseudo wigner ville disstribution, the frequency division Butut obtained utilize Viterbi optimum route search algorithm estimate the instantaneous frequency (InstantaneousFrequency, IF) of gear case rotating shaft.Specifically comprise following 3 steps:
11) gather vibration acceleration signal, low pass is carried out to the vibration acceleration signal gathered and protects phase filtering, extract reflection turning axle and first the low pass guarantor phase filtering that input shaft turns signal is frequently comprised to signal, remove the interference of other orders in signal, the vibration signal (low frequency signal of reflection rotor turns) that to obtain with input shaft vibration signal and low power order thereof be fundamental component.The upper cutoff frequency of lowpass filter is according to the highest running speed setting in machine operation procedure.Protect phase filtering very important, if otherwise have the changing of the relative positions of phase place that IF can be caused to estimate and vibration signal on a timeline asynchronous after filtering simultaneously.
12) time frequency analysis: the timefrequency distributions of the low frequency signal calculated by smoothed pseudo wigner ville disstribution, obtains time frequency distribution map.
The WignerVille distribution of low frequency signal x (t) is defined as,
WignerVille has numerous advantages, particularly timefrequency locality, but can introduce cross term.In order to overcome the deficiency of cross term, can in time domain or time domain, frequency domain windowing simultaneously to overcome cross term, Here it is pseudoWinger distribution (PWVD) and Smoothing Pseudo Winger distribution (SPWVD).WignerVille distribution after timedomain windowed is called pseudo NMalgebra, is defined as follows:
Wherein, in formula, * represents conjugation.H (t) is window function.Conventional window function comprises Gauss window, Hamming window etc.
Then can at time and frequency both direction suppressing crossterms simultaneously in time domain, frequency domain windowing simultaneously.The WignerVille distribution of this both direction windowing process is called smoothed pseudo wigner ville disstribution (SPWVD), and the smoothed pseudo wigner ville disstribution of x (t) is defined as follows:
Wherein, t represents the time, and ω represents frequency, and g (u) is used to do level and smooth window function in frequency axis direction, and τ, u represent the shift factor of time shaft and frequency axis respectively.Namely time frequency distribution map is obtained by formula (1).
In the time frequency distribution map obtained, 13) utilize Viterbi optimum route search algorithm to estimate the instantaneous frequency (InstantaneousFrequency, IF) of gear case rotating shaft.
At 12) step tries to achieve in time frequency distribution map, adopts Viterbi algorithm search spectrum optimum time frequency path, peak.Spectrum blob detection just searches for the peak value in timefrequency distributions successively, and actual IF is as a continually varying physical quantity, and should consider that the change between two adjacent time points can not be excessive, Viterbi algorithm just in time considers this point.Timefrequency distributions proposes during Searching I F following two hypothesis,
(1) if on certain time point, the maximal value of timefrequency distributions is not the IF of signal, then the IF value of signal is that the probability of SPWVD nth (n is less) individual maximal value is very large;
The IF change of (2) 2 continuous time points is not very greatly.
Only seeing that Article 1 is supposed, is in fact spectrum blob detection method.So the crucial part Article 2 hypothesis of Viterbi algorithm, what this hypothesis reflected is the continuity that in rotating machinery variable parameter operation engineering, freerunning frequency changes, and can not occur very large frequency agility within the extremely short time.
Viterbi algorithm ascribes an extremevalue problem to, namely gets the time series [n in SPWVD distribution
_{1}, n
_{2}] between timefrequency distributions analyze, time point n ∈ [n
_{1}, n
_{2}], the frequency values on each time point n is SPW
_{n}, namely
SPW
_{n}＝[SPW(n，ω
_{1})，SPW(n，ω
_{2})，…，SPW(n，ω
_{M})]
M is that the frequency comprised is counted.
SPW
_{n}be SPW ' by the array define after descending sort
_{n}, then f (n, SPW '
_{n}) be defined as the timefrequency distributions sequence after assignment,
f(n,SPW′
_{n})＝j1(j＝0,1,…,M)
The meaning of above formula be exactly the value at all Frequency point places on SPWVD distribution time point n by order from big to small from 0 successively assignment by increment 1 assignment respectively.
Instantaneous frequency is the path making formula (2) minimum,
In formula, f (n, SPW '
_{n}) be construct according to hypothesis (1), g (ω
_{n}, ω
_{n+1}) according to hypothesis (2) structure.G (ω
_{n}, ω
_{n+1}) be continuity in order to meet frequency change and penalty factor that larger frequency change of adjusting the distance adds, can ensure that IF path meets above hypothesis (2).G (ω
_{n}, ω
_{n+1}) be defined as follows:
In formula, the optimal selection of δ is the maximum expected value that 2 continuity point IF change, and the cost function namely changing ( xy≤δ) for less IF is 0, chooses less δ (as δ=3) in an experiment and can obtain good effect.
With δ=1, c=5 is the implementation procedure that example illustrates Viterbi optimum route search algorithm, and get certain local time frequency plane and analyze, time frequency plane has 3 time points, 8 Frequency points, the concrete steps of Viterbi optimum route search algorithm are:
(1) according to formula f (n, SPW '
_{n})=j1 (j=0,1 ..., M) and calculate each point (n, ω
_{i}) corresponding f (n, ω
_{i}) functional value, determine the local optimum path of all frequencies of each time point n and all frequencies of its previous time point n1 and record, choose f (n, SPW '
_{n})+g (ω
_{n}, ω
_{n1}) minimum path as the local optimum path (left arrow line segment) between two time points, as shown in Fig. 1 (a);
(2) the local optimum path of all Frequency points in like manner, determining each time n and all Frequency points of a time point n+1 thereafter is also recorded (to the right arrow line segment), as shown in Fig. 1 (b);
(3) according to formula (3), with the local optimum path between two time points for benchmark, first along the functional value in all paths of local optimum path computing determined by step (1) straight line, the minimum path of record functional value is as shown in arrow left in Fig. 1 (c); Then the local optimum path searched for backward is determined after the same method, as shown in arrow to the right in Fig. 1 (c);
(4) compare to get in Fig. 1 (c) forward, the size of searching route functional value backward, the path selecting minimum function value corresponding is optimal path, as shown in Fig. 1 (d).
In addition, due to the edge effect of timefrequency distributions, there is very big error in the rim value that IF estimates, can remove some data of the right and left in time frequency distribution map.Simultaneously for the feature of multicomponents composition in actual condition middle gear case vibration signal, need to go out on the basis of an IF in abovementioned Viterbi algorithm search, estimate other IF.Suppose now to estimate
it is the signal that respective signal is the strongest.Will
value zero setting in scope, forms new timefrequency representation, and ξ is the constant of a setting.Above Viterbi searching algorithm step is repeated to new timefrequency representation, until estimated all to estimate by the IF of component each in signal.
2) order based on instantaneous Frequency Estimation is composed.Utilize key signal estimation model to carry out pointbypoint integration to instantaneous frequency and estimated key signal accurately.Then carry out calculating order analysis in conjunction with angle domain average technology to vibration acceleration signal, the order obtained based on instantaneous Frequency Estimation is composed.Specifically comprise following four steps:
21) utilizing key signal estimation model to step 13) instantaneous frequency that obtains carries out pointbypoint integration, obtains estimating key signal.On the basis obtaining IF, in order to computation order analysis, need to estimate key signal with IF.Suppose that a rotating shaft works under being in variable rotate speed, rotating speed is ω (t), as shown in Figure 2.The pass of ω (t) and rotation instantaneous frequency f (t) is ω (t)=2 π f (t).Getting the angle that certain time period in Fig. 2 turns over is Δ θ, and Δ θ can be expressed as:
Δθ＝ω(t)dt＝2πf(t)dt
Now suppose to add up from the measurement starting point t of a point as the time, every a bit of Δ θ, then elapsed time t
_{1}, rotor forwards again a point to, and the angle that rotor turns over is 2 π.Obtain following equation thus:
If through t
_{n}, rotating shaft turns over n circle, and rotating shaft turns over angle such as formula shown in (4):
In formula (4), t
_{n}it is the unknown quantity needing to solve.To f (t) pointbypoint integration, record integrated value is the t of integer
_{n}, at t
_{n}the key signal of a step is inserted in position, is namely the estimation key signal obtained by IF.
22) utilize the estimation key signal that obtains to carry out angularly resampling to vibration acceleration signal, obtain the signal after angularly resampling.
Utilize spindle keyway to mark or a pulse can only be collected in each turn of reflecting piece, angularly resampling needs the work done to be exactly the pulse utilizing a series of mark each turn of 2 π radians gathered, utilize corresponding computing method interpolation to be divided into N section, the radian of each section is
for calculating the angularly resampling moment, suppose that rotating shaft rotating shaft within the time period that each is little is made even angular acceleration and run.
The uniform angular velocity rotation equation of shaft rotary corner θ can be described as:
θ(t)＝b
_{0}+b
_{1}t+b
_{2}t
^{2}
In formula: b
_{0}, b
_{1}, b
_{2}for unknown parameter, t is time point.
Get the time point t of three continuous impulses that tachopulse sensor gathers
_{1}, t
_{2}, t
_{3}, with t
_{1}as the start index of angle step, then t
_{1}and t
_{2}, t
_{1}and t
_{3}angle step between time point is respectively
following equation can be listed:
B can be solved by equation
_{0}, b
_{1}, b
_{2}, shown in (6):
Meanwhile, by equation θ (t)=b
_{0}+ b
_{1}t+b
_{2}t
^{2}solve t:
Formula (6) is substituted into timing node t when formula (5) can obtain angularly θ resampling, then carry out angularly resampling, obtain the signal after angularly resampling.The schematic diagram of an angularly resampling is given shown in Fig. 3.In the drawings, in every two key signals collected every
individual insertion resample points (namely often turning sampling 8 points), these time points inserted, from time shaft, are not constant durations, but in angular domain read fortune every with equally spaced
23) angle domain average process is carried out to the signal after angularly resampling, obtain the signal after angle domain average process.
Burst after angularly resampling is defined as
identical with Synchronous time average thought, the algorithmic formula of angular domain synchronized averaging is formula (7), obtains the burst y after angle domain average process by formula (7)
_{n}.
Wherein,
P: the hop count that the signal averaging after angularly resampling is divided into;
N: each section identical sampling number.
Periodic mechanical signal sequences y (n) is had for what obtained by angular domain resampling in rotating machinery or reciprocating machine operational process, suppose that y (n) is made up of periodic signal f (n) and random noise s (n), that is:
x(n)＝f(n)+s(n)
Go to intercept signal x (n) with the cycle T of f (n), intercept altogether M section, be then added, due to the irrelevance of random noise, can obtain:
Again to x (n
_{i}) be averaging, obtain average after output signal y (n
_{i}):
Now random noise is the random noise in original input signal x (t)
thus signal to noise ratio (S/N ratio) improves
doubly.
Therefore, the wave filter of angular domain synchronized sampling is similar with the filter characteristic of time domain synchronized sampling, can stress release treatment and unrelated signal effectively.
24) Fast Fourier Transform (FFT) (FFT) is carried out to the signal after angular domain average treatment, namely obtain reflecting that the order based on instantaneous Frequency Estimation of gearbox fault characteristic information is composed.
Case study on implementation
The validity of this invention of signal authentication is moved with certain planetary reduction gear testing table institute vibration measuring.This testing table comprises motor, speed reduction unit, generator, control desk four part.Planetary reduction gear is a secondary precision speed reduction device, and reduction gear ratio is 32, and the structural parameters of planetary reduction gear are as shown in table 1 below.
Table 1 star reducer structure parameter
According to the parameter of table 1 and the meshing frequency of planetary gear meshing frequency computing formula calculating first order transmission sun gear and planetary gear.Planetary gear meshing frequency computing formula is as follows:
F
_{engagement}=(n
_{1}n
_{h}) × z
_{1}
In formula, n
_{1}: sun gear rotating speed; n
_{h}: planet carrier rotating speed; z
_{1}: the sun gear number of teeth.
F is obtained by above formula
_{engagement}=10.5n
_{1}.10.5 is important parameters, contrasts previously described order concept, J
_{n}=10.5 is planetary gear order parameters when breaking down.
Manual adjustments control box, makes the speed of input motor be a moderating process, in time domain, gathers vibration signal respectively as shown in Figure 4.Sample frequency is 10240Hz, meets sample frequency and arranges requirement.
Setting maximum analysis order according to feature order is 40, and 80 points of namely angularly sampling within the time that axle often circles, the sampling angle of equiangular sampling is 2 π/80.
Carry out to signal the zero phase filtering that cutoff frequency is 50Hz, filtered signal as shown in Figure 5.
SPWVD is utilized to ask for the timefrequency distributions of signal.Getting window function is Hamming window, and length of window is 1/10 of signal length, after obtaining timefrequency distributions, then removes the timefrequency distributions of 1s duration at timefrequency distributions two ends, to reduce end effect.The timefrequency distributions obtained is illustrated in fig. 6 shown below.Can see and timefrequency figure obviously has input shaft to turn 1 times frequently, the timefrequency distributions of 2 times.
Utilize Viterbi algorithm, timefrequency figure searches out 1 times of instantaneous turn of turning frequently respectively to be estimated in as Fig. 7 shown in heavy line frequently, 2 times turn instantaneous turn frequently and frequently estimate that (2 times of optimal paths turning IF frequently of extraction are divided by 2, then obtain 1 times of IF estimation turned frequently) as shown in Fig. 7 dotted line, true instantaneous frequency profile is as shown in Fig. 7 fine line.As can be seen from the figure the result that estimates of instantaneous turn of frequency and true instantaneous frequency profile can extraordinaryly coincide.
First turn IF frequently with 1 times to estimate, to conjunction key phase algorithm for estimating, namely can obtain the key signal estimated, as shown in Figure 8.
Utilize the key signal estimated, interpolation resampling is carried out to original signal and can obtain equal angular resampling signal, as shown in Figure 9.
Angle domain average process is carried out to the signal after angular domain resampling, and carries out FFT conversion., obtain order spectrum as shown in Figure 10.Spectrogram can be seen significantly relevant with the order 10.5 in the transmission of epicyclic gearbox first order spectrum peak, as figure J
_{n}shown in mark, both sides is with obvious fundamental frequency characteristic edge frequency band simultaneously, and namely sideband is rotating speed order 1.At 2 times of meshing frequency 2J
_{n}place, the frequency conversion band being spaced apart fundamental frequency is high and narrow, belongs to the even pitting fault feature of typical gear.Open this epicyclic gearbox, each flank of tooth of discovery first order transmission centre wheel has the even point corrosion of the flank of tooth shown in Figure 11 (in oval circle).
The ANCOT algorithm that this patent proposes, under the strong phase signals of nothing, the Accurate Diagnosis flank of tooth of gear case under the condition under variable rotate speed even pitting fault feature, improves the accuracy of fault diagnosis, ensures effectively reliably for the Fault Diagnosis of Gear Case under variable speed condition provides.
Claims (8)
1. gear distress is without a key phase angular domain average computation order analysis method, it is characterized in that, comprises the following steps:
1) instantaneous Frequency Estimation:
11) gather vibration acceleration signal, low pass is carried out to the vibration acceleration signal gathered and protects phase filtering, obtain low frequency signal x (t) reflecting rotor turns;
The timefrequency distributions of low frequency signal x (t) 12) calculated by smoothed pseudo wigner ville disstribution, obtains time frequency distribution map;
In the time frequency distribution map obtained, 13) utilize Viterbi optimum route search algorithm to estimate the instantaneous frequency of gear case rotating shaft;
2) order based on instantaneous Frequency Estimation is composed:
21) utilizing key signal estimation model to step 13) instantaneous frequency that obtains carries out pointbypoint integration, obtains estimating key signal;
22) utilize the estimation key signal that obtains to carry out angularly resampling to vibration acceleration signal, obtain the signal after angularly resampling;
23) angle domain average process is carried out to the signal after angularly resampling, obtain the signal after angle domain average process;
24) Fast Fourier Transform (FFT) is carried out to the signal after angular domain average treatment, obtain reflecting that the order based on instantaneous Frequency Estimation of gearbox fault characteristic information is composed, gear distress is analyzed.
2. gear distress according to claim 1 is without key phase angular domain average computation order analysis method, it is characterized in that: described step 12) concrete operations be:
Step 11) the smoothed pseudo wigner ville disstribution SPW (t, ω) of low frequency signal x (t) that obtains is by shown in formula (1):
Namely time frequency distribution map is obtained by formula (1); In formula (1), * represents conjugation, and t represents the time, and ω represents frequency, and h (τ) and g (u) is window function, and τ, u represent the shift factor of time shaft and frequency axis respectively.
3. gear distress according to claim 2 is without key phase angular domain average computation order analysis method, it is characterized in that: described window function comprises Gauss window and Hamming window.
4. according to the gear distress in claim 13 described in any one without key phase angular domain average computation order analysis method, it is characterized in that: described step 13) concrete operations be:
Get the time series [n on the smoothed pseudo wigner ville disstribution of x (t)
_{1}, n
_{2}], time point n ∈ [n
_{1}, n
_{2}], the frequency values on each time point n is SPW
_{n}, SPW
_{n}by the array after descending sort for being designated as SPW '
_{n};
From 0 s, press increment 1 successively by order from big to small and distinguish assignment again, the timefrequency distributions sequence after assignment be designated as f (n, SPW '
_{n});
Instantaneous frequency is the path making formula (2) minimum:
Wherein, g (ω
_{n}, ω
_{n+1}) be continuity in order to meet frequency change and the penalty factor that larger frequency change of adjusting the distance adds.
5. gear distress according to claim 4 is without key phase angular domain average computation order analysis method, it is characterized in that: described penalty factor g (ω
_{n}, ω
_{n+1}) such as formula shown in (3):
Wherein δ is the maximum expected value of the Instantaneous frequency variations of 2 continuity points.
6. gear distress according to claim 4 is without key phase angular domain average computation order analysis method, it is characterized in that: described step 21) concrete operations be:
Work under a rotating shaft is in variable rotate speed, rotation instantaneous frequency is f (t); Elapsed time t
_{n}, rotating shaft turns over n circle, and rotating shaft turns over angle such as formula shown in (4):
To f (t) pointbypoint integration, record integrated value is the t of integer
_{n}, at t
_{n}the key signal of a step is inserted in position, is namely the estimation key signal obtained by instantaneous frequency.
7. gear distress according to claim 4 is without key phase angular domain average computation order analysis method, it is characterized in that: described step 22) concrete operations be:
Timing node t when obtaining angularly θ resampling according to formula (5), carries out angularly resampling, obtains the signal after angularly resampling;
B in formula (5)
_{0}, b
_{1}, b
_{2}obtain according to formula (6);
Wherein
8. gear distress according to claim 4 is without key phase angular domain average computation order analysis method, it is characterized in that: described step 23) concrete operations be:
Burst after angularly resampling is designated as
i=1,2,3 ..., N, obtains the burst y after angle domain average process by the formula (7) of angle domain average process
_{n},
Wherein, P is the hop count signal averaging after angularly resampling be divided into, and N is each section of identical sampling number.
Priority Applications (1)
Application Number  Priority Date  Filing Date  Title 

CN201310416504.7A CN103499443B (en)  20130912  20130912  A kind of gear distress is without key phase angular domain average computation order analysis method 
Applications Claiming Priority (1)
Application Number  Priority Date  Filing Date  Title 

CN201310416504.7A CN103499443B (en)  20130912  20130912  A kind of gear distress is without key phase angular domain average computation order analysis method 
Publications (2)
Publication Number  Publication Date 

CN103499443A CN103499443A (en)  20140108 
CN103499443B true CN103499443B (en)  20160120 
Family
ID=49864674
Family Applications (1)
Application Number  Title  Priority Date  Filing Date 

CN201310416504.7A CN103499443B (en)  20130912  20130912  A kind of gear distress is without key phase angular domain average computation order analysis method 
Country Status (1)
Country  Link 

CN (1)  CN103499443B (en) 
Families Citing this family (17)
Publication number  Priority date  Publication date  Assignee  Title 

CN104266747B (en) *  20140609  20171024  中能电力科技开发有限公司  A kind of method for diagnosing faults based on vibration signal order analysis 
CN104077474A (en) *  20140623  20141001  华南理工大学  Meshing frequency and spectrum correction technology based wind power gear box order tracking method 
CN104655380B (en) *  20150316  20171024  北京六合智汇技术有限责任公司  A kind of rotating machinery fault signature extracting method 
CN104819841B (en) *  20150505  20170419  西安交通大学  Builtincodinginformationbased single sensing flexible angledomain averaging method 
CN105277362B (en) *  20151123  20170912  西安交通大学  Gear distress detection method based on encoder multidigit angular signal 
CN105547698B (en) *  20151231  20171114  新疆金风科技股份有限公司  The method for diagnosing faults and device of rolling bearing 
US10094743B2 (en) *  20160314  20181009  Epro Gmbh  Order analysis system 
CN106644467B (en) *  20161227  20190514  华南理工大学  A kind of gearbox nonstationary signal fault signature extracting method 
CN107292067B (en) *  20170817  20210202  湖南纬拓信息科技有限公司  Gear fault diagnosis method based on compressed sensing and bispectrum analysis 
CN107907324A (en) *  20171017  20180413  北京信息科技大学  A kind of Fault Diagnosis of Gear Case method composed based on DTCWT and order 
CN107941510B (en) *  20171019  20190719  西安交通大学  Extracting method based on the angularly Rolling Bearing Fault Character of dual sampling 
CN108225764A (en) *  20171205  20180629  昆明理工大学  It is a kind of based on the highprecision of envelope extraction without key signal Order Tracking and system 
CN108362942B (en) *  20180211  20201120  中国铁道科学研究院集团有限公司  Timefrequency spectrum obtaining method and device for multicomponent signals 
CN108871742A (en) *  20180503  20181123  西安交通大学  A kind of improved no key phase fault feature order extracting method 
CN109813546B (en) *  20190320  20200915  苏州微著设备诊断技术有限公司  Offline detection method for abnormal knocking sound of gear box 
CN110686768A (en) *  20191017  20200114  昆明理工大学  Improved rotating machinery nonstationary vibration signal calculation order ratio analysis method 
CN110907162A (en) *  20191213  20200324  北京天泽智云科技有限公司  Rotating machinery fault feature extraction method without tachometer under variable rotating speed 
Citations (7)
Publication number  Priority date  Publication date  Assignee  Title 

CN102175393A (en) *  20110104  20110907  广东电网公司电力科学研究院  Imbalance phase estimation method based on procession decomposition technology 
US8024120B2 (en) *  20080516  20110920  Turner Larry A  Complex phase locked loop 
FR2967541A1 (en) *  20101117  20120518  Centre Nat Etd Spatiales  Oversampled digital signal i.e. split phaselevel/phase modulation signal, demodulating method, involves calculating frequency estimation from channel estimation, extracting data, and compensating carrier phase with channel estimation 
CN102998110A (en) *  20121129  20130327  西安交通大学  Rotary machine fault characteristic extraction method based on orderholospectrum principle 
CN103196547A (en) *  20130311  20130710  安徽新力电业科技咨询有限责任公司  Method achieving rotary machinery vibration signal synchronization order ratio tracing analysis 
CN103234627A (en) *  20130417  20130807  国家电网公司  Complete alternation synchronous sampling and analyzing method for rotating machinery vibration signals 
CN103278235A (en) *  20130603  20130904  合肥伟博测控科技有限公司  Novel transient oscillation signal angular domain order tracking sampling and analytical method 

2013
 20130912 CN CN201310416504.7A patent/CN103499443B/en active IP Right Grant
Patent Citations (7)
Publication number  Priority date  Publication date  Assignee  Title 

US8024120B2 (en) *  20080516  20110920  Turner Larry A  Complex phase locked loop 
FR2967541A1 (en) *  20101117  20120518  Centre Nat Etd Spatiales  Oversampled digital signal i.e. split phaselevel/phase modulation signal, demodulating method, involves calculating frequency estimation from channel estimation, extracting data, and compensating carrier phase with channel estimation 
CN102175393A (en) *  20110104  20110907  广东电网公司电力科学研究院  Imbalance phase estimation method based on procession decomposition technology 
CN102998110A (en) *  20121129  20130327  西安交通大学  Rotary machine fault characteristic extraction method based on orderholospectrum principle 
CN103196547A (en) *  20130311  20130710  安徽新力电业科技咨询有限责任公司  Method achieving rotary machinery vibration signal synchronization order ratio tracing analysis 
CN103234627A (en) *  20130417  20130807  国家电网公司  Complete alternation synchronous sampling and analyzing method for rotating machinery vibration signals 
CN103278235A (en) *  20130603  20130904  合肥伟博测控科技有限公司  Novel transient oscillation signal angular domain order tracking sampling and analytical method 
NonPatent Citations (1)
Title 

无转速计旋转机械升降速振动信号零相位阶比跟踪滤波;郭瑜等;《机械工程学报》;20040331;第40卷(第3期);5054 * 
Also Published As
Publication number  Publication date 

CN103499443A (en)  20140108 
Similar Documents
Publication  Publication Date  Title 

Feng et al.  Iterative generalized synchrosqueezing transform for fault diagnosis of wind turbine planetary gearbox under nonstationary conditions  
Abboud et al.  The spectral analysis of cyclononstationary signals  
Borghesani et al.  A new procedure for using envelope analysis for rolling element bearing diagnostics in variable operating conditions  
Chow et al.  Induction machine fault diagnostic analysis with wavelet technique  
Chen et al.  Generator bearing fault diagnosis for wind turbine via empirical wavelet transform using measured vibration signals  
US20160033580A1 (en)  Detecting Faults in Turbine Generators  
Wilkinson et al.  Condition monitoring of generators & other subassemblies in wind turbine drive trains  
Li et al.  New procedure for gear fault detection and diagnosis using instantaneous angular speed  
CN102937668B (en)  Electric system lowfrequency oscillation detection method  
CN102901929B (en)  For calculating device and the battery impedance measuring system of cell impedance  
CN101221066B (en)  Engineering nonlinear vibration detecting method  
CN103321854B (en)  A kind of vibration control method for wind generator set tower  
Li et al.  A new noisecontrolled secondorder enhanced stochastic resonance method with its application in wind turbine drivetrain fault diagnosis  
Wang et al.  Rolling element bearing fault diagnosis via fault characteristic order (FCO) analysis  
Zhao et al.  A tacholess order tracking technique for large speed variations  
CN102435844B (en)  Sinusoidal signal phasor calculating method being independent of frequency  
He et al.  Multiscale stochastic resonance spectrogram for fault diagnosis of rolling element bearings  
Tsao et al.  An insight concept to select appropriate IMFs for envelope analysis of bearing fault diagnosis  
Sawalhi et al.  Gear parameter identification in a wind turbine gearbox using vibration signals  
Blunt et al.  Detection of a fatigue crack in a UH60A planet gear carrier using vibration analysis  
CN102483368B (en)  Method for detecting structural defect in mechanical assembly including rotary member  
CN102156043B (en)  Online state monitoring and fault diagnosis system of wind generator set  
Feng et al.  Fault diagnosis of wind turbine planetary gearbox under nonstationary conditions via adaptive optimal kernel time–frequency analysis  
Barszcz et al.  Application of spectral kurtosis for detection of a tooth crack in the planetary gear of a wind turbine  
Feng et al.  Joint envelope and frequency order spectrum analysis based on iterative generalized demodulation for planetary gearbox fault diagnosis under nonstationary conditions 
Legal Events
Date  Code  Title  Description 

C06  Publication  
PB01  Publication  
C10  Entry into substantive examination  
SE01  Entry into force of request for substantive examination  
C14  Grant of patent or utility model  
GR01  Patent grant 