CN105182172A - Vibration signal pattern spectrum-based method for diagnosing winding conditions under sudden short circuit of transformer - Google Patents

Abstract

The invention discloses a vibration signal pattern spectrum-based method for diagnosing winding conditions under sudden short circuit of a transformer. The method includes the following steps that: (1) real-time short circuit monitoring is performed on the transformer, and a section of continuous vibration signals of the box wall of the transformer are acquired when sudden short circuit occurs on the transformer; the steady state vibration signal pattern spectrum entropy of vibration signals which are acquired when short circuit occurs on the transformer every time is calculated through methods mentioned from step (2) to step (6); and (7) the winding conditions are judged according to the pattern spectrum entropy of the steady state vibration signals xs(t): when the change Delta H of the pattern spectrum entropy satisfies an expression described in the descriptions in the invention, it is judged that the winding of the transformer becomes loose or is deformed when Y-th short circuit occurs on the transformer. The method of the invention is accurate and efficient.

Description

Based on winding state diagnostic method under the transformer sudden short circuit of vibration signal morphology spectrum

Technical field

The present invention relates to a kind of signal monitoring method, particularly relate to the diagnostic method of winding state during a kind of transformer generation sudden short circuit.

Background technology

Transformer is one of most important equipment in electric system, and its operation stability on power system security impact greatly.Along with the increase day by day of China's net capacity, capacity of short circuit also constantly increases thereupon, the deformation of transformer winding fault caused because of external short circuit remains high always, is fault comparatively common in transformer operational process, causes very large threat to the safe operation of electric system.

After transformer suffers sudden short circuit, first its winding may occur to loosen or slight deformation, and deformation of transformer winding has cumulative effect, if can not Timeliness coverage and repair Transformer Winding loosen or metaboly, then when transformer loosens or be out of shape to be accumulated to a certain degree the anti-short circuit capability of transformer can be made to decline to a great extent, now, occur even if suffer less short-circuit impact electric current also likely can cause large accident.

Winding deformation can cause Transformer Winding physical strength to decline on the one hand, also coil inside minor insulation distance can be made to change on the other hand, there is minor insulation thin spot, when running into superpotential effect, likely occur between cake or turn-to-turn short circuit causes transformer insulated breakdown accident, or cause shelf depreciation because local field strength increases, insulation harm position can expand gradually, finally cause transformer generation dielectric breakdown accident and cause the expansion of the further state of affairs.

Therefore, in operational process in the routine after transformer experienced by external short circuit accident or after running a period of time is overhauled, how effectively to detect whether Transformer Winding exists to loosen and distortion, and then take effective O&M measure just to seem very important, this is the important means ensureing that transformer safety is run.

At present, it is one of conventional test project of current transformer that winding deformation detects, and the detection method in practical application mainly contains short circuit impedance method and Frequency Response Analysis method.Wherein, by carrying out detection to the short-circuit reactance of transformer, short circuit impedance method judges that whether tested Transformer Winding is qualified.Generally, after operating transformer receives the impact of short-circuit current, or the short-circuit impedance value recorded and original record will be compared when regular routine inspection and judge whether winding there occurs distortion, if short-circuit impedance value changes greatly, such as, be set as in GB that change is more than 3%, then can confirm that winding has remarkable distortion.Up to the present, short circuit impedance method establishes standard in long-term production practices, and criterion is comparatively clear and definite, all clearly gives the criterion of winding deformation degree in international electrical engineering standard IEC60076-5 and GB1095-85.But the sensitivity of this method is very low in a lot of situation, only clearer and more definite reflection can be obtained when coil bulk deformation situation is comparatively serious.

Frequency Response Analysis method is analyzed by the situation of change of the network transfer function curve analyzing Transformer Winding and is judged whether the network electrical quantity of its inside changes, thus infers whether corresponding physical construction there occurs distortion.By a stable sine sweep voltage signal V _ibe applied to one end of tested Transformer Winding, then record this port V simultaneously _iwith the voltage V on other output port _o, thus obtaining one group of Frequency Response curve of this tested winding, its expression formula is

H(jω)＝V _o/V _i

Comparatively short circuit impedance method is high for the measurement sensitivity of Frequency Response Analysis method, but due to the complicacy of its frequency response waveform, needs more experience, the quantitative criteria that more difficult formation is clear and definite, therefore do not form discrimination standard so far to the differentiation of winding situation.

The starting point of short circuit impedance method and these two kinds of methods of Frequency Response Analysis method is all change according to the electrical quantity of Transformer Winding to carry out measurement to it and differentiate, it is comparatively suitable to there is obvious deformation to Transformer Winding in this, but to winding generation slight deformation, especially comparatively clear and definite judgement can not be provided to the relatively loosening and torsional deformation state that Transformer Winding exists.Consider that Transformer Winding loosens or torsional deformation has a great impact its anti-short circuit capability, therefore need to seek more sensitive winding condition diagnosing method.

The ultimate principle of the vibration analysis method that development in recent years is got up is that Transformer Winding is regarded as a physical construction body, then, when winding construction or any change of stressed generation, can be reflected from its mechanical vibration performance change.The vibration of winding is delivered to transformer-cabinet by inside transformer structural connection, so the basket vibration characteristic of the vibration signal that obtains of transformer-cabinet Surface testing and transformer has close relationship.Under considering sudden short circuit, transformer-cabinet surface vibration signals is mainly basket vibration, therefore, and the approach analysis of vibration signal of transformer tank surface can diagnosed as transformer winding fault.Compared with aforementioned electrical measurements, the great advantage of vibration analysis method is by being adsorbed on vibration transducer on transformer box wall to obtain the vibration signal of transformer, the situation of change of winding state is judged by the change analyzing its vibration characteristics, as long as the mechanical property of winding (as malformation, pretightning force loosen) changes, can be reflected from its vibration characteristics change, thus be substantially increased the sensitivity of detection.In addition, the vibration detection be placed in by vibration transducer on tank wall is not directly connected with whole strong power system, can develop into a kind of more accurate, convenient, safe on-line monitoring method.

When transformer suffers sudden short circuit, have huge electrodynamic action in Transformer Winding, make transformer vibrate aggravation, its vibrational waveform present stage that vibration amplitude increases in time domain, vibration amplitude decling phase after vibration amplitude metastable stage and short-circuit current disappear.Therefore, the change indicator of transformer winding state can be sought from the time domain form of vibrational waveform during sudden short circuit, and use it for the efficient diagnosis of transformer winding state.

Summary of the invention

Technical matters to be solved by this invention is: provide a kind of based on winding state diagnostic method under the transformer sudden short circuit of vibration signal morphology spectrum, can be differentiated by transformer box wall vibration signal during on-line monitoring transformer sudden short circuit to winding duty.

Solve the problems of the technologies described above, the technical solution adopted in the present invention is as follows:

Based on a winding state diagnostic method under the transformer sudden short circuit of vibration signal morphology spectrum, comprise the following steps:

(1) real-time short circuit monitoring is carried out to transformer, and, whenever transformer generation sudden short circuit being detected, gather one section of continuous print vibration signal of transformer box wall, and when calculating transformer generation sudden short circuit each time by following steps (2) to the method for step (6) gather vibration signal steady-state vibration signal aspect compose entropy; Wherein, during transformer generation sudden short circuit to transformer box wall gather vibration signal start time should as far as possible close to the initial time of transformer generation sudden short circuit, transformer box wall gathered to end time of vibration signal and judges by existing vibration detection device, with the complete vibration signal collecting transformer box wall and produce due to transformer generation sudden short circuit as far as possible;

The vibration signal collected from transformer box wall during the above-mentioned generation sudden short circuit of transformer is each time expressed as x ₀(t);

(2) according to following step according to vibration signal x ₀(t) structure white noise of the narrowband signal g (t):

2a. is to vibration signal x ₀t () carries out Fourier transform, obtain vibration signal x ₀the spectrum distribution of (t); In this step, Fourier transform is mathematical method conventional in this area, and therefore inventor is no longer described in detail at this;

2b. is according to vibration signal x ₀t the spectrum distribution of () generates white noise signal s ₀(t), wherein white noise signal s ₀t the expression formula of the amplitude A of () is

In formula, K is coefficient, and COEFFICIENT K value is vibration signal x ₀1/10 of (t) amplitude mean value; N _ffor vibration signal x ₀highest frequency component f in (t) spectrum component _hwith the ratio of 50Hz; f _i(i=1,2 ..., N _f) be vibration signal x ₀each frequency component of (t); a _i(i=1,2 ..., N _f) be vibration signal x ₀t amplitude that each frequency component of () is corresponding;

2c. uses Butterworth bandpass filter to white noise signal s ₀t () carries out filtering, obtain white noise of the narrowband signal g (t); The transport function expression formula of described Butterworth bandpass filter is

In formula, ω _cl=K _l2 π f _lfor low-frequency cut-off frequency, f _lfor vibration signal x ₀lowest frequency components in (t) spectrum component, K _lfor low frequency bandwidth factor, usually get 1.5; ω _ch=K _h2 π f _hfor high-frequency cut-off frequency, f _hfor vibration signal x ₀highest frequency component in (t) spectrum component, K _hfor high frequency bandwidth coefficient, usually get 0.8; M is filter order; ω=2 π f is angular frequency, f=50;

(3) by vibration signal x ₀t () is added with white noise of the narrowband signal g (t), obtain superposed signal x (t), according to following step, superposed signal x (t) is decomposed into several intrinsic mode function components (IntrinsicModeFunction, referred to as IMF component):

3a., to superposed signal x (t) differentiate, obtains time series y (t);

The product that 3b. sequences y computing time (t) is adjacent 2

py _i(t)＝y _i(t)×y _i-1(t)

In formula, i=2,3 ..., N-1, wherein, N is that signal is counted;

3c. is according to product py _it () and time series y's (t) is positive and negative, look for all Local modulus maximas eb (t) of superposed signal x (t) and all local minizing points es (t) successively:

Work as py _iduring (t) < 0, if py _i(t) < 0 and y _i-1(t) < 0, then x _i-1t () is local minizing point; If py _i(t) < 0 and y _i-1(t) > 0, then x _i-1t () is Local modulus maxima;

Work as py _iduring (t) > 0, x _i-1t () is non-extreme point;

Work as py _iduring (t)=0, if y _i-1t ()=0, calculates 2 y _i(t) and y _i-2t the product of (), makes py _i(t) '=y _i(t) × y _i-2t (), if py _i(t) ' < 0 and y _i-2(t) < 0, then x _i-1t () is local minizing point; If py _i(t) ' < 0 and y _i-2(t) > 0, then x _i-1t () is Local modulus maxima; If y _i-2(t)=0, then x _i-1t () is non-extreme point;

All Local modulus maximas eb (t) and all local minizing points es (t) couple together with cubic spline functions s (t) and obtain coenvelope line e respectively by 3d. _max(t) and lower envelope line e _mint (), described cubic spline functions s (t) is each minizone [t at superposed signal x (t) _i, t _i+1] (i=1,2 ..., N-1) on be no more than the polynomial expression of three times, its expression formula is

In formula, m _iand m _i+1for cubic spline functions s (t) is at interval [t _i, t _i+1] second derivative values at two-end-point place; In this step, the algorithm of envelope is mathematical method conventional in this area, and therefore inventor is no longer described in detail at this;

3e. is according to the coenvelope line e tried to achieve _max(t) and lower envelope line e _mint () calculates average m (t)=(e of upper and lower envelope _max(t)+e _min(t))/2, superposed signal x (t) is deducted m (t), obtains a new time series y ₁(t);

3f. judges time series y ₁t whether () meet following two conditions simultaneously:

A., in whole signal length, the number of extreme point and zero crossing only must differ one equal or at the most;

B. at any time, the coenvelope line defined by maximum point and the mean value of lower envelope line defined by minimum point are zero;

If meet above-mentioned two conditions, then y simultaneously ₁t () is intrinsic mode function component; If above-mentioned two conditions can not be met, then by y simultaneously ₁t (), as an original component, repeats abovementioned steps 3a ~ 3e, until time series y ₁t () meets above-mentioned two conditions simultaneously, will meet the y of above-mentioned two conditions simultaneously ₁t () is designated as c _i(t), then c _ia t intrinsic mode function component that () is superposed signal x (t), and i=1,2 ..., N _h, N _hfor the intrinsic mode function component total quantity of superposed signal x (t);

3g. is by c _it () separates from superposed signal x (t), obtain difference signal r _i(t) be

r _i(t)＝x(t)-c _i(t)

3h. is by difference signal r _it () carrys out repetition above-mentioned steps 3a ~ 3g as pending signal substituting superposed signal x (t), until meet Stopping criteria, obtain whole N _hindividual intrinsic mode function component, described Stopping criteria is: obtain new time series y _it () is narrow band signal, this narrow band signal refers to that the bandwidth deltaf f of signal is much smaller than centre frequency f _csignal, its concept is that one of ordinary skilled in the art is known;

Through above-mentioned steps, vibration signal x ₀t () has been broken down into several intrinsic mode function components and residual signal sum, its expression formula is

x ₀(t)＝Σc _i(t)+r _i(t)；

(4) vibration signal x is calculated ₀t () and step (3) decompose the related coefficient between each intrinsic mode function component of obtaining, and the intrinsic mode function component chosen corresponding to the maximum related coefficient of numerical value is eigen vibration natural mode of vibration component, is designated as IMF _m;

Wherein, described vibration signal x ₀(t) related coefficient cor of decomposing kth the intrinsic mode function component that obtain middle with step (3) _kcomputing formula be

In formula, x _0it () is vibration signal x ₀t () is at t+i/f ₀the amplitude in moment; for vibration signal x ₀the mean value of (t) amplitude; c _kit () is for a kth intrinsic mode function component is at t+i/f ₀the amplitude in moment; for a kth intrinsic mode function component c _kthe mean value of (t) amplitude; f ₀for vibration signal x ₀the sample frequency of (t);

(5) according to following step according to natural mode of vibration component IMF _mpeak change determination vibration signal x ₀the initial time t of (t) steady-state vibration process _bwith finish time t _e, and remember vibration signal x ₀t () is from initial time t _bto finish time t _etime period in part be steady-state vibration signal x _s(t):

5a. calculates eigen vibration natural mode of vibration component IMF _mspectrum distribution, choose frequency spectrum

In component, the vibration frequency of amplitude maximum is characteristic frequency, is designated as f _m;

5b. is from vibration signal x ₀t the initial time of () starts to calculate vibration signal x successively in chronological order ₀t (), at grid number Ng (t) in each moment, wherein, grid number Ng (t) refers to vibration signal x ₀t () is at time period [t, t+0.5T _m] mean value of interior variable quantity absolute value sum, its computing formula is

T _m＝1/f _m

In formula, N ₀for vibration signal x ₀t () is at time period [t, t+0.5T _m] in data length; x _0i(t) and x _{0 (i+1)}t () is respectively vibration signal at t+i/f ₀moment and t+ (i+1)/f ₀moment amplitude; T _mfor characteristic frequency f _minverse;

5c. is in chronological order successively to vibration signal x ₀t () screens at grid number Ng (t) in each moment, the principle of screening is: if grid number Ng (t) is sometime greater than default threshold value, then retain the grid number in this moment; If Ng (t) is sometime less than default threshold value, then the grid number in this moment is set to 0; Note vibration signal x ₀(t) screen after grid number for screening after grid number Ng'(t); Wherein, the preferred value of threshold value preset is 20% of grid number Ng (t) maximal value and rounds value;

5d. is from grid number Ng'(t after screening) in moment corresponding to first non-zero grid number carry out out Top-hat computing successively, first peak value detected is vibration signal x ₀the initial time t of (t) steady-state vibration process _b; From grid number Ng'(t after screening) in moment corresponding to last non-zero grid number oppositely carry out out Top-hat computing successively, first peak value detected is vibration signal x ₀the finish time t of (t) steady-state vibration process _e; The described expression formula opening Top-hat computing is

In formula, f is grid number Ng'(t); G is structural element; F Θ g carries out erosion operation for using structural element g to signal f; for using structural element g, dilation operation is carried out to signal f Θ g; The fundamental operation of mathematical morphology when described dilation operation and erosion operation, be mathematical method conventional in this area, therefore inventor is no longer described in detail at this;

(6) steady-state vibration signal x is calculated _st the morphology spectrum of () and morphology spectrum entropy, its computing formula is

q(λ)＝PS _f(λ,g)/S(f·λ _ming)

In formula, H (f/g) is morphology spectrum entropy; PS _f(λ, g) is morphology spectrum; F is grid number Ng'(t); G is structural element; λ is yardstick; λ _maxfor yardstick maximal value; λ _minfor yardstick minimum value.

(7) according to steady-state vibration signal x _st the morphology spectrum entropy of () differentiates winding state: when the changes delta H of morphology spectrum entropy meets time, judge that the winding of transformer occurs when the short circuit of Y secondary burst occurs transformer to loosen or distortion, now need to carry out overhaul plan in time, avoid the formation of significant trouble; Wherein, H ₁the steady-state vibration signal aspect gathering vibration signal for institute when sudden short circuit occurs for transformer the 1st time composes entropy, H _ythe steady-state vibration signal aspect gathering vibration signal for institute when sudden short circuit occurs for transformer the Y time composes entropy.

Compared with prior art, the present invention has following beneficial effect:

First, the present invention is by transformer box wall vibration signal during on-line monitoring transformer sudden short circuit, and from transformer sudden short circuit state occurs second time, by when this sudden short circuit state lower gather vibration signal steady-state vibration signal aspect spectrum entropy and the transformer steady-state vibration signal aspect when there is sudden short circuit state for the first time compose entropy and compare differentiation, judge the duty of Transformer Winding, so that can Timeliness coverage problem, transformer is overhauled in time, there is advantage accurately and efficiently.

Second, the present invention constructs white noise of the narrowband signal according to the spectral characteristic of transformer vibration signal in step 2, vibration signal use experience mode decomposition is effectively inhibit to be modal overlap phenomenon in intrinsic mode function decomposable process, improve discomposing effect, improve accuracy winding state under sudden short circuit judged according to steady-state vibration signal aspect spectrum entropy.

Accompanying drawing explanation

Below in conjunction with the drawings and specific embodiments, the present invention is described in further detail:

Fig. 1 be the transformer winding state that monitors in the present embodiment good time vibrational waveform;

Vibrational waveform when Fig. 2 is the transformer winding state deterioration monitored in the present embodiment;

Fig. 3 show transformer winding state in the present embodiment good and loosening time morphology spectrum.

Embodiment

With certain Utilities Electric Co. 35kV transformer for subjects carries out short-circuit impact test.Test mesolow short circuit in winding, high pressure B phase winding is energized, and carries out 3 short-circuit impacts, records the vibration signal in the test of each short-circuit impact, Fig. 1 show transformer winding state good time vibrational waveform, Fig. 2 show Transformer Winding duty worsen time vibrational waveform.

Winding duty when the present invention judges transformer short-circuit according to the following step:

(1) real-time short circuit monitoring is carried out to transformer, and, whenever transformer generation sudden short circuit being detected, gather one section of continuous print vibration signal of transformer box wall, and when calculating transformer generation sudden short circuit each time by following steps (2) to the method for step (6) gather vibration signal steady-state vibration signal aspect compose entropy; Wherein, during transformer generation sudden short circuit to transformer box wall gather vibration signal start time should as far as possible close to the initial time of transformer generation sudden short circuit, transformer box wall gathered to end time of vibration signal and judges by existing vibration detection device, with the complete vibration signal collecting transformer box wall and produce due to transformer generation sudden short circuit as far as possible; In this example, gathering the transformer generation sudden short circuit moment plays vibration signal in 0.45s, with the complete vibration signal collecting transformer box wall and produce due to transformer generation sudden short circuit;

The vibration signal collected from transformer box wall during the above-mentioned generation sudden short circuit of transformer is each time expressed as x ₀(t).

(2) according to following step according to vibration signal x ₀(t) structure white noise of the narrowband signal g (t):

2a. is to vibration signal x ₀t () carries out Fourier transform, obtain vibration signal x ₀the spectrum distribution of (t); In this step, Fourier transform is mathematical method conventional in this area, and therefore inventor is no longer described in detail at this;

2b. is according to vibration signal x ₀t the spectrum distribution of () generates white noise signal s ₀(t), wherein white noise signal s ₀t the expression formula of the amplitude A of () is

In formula, K is coefficient, and COEFFICIENT K value is vibration signal x ₀t 1/10 of () amplitude mean value, gets 0.4; N _ffor vibration signal x ₀highest frequency component f in (t) spectrum component _hwith the ratio of 50Hz; f _i(i=1,2 ..., N _f) be vibration signal x ₀each frequency component of (t); a _i(i=1,2 ..., N _f) be vibration signal x ₀t amplitude that each frequency component of () is corresponding;

2c. uses Butterworth bandpass filter to white noise signal g ₀t () carries out filtering, obtain white noise of the narrowband signal s (t); The transport function expression formula of described Butterworth bandpass filter is

In formula, ω _cl=K _l2 π f _lfor low-frequency cut-off frequency, f _lfor vibration signal x ₀lowest frequency components in (t) spectrum component, K _lfor low frequency bandwidth factor; ω _ch=K _h2 π f _hfor high-frequency cut-off frequency, K _hfor high frequency bandwidth coefficient; M is filter order, gets 4; Wherein, K _l=1.5, K _h=0.8; ω=2 π f is angular frequency, f=50;

(3) by vibration signal x ₀t () is added with white noise of the narrowband signal g (t), obtain superposed signal x (t), according to following step, superposed signal x (t) is decomposed into several IMF components:

3a., to superposed signal x (t) differentiate, obtains time series y (t);

The product that 3b. sequences y computing time (t) is adjacent 2

py _i(t)＝y _i(t)×y _i-1(t)

In formula, i=2,3 ..., N-1, wherein, N is that signal is counted;

3c. is according to product py _it () and time series y's (t) is positive and negative, look for all Local modulus maximas eb (t) of superposed signal x (t) and all local minizing points es (t) successively:

Work as py _iduring (t) < 0, if py _i(t) < 0 and y _i-1(t) < 0, then x _i-1t () is local minizing point; If py _i(t) < 0 and y _i-1(t) > 0, then x _i-1t () is Local modulus maxima;

Work as py _iduring (t) > 0, x _i-1t () is non-extreme point;

Work as py _iduring (t)=0, if y _i-1t ()=0, calculates 2 y _i(t) and y _i-2t the product of (), makes py _i(t) '=y _i(t) × y _i-2t (), if py _i(t) ' < 0 and y _i-2(t) < 0, then x _i-1t () is local minizing point; If py _i(t) ' < 0 and y _i-2(t) > 0, then x _i-1t () is Local modulus maxima; If y _i-2(t)=0, then x _i-1t () is non-extreme point;

All Local modulus maximas eb (t) and all local minizing points es (t) couple together with cubic spline functions s (t) and obtain coenvelope line e respectively by 3d. _max(t) and lower envelope line e _mint (), this cubic spline functions s (t) is each minizone [t at superposed signal x (t) _i, t _i+1] (i=1,2 ..., N-1) on be no more than the polynomial expression of three times, its expression formula is

In formula, m _iand m _i+1for cubic spline functions s (t) is at interval [t _i, t _i+1] second derivative values at two-end-point place;

3e. is according to the coenvelope line e tried to achieve _max(t) and lower envelope line e _mint () calculates average m (t)=(e of upper and lower envelope _max(t)+e _min(t))/2, superposed signal x (t) is deducted m (t), obtains a new time series y ₁(t);

3f. judges time series y ₁t whether () meet following two conditions simultaneously:

A., in whole signal length, the number of extreme point and zero crossing only must differ one equal or at the most;

B. at any time, the coenvelope line defined by maximum point and the mean value of lower envelope line defined by minimum point are zero;

If meet above-mentioned two conditions, then y simultaneously ₁t () is intrinsic mode function component; If above-mentioned two conditions can not be met, then by y simultaneously ₁t (), as an original component, repeats abovementioned steps 3a ~ 3e, until time series y ₁t () meets above-mentioned two conditions simultaneously, will meet the y of above-mentioned two conditions simultaneously ₁t () is designated as c _i(t), then c _ia t intrinsic mode function component that () is superposed signal x (t), and i=1,2 ..., N _h, N _hfor the intrinsic mode function component total quantity of superposed signal x (t);

3g. is by c _it () separates from superposed signal x (t), obtain difference signal r _i(t) be

r _i(t)＝x(t)-c _i(t)

3h. is by difference signal r _it () carrys out repetition above-mentioned steps 3a ~ 3g as pending signal substituting superposed signal x (t), until meet Stopping criteria, obtain whole N _hindividual intrinsic mode function component, Stopping criteria is: obtain new time series y _it () is narrow band signal;

Through above-mentioned steps, initial vibration signal x ₀t () has been broken down into 6 IMF components.

(4) vibration signal x is calculated ₀t () and step (3) decompose the related coefficient between each intrinsic mode function component of obtaining, and the intrinsic mode function component chosen corresponding to the maximum related coefficient of numerical value is eigen vibration natural mode of vibration component, is designated as IMF _m;

Wherein, described vibration signal x ₀(t) related coefficient cor of decomposing kth the intrinsic mode function component that obtain middle with step (3) _kcomputing formula be

In formula, x _0it () is vibration signal x ₀t () is at t+i/f ₀the amplitude in moment; for vibration signal x ₀the mean value of (t) amplitude; c _kit () is for a kth intrinsic mode function component is at t+i/f ₀the amplitude in moment; for a kth intrinsic mode function component c _kthe mean value of (t) amplitude; f ₀for vibration signal x ₀the sample frequency of (t);

(5) according to following step according to natural mode of vibration component IMF _mpeak change determination vibration signal x ₀the initial time t of (t) steady-state vibration process _bwith finish time t _e, and remember vibration signal x ₀t () is from initial time t _bto finish time t _etime period in part be steady-state vibration signal x _s(t):

5a. calculates eigen vibration natural mode of vibration component IMF _mspectrum distribution, the vibration frequency choosing amplitude maximum in spectrum component is characteristic frequency, is designated as f _m, be 100Hz herein;

5b. is from vibration signal x ₀t the initial time of () starts to calculate vibration signal x successively in chronological order ₀t (), at grid number Ng (t) in each moment, wherein, grid number Ng (t) refers to vibration signal x ₀t () is at time period [t, t+0.5T _m] mean value of interior variable quantity absolute value sum, its computing formula is

T _m＝1/f _m

In formula, N ₀for vibration signal x ₀t () is at time period [t, t+0.5T _m] in data length; x _0i(t) and x _{0 (i+1)}t () is respectively vibration signal at t+i/f ₀moment and t+ (i+1)/f ₀moment amplitude; f ₀for vibration signal x ₀t the sample frequency of () is 10240 herein; T _mfor characteristic frequency f _minverse, be 0.01 herein;

5c. is in chronological order successively to vibration signal x ₀t () screens at grid number Ng (t) in each moment, the principle of screening is: if grid number Ng (t) is sometime greater than default threshold value, then retain the grid number in this moment; If Ng (t) is sometime less than default threshold value, then the grid number in this moment is set to 0; Note vibration signal x ₀(t) screen after grid number for screening after grid number Ng'(t); Wherein, the threshold value value preset is 20% of grid number Ng (t) maximal value and rounds value;

5d. is from grid number Ng'(t after screening) in moment corresponding to first non-zero grid number carry out out Top-hat computing successively, first peak value detected is vibration signal x ₀the initial time t of (t) steady-state vibration process _b; From grid number Ng'(t after screening) in moment corresponding to last non-zero grid number oppositely carry out out Top-hat computing successively, first peak value detected is vibration signal x ₀the finish time t of (t) steady-state vibration process _e; The described expression formula opening Top-hat computing is

In formula, f is grid number Ng'(t); G is structural element; F Θ g carries out erosion operation for using structural element g to signal f; for using structural element g, dilation operation is carried out to signal f Θ g; The fundamental operation of mathematical morphology when described dilation operation and erosion operation, be mathematical method conventional in this area, therefore inventor is no longer described in detail at this;

(6) steady-state vibration signal x is calculated _st the morphology spectrum of () and morphology spectrum entropy, its computing formula is

q(λ)＝PS _f(λ,g)/S(f·λ _ming)

In formula, H (f/g) is morphology spectrum entropy; PS _f(λ, g) is morphology spectrum; F is grid number Ng'(t); G is structural element; λ is yardstick; λ _maxfor yardstick maximal value, be 51 herein; λ _minfor yardstick minimum value, be 1 herein.

Fig. 3 be in the present embodiment transformer winding state good and loosening time morphology spectrum.(7) according to steady-state vibration signal x _st the morphology spectrum entropy of () differentiates winding state: when the changes delta H of morphology spectrum entropy meets time, judge that the winding of transformer occurs when the short circuit of Y secondary burst occurs transformer to loosen or distortion, now need to carry out overhaul plan in time, avoid the formation of significant trouble; Wherein, H ₁the steady-state vibration signal aspect gathering vibration signal for institute when sudden short circuit occurs for transformer the 1st time composes entropy, H _ythe steady-state vibration signal aspect gathering vibration signal for institute when sudden short circuit occurs for transformer the Y time composes entropy.

As can be seen from Table 1, in the present embodiment, front 2 short-circuit impacts, the absolute value rate of the change of its morphology spectrum entropy is 2.97%, describes transformer winding state normal.During the 3rd Secondary Shocks test, the absolute value of the change of morphology spectrum entropy is 28.31%, shows that transformer winding state is abnormal, needs to overhaul in time it, avoid the formation of significant trouble.

Table 1

Short-circuit impact number of times	Morphology spectrum entropy	Morphology spectrum progressive increase rate/%
			1	0.70493	-
2	0.68402	2.97
			3	0.90449	28.31

The present invention does not limit to and above-mentioned embodiment; according to foregoing; according to ordinary technical knowledge and the customary means of this area; do not departing under the present invention's above-mentioned basic fundamental thought prerequisite; the present invention can also make the equivalent modifications of other various ways, replacement or change, all drops among protection scope of the present invention.

Claims (6)

1., based on a winding state diagnostic method under the transformer sudden short circuit of vibration signal morphology spectrum, comprise the following steps:

(1) real-time short circuit monitoring is carried out to transformer, and, whenever transformer generation sudden short circuit being detected, gather one section of continuous print vibration signal of transformer box wall, and when calculating transformer generation sudden short circuit each time by following steps (2) to the method for step (6) gather vibration signal steady-state vibration signal aspect compose entropy;

Wherein, the vibration signal collected from transformer box wall during the above-mentioned generation sudden short circuit of transformer is each time expressed as x ₀(t);

(2) vibration signal x is used ₀t () constructs white noise of the narrowband signal g (t);

(3) by vibration signal x ₀t () is added with white noise of the narrowband signal g (t), obtain superposed signal x (t), and superposed signal x (t) is decomposed into several intrinsic mode function components;

(4) vibration signal x is calculated ₀t () and step (3) decompose the related coefficient between each intrinsic mode function component of obtaining, and the intrinsic mode function component chosen corresponding to the maximum related coefficient of numerical value is eigen vibration natural mode of vibration component, is designated as IMF _m;

Wherein, described vibration signal x ₀(t) related coefficient cor of decomposing kth the intrinsic mode function component that obtain middle with step (3) _kcomputing formula be

${\mathrm{cor}}_k=\frac{\underset{i=1}{\overset{N}{\&Sigma;}}\left(x_{0i}\right(t)-{\stackrel{\&OverBar;}{x}}_0(t\left)\right)\left(c_{ki}\right(t)-{\stackrel{\&OverBar;}{c}}_k(t\left)\right)}{\sqrt{\underset{i=1}{\overset{N}{\&Sigma;}}{\left(x_{0i}\right(t)-{\stackrel{\&OverBar;}{x}}_0(t\left)\right)}^2\&CenterDot;\underset{i=1}{\overset{N}{\&Sigma;}}{\left(c_{ki}\right(t)-{\stackrel{\&OverBar;}{c}}_k(t\left)\right)}^2}},k=1,...,N_h$

In formula, x _0it () is vibration signal x ₀t () is at t+i/f ₀the amplitude in moment; for vibration signal x ₀the mean value of (t) amplitude; c _kit () is for a kth intrinsic mode function component is at t+i/f ₀the amplitude in moment; for a kth intrinsic mode function component c _kthe mean value of (t) amplitude; f ₀for vibration signal x ₀the sample frequency of (t);

(5) according to natural mode of vibration component IMF _mpeak change determination vibration signal x ₀the initial time t of (t) steady-state vibration process _bwith finish time t _e, and remember vibration signal x ₀t () is from initial time t _bto finish time t _etime period in part be steady-state vibration signal x _s(t);

(6) steady-state vibration signal x is calculated _st the morphology spectrum of () and morphology spectrum entropy, its computing formula is

$H(f/g)=-\frac{\underset{\mathrm{\&lambda;}={\mathrm{\&lambda;}}_{min}}{\overset{{\mathrm{\&lambda;}}_{\max }}{\&Sigma;}}q\left(\mathrm{\&lambda;}\right)\lg q\left(\mathrm{\&lambda;}\right)}{\lg ({\mathrm{\&lambda;}}_{max}-{\mathrm{\&lambda;}}_{\min }+1)}$

q(λ)＝PS _f(λ,g)/S(f·λ _ming)

In formula, H (f/g) is morphology spectrum entropy; PS _f(λ, g) is morphology spectrum; F is grid number Ng'(t); G is structural element; λ is yardstick; λ _maxfor yardstick maximal value; λ _minfor yardstick minimum value;

(7) according to steady-state vibration signal x _st the morphology spectrum entropy of () differentiates winding state: when the changes delta H of morphology spectrum entropy meets time, judge that the winding of transformer occurs when the short circuit of Y secondary burst occurs transformer to loosen or distortion; Wherein, H ₁the steady-state vibration signal aspect gathering vibration signal for institute when sudden short circuit occurs for transformer the 1st time composes entropy, H _ythe steady-state vibration signal aspect gathering vibration signal for institute when sudden short circuit occurs for transformer the Y time composes entropy.

2. method according to claim 1, is characterized in that: in described step (1a), described transformer box wall is gathered to the initial time of start time close to described transformer generation sudden short circuit of vibration signal.

3. method according to claim 1, is characterized in that: in described step (2), use vibration signal x ₀t method that () constructs white noise of the narrowband signal g (t) comprises the following steps:

2a. is to vibration signal x ₀t () carries out Fourier transform, obtain vibration signal x ₀the spectrum distribution of (t);

2b. is according to vibration signal x ₀t the spectrum distribution of () generates white noise signal s ₀(t), wherein white noise signal s ₀t the expression formula of the amplitude A of () is

$A=K\underset{i=1}{\overset{N_f}{\&Sigma;}}\left(a_i{f}_i\right)/\underset{i=1}{\overset{N_f}{\&Sigma;}}{f}_i,i=1,2,...,N_f$

In formula, K is coefficient, and COEFFICIENT K value is vibration signal x ₀1/10 of (t) amplitude mean value; N _ffor vibration signal x ₀highest frequency component f in (t) spectrum component _hwith the ratio of 50Hz; f _i(i=1,2 ..., N _f) be vibration signal x ₀each frequency component of (t); a _i(i=1,2 ..., N _f) be vibration signal x ₀t amplitude that each frequency component of () is corresponding;

2c. uses Butterworth bandpass filter to white noise signal s ₀t () carries out filtering, obtain white noise of the narrowband signal g (t); The transport function expression formula of described Butterworth bandpass filter is

${\left|H\left(\mathrm{\&omega;}\right)\right|}^2=\frac{1}{\&lsqb;1+{\left(\frac{\mathrm{\&omega;}}{{\mathrm{\&omega;}}_{cl}}\right)}^{2m}\&rsqb;\&lsqb;1+{\left(\frac{{\mathrm{\&omega;}}_{ch}}{\mathrm{\&omega;}}\right)}^{2m}\&rsqb;}$

In formula, ω _cl=K _l2 π f _lfor low-frequency cut-off frequency, f _lfor vibration signal x ₀lowest frequency components in (t) spectrum component, K _lfor low frequency bandwidth factor; ω _ch=K _h2 π f _hfor high-frequency cut-off frequency, f _hfor vibration signal x ₀highest frequency component in (t) spectrum component, K _hfor high frequency bandwidth coefficient; M is filter order; ω=2 π f is angular frequency, f=50.

4. method according to claim 1, is characterized in that: in described step (3), method superposed signal x (t) being decomposed into several intrinsic mode function components comprises the following steps:

3a., to superposed signal x (t) differentiate, obtains time series y (t);

The product that 3b. sequences y computing time (t) is adjacent 2

py _i(t)＝y _i(t)×y _i-1(t)

In formula, i=2,3 ..., N-1, wherein, N is that signal is counted;

3c. is according to product py _it () and time series y's (t) is positive and negative, look for all Local modulus maximas eb (t) of superposed signal x (t) and all local minizing points es (t) successively:

Work as py _iduring (t) < 0, if py _i(t) < 0 and y _i-1(t) < 0, then x _i-1t () is local minizing point; If py _i(t) < 0 and y _i-1(t) > 0, then x _i-1t () is Local modulus maxima;

Work as py _iduring (t) > 0, x _i-1t () is non-extreme point;

Work as py _iduring (t)=0, if y _i-1t ()=0, calculates 2 y _i(t) and y _i-2t the product of (), makes py _i(t) '=y _i(t) × y _i-2t (), if py _i(t) ' < 0 and y _i-2(t) < 0, then x _i-1t () is local minizing point; If py _i(t) ' < 0 and y _i-2(t) > 0, then x _i-1t () is Local modulus maxima; If y _i-2(t)=0, then x _i-1t () is non-extreme point;

All Local modulus maximas eb (t) and all local minizing points es (t) couple together with cubic spline functions s (t) and obtain coenvelope line e respectively by 3d. _max(t) and lower envelope line e _mint (), described cubic spline functions s (t) is each minizone [t at superposed signal x (t) _i, t _i+1] (i=1,2 ..., N-1) on be no more than the polynomial expression of three times, its expression formula is

$\begin{array}{c}s\left(t\right)=m_i\frac{{(t_{i+1}-t)}^3}{6(t_{i+1}-t_i)}+m_{i+1}\frac{{(t-t_i)}^3}{6(t_{i+1}-t_i)}+\frac{x_{i+1}\left(t\right)-x_i\left(t\right)}{t_{i+1}-t_i}\\ -\frac{t_{i+1}-t_i}{6}(m_{i+1}-m_i)+x_i\left(t\right)-m_i\frac{{(t_{i+1}-t_i)}^2}{6}\end{array}$

In formula, m _iand m _i+1for cubic spline functions s (t) is at interval [t _i, t _i+1] second derivative values at two-end-point place;

3e. is according to the coenvelope line e tried to achieve _max(t) and lower envelope line e _mint () calculates average m (t)=(e of upper and lower envelope _max(t)+e _min(t))/2, superposed signal x (t) is deducted m (t), obtains a new time series y ₁(t);

3f. judges time series y ₁t whether () meet following two conditions simultaneously:

A., in whole signal length, the number of extreme point and zero crossing only must differ one equal or at the most;

B. at any time, the coenvelope line defined by maximum point and the mean value of lower envelope line defined by minimum point are zero;

If meet above-mentioned two conditions, then y simultaneously ₁t () is intrinsic mode function component; If above-mentioned two conditions can not be met, then by y simultaneously ₁t (), as an original component, repeats abovementioned steps 3a ~ 3e, until time series y ₁t () meets above-mentioned two conditions simultaneously, will meet the y of above-mentioned two conditions simultaneously ₁t () is designated as c _i(t), then c _ia t intrinsic mode function component that () is superposed signal x (t), and i=1,2 ..., N _h, N _hfor the intrinsic mode function component total quantity of superposed signal x (t);

3g. is by c _it () separates from superposed signal x (t), obtain difference signal r _i(t) be

r _i(t)＝x(t)-c _i(t)

3h. is by difference signal r _it () carrys out repetition above-mentioned steps 3a ~ 3g as pending signal substituting superposed signal x (t), until meet Stopping criteria, obtain whole N _hindividual intrinsic mode function component, described Stopping criteria is: obtain new time series y _it () is narrow band signal;

Through above-mentioned steps, vibration signal x ₀t () has been broken down into several intrinsic mode function components and residual signal sum, its expression formula is

x ₀(t)＝Σc _i(t)+r _i(t)。

5. method according to claim 1, is characterized in that: in described step (5), according to natural mode of vibration component IMF _mpeak change determination vibration signal x ₀the initial time t of (t) steady-state vibration process _b, finish time t _ewith steady-state vibration signal x _st the method for () comprises the following steps:

5a. calculates eigen vibration natural mode of vibration component IMF _mspectrum distribution, the vibration frequency choosing amplitude maximum in spectrum component is characteristic frequency, is designated as f _m;

5b. is from vibration signal x ₀t the initial time of () starts to calculate vibration signal x successively in chronological order ₀t (), at grid number Ng (t) in each moment, wherein, grid number Ng (t) refers to vibration signal x ₀t () is at time period [t, t+0.5T _m] mean value of interior variable quantity absolute value sum, its computing formula is

$Ng\left(t\right)=\frac{1}{0.5T_m}\underset{i=1}{\overset{N_0}{\&Sigma;}}|x_{0i}\left(t\right)-x_{0(i+1)}\left(t\right)|$

T _m＝1/f _m

In formula, N ₀for vibration signal x ₀t () is at time period [t, t+0.5T _m] in data length; x _0i(t) and x _{0 (i+1)}t () is respectively vibration signal at t+i/f ₀moment and t+ (i+1)/f ₀moment amplitude; f ₀for vibration signal x ₀the sample frequency of (t); T _mfor characteristic frequency f _minverse;

5c. is in chronological order successively to vibration signal x ₀t () screens at grid number Ng (t) in each moment, the principle of screening is: if grid number Ng (t) is sometime greater than default threshold value, then retain the grid number in this moment; If Ng (t) is sometime less than default threshold value, then the grid number in this moment is set to 0; Note vibration signal x ₀(t) screen after grid number for screening after grid number Ng'(t);

5d. is from grid number Ng'(t after screening) in moment corresponding to first non-zero grid number carry out out Top-hat computing successively, first peak value detected is vibration signal x ₀the initial time t of (t) steady-state vibration process _b; From grid number Ng'(t after screening) in moment corresponding to last non-zero grid number oppositely carry out out Top-hat computing successively, first peak value detected is vibration signal x ₀the finish time t of (t) steady-state vibration process _e; The described expression formula opening Top-hat computing is

In formula, f is grid number Ng'(t); G is structural element; for using structural element g, erosion operation is carried out to signal f; for using structural element g to signal carry out dilation operation.

6. method according to claim 5, is characterized in that: in described step 5c, and described default threshold value value is 20% of grid number Ng (t) maximal value and rounds value.