Building concealed leakage intelligent monitoring method and system based on pipe network acoustic fingerprint features
By decomposing and analyzing the acoustic signature signal of the pipeline network, an acoustic signature state evolution path diagram and a local topological subgraph are constructed, which solves the problems of timeliness and positioning accuracy in the existing technology of leakage monitoring, and realizes highly sensitive leakage identification and accurate positioning.
Patent Information
- Application Number
- CN202610505777.6
- Authority / Receiving Office
- CN · China
- Patent Type
- Applications(China)
- Current Assignee / Owner
- Filing Date
- 2026-04-16
- Publication Date
- 2026-06-12
AI Technical Summary
Existing building water supply network leakage monitoring technologies have poor timeliness, making it difficult to detect minor leaks in hidden areas in the early stages. Furthermore, existing acoustic fingerprint monitoring schemes have a high false alarm rate and inaccurate leakage source location.
By acquiring the acoustic signature signal sequence collected by the acoustic signature sensing unit at each node of the pipeline network, it is decomposed into steady-state background acoustic signature components and transient impact acoustic signature components, and an acoustic signature state evolution path diagram is constructed. Combined with the local pipeline network topology sub-graph and fingerprint deviation vector field, the leakage source is determined.
It improves the sensitivity and accuracy of hidden leak detection, enabling timely detection of millimeter-level micro-leakage, reducing false alarm rates, and achieving early leak warning.
Smart Images

Figure CN122191468A_ABST
Abstract
Description
Technical Field
[0001] This invention relates to the field of artificial intelligence technology, and more specifically, to a method and system for intelligent monitoring of concealed building leaks based on the acoustic signature characteristics of pipe networks. Background Technology
[0002] Leaks in building water supply networks often occur in concealed areas such as ceilings and manholes. Existing monitoring technologies are mainly divided into two categories: manual inspection and flow / pressure monitoring. Manual inspection relies on personnel experience, has poor timeliness, and a high rate of missed detection, making it impossible to provide early warning of leaks. Flow / pressure monitoring has a high sensing threshold and extremely low sensitivity to millimeter-level leaks, making it difficult to detect early leaks in concealed locations in a timely manner. Some existing voiceprint-based monitoring schemes only extract single-dimensional voiceprint features, failing to fully explore the temporal evolution of voiceprints or combine them with the propagation characteristics of sound signals within the pipe network for source tracing. This results in a high false alarm rate and large deviations in leak source location, failing to meet the needs of monitoring concealed leaks in complex building water supply networks. Summary of the Invention
[0003] The purpose of this invention is to provide a method and system for intelligent monitoring of concealed building leaks based on the acoustic signature characteristics of pipe networks.
[0004] In a first aspect, embodiments of the present invention provide a method for intelligent monitoring of concealed building leaks based on the acoustic signature characteristics of a pipe network, comprising:
[0005] The original voiceprint signal sequence is acquired by the voiceprint sensing units deployed at each node of the water supply network inside the building during a continuous monitoring period. The original voiceprint signal sequence carries the node identifier and the acquisition timestamp identifier corresponding to each voiceprint sensing unit.
[0006] Each of the original acoustic signature signal sequences is decomposed to obtain its corresponding steady-state background acoustic signature component sequence and transient impact acoustic signature component sequence;
[0007] The transient impact acoustic component sequence is grouped according to the node identifier, and the transient impact event time axis of each node is obtained by arranging it in time sequence according to the acquisition timestamp identifier. An event interval fluctuation sequence is constructed based on the time interval between adjacent events.
[0008] The event interval fluctuation sequence is correlated and mapped with the spectral centroid offset trajectory of the steady-state background voiceprint component sequence to generate a voiceprint state evolution path diagram for each node.
[0009] Based on the state trajectory of the voiceprint state evolution path diagram, the voiceprint dynamic behavior fingerprint of each node is constructed, and the voiceprint dynamic behavior fingerprint of each node under current and historical normal working conditions is compared. Abnormal nodes are determined based on the divergence distribution of the fingerprint deviation vector field.
[0010] Based on the abnormal nodes, a local pipeline network topology subgraph is constructed. The direction of abnormal acoustic fingerprint propagation is determined according to the curl direction of the deviation vector field of the fingerprints of adjacent nodes in the local pipeline network topology subgraph. The set of candidate nodes for leakage source is obtained by backtracking the abnormal nodes and the direction of abnormal acoustic fingerprint propagation.
[0011] The node with the largest absolute value of the fingerprint deviation vector field divergence in the candidate node set of the leakage source is selected as the leakage source localization result output.
[0012] In a second aspect, embodiments of the present invention provide a server system, including a server, the server being used to execute the method described in the first aspect.
[0013] Compared to existing technologies, the beneficial effects of this invention include: The method and system for intelligent monitoring of concealed leaks in buildings based on pipeline acoustic signature features, as disclosed in this invention, acquires the original acoustic signature signal sequences carrying node identifiers and timestamps collected by acoustic signature sensing units at each node of the pipeline network. These sequences are decomposed into steady-state background acoustic signature components and transient impact acoustic signature components. The transient impact event timelines for each node are grouped and sorted to construct event interval fluctuation sequences. The spectral centroid offset trajectories of the event interval fluctuation sequences and steady-state background acoustic signature components are correlated to generate acoustic signature state evolution path diagrams for each node. This leads to the construction of acoustic signature dynamic behavior fingerprints. By comparing these fingerprints with those of normal operating conditions, abnormal nodes are identified based on deviation vector field divergence. The direction of abnormal propagation is determined by combining the curl of the vector field of adjacent nodes in the local pipeline network topology. A candidate set of leak sources is obtained by backtracking, and the node with the largest absolute divergence value is selected as the leak source output. This method can improve the sensitivity and accuracy of concealed leak identification. Attached Figure Description
[0014] To more clearly illustrate the technical solutions of the embodiments of the present invention, the accompanying drawings used in the embodiments will be briefly described below. It should be understood that the following drawings only show some embodiments of the present invention and should not be considered as limiting the scope. For those skilled in the art, other related drawings can be obtained based on these drawings without creative effort.
[0015] Figure 1 A flowchart illustrating the steps of the intelligent monitoring method for concealed building leaks based on pipeline acoustic signatures provided in an embodiment of the present invention.
[0016] Figure 2 A schematic block diagram of the structure of a computer device provided in an embodiment of the present invention. Detailed Implementation
[0017] To make the objectives, technical solutions, and advantages of the embodiments of the present invention clearer, the technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some, not all, of the embodiments of the present invention. The components of the embodiments of the present invention described and shown in the accompanying drawings can generally be arranged and designed in various different configurations.
[0018] The specific embodiments of the present invention will now be described in detail with reference to the accompanying drawings.
[0019] In order to solve the technical problems mentioned in the background art Figure 1 This is a flowchart illustrating the intelligent monitoring method for concealed building leaks based on pipeline acoustic signatures provided in this embodiment. The following is a detailed description of the intelligent monitoring method for concealed building leaks based on pipeline acoustic signatures.
[0020] Step S201: Obtain the original voiceprint signal sequence collected by the voiceprint sensing units deployed at each node of the water supply network inside the building during the continuous monitoring period. The original voiceprint signal sequence carries the node identifier and collection timestamp identifier corresponding to each voiceprint sensing unit.
[0021] The server periodically receives monitoring data from acoustic signature sensors deployed at key nodes (such as pipe joints, valves, and branch pipe connections) within the building's water supply network via its communication interface. Each data packet contains a digitized sequence of acoustic signature signals collected within a fixed monitoring period, along with a unique node identifier (e.g., device ID "SN-2024-001") and a precise timestamp identifier (e.g., "2024-05-10 14:30:00.000"). The server organizes these raw data streams according to time windows (e.g., one data packet every 5 minutes) and stores them in a time-series database, categorized by node identifier, forming the raw acoustic signature signal sequence corresponding to each node. For example, the server records the sound pressure signal sequence containing multiple sampling points collected by node "SN-2024-001" between 8:00 AM and 8:05 AM.
[0022] Step S202: Decompose each of the original voiceprint signal sequences to obtain their respective steady-state background voiceprint component sequences and transient impact voiceprint component sequences;
[0023] The server invokes the signal processing module to perform adaptive mode decomposition on the raw acoustic signature signal sequence of each node. Specifically, the server inputs the raw signal sequence of node "SN-2024-001" into a pre-defined separation network. This network iteratively filters the complex raw signal based on an empirical mode decomposition algorithm, decomposing it into eight intrinsic mode function components. Each component represents an oscillation mode at a different time scale. For example, the first component might correspond to high-frequency noise generated by water flow rubbing against the pipe wall, while the last component might correspond to low-frequency vibration caused by the operation of a water pump. Next, the server calculates the variance of the instantaneous frequency of each component over the entire time period. It sets a frequency variance threshold (e.g., 0.1 Hz²). The server categorizes three low-frequency stable components with variances below this threshold into a steady-state background acoustic signature candidate component set, and five high-frequency, bursty components with variances above or equal to this threshold into a transient impact acoustic signature candidate component set. Subsequently, the server linearly superimposes the three components in the steady-state candidate set to reconstruct the steady-state background acoustic signature component sequence of the node, which exhibits a continuous signal with relatively stable amplitude and frequency. Simultaneously, the server superimposes the five components from the transient candidate set to obtain an initial transient impact acoustic component sequence. To filter out minor environmental interference, the server sets an amplitude filtering threshold (e.g., 3 times the root mean square value of the background noise), retaining only waveform pulse segments in the initial sequence whose amplitude exceeds this threshold, ultimately generating a pure transient impact acoustic component sequence, where each pulse segment may correspond to a water hammer effect, valve opening and closing, or potential leakage impact.
[0024] Step S203: Group the transient impact acoustic component sequence according to the node identifier, arrange the transient impact events of each node in time sequence according to the acquisition timestamp identifier, and construct the event interval fluctuation sequence based on the time interval between adjacent events.
[0025] The server initiates the event analysis process. It first extracts all transient impact acoustic component sequences belonging to node "SN-2024-001" from the cache and sorts them in ascending order based on the acquisition timestamp attached to each sequence, forming a coherent transient impact event sequence. The server iterates through this sequence, identifying the starting sampling point of each impact waveform and converting it into an absolute time point (e.g., "14:30:01.235", "14:30:03.112", etc.), thereby generating the transient impact event timeline for that node. Next, the server calculates the difference between the occurrence times of adjacent events on the timeline, obtaining the original event interval sequence (e.g., [1.877s, 0.945s, 2.301s, ...]). To prevent individual abnormal intervals (such as false triggers caused by sudden strong external noise) from interfering with statistical characteristics, the server performs outlier correction: it calculates the difference between each interval value and the mean of its immediate and adjacent interval values. If the absolute value of this difference exceeds a preset deviation threshold (e.g., 50% of the average interval), the outlier value is replaced with the mean of the adjacent intervals. After correction, the server uses a moving average filter with a window size of 3 to smooth the interval sequence to suppress random fluctuations. Finally, the server takes the natural logarithm of each value in the smoothed event interval sequence to generate an event interval fluctuation sequence. The numerical changes in this sequence better conform to the characteristics of a steady-state stochastic process, facilitating subsequent correlation analysis with steady-state signal characteristics.
[0026] Step S204: Associate and map the spectral centroid offset trajectory of the event interval fluctuation sequence with the steady-state background voiceprint component sequence to generate a voiceprint state evolution path diagram for each node.
[0027] The server processes the steady-state background acoustic signature component sequence in parallel. It divides the steady-state component sequence of node "SN-2024-001" into frames of 256 sampling points each, with a frame shift of 128 sampling points, resulting in hundreds of continuous steady-state background acoustic signature analysis frames. For each frame, the server performs a Fast Fourier Transform to obtain its power spectral density. The server calculates the centroid frequency of the spectrum for each frame, i.e., the weighted average frequency of the spectral energy distribution. Arranging the centroid frequencies of all analysis frames in chronological order yields the original trajectory of the spectral centroid. The server calculates the first-order difference of this trajectory to obtain the change in centroid frequency (Δf) between adjacent frames, and generates an offset direction sequence based on the sign of Δf (e.g., +1 indicates a frequency increase, -1 indicates a frequency decrease). By accumulating and summing this direction sequence, the server obtains the cumulative trajectory of the spectral centroid offset, which reflects the overall shift trend of the background acoustic signature spectral centroid over time, possibly related to the slow changes in pipeline pressure. To correlate features across different time scales, the server synchronizes the event interval fluctuation sequence with the cumulative spectral centroid shift trajectory. Using the same analysis time window (e.g., every 10 seconds) as a benchmark, it extracts the average value of the event interval fluctuation sequence within the same window (representing the density of impact events) and the final value of the cumulative spectral centroid shift (representing the cumulative shift of the background spectrum). The server plots the data points (X, Y) corresponding to each time window on a two-dimensional plane, using the cumulative spectral centroid shift as the x-axis and the event interval fluctuation value as the y-axis. These points are then connected by line segments in chronological order, ultimately generating the acoustic signature evolution path diagram for node "SN-2024-001". This diagram vividly illustrates the trajectory of the node's acoustic signature state within the phase plane formed by "background spectrum shift" and "event frequency".
[0028] Step S205: Based on the state trajectory of the voiceprint state evolution path diagram, construct the voiceprint dynamic behavior fingerprint of each node, compare the current and historical normal working conditions of each node with the voiceprint dynamic behavior fingerprint, and determine the abnormal node based on the divergence distribution of the fingerprint deviation vector field.
[0029] The server performs in-depth dynamic analysis on the generated voiceprint state evolution path map to construct a behavioral fingerprint. It first reconstructs the phase space of the two-dimensional coordinate point sequence on the path map. The server automatically determines the time delay parameter τ=5 and the embedding dimension m=3 based on the autocorrelation function method and the spurious nearest neighbor method. Through reconstruction, the two-dimensional trajectory is embedded into a three-dimensional phase space, resulting in a series of reconstructed trajectory points. The server calculates the Euclidean distance between all reconstructed points and their nearest neighbors, taking the upper quartile of their statistical distribution as the global neighborhood radius R. For each reconstructed point in the phase space, the server finds all neighboring points whose Euclidean distance is less than R, forming a local neighborhood point set. The server analyzes the change in the number of points in this set during subsequent evolution time steps, fits its growth curve, and extracts the exponential growth rate as the orbit divergence rate parameter λ for that point. Simultaneously, the server scans the entire phase space, identifies regions repeatedly visited by trajectory points, marks regions visited more than a threshold (e.g., 10 times) as orbit convergence regions, and calculates the geometric center and equivalent radius of these regions. The server combines the divergence rate λ and the convergence region radius r of each node into the original feature vector [λ, r]. After minimax normalization, the server uses the t-SNE nonlinear dimensionality reduction algorithm to project the high-dimensional feature vector onto a two-dimensional feature plane, obtaining a feature projection point. The server calculates the distance from this projection point to the origin of the plane as the dynamic feature intensity parameter S, and calculates the angle between this point and the horizontal axis (first principal axis) of the plane as the dynamic feature direction parameter θ. Finally, the server defines (S, θ) as the acoustic dynamic behavior fingerprint of this node in the current monitoring period and stores it in association with the node ID.
[0030] When identifying anomalous nodes, the server retrieves the behavioral fingerprint set calculated for node "SN-2024-001" under normal operating conditions over the past 30 days from the historical database, and calculates the historical mean of its intensity parameter S_hist and the historical mean of its orientation parameter θ_hist. The server calculates the intensity deviation ΔS = S_curr – S_hist and the orientation deviation Δθ = θ_curr – θ_hist for the fingerprint of the current monitoring period. These two values constitute a two-dimensional deviation vector. The server performs the same calculation for all 50 monitoring nodes within the building, obtaining 50 deviation vectors. Based on the actual (x, y) coordinates of these 50 nodes on the building floor plan, the server constructs a spatial vector field F(x, y) = (ΔS, Δθ). Subsequently, the server calculates the divergence at the coordinates of each node in this vector field. For example, at node "SN-2024-015", the server calculates the divergence value div_F = +0.85 using the central difference method based on its deviation vector from its surrounding neighboring nodes. This positive value indicates that at this node, the fingerprint deviation exhibits an outward divergence characteristic, suggesting that this may be the starting point of an abnormal state or a point significantly affected. The server marks all nodes with absolute divergence values exceeding a threshold (e.g., 0.1) as abnormal nodes, forming a set of abnormal nodes. Positive divergence values are marked as "divergent anomalies," and negative values are marked as "convergent anomalies."
[0031] Step S206: Construct a local pipeline network topology subgraph based on the abnormal nodes, determine the abnormal propagation direction of the acoustic fingerprint according to the curl direction of the deviation vector field of the fingerprint of adjacent nodes in the local pipeline network topology subgraph, and obtain a set of candidate nodes for leakage source by combining the abnormal nodes and the abnormal propagation direction of the acoustic fingerprint.
[0032] The server loads a digital topology map of the building's water supply network, which defines the connections between all sensor nodes using a graph structure. The server selects the node "SN-2024-015" with the largest divergence value from the set of anomalous nodes as the center, and performs a breadth-first search in the global topology map to extract all nodes within 3 hops of it (e.g., including its upstream nodes "SN-2024-010" and "SN-2024-008" and its downstream node "SN-2024-018"), constructing a local sub-graph of the network topology.
[0033] Within the subgraph, the server obtains the two-dimensional deviation vector for each node. For a connecting edge in the subgraph, such as a pipe connecting node A (“SN-2024-010”) and node B (“SN-2024-015”), the server calculates the difference between the deviation vectors of node B and node A, and performs a cross product operation on this difference vector with the pipe direction vector from A to B to estimate the curl component on that edge. By combining the calculations of all edges in the subgraph, the server determines that the curl value at node “SN-2024-015” is negative, showing a clockwise rotation trend, while the curl value of its upstream node “SN-2024-010” is positive. According to preset rules, when adjacent nodes have opposite curl signs, combined with the pipe direction and rotation direction, the server infers that the propagation direction of the acoustic fingerprint anomaly (i.e., a specific dynamic fingerprint deviation pattern) is from node “SN-2024-010” to node “SN-2024-015”. The server performs similar calculations on the connection relationships between all abnormal nodes, generating a mapping table of abnormal voiceprint propagation directions.
[0034] Based on this mapping table, the server initiates a reverse backtracking process. Starting with all anomalous nodes, it reverses the propagation direction and traces upstream. For example, it traces back from "SN-2024-015" to find "SN-2024-010," and then from "SN-2024-010" it might trace back to "SN-2024-005." The server records all tracing paths, forming multiple backtracking trees with the initial anomalous node as the leaf node. It merges nodes on all backtracking paths and counts the frequency of each node appearing in different backtracking paths. The server sorts candidate nodes based on their depth in the backtracking tree (hops from the leaf node), backtracking frequency, and the number of their direct downstream branches (branch factor). Nodes with shallower depth (closer to upstream), higher backtracking frequency, and smaller branch factors are given higher priority. Finally, the server generates an ordered set of candidate nodes for the leakage source, such as {"SN-2024-005", "SN-2024-003", "SN-2024-008"}.
[0035] Step S207: Select the node with the largest absolute value of fingerprint deviation vector field divergence in the candidate node set of leakage sources as the leakage source localization result output.
[0036] The server reads the absolute divergence value of the fingerprint deviation vector field calculated in step five for each candidate node from the set of candidate nodes for leakage sources. For example, the absolute divergence value of node "SN-2024-005" is 1.2, that of node "SN-2024-003" is 0.9, and that of node "SN-2024-008" is 0.7. The server sorts these values in descending order and determines "SN-2024-005" as the primary locating node and "SN-2024-003" as the secondary locating node.
[0037] The server performs cross-validation. It checks that the depth of the primary locating node "SN-2024-005" in the backtracking tree is 2, which is less than the preset depth threshold of 3, indicating that its location is reasonable. Next, the server calculates the shortest pipe path distance between "SN-2024-005" and "SN-2024-003" in the pipeline topology map, which is 4 hops, greater than the preset distance threshold of 2, ruling out the possibility that they are closely related sources of the same leakage event. Subsequently, the server calculates the curl field of the local area where the primary locating node "SN-2024-005" is located, confirming that it is itself the node with the most significant curl value, corroborating the divergence analysis conclusions.
[0038] Therefore, the server ultimately determined that node "SN-2024-005" was the most likely source of the leak. The server formatted a leak source location report, including the node identifier, geographical location description, anomaly intensity score, and confidence level, and sent it to the building management platform through its output interface, triggering a corresponding alarm notification. Simultaneously, the server archived and stored all intermediate data, feature fingerprints, and analysis results from this monitoring period for updating the historical database and optimizing model parameters.
[0039] In this embodiment of the invention, the step of performing a voiceprint signal structure decomposition operation on each original voiceprint signal sequence to decompose the steady-state background voiceprint component sequence and the transient impact voiceprint component sequence in the original voiceprint signal sequence can be implemented through the following example.
[0040] The original voiceprint signal sequence is input into a separation network based on adaptive mode decomposition. The separation network decomposes the original voiceprint signal sequence into multiple intrinsic mode function components through an iterative filtering process. Each intrinsic mode function component corresponds to a voiceprint oscillation mode at a different time scale.
[0041] Calculate the instantaneous frequency variance of each intrinsic mode function component. Intrinsic mode function components with instantaneous frequency variance values lower than a preset frequency variance threshold are classified into a steady-state background soundprint candidate component set, and intrinsic mode function components with instantaneous frequency variance values not lower than a preset frequency variance threshold are classified into a transient impact soundprint candidate component set.
[0042] Linear superposition and reconstruction are performed on each intrinsic mode function component in the steady-state background voiceprint candidate component set to generate a steady-state background voiceprint component sequence. The instantaneous frequency of the steady-state background voiceprint component sequence remains stable over time.
[0043] Linear superposition and reconstruction are performed on each intrinsic mode function component in the transient impact acoustic fingerprint candidate component set to generate an initial transient impact acoustic fingerprint component sequence. The initial transient impact acoustic fingerprint component sequence is then filtered by amplitude threshold, and transient impact waveform segments with amplitudes exceeding multiples of the background noise amplitude are retained to obtain the transient impact acoustic fingerprint component sequence.
[0044] Extract the waveform envelope of each transient impact waveform segment in the transient impact acoustic component sequence, calculate the time required for the waveform envelope to rise from the starting point to the peak point as the rising edge duration parameter, and calculate the time required for the waveform envelope to fall from the peak point to the ending point as the falling edge duration parameter.
[0045] The waveform morphology category of the transient impact waveform segment is determined based on the ratio of the rising edge duration parameter to the falling edge duration parameter, and transient impact waveform segments of the same waveform morphology category are combined in chronological order to generate the transient impact subsequence corresponding to that waveform morphology category.
[0046] Statistical analysis of the inter-pulse time interval is performed on the transient impact subsequences corresponding to each waveform morphology category to generate the pulse interval probability density distribution function for each waveform morphology category, and the pulse interval probability density distribution function is used as the temporal structure feature of the transient impact subsequence corresponding to that waveform morphology category.
[0047] The signal components remaining after subtracting the steady-state background acoustic component sequence from the original acoustic signal sequence are compared with the transient impact acoustic component sequence. The residual signal energy value between the two is calculated. When the residual signal energy value is lower than the energy threshold, the separation operation is confirmed to be complete.
[0048] After confirming the completion of the separation operation, the steady-state background acoustic component sequence and the transient impact acoustic component sequence are stored in their respective data buffers, and a timestamp index and a node identifier index are added to each data buffer.
[0049] A synchronization association index table is established between the steady-state background acoustic component sequence and the transient impact acoustic component sequence. The synchronization association index table records the correspondence between the spectral characteristics of the steady-state background acoustic component sequence and the pulse characteristics of the transient impact acoustic component sequence within the same time window.
[0050] In this embodiment of the invention, for example, when the server performs the acoustic signature signal structure decomposition operation, it first calls its built-in adaptive mode decomposition algorithm module. This module receives a 10-second raw acoustic signature signal sequence with a sampling rate of 16kHz from node "SN-2024-001". The server starts the decomposition process, inputting the raw signal sequence into a separation network based on an improved empirical mode decomposition algorithm. Through an iterative filtering process, the network automatically and adaptively decomposes the non-stationary raw signal into eight intrinsic mode function components. These components are arranged from high to low frequency. For example, the IMF1 component may mainly contain high-frequency hissing sounds generated by the friction between water flow and pipe wall, the IMF2 and IMF3 components may contain mid-frequency equipment vibration harmonics, and the IMF7 and IMF8 components may correspond to low-frequency fundamental waves and their harmonics caused by the operation of the water pump. Each component represents a physical oscillation mode at a different time scale.
[0051] Subsequently, the server performs a Hilbert transform on each decomposed IMF component, calculates its instantaneous frequency sequence, and further calculates the variance of this instantaneous frequency sequence over the entire 10-second time window. The server presets a frequency variance threshold of 0.15 Hz². After calculation, the server determines that the instantaneous frequency variances of IMF6, IMF7, and IMF8 are 0.08 Hz², 0.05 Hz², and 0.03 Hz², respectively, all below the threshold, and therefore classifies them as the steady-state background soundprint candidate component set. The instantaneous frequency variances of IMF1 to IMF5 are between 0.2 Hz² and 1.5 Hz², not below the threshold, and therefore are classified as the transient impact soundprint candidate component set.
[0052] Next, the server linearly superimposes three IMF components (IMF6, IMF7, IMF8) from the steady-state background acoustic fingerprint candidate component set to reconstruct a new time series, namely the steady-state background acoustic fingerprint component sequence. This sequence exhibits a stable signal with slowly changing amplitude and a small fluctuation in instantaneous frequency around a central value (e.g., 50Hz), reflecting the continuous background acoustic environment maintained by the water pump and system pressure at this node. Simultaneously, the server linearly superimposes five IMF components (IMF1 to IMF5) from the transient impact acoustic fingerprint candidate component set to generate the initial transient impact acoustic fingerprint component sequence. This sequence contains many pulses with varying amplitudes. To focus on significant impact events, the server first calculates the root mean square of the amplitude of the steady-state background acoustic fingerprint component sequence as a baseline for background noise amplitude. The server sets an amplitude filtering threshold to three times this baseline value and scans the initial transient impact acoustic fingerprint component sequence, retaining only complete waveform segments with peak values exceeding this threshold, ultimately obtaining a clean transient impact acoustic fingerprint component sequence, where each segment may represent an independent event.
[0053] The server further analyzes each retained transient impact waveform segment. It extracts the waveform envelope of the segment, accurately identifying the start, peak, and end points of the waveform. For example, for a typical pulse, the server calculates the time it takes for its envelope to rise from the start point to the peak point as 2.5 milliseconds, recording this as the rise-edge duration parameter; it calculates the time it takes for the envelope to fall from the peak point to the end point as 8 milliseconds, recording this as the fall-edge duration parameter. The server determines the waveform morphology category based on the ratio of these two parameters (2.5 / 8 = 0.3125 in this example). For example, a ratio less than 0.5 is classified as "fast rise, slow fall," a ratio between 0.5 and 2.0 is classified as "symmetrical," and a ratio greater than 2.0 is classified as "slow rise, fast fall." The server extracts and combines all transient impact waveform segments belonging to the "fast rise, slow fall" type according to their occurrence time sequence, generating a "fast rise, slow fall" transient impact subsequence. The same operation is performed on other identified morphology categories to generate multiple morphology subsequences.
[0054] Next, the server performs temporal structure analysis on each morphological subsequence. Taking the "fast rise, slow fall" subsequence as an example, the server calculates the time interval between the start points of all adjacent pulses in the subsequence, forming a set of interval data. The server performs kernel density estimation on this set of data to generate a pulse interval probability density distribution function describing the distribution pattern of pulse intervals. This function reveals the statistical temporal pattern of such events, for example, exhibiting a bimodal distribution with two peaks around 0.5 seconds and 2.0 seconds. This distribution function is stored as the temporal structure feature of this morphological category.
[0055] To verify the completeness of the decomposition, the server performs residual verification. It subtracts the generated steady-state background acoustic component sequence from the original acoustic signal sequence to obtain an intermediate residual signal. The server then calculates the difference between this intermediate residual signal and the final transient impact acoustic component sequence to obtain the final residual signal, and calculates the energy value of this final residual signal. The calculated energy value is only 0.8% of the total energy of the original signal, lower than the preset 1.5% energy threshold. Therefore, the server confirms that the signal separation operation has been effectively completed.
[0056] After confirmation, the server stores the steady-state background acoustic signature component sequence and the transient impact acoustic signature component sequence (including their various morphological subsequences and corresponding temporal features) into different data buffers. Each buffer contains data with a precise acquisition timestamp index (e.g., "2024-05-10 14:30:00-14:30:10") and a node identifier index ("SN-2024-001"), facilitating subsequent retrieval by time and spatial dimensions. Finally, the server establishes a synchronous association index table. This table uses time windows (e.g., 100 milliseconds per window) as index units, recording the instantaneous average frequency, spectral centroid, and other characteristics of the steady-state background acoustic signature component sequence within each time window, as well as whether a transient impact pulse appears within that window, the pulse's morphology, amplitude, and other characteristics. This establishes a precise temporal correspondence between the two types of components, laying the data foundation for subsequent association mapping analysis.
[0057] In this embodiment of the invention, the step of arranging all transient impact acoustic component sequences corresponding to the same node identifier in chronological order according to the acquisition timestamp identifier to generate a transient impact event timeline for each node, and constructing an event interval fluctuation sequence based on the time interval between adjacent transient impact events on the transient impact event timeline, can be implemented through the following example.
[0058] Extract all transient impact acoustic component sequences corresponding to the same node identifier from the data buffer, sort them according to the order of the acquisition timestamp identifiers associated with each transient impact acoustic component sequence, and generate a transient impact event sequence arranged in ascending order of time.
[0059] Traverse the sequence of transient impact events, extract the start time of each transient impact event as the event occurrence time point, and arrange all event occurrence time points in chronological order to form a transient impact event time axis;
[0060] Calculate the time difference between the occurrence times of adjacent events on the time axis of the transient impact event, generate the original sequence of adjacent event intervals, and perform outlier correction on the original sequence of adjacent event intervals, removing event interval values whose deviation from the intervals of adjacent events exceeds the deviation threshold, to obtain the corrected event interval sequence.
[0061] The corrected event interval sequence is smoothed by a sliding window. The sliding window is formed by selecting the adjacent event interval values before and after the current event interval value as the center. The weighted average of all event interval values in the sliding window is calculated as the smoothed event interval value, and the event interval smoothed sequence is generated.
[0062] Calculate the natural logarithmic transformation value of each event interval value in the event interval smoothing sequence, and take the event interval sequence after natural logarithmic transformation as the event interval fluctuation sequence. The numerical distribution of the event interval fluctuation sequence exhibits the statistical characteristics of a steady-state random process.
[0063] Extract the mean and variance parameters of the event interval fluctuation sequence, construct a theoretical probability distribution model of the event interval fluctuation sequence based on the mean and variance parameters, and calculate the fitting residual sequence between each event interval value in the event interval fluctuation sequence and the theoretical probability distribution model.
[0064] The autocorrelation function is calculated on the fitted residual sequence, the autocorrelation coefficient at the delay order is extracted, and a long-range correlation descriptor for the event interval fluctuation sequence is constructed based on the autocorrelation coefficient.
[0065] The event interval fluctuation sequence and its long-range correlation descriptor are associated and stored with the node identifier of the corresponding node to generate a node event interval feature record.
[0066] Perform the above operations on all nodes deployed within the same building to generate a set of node event interval feature records corresponding to each node;
[0067] Establish a topological association mapping between node event interval feature records, pair and store the node event interval feature records of adjacent nodes according to the pipeline connection relationship, and generate a set of event interval association pairs between nodes.
[0068] In this embodiment of the invention, exemplarily, when the server performs transient impact event timeline construction and event interval analysis, it first extracts all transient impact acoustic component sequences associated with the node identifier "SN-2024-001" from its data cache. These sequences are collected and stored within a continuous 5-minute monitoring period, and each sequence is accompanied by a precise collection timestamp identifier, such as "2024-05-10 14:30:05.200", "14:30:10.150", "14:30:15.100", etc. The server sorts all sequences according to the chronological order of these timestamps, generating a transient impact event sequence arranged in ascending order of time. This sequence logically represents an ordered set of all discrete impact events that occurred at the node within a certain period of time.
[0069] The server then iterates through this ordered sequence of events. For each transient impact acoustic component sequence (i.e., an event record) in the sequence, the server analyzes its waveform to accurately identify the starting sampling point of the impact waveform and converts it into an absolute event occurrence time point. For example, the starting time of the first event is determined to be "14:30:05.235", the second is "14:30:06.112", the third is "14:30:08.413", and so on. The server arranges all the identified event occurrence times in chronological order to form a precise transient impact event timeline. This timeline is essentially a list of time points, clearly recording the timing of the impact events at the "SN-2024-001" node.
[0070] Next, the server calculates the time difference between every two adjacent events on the timeline. For example, the first interval Δt1 = 6.112 - 5.235 = 0.877 seconds, and the second interval Δt2 = 8.413 - 6.112 = 2.301 seconds. After calculating all adjacent intervals, an original sequence of adjacent event intervals is generated, such as [0.877s, 2.301s, 1.502s, 0.956s, 15.221s, 1.023s, ...]. The server immediately performs outlier correction on this original sequence to eliminate occasional interference. It sets the deviation threshold to 100% of the average of adjacent intervals. When the interval value of 15.221 seconds is checked, the server finds that its preceding interval is 0.956 seconds and its following interval is 1.023 seconds, which deviates from the average of the preceding and following intervals by a significant margin. Therefore, the server replaces this outlier with the average of the preceding and following intervals (approximately 0.99 seconds). The corrected sequence of events is more reasonable: [0.877s, 2.301s, 1.502s, 0.956s, 0.990s, 1.023s,...].
[0071] To further smooth random fluctuations and highlight trends, the server performs sliding window smoothing on the corrected sequence. It sets the window size to 3 (containing the current value, the previous value, and the next value) and assigns a weight of 0.5 to the center value and 0.25 to each of the values on either side. For example, for an interval of 1.502 seconds (located in the sequence), its smoothed value is calculated as: 0.25 * 2.301 + 0.5 * 1.502 + 0.25 * 0.956 ≈ 1.565 seconds. The server performs this operation on the entire sequence, generating a smoothed event interval sequence.
[0072] The server then performs a natural logarithmic transformation on each interval value in the smoothed sequence. For example, 1.565 seconds is transformed into ln(1.565)≈0.448. The resulting event interval fluctuation sequence has a more stable numerical distribution, more closely resembling a steady-state random process, which facilitates subsequent statistical modeling and analysis.
[0073] The server further extracts the statistical characteristics of the fluctuation sequence. It calculates the mean μ of the sequence as 0.35 (corresponding to an original interval of approximately 1.42 seconds) and the variance σ² as 0.08. Based on these two parameters, the server constructs a theoretical log-normal distribution model as a reference. Next, the server calculates the difference between each actual value in the fluctuation sequence and the predicted value of the theoretical model, forming a fitted residual sequence. For example, at a certain time point, if the actual fluctuation value is 0.448 and the theoretical model predicts a value of 0.400, then the residual is 0.048.
[0074] To investigate whether there is a temporal dependency in the event intervals, the server calculated the autocorrelation function of the fitted residual sequence. The calculation revealed that the autocorrelation coefficients were 0.15 and 0.05 at delays of 1st and 2nd order, respectively, but rapidly decayed to near zero after higher delays. This indicates a short-term positive correlation between event intervals, meaning that a longer interval is more likely to be followed by another longer interval (or vice versa). The server constructed a vector [0.15, 0.05, 0.01, 0.00, -0.02] from the first five autocorrelation coefficients and used it as a long-range correlation descriptor for the event interval fluctuation sequence.
[0075] Finally, the server packages the event interval fluctuation sequence, its statistical parameters (μ, σ²), the fitted residual sequence, and the long-range correlation descriptor of node "SN-2024-001", associates them with the node identifier, generates a complete node event interval feature record, and stores it in the feature database. The server repeats the above entire set of operations for all 50 monitoring nodes in the building, generating a set containing 50 records.
[0076] Based on this, the server establishes relationships between records according to a pre-stored building water supply network topology map. For example, nodes "SN-2024-001" and "SN-2024-010" are directly connected on the physical pipeline. The server extracts the event interval feature records of these two nodes, pairs them to form an inter-node event interval relationship pair, and stores it. This relationship pair will be used for subsequent analysis of the coordinated changes or propagation delays of event patterns of adjacent nodes, providing topological-level relationship evidence for locating the leakage source.
[0077] In this embodiment of the invention, the step of associating and mapping the spectral centroid offset trajectory of the event interval fluctuation sequence with the steady-state background voiceprint component sequence to generate a voiceprint state evolution path diagram for each node can be implemented through the following example.
[0078] The steady-state background acoustic component sequence is segmented into frames to obtain multiple consecutive steady-state background acoustic analysis frames. A Fourier transform is then performed on each steady-state background acoustic analysis frame to generate the steady-state spectrum distribution of each steady-state background acoustic analysis frame.
[0079] Calculate the centroid frequency of each steady-state spectral distribution, where the centroid frequency is the centroid frequency value of the spectral energy distribution. Arrange the centroid frequencies corresponding to each steady-state background acoustic text analysis frame in chronological order to generate the original centroid trajectory.
[0080] The original trajectory of the spectral centroid is differentially calculated to obtain the change in the frequency of the spectral centroid between adjacent steady-state background acoustic text analysis frames, and a sequence of spectral centroid offset directions is generated based on the positive and negative directions of the change.
[0081] The cumulative summation of the spectral centroid offset direction sequence is calculated to obtain the trajectory of the cumulative spectral centroid offset, which reflects the trend of the total amount of cumulative offset of the spectral centroid over time.
[0082] The event interval values in the event interval fluctuation sequence are synchronized with the cumulative trajectory of the spectral centroid offset according to the time alignment rule, and the event interval values and the cumulative values of the spectral centroid offset within the same time window are extracted to form two-dimensional coordinate points.
[0083] Using the cumulative value of the spectral centroid shift as the horizontal axis and the event interval value of the event interval fluctuation sequence as the vertical axis, the two-dimensional coordinate points extracted in each time window are connected sequentially to form a path diagram of the evolution of the voiceprint state.
[0084] Calculate the Euclidean distance between adjacent two-dimensional coordinate points in the voiceprint state evolution path diagram, determine the overall length parameter of the path diagram based on the cumulative value of the Euclidean distance, and calculate the curvature parameter of the path diagram based on the change of the direction angle of each line segment in the path diagram.
[0085] Extract the trajectory closed loop structure formed by two-dimensional coordinate points in the voiceprint state evolution path diagram, identify the geometric center coordinates and closed loop area parameters of the trajectory closed loop structure, and take the geometric center coordinates of the trajectory closed loop structure as the voiceprint state stability center.
[0086] Using the stable center of the voiceprint state as the origin, the voiceprint state evolution path map is divided into multiple angular sector regions. The distribution density of two-dimensional coordinate points in each angular sector region is statistically analyzed to generate a voiceprint state distribution density map.
[0087] The voiceprint state evolution path map, voiceprint state distribution density map, and corresponding node identifiers are associated and stored to generate a node voiceprint state evolution record.
[0088] In this embodiment of the invention, for example, when the server performs the voiceprint state evolution path graph generation operation, it first processes the steady-state background voiceprint component sequence of node "SN-2024-001". This sequence is a continuous signal with a length of 160,000 sampling points (corresponding to a duration of 10 seconds and a sampling rate of 16kHz). The server calls a framing function, dividing the sequence into frames of 256 sampling points each, shifting the frame by 128 sampling points, resulting in 1249 consecutive steady-state background voiceprint analysis frames. Subsequently, the server applies a Hamming window to each frame signal and performs a Fast Fourier Transform to obtain the power spectrum corresponding to each frame signal. For example, the spectrum of the first frame is distributed in the range of 0-8kHz, and the second, third, and so on are similar.
[0089] Next, the server calculates the centroid frequency of the spectrum for each frame. For the first frame, the server uses the amplitude value corresponding to each frequency point in the spectrum as a weight to calculate the weighted average of all frequency points. The calculated centroid frequency of the first frame is 50.2 Hz. The server then calculates the centroid frequency for the second frame as 50.5 Hz, the third frame as 50.1 Hz, and so on, up to the 1249th frame as 52.3 Hz. The server arranges these centroid frequency values according to their corresponding time sequence (i.e., frame order), generating a raw centroid trajectory of 1249 pixels in length.
[0090] To analyze its changing trend, the server performs a first-order difference calculation on this original trajectory. For example, subtracting the centroid frequency of frame 1 (50.2Hz) from the centroid frequency of frame 2 (50.5Hz) yields a change of +0.3Hz; subtracting frame 2 from frame 3 yields -0.4Hz. The server records all changes and generates an offset direction sequence based on their sign: "+1" represents a frequency increase, "-1" represents a frequency decrease, and "0" represents no change. Then, the server accumulates and sums this direction sequence. Starting from 0, it increments by 1 for each "+1" and decrements by 1 for each "-1". For example, if the first three direction values are [+1, -1, +1], then the first three values of the accumulated trajectory are [1, 0, 1]. This accumulated trajectory ultimately reflects the net offset trend of the spectral centroid relative to the initial state.
[0091] Meanwhile, the server reads the event interval fluctuation sequence generated by the same node "SN-2024-001" within the same 10-second monitoring period from the cache. This sequence contains approximately 100 data points, each representing an event interval value that has been smoothed and logarithmically transformed. To align the two sequences with different sampling rates, the server establishes a time alignment rule. It resamples the spectral centroid offset cumulative trajectory with a basic alignment time window of 10 milliseconds, outputting a value within each alignment window (e.g., taking the average of the cumulative values corresponding to all frames within that window). Similarly, the event interval fluctuation sequence is interpolated to also have a representative value within each alignment window. After synchronization processing, the server obtains paired data for 1000 alignment time windows (10 seconds / 0.01 seconds).
[0092] For each alignment time window (e.g., the i-th window), the server extracts two values: the cumulative spectral centroid offset C_i (e.g., -2.5) within that window and the event interval fluctuation I_i (e.g., 0.448). These two values constitute a two-dimensional coordinate point P_i (X=C_i, Y=I_i). The server processes all 1000 windows sequentially, resulting in 1000 such coordinate points.
[0093] Subsequently, the server plots these points on a two-dimensional plane using the cumulative spectral centroid offset C as the horizontal axis (X-axis) and the event interval fluctuation value I as the vertical axis (Y-axis). It connects points P_1, P_2, P_3, ..., P_1000 sequentially with line segments according to time, thus forming a continuous trajectory line reflecting the evolution of the node's voiceprint state over time, i.e., the voiceprint state evolution path diagram.
[0094] After generating the path map, the server performs geometric feature extraction. It calculates the Euclidean distance between adjacent coordinate points in the map, such as the distance d_1 from point P_1 to P_2, and the distance d_2 from point P_2 to P_3, and sums all 999 distances to obtain the overall length parameter L of the path map (e.g., L = 145.7). Simultaneously, the server calculates the direction angle of each line segment (relative to the horizontal axis) and the absolute value of the difference between the direction angles of adjacent line segments. The average of these differences is calculated as the curvature parameter of the path map (e.g., an average curvature of 12 degrees).
[0095] The server also automatically identifies loop structures formed in the path graph. Through scanning, it discovers that the trajectory between point P_200 and point P_400 approximately forms a closed loop. The server calculates the area enclosed by this closed loop polygon (e.g., 8.2 square units) and determines the coordinates of its geometric center, for example, O_c (X=1.5, Y=0.4). The server marks this geometric center as the stable center of the voiceprint state for that path graph.
[0096] Next, the server divides the entire two-dimensional plane into 12 sector regions with a 30-degree angle, using this stable center O_c as the origin. It counts the number of all 1000 coordinate points in the path graph falling within each sector region, and then calculates the point density (number of points / area) for each region. For example, the 0-30 degree sector region has the highest density, with 120 points distributed within it. Based on these density values, the server generates a polar coordinate voiceprint status distribution density map, visually displaying the clustered areas of status points.
[0097] Finally, the server packages and associates all feature data, including the generated voiceprint state evolution path map (containing the coordinate point sequence), path length parameter L, curvature parameter, identified stable center coordinates O_c, closed loop area, and voiceprint state distribution density map, with the node identifier "SN-2024-001" to form a complete node voiceprint state evolution record, and stores it in the feature database for subsequent use in constructing behavioral fingerprints.
[0098] In this embodiment of the invention, the step of identifying the phase space reconstruction features of the state trajectory in the voiceprint state evolution path diagram, extracting the orbit divergence rate parameters and orbit convergence region boundary parameters of the state trajectory in the phase space, and constructing the voiceprint dynamic behavior fingerprint of each node can be implemented through the following example.
[0099] The two-dimensional coordinate point sequence in the voiceprint state evolution path diagram is used as the original state trajectory. The original state trajectory is reconstructed in phase space. The time delay parameter and the embedding dimension parameter are selected to embed the original state trajectory into the high-dimensional phase space, generating a set of reconstructed trajectory points in the high-dimensional phase space.
[0100] Calculate the local neighborhood radius of each reconstructed trajectory point in the set of reconstructed trajectory points. The local neighborhood radius is the Euclidean distance between the reconstructed trajectory point and its nearest neighbor trajectory point. Determine the global neighborhood radius based on the statistical distribution of the local neighborhood radii of all reconstructed trajectory points.
[0101] For each reconstructed trajectory point in the set of reconstructed trajectory points, select all neighboring trajectory points within a spherical neighborhood centered on the reconstructed trajectory point and with the global neighborhood radius as the radius, and construct a local neighborhood point set for the reconstructed trajectory point;
[0102] Calculate the evolution trend of the number of trajectory points in the local neighborhood set of each reconstructed trajectory point over time, and extract the exponential growth rate of the number of trajectory points in the local neighborhood set during the subsequent evolution process as the trajectory divergence rate parameter.
[0103] Identify the local regions in phase space that are repeatedly visited by the reconstructed trajectory point set, count the number of times each local region is visited by the reconstructed trajectory points, and mark the local regions whose number of visits exceeds the access threshold as trajectory convergence regions;
[0104] Extract the boundary coordinates of the orbit convergence region, calculate the geometric center position and region radius parameter of the orbit convergence region based on the boundary coordinates, and use the region radius parameter as the boundary parameter of the orbit convergence region;
[0105] The original dynamic feature vector is generated by combining the orbit divergence rate parameter and the orbit convergence region boundary parameter, and then the original dynamic feature vector is normalized to generate a normalized dynamic feature vector.
[0106] Nonlinear dimension reduction projection is performed on the normalized dynamic eigenvectors to project the normalized dynamic eigenvectors in the high-dimensional phase space onto the two-dimensional feature plane, generating the coordinates of the feature projection points on the feature plane.
[0107] The Euclidean distance between the coordinates of the feature projection point and the origin of the feature plane is calculated as the dynamic feature intensity parameter, and the orientation angle of the coordinates of the feature projection point relative to the horizontal axis of the feature plane is calculated as the dynamic feature direction parameter.
[0108] The dynamic feature intensity parameter and the dynamic feature direction parameter are combined to generate an acoustic dynamic behavior fingerprint, and the acoustic dynamic behavior fingerprint is associated with and stored with the node identifier of the corresponding node.
[0109] In an embodiment of the invention, exemplarily, when the server performs the voiceprint dynamics behavioral fingerprint construction operation, it first reads the voiceprint state evolution path record of node "SN-2024-001" from the feature database and extracts the two-dimensional coordinate point sequence as the original state trajectory. This trajectory contains 1000 points arranged in chronological order, such as P1(-2.5, 0.448), P2(-2.3, 0.455), ..., P1000(3.1, 0.510). The server calls the phase space reconstruction algorithm to first perform time delay analysis on the original trajectory. By calculating the time required for the autocorrelation function to decay to 1-1 / e, the time delay parameter τ=3 is determined. Next, the server uses the spurious nearest neighbor method to determine the embedding dimension parameter m=3. Subsequently, the server embeds the two-dimensional trajectory into the three-dimensional phase space according to the formula V_i=[P_i,P_{i+τ},P_{i+2τ}], generating a set containing approximately 997 three-dimensional vectors (reconstructed trajectory points), such as the three-dimensional coordinates of V1=[P1,P4,P7], V2=[P2,P5,P8], etc.
[0110] Next, the server calculates the distance between each point in the reconstructed trajectory point set and its nearest neighbor (the one with the smallest Euclidean distance). For example, the distance between point V1 and its nearest neighbor V20 is 0.15, and the distance between point V2 and its nearest neighbor V45 is 0.22. The server calculates all 997 such nearest neighbor distances and calculates their distribution. It finds that the upper quartile of these distances is 0.25, and therefore sets the global neighborhood radius R to this value.
[0111] Then, the server constructs a spherical neighborhood in three-dimensional phase space, centered on each reconstructed trajectory point and with a radius of R=0.25. For example, for the center point V1, the server calculates its Euclidean distance to all other points in the set and includes points with a distance less than 0.25 (e.g., V20, V150, V300) in its local neighborhood set N(V1). The server constructs its own local neighborhood set for each of the 997 center points.
[0112] The server then analyzes the evolutionary dynamics of each local neighborhood point set over time. Taking the neighborhood point set N(V1) = {V1, V20, V150, V300} of point V1 as an example, the server tracks the positions of these four points in subsequent evolutionary steps (e.g., the next five time steps). It finds that after one step, three of these four points are still relatively close to each other (still within a sphere of radius R); after two steps, two remain; and after three steps, one remains. The server fits a curve showing the decay of this number over time steps and calculates its exponential decay rate constant. In fact, the orbital divergence rate parameter λ is defined as this decay rate (or an estimate of the Lyapunov exponent). The server calculates an estimated λ value of 0.12 bits / step for point V1. The server performs this calculation for all center points and takes the median of all λ values (e.g., 0.15 bits / step) as the orbital divergence rate parameter for the entire trajectory.
[0113] Simultaneously, the server scans the entire 3D phase space, identifying frequently visited regions by the reconstructed trajectory points. It divides the phase space into small grid cells and counts the number of times each cell is "hit" by a trajectory point. The server finds a small region centered at coordinates (0.5, -0.2, 1.1) that was visited 25 times, far exceeding the set visit threshold (e.g., 15 times). The server marks this region as a trajectory convergence region. The server extracts the coordinates of all trajectory points that visited this region, calculates the boundaries of the point cloud formed by these points (minimum and maximum values in each dimension), and uses this to calculate the geometric center of the region (i.e., the aforementioned center coordinates) and the average radius of the point cloud (e.g., r = 0.18). The server uses this radius r as the boundary parameter of the convergence region.
[0114] Subsequently, the server combines the two extracted core parameters, the orbital divergence rate parameter λ=0.15 and the orbital convergence region boundary parameter r=0.18, into an original dynamic feature vector [0.15, 0.18]. The server calls the normalization module to perform Z-score normalization on this vector using the mean and standard deviation of λ and r for this node learned from historical normal data, to obtain a normalized dynamic feature vector, for example [0.8, -0.3].
[0115] Next, the server uses a trained t-SNE nonlinear dimensionality reduction model to project this two-dimensional normalized feature vector (which is actually in an abstract feature space) onto a predefined two-dimensional feature plane. The projection yields specific coordinates, for example (1.2, 0.5). The server calculates the Euclidean distance from this projected point to the origin (0,0) of the feature plane, obtaining the dynamic feature intensity parameter S = sqrt(1.2). 2 +0.5 2 )≈1.3. At the same time, the server calculates the counterclockwise angle between the point and the horizontal axis (due east) of the feature plane, and obtains the dynamic feature direction parameter θ=arctan(0.5 / 1.2)≈22.6 degrees.
[0116] Finally, the server defines this pair of parameters (S=1.3, θ=22.6°) as the acoustic dynamics behavioral fingerprint of node "SN-2024-001" in the current monitoring period. The server encrypts and stores this fingerprint together with the node identifier "SN-2024-001" and the timestamp "2024-05-10 14:30" in the behavioral fingerprint database, completing the feature modeling of this node in this monitoring period.
[0117] In this embodiment of the invention, the step of comparing the acoustic dynamic behavior fingerprints of each node in the current monitoring period with the acoustic dynamic behavior fingerprints of each node in the historical normal operation period, calculating the fingerprint deviation vector field between the two, and determining the abnormal nodes with non-zero divergence values in the deviation vector field based on the divergence distribution of the fingerprint deviation vector field can be implemented through the following example.
[0118] Extract the acoustic dynamic behavior fingerprints of each node during the historical normal operation cycle from the database, and construct a historical fingerprint feature library. The historical fingerprint feature library contains the historical values of the dynamic feature intensity parameter and the historical value of the dynamic feature direction parameter of each node during the historical normal operation cycle.
[0119] Obtain the acoustic dynamic behavior fingerprint of each node within the current monitoring period, and extract the current value of the dynamic feature intensity parameter and the current value of the dynamic feature direction parameter of each node within the current monitoring period;
[0120] For each node, the difference between the current value and the historical value of the dynamic characteristic intensity parameter is calculated as the intensity deviation value, and the difference between the current value and the historical value of the dynamic characteristic direction parameter is calculated as the direction deviation value.
[0121] The intensity deviation value and the direction deviation value are combined to generate a two-dimensional deviation vector for each node. The magnitude of the two-dimensional deviation vector is determined by the intensity deviation value, and the direction angle of the two-dimensional deviation vector is determined by the direction deviation value.
[0122] All the two-dimensional deviation vectors of the nodes are spatially arranged according to the physical location coordinates of the nodes inside the building, and a fingerprint deviation vector field is constructed with the physical location coordinates as the independent variable and the two-dimensional deviation vector as the dependent variable.
[0123] Divergence calculation is performed on the fingerprint deviation vector field. At the location coordinates of each node, the divergence value of the two-dimensional deviation vector of that node is calculated. The divergence value reflects the degree of divergence or convergence of the two-dimensional deviation vector at that node.
[0124] Nodes with non-zero divergence values are marked as divergence anomaly nodes. The type of divergence anomaly node is determined according to the positive or negative sign of the divergence value. A positive divergence value indicates that the fingerprint deviation vector field at the node exhibits divergence characteristics, while a negative divergence value indicates that the fingerprint deviation vector field at the node exhibits convergence characteristics.
[0125] Obtain the divergence values of neighboring nodes around the divergence anomaly node, calculate the divergence gradient vector between the divergence anomaly node and its neighboring nodes, and determine the anomaly propagation direction of the fingerprint deviation vector field based on the direction of the divergence gradient vector.
[0126] The divergence outliers are sorted according to the absolute value of their divergence to generate a priority sequence. The larger the absolute value of the divergence, the more significant the fingerprint deviation of the node.
[0127] Output all nodes with non-zero absolute divergence values in the divergence anomaly priority sequence as anomaly node set, and append the corresponding divergence value and divergence gradient vector information to each node in the anomaly node set.
[0128] In this embodiment of the invention, for example, when the server performs the abnormal node identification operation, it first extracts all voiceprint dynamic behavioral fingerprint records stored by node "SN-2024-001" during its normal operation over the past 30 days from the behavioral fingerprint database. The server calculates statistics for the dynamic feature intensity parameter S and orientation parameter θ in these historical records to obtain the historical baseline of the node: historical value of intensity parameter S_hist (mean 1.0, standard deviation 0.1), and historical value of orientation parameter θ_hist (mean 20.0 degrees, standard deviation 5.0 degrees). The server performs the same operation on all 50 monitored nodes to construct a complete historical fingerprint feature database.
[0129] Subsequently, the server retrieves the acoustic dynamic behavior fingerprints of each node that have just been calculated within the current monitoring period (e.g., the 5-minute window of "2024-05-10 14:30"). For node "SN-2024-001", the server reads its current fingerprint: the current value of the intensity parameter S_curr=1.3, and the current value of the orientation parameter θ_curr=22.6 degrees.
[0130] Next, the server calculates the deviation of the current fingerprint from the historical baseline for each node. For node "SN-2024-001", the intensity deviation ΔS = S_curr - S_hist_mean = 1.3 - 1.0 = 0.3. The orientation deviation Δθ = θ_curr - θ_hist_mean = 22.6 - 20.0 = 2.6 degrees. These two values together constitute a deviation in polar coordinates. The server converts this into an deviation vector in a two-dimensional plane: the magnitude (amplitude) of this vector is determined by the intensity deviation ΔS (0.3 units), and its orientation angle is the orientation deviation Δθ (2.6 degrees). In a Cartesian coordinate system, this vector can be represented as (ΔS*cos(Δθ), ΔS*sin(Δθ)).
[0131] The server obtains the physical deployment coordinates of all 50 nodes inside the building (e.g., node "SN-2024-001" is located at coordinates (x1, y1) = (10.5, 25.0) meters on the building floor plan). It uses the physical coordinates of each node as the independent variable and the calculated two-dimensional deviation vector of that node as the dependent variable, spatially arranging them to construct a fingerprint deviation vector field covering the entire pipeline monitoring area. This vector field visually displays the direction and magnitude of the acoustic signature dynamics at each node's deviation from its historical normal state.
[0132] Subsequently, the server performs divergence calculations on the vector field to quantify the source or convergence point of the deviation pattern. At the coordinates (x1, y1) of node "SN-2024-001", the server estimates the divergence using the central difference method. It finds the four neighboring nodes of this node in the network topology (e.g., in the east, south, west, and north directions) and obtains their deviation vectors. Assuming the vector of the east neighbor node is (0.2, 0.1), the vector of the west neighbor node is (-0.1, 0.0), the vector of the south neighbor node is (0.0, -0.2), and the vector of the north neighbor node is (0.1, 0.3), the server calculates the divergence value at this point div=(∂Fx / ∂x+∂Fy / ∂y), approximating the partial derivative by the difference between the vectors of the neighboring nodes. The calculated divergence value at node "SN-2024-001" is approximately +0.85.
[0133] The server iterates through all node coordinates to calculate divergence. Nodes with non-zero divergence values (absolute values exceeding a small threshold, such as 0.05) are marked as divergence anomalies. Based on the sign of the divergence value, the server further distinguishes categories: for example, node "SN-2024-001" has a divergence value of +0.85 (positive) and is marked as a "divergence anomaly," indicating that this node may be the "source point" of the anomalous state or a significantly affected point, and its deviation pattern tends to spread to the surrounding areas; while another node "SN-2024-018" has a divergence value of -0.62 (negative) and is marked as a "convergence anomaly," indicating that this node may be the "sink point" of the anomalous influence.
[0134] To determine the potential direction of anomaly propagation, the server calculates the difference vector between the divergence value of each divergence anomaly node and the divergence values of all its directly adjacent nodes, i.e., the divergence gradient vector. For node "SN-2024-001", its divergence gradient vector points in the direction of the fastest increase in its divergence value, which suggests a deviation from the direction from which the disturbance may have propagated or the direction of its increasing intensity.
[0135] Then, the server sorts all divergence anomaly nodes in descending order based on their absolute divergence values, generating a priority sequence for these nodes. For example, node "SN-2024-015" has an absolute divergence value of 1.2, ranking first; node "SN-2024-001" has a value of 0.85, ranking second. A larger absolute divergence value indicates a more significant spatial deviation at that node, making it more likely to be a core anomaly.
[0136] Finally, the server includes all nodes with non-zero divergence values in the priority sequence (i.e., all divergence anomaly nodes) into an anomaly node set as the output for this monitoring cycle. Each node entry in the set includes its detailed divergence value, divergence type (positive / negative), and calculated divergence gradient vector information. For example, the output set is: {Node ID:SN-2024-015, divergence value: +1.2, type: divergence, gradient vector: (0.1, 0.2); Node ID:SN-2024-001, divergence value: +0.85, type: divergence, gradient vector: (0.05, -0.1);...}. This set provides clear spatial anomaly clues for subsequent leakage source localization.
[0137] In this embodiment of the invention, the step of constructing a local network topology subgraph with the abnormal node as the central node, obtaining the curl value of the fingerprint deviation vector field of the acoustic fingerprint of each adjacent node in the local network topology subgraph, and determining the propagation direction of the acoustic fingerprint anomaly between adjacent nodes based on the positive and negative directions of the curl value can be implemented through the following example.
[0138] Based on the design drawings of the building's internal water supply network, the physical connection relationship between all nodes is extracted, and a global network topology map is constructed. The nodes of the global network topology map correspond to the deployment positions of the voiceprint sensing units, and the edges correspond to the water supply pipeline connection relationship between adjacent nodes.
[0139] Select an abnormal node from the set of abnormal nodes as the center node, and perform a breadth-first traversal in the global pipeline topology graph starting from the center node to extract all nodes whose distance from the center node is within a preset number of hops to form a local pipeline topology subgraph.
[0140] Obtain the two-dimensional deviation vector in the fingerprint deviation vector field of each node in the local pipeline network topology subgraph, and use the two-dimensional deviation vector of each node as the deviation vector field value of that node;
[0141] For each node in the local network topology subgraph, calculate the difference vector of the deviation vector field values between the node and all its neighboring nodes, and determine the curl value at the node based on the cross product of the difference vector and the direction vector of the connecting edge between the node and its neighboring nodes.
[0142] Nodes with non-zero curl values are marked as curl-active nodes, and the rotation direction characteristics of the curl-active nodes are determined according to the positive or negative sign of the curl value. A positive curl value indicates that the deviation vector field of the node exhibits a counterclockwise rotation trend, and a negative curl value indicates that the deviation vector field of the node exhibits a clockwise rotation trend.
[0143] For each connecting edge in the local pipeline topology subgraph, obtain the sign of the curl value of the nodes at both ends of the connecting edge. When the signs of the curl values of the nodes at both ends are the same, determine that the direction of the abnormal propagation of the acoustic pattern on the connecting edge is from the node with the larger absolute value of curl to the node with the smaller absolute value of curl.
[0144] When the signs of the curl values at the two ends are different, the rotation direction of the deviation vector field of the node with the positive curl value is compared with the direction of the connecting edge, and the direction of abnormal propagation of the acoustic pattern is determined based on the angle between the rotation direction of the deviation vector field and the direction of the connecting edge.
[0145] Record the abnormal propagation direction of acoustic fingerprints on each connecting edge in the local pipeline network topology sub-graph, and generate an abnormal propagation direction mapping table of acoustic fingerprints with connecting edges as the recording unit. The abnormal propagation direction mapping table of acoustic fingerprints includes the start node identifier, end node identifier and propagation direction identifier of each connecting edge.
[0146] A directed propagation subgraph is constructed based on the soundprint anomaly propagation direction mapping table. The nodes of the directed propagation subgraph are the same as the nodes of the local pipeline topology subgraph. The directed edges of the directed propagation subgraph are determined by the propagation direction recorded in the soundprint anomaly propagation direction mapping table.
[0147] Identify the directed paths in the directed propagation subgraph, extract the starting nodes of the directed paths as candidate nodes for the source of voiceprint anomalies, and extract the ending nodes of the directed paths as candidate nodes for the convergence of voiceprint anomalies.
[0148] In this embodiment of the invention, for example, when the server performs the abnormal propagation direction analysis of voiceprints, it first loads the digital design drawings of the building's water supply network from the storage system and automatically parses the deployment locations of all voiceprint sensing units and the pipe connection relationships between them. The server constructs a global pipe network topology map, where each graph node corresponds to a physical sensor node (e.g., "SN-2024-001"), and each undirected edge represents a directly connected water supply pipe segment. For example, the drawing shows that nodes A, B, and C are connected in series by pipes.
[0149] Subsequently, the server selects the node "SN-2024-015" with the largest absolute divergence value from the set of abnormal nodes output in the previous stage as the central node. In the global network topology graph, the server executes a breadth-first search algorithm starting from this central node. The algorithm defaults to a hop count of 3, meaning it searches for all nodes whose pipeline distance from the central node is no more than 3 hops. After the traversal, the server extracts a total of 12 nodes, including the central node itself, its direct upstream and downstream neighbors (1 hop), the neighbors' neighbors (2 hops), and nodes at the next higher level (3 hops), along with all the connecting edges between them, forming a local network topology subgraph. This subgraph focuses on the network area surrounding the core abnormal point.
[0150] Next, the server extracts the two-dimensional deviation vector corresponding to each node within the local subgraph from the fingerprint deviation vector field data. For example, the deviation vector of node "SN-2024-015" is (0.8, 0.6), the vector of its upstream neighbor "SN-2024-010" is (0.2, 0.3), and the vector of its downstream neighbor "SN-2024-018" is (0.1, -0.2). The server treats these vectors as the deviation vector field values at each node.
[0151] Then, the server calculates the curl at each node in the subgraph to quantify the rotational characteristics of the deviation field. Taking node "SN-2024-015" as an example, the server calculates the difference in its deviation vector with its upstream node "SN-2024-010": ΔV=(0.8-0.2,0.6-0.3)=(0.6,0.3). Simultaneously, it obtains the unit vector of the pipe direction pointing from node "SN-2024-010" to "SN-2024-015", assumed to be (1,0) (due east). The server calculates the two-dimensional cross product (scalar result) of the difference vector ΔV and the pipe direction vector: Curl_z=ΔV_x*Dir_y-ΔV_y*Dir_x=0.60-0.31=-0.3. This value is the contribution of this connecting edge to the curl of the central node. The server performs a similar calculation and sums all the connecting edges of the central node, resulting in an overall curl value of -0.5 for node "SN-2024-015".
[0152] The server traverses all nodes in the subgraph to calculate the curl. It marks nodes with absolute curl values exceeding a threshold (e.g., 0.1) as curl-active nodes. Based on the sign, the server determines that a negative curl value (e.g., -0.5) indicates that the deviation vector field at that node exhibits a clockwise rotation trend; a positive curl value indicates a counterclockwise rotation trend.
[0153] Next, the server determines the abnormal propagation direction on each connecting edge within the subgraph. Taking the edge connecting nodes "SN-2024-010" (curl value +0.2) and "SN-2024-015" (curl value -0.5) as an example, since the curl values of the two nodes have different signs (one positive and one negative), the server adopts the second rule: the curl of node "SN-2024-010" is positive, and its deviation field rotates counterclockwise. The server compares the rotation direction (counterclockwise) with the pipe direction from "SN-2024-010" to "SN-2024-015". Based on the preset angle relationship model, it determines that the abnormal propagation direction of the voiceprint is from node "SN-2024-010" to node "SN-2024-015". For adjacent nodes with the same curl sign, the propagation direction is determined to be from the one with the larger absolute curl value to the one with the smaller absolute curl value.
[0154] The server records the analysis results of all edges within the subgraph and generates a mapping table of abnormal voiceprint propagation directions. Each record contains the starting node ID, the ending node ID, and the direction identifier (e.g., "SN-2024-010→SN-2024-015").
[0155] Based on this mapping table, the server constructs a directed propagation subgraph. The nodes in this subgraph are the same as the original local topology subgraph, but each edge is assigned a direction according to the mapping table. Finally, the server identifies the critical path in this directed graph. It finds the longest directed path: SN-2024-005→SN-2024-008→SN-2024-010→SN-2024-015. The server extracts the starting node "SN-2024-005" as a candidate node for the source of the voiceprint anomaly, and marks the ending node "SN-2024-015" and other nodes with only inbound edges as candidate nodes for the convergence of voiceprint anomalies. This candidate node information provides a clear topological flow basis for subsequent backtracking of the leakage source.
[0156] In this embodiment of the invention, the step of inputting the abnormal node and the abnormal propagation direction of the acoustic fingerprint into the leakage source reverse tracing model, and tracing back step by step in the local pipeline topology subgraph according to the reverse pointing relationship of the abnormal propagation direction of the acoustic fingerprint to generate a set of candidate nodes for leakage sources can be implemented through the following example.
[0157] All nodes in the abnormal node set are used as the initial node set for reverse tracing. The propagation direction in the abnormal propagation direction mapping table is reversed to generate a reverse propagation direction mapping table. The direction identifier of each record in the reverse propagation direction mapping table is opposite to the original propagation direction.
[0158] Starting from each node in the initial node set, search the reverse propagation direction mapping table for the reverse propagation record with that node as the termination node, and extract the starting node in that reverse propagation record as the previous level backtracking node.
[0159] Add the extracted previous backtracking node to the current backtracking node set, and repeat the operation of finding the backpropagation record for each node in the current backtracking node set until no new previous backtracking node can be found.
[0160] Record the sequence of backtracking nodes obtained in each backtracking operation. Each backtracking node sequence corresponds to a backtracking path that starts from the initial node and backtracks to the source node step by step along the reverse propagation direction.
[0161] Multiple backtracking paths are merged, and nodes appearing in all backtracking paths are extracted to form an initial set of candidate nodes for leakage sources. The number of times each node appears in the backtracking path is recorded as a backtracking frequency parameter for each node in the initial set of candidate nodes for leakage sources.
[0162] A backtracking tree structure is constructed based on the topology of the backtracking path. The root node of the backtracking tree structure is the initial node, and the child nodes are the parent nodes pointed to by the reverse propagation direction. The depth of the backtracking tree structure is determined by the length of the backtracking path.
[0163] Calculate the branch factor of each node in the backtrace tree structure, where the branch factor is the number of child nodes that the node has in the backtrace tree structure, and mark the nodes whose branch factors exceed the branch threshold as multi-path convergence nodes.
[0164] For each node in the initial set of candidate leakage sources, extract the depth value, backtracking frequency parameter, and branching factor of the node in the backtracking tree structure, and select the node with the smallest depth value, the largest backtracking frequency parameter, and the smallest branching factor as the candidate leakage source with the highest confidence.
[0165] Sort all nodes in the initial candidate node set of leakage sources in ascending order of depth value and descending order of backtracking frequency parameter to generate a sorted list of candidate nodes for leakage sources.
[0166] Select the nodes with the highest sorted position from the candidate nodes of leakage sources as primary candidate nodes of leakage sources, and select the nodes with the lowest sorted position as secondary candidate nodes of leakage sources, thereby generating a set of candidate nodes of leakage sources containing primary and secondary candidate nodes of leakage sources.
[0167] In this embodiment of the invention, for example, when the server performs the reverse tracing operation of the leakage source, it first sets all nodes in the abnormal node set {SN-2024-015, SN-2024-001, SN-2024-018} determined in the previous step as the initial node set for reverse tracing. Simultaneously, the server reads the generated voiceprint anomaly propagation direction mapping table, which contains records such as "SN-2024-010->SN-2024-015" and "SN-2024-008->SN-2024-010". The server initiates a direction reversal process, reversing the propagation direction of each record in the mapping table to generate a reverse propagation direction mapping table. The original record “SN-2024-010->SN-2024-015” is reversed to “SN-2024-015->SN-2024-010”, meaning that in the tracing model, the anomaly is considered to be likely to propagate backward from node 015 to its upstream node 010.
[0168] Subsequently, the server initiates a multi-threaded backtracking algorithm. It uses each node in the initial node set as the starting point for independent tracing. For example, for the initial node "SN-2024-015", the server uses it as the current tracing starting point and searches the reverse propagation direction mapping table for all records ending with "SN-2024-015". It finds the record "SN-2024-015->SN-2024-010", thus extracting node "SN-2024-010" as the next-level backtracking node. The server adds node "SN-2024-010" to a temporary "current backtracking node set". Next, the server repeats this operation for the new node "SN-2024-010" in the set: it searches the reverse mapping table for a record ending with "SN-2024-010", finds "SN-2024-010->SN-2024-008", and extracts node "SN-2024-008" as the next-level node. This process continues iterating until a node (e.g., "SN-2024-005") is searched and no record ending with it exists in the reverse mapping table, at which point the process terminates.
[0169] The server fully records the tracing chain starting from "SN-2024-015" and generates a backtracking path: SN-2024-015←SN-2024-010←SN-2024-008←SN-2024-005. The server performs the same tracing operation in parallel on other nodes in the initial set, such as "SN-2024-001" and "SN-2024-018", potentially resulting in other backtracking paths, such as "SN-2024-001←SN-2024-005" and "SN-2024-018←SN-2024-015←…".
[0170] After tracing is complete, the server merges all nodes appearing in all backtracking paths to form an initial set of candidate nodes for the leakage source, such as {SN-2024-015, SN-2024-010, SN-2024-008, SN-2024-005, SN-2024-001, SN-2024-018}. The server also counts the number of times each node appears in all backtracking paths as a backtracking frequency parameter. For example, node “SN-2024-005” appears in two paths, with a frequency of 2; node “SN-2024-015” appears in two paths, with a frequency of 2; and node “SN-2024-010” appears in one path, with a frequency of 1.
[0171] Next, the server constructs a backtracking tree structure based on the topological relationships of all backtracking paths. The leaf nodes of this tree are the initial anomaly nodes, and their parent nodes are their parent backtracking nodes. For example, node "SN-2024-005" is the parent node of nodes "SN-2024-010" and "SN-2024-001". The server calculates the branch factor of each node in the tree, i.e., the number of its child nodes. Node "SN-2024-005" has two child nodes and a branch factor of 2; node "SN-2024-015" has one child node (SN-2024-018) and a branch factor of 1. The server marks nodes with branch factors exceeding a threshold (e.g., 2) as multi-path convergence nodes, indicating that anomalies may converge there from multiple downstream paths.
[0172] Then, the server performs a comprehensive evaluation of each node in the initial candidate set. It extracts three key parameters: the node's depth in the backtracking tree (the number of hops from the leaf node, with the root node having a depth of 0), the backtracking frequency parameter, and the branch factor. According to preset rules, the server identifies the node with the smallest depth (closest to the upstream), the highest backtracking frequency (involved in multiple backtracking attempts), and the smallest branch factor (not a complex confluence point) as the candidate node with the highest confidence in being a leakage source. In this example, node "SN-2024-005" has a depth of 0, a frequency of 2, and a branch factor of 2, resulting in the highest overall confidence.
[0173] Finally, the server sorts all candidate nodes. The primary sorting criterion is ascending depth (upstream priority), and the secondary criterion is descending backtracking frequency (higher frequency priority), generating a sorted list of candidate nodes for leakage sources. The sorted result may be: [SN-2024-005, SN-2024-008, SN-2024-010, SN-2024-015, SN-2024-001, SN-2024-018]. The server identifies the first-ranked node "SN-2024-005" as the primary candidate node for leakage sources, and the second and third-ranked nodes "SN-2024-008" and "SN-2024-010" as secondary candidate nodes for leakage sources. Together, they form the final set of candidate nodes for leakage sources, which is used for final location determination.
[0174] In this embodiment of the invention, the step of selecting the node with the largest absolute value of the divergence of the fingerprint deviation vector field from the candidate node set of leakage sources as the leakage source location node, and outputting the node identifier corresponding to the leakage source location node as the leakage source location result, can be implemented through the following example.
[0175] Obtain the absolute value of the divergence of the fingerprint deviation vector field corresponding to each candidate node in the candidate node set of leakage sources, and use the absolute value of the divergence as the anomaly intensity score of the candidate node.
[0176] The abnormal intensity scores of all candidate nodes in the candidate node set of leakage sources are sorted in descending order to generate an abnormal intensity score sorted list, with the candidate node with the highest abnormal intensity score at the top of the sorted list.
[0177] The candidate node with the highest anomaly intensity score is extracted from the anomaly intensity score ranking list as the primary localization node, and the candidate node with the second highest anomaly intensity score is extracted as the secondary localization node.
[0178] Obtain the depth value of the main location node in the backtrace tree structure. When the depth value of the main location node is greater than the depth threshold, mark the node chain consisting of the main location node and all its ancestor nodes as a suspected leakage propagation path.
[0179] Obtain the shortest path distance between the main positioning node and the auxiliary positioning node in the global pipeline topology map. When the shortest path distance is less than the distance threshold, mark all nodes on the path between the main positioning node and the auxiliary positioning node as the node set of the leakage impact range.
[0180] Calculate the curl value of the fingerprint deviation vector field for each node in the node set of the leakage influence range, extract the node with the largest absolute value of the curl value as the leakage center verification node, and compare the leakage center verification node with the main location node;
[0181] When the leakage center verification node is the same as the main location node, the main location node is determined as the final leakage source location node. When the leakage center verification node is different from the main location node, the leakage center verification node and the main location node are jointly marked as dual-source leakage candidate nodes.
[0182] The node with the larger absolute value of divergence is selected from the candidate nodes of dual-source leakage as the final leakage source location node, and the node identifier corresponding to the final leakage source location node is used as the main leakage source location result.
[0183] The main leakage source location result is associated and encapsulated with the node identifiers of other candidate nodes in the abnormal intensity score ranking list to generate a leakage source location result data packet containing the main leakage source identifier and the list of auxiliary candidate node identifiers.
[0184] The data packet containing the leakage source location result is transmitted to the data receiving interface of the monitoring platform. The monitoring platform triggers the leakage warning display operation of the corresponding node based on the main leakage source identifier in the data packet containing the leakage source location result.
[0185] On the pipeline topology visualization interface of the monitoring platform, the node location corresponding to the final leakage source is highlighted, and the abnormal intensity level of each node in the node set of the leakage impact range is marked with different color depths.
[0186] In this embodiment of the invention, for example, when the server performs the final leak source location and result output operation, it first obtains the absolute value of the divergence of the fingerprint deviation vector field corresponding to each node in the candidate node set of the leak source from the feature database. For example, for candidate node "SN-2024-005", the server reads its absolute divergence value as 1.2; for node "SN-2024-008", the absolute divergence value is 0.9; and for node "SN-2024-010", the absolute divergence value is 0.7. The server directly uses these absolute divergence values as the anomaly intensity score for each candidate node.
[0187] Next, the server invokes a sorting algorithm to sort the candidate node set {SN-2024-005, SN-2024-008, SN-2024-010} in descending order of anomaly intensity score, generating an anomaly intensity score sorted list: [(SN-2024-005, 1.2), (SN-2024-008, 0.9), (SN-2024-010, 0.7)]. From this list, the server extracts the node with the highest score, “SN-2024-005”, as the primary location node, and the node with the second highest score, “SN-2024-008”, as the secondary location node.
[0188] Subsequently, the server performs cross-validation and range confirmation. It first queries the backtracking tree structure to obtain the depth value of the main location node "SN-2024-005". The server finds that the depth value of this node is 0 (i.e., the root of the tree), which is less than the preset depth threshold (e.g., 3), so it determines that its location is reasonable and there is no need to mark a long-distance propagation path.
[0189] Then, the server calculates the shortest pipe path distance between the primary locating node "SN-2024-005" and the secondary locating node "SN-2024-008" in the global pipeline topology map. The calculation shows that they are directly connected, with a distance of 1 hop. This distance is less than a preset distance threshold (e.g., 3 hops). Therefore, the server marks these two nodes and all nodes along their path (in this example, themselves) as the set of nodes within the leakage impact range {SN-2024-005, SN-2024-008}.
[0190] To further verify, the server calculated the curl value of the fingerprint deviation vector field for each node within the affected area set. The query revealed that node "SN-2024-005" had a curl value of +0.8, and node "SN-2024-008" had a curl value of +0.3. The server extracted the node with the largest absolute curl value, "SN-2024-005" (|0.8|>|0.3|), and designated it as the leakage center verification node. The server compared this verification node with the main location node and found that both were "SN-2024-005," indicating complete consistency.
[0191] Based on this consistency verification, the server ultimately determined the primary location node "SN-2024-005" as the final leak source location node. The server uses the identifier "SN-2024-005" of this node as the primary leak source location result.
[0192] Next, the server encapsulates the results. It associates the main result "SN-2024-005" with the identifiers of the remaining candidate nodes ("SN-2024-008", "SN-2024-010") in the anomaly intensity score ranking list, generating a structured leak source location result data packet. This data packet contains the main leak source identifier, a list of auxiliary candidate nodes, the anomaly intensity score of each node, and information on its impact range.
[0193] Then, the server transmits the data packet to the designated data receiving interface of the building monitoring platform via its network communication module. Upon receiving the data packet, the monitoring platform's data parsing service immediately triggers a leakage early warning display operation for that node in the system's alarm management module based on the main leakage source identifier "SN-2024-005" in the packet, generating a level-one alarm log.
[0194] Simultaneously, the monitoring platform utilizes its visualization engine to highlight the sensor node "SN-2024-005" in a striking red flashing pattern on the topology visualization interface of the building's water supply network. For another node in the leakage impact area set, "SN-2024-008," the platform, based on its anomaly intensity score of 0.9, marks it as orange (representing moderate anomaly) according to a preset color mapping rule. A location report window also pops up on the platform interface, clearly displaying the main leakage source, supporting evidence, and impact area, completing a full closed loop from server analysis to terminal alarm.
[0195] This invention provides a computer device 100, which includes a processor and a non-volatile memory storing computer instructions. When the computer instructions are executed by the processor, the computer device 100 executes the aforementioned intelligent monitoring method for concealed building leaks based on pipe network acoustic signatures. Figure 2 As shown, Figure 2 This is a structural block diagram of a computer device 100 provided in an embodiment of the present invention. The computer device 100 includes a memory 111, a processor 112, and a communication unit 113. To enable data transmission or interaction, the memory 111, processor 112, and communication unit 113 are electrically connected to each other directly or indirectly. For example, these components can be electrically connected to each other through one or more communication buses or signal lines.
[0196] For illustrative purposes, the foregoing description has been made with reference to specific embodiments. However, the foregoing illustrative discussions are not intended to be exhaustive or to limit the present disclosure to the precise forms disclosed. Numerous modifications and variations are possible in accordance with the foregoing teachings. These embodiments were chosen and described in order to best illustrate the principles of the present disclosure and its practical application, thereby enabling those skilled in the art to best utilize the disclosure and to employ various embodiments with different modifications to suit a particular intended application.
Claims
1. A method for intelligent monitoring of concealed building leaks based on the acoustic signature characteristics of pipe networks, characterized in that, include: The original voiceprint signal sequence is acquired by the voiceprint sensing units deployed at each node of the water supply network inside the building during a continuous monitoring period. The original voiceprint signal sequence carries the node identifier and the acquisition timestamp identifier corresponding to each voiceprint sensing unit. Each of the original acoustic signature signal sequences is decomposed to obtain its corresponding steady-state background acoustic signature component sequence and transient impact acoustic signature component sequence; The transient impact acoustic component sequence is grouped according to the node identifier, and the transient impact event time axis of each node is obtained by arranging it in time sequence according to the acquisition timestamp identifier. An event interval fluctuation sequence is constructed based on the time interval between adjacent events. The event interval fluctuation sequence is correlated and mapped with the spectral centroid offset trajectory of the steady-state background voiceprint component sequence to generate a voiceprint state evolution path diagram for each node. Based on the state trajectory of the voiceprint state evolution path diagram, the voiceprint dynamic behavior fingerprint of each node is constructed, and the voiceprint dynamic behavior fingerprint of each node under current and historical normal working conditions is compared. Abnormal nodes are determined based on the divergence distribution of the fingerprint deviation vector field. Based on the abnormal nodes, a local pipeline network topology subgraph is constructed. The direction of abnormal acoustic fingerprint propagation is determined according to the curl direction of the deviation vector field of the fingerprints of adjacent nodes in the local pipeline network topology subgraph. The set of candidate nodes for leakage source is obtained by backtracking the abnormal nodes and the direction of abnormal acoustic fingerprint propagation. The node with the largest absolute value of the fingerprint deviation vector field divergence in the candidate node set of the leakage source is selected as the leakage source localization result output.
2. The method according to claim 1, characterized in that, The process of decomposing each of the original acoustic signature signal sequences to obtain their respective steady-state background acoustic signature component sequences and transient impact acoustic signature component sequences includes: The original voiceprint signal sequence is input into a separation network based on adaptive mode decomposition. The separation network decomposes the original voiceprint signal sequence into multiple intrinsic mode function components through an iterative filtering process. Each intrinsic mode function component corresponds to a voiceprint oscillation mode at a different time scale. Calculate the instantaneous frequency variance of each intrinsic mode function component. Intrinsic mode function components with instantaneous frequency variance values lower than a preset frequency variance threshold are classified into a steady-state background soundprint candidate component set, and intrinsic mode function components with instantaneous frequency variance values not lower than a preset frequency variance threshold are classified into a transient impact soundprint candidate component set. Linear superposition and reconstruction are performed on each intrinsic mode function component in the steady-state background voiceprint candidate component set to generate the steady-state background voiceprint component sequence. The instantaneous frequency of the steady-state background voiceprint component sequence remains stable over time. Linear superposition and reconstruction are performed on each intrinsic mode function component in the transient impact acoustic fingerprint candidate component set to generate an initial transient impact acoustic fingerprint component sequence. The initial transient impact acoustic fingerprint component sequence is then subjected to amplitude threshold screening, and transient impact waveform segments with amplitudes exceeding multiples of the background noise amplitude are retained to obtain the transient impact acoustic fingerprint component sequence.
3. The method according to claim 1, characterized in that, The process of grouping the transient impact acoustic component sequences according to the node identifiers, arranging them sequentially according to the acquisition timestamp identifiers to obtain the time axis of transient impact events for each node, and constructing an event interval fluctuation sequence based on the time interval between adjacent events includes: Extract all transient impact acoustic component sequences corresponding to the same node identifier from the data buffer, sort them according to the order of the acquisition timestamp identifiers associated with each transient impact acoustic component sequence, and generate a transient impact event sequence arranged in ascending order of time. Traverse the sequence of transient impact events, extract the start time of each transient impact event as the event occurrence time point, and arrange all event occurrence time points in chronological order to form a transient impact event time axis; Calculate the time difference between the occurrence times of adjacent events on the time axis of the transient impact event, generate the original sequence of adjacent event intervals, and perform outlier correction on the original sequence of adjacent event intervals, removing event interval values whose deviation from the interval of adjacent events exceeds the deviation threshold, and obtain the corrected event interval sequence. The corrected event interval sequence is smoothed by a sliding window. The sliding window is formed by selecting the adjacent event interval values before and after the current event interval value as the center. The weighted average of all event interval values in the sliding window is calculated as the smoothed event interval value, and the event interval smoothed sequence is generated. Calculate the natural logarithmic transformation value of each event interval value in the event interval smoothing sequence, and take the event interval sequence after natural logarithmic transformation as the event interval fluctuation sequence. The numerical distribution of the event interval fluctuation sequence exhibits the statistical characteristics of a steady-state random process.
4. The method according to claim 1, characterized in that, The step of associating and mapping the spectral centroid offset trajectories of the event interval fluctuation sequence and the steady-state background voiceprint component sequence to generate a voiceprint state evolution path diagram for each node includes: The steady-state background acoustic text component sequence is processed by frame segmentation to obtain multiple consecutive steady-state background acoustic text analysis frames. Fourier transform is performed on each steady-state background acoustic text analysis frame to generate the steady-state spectrum distribution of each steady-state background acoustic text analysis frame. Calculate the centroid frequency of each steady-state spectral distribution, where the centroid frequency is the centroid frequency value of the spectral energy distribution. Arrange the centroid frequencies corresponding to each steady-state background acoustic text analysis frame in chronological order to generate the original centroid trajectory. The original trajectory of the spectral centroid is differentially calculated to obtain the change in the frequency of the spectral centroid between adjacent steady-state background acoustic text analysis frames, and a sequence of spectral centroid offset directions is generated based on the positive and negative directions of the change. The cumulative summation of the spectral centroid offset direction sequence is calculated to obtain the trajectory of the cumulative spectral centroid offset, which reflects the trend of the total amount of cumulative offset of the spectral centroid over time. The event interval values in the event interval fluctuation sequence are synchronized with the cumulative trajectory of the spectral centroid offset according to the time alignment rule, and the event interval values and the cumulative values of the spectral centroid offset within the same time window are extracted to form two-dimensional coordinate points. Using the cumulative value of the spectral centroid shift as the horizontal axis and the event interval value of the event interval fluctuation sequence as the vertical axis, the two-dimensional coordinate points extracted within each time window are sequentially connected to form a path diagram of the evolution of the voiceprint state.
5. The method according to claim 1, characterized in that, The process of constructing the voiceprint dynamic behavior fingerprint of each node based on the state trajectory of the voiceprint state evolution path graph includes: The two-dimensional coordinate point sequence in the voiceprint state evolution path diagram is used as the original state trajectory. The original state trajectory is reconstructed in phase space. The time delay parameter and the embedding dimension parameter are selected to embed the original state trajectory into the high-dimensional phase space, generating a set of reconstructed trajectory points in the high-dimensional phase space. Calculate the local neighborhood radius of each reconstructed trajectory point in the set of reconstructed trajectory points. The local neighborhood radius is the Euclidean distance between the reconstructed trajectory point and its nearest neighbor trajectory point. Determine the global neighborhood radius based on the statistical distribution of the local neighborhood radii of all reconstructed trajectory points. For each reconstructed trajectory point in the set of reconstructed trajectory points, select all neighboring trajectory points within a spherical neighborhood centered on the reconstructed trajectory point and with the global neighborhood radius as the radius, and construct a local neighborhood point set for the reconstructed trajectory point; Calculate the evolution trend of the number of trajectory points in the local neighborhood set of each reconstructed trajectory point over time, and extract the exponential growth rate of the number of trajectory points in the local neighborhood set during the subsequent evolution process as the trajectory divergence rate parameter. Identify the local regions in phase space that are repeatedly visited by the reconstructed trajectory point set, count the number of times each local region is visited by the reconstructed trajectory points, and mark the local regions whose number of visits exceeds the access threshold as trajectory convergence regions; Extract the boundary coordinates of the orbit convergence region, calculate the geometric center position and region radius parameter of the orbit convergence region based on the boundary coordinates, and use the region radius parameter as the boundary parameter of the orbit convergence region; The original dynamic feature vector is generated by combining the orbit divergence rate parameter and the orbit convergence region boundary parameter, and then the original dynamic feature vector is normalized to generate a normalized dynamic feature vector. Nonlinear dimension reduction projection is performed on the normalized dynamic eigenvectors to project the normalized dynamic eigenvectors in the high-dimensional phase space onto the two-dimensional feature plane, generating the coordinates of the feature projection points on the feature plane. The Euclidean distance between the coordinates of the feature projection point and the origin of the feature plane is calculated as the dynamic feature intensity parameter, and the orientation angle of the coordinates of the feature projection point relative to the horizontal axis of the feature plane is calculated as the dynamic feature direction parameter. The dynamic feature intensity parameter and the dynamic feature direction parameter are combined to generate an acoustic dynamic behavior fingerprint, and the acoustic dynamic behavior fingerprint is associated with and stored with the node identifier of the corresponding node.
6. The method according to claim 1, characterized in that, The comparison of the acoustic dynamic behavior fingerprints of each node under current and historical normal operating conditions, and the determination of abnormal nodes based on the divergence distribution of the fingerprint deviation vector field, includes: Extract the acoustic dynamic behavior fingerprints of each node during the historical normal operation cycle from the database, and construct a historical fingerprint feature library. The historical fingerprint feature library contains the historical values of the dynamic feature intensity parameter and the historical value of the dynamic feature direction parameter of each node during the historical normal operation cycle. Obtain the acoustic dynamic behavior fingerprint of each node within the current monitoring period, and extract the current value of the dynamic feature intensity parameter and the current value of the dynamic feature direction parameter of each node within the current monitoring period; For each node, the difference between the current value and the historical value of the dynamic characteristic intensity parameter is calculated as the intensity deviation value, and the difference between the current value and the historical value of the dynamic characteristic direction parameter is calculated as the direction deviation value. The intensity deviation value and the direction deviation value are combined to generate a two-dimensional deviation vector for each node. The magnitude of the two-dimensional deviation vector is determined by the intensity deviation value, and the direction angle of the two-dimensional deviation vector is determined by the direction deviation value. All the two-dimensional deviation vectors of the nodes are spatially arranged according to the physical location coordinates of the nodes inside the building, and a fingerprint deviation vector field is constructed with the physical location coordinates as the independent variable and the two-dimensional deviation vector as the dependent variable. Divergence calculation is performed on the fingerprint deviation vector field. At the location coordinates of each node, the divergence value of the two-dimensional deviation vector of that node is calculated. The divergence value reflects the degree of divergence or convergence of the two-dimensional deviation vector at that node. Nodes with non-zero divergence values are marked as divergence anomaly nodes. The type of divergence anomaly node is determined by the positive or negative sign of the divergence value. A positive divergence value indicates that the fingerprint deviation vector field at that node exhibits divergence characteristics, while a negative divergence value indicates that the fingerprint deviation vector field at that node exhibits convergence characteristics.
7. The method according to claim 1, characterized in that, The step of constructing a local network topology subgraph based on the abnormal nodes, and determining the propagation direction of the abnormal acoustic signature according to the curl direction of the deviation vector field of the fingerprints of adjacent nodes within the local network topology subgraph, includes: Based on the design drawings of the building's internal water supply network, the physical connection relationship between all nodes is extracted, and a global network topology map is constructed. The nodes of the global network topology map correspond to the deployment positions of the voiceprint sensing units, and the edges correspond to the water supply pipeline connection relationship between adjacent nodes. Select an abnormal node from the set of abnormal nodes as the center node, and perform a breadth-first traversal in the global pipeline topology graph starting from the center node to extract all nodes whose distance from the center node is within a preset number of hops to form a local pipeline topology subgraph. Obtain the two-dimensional deviation vector in the fingerprint deviation vector field of each node in the local pipeline network topology subgraph, and use the two-dimensional deviation vector of each node as the deviation vector field value of that node; For each node in the local network topology subgraph, calculate the difference vector of the deviation vector field values between the node and all its neighboring nodes, and determine the curl value at the node based on the cross product of the difference vector and the direction vector of the connecting edge between the node and its neighboring nodes. Nodes with non-zero curl values are marked as curl-active nodes, and the rotation direction characteristics of the curl-active nodes are determined according to the positive or negative sign of the curl value. A positive curl value indicates that the deviation vector field of the node exhibits a counterclockwise rotation trend, and a negative curl value indicates that the deviation vector field of the node exhibits a clockwise rotation trend. For each connecting edge in the local pipeline topology subgraph, obtain the sign of the curl value of the nodes at both ends of the connecting edge. When the signs of the curl values of the nodes at both ends are the same, determine that the direction of the abnormal propagation of the acoustic pattern on the connecting edge is from the node with the larger absolute value of curl to the node with the smaller absolute value of curl. When the signs of the curl values of the two nodes are different, the rotation direction of the deviation vector field of the node with the positive curl value is compared with the direction of the connecting edge. The direction of abnormal propagation of the acoustic pattern is determined based on the angle between the rotation direction of the deviation vector field and the direction of the connecting edge.
8. The method according to claim 1, characterized in that, The process of combining the abnormal nodes with the propagation direction of the abnormal acoustic signature to obtain the set of candidate nodes for the leakage source includes: All nodes in the abnormal node set are used as the initial node set for reverse tracing. The propagation direction of the abnormal voiceprint propagation direction mapping table is reversed to generate a reverse propagation direction mapping table. The direction identifier of each record in the reverse propagation direction mapping table is opposite to the original propagation direction. Starting from each node in the initial node set, search for a back propagation record with that node as the termination node in the back propagation direction mapping table, extract the corresponding starting node as the previous level backtracking node, repeat the search operation until there is no new previous level backtracking node, and record all backtracking paths. The nodes of each backtracking path are merged to obtain an initial set of candidate nodes for leakage sources. The backtracking frequency parameter of each node in the backtracking path is recorded. A backtracking tree structure with each node in the initial node set as the root node is constructed based on each backtracking path. The branching factor of each node in the backtracking tree structure is calculated. Based on the depth value of each node in the backtracking tree structure, the backtracking frequency parameter, and the branching factor, the initial set of candidate nodes for leakage sources is sorted in ascending order of depth and descending order of backtracking frequency to generate a set of candidate nodes for leakage sources that includes primary and secondary candidate nodes for leakage sources.
9. The method according to claim 1, characterized in that, The step of selecting the node with the largest absolute value of the fingerprint deviation vector field divergence in the candidate node set of leakage sources as the leakage source localization result output includes: Obtain the absolute value of the divergence of the fingerprint deviation vector field corresponding to each candidate node in the candidate node set of leakage sources, and use the absolute value of the divergence as the anomaly intensity score of the candidate node. The abnormal intensity scores of all candidate nodes in the candidate node set of leakage sources are sorted in descending order to generate an abnormal intensity score sorted list, with the candidate node with the highest abnormal intensity score at the top of the sorted list. The candidate node with the highest anomaly intensity score is extracted from the anomaly intensity score ranking list as the primary localization node, and the candidate node with the second highest anomaly intensity score is extracted as the secondary localization node. Obtain the depth value of the main location node in the backtrace tree structure. When the depth value of the main location node is greater than the depth threshold, mark the node chain consisting of the main location node and all its ancestor nodes as a suspected leakage propagation path. Obtain the shortest path distance between the main positioning node and the auxiliary positioning node in the global pipeline topology map. When the shortest path distance is less than the distance threshold, mark all nodes on the path between the main positioning node and the auxiliary positioning node as the node set of the leakage impact range. Calculate the curl value of the fingerprint deviation vector field for each node in the node set of the leakage influence range, extract the node with the largest absolute value of the curl value as the leakage center verification node, and compare the leakage center verification node with the main location node; When the leakage center verification node is the same as the main location node, the main location node is determined as the final leakage source location node. When the leakage center verification node is different from the main location node, the leakage center verification node and the main location node are jointly marked as dual-source leakage candidate nodes. The node with the larger absolute divergence value among the dual-source leakage candidate nodes is selected as the final leakage source location node, and the node identifier corresponding to the final leakage source location node is used as the leakage source location result.
10. A server system, characterized in that, Includes a server, the server being used to perform the method according to any one of claims 1-9.