Pipe gallery inspection robot environment parameter monitoring method and system and storage medium

Abstract

The application discloses a pipeline gallery inspection robot environment parameter monitoring method and system and a storage medium. The pipeline gallery inspection robot environment parameter monitoring method and system and the storage medium comprise the following steps: acquiring a pipeline gallery field space point cloud and a first gas concentration sequence which are synchronously collected by a laser radar and a gas sensor carried by a pipeline gallery inspection robot during pipeline gallery inspection and which are established in space according to the same sampling time; and extracting an obstacle shielding feature in the pipeline gallery according to the pipeline gallery field space point cloud, and determining a space boundary of an air flow vortex area. The pipeline gallery inspection robot environment parameter monitoring method and system and the storage medium improve the accuracy and reliability of gas environment parameter monitoring and abnormal source positioning in a complex pipeline gallery environment.

Description

Technical Field

[0001] This invention relates to the field of robot environmental monitoring technology, specifically to environmental parameter monitoring methods, systems, and storage media for utility tunnel inspection robots. Background Technology

[0002] During the inspection of utility tunnels, using inspection robots equipped with gas sensors to continuously monitor the gas environment parameters inside the tunnels is an important means of timely detection of abnormalities such as gas leaks. Existing methods usually determine abnormal areas based on the gas concentration and its changing trends collected by the robot at different locations, and some methods further utilize concentration gradients to track the location of gas anomaly sources.

[0003] However, utility tunnels typically contain numerous pipes, supports, and other facilities, with long, narrow spaces and complex obstructions. Gas diffusion within these tunnels is easily affected by both obstacles and airflow. Airflow passing over obstacles may be deflected, backflowed, or vortexed, causing leaked gas to accumulate in localized areas or break into multiple discontinuous concentration clusters. Therefore, the actual gas concentration distribution often does not change continuously and smoothly with distance from the anomaly source, but may exhibit significant fragmentation. When a robot passes through vortices or areas of localized gas accumulation, even if the area is not near the actual anomaly source, it may detect high gas concentrations, interfering with subsequent anomaly location determination. Simultaneously, gas sensors typically have response and recovery times. When the inspection robot continuously passes through areas with significant concentration differences, the sensor output may not promptly reflect changes in the actual environmental concentration. For example, after the robot leaves a high-concentration area, the sensor may still maintain a high reading for a period, resulting in numerical tailing due to response delay. The numerical trailing effect, combined with the concentration abrupt changes inherent in the pipe gallery, makes it difficult for the collected gas concentration to accurately reflect the real gas environment at the robot's current location. Therefore, in pipe gallery environments with complex obstructions and local vortex airflow, monitoring environmental parameters based solely on the original gas concentration or local concentration changes is easily affected by the distribution of concentration fragments and sensor response delays. This can lead to confusion between abnormally high concentration areas and real gas anomaly sources, reducing the accuracy of anomaly source location identification and environmental parameter monitoring results. Summary of the Invention

[0004] The purpose of this invention is to provide a method, system and storage medium for monitoring environmental parameters of a utility tunnel inspection robot, thereby solving the problems existing in the prior art.

[0005] To achieve the above objectives, the present invention provides the following technical solution: a method for monitoring environmental parameters of a utility tunnel inspection robot, comprising: S1. Obtain the spatial point cloud and first gas concentration sequence of the utility tunnel site, which are synchronously collected by the lidar and gas sensor carried by the utility tunnel inspection robot during the utility tunnel inspection process and have a spatial correspondence established according to the same sampling time. S2. Extract obstacle and occlusion features within the utility tunnel based on the spatial point cloud at the site, and determine the spatial boundary of the airflow vortex region. S3. Perform sliding window variance calculation on the first gas concentration sequence within the spatial boundary of the airflow vortex region to obtain the concentration distribution feature matrix. S4. If the largest element in the concentration distribution feature matrix is ​​greater than the preset variance threshold, then it is determined that there is a fragment distribution phenomenon in the spatial boundary of the airflow vortex region. S5. Extract the concentration abrupt difference of the first gas concentration sequence at adjacent sampling times from the spatial boundary of the airflow vortex region where fragmentation exists. S6. If the concentration change difference exceeds the preset change threshold, it is determined that the first gas concentration sequence has a numerical tail caused by response delay. S7. The Kalman filter algorithm is used to perform state estimation and noise reduction on the first gas concentration sequence with numerical tails to obtain the second gas concentration sequence. S8. By extracting the spatial gradient features of the second gas concentration sequence and the geometric center distance of the spatial boundary of the airflow vortex region, a source feature vector is constructed. S9. Support vector machine is used to perform high-dimensional spatial mapping and classification of source feature vectors to determine the true source coordinates of abnormal gas environment parameters in the utility tunnel, and the true source coordinates are used as the environmental parameter monitoring results of the utility tunnel inspection robot.

[0006] Preferably, step S1 includes acquiring the original point cloud and concentration sequence, which are synchronously generated by the lidar and gas sensor carried by the utility tunnel inspection robot during the utility tunnel inspection process; performing timestamp alignment processing on the original point cloud and concentration sequence to obtain an aligned sequence; establishing a spatial correspondence between the aligned sequence to obtain a mapping sequence; determining whether there are data gaps in the mapping sequence; if there are data gaps in the mapping sequence, then performing data gap filling processing through a spatial interpolation fitting algorithm to obtain the on-site spatial point cloud and the first gas concentration sequence of the utility tunnel.

[0007] Preferably, step S2 includes acquiring the spatial point cloud of the utility tunnel and extracting surface normals and curvature features; extracting obstacle occlusion features based on the surface normals and curvature features, wherein the obstacle occlusion features include occlusion contours; extracting the velocity gradient matrix based on the occlusion contours; determining the vortex center if the velocity gradient matrix meets a preset threshold; generating boundary coordinates based on the vortex center; and determining the spatial boundary of the airflow vortex region based on the boundary coordinates.

[0008] Preferably, step S3 includes acquiring sampling nodes within the spatial boundary of the airflow vortex region, extracting a first gas concentration sequence based on the sampling nodes, truncating the first gas concentration sequence through a sliding window to obtain a local concentration subsequence, calculating the variance of the local concentration subsequence to obtain a local dispersion value, and constructing a concentration distribution feature matrix based on the local dispersion value.

[0009] Preferably, step S4 includes extracting the maximum element of the concentration distribution feature matrix; comparing the maximum element with a preset variance threshold to obtain a threshold comparison result; if the threshold comparison result indicates that the maximum element of the matrix is ​​greater than the preset variance threshold, then generating a boundary anomaly marker; and determining the presence of debris distribution within the spatial boundary of the airflow vortex region based on the boundary anomaly marker.

[0010] Preferably, step S5 includes obtaining a first gas concentration sequence corresponding to the spatial boundary of the airflow vortex region; performing a difference operation on the first gas concentration sequence to obtain a concentration gradient vector; performing feature extraction on the concentration gradient vector to obtain a difference feature set; determining whether the feature value of the difference feature set is greater than a preset mutation threshold; if the feature value is greater than the preset mutation threshold, then extracting the concentration mutation difference of the first gas concentration sequence at adjacent sampling times from the spatial boundary.

[0011] Preferably, step S6 includes: if the concentration change difference exceeds a preset change threshold, extracting the corresponding adjacent sampling times to construct a delayed response interval; truncating the first gas concentration sequence according to the delayed response interval to obtain a concentration decay feature sequence; obtaining the tail duration and difference fluctuation amplitude according to the concentration decay feature sequence to construct a tail determination condition; if the concentration decay feature sequence meets the tail determination condition, determining that the first gas concentration sequence has a numerical tail caused by response delay.

[0012] Preferably, step S7 includes obtaining a first gas concentration sequence with numerical tails; constructing a state transition matrix and an observation matrix based on the first gas concentration sequence; generating a priori state estimation sequence using the state transition matrix and the observation matrix; if the system noise covariance calculated using the priori state estimation sequence is greater than a preset threshold, determining the Kalman gain matrix based on the observation noise covariance; and denoising the priori state estimation sequence using the Kalman gain matrix to obtain a second gas concentration sequence.

[0013] Preferably, step S8 includes acquiring a second gas concentration sequence, reconstructing the second gas concentration sequence to obtain a concentration distribution matrix, and extracting the spatial gradient features of the concentration distribution matrix; determining the airflow vortex region based on the spatial gradient features, acquiring the spatial boundary of the airflow vortex region and calculating the centroid position of the boundary; calculating the geometric center distance of the spatial boundary through the centroid position of the boundary; and constructing a source feature vector based on the spatial gradient features and the geometric center distance.

[0014] Preferably, step S9 includes using a support vector machine to perform high-dimensional space mapping on the source feature vector to obtain a high-dimensional feature space; extracting a feature distribution matrix based on the high-dimensional feature space; determining the decision boundary hyperplane through the feature distribution matrix and performing classification processing to obtain anomaly classification labels; if the anomaly classification label indicates an anomaly source, extracting the spatial coordinate grid corresponding to the anomaly source; parsing the spatial coordinate grid to determine the true source coordinates; and using the true source coordinates as the environmental parameter monitoring results of the utility tunnel inspection robot.

[0015] An environmental parameter monitoring system for a utility tunnel inspection robot, used to implement the steps of the aforementioned method for monitoring the environmental parameters of a utility tunnel inspection robot, the system comprising: The data acquisition module acquires the spatial point cloud and the first gas concentration sequence of the utility tunnel site, which are synchronously collected by the lidar and gas sensor carried by the utility tunnel inspection robot during the utility tunnel inspection process and establish spatial correspondence according to the same sampling time. The vortex region determination module extracts obstacle and occlusion features within the utility tunnel based on the spatial point cloud at the site, and determines the spatial boundary of the airflow vortex region. The concentration feature calculation module performs sliding window variance calculation on the first gas concentration sequence within the spatial boundary of the airflow vortex region to obtain the concentration distribution feature matrix. The debris distribution determination module determines that if the largest element in the concentration distribution feature matrix is ​​greater than the preset variance threshold, then the debris distribution phenomenon exists within the spatial boundary of the airflow vortex region. The concentration mutation extraction module extracts the concentration mutation difference of the first gas concentration sequence at adjacent sampling times from the spatial boundary of the airflow vortex region where fragmentation exists. The numerical tailing determination module determines that if the concentration change difference exceeds the preset change threshold, the first gas concentration sequence has numerical tailing caused by response delay. The concentration filtering module uses the Kalman filter algorithm to perform state estimation and noise reduction on the first gas concentration sequence with numerical tails to obtain the second gas concentration sequence. The source feature construction module constructs a source feature vector by extracting the spatial gradient features of the second gas concentration sequence and the geometric center distance of the spatial boundary of the airflow vortex region. The anomaly source localization module uses a support vector machine to perform high-dimensional spatial mapping and classification of the source feature vectors, determine the true source coordinates of the gas environment parameter anomalies in the utility tunnel, and use the true source coordinates as the environmental parameter monitoring results of the utility tunnel inspection robot.

[0016] A computer-readable storage medium storing computer instructions thereon, which execute the environmental parameter monitoring method of the utility tunnel inspection robot when the computer instructions are executed.

[0017] As can be seen from the above technical solution, the present invention has the following beneficial effects: This method, system, and storage medium for environmental parameter monitoring of a utility tunnel inspection robot establishes a spatial correspondence between the spatial point cloud of the utility tunnel collected by lidar and the gas concentration sequence collected by gas sensors at the same sampling time. It then identifies airflow vortex regions by combining the characteristics of obstacles within the utility tunnel and uses concentration distribution features to determine the distribution of gas concentration fragments within the vortex regions. This reduces the interference of abnormally high concentration areas caused by obstacles such as pipes and supports, as well as localized vortex airflow, on the monitoring results. Furthermore, it identifies numerical tailing caused by gas sensor response delay based on concentration abrupt changes between adjacent sampling times. Through state estimation and denoising, it obtains a second gas concentration sequence that better reflects the actual gas concentration state at the current spatial location. Based on this, it constructs a source feature vector by combining the spatial gradient characteristics of the gas concentration and its distance from the geometric center of the vortex region. This source feature vector is then classified to determine the true source coordinates of abnormal gas environmental parameters within the utility tunnel. This reduces the combined impact of concentration fragmentation and sensor response delay on abnormal source identification, improving the accuracy and reliability of gas environmental parameter monitoring and abnormal source location in complex utility tunnel environments. Attached Figure Description

[0018] Figure 1 This is a flowchart of the overall method of the present invention; Figure 2 This is a thermogram of the concentration distribution characteristic matrix of the airflow vortex region in the experiment of this invention; Figure 3 This is a time series curve of gas concentration in the experiment of this invention; Figure 4 This is a scatter plot of the results of locating the anomaly source in the pipe gallery plane during the experiment of this invention; Figure 5 This is a bar chart comparing the average positioning errors of the three methods used in the experiment of this invention. Detailed Implementation

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

[0020] Example 1:

[0021] like Figure 1 As shown, the present invention provides a technical solution: a method for monitoring environmental parameters of a utility tunnel inspection robot, comprising: S1. Obtain the spatial point cloud and first gas concentration sequence of the utility tunnel site, which are synchronously collected by the lidar and gas sensor carried by the utility tunnel inspection robot during the utility tunnel inspection process and have a spatial correspondence established according to the same sampling time. S2. Extract obstacle and occlusion features within the utility tunnel based on the spatial point cloud at the site, and determine the spatial boundary of the airflow vortex region. S3. Perform sliding window variance calculation on the first gas concentration sequence within the spatial boundary of the airflow vortex region to obtain the concentration distribution feature matrix. S4. If the largest element in the concentration distribution feature matrix is ​​greater than the preset variance threshold, then it is determined that there is a fragment distribution phenomenon in the spatial boundary of the airflow vortex region. S5. Extract the concentration abrupt difference of the first gas concentration sequence at adjacent sampling times from the spatial boundary of the airflow vortex region where fragmentation exists. S6. If the concentration change difference exceeds the preset change threshold, it is determined that the first gas concentration sequence has a numerical tail caused by response delay. S7. The Kalman filter algorithm is used to perform state estimation and noise reduction on the first gas concentration sequence with numerical tails to obtain the second gas concentration sequence. S8. By extracting the spatial gradient features of the second gas concentration sequence and the geometric center distance of the spatial boundary of the airflow vortex region, a source feature vector is constructed. S9. Support vector machine is used to perform high-dimensional spatial mapping and classification of source feature vectors to determine the true source coordinates of abnormal gas environment parameters in the utility tunnel, and the true source coordinates are used as the environmental parameter monitoring results of the utility tunnel inspection robot.

[0022] In this embodiment, the utility tunnel inspection robot moves inside the utility tunnel along a pre-set inspection path, the lidar continuously collects spatial point cloud data inside the utility tunnel according to a pre-set scanning cycle, and the gas sensor continuously collects concentration data of the target gas according to a pre-set sampling cycle.

[0023] Specifically, the lidar and gas sensor use the same system clock to record the acquisition time, so that each frame of spatial point cloud and each gas concentration value has corresponding time information. Taking the sampling time of each gas concentration value as a reference, the spatial point cloud with the smallest absolute time difference with the sampling time is found in the lidar acquisition time, and the time correspondence between the gas concentration value and the corresponding spatial point cloud is established. When the sampling time of the gas concentration value is located between the acquisition times of two adjacent frames of spatial point clouds, the robot spatial pose corresponding to the sampling time of the gas concentration value is determined according to the robot pose corresponding to the two adjacent frames of spatial point clouds and the time position of the gas concentration value sampling time between the two point cloud acquisition times.

[0024] Furthermore, the spatial transformation relationship between the lidar coordinate system, the gas sensor coordinate system, and the unified coordinate system of the utility tunnel is calibrated before the inspection.

[0025] Specifically, based on the robot's spatial pose at each gas sampling moment and the installation position of the gas sensor relative to the robot body, the sampling position corresponding to each gas concentration value is transformed to the unified coordinate system of the utility tunnel. After the above processing, each gas concentration value has a definite sampling moment and a definite three-dimensional spatial position, thereby forming a spatiotemporal correspondence between the spatial point cloud of the utility tunnel site and the first gas concentration sequence.

[0026] In addition, outlier deletion and spatial downsampling were performed on the spatial point cloud of the utility tunnel site.

[0027] For each spatial point, a spherical neighborhood is established with a pre-defined neighborhood radius, and the number of spatial points contained within the spherical neighborhood is counted. When the number of spatial points within the neighborhood is less than a pre-defined threshold, the current spatial point is identified as an outlier and deleted from the spatial point cloud of the utility tunnel site. After the outlier deletion is completed, the remaining spatial points are divided into three-dimensional voxels according to the pre-defined voxel side length. Each voxel containing spatial points retains a representative point to form point cloud data with unified spatial distribution.

[0028] Furthermore, based on the downsampled spatial point cloud of the utility tunnel, the tunnel walls and obstacles located in the passage space of the utility tunnel are identified.

[0029] Specifically, a local neighborhood is established for each spatial point. The normal direction and curvature of the local surface are determined based on the positional relationship of the spatial points within the local neighborhood. The spatial points are then divided into regions based on the distance, difference in normal direction, and difference in curvature between adjacent spatial points. The set of points that can form the continuous inner wall surface of the utility tunnel is defined as the utility tunnel wall point set. The set of continuous spatial points located inside the space enclosed by the utility tunnel wall and not assigned to the utility tunnel wall point set is defined as the obstacle point set. The obstacle point set corresponds to pipelines, supports, equipment, and other physical structures occupying the internal space of the utility tunnel.

[0030] Furthermore, the longitudinal extension direction of the utility tunnel is determined based on the spatial point cloud at the site, and this direction is used as a reference direction for airflow propagation.

[0031] For each set of obstacle points, its outer contour is extracted, and an obstruction space is established backward from the obstacle's outer contour along the airflow propagation reference direction. The lateral and vertical extent of the obstruction space is determined by the obstacle's outer contour, and the extension direction of the obstruction space is consistent with the airflow propagation reference direction. For obstruction spaces formed by adjacent obstacles, spatial overlap is determined. If two obstruction spaces have an actual overlapping volume, they are assigned to the same candidate region. Subsequently, spatial connectivity is determined for all candidate regions, and interconnected candidate regions are merged. The outermost spatial contour of the merged region is extracted. This outermost spatial contour is determined as the spatial boundary of the airflow vortex region.

[0032] Furthermore, based on the spatial boundary of the airflow vortex region, the first gas concentration data located inside the spatial boundary is read and arranged in chronological order of sampling time.

[0033] Specifically, the number of consecutive sampling points contained in the sliding window is preset. Starting from the beginning of the first gas concentration sequence, the preset number of consecutive gas concentration values ​​are selected as a calculation window each time. For each calculation window, all gas concentration values ​​inside the window are summed and then divided by the number of gas concentration values ​​contained in the window to obtain the average concentration value corresponding to the window. Then, the difference between each gas concentration value in the window and the average concentration value is calculated. Each difference is squared and then all squared results are summed and finally divided by the number of gas concentration values ​​contained in the window to obtain the concentration variance corresponding to the calculation window.

[0034] Furthermore, after a calculation window completes the concentration variance calculation, the sliding window is moved forward one sampling point along the sampling time sequence, and the average concentration calculation and concentration variance calculation are repeated until the sliding window covers the last position of the first gas concentration sequence that allows complete sampling. The concentration variances obtained from each sliding window are arranged according to the corresponding sampling time and spatial sampling position to form a concentration distribution feature matrix. Each element in the concentration distribution feature matrix corresponds to a discrete concentration state within a certain time range and a certain spatial location range within the airflow vortex region.

[0035] As a preferred embodiment, numerical comparisons are performed one by one from all elements of the concentration distribution feature matrix, and the element with the largest value is retained as the maximum variance value corresponding to the airflow vortex region.

[0036] Specifically, the variance threshold is preset based on the background environment calibration results of the gas sensor before the inspection task begins and stored in the parameter storage unit of the pipe gallery inspection robot. The maximum variance value is compared with the preset variance threshold. When the maximum variance value is greater than the preset variance threshold, it is determined that there is a fragment distribution phenomenon within the spatial boundary of the current airflow vortex region. When the maximum variance value is not greater than the preset variance threshold, the current airflow vortex region does not enter the numerical tail recognition process.

[0037] Furthermore, the fragmentation phenomenon refers to the phenomenon where, after gas enters the airflow vortex region, the return path and circulating flow path formed by obstacles change the original diffusion path of the gas, resulting in the formation of multiple concentration distribution zones with spatial intervals within the same airflow vortex region.

[0038] After quantifying the distribution state using sliding window variance, as long as the concentration variance corresponding to at least one window exceeds the preset variance threshold, the existence of fragmented distribution is determined according to the above judgment rules.

[0039] In addition, for airflow vortex regions where fragmentation is identified, the first gas concentration sequence is read in chronological order.

[0040] For any two consecutive sampling times, the gas concentration value corresponding to the previous sampling time is subtracted from the gas concentration value corresponding to the later sampling time, and the absolute value of the calculation result is taken. The obtained absolute value is determined as the concentration change difference between two adjacent sampling times. All adjacent sampling times in the first gas concentration sequence are processed in sequence to obtain the corresponding concentration change difference sequence.

[0041] As a preferred embodiment, the mutation threshold is determined based on the step response calibration results of the gas sensor before the pipe gallery inspection robot performs the inspection task, and is stored in the parameter storage unit.

[0042] Specifically, the concentration abrupt change difference corresponding to each adjacent sampling time is compared with a preset abrupt change threshold in turn; when at least one concentration abrupt change difference is greater than the preset abrupt change threshold, it is determined that the current first gas concentration sequence has numerical tailing caused by the gas sensor response delay; when all concentration abrupt change differences are not greater than the preset abrupt change threshold, Kalman filtering processing for numerical tailing is not performed.

[0043] The formation process of numerical tail is as follows: After the inspection robot moves from one sampling position to the next, the actual gas concentrations corresponding to the two sampling positions change. The gas sensor needs to go through its own determined response process to switch from the previous concentration state to the next concentration state. Therefore, the sensor output value corresponding to the next sampling position still contains the response component corresponding to the previous sampling position. When the inspection robot moves continuously, this response component is continuously retained in the measurement values ​​corresponding to multiple sampling positions along the robot's movement path. This causes the first gas concentration sequence to form extended data that is inconsistent with the actual spatial concentration change position. This extended data is the numerical tail.

[0044] In addition, the Kalman filter algorithm is used for state estimation and noise reduction of the first gas concentration sequence with numerical tails.

[0045] Specifically, the estimated real gas concentration and the corresponding gas concentration change rate at each sampling time are used together as state information. First, the corrected real gas concentration estimate and concentration change rate at the previous sampling time are read. Based on the time interval between adjacent sampling times, the concentration change rate at the previous sampling time is multiplied by the time interval, and then added to the real gas concentration estimate at the previous sampling time to obtain the predicted gas concentration at the current sampling time. The concentration change rate at the previous sampling time is used as the predicted concentration change rate at the current sampling time.

[0046] Furthermore, after obtaining the predicted gas concentration at the current sampling time, the uncertainty of the current prediction result is updated according to the pre-set process noise parameters, and then the actual measured value in the first gas concentration sequence at the current sampling time is read.

[0047] Specifically, the difference between the actual measured value and the predicted gas concentration is calculated, and the correction weight for the current sampling time is determined based on the uncertainty of the prediction result and the measurement noise parameters pre-calibrated by the gas sensor. The difference is multiplied by the correction weight to obtain the concentration correction amount for the current sampling time, and then the concentration correction amount is added to the predicted gas concentration to obtain the corrected estimate of the true gas concentration for the current sampling time. At the same time, the concentration change rate and state uncertainty are updated based on the current correction result, and the updated result is passed to the next sampling time.

[0048] In addition, according to the sampling time sequence of the first gas concentration sequence, state prediction, measurement comparison, weight correction and state update are performed sequentially for each sampling time.

[0049] The corrected gas concentration estimates at each sampling time are rearranged according to the original sampling order to obtain the second gas concentration sequence. Each gas concentration value in the second gas concentration sequence continues to use the sampling time and spatial coordinates corresponding to the first gas concentration sequence, so that the concentration data after state estimation and the spatial location of the utility tunnel are kept in a one-to-one correspondence.

[0050] In this embodiment, spatial gradient features are extracted from the second gas concentration sequence.

[0051] Specifically, for each spatial sampling point in the second gas concentration sequence, a three-dimensional neighborhood with a pre-defined range is established with the spatial sampling point as the center, and adjacent sampling points with positional differences from the current spatial sampling point in the first coordinate direction, the second coordinate direction, and the third coordinate direction are selected from the neighborhood respectively; for each coordinate direction, the gas concentration value of the current spatial sampling point is subtracted from the gas concentration value of the adjacent sampling point, and the resulting concentration difference is divided by the coordinate difference between the two sampling points in the corresponding coordinate direction to obtain the spatial change of concentration in that coordinate direction.

[0052] Furthermore, after obtaining the spatial changes in concentration corresponding to the three coordinate directions in sequence, the three spatial changes in concentration are squared respectively, the three squared results are added together, and then the square root operation is performed on the sum to obtain the spatial gradient amplitude corresponding to the current spatial sampling point.

[0053] The spatial gradient feature of the current spatial sampling point is formed by the spatial change of concentration in three coordinate directions and the magnitude of the spatial gradient. This spatial gradient feature is used to represent the concentration change of the second gas concentration sequence along the three spatial directions at the current spatial sampling location.

[0054] In this embodiment, the geometric center of the region is determined based on all boundary points on the spatial boundary of the airflow vortex region.

[0055] Specifically, the coordinate values ​​of all boundary points in the first coordinate direction are added together and divided by the total number of boundary points to obtain the coordinate value of the geometric center in the first coordinate direction; the same processing method is used to obtain the coordinate values ​​of the geometric center in the second and third coordinate directions; the geometric center position of the spatial boundary of the airflow vortex region is determined by the center coordinates in the three coordinate directions.

[0056] In addition, for each spatial sampling point corresponding to the second gas concentration sequence, the coordinate differences between the spatial sampling point and the geometric center in three coordinate directions are calculated. The three coordinate differences are squared respectively, the three squared results are added together, and the square root of the sum is calculated to obtain the spatial distance between the spatial sampling point and the geometric center of the airflow vortex region. The spatial change of concentration in the three coordinate directions, the spatial gradient amplitude, and the spatial distance corresponding to the current spatial sampling point are combined according to a predetermined feature arrangement order to form the source feature vector corresponding to the current spatial sampling point.

[0057] Before the support vector machine is put into actual inspection, the source location of the utility tunnel space to be monitored is classified.

[0058] Specifically, the space of the pipe gallery to be monitored is divided into non-overlapping candidate source units according to the pre-set spatial dimensions. Each candidate source unit is assigned a unique location category, and the center spatial coordinates of the candidate source unit are recorded. During the training phase, gas anomaly sources with known coordinates are set at the locations corresponding to each candidate source unit in sequence. The inspection robot is controlled to obtain lidar point cloud and gas concentration data according to the same data acquisition process as the actual inspection. The vortex region determination, fragment distribution judgment, numerical tail judgment, Kalman filtering processing, spatial gradient extraction, and geometric center distance calculation are completed according to the aforementioned steps.

[0059] Furthermore, a correspondence is established between each source feature vector obtained during training and the location category of the candidate source unit where the actual gas anomaly source is located, and this correspondence is used as training data for the support vector machine.

[0060] Specifically, the support vector machine uses a pre-defined kernel function to perform high-dimensional space mapping on the source feature vectors and establishes classification boundaries between categories based on the training data corresponding to different position categories. After training is completed, each source feature vector input to the support vector machine obtains its corresponding position category through the established classification boundaries.

[0061] During actual utility tunnel inspection, the source feature vector corresponding to the current airflow vortex region is input into a trained support vector machine. After performing high-dimensional spatial mapping and classification, the support vector machine outputs a location category. Based on this location category, the pre-established correspondence between the location category and the candidate source unit is queried to determine the candidate source unit corresponding to the location category, and the center spatial coordinates of the candidate source unit are read. The center spatial coordinates are determined as the true source coordinates of the abnormal gas environment parameters in the utility tunnel, and the true source coordinates are output and stored as the environmental parameter monitoring results of the utility tunnel inspection robot.

[0062] Example 2:

[0063] S1 includes acquiring the original point cloud and concentration sequence, which are synchronously collected and generated by the lidar and gas sensor carried by the utility tunnel inspection robot during the utility tunnel inspection process; performing timestamp alignment processing on the original point cloud and concentration sequence to obtain an aligned sequence; establishing a spatial correspondence between the aligned sequence to obtain a mapping sequence; determining whether there are data gaps in the mapping sequence; if there are data gaps in the mapping sequence, performing data gap filling processing through a spatial interpolation fitting algorithm to obtain the on-site spatial point cloud and the first gas concentration sequence of the utility tunnel.

[0064] In this embodiment, after the utility tunnel inspection robot starts its inspection, the lidar continuously outputs raw point cloud frames according to the scanning cycle set by the device, and the gas sensor continuously outputs gas concentration values ​​according to the sampling cycle set by the device.

[0065] The lidar and gas sensor use the same clock source provided by the inspection robot controller. Each frame of raw point cloud is written with the corresponding point cloud acquisition timestamp, and each gas concentration value is written with the corresponding concentration acquisition timestamp. The raw point cloud is arranged in ascending order according to the point cloud acquisition timestamp, and the gas concentration values ​​are arranged in ascending order according to the concentration acquisition timestamp, thus forming the raw point cloud and concentration sequence. The lidar scanning cycle, gas sensor sampling cycle, and timestamp timing accuracy are all read from the configuration parameters of the corresponding devices and remain unchanged during this inspection task.

[0066] Furthermore, timestamp alignment uses the sampling time of the gas sensor as a unified reference time.

[0067] Specifically, for each gas concentration sampling time, firstly, the previous frame of the point cloud with the smallest time difference and a timestamp before the sampling time is retrieved from the original point cloud, and then the next frame of the point cloud with the smallest time difference and a timestamp after the sampling time is retrieved. When there is a point cloud frame in the original point cloud with the same timestamp as the gas concentration sampling timestamp, the point cloud frame is directly associated with the corresponding gas concentration value. When there is no matching timestamp, the time interval between the gas concentration sampling time and the sampling time of the previous frame of the point cloud is calculated, and then the time interval between the previous frame of the point cloud and the next frame of the point cloud is calculated. The time interpolation ratio is obtained by dividing the two time intervals mentioned above.

[0068] In addition, the point cloud of the previous frame and the point cloud of the next frame are converted to the unified coordinate system of the utility tunnel before the time interpolation is performed.

[0069] In specific processing, the robot pose corresponding to the two point cloud acquisition times is read. The point coordinates in the LiDAR's own coordinate system are first transformed to the robot's body coordinate system according to the pre-calibrated installation extrinsic parameters between the LiDAR and the robot body, and then transformed to the unified coordinate system of the pipe gallery according to the robot pose at the corresponding acquisition time. The robot translation pose is processed separately in three coordinate directions. The translation position corresponding to the gas concentration sampling time is obtained by linear interpolation of the two point cloud acquisition positions according to the time interpolation ratio. The robot rotation pose is obtained by spherical linear interpolation of the attitude quaternions of the two acquisition times, and the interpolation ratio is consistent with the aforementioned time interpolation ratio. Thus, the complete robot pose corresponding to the gas concentration sampling time is determined.

[0070] Furthermore, for the point cloud frames located on both sides of the gas concentration sampling time, the point cloud after being converted to the unified coordinate system of the pipe gallery is spatially divided according to the same voxel size as the point cloud downsampling in the embodiment.

[0071] The average coordinates of points within the same spatial voxel are used to form representative voxel points. A correspondence is established between representative points at the same voxel position in the point cloud of the previous frame and the point cloud of the next frame. For two representative points with established correspondence, the difference between their three spatial coordinate components is calculated. The difference between each coordinate component is multiplied by the aforementioned time interpolation ratio and added to the corresponding coordinate component of the previous frame. The resulting three coordinate components form the interpolation point corresponding to the gas concentration sampling time. All interpolation points are combined according to their original voxel positions to form an aligned point cloud for the gas concentration sampling time. The above process is repeated according to the gas concentration sampling time order to obtain an aligned sequence with the gas concentration sampling time as a unified time reference.

[0072] In addition, when establishing spatial correspondences for alignment sequences, the installation position of the gas sensor is represented by the sensor extrinsic parameters calibrated before the inspection.

[0073] Specifically, for each alignment moment, the gas sensor detection position is first transformed to the robot body coordinate system based on the fixed installation position of the gas sensor relative to the robot body. Then, the robot pose determined at that moment is transformed to the unified coordinate system of the pipe gallery to obtain the three-dimensional sampling coordinates of the current gas concentration value. Subsequently, the current gas concentration value, three-dimensional sampling coordinates, corresponding alignment point cloud and sampling timestamp are written into the same mapping data unit. Each mapping data unit is arranged in ascending order according to the sampling timestamp to form a mapping sequence.

[0074] Therefore, each effective gas concentration value in the mapping sequence corresponds to a specific spatial location of the pipe gallery and is associated with point cloud data at the same reference time.

[0075] As a preferred embodiment, data holes in the mapping sequence are determined directly according to the device sampling period.

[0076] Specifically, the first effective gas concentration sampling time is taken as the starting time, and reference sampling times that should exist in this inspection process are generated sequentially according to the fixed sampling period set by the gas sensor. Each reference sampling time is checked against the mapping sequence one by one. When a reference sampling time does not have a corresponding mapping data unit, or the corresponding mapping data unit is missing any of the following: effective gas concentration value, three-dimensional sampling coordinates, or alignment point cloud, the position corresponding to the reference sampling time is determined as a data hole. This determination is not based on the concentration change amplitude, thereby avoiding misjudging real concentration changes as missing data.

[0077] Furthermore, after identifying the data hole, the spatial sampling position corresponding to the data hole is first restored.

[0078] Specifically, the nearest effective robot pose before and after the data hole is read. The time interval between the sampling time corresponding to the data hole and the previous effective pose is calculated, and this time interval is divided by the time interval between the two effective poses to obtain the position interpolation ratio. The difference between the subsequent effective pose and the previous effective pose in the three translational coordinate directions is multiplied by the position interpolation ratio, and then added to the corresponding coordinate value of the previous effective pose to obtain the robot translational position at the time of the data hole. The rotational attitude is subjected to quaternion spherical linear interpolation according to the same time ratio. Then, coordinate transformation is performed according to the fixed installation extrinsic parameters of the gas sensor to obtain the spatial sampling coordinates of the gas sensor corresponding to the data hole.

[0079] In this embodiment, the gas concentration data is completed using interpolation fitting based on spatial distance.

[0080] Specifically, using the spatial sampling coordinates of the data hole as the retrieval center, four effective gas concentration sampling points with the smallest spatial distance to the location are selected from the mapping sequence. The spatial distance between each effective sampling point and the data hole is calculated by squaring the coordinate differences in the three coordinate directions, adding the three squared results, and taking the square root of the sum. Then, the reciprocals of the four spatial distances are calculated and added together. The interpolation weight of the effective sampling point is obtained by dividing each reciprocal distance by the sum. The gas concentration values ​​of the four effective sampling points are multiplied by their respective interpolation weights, and the four products are added together. The result is used as the complete gas concentration value of the data hole location. When the retrieved effective sampling point coincides with the spatial coordinates of the data hole, the gas concentration value of the effective sampling point is directly used, and the reciprocal distance calculation is not performed.

[0081] When there are missing alignment point clouds, the effective alignment point cloud with the closest time distance before and after the data hole is used as the fitting input.

[0082] Specifically, since both point clouds are already in the unified coordinate system of the utility tunnel, they are divided into voxel grids in the same way as the aforementioned time alignment stage, and the average value of the point coordinates within the corresponding voxels is calculated to form representative points. A correspondence relationship of representative points is established for voxels with the same spatial position in the two point clouds. Based on the time ratio between the sampling time of the data hole and the two valid point cloud times, linear interpolation is performed on the three coordinate components of the corresponding representative points to form the complete point cloud at the time of the data hole. This processing only performs coordinate fitting on the corresponding voxels that have valid point data in both frames. Voxels that do not form a correspondence do not generate interpolation points, so as not to construct new structural points in the original point cloud where there is no spatial constraint.

[0083] Furthermore, after each data hole is filled, the filled gas concentration value, filled spatial sampling coordinates, and filled point cloud are written back to the position corresponding to its original reference sampling time.

[0084] After all data voids are processed, the mapping sequence is rechecked according to the gas sensor sampling cycle to confirm that each reference sampling time has gas concentration value, three-dimensional sampling coordinates and corresponding point cloud data. The point clouds in the mapping sequence are collected according to the unified coordinate system of the utility tunnel, and the points with overlapping coordinates are merged according to the predetermined voxel size to form the spatial point cloud of the utility tunnel site. The gas concentration values ​​corresponding to each reference sampling time are arranged in chronological order to form the first gas concentration sequence.

[0085] Example 3:

[0086] S2 includes acquiring the spatial point cloud of the utility tunnel and extracting surface normals and curvature features; extracting obstacle occlusion features based on the surface normals and curvature features, wherein the obstacle occlusion features include occlusion contours; extracting the velocity gradient matrix based on the occlusion contours; determining the vortex center if the velocity gradient matrix meets a preset threshold; generating boundary coordinates based on the vortex center; and determining the spatial boundary of the airflow vortex region based on the boundary coordinates.

[0087] In this embodiment, the spatial point cloud of the utility tunnel obtained in Embodiment 1 is first read, and the surface normal and curvature features of each point in the spatial point cloud are extracted.

[0088] Specifically, taking the current point cloud point as the center, neighboring points are retrieved within a preset neighborhood radius. The coordinates of the current point and each neighboring point in three coordinate directions are averaged to obtain the center position of the neighborhood point set. Then, the coordinate deviations of each neighboring point relative to this center position in three coordinate directions are calculated, and a local spatial discrete relationship is formed based on all coordinate deviations. The local spatial discrete relationship is decomposed into three mutually perpendicular feature directions. The feature direction corresponding to the minimum feature value is determined as the surface normal direction of the current point cloud point. To ensure that the normal directions are consistent, the obtained surface normals are uniformly oriented according to the direction of the free space inside the pipe gallery, so that the normals of neighboring points all face the internal space of the pipe gallery.

[0089] In this embodiment, the curvature feature is calculated using the same neighborhood point set.

[0090] Specifically, the three feature values ​​are added together to obtain the total discrete value of the neighborhood point set. The smallest feature value is then divided by the total discrete value, and the resulting ratio is used as the curvature feature value of the current point cloud point. All valid points in the spatial point cloud of the utility tunnel are processed in the same way to obtain the surface normal and curvature feature corresponding to each point cloud point. The neighborhood radius is preset according to the resolution of the lidar point cloud. The principle for setting the value is: the neighborhood diameter is greater than 3 times the average distance between adjacent points and less than the width of the smallest obstacle structure to be identified in the utility tunnel. This ensures that each neighborhood contains points that form local surfaces, while avoiding merging independent structures into the same neighborhood.

[0091] When extracting obstacle occlusion features based on surface normals and curvature characteristics, the longitudinal centerline of the utility tunnel is first obtained by fitting the point cloud of the tunnel wall, and the tangential direction of the centerline in the current area is used as the reference direction for airflow propagation.

[0092] The angle between the surface normal of each point cloud point and the reference direction of airflow propagation is calculated, and structural segmentation is performed in conjunction with curvature features. The curvature threshold is determined by the standard pipe gallery wall point cloud collected before inspection. Specifically, for the pipe gallery wall point cloud in the barrier-free area, all curvature feature values ​​are calculated using the aforementioned method. The arithmetic mean of all curvature feature values ​​is obtained. Then, the absolute value of the difference between each curvature feature value and the arithmetic mean is calculated. All absolute values ​​are averaged again. Finally, the sum of the arithmetic mean of the curvature feature values ​​and the average absolute deviation is determined as the curvature threshold. Points with curvature feature values ​​greater than the curvature threshold are included in the obstacle candidate point set.

[0093] Perform spatial connectivity processing on the set of obstacle candidate points.

[0094] The three-dimensional Euclidean distance between any two points is compared with a preset connectivity distance threshold. When the distance between two points is not greater than the connectivity distance threshold, the two points are assigned to the same connected region. The connectivity distance threshold is twice the side length of the downsampled voxel of the lidar. After the connectivity processing is completed, the set of points that are continuously connected to the fitting surface of the inner wall of the pipe gallery is deleted, and the remaining connected point set is determined as the set of obstacle points corresponding to pipelines, supports, equipment shells or protruding components.

[0095] For each set of obstacle points, the obstacle points are projected onto a cross section perpendicular to the airflow propagation reference direction. Edge points located on the periphery of the point set are extracted from the projected cross section, and connections are established from near to far according to the distance between adjacent edge points, so that all edge points form a closed profile.

[0096] Then, restore each edge point in the closed profile to its original three-dimensional coordinate position. The resulting three-dimensional closed profile is the occlusion profile. The occlusion profile defines the cross-sectional area occupied by the obstacle structure when the airflow propagates longitudinally along the pipe gallery.

[0097] When extracting the velocity gradient matrix based on the occlusion profile, continuous analysis sections are established along the longitudinal direction of the pipe gallery in front of the occlusion profile, at the location of the occlusion profile, and behind the occlusion profile. The spacing between adjacent analysis sections is taken as the side length of the point cloud voxels. Each analysis section is divided into regular grids according to the side length consistent with the downsampling of the point cloud. Grids whose center falls inside the occlusion profile are marked as obstacle grids, and the remaining grids located inside the pipe gallery are marked as flow grids.

[0098] Furthermore, the reference flow velocity is obtained through calibration data from the airflow sensor set up before the inspection robot enters the current pipe section.

[0099] Specifically, airflow velocity is continuously collected for 30 sampling cycles within a straight pipe section without obstruction. The 30 sampled values ​​are added together and divided by 30, and the resulting value is used as the reference velocity for that pipe section. For each analysis section, the total number of flow grids is counted, and the effective flow area of ​​the analysis section is obtained by multiplying the number of flow grids by the area of ​​a single grid. The ratio of the effective flow area of ​​the reference section to the effective flow area of ​​the current analysis section is then calculated, and this ratio is multiplied by the reference velocity. The result is used as the cross-sectional velocity of the current analysis section.

[0100] To obtain the local flow velocity of each grid inside the cross section, the spatial distance from the center of each flow grid to the nearest edge point of the occlusion profile is calculated.

[0101] Specifically, the spatial distances of all circulation grids are summed, and the spatial distance of each circulation grid is divided by the sum of these distances to obtain the distance weight of the corresponding grid. The cross-sectional velocity is multiplied by the number of circulation grids, and then multiplied by the distance weight of each grid. The result is then distributed to the corresponding circulation grids to obtain the local velocity of each circulation grid. For grids directly adjacent to the occlusion profile, their local velocity is calculated according to the aforementioned distance weights, without assigning an empirical velocity value, thus ensuring that the local velocity is completely determined by the occlusion profile, the effective circulation area, and the spatial distance.

[0102] Furthermore, after obtaining the local flow velocity of each flow grid, the flow velocity gradient in the three coordinate directions is calculated respectively.

[0103] Specifically, along the longitudinal direction of the utility tunnel, the local velocity of the two adjacent grids before and after the current grid is read. The local velocity of the previous grid is subtracted from the local velocity of the next grid, and then divided by the actual distance between the centers of the two grids to obtain the longitudinal velocity gradient of the current grid. The same method is used to calculate the corresponding velocity gradients along the transverse and vertical directions of the utility tunnel. When the grid is located at the boundary and there are only two adjacent grids on one side, the velocity difference between the current grid and its adjacent grid is divided by the distance between the centers of the two grids. The velocity gradients of each grid in the three directions are written into the corresponding matrix cells according to the grid spatial position to form a velocity gradient matrix.

[0104] In addition, the flow velocity gradient threshold is obtained through calibration of unobstructed pipe sections.

[0105] Specifically, in the same unobstructed pipe section used to obtain the reference flow velocity, the absolute value of the flow velocity gradient of all grids is calculated according to the aforementioned grid division and gradient calculation method; the absolute values ​​of all flow velocity gradients are added together and divided by the number of effective grids to obtain the average gradient value; then the absolute value of the difference between each flow velocity gradient absolute value and the average gradient value is calculated, and the absolute values ​​of all differences are added together and divided by the number of effective grids to obtain the average absolute gradient deviation; the sum of the average gradient value and twice the average absolute gradient deviation is used as the preset flow velocity gradient threshold; this threshold is written into the parameter storage unit of the inspection robot after the calibration of the same pipe section is completed, and remains unchanged during this inspection process.

[0106] In addition, each matrix cell in the velocity gradient matrix is ​​read one by one, and the absolute values ​​of the velocity gradients in the three directions of longitudinal, transverse and vertical are taken respectively. If the absolute value of the velocity gradient in any one of the directions is greater than the preset velocity gradient threshold, the corresponding grid is determined as a gradient anomaly grid.

[0107] Perform 3D connectivity processing on all gradient anomaly meshes, and classify gradient anomaly meshes with adjacent faces into the same connected region.

[0108] In this embodiment, for each connected region, it is further determined whether it meets the conditions for determining the vortex center.

[0109] Specifically, each grid in the connected region is taken as the grid to be determined. The grids adjacent to it in the same cross section are read, and the velocity gradient directions of each adjacent grid are compared in a clockwise order. When the gradient direction of the adjacent grids changes continuously around the grid to be determined, and there are at least two sets of relative position grids whose gradient directions point to both sides of the grid to be determined, the grid to be determined is determined as a candidate grid for the vortex center. If there are more than two candidate grids for the vortex center in the same connected region, the absolute values ​​of the velocity gradients in the three directions of the adjacent grids around each candidate grid are added together. The candidate grid with the largest sum is determined as the grid where the vortex center is located, and the coordinates of the center of this grid are used as the coordinates of the vortex center.

[0110] Furthermore, after determining the vortex center, a layer-by-layer expansion is performed starting from the grid where the vortex center is located.

[0111] Specifically, for each additional layer, all grids adjacent to the outer edge of the current region are read, and it is determined whether at least one of the absolute values ​​of the velocity gradients in the three directions of the adjacent grids is greater than a preset velocity gradient threshold. Grids that meet this condition are merged into the current vortex region; grids whose absolute values ​​of the velocity gradients in the three directions are not greater than the preset velocity gradient threshold are not merged into the vortex region. The above expansion process is repeated until all adjacent grids on the outer periphery of the vortex region no longer meet the merging condition.

[0112] The outermost mesh that is ultimately incorporated into the vortex region is defined as the boundary mesh. The vertex coordinates of each boundary mesh facing outward from the vortex region are extracted, and the vertices are connected sequentially according to the connection relationship of adjacent meshes to form a closed three-dimensional boundary coordinate set.

[0113] The space enclosed by this set of boundary coordinates is defined as the airflow vortex region, and this set of boundary coordinates constitutes the spatial boundary of the airflow vortex region described in Embodiment 1. The first gas concentration sequence corresponding to this spatial boundary is then selected, and the process proceeds to the sliding window variance calculation step described in Embodiment 1.

[0114] Example 4:

[0115] S3 includes acquiring sampling nodes within the spatial boundary of the airflow vortex region, extracting a first gas concentration sequence based on the sampling nodes, truncating the first gas concentration sequence through a sliding window to obtain a local concentration subsequence, calculating the variance of the local concentration subsequence to obtain the local dispersion value, and constructing a concentration distribution feature matrix based on the local dispersion value.

[0116] In this embodiment, the spatial boundary of the airflow vortex region determined in Embodiment 3 and the first gas concentration sequence formed in Embodiment 1 are read first. Each gas concentration value in the first gas concentration sequence retains the corresponding sampling timestamp and three-dimensional sampling coordinates, so the concentration data located in the airflow vortex region can be filtered according to the spatial location.

[0117] For each sampling node in the first gas concentration sequence, the inclusion relationship between the three-dimensional sampling coordinates of the sampling node and the closed spatial boundary of the airflow vortex region is determined.

[0118] In specific processing, starting from the sampling node, a ray is drawn along the positive direction of the first coordinate of the unified coordinate system of the pipe gallery, and the number of valid intersections between the ray and the boundary surface of the airflow vortex region is counted. When the number of intersections is odd, the sampling node is determined to be located inside the airflow vortex region; when the number of intersections is even, the sampling node is determined to be located outside the airflow vortex region. When the sampling node falls on the boundary surface, it is directly regarded as a valid sampling node within the boundary. After completing the judgment of all sampling nodes, only the sampling nodes located inside the airflow vortex region and on the boundary surface are retained.

[0119] Furthermore, the retained sampling nodes are arranged in order of their timestamps, and the gas concentration values ​​associated with each sampling node are read, thereby obtaining the first gas concentration sequence within the spatial boundary of the airflow vortex region.

[0120] This sequence uses the original sampling time and spatial coordinates from Example 1, without interpolating the gas concentration values ​​again or changing the time order between the sampling nodes.

[0121] A sliding window is used to extract the first gas concentration sequence.

[0122] In this embodiment, each sliding window contains 5 consecutive sampling nodes, and the window moves in steps of 1 sampling node. The first sliding window reads the 1st to 5th gas concentration values ​​in chronological order and forms the first local concentration subsequence. After completing this truncation, the starting position of the window is moved back by 1 sampling node, and the 2nd to 6th gas concentration values ​​are read to form the next local concentration subsequence. The window moves continuously in the above order until there are fewer than 5 remaining gas concentration values, at which point the window truncation stops.

[0123] In this embodiment, the sliding window contains 5 sampling nodes, the values ​​of which are predetermined based on the gas sensor sampling period and the inspection speed of the inspection robot.

[0124] Specifically, before the inspection task begins, the total time corresponding to 5 consecutive sampling cycles of the gas sensor is multiplied by the robot's set inspection speed to obtain the inspection distance covered by the 5 sampling nodes. When the inspection distance is less than the length of the airflow vortex region along the robot's travel direction determined in Example 3, the 5 sampling nodes are used as the sliding window length for this inspection. After the parameter setting is completed, the window length and window movement step size remain unchanged in the same inspection task.

[0125] For each local concentration subsequence, calculate the local dispersion value.

[0126] First, the five gas concentration values ​​in the local concentration subsequence are summed, and the sum is divided by 5 to obtain the average concentration value of the local concentration subsequence. Then, the average concentration value is subtracted from each of the five gas concentration values ​​to obtain five concentration deviation values. Each of the five concentration deviation values ​​is squared, and the five squared results are summed and divided by 5. The resulting value is used as the local dispersion value of the local concentration subsequence. The local dispersion value is determined by all concentration data within the corresponding window, without using a single concentration peak as a substitute.

[0127] Each local discreteness value is spatially correlated with the sliding window that generated that value.

[0128] For a sliding window containing 5 sampling nodes, the sampling node that is 3rd in time sequence is taken as the center node of the window, and the three-dimensional sampling coordinates of the center node are determined as the spatial position corresponding to the current local dispersion value. Thus, each local dispersion value has a definite spatial coordinate, and the same value does not correspond to multiple spatial positions.

[0129] In constructing the concentration distribution feature matrix based on the local dispersion value, the spatial grid used in Example 3 to divide the airflow vortex region is used as the matrix position reference; the spatial grids are numbered sequentially according to the longitudinal position of the pipe gallery, and then the corresponding position of each grid in the matrix is ​​determined according to the horizontal and vertical positions; the three-dimensional sampling coordinates of the center node of each window are read to determine its spatial grid, and the corresponding local dispersion value is written into the matrix position corresponding to the spatial grid.

[0130] When the same spatial grid corresponds to only one local dispersion value, the value is directly written into the corresponding matrix cell; when the same spatial grid corresponds to two or more local dispersion values, all local dispersion values ​​corresponding to the spatial grid are added together and then divided by the actual number of local dispersion values, and the resulting arithmetic mean is used as the value of the matrix cell; spatial grids without a window center node are filled with the value 0; after all local dispersion values ​​are written into the space, a concentration distribution feature matrix is ​​formed.

[0131] Each non-zero matrix element in the concentration distribution feature matrix corresponds to a defined spatial grid within the airflow vortex region, and characterizes the degree of gas concentration dispersion within the range of five consecutive sampling nodes corresponding to that spatial grid.

[0132] Example 5:

[0133] S4 includes extracting the concentration distribution feature matrix to obtain the largest element of the matrix; comparing the largest element of the matrix with a preset variance threshold to obtain a threshold comparison result; if the threshold comparison result indicates that the largest element of the matrix is ​​greater than the preset variance threshold, a boundary anomaly marker is generated; and based on the boundary anomaly marker, it is determined that there is a fragment distribution phenomenon within the spatial boundary of the airflow vortex region.

[0134] In this embodiment, the concentration distribution feature matrix obtained in Example 4 is read; each effective matrix element in the matrix corresponds to a specific spatial grid within the airflow vortex region, and its value is the local dispersion value calculated according to 5 consecutive sampling nodes for the corresponding local concentration subsequence.

[0135] First, read the matrix elements one by one according to the row and column order, and take the first valid matrix element as the current maximum value. Then, compare the next matrix element with the current maximum value. If the next matrix element is greater than the current maximum value, replace the current maximum value with the next matrix element and record the corresponding spatial grid position of the matrix element. If the next matrix element is not greater than the current maximum value, keep the current maximum value unchanged. After comparing all valid matrix elements in sequence, the final value retained is the maximum element of the matrix, and its corresponding spatial grid position is saved at the same time.

[0136] As a preferred embodiment, the variance threshold is determined using background concentration data before the inspection task begins.

[0137] During the calibration process, the pipe gallery section corresponding to the actual inspection is selected. Without any abnormal gas sources to be detected, the pipe gallery inspection robot is made to move at the same speed as the formal inspection. The gas sensor continuously collects the background gas concentration according to the sampling period used in the formal inspection. The collected background data is processed in the same way as described in Examples 2 to 4 to complete timestamp alignment, spatial correspondence, vortex region sampling node screening, 5-point sliding window interception, and local dispersion calculation to obtain a set of background local dispersion values.

[0138] Furthermore, the sum of all the background local dispersion values ​​is added together and then divided by the number of background local dispersion values ​​to obtain the average background dispersion value.

[0139] Then, the difference between each local background dispersion value and the average background dispersion value is calculated, and the absolute value of the difference is taken. All absolute values ​​are added together and then divided by the number of local background dispersion values ​​to obtain the average absolute deviation of the background dispersion. The average background dispersion value is added to twice the average absolute deviation of the background dispersion, and the resulting value is set as the variance threshold and written into the parameter storage unit of the utility tunnel inspection robot. The variance threshold remains unchanged in the same inspection task.

[0140] Furthermore, after obtaining the maximum element of the matrix, the maximum element of the matrix is ​​directly compared with the variance threshold.

[0141] Specifically, when the largest element of the matrix is ​​greater than the variance threshold, the threshold comparison result is written as an abnormal state; when the largest element of the matrix is ​​equal to or less than the variance threshold, the threshold comparison result is written as a normal state. This comparison process only uses the largest element of the matrix formed in Example 4 and the pre-calibrated variance threshold, without re-correcting the concentration of the first gas concentration sequence or changing the spatial position corresponding to the matrix element.

[0142] When the threshold comparison result is abnormal, write the boundary anomaly identifier into the data record corresponding to the current airflow vortex region and set the status value of the boundary anomaly identifier to 1.

[0143] Simultaneously, the boundary number of the current airflow vortex region, the maximum element of the matrix, the spatial grid coordinates corresponding to the maximum element of the matrix, and the corresponding sampling time are associated and stored with the boundary anomaly identifier; when the threshold comparison result is in a normal state, the state value of the boundary anomaly identifier is set to 0, and it does not proceed to the subsequent concentration mutation difference extraction step.

[0144] Among them, the boundary anomaly marker is used to make the final determination of the debris distribution phenomenon within the airflow vortex region.

[0145] When the state value of the boundary anomaly indicator is 1, it indicates that there is at least one spatial grid in the current airflow vortex region, and the dispersion of its corresponding local concentration subsequence has exceeded the variance threshold determined based on the background concentration data. Therefore, the airflow vortex region is judged to have a fragmentation distribution phenomenon. When the state value of the boundary anomaly indicator is 0, it is determined that the current airflow vortex region does not meet the fragmentation distribution judgment condition.

[0146] After determining that there is a fragmentation distribution phenomenon, the spatial grid position corresponding to the largest element of the matrix and the spatial boundary of the current airflow vortex region are retained, and the first gas concentration sequence arranged according to the sampling time within the spatial boundary is input into S5 as described in Example 1 to extract the concentration abrupt difference between adjacent sampling times.

[0147] Example 6:

[0148] S5 includes obtaining the first gas concentration sequence corresponding to the spatial boundary of the airflow vortex region; performing a difference operation on the first gas concentration sequence to obtain a concentration gradient vector; performing feature extraction on the concentration gradient vector to obtain a difference feature set; determining whether the feature value of the difference feature set is greater than a preset mutation threshold; if the feature value is greater than the preset mutation threshold, then extracting the concentration mutation difference of the first gas concentration sequence at adjacent sampling times from the spatial boundary.

[0149] In this embodiment, the airflow vortex region with a boundary anomaly flag value of 1 in Embodiment 5 is read, and the spatial boundary coordinates of the airflow vortex region are retrieved.

[0150] Specifically, based on the three-dimensional sampling coordinates retained by each sampling node in the first gas concentration sequence, the sampling nodes are matched with the spatial boundary, and the sampling nodes located inside the spatial boundary and on the boundary surface are retained and arranged in order of sampling timestamp; the gas concentration values ​​corresponding to each sampling node are extracted sequentially to form the first gas concentration sequence corresponding to the current airflow vortex region; the first gas concentration sequence uses the concentration values, sampling times and three-dimensional coordinates after spatial mapping in Example 2, without further interpolation or numerical correction.

[0151] In the process of performing differential operations on the first gas concentration sequence, starting from the second sampling node in the time sequence, the gas concentration value of the current sampling node is subtracted from the gas concentration value of the previous sampling node to obtain a concentration difference value. Then, the gas concentration value of the third sampling node is subtracted from the gas concentration value of the second sampling node to obtain the next concentration difference value. This process is repeated until the last sampling node. All concentration difference values ​​are arranged in the order of their corresponding sampling times to form a concentration gradient vector. Each concentration difference value is associated with the two sampling nodes before and after it, thus preserving the sampling time and spatial location of the concentration change.

[0152] When extracting features from the concentration gradient vector, the absolute value of each concentration difference is taken, and the resulting value is used as the concentration change feature value between the corresponding adjacent sampling nodes. After taking the absolute value, only the concentration change amplitude is retained, and the direction of concentration increase or decrease is not used as the basis for judging abrupt changes. Each concentration change feature value is associated with the corresponding previous sampling time, the next sampling time, the coordinates of the previous sampling node, and the coordinates of the next sampling node, and arranged in chronological order to obtain the difference feature set.

[0153] As a preferred embodiment, the preset mutation threshold is calibrated using background response data and step response data of the gas sensor.

[0154] During calibration, the gas sensor is first placed in a standard background gas environment with a constant concentration, and 30 gas concentration values ​​are continuously collected according to the sampling cycle used in the actual inspection. Using the aforementioned differential method, the concentration value of the previous sampling time is subtracted from the concentration value of the next sampling time in turn, and the absolute value of each differential result is taken to obtain a background difference sequence. All values ​​in the background difference sequence are added together and divided by the number of background differences to obtain the average background difference. Then, the absolute value of the difference between each background difference and the average background difference is calculated, and all absolute values ​​are added together and divided by the number of background differences to obtain the average absolute deviation of the background difference.

[0155] Furthermore, after completing the background response calibration, the gas sensor is directly transferred from the standard background gas environment to a standard target gas environment with a known concentration, and 30 gas concentration values ​​are continuously collected at the same sampling period as the actual inspection.

[0156] Specifically, the concentration difference between adjacent sampling times is calculated sequentially and its absolute value is taken to obtain a step response difference sequence. The difference data formed from the first sampling time to the point before the sensor output reaches the standard target gas concentration stability value is extracted from the step response difference sequence, and the minimum difference is determined. The average background difference is added to three times the average absolute deviation of the background difference, and the result is compared with the minimum step response difference. The smaller of the two values ​​is taken as the mutation threshold. This ensures that the mutation threshold is higher than the difference range corresponding to the natural fluctuation of the background concentration, while falling within the concentration change range that the gas sensor can produce when it has a step response. After the mutation threshold is determined, it is written into the parameter storage unit of the inspection robot and remains unchanged in the same inspection task.

[0157] Furthermore, after completing the difference feature set, the feature values ​​are read one by one according to the sampling time sequence and compared with the preset mutation threshold respectively.

[0158] When a certain feature value is greater than the preset mutation threshold, the two sampling nodes before and after generating the feature value are determined as a mutation sampling pair; when the feature value is equal to or less than the preset mutation threshold, the corresponding sampling node is not included in the mutation sampling pair; after comparing all feature values ​​in the difference feature set, all mutation sampling pairs that meet the mutation conditions in the current airflow vortex region are obtained.

[0159] In addition, for each mutation sampling pair, the original gas concentration values ​​corresponding to the previous sampling time and the next sampling time are reread from the first gas concentration sequence.

[0160] The gas concentration value at the previous sampling time is subtracted from the gas concentration value at the next sampling time, and the absolute value of the calculation result is taken. The resulting value is determined as the concentration mutation difference corresponding to the mutation sampling pair. The original concentration value in the first gas concentration sequence is reread and calculated so that the concentration mutation difference comes directly from the first gas concentration sequence in Example 1, without using other concentration values ​​formed during the intermediate processing.

[0161] Each concentration abrupt change difference is associated with and stored in relation to its corresponding two previous and two subsequent sampling times, the three-dimensional coordinates of the two previous and two subsequent sampling nodes, and the boundary number of the current airflow vortex region, and then output in the order of sampling time.

[0162] Example 7:

[0163] S6 includes: if the concentration change difference exceeds a preset change threshold, extracting the corresponding adjacent sampling time to construct a delayed response interval; truncating the first gas concentration sequence according to the delayed response interval to obtain a concentration decay feature sequence; obtaining the tail duration and difference fluctuation amplitude according to the concentration decay feature sequence to construct a tail determination condition; if the concentration decay feature sequence meets the tail determination condition, determining that the first gas concentration sequence has a numerical tail caused by response delay.

[0164] In this embodiment, the concentration mutation difference obtained in Example 6 and its corresponding two adjacent sampling times are read, and the mutation threshold determined in Example 6 is used for judgment.

[0165] When the concentration change difference is greater than the change threshold, the previous sampling time and the next sampling time are used as the time reference for this delayed response analysis. The previous sampling time corresponds to the concentration state before the concentration change occurs, and the next sampling time corresponds to the first concentration state after the concentration change occurs. The change threshold continues to be the value obtained by background response calibration and step response calibration in Example 6.

[0166] The length of the delay response interval is determined based on the step response calibration data of the gas sensor.

[0167] Before inspection, the gas sensor is placed in the first standard concentration environment. After the absolute value of the difference between five consecutive sampling values ​​and the first standard concentration value is not greater than the sensor's allowable indication error, the gas sensor is transferred to the second standard concentration environment within one sampling cycle. The concentration value is continuously recorded starting from the first sampling time after entering the second standard concentration environment. When the absolute value of the difference between five consecutive sampling values ​​and the second standard concentration value is again not greater than the allowable indication error, the first sampling time among these five consecutive sampling values ​​is taken as the stable time. The time difference between the first sampling time and the stable time in the second standard concentration environment is recorded as the duration of one response. Ten calibrations are completed under the same conditions. The durations of the ten responses are added together and divided by 10 to obtain the delayed response calibration duration.

[0168] During actual inspection, when a concentration change exceeds the change threshold, the first sampling time corresponding to the concentration change is taken as the start time of the delay response interval, and the time obtained by extending the delay response calibration time backward from the next sampling time is taken as the end time of the delay response interval. If the delay response intervals formed by two consecutive concentration change differences overlap, the earliest start time and the latest end time are retained, and the overlapping part is merged into one delay response interval to avoid repeated processing of the same segment of sensor response data.

[0169] Furthermore, based on the determined delay response interval, all concentration values ​​with timestamps within that interval are extracted from the first gas concentration sequence, while maintaining the original sampling order of each concentration value.

[0170] Specifically, five consecutive gas concentration values ​​are taken from the end of the delayed response interval. The sum of these five gas concentration values ​​is then divided by 5 to obtain the stable reference concentration corresponding to that interval. Subsequently, the stable reference concentration is subtracted from each gas concentration value within the interval, and the absolute value of the difference is taken. The absolute values ​​are arranged in order of sampling time to form a concentration decay characteristic sequence. After adopting this processing method, regardless of whether the actual concentration relationship before and after the abrupt change is from large to small or from small to large, the concentration decay characteristic sequence represents the remaining deviation between the sensor output value and the stable concentration after the abrupt change.

[0171] As a preferred embodiment, the duration of the tail is determined based on the concentration decay characteristic sequence and the stability determination threshold.

[0172] Specifically, the stability threshold is determined using the permissible indication error obtained from the gas sensor's calibration or verification. Starting from the next sampling time, the concentration decay characteristic sequence is read point by point. When five consecutive concentration deviations are all below the permissible indication error for the first time, the sampling time corresponding to the first of these five concentration deviations is determined as the tailing end time. The time difference between the next sampling time and the tailing end time is taken as the tailing duration. If no five consecutive concentration deviations within the delay response interval meet the above conditions, the end time of the delay response interval is taken as the tailing end time.

[0173] Furthermore, the amplitude of the difference fluctuation is calculated from adjacent data in the concentration decay characteristic sequence.

[0174] The concentration deviation is subtracted from the previous concentration deviation in turn, and the absolute value of the difference is taken to form an attenuation difference sequence. All values ​​in the attenuation difference sequence are compared, the maximum and minimum values ​​are extracted, and the minimum value is subtracted from the maximum value. The result is used as the difference fluctuation amplitude of the delay response interval. This value is used to characterize the degree of dispersion of the change amplitude between adjacent sampling times during the process of the sensor output converging towards the stable reference concentration.

[0175] As a preferred embodiment, the duration threshold in the tail determination condition is determined using the aforementioned 10 step response calibration data.

[0176] Specifically, the duration of each of the 10 responses was recorded, and the shortest response duration was used as the duration threshold. The difference fluctuation threshold was determined using gas sensor sampling data under a stable background environment: 30 concentration values ​​were continuously collected under the condition that the background gas concentration remained constant, and the absolute value of the difference between adjacent concentration values ​​was calculated sequentially. All absolute values ​​were added together and divided by the number of absolute values ​​to obtain the average background difference. The absolute value of the difference between each absolute value and the average background difference was then calculated, and the values ​​were added together and divided by the number to obtain the average absolute deviation of the background difference. The average background difference was added to three times the average absolute deviation of the background difference, and the resulting value was used as the difference fluctuation threshold.

[0177] Therefore, the tail duration is not less than the duration threshold and the difference fluctuation amplitude is greater than the difference fluctuation threshold as the tail determination conditions; the current concentration decay feature sequence is judged, and when the above two conditions are met at the same time, it is determined that the corresponding first gas concentration sequence has a numerical tail caused by the gas sensor response delay; if either condition is not met, no numerical tail identifier is generated.

[0178] Example 8:

[0179] S7 includes obtaining a first gas concentration sequence with numerical tails, constructing a state transition matrix and an observation matrix based on the first gas concentration sequence, generating a priori state estimation sequence using the state transition matrix and the observation matrix, determining a Kalman gain matrix based on the observation noise covariance if the system noise covariance calculated from the priori state estimation sequence is greater than a preset threshold, and denoising the priori state estimation sequence using the Kalman gain matrix to obtain a second gas concentration sequence.

[0180] In this embodiment, the first gas concentration sequence with numerical tails determined in Example 7 is obtained and arranged in order of sampling timestamps.

[0181] Each sampling node in the first gas concentration sequence retains the corresponding gas concentration value, sampling time, and three-dimensional sampling coordinates. The Kalman filter processes the numerical tailing interval. Its input data is the continuous gas concentration value with numerical tailing, and its output data is the second gas concentration sequence that maintains a one-to-one correspondence with the original sampling time and the original three-dimensional sampling coordinates.

[0182] When constructing the filter state, two state components are set for each sampling time. The first state component represents the gas concentration estimate corresponding to that sampling time, and the second state component represents the gas concentration change per unit time.

[0183] The gas concentration estimate at the first sampling time is the first gas concentration value in the first gas concentration sequence; the first concentration change is obtained by subtracting the first gas concentration value from the second gas concentration value and then dividing by the time interval between the two sampling times; thus, subsequent state predictions can simultaneously utilize the current concentration level and the concentration change trend over time.

[0184] The state transition matrix is ​​set to 2 rows and 2 columns, and its connection relationship is determined according to the concentration continuous change model.

[0185] Specifically, the first line is used to predict the gas concentration at the current sampling time from the gas concentration and concentration change at the previous sampling time, wherein the gas concentration at the previous sampling time remains at its original value, and the concentration change at the previous sampling time is multiplied by the time interval between adjacent sampling times and then added to the gas concentration; the second line is used to transmit the concentration change, so that the concentration change at the previous sampling time is used as the predicted value of the concentration change at the current sampling time; according to this connection relationship, the state transition matrix can convert the previous sampling state into the predicted state at the current sampling time.

[0186] In addition, the observation matrix is ​​set to 1 row and 2 columns, and only a direct connection is established between the actual measured values ​​of the gas sensor and the gas concentration state.

[0187] In the observation matrix, the position corresponding to the gas concentration state is set to 1, and the position corresponding to the concentration change is set to 0. This ensures that the actual concentration value collected by the gas sensor is used to correct the predicted concentration, while the concentration change is updated based on the continuous sampling status and the filtering correction results. This setting is consistent with the data format of the gas sensor's actual output, which only contains the concentration measurement value.

[0188] When generating the prior state estimation sequence using the state transition matrix, the processing starts from the second sampling node of the first gas concentration sequence.

[0189] Specifically, the gas concentration estimate and concentration change after filtering correction at the previous sampling time are read. The concentration change at the previous sampling time is multiplied by the time interval between the current sampling time and the previous sampling time, and then added to the gas concentration estimate at the previous sampling time. The result is used as the prior concentration estimate at the current sampling time. The concentration change at the previous sampling time is used as the prior concentration change at the current sampling time. All sampling nodes within the numerical tail interval are processed in sequence to form the prior state estimation sequence.

[0190] In this embodiment, the system noise covariance is used to describe the deviation between the actual gas concentration change and the above-mentioned state transition prediction; its parameters are calibrated before the inspection according to the same gas sampling period as the actual inspection.

[0191] During calibration, 30 concentration values ​​are continuously collected in a standard gas environment where the gas concentration remains constant. The predicted concentration values ​​at adjacent sampling times are calculated according to the aforementioned state transition relationship. The actual concentration value at each sampling time is then subtracted from the corresponding predicted concentration value to obtain a prediction residual sequence formed by the 30 sampling data. All prediction residuals are summed and divided by the number of prediction residuals to obtain the average residual. The difference between each prediction residual and the average residual is calculated, and the obtained difference is squared. All squared results are summed and divided by the number of prediction residuals. The result is used as the system noise variance corresponding to the concentration state.

[0192] As a preferred embodiment, the system noise parameters corresponding to the concentration change are determined in the same manner.

[0193] Specifically, the concentration change between adjacent sampling times is first calculated, and then the difference between two consecutive concentration changes is calculated, thus forming a residual sequence of changes. After obtaining the average value of the residual sequence of changes, the square of the difference between each residual and the average value is calculated. All squared results are added together and divided by the number of residuals of changes to obtain the system noise variance corresponding to the concentration change. At the same time, the concentration prediction residual and the residual of change at the same sampling location are subtracted from their respective average values ​​and multiplied. All products are then added together and divided by the number of valid data to obtain the cross covariance between the two states. The concentration state variance, concentration change variance, and cross covariance are written into a 2x2 matrix according to the correspondence between the two state components to form the system noise covariance.

[0194] To perform a comparison between the system noise covariance and a preset threshold, the variance value directly corresponding to the gas concentration state is extracted from the system noise covariance as the covariance judgment value.

[0195] Specifically, the preset threshold is determined using 10 sets of calibration data under a stable gas environment. For each set, 30 concentration values ​​are collected and the system noise variance for one concentration state is calculated using the method described above. The 10 system noise variances are added together and divided by 10 to obtain the average system noise variance. The absolute value of the difference between each system noise variance and the average system noise variance is then calculated. The absolute values ​​are added together and divided by 10 to obtain the mean absolute deviation of the system noise variance. The average system noise variance is added to three times the mean absolute deviation of the system noise variance, and the resulting value is set as the preset threshold. This threshold is written into the parameter storage unit before the start of this inspection task and remains unchanged during the inspection.

[0196] As a preferred embodiment, the observation noise covariance is determined based on the repeatability error of the gas sensor itself.

[0197] Specifically, the gas sensor is placed in a standard gas environment with a constant concentration, and 30 gas concentration values ​​are continuously collected according to the sampling cycle used in the formal inspection. The 30 concentration values ​​are added together and divided by 30 to obtain the calibrated average concentration. The calibrated average concentration is subtracted from each concentration value, the difference is squared, and the sum of all squared results is divided by 30. The obtained value is used as the observation noise covariance. The observation noise covariance corresponds to the gas sensor used and remains unchanged in the same inspection task.

[0198] When the system noise covariance at the current sampling time is greater than the preset threshold, the Kalman gain calculation is initiated.

[0199] First, read the updated state error covariance from the previous sampling time. Apply the state transition matrix to this state error covariance to propagate the state error from the previous time to the current sampling time. Then, add the current system noise covariance to obtain the prior state error covariance for the current sampling time. This prior state error covariance is also 2 rows and 2 columns, which includes the predicted concentration error, the predicted concentration change error, and the correlation error between the two.

[0200] Subsequently, the error component corresponding to the gas concentration measurement is extracted from the prior state error covariance using the observation matrix, and then added to the observation noise covariance to obtain the total error corresponding to the current gas concentration observation.

[0201] Then, the two error components related to gas concentration observation in the prior state error covariance are divided by the total error amount to form a 2-row, 1-column Kalman gain matrix, where the first gain is used to correct the gas concentration state and the second gain is used to correct the concentration change. Thus, the Kalman gain is controlled by both the state prediction error and the sensor observation error.

[0202] When using the Kalman gain matrix for denoising, the actual concentration value in the first gas concentration sequence at the current sampling time is read, and the prior concentration estimate at the current sampling time is subtracted from the actual concentration value to obtain the observation residual.

[0203] The observation residuals are multiplied by the two gains in the Kalman gain matrix to obtain the concentration correction and the concentration change correction. The concentration correction is added to the prior concentration estimate to obtain the posterior concentration estimate at the current sampling time. The concentration change correction is added to the prior concentration change to obtain the corrected concentration change at the current sampling time.

[0204] Furthermore, after completing the state correction, the error portion that has been corrected by actual observations in the current prior state error covariance is reduced according to the Kalman gain matrix and the observation matrix to obtain the posterior state error covariance at the current sampling time. The posterior concentration estimate, the corrected concentration change, and the posterior state error covariance are used together as the input for the state prediction at the next sampling time, thereby continuously completing the state estimation within the interval with numerical tails in the order of sampling time.

[0205] For sampling times where the system noise covariance judgment value is not greater than the preset threshold, the prior concentration estimate obtained according to the state transition matrix is ​​maintained, and gain correction based on abnormal system noise is not initiated; when the judgment value is greater than the preset threshold again, the Kalman gain matrix of the current sampling time is recalculated according to the aforementioned process and observation correction is implemented.

[0206] This makes the filtering process consistent with the system noise covariance determination relationship defined in Example 8.

[0207] After all sampling nodes have been processed, the gas concentration estimates corrected by Kalman filtering at each sampling time are extracted sequentially and arranged in the original sampling time order of the first gas concentration sequence. At the same time, the three-dimensional spatial coordinates of the corresponding sampling nodes are inherited to form the second gas concentration sequence.

[0208] Example 9:

[0209] S8 includes obtaining a second gas concentration sequence, reconstructing the second gas concentration sequence to obtain a concentration distribution matrix, and extracting the spatial gradient features of the concentration distribution matrix; determining the airflow vortex region based on the spatial gradient features, obtaining the spatial boundary of the airflow vortex region and calculating the centroid position of the boundary; calculating the geometric center distance of the spatial boundary through the centroid position of the boundary; and constructing the source feature vector based on the spatial gradient features and the geometric center distance.

[0210] In this embodiment, the second gas concentration sequence obtained in Example 8 is read; each gas concentration value in the second gas concentration sequence retains the corresponding sampling timestamp and three-dimensional sampling coordinates; the spatial division of the concentration distribution matrix follows the spatial grid used in Example 3 when determining the airflow vortex region, so that the second gas concentration data and the already determined airflow vortex region spatial boundary adopt the same pipe gallery coordinate reference.

[0211] During the concentration distribution reconstruction, the sampling nodes in the second gas concentration sequence are read one by one, the spatial grid where the sampling node is located is determined based on the three coordinate values ​​of the sampling node, and the corresponding second gas concentration value is written into the spatial grid.

[0212] Specifically, when there is only one valid sampling node in a spatial grid, the gas concentration value of that sampling node is used as the concentration value of the corresponding matrix unit; when there are two or more valid sampling nodes in a spatial grid, all gas concentration values ​​in the spatial grid are added together and then divided by the actual number of valid sampling nodes in the spatial grid, and the resulting arithmetic mean is used as the concentration value of the corresponding matrix unit; all concentration data are written according to the arrangement order of the spatial grid in the longitudinal, transverse and vertical directions of the pipe gallery to form a concentration distribution matrix.

[0213] Furthermore, after the concentration distribution matrix is ​​formed, each spatial grid with an effective concentration value is used as the grid to be calculated, and the concentration gradients in the longitudinal, transverse and vertical directions of the pipe gallery are calculated respectively.

[0214] When there are two valid adjacent grids in the current direction, the concentration gradient of the current grid in that direction is obtained by subtracting the gas concentration value of the front grid from the gas concentration value of the rear grid and then dividing by the actual distance between the centers of the two grids. When there is only one valid adjacent grid in the current direction, the concentration gradient of the current grid in that direction is obtained by subtracting the gas concentration value of the current grid from the gas concentration value of the adjacent grid and then dividing by the actual distance between the centers of the two grids. When there are no valid adjacent grids in the current direction, that direction is not included in the gradient calculation of this sampling node.

[0215] After obtaining the concentration gradients in the three coordinate directions, the three concentration gradients are squared respectively, the three squared results are added together, and the square root of the sum is taken to obtain the gradient magnitude of the current spatial grid.

[0216] The spatial gradient features of the current spatial grid are composed of the longitudinal, transverse, and vertical concentration gradients and their magnitudes. All valid spatial grids in the concentration distribution matrix are processed according to the same calculation process to obtain the spatial gradient features corresponding to each sampling spatial location.

[0217] In determining the airflow vortex region based on spatial gradient characteristics, the airflow vortex region and its spatial boundary obtained in Example 3 are used as candidate regions without changing the boundary coordinates determined by the obstacle occlusion characteristics. The effective spatial grid inside each candidate airflow vortex region is read, and the gradient magnitude of each grid is compared with the gradient determination threshold. Grids with gradient magnitudes greater than the gradient determination threshold are marked as effective grids of the concentration gradient, and then the effective grids of the concentration gradient are checked for connectivity according to the surface adjacency relationship between the grids.

[0218] As a preferred embodiment, the gradient determination threshold is determined using background concentration data under conditions of no gas anomaly source before the formal inspection.

[0219] Specifically, under the same sampling period and robot movement speed as the formal inspection, 30 background gas concentration samples are continuously obtained. After Kalman filtering as described in Example 8, the gradient magnitudes of all effective grids are obtained according to the aforementioned concentration distribution reconstruction and spatial gradient calculation process. The average background gradient is obtained by summing all gradient magnitudes and dividing by the number of gradient magnitudes. The absolute value of the difference between each gradient magnitude and the average background gradient is then calculated. The absolute values ​​are summed and divided by the number of gradient magnitudes to obtain the average absolute deviation of the background gradient. The average background gradient is added to three times the average absolute deviation of the background gradient, and the resulting value is used as the gradient determination threshold, which remains unchanged in the same inspection task.

[0220] When a continuous region consisting of an effective concentration gradient grid exists within a certain airflow vortex region determined in Example 3, and this continuous region contains the spatial grid corresponding to the largest element of the matrix recorded in Example 5, the airflow vortex region is determined as the airflow vortex region corresponding to the current second gas concentration sequence, and the original spatial boundary coordinates of the region are read.

[0221] Therefore, the spatial gradient feature is used to confirm the correspondence between the current concentration data and the existing airflow vortex region, without generating a vortex spatial boundary different from that in Example 3.

[0222] As a preferred embodiment, the location of the boundary centroid is determined based on the actual space enclosed by the spatial boundary of the airflow vortex region.

[0223] The boundary is filled with the spatial grid used in Example 3, retaining all valid grids whose centers are located inside or on the boundary surface. Since each spatial grid has the same volume, the coordinates of all valid grid centers in the longitudinal direction of the pipe gallery are added together and divided by the number of valid grids to obtain the longitudinal coordinates of the boundary centroid. The transverse and vertical coordinates are calculated using the same method. The three coordinates together determine the location of the boundary centroid of the airflow vortex region. The centroid is calculated using a uniform-volume grid inside the boundary, so that the location of the centroid is determined by the actual volume range enclosed by the spatial boundary and is not affected by changes in the sampling density of the boundary points.

[0224] Furthermore, after obtaining the location of the boundary centroid, the geometric center distance from each effective sampling node in the second gas concentration sequence to the boundary centroid is calculated.

[0225] For any sampling node, calculate the coordinate differences between the sampling node and the boundary centroid in the longitudinal, lateral, and vertical directions respectively. Then, square the three coordinate differences respectively, add the three squared results, and take the square root of the sum. The resulting value is the geometric center distance corresponding to the sampling node. Process all valid sampling nodes in the current airflow vortex region in the same way so that each sampling node obtains the geometric center distance corresponding to its three-dimensional spatial position.

[0226] When constructing the source feature vector, a valid sampling node is used as a feature generation unit. The vertical concentration gradient, horizontal concentration gradient, vertical concentration gradient, gradient magnitude, and geometric center distance from the sampling node to the boundary centroid are written in a fixed order.

[0227] The above five values ​​constitute the source feature vector of the current sampling node; the three-dimensional coordinates of the sampling node, the sampling timestamp, and the boundary number of the corresponding airflow vortex region are stored as the position index of the source feature vector, and are not used as replacement data for the aforementioned five feature values.

[0228] After processing all valid sampling nodes according to the sampling order of the second gas concentration sequence, a source feature vector corresponding to each sampling node is obtained. The source feature vector simultaneously includes the direction and magnitude of gas concentration change in the pipe gallery space, as well as the positional relationship of the sampling node relative to the centroid of the boundary of the airflow vortex region.

[0229] Example 10:

[0230] S9 includes using a support vector machine to perform high-dimensional spatial mapping on the source feature vectors to obtain a high-dimensional feature space; extracting a feature distribution matrix based on the high-dimensional feature space; determining the decision boundary hyperplane through the feature distribution matrix and performing classification processing to obtain anomaly classification labels; if the anomaly classification label indicates an anomaly source, extracting the spatial coordinate grid corresponding to the anomaly source; parsing the spatial coordinate grid to determine the true source coordinates; and using the true source coordinates as the environmental parameter monitoring results of the utility tunnel inspection robot.

[0231] In this embodiment, the input to the support vector machine is the source feature vector obtained in Embodiment Nine.

[0232] Specifically, each source feature vector sequentially includes the longitudinal concentration gradient, lateral concentration gradient, vertical concentration gradient, gradient magnitude, and geometric center distance from the sampling node to the centroid of the boundary of the airflow vortex region, while retaining the spatial grid number corresponding to the source feature vector; the support vector machine is trained before the inspection task is executed, and only feature mapping, classification, and coordinate analysis are performed during the actual inspection process.

[0233] The training samples were obtained through calibration experiments at known anomaly source locations.

[0234] Specifically, using the same spatial grid division method as in Example 9, the pipe gallery to be monitored is divided into spatial coordinate grids of fixed size. Standard gas release points are sequentially set at the center of each spatial grid to be calibrated, so that the inspection robot completes data collection according to the actual inspection speed and sampling cycle, and sequentially executes the processing steps described in S1 to S8 to obtain the source feature vectors corresponding to each sampling node. The source feature vectors corresponding to the spatial grid where the standard gas release point is located are marked as abnormal source categories, and the source feature vectors corresponding to the other spatial grids in the same airflow vortex area are marked as non-abnormal source categories, thereby forming support vector machine training samples.

[0235] In addition, standardization was performed on the five features before training.

[0236] Taking one feature as an example, the feature values ​​of all training samples are summed and divided by the number of training samples to obtain the average value of the feature. The difference between each feature value and the average value is calculated, and each difference is squared. All squared results are summed and divided by the number of training samples. The square root of the result is then taken to obtain the standard deviation of the feature. Each original feature value is subtracted from its corresponding average value and then divided by the standard deviation of the feature to obtain the standardized feature value. All five features are processed in this way. The five sets of average values ​​and standard deviations obtained from training are saved. During actual inspection, the same parameters are used to process the source feature vector to be classified.

[0237] In this embodiment, the support vector machine uses a radial basis kernel function for high-dimensional space mapping.

[0238] Specifically, for any two standardized source feature vectors, the differences between the five corresponding features are calculated sequentially. The five differences are squared and summed to obtain the squared feature distance between the two source feature vectors. This value is then multiplied by the kernel width parameter and the negative value is taken. The resulting value is then exponentially calculated to obtain the kernel similarity between the two source feature vectors. According to the order of the training samples, the kernel similarity is calculated for any two training samples one by one, and the calculation results are written into the corresponding row and column positions to form a feature distribution matrix. Each element in this matrix represents the degree of similarity between the corresponding two training samples after mapping to the high-dimensional feature space.

[0239] In this embodiment, the kernel width parameter and penalty parameter are determined by 5-fold cross-validation.

[0240] Specifically, the training samples for the anomaly source category and the training samples for the non-anomaly source category are divided into 5 groups according to their sample numbers. These groups are then combined to form 5 sets of validation data. Four sets are used as training data each time, and the remaining set is used as validation data. This process is repeated 5 times to ensure each set of data completes one validation. The penalty parameter is set to 1, 10, 100, and 1000, and the kernel width parameter is set to 0.01, 0.1, 1, and 10, respectively. Five-fold cross-validation is performed on each of the 16 parameter combinations. After each parameter validation, the number of samples whose classification results match the known categories is counted, and this number is divided by the total number of validation samples to obtain the classification accuracy of that parameter combination. The parameter combination with the highest classification accuracy is selected as the training parameter. When there are samples with the same classification accuracy, the combination with the smaller penalty parameter value is preferred. If the penalty parameters are still the same, the combination with the smaller kernel width parameter value is used.

[0241] Furthermore, after determining the training parameters, a classification model is obtained based on the category labels and feature distribution matrix of the training samples.

[0242] During training, the objective is to maximize the margin between the anomaly and non-anomaly categories, and a penalty parameter limits the extent to which samples cross the classification boundary. After calculation, each training sample is assigned a corresponding weight coefficient, with training samples whose weight coefficients are not zero being identified as support vectors. The weight coefficients, class labels, and kernel similarities between each support vector and other support vectors are substituted into the classification relation to determine the classification bias. The support vectors, weight coefficients, and classification bias together determine the decision boundary hyperplane in the high-dimensional feature space. After training, the mean, standard deviation, kernel width parameter, support vectors, support vector weight coefficients, and classification bias of the five features are saved as classification model parameters for the actual inspection phase.

[0243] During the actual inspection process, the source feature vectors formed in Example 9 were read one by one.

[0244] First, the mean and standard deviation saved during the training phase are used to perform the same standardization process on the five features. Then, the kernel similarity between the current source feature vector and all support vectors is calculated. Each kernel similarity is multiplied by the weight coefficient and class label of the corresponding support vector. All products are then added together and the classification bias is added. The resulting value is used as the classification judgment value of the current source feature vector relative to the decision boundary hyperplane.

[0245] When the classification judgment value is greater than 0, the abnormal classification label corresponding to the current source feature vector is determined as an abnormal source; when the classification judgment value is equal to or less than 0, the abnormal classification label is determined as a non-abnormal source.

[0246] Therefore, the anomaly classification label is directly obtained by the support vector machine after training the five source features formed in Example 9, without introducing additional judgment conditions unrelated to the source features.

[0247] When the anomaly classification label indicates an anomaly source, the spatial grid number of the source feature vector saved in Example 9 is read, and the grid boundary coordinates in the unified coordinate system of the utility tunnel are queried according to the spatial grid number to obtain the spatial coordinate grid corresponding to the anomaly source.

[0248] If only one spatial grid in the same airflow vortex region is identified as an anomaly source, that spatial grid is directly identified as the grid where the anomaly source is located.

[0249] In addition, when there are more than two anomalous source spatial grids in the same airflow vortex region, the classification judgment value corresponding to each spatial grid is read. The values ​​are compared from largest to smallest, and the spatial grid with the largest classification judgment value is determined as the grid where the anomalous source is located. The classification judgment value is used to reflect the classification position of the current source feature vector relative to the decision boundary hyperplane. Therefore, this screening process continues to use the classification result of the same support vector machine.

[0250] Furthermore, after determining the spatial coordinate grid where the anomaly source is located, the two boundary coordinates of the grid along the longitudinal direction of the pipe gallery are read respectively. The two coordinate values ​​are added together and divided by 2 to obtain the longitudinal center coordinates. The transverse center coordinates and vertical center coordinates are calculated in the same way.

[0251] The center position of this spatial coordinate grid is composed of the center coordinates in three directions, and this center position is determined as the true source coordinate of the gas environment parameter anomaly.

[0252] Finally, the true source coordinates are associated with the corresponding airflow vortex region number, spatial grid number, anomaly classification label, detection time, and second gas concentration data, and output as the environmental parameter monitoring results of the utility tunnel inspection robot. Thus, the source feature vector formed in Example 9, after high-dimensional spatial mapping, decision boundary classification, and spatial grid backtracking, can obtain the true source coordinates directly corresponding to the unified coordinate system of the utility tunnel.

[0253] Example 11: An environmental parameter monitoring system for a utility tunnel inspection robot, used to implement the steps of the aforementioned method for monitoring the environmental parameters of a utility tunnel inspection robot, the system comprising: The data acquisition module acquires the spatial point cloud and the first gas concentration sequence of the utility tunnel site, which are synchronously collected by the lidar and gas sensor carried by the utility tunnel inspection robot during the utility tunnel inspection process and establish spatial correspondence according to the same sampling time. The vortex region determination module extracts obstacle and occlusion features within the utility tunnel based on the spatial point cloud at the site, and determines the spatial boundary of the airflow vortex region. The concentration feature calculation module performs sliding window variance calculation on the first gas concentration sequence within the spatial boundary of the airflow vortex region to obtain the concentration distribution feature matrix. The debris distribution determination module determines that if the largest element in the concentration distribution feature matrix is ​​greater than the preset variance threshold, then the debris distribution phenomenon exists within the spatial boundary of the airflow vortex region. The concentration mutation extraction module extracts the concentration mutation difference of the first gas concentration sequence at adjacent sampling times from the spatial boundary of the airflow vortex region where fragmentation exists. The numerical tailing determination module determines that if the concentration change difference exceeds the preset change threshold, the first gas concentration sequence has numerical tailing caused by response delay. The concentration filtering module uses the Kalman filter algorithm to perform state estimation and noise reduction on the first gas concentration sequence with numerical tails to obtain the second gas concentration sequence. The source feature construction module constructs a source feature vector by extracting the spatial gradient features of the second gas concentration sequence and the geometric center distance of the spatial boundary of the airflow vortex region. The anomaly source localization module uses a support vector machine to perform high-dimensional spatial mapping and classification of the source feature vectors, determine the true source coordinates of the gas environment parameter anomalies in the utility tunnel, and use the true source coordinates as the environmental parameter monitoring results of the utility tunnel inspection robot.

[0254] In this embodiment, the data acquisition module receives the spatial point cloud of the pipe gallery collected by the lidar and the gas concentration data collected by the gas sensor, aligns the two according to the sampling time, and then establishes the correspondence between the gas concentration value and the spatial position of the pipe gallery by combining the spatial pose of the inspection robot, thus forming the spatial point cloud of the pipe gallery and the first gas concentration sequence.

[0255] Furthermore, the vortex region determination module reads the spatial point cloud of the pipe gallery, extracts the surface normal, curvature, and obstacle contours of the point cloud, determines the velocity change region based on the occlusion relationship of the obstacles on the flow section of the pipe gallery, and forms the spatial boundary of the airflow vortex region accordingly; the concentration feature calculation module filters the first gas concentration sequence within the spatial boundary, continuously extracts concentration data using a sliding window, calculates the variance of the concentration values ​​in each window, and arranges them according to their corresponding spatial positions to form a concentration distribution feature matrix.

[0256] The fragment distribution determination module extracts the maximum element in the concentration distribution feature matrix and compares it with a pre-calibrated variance threshold. When the maximum element is greater than the variance threshold, the current airflow vortex region is determined to have fragment distribution. The concentration mutation extraction module then calculates the absolute value of the difference between adjacent concentration values ​​according to the sampling time sequence and extracts the concentration mutation difference that meets the mutation screening conditions.

[0257] The numerical tailing determination module compares the concentration mutation difference with the preset mutation threshold, and combines the decay duration of continuous concentration data after the mutation and the fluctuation state of the difference to determine whether the gas sensor has a response delay; when the tailing determination condition is met, it is determined that there is numerical tailing in the corresponding first gas concentration sequence.

[0258] The concentration filtering module establishes the concentration state and concentration change state for the first gas concentration sequence with numerical tails, and uses Kalman filtering to perform state prediction, observation correction and error update to obtain the second gas concentration sequence corresponding to the original sampling time and spatial coordinates.

[0259] The source feature construction module calculates the longitudinal, lateral and vertical concentration gradients based on the distribution of the second gas concentration sequence in the pipe gallery space, and calculates the spatial distance from each sampling node to the geometric center of the airflow vortex region. The concentration gradient, gradient amplitude and geometric center distance are combined in a fixed order to form the source feature vector.

[0260] The anomaly source localization module inputs the source feature vector into the trained support vector machine, completes the high-dimensional space mapping through the kernel function, and outputs anomaly classification labels based on the decision boundary. When the classification result indicates an anomaly source, the corresponding grid coordinates are read according to the spatial grid number associated with the source feature vector, and the true source coordinates of the gas environment parameter anomaly source are determined by the center coordinates of the spatial grid, which serve as the environmental parameter monitoring results of the pipeline inspection robot.

[0261] Example 12: A computer-readable storage medium storing computer instructions thereon, which execute the environmental parameter monitoring method of the utility tunnel inspection robot when the computer instructions are executed.

[0262] In this embodiment, a computer-readable storage medium is used to store computer instructions for implementing the environmental parameter monitoring method of the utility tunnel inspection robot.

[0263] The computer-readable storage medium includes any one of solid-state memory, flash memory, read-only memory, or non-volatile memory.

[0264] When the processor reads and executes the computer instructions, it performs data acquisition, airflow vortex region determination, concentration feature calculation, debris distribution determination, concentration mutation extraction, numerical tail determination, Kalman filtering, source feature construction, and anomaly source localization processes according to a pre-set program.

[0265] During execution, the processor receives data collected by the lidar and gas sensors, completes the spatiotemporal correspondence between point cloud and gas concentration data according to computer instructions, and sequentially calls the corresponding data processing programs to complete the identification of airflow vortex regions, concentration sequence analysis, and anomaly source location calculation. After the support vector machine outputs the spatial grid corresponding to the anomaly source, the processor parses the coordinate information of the spatial grid to obtain the true source coordinates of the gas environmental parameter anomaly source, and outputs or stores the coordinates as the environmental parameter monitoring result.

[0266] Example 13: like Figure 2 and Figure 5 Verification experiments were conducted for this application. The experimental platform and conditions are shown in Table 1. Table 1 Experimental Platform and Conditions

[0267] Note: This experiment uses a standard gas mixing device to simulate a fixed-point leak, and the leak flow rate and concentration can be repeatedly set; the inspection robot reciprocates at a set speed, and the lidar and gas sensor collect data synchronously; all data can be repeatedly acquired under the above conditions.

[0268] Experiment 1: Detection of Concentration Debris Distribution in Airflow Vortex Region Near the shielded area of ​​the support (x≈102-115m), a local airflow vortex is formed on the leeward side of the obstacle, and the gas concentration shows a fragmented high value distribution; according to S1-S4: the lidar point cloud and the gas sensor concentration sequence are aligned according to the sampling time and a spatial correspondence is established; the shielding contour is extracted from the surface normal and curvature of the point cloud, the velocity gradient matrix is ​​calculated and the vortex center is determined, and the spatial boundary of the vortex region is generated; the variance of the concentration sequence of 20 sampling nodes within the boundary is calculated by a sliding window with a window width of 5s, and a concentration distribution feature matrix (20 nodes × 15 windows) is constructed.

[0269] The variance threshold for determining fragment distribution is set to 1.5 ppm. 2 If the largest element of the feature matrix is ​​greater than a threshold, then fragmentation is determined to exist within the vortex region. Experimental results are shown below. Figure 2 (i.e., the heat map of the concentration distribution characteristic matrix).

[0270] Experiment 2: Response Delay Numerical Tail Recognition and Kalman Filtering for Noise Reduction The inspection robot patrols the channel from upstream to downstream of the leak source, and the sensor-measured concentration sequence is as follows: Figure 3 As shown; the leak occurred and increased from t=40s, and the concentration linearly decreased from t=55-90s as the robot moved away from the leak source; due to the influence of the sensor's first-order inertial response, the original sequence had obvious peak undershoot (measured peak 3.92ppm, true peak 5.20ppm) and fallback tail; according to S5-S7: extract the concentration change difference between adjacent sampling times and compare it with the change threshold to identify the numerical tail caused by the response delay; use Kalman filtering (the state vector contains the true concentration and sensor output, and the observation is the sensor reading) to estimate the state, and perform a moving average with a window width of 7 to denoise the filtered output to obtain the second gas concentration sequence.

[0271] Experiment 3: Comparison Experiment on Anomaly Source Localization Accuracy Under the same utility tunnel scenario and leakage conditions, five repeated location experiments were conducted using three different methods: Method A (traditional peak localization): directly take the peak position of the original concentration sequence as the coordinates of the anomaly source; Method B (removing only vortex interference): First, the interference of debris distribution in the vortex region is removed, but the sensor response delay tail is not processed, and the location is directly based on the processed peak value; Method C (this method): The fragment distribution identification and removal in the vortex region, response delay tail identification and Kalman filtering noise removal are completed in sequence. Then, the source feature vector (spatial gradient feature + distance from the geometric center of the vortex region) is constructed according to S8-S9. The true source coordinates are determined by support vector machine classification.

[0272] Each experiment used the same standard leak source (x=145m, y=1.5m) to generate a leak. The location coordinates of the three methods were recorded, and the Euclidean distance from the actual source coordinates was calculated as the location error. The experimental results are shown in […]. Figure 4 and Figure 5 .

[0273] Experimental results and comparison table: Table 2 Comparison of response delay tail suppression effects (Experiment 2, single measurement)

[0274] Note: Peak deviation is the difference between the measured peak value and the true peak value of 5.20 ppm; Tail duration is the portion of the time required for the concentration to drop below 1.30 ppm after the leakage process ends (t≥55s) that exceeds the ideal profile drop time (approximately 87s); Kalman filter state estimation compensates for the sensor's first-order inertial delay, making the sequence closer to the gas concentration state at the real spatial location.

[0275] Table 3. Detection results of debris distribution in the airflow vortex region (Experiment 1)

[0276] Note: The maximum element concentration in the three high-variance clusters within the vortex region is greater than 1.5 ppm. 2 The threshold value shows that the maximum variance in the non-fragmented region is only 0.43 ppm. 2 This demonstrates that the method can clearly distinguish between the fragment distribution in the vortex region and normal concentration fluctuations, providing a basis for subsequent fragment interference removal.

[0277] Table 4 Comparison of anomaly source localization errors among the three methods (Experiment 3, 5 repetitions)

[0278] Figure 2 In the diagram, the three high-variance clusters marked with white boxes correspond to the debris distribution areas at nodes 4, 11, and 16 within the vortex region, respectively, with maximum variances of 3.26, 2.88, and 3.07 ppm. 2 All significantly exceeded 1.5 ppm. 2 Judgment threshold; the maximum variance in other regions does not exceed 0.43 ppm. 2 This figure verifies the ability of this method to detect the concentration fragment distribution phenomenon in the vortex region.

[0279] Figure 3In the figure, the gray dashed line represents the ideal concentration profile, the orange curve represents the original sensor sequence, and the blue curve represents the sequence after filtering using this method. It can be seen that the original sequence has obvious peak undershoot (the measured peak value is only 3.92 ppm) and a tail of about 23 seconds. After state estimation and denoising using Kalman filtering using this method, the peak value is restored to 5.24 ppm (close to the true value of 5.20 ppm), the tail is shortened to about 4 seconds, the overall sequence is closer to the ideal profile, and the RMSE is reduced from 0.805 ppm to 0.315 ppm.

[0280] Figure 4 The planar distribution of the localization results of the three methods in 5 repeated experiments is shown; the red pentagram represents the true source (x=145m, y=1.5m); the green dots represent the localization results of method C (this method), which are closely clustered around the true source; the orange squares represent the localization results of method B (only vortex interference is removed), which are distributed about 2m downstream of the source; the red triangles represent the localization results of method A (traditional peak localization), with a deviation of more than 3m; the localization results of this method in all 5 experiments all fall within a radius of 0.8m.

[0281] Figure 5 The mean error bar represents the standard deviation of five repeated experiments. Method C (this method) has an average positioning error of 0.53m, which is about 83.5% lower than Method A (3.21m) and about 74.1% lower than Method B (2.05m), and has the smallest standard deviation (±0.19m), indicating the best positioning stability.

[0282] From Table 3 and Figure 2 As can be seen, this method constructs a concentration distribution feature matrix using a sliding window variance, which can significantly distinguish fragmented high-concentration regions within the vortex region from normal concentration fluctuations: the variance of the fragmented region (2.88–3.26 ppm) 2 The concentration is approximately in the non-fragmented region (≤0.43ppm). 2 The value is more than 7 times that of the original source; this identification result provides a basis for eliminating abnormal sources in the subsequent identification, and avoids misjudging the high value of debris in the vortex area as the real leakage source.

[0283] From Table 2 and Figure 3 As can be seen, after Kalman filtering state estimation and denoising, the RMSE of the concentration sequence decreased from 0.805 ppm to 0.315 ppm (a reduction of 60.9%), the measured peak value recovered from 3.92 ppm to 5.24 ppm (the peak deviation decreased from −1.28 ppm to +0.04 ppm), and the fallback tail duration was shortened from 23.0 s to 4.0 s (a reduction of 82.6%). This indicates that the method effectively compensates for the numerical tail caused by the sensor's first-order inertial response, making the second gas concentration sequence more realistically reflect the actual gas concentration state at the current spatial location.

[0284] From Table 4, Figure 4 and Figure 5It is evident that, under actual working conditions where vortex debris interference and sensor response delay coexist, the traditional peak positioning method (Method A) has an average error of 3.21m, the method that only eliminates vortex interference (Method B) has an average error of 2.05m, while the average error of our method (Method C) is reduced to 0.53m, which is about 83.5% lower than Method A and about 74.1% lower than Method B; at the same time, the standard deviation of the five repeated experiments is the smallest (±0.19m), indicating the best positioning consistency and stability.

[0285] Although embodiments of the invention have been shown and described, it will be understood by those skilled in the art that various changes, modifications, substitutions and alterations can be made to these embodiments without departing from the principles and spirit of the invention, the scope of which is defined by the appended claims and their equivalents.

Claims

1. A method for monitoring environmental parameters of a utility tunnel inspection robot, characterized in that, include: S1. Obtain the spatial point cloud and first gas concentration sequence of the utility tunnel site, which are synchronously collected by the lidar and gas sensor carried by the utility tunnel inspection robot during the utility tunnel inspection process and have a spatial correspondence established according to the same sampling time. S2. Extract obstacle and occlusion features within the utility tunnel based on the spatial point cloud at the site, and determine the spatial boundary of the airflow vortex region. S3. Perform sliding window variance calculation on the first gas concentration sequence within the spatial boundary of the airflow vortex region to obtain the concentration distribution feature matrix. S4. If the largest element in the concentration distribution feature matrix is ​​greater than the preset variance threshold, then it is determined that there is a fragment distribution phenomenon in the spatial boundary of the airflow vortex region. S5. Extract the concentration abrupt difference of the first gas concentration sequence at adjacent sampling times from the spatial boundary of the airflow vortex region where fragmentation exists. S6. If the concentration change difference exceeds the preset change threshold, it is determined that the first gas concentration sequence has a numerical tail caused by response delay. S7. The Kalman filter algorithm is used to perform state estimation and noise reduction on the first gas concentration sequence with numerical tails to obtain the second gas concentration sequence. S8. By extracting the spatial gradient features of the second gas concentration sequence and the geometric center distance of the spatial boundary of the airflow vortex region, a source feature vector is constructed. S9. Support vector machine is used to perform high-dimensional spatial mapping and classification of source feature vectors to determine the true source coordinates of abnormal gas environment parameters in the utility tunnel, and the true source coordinates are used as the environmental parameter monitoring results of the utility tunnel inspection robot.

2. The environmental parameter monitoring method for the utility tunnel inspection robot according to claim 1, characterized in that: S1 includes: The original point cloud and concentration sequence are obtained. The original point cloud and concentration sequence are generated synchronously by the lidar and gas sensor carried by the pipe gallery inspection robot during the pipe gallery inspection process. The aligned sequence is obtained by aligning the original point cloud with the concentration sequence using timestamps. A spatial correspondence is established between the aligned sequences to obtain the mapping sequence; Determine if there are data gaps in the mapping sequence; If there are gaps in the mapping sequence, the gaps are filled by a spatial interpolation fitting algorithm to obtain the spatial point cloud of the utility tunnel and the first gas concentration sequence.

3. The environmental parameter monitoring method for the utility tunnel inspection robot according to claim 1, characterized in that: S2 includes: Acquire spatial point cloud data of the utility tunnel and extract surface normals and curvature features; Obstacle occlusion features are extracted based on surface normals and curvature features, wherein the obstacle occlusion features include occlusion contours; Extract the flow velocity gradient matrix based on the occlusion contour; If the velocity gradient matrix satisfies a preset threshold, then the vortex center is determined. The boundary coordinates are generated based on the vortex center, and the spatial boundary of the airflow vortex region is determined based on the boundary coordinates.

4. The environmental parameter monitoring method for the utility tunnel inspection robot according to claim 1, characterized in that: S3 includes: Obtain sampling nodes within the spatial boundary of the airflow vortex region, and extract the first gas concentration sequence based on the sampling nodes; By extracting the first gas concentration sequence through a sliding window, a local concentration subsequence is obtained; Variance is calculated for local concentration subsequences to obtain the numerical value of local dispersion. Concentration distribution feature matrix is ​​constructed based on the numerical values ​​of local discreteness.

5. The environmental parameter monitoring method for the utility tunnel inspection robot according to claim 1, characterized in that: S4 includes: Extract the maximum element from the concentration distribution feature matrix; The maximum element of the matrix is ​​compared with a preset variance threshold to obtain the threshold comparison result; If the threshold comparison result indicates that the largest element of the matrix is ​​greater than the preset variance threshold, a boundary anomaly marker is generated. Based on the boundary anomaly markers, it was determined that there was a fragment distribution phenomenon within the spatial boundary of the airflow vortex region.

6. The environmental parameter monitoring method for the utility tunnel inspection robot according to claim 1, characterized in that: S5 includes: Obtain the first gas concentration sequence corresponding to the spatial boundary of the airflow vortex region; The concentration gradient vector is obtained by performing a difference operation on the first gas concentration sequence; A set of difference features is obtained by extracting features from the concentration gradient vector; Determine whether the feature values ​​of the difference feature set are greater than a preset mutation threshold; If the feature value is greater than the preset mutation threshold, the concentration mutation difference of the first gas concentration sequence at adjacent sampling times is extracted from the spatial boundary.

7. The environmental parameter monitoring method for the utility tunnel inspection robot according to claim 1, characterized in that: S6 includes: If the concentration change difference exceeds the preset change threshold, the corresponding adjacent sampling times are extracted to construct the delayed response interval; The concentration decay characteristic sequence is obtained by truncating the first gas concentration sequence based on the delayed response interval; The duration of the trailing effect and the amplitude of the difference fluctuation are obtained from the concentration decay characteristic sequence to construct the criteria for determining the trailing effect. If the concentration decay feature sequence meets the tailing determination condition, then it is determined that the first gas concentration sequence has a numerical tailing caused by response delay.

8. The environmental parameter monitoring method for the utility tunnel inspection robot according to claim 1, characterized in that: S7 includes: Obtain the first gas concentration sequence with numerical tails, and construct the state transition matrix and observation matrix based on the first gas concentration sequence; A priori state estimation sequence is generated using the state transition matrix and the observation matrix; If the system noise covariance calculated from the prior state estimation sequence is greater than a preset threshold, then the Kalman gain matrix is ​​determined based on the observed noise covariance. The Kalman gain matrix is ​​used to denoise the prior state estimation sequence to obtain the second gas concentration sequence.

9. The environmental parameter monitoring method for the utility tunnel inspection robot according to claim 1, characterized in that: S8 includes: Obtain the second gas concentration sequence, reconstruct the second gas concentration sequence to obtain the concentration distribution matrix, and extract the spatial gradient features of the concentration distribution matrix; Based on the spatial gradient characteristics, the airflow vortex region is determined, the spatial boundary of the airflow vortex region is obtained, and the position of the boundary centroid is calculated. Calculate the distance between the geometric centers of the spatial boundary by using the location of the boundary centroid; The source feature vector is constructed based on the spatial gradient characteristics and the distance to the geometric center.

10. The environmental parameter monitoring method for the utility tunnel inspection robot according to claim 1, characterized in that: S9 includes: A support vector machine is used to perform high-dimensional space mapping on the source feature vectors to obtain a high-dimensional feature space; The feature distribution matrix is ​​extracted from the high-dimensional feature space. The decision boundary hyperplane is determined by the feature distribution matrix and then classified to obtain anomaly classification labels. If the anomaly classification label indicates an anomaly source, then extract the spatial coordinate grid corresponding to the anomaly source; The spatial coordinate grid is analyzed to determine the true source coordinates, and these true source coordinates are used as the environmental parameter monitoring results for the utility tunnel inspection robot.

11. An environmental parameter monitoring system for a utility tunnel inspection robot, used to implement the steps of the environmental parameter monitoring method for a utility tunnel inspection robot as described in any one of claims 1-9, characterized in that, The system includes: The data acquisition module acquires the spatial point cloud and the first gas concentration sequence of the utility tunnel site, which are synchronously collected by the lidar and gas sensor carried by the utility tunnel inspection robot during the utility tunnel inspection process and establish spatial correspondence according to the same sampling time. The vortex region determination module extracts obstacle and occlusion features within the utility tunnel based on the spatial point cloud at the site, and determines the spatial boundary of the airflow vortex region. The concentration feature calculation module performs sliding window variance calculation on the first gas concentration sequence within the spatial boundary of the airflow vortex region to obtain the concentration distribution feature matrix. The debris distribution determination module determines that if the largest element in the concentration distribution feature matrix is ​​greater than the preset variance threshold, then the debris distribution phenomenon exists within the spatial boundary of the airflow vortex region. The concentration mutation extraction module extracts the concentration mutation difference of the first gas concentration sequence at adjacent sampling times from the spatial boundary of the airflow vortex region where fragmentation exists. The numerical tailing determination module determines that if the concentration change difference exceeds the preset change threshold, the first gas concentration sequence has numerical tailing caused by response delay. The concentration filtering module uses the Kalman filter algorithm to perform state estimation and noise reduction on the first gas concentration sequence with numerical tails to obtain the second gas concentration sequence. The source feature construction module constructs a source feature vector by extracting the spatial gradient features of the second gas concentration sequence and the geometric center distance of the spatial boundary of the airflow vortex region. The anomaly source localization module uses a support vector machine to perform high-dimensional spatial mapping and classification of the source feature vectors, determine the true source coordinates of the gas environment parameter anomalies in the utility tunnel, and use the true source coordinates as the environmental parameter monitoring results of the utility tunnel inspection robot.

12. A computer-readable storage medium, characterized in that, The computer instructions stored thereon execute the environmental parameter monitoring method of the utility tunnel inspection robot according to any one of claims 1-10 when the computer instructions are executed.

Images