CN104568444A - Method for extracting fault characteristic frequencies of train rolling bearings with variable rotational speeds - Google Patents

Abstract

The invention discloses a method for extracting fault characteristic frequencies of train rolling bearings with variable rotational speeds, and belongs to the field of technologies for diagnosing faults and processing signals. The method includes steps of analyzing vibration signals of the bearings with the variable rotational speeds in time and frequency domains, searching local peak values of the vibration signals and extracting instantaneous frequency values corresponding to the rotational speeds at different moments; fitting instantaneous frequencies by the aid of neural networks, acquiring rotational speed curves of reference spindles, re-sampling original signals at uniform angles on the basis of the rotational speed curves and analyzing order ratios of the original signals on the basis of the rotational speed curves; separating signals with mixed order ratio signals by the aid of fixed-point independent component analysis and spectrum peak search technologies to acquire order ratio component characteristics of fault components of the bearings. The method has the advantages that the method is used for estimating the rotational speeds of the train bearings without tachometers in real time, the instable fault bearing signals can be converted into the stable signals in uniform-angle domains, independent order ratio components can be effectively separated from the signals, and the method is favorable for extracting the fault characteristic frequencies of the train bearings and detecting the fault characteristic frequencies of the train bearings in an online manner.

Description

Variable speed train Rolling Bearing Fault Character frequency extraction method

Technical field

The present invention relates to plant equipment Acoustic Based Diagnosis technical field, particularly variable speed state to detrain Rolling Bearing Fault Character frequency extraction method.

Background technology

Safety is the eternal theme of transportation by railroad.China now has railway freight-car more than 60 ten thousand, rolling bearing is as important component of machine wherein, due to long-term heavy duty, touch mill operation, very easily there is the faults such as sur-face peeling, light then cause withdrawal of train late, heavy then cause hot box, fire axle, cut axle, the serious accident such as derailing is overturned, on-line fault diagnosis is implemented to it most important.According to rolling stock section's finding, bearing maintenance process passes through hand-turning bearing outer ring primarily of experience workman, the mode of listening with people's ear judges whether it has different sound or fault, the method affects greatly by subjective factor, be unfavorable for that industry standardization manages, and workman, in order to reduce the liability exposure of oneself, takes to judge by accident and mode of not failing to judge sometimes, cause that false alarm rate is too high, overhaul efficiency is low.

For rolling bearing fault diagnosis problem, develop multiple diagnostic techniques.Divide by the physical property of signal and mainly comprise IR thermometry, oil analysis method, oil film thickness analytic approach, vibration and acoustic method.IR thermometry has simply, be easy to the feature that realizes, be applied to Truck Train Inspection, but temperature rise belongs to rolling bearing fault iate feature, and for early stage spot corrosion, the minor failure such as to peel off and then cannot detect, there is larger detection risk, pre-alerting ability is weak.Oil analyzing technology judges bearing working state by the physicochemical property and contained metal worn particle size, pattern and concentration analyzing lubricant grease itself, can be used for early diagnosis, but there is sampling inconvenience, poor real shortcoming.Oil film thickness analysis judges lubricating status by measuring oil film resistance, and the detectability of the faults such as effects on surface peels off, crackle is more weak, and is not suitable for the situations such as low speed, turning axle do not expose.Vibration or acoustics diagnose method relatively ripe, obtain extensive investigation and application, axis of rolling support rail limit acoustics diagnose system (the Trackside Acoustic Detection System that ripe commercial system has the U.S. to roll, TADS), this system has higher Detection accuracy, but it can only detect the bearing that fault is serious, and it is expensive, single measuring point needs about 600,000 dollars (disposing 60 measuring points altogether), renewal of the equipment expense after also not comprising later maintenance and several years.

Characteristic frequency is the important evidence of rolling bearing acoustic fault recognition and classification.Current main stream approach utilizes resonance and demodulation (Resonance Demodulation) technology to extract fault characteristic frequency, namely adopt bandpass filter by the resonant frequency of the weak impact signal madulation of low frequency to high frequency, the fault characteristic frequency of different parts can be obtained again by envelope demodulation and spectrum analysis, orient inner ring, outer ring or roller part fault.Because traditional resonance and demodulation method adopts the envelope-demodulation method of Hilbert conversion, there is bandpass filter centre frequency and bandwidth parameter is difficult to problem identificatioin.For this reason, some scholars propose and Xi Er baud-Huang (Hilbert – Huang Transform, HHT), improving one's methods of composing that kurtosis (Spectral Kurtosis, SR) etc. combines.Patent (CN201210346147) discloses bearing vibration signal characteristic abstraction and analytical approach under a kind of initial failure state, can give prominence to the modulation signature in vibration signal under reflection initial failure state; Patent CN201410140890 discloses a kind of Rolling Bearing Fault Character extracting method and system, is optimized chooses resonance and demodulation frequency band, improves the quality that fault signature extracts to a certain extent, etc.

Although achieve some achievements in research in bearing fault characteristics frequency abstraction, directly apply to variable motion detrain bearing state detect, poor effect.For this reason, the present invention organically blends the multiple technologies such as Short Time Fourier Transform, order ratio analysis, independent component analysis, neural network, provide a kind of new train Rolling Bearing Fault Character frequency extraction method, the method can be used for variable motion lower bearing on-line checkingi, has good application prospect.

Summary of the invention

The object of this invention is to provide a kind of variable speed train Rolling Bearing Fault Character frequency extraction method, by carrying out Time-Frequency Analysis to bearing vibration signal, in conjunction with local peaking's adaptable search, extracting not the instantaneous frequency values of rotating speed in the same time; Adopt neural network to carry out matching to instantaneous frequency, obtain the speed curves with reference to main shaft, on this basis angularly resampling and order ratio analysis are carried out to original signal; Adopt point of fixity independent component analysis and spectrum peak search method, than signal, independent source separating treatment is carried out to mixing rank, feature is compared on the rank obtaining bearing fault parts, the method is applicable to bearing variable motion and fault characteristic frequency extracts under not having hardware velocity gauge situation, can meet the inline diagnosis demand of actual trains rolling bearing.

To achieve these goals, the present invention proposes following technical scheme:

Step one, install acceleration transducer on housing washer surface, number of sensors is greater than and equals to analyze the quantity of independent source here, gathers vibration signal during bearing movable;

Step 2, utilize short time discrete Fourier transform method to carry out time frequency analysis to variable speed bearing vibration signal, the Frequency point that in Local Search time-frequency energy distribution, energy is maximum, extracts not the instantaneous frequency that the speed of mainshaft is corresponding in the same time;

Step 3, design BP neural network model, carry out matching to discrete instantaneous frequency point, according to the theory relation between rotating speed and instantaneous frequency, calculate the real-time rotate speed of train rolling bearing main shaft;

The real-time rotate speed information of step 4, foundation bearing, is carried out equal angular resampling to vibration signal, calculates equal angular sampling markers, by interpolating function, original signal is transformed to the stationary signal of angle domain;

Step 5, carry out equalization, Whitening to mixed signal, introduce point of fixity independent component analysis and spectrum peak search technology, be separated than signal mixing rank, feature is compared on the rank obtaining bearing fault parts.

The advantage of the method is:

(1) adopt software approach completely, without the need to hardware units such as velocity gauges, variable speed motion can be estimated and to detrain the real-time rotate speed of rolling bearing, there is stronger use dirigibility;

(2) by carrying out angularly resampling to original vibration signal, unstable state acoustical signal being transformed to the steady-state signal in angularly territory, the noise irrelevant with rotating speed can be removed, advantageously in the identification of consequent malfunction feature extraction and classifying;

(3) by introducing point of fixity independent component analysis and spectrum peak search technology, effectively can isolate mixing rank than component, extracting the characteristic frequency of be concerned about parts, be applicable to the on-line checkingi of train rolling bearing.

Accompanying drawing explanation

Fig. 1 is igneous rock cracks simulated experiment platform and sensor deployment figure in the invention process, in figure:

Bearing on the left of bearing 5-on the right side of 1-left side wheel 2-right side wheels 3-rotating shaft 4-

6-motor (containing frequency converter) 7-hydraulic pump 8-friction pulley 9-supporting seat 10 (11,12,13)-acceleration transducer

Fig. 2 is the structural drawing of train rolling bearing, in figure:

1-outer ring 2-inner ring 3-roller 4-retainer

Fig. 3 is method flow diagram of the present invention;

Embodiment

Below in conjunction with accompanying drawing, the present invention is described in further detail:

Fig. 1 gives the structure of igneous rock cracks simulator stand, mainly comprise wheel (1,2), rotating shaft (3), bearing (4,5), motor (6), hydraulic pump (7), friction pulley (8), supporting seat (9); Lance 2052 acceleration transducer (10-13) is disposed outside each rolling bearing; Bearing drives friction pulley to rotate by motor, and outer ring geo-stationary is motionless; Fig. 2 gives the outer ring stripping of experiment employing and the double-row conical bearing bearing (197726 model) of inner ring stripping fault, be made up of a few part such as outer ring (1), inner ring (2), roller (3), retainer (4), wherein, roller number N=20, roller diameter d=23.7mm, bearing pitch diameter D=180mm, contact angle a=10 degree, by regulating Frequency Converter Control bearing from 50rpm accelerated motion to 240rpm, systematic sampling rate 25.6kHz;

The present invention is a kind of variable speed train Rolling Bearing Fault Character frequency extraction method, and algorithm flow as shown in Figure 3, specifically comprises the following steps:

1, outside every side roll bearing, L acceleration transducer being installed, here meeting in the quantity situation that observation sensor quantity is more than or equal to independent failure source, choosing L=2, pickup L road vibration signal x (n);

2, Short-time Fourier analysis is done to single channel vibration signal, peak point in Local Search time-frequency energy distribution, extract the instantaneous frequency that speed of mainshaft is not corresponding in the same time;

1. discrete Short Time Fourier Analysis is done to bearing vibration signal x (n), when obtaining-frequency division cloth function F (k, f):

$F(k,f)=\underset{n=-\∞}{\overset{n=\∞}{\mathrm{\Σ}}}x\left(n\right)w^{*}(n-k)e^{-j2\mathrm{\πnf}}---\left(1\right)$

Wherein, w ^*() represents the conjugate complex number of window function, and adopt Hanning window here, its window width chooses 512 points;

2. time by asking-mould square of frequently distribution function, can obtain vibration signal time-frequency energy distribution P (k, f)=| F (k, f) | ²;

3. instantaneous frequency corresponding to rotating speed cannot accurately be followed the tracks of in order to avoid noise is excessive, adopt local auto-adaptive search peak mode, obtain the instantaneous frequency of main shaft, concrete grammar is: intercept relatively stable time window in the actual changed speed operation process of bearing and analyze, estimate theoretical turn of frequency value in each time window, suppose t ₀it is f that moment turns estimated value frequently ₀, evaluated error Δ f, then with interval [f ₀-Δ f, f ₀+ Δ f] maximal value turn frequently pursuit gain, simultaneously as the updated value at the Frequency Estimation center of subsequent time as this moment; It should be noted that, algorithm does not limit and turns fundamental frequency feature frequently, if turn frequently, multicomponent energy is less, then the higher hamonic wave frequently that turns of searching for energy higher is followed the tracks of, and is then scaled and turns fundamental frequency frequently; Method search one by one like this, until searched for free, obtain the real-time instantaneous frequency of main shaft;

3, design the instantaneous frequency model of fit of BP neural network, according to the theory relation between rotating speed and instantaneous frequency, calculate the real-time rotate speed of train rolling bearing main shaft;

1. set up the data fitting model of neural network, wherein, the neuron number p of input layer gets 1, using discrete time point as input vector; Output layer neuron number q is 1, the object vector using discrete instantaneous frequency as training function; The neuron number of hidden layer, can obtain according to experimental formula below:

$n=\sqrt{p+q}+\mathrm{\β}---\left(2\right)$

Wherein, β is nondimensional corrected parameter, β=1 ~ 10;

The span of hidden layer neuron number is [1 12], designs variable BP neural network model, considers time and the error precision of training, optimizes and determines that hidden layer neuron number is 5;

2. choose mean square error function as the deviation between neural computing result and desired output, computing method are as follows:

$J=\frac{1}{2}\underset{k=1}{\overset{q}{\mathrm{\Σ}}}{(y_k^o-{\stackrel{\‾}{y}}_k)}^2---\left(3\right)$

Wherein, represent that k moment neural computing obtains output frequency value, represent the desirable output frequency value in k moment, if error J is less than given minimum positive number, then neural metwork training terminates, thus m-instantaneous frequency profile when obtaining;

3. according to relation: r=60*f between bearing spindle rotating speed r and instantaneous frequency f, and T/F curve, calculate the real-time rotate speed of bearing;

4, according to the real-time rotate speed information of bearing, carry out angularly resampling to original vibration signal, non-stationary signal is transformed to the stationary signal of angle domain, concrete grammar is:

1. calculate the sampling markers in angularly territory, can regard approximate uniform variable motion as each window inner bearing motion analysis time, its angle θ turned over can be expressed as:

θ(t)＝b ₀+b ₁t+b ₂t ²(4)

B in formula ₀, b ₁, b ₂for undetermined coefficient, method for solving is: the value θ (t supposing the continuous sampling moment ₁)=0, θ (t ₂)=Δ θ, θ (t ₃)=2 Δ θ, then simultaneous can try to achieve coefficient b ₀, b ₁, b ₂, substitute into formula (4), can try to achieve the moment t in this time window corresponding to any rotation angle, then the fundamental formular of resampling moment calculating is:

$t_n=\frac{1}{2b_2}[\sqrt{b_1^2+4b_2({\mathrm{\θ}}_k-b_0)}-b_1]---\left(5\right)$

θ _kfor elapsed time t _nthe angle turned over; Method is until searched for all time windows like this, obtains bearing rotational angle and time relationship;

2. according to calculating the angularly resampling moment t tried to achieve _n, by interpolating function, angularly resampling is carried out to raw data, thus obtains angle domain sampled signal:

$x\left(t_n\right)=\underset{k}{\mathrm{\Σ}}x(k\·\mathrm{\Δt})*h(k\·\mathrm{\Δt}-t_n)---\left(6\right)$

Wherein, Δ t is time-domain sampling interval, and h () is interpolating function, gets

5, introduce point of fixity Independent Component Analysis Technology, be separated than signal mixing rank, concrete grammar is:

1. carry out centralization process to mixing rank than signal x, obtain new data x=x-E (x) that average is 0, wherein E (x) represents the average of x;

2. to going the signal data after average to carry out whiten process, detailed process is: carry out svd x=U Σ V ' to x, wherein, Σ is pre-diagonal matrix, U and V is the square formation of orthonomalization, and V ' is the transposed matrix of V, does following decomposition to the covariance matrix C of x:

C＝xx′＝[UΣV′][VΣ′U′]＝UΣ ²U′＝UΛU′ (7)

Wherein, representation feature value matrix, λ _i(i=1,2 ..., m) be the characteristic root of C, and λ _i> 0, U=[u ₁, u ₂..., u _m] be feature matrix, be the left singular matrix in svd;

With transformation operator s=Λ ^-1/2u ' multiplier, according to x, obtains new data Z:

Z＝sx＝Λ ^-1/2U′x (8)

Because ZZ '=[Λ ^-1/2u ' x] [x ' U Λ ^-1/2]=Λ ^-1/2u ' [U Λ U '] U Λ ^-1/2=Λ ^-1/2Λ Λ ^-1/2=I, so, through Whitening, the new data Z that average is 0, uncorrelated, unit variance is 1 can be obtained;

3. adopt negentropy maximization as being separated the evaluation index of mixing rank than signal, because according to central limit theorem, if stochastic variable Z is by many mutually independent random variables S _i(i=1,2,3 ... n) sum composition, as long as S _ihave limited average and variance, no matter then which kind of distribution it is, stochastic variable Z is S comparatively _icloser to Gaussian distribution, i.e. S _ithere is stronger non-Gaussian system, therefore, the non-Gaussian system tolerance by separating resulting represents the mutual independence between separating resulting, when non-Gaussian system reaches maximum, then show that the detachment process of isolated component completes, adopt the correction form negentropy of entropy to measure non-Gaussian system here:

N _g(Z)＝H(Z _Gauss)-H(Z) (9)

Wherein, Z _gaussbe have mutually homoscedastic Gaussian random variable with Z, the differential entropy that H () is stochastic variable, have in mutually homoscedastic stochastic variable, the stochastic variable of Gaussian distribution has maximum differential entropy, and non-Gaussian system is stronger, and its differential entropy is less, N _g(Z) value is larger, considers that the probability density function of Z is unknown, adopts following approximate formula:

N _g(Z)＝E(g(z))-E(g(Z _Gauss)) ²(10)

Wherein, E () is mean operation, and g () is nonlinear function, and choosing of g is not only, selects here

Learning rules are that searching direction makes rank compared estimate value Y=W'X have maximum non-Gaussian system, and W refers to separation matrix here, adopt Newton iteration method to solve, for inverting of simplification matrix, and because data are by albefaction, the norm equaling to retrain W is 1, adopts the following iterative formula simplified:

$W(k+1)=E[Z\·g(W^{\′}\left(k\right)\·Z)]-E\left[\stackrel{\‾}{g}(W^{\′}\left(k\right)\·Z)\right]\·W\left(k\right)---\left(11\right)$

the derivative of representative function g (), E [] expression is averaging, and k represents the number of times of iterative computation;

4. select initially-separate matrix W (k) arbitrarily, require the 2-norm of W (k) || W (k) || ₂=1, according to formula (11), by the calculating that iterates, until convergence, obtain the signal of initial gross separation;

6, adopt spectrum peak search method, extract the rank of bearing fault parts than component characterization, concrete grammar is:

Utilize the further Analyze & separate signal of Short Time Fourier Transform, the maximum component of energy in separation signal is found out by spectrum peak search, record its instantaneous frequency, then according to frequency information, signal is reconstructed, remove the component of signal that energy is low, extract desired signal component, experimental calculation obtains bearing inner race, the rank of outer ring trouble unit are respectively 11.2 and 8.6 than characteristic component; In theory, the rank of inner ring and outer ring are calculated as follows than feature:

$\left\{\begin{array}{c}{I}_{\mathrm{outter}}=\frac{f_{\mathrm{Outer}}}{f_r}\frac{\frac{N}{2}{f}_r(1-\frac{d}{D}\cos \mathrm{\α})}{f_r}=\frac{N}{2}(1-\frac{d}{D}\cos \mathrm{\α})=\frac{20}{2}(1-\frac{23.7}{180}\cos 10)=8.7\\ f_{\mathrm{Inner}}=\frac{f_{\mathrm{Inner}}}{f_r}=\frac{\frac{N}{2}{f}_r(1+\frac{d}{D}\cos \mathrm{\α})}{f_r}=\frac{N}{2}(1-\frac{d}{D}\cos \mathrm{\α})=\frac{20}{2}(1+\frac{23.7}{280}\cos 10)=11.296\end{array}\right.$

Wherein, I _outer, I _innerrepresent fault signature rank, outer ring ratio and inner ring fault signature rank ratio respectively, f _outer, f _innerrepresent outer ring fault characteristic frequency and inner ring fault characteristic frequency respectively, f _rrepresent turning frequently of bearing, equal rotating speed divided by 60; From result, substantially identical between actual computation value and theoretical value.

Finally it should be noted that above embodiment only in order to technical scheme of the present invention and unrestricted to be described, can refer to technical scheme of the present invention and modify or equivalently to replace.

Claims (3)

1. variable speed train Rolling Bearing Fault Character frequency extraction method, is characterized in that comprising the following steps:

Step one, utilize short time discrete Fourier transform method to carry out time frequency analysis to variable speed bearing vibration signal, the Frequency point that in Local Search time-frequency energy distribution, energy is maximum, extracts not the instantaneous frequency that the speed of mainshaft is corresponding in the same time;

(1) vibration signal x (n) of given bearing, do discrete Short Time Fourier Transform to it, circular is:

Wherein, w ^*() represents the conjugate complex number of window function, and F (k, f) represents the frequency spectrum of signal at k frequency f analysis time place;

(2) time by asking-mould square of frequently distribution function, can obtain vibration signal time-frequency energy distribution P (k, f)=| F (k, f) | ²;

(3) adopt local peaking's way of search, the instantaneous frequency f of acquisition main shaft, concrete mode is: in the actual changed speed operation process of intercepting bearing, time window is analyzed relatively stably, estimates theoretical turn of frequency value in each time window, supposes t ₀it is f that moment turns estimated value frequently ₀, evaluated error Δ f, then with interval [f ₀-Δ f, f ₀+ Δ f] maximal value turn frequently pursuit gain, simultaneously as the updated value at the Frequency Estimation center of subsequent time as this moment; It should be noted that, algorithm does not limit and turns fundamental frequency feature frequently, if turn frequently, multicomponent energy is less, then the higher hamonic wave frequently that turns of searching for energy higher is followed the tracks of, and is then scaled and turns fundamental frequency frequently; Method search one by one like this, until searched for free, obtain the real-time instantaneous frequency of main shaft;

The instantaneous frequency model of fit of step 2, design BP (Back Propagation) neural network, calculates the real-time rotate speed of train rolling bearing main shaft;

(1) set up three-layer neural network model, wherein, the neuron number of input layer is 1, using discrete time point as input vector; Output layer neuron number is 1, the object vector using instantaneous frequency corresponding to discrete time point as training function; The neuron number of hidden layer, can obtain according to experimental formula below:

Wherein, β is nondimensional corrected parameter, β=1 ~ 10, and p is input neuron number, and q is output neuron number; According to the span of hidden layer neuron number, design variable BP neural network model, finally determined the quantity of hidden layer neuron by application condition;

(2) using mean square error function as the deviation between neural computing result and desired output, computing method are as follows:

Wherein, represent that k moment neural computing obtains output frequency value and desirable output frequency value respectively; If error J is less than given minimum positive number, then neural metwork training terminates, thus m-instantaneous frequency profile when obtaining;

(3) according to relation: r=60*f between bearing spindle rotating speed r and instantaneous frequency f, and T/F curve, the real-time rotate speed of bearing can be calculated;

The real-time rotate speed information of step 3, foundation bearing, carries out equal angular resampling to vibration signal, original signal is transformed to the stationary signal of angle domain;

Step 4, employing point of fixity independent component analysis method, be separated than signal mixing rank;

Step 5, employing spectrum peak search process separation signal, obtain trouble unit rank than component characterization, concrete grammar is:

Utilize the further Analyze & separate signal of Short Time Fourier Transform, the maximum component of energy in separation signal is found out by spectrum peak search, record its instantaneous frequency, then according to frequency information, signal is reconstructed, remove other component of signals that energy is low, extract desired signal component, component characterization is compared on the rank obtaining trouble unit.

2. variable speed train Rolling Bearing Fault Character frequency extraction method according to claim 1, it is characterized in that: according to the real-time rotate speed information of bearing, carry out equal angular resampling to vibration signal, original signal is transformed to the stationary signal of angle domain, concrete grammar is:

(1) calculate the sampling markers in angularly territory, can regard approximate uniform variable motion as each window inner bearing motion analysis time, its angle θ turned over can be expressed as:

θ(t)＝b ₀+b ₁t+b ₂t ²(4)

B in formula ₀, b ₁, b ₂for undetermined coefficient, by method for solving in this time window be: the value θ (t supposing the continuous sampling moment ₁)=0, θ (t ₂)=Δ θ, θ (t ₃)=2 Δ θ, then simultaneous can try to achieve coefficient b ₀, b ₁, b ₂, substitute into formula (4), can try to achieve the moment t in this time window corresponding to any rotation angle, then the fundamental formular of resampling moment calculating is

θ _kfor elapsed time t _nthe angle turned over; Method is until searched for all time windows like this, obtains bearing and forwards angle and time relationship;

(2) according to calculating the angularly resampling moment t tried to achieve _n, by interpolating function, angularly resampling is carried out to raw data, thus obtains angle domain sampled signal:

Wherein, Δ t is time-domain sampling interval, and h () is interpolating function, gets .

3. variable speed train Rolling Bearing Fault Character frequency extraction method according to claim 1, is characterized in that: adopt point of fixity independent component analysis method, and be separated than signal mixing rank, concrete grammar is:

(1) carry out centralization process to mixing rank than signal x, obtain new data x=x-E (x) that average is 0, wherein E (x) represents the average of x;

(2) to going the signal data after average to carry out whiten process, detailed process is: carry out svd x=U Σ V ' to x, wherein, Σ is pre-diagonal matrix, U and V is the square formation of orthonomalization, and V ' is the transposed matrix of V, does following decomposition to the covariance matrix C of x:

C＝xx′＝[UΣV′][VΣ'U′]＝UΣ ²U′＝UΛU′ (7)

Wherein, representation feature value matrix, λ _i(i=1,2 ..., m) be the characteristic root of C, and λ _i> 0, U=[u ₁, u ₂..., u _m] be feature matrix, be the left singular matrix in svd;

With transformation operator s=Λ ^-1/2u ' multiplier, according to x, obtains the new data Z that average is 0, uncorrelated, unit variance is 1:

Z＝sx＝Λ ^-1/2U′x (8)

(3) mutual independence between separating resulting is represented by the non-Gaussian system tolerance of separating resulting, here adopt the correction form negentropy of entropy to measure non-Gaussian system, also namely using negentropy maximization as being separated the evaluation index of mixing rank than signal, the concrete grammar of negentropy is:

N _g(Z)＝H(Z _Gauss)-H(Z) (9)

Wherein, Z _gaussbe have mutually homoscedastic Gaussian random variable with Z, the differential entropy that H () is stochastic variable, have in mutually homoscedastic stochastic variable, the stochastic variable of Gaussian distribution has maximum differential entropy, and non-Gaussian system is stronger, and its differential entropy is less, N _g(Z) value is larger, considers that the probability density function of Z is unknown, adopts following approximate formula:

N _g(Z)＝E(g(Z))-E(g(Z _Gauss)) ²(10)

Wherein, E () is mean operation, and g () is nonlinear function, and choosing of g is not only, selects here

(4) finding a direction makes rank compared estimate value Y=W'X have maximum non-Gaussian system, and W refers to separation matrix here, adopts Newton iteration method to solve, for inverting of simplification matrix, and because data are by albefaction, the norm equaling to retrain W is 1, adopts the following iterative formula simplified:

the derivative of representative function g (), E [] expression is averaging;

(5) select initially-separate matrix W (k) arbitrarily, require the 2-norm of W (k) || W (k) || ₂=1, according to formula (11), by the calculating that iterates, until convergence, obtain the signal of initial gross separation.