Rock fissure three-dimensional trajectory inversion method and system based on multi-source sensor fusion

Abstract

The application provides a rock fissure three-dimensional trajectory inversion method and system based on multi-source sensor fusion, and the system comprises: a plurality of packaged detectors, which are internally integrated with a motion measurement unit, a magnetic field measurement unit and a width measurement unit; the motion measurement unit is used for collecting IMU data of the packaged detector, the environment measurement unit is used for collecting magnetic field characteristic data of the rock fissure, and the width measurement unit is used for collecting width data of the rock fissure; a plurality of magnetic markers are arranged at intervals on the surface or inside of the rock mass; an injection device is used for injecting fluid or gas into the rock fissure, sending the plurality of packaged detectors into the rock fissure, and driving the plurality of packaged detectors to move in the rock fissure; a relay network is arranged on the surface of the rock mass and is used for transmitting motion data, magnetic field characteristic data and width data collected by the packaged detector; and a mapping module is used for inversely obtaining the three-dimensional topography of the rock fissure according to the IMU data, the magnetic field characteristic data and the width data.

Description

Technical Field

[0001] This application belongs to the field of rock mass engineering exploration technology, and in particular relates to a method and system for inverting the three-dimensional trajectory of rock fractures based on multi-source sensor fusion. Background Technology

[0002] Rock fissures, as a major defect in rock mass structure, rely heavily on their three-dimensional spatial morphological parameters, including strike, width, depth, branch distribution, and internal filling characteristics, as core evidence for geotechnical stability evaluation and geological hazard early warning. Accurately obtaining the three-dimensional contour information of rock fissures is of paramount scientific and engineering value for geotechnical engineering fields such as tunnel excavation, slope stabilization, mining, and oil and gas fracturing. However, due to the heterogeneity of rock media and the concealment and complexity of rock fissure spaces, achieving high-precision, non-destructive mapping of the three-dimensional morphology of rock fissures has long been a challenge in the field of geotechnical engineering exploration.

[0003] The current rock fracture detection technology mainly includes the following four technical routes, and each method has its inherent limitations: (1) The first type is the borehole exploration method. A representative technology is Chinese application CN103114848A, which reconstructs the spatial distribution of stratigraphic fractures by drilling rock cores and conducting laboratory analysis. It requires drilling a large number of exploration holes in the rock mass, extracting rock core samples or lowering optical probes through the exploration holes, and then obtaining fracture information at the exploration hole wall. However, this type of method will cause irreversible damage to the rock mass structure. In addition, the drilling cost is extremely high, with the comprehensive cost of a single hole reaching hundreds of yuan per linear meter. Moreover, the information obtained through the exploration holes is local data near the exploration holes, and it cannot accurately reconstruct the three-dimensional morphology of rock fractures for non-linear or branched complex fracture networks. (2) The second type is the geophysical exploration method. It uses ground-penetrating radar to transmit electromagnetic waves to invert the rock mass structure. Its principle is to invert the underground structure by transmitting and receiving physical wave fields and using the difference in physical properties between the rock mass and the fractures. However, it is difficult to balance detection depth and resolution. For fine rock fissures with a width of less than centimeters, the identification error often exceeds 1 meter. The detection effect in rock fissures filled with water or mud is extremely poor, and its three-dimensional imaging results rely on complex inversion algorithms, resulting in serious ambiguity and strong subjectivity in interpretation. (3) The third category is magnetic marker detection method. For example, Chinese application CN108051497A, by injecting slurry containing magnetic materials into rock fissures, and after the slurry solidifies, a three-dimensional magnetic field scanning system is set up outside the rock mass to measure the spatial magnetic field distribution. The distribution pattern of the fissure is inverted through the magnetic field strength. Although this method introduces magnetic field detection means, it still has inherent defects such as relying on external equipment to scan and obtain the magnetic field distribution, being unable to obtain continuous trajectory data of sensors inside the fissure, the magnetic material being fixed in the fissure after grouting, making it difficult to accurately restore the three-dimensional orientation and spatial distribution of the fissure, and the slurry not being able to completely fill narrow or curved fissures, resulting in detection blind zones. (4) The fourth category is micro-robot detection method. The principle is to use a miniature mobile robot to directly probe inside rock fissures. However, the uneven surface inside rock fissures can easily cause the mobile robot's movement mechanism to jam or tip over. In addition, due to the minimum size of the power system and transmission mechanism, the diameter of current miniature mobile robots is generally greater than 3cm, making it difficult to enter narrower cracks, and their large turning radius makes it impossible to navigate sharp bends.

[0004] In summary, it is necessary to provide a new type of rock fissure detection technology to achieve high-precision detection of rock fissures. Summary of the Invention

[0005] In view of this, this application provides a method and system for inverting the three-dimensional trajectory of rock fissures based on multi-source sensor fusion, which aims to achieve high-precision detection of the three-dimensional morphology of rock fissures without damaging the rock mass or with minimal damage to the rock mass.

[0006] In a first aspect, this application provides a rock fracture three-dimensional trajectory inversion system based on multi-source sensor fusion, comprising: Several encapsulated detectors are provided, each encapsulating detector integrating a motion measurement unit, a magnetic field measurement unit, and a width measurement unit. The motion measurement unit is used to collect IMU data from the encapsulated detectors, the environmental measurement unit is used to collect magnetic field characteristic data of rock fissures, and the width measurement unit is used to collect width data of rock fissures. Several magnetic markers are placed at intervals on the surface or inside the rock mass; An injection device is used to inject fluid or gas into rock fissures, send several encapsulated detectors into the rock fissures, and drive several encapsulated detectors to move within the rock fissures. A relay network, set up on the rock surface, is used to transmit motion data, magnetic field characteristic data, and width data collected by the encapsulated detector; The mapping module is used to invert the three-dimensional morphology of the rock fissure based on IMU data, magnetic field characteristic data, and width data.

[0007] Optionally, the encapsulated detector also integrates a pressure measurement unit for collecting the pressure generated when the encapsulated detector comes into contact with the rock mass.

[0008] Secondly, this application provides a method for inverting three-dimensional rock fracture trajectories based on multi-source sensor fusion, which is applied to the aforementioned three-dimensional rock fracture trajectory inversion system based on multi-source sensor fusion, including: Several encapsulated detectors were sent into the rock crevice, and the detectors were driven to move within the crevice. The IMU data of the encapsulated detectors, the magnetic field characteristics of the rock crevice, and the width of the rock crevice were collected. The initial trajectory is determined based on IMU data collected by several encapsulated detectors; Based on the magnetic field characteristic data collected by several encapsulated detectors, the initial trajectory is optimized to obtain the final trajectory; Based on the final trajectory and the width data of the rock fissure, the three-dimensional morphology of the rock fissure is determined.

[0009] Optionally, the step of determining several initial trajectories based on IMU data acquired by several encapsulated detectors includes: Savitzky-Golay filtering was applied to the IMU data. The attitude of the encapsulated detector as it moves within the rock crevice is determined using the quaternion method based on the angular velocity data in the IMU data. Based on the attitude of the encapsulated detector and the acceleration data in the IMU data, the trajectory of the encapsulated detector in the rock crevice is determined as the initial trajectory.

[0010] Optionally, when determining the movement trajectory of the encapsulated detector in the rock crevice, the zero-velocity of the encapsulated detector can be determined based on IMU data, and the speed of the encapsulated detector can be calibrated based on the zero-velocity determination.

[0011] Optionally, the process of optimizing the initial trajectory based on magnetic field characteristic data collected by several encapsulated detectors includes: Obtain the magnetic field characteristic data of the magnetic beacon, and determine several position constraint matching pairs based on the magnetic field characteristic data of the magnetic beacon and the magnetic field characteristic data of the encapsulated detector. Based on several position constraint matching pairs, strong constraints are applied to the position points in the initial trajectory that correspond to the target position matching pairs. A cost function is constructed, and an optimization method is used to optimize the non-strongly constrained position points in the initial trajectory based on the cost function to obtain the final trajectory.

[0012] Optionally, the step of determining several position constraint matching pairs based on the magnetic field gradient sequence, according to the magnetic field characteristic data of the magnetic marker and the magnetic field characteristic data of the encapsulated detector, includes: The magnetic field feature data of the magnetic beacon and the magnetic field feature data collected by the encapsulated detector are divided into several key frames according to a fixed number of frames. The time domain features and frequency domain features of each key frame are extracted, and the time domain features and frequency domain features are fused into a feature vector. The gradient of the feature vector of several frames is extracted to obtain the magnetic field gradient sequence of the magnetic beacon and the encapsulated detector. Based on the magnetic field gradient sequence, the magnetic field feature data of the magnetic marker is matched with the magnetic field feature data collected by the encapsulated detector to obtain several candidate position matching pairs. Based on the pressure feature data collected by the encapsulated detector, several tactile events are determined to form a tactile signature sequence; the sampling time point of the magnetic field feature in the candidate position matching pair is determined, and tactile events under a preset time window before and after the sampling time point of the magnetic field feature in the candidate position matching pair are extracted to form a first local tactile subsequence and a second local tactile subsequence; based on the similarity between the first local tactile subsequence and the second local tactile subsequence, position constraint matching pairs are selected from the candidate position matching pairs.

[0013] Optionally, the method, based on a magnetic field gradient sequence, involves matching the magnetic field feature data of a magnetic marker with the magnetic field feature data acquired by an encapsulated detector, including: Calculate the cosine similarity between the frequency domain characteristics of the magnetic marker and the frequency domain characteristics of the encapsulated detector. Based on the cosine phase velocity, perform a coarse matching between the trajectory points of the encapsulated detector and the position of the magnetic marker to obtain a coarse matching pair. The magnetic field gradient sequences of several coarse matching pairs are compared using the dynamic time warping method or the longest common subsequence method to obtain candidate position matching pairs from the coarse matching pairs.

[0014] Optionally, the cost function includes:

[0015] In the formula: This represents the attitude constraint loss. This represents the trajectory smoothing loss. Indicates rock wall penetration loss, , , These are the weighting coefficients.

[0016] Optionally, the optimization method used may include particle swarm optimization, genetic algorithm or L-BFGS-B algorithm.

[0017] The beneficial effects of the technical solution provided in this application include: (1) This application uses a packaged detector plus a fluid or gas injection device to detect the morphology of rock fissures. Only a single centimeter-sized injection hole needs to be formed on the rock surface. The packaged detector is then sent into the rock fissure through the injection hole. There is no need for dense drilling, and the damage to the rock structure is minimal. It is particularly suitable for safety inspection and long-term monitoring of existing projects. In addition, the packaged detector moves passively, which eliminates the risk of jamming in the moving mechanism, as is the case with micro-robots. The packaged detector can pass smoothly through narrow cracks, sharp bends and water accumulation areas with a width of not less than 0.8 cm, covering the detection blind spots that traditional micro-robots cannot reach.

[0018] (2) As a further improvement, the initial trajectory of the rock fissure is formed by IMU data, and the initial trajectory is optimized by combining magnetic field data to limit the trajectory drift of the initial trajectory and improve the accuracy of the final trajectory. Attached Figure Description

[0019] To more clearly illustrate the technical solutions in this invention or the prior art, the drawings used in the description of the embodiments or the prior art will be briefly introduced below. Obviously, the drawings described below are some embodiments of this invention. For those skilled in the art, other drawings can be obtained from these drawings without creative effort.

[0020] Figure 1 This is a structural block diagram of a rock fracture three-dimensional trajectory inversion system based on multi-source sensor fusion provided in an embodiment of this application; Figure 2 A schematic diagram of the structure of a rock fissure three-dimensional trajectory inversion system based on multi-source sensor fusion provided in an embodiment of this application; Figure 3 This is a schematic diagram of the structure of a packaged detector provided in one embodiment of this application; Figure 4A flowchart of a rock fracture three-dimensional trajectory inversion method based on multi-source sensor fusion provided in an embodiment of this application; Figure 5 A schematic diagram illustrating the initial movement trajectory generation process provided in an embodiment of this application; Figure 6 A schematic diagram of an initial movement trajectory optimization process provided in an embodiment of this application; Figure 7 A trajectory optimization comparison diagram provided for an embodiment of this application. Detailed Implementation

[0021] To make the objectives, technical solutions, and advantages of this invention clearer, the technical solutions of this invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some, not all, of the embodiments of this invention. All other embodiments obtained by those skilled in the art based on the embodiments of this invention without creative effort are within the scope of protection of this invention.

[0022] See Figures 1 to 3 This application provides a three-dimensional trajectory inversion system for rock fissures based on multi-source sensor fusion, comprising: Several encapsulated detectors 1 are provided, each encapsulated detector integrating a motion measurement unit 11, an environmental measurement unit 12, and a width measurement unit 13. The motion measurement unit 11 is used to collect motion data of the encapsulated detector, the environmental measurement unit 12 is used to collect magnetic field characteristic data of the rock fissure 7, and the width measurement unit 13 is used to collect width data of the rock fissure 7. Several magnetic markers 5 are spaced out on the surface or inside the rock mass; The injection device 2 is used to inject fluid or gas into the rock fissure, send several encapsulated detectors 1 into the rock fissure, and drive several encapsulated detectors 1 to move within the rock fissure 7. Relay network 3 is installed on the surface of the rock mass 6 to receive and forward data measured by the encapsulated detector 1; The mapping module 4 is used to invert the three-dimensional morphology of rock fissure 7 based on motion data, magnetic field characteristic data, and width data.

[0023] This application provides a three-dimensional trajectory inversion system for rock fissures based on multi-source sensor fusion. The system includes several encapsulated detectors, each integrating a motion measurement unit, an environmental measurement unit, and a width measurement unit. These detectors can collect motion data, magnetic field characteristic data of the rock fissure, and width data. By employing an injection device, several encapsulated detectors are injected into the rock fissure. The detectors move along the fissure under the influence of the fluid or gas injected. These detectors can flow to different branches of the fissure. The motion, environmental, and width data collected by these detectors effectively reflect the main body and branching directions of the rock fissure, allowing for a more accurate inversion of the three-dimensional contour of the fissure. A relay network along the rock mass receives data from the encapsulated detectors within the fissure, reducing the data transmission path length and minimizing the impact of the rock mass on data transmission.

[0024] See Figure 3 In some embodiments, as shown in (A), the encapsulated detector 1 has a streamlined shape with a maximum outer diameter less than or equal to 2 cm. The housing 15 of the encapsulated detector 1 has a streamlined surface to reduce movement resistance, and integrates an ultra-thin flexible battery 18 and / or an energy harvester (not shown in the figure), providing a battery life of no less than 6 hours. In some embodiments, in rock fissures with seepage, the detector can integrate a micro hydroelectric turbine as an energy harvester, using water flow for charging to achieve unlimited battery life; it can also integrate a radio frequency energy harvesting module as an energy harvester, obtaining energy from the electromagnetic waves emitted by the repeater, with an operating frequency matching the repeater's 100-300 kHz.

[0025] See Figure 3 (A) and Figure 3 (D) in the diagram represents the axial and radial cross-sectional views of the encapsulated detector, respectively.

[0026] In some embodiments, a superhydrophobic coating 16 may be formed on the surface of the housing 15 of the encapsulated detector 1, and the interior of the housing 15 may be waterproofed by waterproof sealant 17.

[0027] See Figure 3 (B) in the figure specifically shows the environmental measurement unit 12, which is a magnetic field sensor (Magn) array.

[0028] In some embodiments, the motion measurement unit 11 includes an inertial measurement unit.

[0029] In some embodiments, the environmental measurement unit 12 includes a magnetic field sensor array.

[0030] In other embodiments, the magnetic field sensor array can be replaced with a quantum magnetometer, which improves the resolution to the pT level for use in weak magnetic environments.

[0031] In some embodiments, the width measuring unit 13 includes a miniature ultrasonic ranging sensor.

[0032] In some embodiments, the inertial measurement unit may be an MPU-605, MPU-9250, or a fiber optic gyroscope, suitable for ultra-long-distance (e.g., greater than 50m) rock fissure detection.

[0033] See you again Figure 3 In some embodiments, as shown in (A), the encapsulated detector 1 may also integrate a pressure measurement unit 14 for collecting the pressure generated when the encapsulated detector comes into contact with the rock mass. The pressure measurement unit 14 may be a flexible pressure tactile sensor.

[0034] See Figure 3 In (C), the flexible pressure tactile sensor generally includes a stacked flexible substrate layer 141, a sensing layer 142, and a connector 143.

[0035] In some embodiments, if the injection device is a pneumatic injection device, the injection device 2 includes a high-pressure air pump and an injection conduit 21. If the injection device is a fluid injection device, the injection device includes an inert fluid storage tank and an injection conduit 21.

[0036] In some embodiments, the injection conduit 21, with an inner diameter not exceeding 2 cm, is used to generate air pressure or liquid flow entrainment force to send the encapsulated detector into the rock crevice through the injection hole and to make the encapsulated detector move passively. The moving speed can be controlled by pressure regulation to be between 0.5 and 2 cm per second.

[0037] In some embodiments, if the rock fissure is a vertical or steeply inclined rock fissure, gravity can be used for auxiliary injection, such as increasing the weight of the encapsulated detector so that the encapsulated detector slides down by gravity, in conjunction with flushing by top fluid; for extremely narrow fissures (<0.5cm), the detector can be further miniaturized to a diameter of 0.5cm, and high-pressure gas pulse injection can be used with a pulse frequency of 1Hz and a duty cycle of 30%.

[0038] In some embodiments, the relay network 3 includes a low-frequency electromagnetic relay network, consisting of multiple miniature low-frequency electromagnetic repeaters 31, with an operating frequency between 100 and 300 kHz, a single unit transmission distance of not less than 5 meters, a diameter of not more than 8 millimeters, and is deployed along the rock fissure path to receive and forward sensor data emitted by the detector.

[0039] In some embodiments, the relay network 3 can use fiber optic microcables instead of electromagnetic relays to achieve ultra-low loss transmission; the detector can use backscatter communication technology to reflect external carrier waves to transmit data, and the power consumption at the detector end can be reduced to the μW level.

[0040] In some embodiments, the mapping module 4 includes: The embedded computing module, based on a low-power embedded platform, collects and preprocesses sensor data in real time, and runs trajectory calculation and optimization algorithms; and the 3D drawing terminal, equipped with point cloud processing software and custom algorithm plugins, is used to receive optimized trajectory data and generate a 3D contour model of the rock fissure.

[0041] In some embodiments, the magnet used is a nanoscale neodymium iron boron magnet.

[0042] In some embodiments, a rock fissure 7 is formed in the rock mass 6. The rock fissure 7 includes a main part and a branch part 71. This application can detect the main part of the rock fissure and can also detect the branch part of the rock fissure to a certain extent.

[0043] Figure 4 This is a flowchart illustrating a method for inverting three-dimensional rock fracture trajectories based on multi-source sensor fusion, provided in one embodiment of this application. See also... Figure 4 ,include: S101, drive the encapsulated detector to move along the rock fissure, and acquire the motion data of the encapsulated detector, the magnetic field characteristic data of the rock fissure, and the width data of the rock fissure collected during the movement.

[0044] First, injection preparation and relay deployment are carried out: an injection hole with a diameter of no more than 2 cm is drilled at the entrance of the rock fissure, and miniature low-frequency electromagnetic repeaters are deployed on the surface or inside the rock mass every 3 to 5 meters along the estimated path to build an uninterrupted data transmission network.

[0045] Choose either pneumatic or fluid injection mode based on the rock fracture environment; additionally, several nanoscale neodymium iron boron magnets can be installed on or within the rock mass to construct artificial feature anchor points.

[0046] Next, sensor calibration and injection are performed: the detector is placed at the front end of the injection catheter, static IMU data is collected, and the zero bias values ​​of acceleration and angular velocity of the IMU data are calculated for zero bias correction of the IMU data.

[0047] The injection device is activated, and the injection pressure is adjusted to allow the encapsulated detector to inject into the rock fissure at a set speed. The detector moves passively under pneumatic thrust or fluid entrainment, simultaneously acquiring IMU data, magnetic field characteristic data of the rock fissure, width data of the rock fissure, and optional tactile pressure data (if a pressure measurement unit is integrated, tactile pressure data is further acquired). During the movement, the ultrasonic detector monitors the width of the rock fissure in real time, and automatically reduces the injection speed when the width approaches the diameter of the detector to avoid jamming.

[0048] S102. Based on the IMU data collected by several encapsulated detectors, determine several initial trajectories.

[0049] One encapsulated detector corresponds to one initial trajectory. In the process of determining the initial trajectory and the final trajectory, several initial trajectories can be optimized separately. After optimization, the several initial trajectories can be combined into one final trajectory. Alternatively, several initial trajectories can be combined into one and then optimized to obtain one final trajectory.

[0050] In some embodiments, step S102 includes: Step 1: Perform Savitzky-Golay filtering on the IMU data.

[0051] Specifically, a Savitzky-Golay filter is used to smooth the IMU data, thereby suppressing high-frequency fluid disturbance noise.

[0052] The window length and order used in the Savitzky-Golay filter can be set according to specific circumstances.

[0053] Step 2: Use the quaternion method to determine the attitude of the encapsulated detector as it moves within the rock crevice based on the angular velocity data in the IMU data.

[0054] In some embodiments, quaternion attitude update is performed based on IMU data to obtain the attitude of the encapsulated detector as it moves in the crack. The specific process is as follows: Step 1: Angular velocity of the encapsulated detector Zero bias calibration is performed as follows:

[0055] In the formula: This represents the angular velocity of the calibrated encapsulated detector. , , This indicates the angular velocity of the encapsulated detector in , , Zero bias in three directions.

[0056] The angular velocity increment of the packaged detector in quaternion form Represented as: .

[0057] In the formula: This indicates the sampling time step, corresponding to the IMU's sampling frequency. For example, (Corresponding to a sampling frequency of 200Hz).

[0058] The second step is to perform quaternion pose update. The rules for quaternion pose update are as follows:

[0059] In the formula: The quaternion in quaternion form The attitude of the encapsulated detector after the next attitude update. The quaternion form of the first The attitude of the encapsulated detector after the next attitude update. This represents quaternion multiplication.

[0060] After obtaining the attitude of the encapsulated detector in the rock crevice through the above quaternion attitude update process, the attitude is normalized as follows:

[0061] In the formula: This indicates the orientation of the packaged detector before normalization. This indicates the normalized orientation of the encapsulated detector.

[0062] Step 3: Based on the attitude of the encapsulated detector and the acceleration data in the IMU data, determine the movement trajectory of the encapsulated detector in the rock crevice as the initial trajectory.

[0063] The transformation matrix is ​​determined based on the attitude of the encapsulated detector. The acceleration of the encapsulated detector is transformed from the carrier coordinate system to the world coordinate system, and the trajectory of the encapsulated detector in the rock crevice is determined based on the transformed acceleration, which serves as the initial trajectory.

[0064] In some embodiments, the process of determining the initial trajectory is as follows:

[0065] In the formula: Indicates the first The trajectory points of the encapsulated detector obtained from the second trajectory update. Indicates the first The trajectory points of the encapsulated detector obtained from the second trajectory update. Indicates the first The velocity of the encapsulated detector after the next trajectory update. Indicates the sampling time step.

[0066] The speed of the encapsulated detector during each update process The calculation process is as follows:

[0067] In the formula: Represents the first in world coordinates The acceleration of the encapsulated detector after the next trajectory update.

[0068] During each update, the coordinate system transformation process for acceleration is as follows:

[0069] In the formula: Indicates according to the first The attitude determination transformation matrix of the encapsulated detector after each attitude update, with the same time step used for each attitude update as for the position update.

[0070] In some optional embodiments, during the process of determining the initial trajectory, the encapsulated detector can also be zero-velocity determined based on IMU data, and the velocity of the encapsulated detector can be calibrated based on the zero-velocity determination. The process is as follows: Calculate the acceleration vector of the denoised encapsulated detector Length of the module and angular velocity vector Length of the module .

[0071] If satisfied and If the encapsulated detector is determined to be stationary, its velocity will be set to zero to reduce trajectory drift of the initial trajectory.

[0072] Finally, by fitting the obtained position points of the encapsulated detector, the initial trajectory can be obtained.

[0073] See Figure 5 The process of forming the initial trajectory described in this application is specifically illustrated.

[0074] S103. Based on the magnetic field characteristic data collected by several encapsulated detectors, the initial trajectory is optimized to obtain the final trajectory.

[0075] In some embodiments, step S103 includes: S1031. Obtain the magnetic field characteristic data of the magnetic beacon, and determine several position constraint matching pairs based on the magnetic field characteristic data of the magnetic beacon and the magnetic field characteristic data of the encapsulated detector.

[0076] S1032. Based on several position constraint matching pairs, strong constraints are applied to the position points in the initial trajectory that correspond to the target position matching pairs, and a cost function is constructed. Based on the cost function, an optimization method is used to optimize the non-strongly constrained position points in the initial trajectory to obtain the final trajectory.

[0077] In some embodiments, step S1031 includes: Step 1: Obtain the magnetic field feature data of the magnetic beacon. Divide the magnetic field feature data of the magnetic beacon and the magnetic field feature data collected by the encapsulated detector into several key frames according to a fixed number of frames. Extract the time domain features and frequency domain features of each key frame. Fuse the time domain features and frequency domain features into a feature vector. Extract the gradient of the feature vector of several frames to obtain the magnetic field gradient sequence of the magnetic beacon and the encapsulated detector.

[0078] Among them, an appropriate overlap rate can be set to generate hundreds of keyframes.

[0079] In some embodiments, the extracted time-domain features include the mean, variance, peak value, and root mean square of the magnetic field feature data.

[0080] In some embodiments, the extracted frequency domain features include: dominant frequencies at different amplitudes extracted by Fourier transform.

[0081] Subsequently, the time-domain features and frequency-domain features are fused into a high-dimensional feature vector, and the feature vector is then normalized.

[0082] After obtaining the feature vectors corresponding to several frames of magnetic field feature data, the magnetic field gradient sequence of the feature vectors is obtained by calculating the gradient of the feature vectors. This magnetic field gradient sequence can be used as a topological feature.

[0083] Step 2: Based on the magnetic field gradient sequence, match the magnetic field feature data of the magnetic marker with the magnetic field feature data collected by the encapsulated detector to obtain several candidate position matching pairs.

[0084] In some embodiments, step 2 includes: The first step is to calculate the cosine similarity between the frequency domain characteristics of the magnetic marker and the frequency domain characteristics of the encapsulated detector. Based on the cosine phase velocity, the trajectory points of the encapsulated detector and the position of the magnetic marker are coarsely matched to obtain a coarse matching pair.

[0085] The second step involves comparing the magnetic field gradient sequences of several candidate position matching pairs using the dynamic time warping method or the longest common subsequence method, and obtaining candidate position matching pairs from several coarse matching pairs.

[0086] One approach is to use the longest common subsequence method to compare the magnetic field gradient sequences of coarse matching pairs. The matching condition is that the Euclidean distance is less than a threshold, and the matching pairs with normalized similarity greater than the threshold are taken as the target matching pairs.

[0087] In some embodiments, the calculation process of the Longest Common Subsequence (LCSS) method is as follows:

[0088] In the formula: Represents a sequence with sequence The similarity value. Represents a sequence with sequence The length of the longest common subsequence. Represents a sequence The total length. Represents a sequence The total length.

[0089]

[0090] In the formula: Indicates removing the sequence with sequence The length of the longest common subsequence of the remaining parts after the last element. Represents a sequence The first in Each feature value. Represents a sequence The first in Each feature value. This represents the distance threshold.

[0091] In some embodiments, magnetic field feature matching can also use Dynamic Time Warping (DTW) instead of LCSS, and the specific calculation process is as follows:

[0092] In the formula: Represents a sequence with sequence The similarity value. Indicates the alignment path.

[0093] In some embodiments, to further improve the reliability and robustness of feature matching, this method, after completing the coarse matching based on frequency domain features, also performs a mismatch removal step based on the magnetic field gradient topology map. Specifically: In the keyframe segmentation stage, for the magnetic field triaxial data of each frame, not only are its temporal and frequency domain features extracted, but the magnetic field gradient changes between consecutive frames are also calculated. Let the first frame be... The magnetic field vector of the frame is Then the first Frame magnetic field gradient It can be approximated as Take continuous Gradient sequence of frames The construction dimension is The topological feature vector of the magnetic field gradient. This topological feature describes the spatial variation trend of the magnetic field in a local region, reflects the spatial variation pattern of the magnetic field caused by the distribution of magnetic minerals in rock fissures or artificial magnetic markers, and has local spatial uniqueness.

[0094] In the feature matching stage, for coarse matching pairs Extract the corresponding magnetic field gradient topological feature vectors respectively. and Calculate the normalized Euclidean distance between the two. ,in It is a very small positive number to prevent division by zero. If Greater than the preset threshold (Preferred, = 0.15), then the coarse match is determined to be a false match and is discarded; if Less than or equal to If the match is not found, the pair will be retained as a candidate and proceed to the subsequent fine-matching stage.

[0095] The aforementioned mismatch elimination strategy based on the magnetic field gradient topology map effectively filters out mismatched pairs caused by similar frequency domain features but mismatched spatial locations, significantly improving the accuracy of feature matching. Experiments show that this strategy can increase the matching accuracy from approximately 96% to over 98% without significantly increasing computational overhead.

[0096] Step 3: Based on the pressure feature data collected by the encapsulated detector, determine several tactile events to form a tactile signature sequence; determine the sampling time point of the magnetic field feature in the candidate position matching pair, extract the tactile events before and after the sampling time point of the magnetic field feature in the candidate position matching pair under a preset time window to form a first local tactile subsequence and a second local tactile subsequence; based on the similarity between the first local tactile subsequence and the second local tactile subsequence, select position constraint matching pairs from the candidate position matching pairs.

[0097] Step 3 is an optional step. In some other embodiments, if the encapsulated detector integrates a flexible pressure sensor, the pressure data generated by the detector's contact with the rock wall during its movement can be used for auxiliary positioning and matching verification.

[0098] Specifically, when the flexible pressure sensor passes through a curved or narrow section of a rock crevice, it generates pressure pulse signals due to compression against the rock wall. This method performs real-time peak detection on the time series output by the pressure sensor and sets a pressure threshold. (Preferred, =0.5N), when the continuously detected pressure value exceeds A tactile event is recorded when at least three sampling points are collected. All tactile events during the detection process are arranged chronologically to construct a tactile signature sequence. ,in Indicates the first The moment when a tactile event occurs.

[0099] This tactile signature reflects the geometric changes along the rock fissure—characteristic tactile response patterns are generated at locations such as bends, narrow sections, and wall protrusions. Since these geometric features are inherent properties of the rock fissure, the tactile signature, similar to magnetic field characteristics, also exhibits spatial repeatability. Therefore, the tactile signature can be used as a high-confidence local feature in the matching process, specifically as follows: Once the software completes the matching based on magnetic field characteristics and obtains candidate matching pairs, for any candidate matching pair... Determine the magnetic field matching frame respectively and The corresponding time, and retrieve the time window before and after this frame. (Preferred, Construct local tactile subsequences from tactile events within a time interval of 1.0 second. and .

[0100] Subsequently, the temporal alignment similarity between the two tactile subsequences was calculated. ,in and These are the temporal distribution characteristics of tactile events in the two subsequences (such as the average interval between events or the arrival time sequence). To adjust the bandwidth parameters.

[0101] like Greater than the preset threshold (Preferred, If the confidence level is 0.7, then the candidate pair is considered to still satisfy the tactile consistency constraint, and its confidence level is further increased; if If the magnetic field similarity is below a preset threshold and the confidence level of the magnetic field match is not high (e.g., the cosine similarity is between 0.70 and 0.85), the candidate matching pair is marked as a suspicious pair, and the subsequent optimization algorithm will automatically decide whether to retain it based on global consistency.

[0102] In a preferred embodiment, the tactile signature and magnetic field features employ a hierarchical joint matching strategy. First, the magnetic field features undergo coarse and fine matching to obtain an initial matching set. Then, for For each pair of matches, calculate its tactile signature consistency score.

[0103] If a matched pair has a magnetic field similarity ranking in the top 10% and a tactile consistency score above a threshold, the match is marked as a "high-confidence match"; if the magnetic field similarity is low but the tactile consistency score is extremely high (e.g., ...), the match is marked as a "high-confidence match". If the value is greater than 0.85, the match will be reserved separately as a "tactile enhancement match" for use as a strong constraint in subsequent trajectory optimization.

[0104] Finally, all verified matching pairs are used as attitude constraints and input into the trajectory optimization step to participate in the construction and solution of the cost function. Through the joint strategy of magnetic field gradient topology map mismatch elimination and tactile signature-assisted verification, this method can still maintain a high feature matching accuracy in complex rock fissure environments with sparse or highly repetitive magnetic field features, providing a reliable data association foundation for subsequent multi-source data fusion and trajectory optimization.

[0105] See Figure 6 This paper specifically illustrates the process by which this application obtains high-confidence feature matching pairs based on magnetic field feature data.

[0106] In some embodiments, step S1032 includes: Step 1: Construct a cost function based on affine transformation, taking into account attitude constraints, trajectory smoothing constraints, and rock wall penetration penalty terms.

[0107] In some embodiments, the cost function is specifically expanded as follows.

[0108] Suppose that a total of [number] features are obtained after feature matching. For high-confidence matching keyframes, denote the first... The spatial positions of the matching keyframes in the two tracks are respectively and To eliminate potential rigid body transformation deviations between different detection trajectories, an affine transformation is introduced. ,in This is the scaling factor. For rotation matrix, It is a translation vector.

[0109] The cost function consists of the weighted sum of the following three terms: The first term is the pose constraint, which physically means aligning matched keyframe pairs as much as possible in space. Since matched keyframes correspond to the same physical location within the crevice, the spatial distance between them should be minimized after optimization. This term is defined as the sum of the squares of the Euclidean distances of all matched keyframe pairs after affine transformation:

[0110] In the formula: This represents the attitude constraint loss. and The first For matching the position vectors of keyframes in the two trajectories. This represents a transformation matrix, used to transform the coordinate system of one trajectory to another.

[0111] The second term is the trajectory smoothing constraint, which physically aims to suppress abrupt changes in the local trajectory caused by fluid disturbances, ensuring the continuity and naturalness of the trajectory. This term is defined as the sum of the squares of the distances between adjacent trajectory points, specifically expressed as:

[0112] In the formula: This represents the trajectory smoothing constraint loss. Indicates the current trajectory's first... The position of each trajectory point. Indicates the current trajectory's first... The position of each trajectory point. This represents the total number of points on a single trajectory. In multi-trip detection scenarios, the smoothing term of all trajectories can be summed.

[0113] The third term is the rock wall penetration penalty, which is one of the key features that distinguishes this invention from traditional trajectory optimization methods. Its physical meaning is: the rock fissure wall is a solid medium, and trajectory points at different times are all located inside the rock fissure. The distance between any two trajectory points cannot be less than half the sum of the rock fissure widths at the corresponding times; otherwise, it means that the two trajectory points have "penetrated" the rock wall. Based on this physical fact, the penetration penalty term is defined as:

[0114] In the formula: Indicates rock wall penetration loss, and The first The and the first The width of the rock fissure measured by ultrasound at each trajectory point. Indicates the current trajectory's first... The position of each trajectory point. Indicates the current trajectory's first... The position of each trajectory point. This represents the total number of points on a single trajectory. When the distance between two trajectory points is less than half the sum of their corresponding widths, this term applies a penalty proportional to the square of the penetration depth; when the distance between two trajectory points is greater than or equal to this threshold, the penalty is zero.

[0115] The total cost function is a weighted sum of three terms:

[0116] In the formula: Represents the cost function, , , These are weighting coefficients used to balance the relative importance of various constraints. The attitude constraint has the highest weight to ensure proper alignment of matching keyframes; the rock wall penetration penalty is next to it to ensure that the optimization result conforms to physical laws; the smoothness constraint has a relatively low weight to ensure trajectory smoothness while avoiding over-smoothing that could lead to loss of detail.

[0117] Step 2: Based on the cost function, use an optimization method to optimize the non-strongly constrained position points in the initial trajectory to obtain the final trajectory.

[0118] In some embodiments, the optimization method employed includes particle swarm optimization, genetic algorithm, or L-BFGS-B algorithm. For example, the L-BFGS-B algorithm is used to minimize the aforementioned cost function to obtain the globally optimized trajectory.

[0119] The L-BFGS-B algorithm is a quasi-Newton method suitable for large-scale bounded constraint optimization problems. In each iteration, an approximate Hessian inverse matrix is ​​constructed using the gradient information of the objective function to determine the search direction, and then the optimal step size is determined through line search to ensure that the objective function decreases sufficiently. At the same time, the algorithm has a built-in variable boundary constraint handling mechanism, which can directly incorporate the boundary constraints of the trajectory point coordinates into the optimization process.

[0120] In this invention, the optimization variable is the three-dimensional spatial coordinates of all trajectory points. Let the total number of trajectory points be... Then the total number of optimization variables For a typical detection (detection time approximately 14 seconds, sampling rate 200Hz), Approximately 2800 points. It has approximately 8400 dimensions. The boundary constraints for each variable are: ,in The coordinates of the initial trajectory points, For the maximum allowable offset (preferably, =50cm), to prevent optimization from diverging to solutions that have no physical meaning.

[0121] The solution process of the L-BFGS-B algorithm is as follows: First, the initial relative trajectory obtained in step S102 is used as the initial value for iteration. In each iteration, the cost function value and its gradient vector at the current trajectory point coordinates are calculated; the finite memory BFGS update formula is used to construct an approximate Hessian inverse matrix from the current gradient, historical gradient difference, and variable difference information, and then the search direction is calculated; a line search is performed along the search direction to find the step size that satisfies the Wolfe condition; the variables are updated and coordinates exceeding the boundary constraints are clamped. The above process is iterated repeatedly until the convergence condition is met (such as the gradient norm being less than the preset tolerance, the variable change being less than the preset tolerance, or the maximum number of iterations of 200 steps being reached), and the optimized trajectory point coordinates are output.

[0122] In the above iteration process, the three constraints in the cost function each exert different forces: the attitude constraint creates an attractive force between matched keyframe pairs, pulling the trajectory towards the matching position and correcting cumulative drift; the smoothness constraint creates a "straightening" force between adjacent trajectory points, suppressing local abrupt changes caused by fluid disturbances and maintaining the overall smoothness of the trajectory; and the penetration penalty creates a repulsive force between trajectory point pairs that are too close together, pushing them away to a position that satisfies the rock wall intrusion constraint. These three forces work synergistically during the optimization process, through weighting coefficients. , , By adjusting the relative strength, the trajectory can ultimately achieve an optimal balance between matching accuracy, smoothness, and physical realism.

[0123] In other embodiments, if the cost function exhibits multimodal characteristics or poses a risk of local minima, Particle Swarm Optimization (PSO) or Genetic Algorithm (GA) can be used instead of the L-BFGS-B algorithm. PSO maintains a set of candidate solutions (particles), with each particle updating its velocity and position based on its own historical best position and the population's historical best position, searching in parallel within the solution space, thus possessing strong global search capabilities. Genetic Algorithm iteratively evolves the population through selection, crossover, and mutation operations, and is also applicable to non-convex optimization problems. Although both have higher computational costs (typically 5-10 times that of L-BFGS-B), they are suitable for scenarios where real-time performance is less critical and global optimality is more important.

[0124] See Figure 7 , Figure 7 The process of optimizing the movement trajectory is illustrated in detail.

[0125] S104. Determine the three-dimensional morphology of the rock fissure based on the final trajectory and the width data of the rock fissure.

[0126] The optimized final trajectory and width data are input into the 3D drawing terminal. The corresponding width is expanded to both sides based on the center line to generate point clouds of the rock fissure wall. After surface reconstruction, the 3D contour of the rock fissure is obtained.

[0127] The optimized trajectory represents the three-dimensional spatial position of the rock fissure centerline. Combined with the width data measured by the ultrasonic sensor, the complete rock fissure wall morphology can be reconstructed.

[0128] The optimized trajectory undergoes preprocessing, including removing obvious outliers and using median filtering for local smoothing to ensure trajectory continuity. For cases where the sampling point spacing in the trajectory is uneven, spline interpolation is used for uniform resampling, providing a regularized input for subsequent point cloud generation.

[0129] Calculate the local normal plane at each trajectory point. The normal plane is constructed as follows: First, calculate the tangent direction at that point, usually obtained by the difference between adjacent points; at the trajectory endpoints, one-sided difference is used. Then, construct two orthogonal basis vectors based on the tangent direction. One is obtained by the cross product of the tangent vector and a reference vector (e.g., the perpendicular direction), and the other is obtained by the cross product of the tangent vector and the first basis vector, thus forming a local coordinate system. During the calculation, it is necessary to handle degenerate cases such as the reference vector being parallel to the tangent direction; in this case, the reference vector needs to be switched.

[0130] In the local coordinate system, multiple directions are uniformly extended in a plane perpendicular to the tangent, centered on the trajectory point. The extension distance in each direction is half the width of the rock fissure at that point. The number of extension directions can be adjusted according to the required point cloud density; typically, extending in 8 to 12 directions per point is sufficient to achieve good wall coverage. For areas with drastic width variations, the number of extension directions can be appropriately increased.

[0131] The wall points generated from all trajectory points are merged to obtain the original point cloud data. The original point cloud may contain a small number of noisy points and outliers, which are cleaned using statistical filtering or radius filtering. For regions with uneven point cloud density, downsampling or upsampling can be performed to provide a uniformly distributed point set for subsequent reconstruction.

[0132] Point cloud reconstruction employs the Poisson surface reconstruction algorithm. This algorithm transforms the point cloud and its normal information into a continuous mesh model by solving the Poisson equation, effectively handling complex topologies and closed surfaces. Before reconstruction, the normal of each point needs to be estimated, typically using Principal Component Analysis (PCA) based on a local neighborhood-fitted plane. For special regions such as wall edges and branch cracks, manual adjustment of the normal direction may be necessary to ensure consistency.

[0133] The reconstructed 3D mesh model can undergo further post-processing, including hole filling, smoothing filtering, and simplification / reduction of surface areas. The final model can be exported to common 3D file formats (such as STL, OBJ, and PLY), facilitating import into other analysis software for geometric measurement and visualization.

[0134] Example 1: Three-dimensional mapping of rock fissures in the foundation of a hydropower station dam.

[0135] (1) Detection scenario and target: A through-rock fissure has developed in the foundation of a hydropower station dam. The fissure is about 15 meters long and its width varies from 1.2 to 18 centimeters. There are 3 bends (the smallest bend radius is about 8 centimeters). Some areas are filled with water due to seepage (the water depth is 2 to 5 centimeters). It is necessary to accurately map the three-dimensional outline of the fissure without damaging the dam structure to provide basic data for the seepage stability analysis of the dam.

[0136] (2) System configuration: Encapsulated detector: Utilizing a waterproof, pressure-resistant, and lightweight housing, approximately 1cm in diameter, several centimeters in length, and weighing only a few grams, it boasts IP68 waterproof rating and sufficient pressure resistance. Internally integrated: a miniature IMU (containing a three-axis accelerometer and a three-axis gyroscope, with a sampling rate of at least 100Hz), a high-sensitivity magnetic field sensor array (resolution up to nT levels), a miniature ultrasonic ranging sensor (accuracy better than 0.1mm), an optional flexible pressure sensor for sensing wall contact, and an ultra-thin flexible battery (providing several hours of battery life). The housing surface is designed with a streamlined shape to reduce drag.

[0137] Injection device: includes a high-pressure gas source and an inert fluid storage tank, as well as an injection conduit with an inner diameter of about 2 cm, and a pressure and flow regulating component at the front end.

[0138] Low-frequency electromagnetic repeater network: It consists of multiple miniature repeaters, with operating frequencies in the hundreds of kHz range. Each unit has an effective transmission distance of several meters. It is deployed at appropriate intervals along the rock fissure path, and the end is connected to the ground receiving equipment via wired connection.

[0139] Embedded computing module: Employs a low-power microcontroller with sufficient processing power and storage space to acquire and preprocess sensor data in real time.

[0140] (3) Experimental steps Step 1: Drive the encapsulated detector to move along the rock fissure and acquire motion data of the encapsulated detector, rock fissure magnetic field characteristic data and rock fissure width data collected during the movement.

[0141] Drill an injection hole approximately 2 cm in diameter at the entrance to the rock fissure, ensuring the injection conduit can be inserted tightly. Push multiple repeaters sequentially into the fissure through the injection conduit, spacing them approximately several meters apart. Activate the repeaters and test the communication link with the ground station to confirm that each node is functioning correctly.

[0142] Depending on the water accumulation within the rock fissures, a fluid injection mode is selected, and an appropriate amount of inert fluid is injected into the storage tank. Additionally, several nano-sized neodymium iron boron magnetic markers can be spaced out on or within the rock mass to construct artificial feature anchors.

[0143] The encapsulated detector was placed inside the sealed chamber at the front end of the injection catheter, kept stationary, and 1000 frames of static IMU data were acquired. The IMU data was then calibrated. The calibration method is as follows: Accelerometer zero bias calculation: With the sensor remaining absolutely stationary, acceleration data is collected over a period of time, and its zero bias is calculated. , , The average value on the three axes, with the average value divided by the mean. After calculating the gravitational acceleration in the direction of gravity, the accelerometer zero bias is obtained, and the calculation result is: m / s 2 .

[0144] Angular velocity zero bias calculation: With the sensor remaining absolutely stationary, angular velocity data is collected over a period of time, and its average value on the X, Y, and Z axes is calculated to obtain the gyroscope zero bias. The calculation result is as follows: rad / s.

[0145] After calculating the zero-bias values ​​of angular velocity and angular velocity, the zero-bias values ​​of angular velocity and angular velocity measured by the sensor are removed to complete the calibration.

[0146] Start the high-pressure air pump and adjust the output pressure to an appropriate value so that the detector begins to inject at a speed of approximately 1 cm / s. When the detector enters the water accumulation area, switch to fluid injection mode and fine-tune the pressure to maintain a stable speed.

[0147] During the detection process, the ultrasonic sensor monitored the width in front in real time. When the detected width was ≤1.2cm (i.e., 1.2 times the detector diameter), the injection device automatically reduced the pressure and the moving speed to 0.6cm / s. The encapsulated detector successfully passed through three curved sections without any jamming. When the encapsulated detector reached the end of the rock fissure (where the injection resistance suddenly increased), the injection automatically stopped. A total of approximately 2800 frames of IMU data, 2800 frames of magnetic field data, 140 frames of ultrasonic data, and 2800 frames of tactile data were collected.

[0148] Step 2: Based on the motion data of the encapsulated detector, calculate the initial movement trajectory of the encapsulated detector. The embedded computing module receives sensor data in real time and processes it according to the following algorithm: 1. Savitzky-Golay filtering of IMU data: Savitzky-Golay filter is used to smooth the IMU acceleration and angular velocity, effectively suppressing high-frequency noise introduced by fluid disturbance.

[0149] 2. The attitude of the encapsulated detector as it moves within the rock fissure is determined using the quaternion method based on the angular velocity data in the IMU data.

[0150] 3. Based on the attitude of the encapsulated detector and the acceleration data in the IMU data, determine the trajectory of the encapsulated detector in the rock crevice as the initial trajectory; zero-velocity detection is added during the initial trajectory calculation; calculate the acceleration modulus after denoising. and angular velocity modulus If satisfied and If the velocity is zero, it is determined to be in a stationary state, and the velocity is set to zero to reduce trajectory drift.

[0151] Step 3: Based on the magnetic field characteristic data collected by several encapsulated detectors, optimize several initial trajectories to obtain the final trajectory.

[0152] The magnetic field data was divided into keyframes at a fixed number of frames (corresponding to a movement distance of about 0.5cm), and an appropriate overlap rate was set to generate hundreds of keyframes.

[0153] The first step is to acquire the magnetic field feature data of the magnetic beacon. The magnetic field feature data of the magnetic beacon and the magnetic field feature data acquired by the encapsulated detector are divided into several key frames according to a fixed number of frames. The time domain features and frequency domain features of each key frame are extracted and fused into a feature vector. The gradient of the feature vector of several frames is then extracted to obtain the magnetic field gradient sequence of the magnetic beacon and the encapsulated detector.

[0154] The second step involves matching the magnetic field feature data of the magnetic marker with the magnetic field feature data collected by the encapsulated detector based on the magnetic field gradient sequence, thereby obtaining several target position matching pairs.

[0155] The third step is to construct and minimize a cost function that includes trajectory alignment error and trajectory smoothness, using the feature matching pairs as constraints, to achieve global optimization of the initial movement trajectory.

[0156] After optimization, the trajectory's endpoint drift was reduced from over ten centimeters initially to several centimeters, and the average relative error was also significantly reduced. The optimized trajectory is smoother in space, and the alignment accuracy with feature matching points is greatly improved.

[0157] Step 4: Based on the movement trajectory and the rock fissure width data, invert the three-dimensional contour of the rock fissure.

[0158] The optimized trajectory represents the three-dimensional spatial position of the rock fissure centerline. Combined with the width data measured by the ultrasonic sensor, the complete rock fissure wall morphology can be reconstructed.

[0159] The optimized trajectory undergoes preprocessing, including removing obvious outliers and using median filtering for local smoothing to ensure trajectory continuity. For cases where the sampling point spacing in the trajectory is uneven, spline interpolation is used for uniform resampling, providing a regularized input for subsequent point cloud generation.

[0160] Calculate the local normal plane at each trajectory point. The normal plane is constructed as follows: First, calculate the tangent direction at that point, usually obtained by the difference between adjacent points; at the trajectory endpoints, one-sided difference is used. Then, construct two orthogonal basis vectors based on the tangent direction. One is obtained by the cross product of the tangent vector and a reference vector (e.g., the perpendicular direction), and the other is obtained by the cross product of the tangent vector and the first basis vector, thus forming a local coordinate system. During the calculation, it is necessary to handle degenerate cases such as the reference vector being parallel to the tangent direction; in this case, the reference vector needs to be switched.

[0161] In the local coordinate system, multiple directions are uniformly extended in a plane perpendicular to the tangent, centered on the trajectory point. The extension distance in each direction is half the width of the rock fissure at that point. The number of extension directions can be adjusted according to the required point cloud density; typically, extending in 8 to 12 directions per point is sufficient to achieve good wall coverage. For areas with drastic width variations, the number of extension directions can be appropriately increased.

[0162] The wall points generated from all trajectory points are merged to obtain the original point cloud data. The original point cloud may contain a small number of noisy points and outliers, which are cleaned using statistical filtering or radius filtering. For regions with uneven point cloud density, downsampling or upsampling can be performed to provide a uniformly distributed point set for subsequent reconstruction.

[0163] Point cloud reconstruction employs the Poisson surface reconstruction algorithm. This algorithm transforms the point cloud and its normal information into a continuous mesh model by solving the Poisson equation, effectively handling complex topologies and closed surfaces. Before reconstruction, the normal of each point needs to be estimated, typically using Principal Component Analysis (PCA) based on a local neighborhood-fitted plane. For special regions such as wall edges and branch cracks, manual adjustment of the normal direction may be necessary to ensure consistency.

[0164] The reconstructed 3D mesh model can undergo further post-processing, including hole filling, smoothing filtering, and simplification / reduction of surface areas. The final model can be exported to common 3D file formats (such as STL, OBJ, and PLY), facilitating import into other analysis software for geometric measurement and visualization.

[0165] The model generated in this embodiment clearly shows the main orientation of the rock fissure (roughly north-northeast), with its width gradually increasing from about 1 cm at the entrance to nearly 20 cm at the end. The angles of the three curved sections are approximately 20 degrees, 30 degrees, and nearly 20 degrees, respectively. At a distance of about 8 meters from the entrance, the model clearly shows a branching fissure, about 2 meters long and with an average width of about 3.5 cm. This three-dimensional model matches well with the previous ground-penetrating radar detection results, verifying the effectiveness of the method of this invention.

[0166] (4) Effect verification To verify the mapping accuracy of this invention, several verification holes (3cm in diameter) were drilled at key locations after the exploration was completed, and the width and location of the rock fissures at the corresponding points were measured using borehole television.

[0167] Table 1: Results of Surveying Accuracy Verification

[0168] Table 1 shows the comparison results. The orientation error at the verification hole is less than 3°, the width error is less than 0.2cm, and the maximum error in the depth direction (calculated from the position of the injection hole) is 16cm, which is better than the design specifications.

[0169] Example 2: Long-term monitoring of rock fissures on highway slopes.

[0170] The strategies and methods provided in this application can also be applied to other fields besides the embodiments described above: in monitoring shale gas fracturing fractures in oil and gas fields, the detector is equipped with pressure and temperature sensors, and as fracturing fluid is injected into the fracture, the pressure and temperature changes during the fracture propagation process are monitored in real time; in geothermal resource exploration, the detector can measure the temperature gradient within the rock fracture; in urban underground space exploration, it can be used to detect the development of soil fractures around underground pipelines.

[0171] Finally, it should be noted that the above embodiments are only used to illustrate the technical solutions of the present invention, and not to limit them; although the present invention has been described in detail with reference to the foregoing embodiments, those skilled in the art should understand that modifications can still be made to the technical solutions described in the foregoing embodiments, or equivalent substitutions can be made to some of the technical features; and these modifications or substitutions do not cause the essence of the corresponding technical solutions to deviate from the spirit and scope of the technical solutions of the embodiments of the present invention.

Claims

1. A rock fracture three-dimensional trajectory inversion system based on multi-source sensor fusion, characterized in that, include: Several encapsulated detectors are provided, each encapsulating detector integrating a motion measurement unit, a magnetic field measurement unit, and a width measurement unit. The motion measurement unit is used to collect IMU data from the encapsulated detectors, the environmental measurement unit is used to collect magnetic field characteristic data of rock fissures, and the width measurement unit is used to collect width data of rock fissures. Several magnetic markers are placed at intervals on the surface or inside the rock mass; An injection device is used to inject fluid or gas into rock fissures, send several encapsulated detectors into the rock fissures, and drive several encapsulated detectors to move within the rock fissures. A relay network, set up on the rock surface, is used to transmit motion data, magnetic field characteristic data, and width data collected by the encapsulated detector; The mapping module is used to invert the three-dimensional morphology of the rock fissure based on IMU data, magnetic field characteristic data, and width data.

2. The rock fracture three-dimensional trajectory inversion system based on multi-source sensor fusion according to claim 1, characterized in that, The encapsulated detector also integrates a pressure measurement unit for collecting the pressure generated when the encapsulated detector comes into contact with the rock mass.

3. A method for inverting the three-dimensional trajectory of rock fissures based on multi-source sensor fusion, characterized in that, The application to the rock fracture three-dimensional trajectory inversion system based on multi-source sensor fusion as described in claim 1 or 2 includes: Several encapsulated detectors were sent into the rock crevice, and the detectors were driven to move within the crevice. The IMU data of the encapsulated detectors, the magnetic field characteristics of the rock crevice, and the width of the rock crevice were collected. The initial trajectory is determined based on IMU data collected by several encapsulated detectors; Based on the magnetic field characteristic data collected by several encapsulated detectors, the initial trajectory is optimized to obtain the final trajectory; Based on the final trajectory and the width data of the rock fissure, the three-dimensional morphology of the rock fissure is determined.

4. The method for inverting three-dimensional rock fracture trajectories based on multi-source sensor fusion according to claim 3, characterized in that, The steps for determining several initial trajectories based on IMU data acquired from several encapsulated detectors include: Savitzky-Golay filtering was applied to the IMU data. The attitude of the encapsulated detector as it moves within the rock crevice is determined using the quaternion method based on the angular velocity data in the IMU data. Based on the attitude of the encapsulated detector and the acceleration data in the IMU data, the trajectory of the encapsulated detector in the rock crevice is determined as the initial trajectory.

5. The method for inverting three-dimensional rock fracture trajectories based on multi-source sensor fusion according to claim 4, characterized in that, When determining the movement trajectory of the encapsulated detector in the rock crevice, the zero-velocity of the encapsulated detector can also be determined based on IMU data, and the speed of the encapsulated detector can be calibrated based on the zero-velocity determination.

6. The method for inverting three-dimensional rock fracture trajectories based on multi-source sensor fusion according to any one of claims 3 to 5, characterized in that, The process of optimizing the initial trajectory based on magnetic field characteristic data collected by several encapsulated detectors includes: Obtain the magnetic field characteristic data of the magnetic beacon, and determine several position constraint matching pairs based on the magnetic field characteristic data of the magnetic beacon and the magnetic field characteristic data of the encapsulated detector. Based on several position constraint matching pairs, strong constraints are applied to the position points in the initial trajectory that correspond to the target position matching pairs. A cost function is constructed, and an optimization method is used to optimize the non-strongly constrained position points in the initial trajectory based on the cost function to obtain the final trajectory.

7. The method for inverting three-dimensional rock fracture trajectories based on multi-source sensor fusion according to claim 6, characterized in that, Based on the magnetic field gradient sequence, the steps for determining several position constraint matching pairs according to the magnetic field characteristic data of the magnetic marker and the magnetic field characteristic data of the encapsulated detector include: The magnetic field feature data of the magnetic beacon and the magnetic field feature data collected by the encapsulated detector are divided into several key frames according to a fixed number of frames. The time domain features and frequency domain features of each key frame are extracted, and the time domain features and frequency domain features are fused into a feature vector. The gradient of the feature vector of several frames is extracted to obtain the magnetic field gradient sequence of the magnetic beacon and the encapsulated detector. Based on the magnetic field gradient sequence, the magnetic field feature data of the magnetic marker is matched with the magnetic field feature data collected by the encapsulated detector to obtain several candidate position matching pairs. Based on the pressure feature data collected by the encapsulated detector, several tactile events are determined to form a tactile signature sequence; the sampling time point of the magnetic field feature in the candidate position matching pair is determined, and tactile events under a preset time window before and after the sampling time point of the magnetic field feature in the candidate position matching pair are extracted to form a first local tactile subsequence and a second local tactile subsequence; based on the similarity between the first local tactile subsequence and the second local tactile subsequence, position constraint matching pairs are selected from the candidate position matching pairs.

8. The method for inverting three-dimensional rock fracture trajectories based on multi-source sensor fusion according to claim 7, characterized in that, The method, based on a magnetic field gradient sequence, involves matching the magnetic field feature data of a magnetic marker with the magnetic field feature data acquired by an encapsulated detector. The process includes: Calculate the cosine similarity between the frequency domain characteristics of the magnetic marker and the frequency domain characteristics of the encapsulated detector. Based on the cosine phase velocity, perform a coarse matching between the trajectory points of the encapsulated detector and the position of the magnetic marker to obtain a coarse matching pair. The magnetic field gradient sequences of several coarse matching pairs are compared using the dynamic time warping method or the longest common subsequence method to obtain candidate position matching pairs from the coarse matching pairs.

9. The method for inverting three-dimensional rock fracture trajectories based on multi-source sensor fusion according to claim 6, characterized in that, The cost function includes: In the formula: This represents the attitude constraint loss. This represents the trajectory smoothing loss. Indicates rock wall penetration loss, , , These are the weighting coefficients.

10. The method for inverting three-dimensional rock fracture trajectories based on multi-source sensor fusion according to claim 6, characterized in that, The optimization methods used include particle swarm optimization, genetic algorithm, or L-BFGS-B algorithm.

Images