A filter performance test method and system based on intelligent sensors

Abstract

The application relates to the technical field of filter detection, and discloses a filter performance test method and system based on an intelligent sensor, which comprises the following steps: acquiring a standard sampling sequence and a high-frequency sampling sequence containing inlet pressure, outlet pressure, flow and vibration acceleration, performing resampling processing to obtain a structured signal set; performing phase compensation according to the structured signal set to obtain a smooth zero-phase sequence; performing pressure difference index calculation according to the smooth zero-phase sequence to obtain a comprehensive pressure difference index; performing clogging point positioning according to the comprehensive pressure difference index to obtain a clogging starting point; performing abnormal feature elimination according to the clogging starting point to obtain a final vector set; performing risk probability mapping according to the final vector set to obtain a risk probability curve; and performing failure time prediction according to the risk probability curve to obtain a critical failure time and send a maintenance warning signal. The method can accurately identify clogging risks and give early warnings.

Description

Technical Field

[0001] This invention relates to the field of filter testing technology, and in particular to a filter performance testing method and system based on intelligent sensors. Background Technology

[0002] Currently, even minor fluctuations in pressure and flow during filter operation can gradually accumulate and cause localized blockages. These blockages, in turn, exacerbate abnormal vibrations and accelerate component wear. This chain reaction from data fluctuations to physical damage is difficult to capture using conventional methods. Accurately capturing the dynamic changes of multi-point data throughout the entire operation process and predicting performance degradation and blockage risks based on these changes has become a challenge for improving testing accuracy. This requires the use of intelligent sensors to collect pressure, flow, and vibration signals in real time, and through signal processing and trend analysis, to achieve a precise assessment of the filter's condition.

[0003] In one existing technology, pressure sensors deployed at the filter inlet and outlet collect pressure signals at a fixed sampling rate, calculate the pressure difference between the inlet and outlet, and compare it with a preset fixed pressure difference threshold. When the pressure difference exceeds the threshold, the system issues a blockage alarm. Simultaneously, a flow sensor collects flow signals at a lower sampling rate to help determine the filter's flow capacity. This method relies solely on the instantaneous value of the pressure difference for judgment, excluding vibration signals. The pressure difference calculation is based on a single sensor pair, without redundancy verification of multi-channel data. When the pressure difference fluctuates around the threshold, the system cannot distinguish whether it is caused by localized blockage or changes in operating conditions, and the alarm timing often lags behind the actual occurrence of blockage. Existing technologies only perform threshold comparisons on fixed parameters, ignoring the dynamic correlation between pressure, flow rate, and vibration, resulting in the inability to detect subtle early signs of blockage.

[0004] Therefore, existing technologies lack sufficient accuracy in identifying congestion risks and suffer from delayed early warning responses. Summary of the Invention

[0005] This invention provides a filter performance testing method and system based on intelligent sensors to accurately identify clogging risks and provide early warnings.

[0006] In a first aspect, to address the aforementioned technical problems, the present invention provides a filter performance testing method based on intelligent sensors, comprising: Acquire standard and high-frequency sampling sequences containing inlet pressure, outlet pressure, flow rate, and vibration acceleration; Based on the standard sampling sequence, the sampling frequency is adjusted to obtain a structured signal set; Based on the structured signal set, cross-correlation calculation is performed to obtain the signal delay. Based on the signal delay, phase compensation is performed to obtain a smooth zero-phase sequence. Based on the smoothed zero-phase sequence, pressure difference constants are eliminated to obtain a comprehensive pressure difference index; Based on the comprehensive pressure difference index, the periodic interference frequency is obtained by performing spectral analysis on the high-frequency sampling sequence, and the noise type is obtained by noise feature matching on the high-frequency sampling sequence. The blockage initiation point is determined by combining the variation characteristics of the comprehensive pressure difference index, the periodic interference frequency, and the noise type. Based on the blockage initiation point and the smoothed zero-phase sequence, abnormal features are removed to obtain the final vector set; Curve fitting is performed on the final vector set to obtain the predicted pressure difference curve. The deviation of the predicted pressure difference curve is corrected using the comprehensive pressure difference index to obtain the corrected pressure difference curve. The risk probability curve is obtained by combining the corrected pressure difference curve with the pre-constructed risk mapping function. Based on the risk probability curve, failure time is predicted to obtain the critical failure time. Based on the critical failure time, a maintenance alarm signal is issued and the replacement time is displayed.

[0007] Secondly, the present invention provides a filter performance testing system based on intelligent sensors, comprising: The data acquisition module is used to acquire standard sampling sequences and high-frequency sampling sequences containing inlet pressure, outlet pressure, flow rate, and vibration acceleration. The frequency adjustment module is used to adjust the sampling frequency according to the standard sampling sequence to obtain a structured signal set; The phase compensation module is used to perform cross-correlation calculation based on the structured signal set to obtain the signal delay, and to perform phase compensation based on the signal delay to obtain a smooth zero-phase sequence. The index calculation module is used to remove pressure difference anomalies based on the smoothed zero-phase sequence to obtain a comprehensive pressure difference index. The blockage location module is used to perform spectral analysis on the high-frequency sampling sequence to obtain the periodic interference frequency based on the comprehensive pressure difference index, and to perform noise feature matching on the high-frequency sampling sequence to obtain the noise type. The blockage initiation point is determined by combining the variation characteristics of the comprehensive pressure difference index, the periodic interference frequency, and the noise type. The anomaly removal module is used to remove anomaly features based on the blockage initiation point and the smoothed zero-phase sequence to obtain the final vector set; The risk mapping module is used to perform curve fitting on the final vector set to obtain the predicted pressure difference curve, use the comprehensive pressure difference index to correct the deviation of the predicted pressure difference curve to obtain the corrected pressure difference curve, and combine the pre-constructed risk mapping function to perform risk probability mapping on the corrected pressure difference curve to obtain the risk probability curve. The output module is used to predict the failure time based on the risk probability curve, obtain the critical failure time, and issue a maintenance alarm signal and display the replacement time based on the critical failure time.

[0008] Compared with the prior art, the present invention has the following beneficial effects: (1) This invention collects standard sampling sequences and high-frequency sampling sequences of inlet pressure, outlet pressure, and vibration acceleration. It performs spectral analysis on the standard sampling sequences and adaptively increases the sampling frequency according to the dominant frequency. Then, it obtains a structured signal set through filtering and splicing. Since the frequencies of different sensors can be dynamically adjusted according to the dominant frequency of the signal, it not only preserves high-frequency details but also avoids data redundancy. This solves the problem that the fixed sampling rate of the prior art is difficult to capture dynamic changes, and improves the adaptability of data acquisition and signal fidelity.

[0009] (2) When the pressure difference exceeds the threshold, this invention calls the vibration acceleration in the high-frequency sampling sequence for spectrum analysis, extracts the periodic interference frequency, and constructs a Chebyshev filter to filter out equipment vibration noise in the pressure signal. Then, it identifies fluid impact noise through peak factor and skewness, and calculates the envelope energy sequence and screens energy mutation points to locate the blockage initiation time. Because the vibration interference is accurately suppressed, the characteristics of fluid impact noise are highlighted, solving the problem of misjudgment or delayed location of blockage caused by vibration interference in the prior art, and realizing accurate identification of the blockage initiation point; verified by bench simulation test, such as Figure 2 As shown, the present invention successfully identified the blockage initiation point (t=5.20s) under simulated blockage conditions, which is about 3.3 seconds earlier than the traditional fixed threshold method (alarm time about t=8.50s), with a relative advance of 38%. During the entire test cycle, no false alarms were caused by non-blockage factors such as pipeline eddies and bubbles, and the false alarm rate was 0%.

[0010] (3) Based on the extracted multidimensional statistical features, this invention uses the least squares method to fit and predict the pressure difference curve. It calculates the deviation compensation amount using the measured comprehensive pressure difference index and superimposes the correction curve. Then, it is converted into a risk probability curve through a risk mapping function. Finally, the critical failure time is obtained by extrapolation based on the risk probability curve, and an early warning is issued. Because the dynamic correction mechanism makes the prediction curve fit the actual working conditions, the risk probability mapping intuitively transforms the pressure difference change into the degree of performance degradation, which solves the problem of the delayed early warning time in the prior art and realizes the early warning of blockage risk and the scientific decision-making of maintenance timing. Attached Figure Description

[0011] Figure 1 This is a schematic diagram of the filter performance testing method based on intelligent sensors provided in the first embodiment of the present invention; Figure 2 This is a schematic diagram of the blockage initiation point location provided in the first embodiment of the present invention; Figure 3 This is a schematic diagram of the filter performance testing system based on intelligent sensors provided in the second embodiment of the present invention. Detailed Implementation

[0012] The technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of the present invention, and not all embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of the present invention.

[0013] Reference Figure 1 The first embodiment of the present invention provides a filter performance testing method based on intelligent sensors, including the following steps: S11, acquire standard sampling sequences and high-frequency sampling sequences containing inlet pressure, outlet pressure, flow rate and vibration acceleration; S12, adjust the sampling frequency according to the standard sampling sequence to obtain a structured signal set; S13, perform cross-correlation calculation based on the structured signal set to obtain the signal delay, and perform phase compensation based on the signal delay to obtain a smooth zero-phase sequence; S14. Based on the smoothed zero-phase sequence, pressure difference constants are eliminated to obtain the comprehensive pressure difference index. S15, based on the comprehensive pressure difference index, perform spectral analysis on the high-frequency sampling sequence to obtain the periodic interference frequency, and perform noise feature matching on the high-frequency sampling sequence to obtain the noise type. Combine the variation characteristics of the comprehensive pressure difference index, the periodic interference frequency, and the noise type to determine the blockage initiation point. S16, based on the blockage initiation point and the smoothed zero-phase sequence, perform abnormal feature removal to obtain the final vector set; S17, curve fitting is performed on the final vector set to obtain the predicted pressure difference curve, the deviation of the predicted pressure difference curve is corrected using the comprehensive pressure difference index to obtain the corrected pressure difference curve, and the risk probability is mapped to the corrected pressure difference curve by combining the pre-constructed risk mapping function to obtain the risk probability curve. S18. Based on the risk probability curve, predict the failure time to obtain the critical failure time, and issue a maintenance alarm signal and display the replacement time based on the critical failure time.

[0014] In step S11, a standard sampling sequence and a high-frequency sampling sequence containing inlet pressure, outlet pressure, flow rate, and vibration acceleration are obtained.

[0015] Specifically, multiple pressure sensors, flow sensors, and vibration acceleration sensors are deployed in the inlet and outlet pipes of the filter, respectively. For example, three pressure sensors and two flow sensors are deployed at the inlet, and two pressure sensors and three flow sensors are deployed at the outlet. An acceleration sensor is also installed on the outer wall of both the inlet and outlet pipes to collect vibration acceleration signals. All sensors acquire data at an initial sampling rate of 10 Hz, and the acquisition time of each sensor is synchronized through a network time protocol, ensuring that the time error between different sensors is less than one millisecond. Each sensor's data acquisition is accompanied by a precise timestamp, forming a raw sampling sequence. The system simultaneously acquires data at two sampling rates: a standard sampling sequence at 10 Hz and a high-frequency sampling sequence at 10 kHz. The high-frequency sampling sequence only includes inlet pressure, outlet pressure, and vibration acceleration signals; flow signals are not acquired at high frequency. All acquired data is transmitted to the central processing unit via industrial Ethernet and stored in timestamp order, forming a raw sampling sequence set containing timestamps. The standard sampling sequence is used for routine monitoring and trend analysis, while the high-frequency sampling sequence is used for transient feature capture when the pressure difference exceeds a threshold. This step provides a complete data foundation for subsequent adaptive sampling frequency adjustment and blockage initiation point location.

[0016] In step S12, the sampling frequency is adjusted according to the standard sampling sequence to obtain a structured signal set, including: Based on the standard sampling sequence, timestamp alignment is performed to obtain a synchronous sampling sequence, and based on the synchronous sampling sequence, spectrum analysis is performed using the Fast Fourier Transform algorithm to obtain the initial main frequency; When the initial main frequency exceeds a preset main frequency threshold, the preset sampling frequency is increased to twice the initial main frequency to obtain a high-density sampling sequence; Based on the high-density sampling sequence and the synchronous sampling sequence, interpolation calculation is performed using the Lagrange interpolation algorithm to obtain the resampling sequence; When the initial clock frequency does not exceed the clock frequency threshold, the synchronous sampling sequence is directly used as the resampling sequence; Based on the resampled sequence, the data is spliced ​​together according to the time dimension to obtain a structured signal set.

[0017] Specifically, a timestamp alignment operation is performed on the acquired standard sampling sequence. Data collected by different sensors in the standard sampling sequence carries its own timestamp. Due to slight differences in sensor startup times, timestamp alignment is needed to match the data points collected by each sensor at the same physical moment. Specifically, using the sampling time corresponding to the highest standard sampling frequency in the system as a benchmark, linear interpolation is performed on the data sequences of other sensors to ensure that all sensors have the same time point sequence on the time axis, thus obtaining a synchronized sampling sequence. Each time point in this synchronized sampling sequence contains synchronized measurements of inlet pressure, outlet pressure, flow rate, and vibration acceleration.

[0018] Spectral analysis was performed on the pressure signal in the synchronous sampling sequence. A five-minute segment of the pressure signal from the synchronous sampling sequence was taken, and a Fast Fourier Transform (FFT) was performed on this segment to convert the time-domain signal into a frequency-domain distribution. The frequency component with the largest amplitude in the frequency domain was calculated to obtain the initial dominant frequency. The preset dominant frequency threshold was determined by statistical spectral analysis of pressure signals collected under historical normal operating conditions. Specifically, one thousand pressure signal segments under normal operating conditions were collected, their dominant frequencies were calculated, and these dominant frequency values ​​were sorted in ascending order. The 95th percentile was taken as the dominant frequency threshold. If the calculated initial dominant frequency is greater than this threshold, it indicates that there are high-frequency variation components in the signal, and the sampling frequency needs to be increased to avoid aliasing; if the initial dominant frequency is less than or equal to this threshold, no increase is needed.

[0019] When the initial master frequency exceeds a preset master frequency threshold, the sampling frequency needs to be increased to twice the initial master frequency. Since the original standard sampling sequence has a sampling frequency of 10 Hz, it is impossible to directly obtain data points with higher frequencies. Therefore, a Lagrange interpolation algorithm is used to insert new data points between existing sampling points. The specific process of the Lagrange interpolation algorithm is as follows: for the time point requiring interpolation, one known sampling point is selected before and after that point, and Lagrange basis functions are constructed. Each basis function is a polynomial in the form of a product of its numerator and denominator. Each basis function is multiplied by the amplitude of the corresponding sampling point, and then summed to obtain the amplitude of the interpolation point. In this way, the data of each sensor channel in the synchronous sampling sequence is interpolated from the preset sampling frequency to a sampling frequency twice the initial master frequency, resulting in a high-density sampling sequence. This high-density sampling sequence only contains the signal channels whose frequency needs to be increased; for signal channels that do not exceed the threshold, their data are still retained in the synchronous sampling sequence. The preset sampling frequency is the standard sampling frequency at the initial acquisition, which is 10 Hz in this embodiment.

[0020] In this embodiment, the Lagrange interpolation is fixed to three-point interpolation, that is, for each time point to be interpolated, one known sampling point before and one before it are selected, and a quadratic Lagrange polynomial is constructed from three known points. Three-point interpolation can avoid unstable oscillations of higher-order polynomials, and at the same time, the amount of computation is small, making it suitable for embedded real-time processing.

[0021] The high-density sampling sequence is fused with the un-upgraded signal channels in the synchronous sampling sequence, and all signals are unified to the same sampling frequency using the Lagrange interpolation algorithm. The highest sampling frequency among all signals is selected as the target frequency, and the Lagrange interpolation upsampling operation is performed on signal channels with frequencies lower than this target frequency, so that all signal channels have the same sampling density on the time axis, resulting in a resampled sequence.

[0022] The resampled sequences are concatenated along the time dimension. Each signal channel in the resampled sequence is still stored independently, with each channel corresponding to a single column of data. The data concatenation operation merges inlet pressure, outlet pressure, flow rate, and vibration acceleration at the same timestamp into a single record. All timestamp records are arranged chronologically, forming a two-dimensional table structure. Each row corresponds to a time point, and each column corresponds to a type of sensor signal. This two-dimensional table is the structured signal set, used in subsequent steps.

[0023] It should be noted that the initial sampling rate of the standard sampling sequence is set according to the highest effective frequency of the filter pressure signal, ensuring that this highest effective frequency is less than half of the sampling rate, thus satisfying the Nyquist theorem. The initial dominant frequency is the true frequency component extracted based on the non-aliased sampling sequence. When the initial dominant frequency exceeds the preset dominant frequency threshold, it indicates that there are noteworthy high-frequency changes in the signal, but this frequency is still less than half of the sampling rate. At this time, the Lagrange interpolation is only used to refine the discrete time points to improve the time resolution of subsequent cross-correlation calculations, without changing the frequency domain characteristics of the signal. The interpolation adopts three-point Lagrange interpolation, that is, each interpolation point is determined by one original sampling point before and after it, in order to suppress the oscillation phenomenon of higher-order interpolation. Through the above processing, the resampled sequence obtains a finer time grid while retaining the authenticity of the original spectrum.

[0024] In step S13, cross-correlation calculation is performed based on the structured signal set to obtain the signal delay. Phase compensation is then performed based on the signal delay to obtain a smoothed zero-phase sequence, including: Based on the structured signal set, the cutoff frequency and attenuation amplitude are extracted by fast Fourier transform to obtain filtering parameters. Based on the filtering parameters, the structured signal set is convolved by the Butterworth filtering algorithm to obtain a preliminary filtering sequence. Based on the structured signal set and the preliminary filtering sequence, cross-correlation function is calculated to obtain a cross-correlation function sequence and a time shift sequence, and the time shift in the time shift sequence corresponding to the peak value in the cross-correlation function sequence is used as the signal delay. When the signal delay is greater than the preset signal delay threshold, phase compensation is performed using a zero-phase filtering algorithm based on the preliminary filtering sequence to obtain an initial zero-phase sequence. When the signal delay is not greater than the signal delay threshold, the preliminary filtering sequence is directly used as the initial zero-phase sequence; Based on the initial zero-phase sequence, a smoothed zero-phase sequence is obtained by performing boundary weighted smoothing using a preset Hanning window function.

[0025] Specifically, the cutoff frequency and attenuation amplitude are extracted from the structured signal set. The structured signal set is a two-dimensional table containing inlet pressure, outlet pressure, flow rate, and vibration acceleration. Each row corresponds to a time point, and each column corresponds to a sensor signal. The spectral distribution characteristics of the pressure signal column in this signal set are extracted. A continuous ten-minute pressure signal segment is taken, and a Fast Fourier Transform is performed to obtain the frequency domain amplitude distribution. In the frequency domain, a search is performed from low to high frequencies to find the frequency corresponding to the amplitude decreasing to 0.7 times the maximum amplitude. This frequency is taken as the cutoff frequency; simultaneously, the ratio of the amplitude at this frequency to the maximum amplitude is recorded as the attenuation amplitude. This yields a set of filtering parameters, including the cutoff frequency and attenuation amplitude. The preset adjustment range of the filtering parameters is determined through statistical analysis of historical normal operation data. Specifically, one thousand sets of pressure signals under normal operating conditions are collected, their cutoff frequencies and attenuation amplitudes are calculated for each, and the minimum value is taken as the lower limit, and the maximum value as the upper limit. Subsequent filtering parameters are dynamically determined within this range based on the actual signal.

[0026] Based on the filtering parameters, Butterworth filtering is performed on the structured signal set. The Butterworth filter is fixed at order four. First, the normalized cutoff frequency is calculated based on the cutoff frequency and sampling frequency. Then, the numerator and denominator coefficients of the fourth-order Butterworth filter are obtained through bilinear transformation; these coefficients are a set of fixed values. The filtering process uses convolution operations. For each column of signal in the structured signal set, sixteen points before and after each data point are taken, and a weighted sum is performed with the filter coefficients. The weighted summation method is to multiply the corresponding data point with the filter coefficient and then accumulate the sum to obtain the filtered output value at that point. After calculating this for all data points, a preliminary filtering sequence is obtained. This preliminary filtering sequence has the same time and signal dimensions as the structured signal set.

[0027] Calculate the signal delay. Take the pressure signal column from the structured signal set and the corresponding pressure signal column from the preliminary filtered sequence. Use the pressure signal column from the structured signal set as the reference sequence and the pressure signal column from the preliminary filtered sequence as the test sequence. Calculate the cross-correlation function between the reference sequence and the test sequence. Specifically, shift the test sequence relative to the reference sequence sequentially to the left and right, each shift being one sampling interval, covering a range from -100 to +100 sampling points. For each shift, calculate the sum of the dot products of the reference sequence and the shifted test sequence, then divide by the square root of the sequence length to obtain the cross-correlation value for that shift. The cross-correlation values ​​corresponding to all shifts constitute a cross-correlation function sequence, and all shifts constitute a time shift sequence. Search for the maximum value in the cross-correlation function sequence and record the shift in the time shift sequence corresponding to that maximum value; this shift is the signal delay. If the signal delay is positive, it indicates that the preliminary filtered sequence lags behind the structured signal set; if it is negative, it indicates that it leads.

[0028] The preset signal delay threshold is determined by the theoretical group delay of the filter. For a fourth-order Butterworth filter, the theoretical group delay is approximately twice the filter order divided by the cutoff frequency; 120% of this theoretical value is taken as the delay threshold. When the absolute value of the calculated signal delay exceeds this threshold, phase compensation is required; otherwise, the preliminary filtering sequence is used directly.

[0029] When the signal delay exceeds a preset signal delay threshold, a zero-phase filtering algorithm is applied to the initial filtered sequence for phase compensation. The specific process of the zero-phase filtering algorithm is as follows: the initial filtered sequence is input forward into a Butterworth filter to obtain the first filtering result; then, the time order of the first filtering result is reversed, and it is input into the same Butterworth filter again; finally, the second filtering result is reversed again to obtain the initial zero-phase sequence. After two forward and reverse filtering operations, the phase delays cancel each other out, achieving a zero-phase offset between the output and input sequences. When the signal delay is not greater than the threshold, the initial filtered sequence is directly used as the initial zero-phase sequence.

[0030] Boundary-weighted smoothing is performed on the initial zero-phase sequence. Since the filtering process introduces boundary abrupt changes at both ends of the sequence, weighted smoothing is required for the first 5% and last 5% of data points. A Hanning window function is used to generate weight coefficients. The Hanning window function is calculated as follows: for the nth point within the window, with a total window length of N, the weight coefficient equals 0.5 multiplied by the number within parentheses minus the cosine function, where the independent variable of the cosine function is twice pi multiplied by n divided by N minus one. The original amplitude of each point in the boundary region is multiplied by the corresponding Hanning window weight coefficient, and then weighted and averaged with the amplitude of adjacent internal points to gradually transition the boundary values ​​to a stable state. After boundary-weighted smoothing, a smoothed zero-phase sequence is obtained. This sequence preserves the phase characteristics of the original signal, eliminates boundary abrupt changes, and provides a high-quality signal foundation for subsequent pressure difference calculations.

[0031] In step S14, pressure difference constants are eliminated based on the smoothed zero-phase sequence to obtain a comprehensive pressure difference index, including: Extract the inlet pressure and the outlet pressure from the smoothed zero-phase sequence to obtain the inlet pressure sequence and the outlet pressure sequence; Based on the inlet pressure sequence and the outlet pressure sequence, calculate the pressure difference of all sensor combinations at each time point to obtain the initial pressure difference set; Calculate the Euclidean distance between the pressure difference vector and the mean vector at each moment in the initial pressure difference set to obtain the deviation Euclidean distance; By removing the pressure differences in the initial pressure difference set whose deviation Euclidean distance is greater than a preset Euclidean distance threshold, a consistency difference sequence is obtained. The arithmetic mean of the consistency difference sequence is calculated to obtain the comprehensive pressure difference index.

[0032] Specifically, inlet and outlet pressure data are extracted from the smoothed zero-phase sequence. The smoothed zero-phase sequence is a two-dimensional table containing multiple columns of signals, including inlet pressure, outlet pressure, flow rate, and vibration acceleration. All columns corresponding to the inlet pressure are extracted, including measurements from multiple pressure sensors at the inlet, to obtain the inlet pressure sequence; all columns corresponding to the outlet pressure are extracted to obtain the outlet pressure sequence. Each row of each sequence corresponds to a time point, and each column corresponds to a sensor.

[0033] Based on the inlet and outlet pressure sequences, calculate the pressure difference for all sensor combinations at each time point. For each time point, take all inlet pressure values ​​from the inlet pressure sequence and all outlet pressure values ​​from the outlet pressure sequence. Subtract each outlet pressure value from each inlet pressure value to obtain a set of pressure differences for that time point. For example, if there are three pressure sensors at the inlet and two at the outlet, six pressure differences can be calculated for each time point. Organize the pressure differences for all time points in chronological order, with each pressure difference forming a vector. The set of vectors from all time points forms the initial pressure difference set.

[0034] Calculate the Euclidean distance between the pressure difference vector and the mean vector at each moment in the initial pressure difference set. First, calculate the arithmetic mean of the pressure difference vectors at all moments in the initial pressure difference set to obtain the mean vector. Each component of this mean vector corresponds to the average pressure difference value of the same sensor combination at all moments. Then, for the pressure difference vector at each moment, calculate its Euclidean distance to the mean vector. Subtract each corresponding component of the vector, square the result, sum all the squares, and take the square root to obtain the deviation Euclidean distance at that moment. The deviation Euclidean distances at all moments form a sequence.

[0035] The preset Euclidean distance threshold is determined through the statistical distribution of historical normal operating data. Pressure difference vectors at 1000 time points during normal filter operation are collected, and the deviation Euclidean distance at each time point is calculated, resulting in 1000 deviation Euclidean distance values. These values ​​are sorted in ascending order, and the 95th percentile is used as the Euclidean distance threshold. When the deviation Euclidean distance at a given time point exceeds this threshold, the pressure difference data at that time is considered an outlier and is removed from the initial pressure difference set; otherwise, it is retained. After removing outliers, the remaining pressure difference data are arranged in chronological order to form a consistency difference sequence.

[0036] Calculate the arithmetic mean of the consistency difference sequence. Sum all pressure differences in the consistency difference sequence and divide by the total number of pressure differences to obtain the comprehensive pressure difference index. This index integrates reliable pressure difference information from multiple sensor combinations and multiple time points, representing the average pressure difference of the filter under the current operating conditions, and is used for threshold comparison and trend prediction in subsequent steps.

[0037] In step S15, based on the comprehensive pressure difference index, spectral analysis is performed on the high-frequency sampling sequence to obtain the periodic interference frequency, and noise feature matching is performed on the high-frequency sampling sequence to obtain the noise type. Combining the variation characteristics of the comprehensive pressure difference index, the periodic interference frequency, and the noise type, the blockage initiation point is determined, including: When the comprehensive pressure difference index is greater than the preset pressure difference threshold, the vibration acceleration in the high-frequency sampling sequence is analyzed by the fast Fourier transform algorithm to obtain the vibration spectrum sequence, and the frequency corresponding to the peak value in the vibration spectrum sequence is extracted to obtain the periodic interference frequency. Based on the periodic interference frequency, the inlet pressure and outlet pressure of the high-frequency sampling sequence are filtered using the Chebyshev filtering algorithm to obtain a filtered pressure sequence. The peak factor is obtained by calculating the ratio of the absolute value of the peak to the root mean square of the filtered pressure sequence, and the pressure skewness is obtained by calculating the ratio of the third central moment to the cube of the standard deviation of the filtered pressure sequence. Based on the peak factor and the pressure skewness, feature matching is performed using a preset noise feature library to obtain the noise type, and the rate of change of the comprehensive pressure difference index per unit time is calculated to obtain the pressure difference change rate. When the noise type is a preset fluid impact noise and the differential pressure change rate is greater than a preset trend threshold, the envelope energy is calculated using the Hilbert transform algorithm based on the filtered pressure sequence to obtain the envelope energy sequence. The time nodes corresponding to the envelope energy exceeding the preset envelope energy threshold in the envelope energy sequence are then selected to obtain the blockage initiation point. When the noise type is not the preset fluid impact noise or the pressure difference change rate is not greater than the trend threshold, the standard sampling sequence and the high-frequency sampling sequence are reacquired, and the comprehensive pressure difference index is continued to be detected.

[0038] Specifically, the system determines whether the overall pressure difference index exceeds a preset pressure difference threshold. This threshold is determined through statistical analysis of historical fault data. One hundred overall pressure difference index values ​​are collected when the filter experiences early blockage, and the minimum value among them is taken as the pressure difference threshold. If the current overall pressure difference index is less than or equal to this threshold, the subsequent steps are skipped, and monitoring continues. If it exceeds the threshold, a high-frequency analysis process is initiated.

[0039] Vibration acceleration data is extracted from the high-frequency sampling sequence. The high-frequency sampling sequence, acquired in step S11 at a sampling rate of 10 kHz, includes inlet pressure, outlet pressure, and vibration acceleration. Vibration acceleration signal segments one second before and after the current moment are taken, and a Fast Fourier Transform (FFT) is performed on these segments. The FFT converts the time-domain signal into a frequency-domain distribution, outputting the amplitude corresponding to each frequency. All frequency points are iterated to find the frequency with the largest amplitude; this frequency is the periodic interference frequency. For example, if there is a motor operating at 50 Hz near the equipment, a peak will appear at 50 Hz in the spectrum; this frequency is the periodic interference frequency.

[0040] Based on the periodic interference frequency, Chebyshev filtering is applied to the inlet and outlet pressures in the high-frequency sampling sequence. The design parameters of the Chebyshev filter are: stopband bandwidth set to ±5 Hz of the periodic interference frequency, stopband attenuation set to 40 dB, and passband ripple set to 0.1 dB. First, the numerator and denominator coefficients of the Chebyshev filter are calculated based on these parameters, obtained through the recursive relationship of the Chebyshev polynomial. Then, filtering is performed on the high-frequency inlet and outlet pressure sequences; that is, for each data point, the filter coefficients are weighted and summed with the 32 points before and after that point to obtain the new filtered value. After filtering, the periodic interference frequency component is attenuated, resulting in the filtered pressure sequence.

[0041] Calculate the peak factor and pressure skewness of the filtered pressure sequence. The peak factor is calculated as follows: take the absolute values ​​of all pressure values ​​in the filtered pressure sequence and find the maximum value as the peak; calculate the square of each pressure value in the sequence, sum them, divide by the sequence length, and then take the square root to obtain the root mean square value; divide the peak value by the root mean square value to obtain the peak factor. The pressure skewness is calculated as follows: calculate the arithmetic mean of all pressure values ​​in the filtered pressure sequence as the mean; calculate the cube of the difference between each pressure value and the mean, sum them, and divide by the sequence length to obtain the third central moment; calculate the square of the difference between each pressure value and the mean, sum them, and divide by the sequence length to obtain the variance; take the square root of the variance to obtain the standard deviation, and then take the cube root of the variance to obtain the cube of the standard deviation; divide the third central moment by the cube of the standard deviation to obtain the pressure skewness.

[0042] The noise feature library is pre-built to identify noise types based on peak factor and pressure skewness. The construction process involves collecting one thousand signal samples of known noise types, each containing a pressure signal segment and its corresponding noise type label. The peak factor and pressure skewness are calculated for each sample to obtain two-dimensional feature points. These feature points of all fluid impact noise samples are plotted on a plane, and the convex hull boundary of these feature points is calculated to obtain the discrimination region for fluid impact noise. The pre-built noise feature library contains the boundary equations for this discrimination region. During use, the currently calculated peak factor and pressure skewness are used as coordinate points to determine whether the point falls within the discrimination region for fluid impact noise. If it falls within the region, the noise type is determined to be fluid impact noise; otherwise, it is another type of noise.

[0043] It should be further explained that the noise feature database contains discrimination regions for various noise types. Each noise type corresponds to a two-dimensional region formed by the convex hull of sample feature points. During matching, the currently calculated peak factor and pressure skewness are used as a two-dimensional coordinate point, and it is sequentially determined whether the point is located within the convex hull region of each noise type. If it is located within a certain region, it is determined to be the corresponding noise type. If it is located within multiple regions, the type corresponding to the region with the smallest convex hull area is selected. If it is not located within any region, it is determined to be an unknown noise type. In this case, the blockage initiation point location is not triggered, and monitoring continues.

[0044] After determining the noise type, the rate of change of the comprehensive pressure difference index over a unit of time is further calculated to obtain the pressure difference change rate. The pressure difference change rate is the ratio of the difference in the comprehensive pressure difference index between adjacent unit time intervals to the unit time itself. Specifically, the comprehensive pressure difference index at the end of the current unit time interval is subtracted from the comprehensive pressure difference index at the beginning of the current unit time interval to obtain the pressure difference increment. This increment is then divided by the unit time length to obtain the pressure difference change rate. The pressure difference change rate is used to characterize the upward trend of the comprehensive pressure difference index over time.

[0045] The preset trend threshold is determined through statistical analysis of historical operating data. Specifically, the comprehensive pressure difference index is collected under normal operating conditions, pipeline vortex conditions, media containing air bubbles conditions, valve disturbance conditions, and fluid pressure pulsation conditions. The rate of change of pressure difference per unit time under each condition is calculated to obtain a sample set of non-clogging pressure difference change rates. The sample set of non-clogging pressure difference change rates is arranged in ascending order of values, and the 95th percentile is taken as the trend threshold. Therefore, when the current rate of change of pressure difference is greater than the trend threshold, it indicates that the increase in the comprehensive pressure difference index has exceeded the normal fluctuation range under non-clogging disturbance conditions.

[0046] The causes of fluid impact noise are not limited to filter element blockage. Pipeline eddies, media bubbles, valve disturbances, and fluid pressure pulsations can also produce impact characteristics. If the noise type is determined solely based on the peak factor and pressure skewness, the aforementioned non-clogging disturbances may be misjudged as filter element blockage. Therefore, this embodiment further introduces the differential pressure change rate in addition to the noise type determination. When the noise type is the preset fluid impact noise and the differential pressure change rate is greater than the trend threshold, it indicates that the pressure difference between the filter inlet and outlet is continuously increasing per unit time, which is consistent with the differential pressure change pattern when filter element blockage gradually worsens. When the differential pressure change rate is not greater than the trend threshold, it indicates that the pressure difference does not show a continuous upward trend, and the current fluid impact noise is more likely caused by pipeline eddies, media bubbles, valve disturbances, or fluid pressure pulsations.

[0047] When the noise type is fluid impact noise and the rate of change of pressure difference is greater than the trend threshold, a Hilbert transform is performed on the filtered pressure sequence to calculate the envelope energy sequence. The specific process of the Hilbert transform is as follows: a Fast Fourier Transform (FFT) is performed on the filtered pressure sequence, the amplitude of the positive frequency portion is multiplied by two, the negative frequency portion is set to zero, and then an Inverse Fast Fourier Transform (IFFT) is performed to obtain the analytic signal. The real part of the analytic signal is the original signal, and the imaginary part is the result of the Hilbert transform. The sum of the squares of the real and imaginary parts of the analytic signal at each time point is calculated, and then the square root is taken to obtain the envelope amplitude sequence. The amplitude at each point in the envelope amplitude sequence is squared to obtain the envelope energy sequence.

[0048] The preset envelope energy threshold is determined by statistically analyzing the envelope energy baseline under normal operating conditions. Envelope energy data for one minute continuously under normal operating conditions is collected from the filter. The arithmetic mean of these data is calculated as the baseline value, and this baseline value is multiplied by three to obtain the envelope energy threshold. Starting from the beginning of the envelope energy sequence, the envelope energy value at each time point is compared with this threshold in chronological order. The time point at which the first envelope energy value continuously exceeds the threshold is identified. This time point is the clogging initiation point. If no instance of exceeding the threshold occurs throughout the entire sequence, the clogging initiation point is not output, and monitoring continues.

[0049] When the noise type is not the preset fluid impact noise or the differential pressure change rate is not greater than the trend threshold, the envelope energy calculation and blockage start point location are not performed. The process returns to step S11 to continue acquiring real-time data, and the comprehensive differential pressure index is updated in step S14. That is, the high-frequency sampling sequence segment currently used for judgment is discarded, and the process returns to step S11 to continue acquiring new standard sampling sequences and high-frequency sampling sequences, and the process starts again from step S12.

[0050] like Figure 2 As shown in the figure, this graph displays the characteristic curve of the filter inlet pressure difference changing over time, showing the complete process from the initial stable stage to the later nonlinear accelerated increase in pressure difference. The graph also marks the "blockage surge initiation point" identified based on this invention. Figure 2 This intuitively demonstrates the physical characteristics of the pressure difference rising slowly in the early stage of filter blockage and accelerating in the later stage, thus effectively supporting the necessity and rationality of the second-order dynamic extrapolation strategy in step S18 of this application, and verifying that the nonlinear prediction model can accurately capture the real pressure difference change trend.

[0051] It is worth noting that, such as Figure 2As shown, under simulated conditions, the inlet pressure difference exhibits a slow, fluctuating upward trend in the initial stage (0 to 5.2s), followed by an accelerated upward phase. The algorithm of this embodiment successfully identified the "surge start point" as t=5.2s. If a traditional fixed threshold method is used, for example, using a pressure difference reaching 1.00 MPa as the alarm condition, the alarm time would typically lag to around t=8.5s under this simulated condition. In contrast, this embodiment advances the warning time by approximately 3.3 seconds, an advance of about 38%. Furthermore, throughout the entire monitoring period, this embodiment did not misjudge non-blocking disturbances such as water hammer in the pipeline as blockages, with a false alarm rate of 0%. The above experimental data fully demonstrates that this invention, by introducing a nonlinear extrapolation strategy and a dual constraint mechanism, can effectively capture the weak characteristics of the initial stage of blockage, significantly improving the timeliness and reliability of the warning.

[0052] In step S16, based on the blockage initiation point and the smoothed zero-phase sequence, abnormal feature removal is performed to obtain the final vector set, including: Centered on the blockage initiation point, the smooth zero-phase sequence is segmented through a preset time truncation window to obtain the target signal segment; The mean, variance, standard deviation, peak value, and kurtosis of the target signal segment are calculated to obtain the initial feature vector; Based on the initial feature vector, abnormal feature vectors are removed using the Laida criterion to obtain a cleaned feature vector; The cleaned feature vectors are normalized to obtain the final vector set.

[0053] Specifically, the target signal segment is extracted from the smoothed zero-phase sequence, centered on the blockage initiation point. The smoothed zero-phase sequence is a two-dimensional table output in step S13, containing signals such as inlet pressure, outlet pressure, flow rate, and vibration acceleration. The preset time truncation window length is determined through historical data statistics. Specifically, one hundred blockage events are collected. For each blockage event, the duration range in which the change in signal characteristic values ​​before and after the blockage initiation point exceeds 20% of the baseline value is calculated. The 90th percentile of these duration ranges is taken as the window length, for example, two seconds before and after. Centered on the timestamp corresponding to the blockage initiation point, half the window length is extracted forward and backward to obtain a continuous time segment. This segment contains data from all signal channels within this time interval and is called the target signal segment.

[0054] Calculate the mean, variance, standard deviation, peak value, and kurtosis for each signal channel in the target signal segment. For a data sequence of a given signal channel, the mean is calculated by summing all data points and then dividing by the number of data points. The variance is calculated by first calculating the difference between each data point and the mean, squared the differences, summing them, and then dividing by the number of data points minus one. The standard deviation is the square root of the variance. The peak value is calculated by iterating through all data points, finding the data point with the largest absolute value, and taking that absolute value as the peak value. The kurtosis is calculated by first calculating the fourth power of the difference between each data point and the mean, summing them, and then dividing by the number of data points to obtain the fourth central moment; then calculating the square of the variance, dividing the fourth central moment by the square of the variance, and obtaining the kurtosis. Calculate the above five statistics for each signal channel in the target signal segment, and concatenate the calculation results of all channels in the order of inlet pressure, outlet pressure, flow rate, and vibration acceleration into a one-dimensional vector, which is called the initial feature vector.

[0055] The initial feature vectors are subjected to anomaly removal according to the Raida criterion. The Raida criterion is based on the normal distribution assumption, which states that under normal conditions, each feature value in the feature vector should fall within a range of plus or minus three standard deviations of the mean. First, one thousand initial feature vectors are collected under normal filter operation, each containing multiple feature values. For each feature dimension, the arithmetic mean and standard deviation of these one thousand feature values ​​are calculated. The mean minus three standard deviations is used as the lower threshold, and the mean plus three standard deviations is used as the upper threshold. For the current initial feature vector to be detected, its feature values ​​are checked dimension by dimension to see if they fall within the threshold range of the corresponding dimension. If any feature value exceeds the range, the entire initial feature vector is determined to be anomaly and is removed; otherwise, it is retained. The feature vectors retained after removal are called cleaned feature vectors. If the initial feature vector is removed, no further processing is performed, and the process returns to step S14 for re-monitoring.

[0056] It is worth noting that for multi-dimensional feature vectors, this embodiment adopts a method of independent testing of each dimension. First, the mean and standard deviation of each feature dimension are calculated from the feature vectors under historical normal operating conditions. For the feature vector to be detected, if the value of any feature dimension exceeds the range of plus or minus three times the standard deviation of the corresponding mean, the feature vector is determined to be abnormal and is removed; if all feature dimensions are within the range, it is retained. This method is simple to calculate and suitable for real-time online processing. Moreover, since the physical meaning of each feature dimension is independent, such as the mean reflecting the amplitude level and kurtosis reflecting the impact characteristics, independent judgment of each dimension can effectively identify abnormal deviations in a single dimension.

[0057] The cleaned feature vectors are normalized to obtain the final vector set. Normalization uses a minimum-maximum method, scaling the values ​​of each feature dimension to between zero and one. First, from the cleaned feature vectors collected during historical normal operation, the minimum and maximum values ​​for each feature dimension are calculated. Then, for each feature value in the current cleaned feature vector, the minimum value of that dimension is subtracted, and then divided by the difference between the maximum and minimum values ​​of that dimension to obtain the normalized feature value. All normalized feature values ​​constitute the final vector set. This set eliminates the differences in units and numerical ranges between different features, providing input data with a uniform scale for subsequent curve fitting and risk prediction.

[0058] In step S17, curve fitting is performed on the final vector set to obtain the predicted pressure difference curve. The predicted pressure difference curve is then corrected for deviation using the comprehensive pressure difference index to obtain the corrected pressure difference curve. Finally, a risk probability curve is obtained by combining the corrected pressure difference curve with a pre-constructed risk mapping function, including: Based on the final vector set, curve fitting is performed using the least squares method to obtain the predicted pressure difference curve; The difference between the predicted pressure difference curve at the current moment and the comprehensive pressure difference index is calculated to obtain the prediction deviation; The predicted deviation is weighted and calculated according to the preset deviation compensation coefficient to obtain the deviation compensation amount; The deviation compensation amount is superimposed on the predicted pressure difference curve to obtain the corrected pressure difference curve.

[0059] The process of constructing the risk mapping function includes: Obtain historical pressure difference sequences and historical risk sequences; Calculate the mean, peak value, and rate of change in the historical pressure difference sequence to obtain the pressure difference feature set; Based on the pressure difference feature set and the historical risk sequence, a risk mapping function is obtained by performing polynomial regression fitting using the least squares method.

[0060] Specifically, the predicted pressure difference curve is obtained by curve fitting based on the final vector set. The final vector set is a two-dimensional table, with each row corresponding to a time window and each column corresponding to a normalized feature value, which contains the statistical characteristics of the pressure difference. The mean pressure difference feature in the final vector set is selected as the fitting object and arranged into a sequence in chronological order. The least squares method is used to perform polynomial curve fitting on this sequence, with the polynomial order set to third. The specific process of the least squares method is as follows: assuming the fitted curve is a cubic polynomial with unknown coefficients, the sum of the squares of the differences between the mean actual pressure difference at each time point and the predicted value at that time point on the polynomial is obtained. The coefficients that minimize the sum of squares of errors are obtained by solving the system of linear equations, i.e., by taking the partial derivative of each coefficient and setting it to zero. The time variable is then substituted into the polynomial to calculate the predicted pressure difference value at each time point in the future. These predicted values ​​are then connected in chronological order to form a curve, thus obtaining the predicted pressure difference curve.

[0061] Calculate the prediction deviation. Take the predicted value of the pressure difference curve at the current moment, i.e., the ordinate value of the curve corresponding to the current time point. Obtain the comprehensive pressure difference index from step S14, which is the actual measured pressure difference value at the current moment. Subtract the predicted pressure difference of the predicted pressure difference curve at the current moment from the comprehensive pressure difference index to obtain the prediction deviation. If the predicted value is greater than the measured value, the prediction deviation is negative; if the predicted value is less than the measured value, the prediction deviation is positive.

[0062] The deviation compensation amount is obtained by weighting the predicted deviation according to a preset deviation compensation coefficient. The deviation compensation coefficient is determined through historical data statistics. Specifically, one hundred sets of data on predicted deviation and actual error are collected under normal operating conditions of the filter. The linear regression coefficient between the predicted deviation and actual error is calculated, and this coefficient is used as the deviation compensation coefficient, with a value between 0.5 and 1.5. The predicted deviation is multiplied by the deviation compensation coefficient to obtain the deviation compensation amount. For example, if the deviation compensation coefficient is 1.2 and the predicted deviation is -0.1 MPa, then the deviation compensation amount is -0.12 MPa.

[0063] The deviation compensation is superimposed onto the predicted pressure difference curve to obtain the corrected pressure difference curve. The superposition operation involves adding the deviation compensation to the predicted pressure difference value at each future time point on the predicted pressure difference curve to obtain the corrected pressure difference value at that time point. All corrected pressure difference values ​​are connected in chronological order to form the corrected pressure difference curve. Because the deviation compensation is calculated based on the difference between the measured and predicted values ​​at the current moment, and the same compensation value is applied to the entire future curve, the corrected curve more closely reflects actual operating conditions.

[0064] Based on the corrected pressure difference curve, a risk probability curve is obtained by mapping the risk probability using a pre-constructed risk mapping function. The risk mapping function converts the pressure difference value into a performance degradation ratio. The construction process of this function is as follows: Historical pressure difference sequences and historical risk sequences are obtained. The historical pressure difference sequence consists of one thousand sets of pressure difference data collected from the filter's historical operating data. Each set of data includes the mean, peak, and rate of change of the pressure difference within a time window. The historical risk sequence represents the actual performance degradation ratio of the filter at the corresponding time point, obtained through shutdown testing or laboratory testing. The mean, peak, and rate of change of each set of data in the historical pressure difference sequence are calculated. The mean is calculated as the arithmetic mean of all pressure difference values; the peak is the pressure difference value with the largest absolute value in the sequence; the rate of change is the ratio of the difference between pressure difference values ​​at adjacent times to the time interval, and the average of the absolute values ​​is taken. These three features constitute the pressure difference feature set, with each sample being a three-dimensional vector. Using the pressure difference feature set as the independent variable and the historical risk sequence as the dependent variable, a third-order polynomial regression is performed using the least squares method to obtain the risk mapping function. This function accepts three inputs: the average pressure difference, the peak pressure difference, and the rate of change, and outputs the corresponding risk probability value. When using this function in step S17, the average pressure difference, peak pressure difference, and rate of change for each time point are extracted from the corrected pressure difference curve, input into the risk mapping function, and the risk probability value for that time point is obtained. The risk probability values ​​for all time points are then concatenated sequentially to obtain the risk probability curve. This curve reflects the trend of filter performance degradation over time.

[0065] In step S18, based on the risk probability curve, failure time is predicted to obtain the critical failure time, and based on the critical failure time, a maintenance alarm signal is issued and the replacement time is displayed, including: The warning start point is obtained by filtering out the time node corresponding to the first time that the risk probability in the risk probability curve exceeds the preset risk probability threshold in chronological order. Calculate the time elapsed from the warning start point to the risk probability reaching the preset failure threshold to obtain the prediction cycle length; The risk probability curve is extended to the future time point corresponding to the prediction period length to obtain the complete risk curve; Extract the time nodes corresponding to the failure thresholds from the complete risk curve to obtain the critical failure time.

[0066] Specifically, the warning start point is obtained by filtering the time node corresponding to the first time the risk probability in the risk probability curve exceeds the preset risk probability threshold in chronological order. The risk probability curve is a sequence of performance degradation ratios over time output in step S17, with time on the horizontal axis and risk probability value on the vertical axis. The preset risk probability threshold is determined by statistical analysis of historical fault data. One hundred risk probability values ​​are collected when the filter experiences irreversible blockage, and the minimum value among these values ​​is taken as the risk probability threshold, for example, 60%. Starting from the starting time point of the risk probability curve, the risk probability value at each time point is compared with the threshold. When the first time point with a risk probability value greater than the threshold is encountered, that time point is recorded as the warning start point. If the entire curve does not exceed the threshold, the subsequent warning process is not triggered, and the process returns to step S14 to continue monitoring.

[0067] The prediction period length is calculated by determining the time elapsed from the warning start point to the risk probability reaching a preset failure threshold. The failure threshold, set to 100%, represents the risk probability value corresponding to complete filter failure. Starting from the warning start point, the risk probability curve is traversed along the time direction to find the point where the risk probability value first reaches or exceeds 100%, which is recorded as the failure point. The timestamp of the failure point is subtracted from the timestamp of the warning start point to obtain the time difference, which is the prediction period length. For example, if the warning start point is at the 10th minute and the failure point is at the 30th minute, the prediction period length is 20 minutes. If the risk probability curve does not reach the failure threshold within the predictable range, an extrapolation estimate is made based on the trend at the end of the risk probability curve.

[0068] The risk probability curve is extended to future time points corresponding to the prediction period length to obtain a complete risk curve. The extension operation utilizes the risk mapping function corresponding to the corrected pressure difference curve fitted in step S17. Specifically, starting from the end of the current risk probability curve, multiple consecutive time points at the end of the curve are selected as the calculation objects. The difference between the risk probability values ​​at adjacent time points is calculated to obtain the first-order change, and the difference between adjacent first-order changes is calculated to obtain the second-order change. When the second-order change is greater than zero, it indicates that the rate of increase in risk probability is increasing. Subsequent risk probability values ​​are then calculated point-by-point according to the current risk probability value, the first-order change, and the second-order change. When the second-order change is not greater than zero, it indicates that the risk probability has not shown an accelerating upward trend. Subsequent risk probability values ​​are then calculated according to the first-order change. The calculated risk probability values ​​are appended to the original curve in chronological order to form a complete risk probability trajectory from the warning start point to the failure threshold, called the complete risk curve.

[0069] It should be noted that the extension of the risk probability curve is not always a linear extrapolation. In the early stages of filter clogging, the pressure difference typically increases slowly; as the clogging worsens, the pressure difference may increase non-linearly and at an accelerated pace. Since the risk probability curve is mapped from the modified pressure difference curve, it may also exhibit an accelerated upward trend in the later stages. Therefore, in this embodiment, when extending the risk probability curve, the first-order change represents the current rate of increase in risk probability, and the second-order change represents the change in the rate of increase. When the second-order change is greater than zero, an acceleration term is added to the subsequent calculation to reflect the non-linear increase in risk probability in the later stages of clogging; when the second-order change is not greater than zero, extrapolation is performed according to the first-order change. If the risk probability value exceeds 100% during the calculation process, it is truncated to 100%.

[0070] The critical failure time is obtained by extracting the time points corresponding to the failure thresholds from the complete risk curve. Specifically, starting from the beginning of the complete risk curve, the risk probability value is checked point by point in chronological order to determine whether it reaches the failure threshold. When the risk probability value first reaches or exceeds the failure threshold, that time point is taken as the critical failure time. If there are multiple time points in the complete risk curve where the risk probability value reaches the failure threshold, the earliest time point in chronological order is selected as the critical failure time.

[0071] Based on the critical failure time, a maintenance alarm signal is issued and the replacement time is displayed. Specifically, the current time is compared with the critical failure time to calculate the remaining operational time. When the remaining operational time is less than the preset maintenance lead time, the system issues a maintenance alarm signal through an audible and visual alarm, a display terminal, or a remote communication module, and displays the recommended replacement time, remaining operational time, and current risk probability value on the display interface. Through this process, the system can provide maintenance prompts before the filter reaches a critical failure state, avoiding decreased filtration performance or abnormal equipment shutdown due to increased clogging.

[0072] Reference Figure 3 The second embodiment of the present invention provides a filter performance testing system based on intelligent sensors, comprising: The data acquisition module is used to acquire standard sampling sequences and high-frequency sampling sequences containing inlet pressure, outlet pressure, flow rate, and vibration acceleration. The frequency adjustment module is used to adjust the sampling frequency according to the standard sampling sequence to obtain a structured signal set; The phase compensation module is used to perform cross-correlation calculation based on the structured signal set to obtain the signal delay, and to perform phase compensation based on the signal delay to obtain a smooth zero-phase sequence. The index calculation module is used to remove pressure difference anomalies based on the smoothed zero-phase sequence to obtain a comprehensive pressure difference index. The blockage location module is used to perform spectral analysis on the high-frequency sampling sequence to obtain the periodic interference frequency based on the comprehensive pressure difference index, and to perform noise feature matching on the high-frequency sampling sequence to obtain the noise type. The blockage initiation point is determined by combining the variation characteristics of the comprehensive pressure difference index, the periodic interference frequency, and the noise type. The anomaly removal module is used to remove anomaly features based on the blockage initiation point and the smoothed zero-phase sequence to obtain the final vector set; The risk mapping module is used to perform curve fitting on the final vector set to obtain the predicted pressure difference curve, use the comprehensive pressure difference index to correct the deviation of the predicted pressure difference curve to obtain the corrected pressure difference curve, and combine the pre-constructed risk mapping function to perform risk probability mapping on the corrected pressure difference curve to obtain the risk probability curve. The output module is used to predict the failure time based on the risk probability curve, obtain the critical failure time, and issue a maintenance alarm signal and display the replacement time based on the critical failure time.

[0073] It should be noted that the filter performance testing system based on intelligent sensors provided in this embodiment of the invention is used to execute all the process steps of the filter performance testing method based on intelligent sensors in the above embodiment. The working principles and beneficial effects of the two are one-to-one, so they will not be described again.

[0074] It should be noted that the system embodiments described above are merely illustrative. The units described as separate components may or may not be physically separate, and the components shown as units may or may not be physical units; that is, they may be located in one place or distributed across multiple network units. Some or all of the modules can be selected to achieve the purpose of this embodiment according to actual needs. Furthermore, in the accompanying drawings of the system embodiments provided by this invention, the connection relationships between modules indicate that they have communication connections, which can be specifically implemented as one or more communication buses or signal lines. Those skilled in the art can understand and implement this without any creative effort.

[0075] The specific embodiments described above further illustrate the purpose, technical solution, and beneficial effects of the present invention. It should be understood that the above descriptions are merely specific embodiments of the present invention and are not intended to limit the scope of protection of the present invention. In particular, it should be noted that any modifications, equivalent substitutions, improvements, etc., made within the spirit and principles of the present invention should be included within the scope of protection of the present invention for those skilled in the art.

Claims

1. A filter performance testing method based on intelligent sensors, characterized in that, include: Acquire standard and high-frequency sampling sequences containing inlet pressure, outlet pressure, flow rate, and vibration acceleration; Based on the standard sampling sequence, the sampling frequency is adjusted to obtain a structured signal set; Based on the structured signal set, cross-correlation calculation is performed to obtain the signal delay. Based on the signal delay, phase compensation is performed to obtain a smooth zero-phase sequence. Based on the smoothed zero-phase sequence, pressure difference constants are eliminated to obtain a comprehensive pressure difference index; Based on the comprehensive pressure difference index, the periodic interference frequency is obtained by performing spectral analysis on the high-frequency sampling sequence, and the noise type is obtained by noise feature matching on the high-frequency sampling sequence. The blockage initiation point is determined by combining the variation characteristics of the comprehensive pressure difference index, the periodic interference frequency, and the noise type. Based on the blockage initiation point and the smoothed zero-phase sequence, abnormal features are removed to obtain the final vector set; Curve fitting is performed on the final vector set to obtain the predicted pressure difference curve. The deviation of the predicted pressure difference curve is corrected using the comprehensive pressure difference index to obtain the corrected pressure difference curve. The risk probability curve is obtained by combining the corrected pressure difference curve with the pre-constructed risk mapping function. Based on the risk probability curve, failure time is predicted to obtain the critical failure time. Based on the critical failure time, a maintenance alarm signal is issued and the replacement time is displayed.

2. The filter performance testing method based on intelligent sensors according to claim 1, characterized in that, The step of adjusting the sampling frequency according to the standard sampling sequence to obtain a structured signal set includes: Based on the standard sampling sequence, timestamp alignment is performed to obtain a synchronous sampling sequence, and based on the synchronous sampling sequence, spectrum analysis is performed using the Fast Fourier Transform algorithm to obtain the initial main frequency; When the initial main frequency exceeds a preset main frequency threshold, the preset sampling frequency is increased to twice the initial main frequency to obtain a high-density sampling sequence; Based on the high-density sampling sequence and the synchronous sampling sequence, interpolation calculation is performed using the Lagrange interpolation algorithm to obtain the resampling sequence; When the initial clock frequency does not exceed the clock frequency threshold, the synchronous sampling sequence is directly used as the resampling sequence; Based on the resampled sequence, the data is spliced ​​together according to the time dimension to obtain a structured signal set.

3. The filter performance testing method based on intelligent sensors according to claim 1, characterized in that, The step of performing cross-correlation calculations based on the structured signal set to obtain the signal delay, and then performing phase compensation based on the signal delay to obtain a smoothed zero-phase sequence includes: Based on the structured signal set, the cutoff frequency and attenuation amplitude are extracted by fast Fourier transform to obtain filtering parameters. Based on the filtering parameters, the structured signal set is convolved by the Butterworth filtering algorithm to obtain a preliminary filtering sequence. Based on the structured signal set and the preliminary filtering sequence, cross-correlation function is calculated to obtain a cross-correlation function sequence and a time shift sequence, and the time shift in the time shift sequence corresponding to the peak value in the cross-correlation function sequence is used as the signal delay. When the signal delay is greater than the preset signal delay threshold, phase compensation is performed using a zero-phase filtering algorithm based on the preliminary filtering sequence to obtain an initial zero-phase sequence. When the signal delay is not greater than the signal delay threshold, the preliminary filtering sequence is directly used as the initial zero-phase sequence; Based on the initial zero-phase sequence, a smoothed zero-phase sequence is obtained by performing boundary weighted smoothing using a preset Hanning window function.

4. The filter performance testing method based on intelligent sensors according to claim 1, characterized in that, The step of removing pressure difference constants based on the smoothed zero-phase sequence to obtain a comprehensive pressure difference index includes: Extract the inlet pressure and the outlet pressure from the smoothed zero-phase sequence to obtain the inlet pressure sequence and the outlet pressure sequence; Based on the inlet pressure sequence and the outlet pressure sequence, calculate the pressure difference of all sensor combinations at each time point to obtain the initial pressure difference set; Calculate the Euclidean distance between the pressure difference vector and the mean vector at each moment in the initial pressure difference set to obtain the deviation Euclidean distance; By removing the pressure differences in the initial pressure difference set whose deviation Euclidean distance is greater than a preset Euclidean distance threshold, a consistency difference sequence is obtained. The arithmetic mean of the consistency difference sequence is calculated to obtain the comprehensive pressure difference index.

5. The filter performance testing method based on intelligent sensors according to claim 1, characterized in that, The step of performing spectral analysis on the high-frequency sampling sequence based on the comprehensive pressure difference index to obtain the periodic interference frequency, and performing noise feature matching on the high-frequency sampling sequence to obtain the noise type, and combining the variation characteristics of the comprehensive pressure difference index, the periodic interference frequency, and the noise type to determine the blockage initiation point, includes: When the comprehensive pressure difference index is greater than the preset pressure difference threshold, the vibration acceleration in the high-frequency sampling sequence is analyzed by the fast Fourier transform algorithm to obtain the vibration spectrum sequence, and the frequency corresponding to the peak value in the vibration spectrum sequence is extracted to obtain the periodic interference frequency. Based on the periodic interference frequency, the inlet pressure and outlet pressure of the high-frequency sampling sequence are filtered using the Chebyshev filtering algorithm to obtain a filtered pressure sequence. The peak factor is obtained by calculating the ratio of the absolute value of the peak to the root mean square of the filtered pressure sequence, and the pressure skewness is obtained by calculating the ratio of the third central moment to the cube of the standard deviation of the filtered pressure sequence. Based on the peak factor and the pressure skewness, feature matching is performed using a preset noise feature library to obtain the noise type, and the rate of change of the comprehensive pressure difference index per unit time is calculated to obtain the pressure difference change rate. When the noise type is a preset fluid impact noise and the differential pressure change rate is greater than a preset trend threshold, the envelope energy is calculated using the Hilbert transform algorithm based on the filtered pressure sequence to obtain the envelope energy sequence. The time nodes corresponding to the envelope energy exceeding the preset envelope energy threshold in the envelope energy sequence are then selected to obtain the blockage initiation point. When the noise type is not the preset fluid impact noise or the pressure difference change rate is not greater than the trend threshold, the standard sampling sequence and the high-frequency sampling sequence are reacquired, and the comprehensive pressure difference index is continued to be detected.

6. The filter performance testing method based on intelligent sensors according to claim 1, characterized in that, The process involves removing abnormal features based on the blockage initiation point and the smoothed zero-phase sequence to obtain a final vector set, including: Centered on the blockage initiation point, the smooth zero-phase sequence is segmented through a preset time truncation window to obtain the target signal segment; The mean, variance, standard deviation, peak value, and kurtosis of the target signal segment are calculated to obtain the initial feature vector; Based on the initial feature vector, abnormal feature vectors are removed using the Laida criterion to obtain a cleaned feature vector; The cleaned feature vectors are normalized to obtain the final vector set.

7. The filter performance testing method based on intelligent sensors according to claim 1, characterized in that, The process of obtaining a predicted pressure difference curve by curve fitting the final vector set, and then correcting the predicted pressure difference curve using the comprehensive pressure difference index to obtain a corrected pressure difference curve, includes: Based on the final vector set, curve fitting is performed using the least squares method to obtain the predicted pressure difference curve; The difference between the predicted pressure difference curve at the current moment and the comprehensive pressure difference index is calculated to obtain the prediction deviation; The predicted deviation is weighted and calculated according to the preset deviation compensation coefficient to obtain the deviation compensation amount; The deviation compensation amount is superimposed on the predicted pressure difference curve to obtain the corrected pressure difference curve.

8. The filter performance testing method based on intelligent sensors according to claim 1, characterized in that, The process of constructing the risk mapping function includes: Obtain historical pressure difference sequences and historical risk sequences; Calculate the mean, peak value, and rate of change in the historical pressure difference sequence to obtain the pressure difference feature set; Based on the pressure difference feature set and the historical risk sequence, a risk mapping function is obtained by performing polynomial regression fitting using the least squares method.

9. The filter performance testing method based on intelligent sensors according to claim 1, characterized in that, The step of predicting the failure time based on the risk probability curve to obtain the critical failure time includes: The warning start point is obtained by filtering out the time node corresponding to the first time that the risk probability in the risk probability curve exceeds the preset risk probability threshold in chronological order. Calculate the time elapsed from the warning start point to the risk probability reaching the preset failure threshold to obtain the prediction cycle length; The risk probability curve is extended to the future time point corresponding to the prediction period length to obtain the complete risk curve; Extract the time nodes corresponding to the failure thresholds from the complete risk curve to obtain the critical failure time.

10. A filter performance testing system based on intelligent sensors, characterized in that, include: The data acquisition module is used to acquire standard sampling sequences and high-frequency sampling sequences containing inlet pressure, outlet pressure, flow rate, and vibration acceleration. The frequency adjustment module is used to adjust the sampling frequency according to the standard sampling sequence to obtain a structured signal set; The phase compensation module is used to perform cross-correlation calculation based on the structured signal set to obtain the signal delay, and to perform phase compensation based on the signal delay to obtain a smooth zero-phase sequence. The index calculation module is used to remove pressure difference anomalies based on the smoothed zero-phase sequence to obtain a comprehensive pressure difference index. The blockage location module is used to perform spectral analysis on the high-frequency sampling sequence to obtain the periodic interference frequency based on the comprehensive pressure difference index, and to perform noise feature matching on the high-frequency sampling sequence to obtain the noise type. The blockage initiation point is determined by combining the variation characteristics of the comprehensive pressure difference index, the periodic interference frequency, and the noise type. The anomaly removal module is used to remove anomaly features based on the blockage initiation point and the smoothed zero-phase sequence to obtain the final vector set; The risk mapping module is used to perform curve fitting on the final vector set to obtain the predicted pressure difference curve, use the comprehensive pressure difference index to correct the deviation of the predicted pressure difference curve to obtain the corrected pressure difference curve, and combine the pre-constructed risk mapping function to perform risk probability mapping on the corrected pressure difference curve to obtain the risk probability curve. The output module is used to predict the failure time based on the risk probability curve, obtain the critical failure time, and issue a maintenance alarm signal and display the replacement time based on the critical failure time.

Images