Method for detecting pipe burst in water distribution systems based on pressure disturbance extraction

Abstract

The present invention introduced a method for detecting pipe bursts in water distribution systems based on pressure disturbance extraction, including following steps: (1) collecting and pre-processing monitored pressure data, and establishing a matrix of monitored pressures; (2) analyzing time-domain and frequency-domain features of the matrix of monitored pressures by Fourier transform, extracting disturbances, and generating pressure disturbances matrices; (3) identifying outliers in the pressure disturbances matrices by isolation forest algorithm; and (4) further identifying detected outliers by calculating and qualitative index A and quantitative index B and outputting a result of pipe burst detection. Compared with previous methods, the pipe burst detection method introduced in the present invention is accurate and reliable, and is more applicable to large-scale complex water distribution systems.

Claims

1. A method for detecting pipe bursts in water distribution systems based on pressure disturbance extraction, the method comprising: (1) collecting and pre-processing monitored pressure data, and establishing pressure matrices; (11) numbering date and time of historical data; (12) for a moment t=(x,y) to be detected, extracting SCADA monitored pressure data in m days before the moment t, wherein (x,y) represents the moment y in the x.sup.th day; (13) pre-processing the extracted monitored pressure data; (14) dividing the pre-processed pressure data at a sensor into m row vectors, each vector having a length of N, N being the number of samples at one sensor in a day, and storing the row vectors in pre-processed pressure matrices P, wherein, for the i.sup.th sensor, the corresponding matrix P.sub.1 is: $P_i=\left(\begin{array}{cccccc}{P}_{x-m,y+1}& \mathrm{.Math.}& p_{x-m,N}& P_{x-m+1,1}& \mathrm{.Math.}& p_{x-m+1,y}\\ p_{x-m,y+1}& \mathrm{.Math.}& P_{x-m-1}& p_{x-m+2,1}& \mathrm{.Math.}& P_{x-m+2,y}\\ \mathrm{.Math.}& \ddots & \mathrm{.Math.}& \mathrm{.Math.}& \ddots & \mathrm{.Math.}\\ P_{x-1,y+1}& \mathrm{.Math.}& P_{x-1,N}& P_{x,1}& \mathrm{.Math.}& P_{x,y}\end{array}\right),$ where element p.sub.j,k in P.sub.i represents the monitored pressure at the sensor i at the moment k in the j.sup.th day, the m.sup.th row in P.sub.1 represents the monitored pressure at the sensor i in a day before the moment to be detected, and for n sensors, total n matrices of monitored pressures are generated; (2) analyzing time-domain and frequency-domain features of the matrix of monitored pressures by Fourier transform, extracting disturbances, and generating pressure disturbances matrices; (3) identifying outliers in the pressure disturbances matrices by isolation forest algorithm; and (4) identifying detected outliers by calculating a qualitative index A and a quantitative index B, and outputting a result of pipe burst detection.

2. The method of claim 1, wherein the step (13) comprises: (131) considering a sensor with plenty of blank data as an invalid sensor, and abandoning the monitored pressure data at the sensor; (132) linearly interpolating the remained monitored pressure data at sensors and filling missing values; and (133) setting an upper threshold and a lower threshold to remove the monitored pressure data, which is abnormal in remained data.

3. A method for detecting pipe bursts in water distribution systems based on pressure disturbance extraction, the method comprising: (1) collecting and pre-processing monitored pressure data, and establishing pressure matrices; (2) analyzing time-domain and frequency-domain features of the matrix of monitored pressures by Fourier transform, extracting disturbances, and generating pressure disturbances matrices; (21) processing row vectors in pre-processed pressure matrices by Fourier transform, wherein, given that each row in pre-processed pressure matrices is a discrete signal x[n] having a length of N, then the Fourier transform X[k] is: $X\left[n\right]=\underset{n=0}{\overset{N-1}{.Math.}}x\left[n\right]e^{-j\frac{2\pi }{N}\mathrm{kn}},0\le k\le N-1,$ where X[k] is the Fourier transform of x[n], and $e^{-j\frac{2\pi }{N}\mathrm{kn}}$ is a periodic complex-exponential function of $\frac{2\pi }{N};$ (22) extracting high-frequency components (disturbances)in the monitored signal, the high-frequency component X′[k] in the frequency domain being extracted as follows: $X^{\prime }\left[k\right]=\{\begin{array}{cc}X\left[k\right],& \mu \le k\le N-\mu \\ 0,& \mathrm{otherwise}\end{array},0\le k\le N-1$ where, μ is a parameter which controls the number of selected high-frequency components; (23) reconstructing the extracted high-frequency component X′[k] in the frequency domain into a time-domain signal x′[n] by inverse Fourier transform: $x^{\prime }\left[n\right]=\frac{1}{N}\underset{k=0}{\overset{N-1}{.Math.}}{X}^{\prime }\left[k\right]e^{j\frac{2\pi }{N}\mathrm{kn}},0\le n\le N-1,$ (24) generating a matrix H of disturbed pressures, wherein, for the i.sup.th sensor, the corresponding disturbance matrix H.sub.i is: $H_i=\left(\begin{array}{cccccc}{h}_{x-m,y+1}& \mathrm{.Math.}& h_{x-m,N}& h_{x-m+1,1}& \mathrm{.Math.}& h_{x-m+1,y}\\ h_{x-m+1,y+1}& \mathrm{.Math.}& h_{x-m+1,N}& h_{x-m+2,1}& \mathrm{.Math.}& h_{x-m+2,y}\\ \mathrm{.Math.}& \ddots & \mathrm{.Math.}& \mathrm{.Math.}& \ddots & \mathrm{.Math.}\\ h_{x-1,y+1}& \mathrm{.Math.}& h_{x-1,N}& h_{x,1}& \mathrm{.Math.}& h_{x,y}\end{array}\right),$ where element h.sub.j,k in H.sub.i represents the high-frequency component extracted from the monitored pressure at the sensor i at the moment k in the j.sup.th day. For n sensors, total n matrices of disturbed pressures are generated; (3) identifying outliers in the pressure disturbances matrices by isolation forest algorithm; and (4) identifying detected outliers by calculating a qualitative index A and a quantitative index B, and outputting a result of pipe burst detection.

4. A method for detecting pipe bursts in water distribution systems based on pressure disturbance extraction, the method comprising: (1) collecting and pre-processing monitored pressure data, and establishing pressure matrices; (2) analyzing time-domain and frequency-domain features of the matrix of monitored pressures by Fourier transform, extracting disturbances, and generating pressure disturbances matrices; (3) identifying outliers in the pressure disturbances matrices by isolation forest algorithm; (31) establishing a detection matrix D using pre-processed pressure matrices: $D=\left(\begin{array}{ccc}{d}_{x-m+1,1}& \mathrm{.Math.}& d_{x-m+1,n}\\ \mathrm{.Math.}& \ddots & \mathrm{.Math.}\\ d_{x,1}& \mathrm{.Math.}& d_{x,n}\end{array}\right),$ where the i.sup.th column in the detection matrix D corresponds to the last column in the matrix H.sub.i of high-frequency components, each row in the detection matrix D is a row vector having a length of n, and the last row in the detection matrix D is the high-frequency components extracted at sensors at the moment to be detected; (32) inputting the detection matrix D into the isolation forest algorithm, each row in the detection matrix D being a detection sample, establishing t isolation trees, randomly selecting one dimension x as a study object for each isolation tree, 1≤x ≤n, randomly selecting a boundary value in this dimension to divide samples into two parts, and repeating this process until all samples are separated from others; (33) calculating an average path length of each sample in all isolation trees, wherein, for each sample in each isolation tree, the number of separations for this detection sample from other detection samples is the path length of this detection sample in this isolation tree, selecting k detection samples with minimal average path length as abnormal samples; (34) if the last row in the detection matrix D is selected as an abnormal sample, then identifying the moment to be detected as an abnormal state and passing the results to step (4), otherwise outputting a result indicating that no pipe burst occurs; and (4) identifying detected outliers by calculating a qualitative index A and a quantitative index B, and outputting a result of pipe burst detection.

5. The method of claim 4, wherein the step (4) comprises: (41) calculating the qualitative index A: $A=\{\begin{array}{c}0,{\mathrm{Sum}}_m\ge \stackrel{\_}{\mathrm{Sum}}\\ 1,{\mathrm{Sum}}_m<\stackrel{\_}{\mathrm{Sum}}\end{array}\mathrm{where},\stackrel{\_}{\mathrm{Sum}}=\frac{1}{m-1}\underset{i=1}{\overset{m-1}{.Math.}}{\mathrm{Sum}}_i,{\mathrm{Sum}}_i=\underset{j=1}{\overset{n}{.Math.}}D\left(i,j\right),i=1,2,L,m,$ D(i,j) represents the element in the i.sup.th row and the j.sup.th column in the detection matrix D, m represents the total number of rows in the detection matrix D, n represents the total number of columns in the detection matrix D, and Sum.sub.m represents the sum of the last row in the detection matrix D; (42) calculating the quantitative index B: $B=\{\begin{array}{c}0,{\mathrm{Sum}}_m\in \left(\stackrel{\_}{\mathrm{Sum}}-\mathrm{\lambda \delta },\stackrel{\_}{\mathrm{Sum}}+\mathrm{\lambda \delta }\right)\\ 1,\mathrm{otherwise}\end{array}$ where, λ is a parameter which controls the number of selected high-frequency components, and δ represents the standard deviation of Sum.sub.i=1, 2, . . . , m−1; and (43) outputting a result indicating that pipe burst occurs if both the qualitative index A and the quantitative index B are 1, or otherwise outputting a result indicating that no pipe burst occurs.

Description

BRIEF DESCRIPTION OF THE DRAWINGS

(1) FIG. 1 is a flowchart of the method for detecting pipe burst in water distribution systems based on pressure disturbance extraction.

(2) FIG. 2 is a graph of the monitored pressure data and the extracted high-frequency components at a certain sensor during a real-life pipe burst.

DETAILED DESCRIPTION OF THE PRESENT INVENTION

(3) The present invention will be described below in detail with reference to the accompanying drawings by specific embodiments. It is to be noted that the following description of the implementations is merely illustrative and the present invention is not intended to limit the applications or uses thereof. Furthermore, the present invention is not limited to the following implementations.

Embodiment

(4) As shown in FIG. 1, a method for detecting pipe bursts in water distribution systems based on pressure disturbance extraction is provided, comprising following steps:

(5) (1) collecting and pre-processing monitored pressure data, and establishing pressure matrices;

(6) (2) analyzing time-domain and frequency-domain features of the matrix of monitored pressures by Fourier transform, extracting disturbances, and generating pressure disturbances matrices;

(7) (3) identifying outliers in the pressure disturbances matrices by the isolation forest algorithm;

(8) (4) further identifying detected outliers by calculating and qualitative index A and quantitative index B and outputting a result of pipe burst detection.

(9) The step (1) can be implemented by following steps:

(10) (11) numbering the date and time of historical data. For example, 1 corresponds to January 1 and 2 corresponds to January 2; when the sampling frequency at a sensor is T min, 1 corresponds to 0:00 and 2 corresponds to 0: T;

(11) (12) for a moment t=(x,y) to be detected, extracting SCADA monitored pressure data in m days before the moment t, where (x,y) represents the moment y in the x.sup.th day;

(12) (13) pre-processing the extracted monitored pressure data;

(13) (14) dividing the pre-processed pressure data at a sensor into m row vectors, each vector having a length of N, N being the number of samples at one sensor in a day, and storing the row vectors in pre-processed pressure matrices P. For the i.sup.th sensor, the corresponding matrix P.sub.i is:

(14) $P_i=\left(\begin{array}{cccccc}{p}_{x-m,y+1}& L& p_{x-m,N}& p_{x-m+1,1}& L& p_{x-m+1,y}\\ p_{x-m+1,y+1}& L& p_{x-m+1,N}& p_{x-m+2,1}& L& p_{x-m+2,y}\\ M& O& M& M& O& M\\ p_{x-1,y+1}& L& p_{x-1,N}& p_{x,1}& L& p_{x,y}\end{array}\right),$
where the element p.sub.j,k in P.sub.i represents the monitored pressure at the sensor i at the moment k in the j.sup.th day, the m.sup.th row in P.sub.i represents the monitored pressure at the sensor i in a day before the moment to be detected. For n sensors, total n matrices of monitored pressures are generated.

(15) The step (13) can be implemented by following steps:

(16) (131) considering a sensor with plenty of blank data as an invalid sensor, and abandoning the monitored pressure data at the sensor;

(17) (132) linearly interpolating the remained monitored pressure data at sensors and filling missing values; and

(18) (133) setting an upper threshold and a lower threshold to remove the monitored pressure data which is abnormal in remained data.

(19) At the end of data pre-processing, time-domain and frequency-domain features of the monitored data are analyzed by Fourier transform. The disturbance extraction is implemented by associating high-frequency components in the monitored data with the abrupt pressure drop caused by pipe burst. Therefore, the step (2) can be implemented by following steps:

(20) (21) considering a row in a matrix of monitored pressures at a sensor as a discrete aperiodic signal, decomposing this signal into a series of linear combinations of complex-exponential signals by Fourier transform. Given that each row in pre-processed pressure matrices is a discrete signal x[n] having a length of N, then the Fourier transform X[k] is:

(21) $X\left[k\right]=\underset{n=0}{\overset{N-1}{.Math.}}x\left[n\right]e^{-j\frac{2\pi }{N}\mathrm{kn}},0\le k\le N-1,$
where, X[k] is the Fourier transform of x[n], and

(22) $e^{-j\frac{2\pi }{N}\mathrm{kn}}$
is a periodic complex-exponential function of

(23) $\frac{2\pi }{N};$

(24) (22) extracting high-frequency components (disturbances) in the monitored signal. The high-frequency component X′[k] in the frequency domain is extracted as follows:

(25) $X^{\prime }\left[k\right]=\{\begin{array}{c}X\left[k\right],\mu \le k\le N-\mu \\ 0,\text{otherwise}\end{array},0\le k\le N-1$
where, μ is a parameter for controlling the selection of the number of high-frequency terms;

(26) (23) reconstructing the extracted high-frequency component X′[k] in the frequency domain into a time-domain signal x′[n] by inverse Fourier transform:

(27) $x^{\prime }\left[n\right]=\frac{1}{N}\underset{k=0}{\overset{N-1}{.Math.}}{X}^{\prime }\left[k\right]e^{j\frac{2\pi }{N}\mathrm{kn}},0\le n\le N-1$

(28) (24) generating a matrix H of disturbed pressures, wherein, for the i.sup.th sensor, the corresponding matrix H.sub.i of disturbed pressures is:

(29) $H_i=\left(\begin{array}{cccccc}{h}_{x-m,y+1}& L& h_{x-m,N}& h_{x-m+1,1}& L& h_{x-m+1,y}\\ h_{x-m+1,y+1}& L& h_{x-m+1,N}& h_{x-m+2,1}& L& h_{x-m+2,y}\\ M& O& M& M& O& M\\ h_{x-1,y+1}& L& h_{x-1,N}& h_{x,1}& L& h_{x,y}\end{array}\right),$
where element h.sub.j,k in H.sub.i represents the high-frequency component extracted from the monitored pressure at the sensor i at the moment k in the j.sup.th day. For n sensors, total n matrices of disturbed pressures are generated.

(30) The step (3) can be implemented by following steps:

(31) (31) establishing a detection matrix D using pre-processed pressure matrices:

(32) $D=\left(\begin{array}{ccc}{d}_{x-m+1,1}& L& d_{x-m+1,n}\\ M& O& M\\ d_{x,1}& L& d_{x,n}\end{array}\right),$
where, the i.sup.th column in the detection matrix D corresponds to the last column in the matrix H.sub.i of high-frequency components, each row in the detection matrix D is a row vector having a length of n, and the last row in the detection matrix D is the high-frequency components extracted at sensors at the moment to be detected;

(33) (32) inputting the detection matrix D into the isolation forest algorithm with each row in the detection matrix D being a detection sample, establishing t isolation trees. randomly selecting one dimension x as a study object for each isolation tree, 1≤x≤n, randomly selecting a boundary value in this dimension to divide samples into two parts, and repeating this process until all samples are separated from others;

(34) (33) calculating an average path length of each sample in all isolation trees, wherein, for each sample in each isolation tree, the number of separations for this detection sample from other detection samples is the path length of this detection sample in this isolation tree, selecting k detection samples with minimal average path length as abnormal samples;

(35) (34) if the last row in the detection matrix D is selected as an abnormal sample, then identifying the moment to be detected as an abnormal state and passing the results to step (35), otherwise outputting a result indicating that no pipe burst occurs.

(36) To avoid the impact of pressure fluctuations caused by normal demand deviations in water distribution systems, the detected abnormal states should be further identified based on the system pressure response features during pipe bursts. That is, the step (4) is to further identify the detected abnormal state, which is implement by calculating qualitative index A and quantitative index B. The step (4) can be implemented by following steps:

(37) (41) calculating a qualitative index A:

(38) 0$A=\{\begin{array}{c}0,{\mathrm{Sum}}_m\ge \stackrel{\_}{\mathrm{Sum}}\\ 1,{\mathrm{Sum}}_m<\stackrel{\_}{\mathrm{Sum}}\end{array}\mathrm{where},\stackrel{\_}{\mathrm{Sum}}=\frac{1}{m-1}\underset{i=1}{\overset{m-1}{.Math.}}{\mathrm{Sum}}_i,{\mathrm{Sum}}_i=\underset{j=1}{\overset{n}{.Math.}}D\left(i,j\right),i=1,2,L,m,$
D(i,j) represents the element in the i.sup.th row and the j.sup.th column in the detection matrix D, m represents the total number of rows in the detection matrix D, n represents the total number of columns in the detection matrix D, and Sum.sub.m represents the sum of the last row in the detection matrix D;

(39) (42) calculating a quantitative index B:

(40) $B=\{\begin{array}{c}0,{\mathrm{Sum}}_m\in \left(\stackrel{\_}{\mathrm{Sum}}-\mathrm{\lambda \delta },\stackrel{\_}{\mathrm{Sum}}+\mathrm{\lambda \delta }\right)\\ 1,\mathrm{otherwise}\end{array}$
where, λ is a parameter which controls the number of selected high-frequency components, and δ represents the standard deviation of Sum.sub.i, i=1, 2, . . . , m−1; and

(41) (43) outputting a result indicating that pipe burst occurs if both qualitative index A and quantitative index B are 1, otherwise outputting a result indicating that no pipe burst occurs.

(42) To further demonstrate the implement of the method, a real-life pipe burst which occurred in a water distribution system at 3:15 on Mar. 28, 2016 is taken as an example.

(43) (1) Pressure data is collected and pre-processed, and pre-processed pressure matrices are established.

(44) There are total 23 sensors in the water distribution system, and the sampling interval at the sensors is 15 min. The monitored data at all sensors in 15 days before this moment is collected. The monitored values above 50 m and below 0 m are deleted. After the pre-processing, there are total 18 valid sensors. The pre-processed data is stored in pre-processed pressure matrices P. 18 pre-processed pressure matrices are generated. Since the sampling interval at the sensors is 15 min, 24 h/15 min=96 data are sampled every day at each sensor. Each matrix of monitored pressures has a dimension of 15 rows×96 columns.

(45) (2) Time-domain and frequency-domain features of each row in pre-processed pressure matrices are analyzed by Fourier transform, disturbances are extracted and disturbances matrices are generated.

(46) Each row vector (having a length of 96) in each matrix of monitored pressures is considered as one discrete signal, and the spectral composition of this signal is analyzed by Fourier transform. High-frequency components are extracted, where the value of the parameter μ is 10. The extracted high-frequency components are then converted to a time-domain signal by inverse Fourier transform. By taking one sensor as an example, the measured data and the extracted high-frequency value data at this sensor are shown in FIG. 2. The extracted disturbances are stored in H.

(47) (3) Outliers in the matrix of disturbed pressures are detected by using the isolation forest algorithm.

(48) The last columns in the 18 high-frequency components matrices are extracted to the detection matrix D. The detection matrix D is input to the isolation forest algorithm, and each row in the detection matrix D is considered as a sample. 100 isolation trees are generated. An average path length of each sample in all isolation trees is calculated, and k=1 abnormal sample are selected. The moment to be detected is selected as an abnormal sample because. It has minimal average path length. Therefore, the moment to be detected is determined as an abnormal state.

(49) (4) The detected outliers are further identified by qualitative index A and quantitative index B, and a result of pipe burst detection is output.

(50) The qualitative index A and the quantitative index B are calculated. Both of them are 1. Therefore, a result indicating that pipe burst occurs is output. There is a pipe burst repair record indicating that a pipe burst was reported by consumers at 9:55 in that day. Workers found that a pipe with a diameter of 300 mm failed after they arrived on 14:41. Therefore, the detection result is deemed as correct according to the pipe burst repair record. In this instance, the alarm generated by the present invention is approximately 7 hours earlier than the received consumer complaint. The time from the occurrence of pipe burst to the detection of pipe burst is successfully shortened.

(51) To verify the rationality of the detection result, the detection result is compared with the actual pipe burst repair record. This verifies the detection accuracy of this method.

(52) The above embodiments are merely illustrative and are not intended to limit the scope of the present invention. Those embodiments can be implemented in various other ways, and various omissions, substitutions and changes can be made without departing from the scope of the present invention.