A Distributed Sound Sensing Method for Detecting Road Voids Based on Active Excitation Sources

Abstract

This invention discloses a distributed acoustic sensing method for detecting road cavities based on an active excitation source. This method, based on Distributed Acoustic Sensing (DAS), first acquires the DAS signal sampling data sequence to be processed; preprocesses the data from each DAS channel; calculates the short-time zero-crossing rate to remove channels with abnormal data records; then, uses sliding window detection to perform spatiotemporal signal localization, extracts the signal envelope within the signal region, and obtains the spatiotemporal distribution information of the signal; finally, it extracts the spatial correlation features and energy distribution features of the signal, and achieves automatic detection of road cavities based on feature-based decision-making. This method fully utilizes the all-weather sensing and high-resolution accuracy features of the DAS array, achieving accurate estimation of road cavity detection with relatively low computational load, enabling large-area sensing, and is suitable for real-time engineering applications.

Description

Technical Field

[0001] This invention belongs to the field of distributed acoustic sensing technology, and particularly relates to a distributed acoustic sensing method for detecting road cavities based on an active excitation source. Background Technology

[0002] With the accelerated pace of urbanization, the monitoring and management of urban road cavities face severe challenges. In the current process of urban space development, the intertwined influence of various municipal engineering projects and the complexity of underground pipeline networks have led to the gradual emergence of defects in underground road structures. These highly concealed, sudden, and destructive potential hazards have become significant sources of risk threatening urban public safety. Existing monitoring methods are mainly implemented in two scenarios: firstly, the emergency response phase after a collapse accident, where the risk has already translated into actual losses, causing not only economic losses but also attracting public attention; secondly, periodic road maintenance and inspection, which, while providing early warnings, suffers from high technical costs and intensive resource investment. This coexistence of passive response and high-cost prevention highlights the urgency and systemic challenges of urban underground space safety management. Due to the strong resistance to electromagnetic interference, good electrical insulation, and light transmission characteristics of optical fibers, distributed sensing systems using optical fibers as sensitive elements and signal transmission media have received increasing attention in recent years for applications such as urban road monitoring and intrusion detection. Therefore, applications such as DAS-based urban road underground space safety monitoring are of great significance.

[0003] The existing methods for detecting cavities in urban roads mainly include: (1) Ground penetrating radar method, which uses radar echoes to image and reflect the underground space structure. It can perceive the underground space structure with high precision, but this method cannot monitor the underground space structure in real time and has a limited detection range. At the same time, the deployment of ground penetrating radar takes a long time and cannot meet the requirements of monitoring the underground space structure of a large number of urban roads; (2) Passive noise imaging method based on DAS, which uses DAS to receive background noise signals such as vehicle vibration and subway driving signals. After long-term accumulation and superposition processing, it can reflect the structure of underground space to a certain extent. However, this method requires long-term recording of background noise signals and is only of reference significance for monitoring underground structures under real-time conditions; (3) Transient surface wave imaging method based on active excitation source of DAS, which is widely used because the signal characteristics of active excitation source are more obvious and can be controlled manually. However, the surface wave imaging method deteriorates sharply under low signal-to-noise ratio conditions and cannot truly reflect the underground space structure information. At the same time, the optical fibers under urban roads are irregularly arranged, and the signals received by each DAS sensor cannot be effectively combined, and the imaging results often have large deviations. Summary of the Invention

[0004] The purpose of this invention is to provide a distributed acoustic sensing method for detecting road cavities based on active excitation sources, so as to solve the technical problems of poor real-time performance and limited detection range in underground space structure monitoring.

[0005] To solve the above-mentioned technical problems, the specific technical solution of the present invention is as follows:

[0006] A distributed acoustic sensing method for detecting road voids based on active excitation sources includes the following steps:

[0007] Step 1: Parameter initialization. Obtain the DAS signal sampling data sequence x(l,n), l=1,2,...,L,n=1,2,...,N, where l is the discrete time index, n is the DAS channel number, L is the number of sampling points, and N is the number of sensor channels.

[0008] Step 2: Preprocess the DAS signal sampled data sequence x(l,n) to be processed;

[0009] Step 3: Calculate the short-time zero-crossing rate Z(n) for all channels and remove DAS channels with abnormal short-time zero-crossing rates;

[0010] Step 4: Extract signal energy from different regions using a sliding window to locate all active source signal regions D. m m = 1, 2, ..., M, where m is the index of the active source signal region and M is the total number of active source signal regions;

[0011] Step 5: For the current active source signal region D m Extract the signal envelope, estimate the signal start time t1 and end time t2, as well as the signal start channel C1 and end channel C2, and determine whether it is a valid active source signal region. If so, estimate the signal source channel n. s Otherwise, proceed directly to step 7;

[0012] Step 6: For the effective active source signal region output in Step 5, calculate the spatial correlation feature r(n), energy distribution kurtosis feature K, and energy symmetry feature S of the active source signal, and determine whether there are holes in each channel based on the joint features.

[0013] Step 7: Determine whether hole detection in all areas is complete. If complete, output all channels with detected holes; otherwise, set m = m + 1 and return to step 5.

[0014] Furthermore, in step 1, the parameters are initialized using the following method to obtain the DAS signal sampling data sequence x(l,n), l=1,2,...,L,n=1,2,...,N, specifically including the following steps:

[0015] Step 1.1: Initialize the sensing parameters for cavity detection, specifically including the initialization of the following parameters:

[0016] sampling frequency f s Initialize to: 2 < f s <8kHz;

[0017] upper limit frequency f of bandpass filter H Initialize to 50 < f H <80Hz;

[0018] lower cutoff frequency f of bandpass filter L Initialize to 1 < f L <15Hz;

[0019] Minimum value of signal amplitude ε x Initialize to 10 -6 <ε x <10 -4 real numbers;

[0020] Short-time zero-crossing rate abnormal lower limit threshold z L Initialize to 0.05 < z L Real numbers less than 0.5;

[0021] Signal region detection sliding window channel window length W c Initialize to 40≤W c Integers ≤60;

[0022] Signal region detection sliding window time window length W t Initialized to: 2f s ≤W t ≤5f s Integers;

[0023] Signal region detection sliding window channel step S c Initialize to 8≤S c Integers ≤ 12;

[0024] Signal region detection sliding window time step S t Initialized to 0.25f s ≤S t ≤f s Integers;

[0025] Signal region energy threshold E S Initialized to 0.01 < E S Real numbers less than 5;

[0026] Spatial correlation feature hard threshold th hard Initialize 0 < th hard Real numbers less than 0.5;

[0027] Energy kurtosis characteristic determination threshold th K Initialize to 0.5 < th K Real numbers less than 0.9;

[0028] Energy symmetry characteristic determination threshold th S Initialize to 0.1 < th S Real numbers less than 0.45.

[0029] Step 1.2: Obtain the DAS signal sampling data sequence x(l,n), l=1,2,...,L,n=1,2,...,N: The real-time acquisition data of L sampling points simultaneously received by N sensors in an array is used as the data sequence x(l,n) to be processed, or the DAS signal data containing L sampling points acquired by N sensors is extracted from the memory as the data sequence to be processed.

[0030] Furthermore, in step 2, the DAS signal sampling data sequence x(l,n) is preprocessed using the following method, specifically including the following steps:

[0031] Step 2.1, Design frequency range is [f L ,f H The bandpass filter is used to filter the DAS signal sampled data sequence x(l,n) to obtain the filtered data x. b (l,n);

[0032] Step 2.2: Filter the data x b (l,n) is normalized to obtain the preprocessed data x. p (l,n):

[0033] x p (l,n)=x b (l,n) / max(x max (n),ε x )

[0034] Where, x max (n) represents all time frames x of channel n. b (l,n), l=1,2,...,L, the maximum absolute value,ε x The minimum value of the signal amplitude initialized in step 1, max(·) means taking the larger of the two values.

[0035] 4. The distributed acoustic sensing method for detecting road voids based on an active excitation source according to claim 1, characterized in that, in step 3, the short-time zero-crossing rate Z(n) of all channels is calculated using the following method, and DAS channels with abnormal data recordings are removed, specifically including the following steps:

[0036] Step 3.1: Calculate the preprocessed data x for each channel n = 1, 2, ..., N. p The short-time zero-crossing rate Z(n) of (l,n):

[0037]

[0038] Where sign(·) is the sign function, defined as:

[0039]

[0040] Where |·| represents the modulo operation;

[0041] Step 3.2: Calculate the upper limit of the short-time zero-crossing rate anomaly Z. H and abnormal lower limit Z L :

[0042]

[0043] Among them, z H and z L These are the upper and lower thresholds for short-time zero-crossing rate anomalies initialized in step 1, respectively.

[0044] Step 3.3: For each channel in n = 1, 2, ..., N, determine:

[0045]

[0046] If P(n) = 1, then retain the data x from that channel. p If (l,n), otherwise discard the data for that channel and obtain the data x after removing the abnormal channels. d (l,n).

[0047] Furthermore, in step 4, the signal energy in each region of the DAS data is extracted using the following method to locate the active source signal region:

[0048] Step 4.1: The signal region detection sliding window traverses the entire data segment to obtain the average energy E(p,q) within each window;

[0049]

[0050] in, For time-domain window index, For channel window indexing, W is the floor function. t The time window length W for detecting the signal region initialized in step 1. c Detect the sliding window channel length for the signal region initialized in step 1; S tThe time step of the sliding window for detecting the signal region initialized in step 1; S c For the signal region initialized in step 1, detect the sliding window channel stepping;

[0051] Step 4.2: Iterate through the energy of each window and select an energy E(p,q) > E. S M windows are used as the active source signal detection area D m m = 1, 2, ..., M, E S The signal region energy threshold initialized in step 1;

[0052] Step 4.3: Initialize the active source signal detection region index m = 1.

[0053] Furthermore, in step 5, the signal envelope within the current region is extracted using the following method, and the signal start time t1 and end time t2, as well as the signal start channel C1 and end channel C2, and the signal source channel n are estimated. s Specifically, it includes the following steps:

[0054] Step 5.1: Extract the current signal detection region D using smoothing filtering. m Data x after internal preprocessing d The signal envelopes of each channel in (l,n);

[0055] Step 5.2: Perform peak detection on the extracted signal envelopes of each channel, estimate the signal start time t1, end time t2, start channel C1, and end channel C2, and calculate the signal start time frame. Termination Frame Store the signal time range in Ω t , i.e. Ω t = [L1, L2], the signal channel range is stored in Ω c , i.e. Ω c = [C1, C2];

[0056] Step 5.3: Determine whether L1≠L2 and C1≠C2 are true. If they are true, proceed to step 5.4; otherwise, proceed to step 7.

[0057] Step 5.4: Based on the signal spatiotemporal detection results, estimate the signal source channel n. s :

[0058]

[0059] Furthermore, in step 6, based on the spatiotemporal location of the signal extracted in step 5, the spatial correlation feature r(n), energy distribution kurtosis feature K, and energy symmetry feature S of the active source signal are extracted using the following method, specifically including the following steps:

[0060] Step 6.1: Calculate the spatial correlation characteristic r(n) of the active source signal:

[0061]

[0062] in, Indicates that channel n is in Ω t The time-domain average amplitude within the range is calculated as follows:

[0063]

[0064] Step 6.2, Extract Ω c Kurtosis characteristics K of energy distribution within the range:

[0065]

[0066] in, The kurtosis of the energy distribution is calculated as follows:

[0067]

[0068] Among them, E ave (n) represents the time-domain average energy of the active source signal in channel n. The average energy of the active source signal is calculated as follows:

[0069]

[0070] Step 6.3: Extract signal source channel n s The energy distribution symmetry characteristic S in the vicinity:

[0071]

[0072] Where ξ is the mean square error about the symmetry of the energy distribution, calculated as:

[0073]

[0074] Step 6.4, for the signal channel range Ω c For each channel in the array, determine:

[0075]

[0076] Among them, th soft The soft threshold for spatial correlation features is calculated as Ω. c The average value of r(n) for all spatial correlation characteristic valley channels, th hard The hard threshold for the spatial correlation feature initialized in step 1, th K The threshold for determining the energy kurtosis feature initialized in step 1, th SThe threshold for determining the energy symmetry feature initialized in step 1 is min(·), which means taking the smaller of the two values. If H(n) = 1, then it is determined that there is a hole in channel n and it is recorded; otherwise, it is not recorded.

[0077] Furthermore, in step 7, the following method is used to determine whether the data processing is complete, specifically including the following steps: Step 7.1, determine whether the active source signal region index m = M is true. If it is true, the DAS data processing ends and all channels with holes are output. Otherwise, set the active source signal region index m = m + 1 and return to step 5.

[0078] The distributed fiber optic acoustic sensing method for detecting road cavities based on an active excitation source, as proposed in this invention, has the following advantages:

[0079] 1. This invention utilizes the characteristic that the main frequency component of the signal deviates when the DAS channel is abnormally coupled or experiences continuous abnormal fluctuations. Through short-time zero-crossing rate analysis, abnormal DAS channels are effectively identified and eliminated, as shown in step 3, thereby improving the robustness of road cavity detection.

[0080] 2. This invention uses sliding window energy detection to initially locate the signal region, and performs spatiotemporal localization of the active source signal only within the selected target region, as shown in steps 4 and 5. This avoids the extraction of signal parameters across the entire spatiotemporal range, significantly reduces the computational complexity of the method, and improves the application efficiency of the method in large-scale urban road cavity detection.

[0081] 3. This invention makes full use of the propagation characteristics of active source signals in underground space. As shown in the corresponding processing steps in step 6, the extracted feature parameters have good distinguishability between road cavity structures and non-cavity structures. They can still reflect the structural characteristics of underground space well even at low signal-to-noise ratios. Furthermore, by combining spatial correlation features and energy distribution features, the robustness and accuracy of road cavity detection are further improved. Attached Figure Description

[0082] Figure 1 This is a schematic flowchart of the method of the present invention;

[0083] Figure 2 This is a partial channel time-domain waveform diagram of the DAS signal sampling data sequence in Example 1;

[0084] Figure 3 This is a local channel time-domain waveform diagram of the preprocessed signal in Example 1;

[0085] Figure 4 This is a short-time zero-crossing rate detection graph for the entire channel in Example 1;

[0086] Figure 5This is a local channel signal energy distribution diagram of Example 1;

[0087] Figure 6 The active source signal time-domain waveform diagram of Example 1;

[0088] Figure 7 This is a spatial correlation feature map of the active source signal in Example 1;

[0089] Figure 8 This is a partial channel time-domain waveform diagram of the DAS signal sampling data sequence in Example 2;

[0090] Figure 9 This is a partial channel time-domain waveform diagram of the preprocessed signal in Example 2;

[0091] Figure 10 This is a short-time zero-crossing rate detection diagram for the entire channel in Example 2;

[0092] Figure 11 This is a local channel signal energy distribution diagram of Example 2;

[0093] Figure 12 The active source signal time-domain waveform diagram for Example 2;

[0094] Figure 13 This is a spatial correlation feature map of the active source signal in Example 2. Detailed Implementation

[0095] To better understand the purpose, structure, and function of this invention, the following detailed description of a distributed acoustic sensing method for detecting road cavities based on an active excitation source is provided in conjunction with the accompanying drawings.

[0096] like Figure 1 As shown, this invention proposes a distributed acoustic sensing method for detecting road cavities based on active excitation sources, comprising the following steps:

[0097] Step 1: Parameter initialization. Obtain the DAS signal sampling data sequence x(l,n), l=1,2,...,L,n=1,2,...,N, where l is the discrete time index, n is the DAS channel number, L is the number of sampling points, and N is the number of sensor channels.

[0098] Step 2: Preprocess the DAS signal sampled data sequence x(l,n) to be processed;

[0099] Step 3: Calculate the short-time zero-crossing rate Z(n) for all channels and remove DAS channels with abnormal short-time zero-crossing rates;

[0100] Step 4: Extract signal energy from different regions using a sliding window to locate all active source signal regions D. mm = 1, 2, ..., M, where m is the index of the active source signal region and M is the total number of active source signal regions;

[0101] Step 5: For the current active source signal region D m Extract the signal envelope, estimate the signal start time t1 and end time t2, as well as the signal start channel C1 and end channel C2, and determine whether it is a valid active source signal region. If so, estimate the signal source channel n. s Otherwise, proceed directly to step 7;

[0102] Step 6: For the effective active source signal region output in Step 5, calculate the spatial correlation feature r(n), energy distribution kurtosis feature K, and energy symmetry feature S of the active source signal, and determine whether there are holes in each channel based on the joint features.

[0103] Step 7: Determine whether hole detection in all areas is complete. If complete, output all channels with detected holes; otherwise, set m = m + 1 and return to step 5.

[0104] Furthermore, in step 1, the parameters are initialized using the following method to obtain the DAS signal sampling data sequence x(l,n), l=1,2,...,L,n=1,2,...,N, specifically including the following steps:

[0105] Step 1.1: Initialize the sensing parameters for cavity detection, specifically including the initialization of the following parameters:

[0106] sampling frequency f s Initialize to: 2 < f s <8kHz;

[0107] upper limit frequency f of bandpass filter H Initialize to 50 < f H <80Hz;

[0108] lower cutoff frequency f of bandpass filter L Initialize to 1 < f L <15Hz;

[0109] Minimum value of signal amplitude ε x Initialize to 10 -6 <ε x <10 -4 real numbers;

[0110] Short-time zero-crossing rate abnormal lower limit threshold z L Initialize to 0.05 < z L Real numbers less than 0.5;

[0111] Signal region detection sliding window channel window length W cInitialize to 40≤W c Integers ≤60;

[0112] Signal region detection sliding window time window length W t Initialized to: 2f s ≤W t ≤5f s Integers;

[0113] Signal region detection sliding window channel step S c Initialize to 8≤S c Integers ≤ 12;

[0114] Signal region detection sliding window time step S t Initialized to 0.25f s ≤S t ≤f s Integers;

[0115] Signal region energy threshold E S Initialized to 0.01 < E S Real numbers less than 5;

[0116] Spatial correlation feature hard threshold th hard Initialize 0 < th hard Real numbers less than 0.5;

[0117] Energy kurtosis characteristic determination threshold th K Initialize to 0.5 < th K Real numbers less than 0.9;

[0118] Energy symmetry characteristic determination threshold th S Initialize to 0.1 < th S Real numbers less than 0.45.

[0119] Step 1.2: Obtain the DAS signal sampling data sequence x(l,n), l=1,2,...,L,n=1,2,...,N: The real-time acquisition data of L sampling points simultaneously received by N sensors in an array is used as the data sequence x(l,n) to be processed, or the DAS signal data containing L sampling points acquired by N sensors is extracted from the memory as the data sequence to be processed.

[0120] In step 1, to balance the computational complexity and estimation accuracy of this invention, the sampling frequency f is... s The preferred value is 4kHz, and the upper limit frequency f of the bandpass filter is... H The preferred value is 60Hz, and the lower cutoff frequency f of the bandpass filter is... L The preferred value is 5Hz, and the minimum signal amplitude ε x The preferred value is 10 -5Short-term zero-crossing rate anomaly upper limit threshold z H The preferred value is 2, and the short-time zero-crossing rate anomaly lower limit threshold z L The preferred value is 0.5, and the channel size W of the signal region detection sliding window is... c The preferred value is 50, and the sliding window time size W for signal region detection is... t The preferred value is 12, and the step window length S of the sliding window channel for signal region detection is... c The preferred value is 20, and the time step S of the sliding window for signal region detection is... t The preferred value is 3, and the signal region energy threshold E S The preferred value is 0.5, and the hard threshold for spatial correlation features is th. hard The preferred value is 0.3, and the threshold for determining the energy kurtosis feature is th. K The preferred value is 0.6, and the energy symmetry characteristic determination threshold is th. S The preferred value is 0.3.

[0121] Furthermore, in step 2, the DAS signal sampling data sequence x(l,n) is preprocessed using the following method, specifically including the following steps:

[0122] Step 2.1, Design frequency range is [f L ,f H The bandpass filter is used to filter the DAS signal sampled data sequence x(l,n) to obtain the filtered data x. b (l,n);

[0123] Step 2.2: Filter the data x b (l,n) is normalized to obtain the preprocessed data x. p (l,n):

[0124] x p (l,n)=x b (l,n) / max(x max (n),ε x )

[0125] Where, x max (n) represents all time frames x of channel n. b (l,n), l=1,2,...,L, the maximum absolute value,ε x The minimum value of the signal amplitude initialized in step 1, max(·) means taking the larger of the two values.

[0126] 4. The distributed acoustic sensing method for detecting road voids based on an active excitation source according to claim 1, characterized in that, in step 3, the short-time zero-crossing rate Z(n) of all channels is calculated using the following method, and DAS channels with abnormal data recordings are removed, specifically including the following steps:

[0127] Step 3.1: Calculate the preprocessed data x for each channel n = 1, 2, ..., N. p The short-time zero-crossing rate Z(n) of (l,n):

[0128]

[0129] Where sign(·) is the sign function, defined as:

[0130]

[0131] Where |·| represents the modulo operation;

[0132] Step 3.2: Calculate the upper limit of the short-time zero-crossing rate anomaly Z. H and abnormal lower limit Z L :

[0133]

[0134] Among them, z H and z L These are the upper and lower thresholds for short-time zero-crossing rate anomalies initialized in step 1, respectively.

[0135] Step 3.3: For each channel in n = 1, 2, ..., N, determine:

[0136]

[0137] If P(n) = 1, then retain the data x from that channel. p If (l,n), otherwise discard the data for that channel and obtain the data x after removing the abnormal channels. d (l,n).

[0138] Furthermore, in step 4, the signal energy in each region of the DAS data is extracted using the following method to locate the active source signal region:

[0139] Step 4.1: The signal region detection sliding window traverses the entire data segment to obtain the average energy E(p,q) within each window;

[0140]

[0141] in, For time-domain window index, For channel window indexing, W is the floor function. t The time window length W for detecting the signal region initialized in step 1. c Detect the sliding window channel length for the signal region initialized in step 1; S t The time step of the sliding window for detecting the signal region initialized in step 1; S c For the signal region initialized in step 1, detect the sliding window channel stepping;

[0142] Step 4.2: Iterate through the energy of each window and select an energy E(p,q) > E. S M windows are used as the active source signal detection area D m m = 1, 2, ..., M, E S The signal region energy threshold initialized in step 1;

[0143] Step 4.3: Initialize the active source signal detection region index m = 1.

[0144] Furthermore, in step 5, the signal envelope within the current region is extracted using the following method, and the signal start time t1 and end time t2, as well as the signal start channel C1 and end channel C2, and the signal source channel n are estimated. s Specifically, it includes the following steps:

[0145] Step 5.1: Extract the current signal detection region D using smoothing filtering. m Data x after internal preprocessing d The signal envelopes of each channel in (l,n);

[0146] Step 5.2: Perform peak detection on the extracted signal envelopes of each channel, estimate the signal start time t1, end time t2, start channel C1, and end channel C2, and calculate the signal start time frame. Termination Frame Store the signal time range in Ω t , i.e. Ω t = [L1, L2], the signal channel range is stored in Ω c , i.e. Ω c = [C1, C2];

[0147] Step 5.3: Determine whether L1≠L2 and C1≠C2 are true. If they are true, proceed to step 5.4; otherwise, proceed to step 7.

[0148] Step 5.4: Based on the signal spatiotemporal detection results, estimate the signal source channel n. s :

[0149]

[0150] Furthermore, in step 6, based on the spatiotemporal location of the signal extracted in step 5, the spatial correlation feature r(n), energy distribution kurtosis feature K, and energy symmetry feature S of the active source signal are extracted using the following method, specifically including the following steps:

[0151] Step 6.1: Calculate the spatial correlation characteristic r(n) of the active source signal:

[0152]

[0153] in, Indicates that channel n is in Ω t The time-domain average amplitude within the range is calculated as follows:

[0154]

[0155] Step 6.2, Extract Ω c Kurtosis characteristics K of energy distribution within the range:

[0156]

[0157] in, The kurtosis of the energy distribution is calculated as follows:

[0158]

[0159] Among them, E ave (n) represents the time-domain average energy of the active source signal in channel n. The average energy of the active source signal is calculated as follows:

[0160]

[0161] Step 6.3: Extract signal source channel n s The energy distribution symmetry characteristic S in the vicinity:

[0162]

[0163] Where ξ is the mean square error about the symmetry of the energy distribution, calculated as:

[0164]

[0165] Step 6.4, for the signal channel range Ω c For each channel in the array, determine:

[0166]

[0167] Among them, th soft The soft threshold for spatial correlation features is calculated as Ω. cThe average value of r(n) for all spatial correlation characteristic valley channels, th hard The hard threshold for the spatial correlation feature initialized in step 1, th K The threshold for determining the energy kurtosis feature initialized in step 1, th S The threshold for determining the energy symmetry feature initialized in step 1 is min(·), which means taking the smaller of the two values. If H(n) = 1, then it is determined that there is a hole in channel n and it is recorded; otherwise, it is not recorded.

[0168] Furthermore, in step 7, the following method is used to determine whether the data processing is complete, specifically including the following steps: Step 7.1, determine whether the active source signal region index m = M is true. If it is true, the DAS data processing ends and all channels with holes are output. Otherwise, set the active source signal region index m = m + 1 and return to step 5.

[0169] Example 1:

[0170] Field verification of DAS data collection and void detection was conducted on Youhu Road, Jiangning District, Nanjing City.

[0171] The following section of DAS data, consisting of 1200 DAS channels and lasting 18 seconds, is selected for hole detection:

[0172] Based on step 1, set the sampling frequency f. s =4kHz, upper limit frequency f of bandpass filter H =60Hz, lower cutoff frequency of bandpass filter f L =5Hz, minimum signal amplitude ε x =10 -5 Short-term zero-crossing rate anomaly upper limit threshold z H =2, short-time zero-crossing rate abnormal lower limit threshold z L =0.5, the size of the sliding window channel W for signal region detection c =50, the sliding window time size W for signal region detection t =2f s Signal region detection sliding window channel step window length S c =20, signal region detection sliding window time step S t =0.5f s Signal region energy threshold E S =0.1, hard threshold for spatial correlation features th hard =0.4, threshold for determining energy kurtosis characteristics. K =0.6, threshold for determining energy symmetry characteristics. S =0.3, the local channel time-domain waveform of the read DAS signal sampling data sequence is as follows Figure 2 As shown;

[0173] Based on step 2, the designed frequency range is [f L ,f H The bandpass filter is used to filter the DAS signal sampled data sequence x(l,n), and the filtered data x is then processed. b (l,n) is normalized to obtain the preprocessed data x. p (l,n) such as Figure 3 As shown;

[0174] Based on step 3, calculate the short-time zero-crossing rate Z(n) for all channels. The calculation results are as follows: Figure 4 As shown, removing short-term zero-crossing rates greater than the upper limit Z of anomalies H And less than the abnormal lower limit Z L From the channel data, obtain the updated data x. d (l,n);

[0175] Based on step 4, the signal energy of each region of the DAS data is extracted using a sliding window. The extracted local signal energy results are as follows: Figure 5 As shown, one active source signal region D1 is located. Figure 6 As shown, its channel range is [770, 820], and its time range is [1s, 3s].

[0176] Based on step 5, extract the signal envelope of the signal detection region, and estimate the signal start time t1 = 1.7s, the end time t2 = 2.0s, the signal start channel C1 = 720, the end channel C2 = 770, and the signal source channel n. s =745;

[0177] Based on step 6, the spatial correlation feature r(n), energy distribution kurtosis feature K, and energy symmetry feature S of the active source signal are extracted, and it is determined whether there are holes in the channel. The extracted spatial correlation feature r(n) of the active source signal is as follows: Figure 7 As shown, the soft threshold th for spatial correlation features is calculated. soft =0.53, extract Ω c The energy distribution kurtosis characteristic within the range is K = 0.23, and the signal source channel n s The energy distribution in the vicinity is symmetrical, with a characteristic S = 0.08, for the signal channel range Ω. c Perform a void check on each channel in the process;

[0178] Based on step 7, determine whether the data processing is complete. If the processing is complete, output all channels with detected voids; otherwise, set the active source signal region index m = m + 1 and return to step 5.

[0179] The data segment contains one road cavity in the experimental section, located in lane 799. This method successfully detected the lane where the cavity is located.

[0180] Example 2:

[0181] DAS data collection and field verification of the absence of void features were conducted on Daishan Road, Yuhuatai District, Nanjing City.

[0182] The following section of DAS data, consisting of 400 DAS channels and lasting 18 seconds, is selected for hole detection:

[0183] Based on step 1, set the sampling frequency f. s =4kHz, upper limit frequency f of bandpass filter H =60Hz, lower cutoff frequency of bandpass filter f L =5Hz, minimum signal amplitude ε x =10 -5 Short-term zero-crossing rate anomaly upper limit threshold z H =2, short-time zero-crossing rate abnormal lower limit threshold z L =0.5, the size of the sliding window channel W for signal region detection c =30, the sliding window time size W for signal region detection t =2f s Signal region detection sliding window channel step window length S c =10, signal region detection sliding window time step S t =0.5f s Signal region energy threshold E S =2.5, hard threshold for spatial correlation features th hard =0.5, threshold for determining energy kurtosis characteristics. K =0.6, threshold for determining energy symmetry characteristics. S =0.3, the local channel time-domain waveform of the read DAS signal sampling data sequence is as follows Figure 8 As shown;

[0184] Based on step 2, the designed frequency range is [f L ,f H The bandpass filter is used to filter the DAS signal sampled data sequence x(l,n), and the filtered data x is then processed. b (l,n) is normalized to obtain the preprocessed data x. p (l,n) such as Figure 9 As shown;

[0185] Based on step 3, calculate the short-time zero-crossing rate Z(n) for all channels. The calculation results are as follows: Figure 4 As shown, removing short-term zero-crossing rates greater than the upper limit Z of anomaliesH And less than the abnormal lower limit Z L From the channel data, obtain the updated data x. d (l,n);

[0186] Based on step 4, the signal energy of each region of the DAS data is extracted using a sliding window. The extracted local signal energy is as follows: Figure 11 As shown, one active source signal region D1 is located. Figure 12 As shown, its channel range is [240, 270], and its time range is [8s, 10s].

[0187] Based on step 5, extract the signal envelope of the signal detection region, and estimate the signal start time t1 = 8.5s, end time t2 = 8.7s, signal start channel C1 = 240, end channel C2 = 270, and signal source channel n. s =255;

[0188] Based on step 6, the spatial correlation feature r(n), energy distribution kurtosis feature K, and energy symmetry feature S of the active source signal are extracted, and it is determined whether there are holes in the channel. The extracted spatial correlation feature r(n) of the active source signal is as follows: Figure 7 As shown, the soft threshold th for spatial correlation features is calculated. soft =0.64, extract Ω c The energy distribution kurtosis characteristic within the range is K = 0.89, and the signal source channel n s The energy distribution in the vicinity is symmetrical, with a characteristic value of S = 0.77, for the signal channel range Ω. c Perform a void check on each channel in the process;

[0189] Based on step 7, determine whether the data processing is complete. If the processing is complete, output all channels with detected voids; otherwise, set the active source signal region index m = m + 1 and return to step 5.

[0190] This experiment was conducted on a closed road and selected a section of road without cavities for feature verification. The method did not detect any cavities in this area.

[0191] It is understood that the present invention has been described through some embodiments, and those skilled in the art will recognize that various changes or equivalent substitutions can be made to these features and embodiments without departing from the spirit and scope of the invention. Furthermore, under the teachings of the present invention, these features and embodiments can be modified to adapt to specific situations and materials without departing from the spirit and scope of the invention. Therefore, the present invention is not limited to the specific embodiments disclosed herein, and all embodiments falling within the scope of the claims of this application are within the protection scope of the present invention.

Claims

1. A distributed acoustic sensing method for detecting road voids based on active excitation sources, characterized in that, Includes the following steps: Step 1: Parameter initialization, obtain the DAS signal sampling data sequence to be processed. l is the discrete-time index, n is the DAS channel number, L is the number of sampling points, and N is the number of sensor channels; Step 2: Sample data sequence of the DAS signal to be processed Preprocessing is required; Step 3: Calculate the short-time zero-crossing rate for all channels. Remove DAS channels with short-term zero-crossing rate anomalies; Step 4: Extract signal energy from different regions using a sliding window to locate all active source signal regions. , m is the active source signal region index, and M is the total number of active source signal regions; Step 5: For the current active source signal region Extract the signal envelope and estimate the signal start time. and termination time and signal start channel and termination channel Determine whether it is a valid active source signal region; if so, estimate the signal source channel. Otherwise, proceed directly to step 7; Step 6: Calculate the spatial correlation characteristics of the active source signals for the effective active source signal region output in Step 5. Energy distribution kurtosis characteristics K and energy symmetry characteristics Based on joint features, determine whether there are voids in each channel; Step 7: Determine if hole detection in all regions is complete. If complete, output all channels with detected holes; otherwise, let... Return to step 5; In step 4, the signal energy in each region of the DAS data is extracted using the following method to locate the active source signal region: Step 4.1: The signal region detection sliding window traverses the entire data segment to obtain the average energy within each window. ; ； in, For time-domain window index, For channel window indexing, This is the floor function. The time window length for detecting the signal region initialized in step 1. Detect the sliding window channel length for the signal region initialized in step 1; The time step of the sliding window is detected for the signal region initialized in step 1; For the signal region detection sliding window channel stepping initialized in step 1, This represents the data after removing abnormal channels; Step 4.2: Traverse the energy of each window and select the energy. of Each window serves as the active source signal detection area. , The signal region energy threshold initialized in step 1; Step 4.3: Initialize the active source signal detection region index m = 1; In step 5, the signal envelope within the current region is extracted using the following method, and the signal start time is estimated. and termination time and signal start channel and termination channel Signal source channel Specifically, it includes the following steps: Step 5.1: Extract the current signal detection region using smoothing filtering. Data after internal preprocessing The signal envelopes of each channel; Step 5.2: Perform peak detection on the extracted signal envelopes of each channel to estimate the signal start time. End time Starting Channel and termination channel Calculate the start time frame of the signal Termination frame And store the signal time range ,Right now Signal channel range stored ,Right now ; Step 5.3, Judgment Check if the condition is met. If it is met, proceed to step 5.4; otherwise, proceed to step 7. Step 5.4: Estimate the signal source channel based on the signal spatiotemporal detection results. : 。 2. The distributed acoustic sensing method for detecting road cavities based on active excitation sources according to claim 1, characterized in that, In step 1, the parameters are initialized using the following method to obtain the DAS signal sampling data sequence to be processed. Specifically, it includes the following steps: Step 1.1: Initialize the sensing parameters for cavity detection, specifically including the initialization of the following parameters: Sampling frequency Initialize to: kHz; upper limit frequency of bandpass filter Initialize to Hz; lower cutoff frequency of bandpass filter Initialize to Hz; Minimum value of signal amplitude Initialize to real numbers; Short-time zero-crossing rate anomaly upper limit threshold Initialize to real numbers; Short-time zero-crossing rate abnormal lower limit threshold Initialize to real numbers; Signal region detection sliding window channel window length Initialize to Integers; Signal region detection sliding window time window length Initialize to: Integers; Signal region detection sliding window channel stepping Initialize to Integers; Signal Region Detection Sliding Window Time Step Initialize to Integers; Signal region energy threshold Initialize to real numbers; Spatial correlation feature hard threshold initialization real numbers; Energy kurtosis characteristic threshold Initialize to real numbers; Energy symmetry characteristic determination threshold Initialize to real numbers; Step 1.2: Obtain the DAS signal sampling data sequence to be processed. The real-time data acquired simultaneously from L sampling points by N sensors in an array is used as the data sequence to be processed. Alternatively, DAS signal data containing L sampling points collected by N sensors can be extracted from the memory as a data sequence to be processed.

3. The distributed acoustic sensing method for detecting road cavities based on active excitation sources according to claim 2, characterized in that, In step 2, the DAS signal data sequence is sampled using the following method. Preprocessing includes the following steps: Step 2.1, Design frequency range is A bandpass filter for sampling DAS signal data sequences Perform filtering to obtain filtered data. ; Step 2.2: Filter the data. Amplitude normalization is performed to obtain the preprocessed data. : ； in, Represents all time frames of channel n. Maximum absolute value The minimum value of the signal amplitude initialized in step 1, This indicates taking the larger of the two values.

4. The distributed acoustic sensing method for detecting road cavities based on active excitation sources according to claim 1, characterized in that, In step 3, the short-time zero-crossing rate of all channels is calculated using the following method. And remove DAS channels with abnormal data records, specifically including the following steps: Step 3.1: Calculate the preprocessed data for each channel n = 1, 2, ..., N. short-time zero crossing rate : ； Among them, sign For symbolic functions, it is defined as: ； in, For modulo operations; Step 3.2: Calculate the upper limit of the short-time zero-crossing rate anomaly. and abnormal lower limit : ， ； in, and These are the upper and lower thresholds for short-time zero-crossing rate anomalies initialized in step 1, respectively. Step 3.3, for For each channel in the array, determine: ； like If so, then retain the data from that channel. Otherwise, discard the data from that channel and obtain the data after removing the abnormal channels. .

5. The distributed acoustic sensing method for detecting road cavities based on active excitation sources according to claim 1, characterized in that, In step 6, based on the spatiotemporal location of the signal extracted in step 5, the spatial correlation features of the active source signal are extracted using the following method. Energy distribution kurtosis characteristics K and energy symmetry characteristics Specifically, it includes the following steps: Step 6.1: Calculate the spatial correlation characteristics of the active source signal. : ； in, Indicates channel exist The time-domain average amplitude within the range is calculated as follows: ； Step 6.2, Extraction Kurtosis characteristics K of energy distribution within the range: ； in, The kurtosis of the energy distribution is calculated as follows: ； in, Indicates channel The time-domain average energy of the active source signal, The average energy of the active source signal is calculated as follows: ， ； Step 6.3: Extract the signal source channel Symmetrical characteristics of energy distribution in the vicinity : ； Where ξ is the mean square error about the symmetry of the energy distribution, calculated as: ； Step 6.4, for the signal channel range For each channel in the array, determine: ； in, The soft threshold for spatial correlation features is calculated as follows: Valley channels of all spatial correlation features The average value, The hard threshold for the spatial correlation feature initialized in step 1, The threshold for determining the energy kurtosis feature initialized in step 1. The threshold for determining the energy symmetry feature initialized in step 1, min This indicates taking the smaller of the two values; if H(n) = 1, then the channel is determined. If a void exists, it is recorded; otherwise, it is not recorded.

6. The distributed acoustic sensing method for detecting road cavities based on active excitation sources according to claim 1, characterized in that, In step 7, the following method is used to determine whether the data processing is complete, specifically including the following steps: Step 7.1: Determine the active source signal region index If the condition is met, the DAS data processing ends, and all channels with detected holes are output; otherwise, the active source signal region index is set. Return to step 5.

Images