CN105606894B - Instantaneous Frequency Estimation method based on simulated annealing - Google Patents
Instantaneous Frequency Estimation method based on simulated annealing Download PDFInfo
- Publication number
- CN105606894B CN105606894B CN201610060143.0A CN201610060143A CN105606894B CN 105606894 B CN105606894 B CN 105606894B CN 201610060143 A CN201610060143 A CN 201610060143A CN 105606894 B CN105606894 B CN 105606894B
- Authority
- CN
- China
- Prior art keywords
- frequency
- moment
- time
- instantaneous frequency
- simulated annealing
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Expired - Fee Related
Links
Classifications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01R—MEASURING ELECTRIC VARIABLES; MEASURING MAGNETIC VARIABLES
- G01R23/00—Arrangements for measuring frequencies; Arrangements for analysing frequency spectra
- G01R23/16—Spectrum analysis; Fourier analysis
Landscapes
- Physics & Mathematics (AREA)
- Mathematical Physics (AREA)
- General Physics & Mathematics (AREA)
- Measurement Of Mechanical Vibrations Or Ultrasonic Waves (AREA)
Abstract
The invention discloses the instantaneous Frequency Estimation methods based on simulated annealing, introduce simulated annealing (SA), it are combined with Short Time Fourier Transform (STFT), it is proposed that STFT-SA instantaneous frequency estimation algorithms.STFT-SA algorithms are based on signal time-frequency spectrum, seek the thought of optimal path using simulated annealing, in conjunction with the characteristics of rotating machinery non-stationary signal, the single order instantaneous frequency distilling for realizing characteristic of rotating machines vibration signal all has extraordinary effect for vibration signals such as strong noise, neighbouring rank ratios.
Description
Technical field
The invention belongs to rotating machinery, fault diagnosis technology field, more particularly to based on the instantaneous of simulated annealing
Frequency estimating methods.
Background technology
Rotating machinery has been part indispensable in industrial production, in order to ensure the normal operation of rotating machinery, to rotation
Tool of making a connection carries out state-detection and fault diagnosis is particularly important.The vibration signal in rotating machinery lifting speed stage contains abundant
Status information is studied the system defect for more easily finding to be generally difficult to find to this stage.Therefore it can effectively carry
The fault signature for taking this stage is the important link of rotary machinery fault diagnosis, has weight to the normal operation of rotating machinery
Want meaning.
In order to be best understood from the fault message of rotating machinery, the various methods for instantaneous Frequency Estimation are suggested.Peak
Value search method, hidden Markov combination peak searching algorithm and improvement peak searching algorithm etc. are provided to solve noise to instantaneous frequency
The problem of rate estimation influences.But rotating machinery is in certain circumstances, such as neighbouring rank than generation, instantaneous Frequency Estimation
Difficulty increase.STFT-VF(Short-time Fourier transform_Viterbi algorithm fit,STFT-
VF) there is degree of precision in strong noise and neighbouring rank ratio.The neighbouring rank score of quadrature compensation matrix combination discrete Fourier
Analysis method and the neighbouring order ratio analysis method of independent component analysis combination peak value searching can efficiently extract out neighbouring rank ratio
Instantaneous frequency, but in place of both methods all Shortcomings.Such as the neighbouring rank ratio of quadrature compensation matrix combination discrete Fourier
Annoyance level is affected between variation and neighbouring rank ratio of the analysis method by rank than amplitude.And independent component analysis separation letter
Number when required source signal number must be no less than the number of isolated component, this makes the neighbour of independent component analysis combination peak value searching
There is limitation in actual use in nearly order ratio analysis method.
Invention content
In order to solve the technical issues of above-mentioned background technology proposes, the present invention is intended to provide the wink based on simulated annealing
When frequency estimating methods, overcome the prior art neighbouring rank than instantaneous Frequency Estimation the problem of.
In order to achieve the above technical purposes, the technical scheme is that：
Instantaneous Frequency Estimation method based on simulated annealing, includes the following steps：
(1) STFT transformation is carried out to the vibration signal of acquisition, and obtains the time-frequency spectrum of signal；
(2) original state is determined on the time-frequency spectrum that step (1) obtains；
(3) time-frequency spectrum obtained according to step (1) finds out all candidates' near t moment frequency values point
Point where t-1 moment frequency values, when the energy value at these candidate t-1 moment frequency values points is replaced t-1 successively
The value for carving original state, generates new state；
(4) energy for the new state for generating step (3) obtains the increment of state change compared with original state, selects
Acceptance criterion receives new state, finds out the maximum new state of energy value and the t-1 moment frequency of the new state corresponding candidate
Point where rate value, using the point as the estimated value of t-1 moment frequency values points；
(5) return to step (3) find out t-2 moment frequency values points according to the estimated value of t-1 moment frequency values points
Estimated value, regular cycles step (3)-(4) according to this, to estimate the instantaneous frequency at all moment；
(6) instantaneous frequency that step (5) estimates is fitted, eliminates local error.
A kind of preferred embodiment based on the above-mentioned technical proposal, the detailed process of step (1)：
The very short window function of a time width is given, window function changes over time intercept signal, then does such as formula (1) institute
The Fourier transform shown finds out STFT frequency spectrums further according to formula (2), and STFT frequency spectrums are the time-frequency energy distribution of STFT：
P (t, f)=| Sx(t,f)|2 (2)
Wherein, x (t) is vibration signal, and γ (t) is window function, and * represents complex conjugate, and f is time variable, and t and τ are equivalent
Time variable.
A kind of preferred embodiment based on the above-mentioned technical proposal, in step (2), on the time-frequency spectrum that step (1) obtains
The summation for the energy value that single order rotating speed final moment corresponding Frequency point is expert at is found as original state.
A kind of preferred embodiment based on the above-mentioned technical proposal is connect in step (4) using Metropolis acceptance criterions
By new state.
A kind of preferred embodiment based on the above-mentioned technical proposal estimates step (5) using least square method in step (6)
The instantaneous frequency counted out is fitted.
The advantageous effect brought using above-mentioned technical proposal：
The present invention introduces the thought that optimal path is sought in simulated annealing based on the time-frequency spectrum of vibration signal,
In conjunction with the characteristics of rotating machinery lifting speed stage non-stationary signal, find from initial time to the local energy at final moment
It is worth maximum path, the instantaneous frequency of reference axis is obtained finally by fitting, small, the fast and with high accuracy spy of speed with calculation amount
Point.
Description of the drawings
Fig. 1 is flow chart of the method for the present invention；
Fig. 2 is the location estimation schematic diagram of t-1 moment Frequency points；
Fig. 3 is the time-frequency figure for emulating signal；
Fig. 4 is the single order instantaneous frequency figure of STFT_SA algorithms estimation；
Fig. 5 is the first order component instantaneous frequency figure for emulating signal；
Fig. 6 is that big end vertical direction axis shakes ramp-up stage time-domain signal figure；
Fig. 7 is that big end vertical direction axis shakes ramp-up stage vibration signal time-frequency spectrum；
Fig. 8 is STFT-SA instantaneous Frequency Estimation figures；
Fig. 9 is the instantaneous frequency figure of the first order component acquired by tach signal.
Specific implementation mode
Below with reference to attached drawing, technical scheme of the present invention is described in detail.
Flow chart of the method for the present invention as shown in Figure 1, the instantaneous Frequency Estimation method based on simulated annealing, including with
Lower step：
Step 1 carries out STFT transformation to the vibration signal of acquisition, and obtains the time-frequency spectrum of signal.
The very short window function of a time width is given, window function changes over time intercept signal, then does such as formula (1) institute
The Fourier transform shown finds out STFT frequency spectrums further according to formula (2), and STFT frequency spectrums are the time-frequency energy distribution of STFT：
P (t, f)=| Sx(t,f)|2 (2)
Wherein, x (t) is vibration signal, and γ (t) is window function, and * represents complex conjugate, and f is time variable, and t and τ are equivalent
Time variable.
Step 2 finds the energy that single order rotating speed final moment corresponding Frequency point is expert on obtained time-frequency spectrum
The summation of value is as original state.
It is to refer to the corresponding instantaneous frequency of rotating speed for need to extract in rotating machinery order ratio analysis, this process is total
There is the stage of a stabilized (steady-state) speed.Such as it to the rotating machinery that a rated speed is 3300RPM, is carried out by taking its ramp-up stage as an example
Analysis, the single order rotating speed final moment, corresponding frequency should be near 55Hz on time-frequency spectrum.When therefore can select final
It is original state to carve the energy value summation that the corresponding points of 55Hz are expert at.
Step 3, according to time-frequency spectrum, the t-1 moment frequencies of all candidates are found out near t moment frequency values point
Energy value at these candidate t-1 moment frequency values points is replaced t-1 moment original states by the point where being worth successively
Value, generates new state.
When seeking the spectrogram of signal with STFT, a segment signal is intercepted by window function translation and its Fourier is asked to become
The signal level of coverage for changing, and being intercepted before and after translation reaches 80% or so, so that the same single order of t → t-1 is than signal
Frequency values show very close on time-frequency figure, that is to say, that the point where the frequency values at t-1 moment is very close to t moment frequency
The point of value, that is, the frequency values at t-1 moment point in the left side of point of t moment frequency values, the upper left corner and lower-left angular direction
One in point, if as shown in Fig. 2, the instantaneous frequency of t moment is putting 1, then with single order than the instantaneous frequency at middle t-1 moment
Certainly some in 2,3,4 three points.
Step 4, the new state for generating step 3 energy compared with original state, obtain the increment of state change, select
Acceptance criterion receives new state, finds out the maximum new state of energy value and the t-1 moment frequency of the new state corresponding candidate
Point where rate value, using the point as the estimated value of t-1 moment frequency values points.
Step 5, return to step 3 find out t-2 moment frequency values place according to the estimated value of t-1 moment frequency values points
The estimated value of point, regular cycles step 3-4 according to this, to estimate the instantaneous frequency at all moment.
First, for emulating signal.
Emulate a signal, sample frequency 5kHZ, sampling time 5s, reference frequency linear change from 20Hz to 100Hz.
The frequency representation of n-th of sampling instant is f (n)=20+80*n/N, and N is sampling number, the angular velocimeter of n-th of sampling instant
The π f (n) of w (n)=2 are shown as,For the phase angle that n-th of sampling instant turns over, Δ t is between the sampling time
Every.The wherein amplitude of fundamental frequency signal saltus step at 500, as shown in formula (3).
The amplitude of second order frequency signal is sinusoidal variations, as shown in formula (4).
B (n)=0.8+sin (2 π n/N) (4)
The amplitude of three order frequency signals is in 0 to 1 linear change.The multi -components emulation signal that then three frequency contents are set up is such as
Shown in formula (5).
η (n) is 30% white Gaussian noise in formula.
Short time discrete Fourier transform (STFT) is done to signal, it is as shown in Figure 3 to obtain time-frequency spectrum.It is used after determining original state
Simulated annealing carries out instantaneous Frequency Estimation to the first order component of signal, and the results are shown in Figure 4.By the first order component of estimation
The instantaneous frequency of instantaneous frequency and the first order component of emulation compares, i.e. comparison diagram 4 and Fig. 5, can obtain the degree of fitting of the two
It is 99.4%.0.1s is only consumed to the time spent in instantaneous Frequency Estimation of first order component.
Then, by taking actual vibration signal as an example.
It is being worked horizontal spiral centrifuge (abbreviation decanter centrifuge) using the dynamic signal analyzer of OROS R3X systems
The vibration signal of rotating speed operation phase carries out test experiments.Sample frequency is 12.8kHz, sampling time 20s.Take big end vertical
Axis of orientation shakes for ramp-up stage vibration signal, and Fig. 6 is its time-domain signal；Spectrum analysis, spectrum analysis are carried out to vibration signal
The results are shown in Figure 7.Single order instantaneous frequency is carried out to decanter centrifuge vibration signal using STFT-SA instantaneous frequency estimation algorithms
Estimation, the results are shown in Figure 8.The instantaneous frequency for acquiring 1 order component according to tach signal simultaneously, as shown in figure 9, by itself and Fig. 8
It is compared, acquires the degree of fitting of the instantaneous frequency that STFT-SA estimates and the reference axis instantaneous frequency that pulse signal is converted to
It is 99.1%.And the time of 1s is only spent to the instantaneous Frequency Estimation of first order component.
Above example is merely illustrative of the invention's technical idea, and protection scope of the present invention cannot be limited with this, every
According to technological thought proposed by the present invention, any change done on the basis of technical solution each falls within the scope of the present invention
Within.
Claims (5)
1. the instantaneous Frequency Estimation method based on simulated annealing, which is characterized in that include the following steps：
(1) STFT transformation is carried out to the vibration signal of acquisition, and obtains the time-frequency spectrum of signal；
(2) original state is determined on the time-frequency spectrum that step (1) obtains；
(3) time-frequency spectrum obtained according to step (1), when finding out the t-1 of all candidates near t moment frequency values point
The point where frequency values is carved, replaces the t-1 moment initial successively the energy value at these candidate t-1 moment frequency values points
The value of state generates new state；
(4) energy for the new state for generating step (3) obtains the increment of state change compared with original state, and selection receives
Criterion receives new state, finds out the t-1 moment frequency values of the maximum new state of energy value and the corresponding candidate of the new state
The point at place, using the point as the estimated value of t-1 moment frequency values points；
(5) return to step (3) find out estimating for t-2 moment frequency values points according to the estimated value of t-1 moment frequency values points
Evaluation, regular cycles step (3)-(4) according to this, to estimate the instantaneous frequency at all moment；
(6) instantaneous frequency that step (5) estimates is fitted, eliminates local error.
2. the instantaneous Frequency Estimation method based on simulated annealing according to claim 1, which is characterized in that step (1)
Detailed process：
The very short window function of a time width is given, window function changes over time intercept signal, then does as shown in formula (1)
Fourier transform finds out STFT frequency spectrums further according to formula (2), and STFT frequency spectrums are the time-frequency energy distribution of STFT：
P (t, f)=| Sx(t,f)|2 (2)
Wherein, x (t) is vibration signal, and γ (t) is window function, and * represents complex conjugate, and f is time variable, when t and τ are equivalent
Between variable.
3. the instantaneous Frequency Estimation method based on simulated annealing according to claim 1, it is characterised in that：In step
(2) in, the energy value that single order rotating speed final moment corresponding Frequency point is expert at is found on the time-frequency spectrum that step (1) obtains
Summation as original state.
4. the instantaneous Frequency Estimation method based on simulated annealing according to claim 1, it is characterised in that：In step
(4) in, new state is received using Metropolis acceptance criterions.
5. the instantaneous Frequency Estimation method based on simulated annealing according to claim 1, it is characterised in that：In step
(6) in, the instantaneous frequency estimated to step (5) using least square method is fitted.
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201610060143.0A CN105606894B (en) | 2016-01-28 | 2016-01-28 | Instantaneous Frequency Estimation method based on simulated annealing |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN201610060143.0A CN105606894B (en) | 2016-01-28 | 2016-01-28 | Instantaneous Frequency Estimation method based on simulated annealing |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN105606894A CN105606894A (en) | 2016-05-25 |
| CN105606894B true CN105606894B (en) | 2018-11-02 |
Family
ID=55986982
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN201610060143.0A Expired - Fee Related CN105606894B (en) | 2016-01-28 | 2016-01-28 | Instantaneous Frequency Estimation method based on simulated annealing |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN105606894B (en) |
Families Citing this family (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN106370403A (en) * | 2016-08-22 | 2017-02-01 | 南京信息工程大学 | Instant frequency estimation method based on edge detection |
| CN107622035B (en) * | 2017-09-30 | 2020-07-17 | 中国人民解放军战略支援部队航天工程大学 | Polynomial phase signal self-adaptive time-frequency transformation method based on simulated annealing |
| CN108680824A (en) * | 2018-05-16 | 2018-10-19 | 南方电网科学研究院有限责任公司 | Distributed Wave Recording Synchronization Method, Device, Equipment and Medium |
| CN115225450A (en) * | 2022-09-20 | 2022-10-21 | 南京艾泰克物联网科技有限公司 | Multi-data machine room virtualization cluster management system based on edge computing |
Citations (6)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN101603854A (en) * | 2009-07-15 | 2009-12-16 | 南京信息工程大学 | The rotating machinery non-stationery vibration signal instantaneous frequency estimation algorithm in start and stop period |
| US8127138B1 (en) * | 2008-09-29 | 2012-02-28 | The United States Of America As Represented By The Secretary Of The Navy | Method for embedding information in sonar |
| CN102650658A (en) * | 2012-03-31 | 2012-08-29 | 机械工业第三设计研究院 | Time-varying non-stable-signal time-frequency analyzing method |
| CN104407393A (en) * | 2014-12-08 | 2015-03-11 | 中国石油天然气集团公司 | Time frequency electromagnetic-based adaptive genetic simulated annealing inversion method and system |
| CN104749432A (en) * | 2015-03-12 | 2015-07-01 | 西安电子科技大学 | Estimation method of multi-component non-stationary signal instantaneous frequency based on focusing S-transform |
| CN105072440A (en) * | 2015-08-27 | 2015-11-18 | 中国电子科技集团公司第二十六研究所 | Method for extracting transient signal parameters by acousto-optical spectrum analyzer |
-
2016
- 2016-01-28 CN CN201610060143.0A patent/CN105606894B/en not_active Expired - Fee Related
Patent Citations (6)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US8127138B1 (en) * | 2008-09-29 | 2012-02-28 | The United States Of America As Represented By The Secretary Of The Navy | Method for embedding information in sonar |
| CN101603854A (en) * | 2009-07-15 | 2009-12-16 | 南京信息工程大学 | The rotating machinery non-stationery vibration signal instantaneous frequency estimation algorithm in start and stop period |
| CN102650658A (en) * | 2012-03-31 | 2012-08-29 | 机械工业第三设计研究院 | Time-varying non-stable-signal time-frequency analyzing method |
| CN104407393A (en) * | 2014-12-08 | 2015-03-11 | 中国石油天然气集团公司 | Time frequency electromagnetic-based adaptive genetic simulated annealing inversion method and system |
| CN104749432A (en) * | 2015-03-12 | 2015-07-01 | 西安电子科技大学 | Estimation method of multi-component non-stationary signal instantaneous frequency based on focusing S-transform |
| CN105072440A (en) * | 2015-08-27 | 2015-11-18 | 中国电子科技集团公司第二十六研究所 | Method for extracting transient signal parameters by acousto-optical spectrum analyzer |
Non-Patent Citations (1)
| Title |
|---|
| Crazy Climber算法与重采样技术在消除多普勒效应及列车轴承诊断中的应用;王超 等;《振动工程学报》;20141031;第27卷(第5期);770-779 * |
Also Published As
| Publication number | Publication date |
|---|---|
| CN105606894A (en) | 2016-05-25 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN105606894B (en) | Instantaneous Frequency Estimation method based on simulated annealing | |
| CN104266747B (en) | A kind of method for diagnosing faults based on vibration signal order analysis | |
| CN110763462B (en) | Time-varying vibration signal fault diagnosis method based on synchronous compression operator | |
| CN105866250B (en) | Ventilating vane method for crack based on vibration | |
| CN106370403A (en) | Instant frequency estimation method based on edge detection | |
| CN110243605A (en) | Multi-source time-frequency crestal line extracting method | |
| CN106248356B (en) | A kind of rotary machinery fault diagnosis method based on kurtosis index | |
| CN103259485B (en) | Method of improving identification precision of speedless sensor under condition of unbalanced network voltage | |
| CN102967414B (en) | Method for extracting imbalanced components of micro-speed-difference double-rotor system based on frequency spectrum correction | |
| CN108871742A (en) | A kind of improved no key phase fault feature order extracting method | |
| Chen et al. | Nonstationary signal denoising using an envelope-tracking filter | |
| He et al. | Feature extraction of acoustic signals based on complex Morlet wavelet | |
| CN114061678B (en) | Digital driving method for Coriolis flowmeter | |
| CN105510066B (en) | A kind of rotatory mechanical system method for diagnosing faults based on adaptive noise reduction algorithm | |
| CN105865793A (en) | Method for improving vibration monitoring precision of rotor aeroengine | |
| Kabla et al. | Bearing fault diagnosis using Hilbert-Huang transform (HHT) and support vector machine (SVM) | |
| CN111738068A (en) | Transmission shaft fault diagnosis method and system under working condition of rotating speed fluctuation | |
| Ding et al. | Multiple instantaneous frequency ridge based integration strategy for bearing fault diagnosis under variable speed operations | |
| Geng et al. | Fault identification of rolling bearing with variable speed based on generalized broadband mode decomposition and distance evaluation technique | |
| WO2015083397A1 (en) | Calculation device | |
| CN113586177B (en) | Blade natural frequency identification method based on single-blade-end timing sensor | |
| Song et al. | Multispectral balanced automatic fault diagnosis for rolling bearings under variable speed conditions | |
| CN104076203B (en) | A kind of intrasonic harmonic detection method considering that negative frequency affects | |
| CN109245655A (en) | The mover initial position detection method of the planar motor of position-sensor-free | |
| Zarei et al. | Broken rotor bars detection via Park's vector approach based on ANFIS |
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: 20181102 Termination date: 20220128 |
|
| CF01 | Termination of patent right due to non-payment of annual fee |