A multi-source heterogeneous sensor time-frequency feature extraction and fault diagnosis method and system

Abstract

The application discloses a multi-source heterogeneous sensor time-frequency feature extraction and fault diagnosis method and system, wherein the method comprises the following steps: S1, multi-source heterogeneous sensor data acquisition and high-precision space-time alignment, implementing MAC layer hardware level time synchronization and dynamic boundary constraint frequency mapping; S2, multi-dimensional time-frequency feature space construction and fusion, extracting time domain high-order dimensionless statistics and frequency domain optimal base wavelet energy, and constructing multi-order time-frequency collaborative damage sensitivity index; S3, space topology clustering and self-diagnosis stripping logic, constructing a single-point-space joint heterogeneous source discriminant function, three-branch joint analysis, and outputting an analysis result; and S4, flow data evolution anomaly qualitative classification based on a robust random split forest. The application also discloses a multi-source heterogeneous sensor time-frequency feature extraction and fault diagnosis system. The application can capture early hidden damage in advance, and accurately distinguish sensor hardware faults, extreme environmental events and real structure damage.

Description

Technical Field

[0001] This invention relates to the field of predictive maintenance technology for industrial internet and intelligent equipment, and in particular to a method and system for extracting time-frequency features and diagnosing faults in multi-source heterogeneous sensors. Background Technology

[0002] In recent years, with the deepening of the national strategies of building a strong transportation network and the energy internet, the scale and number of large-scale infrastructure projects (such as long-span bridges, towering cable towers, and power transmission towers) have grown rapidly. These infrastructure projects are key nodes in the national economic operation, and their structural safety directly affects public safety and the stability of energy supply. To monitor the health status of these critical structures in real time, dense sensor networks are typically deployed in key structural areas (such as cable anchorage zones, piers, and tower foundations), encompassing vibration sensors (to acquire modal parameters), strain gauges (to sense loads), displacement gauges (to monitor deformation), and meteorological and temperature sensors (to monitor environmental impacts). However, for the long-term operational monitoring of large-scale infrastructure, existing technical solutions face significant technical challenges when processing cross-regional, large-scale, multi-source, heterogeneous data. First, there is the problem of "spatiotemporal asynchrony" across regions. Infrastructure monitoring points are widely distributed, making cabling extremely difficult. Distributed wireless data acquisition is often employed, but due to limitations in network transmission and the processing capabilities of individual terminals, data streams are frequently accompanied by significant asynchrony. In large-scale structural analysis, without nanosecond-level time alignment, minute modal parameter shifts (such as changes in natural frequency caused by cable tension variations) can be masked by data synchronization errors, leading to "false conclusions" in structural safety assessments.

[0003] Second, there's the "false alarm storm" caused by environmental interference. Large infrastructures are situated in completely open natural environments, making them highly susceptible to extreme weather (strong winds, heavy rain), ambient temperature drift, and unrelated transient interference (passing heavy vehicles, pedestrians, bird strikes). Traditional threshold alarm systems struggle to distinguish these environmental responses from minor structural damage. Especially in strong winds, sensors often experience signal anomalies due to probe installation vibrations or loose cables. If the system cannot effectively differentiate between "sensor malfunctions" and "structural damage," it will frequently generate invalid warnings, reducing maintenance personnel's trust in the monitoring system.

[0004] Third, the insidious nature of structural damage evolution and the difficulty of long-term monitoring. The damage process of large structures during their service life is often long-term and gradual. Traditional monitoring models mostly rely on short-term static monitoring and lack the ability to learn dynamically from convective data over a long period of time. How to accurately identify the characteristic deviations caused by early structural damage in complex noise backgrounds and distinguish them from random environmental disturbances is a key challenge in ensuring the long-term safety of infrastructure.

[0005] Therefore, this invention urgently needs to design a method for extracting time-frequency features and diagnosing faults in multi-source heterogeneous sensors to solve the above-mentioned technical problems. Summary of the Invention

[0006] The purpose of this invention is to address the technical problems existing in monitoring networks for large-scale infrastructure (such as long-span bridges and tall towers) including spatiotemporal asynchrony of heterogeneous data, coupling and aliasing of strong environmental noise and physical micro-damage features, and the difficulty in accurately separating sensor hardware failures from actual structural damage. This invention provides a method and system for extracting time-frequency features and diagnosing faults from multi-source heterogeneous sensors. It achieves absolute spatiotemporal benchmark alignment of multimodal data by introducing nanosecond-level clock dynamic drift compensation and zero-phase polynomial resampling technology; constructs a noise-resistant high-dimensional feature space that integrates high-order dimensionless statistics and optimal base wavelet packet frequency band energy; and innovatively proposes a multi-branch self-diagnostic logic structure based on nonlinear spatial topological clustering (Ward-Spearman) to accurately isolate false alarms caused by extreme environments and hardware failures; finally, it combines a robust flow analysis model to achieve long-term evolution assessment and qualitative classification of structural states.

[0007] The technical solution of this invention is: The present invention provides a method for time-frequency feature extraction and fault diagnosis of multi-source heterogeneous sensors, comprising: S1: Acquire multi-source heterogeneous sensor signals from infrastructure deployment; perform hardware-level time synchronization and unified timing on the signals, and dynamically compensate for clock drift rate; then fuse signals with different sampling rates by using zero-phase frequency domain digital filtering downsampling and spline interpolation upsampling with natural boundary constraints, and use adaptive empirical mode decomposition with extreme value mirror extension to remove slowly varying trend terms; output the reconstructed multi-source clean signal sequence. S2: Based on the multi-source clean signal sequence, the signal sequence is synchronously truncated through a non-overlapping time sliding window, and damage-sensitive features are extracted from both the time and frequency domains: in the time domain, high-order dimensionless statistics sensitive to mutations are extracted, and in the frequency domain, the spectral centroid frequency and wavelet packet energy features are extracted; the time-domain and frequency-domain features are spliced ​​and fused to construct a multi-dimensional fused feature space; based on the abnormal quantities in each dimension of the feature space, a multi-order time-frequency collaborative damage sensitivity index with a multiplicative collaborative excitation mechanism is constructed, and the multi-dimensional abnormal signals are fused into a single scalar through multiplicative coupling; the current node status is determined based on the comparison result of the sensitivity index and a preset threshold, and a normal node or suspicious node label is output; S3: Based on nodes marked as suspicious, construct a single-point-space joint heterogeneous source discrimination function containing three components: single-point anomaly exceeding the standard multiple, local spatial coupling collapse degree, and global spatial synchronization degree. Execute three-branch hard judgment logic according to the combination mode of the three component values, and output independent hardware fault isolation instructions, local real damage location warnings, or global extreme value environmental event silent monitoring instructions respectively.

[0008] Furthermore, it also includes step S4: using the high-dimensional fusion feature vector generated in step S2 as real-time feature points, inserting them into each tree structure respectively, calculating the collaborative offset score of each tree, and combining the collaborative offset scores of each tree to obtain the anomaly score at the current time; maintaining a time evaluation sliding window, calculating the moving average of the anomaly score, and obtaining the first and second derivatives of the anomaly score; executing the judgment logic based on the anomaly score, the moving average, and the values ​​of the first and second derivatives, and outputting the judgment result.

[0009] Furthermore, step S1 specifically includes the following sub-steps: S11: Real-time acquisition of multi-source heterogeneous signals output by multi-source heterogeneous sensors deployed at key stress nodes of infrastructure; S12: Through the interaction of the timing protocol between the master node and the terminal sensor node, calculate the one-way transmission delay and instantaneous clock offset under the asymmetric network path, and estimate the relative drift rate using the clock offset sequence of multiple consecutive timing cycles. Based on this, dynamically adjust the local clock and embed a unified global absolute timestamp for each sensor data packet. S13: For high-frequency dynamic signals with sampling rates higher than the reference frequency, a multi-order Butterworth low-pass anti-aliasing filter is constructed and its amplitude-frequency response function is calculated. The signal is then downsampled after being filtered using zero-phase filtering technology, i.e., the signal passes through the filter in the forward direction, is time-reversed, passes through the same filter again, and is reversed back to the original timing sequence. For low-frequency static signals with sampling rates lower than the reference frequency, upsampling is performed using multiple spline interpolation with natural boundary constraints. A tridiagonal matrix is ​​constructed and the polynomial coefficients are obtained by solving the Thomas algorithm. The known low-frequency data points are substituted into the polynomials corresponding to the polynomial coefficients to calculate the interpolation results at the reference sampling time. S14: For the signal after frequency unification, the extreme value mirror extension method is used to virtually construct the extreme point closed loop at both ends of the signal. The signal is used as the current residual term. Iterative sieving is performed by adaptive empirical mode decomposition. The stopping criterion is that the standard deviation of at least two consecutive iterations is lower than the preset dynamic threshold. The intrinsic mode functions are extracted from the residual term in sequence. The above sieving process is repeated until the residual term becomes a monotonic function or the number of poles is insufficient to construct the envelope. The decomposition is terminated when the residual term becomes a monotonic function or the number of poles is insufficient to construct the envelope. The residual term is the slowly varying trend term. After removing the slowly varying trend term, the reconstructed multi-source pure signal sequence is output.

[0010] Furthermore, step S2 specifically includes the following sub-steps: S21: Synchronously extract the multi-source pure signal sequence through a non-overlapping time sliding window; S22: In the time domain dimension, extract kurtosis, margin factor, skewness, impulse factor and waveform factor from each window sequence to form a high-order dimensionless parameter group in the time domain; S23: In the frequency domain, perform a fast Fourier transform on each window sequence to extract the spectral centroid frequency, and use a bottom-up pruning strategy driven by orthogonal wavelet packet decomposition and Shannon entropy cost function to obtain the optimal basis. Calculate the energy integral and energy ratio of the reconstruction coefficients of each frequency band under the optimal basis. After splicing and normalizing the time domain parameter group, spectral centroid frequency and frequency band energy ratio, generate a high-dimensional fusion feature vector to complete the construction of the multi-dimensional fusion feature space. S24: On the multi-dimensional fusion feature space, construct a multi-order time-frequency collaborative damage sensitivity index, fuse the abnormal quantities of each dimension in the feature space into a single scalar through multiplicative coupling, and compare the index with a preset multi-level early warning threshold. If the sensitivity index exceeds the limit, mark the node as a suspicious node and output it to step 3; otherwise, output the normal state.

[0011] Furthermore, step S3 specifically includes the following sub-steps: S31: Based on the historical pure signal sequence of all nodes, the median of each dimension of each node is used to form a benchmark representative vector; hierarchical clustering is performed according to the rank correlation coefficient between each node to obtain the spatial monitoring cluster. S32: Based on the current high-dimensional fusion feature vector of each suspicious node, calculate the increment of the current feature vector relative to the historical baseline level to obtain the single-point anomaly exceeding the standard multiple of each suspicious node. S33: Calculate the degree of local spatial coupling collapse between each suspicious node and its preset number of spatial nearest neighbors within its cluster; S34: Calculate the global spatial synchronization degree of the cluster to which each suspicious node belongs; S35: The single-point anomaly exceeding the standard multiple, the degree of local spatial coupling collapse, and the degree of global spatial synchronization are fused by nonlinear weighting with adaptive dynamic balance coefficients to construct a single-point-space joint heterogeneous source discrimination function; wherein the dynamic balance coefficient is adaptively selected according to the proportion of nodes exceeding the limit within the cluster; the three-component value forms of the discrimination function are jointly analyzed, and a three-branch hard decision logic is executed according to the combination mode of the three-component value forms, and an independent hardware fault isolation instruction, a local real damage location warning, or a global extreme value environmental event silent monitoring instruction is output according to the decision result.

[0012] Furthermore, the three-branch hard decision logic includes: The first branch, when the single-point anomaly exceeds the limit multiple... The degree of local spatial coupling collapse exceeds Furthermore, the global space synchronization level is lower than that of other systems. When the node malfunction is determined to be due to sensor hardware failure, soft isolation is performed and an operation and maintenance inspection work order is generated; The second branch, when the single-point anomaly exceeds the limit multiple... The degree of local spatial coupling collapse exceeds And the global space synchronization level exceeds When the node is located, it is determined that there is actual structural damage in the local area where the node is located, triggering a regional structural anomaly warning and outputting the spatial location of the node. The third branch considers nodes within a cluster that exceed the threshold for the proportion of nodes exceeding the limit, where the single-point anomaly exceeds the limit multiple. The degree of local spatial coupling collapse of all nodes is lower than Furthermore, the global space synchronization level is lower than that of other systems. When it is determined that the infrastructure is experiencing a global extreme environmental load, all alarms are automatically filtered and the monitoring level is upgraded; in, and These are the strong and moderate thresholds for the multiple of a single-point anomaly exceeding the standard, respectively. , , These represent high, medium, and low thresholds for the degree of local spatial coupling collapse. and These represent the high and low thresholds for global spatial synchronization, respectively.

[0013] Furthermore, step S4 specifically includes the following sub-steps: S41: Construct an isolated forest composed of multiple robust random splitting trees. When constructing each tree, calculate the range of each dimension of the current feature space at each node, extract the cutting dimension with the range as the weight, and randomly select the cutting point uniformly within the value range of the dimension to split the node. S42: Insert real-time feature points into each tree respectively, calculate the collaborative offset score of each tree, and combine the collaborative offset scores of each tree to obtain the anomaly score at the current time; maintain the time evaluation sliding window, calculate the moving average of the anomaly score, and use the finite difference method to calculate the first and second derivatives of the anomaly score. S43: When the abnormal score shows a steep unipolar pulse, the first and second derivatives alternate between positive and negative, and the moving average does not show a significant jump, it is judged as an instantaneous point anomaly and filtered; when the abnormal score exceeds the preset threshold for multiple consecutive cycles, the first derivative is continuously positive, and the moving average shows a continuous step divergence trend, it is judged that the structural dynamic characteristics have undergone irreversible degradation, triggering a structural degradation alarm.

[0014] Furthermore, in step S2, the multi-order time-frequency collaborative damage sensitivity index It can be obtained through the following formula: ; In the formula, Indicates the kurtosis rate in the time domain. Indicates the rate of increase of the time-domain margin factor. This represents the absolute offset rate of time-domain skewness. This represents the frequency domain spectral centroid offset rate. This represents the energy transition rate of the high-frequency wavelet packet band; This is the synergistic incentive index for each indicator.

[0015] Furthermore, in step S3, the single-point anomaly exceedance multiple is obtained by the following formula: ; The degree of local spatial coupling collapse is obtained by the following formula: ; The global spatial synchronization level is obtained using the following formula: ; The single-point-spatial joint heterogeneous source discrimination function is obtained by the following formula: ; In the above formulas, Suspicious node The single-point anomaly exceeding the limit at the current time t. The Euclidean distance between the feature vectors of the current window and the previous window; An adaptive threshold that is dynamically updated; Suspicious node The degree of local spatial coupling collapse with its K nearest neighbors within its cluster. Suspicious node The historical Spearman correlation coefficient benchmark established with node j during the stable period. This represents the real-time Spearman correlation coefficient at the current moment. Suspicious node Cluster The degree of global space synchronization, The historical Spearman correlation coefficient between any nodes u and v within the cluster during the stable period serves as a benchmark. This represents the real-time Spearman correlation coefficient at the current moment. , , This is the dynamic balance coefficient.

[0016] The present invention discloses a multi-source heterogeneous sensor time-frequency feature extraction and fault diagnosis system, comprising a memory, a processor, and a computer program stored in the memory and executable on the processor; the processor executes the computer program to implement the steps of the method described above.

[0017] The beneficial effects of this invention are: 1. In response to the problem of "severe spatiotemporal asynchrony and modal phase distortion" in cross-regional monitoring of large structures, this invention realizes absolute spatiotemporal nanosecond-level reconstruction based on the continuity of physical derivatives.

[0018] Traditional monitoring systems typically rely on software-level NTP protocol alignment and employ crude decimation downsampling or step-hold methods during heterogeneous frequency mapping. This not only introduces millisecond-level time misalignment but also causes severe phase delay and spectral aliasing, directly leading to the failure of bridge modal analysis (such as phase difference calculation). This invention not only introduces PTP nanosecond-level timing with dynamic drift rate compensation at the MAC layer but also innovatively integrates zero-phase frequency domain digital filtering (Filtfilt) with a cubic spline interpolation algorithm with natural boundary constraints. This mechanism, while mapping multi-source heterogeneous data (such as 10kHz vibration and 1Hz temperature) to a unified reference surface, strictly ensures zero total phase shift, perfectly preserving the physical continuity of the original signal's amplitude, first-order rate (velocity), and second-order derivative (acceleration), fundamentally eliminating the spatiotemporal error barrier during heterogeneous data fusion.

[0019] 2. To address the shortcomings of strong background noise and nonlinear environmental interference in easily masking minor structural damage, this invention constructs a high-dimensional time-frequency composite feature space that is exponentially sensitive to early hidden minor damage.

[0020] In complex open environments, the transient impacts of minor damage to bridges and towers (such as microfractures in the main cable's internal steel wires or the propagation of deep microcracks in concrete) are easily drowned out by strong background noise from traffic flow or wind loads. Traditional monitoring indicators based on root mean square (RMS) or fundamental FFT frequency domain amplitudes are essentially unresponsive in the early stages of minor damage. This invention, by forcibly introducing higher-order dimensionless statistics such as kurtosis and margin factor, and combining them with the optimal base wavelet packet (WPT) frequency band energy based on bottom-up pruning of the Shannon entropy cost function, enables the system to possess extremely high sensitivity to non-periodic transient microphysical impacts. Even before significant anomalies are observed in the RMS index, the high-dimensional features constructed in this invention can capture early, hidden structural damage through exponential synergistic leaps.

[0021] 3. To address the maintenance pain point that harsh environments can easily cause sensors to fail and trigger an "alarm storm," this invention achieves precise logical separation between "false sensor failures," "legitimate environmental responses," and "real structural damage."

[0022] Traditional single-point threshold alarm systems cannot distinguish between "loose sensor probes / zero-point drift," "extreme winds / heavy-loaded convoys passing through," and "structural physical degradation," resulting in a high rate of invalid alarms at maintenance sites and severely consuming management resources. This invention breaks through the limitations of isolated single-point data. Based on the inherent spatial mechanical coupling laws of large structures, it innovatively proposes a Ward-Spearman nonlinear spatial topological clustering algorithm based on the Lance-Williams update formula. In real-time monitoring, the system's embedded three-branch Boolean logic decision tree, through the step-cliff-like decay of correlation coefficients between nodes, can accurately "isolate and shield" single hardware failures; simultaneously, by identifying highly synchronized global resonances within clusters, it automatically filters out legitimate extreme environment intrusions. This mechanism completely solves the "alarm fatigue" problem that has plagued the industry, ensuring that the system only triggers high-level warnings for genuine structural load-bearing capacity degradation.

[0023] 4. To address the problem of static model failure caused by slow baseline drift during the service life of large facilities, this invention realizes non-parametric dynamic evolution assessment and qualitative classification of streaming data.

[0024] The dynamic baseline of long-life infrastructure naturally drifts with service life (e.g., material creep, slow prestress decay), and traditional static machine learning models that rely on fixed thresholds or require offline retraining are prone to generating numerous misjudgments. This invention introduces a Robust Randomized Partition Forest (RRCF) model, which achieves unsupervised streaming anomaly capture without retraining by solving for bounding box dimension segmentation probabilities and Co-Disp expectations. More significantly, this invention combines finite difference first derivatives and moving mean tracking within a time-evaluation sliding window, endowing the system with the ability to perform qualitative classification from a temporal morphological perspective: it can accurately distinguish and filter "point anomalies" caused by a single gust of wind or bird strike, and firmly lock onto "trend shifts" caused by structural tearing or yielding, providing solid, reliable, and adaptive mathematical decision support for predictive maintenance of bridges and towers. Attached Figure Description

[0025] Figure 1 This is a schematic diagram of the method flow according to an embodiment of the present invention; Figure 2 This is a schematic diagram of the spatial monitoring topology clustering and three-branch self-diagnosis logic in an embodiment of the present invention; Figure 3 yes Figure 2 The diagram shown is an enlarged view of the spatial topological clustering on the left side of the embodiment. Figure 4 yes Figure 2 An enlarged schematic diagram of the three-branch self-diagnostic logic on the right side of the embodiment shown. Detailed Implementation

[0026] The present invention will be further described in detail below with reference to the accompanying drawings and specific embodiments.

[0027] This embodiment, based on the electrical and physical characteristics of multi-source heterogeneous sensors, achieves long-term monitoring and intelligent self-diagnosis of infrastructure status through a four-step core process: "multi-source spatiotemporal alignment denoising—high-dimensional time-frequency feature extraction—spatial topology clustering and fault stripping—qualitative classification of streaming data anomalies." Figure 1 As shown, the main steps and extremely detailed calculation logic are as follows: S101: Acquire signals from multi-source heterogeneous sensors deployed at key stress nodes of infrastructure; perform hardware-level time synchronization and unified timing on the signals, and dynamically compensate for clock drift rate; then fuse signals with different sampling rates by using zero-phase frequency domain digital filtering downsampling and spline interpolation upsampling with natural boundary constraints, and use adaptive empirical mode decomposition with extreme value mirror extension to remove slowly varying trend terms; output the reconstructed multi-source clean signal sequence.

[0028] This embodiment acquires multi-source heterogeneous sensor signals from key stress-bearing nodes of infrastructure (such as bridge piers, main cables, and tower bases) in real time. For asynchronous continuous data streams with extremely wide sampling frequency ranges, such as vibration (e.g., 2000Hz), strain (e.g., 100Hz), static displacement (e.g., 1Hz), and meteorological temperature, hardware-level time synchronization and frequency mapping with dynamic boundary constraints are implemented at the MAC layer. Specific sub-steps are as follows: A01: Data acquisition from multi-source heterogeneous sensors.

[0029] This embodiment selects a major bridge in China (a double-tower, double-cable-stayed steel box girder bridge with a main span of 800m, tower height of 150m, average annual wind speed at the bridge site of 6.5m / s, and extreme wind speeds reaching 42m / s) as the application object. A heterogeneous sensor network is deployed in key locations such as the bridge tower anchorage zone, the mid-span of the main girder (1 / 4 section, 1 / 2 section, 3 / 4 section), the critical sections of the stay cables (1m from the anchorage end), and the pier foundation caps. The specific configuration is shown in Table 1.

[0030] Table 1 Sensor Node Configuration and Sampling Parameters ; The aforementioned 94 sensor nodes are aggregated to the edge computing node via a fiber optic ring network (redundant topology, 1Gbps bandwidth). All sensor data packets are transmitted using the UDP protocol, and each data packet header includes the node ID and the original sampling time. (Local crystal oscillator count), data type identifier, and checksum field.

[0031] For ease of description, this embodiment uses a unified numbering system for each sensor node: nodes named with the prefix "ACC" are accelerometer nodes (e.g., ACC-08 represents the 8th accelerometer in the mid-span area of ​​the main beam), nodes named with the prefix "STRAIN" are strain sensor nodes, nodes named with the prefix "DISP" are displacement meter nodes, and nodes named with the prefix "WIND" are anemometer nodes.

[0032] This step aims to address the "spatiotemporal asynchrony" problem caused by the extremely wide sampling frequency range (from 1Hz static displacement to 2560Hz high-frequency vibration) and highly nonlinear network transmission delay in cross-regional distributed data acquisition networks. Traditional NTP protocols can only provide millisecond-level synchronization accuracy. For vibration signals with frequencies as high as 2560Hz, millisecond-level errors will result in phase ambiguity of approximately 0.001×2560×360 ≈ 921.6, rendering modal analysis completely ineffective. Therefore, this embodiment introduces PTP nanosecond-level timing with dynamic drift rate compensation at the MAC layer, and combines zero-phase filtering and natural boundary spline interpolation to achieve a unified reference mapping for heterogeneous data.

[0033] A02: High-precision time synchronization protocol interaction and dynamic clock drift rate compensation.

[0034] Due to the kilometer-scale span of large-scale facilities, the transmission delay of distributed wireless or wired data acquisition networks exhibits highly nonlinearity. Traditional systems rely solely on the software application layer's time protocol (NTP), with errors on the order of milliseconds, leading to severe phase misalignment in high-frequency vibration signals. This system employs a precise time protocol (IEEE 1588 PTP) at the underlying layer and utilizes a hardware timestamp unit (TSU) at the MAC layer to record transmission and reception times. The master node and terminal sensor nodes interact through four steps: Sync, Follow_Up, Delay_Req, and Delay_Resp messages. Let the master clock transmission time be t1, the terminal reception time be t2; the terminal transmission time be t3, and the master clock reception time be t4. The terminal side captures these four moments (t1-t4) using hardware to accurately calculate clock deviation and link delay.

[0035] Specific example: In this case, the 94 sensor nodes are distributed across the entire 2.4km length of the bridge (including the north and south approach bridges). The physical fiber optic asymmetric delay between each distributed acquisition terminal and the central clock source is significantly affected by ambient temperature (the temperature coefficient of the fiber optic refractive index is approximately 1.1 × 10⁻⁶). -5 At a temperature difference of 40℃, the unidirectional time delay variation of a 1.2km optical fiber is approximately 1.2 × 10⁻⁶℃. 3 ×1.1×10 −5 ×40×5×10 −9≈ 26.4ns, cumulative effect cannot be ignored). The system's underlying layer adopts the IEEE 1588v2 Precision Time Protocol (PTP). At the MAC layer, the hardware timestamp unit (TSU) of the DP83640 physical layer chip records the transmission and reception times, with a timestamp resolution of 8ns. The master clock (Grandmaster) is deployed in the south tower rack, locking to UTC absolute time via a GPS / BeiDou dual-mode time synchronization module, achieving a timekeeping accuracy of ±15ns and a discipline timekeeping period of 1s. Taking the IEPE accelerometer node at the mid-span of the main beam (node ​​ID: ACC-08, approximately 1.2km fiber optic link from the master clock) as an example, for a certain sampling period (T... s The original timestamps of the four interactions within 1 second are recorded as follows: t1 (Master clock sends Sync message): 1, 000, 000, 000, 000, 000.000 ns; t2 (Terminal receives Sync message): 1000, 000, 000, 125, 400.320 ns; t3 (Terminal sends Delay_Req message): 1, 000, 000, 000, 200, 000.000 ns; t4 (Master clock receives Delay_Req messages): 1, 000, 000, 000, 300, 600.680 ns; The terminal calculates the one-way transmission delay d and instantaneous clock offset in real time under the asymmetric network path: ; ; Where t2-t1 is the total downlink delay of the Sync message from the master clock to the terminal (including the sum of transmission delay and the positive offset of the terminal relative to the master clock), and t4-t3 is the total uplink delay of the Delay_Req message from the terminal to the master clock (including the sum of transmission delay and the negative offset of the terminal relative to the master clock). Adding the two together and dividing by 2 exactly cancels out the terminal clock offset, yielding the pure one-way transmission delay; subtracting the two together and dividing by 2 cancels out the transmission delay, yielding the pure clock offset. This symmetrical path assumption holds largely in fiber optic Ethernet (measured uplink and downlink delay difference is less than 5ns).

[0036] Furthermore, to address the changes in the physical characteristics of the terminal crystal oscillator under extreme temperature differences, this embodiment introduces clock drift rate for dynamic compensation. The relative drift rate of the local crystal oscillator is calculated using a sequence of offsets for M consecutive periods. The system utilizes Adjust the accumulation step size of the local clock to embed a 64-bit nanosecond-level absolute timestamp T, based on a unified global reference, into every sensor data packet entering the system. abs .

[0037] Specific example: Actual measurements revealed that due to the extreme temperature difference between the bridge tower during the day (temperature inside the steel box can reach 40℃) and at night (temperature can reach -5℃ in winter), the physical characteristics of the terminal crystal oscillator (TCXO, nominal frequency 25MHz, frequency temperature stability ±1ppm, aging rate ±1ppm / year) change, causing the offset to exhibit a deterministic trend drift rather than white noise. Under a cyclic sliding window with a system maintenance length of M = 120 (corresponding to 120 seconds, sufficiently smoothing instantaneous network jitter while timely tracking temperature changes), the offset sequence for M consecutive cycles was recorded. The drift rate was estimated using a least-squares linear fitting method. Its normal equation can be written as: ; in, The physical meaning is the regression slope of the clock offset over time, that is, the average rate of change of the clock offset per unit time. Substituting the measured offset data of the ACC-08 node for 120 consecutive cycles... The calculation formula is obtained during the daytime (ambient temperature approximately 35°C). (Positive drift indicates that the terminal crystal oscillator frequency is too fast, and the local clock cumulatively leads the master clock.) Calculated at night (ambient temperature approximately 0°C). (Negative drift indicates that the terminal crystal oscillator frequency is too slow). The system recalculates every minute. and with the current The cumulative step size of the local clock is adjusted in real time, and its frequency-locked loop update formula is: ; in, =1 / 2560s is the nominal sampling interval of the ACC-08 node. The global absolute timestamp embedded for the nth sampling point. This is the instantaneous offset obtained from the most recent synchronization via the PTP protocol (synchronization occurs hourly to reduce network overhead). (Second item on the right) This is the drift rate compensation term, and its physical meaning is: due to At nominal interval The actual accumulated time error is Seconds, the length of this sampling interval needs to be changed from Adjusted to This is used to offset the error. After dynamic compensation using the above formula, the synchronization accuracy of the ACC-08 node stabilized within ±32 ns during continuous 24-hour operation, an improvement of approximately 37,500 times compared to the traditional NTP protocol (error ±1.2 ms). Each compensated data packet header additionally records the estimated cumulative error e after compensation. est (Unit: ns) is used as the confidence label; if e est If the time exceeds 100 ns, the data packet is marked as "low confidence" for reference in subsequent steps.

[0038] After completing the above processing, the system embeds a 64-bit global reference absolute timestamp into each sensor data packet of the ACC-08 node (each data packet contains 2560 sampling points, corresponding to a 1-second data block, and is packaged into data blocks to reduce protocol overhead). The embedding format is as follows: the seconds portion occupies 32 bits (calculated from GPS epoch January 1, 2000, covering approximately 136 years), and the nanosecond portion occupies 32 bits (1 ns resolution), with a total range sufficient to cover the entire lifespan of the bridge. Simultaneously, the system adds a compensated cumulative error estimate e to the data packet header. est If e est Values ​​greater than 100 ns are marked as "low confidence" for reference in subsequent steps.

[0039] A03: Zero-phase frequency domain digital filtering downsampling and natural boundary constraint upsampling.

[0040] Specifically, the system-level reference mapping frequency f is set. base .

[0041] For sampling rate f H > f base For high-frequency dynamic signals (such as acceleration), traditional direct decimation downsampling can cause spectral aliasing, while conventional low-pass filtering introduces severe phase delay, disrupting the phase difference of mode shapes. This system constructs an Nth-order Butterworth low-pass anti-aliasing filter and calculates its amplitude-frequency response function. The system employs zero-phase forward and backward filtering (Filtfilt), where the signal first passes forward through a difference equation filter, then the sequence is reversed and passed backward through the same filter, ensuring a strictly zero total phase shift and guaranteeing phase synchronization of the downsampled signal across nodes.

[0042] For sampling rate f H < f base For low-frequency static signals (such as temperature and cable tension), the traditional step-preservation method results in discontinuities in the first derivative (rate of change) of the signal. This embodiment employs cubic spline interpolation with natural boundary conditions. That is, given the low-frequency data points... In the interval The piecewise function on is constructed as follows: ; The boundary conditions are set such that the second derivative at both endpoints is zero (i.e., and System Construction The tridiagonal matrix is ​​used to quickly solve the system of equations using the Thomas Algorithm to obtain the polynomial coefficients. This mixed... The alignment mechanism perfectly preserves the physical continuity of the original signal amplitude, first-order rate, and second-order acceleration.

[0043] Specific example: In this case, the system-level reference mapping frequency f is set. base = 128 Hz. The technical basis for choosing 128 Hz as the reference frequency is as follows: the first 20 modal analyses of the bridge show that the highest effective modal frequency is approximately 56 Hz (the 18th torsional mode), according to the Nyquist sampling theorem f base ≥ 2 × 56 = 112 Hz. Selecting 128 Hz provides approximately 14% engineering margin while also considering computational efficiency (each 4-second window contains 512 sampling points, making the FFT computation moderate). All sensor data, regardless of the original sampling rate, are uniformly mapped to this reference frequency before subsequent fusion analysis. Sampling processing methods include: (1) Downsampling of high-frequency signals: Taking ACC-08 (2560 Hz → 128 Hz) as an example: For the original sampling rate A high-frequency acceleration signal of 2560 Hz (including 24 IEPE accelerometers) will, if downsampled directly by an integer multiple of 20:1, suffer from irreversible spectral aliasing distortion according to the Nyquist-Shannon sampling theorem. This occurs because frequency components above 64 Hz will be mirrored and folded into the 0-64 Hz passband. For example, the 78.2 Hz structural mode component will be folded into the passband of 0-64 Hz. The spurious Hz frequency appears in the downsampled signal, causing complete errors in mode identification. Although conventional low-pass filtering can suppress high-frequency components, IIR filters (such as Chebyshev type) introduce significant nonlinear phase at the passband edge. The group delay can reach ±12 sampling points near the cutoff frequency (corresponding to a phase error of ±12 / 2560×360°×78.2≈±131°), which destroys the mode shape phase difference across nodes.

[0044] This system constructs a 6th-order Butterworth low-pass anti-aliasing filter, whose squared amplitude-frequency response function is: ; Where N = 6 is the filter order. This is the cutoff angular frequency. Here, 0.45 × f is chosen. base Instead of 0.5×f base This is to reserve a 10% transition band (57.6 ~ 64Hz), making the filter's amplitude-frequency response flatter in the passband and resulting in faster stopband attenuation. The stopband attenuation slope of a 6th-order Butterworth filter is... At 80Hz, the theoretical attenuation is approximately It can effectively suppress the 78.2Hz modal component (attenuating it to about 6.3% of the original amplitude), while the amplitude fluctuation within the passband (<57.6Hz) is less than 0.5dB.

[0045] The corresponding digital filter difference equation is: ; in, The feedback coefficients are the coefficients of the denominator polynomial, normalized by a0=1. These are the feedforward coefficients (numerator polynomial coefficients). These coefficients are obtained by mapping the simulated Butterworth prototype H(s) using the bilinear transform method, ensuring that the frequency response of the digital filter is highly consistent with the simulated prototype in the low-frequency range (pre-distortion processing ensures accurate cutoff frequency). The specific values ​​are (stored as double-precision floating-point numbers): , .

[0046] The system employs zero-phase forward-backward digital filtering technology (Filtfilt). The core of this technology lies in using two filtering operations in opposite directions to cancel phase distortion. The specific mathematical process is as follows: ; ; ; in, This is the result of the first forward filtering, where x[n] is the original 2560Hz signal sequence. This indicates that one-dimensional digital filtering is performed according to the difference equation of the digital filter; This is the result of inverse filtering; This results in a zero-phase output signal. The first forward filter introduces a phase delay. The inverted filter then introduces the same phase delay again. However, due to the signal reversal, the two delays cancel each other out, and the phase frequency response of the overall transfer function is always zero. The overall transfer function of this operation... Phase response (The rounding error in the numerical calculation is approximately ±0.02°). After filtering, the result is expressed as f. base Resampling, which involves sampling 1 point out of every 20 points: ; in, is the downsampled output signal, and m is the downsampling index.

[0047] Taking a 1024-point window of data from ACC-08 (corresponding to a 0.4-second time-domain waveform at the original 2560Hz) as an example, this window contains the fourth torsional mode of the structure (frequency 78.2Hz, corresponding to the coupled vibration of the stay cable and the main beam). After a traditional single Butterworth filter, this 78.2Hz component has a phase lag. (Group delay approximately 1.52ms), amplitude attenuation 0.7dB; however, after processing with Filtfilt in this embodiment, the phase shift is strictly 0.00° ± 0.02°, and the amplitude fidelity reaches 99.8% (attenuation only 0.009dB). Then, it is downsampled to 128Hz according to the downsampling formula, which completely preserves the modal shape phase difference information between nodes, providing high-quality phase synchronization data for subsequent modal analysis.

[0048] (2) Upsampling of low-frequency signals: Taking the displacement meter DISP-03 (1 Hz → 128 Hz) as an example: For the original sampling rate f L For a low-frequency static signal of 1Hz (including 12 magnetostrictive displacement gauges), traditional zero-order hold (ZOH) interpolation, while simple, has two major drawbacks: First, a (sin x) / x roll-off effect appears in the spectrum. At the new sampling rate of 128Hz, holding the original 1Hz signal will generate a large number of high-frequency image components. Second, the stepped output mathematically results in huge spikes (theoretically infinite) in the first derivative (velocity), and the second derivative (acceleration) is nonexistent, which is physically unacceptable. For example, when a bridge experiences continuous and smooth displacement changes as a train passes, ZOH interpolation will produce "step-like" artifacts, causing severe distortion in the subsequently calculated velocity and acceleration waveforms.

[0049] This embodiment employs natural cubic spline interpolation constrained by natural boundary conditions. The mathematical essence of this method is to construct a cubic polynomial between adjacent original data points, ensuring that the displacement, velocity, and acceleration of the entire interpolation curve are continuous, and that it bends "freely" at the endpoints (with the second derivative being zero), thereby guaranteeing physical continuity.

[0050] Take 24 consecutive hours of raw data from DISP-03 (a total of 86,400 points, time intervals) Equal intervals, = 0, 1, ..., n-1, n = 86400). In the interval... Define a piecewise cubic polynomial: ; in, This represents the displacement value of the left endpoint of the interval. The first derivative (velocity) is the left endpoint. It is half of the second derivative (acceleration) at the left endpoint. It is one-sixth of the third derivative (jerk). The four undetermined coefficients require four constraints to determine. The interpolation function satisfies the following continuity constraints: ; ; ; The above constraints ensure that the interpolation curve passes through all the original data points precisely without introducing any amplitude deviation; the first derivative of displacement with respect to time (i.e., velocity) is continuous along the entire time axis without any abrupt changes; the second derivative of displacement with respect to time (i.e., acceleration) is also continuous, which conforms to the smooth change of displacement with time in actual physical processes.

[0051] make For the first The second derivative values ​​at each of the original data points Based on the theoretical derivation of cubic splines, Satisfy the following three bending moment equations: ; Where, the left side of the equal sign is , , A weighted linear combination, where the weights are determined by the lengths of adjacent intervals. , The longer the interval, the greater the contribution of M at that end to the center point, which aligns with the intuitive meaning of torque distribution in mechanics. The right-hand side of the equation is the difference in average slopes over adjacent intervals divided by the sum of the interval lengths; its physical meaning is a discrete measure of curvature change. When = When = 1, the equation simplifies to: ; The natural boundary condition is that the second derivative at both endpoints is zero, i.e., M0 = 0, M n = 0, which means that at the two endpoints of the interpolation interval, the curve is not subject to external torque and can bend freely, ensuring that the interpolation will not warp abnormally at the endpoints. This property is crucial for the physical rationality of the displacement signal boundary.

[0052] For this case, the sampling is at a medium interval ( In the case of (=1s), the three bending moment equations construct a... The system of linear equations is a strictly tridiagonal matrix, diagonally dominant (4 elements on the main diagonal and 1 element on the secondary diagonal). This system can be solved efficiently in O(n) time using the Thomas algorithm. The Thomas algorithm consists of two stages: forward elimination and back substitution. Forward elimination : ; Back-substitution solution : ; in, , , , The subdiagonal, main diagonal, upper subdiagonal, and right-hand vector of a tridiagonal matrix are respectively located in the rows. The element at that position. For the case of equal intervals, =1, =4, =1, = Forward elimination transforms the coefficient matrix into an upper triangular matrix, while back substitution directly yields all the coefficients. Seeking Then, back-substitute to calculate the four coefficients of the polynomial for each interval: ; in, It is directly equal to the displacement value of the left endpoint, ensuring that the interpolation curve passes precisely through the original data points; This is the first derivative correction value after taking into account the curvature contribution; It is half of the second derivative at the left endpoint; It is one-sixth of the rate of change of the second derivative within the interval. These four coefficients completely determine the displacement value and its derivatives at any time within the interval.

[0053] During interpolation, in each original interval Within a 1-second interval, 127 new points are inserted at equal intervals with a step size of 1 / 128 seconds, resulting in a total of 86400 × 128 = 11059200 interpolation points, fully covering the 24-hour duration. For any interpolation time... The displacement, velocity, and acceleration values ​​are calculated using a cubic spline piecewise polynomial and its first and second derivatives, ensuring strict physical continuity. Experimental results show that the correlation coefficient between the interpolated displacement signal and the velocity waveform during train passage is 0.97, and the correlation coefficient between the acceleration waveform and the original accelerometer waveform is 0.94, significantly better than the 0.51 and 0.38 of ZOH interpolation.

[0054] A04: Extreme value mirror extension and adaptive environment strong noise stripping.

[0055] Specifically, bridge towers undergo significant, slow, nonlinear deformation under the influence of solar radiation and temperature differences. This deformation, superimposed on the strain signal, severely masks the subtle high-frequency damage characteristics of the structure itself. The system utilizes an adaptive empirical mode decomposition (EMD) algorithm for trend stripping. To address the data divergence problem caused by the endpoint effect in traditional EMD, this embodiment introduces an extreme value mirror extension method, virtually constructing a closed loop at the extreme point at both ends of the signal. The original stationary signal is then extracted. Mean of upper and lower envelopes Iterative calculation Then, the standard deviation (SD) is introduced as a strict adaptive iteration stopping criterion.

[0056] During the iterative screening process, when SD is less than the preset dynamic threshold twice consecutively... At this point, the current round of screening is stopped, and an intrinsic mode function (IMF) is extracted from the remaining signal. The final signal is orthogonally decomposed into: ; The system removes low-frequency residual terms that represent gradual changes in temperature / solar radiation. The output is the reconstructed multi-source clean signal sequence.

[0057] Specific example: The strain drift caused by thermal expansion and contraction of the steel box girder of a large cable-stayed bridge under daily temperature difference (ΔT ≈ 25℃) can reach 180°. (Based on the coefficient of linear expansion of steel) The beam is 800m long, with a temperature difference of 25℃, and the strain changes... In reality, due to the partial release of structural constraints, the measured value is approximately 180. The amplitude of this low-frequency trend term is much larger than that of structural micro-damage (5~20). The dynamic responses generated by these two methods differ by more than an order of magnitude. Traditional polynomial fitting detrending methods cannot adapt to non-stationary temperature changes (such as abrupt changes in the rate of temperature change caused by alternating sunny and cloudy weather or the passage of cold air), and the selection of the polynomial order lacks objective basis. This system adopts adaptive empirical mode decomposition (EMD) combined with extreme value mirror continuation to achieve data-driven adaptive trend stripping.

[0058] Taking the continuous 10-hour strain signal from the STRAIN-12 fiber optic strain gauge on the bottom plate of the main beam as an example. The original sampling rate is 200Hz. To reduce the computational complexity of EMD (the iterative sieving computational complexity of EMD is O(N...),... 2(Direct processing of a 200Hz signal is time-consuming). First, downsample to 5Hz (downsampling is performed according to step A03), resulting in N = 10 × 3600 × 5 = 180,000 data points (for clarity, this case selects a subsequence of T = 1,800 seconds, or 30 minutes, for detailed demonstration; the processing flow for the complete 10 hours of data is exactly the same). This subsequence is denoted as X(t) (t = 0, 1, …, 1799, corresponding to a time resolution of 0.2 seconds).

[0059] The key challenge of EMD decomposition lies in the "endpoint effect": due to the lack of data at both ends of the signal, the upper and lower envelopes diverge at the endpoints, and this divergence "contaminates" the entire decomposition process. To address this, this embodiment introduces the extreme value mirror extension method before each sieving iteration, extrapolating three extreme points at each end of the signal to form a closed-loop extension sequence before constructing the envelope.

[0060] Taking the extension at the left endpoint t0 of the signal as an example: First, extract the first 5 local maxima points (in chronological order) at the left end of the signal, denoted as... (In this case subsequence, the time coordinates of the first 5 local maxima are as follows) , , , , ). Mirror the left endpoint t0 to construct a virtual extreme point: ; The physical meaning of this formula is: using the left endpoint as a mirror, it "reflects" the extreme points inside the signal to the outside of the endpoint, causing the envelope to extend naturally at the endpoint. For example, the first maximum point (2.0s, 12.3) ) mapped to The second mapping is Similarly, the minimum points are handled. After extension, three virtual maxima and three virtual minima are added to the left end, so that the envelope fit no longer diverges at the endpoints.

[0061] After the extension is completed, cubic Hermite interpolation is used to fit all the maxima and minima respectively to obtain the upper envelope e. max (t) and lower envelope e min (t). Calculate the mean of the envelope. ,make The goal of EMD screening is to make h k (t) becomes a signal that satisfies two conditions of the intrinsic mode function (IMF): ① the number of extreme points and the number of zero crossings are equal or differ by at most 1 throughout the entire data segment; ② at any time, the upper and lower envelopes defined by the local maxima and local minima have a mean of zero.

[0062] Introducing the strict standard deviation (SD) stopping criterion to determine whether the IMF conditions are met: ; Among them, h k-1 (t) represents the residue from the previous iteration, h k (t) represents the residue of the current iteration. This is the disturbance constant (its order of magnitude is much lower than the strain gauge A / D quantization noise of 0.1). The square of 0.01 has no effect on the result. The numerator of this criterion is the sum of the squares of the differences between the residuals of the two iterations, and the denominator is a normalization factor to prevent division by zero. Its physical meaning is: when the change between two iterations is sufficiently small relative to the signal energy... It is assumed that the convergence has reached a stable IMF component. A dynamic threshold is set. (Empirical value, determined based on the signal-to-noise ratio characteristics of the strain signal in this case, when...) When the concentration is 0.03, it ensures sufficient decomposition without over-sieving.

[0063] Iterative calculation record: First iteration: The calculated SD = 0.86 (due to the inaccurate initial estimate of the envelope mean, the value varies greatly), which is much greater than 0.03, and therefore does not meet the requirement. Fourth iteration: The envelope gradually converges, and the calculated SD = 0.112, which is still greater than 0.03; 7th iteration: SD = 0.041, which is still greater than 0.03 but close to the threshold; 8th iteration: SD = 0.028 is calculated, which satisfies the condition < 0.03 for the first time; 9th iteration: To ensure the stability of convergence (to prevent premature termination due to accidental numerical fluctuations), the iteration continues, and SD is calculated to be 0.031, which is slightly greater than 0.03 (due to the slight fluctuation caused by the gradual change in temperature drift). 10th iteration: SD = 0.027 was calculated, again satisfying < 0.03. The threshold condition was met twice consecutively, the iteration stopped, and the first intrinsic mode function IMF1 = h was successfully extracted. 10 .

[0064] Repeat the above screening process to extract IMF2, IMF3, ... from the residual signal in sequence until the residual term r is reached. n (t) becomes a monotonic function (fewer than 2 poles) or the number of poles is insufficient to construct the envelope. Ultimately, the strain signal in this case is decomposed into: ; Among them, IMF1~IMF4 correspond to the elastic vibration modes of the bridge (frequency range 0.5~8Hz, such as the first vertical bending mode of the main beam at 0.68Hz and the second torsion mode at 3.2Hz, etc.), IMF5~IMF8 correspond to the transient response caused by environmental noise and traffic load (frequency range 8~25Hz), and r8(t) is the low-frequency residual trend term (frequency < 0.1Hz).

[0065] The system calculates the normalized cross-correlation function between the residual term r8(t) and the temperature T(t) of the steel box girder recorded by the synchronous PT100 temperature sensor: ; Where τ is the time delay (in seconds) and N is the sequence length (N = 1,800). The normalized correlation coefficient is... ,in and These are the variances of the residual term and the temperature series, respectively. Substituting the measured data, we obtain... This verifies that the energy of r8(t) is mainly caused by temperature changes (the residual term waveform is almost synchronized with the temperature waveform, and the correlation coefficient of 0.93 indicates that the temperature factor explains about 86% of the residual term variance). The system automatically removes r8(t) and reconstructs a pure strain signal. ; After reconstruction, the daily average value of the strain signal increased from the original 119.6. (Almost entirely due to temperature) dropped to 2.3 The peak value was reduced from the original 140 Dropped to 8.5 The signal-to-noise ratio of micro-damage is reduced from the original The gain is increased to +16dB, reaching 24dB, enabling subsequent higher-order statistics to effectively detect early micro-damage. The final output is a reconstructed, clean signal sequence.

[0066] Step S102: Based on the multi-source clean signal sequence, the signal sequence is synchronously truncated through a non-overlapping time sliding window, and damage-sensitive features are extracted from both the time and frequency domains: in the time domain, high-order dimensionless statistics sensitive to mutations are extracted, and in the frequency domain, the spectral centroid frequency and wavelet packet energy features are extracted; the time-domain and frequency-domain features are spliced ​​and fused to construct a multi-dimensional fused feature space; based on the abnormal quantities in each dimension of the feature space, a multi-order time-frequency collaborative damage sensitivity index with a multiplicative collaborative excitation mechanism is constructed, and the multi-dimensional abnormal signals are fused into a single scalar through multiplicative coupling; the current node status is determined based on the comparison result of the sensitivity index and a preset threshold, and a normal node or suspicious node label is output.

[0067] Specifically, for the multi-source clean signal sequence after alignment and denoising in step S101, a multi-dimensional fusion feature space is constructed by simultaneously extracting time-domain high-order dimensionless statistics and frequency-domain optimal basis wavelet energy, which are highly sensitive to structural micro-damage, through a non-overlapping time sliding window.

[0068] In this case, after step S101, all data have been mapped to a 128Hz reference frequency and environmentally variable components have been removed. This step uses a sliding window (length L = 2048 points / 16 seconds, step size 512 points / 4 seconds, overlap rate 75%) to simultaneously extract time-domain higher-order statistics and frequency-domain wavelet energy, constructing a 16-dimensional fused feature space. Specific sub-steps are as follows: B01: Extraction of high-order dimensionless sensitive features in the time domain.

[0069] Early-stage localized minor damage to infrastructure (such as micro-fractures in prestressed steel strands and fretting wear in anchorage zones) manifests as extremely short, transient, non-periodic physical impacts. Within a sliding time window, the system not only calculates the root mean square (RMS) but also forcibly extracts a series of high-order dimensionless parameter groups. For example, it calculates quinary high-order statistics in the time domain within the window (2048 points). Let the pure signal sequence be... mean Standard deviation The definitions of each statistic are as follows: ; ; ; ; ; The above formulas are respectively: kurtosis (The fourth moment, which is most sensitive to transient shocks), margin factor (Ratio of extreme values ​​to the square root mean that are highly sensitive to weak shocks), skewness (Reflecting signal distribution asymmetry), pulse factor (Ratio of extreme values ​​to absolute mean), waveform factor (The ratio of the root mean square to the absolute mean).

[0070] By using the aforementioned high-order dimensionless parameter group, the system can accurately separate and sensitively amplify the weak non-stationary transient impacts caused by early micro-losses from strong background noise and steady vibrations, fundamentally making up for the shortcomings of traditional root mean square (RMS) which easily misses transient spike signals.

[0071] Taking ACC-08 as an example, the normal operating window (day 120): μ = 0.005 m / s², σ = 0.092 m / s², calculate K = 3.12, Cf = 5.34, S k = 0.08, I f = 6.12, S f = 0.94. Early morning window on day 213: RMS only rose from 0.092 to 0.101 (+9.8%, completely missed by the traditional threshold method), but calculated K=9.09 (a jump of 191%), C f = 28.93 (a jump of 442%) (From positive to negative), I f = 7.52. Continuous monitoring revealed a similar pulse occurring approximately every 47 seconds. The bridge inspection vehicle confirmed that the transient stress wave was caused by friction between the steel strands due to damage to the PE sheath of the LZ12 cable on the south side. This method provides an early warning 63 days earlier than the traditional RMS method.

[0072] B02: Frequency domain energy centroid and optimal basic wavelet packet frequency band feature extraction.

[0073] Perform a Fast Fourier Transform (FFT) on the window sequence to extract the spectral centroid frequency (FC) and quantize the overall macroscopic shift of the structural mode natural frequencies: ; in, It is the power spectral density. The highest cutoff frequency is the Nyquist frequency of the signal at the current sampling rate.

[0074] Furthermore, to finely characterize the energy distribution under non-stationary wind loads, orthogonal Daubechies family wavelets are used for j-level wavelet packet decomposition (WPT). Recursion is then performed using a dual-channel orthogonal mirror filter bank (QMF). ; ; Where h(n) and g(n) are the coefficients of the low-pass and high-pass filters, respectively. The system uses Shannon entropy as the cost function to perform bottom-up pruning on the wavelet packet tree to obtain the optimal tree structure. The energy integral E of the reconstruction coefficients in the m-th frequency band under the optimal basis is calculated. m And calculate the energy percentage. The system concatenates abrupt change-sensitive time-domain statistics with frequency-domain macroscopic offsets and fine wavelet energy proportions, and then applies Z-score normalization to generate a high-dimensional fused feature vector F. t .

[0075] Specific example: Performing a 1024-point Fast Fourier Transform (FFT) on a window sequence (from zero-padding to 2048 points, frequency resolution 0.0625Hz): ; Calculate the spectral centroid frequency FC and quantize the overall macroscopic shift of the structural mode natural frequencies: ; In the 213-day window, the FC drifted from the baseline of 1.03 Hz to 1.41 Hz (+36.9%), indicating a significant nonlinear response in the structure. Furthermore, to finely characterize the energy distribution under non-stationary wind loads, orthogonal Daubechies 6 (db6) wavelets were used for 5-level wavelet packet decomposition, generating 32 equal-width frequency bands (each 2 Hz). Recursion was then performed using a dual-channel orthogonal mirror filter bank (QMF). ; ; Where m is the number of frequency bands, and h0 is the low-pass filter coefficient. The filter is a conjugate orthogonal mirror high-pass filter that satisfies the complete reconstruction condition. The system uses normalized logarithmic energy entropy as the cost function for bottom-up pruning. ; in, These are the wavelet packet coefficients at the current wavelet packet node. It is the cost function value of the current wavelet packet node, that is, the normalized logarithmic energy entropy calculated on the wavelet packet coefficients at this node.

[0076] The pruning decision condition is: if the entropy value of the parent node is... If the energy is already concentrated, retain the parent node; otherwise, retain the child node. The final set of optimal base coverage frequency bands is: .

[0077] Calculate the energy integral of the reconstruction coefficients of the m-th frequency band under the optimal basis. and energy percentage: ; in, It is the single-branch reconstructed signal (time series) of the m-th optimal fundamental frequency band. It is the first Energy integral of each frequency band.

[0078] During the 213-day window, the energy share of the B3 (4-6Hz) band jumped from 0.7% at the baseline to 8.31% (approximately 11.9 times), which correlated with kurtosis K and margin factor C. f The abrupt changes corroborate each other, further confirming the existence of high-frequency nonlinear response from a frequency domain perspective.

[0079] The time-domain quintuple (S) k , K, S f , If C f The frequency domain macroscopic quantities (FC and modal frequency shift) and the frequency band energy proportions of the nine optimal basis wavelet packets are concatenated to form a 16-dimensional original feature vector X. t Online Z-score standardization using a sliding window: ; in, The mean of the sliding window at the current moment. Let be the variance of the sliding window at the current time; α = 0.01 (extremely slow update rate ensures long-period baseline stability), μ0 = 0, σ0 = 1, ultimately generating a 16-dimensional fused feature vector F. t .

[0080] B03: Construct a multi-order time-frequency synergistic damage sensitivity index.

[0081] To address the shortcomings of traditional methods that use time-domain and frequency-domain features independently and lack synergistic amplification effects, making it difficult to effectively identify early-stage micro-damage signals, a dedicated damage sensitivity index with a nonlinear synergistic excitation mechanism is constructed. This index integrates multidimensional anomalous signals such as kurtosis jumps, margin factor surges, spectral centroid shifts, and energy proportion transitions in specific high-frequency wavelet packet bands into a single scalar through multiplicative coupling. This enables exponentially sensitive capture of early-stage, concealed micro-damage, providing a high-confidence pre-trigger signal for subsequent early warning judgments.

[0082] Based on the time-domain quinary higher-order statistics, spectral centroid frequency, and optimal base wavelet packet band energy proportion extracted in steps B01 and B02, a specific damage sensitivity index with a multiplicative co-excitation mechanism is constructed. Let the rate of change of each index relative to the stable baseline in the current window be: ; ; ; ; ; These rates of change represent the kurtosis rate in the time domain, respectively. The rate of increase in time margin factor absolute offset of time-domain skewness Frequency domain spectral centroid offset and high-frequency wavelet packet band energy transition rate (Select the optimal energy percentage of the intermediate frequency band B3) (as a high-frequency damage-sensitive indicator band). To prevent division by zero, we take 10. -6 ; The base is the skewness reference damping constant (to prevent excessive amplification of the rate of change when the reference skewness is close to zero), and is set to 0.1; the subscript "base" indicates the median reference value of the corresponding index during the stable period.

[0083] The physical meanings of the above rates of change are as follows: It reflects the relative change in the sharpness of the signal waveform and is most sensitive to transient impacts; It reflects the relative change of signal extrema with respect to the overall energy distribution and has extremely high recognition of weak and short pulses; The degree of disruption to the symmetry of signal distribution is reflected in structural damage, which is usually accompanied by the release of asymmetric stress waves. It reflects the overall macroscopic shift of the structural modal natural frequencies and characterizes the relative change in structural stiffness; The energy percentage transitions reflect the high-frequency damage-related frequency bands (corresponding to the higher-order vibration modes in the cable-stayed cable surface and the nonlinear responses generated by the opening and closing of microcracks).

[0084] Taking the LZ12 case failure incident on day 213 as an example: Normal operating window: , , , , Measured value of the damage window on day 213: , , , , .

[0085] Calculate each rate of change: ; ; ; ; ; Traditional methods typically treat these features as independent variables, using linear weighting or simple concatenation before inputting them into a classifier. The subtle changes in each feature, under linear superposition, are insufficient to form an effective detection signal. This embodiment constructs a multi-order time-frequency collaborative damage sensitivity index with a multiplicative collaborative excitation mechanism. : ; In the formula, The synergistic incentive index for each indicator has the physical meaning of being the amplification factor of the rate of change of each dimension in multiplicative fusion. The larger the value, the more important that dimension is. Its contribution shows superlinear growth. The physical meaning is: the numerator is the multiplicative combination of the five characteristic rates of change (amplification effect), and the denominator is the arithmetic sum of the five characteristic rates of change (normalization benchmark). When each rate of change is small, When multiple rates of change increase significantly at the same time, the numerator grows exponentially while the denominator grows only linearly. Significantly greater than 1, achieving coordinated amplification of multi-physics field anomalous signals.

[0086] After optimization by grid search, the following was obtained: , (The margin factor is most sensitive to minor shocks and is given the highest incentive index.) , , Substituting the above rate of change: ; To suppress false alarms caused by single transient disturbances, a time-domain smoothing mechanism is introduced to calculate the cumulative weighted sensitivity index: ; in, Set the smoothing coefficient to 0.3; Among them, the historical base value : ; The system sets four warning thresholds: when Output "Normal" when it is active; when it is not active. When triggering "Blue Alert" (initial anomaly detected, encrypted monitoring recommended); when When an "orange alert" is triggered (multi-feature coordinated anomaly, early structural damage may occur); when A "red alert" is triggered when multiple features exhibit exponential and coordinated leaps, indicating significant physical degradation of the structure, requiring immediate action. The thresholds mentioned above are set based on the following condition: under normal operating conditions, the rate of change of each feature is less than 0.1. Not exceeding 1.5; when the single characteristic exhibits moderate variation (rate of change approximately 0.3~0.5). Not exceeding 2.0; when a certain characteristic undergoes extreme changes (such as a significant shift in skewness due to an impact pulse), It will enter the blue attention zone of 2.0~20; when multiple features change significantly and simultaneously, It will enter the orange alert range of 20-200; when multiple features simultaneously exhibit extreme changes and produce a multiplicative amplification effect, If the value exceeds 200, a red alarm will be triggered. In this embodiment, if... If the node is marked as a "suspicious node", proceed to step S103 for spatial positioning diagnosis; otherwise, output "normal" status directly.

[0087] In contrast, if a linear weighted fusion method is used (with all weights set to 0.2): ; Linear fusion can only reach the "blue concern" level ( (This indicates that) no orange or red alerts can be triggered. The multiplicative ratio mechanism achieves a synergistic amplification effect approximately 129 times higher than linear fusion, demonstrating the invention's hypersensitive ability to capture multi-feature synergistic anomalies.

[0088] S103: Based on nodes marked as suspicious, construct a single-point-space joint heterogeneous source discrimination function containing three components: single-point anomaly exceeding the standard multiple, local spatial coupling collapse degree, and global spatial synchronization degree. Execute three-branch hard judgment logic according to the combination mode of the three component values, and output independent hardware fault isolation instructions, local real damage location warnings, or global extreme value environmental event silent monitoring instructions, respectively.

[0089] This step breaks down the "island" judgment of single-point data, utilizing the inherent spatial physical coupling law between measurement points on large structures to achieve precise differentiation between "sensor failure," "local damage," and "global environmental response" by quantifying spatial correlation. For example... Figures 2-4 As shown, the specific sub-steps are as follows: C01: Construct a single-point-space joint heterogeneous source discrimination function.

[0090] A 30-day stable baseline dataset, i.e., a historical clean signal sequence, is extracted. A window feature vector is generated every 4 seconds, totaling 648,000 samples. The median of each dimension is taken for each node to obtain the baseline representative vector. Calculate the Spearman rank correlation coefficient between any two nodes. If equal eigenvalues ​​exist, introduce a same-rank correction factor: ; in, Let t be the number of groups with equal feature values ​​in node X. X, g Let G be the number of equal values ​​in the g-th group. Y and t Y, g Similarly, the corrected Spearman correlation coefficient is: ; In the formula, As a level, The average grade. The ρ values ​​of ACC-08 and ACC-11. s = 0.92 (strong spatial coupling), with ρ of WIND-01 s= 0.23 (weak correlation). After calculating the full 94×94 correlation matrix, Ward's minimum variance hierarchical clustering is applied. For ease of explanation, clusters A to C are used as examples to illustrate the clustering details; the clustering principles for the remaining clusters are completely consistent with the three clusters mentioned above. The following detailed calculation explanation uses the Ward clustering process of clusters A to C as an example.

[0091] The sum of squared deviations within cluster A is: ; The ESS increment when merging clusters A and B is: ; New cluster after merging The distance update to the existing cluster C follows the Lance-Williams recurrence relation: ; The weighting coefficients are: ; In the formula, , , These represent the number of nodes (sensors) contained in clusters A, B, and C, respectively.

[0092] Clustering stops when the merging distance jumps abruptly from 0.41 to 2.13 (the jump is 34 times that of the previous step), and the 94 nodes are automatically divided into 9 spatial monitoring clusters. .

[0093] For each "suspicious node" output in step S102 ( Construct the following three discriminant components: Component 1: Single-point anomaly exceeding the standard multiple .

[0094] Let node The relative increment of the feature (exceedance factor) at the current time t is: ; in, This is the Euclidean distance (relative feature increment) between the feature vectors of the current window and the previous window. This is an adaptive threshold that is dynamically updated using the EWMA algorithm. The physical meaning of a node The intensity multiple of the current abnormal impact; the larger the value, the more serious the deviation of the node's characteristics from the historical baseline.

[0095] The following is a real-world example demonstrating the complete process for identifying fault types: Case 1: Sensor independent hardware failure identification - ACC-09 probe loose.

[0096] cluster Inner node ACC-09 Current threshold Exceeding the standard by multiples: ; Component 2: Degree of Local Spatial Coupling Collapse .

[0097] Define nodes The degree of local spatial coupling collapse with its top K = 5 spatial nearest neighbors (sorted by historical correlation coefficients) within its cluster: ; in, For nodes The historical Spearman correlation coefficient benchmark established with node j during the post-completion stabilization period. This represents the real-time Spearman correlation coefficient at the current moment. Represents a node The set of the top K spatial nearest neighbors within the spatial monitoring cluster. The physical meaning of a node The degree of spatial correlation collapse of nodes coupled with surrounding mechanical forces; the larger the value, the more significant the collapse. The more disconnected the response patterns of surrounding nodes are in space, the more crucial the distinction between "single-point hardware failure" and "real structural damage" becomes.

[0098] Case 1: Comparison of historical and real-time correlation coefficients between ACC-09 and its 5 nearest neighbors: Historical baseline: 0.91, 0.88, 0.85, 0.82, 0.79; Real-time correlation coefficient: 0.12, 0.09, 0.08, 0.07, 0.06; Subtraction yields: 0.79, 0.79, 0.77, 0.75, 0.73.

[0099] ; Component 3: Global Spatial Synchronization Degree .

[0100] Define the cluster Global space synchronization level: ; This index is the global average of the correlation coefficient residuals of all nodes within the cluster. Its physical meaning is the consistency of the synchronous response of the entire space monitoring cluster: if all nodes exceed the synchronization limit and... This indicates that the overall spatial coupling relationship of the structure has not fundamentally changed, and it belongs to global environmental intrusion; if A significant increase indicates a regional spatial correlation break within the cluster, corresponding to local real damage.

[0101] Cluster in Case 1 There are 16 nodes in total. The residuals of the correlation coefficients between all nodes are less than 0.03. The overall mean is: ; By fusing the above three physical quantities through nonlinear weighting with adaptive balance coefficients, a single-point-space joint heterogeneous source discrimination function is constructed: ; Among them, the dynamic balance coefficient , , satisfy . The physical meaning is: when there are more nodes that exceed the limit (p(t) increases), the probability of a single point anomaly being regarded as a global event decreases, and the weight is adjusted down. The physical meaning is: when the global synchronous collapse degree The larger the value, the higher the confidence level of local spatial coupling collapse, and the higher the weight is accordingly. A fallback weight is provided for the global synchronization term to ensure the normalization of the three terms. This adaptive mechanism ensures that the discriminant function can provide highly reliable heterogeneous source stripping results under different spatial response modes.

[0102] Let the percentage of nodes exceeding the limit in cluster Ω at the current time be: ; in, This is an indicator function. , , The value of is adaptively adjusted according to the overall proportion of nodes exceeding the limit within the cluster, p(t).

[0103] In this case, the proportion of nodes exceeding the limit within the cluster ; The adaptive update rule for the dynamic balance coefficient is as follows: ; ; ; ; C02: Self-diagnostic logic tree based on the three-component form of the discriminant function.

[0104] Real-time calculation of each node And for the three components , and The value forms are analyzed jointly to form the following three-branch hard decision rule: Branch 1 (Independent hardware fault isolation of sensors): If a certain node satisfy (Single-point abnormally strong), and (Local spatial correlation cliff-like collapse), but (If global synchronization is intact and the correlation coefficients between other nodes do not decrease significantly relative to the historical baseline), then the node is determined to be... The anomaly stems from pure hardware failures such as zero-point drift of the sensor itself, probe detachment, or communication packet loss. The system soft-isolates the node, preventing it from triggering any structural alarms, and only records data quality logs and generates maintenance inspection work orders.

[0105] Branch 2 (Local Real Physical Damage Qualitative Analysis): If a certain node (or a subset of local nodes) has... Continuously increasing (multiple consecutive windows) ),and (The strong historical correlation between this node and its spatial nearest neighbor has significantly loosened), while (If the global spatial correlation of the cluster shows a regional break, i.e., the correlation coefficient between this node and other distant nodes decreases synchronously), then the node is determined to be... If the local area experiences a degradation in the actual structural bearing capacity or material fatigue yielding (such as concrete cracking or steel cable tension instability), the system triggers a regional structural anomaly warning and precisely locates the physical spatial position represented by that node.

[0106] Branch 3 (Qualitative Analysis of Global Extreme Environment Response): If the cluster exceeds p th proportional nodes (Global synchronization limit exceeded), and it applies to all nodes. (The spatial correlation between any node and its neighborhood did not substantially collapse), and at the same time If the difference between the global correlation coefficient matrix and the Frobenius norm of the historical benchmark is less than the threshold, indicating that the overall spatial coupling structure is intact, then it is determined that the infrastructure is being subjected to a global extreme environmental load intrusion within a legal range (such as a strong typhoon, extreme temperature changes, or a convoy of heavy vehicles crossing a bridge simultaneously). The system automatically filters all alarms, silently upgrades the monitoring level (increasing the sampling rate and storage redundancy) without generating any audible or visual alarms, and automatically downgrades back to the normal monitoring state after the environmental event subsides.

[0107] In this case, the threshold settings are as follows: , , , , , , pth = 0.8.

[0108] The meanings of each threshold are as follows: and These are the strong and moderate thresholds for the multiple of a single point of anomaly exceeding the standard, used to distinguish different degrees of severity of a single point of anomaly; , , These are high, medium, and low thresholds for the degree of local spatial coupling collapse, used to determine whether the degradation of the correlation between a node and its neighboring space has reached a state of fracture, loosening, or integrity. and These are the high and low thresholds for global spatial synchronization, used to distinguish between regional damage and global environmental response; This is the threshold for the proportion of nodes exceeding the limit within the cluster, used to determine whether an anomaly is a global event.

[0109] Analysis of the three-component form in Case 1: (Indicates a strong single-point violation of the limit) (This indicates a collapse of local correlation). (Indicates global synchronization is complete). All conditions for branch one are met. , , The system determined that ACC-09 was a sensor hardware failure. It then soft-isolated the sensor, blocking the alarm and generating a maintenance work order recommending checking the ACC-09 probe's installation torque and grounding. On-site inspection confirmed that the probe's fixing bolts were loosened by approximately 1.5 turns; after retightening, the fault was resolved, preventing another high-level false alarm.

[0110] Case 2: Damage to the PE sheath of the LZ12 cable stay.

[0111] ACC-08 was marked as a suspicious node. (and has been in 7 consecutive windows) (and continues to increase).

[0112] calculate The historical correlation coefficients with the 5 nearest neighbors were 0.94, 0.92, 0.90, 0.88, and 0.85, respectively; the real-time correlation coefficients were 0.42, 0.38, 0.35, 0.32, and 0.30, respectively; subtracting these values ​​yielded 0.52, 0.54, 0.55, 0.56, and 0.55.

[0113] ; calculate :cluster The correlation coefficients of the internal node pairs have decreased significantly, and the full matrix calculation shows that: ; Satisfy the condition of branch two ( And it continues to increase. , The system determined that actual structural degradation had occurred in the local area (1 / 4 section of the main beam) where ACC-08 was located. The system triggered an orange alert and pinpointed the location to the south-side cable-stayed bridge LZ12. A bridge inspection vehicle confirmed that the PE sheath was damaged and that transient stress waves were generated by friction between the steel strands, providing a warning 63 days in advance.

[0114] Case 3: A powerful typhoon passes through.

[0115] During the typhoon's passage, all 94 nodes ,part Taking ACC-08 as an example, the real-time correlation coefficient is almost identical to the historical benchmark: ; ; All nodes All less than 0.1. Cluster Global space synchronization level: ; The percentage of nodes exceeding the limit within the cluster, p(t), is 1.0 (100% exceeding the limit). This satisfies all three branch conditions (percentage of nodes exceeding the limit > 80%). , The system was identified as experiencing a global extreme environmental load intrusion. The system silently upgraded its monitoring level (temporarily increasing the sampling rate to 256Hz), blocking approximately 380 potentially invalid alarms. It automatically downgraded back to normal status 36 hours after the typhoon passed.

[0116] S104: The high-dimensional fusion feature vector generated in step S102 is used as real-time feature points and inserted into each tree structure respectively. The collaborative offset score of each tree is calculated, and the anomaly score at the current time is obtained by combining the collaborative offset scores of each tree. The time evaluation sliding window is maintained, the moving average of the anomaly score is calculated, and the first and second derivatives of the anomaly score are obtained. The judgment logic is executed according to the anomaly score, the moving average, and the value form of the first and second derivatives, and the judgment result is output.

[0117] Traditional fixed-parameter models are prone to failure when faced with the aforementioned legitimate extreme value environments. This step is specifically designed for endless streaming data, introducing a non-parametric anomaly detection model that requires no fixed threshold or offline training and supports online incremental updates. The specific sub-steps are as follows: D01: Dynamic Boundaries and Unsupervised Flow Modeling.

[0118] The high-dimensional fusion feature vector is used as F t Input, construct from N treeAn isolated forest composed of robust random splitting trees (RRCF). When constructing tree nodes, the range of each dimension of the current feature space is calculated to form bounding boxes. This is done using probability... ( The first node in the current tree node sample set Randomly select the cutting dimension (the width of the value range of each feature dimension) and in... Cutting points are selected using a uniform distribution within the space. This mechanism exhibits extremely high sensitivity to outliers that are free in multidimensional space.

[0119] Specific example: Maintaining an online sample buffer (reservoir sampling) of size N = 10,000. Constructing N tree =200 Robust Randomized Cross-Sectional Trees (RRCF). Each tree's root node contains all samples, and the range of each of the 16 dimensions is calculated: ; In the formula, This is the set of all samples contained in the root node v.

[0120] With probability (+0.01 to prevent division by zero errors) Extract the cutting dimension (prioritize dimensions with large ranges), in [L i , R i Uniformly and randomly select cutting points within: ; in, These represent the current node. within, no. The minimum and maximum values ​​of each dimension.

[0121] Based on the sample Divide into left and right child nodes: ; Termination condition: The recursive construction reaches a leaf node sample count ≤ 30 or a tree height ≥ 20.

[0122] D02: Expected solution for CoDisp topology depth fraction.

[0123] When real-time feature points When data flows in, there is no need to retrain on the entire historical data; only the calculation of the data will be performed. The amount of topological change (i.e., collaborative offset fraction) caused by inserting into the current tree structure: ; Where |S| represents the total number of root samples (fixed at 10,000), |S Y | represents the number of samples in the subtree rooted at node y. The physical meaning of this formula: If Falling into a rare partition region (denominator | S) Y |Small, and| S | - | S Y The ratio increases dramatically when the area is large; however, it decreases when the area is dense. (System calculation) Expected value of co-offset fraction across 200 trees After removing the extreme values ​​from the top 10% and bottom 10%, we get: ; D03: Differential Evolution and Temporal Morphological Qualitative Classification.

[0124] System maintenance time assessment sliding window (width W = 500, corresponding to 2000 seconds), calculate the sliding average: ; And its first and second derivatives (using finite differences) are used to determine the trend direction and pulse width, respectively: ; ; In the above formulas, The sampling interval is denoted as .

[0125] If the real-time CoDisp exhibits a steep, unipolar, high-amplitude pulse with extremely rapid alternation of positive and negative derivatives and no significant jump in the moving average MA, the system is characterized as an instantaneous "point anomaly" (such as a single anomalous gust of wind) and is automatically filtered out.

[0126] If the real-time CoDisp crosses the extreme value red line for multiple consecutive periods, the first derivative remains positive, and the moving average MA shows a continuous step divergence without falling back, then the system is characterized as a dangerous "trend shift." This indicates that the structural dynamics have undergone irreversible and severe physical degradation, at which point the system triggers the highest priority structural degradation alarm work order.

[0127] The following is a specific example: Type A – Instantaneous Point Anomaly (e.g., heavy-load vehicle fleet): E[CoDisp] pulse peak value 5.23 (baseline 0.85), lasting only one window (4 seconds), first derivative alternation frequency 8 times / second, second derivative alternation frequency 12 times / second, MA increment only 0.06 < 0.1. The system determines this as an instantaneous point anomaly and automatically filters it.

[0128] Type B – Trend Drift (e.g., filament erosion): E[CoDisp] continuously increases from 2.87 to 7.84, MA increases from 1.87 to 5.94, with a cumulative increment of 4.07 > 0.5. First derivative least squares fitting slope: ; Where N d = 45 (45 days of data points), t i For days, MA i This is the daily moving average. The system triggered a trend drift warning; after the index was changed, E[CoDisp] fell back to the baseline.

[0129] In the above determination, 0.1 and 0.5 are the single-step increment threshold and cumulative increment threshold for MA, respectively. Specifically, 0.1 (corresponding to approximately 10% of the baseline MA) is used to determine whether a single event causes a significant jump in MA: if the MA increment is less than this threshold, it indicates that the anomaly is a transient pulse and does not change the long-term trend; 0.5 (corresponding to approximately 50% of the baseline MA) is used to determine whether an irreversible step drift has occurred in MA: if the cumulative increment of MA exceeds this threshold, it indicates that the anomaly has been accumulating and causing baseline shift, corresponding to physical degradation of the structure. These thresholds are empirical values, all set with the baseline MA value as a reference. The range is 5%~15% of the baseline (single-step increment threshold) and 2~4 times the baseline standard deviation (cumulative increment threshold). The specific values ​​can be adaptively configured according to the dynamic characteristics of the monitored structure and the fluctuation range of historical data.

[0130] Type C - Periodic environmental disturbances (e.g., daily afternoon thermal deflection in summer): E[CoDisp] exhibits a 24-hour periodic pattern, with daily peak values ​​not accumulating and increasing. After the system automatically learns this pattern, it only increases the monitoring level and suppresses alarm output.

[0131] This embodiment also provides a computer system, including a memory, a processor, and a computer program stored in the memory and executable on the processor, wherein the processor executes the computer program to implement the method steps of the above embodiment.

[0132] This embodiment also provides a computer storage medium on which a computer program is stored, which, when executed by a processor, implements the method steps in the above embodiment.

[0133] Obviously, the embodiments described above are only some embodiments of the present invention, and not all embodiments. All other embodiments obtained by those skilled in the art based on the embodiments of the present invention without inventive effort are within the scope of protection of the present invention.

Claims

1. A method for time-frequency feature extraction and fault diagnosis of multi-source heterogeneous sensors, characterized in that, include: S1: Acquire multi-source heterogeneous sensor signals from infrastructure deployment; perform hardware-level time synchronization and unified timing on the signals, and dynamically compensate for clock drift rate; then fuse signals with different sampling rates by using zero-phase frequency domain digital filtering downsampling and spline interpolation upsampling with natural boundary constraints, and use adaptive empirical mode decomposition with extreme value mirror extension to remove slowly varying trend terms; output the reconstructed multi-source clean signal sequence. S2: Based on the multi-source clean signal sequence, the signal sequence is synchronously truncated through a non-overlapping time sliding window, and damage-sensitive features are extracted from both the time domain and the frequency domain: in the time domain, high-order dimensionless statistics sensitive to sudden changes are extracted, and in the frequency domain, the spectral centroid frequency and wavelet packet energy features are extracted. A multi-dimensional fused feature space is constructed by splicing and fusing time-domain and frequency-domain features; based on the abnormal quantities of each dimension in the feature space, a multi-order time-frequency collaborative damage sensitivity index with a multiplicative collaborative excitation mechanism is constructed, and the multi-dimensional abnormal signals are fused into a single scalar through multiplicative coupling; the current node status is determined based on the comparison result of the sensitivity index and a preset threshold, and a normal node or suspicious node label is output. S3: Based on nodes marked as suspicious, construct a single-point-space joint heterogeneous source discrimination function containing three components: single-point anomaly exceeding the standard multiple, local spatial coupling collapse degree, and global spatial synchronization degree. Execute three-branch hard judgment logic according to the combination mode of the three component values, and output independent hardware fault isolation instructions, local real damage location warnings, or global extreme value environmental event silent monitoring instructions respectively.

2. The method for time-frequency feature extraction and fault diagnosis of multi-source heterogeneous sensors according to claim 1, characterized in that, The process also includes step S4: using the high-dimensional fusion feature vector generated in step S2 as real-time feature points, inserting them into each tree structure, calculating the collaborative offset score of each tree, and combining the collaborative offset scores of each tree to obtain the anomaly score at the current time; maintaining a time evaluation sliding window, calculating the moving average of the anomaly score, and obtaining the first and second derivatives of the anomaly score; executing the judgment logic based on the anomaly score, the moving average, and the values ​​of the first and second derivatives, and outputting the judgment result.

3. The method for time-frequency feature extraction and fault diagnosis of multi-source heterogeneous sensors according to claim 1, characterized in that, Step S1 specifically includes the following sub-steps: S11: Real-time acquisition of multi-source heterogeneous signals output by multi-source heterogeneous sensors deployed at key stress nodes of infrastructure; S12: Through the interaction of the timing protocol between the master node and the terminal sensor node, calculate the one-way transmission delay and instantaneous clock offset under the asymmetric network path, and estimate the relative drift rate using the clock offset sequence of multiple consecutive timing cycles. Based on this, dynamically adjust the local clock and embed a unified global absolute timestamp for each sensor data packet. S13: For high-frequency dynamic signals with sampling rates higher than the reference frequency, a multi-order Butterworth low-pass anti-aliasing filter is constructed and its amplitude-frequency response function is calculated. The signal is then downsampled after being filtered using zero-phase filtering technology, i.e., the signal passes through the filter in the forward direction, is time-reversed, passes through the same filter again, and is reversed back to the original timing sequence. For low-frequency static signals with sampling rates lower than the reference frequency, upsampling is performed using multiple spline interpolation with natural boundary constraints. A tridiagonal matrix is ​​constructed and the polynomial coefficients are obtained by solving the Thomas algorithm. The known low-frequency data points are substituted into the polynomials corresponding to the polynomial coefficients to calculate the interpolation results at the reference sampling time. S14: For the signal after frequency unification, the extreme value mirror extension method is used to virtually construct the extreme point closed loop at both ends of the signal. The signal is used as the current residual term. Iterative sieving is performed by adaptive empirical mode decomposition. The stopping criterion is that the standard deviation of at least two consecutive iterations is lower than the preset dynamic threshold. The intrinsic mode functions are extracted from the residual term in sequence. The above sieving process is repeated until the residual term becomes a monotonic function or the number of poles is insufficient to construct the envelope. The decomposition is terminated when the residual term becomes a monotonic function or the number of poles is insufficient to construct the envelope. The residual term is the slowly varying trend term. After removing the slowly varying trend term, the reconstructed multi-source pure signal sequence is output.

4. The method for time-frequency feature extraction and fault diagnosis of multi-source heterogeneous sensors according to claim 1, characterized in that, Step S2 specifically includes the following sub-steps: S21: Synchronously extract the multi-source pure signal sequence through a non-overlapping time sliding window; S22: In the time domain dimension, extract kurtosis, margin factor, skewness, impulse factor and waveform factor from each window sequence to form a high-order dimensionless parameter group in the time domain; S23: In the frequency domain, perform a fast Fourier transform on each window sequence to extract the spectral centroid frequency, and use a bottom-up pruning strategy driven by orthogonal wavelet packet decomposition and Shannon entropy cost function to obtain the optimal basis. Calculate the energy integral and energy ratio of the reconstruction coefficients of each frequency band under the optimal basis. After splicing and normalizing the time domain parameter group, spectral centroid frequency and frequency band energy ratio, generate a high-dimensional fusion feature vector to complete the construction of the multi-dimensional fusion feature space. S24: On the multi-dimensional fusion feature space, construct a multi-order time-frequency collaborative damage sensitivity index, fuse the abnormal quantities of each dimension in the feature space into a single scalar through multiplicative coupling, and compare the index with a preset multi-level early warning threshold. If the sensitivity index exceeds the limit, mark the node as a suspicious node and output it to step 3; otherwise, output the normal state.

5. The method for time-frequency feature extraction and fault diagnosis of multi-source heterogeneous sensors according to claim 1, characterized in that, Step S3 specifically includes the following sub-steps: S31: Based on the historical pure signal sequence of all nodes, the median of each dimension of each node is taken to form a baseline representative vector; hierarchical clustering is performed according to the rank correlation coefficient between each node to obtain spatial monitoring clusters; S32: Based on the current high-dimensional fusion feature vector of each suspicious node, the increment of the current feature vector relative to the historical baseline level is calculated to obtain the single-point anomaly exceeding multiple of each suspicious node. S33: Calculate the degree of local spatial coupling collapse between each suspicious node and its preset number of spatial nearest neighbors within its cluster; S34: Calculate the global spatial synchronization degree of the cluster to which each suspicious node belongs; S35: The single-point anomaly exceeding the standard multiple, the degree of local spatial coupling collapse, and the degree of global spatial synchronization are fused by nonlinear weighting with adaptive dynamic balance coefficients to construct a single-point-space joint heterogeneous source discrimination function; wherein the dynamic balance coefficient is adaptively selected according to the proportion of nodes exceeding the limit within the cluster; the three-component value forms of the discrimination function are jointly analyzed, and a three-branch hard decision logic is executed according to the combination mode of the three-component value forms, and an independent hardware fault isolation instruction, a local real damage location warning, or a global extreme value environmental event silent monitoring instruction is output according to the decision result.

6. The method for time-frequency feature extraction and fault diagnosis of multi-source heterogeneous sensors according to claim 5, characterized in that, The three-branch hard decision logic includes: The first branch, when the single-point anomaly exceeds the limit multiple... The degree of local spatial coupling collapse exceeds Furthermore, the global space synchronization level is lower than that of other systems. When the node malfunction is determined to be due to sensor hardware failure, soft isolation is performed and an operation and maintenance inspection work order is generated; The second branch, when the single-point anomaly exceeds the limit multiple... The degree of local spatial coupling collapse exceeds And the global space synchronization level exceeds When the node is located, it is determined that there is actual structural damage in the local area where the node is located, triggering a regional structural anomaly warning and outputting the spatial location of the node. The third branch considers nodes within a cluster that exceed the threshold for the proportion of nodes exceeding the limit, where the single-point anomaly exceeds the limit multiple. The degree of local spatial coupling collapse of all nodes is lower than Furthermore, the global space synchronization level is lower than that of other systems. When it is determined that the infrastructure is experiencing a global extreme environmental load, all alarms are automatically filtered and the monitoring level is upgraded; in, and These are the strong and moderate thresholds for the multiple of a single-point anomaly exceeding the standard, respectively. , , These represent high, medium, and low thresholds for the degree of local spatial coupling collapse. and These represent the high and low thresholds for global spatial synchronization, respectively.

7. The method for time-frequency feature extraction and fault diagnosis of multi-source heterogeneous sensors according to claim 2, characterized in that, Step S4 specifically includes the following sub-steps: S41: Construct an isolated forest composed of multiple robust random splitting trees. When constructing each tree, calculate the range of each dimension of the current feature space at each node, extract the cutting dimension with the range as the weight, and randomly select the cutting point uniformly within the value range of the dimension to split the node. S42: Insert the real-time feature points into each tree, calculate the collaborative offset score of each tree, and combine the collaborative offset scores of all trees to obtain the anomaly score at the current time. A maintenance time assessment sliding window is used to calculate the moving average of the anomaly scores, and the first and second derivatives of the anomaly scores are calculated using the finite difference method. S43: When the abnormal score shows a steep unipolar pulse, the first and second derivatives alternate between positive and negative, and the moving average does not show a significant jump, it is judged as an instantaneous point anomaly and filtered; when the abnormal score exceeds the preset threshold for multiple consecutive cycles, the first derivative is continuously positive, and the moving average shows a continuous step divergence trend, it is judged that the structural dynamic characteristics have undergone irreversible degradation, triggering a structural degradation alarm.

8. The method for time-frequency feature extraction and fault diagnosis of multi-source heterogeneous sensors according to claim 1 or 4, characterized in that, In step S2, the multi-order time-frequency coordinated damage sensitivity index It can be obtained through the following formula: ； In the formula, Indicates the kurtosis rate in the time domain. Indicates the rate of increase of the time-domain margin factor. This represents the absolute offset rate of time-domain skewness. This represents the frequency domain spectral centroid offset rate. This represents the energy transition rate of the high-frequency wavelet packet band; This is the synergistic incentive index for each indicator.

9. The method for time-frequency feature extraction and fault diagnosis of multi-source heterogeneous sensors according to claim 1 or 5, characterized in that, In step S3, the single-point anomaly exceedance multiple is obtained by the following formula: ； The degree of local spatial coupling collapse is obtained by the following formula: ； The global spatial synchronization level is obtained using the following formula: ； The single-point-spatial joint heterogeneous source discrimination function is obtained by the following formula: ； In the above formulas, Suspicious node The single-point anomaly exceeding the limit at the current time t. The Euclidean distance between the feature vectors of the current window and the previous window; An adaptive threshold that is dynamically updated; Suspicious node The degree of local spatial coupling collapse with its K nearest neighbors within its cluster. Suspicious node The historical Spearman correlation coefficient benchmark established with its nearest neighbor node j during the stable period. This represents the real-time Spearman correlation coefficient at the current moment. Suspicious node Cluster The degree of global space synchronization, The historical Spearman correlation coefficient between any nodes u and v within the cluster during the stable period serves as a benchmark. This represents the real-time Spearman correlation coefficient at the current moment. , , This is the dynamic balance coefficient.

10. A system for extracting time-frequency features and diagnosing faults in multi-source heterogeneous sensors, comprising a memory, a processor, and a computer program stored in the memory and executable on the processor; characterized in that, When the processor executes the computer program, it implements the steps of the method according to any one of claims 1 to 9.

Images