CN104964736A - Optical fiber invasion vibration source identification method based on time-frequency characteristic EM classification - Google Patents

Abstract

The invention provides a fiber intrusion identification method based on a time-frequency two-dimensional characteristic EM (maximum expectation) classification algorithm aiming at a fiber intrusion early warning system, which comprises the following steps: calculating the duty ratio of the vibration data obtained through detection, and extracting the time domain characteristics of the signal; carrying out Fourier transform on the signal obtained after wavelet denoising, calculating a frequency center, and extracting signal frequency domain characteristics; and finally, the time domain and frequency domain characteristics obtained by the method are used as the input of a two-dimensional EM classifier to identify the vibration source of the optical fiber vibration signal. The method identifies the intrusion signal in time-frequency two-dimension, can effectively distinguish the manual signal from the mechanical signal, and has higher accuracy.

Description

Based on the fiber optic intrusion Recognition of Vibration Sources method that time-frequency characteristic EM classifies

Technical field

The present invention relates to the fiber optic intrusion Recognition of Vibration Sources method based on time-frequency two-dimensional feature EM (greatest hope) sorting algorithm of the Recognition of Vibration Sources of fiber optic intrusion system.

Background technology

Along with the fast development of global economy, society increases increasingly to the demand of the energy especially hydrocarbon resources.In national energy strategy, the construction and development of Oil & Gas Storage be related to for the development of the national economy and social development provide for a long time, stablize, the strategy overall situation of economic, safe energy safeguard.Nowadays, underground oil and gas conveyance conduit has become the main artery of energy transport, and the problem of pipe safety protection is put in face of people with also becoming increasingly conspicuous.Pipeline, once leakage is very easily blasted, not only affects the normal transport of the energy, also brings about great losses to the life of the country and people masses, property.Therefore, the early warning system of pipe safety is widely used background.

Carrying out using optical fiber sensing system, on the basis of monitoring in real time, classifying to detected vibration signal, identifying the external event source causing vibration, being convenient to adopt an effective measure in time, preventing harmful generation invaded.By the vibration event of optical fiber sensing system detecting optical cable periphery, gather the various vibration signals of petroleum pipe line periphery, extract signal characteristic parameter, the classification of realize target and identification.In the face of the vibration signal of large amount of complex, how accurately to identify that target vibration source is the difficult point of safety pre-warning system research.Recognition of Vibration Sources is behavior based on vibration source and attributive character thereof, take computing machine as instrument, adopts pattern recognition theory, sets up a special kind of skill of vibration signal and vibration source corresponding relation.System carries out pre-service, feature extraction and identification to the vibration signal that Fiber Duct collects, and carries out safe early warning according to the type of its feature determination destructive insident, thus realizes ensureing oil-gas pipeline safety, the object prevented trouble before it happens.

The subject matter that existing research exists lacks suitable model, particularly do not set up suitable Signal analysis model, therefore, needs to set up a kind of effective model to realize the identification of vibration signal, to reduce the error rate of Recognition of Vibration Sources.

Summary of the invention

The present invention relates to the Recognition of Vibration Sources of fiber optic intrusion system, it is identified and carries out Intrusion Signatures extraction to system by time domain, frequency domain two dimension, thus obtained vibration source attribute and specifying information.

The fiber optic intrusion Recognition of Vibration Sources method that the present invention is based on time-frequency two-dimensional feature EM sorting algorithm, in order to solve the identification to fiber optic intrusion signal, is determined to invade signal type.

Based on the fiber optic intrusion Recognition of Vibration Sources method of time-frequency two-dimensional feature EM sorting algorithm, it is characterized in that comprising: detect vibration signal, the data of the vibration position detected are put 1, other position datas set to 0, and calculate dutycycle and are:

$\mathrm{\η}=\frac{m}{m+n}---\left(1\right)$

Wherein, m is the number of in every frame data 1, and n is the number of in every frame data 0.

Afterwards Fourier transform is carried out to the signal after Wavelet Denoising Method process, calculated rate center, to extract signal frequency domain feature.First energy is asked to vibration signal

$E=\left|\right|s\left(t\right)|{|}^2={\&Integral;\left|s\left(t\right)\right|}^{2}dt=\frac{1}{2\mathrm{\&pi;}}{\&Integral;\left|S\left(j\mathrm{\&Omega;}\right)\right|}^{2}d\mathrm{\&Omega;}<\mathrm{\&infin;}---\left(2\right)$

Wherein, t is the time, and E is signal energy, and s (t) is vibration signal, the Fourier transform that S (j Ω) is s (t).The center frequency calculating signal s (t) is afterwards

$\mathrm{\&zeta;}\left(\mathrm{\&Omega;}\right)=\frac{1}{2\mathrm{\&pi;}E}\&Integral;\mathrm{\&Omega;}|S\left(j\mathrm{\&Omega;}\right){|}^{2}d\mathrm{\&Omega;}---\left(3\right)$

Wherein, Ω is frequency, the center frequency that ζ (Ω) is s (t), and E is signal energy, the Fourier transform that S (j Ω) is s (t).

The time domain obtained using said method, frequency domain character are carried out fiber-optic vibration signal Recognition of Vibration Sources as the input of two-dimentional EM sorter.

Accompanying drawing explanation

The enforcement of Fig. 1 the inventive method and proof procedure;

Fig. 2 temporal signatures extracts process flow diagram;

Fig. 3 frequency domain character extracts process flow diagram;

Fig. 4 frequency domain character extracts result figure;

Fig. 5 EM classifier design block diagram;

Fig. 6 EM sorter is for the recognition result of manual signal and mechanical signal.

Specific embodiments

Below in conjunction with accompanying drawing, technical scheme is according to an embodiment of the invention described.

Fig. 1 is the overall procedure of time-frequency two-dimensional recognition methods according to an embodiment of the invention.In this embodiment, say that the object identified comprises: manual signal, it is for owing to using on-electric class instrument and the vibration signal that produces, as pick; Mechanical signal, its vibration signal for producing owing to using electronic class instrument, as electric drill, electric pick.

The time-frequency two-dimensional recognizer of embodiment as shown in Figure 1 comprises:

S101: extract signal temporal signatures, calculates vibration data dutyfactor value;

S102: extract signal frequency domain feature, vibration signal is carried out Fourier transform and calculated rate center;

S103: the time domain extracted, frequency domain character are carried out fiber-optic vibration signal Recognition of Vibration Sources as the input of two-dimentional EM sorter.

According to an embodiment of the inventionly carry out the process of temporal signatures extraction as shown in Figure 2 to signal, it comprises:

The data of the vibration position detected are put 1 by S201: detect vibration signal, and the data of other positions set to 0;

S202: add up the number m of the 1 and number n of 0 in every segment data;

S203: computed duty cycle

S204: using vectorial as temporal signatures stored in matrix for the dutycycle numerical value calculated.

Use wavelet analysis denoising is carried out to vibration data, afterwards frequency domain character extraction is carried out to the signal after de-noising, frequency domain character leaching process according to an embodiment of the invention as shown in Figure 3:

S301: adopt the conversion of N=1024 point quick Fourier to carry out analysis of spectrum to the vibration signal after de-noising, signal length is 5s, and sample frequency is 1kHz, mechanical signal and manual signal time-domain diagram and spectrogram are as shown in Figure 4;

S302: the center frequency calculating vibration signal under all kinds of label, first asks energy to vibration signal

$E=\left|\right|s\left(t\right)|{|}^2={\&Integral;\left|s\left(t\right)\right|}^{2}dt=\frac{1}{2\mathrm{\&pi;}}{\&Integral;\left|S\left(j\mathrm{\&Omega;}\right)\right|}^{2}d\mathrm{\&Omega;}<\mathrm{\&infin;}---\left(4\right)$

Wherein, t is the time, and E is signal energy, and s (t) is vibration signal, the Fourier transform that S (j Ω) is s (t);

Then the center frequency calculating signal x (t) is

$\mathrm{\&zeta;}\left(\mathrm{\&Omega;}\right)=\frac{1}{2\mathrm{\&pi;}E}\&Integral;\mathrm{\&Omega;}|S\left(j\mathrm{\&Omega;}\right){|}^{2}d\mathrm{\&Omega;}---\left(5\right)$

Wherein, Ω is frequency, the center frequency that ζ (Ω) is s (t), and E is signal energy, the Fourier transform that S (j Ω) is s (t);

S303: using vectorial as frequency domain character stored in matrix for the dutycycle calculated.

Time domain obtained above, frequency domain character vector are classified as the input of two-dimentional EM sorter.As shown in Figure 5, it comprises classification process according to an embodiment of the invention:

First, dutycycle, center frequency two feature are generated two-dimensional feature vector also as the sample to be sorted of input sorter, namely $x=\left[\begin{array}{c}\mathrm{\&eta;}\\ \mathrm{\&zeta;}\end{array}\right];$

Then, initialization distribution parameter θ and two class data accounting a, wherein parameter θ comprises average μ and covariance cov, namely $\mathrm{\&mu;}=\left[\begin{array}{c}E\left(\mathrm{\&eta;}\right)\\ E\left(\mathrm{\&zeta;}\right)\end{array}\right]$ With $\mathrm{cov}=\left[\begin{array}{cc}\mathrm{cov}(\mathrm{\&eta;},\mathrm{\&eta;})& \mathrm{cov}(\mathrm{\&eta;},\mathrm{\&zeta;})\\ \mathrm{cov}(\mathrm{\&zeta;},\mathrm{\&eta;})& \mathrm{cov}(\mathrm{\&zeta;},\mathrm{\&zeta;})\end{array}\right];$ If the distribution parameter θ of first kind vibration ₁for ${\mathrm{\&mu;}}_1=\left[\begin{array}{c}{\mathrm{\&mu;}}_{11}\\ {\mathrm{\&mu;}}_{12}\end{array}\right]$ With ${\mathrm{cov}}_1=\left[\begin{array}{cc}{c}_{11}& c_{12}\\ c_{21}& c_{22}\end{array}\right],$ The distribution parameter θ of Equations of The Second Kind vibration ₂for ${\mathrm{\&mu;}}_2=\left[\begin{array}{c}{\mathrm{\&mu;}}_{21}\\ {\mathrm{\&mu;}}_{22}\end{array}\right]$ With ${\mathrm{cov}}_2=\left[\begin{array}{cc}{d}_{11}& d_{12}\\ d_{21}& d_{22}\end{array}\right],$ And set two class vibration data accountings as $a=\left[\begin{array}{c}{k}_1\\ k_2\end{array}\right],$ Wherein k ₁+ k ₂=1;

Finally, carry out calculation expectation (Expectation) and maximization (Maximization) two step, this flow process according to an embodiment of the invention comprises:

S501: calculation expectation step (E step): the posterior probability (i.e. the expectation of recessive variable) calculating recessive variable according to the model parameter of initial parameter value or last iteration.Current estimated value as recessive variable:

Q _i(z ⁽ⁱ⁾):＝p(z ⁽ⁱ⁾|x ⁽ⁱ⁾；θ) (6)

Wherein, x ⁽ⁱ⁾for sample data to be sorted, z ⁽ⁱ⁾for each sample x ⁽ⁱ⁾corresponding classification, Q _i(z ⁽ⁱ⁾) be the posterior probability of recessive variable;

S502: maximization steps (M step): likelihood function is maximized to obtain new parameter value:

$\mathrm{\&theta;}:=\arg \underset{\mathrm{\&theta;}}{max}\underset{i}{\&Sigma;}\underset{z^{\left(i\right)}}{\&Sigma;}{Q}_i\left(z^{\left(i\right)}\right)log\frac{p(x^{\left(i\right)},z^{\left(i\right)};\mathrm{\&theta;})}{Q_i\left(z^{\left(i\right)}\right)}---\left(7\right)$

Now, judge whether parameter (comprising average μ and covariance cov) restrains, if do not restrain, then return S501 and repeat this process, if restrain, algorithm terminates, output estimation parameter value and classification results.

The present inventor, for above-mentioned time-frequency two-dimensional recognition methods, carries out classification to the manual signal of actual measurement and mechanical signal and emulates; Wherein manual signal is the vibration signal because people uses on-electric class instrument to produce, and as pick, mechanical signal is the vibration signal because people uses electronic class instrument to produce, as electric drill, electric pick.Wherein, first respectively time domain is carried out to it, frequency domain character extracts and obtain the two-dimensional feature vector of manual signal and mechanical signal, then proper vector is inputted EM algorithm classification device it is identified.Recognition result as shown in Figure 6.In the figure, horizontal ordinate represents frequency domain character, and ordinate represents temporal signatures.As can be seen from this simulation result, effectively manual signal and mechanical signal can be distinguished by time-frequency two-dimensional recognition methods, indicate that the present invention has significant effect.

Compared with existing detection method, advantage of the present invention comprises:

(1) method of the present invention can effectively realize fiber optic intrusion identification;

(2) method of the present invention can remove signal major part noise by Wavelet-based Denoising effect, and more accurately identifying signal provides convenient;

(3) manual signal and mechanical signal can be differentiated by the identification of time-frequency two-dimensional feature by method of the present invention effectively, and accuracy is higher.

Claims (5)

1., based on the fiber optic intrusion Recognition of Vibration Sources method of time-frequency two-dimensional feature greatest hope sorting algorithm, it is characterized in that comprising:

A) calculate the dutycycle of the vibration data obtained afterwards after testing, extract signal temporal signatures;

B) Fourier transform is carried out to vibration signal, calculated rate center, extract signal frequency domain feature;

C) time domain obtained using said method, frequency domain character are carried out fiber-optic vibration signal Recognition of Vibration Sources as the input of two-dimentional EM sorter.

2. method according to claim 1, is characterized in that described steps A) comprise further:

Detect vibration signal, the data of the vibration position detected are put 1, other position datas set to 0, and calculate dutycycle:

$\mathrm{\&eta;}=\frac{m}{m+n}---\left(1\right)$

Wherein, m is the number of in every frame data 1, and n is the number of in every frame data 0.

3. method according to claim 1, is characterized in that described step B) comprise further:

Wavelet Denoising Method process is carried out to vibration signal, carries out fast fourier transform afterwards, and obtain center frequency, wherein saidly ask the operation of center frequency to comprise:

Energy is asked to vibration signal

$E=\left|\right|s\left(t\right)|{|}^2={\&Integral;\left|s\left(t\right)\right|}^{2}dt=\frac{1}{2\mathrm{\&pi;}}{\&Integral;\left|S\left(j\mathrm{\&Omega;}\right)\right|}^{2}d\mathrm{\&Omega;}<\mathrm{\&infin;}---\left(2\right)$

Wherein, t is the time, and E is signal energy, and s (t) is vibration signal, the Fourier transform that S (j Ω) is s (t);

Calculate the center frequency of signal s (t)

$\mathrm{\&zeta;}\left(\mathrm{\&Omega;}\right)=\frac{1}{2\mathrm{\&pi;}E}\&Integral;\mathrm{\&Omega;}|S\left(j\mathrm{\&Omega;}\right){|}^{2}d\mathrm{\&Omega;}---\left(3\right)$

Wherein, Ω is frequency, the center frequency that ζ (Ω) is s (t), and E is signal energy, the Fourier transform that S (j Ω) is s (t).

4., according to the method for Claims 2 or 3, it is characterized in that described step C) comprise further:

The dutycycle extracted from different types of signal and center frequency two feature are generated two-dimensional feature vector also as the sample to be sorted of input sorter, namely $x=\left[\begin{array}{c}\mathrm{\&eta;}\\ \mathrm{\&zeta;}\end{array}\right];$

Adopt EM algorithm classification device to identify proper vector, wherein said EM algorithm comprises:

D1) initialization distribution parameter θ and two class data accounting a, wherein parameter θ comprises average μ and covariance cov, namely $\mathrm{\&mu;}=\left[\begin{array}{c}E\left(\mathrm{\&eta;}\right)\\ E\left(\mathrm{\&zeta;}\right)\end{array}\right]$ With $\mathrm{cov}=\left[\begin{array}{cc}\mathrm{cov}(\mathrm{\&eta;},\mathrm{\&eta;})& \mathrm{cov}(\mathrm{\&eta;},\mathrm{\&zeta;})\\ \mathrm{cov}(\mathrm{\&zeta;},\mathrm{\&eta;})& \mathrm{cov}(\mathrm{\&zeta;},\mathrm{\&zeta;})\end{array}\right];$ If the distribution parameter θ of first kind vibration ₁for ${\mathrm{\&mu;}}_1=\left[\begin{array}{c}{\mathrm{\&mu;}}_{11}\\ {\mathrm{\&mu;}}_{12}\end{array}\right]$ With ${\mathrm{cov}}_1=\left[\begin{array}{cc}{c}_{11}& c_{12}\\ c_{21}& c_{22}\end{array}\right],$ The distribution parameter θ of Equations of The Second Kind vibration ₂for ${\mathrm{\&mu;}}_2=\left[\begin{array}{c}{\mathrm{\&mu;}}_{21}\\ {\mathrm{\&mu;}}_{22}\end{array}\right]$ With ${\mathrm{cov}}_2=\left[\begin{array}{cc}{d}_{11}& d_{12}\\ d_{21}& d_{22}\end{array}\right],$ And set two class vibration data accountings as $a=\left[\begin{array}{c}{k}_1\\ k_2\end{array}\right],$ Wherein k ₁+ k ₂=1;

D2) posterior probability of recessive variable is calculated according to the model parameter of initial parameter value or last iteration, i.e. the expectation of recessive variable, the current estimated value as hidden variable:

Q _i(z ⁽ⁱ⁾):＝p(z ⁽ⁱ⁾|x ⁽ⁱ⁾；θ) (4)

Wherein, x ⁽ⁱ⁾for sample data to be sorted, z ⁽ⁱ⁾for each sample x ⁽ⁱ⁾corresponding classification, Q _i(z ⁽ⁱ⁾) be the posterior probability of recessive variable;

D3) likelihood function is maximized to obtain that comprise average μ and covariance cov, new parameter value:

$\mathrm{\&theta;}:=\arg \underset{\mathrm{\&theta;}}{\max }\underset{i}{\mathrm{\&Sigma;}}\underset{z^{\left(i\right)}}{\mathrm{\&Sigma;}}{Q}_i\left(z^{\left(i\right)}\right)\log \frac{p(x^{\left(i\right)},z^{\left(i\right)};\mathrm{\&theta;})}{Q_i\left(z^{\left(i\right)}\right)}---\left(5\right)$

D4) judge whether this new parameter value comprising average μ and covariance cov restrains, if do not restrain, returns D2) step, if restrain, algorithm terminates.

5. method according to claim 4, is characterized in that described different types of signal comprises mechanical signal and manual signal two class signal.