Karst aquifer monitoring method and system based on internet of things

Abstract

The application relates to the field of hydrogeological monitoring, and discloses a karst aquifer monitoring method and system based on the Internet of Things. The method comprises the following steps: collecting medium fluctuation signals and extracting propagation difference characteristics, and dividing high and low permeability regions; formulating a differentiated sampling strategy according to the region type; dynamically configuring sensor node parameters and collecting signals, determining a fluctuation source azimuth angle estimation value based on time series analysis; fusing multi-sensor node data to correct the azimuth angle estimation deviation, and obtaining an accurate azimuth angle; evaluating signal loss risk and generating an energy consumption optimization transmission scheme; extracting resource consumption indicators and feeding back to the sampling strategy, so as to realize dynamic balance of monitoring accuracy and energy consumption. The application combines adaptive sampling, multi-sensor node fusion positioning and feedback energy consumption regulation, solves the technical problems of poor adaptive capacity, low positioning accuracy and high energy consumption in karst aquifer monitoring, and improves monitoring accuracy and system endurance.

Description

TECHNICAL FIELD

[0001] The present application relates to the field of hydrogeological monitoring, and particularly relates to a karst aquifer monitoring method and system based on Internet of Things. BACKGROUND

[0002] As a key groundwater resource reservoir, the karst aquifer has a highly heterogeneous geological structure, with complex characteristics of interweaving fissures, conduits and low-permeability matrix, and is easily disturbed by mining, engineering and other activities, causing rapid changes in water flow and pressure field. Real-time and accurate monitoring of the dynamics of such aquifers and locating the disturbance source are of great significance to the safety of water resources and the prevention of geological disasters.

[0003] The existing technology mainly relies on deploying sensor networks with fixed parameters for monitoring, which has significant defects in dealing with the strong heterogeneity of karst aquifers. First, the system generally uses a uniform fixed sampling frequency, which cannot adapt to the dynamic characteristics of different medium regions. In high-permeability conduit or fissure regions, the fixed low-frequency sampling cannot capture the key transient signals, resulting in distorted monitoring data; while in low-permeability matrix regions, fixed high-frequency sampling produces a large amount of redundant data, unnecessarily occupying communication bandwidth and storage resources. Second, in the positioning of fluctuation sources (such as water inrush points), the existing methods are mostly based on simplified homogeneous medium models for triangulation, ignoring the huge differences in signal propagation speed and attenuation degree in different permeability media, resulting in large deviations in positioning results and making it difficult to meet the demand for accurate tracing. Finally, the energy management strategy of the system is extensive, and the sensor nodes usually collect and communicate continuously at a fixed power, lacking a dynamic scheduling mechanism according to geological features and the urgency of monitoring tasks, causing a great waste of energy and seriously restricting the long-term online service capability of Internet of Things nodes in the field without city power supply. These defects together result in the lack of precision, adaptability and energy efficiency of existing monitoring systems.

[0004] In view of the above defects, the present application combines adaptive sampling driven by geological features, accurate positioning through multi-sensor node data fusion, and dynamic energy regulation based on risk assessment, solving the technical problems of poor adaptability, low positioning accuracy and high energy consumption of existing monitoring systems in the strong heterogeneity environment of karst aquifers, improving the positioning accuracy of fluctuation sources, and significantly prolonging the overall endurance of the monitoring system. SUMMARY

[0005] The present application provides a karst aquifer monitoring method and system based on Internet of Things, solving the technical problems of poor adaptability, low positioning accuracy and high energy consumption of existing monitoring systems in the strong heterogeneity environment of karst aquifers, improving the positioning accuracy of fluctuation sources, and significantly prolonging the overall endurance of the monitoring system.

[0006] In a first aspect, the application provides a karst aquifer monitoring method based on Internet of Things, the method comprising:

[0007] In step S101, medium fluctuation signals are collected from a plurality of sensor nodes through an Internet of Things communication network, propagation difference characteristics of the medium fluctuation signals in different medium attribute regions are extracted, and a medium attribute region type distribution is determined according to a comparison result of the propagation difference characteristics and a preset region threshold value;

[0008] In step S102, a high-frequency sampling strategy is adopted for a high-permeability medium attribute region, a low-frequency sampling strategy is adopted for a low-permeability medium attribute region, and a sampling frequency adjustment scheme is generated according to the medium attribute region type distribution;

[0009] In step S103, sensor node parameters are dynamically configured according to the sampling frequency adjustment scheme, medium fluctuation signals are collected by using the configured sensor node parameters, and a wave source azimuth angle estimation value is determined based on time sequence analysis and signal time delay calculation;

[0010] In step S104, an accurate wave source azimuth angle is obtained by fusing data of adjacent sensor nodes in the communication network, correcting azimuth angle estimation deviation through an optimization algorithm, and according to the wave source azimuth angle estimation value;

[0011] In step S105, a signal loss risk is evaluated based on the accurate wave source azimuth angle and the medium attribute region type distribution, and an energy consumption optimization transmission scheme is generated according to a loss risk evaluation result;

[0012] In step S106, a resource consumption index is extracted from the energy consumption optimization transmission scheme, the resource consumption index is fed back to the sampling frequency adjustment scheme through a feedback loop mechanism, and a dynamic balance between monitoring system precision and energy consumption is achieved.

[0013] In a second aspect, the application provides a karst aquifer monitoring system based on Internet of Things, which is used to implement the karst aquifer monitoring method based on Internet of Things, and the system comprises:

[0014] A signal collection module is configured to collect medium fluctuation signals from a plurality of sensor nodes through an Internet of Things communication network, extract propagation difference characteristics of the medium fluctuation signals in different medium attribute regions, and determine a medium attribute region type distribution according to a comparison result of the propagation difference characteristics and a preset region threshold value;

[0015] A sampling strategy module is configured to generate a sampling frequency adjustment scheme according to the medium attribute region type distribution, adopt a high-frequency sampling strategy for a high-permeability medium attribute region, and adopt a low-frequency sampling strategy for a low-permeability medium attribute region;

[0016] An azimuth estimation module is configured to dynamically configure parameters of each sensor node according to the sampling frequency adjustment scheme, collect medium fluctuation signals by using the configured sensor node parameters, and determine an azimuth angle estimation value of a fluctuation source based on time sequence analysis and calculation of signal time delay;

[0017] An azimuth correction module is configured to fuse data of adjacent sensor nodes in a communication network according to the azimuth angle estimation value of the fluctuation source, correct an azimuth angle estimation deviation by using an optimization algorithm, and obtain an accurate azimuth angle of the fluctuation source.

[0018] An energy consumption optimization module is configured to evaluate a signal loss risk based on the accurate azimuth angle of the fluctuation source and a type distribution of the medium attribute region, and generate an energy consumption optimization transmission scheme according to an evaluation result of the loss risk.

[0019] A feedback control module is configured to extract a resource consumption index from the energy consumption optimization transmission scheme, feed back the resource consumption index to the sampling frequency adjustment scheme through a feedback loop mechanism, and realize dynamic balance between monitoring system accuracy and energy consumption.

[0020] The present application provides a karst aquifer monitoring method and system based on Internet of Things, which solves the technical problems of poor self-adaptability, low positioning accuracy, and high energy consumption of the existing monitoring system in a strong heterogeneous environment of a karst aquifer, improves the positioning accuracy of a fluctuation source, and significantly prolongs the overall endurance of the monitoring system. Compared with the prior art, the technical scheme of the present application has at least the following advantages:

[0021] First, by analyzing the propagation speed and attenuation characteristics of medium fluctuation signals, dynamic identification and division of high and low permeability medium attribute regions are realized, providing accurate geological basis for subsequent differentiated monitoring strategies.

[0022] Second, based on the type distribution of the medium attribute region, a sampling strategy combining high and low frequencies is adaptively implemented, which not only ensures the capture ability of rapidly changing signals in high permeability regions, but also avoids data redundancy and resource waste in low permeability regions.

[0023] Third, by fusing time sequence data of multiple sensor nodes and using an optimization algorithm to correct estimation deviation, the uncertainty of signal propagation in heterogeneous media is overcome, and the positioning accuracy of the azimuth angle of the fluctuation source is improved.

[0024] Fourth, by establishing a signal loss risk evaluation mechanism and introducing a feedback loop of resource consumption index, the data transmission priority and system power consumption are dynamically adjusted, the continuous optimization and balance of monitoring accuracy and energy consumption are realized, and the long-term working ability of the system in the field environment is enhanced. BRIEF DESCRIPTION OF DRAWINGS

[0025] In order to more clearly illustrate the technical solutions of the embodiments of the present application, the drawings needed to be used in the embodiment description will be briefly introduced. Obviously, the drawings in the following description are some embodiments of the present application, and other drawings can be obtained by those skilled in the art without creative labor.

[0026] Figure 1 The flowchart of the karst aquifer monitoring method based on the Internet of Things in the present application;

[0027] Figure 2 The flowchart of generating a sampling frequency adjustment scheme in the karst aquifer monitoring process based on the Internet of Things in the present application;

[0028] Figure 3 The core performance index comparison in the present application;

[0029] Figure 4 The data acquisition efficiency comparison in the present application;

[0030] Figure 5 The system endurance comparison in the present application;

[0031] Figure 6 The system comprehensive performance score in the present application;

[0032] Figure 7 The structure diagram of the karst aquifer monitoring system based on the Internet of Things in the present application. DETAILED DESCRIPTION

[0033] The embodiments of the present application provide a karst aquifer monitoring method and system based on the Internet of Things. The terms "first", "second", "third", "fourth" and the like (if any) in the specification and claims of the present application and the above drawings are used to distinguish similar objects, and do not necessarily have to be used to describe a specific order or sequence. It should be understood that the data used in this way can be interchanged under appropriate circumstances, so that the embodiments described herein can be implemented in an order other than that illustrated or described herein. In addition, the term "includes" or "has" and any variations thereof are intended to cover non-exclusive inclusion, for example, a process, method, system, product or device including a series of steps or units does not have to be limited to those steps or units clearly listed, but can include other steps or units that are not clearly listed or inherent to these processes, methods, products or devices.

[0034] For ease of understanding, the specific flow of the embodiments of the present application will be described below. Please refer to Figure 1 One embodiment of the karst aquifer monitoring method based on the Internet of Things in the present application includes:

[0035] In step S101, medium fluctuation signals are collected from a plurality of sensor nodes through an Internet of Things communication network, propagation difference characteristics of the medium fluctuation signals in different medium attribute regions are extracted, and a medium attribute region type distribution is determined according to a comparison result of the propagation difference characteristics and a preset region threshold.

[0036] In an embodiment, the step of extracting the propagation difference characteristics of the medium fluctuation signals in different medium attribute regions can specifically include the following steps:

[0037] Medium fluctuation signal data are collected from a plurality of sensor nodes deployed in a karst aquifer monitoring region through an Internet of Things communication network, and first signal data sets are obtained after data preprocessing;

[0038] Frequency domain characteristics of the first signal data sets are extracted by using a fast Fourier transform algorithm, and second signal data sets are obtained;

[0039] It is determined whether the frequency domain characteristics of the second signal data sets satisfy a preset signal-to-noise ratio threshold, and if yes, the second signal data sets are decomposed by using a wavelet transform algorithm to obtain a propagation difference characteristic set;

[0040] According to the propagation difference characteristic set, a propagation time difference of medium fluctuation signals between sensor nodes is calculated to obtain a propagation speed distribution;

[0041] The propagation speed distribution is matched with a preset medium region division model to determine boundaries of each medium region, and a medium attribute region division result is obtained;

[0042] The propagation difference characteristic set is analyzed by applying a preset attenuation model to calculate attenuation coefficients of medium fluctuation signals in each divided medium attribute region, and an attenuation degree distribution is obtained;

[0043] Based on the propagation speed distribution and the attenuation degree distribution, propagation difference characteristics of medium fluctuation signals in different medium attribute regions of a karst aquifer are generated.

[0044] Specifically, in the karst aquifer monitoring scene based on Internet of Things, the Internet of Things communication network is built using the LoRaWAN protocol, which is adapted to the terrain shielding problem that may exist in the karst aquifer monitoring area. The deployed sensor nodes use ultrasonic sensors and are distributed in the 10km×10km monitoring area with a grid density of 500m×500m. Each sensor node has a built-in GPS module to record the spatial coordinates (accuracy 1m). The initial frequency of 1kHz is used to collect medium fluctuation signal data. The collected data includes timestamp (accurate to milliseconds), signal amplitude and sensor node ID. The collected raw data is preprocessed, high-frequency interference signals are filtered using a 5th order Butterworth low-pass filter (cutoff frequency 500Hz), and abnormal data points exceeding 3 times the standard deviation of the normal fluctuation range are removed by the 3σ criterion. The first signal data set is obtained after preprocessing. The fast Fourier transform algorithm is applied to the first signal data set, the calculation window length is set to 1024 points, and the overlap rate is 50%. The time domain signal is converted to the frequency domain signal, the signal main frequency, frequency bandwidth and each frequency component amplitude ratio are extracted, and the second signal data set is formed. A preset signal-to-noise ratio threshold (20dB based on historical monitoring data statistics of karst aquifer) is introduced to calculate the ratio of effective signal amplitude to noise signal amplitude in the second signal data set. If it meets the threshold requirement, the Daubechies D4 wavelet basis function is used to perform 5-layer wavelet decomposition on the second signal data set. This decomposition process can capture the local characteristics of the signal, and the propagation difference feature set containing the signal mutation point position, different frequency component propagation time delay and other information is obtained after decomposition.

[0045] According to the propagation difference feature set, the key feature points such as peak points and mutation points of the signals collected by each sensor node are extracted, the collection time difference of the same feature point between different sensor nodes is calculated, and the spatial distance between the sensor nodes is obtained by the GPS coordinates and the distance formula between two points. According to the relationship of “propagation speed = distance between sensor nodes / propagation time difference”, the propagation speed distribution of different regions is derived, for example, the distance between sensor nodes in a certain grid unit is 500m, the propagation time difference is 1.47ms, and the corresponding propagation speed is about 340m / s. The distance between sensor nodes in another grid unit is 500m, the propagation time difference is 1.72ms, and the corresponding propagation speed is about 290m / s. The propagation speed distribution is matched with the preset medium region division model (based on known geological data of karst aquifer, including limestone fissure area 320-360m / s, cave filling area 290-320m / s, and clay layer 260-290m / s). By comparing the propagation speed of each grid unit with the speed range of different medium properties in the model, the boundary of each medium region is determined, and the boundary is fitted by the coordinates of adjacent different types of grid units to generate the final medium property region division result. At the same time, the propagation difference feature set is analyzed using a preset attenuation model. The attenuation model uses the formula ,in β This is the attenuation coefficient (unit: dB / m). Let be the initial amplitude of the signal at a certain sensor node. Let L be the amplitude of the signal after propagation to adjacent sensor nodes, and L be the spatial distance between the two sensor nodes. Select the amplitude corresponding to the same wave event from the set of propagation difference features. , The data, combined with the distance L between sensor nodes, is used to calculate the attenuation coefficient of each divided region. For example, a group of sensor nodes in a limestone fracture zone... =5V =4.8V, L=500m, substituting into the formula, we get... A group of sensor nodes in a clay layer =5V =4.2V, L=500m, calculated as follows This leads to the distribution of attenuation levels, which corresponds one-to-one with the results of the media property region division. Each media region is labeled with its average attenuation coefficient and coefficient fluctuation range.

[0046] Based on the distribution of propagation velocity and attenuation degree, the spatial correspondence between the two types of data is integrated: the limestone fracture zone corresponds to a propagation velocity range of 320–360 m / s and an attenuation coefficient range of 0.003–0.004 dB / m; the karst cave filling zone corresponds to a propagation velocity range of 290–320 m / s and an attenuation coefficient range of 0.005–0.007 dB / m; and the clay layer corresponds to a propagation velocity range of 260–290 m / s and an attenuation coefficient range of 0.007–0.009 dB / m. The two types of data are superimposed using a geographic information system to generate the propagation difference characteristics of the medium wave signal in different medium property regions of the karst aquifer. The characteristic data includes the spatial location of each medium region, the statistical value of propagation velocity, the statistical value of attenuation coefficient, and the corresponding original signal segment number.

[0047] In one specific embodiment, determining the distribution of medium attribute region types based on the comparison result of the propagation difference characteristics and a preset region threshold may specifically include the following steps:

[0048] Load a preset region threshold, which includes a preset propagation speed threshold and a preset attenuation coefficient threshold;

[0049] Determine whether the propagation speed in the propagation difference characteristics of any medium attribute region exceeds the preset propagation speed threshold. If yes, mark the corresponding region as a high-permeability medium attribute region; otherwise, mark the corresponding region as a low-permeability medium attribute region.

[0050] Integrate the labeling results of all regions and the propagation difference characteristics, associate the spatial coordinate information of the corresponding region, and generate the medium attribute region type distribution.

[0051] Specifically, in the karst aquifer monitoring scene based on the Internet of Things, when the medium wave signal propagation difference characteristic extraction is completed, the step of determining the medium attribute region type distribution is entered. First, a preset region threshold is loaded, which is based on the historical geological exploration data of the karst aquifer and a large number of monitoring sample statistics, stored in the system database, and includes two types of core parameters, namely, a preset propagation speed threshold and a preset attenuation coefficient threshold. The preset propagation speed threshold is set to 310 m / s, and the preset attenuation coefficient threshold is set to 0.006 dB / m. The two thresholds correspond to the division standards of medium permeability characteristics and signal attenuation characteristics, and are directly associated with the subsequent medium attribute region labeling logic.

[0052] After loading the threshold, the propagation difference characteristic data generated in the early stage is called, which includes the propagation speed value, attenuation coefficient value and corresponding spatial coordinate information of each monitoring sub-region. Each sub-region is divided into a 500m×500m grid, covering a 10km×10km monitoring range. Each sub-region corresponds to a unique spatial coordinate code (such as X01Y01, X01Y02), and each coordinate code is associated with a set of propagation speed values and attenuation coefficient values. The system extracts the propagation speed value of each sub-region in order of spatial coordinate code, and compares it with the preset propagation speed threshold (310 m / s). The comparison logic uses direct comparison of numerical values: if the propagation speed value of a sub-region is greater than 310 m / s, for example, the propagation speed of sub-region X02Y03 is 335 m / s, which satisfies the condition of "exceeding the preset propagation speed threshold", then the sub-region is marked as a high-permeability medium attribute region; if the propagation speed value of a sub-region is less than or equal to 310 m / s, for example, the propagation speed of sub-region X05Y04 is 290 m / s, which does not satisfy the condition of "exceeding the preset propagation speed threshold", then the sub-region is marked as a low-permeability medium attribute region. In the labeling process, the attenuation coefficient values of each sub-region are recorded synchronously, although the attenuation coefficient values do not directly participate in attribute labeling, but need to be stored as auxiliary data, for example, the attenuation coefficient of high-permeability medium attribute region X02Y03 is 0.004 dB / m, and the attenuation coefficient of low-permeability medium attribute region X05Y04 is 0.008 dB / m, to ensure the complete correspondence between the propagation difference characteristics and the region labeling results.

[0053] After all sub-regional labels are completed, the data integration process is started. The system first establishes a four-dimensional correlation table of "spatial coordinate code - regional label - propagation speed - attenuation coefficient", matches the spatial coordinate code of each sub-region with its corresponding regional label result (high permeability / low permeability), propagation speed value, and attenuation coefficient value one by one, and forms a structured data set. Then, the geographic information system (GIS) interface is called to map the spatial coordinate code in the four-dimensional correlation table to the specific location on the monitoring area map, while the regional label result is displayed in different color layers (such as high permeability area in red and low permeability area in blue), and the propagation speed value and the attenuation coefficient value are bound to the corresponding area in the form of a floating data box. Clicking on any sub-region on the map can view the complete propagation difference characteristic data. The final generated medium property regional type distribution is stored in a vector map file format, which supports direct calling of subsequent sampling frequency adjustment, transmission strategy optimization, etc. The file contains data version number and generation timestamp, which is convenient for tracing data source and updating record.

[0054] In step S102, according to the medium property regional type distribution, a high-frequency sampling strategy is adopted for the high-permeability medium property region, and a low-frequency sampling strategy is adopted for the low-permeability medium property region, to generate a sampling frequency adjustment scheme.

[0055] In a specific embodiment, the step S102 can specifically include the following steps:

[0056] For the medium property regional type distribution, a hierarchical analysis method is used to classify the medium property regional type to obtain a classified region set;

[0057] Based on the classified region set, a support vector machine algorithm is used to distinguish the boundaries of the high-permeability medium property region and the low-permeability medium property region, to determine the specific range of each region;

[0058] For the specific range of each region, the historical data of the medium fluctuation signal in each region is used to calculate the change rate of the medium fluctuation signal in each region, to obtain a signal change rate distribution;

[0059] It is judged whether the signal change rate of a region in the signal change rate distribution is greater than a preset rate threshold;

[0060] If yes, a high-frequency sampling strategy is assigned to the high-permeability medium property region, and a high-frequency sampling parameter set is determined, and if no, a low-frequency sampling strategy is assigned to the low-permeability medium property region, and a low-frequency sampling parameter set is determined, wherein the sampling parameter set includes sampling frequency, sampling period, and sampling duration each time;

[0061] According to the high-frequency sampling parameter set and the low-frequency sampling parameter set, an initial scheme of the sampling frequency corresponding to each region is generated;

[0062] The sampling frequency of adjacent regions in the initial sampling frequency scheme is smoothed by using a linear interpolation method to generate a sampling frequency adjustment scheme.

[0063] Specifically, please refer to Figure 2 In the karst aquifer monitoring scene based on the Internet of Things, after determining the distribution of the medium attribute region type, the layered analysis method is used to classify the medium attribute region type according to the distribution. The layered analysis method takes the medium attribute label as the layering basis, splits the 400 sub-region data containing a 10km*10km monitoring range and divided by 500m*500m into two attribute labels of “high permeability” and “low permeability”, forms a high permeability region layer and a low permeability region layer, and each layer of data is associated with the spatial coordinates (such as X01Y01 to X20Y20) of the sub-region, the propagation speed and the attenuation coefficient. Thus, a classified region set is obtained, and the two types of region data in the classified region set are not crossed and collectively cover all monitoring sub-regions. Based on the classified region set, a support vector machine algorithm is used to distinguish the boundaries of the high permeability medium attribute region and the low permeability medium attribute region. The input features of the support vector machine algorithm are the spatial coordinates (x, y) of each sub-region, and the output label is “1” (high permeability) or “0” (low permeability). The kernel function of the algorithm is selected as the radial basis function (RBF), the penalty parameter C is set to 10, and the gamma value is set to 0.1. Through training on 80% of the sub-region data in the classified region set, a classification model is established, and the remaining 20% data is used to verify the accuracy of the model (requiring ≥95%). After the model training is completed, the spatial coordinates of all sub-regions in the monitoring region are input, and the attribute prediction label corresponding to each coordinate point is output. The boundary candidate points are determined according to the coordinate points with varying prediction labels, and the boundary candidate points are linearly fitted to obtain a continuous boundary line coordinate set (such as X05Y03 to X05Y18 to X12Y18). The coordinate points are connected to form the region boundary, and the specific range of the high permeability and low permeability regions is determined. The attribute labels of the sub-regions on both sides of the boundary line are completely distinguished, and there is no attribute ambiguous region.

[0064] According to the specific range of each region, the signal change rate is calculated using the historical data of the medium wave signal in each region. The historical data is the signal data collected by each sensor node at an initial sampling frequency of 1kHz in the past 60 days. Each sub-region corresponds to one sensor node, and the amount of historical data is 60*86400*1000 groups for each sensor node (calculated according to 86400 seconds per day and 1kHz sampling frequency). The historical data of all sensor nodes in the corresponding region is selected according to the region range. For the historical data of a single sensor node, 100 consecutive sampling values are taken in time sequence, and the signal amplitude difference between the 100th group and the 1st group is calculated , and then divided by the time interval (99 ms, 1 ms interval for 1 kHz sampling, 100 groups of data interval 99 ms), get the signal change rate of the sensor node , take the average value of V of all sensor nodes in a certain area as the signal change rate of the area, and store the signal change rates of all areas according to the spatial coordinates to obtain the signal change rate distribution, for example, the signal change rate of the high permeability area X03Y06 is 0.9V / s, and the signal change rate of the low permeability area X15Y12 is 0.4V / s.

[0065] Call the preset rate threshold (based on the statistical fluctuation characteristics of different media of karst aquifer, set to 0.6V / s), judge whether the signal change rate of each area is greater than the threshold. If the signal change rate of a certain area is greater than 0.6V / s, such as high permeability area X03Y06 (0.9V / s), assign high frequency sampling strategy, determine high frequency sampling parameter set: sampling frequency is set to 60Hz (sampling interval ≈16.67ms), sampling period is set to 5 minutes (start sampling once every 5 minutes), and each collection duration is set to 20 seconds (1200 groups of data are collected each time); if the signal change rate of a certain area is less than or equal to 0.6V / s, such as low permeability area X15Y12 (0.4V / s), assign low frequency sampling strategy, determine low frequency sampling parameter set: sampling frequency is set to 25Hz (sampling interval 40ms), sampling period is set to 20 minutes, and each collection duration is set to 10 seconds (250 groups of data are collected each time), sampling parameter set is associated with area spatial coordinates one by one to ensure that parameters match corresponding areas. According to the high and low frequency sampling parameter set, the parameters are distributed to each sub-area according to the spatial range of the area, and each sub-area corresponds to a group of sampling parameters to generate an initial sampling frequency scheme. In the initial scheme, the sampling frequency at the boundary of high permeability and low permeability areas jumps (such as adjacent sub-areas X05Y03 is 60Hz, X05Y04 is 25Hz), which needs to be smoothed by linear interpolation method. Select 2 adjacent sub-areas (set as sub-area P and sub-area Q, interval 500m) at the boundary, the sampling frequency of P , the sampling frequency of Q Hz, 10 interpolation points are divided between P and Q at an interval of 50m. The sampling frequency of the nth interpolation point is , for example, n=1 Hz, n=5 Hz. The sampling frequency corresponding to the interpolation point is supplemented to the initial scheme, so that the sampling frequency of the boundary area is continuously and gradually changed from 60Hz to 25Hz, and finally the sampling frequency adjustment scheme is generated, which contains the sampling frequency, period and collection duration of all sub-areas, and is associated with spatial coordinates and parameter effective time.

[0066] Step S103, dynamically configuring parameters of each sensor node according to the sampling frequency adjustment scheme, collecting medium fluctuation signals by using the configured sensor node parameters, calculating signal time delay based on time series analysis, and determining the azimuth angle estimation value of the fluctuation source.

[0067] In a specific embodiment, the step S103 can specifically include the following steps:

[0068] Extracting the dynamic parameter configuration corresponding to each sensor node from the sampling frequency adjustment scheme, the dynamic parameter configuration including sampling frequency, signal gain, data acquisition period, and each acquisition duration;

[0069] Using a data compression algorithm to encode the dynamic parameter configuration to obtain a compressed parameter set, and then transmitting the compressed parameter set to each sensor node in a broadcast manner;

[0070] Each sensor node decodes the compressed parameter set, updates its own operating parameters according to the decoding result, and obtains a parameter synchronization state;

[0071] According to the parameter synchronization state, activating each sensor node, and making each sensor node collect the medium fluctuation signals of the karst aquifer according to the updated parameters to generate a multi-channel medium fluctuation signal data set;

[0072] Based on the multi-channel medium fluctuation signal data set, using a time series analysis method to calculate the time delay between the signal data collected by any two sensor nodes to obtain a time delay distribution set;

[0073] Judging whether each time delay value in the time delay distribution set is greater than a preset delay threshold, if yes, using a triangular positioning algorithm to calculate the azimuth angle of the fluctuation source to obtain a preliminary estimation value of the azimuth angle of the fluctuation source;

[0074] Smoothing the preliminary estimation value of the azimuth angle of the fluctuation source using a Kalman filtering algorithm to determine the azimuth angle estimation value of the fluctuation source.

[0075] Specifically, the dynamic parameter configuration is encoded by using a data compression algorithm. The LZ77 compression algorithm is selected, the sliding window size is set to 4096 bytes, the lookup buffer size is set to 1024 bytes, the "sensor node ID-dynamic parameter configuration" mapping table is converted into a binary data stream and then compressed, the data amount after compression is reduced from the original 12 KB to 3.2 KB, and the compressed parameter set is obtained. The compressed parameter set is transmitted to each sensor node in a broadcast manner through an Internet of Things communication network (using the LoRaWAN protocol, the transmission rate is 250 kbps, and the transmission power is 17 dBm). The broadcast period is set to 10 seconds, and it is ensured that all sensor nodes receive the complete compressed parameter set within 1 minute. After each sensor node receives the compressed parameter set, the built-in LZ77 decoding module decodes according to the same algorithm parameters (sliding window 4096 bytes, lookup buffer 1024 bytes) as the encoding end, restores the "sensor node ID-dynamic parameter configuration" mapping table, and the sensor node matches the corresponding dynamic parameters according to its own ID, updates the sampling frequency, signal gain, sampling period and sampling duration parameters in the internal register, and sends a "parameter update completion" confirmation signal to the control center after the update is completed. The control center counts the number of confirmation signals, and when the number of sensor nodes receiving the confirmation signals accounts for ≥98% of the total number of sensor nodes, it is determined that the parameter synchronization is completed, the parameter synchronization state is obtained, and the unsynchronized sensor nodes will trigger a secondary broadcast to supplement the transmission.

[0076] According to the parameter synchronization state, the control center sends an "activation collection" instruction to all sensor nodes that have completed synchronization. After receiving the instruction, the sensor nodes collect the karst aquifer medium wave signal according to the updated parameters: the high-permeability area sensor nodes collect 20 seconds of signal every 5 minutes (1200 groups of data each time) at a sampling frequency of 60 Hz and a signal gain of 1.2 dB, and the low-permeability area sensor nodes collect 10 seconds of signal every 20 minutes (250 groups of data each time) at a sampling frequency of 25 Hz and a signal gain of 1.2 dB. The collected signal data is stored in the format of "sensor node ID-collection timestamp-signal amplitude", uploaded to the control center through the Internet of Things network, and the control center aggregates all sensor node data according to the collection timestamp to generate a multi-channel medium wave signal data set. The number of channels in the data set is consistent with the number of sensor nodes (400 channels). Based on the multi-channel medium wave signal data set, a time delay between the signal data collected by any two sensor nodes is calculated using a time series analysis method. The cross-correlation algorithm is selected, the window length is set to 512 points, and the overlap rate is set to 50%. The signal data of any two sensor nodes (such as sensor node A and sensor node B) is cross-correlated, and the time point corresponding to the peak value of the cross-correlation function is found. This time point is the time difference between the signals of the two sensor nodes, which is the time delay value. The time delay value is calculated for all sensor node pairs, stored in the format of "sensor node pair-time delay value", and a time delay distribution set is obtained. For example, the time delay value between sensor node A and sensor node B is 0.023 seconds, and the time delay value between sensor node A and sensor node C is 0.031 seconds.

[0077] A preset delay threshold (set to 0.01 seconds based on the signal propagation characteristics of karst aquifers) is invoked to determine whether each time delay value in the time delay distribution set is greater than the threshold. If the time delay value of a certain sensor node pair is greater than 0.01 seconds, such as sensor node A and sensor node B (0.023 seconds) or sensor node A and sensor node C (0.031 seconds), then at least 3 such sensor node pairs (such as AB, AC, BC) are selected, and the azimuth angle of the wave source is calculated using a triangulation algorithm. The triangulation algorithm establishes a planar coordinate system with the spatial coordinates of three sensor nodes (e.g., sensor node A: X1=1000m, Y1=1000m; sensor node B: X2=1500m, Y2=1000m; sensor node C: X3=1000m, Y3=1500m) as vertices. The distance difference from the wave source to each sensor node is calculated based on the time delay (distance difference = signal propagation speed × time delay; the signal propagation speed in karst aquifers is taken as 3000m / s, e.g., distance difference between A and B = 3000 × 0.023 = 69m). A hyperbola equation is established based on the distance difference, and the coordinates of the intersection point (wave source location) are solved. Then, the azimuth angle is calculated using the wave source location and the coordinates of a reference sensor node (e.g., sensor node A) to obtain a preliminary estimate of the wave source azimuth angle, for example, 42.5°. The preliminary estimate of the wave source azimuth angle is smoothed using a Kalman filter algorithm, with the azimuth angle set as the state variable; the state equation is... (in, For the first k The true azimuth angle at the time of the second observation. For the first The estimated azimuth angle output from the second filter. The process noise has a mean of 0 and a standard deviation of . Gaussian distribution); the observation equation is set as follows: ,in For the first k The initial estimate of the azimuth angle obtained this time (such as the result calculated by the triangulation algorithm). For the observation matrix, To observe the noise, it follows a mean of 0 and a standard deviation of . The Gaussian distribution is used. The preliminary estimate obtained each time is used as the observation value and input into the Kalman filter algorithm. The prediction and update steps are performed in sequence to calculate the filtered azimuth value. After processing 10 sets of preliminary estimates, the output stable value is the azimuth estimate of the wave source. For example, the final estimate is 42.3°, and the error range is controlled within ±0.5°.

[0078] Step S104, according to the fluctuation source azimuth angle estimation value, fusing the data of adjacent sensor nodes in the communication network, correcting the azimuth angle estimation deviation through an optimization algorithm, and obtaining the accurate fluctuation source azimuth angle.

[0079] In a specific embodiment, the step S104 can specifically include the following steps:

[0080] Collecting multi-channel medium fluctuation signal data of a monitoring area corresponding to the fluctuation source azimuth angle estimation value, removing noise from the multi-channel medium fluctuation signal data, and extracting effective time sequence signal data to obtain a time sequence signal data set;

[0081] Using a fast Fourier transform method to decompose the frequency components of the time sequence signal data set to obtain a frequency distribution set;

[0082] If there is a frequency peak value in the frequency distribution set that meets a preset significance condition, then a phase difference between each channel is calculated through a frequency domain analysis method to obtain a phase difference set;

[0083] Using a least square method to fit the corresponding relationship between the phase difference set and the sensor node position to obtain a fluctuation source preliminary position estimation;

[0084] Using a weighted average algorithm to correct the deviation of the fluctuation source preliminary position estimation to obtain the accurate position coordinates of the fluctuation source;

[0085] According to the accurate position coordinates of the fluctuation source, the relative angles of each sensor node and the fluctuation source are calculated to generate an optimized azimuth angle set;

[0086] Using a Kalman filter algorithm to smooth the optimized azimuth angle set to obtain the accurate fluctuation source azimuth angle.

[0087] Specifically, according to the previously obtained wave source azimuth estimation value, multi-channel medium wave signal data corresponding to the monitoring area is collected. The monitoring area is deployed with 400 sensor nodes in a 500m x 500m grid, and the collection range is a circular area with a radius of 2km centered on the direction corresponding to the azimuth estimation value, involving 36 sensor nodes. These sensor nodes operate according to the updated sampling parameters: sensor nodes in the high permeability area sample at a frequency of 60Hz, and sensor nodes in the low permeability area sample at a frequency of 25Hz. The single collection time is set to 30 seconds, and finally 36 channels of medium wave signal data are generated, in which the high permeability area sensor nodes generate 1800 groups of data per channel, and the low permeability area sensor nodes generate 750 groups of data per channel. The collected multi-channel data is denoised, a 5th order Butterworth low-pass filter is used, the cutoff frequency is set to 100Hz, and high-frequency interference signals above 100Hz are filtered out; then through the 3σ criterion, abnormal values with amplitude exceeding ±3 times the standard deviation of all data are removed, and effective signals directly related to medium wave are retained, arranged in the format of "channel number-collection timestamp-signal amplitude", and formed into a time series signal data set. The data set corresponds one-to-one to the 36 sensor nodes, and each channel data contains complete time series and signal amplitude information.

[0088] The Fast Fourier Transform (FFT) method is used to decompose the frequency components of the time-series signal dataset. The calculation window length is set to 1024 points, with an overlap rate of 50%. Transformation parameters are adjusted according to the actual sampling frequency of each channel to convert the time-domain signal to a frequency-domain signal. The signal amplitude corresponding to each frequency point is extracted and stored in the format of "channel number-frequency-amplitude," resulting in a frequency distribution set that reflects the signal strength distribution at different frequencies. A pre-defined significance condition is established: the amplitude at a certain frequency point must be greater than twice the average amplitude of all amplitudes in the frequency distribution set for that channel, and the amplitudes of the three consecutive adjacent frequency points before and after that frequency point must show a decreasing trend. This condition is used to determine whether a frequency peak exists in the frequency distribution set. For example, if the amplitude of channel 1 at 8Hz is 1.2V, exceeding twice the average amplitude of 0.5V for that channel, and the amplitudes at adjacent 5Hz, 6Hz, and 7Hz are 0.4V, 0.6V, and 0.9V respectively, the significance condition is met, and a valid frequency peak is determined to exist. This peak frequency will be used as the benchmark for subsequent phase difference calculations. If a frequency peak exists that meets the conditions, the phase difference between each channel is calculated using the frequency corresponding to that peak (e.g., 8Hz) as the reference. Before the calculation, the timestamps of each sensor node are calibrated using the NTP protocol of the IoT network to ensure that the time synchronization accuracy is controlled within ±1ms, thus avoiding the impact of time deviation on the phase difference calculation. Then, the phase information of the frequency domain signal of each channel at the reference frequency is extracted, and the phase difference value of any two channels is calculated. For example, the phase difference between channel 1 and channel 2 is 15°, and the phase difference between channel 1 and channel 3 is 22°. The phase difference is recorded in the format of "channel pair - phase difference" to form a phase difference set. The phase difference value of each channel pair in this set is associated with the time difference of the two sensor nodes receiving the same wave signal.

[0089] The least squares method is used to fit the correspondence between the phase difference set and the sensor node position. The fitting process is based on the physical relationship between phase difference and distance difference: let the coordinates of the wave source position be... The coordinates of a certain sensor node are The distance between the sensor node and the wave source is The phase difference between two sensor nodes is positively correlated with the distance difference, combined with the signal wavelength. By constructing a correlation model and minimizing the sum of squared errors between the calculated and actual phase difference values ​​using the least squares method, the optimal coordinates of the wave source location are found. This is the initial location estimate of the wave source. A weighted average algorithm is used to correct the deviation of the initial location estimate. The weights are determined based on the signal-to-noise ratio (SNR) of the signals acquired by each sensor node. ( The average amplitude of the sensor node signal. (mean of noise amplitude) The larger the value, the higher the weight. The larger the signal-to-noise ratio (SNR) of the sensor nodes, the less interference from low SNR sensor node data on the correction results. Ten adjacent sensor nodes surrounding the initial position estimate are selected, and their signal data are reacquired and the SNR calculated. For example, sensor node 1 has an SNR of 35 dB, and sensor node 2 has an SNR of 32 dB. Then, based on the deviation between the calculated position of the wave source and the initial position from each sensor node, a weighted average is performed using corresponding weights to obtain the precise position coordinates of the wave source. For example, the precise position coordinates are (1502 m, 2001 m). These coordinates, by fusing high-confidence data from multiple sensor nodes, correct the deviation in the initial position estimate.

[0090] Based on the precise location coordinates, the relative angles between each sensor node and the wave source are calculated. Ten key sensor nodes are selected and evenly distributed around the wave source. The location coordinates of each sensor node are obtained; for example, the coordinates of sensor node A are (1000m, 2001m), and the coordinates of sensor node B are (1502m, 1500m). The coordinates are then calculated using the formula... ( For sensor nodes i The relative angle with respect to the wave source, To calculate the precise coordinates of the wave source, the angle between each sensor node and the wave source is calculated. For example, the relative angle between sensor node A and the wave source is 180°, and that of sensor node B is 270°. The relative angles of all sensor nodes are organized into "sensor node number - relative angle" to generate an optimized azimuth angle set. Each angle value in this set directly corresponds to the spatial position of the sensor node and the wave source. A Kalman filter algorithm is used to smooth the optimized azimuth angle set. The state variable for the next iteration is the relative angle. The state equation is ( The noise is the process noise, which follows a mean of 0 and a standard deviation of . (Gaussian distribution), the observation equation is: ( To optimize the angle values ​​in the azimuth set, For the observation matrix, To observe the noise, it follows a mean of 0 and a standard deviation of . (Gaussian distribution). The 10 angle values ​​in the optimized azimuth angle set are continuously filtered. The mean of the optimized azimuth angle set is used as the initial estimate. The final stable angle value is the accurate azimuth angle of the fluctuation source, for example, 42.1°, with the error range controlled within ±0.3°.

[0091] Step S105: Based on the precise azimuth angle of the fluctuation source and the distribution of the medium property region type, assess the signal loss risk and generate an energy-optimized transmission scheme based on the loss risk assessment results.

[0092] In a specific embodiment, the step S105 can specifically include the following steps:

[0093] Collecting medium fluctuation signal data of the multi-sensor nodes in the monitoring range corresponding to the accurate fluctuation source azimuth, and extracting the propagation time delay of the medium fluctuation signal by using a time domain analysis method to obtain a time delay distribution set;

[0094] According to the time delay distribution set, combining the medium attribute parameters of each region, and using a path loss model to calculate the attenuation degree of the medium fluctuation signal in each region to obtain a propagation loss distribution;

[0095] Based on the propagation loss distribution and the time delay distribution set, using a Monte Carlo method to simulate the propagation path of the medium fluctuation signal to obtain a loss risk assessment value;

[0096] Judging whether the loss risk assessment value is higher than a preset risk threshold value;

[0097] If yes, adjusting the transmission protocol parameters of the communication network, using a dynamic scheduling method to improve the data transmission priority of the corresponding sensor nodes, and determining a high-priority transmission scheme;

[0098] If no, analyzing the transmission load of the current communication network, using a fixed scheduling method to maintain the standard transmission priority of each sensor node, and determining a standard transmission scheme;

[0099] Using a linear programming method to adjust the transmission power of the sensor nodes, optimizing the energy consumption configuration of the high-priority transmission scheme and the standard transmission scheme, and obtaining an optimized energy consumption configuration;

[0100] According to the optimized energy consumption configuration, generating scheduling instructions for data transmission between each sensor node, integrating the scheduling instructions with the corresponding transmission scheme, and generating an energy consumption optimized transmission scheme.

[0101] Specifically, based on the previously obtained precise azimuth angle of the fluctuation source and the distribution of media property types, media fluctuation signal data from multiple sensor nodes within the monitoring range corresponding to that azimuth angle were collected. The monitoring range covers a circular area with a radius of 3 km centered on the fluctuation source, involving 52 sensor nodes. The sensor nodes are deployed in a 500m×500m grid, including 28 sensor nodes in high-permeability media property areas and 24 sensor nodes in low-permeability media property areas. All sensor nodes collect signals according to the updated sampling parameters (60Hz in high-permeability areas and 25Hz in low-permeability areas), with a single acquisition duration of 20 seconds, generating 52 channels of media fluctuation signal data. The propagation delay of signals from each sensor node is extracted using time-domain analysis. Based on the location of the wave source, the time difference between the received signal of each sensor node and the signal generated by the wave source is calculated. For example, the delay of sensor node 1 is 0.0021 seconds and the delay of sensor node 2 is 0.0035 seconds. The delay is stored in the format of "sensor node number - propagation delay" to obtain a set of delay distributions. The delay values ​​in this set are directly related to the distance from the sensor node to the wave source and the properties of the medium in the area. Sensor nodes in high-permeability areas generally have lower delay values ​​than those in low-permeability areas because the signal propagation speed is faster.

[0102] Based on the time delay distribution set and combined with the medium property parameters of each region (medium density 1.2 g / cm³ and porosity 25% in high-permeability regions; medium density 1.8 g / cm³ and porosity 12% in low-permeability regions), the attenuation degree of the medium fluctuation signal in each region is calculated using a path loss model. The logarithmic distance model is selected as the path loss model, and the formula is as follows: (L is the path loss value,) For reference distance Loss value at 100m, high permeability zone =2dB, low permeability area =5dB; The path loss index represents high-penetration areas. =1.8, low permeability area =2.5; The distance from the sensor node to the wave source is calculated from the time delay and the signal propagation speed. The distances from each sensor node to the wave source are then calculated. Value, corresponding region and Substituting into the formula, the path loss value of each sensor node is calculated. For example, sensor node 1 in a high-permeability area. =6.3km dB; Sensor node 2 in low-permeability region =8.75km dB, and a propagation loss distribution is obtained. Based on the propagation loss distribution and the set of delay distributions, the propagation path of the medium fluctuation signal is simulated by using a Monte Carlo method. The simulation parameters are set as follows: 1000 iterations, 100 random propagation paths generated in each iteration, the medium attribute parameters of each region are referred to in the path generation (the path generation probability of the high permeability region is 60%, and the path generation probability of the low permeability region is 40% because the signal propagation resistance is small in the high permeability region), and the path weight is adjusted in combination with the propagation loss distribution (the smaller the loss value, the greater the path weight of the region). The total loss value of each path is calculated in each iteration, and the proportion of the paths with a total loss value exceeding a preset loss threshold (45 dB) in all iterations is calculated. The proportion is the loss risk assessment value. For example, if the proportion of the paths with a total loss exceeding 45 dB in the simulation result is 18%, the loss risk assessment value is 18%, which reflects the probability of signal loss caused by loss in the propagation process.

[0103] A preset risk threshold is called, which is set to 15% based on the requirement of karst aquifer monitoring for signal integrity. It is judged whether the loss risk assessment value is higher than the threshold. If the loss risk assessment value is higher than 15% (such as 18%), the transmission protocol parameters of the communication network are adjusted, and the communication network adopts the LoRaWAN protocol. The adjustments include: the spreading factor is increased from 12 to 10 (to improve the transmission rate), the coding rate is adjusted from 4 / 5 to 4 / 3 (to reduce the coding redundancy to reduce the transmission time), and the receiving sensitivity is optimized from -148 dBm to -145 dBm (to enhance the signal receiving capability). At the same time, a dynamic scheduling method is used to improve the data transmission priority of the corresponding sensor nodes, the scheduling period is set to 10 seconds, the sensor nodes are sorted according to the loss risk value, the top 30% of the sensor nodes are assigned a high priority (priority level 1), and the remaining sensor nodes are assigned a medium priority (level 2), and a high priority transmission scheme is determined. If the loss risk assessment value is lower than or equal to 15%, the transmission load of the current communication network is analyzed, the data packet transmission amount of the network in the past 5 minutes is counted (such as an average of 120 data packets per second, and the maximum load of the network is 200 data packets per second), and it is determined that the current load rate is 60% and is in a low load state. A fixed scheduling method is used to maintain the standard transmission priority of each sensor node, all sensor nodes are uniformly assigned a standard priority (level 3), and the scheduling period is set to 20 seconds to determine a standard transmission scheme.

[0104] Regardless of the high priority transmission scheme or the standard transmission scheme, the sensor node transmission power needs to be adjusted by using a linear programming method to optimize the energy consumption configuration. The objective function of the linear programming is set as ( is the transmission power of the sensor node i , and is the transmission power of the sensor node itransmission length of a single day), constraints include: the sensor node transmission power the value range of the sensor node transmission power (5mW≤ ≤50mW, based on the hardware specifications of the sensor node), the signal-to-noise ratio of the sensor node transmission signal needs to be ≥20dB, and the total energy consumption of the network needs to be ≤5000mWh / day (based on the battery capacity). The transmission length of each sensor node (high priority sensor node single day transmission length 1800 seconds, standard priority sensor node 1200 seconds), the signal-to-noise ratio requirement is substituted into the constraint condition, and the optimal transmission power of each sensor node is obtained by solving the objective function, for example, the optimal transmission power of high priority sensor node 1 is 35mW, and the optimal transmission power of standard priority sensor node 2 is 20mW. Obtain the optimized energy consumption configuration, which minimizes the product of the sensor node transmission power and the transmission length under the premise of meeting the signal transmission quality, and reduces the total energy consumption. According to the optimized energy consumption configuration, generate the scheduling instructions for data transmission between each sensor node, and the instruction content includes sensor node ID, transmission priority, transmission power, scheduling period and transmission period (such as sensor node 1: priority 1, 35mW, 10 seconds period, 00:00~24:00 every 1 second transmission), integrate the scheduling instructions with the corresponding high priority transmission scheme or standard transmission scheme, and arrange them in the format of "scheme type-scheduling instruction set-energy consumption configuration parameter", and generate the energy consumption optimization transmission scheme. The scheduling instructions of each sensor node in the scheme are directly related to the medium attribute and loss risk assessment value of the area where it is located, ensuring the pertinence and energy consumption optimization effect of the scheme.

[0105] Step S106, extract the resource consumption index from the energy consumption optimization transmission scheme, and feed back the resource consumption index to the sampling frequency adjustment scheme through a feedback loop mechanism to realize the dynamic balance between the monitoring system accuracy and the energy consumption.

[0106] In a specific embodiment, the step S106 can specifically include the following steps:

[0107] extract the resource consumption index from the energy consumption optimization transmission scheme, the resource consumption index including transmission power, total energy consumption of data transmission, single transmission energy consumption and data transmission length;

[0108] based on the resource consumption index, use time series analysis to analyze the change trend of the resource consumption of each sensor node, and combine the sampling frequency parameter to obtain the correlation between the sampling frequency and the resource consumption;

[0109] According to the correlation, a resource consumption prediction model is constructed using a linear regression algorithm, and a resource consumption prediction value for a future preset period under the current sampling frequency is output;

[0110] determining whether the resource consumption prediction value is higher than a preset energy consumption threshold;

[0111] If yes, dynamically adjusting a sampling frequency parameter based on the correlation relationship to generate a high-frequency sampling optimization sub-scheme;

[0112] If no, maintaining the current sampling frequency parameter to generate a standard sampling maintenance sub-scheme;

[0113] For the high-frequency sampling optimization sub-scheme or the standard sampling maintenance sub-scheme, adopting a least square method to optimize data collection power of each sensor node to determine a power optimization configuration;

[0114] Feeding back the power optimization configuration and the sampling frequency parameter obtained after judgment to the sampling frequency adjustment scheme to update parameter configuration of the sampling frequency adjustment scheme;

[0115] Based on the updated sampling frequency adjustment scheme, re-collecting medium wave signal of karst aquifer and verifying to ensure that monitoring accuracy meets a preset accuracy threshold and resource consumption is lower than a preset energy consumption threshold, thereby realizing dynamic balance between monitoring system accuracy and energy consumption.

[0116] Specifically, based on the resource consumption index, time series analysis is adopted to analyze change trend of resource consumption of each sensor node, resource consumption data of the past 7 days is selected, time series is constructed per day, data is smoothed by using a moving average method (window size 3 days), and it is analyzed that resource consumption of sensor nodes in a high permeability area increases by 5% per day due to a sampling frequency of 60 Hz, and resource consumption of sensor nodes in a low permeability area increases by 2% per day due to a sampling frequency of 25 Hz. Meanwhile, in combination with current sampling frequency parameters of each sensor node, average resource consumption values under different sampling frequencies are counted, such as 92 mWh per sensor node per day corresponding to 60 Hz and 45 mWh per sensor node per day corresponding to 25 Hz, a positive correlation relationship between sampling frequency and resource consumption is obtained, that is, the higher the sampling frequency, the greater the resource consumption. According to the correlation relationship, a linear regression algorithm is adopted to construct a resource consumption prediction model, taking sampling frequency as an independent variable and resource consumption as a dependent variable, the model form is y = bx + a, b is a regression coefficient, and a is an intercept. Substituting statistical data (92 mWh corresponding to 60 Hz and 45 mWh corresponding to 25 Hz) into the model, it is calculated that b = 1.34, a = 11.6, and the model is y = 1.34x + 11.6. f y ​​​​​​​The preset time period is set to 7 days, the current sampling frequency of each sensor node is substituted into the model, and the resource consumption prediction value of the future 7 days is output. For example, the prediction value of the sensor node in the high permeability area is 1.34*60+11.6=92 mWh / day, and the total prediction value of 52 sensor nodes is 5024 mWh / day.

[0117] A preset energy consumption threshold is called, which is set to 5500 mWh per day based on the battery capacity of the sensor node (2000 mWh per sensor node, 104000 mWh per month for 52 sensor nodes), and it is judged whether the resource consumption prediction value (5024 mWh / day) is higher than the threshold. If the prediction value is higher than the threshold, for example, the sampling frequency is temporarily increased to 70 Hz due to the activity of the fluctuation source in a certain period, and the prediction value is 5600 mWh / day, then based on the correlation, the sampling frequency parameter is dynamically adjusted, and the sampling frequency of the sensor node in the high permeability area is reduced from 70 Hz to 65 Hz, and the sampling frequency of the sensor node in the low permeability area is maintained at 25 Hz, to generate a high-frequency sampling optimization sub-scheme; if the prediction value is lower than the threshold, for example, the current 5024 mWh / day, the current sampling frequency parameter (60 Hz for high permeability, 25 Hz for low permeability) is maintained, and a standard sampling maintenance sub-scheme is generated. For the high-frequency sampling optimization sub-scheme or the standard sampling maintenance sub-scheme, the least square method is used to optimize the data acquisition power of each sensor node. Taking the sensor node acquisition power as the optimization variable, the objective function is min∑(acquisition power*acquisition time length), and the constraint condition is that the signal-to-noise ratio of the acquisition signal is greater than or equal to 20 dB. Substituting the acquisition time length (20 seconds per time for high permeability, 10 seconds per time for low permeability) and the signal-to-noise ratio requirement of each sensor node, the optimized acquisition power of the sensor node in the high permeability area is 18 mW, and the optimized acquisition power of the sensor node in the low permeability area is 12 mW, and the power optimization configuration is determined. The power optimization configuration and the sampling frequency parameter obtained after judgment are fed back to the sampling frequency adjustment scheme, and the sampling frequency and acquisition power parameters in the scheme are updated, such as the sensor node parameter in the high permeability area is updated to 60 Hz and 18 mW, and the sensor node parameter in the low permeability area is updated to 25 Hz and 12 mW. Based on the updated scheme, the medium fluctuation signal is re-acquired for 30 minutes, the monitoring accuracy (such as fluctuation source positioning error) of the acquired data is calculated, and the error is ensured to be less than or equal to ±0.5° (preset accuracy threshold), and the actual resource consumption is also ensured to be less than 5500 mWh / day, so as to realize the dynamic balance of the monitoring system accuracy and energy consumption.

[0118] Please refer to Figure 3 , Figure 3 For the comprehensive comparison of the core performance indicators, the results show that the method is superior to the traditional technology in the three core performance indicators of positioning accuracy, azimuth accuracy and energy consumption efficiency, which shows that the method significantly improves the accuracy of the fluctuation source positioning and the system energy consumption control ability.

[0119] Please refer to Figure 4 , Figure 4For data acquisition efficiency comparison, the results show that the effective data amount of the method is 2.9M, and under the premise of ensuring the effective data amount, the redundant data is reduced by (8.5-4.2) / 8.5*100%≈50.6% compared with the traditional technology, which shows that the sampling strategy of the method can efficiently filter invalid data and improve the effectiveness of data acquisition.

[0120] Please refer to Figure 5 , Figure 5 For system endurance comparison, the results show that the system endurance time of the method is 7.1 days, which is improved by (7.1-4.0) / 4.0*100%=77.5% compared with the traditional technology, which shows that the energy consumption optimization transmission strategy of the method significantly prolongs the continuous working time of the system.

[0121] Please refer to Figure 6 , Figure 6 For system comprehensive performance score, the results show that the comprehensive performance score of the method is 86.6 points, which is improved by 27.2 points compared with the traditional technology, which shows that the method realizes all-round breakthrough in the overall performance of karst aquifer monitoring, and the technical advantage is obvious.

[0122] Please refer to Figure 7 , the karst aquifer monitoring system based on the Internet of Things in the embodiments of the present application is described, which comprises:

[0123] The signal acquisition module is configured to acquire medium fluctuation signals from a plurality of sensor nodes through an Internet of Things communication network, extract the propagation difference characteristics of the medium fluctuation signals in different medium attribute regions, and determine the medium attribute region type distribution according to the comparison results of the propagation difference characteristics and the preset region threshold.

[0124] The sampling strategy module is configured to generate a sampling frequency adjustment scheme according to the medium attribute region type distribution, adopt a high-frequency sampling strategy for high-permeability medium attribute regions, and adopt a low-frequency sampling strategy for low-permeability medium attribute regions.

[0125] The orientation estimation module is configured to dynamically configure sensor node parameters according to the sampling frequency adjustment scheme, acquire medium fluctuation signals using the configured sensor node parameters, calculate signal time delay based on time series analysis, and determine the fluctuation source azimuth angle estimation value.

[0126] The orientation correction module is configured to fuse the data of adjacent sensor nodes in the communication network according to the fluctuation source azimuth angle estimation value, correct the azimuth angle estimation deviation by an optimization algorithm, and obtain the accurate fluctuation source azimuth angle.

[0127] An energy consumption optimization module is configured to evaluate a signal loss risk based on the accurate wave source azimuth and the medium attribute region type distribution, and generate an energy consumption optimization transmission scheme according to a loss risk evaluation result.

[0128] A feedback control module is configured to extract a resource consumption index from the energy consumption optimization transmission scheme, and feed back the resource consumption index to the sampling frequency adjustment scheme through a feedback loop mechanism to realize a dynamic balance between a monitoring system precision and energy consumption.

[0129] The above embodiments are only used to illustrate the technical solutions of the present application, but not limit the same; although the present application has been described in detail with reference to the foregoing embodiments, those skilled in the art should understand that the technical solutions recorded in the foregoing embodiments can be modified, or some technical features can be replaced by equivalents; and these modifications or replacements do not make the essence of the corresponding technical solutions deviate from the spirit and scope of the technical solutions of the embodiments of the present application.

Claims

1. A karst aquifer monitoring method based on the Internet of Things, characterized in that, The method comprises the following steps: Step S101, collecting medium fluctuation signals from a plurality of sensor nodes through an Internet of Things communication network, extracting propagation difference characteristics of the medium fluctuation signals in different medium attribute regions, and determining a medium attribute region type distribution according to a comparison result of the propagation difference characteristics and a preset region threshold value; Step S102, generating a sampling frequency adjustment scheme according to the medium attribute region type distribution, using a high-frequency sampling strategy for a high-permeability medium attribute region, and using a low-frequency sampling strategy for a low-permeability medium attribute region; Step S103, dynamically configuring parameters of each sensor node according to the sampling frequency adjustment scheme, collecting medium fluctuation signals by using the configured sensor node parameters, and determining a fluctuation source azimuth angle estimation value based on time sequence analysis and signal time delay calculation; Step S104, fusing data of adjacent sensor nodes in the communication network according to the fluctuation source azimuth angle estimation value, correcting azimuth angle estimation deviation by using an optimization algorithm, and obtaining an accurate fluctuation source azimuth angle; Step S105, evaluating a signal loss risk based on the accurate fluctuation source azimuth angle and the medium attribute region type distribution, and generating an energy consumption optimization transmission scheme according to a loss risk evaluation result; Step S106, extracting a resource consumption index from the energy consumption optimization transmission scheme, and feeding back the resource consumption index to the sampling frequency adjustment scheme through a feedback loop mechanism to realize dynamic balance between monitoring system accuracy and energy consumption.

2. The method of claim 1, wherein, The extraction of the propagation difference characteristics of the medium fluctuation signals in different medium attribute regions comprises: Collecting medium fluctuation signal data from a plurality of sensor nodes deployed in a karst aquifer monitoring region through an Internet of Things communication network, and obtaining a first signal data set after preprocessing the data; Extracting frequency domain characteristics of the first signal data set by using a fast Fourier transform algorithm, and obtaining a second signal data set; Determining whether the frequency domain characteristics of the second signal data set satisfy a preset signal-to-noise ratio threshold value, and if so, decomposing the second signal data set by using a wavelet transform algorithm to obtain a propagation difference characteristic set; According to the propagation difference characteristic set, calculating a propagation time difference of medium fluctuation signals between each sensor node to obtain a propagation speed distribution; Matching the propagation speed distribution with a preset medium region division model to determine the boundaries of each medium region and obtain a medium attribute region division result; Applying a preset attenuation model to analyze the propagation difference characteristic set to calculate attenuation coefficients of medium fluctuation signals in each divided medium attribute region and obtain an attenuation degree distribution; Based on the propagation speed distribution and the attenuation degree distribution, generating propagation difference characteristics of medium fluctuation signals in different medium attribute regions of the karst aquifer.

3. The method of claim 2, wherein, The determination of the medium attribute region type distribution according to the comparison result of the propagation difference characteristics and the preset region threshold value comprises: Loading a preset region threshold value, wherein the region threshold value comprises a preset propagation speed threshold value and a preset attenuation coefficient threshold value; judging whether the propagation speed in the propagation difference characteristic of any medium attribute region exceeds the preset propagation speed threshold value, if yes, marking the corresponding region as a high-permeability medium attribute region, and if not, marking the corresponding region as a low-permeability medium attribute region; integrating the marking results of all regions and the propagation difference characteristic, associating the spatial coordinate information of the corresponding region, and generating a medium attribute region type distribution.

4. The method of claim 3, wherein, The step S102 comprises: For the medium attribute region type distribution, a hierarchical analysis method is used to classify the medium attribute region type to obtain a classified region set; Based on the classified region set, a support vector machine algorithm is used to distinguish the boundaries of the high-permeability medium attribute region and the low-permeability medium attribute region to determine the specific range of each region; For the specific range of each region, the historical data of the medium wave signal in each region is used to calculate the change rate of the medium wave signal in each region to obtain a signal change rate distribution; judging whether the signal change rate of a region in the signal change rate distribution is greater than a preset rate threshold value; If yes, a high-frequency sampling strategy is assigned to the high-permeability medium attribute region to determine a high-frequency sampling parameter set, and if not, a low-frequency sampling strategy is assigned to the low-permeability medium attribute region to determine a low-frequency sampling parameter set, wherein the sampling parameter set includes sampling frequency, sampling period, and each collection duration; According to the high-frequency sampling parameter set and the low-frequency sampling parameter set, a sampling frequency initial scheme corresponding to each region is generated; The sampling frequencies of adjacent regions in the sampling frequency initial scheme are smoothed by using a linear interpolation method to generate a sampling frequency adjustment scheme.

5. The method of claim 1, wherein, The step S103 comprises: Extracting the dynamic parameter configuration corresponding to each sensor node from the sampling frequency adjustment scheme, wherein the dynamic parameter configuration includes sampling frequency, signal gain, data acquisition period, and each collection duration; encoding the dynamic parameter configuration by using a data compression algorithm to obtain a compressed parameter set, and then transmitting the compressed parameter set to each sensor node in a broadcast manner; Each sensor node decodes the compressed parameter set, updates its own operating parameters according to the decoding result, and obtains a parameter synchronization state; According to the parameter synchronization state, each sensor node is activated, and each sensor node collects the medium wave signal of the karst aquifer according to the updated parameters to generate a multi-channel medium wave signal data set; Based on the multi-channel medium wave signal data set, a time delay distribution set is obtained by using a time sequence analysis method to calculate the time delay between the signal data collected by any two sensor nodes; judging whether each time delay value in the time delay distribution set is greater than a preset delay threshold value, if yes, using a triangular positioning algorithm to calculate the azimuth of the wave source to obtain a preliminary estimated value of the azimuth of the wave source; The preliminary estimated value of the azimuth of the wave source is smoothed by using a Kalman filtering algorithm to determine the estimated value of the azimuth of the wave source.

6. The method of claim 1, wherein, The step S104 comprises: Collecting multi-channel medium fluctuation signal data of a monitoring area corresponding to the fluctuation source azimuth angle estimation value, removing noise from the multi-channel medium fluctuation signal data, and extracting effective time sequence signal data to obtain a time sequence signal data set; Using a fast Fourier transform method to decompose the frequency components of the time sequence signal data set to obtain a frequency distribution set; If there is a frequency peak value in the frequency distribution set that meets a preset significance condition, then the phase difference between each channel is calculated by a frequency domain analysis method to obtain a phase difference set; Using a least squares method to fit the correspondence between the phase difference set and the sensor node position to obtain a preliminary fluctuation source position estimation; Using a weighted average algorithm to correct the deviation of the preliminary fluctuation source position estimation to obtain accurate position coordinates of the fluctuation source; According to the accurate position coordinates of the fluctuation source, the relative angle between each sensor node and the fluctuation source is calculated to generate an optimized azimuth angle set; Using a Kalman filtering algorithm to smooth the optimized azimuth angle set to obtain an accurate fluctuation source azimuth angle.

7. The method of claim 1, wherein, The step S105 comprises: Collecting medium fluctuation signal data of a plurality of sensor nodes in a monitoring range corresponding to the accurate fluctuation source azimuth angle, and extracting medium fluctuation signal propagation time delay using a time domain analysis method to obtain a time delay distribution set; According to the time delay distribution set, combining the medium attribute parameters of each region, using a path loss model to calculate the attenuation degree of the medium fluctuation signal in each region to obtain a propagation loss distribution; Based on the propagation loss distribution and the time delay distribution set, using a Monte Carlo method to simulate the propagation path of the medium fluctuation signal to obtain a loss risk assessment value; Determining whether the loss risk assessment value is higher than a preset risk threshold value; If yes, adjusting the transmission protocol parameters of the communication network, using a dynamic scheduling method to improve the data transmission priority of the corresponding sensor nodes to determine a high-priority transmission scheme; If no, analyzing the transmission load of the current communication network, using a fixed scheduling method to maintain the standard transmission priority of each sensor node to determine a standard transmission scheme; Using a linear programming method to adjust the transmission power of the sensor nodes to optimize the energy consumption configuration of the high-priority transmission scheme and the standard transmission scheme to obtain an optimized energy consumption configuration; According to the optimized energy consumption configuration, generating scheduling instructions for data transmission between sensor nodes, integrating the scheduling instructions with the corresponding transmission scheme to generate an energy-optimized transmission scheme.

8. The method of claim 1, wherein, The step S106 comprises: Extracting a resource consumption index from the energy-optimized transmission scheme, the resource consumption index including transmission power, total data transmission energy consumption, single transmission energy consumption, and data transmission time length; Based on the resource consumption index, using time series analysis to analyze the change trend of resource consumption of each sensor node, and simultaneously combining a sampling frequency parameter to obtain the correlation between the sampling frequency and the resource consumption; According to the correlation, using a linear regression algorithm to construct a resource consumption prediction model to output a resource consumption prediction value for a future preset time period under the current sampling frequency; Determining whether the resource consumption prediction value is higher than a preset energy consumption threshold value; If yes, a high-frequency sampling optimization sub-scheme is generated based on the correlation relationship for dynamically adjusting the sampling frequency parameter; If no, a standard sampling maintenance sub-scheme is generated by maintaining the current sampling frequency parameter; For the high-frequency sampling optimization sub-scheme or the standard sampling maintenance sub-scheme, a least square method is used to optimize the data collection power of each sensor node to determine a power optimization configuration; The power optimization configuration and the sampling frequency parameter obtained after judgment are fed back to the sampling frequency adjustment scheme to update the parameter configuration of the sampling frequency adjustment scheme; Based on the updated sampling frequency adjustment scheme, the medium fluctuation signal of the karst aquifer is re-collected and verified to ensure that the monitoring accuracy meets the preset accuracy threshold and the resource consumption is lower than the preset energy consumption threshold, thereby achieving dynamic balance between the monitoring system accuracy and energy consumption.

9. A karst aquifer monitoring system based on Internet of Things, used to implement the karst aquifer monitoring method based on Internet of Things according to any one of claims 1 to 8, characterized in that, The karst aquifer monitoring system based on the Internet of Things comprises: A signal collection module is configured to collect medium fluctuation signals from multiple sensor nodes through an Internet of Things communication network, extract the propagation difference characteristics of the medium fluctuation signals in different medium attribute regions, and determine the medium attribute region type distribution according to the comparison result of the propagation difference characteristics and a preset region threshold. A sampling strategy module is configured to generate a sampling frequency adjustment scheme by adopting a high-frequency sampling strategy for high-permeability medium attribute regions and a low-frequency sampling strategy for low-permeability medium attribute regions according to the medium attribute region type distribution. An azimuth estimation module is configured to dynamically configure parameters of each sensor node according to the sampling frequency adjustment scheme, collect medium fluctuation signals using the configured sensor node parameters, and determine a fluctuation source azimuth angle estimation value based on time series analysis and signal time delay calculation. An azimuth correction module is configured to fuse data of adjacent sensor nodes in the communication network according to the fluctuation source azimuth angle estimation value, correct the azimuth angle estimation deviation through an optimization algorithm, and obtain an accurate fluctuation source azimuth angle. An energy consumption optimization module is configured to evaluate signal loss risk based on the accurate fluctuation source azimuth angle and the medium attribute region type distribution, and generate an energy consumption optimization transmission scheme according to the loss risk evaluation result. A feedback control module is configured to extract a resource consumption index from the energy consumption optimization transmission scheme, feed back the resource consumption index to the sampling frequency adjustment scheme through a feedback loop mechanism, and achieve dynamic balance between the monitoring system accuracy and energy consumption.

Images