An intelligent sensor integrated hydrological water gauge system and method
Patent Information
- Application Number
- CN202611282051.7
- Authority / Receiving Office
- CN · China
- Patent Type
- Applications(China)
- Current Assignee / Owner
- Filing Date
- 2026-08-24
- Publication Date
- 2026-10-23
AI Technical Summary
[0007]本发明要解决的技术问题是现有水文水尺在汛期强光、无光、高含沙量与水面破碎条件下水位判读易失效,静压测量因含沙量与流速变化产生随工况改变的系统偏差,多通道按固定权重融合且以含被检通道自身的融合值构造判据使高权重通道的偏置难以检出,预警以单一阈值触发而虚警偏多
[0020]本发明设计了相互正交的双极性驱动码按时隙控制各发光刻度单元的近红外发光二极管亮灭,并对各区域的平均灰度与其自身驱动码作相关解调,未被该驱动码调制的环境光分量与相邻发光刻度单元的漏光分量在相关求和中相互抵消,水位判读所依赖的光信号被从环境照度中分离出来。本发明以水面以上各发光刻度单元调制幅度的中位数为参考值对各调制幅度作归一化,镜头表面附着物与空气中雨雾造成的等比例衰减在比值中被约去,判读阈值按已淹没发光刻度单元的调制深度中位数自动确定,判读结果不随水体透过率与镜头透过率的变化而改变。水位判读摆脱了对刻度与背景之间灰度对比度的依赖,白天强光、夜间无光与高含沙量三种条件下的判读条件趋于一致,长期运行中的镜头污染在较长时期内不需要人工清洗即可维持判读有效。相关解调在给出调制幅度的同时给出各调制深度的离散程度,光学通道在输出水位的同时输出与当前成像质量相对应的观测方差,后续的加权融合据此获得了随工况变化的权重依据。
Smart Images

Figure CN122961521A_ABST
Abstract
Description
Technical Field
[0001] This invention belongs to the field of liquid level measurement technology in hydrological monitoring, and particularly relates to an intelligent sensor integrated hydrological water gauge system and method. Background Technology
[0002] River water levels are fundamental observational data for flood control scheduling and water resource management. The process of water level changes and the rate of change of water level over time are direct indicators of whether a flood has formed and whether an early warning needs to be issued. Currently, traditional water level observation methods mainly include manual water gauge readings, electronic water gauges, float-type water level gauges, hydrostatic water level gauges, and radar and ultrasonic non-contact water level gauges. Manual water gauges rely on observers reading the scale; electronic water gauges use electrodes arranged along the elevation to determine the water surface location; float-type water level gauges rely on a float in the well that rises and falls with the water surface; hydrostatic water level gauges use pressure measured by an underwater pressure sensor to calculate water depth; and non-contact water level gauges use electromagnetic or sound waves emitted towards the water surface and the round-trip time of the echo to calculate distance. In recent years, video water gauges combining camera units and image processing have been used to replace manual visual readings. Their operation involves imaging the water gauge scale and identifying the boundary between the water surface and the scale. To improve the reliability of observations, hydrological stations often deploy two or more water level measuring devices simultaneously on the same cross section, and output the average or weighted average of the measurements from each device.
[0003] Traditional observation methods reveal several technical problems that need to be addressed during flood season. Video water level gauges rely on passive imaging and the grayscale contrast between the scale and the background. Nighttime supplemental lighting produces specular reflections on the water surface, while low-angle sunlight during the day causes glare. Increased sediment content during flood season reduces the contrast of the scale below the water surface. Water film, mud, and mist adhering to the lens surface cause overall image attenuation, making water level readings most prone to failure during the periods when readings are most needed. Static pressure level gauges typically use a constant based on the density of clear water and assume zero flow velocity at the pressure gauge when converting pressure to water depth. Increased sediment content during flood season causes the actual water density to deviate from the density of clear water, and increased river flow velocity causes pressure measurement errors at the pressure gauge. Both of these factors change with operating conditions, and the resulting deviations cannot be eliminated by fixed zero-point correction. The signs of the two deviations may also be opposite and cancel each other out under certain operating conditions. Relying solely on the comparison results under that operating condition can easily lead to erroneous judgments that do not require compensation. Non-contact water level gauges experience increased echo fluctuations during flood season when the water surface breaks up, churns, and carries foam and floating objects, leading to jumps in measurement values or even loss of tracking. Existing devices typically only output the measurement value itself and do not output real-time quality indicators that can be used for weighting.
[0004] When multiple water level measuring devices are deployed simultaneously, existing devices typically use fixed weights for averaging or weighted averaging. Even when a channel experiences short-term degradation, its error still enters the output result at a fixed proportion. The rate of change of water level with respect to time, as the derivative of water level, is more sensitive to observation noise, and its estimation is prone to distortion during periods of channel degradation in the flood season. Existing devices often construct a criterion for determining whether a channel is abnormal based on the difference between the channel's measurement and the fused value including itself. The channel under test appears simultaneously in both the minuend and subtrahend, compressing its residuals. Channels with higher weights experience stronger compression, making the most reliable measurement the most difficult to detect if it shows a bias. Furthermore, the sampling rates, exposure times, and communication delays of each sensor are different. Without a unified time reference and alignment rules, the actual times corresponding to each channel are not consistent. This time mismatch enters the result as water level deviation during rapid rises and falls in water level.
[0005] Existing flood warnings are mostly triggered by a single water level threshold or a fixed rate of rise threshold. Sensor malfunctions, piezometer blockages, and obstructions from floating debris can all cause water level readings to rise. Threshold-based criteria cannot distinguish these false rises from actual flood events, resulting in a high false alarm rate. Equipment alarms are often mixed with hydrological criteria in the same logic, making equipment malfunctions easily masked when hydrological elements change synchronously. Ice jams and ice dams causing ice flooding and flood overflow risks differ in their causes and treatment measures, and existing devices lack criteria to distinguish between the two. Furthermore, recovery actions after a channel is deemed abnormal often lack corresponding execution structures and recovery confirmation steps, requiring manual on-site handling at monitoring stations and resulting in short maintenance-free periods.
[0006] Therefore, there is a need for a hydrological gauge system and method that can provide reliable water levels and water level rise rates under complex flood conditions, and can make low false alarm level warnings based on these values. Summary of the Invention
[0007] The technical problems to be solved by the present invention are: the water level reading of existing hydrological gauges is prone to failure under conditions of strong light, no light, high sediment content and broken water surface during the flood season; static pressure measurement has a system deviation that changes with the operating conditions due to changes in sediment content and flow velocity; multi-channel fusion with fixed weights and the construction of criteria with the fusion value of the channel being tested makes it difficult to detect the bias of high-weight channels; and early warning is triggered by a single threshold, resulting in more false alarms.
[0008] This invention provides an operation method for an intelligent sensor-integrated hydrological gauge system, comprising the following steps: Step S1: Multiple luminous scale units and corresponding camera units are evenly spaced on the water-facing side of the gauge column. Pressure measurement holes are opened in the low-disturbance zone of the column sidewall and connected to pressure transmitters. Layered flow velocity units, water temperature chains, and turbidity units are installed along the elevation. A radar unit is installed at the top of the column. Data from each unit is timestamped and aligned according to a common clock. Step S2: The near-infrared LEDs of each luminous scale unit are controlled to turn on and off according to time slots using mutually orthogonal bipolar drive codes. The camera unit synchronously acquires multiple frames of images. The average grayscale of the region of interest of each unit is correlated and demodulated with its drive code to obtain the modulation amplitude. Step S3: The modulation depth is obtained by normalizing the median modulation amplitude of each unit above the water surface. The reading threshold is determined according to the median modulation depth of the submerged units. The optical channel water level and its observation variance are interpolated between adjacent units that cross the threshold. Step S4: Calculate sediment concentration from turbidity and turbidity density from water temperature. Use this density to convert the gauge pressure of the piezometer into pressure head, and correct the dynamic pressure deficit using the flow velocity at that elevation to obtain the static pressure channel water level and its observation variance. Step S5: Obtain the radar channel water level and its observation variance from the radar unit. For each channel, use the fused values of the other channels for a leave-one check and isolate abnormal channels. For healthy channels, use inverse variance weighting to obtain the fused water level and fused variance. Step S6: Using the fused water level as the observation and the fused variance as the time-varying observation noise variance, recursively filter the model with water level and rise rate as state variables, outputting the filtered water level and rise rate. Step S7: Calculate the cross-sectional flow rate from the stratified flow velocity and the stratified cross-sectional area varying with the filtered water level. Calculate the correlation coefficient between the filtered water level and the flow rate within a sliding window. Step S8: When the filtered water level is lower than the warning water level and the rise rate is greater than the starting value, calculate the alarm arrival time. When the correlation coefficient is not lower than the threshold, output flood warnings according to the alarm arrival time. When the coefficient is negative, the rise rate is positive, and the temperature of the submerged water temperature measuring point is not higher than the set value, output ice jam warnings. When the channel is continuously isolated, output equipment alarms.
[0009] Preferably, step S1 is further divided into steps S11 to S14. Each light-emitting scale unit includes one near-infrared and one white light-emitting diode that share a light guide plate and the center elevation is determined according to the segment height. The camera unit is equipped with a narrow-band filter and is triggered synchronously by the frame. The pressure measuring hole is located between the two flow guides.
[0010] Preferably, step S2 is further subdivided into steps S21 to S25. The driving code is taken from each row of the Walsh-Hadamard matrix except for the row with all values of +1. During the calibration stage, each region of interest is determined and the elevation correction normal is calibrated. After demodulation, crosstalk anomalous units are removed according to cross-correlation.
[0011] Preferably, step S3 is further divided into steps S31 to S34. The threshold is the average of the theoretical value of the modulation depth above the water surface and the median of the modulation depth of the submerged unit. The interpolation result is superimposed with the elevation-corrected normal to obtain the water level of the optical channel. Its observation variance is given by the segment height, the degree of dispersion of the modulation depth and the slope at the crossing point.
[0012] Preferably, step S4 is further divided into steps S41 to S44, where the sediment content is linearly converted from turbidity, the turbid water density is converted according to the two-phase mixing relationship, and the static pressure channel water level is calculated by subtracting the dynamic pressure correction term, which is expressed by a signed pressure coefficient and the flow velocity at that elevation, from the pressure head.
[0013] Preferably, step S6 is further divided into steps S61 to S64. The inspected channel does not participate in the formation of its reference value. The standardized residual uses the square root of the sum of the variance of the inspected channel and the left-one fusion variance as the denominator. Isolation and acceptance adopt hysteresis thresholds and zero-point deviation is corrected by exponential smoothing deviation state constrained by the upper limit.
[0014] Preferably, step S7 is further divided into steps S71 to S73, using a constant velocity model and continuous white noise acceleration process noise recursively according to Kalman filtering, and only making predictions when there is no fused water level.
[0015] Preferably, step S8 is further subdivided into steps S81 to S87. The cross-sectional flow rate is summed according to the conversion coefficient of the flow velocity of each layer and the product of the flow velocity and the cross-sectional area of the layer. The warning is divided into four levels: blue, yellow, orange and red, according to the time of arrival and is subject to correlation coefficient gating. The ice jam judgment criteria are composed of negative correlation, positive expansion rate and low temperature.
[0016] The present invention also provides a system for implementing the above method, including a water gauge column, an array of luminous scale units arranged at equal intervals along the water-facing surface, an encoding drive module, a camera unit facing the water-facing surface, an optical calculation unit, a static pressure measurement unit, a turbidity unit, a water temperature chain, a density dynamic pressure compensation unit, a stratified flow velocity unit, a radar unit located on the top of the column, a fusion arbitration unit, a flow calculation unit, a rise rate filtering unit, an early warning classification unit, a communication unit, and an edge computing master control unit with a built-in common clock; the encoding drive module is connected to each drive loop and triggers the camera unit with a frame synchronization output terminal; the camera unit is connected to the optical calculation unit; the static pressure measurement unit, turbidity unit, water temperature chain, and stratified flow velocity unit are connected to the density dynamic pressure compensation unit via a data bus; the optical calculation unit, density dynamic pressure compensation unit, and radar unit are respectively connected to the fusion arbitration unit; the fusion arbitration unit is respectively connected to the rise rate filtering unit and the early warning classification unit; the stratified flow velocity unit and the rise rate filtering unit are respectively connected to the flow calculation unit; the flow calculation unit is connected to the early warning classification unit; and the early warning classification unit is connected to the communication unit.
[0017] Preferably, the system further includes a self-diagnostic calibration unit, a backwashing assembly, and a complementary power supply unit. The self-diagnostic calibration unit is connected to the health status output terminal of the fusion arbitration unit, and its output terminal is connected to the control terminals of the coding drive module, the backwashing assembly, and the radar unit, respectively. The backwashing assembly includes a flushing pump, a solenoid valve, a check valve, and a flushing pipeline. The complementary power supply unit includes a photovoltaic module and a lithium battery pack.
[0018] This invention enables water level interpretation to be independent of ambient illumination and grayscale contrast, eliminates the systemic bias of the static pressure channel as operating conditions change, automatically adjusts the weight of each channel according to the observation variance, and detects the bias of high-weight channels, thus providing graded early warnings with low false alarms during the flood season.
[0019] The beneficial effects of this invention are as follows:
[0020] This invention designs orthogonal bipolar driving codes to control the on / off state of near-infrared LEDs in each luminous scale unit according to time slots. It also performs correlation demodulation on the average grayscale of each region and its own driving code. The ambient light component not modulated by this driving code cancels out the leakage light component of adjacent luminous scale units in the correlation summation, thus separating the light signal upon which water level reading depends from the ambient illuminance. This invention normalizes each modulation amplitude using the median of the modulation amplitude of each luminous scale unit above the water surface as a reference value. The proportional attenuation caused by lens surface deposits and rain / fog in the air is canceled out in the ratio. The reading threshold is automatically determined based on the median modulation depth of the submerged luminous scale units, and the reading result does not change with variations in water transmittance and lens transmittance. Water level reading eliminates the dependence on grayscale contrast between the scale and the background. The reading conditions tend to be consistent under three conditions: strong daylight, no light at night, and high sand content. Lens contamination from long-term operation can maintain effective reading for a considerable period without manual cleaning. Correlation demodulation provides the modulation amplitude and the degree of dispersion of each modulation depth. The optical channel outputs the observation variance corresponding to the current imaging quality while outputting the water level. The subsequent weighted fusion obtains the weight basis that varies with the working conditions based on this.
[0021] This invention calculates sediment concentration using turbidity from the same water gauge and turbidity density using water temperature. Pressure at the piezometer is converted to pressure head based on actual water density, eliminating conversion errors caused by increased sediment concentration. The invention compensates for pressure tapping errors caused by flow around the flow using the flow velocity at the piezometer elevation and the signed pressure tapping coefficient. The compensation amount varies with flow velocity. The piezometer is positioned in the flow-unseparated region between two flow guides, ensuring the pressure tapping coefficient remains stable after initial calibration. The two compensations target two independent and oppositely signed systematic deviations, one varying with sediment concentration and the other with the square of the flow velocity. They cancel each other out under specific operating conditions, bringing the total deviation close to zero. Simultaneous application of both compensations reduces the deviation of the static pressure channel under high flow velocity and high sediment concentration conditions to the level of random noise. The accuracy of the static pressure channel becomes comparable to that of optical and radar channels, satisfying the premise for multi-channel fusion.
[0022] In this invention, each water level channel provides its observation variance in real time based on its current physical quality indicators. The observation variance of the optical channel is derived from the steepness of the modulation depth step and the dispersion of each modulation depth above the water surface. The observation variances of the hydrostatic and radar channels are derived from the variance of their respective detrended window samples. Channel degradation is reflected in the observation variance in real time. This invention performs a leave-one consistency test on each channel using the fusion values of the other channels. The channel under test does not participate in the formation of the reference value, and the residual has the same statistical meaning as its judgment threshold. The channel with the highest weight can also be detected if it has a bias. The situation where the criterion is diluted by itself when constructing the residual using the fusion value including the channel under test is avoided. This invention obtains the fused water level and fused variance by weighting the observation variance of each healthy channel. The weight of the degraded channel automatically decreases, and the slowly varying zero-point deviation of each channel is corrected exponentially using a deviation state constrained by the upper limit. The correction amount does not increase cumulatively with the number of updates. This invention uses the fusion variance as the time-varying observation noise variance for recursive filtering. The water level rise rate, as the derivative of water level with respect to time, is more sensitive to observation noise. The filter automatically reduces the gain and relies more on model prediction during channel deterioration periods, and automatically increases the gain and tracks the actual rise and fall faster during channel recovery periods. The estimation of the water level rise rate remains usable during channel deterioration periods in the flood season.
[0023] This invention uses the arrival time of the alert as the grading metric for early warning, with the warning level directly corresponding to the available time margin for action. The same water level yields different levels under both gradual and rapid rises. This invention uses the correlation coefficient between the filtered water level and the cross-sectional flow rate within a sliding window as the confidence gate for early warning. In real flood events, the flood wave passing through the monitoring station causes both the water level and the flow rate to increase synchronously, and the two are highly positively correlated. False rises caused by pressure gauge blockage, floating debris obstruction, and channel zero-point drift are not accompanied by a synchronous increase in flow rate; such events are largely eliminated under the same hit rate conditions. This invention uses a negative correlation coefficient, a positive water level rise rate, and a temperature at the first water temperature measuring point below the water surface not exceeding a set value as criteria for ice blockage. Ice blockage reduces the cross-sectional area of water passage and causes the upstream water level to rise while the flow rate decreases instead of increasing. This is the opposite of the characteristic of both increasing in the same direction during a flood. Two types of dangerous situations with different causes and different treatment measures are distinguished on the same statistical measure. Equipment alarms are triggered by three factors: continuous abnormality in the health status of the channel, failure of sensor self-test, and failure to be re-accepted after the recovery action is executed. The correlation between the equipment and hydrological elements is decoupled, and equipment failure will not be masked when hydrological elements change synchronously. Attached Figure Description
[0024] Figure 1 The following are comparison diagrams of the optical channel environmental adaptability of Examples 1 to 3 and Comparative Examples 1 and 2. Among them, (a) is a semi-logarithmic curve of the relative error of modulation amplitude as a function of the DC component of ambient background light, (b) is a scatter plot and theoretical curve of the relative error of modulation amplitude as a function of code length, and (c) is a semi-logarithmic curve of the residual error of water level reading as a function of the lens attenuation factor.
[0025] Figure 2 The graph shows the comparison of water level deviation in the static pressure channel between Example 1 and Comparative Examples 3, 4, and 5. (a) is the calculated value of the deviation as a function of flow velocity and head, and the least squares fitted line graph. (b) is the calculated value of the deviation as a function of sediment concentration, and the least squares fitted line graph.
[0026] Figure 3 The following is a comparison chart of fault isolation and inflation rate estimation between Examples 1 to 3 and Comparative Examples 6, 7 and 9. Among them, (a) is a dual ordinate curve of the fault detection rate of Examples 1 and 9 as a function of the injection bias amplitude, (b) is a semi-logarithmic time history curve of the moving root mean square error of water level inflation rate, (c) is a grouped histogram of the root mean square error of water level inflation rate, and (d) is a semi-logarithmic histogram of the root mean square error of water level.
[0027] Figure 4 The diagram shows a comparison of the early warning performance and hazard identification of Examples 1, 4 and Comparative Example 8. (a) is a curve of the working characteristics of the early warning subjects, and (b) is a scatter plot of the discrimination plane for three types of working conditions: flood, ice blockage and equipment failure.
[0028] Figure 5 This is a schematic diagram of the composition structure of an intelligent sensor integrated hydrological water level gauge system according to the present invention. Detailed Implementation
[0029] The technical solutions of this invention will be clearly and completely described below. In addition, the forms of the various structures described in the following embodiments are merely illustrative. This invention is not limited to the structures described in the following embodiments. All other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of this invention.
[0030] An operation method for an intelligent sensor integrated hydrological gauge system includes steps S1-S8, wherein step S1 completes the arrangement of sensing components and time alignment of multi-channel data, including steps S11 to S14.
[0031] In step S11, N luminous scale units are arranged at equal intervals from bottom to top along the water-facing side of the water gauge column. Each unit consists of a light guide plate, a near-infrared LED with a center wavelength of 850nm, and a white LED. The two LEDs couple light into the light from the lower end face of the light guide plate. Each unit is independently driven by its own constant current branch. The center elevation of the k-th luminous scale unit is determined by the following formula:
[0032] ;
[0033] In the formula, The center elevation (m) of the kth illuminated scale unit; The zero point elevation of the water level gauge (m) is determined by fourth-order leveling from the basic leveling point; k is the serial number of the luminous scale unit, which is an integer from 1 to N from bottom to top; N is the total number of luminous scale units. The height (m) of the illuminated scale unit.
[0034] In step S12, a camera unit is installed on a bracket facing the water surface. Its photosensitive device does not have a near-infrared cutoff filter. A narrow-band filter with a center wavelength of 850nm is installed in front of the lens to filter out visible light components. The frame rate is determined according to the code length and the reading period. Its external trigger input terminal is connected to the frame synchronization output terminal of the encoding drive module. The image frame and the encoding time slot correspond one-to-one.
[0035] In step S13, a pressure measuring hole is made flush with the column surface on the side wall between the flow guides on both sides of the water gauge column. The axis of the hole is perpendicular to the side wall and orthogonal to the incoming flow direction. This area is the unseparated flow zone, and the pressure tapping coefficient is valid for a long time after the initial water tank calibration. The pressure measuring hole is connected to a gauge pressure transmitter with an atmospheric pressure reference tube via a pressure guide pipe.
[0036] In step S14, velocity sensors for each layer of the stratified velocity unit, measuring points for the water temperature chain, and turbidity units are arranged along the elevation. The velocity sensors are mounted on the extended supports on the side of the column, offset from the column itself. The water temperature chain covers the elevation of the pressure measuring hole and the range of water level fluctuations. A radar unit is installed on the top of the column when viewed from above. The sampling data from each unit is timestamped by a common clock within the edge computing main control unit and aligned to the same moment using zero-order hold or linear interpolation according to the fusion period. This avoids inconsistencies in the actual times corresponding to different channels, which could lead to water level deviations during rises and falls in water levels.
[0037] Step S2 completes the near-infrared orthogonal coding drive and correlation demodulation, including steps S21 to S25.
[0038] In step S21, take N rows from the L-order Walsh-Hadamard matrix except for the row where all elements are +1, and assign the k-th row to the k-th luminous scale unit as the driving code.
[0039] The Walsh-Hadamard matrix is a square matrix with elements of +1 or -1, and each row is pairwise orthogonal. An L-order Walsh-Hadamard matrix is a matrix with L rows and L columns, where the L row vectors form an orthogonal binary sequence of length L. The sum of the products of corresponding elements in any two rows equals L when the two rows are the same, and equals 0 when the two rows are different. When L is a power of 2, it is recursively generated using the Sylvester construction method, i.e., the top left, top right, and bottom left sections of the 2x matrix are taken from the original matrix, while the bottom right section is taken as the inverse of the original matrix. After removing the row with all elements +1, the sum of the elements in the remaining rows is 0. The code length is a power of 2 not less than the sum of N and 1, and is temporarily doubled by the self-diagnostic calibration unit when the optical channel is disturbed.
[0040] In step S22, within the nth time slot, the near-infrared LED of the kth luminous scale unit is lit when the drive code value is +1, and turned off when it is -1. n is the time slot number, an integer from 1 to L. During the image acquisition period, all white LEDs are turned off, and the exposure time and gain of the camera unit are locked, ensuring that the grayscale of each frame maintains the same linear relationship with the incident irradiance. Outside this period, the white LEDs are lit as needed for nighttime visual readings and inspections.
[0041] In step S23, during the installation and calibration phase, the near-infrared LEDs of each luminous scale unit are individually illuminated segment by segment. The region of interest corresponding to each unit is determined by the position of the bright area in the image and stored in the edge computing main control unit. During the same calibration process, the optical channel elevation constant is determined. This is done by first independently measuring the water surface elevation using a level or water level probe under conditions of a clean lens and stable water surface. Then, following steps S31 to S33, the elevation at which the modulation depth crosses the interpretation threshold is calculated. The difference between the former and the latter is taken as the constant and stored in the edge computing main control unit. The emitted light from the submerged unit is absorbed along the optical path. The position where the interpretation threshold is crossed is located at a certain depth below the actual water surface. This depth is determined by the water attenuation coefficient, the observation angle, and the interpretation threshold, and remains essentially unchanged after installation. This constant is repeatedly calibrated at different water levels, and the median is taken.
[0042] In step S24, the camera unit continuously acquires L frames of images according to the frame synchronization signal, calculates the average gray value of the k-th region of interest frame by frame, discards the current image and re-acquires when the frame number is missing or the gray value of the region is saturated, and then performs correlation demodulation according to the following formula to obtain the modulation amplitude:
[0043] ;
[0044] In the formula, L represents the modulation amplitude (grayscale level) of the k-th luminous scale unit; L is the code length (frame). This represents the average grayscale value (grayscale level) of the k-th region of interest in the n-th frame of the image. The value of the driving code for the k-th luminous scale unit in the n-th time slot; To sum the time slot numbers n from 1 to L; To obtain the absolute value, the grayscale of the region of interest is superimposed from the ambient background component, the illuminated component of this unit, and the light leakage component of adjacent units. The three components are multiplied by the driving code of this unit and summed over the time slot. The background component tends to 0 because the sum of the elements of the driving code is 0, the light leakage component is 0 because the driving codes are orthogonal, and only the component of this unit itself is accumulated because the square of the driving code value is always 1.
[0045] In step S25, the frame-by-frame average gray value of the kth region of interest is cross-correlated with the driving codes of other light-emitting scale units. The maximum value of the cross-correlation values is taken. If the ratio of the maximum value to the modulation amplitude of the unit is greater than the elimination determination ratio, the unit is determined to be a crosstalk abnormal unit. Its modulation amplitude is eliminated and linearly interpolated by the modulation amplitudes of the two adjacent non-eliminated units. The elimination determination ratio is taken from the upper envelope of the ratio of the cross-correlation value to the modulation amplitude of each unit during the installation and calibration stage.
[0046] Step S3 completes the calculation of the optical channel water level and observation variance, including steps S31 to S34.
[0047] In step S31, the light emission scale unit above the water surface is determined by adding the optical channel water level output from the previous fusion cycle to the two segment heights. For the first run or when there is no effective optical channel water level in the previous cycle, the radar channel water level or static pressure channel water level at the same moment is used instead. The median of the modulation amplitude of each unit above the water surface is taken as the reference value, and the modulation depth is obtained by normalization using the following formula:
[0048] ;
[0049] In the formula, The modulation depth (dimensionless) of the k-th luminous scale unit. The meaning is the same as in step S24; The reference modulation amplitude (grayscale level) is used. Lens attachments and rain / fog cause approximately proportional attenuation to each element, which is canceled out in the ratio.
[0050] In step S32, the median of the modulation depth of each unit whose serial number is more than two times lower than the serial number of the light-emitting scale unit above the water surface is taken, and the reading threshold is determined by the following formula:
[0051] ;
[0052] In the formula, The threshold for interpretation (dimensionless); The value is the median (dimensionless) of the modulation depth of all units whose serial numbers are more than two below the serial number of the luminous scale unit above the water surface. The 1 in the formula represents the theoretical value of the normalized modulation depth of the units above the water surface. This threshold is taken as the midpoint between the segment above and below the water surface and is automatically adjusted according to the water penetration conditions. The modulation depth of each unit above the water surface tends to be constant, while the modulation depth of each unit below the water surface decreases with increasing submersion depth and quickly stabilizes. A steep transition zone is formed near the water surface between the two. This is because water absorbs near-infrared radiation much more strongly than visible light, and the scattering of suspended particles further amplifies this attenuation when the sand content increases.
[0053] In step S33, a bottom-up search is performed for adjacent unit pairs where the modulation depth of the current unit is lower than the readout threshold and the modulation depth of the adjacent unit above it is not lower than the readout threshold. If neither unit is eliminated, the unit pair is considered valid. If there are more than two valid unit pairs, the smaller absolute value of the difference between its interpolation result and the fusion level of the previous fusion cycle is selected. The index of the lower unit in the selected unit pair is recorded as k. The elevation where the modulation depth crosses the readout threshold is obtained by linear interpolation using the following formula:
[0054] ;
[0055] In the formula, The elevation (m) at which the modulation depth crosses the interpretation threshold; The modulation depth (dimensionless) of the (k+1)th luminous scale unit. , , and The meaning is the same as in steps S11, S31, and S32. The two modulation depths are located on either side of the interpretation threshold to prevent extrapolation of the results; a multiple-solution selection rule is set because wave rise, surface foam, and splashing water will create additional drops beyond the actual water surface. Then, the optical channel water level is obtained by superimposing the elevation and converting it to a positive constant using the following formula:
[0056] ;
[0057] In the formula, The optical channel water level (m); The optical channel elevation is corrected to a constant (m), which is obtained by calibration in step S23 and remains unchanged during operation; The meaning is the same as the previous formula. This formula subtracts the inherent depth of the crossing point relative to the actual water surface, so that the three water level channels are located on the same elevation datum.
[0058] In step S34, the observation variance of the optical channel is given by the following formula:
[0059] ;
[0060] In the formula, The observation variance of the optical channel (m 2 ); The sample standard deviation (dimensionless) of the modulation depth of each unit above the water surface. , and The meaning is the same as in step S33. This formula is obtained from the interpolation formula in step S33 through first-order error propagation, with the elevation constant replaced by a normal and not included in the propagation. The observation variance changes in real time with the imaging quality rather than being a preset constant, and its value is not less than one-twelfth of the square of the segment height. If there is no effective cell pair or the difference between the two modulation depths in the effective cell pair is lower than the set lower limit, the optical channel level will not be output in this fusion cycle, and the channel will be treated as a missing measurement.
[0061] Step S4 completes the dual compensation of the density of the sediment-laden water and the pressure measurement error of the piezometer, including steps S41 to S44.
[0062] In step S41, the sediment content is calculated from the turbidity using the following formula:
[0063] ;
[0064] In the formula, Sand content (kg / m³) 3 ); a is the turbidity conversion slope (kg·m -3·NTU -1 U represents the turbidity measured in the turbidity unit (NTU); b represents the turbidity conversion intercept (kg / m³). 3 The slope and intercept are obtained by least-squares linear fitting of the sediment concentration measured by the drying and weighing method on water samples collected synchronously at this station with the synchronous turbidity value. A comparative measurement is performed before and after the flood season to update the values. If the turbidity exceeds the fitted range, the sediment concentration is invalidated and the most recent valid value is used.
[0065] In step S42, the density of turbid water is calculated from the sand content and water temperature using the following formula:
[0066] ;
[0067] In the formula, Density of turbid water (kg / m³) 3 ); Density of clean water (kg / m³) 3 The value is obtained by interpolating the pure water density table from the actual measured water temperature at the elevation of the pressure measuring hole using the water temperature chain. The density of sediment particles (kg / m³) 3 ), take the measured value of sediment particle density at this station; The meaning is the same as in step S41. This formula is obtained from the mass balance of a unit volume of turbid water, that is, the mass of sediment is added to the mass of the part of the volume occupied by water; the increase of sediment content increases the density of turbid water, and the water level is overestimated when converted according to the density of clear water.
[0068] In step S43, the interpolated velocity of the stratified velocity unit at the pressure measuring orifice elevation is taken, and the static pressure channel water level is obtained by the following formula:
[0069] ;
[0070] In the formula, The water level in the static pressure channel (m); The elevation of the pressure measuring hole (m); p is the gauge pressure output by the pressure transmitter (Pa); g is the acceleration due to gravity (m / s²). 2 ); The pressure tapping coefficient (dimensionless) is the pressure tapping coefficient of the pressure measuring hole, which is a signed coefficient that can be positive or negative; The water flow velocity (m / s) at the elevation of the pressure gauge hole. The meaning is the same as in step S42. The basis for the third item is the general definition of the pressure coefficient, which is the difference between the pressure at a certain point on the surface and the static pressure of the undisturbed incoming flow divided by the incoming flow pressure. The pressure coefficient is obtained by keeping the water surface elevation constant in the tank and changing the incoming flow velocity in stages, and by performing a least-squares linear fit on the velocity head using the difference between the static pressure channel water level without dynamic pressure correction and the actual water surface elevation, and taking its slope.
[0071] In step S44, take a length of The sample variance of the static pressure channel water level after de-linearizing within the sliding window is the observation variance of the static pressure channel. The observation variance is not less than one-twelfth of the square of the ratio of the pressure transmitter's pressure resolution to the product of the turbid water density and gravitational acceleration. Delinearization removes the fluctuations in water level itself.
[0072] Step S5 completes radar channel calculation, channel health determination and inverse variance fusion, including steps S61 to S64.
[0073] In step S61, the water level of the radar channel is obtained by subtracting the measured distance to the water surface from the installation elevation of the radar unit, and the length is taken as... The sample variance of the radar channel water level after detrending within the sliding window is the observation variance of the radar channel, which is not less than one-twelfth of the square of the radar element resolution. During the flood season, when the water surface breaks up, the echo fluctuations intensify, the detrended sample variance increases accordingly, and the channel weight decreases.
[0074] In step S62, for the i-th channel, perform leave-one fusion with the remaining channels as follows and calculate the standardized residual:
[0075] ;
[0076] ;
[0077] ;
[0078] In the formula, Reserve a fusion water level (m) for the i-th channel; Leave-one fusion variance (m) for the i-th channel 2 ); Let be the standardized residual (dimensionless) of the i-th channel; i and m are the channel numbers, which are taken as optical channel, hydrostatic channel and radar channel respectively; and Let be the water levels (m) of the i-th channel and the m-th channel, respectively. and The observation variances (m) of the i-th channel and the m-th channel are respectively. 2 ); This is to sum over all channels except the i-th channel. The weights of the first two equations are proportional to the reciprocal of their respective variances, assuming the combined variance is minimized. The denominator of the third equation is the square root of the sum of the variance of the tested channel and the variance of the leave-one fusion, since the tested channel does not participate in leave-one fusion and the two are independent. If the residual is constructed from the difference between the tested channel and the fusion value containing itself, this difference is equal to 1 minus the channel's weight multiplied by the difference between it and the leave-one fusion value. The larger the weight, the stronger the compression. The channel with the smallest variance is the most difficult to detect if it shows bias.
[0079] In step S63, when the standardized residual of the i-th channel continuously exceeds the isolation threshold for a set number of cycles, its health status is set to abnormal and it is isolated. When the standardized residual of the isolated channel continuously falls below the acceptance threshold for a set number of cycles, its health status is restored. The acceptance threshold is lower than the isolation threshold, and the two constitute hysteresis. When the number of healthy channels is not less than 2, the fused water level is obtained by weighting the reciprocals of their observed variances, and the fused variance is taken as the reciprocal of the sum of the reciprocals of the observed variances of each healthy channel. When the number of healthy channels is 1, the water level of that channel is used as the fused water level, its observed variance is used as the fused variance, and a downgrade flag is set. When the number of healthy channels is 0, the fused water level is not output in this fusion cycle, and the data does not enter steps S6 to S8.
[0080] In step S64, during the period when the absolute value of the water level rise rate is lower than the set value and the number of healthy channels is not less than 2, the deviation status of the i-th channel is updated according to the following formula.
[0081] ;
[0082] In the formula, and Let m be the deviation state of the i-th channel after the current update and the previous update, respectively. The update coefficient (dimensionless) takes a value greater than 0 and less than 1; To obtain the median of all samples within the window; and The meaning is the same as in step S62. This formula is a first-order exponential smoothing, where the deviation state gradually approaches the median of the deviation within the window without accumulating and increasing with the number of updates. The update coefficient is set to a small value to only track slow-changing processes such as pressure gauge hole siltation, sensor zero-point drift, and slow column settlement. The absolute value of the deviation state is constrained by the upper limit of the deviation. When it exceeds the upper limit, the health status of the channel is set to abnormal. Thereafter, the sum of the channel water level and its deviation state is used in steps S62 and S63.
[0083] Step S6 completes the joint recursive estimation of water level and water level rise rate, including steps S71 to S73.
[0084] In step S71, a discrete state equation is established using water level and water level rise rate as state variables. The state transition matrix and observation matrix are respectively taken as follows:
[0085] ;
[0086] ;
[0087] In the formula, F is the state transition matrix; H is the observation matrix; The filtering period is denoted by s. This model is a constant velocity model, and the rise rate is indirectly estimated by a recursive process. The initial values of the state variables are taken as the first effective fused water level and zero, and the diagonal elements of the initial covariance matrix are taken as the first fused variance and the square of the water level rise rate range, respectively.
[0088] In step S72, the process noise covariance matrix is obtained using the following formula:
[0089] ;
[0090] In the formula, Here, is the process noise covariance matrix; q is the process noise intensity (m) driving the water level rise rate change. 2 / s 3 ); The meaning is the same as in step S71. This matrix is the discretization result of the continuous white noise acceleration model. Its derivation is to apply the process noise only to the rate of rise component and integrate it within one filtering period. The process noise intensity is taken from the continuous water level sequence during the non-flood season of this station, and the value is selected so that the average value of the normalized squared innovation of the one-step prediction tends to 1.
[0091] In step S73, the Kalman filter is used to recursively calculate the observed water level and the observed noise variance, using the fused water level as the observation measurement and the fused variance as the observed noise variance. The prediction equation is: the state prior estimate equals the product of the state transition matrix and the state posterior estimate of the previous cycle; the prior covariance matrix equals the product of the state transition matrix, the previous cycle's posterior covariance matrix and its transpose, plus the process noise covariance matrix. The update equation is: the gain matrix equals the product of the prior covariance matrix and the transpose of the observation matrix divided by the observation matrix, the product of the prior covariance matrix and the transpose of the observation matrix, and the observed noise variance; the innovation equals the observation measurement minus the product of the observation matrix and the state prior estimate; the state posterior estimate equals the prior estimate plus the product of the gain matrix and the innovation; the posterior covariance matrix equals the identity matrix minus the product of the gain matrix and the observation matrix, multiplied by the prior covariance matrix. The filtered water level and water level rise rate are recursively output. The observation noise variance is not taken as a constant, but as the fusion variance given in step S63. The filter automatically reduces the gain during channel degradation periods and automatically increases the gain during recovery periods. If there is no fusion water level in this fusion cycle, only the prediction equation is executed. If the number of consecutive missing measurements exceeds the set number of cycles, the filter is reinitialized according to step S71.
[0092] Steps S7 and S8 complete the cross-sectional flow calculation, hydrological consistency judgment and graded early warning, including steps S81 to S87.
[0093] In step S81, the cross-sectional flow rate is calculated using the following formula:
[0094] ;
[0095] In the formula, Q is the cross-sectional flow rate (m³ / s). 3 / s); To sum the values of layer numbers j from 1 to J; is the velocity conversion factor (dimensionless) for the j-th layer. The velocity (m / s) measured by the velocity sensor at the j-th layer. The cross-sectional area (m²) of the j-th layer 2 J represents the number of velocity sensor layers below the water surface, determined periodically according to the filtered water level. The layered cross-sectional area is obtained by dividing the measured cross-section at this station according to the filtered water level periodically. The midpoint between the elevations of two adjacent velocity sensor layers is the layer interface, and the upper boundary of the top layer is taken as the filtered water level. The velocity conversion factor converts the point velocity at the column into the average cross-sectional velocity of the corresponding layer. Calibration is performed using the velocity area method or the synchronous measured flow rate of the acoustic Doppler current profiler as a reference. The calibration covers low, medium, and high water levels. When the cross-section is eroded or silted up, the layered area is updated according to the measured cross-section after the flood season.
[0096] In step S82, in a length of Within a sliding window, the Pearson correlation coefficient between the filtered water level series and the cross-sectional flow series is calculated. This coefficient is equal to the product of the covariance of the two series and the standard deviation of the two series. Subsequently, within a set time lag range, the time-lag cross-correlation coefficient between the filtered water level series and the sediment concentration series is calculated, and the larger absolute value and its corresponding time lag are taken as the auxiliary confidence level. The response of sediment concentration to water level has a lag, so it is only used as an auxiliary confidence level and not as a necessary condition for flood identification.
[0097] In step S83, when the filtered water level is lower than the warning water level and the water level rise rate is greater than the starting value, the alarm time is calculated according to the following formula:
[0098] ;
[0099] In the formula, Time to police arrival (s); This is the warning water level for this station (m); The filtered water level (m); The water level rise rate (m / s) is the alarm arrival time. The alarm arrival time represents the time required for the water level to reach the warning level at the current rise rate. The alarm arrival time is not output when the filtered water level is not lower than the warning level or the rise rate is not greater than the calculated value. The calculated value is slightly higher than the random fluctuation range of the rise rate estimate when the water level is stable.
[0100] In step S84, a rising water level is determined when the correlation coefficient between the filtered water level and the cross-sectional flow is not lower than the first threshold. A yellow warning is issued when the alarm arrival time is lower than the first time limit, an orange warning when it is lower than the second time limit, and a red warning when it is lower than the third time limit. A blue warning is issued in all other cases. A red warning is issued when the first time limit is greater than the second time limit and the second time limit is greater than the third time limit, and the filtered water level is not lower than the warning level. When the auxiliary confidence level is not lower than the second threshold, the message indicates an increased confidence level without changing the warning level. In a real flood, the water level and flow rate increase synchronously and are highly positively correlated. False rises caused by piezometer blockage, floating debris obstruction, and zero-point drift are not accompanied by an increase in flow rate. Each time limit is determined according to the preparation time required for the station's flood control response, and the two thresholds are selected based on the station's historical flood and false alarm event sequences.
[0101] In step S85, if the correlation coefficient between the filtered water level and the cross-sectional flow rate is not higher than the negative of the first threshold, the water level rise rate is greater than the ice jam starting value, and the temperature of the first water temperature measuring point below the filtered water level from top to bottom is not higher than the freezing judgment temperature, and all three conditions are met simultaneously and the window number is continuously set, it is determined to be an ice jam hazard and an ice jam warning is output. Taking the first water temperature measuring point below the filtered water level is to avoid mistaking air temperature for water temperature. Ice jam reduces the cross-sectional area of water flow while the upstream water level rises, and the flow rate decreases instead of increasing, which is the opposite of the characteristic of both increasing in the same direction during a flood; requiring the rise rate to be positive is to exclude the reverse change during the receding period.
[0102] In step S86, if the channel health status is abnormal and the set number of cycles continues, the sensor self-test fails, or the channel is not re-accepted after the recovery action is executed, the output device will issue an alarm and mark the corresponding channel. This output is unrelated to the correlation coefficient in step S82. The recovery action is executed according to the type of isolated channel. When an optical channel is isolated, the drive current of the near-infrared LED is increased and the code length is doubled before re-measurement. When a static pressure channel is isolated, the backwashing component is activated to remove the sediment in the pressure measuring hole and the pressure guiding tube, and the deviation status update in step S64 is executed again after the backwashing is completed. When a radar channel is isolated, the backup beam is switched and the measurement cycle is extended. After the recovery action is executed, the channel will be automatically re-accepted if the standardized residual meets the acceptance conditions; otherwise, maintenance will be reported.
[0103] In step S87, the warning level, hazard type, filtered water level, water level rise rate, alarm arrival time, and channel health status are encoded into a message and sent in parallel through each channel of the communication unit. If any channel receives a successful acknowledgment, the message is considered delivered. If no acknowledgment is received, the message is retried at a set interval.
[0104] This invention also designs an intelligent sensor integrated hydrological water gauge system, such as... Figure 5As shown, it includes a water gauge column, an array of luminous scale units, an encoding drive module, a camera unit, an optical calculation unit, a static pressure measurement unit, a turbidity unit, a water temperature chain, a density dynamic pressure compensation unit, a stratified flow velocity unit, a radar unit, a fusion arbitration unit, a flow calculation unit, a rate filtering unit, an early warning classification unit, a communication unit, an edge computing main control unit, a self-diagnostic calibration unit, a backwashing component, and a complementary power supply unit.
[0105] The water gauge column is vertically fixed to the bank or a dedicated foundation at the test section. The water-facing side is flat and has flow guides on both sides. The pressure measuring hole of the static pressure measuring unit is located on the side wall between the two flow guides. The luminous scale unit array is arranged at equal intervals from bottom to top along the water-facing side. Each unit contains one near-infrared light-emitting diode and one white light-emitting diode sharing a light guide plate. The controlled ends of the two are connected to the encoding output and illumination output of the encoding drive module, respectively. The encoding output of the encoding drive module is connected to each drive circuit, and its frame synchronization output is connected to the external trigger input of the camera unit. The illumination output is turned off during the image acquisition period.
[0106] The camera unit is mounted on a bracket facing the water-facing side, and its image output is connected to the optical calculation unit. The optical calculation unit stores the positions of each region of interest and the elevation correction constants of the optical channels. The static pressure measurement unit includes a pressure measuring hole, a pressure guide tube, and a gauge pressure transmitter with an atmospheric pressure reference tube. The turbidity unit and water temperature chain are located on the water-facing side. The velocity sensors of each layer of the stratified flow velocity unit are located on the extended support arm. All four of these are connected to the density dynamic pressure compensation unit via a data bus.
[0107] The radar unit, viewed from above, is located at the top of the column and is a radar sensor with two or more switchable beams. Its beam control terminal is connected to the self-diagnostic calibration unit. The output terminals of the optical calculation unit, density dynamic pressure compensation unit, and radar unit are respectively connected to the fusion arbitration unit. The water level output terminal of the fusion arbitration unit is connected to the rise rate filter unit, and its health status output terminal is connected to the early warning classification unit and the self-diagnostic calibration unit. The output terminals of the stratified flow velocity unit and the rise rate filter unit are respectively connected to the flow calculation unit. The output terminals of the flow calculation unit, turbidity unit, water temperature chain, and rise rate filter unit are respectively connected to the early warning classification unit. The output terminal of the early warning classification unit is connected to the communication unit, which includes three types of channels: narrowband spread spectrum, cellular mobile communication, and satellite short message.
[0108] The optical calculation unit, density dynamic pressure compensation unit, fusion arbitration unit, expansion rate filtering unit, and early warning classification unit are integrated into the edge computing main control unit. The main control unit has a built-in common clock, and the sampling data from each unit is timestamped by the common clock and transmitted via the data bus. The output of the self-diagnostic calibration unit is connected to the control terminals of the encoding drive module, backflushing assembly, and radar unit, respectively, outputting three types of commands according to the type of isolated channel: increasing drive current and doubling code length, initiating flushing, and switching beams. The backflushing assembly includes a flushing pump, solenoid valve, check valve, and flushing pipeline connected to the pressure measuring port. The flushing pump inlet forms a loop with the flushing pipeline via a storage tank. The complementary power supply unit includes photovoltaic modules, a charging control circuit, and a lithium battery pack, with its output connected to the power supply input terminals of each unit. All the above structures and units are customized according to the method of this invention.
[0109] Example 1: This example uses a river hydrological station with an approximately trapezoidal test section, representing a width of 45m. The annual water level variation is 4.5m, and the warning water level is taken as 4.50m relative to the zero point of the water gauge.
[0110] Arrange the sensing components according to step S1. According to step S11, arrange them at equal intervals from bottom to top along the water-facing side of the water gauge column. =250 illuminated scale units, segment height =0.020m, the array covers 5.0m upwards from the zero point of the water gauge, exceeding the historical water level fluctuation range and leaving a margin. Each luminous scale unit contains one near-infrared light-emitting diode with a center wavelength of 850nm and one white light-emitting diode, sharing a section of light guide plate. The center elevation of each unit is based on Confirm the zero point elevation of the water draft gauge. The measurement is performed using fourth-order leveling. Following step S12, a camera unit is installed directly facing the water-facing surface, with a frame rate of 50 frames / s. A narrow-band filter with a center wavelength of 850nm and a half-width of 30nm is placed in front of the lens, and its external trigger input is connected to the frame synchronization output of the encoding drive module. Following step S13, a pressure measuring hole is drilled flush with the cylindrical surface on the sidewall between the two guide fairings. The elevation of the pressure measuring hole is... =1.00m, a gauge pressure transmitter with an atmospheric pressure reference tube connected via a pressure guide pipe, measuring range 100kPa, pressure resolution 10Pa. Following step S14, arrange 6 layers of flow velocity sensors, 8 water temperature measuring points, and 1 turbidity unit along the elevation. A radar unit is installed on the top of the column (viewed from above), at an elevation of 8.00m, with a resolution of 1mm. Data from each unit is timestamped using a common clock and fused according to the fusion period. =1s alignment.
[0111] Perform encoding and demodulation in step S2. Obtain the code length in step S21. =256, the driving code is assigned to the first 250 rows of the 256th order Walsh-Hadamard matrix, excluding the row where all elements are +1. In step S22, the near-infrared LEDs are controlled to turn on and off in time slots, while white light is turned off and exposure time and gain are locked during image acquisition. In step S23, during the installation and calibration phase, the LEDs are illuminated segment by segment to determine each region of interest, and the elevation correction constant of the optical channel is simultaneously determined under the condition of a clean lens with known water surface elevation. In step S24, 256 frames are continuously acquired and correlation demodulation is performed to obtain the modulation amplitude. In step S25, crosstalk aberration units are eliminated using the cross-correlation value, with an elimination ratio of 0.30.
[0112] The optical channel is calculated in step S3. In step S31, the median of the modulation amplitude of the units above the water surface is taken as the reference value for normalization. In step S32, the interpretation threshold is determined. In step S33, the optical channel water level is obtained by interpolation between two adjacent units that cross the threshold. In step S34, the observation variance of the optical channel is given, and the lower limit of the variance is taken as one-twelfth of the square of the segment height.
[0113] Solve the static pressure channel according to step S4. Turbidity conversion slope. =0.85kg·m -3 ·NTU -1 Turbidity conversion intercept =−1.2kg / m 3 Both were obtained by fitting the results of sand content measured by the drying and weighing method with the results of turbidity measurement; silt particle density =2650kg / m 3 Pressure coefficient =−0.35, obtained from the water tank calibration; gravitational acceleration =9.807m / s 2 .
[0114] Merge according to step S5. Slide window length. =600s, Isolation threshold 3.0, Admission threshold 1.5, Number of consecutive cycles 5, Deviation state update coefficient =0.05, upper limit of deviation 0.05m.
[0115] Filter recursively according to step S6. Process noise intensity. =1×10 -11 m 2 / s 3 .
[0116] Calculate the flow rate and output the warning according to steps S7 and S8. The stratified cross-sectional area is divided by the measured cross-section according to the filtered water level, and the flow velocity conversion coefficient is obtained by the flow velocity area method. The starting value of the flood rate is 0.05 m / h, the starting value of the ice jam is 0.10 m / h, the first threshold is 0.60, the second threshold is 0.50, and the three time limits are 3h, 2h and 1h respectively.
[0117] Example 2: In this example, the segment height is 0.010m, the total number of units is 500, the code length is 512, the frame rate is 100 frames / s, the fusion period is 0.5s, the sliding window length is 300s, the deviation state update coefficient is 0.02, and the rest is the same as in Example 1.
[0118] Example 3: In this example, the segment height is 0.050m, the total number of units is 100, the code length is 128, the frame rate is 25 frames / s, the fusion period is 2s, the sliding window length is 1800s, the deviation state update coefficient is 0.10, and the rest is the same as in Example 1.
[0119] Example 4 uses an ice-season configuration. The freezing temperature is set to 0.5℃, the duration window for ice jam criterion is set to 3, and the rest is the same as in Example 1.
[0120] Comparative Example 1: Steps S2 and S25 were not used in this comparative example. Each light-emitting scale unit was always lit. The camera unit collected the same number of frames and then took the average. The intensity estimate of the unit was obtained by subtracting the average gray value of the neighboring background area from the average gray value of the region of interest of the unit. The rest was the same as in Example 1.
[0121] Comparative Example 2: This comparative example does not use the self-reference normalization in step S31 and the adaptive interpretation threshold in step S32. Instead, it directly uses the factory nominal brightness for normalization and takes a fixed interpretation threshold of 0.5. The rest is the same as in Example 1.
[0122] Comparative Example 3: This comparative example did not use the turbid water density conversion in step S42 or the pressure compensation in step S43. The static pressure channel was based on a clean water density of 1000 kg / m³. 3 The conversion is the same as in Example 1.
[0123] Comparative Example 4 retains the turbid water density conversion in step S42, but cancels the pressure compensation in step S43. The rest is the same as in Example 1.
[0124] Comparative Example 5: This comparative example does not use the turbid water density conversion in step S42, but retains the pressure compensation in step S43. The rest is the same as in Example 1.
[0125] Comparative Example 6: This comparative example does not use the leave-one-out fusion in step S62. The three channels are weighted arithmetic averages and anomaly detection is performed using the standardized residuals containing their own fusion values. The variance of the observation noise of the recursive filter is taken as a constant. The rest is the same as in Example 1.
[0126] Comparative Example 7: This comparative example does not use the real-time estimation of observation variance in steps S34, S44 and S61. The observation variance of each channel is taken as the factory nominal constant, the fusion weight is fixed, and the observation noise variance of the recursive filter is taken as a constant. The rest is the same as in Example 1.
[0127] Comparative Example 8: This comparative example does not use the consistency gating of step S84, but only triggers the warning with the water level threshold and the fixed rise rate threshold. The rest is the same as in Example 1.
[0128] Comparative Example 9 retains the inverse variance weighted fusion of step S62 but eliminates the leave-one-out process. That is, it constructs a standardized residual using the fusion value that includes the channel under test and makes anomaly detection based on it. The rest is the same as in Example 1. The only difference between this comparative example and Example 1 is the leave-one-out technical feature.
[0129] Experimental Example 1, using the methods of Examples 1, 2, and 3 and Comparative Examples 1 and 2, verifies the adaptability of the optical channel to ambient illumination, code length, and lens attenuation. The results are as follows: Figure 1 As shown.
[0130] The experimental procedure is as follows: following the relevant demodulation formula in step S24, the ambient background light DC component is demodulated 60 times within the range from gray level 1 to gray level 120. The deviation of the modulation amplitude estimate from the true value is used as the relative error. The upper limit of the background light is set so that the peak gray level of the bright frame is lower than the full scale of the photosensitive device to avoid saturation. Then, the above demodulation is repeated with code lengths of 128, 256, and 512 respectively. Finally, following steps S31 to S33, the reading is repeated 40 times under the condition that the lens attenuation factor is reduced from 1.00 to 0.40. The change in the difference between the reading water level and the actual water level relative to the elevation change normal value under the clean lens condition is used as the residual error.
[0131] Figure 1 In (a), the horizontal axis represents the DC component of ambient background light, expressed in gray levels, and the vertical axis represents the relative error of modulation amplitude, expressed as a percentage. Figure 1 In (b), the horizontal axis represents the code length, and the vertical axis is the same as before. Figure 1 In the middle (c), the horizontal axis represents the lens attenuation factor and the vertical axis represents the residual error of water level interpretation, with the unit being mm.
[0132] Depend on Figure 1 As can be seen in (a), the relative error of Example 1 remains stable between 0.11% and 0.14% within a two-order-of-magnitude range of background light, and its value shows no visible correlation with the background light intensity. The relative error of Comparative Example 1 monotonically increases from 2.91% to 10.49%, differing by one to two orders of magnitude, and the difference widens further with increasing background light. This indicates that the encoding-driven and correlation-demodulated methods separate the light signal upon which water level determination depends from the ambient illuminance, making the determination conditions consistent under both strong daylight and no-light conditions.
[0133] In step S21 of this invention, all rows in the Walsh-Hadamard matrix except the row where all elements are +1 are taken as driving codes. The remaining rows have an equal number of elements with values of +1 and -1, meaning the sum of the elements in each driving code within one code period is zero. In the correlation demodulation of step S24, the grayscale of the region of interest is composed of three parts: the ambient background component, the local unit's illumination component, and the leakage light component from adjacent units. These three components are multiplied by the driving code of the local unit and summed over the time slot. The background component remains approximately constant within one code period, while the sum of the driving code elements is zero, and the sum tends towards zero. The leakage light component of adjacent units contains the driving code of the adjacent units and is orthogonal to the driving code of the local unit; its sum is also zero. Only the component of the local unit itself accumulates because the square of the driving code value is always 1. Comparative Example 1 uses the grayscale of the neighboring background area to replace the background of the target area. The spatial non-uniformity of both and the leakage light from adjacent units are directly included in the estimation value, so its error increases linearly with the background light intensity.
[0134] Depend on Figure 1 As can be seen in (b), the relative error of Example 3 is 0.19% at a code length of 128, Example 1 is 0.12% at a code length of 256, and Example 2 is 0.09% at a code length of 512. These three values coincide with the theoretical curve plotted using the reciprocal of the square root of the code length. In contrast, Comparative Example 1 has an average relative error of 4.68% in the same background light region, which is more than an order of magnitude higher than the examples. This indicates that the residual error after demodulation is determined solely by the readout noise of the photosensitive device. Increasing the code length can reduce this error according to a predictable pattern. However, the main error in Comparative Example 1 is the background estimation bias, which cannot be eliminated by increasing the number of acquisition frames. This conclusion is also derived from steps S21 and S24 because the readout noise is unrelated to the driving code. It accumulates in the correlation summation using a random walk and decays after being divided by the code length using the reciprocal of the square root of the code length.
[0135] Depend on Figure 1 As can be seen in (c), the residual error of Example 1 does not exceed 0.20 mm throughout the entire attenuation range, and that of Example 3 does not exceed 1.24 mm. Neither of these errors shows a visible correlation with the lens attenuation factor. In Comparative Example 2, the residual error increases from 1.38 mm to 22.99 mm, and after the attenuation factor drops below 0.50, there are no longer any valid intersection points, resulting in reading failure. This indicates that self-reference normalization and adaptive reading thresholds ensure that the water level reading does not change with lens transmittance, and lens contamination during long-term operation can maintain effective reading without frequent manual cleaning.
[0136] The water film and mud adhering to the lens surface, rain and fog in the air, and the slow changes in exposure parameters cause an approximately proportional attenuation to all light-emitting scale units. Step S31 normalizes the modulation amplitude of each unit above the water surface using the median value as a reference. This proportional attenuation is canceled out in the ratio, so the modulation depth only reflects the relative transmittance of each unit's own optical path. The judgment threshold in step S32 is taken as the midpoint between the two segments above and below the water surface, and is automatically adjusted according to the water transmittance. Comparative Example 2 normalizes the nominal brightness at the factory, and the attenuation factor is directly incorporated into the modulation depth. Its fixed threshold of 0.5 is no longer located between the two segments above and below the water surface when the transmittance drops below this threshold, and the judgment becomes invalid. The residual error of Example 3 is slightly higher than that of Example 1. This is because the modulation depth difference corresponding to the interpolation denominator is smaller when the segment height is larger, and the water level error caused by the same modulation depth noise increases accordingly according to the variance formula in step S34.
[0137] Experimental Example 2, using the methods of Example 1 and Comparative Examples 3, 4, and 5, verifies the effect of turbid water density conversion and pressure compensation on eliminating deviations in the static pressure channel system. The results are as follows: Figure 2 As shown.
[0138] The experiment proceeded as follows: under conditions of a water depth of 4.00m and a water temperature of 18℃, the sand content was initially fixed at 20kg / m³. 3 The flow velocity was increased from 0 to 2.5 m / s, and then the flow velocity was fixed at 1.5 m / s while the sand content was increased from 0 to 50 kg / m³. 3 Calculate the static pressure channel water level under the four compensation configurations according to steps S42 and S43 respectively, subtract the true value from the actual value, and perform least squares linear fitting on the obtained deviation.
[0139] Figure 2 In (a), the horizontal axis represents the ratio of the square of the flow velocity to twice the gravitational acceleration, in meters (m), and the vertical axis represents the water level deviation in the static pressure channel, in millimeters (mm). Figure 2 The horizontal axis of (b) represents the sediment content, in kg / m³. 3 The vertical axis is the same as before. The markers in both figures are the calculated values for each working point, and the straight line of the same color as the marker is the least squares fitted line for that set of calculated values. Each entry in the legend contains both a marker and a straight line, representing the calculated value and its fitted line for that configuration, respectively.
[0140] Depend on Figure 2 As can be seen in (a), the deviation of Comparative Example 4 increases linearly with the flow velocity and head, reaching −114.1 mm at a flow velocity of 2.5 m / s, with a fitting slope of −0.36, consistent with the calibrated pressure coefficient of −0.35; the deviation of Comparative Example 3 also shows a linear relationship but shifts upwards overall, crossing zero at a flow velocity close to 1.6 m / s; the calculated values of Example 1 are randomly distributed near the zero line, with residuals not exceeding 5 mm. Figure 2As can be seen in (b), the deviation of Comparative Example 5 increases linearly with the sand content, with a sand content of 50 kg / m³. 3 The thickness reached 118.5 mm, and the fitting slope was 2.43 mm·m. 3 / kg, multiplied by the water depth by 1 and subtracted from the ratio of the density of clear water to the density of sediment particles, yields a theoretical value of 2.47 mm·m. 3 The residual for Example 1 was consistent with the standard weight ( / kg); the deviation for Comparative Example 3 was 78.3 mm; and the residual for Example 1 was no more than 3 mm.
[0141] The conclusion drawn from these two sets of results is that there are two independent systematic deviations in the static pressure channel. One is proportional to the sediment concentration, and the other is proportional to the flow velocity head. Only by applying both compensations at the same time can the deviation be reduced to the level of random noise. Furthermore, the two deviations have opposite signs and cancel each other out at a certain flow velocity. Therefore, it cannot be concluded that no compensation is needed just because the deviation is close to zero under a certain operating condition.
[0142] In step S42 of this invention, increasing the sediment content increases the density of the turbid water, thus increasing the pressure at the pressure measuring hole at the same water depth. If the water level is overestimated based on the density of clear water, step S42 calculates the turbid water density from the sediment content and water temperature, and then calculates the pressure head accordingly, eliminating this deviation. The fitting slope in the figure is consistent with the theoretical value. The pressure measuring hole is located in the flow acceleration zone between the two guide shields. The measured static pressure at this location is lower than the undisturbed static pressure; the difference is equal to the pressure coefficient multiplied by the incoming flow pressure. Step S43 compensates for this difference by using the product of the signed pressure coefficient and the flow head. When the pressure coefficient is negative, the compensation term increases the calculated water level. Therefore, the fitting slope of Comparative Example 4 is exactly equal to the calibrated pressure coefficient.
[0143] Experimental Example 3, using the methods of Examples 1, 2, and 3 and Comparative Examples 6, 7, and 9, verifies the effect of the leave-one-out consistency test and the variance of time-varying observation noise. The results are as follows: Figure 3 As shown.
[0144] The experiment was divided into two parts. The first part involved a Monte Carlo experiment with three independent, unbiased channels and observation standard deviations of 6mm, 10mm, and 20mm, respectively. Bias increments from 0 to 100mm were injected into designated channels, and each frame was repeated 4000 times. Standardized residuals were calculated using the leave-one-out criterion of Example 1 and the self-inclusion criterion of Comparative Example 9, and the detection rate within a single fusion cycle was statistically analyzed. The decision thresholds for both criteria were determined based on the principle of the same false alarm rate. The leave-one-out criterion was set to 3.0, whose false alarm rate is consistent with the tail probability of three times the two-sided standard deviation of a standard normal distribution. The threshold for the self-inclusion criterion was obtained by searching for the same false alarm rate using a bisection method; its value was lower than 3.0 and was output by the program in the command window. This part statistically analyzed the detection rate within a single fusion cycle. In actual operation, step S63 also required that the standardized residuals continuously exceed the isolation threshold for a set number of cycles before isolation was performed. Therefore, the actual isolation was more robust, but required more cycles. The second part generates a 4-hour true value sequence of single-peak fluctuation water levels. From 2880s to 4320s, the standard deviation of the optical channel observations is increased by 5 times to represent rain, fog, and lens attachments. From 7200s to 8280s, the radar channel shows a 0.15m offset to represent a lock-out jump. From 10800s to 11880s, the standard deviation of the hydrostatic channel observations is increased by 3 times to represent intensified turbulence. The process is then repeated from steps S71 to S73 and subtracted from the true value.
[0145] Figure 3 In (a), the horizontal axis represents the injection bias amplitude in mm, the left vertical axis represents the fault detection rate on the optical channel, and the right vertical axis represents the difference in detection rates of the two criteria on the radar channel. Figure 3 In (b), the horizontal axis represents time in seconds, and the vertical axis represents the moving root mean square of the water level rise rate error on a logarithmic scale in m / h. The three light-colored backgrounds represent the three fault injection periods, as explained in the legend. Figure 3 The vertical axis of (c) represents the root mean square error of the water level rise rate, in m / h, and is divided into two groups: the full-time period and the fault period. Figure 3 The vertical axis of (d) represents the root mean square error of the water level, using a logarithmic scale, with the unit being mm.
[0146] Depend on Figure 3As shown in (a), for the optical channel with the smallest observation variance, Example 1 has a detection rate of 0.42 at an offset of 30 mm, which increases to 0.95 at 50 mm; Comparative Example 9 has a detection rate of 0 at 30 mm, which only reaches 0.47 at 70 mm, meaning that approximately 2.3 times the offset amplitude is required to achieve the same detection level as Example 1. For the radar channel with the largest observation variance, the absolute value of the difference in detection rates between the two criteria does not exceed 0.009 over the entire amplitude range. This indicates that the defect of the criterion containing its own residuals is not universally present in all channels, but is concentrated in the channel with the highest weight, that is, concentrated in the measurement value that the system trusts the most, which is the most dangerous situation; the leave-one consistency test eliminates this defect, making the criteria of each channel have a consistent statistical meaning.
[0147] This invention sets the weight of a channel to be equal to the reciprocal of its observation variance divided by the sum of the reciprocals of the observation variances of all channels. The difference between the channel's measurement and the fused value including itself is exactly equal to the proportion by which the numerator of the residual is compressed to 1 minus the weight. Channels with smaller variances have larger weights and stronger compression. Simultaneously, the denominator used for the self-inclusion criterion is the square root of the sum of the observed variance of the tested channel and the fused variance. Since the fused variance is smaller than the leave-one fusion variance, this denominator is only slightly smaller than the denominator of the leave-one criterion, insufficient to compensate for the numerator compression. Based on the three-channel observation variances in this experimental example, the weight of the optical channel is 0.69, and its test statistic including itself is only 0.43 times that of the leave-one criterion, while the weight of the radar channel is 0.06, with a ratio of 0.94. Under the condition of the same false alarm rate, the optical channel requires an offset amplitude of about 2.2 times to achieve the same detection rate with its own criterion. The two are basically the same in the radar channel. This corresponds perfectly to the two detection rate curves in the figure and the difference in detection rate that is close to zero. The true standard deviation with its own residual is equal to 1 minus the weight and then multiplied by the square root of the sum of the observed variance of the tested channel and the left-one fusion variance. However, the denominator used for the criterion is not equal to the true standard deviation. The larger the weight of the tested channel, the further the two deviate, and the criterion is thus diluted by itself. The left-one fusion makes the tested channel not participate in the formation of the reference value. The two are independent of each other, and the variance of their difference is exactly equal to the sum of the observed variance of the tested channel and the left-one fusion variance. Therefore, the criterion has the same statistical meaning in each channel.
[0148] Depend on Figure 3As shown in (b), the moving mean square error of the inflation rate in Example 1 was basically below that of Comparative Example 7 during the three fault injection periods, and close during the fault-free periods, with medians of 0.048 m / h and 0.050 m / h, respectively. The curve of Comparative Example 6 rose to 0.97 m / h during the period when the radar channel showed bias, more than ten times that of Example 1 during the same period; however, during the period when the observation variance of the optical channel increased, Comparative Example 6 was lower than Comparative Example 7. This is because Comparative Example 7 still allocated the fixed highest weight to the already deteriorated optical channel, while the equal-weighted average of Comparative Example 6 limited the weight of this channel to one-third. This phenomenon further illustrates the defect of fixed weights when the channel deteriorates. All three curves showed a peak near the flood peak, which originated from the lag of the constant velocity model at the point of rapid inflation rate change. This peak was common to all three methods and did not affect the comparison. Figure 3 As can be seen in (c), the root mean square error of the inflation rate during the fault period for Examples 1, 2, and 3 are 0.051 m / h, 0.050 m / h, and 0.060 m / h, respectively, while for Comparative Example 7 it is 0.106 m / h and for Comparative Example 6 it is 0.301 m / h; the corresponding values for the entire period are 0.176 m / h, 0.204 m / h, 0.146 m / h, 0.183 m / h, and 0.328 m / h, respectively. Figure 3 As can be seen from Figure (d), the root mean square errors of water level in Examples 1, 2, and 3 are 1.01 mm, 0.86 mm, and 1.35 mm, respectively, while those in Comparative Example 7 are 1.55 mm and Comparative Example 6 are 13.08 mm. This indicates that channel deterioration only occurs in local periods, and its impact on the water level itself is limited, while its impact on the rate of water level rise is significantly amplified. The real-time observation variance suppresses the error precisely during these periods.
[0149] In this invention, steps S34, S44, S61, and S73 provide real-time observation variance based on the current physical quality indicators of each channel. During rain and fog periods, the modulation depth step of the optical channel becomes gentler, and the sample standard deviation increases, leading to an increase in its observation variance and an automatic decrease in its weight. This ensures consistency between the fusion level and the fusion variance. Step S73 then uses this fusion variance as the observation noise variance for recursive filtering. During deterioration periods, the filter automatically reduces its gain, relying more on model predictions; during recovery periods, it automatically increases its gain, tracking real fluctuations more quickly. In Comparative Example 7, both the weights and the observation noise variance are constants. Deteriorated channels still participate in fusion according to their nominal weights. The filter's estimation of its own prediction uncertainty is inconsistent with reality, leading to overconfidence. The observation noise is amplified in the inflation rate after differential processing. The equal-weighted average in Comparative Example 6 completely fails to reflect the accuracy differences between channels, and its level error is an order of magnitude higher than in Example 1. In Example 2, the segment height and fusion period are taken at the lower limit, resulting in the smallest level error but a slightly higher inflation rate error. The values in Example 3 are the opposite.
[0150] Experimental Example 4, using the methods of Examples 1 and 4 and Comparative Example 8, verifies the early warning performance of the alarm time classification and water level-flow consistency gating, as well as the separability of the ice blockage criterion. The results are as follows: Figure 4 As shown.
[0151] The experiment involved generating 300 real flood events, 300 false flood events (due to malfunctions), and 120 ice blockage events. Each event lasted 1 hour, with a sampling interval of 60 seconds. In real flood events, the water level and cross-sectional flow rate rose synchronously according to the water level-flow rate relationship. In false flood events, the water level rose while the flow rate remained constant, with 10% accompanied by slight natural flow fluctuations. In ice blockage events, the water level rose while the flow rate decreased. The Pearson correlation coefficient between the filtered water level sequence and the cross-sectional flow rate sequence was calculated in step S82. The alarm arrival time was calculated in step S83, and the receiver operating characteristic (ROC) curve was obtained by scanning the alarm arrival time threshold. During actual operation, step S84 selected the first of the three time limits. The auxiliary confidence score obtained by cross-correlation of sediment concentration with time delay in step S82 was only marked in the report and did not change the warning level; therefore, it was not plotted separately.
[0152] Figure 4 In (a), the horizontal axis represents the false alarm rate, and the vertical axis represents the hit rate. Figure 4 In (b), the horizontal axis represents the correlation coefficient between the filtered water level and the cross-sectional flow rate, and the vertical axis represents the water level rise rate, with the unit being m / h. The three types of working conditions are represented by different colors and different mark shapes. The two vertical reference lines are the first threshold and its opposite, respectively, and the horizontal reference line is the starting value for ice jam flood.
[0153] Depend on Figure 4 As can be seen in (a), the receiver operating characteristic curve (ROC) of Comparative Example 8 is close to the diagonal, with an area under the curve of 0.663; the curve of Example 1 is significantly closer to the upper left corner. When the hit rate is 0.90, the false alarm rate of Comparative Example 8 is 0.757 while that of Example 1 is 0.073; when the hit rate is 0.95, the false alarm rate of Comparative Example 8 is 0.867 while that of Example 1 is 0.080. This indicates that when only water level and flood rate are used as criteria, false flooding caused by faults is statistically almost indistinguishable from real floods. After adding a consistency gating of water level and flow rate, the false alarm rate decreases by about an order of magnitude at the same hit rate.
[0154] In a real flood, the flood wave passing through the monitoring station causes the water level and flow rate to increase synchronously, and the two are highly positively correlated. False rises caused by piezometer blockage, floating debris obstruction, or channel zero-point drift do not accompany the synchronous increase in flow, resulting in a low correlation coefficient. The gating system excludes such events accordingly. Residual false alarms arise from a few events where natural flow fluctuations occur during the fault period. In these cases, the two factors do indeed change synchronously, and cannot be distinguished solely by the station's data. This also demonstrates that the gating system reduces false alarms rather than completely eliminating them. Step S84 uses the alarm arrival time rather than the water level itself as the classification metric, allowing the warning level to directly correspond to the available time margin for action. Different levels are given for the same water level under the two scenarios of gradual and rapid rises. Using the correlation coefficient between water level and cross-sectional flow as the basis for confidence gating is because, in a real flood, the water level and flow rate increase synchronously when the flood wave passes through the monitoring station, and the two are highly positively correlated. False rises caused by piezometer blockage, floating debris obstructing the radar beam, or zero-point drift in a certain channel do not accompany the synchronous increase in flow, resulting in a low correlation coefficient. False rises do not trigger flood warnings. The values of each time limit should be determined according to the preparation time required for each level of response in the flood control response process of this station, and the first threshold should be selected according to the statistical results of the historical flood sequence and historical false alarm event sequence of this station.
[0155] Depend on Figure 4 As shown in (b), the three types of operating conditions form three non-overlapping clusters on the discrimination plane. The correlation coefficient of the real flood events is 0.926 with a standard deviation of 0.070, the correlation coefficient of the ice jam events is -0.833 with a standard deviation of 0.109, and the correlation coefficient of the fault-related false flood events is 0.071 with a standard deviation of 0.251. All ice jam events fall to the left of the opposite number of the first threshold, and the water level rise rate is higher than the starting value of the ice jam. The discrimination hit rate according to step S85 is 0.95. This indicates that the same statistic can simultaneously distinguish the three causes of flood, ice jam, and equipment failure on the same plane, but the required treatment measures for these three situations are not the same.
[0156] When ice flowers and floating ice accumulate in bends, shallows, or upstream of buildings, forming ice jams, the cross-sectional area of the water passage is reduced by the ice, increasing water flow resistance. This causes the upstream water level to rise, while the flow rate through the monitoring station decreases instead of increasing. Therefore, water level and flow rate show a negative correlation, the opposite of the characteristic of both increasing in the same direction during floods. Water level drift caused by equipment failure does not change the flow rate, and the correlation coefficient is near zero. Step S85 also requires a positive water level rise rate and requires the above conditions to be continuously set for a certain number of windows to exclude inverse changes in water level and flow rate during the receding period and short-term disturbances. The causes and handling measures for ice jam hazards and flood hazards differ; distinguishing the hazard type in the report helps flood control departments take appropriate measures.
[0157] Although embodiments of this application have been shown and described above, it is understood that the above embodiments are exemplary and should not be construed as limiting this application. Those skilled in the art can make changes, modifications, substitutions and variations to the above embodiments within the scope of this application.
Claims
1. A method for operating a smart sensor integrated hydrological gauge system, characterized in that, Includes the following steps: S1. Multiple luminous scale units and a camera unit facing each other are set at equal intervals on the water-facing side of the water gauge column. Pressure measuring holes are opened in the low-disturbance area of the column side wall and connected to pressure transmitters. Layered flow velocity units, water temperature chains and turbidity units are set along the elevation. A radar unit is set on the top of the column. Data from each unit is timestamped by a common clock and aligned according to the fusion cycle before being sent to the edge computing main control unit. S2. The encoding driving module controls the near-infrared light-emitting diodes of each light-emitting scale unit to turn on and off according to the time slot using mutually orthogonal bipolar driving codes. The camera unit synchronously acquires multiple frames of images, and the average gray level of the region of interest of each unit is correlated with its driving code to obtain the modulation amplitude. S3. Normalize the modulation depth by the median of the modulation amplitude of each unit above the water surface, determine the reading threshold by the median of the modulation depth of the submerged units, and interpolate between two adjacent units that cross the threshold to obtain the optical channel water level and its observation variance. S4. Calculate the sediment content from the turbidity and the density of the turbid water from the water temperature. Use this density to convert the gauge pressure of the pressure measuring hole into a pressure head and use the flow velocity at this elevation to correct the dynamic pressure deficit, so as to obtain the water level of the static pressure channel and its observation variance. S5. Obtain the radar channel water level and its observation variance from the radar unit. Perform a leave-one-out test on each channel using the fusion value of the other channels to determine the health status and isolate abnormal channels. Obtain the fusion water level and fusion variance for healthy channels by weighting the variance inversely. S6. Using the fused water level as the observation and the fused variance as the time-varying observation noise variance, recursively filter the model with water level and rise rate as state variables, and output the filtered water level and water level rise rate. S7. Calculate the cross-sectional flow rate from the stratified flow velocity and the stratified cross-sectional area that varies with the filtered water level, and calculate the correlation coefficient between the filtered water level and the cross-sectional flow rate within the sliding window. S8. When the filtered water level is lower than the warning water level and the rise rate is greater than the starting value, the alarm time is calculated. When the correlation coefficient is not lower than the threshold, the flood warning is output in stages according to the alarm time. When the coefficient is negative, the rise rate is positive and the temperature of the submerged water temperature measuring point is not higher than the set value, the ice jam warning is output. When the channel isolation continues for a set number of cycles, the equipment alarm is output and the recovery action is triggered.
2. The method according to claim 1, characterized in that, Step S1 includes: S11. Arrange N luminous scale units at equal intervals from bottom to top along the water-facing side of the water gauge column. Each luminous scale unit includes a near-infrared light-emitting diode with a center wavelength of 850nm and a white light-emitting diode, all sharing a section of light guide plate. The driving circuits of each luminous scale unit are independent of each other. The center elevation of the k-th luminous scale unit is determined by the following formula: ; In the formula, The center elevation of the k-th illuminated scale unit; is the zero elevation of the water gauge; k is the serial number of the luminous scale unit, which is an integer from 1 to N; N is the total number of luminous scale units; The segment height of the illuminated scale unit; S12. Install a camera unit on the bracket facing the water-facing side of the water gauge column. The photosensitive device of the camera unit does not have a near-infrared cutoff filter. A narrow-band filter with a center wavelength of 850nm is installed in front of its lens. The external trigger input terminal of the camera unit is connected to the frame synchronization output terminal of the encoding drive module. S13. A pressure measuring hole is opened on the side wall between the two guide shrouds on both sides of the water gauge column, flush with the column surface. The axis of the pressure measuring hole is perpendicular to the side wall and orthogonal to the direction of the incoming flow. The pressure measuring hole is connected to a gauge pressure transmitter with an atmospheric pressure reference tube through a pressure guide tube. This arrangement makes the pressure measuring hole located in the low disturbance zone where the flow is not separated, and the pressure tapping coefficient remains stable after calibration. S14. Arrange the velocity sensors of each layer of the layered velocity unit, the measuring points of the water temperature chain, and the turbidity unit along the elevation of the water gauge column. Install the radar unit on the top of the water gauge column from above. The sampling data of each unit is timestamped by a common clock in the edge computing main control unit and aligned to the same time according to the fusion period with zero-order hold or linear interpolation to generate time-aligned multi-channel data.
3. The method according to claim 1, characterized in that, Step S2 includes: S21. Take N rows from the L-order Walsh-Hadamard matrix except for the row where all elements are +1, and assign the k-th row to the k-th luminous scale unit as its driving code. L is the code length, which is an integer power of 2 that is not less than the sum of N and 1. S22. In the nth time slot, when the driving code of the kth light-emitting scale unit is +1, the near-infrared light-emitting diode of that light-emitting scale unit is lit; when it is -1, the near-infrared light-emitting diode is turned off. n is the time slot number, which is an integer from 1 to L. During the period of acquiring multiple frames of images, each white light-emitting diode is turned off and the exposure time and gain of the camera unit are locked. S23. During the installation and calibration phase, the near-infrared light-emitting diodes of each light-emitting scale unit are lit up one by one. The region of interest corresponding to each light-emitting scale unit is determined by the position of the bright area in the image and stored in the edge computing main control unit. In the same calibration process, under the condition that the lens is clean and the water surface is stable, the water surface elevation is independently measured by leveling or water level probe. The elevation at the point where the modulation depth crosses the interpretation threshold is calculated according to steps S31 to S33. The difference between the independently measured water surface elevation and the elevation is taken as the optical channel elevation correction constant and stored in the edge computing main control unit. S24. The camera unit continuously acquires L frames of images according to the frame synchronization signal. For the k-th region of interest, the average gray value is calculated frame by frame. If the frame number is missing or the region is gray-saturated, the current image is discarded and re-acquired. Then, correlation demodulation is performed according to the following formula to obtain the modulation amplitude: ; In the formula, is the modulation amplitude of the k-th luminous scale unit; L is the code length; The average gray value of the k-th region of interest in the n-th frame image; The value of the driving code for the k-th luminous scale unit in the n-th time slot; To sum the time slot numbers n from 1 to L; To obtain the absolute value; this step suppresses the influence of the ambient light DC component on the modulation amplitude; S25. Cross-correlate the frame-by-frame average gray value of the kth region of interest with the driving codes of other light-emitting scale units respectively. Take the maximum value among the obtained cross-correlation values. When the ratio of the maximum value to the modulation amplitude of the kth light-emitting scale unit is greater than a set ratio, the kth light-emitting scale unit is determined to be a crosstalk abnormal unit. Its modulation amplitude is removed and linearly interpolated by the modulation amplitudes of the two adjacent light-emitting scale units that have not been removed, so as to obtain the modulation amplitude of each light-emitting scale unit after crosstalk verification.
4. The method according to claim 1, characterized in that, Step S3 includes: S31. The optical channel water level output from the previous fusion cycle is added to the two segment heights to determine the luminous scale units above the water surface. For the first run or when there is no effective optical channel water level in the previous fusion cycle, the radar channel water level or hydrostatic channel water level at the same moment is used instead. The median of the modulation amplitude of each luminous scale unit above the water surface is taken as a reference value. The modulation depth is obtained by normalizing the modulation amplitude of each luminous scale unit using the following formula: ; In the formula, The modulation depth of the k-th luminous scale unit; The modulation amplitude of the k-th luminous scale unit; For reference modulation amplitude; S32. Take the median of the modulation depth of each light-emitting scale unit whose serial number is more than two times lower than the water surface, and determine the reading threshold using the following formula: ; In the formula, The threshold for judgment; The median of the modulation depth of each light-emitting scale unit whose serial number is more than two times the serial number of the light-emitting scale unit above the water surface; S33. Search from bottom to top for adjacent unit pairs where the modulation depth of this unit is lower than the judgment threshold and the modulation depth of the adjacent unit above it is not lower than the judgment threshold. If neither unit is eliminated, the unit pair is considered valid. If there are more than two valid unit pairs, take the smaller absolute value of the difference between its interpolation result and the fusion level of the previous fusion cycle. Record the index of the lower unit in the selected unit pair as k. Obtain the elevation where the modulation depth crosses the judgment threshold by linear interpolation according to the following formula: ; In the formula, To modulate the elevation at the point where the depth of modulation crosses the interpretation threshold; The center elevation of the k-th illuminated scale unit; The segment height of the illuminated scale unit; The threshold for judgment; The modulation depth of the k-th luminous scale unit; The modulation depth of the (k+1)th luminous scale unit is then superimposed with the positive constant of the optical channel elevation calibrated in step S23 according to the following formula to obtain the optical channel water level: ; In the formula, For optical channel water level; To modulate the elevation at the point where the depth of modulation crosses the interpretation threshold; The elevation of the optical channel is corrected to a constant value, which is obtained by calibration in step S23 and remains unchanged during operation; S34. The observation variance of the optical channel is given by the following formula: ; In the formula, For the observation variance of the optical channel; The segment height of the illuminated scale unit; The standard deviation of the modulation depth of each luminous scale unit above the water surface; The modulation depth of the (k+1)th luminous scale unit; The modulation depth of the k-th luminous scale unit; the observation variance is not less than one-twelfth of the square of the segment height; if there is no effective unit pair or the difference between the two modulation depths in the effective unit pair is less than the set lower limit, the optical channel level will not be output in this fusion cycle, and the channel will be treated as missing measurement.
5. The method according to claim 1, characterized in that, Step S4 includes: S41. Calculate sediment content from turbidity using the following formula: ; In the formula, denoted as ρ, where α is the turbidity conversion slope; U is the turbidity measured by the turbidity unit; b is the turbidity conversion intercept; α and b are obtained by least-squares linear fitting of the turbidity value obtained by the drying and weighing method from water samples collected synchronously at this station section. When the turbidity exceeds the turbidity range used for fitting, the turbidity of this fusion cycle is set to invalid and the most recent valid value is used. S42. Calculate the density of turbid water from the sediment content and water temperature using the following formula: ; In the formula, The density of turbid water; The density of pure water is obtained by interpolating from the pure water density table by measuring the water temperature at the elevation of the pressure test hole using the water temperature chain. It refers to the sand content; For the density of sediment particles, the measured value of sediment particle density at this station is used. S43. Take the interpolated velocity of the stratified velocity unit at the elevation of the pressure measuring hole, and obtain the water level in the static pressure channel using the following formula: ; In the formula, This refers to the water level in the static pressure channel. ρ is the elevation of the pressure measuring hole; p is the gauge pressure output by the pressure transmitter. Where is the density of the turbid water; g is the acceleration due to gravity; The pressure measurement coefficient is a signed coefficient that can be positive or negative. The water level in the tank is kept constant and the inflow velocity is changed in stages. The difference between the static pressure channel water level without dynamic pressure correction and the actual water level in the tank is fitted with a least squares straight line to the ratio of the square of the velocity to twice the gravitational acceleration, and the slope is obtained. The water flow velocity at the elevation of the pressure measuring hole; S44, take the length as The sample variance of the static pressure channel water level within the sliding window after delinearization is the static pressure channel observation variance. This observation variance is not less than one-twelfth of the square of the ratio of the pressure transmitter's pressure resolution to the product of the turbid water density and gravitational acceleration. This represents the length of the sliding window.
6. The method according to claim 1, characterized in that, Step S5 includes: S61. The water level of the radar channel is obtained by subtracting the measured distance to the water surface from the installation elevation of the radar unit. The length is taken as... The sample variance of the radar channel water level after delinearization within the sliding window is the radar channel observation variance, which is not less than one-twelfth of the square of the radar element resolution. S62. For the i-th channel, perform leave-one-out fusion on the remaining channels and calculate the standardized residuals using the following formula: ; ; ; In the formula, Reserve one fusion water level for the i-th channel; Let the leave-one fusion variance be the i-th channel; Let be the standardized residual of the i-th channel; i and m are the channel numbers, which are taken as optical channel, hydrostatic channel and radar channel respectively; and These are the water levels of the i-th channel and the m-th channel, respectively. and These are the observation variances of the i-th channel and the m-th channel, respectively; To sum over all channels except the i-th channel; S63. When the standardized residual of the i-th channel continuously exceeds the isolation threshold for a set number of cycles, the health status of the channel is set to abnormal and isolated. When the standardized residual of the isolated channel continuously falls below the acceptance threshold for a set number of cycles, its health status is restored. The acceptance threshold is lower than the isolation threshold. When the number of healthy channels is not less than 2, the fusion water level is obtained by weighting the reciprocals of their observation variances, and the reciprocal of the sum of the reciprocals of the observation variances of each healthy channel is taken as the fusion variance. When the number of healthy channels is 1, the water level of that channel is taken as the fusion water level and its observation variance is taken as the fusion variance, and a downgrade flag is set. When the number of healthy channels is 0, the fusion water level is not output in this fusion cycle, and the data is not entered into steps S6 to S8. S64. During the period when the absolute value of the water level rise rate is lower than the set value and the number of healthy channels is not less than 2, the deviation status of the i-th channel is updated according to the following formula: ; In the formula, and These represent the deviation states of the i-th channel after the current update and the previous update, respectively. The updated coefficients take values greater than 0 and less than 1. To obtain the median of all samples within the window; Reserve one fusion water level for the i-th channel; Let be the water level of the i-th channel; the absolute value of the deviation state is constrained by a set upper limit, and when it exceeds the upper limit, the health status of the channel is set to abnormal; thereafter, the sum of the water level of the channel and its deviation state is used in steps S62 and S63.
7. The method according to claim 1, characterized in that, Step S6 includes: S71. Establish discrete state equations using water level and water level rise rate as state variables. The state transition matrix and observation matrix are respectively taken as follows: ; ; In the formula, F is the state transition matrix; H is the observation matrix; The filtering period is defined as follows: the initial values of the state variables are taken as the first effective fused water level and zero, and the diagonal elements of the initial covariance matrix are taken as the first fused variance and the square of the water level rise range, respectively. S72. Obtain the process noise covariance matrix using the following formula: ; In the formula, is the process noise covariance matrix; q is the process noise intensity driving the change in water level rise rate. The continuous water level sequence during the non-flood season of this station is taken, and the value is selected within the set upper and lower limits to make the average value of the normalized squared innovation of the one-step prediction tend to 1. The filter period; S73. Using the merged water level as the observation and the merged variance as the observation noise variance, the filter water level and water level rise rate are output by recursively applying the prediction equation and update equation of the Kalman filter cycle by cycle. When there is no merged water level in the current fusion cycle, only the prediction equation is executed and the update equation is not executed. When the continuous missing measurement exceeds the set number of cycles, it is re-initialized according to S71 to obtain the filter water level and water level rise rate that are still continuously usable during the channel deterioration period.
8. The method according to claim 1, characterized in that, Steps S7 and S8 include: S81. Calculate the cross-sectional flow rate using the following formula: ; In the formula, Q is the cross-sectional flow rate; To sum the values of layer numbers j from 1 to J; The velocity conversion factor for the j-th layer is determined by taking the synchronous measured flow rate of the standard velocity-area method or acoustic Doppler velocity profiler as a reference, and according to the ratio of the measured velocity at each layer point to the average velocity of the corresponding layer cross section. The velocity is the velocity measured by the velocity sensor at the j-th layer. is the cross-sectional area of the j-th layer, which is obtained by dividing the measured cross-section of the station according to the filtered water level cycle by cycle; J is the number of flow velocity sensor layers below the water surface, which is determined cycle by cycle according to the filtered water level; S82, in a length of The Pearson correlation coefficient between the filtered water level sequence and the cross-sectional flow sequence is calculated within the sliding window, and the time-delay cross-correlation coefficient between the filtered water level sequence and the sediment concentration sequence is calculated within the set time delay range. The larger absolute value and the corresponding time delay are taken as the auxiliary confidence level. S83. When the filtered water level is lower than the warning water level and the water level rise rate is greater than the starting value, the alarm arrival time shall be calculated according to the following formula. When the filtered water level is not lower than the warning water level or the water level rise rate is not greater than the starting value, the alarm arrival time shall not be output: ; In the formula, Time of arrival at the police station; This is the warning water level for this station; For filtering water level; The rate of increase in water level; S84. When the correlation coefficient between the filtered water level and the cross-sectional flow is not lower than the first threshold, it is determined to be a rising water process. If the alarm arrival time is lower than the first time limit, a yellow warning is issued; if it is lower than the second time limit, an orange warning is issued; if it is lower than the third time limit, a red warning is issued; otherwise, a blue warning is issued. The first time limit is greater than the second time limit, and the second time limit is greater than the third time limit. When the filtered water level is not lower than the warning water level, a red warning is issued. When the auxiliary confidence level is not lower than the second threshold, the confidence level is marked as increased in the message; when the auxiliary confidence level is lower than the second threshold, the warning level is not changed. S85. When the correlation coefficient between the filtered water level and the cross-sectional flow rate is not higher than the negative number of the first threshold, the water level rise rate is greater than the ice jam start value, and the temperature of the first water temperature measuring point below the filtered water level from top to bottom is not higher than the freezing judgment temperature, and the above three conditions are met simultaneously and the window number is continuously set, it is judged as an ice jam hazard and an ice jam warning is output. S86. When the channel health status is abnormal and the set number of cycles continues, when the sensor self-test fails, or when the channel is not re-accepted after the recovery action is executed, the output device alarms and marks the corresponding channel. This output is unrelated to the correlation coefficient of step S82. S87. The warning level, hazard type, filtered water level, water level rise rate, alarm arrival time, and channel health status are compiled into a message and sent in parallel through each channel of the communication unit. If any channel receives a successful acknowledgment, the message is considered delivered. If none of the channels receive a acknowledgment, the message will be retried at the set interval.
9. A smart sensor integrated hydrological gauge system, used to implement the method described in any one of claims 1-8, characterized in that, The system includes a water gauge column, an array of luminous scale units evenly spaced along the water-facing side of the water gauge column, an encoding drive module, a camera unit facing the water-facing side, an optical calculation unit, a static pressure measurement unit with pressure measuring holes on the side wall of the column, a turbidity unit, a water temperature chain, a density dynamic pressure compensation unit, a stratified flow velocity unit, a radar unit located on top of the water gauge column, a fusion arbitration unit, a flow calculation unit, a rise rate filtering unit, an early warning classification unit, a communication unit, and an edge computing master control unit with an internal common clock. The output of the encoding drive module is connected to the drive circuit of each luminous scale unit, its frame synchronization output is connected to the external trigger input of the camera unit, and the image output of the camera unit is connected to the optical calculation unit. The static pressure measurement unit, turbidity unit, water temperature chain, and stratified flow velocity unit are connected to the density dynamic pressure compensation unit via a data bus. The outputs of the optical calculation unit, density dynamic pressure compensation unit, and radar unit are respectively connected to the fusion arbitration unit. The water level output of the fusion arbitration unit is connected to the rise rate filtering unit, and its health status output is connected to the early warning classification unit. The outputs of the stratified flow velocity unit and the expansion rate filtering unit are respectively connected to the flow calculation unit; the outputs of the flow calculation unit, the turbidity unit, the water temperature chain and the expansion rate filtering unit are respectively connected to the early warning classification unit, and the output of the early warning classification unit is connected to the communication unit; the optical calculation unit, the density dynamic pressure compensation unit, the fusion arbitration unit, the expansion rate filtering unit and the early warning classification unit are integrated into the edge computing main control unit.
10. The system according to claim 9, characterized in that, The water level gauge column has a flat surface facing the water, with flow guides on both sides. The pressure measuring hole of the static pressure measurement unit is located on the side wall between the two flow guides, with the hole axis perpendicular to the side wall and orthogonal to the incoming flow direction. Each light-emitting scale unit in the light-emitting scale unit array contains one near-infrared LED and one white LED sharing a light guide plate. The controlled ends of both are connected to the encoding output and illumination output of the encoding drive module, respectively. The encoding drive module turns off the illumination output during image acquisition. The photosensitive device of the camera unit does not have a near-infrared cutoff filter; its lens has a narrow-band filter with the same center wavelength as the near-infrared LED. The water temperature chain's measuring points are distributed along the elevation and cover the elevation of the pressure measuring hole and the water level variation range, with stratified flow velocity... The flow velocity sensors for each layer of the system are mounted on the extended support arms on the side of the water gauge column; the radar unit is a radar sensor with two or more switchable beams, and its beam control terminal is connected to the self-diagnostic calibration unit; the system also includes a self-diagnostic calibration unit, a backwashing assembly, and a complementary power supply unit. The input terminal of the self-diagnostic calibration unit is connected to the health status output terminal of the fusion arbitration unit, and its output terminal is connected to the control terminals of the coding drive module, the backwashing assembly, and the radar unit, respectively. The backwashing assembly includes a flushing pump, a solenoid valve, a check valve, and a flushing pipeline connected to the pressure measuring hole. The inlet of the flushing pump forms a loop with the flushing pipeline via a storage tank; the complementary power supply unit includes a photovoltaic module, a charging control circuit, and a lithium battery pack, and its output terminal is connected to the power supply input terminal of each unit in the system.
