Ice layer internal gap identification and positioning method based on airborne shallow ice penetrating radar

Abstract

The invention relates to the technical field of three-dimensional identification, in particular to an ice layer internal gap identification and positioning method based on an airborne shallow ice penetrating radar. The method comprises the following steps: carrying out coverage type flight detection on a target shallow layer area through a fixed-wing aircraft carrying an ice penetrating radar, and carrying out time-frequency enhancement and ice layer layering identification at the same time to obtain ice layer layering structure data; carrying out interlayer interface analysis and gap suspected area detection on the ice layer layered structure data to obtain ice layer internal gap suspected area data; performing three-dimensional shape geometric analysis and gap type classification on the ice layer internal gap suspected region data to obtain ice layer gap type classification data; based on the ice layer gap type classification data, gap distribution evolution fine reconstruction and visual positioning expression are carried out on the ice layer layered structure data, an ice layer internal gap identification and positioning report is generated, and the ice layer internal gap identification and positioning report is sent to the terminal. According to the invention, accurate identification and positioning of the internal gap of the ice layer can be realized.

Description

Technical Field

[0001] This invention relates to the field of three-dimensional recognition technology, and in particular to a method for identifying and locating internal gaps in ice layers based on airborne shallow ice-detecting radar. Background Technology

[0002] With the deepening of polar scientific research, ice layer detection technology has become an important tool for studying climate change, glacier dynamics, and their impact on sea-level changes. Airborne radar systems, especially shallow ice-penetrating radar, are widely used in ice layer detection because they can penetrate ice layers and obtain information about their internal structure. This technology mainly detects the thickness, density, porosity, and other characteristics of ice layers through electromagnetic wave reflection, providing a detailed description of the ice layer structure, which is particularly important for identifying anomalous structures within the ice layer. However, existing airborne shallow ice-penetrating radar technologies mostly focus on detecting the thickness, surface morphology, and larger internal structures of ice layers, such as the glacier base, the bottom of the ice layer, and geological layers. Detecting and locating tiny fissures or cracks within the ice layer remains a significant challenge. Because these fissures are often small, unevenly distributed, and have strong scattering characteristics, traditional radar systems are often affected by signal attenuation, noise interference, and insufficient resolution when acquiring reflected signals from these fissures, making it difficult to effectively identify and locate these tiny fissures. Summary of the Invention

[0003] Therefore, it is necessary for the present invention to provide a method for identifying and locating gaps inside ice layers based on airborne shallow ice-detecting radar, in order to solve at least one of the above-mentioned technical problems.

[0004] To achieve the above objectives, a method for identifying and locating internal ice gaps based on airborne shallow-layer ice-penetrating radar includes the following steps: Step S1: Use a fixed-wing aircraft equipped with ice-detecting radar to conduct a comprehensive flight detection of the shallow area of ​​the target to obtain multi-channel ice layer echo signal data; perform time-frequency enhancement and ice layer layer identification on the multi-channel ice layer echo signal data to obtain ice layer layer structure data; Step S2: Perform interlayer interface analysis on the ice layer layer structure data to obtain the internal interface feature data of the ice layer; based on the internal interface feature data of the ice layer, perform crack detection on the ice layer layer structure data to obtain crack suspected area data inside the ice layer. Step S3: Perform three-dimensional morphological geometric analysis on the suspected area of ​​the ice layer's internal fissures to obtain the three-dimensional morphological geometric parameter data of the ice layer fissures; classify the fissure types based on the three-dimensional morphological geometric parameter data of the ice layer fissures to obtain the ice layer fissure type classification data; Step S4: Based on the classification data of ice layer gap types, perform fine reconstruction of the gap distribution evolution of ice layer layer structure data to generate a three-dimensional fine model of the gaps inside the ice layer; based on the three-dimensional fine model of the gaps inside the ice layer, perform visualization and localization expression and generate a report on the identification and localization of the gaps inside the ice layer, and send the report on the identification and localization of the gaps inside the ice layer to the terminal.

[0005] Furthermore, step S1 includes the following steps: Step S11: Use a fixed-wing aircraft equipped with ice-detecting radar to conduct a comprehensive flight detection of the shallow area of ​​the target to obtain multi-channel ice layer echo signal data; Step S12: Perform spatiotemporal synchronization and noise suppression processing on the multi-channel ice layer echo signal data to generate noise-suppressed ice layer echo signal data; Step S13: Perform time-frequency feature enhancement processing on the noise-suppressed ice echo signal data to obtain enhanced ice echo signal data; Step S14: Perform wavelet transform multi-scale decomposition on the enhanced ice layer echo signal data to obtain the ice layer echo coefficient data at different scales; perform multi-scale feature extraction on the ice layer echo coefficient data at different scales to extract time-domain features, frequency-domain features, and time-frequency joint features respectively, and perform feature dimensionality reduction and fusion through principal component analysis to obtain multi-scale feature data corresponding to the internal structure of the ice layer. Step S15: Based on multi-scale feature data, a deep belief network combined with a conditional random field model is used to construct an ice layer layer identification model, and the ice layer layer identification model is used to identify the ice layer layer structure data of the multi-scale feature data.

[0006] Furthermore, step S2 includes the following steps: Step S21: Identify ice layer interfaces in the ice layer layer structure data to obtain a set of ice layer layer structure interfaces; Step S22: Perform interlayer interface identification analysis on the different interlayer interfaces in the ice layer within the ice layer layer structure profile interface set to analyze the physical, chemical and structural characteristics between different interlayer interfaces, including interface thickness and particle distribution, and obtain the ice layer interlayer interface feature identification set. Step S23: Based on the ice layer interlayer interface feature identification set, perform spatiotemporal evolution analysis of ice layer interfaces between different interlayer interfaces in the ice layer layer structure profile interface set, measure the ice layer depth annual rings between different interlayer interfaces to infer the generation, development and degradation laws of different ice layer interfaces over time, and generate spatiotemporal evolution simulation data corresponding to different interlayer interfaces inside the ice layer. Step S24: Based on the spatiotemporal evolution simulation data corresponding to different interlayer profiles within the ice layer, perform hierarchical clustering of the interfaces between different interlayer profiles in the ice layer within the ice layer layer structure profile set to generate an internal interface hierarchical classification dataset of the ice layer; based on the internal interface hierarchical classification dataset of the ice layer, perform interlayer interface feature coupling analysis on the interface features between different interlayer profiles within the ice layer interlayer interface feature recognition set to obtain internal interface feature data of the ice layer. Step S25: Detect suspected gap regions in the ice layer layer structure data based on the internal interface feature data of the ice layer to obtain suspected gap region data inside the ice layer.

[0007] Furthermore, step S21 includes the following steps: The layer interface features of the ice layer structure data are analyzed to obtain the density change, temperature gradient and pressure change values ​​corresponding to each layer interface in the ice layer. Based on the density change, temperature gradient and pressure change values ​​corresponding to each layer interface in the ice layer, the abrupt boundary points between the layers are identified. The boundary point set between each layer interface in the ice layer is identified by setting the abrupt threshold corresponding to density, temperature gradient and pressure. Based on the set of boundary demarcation points between the interfaces of each layer within the ice layer, the ice layer layer structure data is used to identify and calibrate the ice layer interface to obtain the ice layer layer structure interface set.

[0008] Furthermore, step S25 includes the following steps: Step S251: Divide the ice layer layer structure data into small units inside the ice layer to obtain the small structural units inside each ice layer; Step S252: Based on the internal interface feature data of the ice layer, perform feature quantization between the micro-structural units of each ice layer to quantify and calculate the interface curvature, interface roughness and interface reflectivity between the micro-units of each ice layer. Step S253: Based on the interface curvature, interface roughness and interface reflectivity between the micro units inside each ice layer, evaluate the interface jump of the micro structural units inside each ice layer to obtain the order of the interface jump between the micro units inside each ice layer. Step S254: Based on the transition order of the internal interface between the micro-units inside each ice layer, perform crack detection on the ice layer layer structure data to obtain crack suspected region data inside the ice layer.

[0009] Furthermore, step S253 includes the following steps: Interface distribution coupling analysis is performed based on the interface curvature, interface roughness, and interface reflectivity between the micro-units inside each ice layer to generate the ice layer interface coupling distribution field between the micro-units inside each ice layer. The interface jump barrier is calculated for the coupling distribution field of the ice layer interface between the micro units inside each ice layer, so as to calculate the height of the interface jump barrier at the interface of the micro units inside different ice layers and obtain the interface jump barrier matrix inside the ice layer. By identifying the abrupt change singularity of the interface jump barrier matrix inside the ice layer, the abrupt change singularity of the ice layer interface jump barrier between each small unit inside the ice layer is obtained. By obtaining the intensity and distribution density of the ice layer interface jump singularity corresponding to the abrupt change of the potential barrier between the small units inside each ice layer, and by evaluating the interface jump of the small structural units inside each ice layer based on the intensity and distribution density of the ice layer interface jump singularity, the order of the ice layer interface jump between the small units inside each ice layer can be obtained.

[0010] Furthermore, step S254 includes the following steps: The probability density of ice gaps is estimated based on the order of the interface transition between the micro-units within each ice layer, so as to obtain the probability density of ice gaps between the micro-units within each ice layer. Based on the probability density of gaps in the ice layer between the micro-units inside each ice layer, the ice layer layer structure data is used to detect suspected gap regions. Then, the suspected gap regions in the ice layer between the micro-units inside each ice layer are divided by adaptive threshold segmentation to obtain the initial suspected gap region data inside the ice layer. Morphological processing was performed on the initial suspected crack region data inside the ice layer to eliminate false crack regions using opening and closing operations, and region filtering was performed in combination with geometric constraints to obtain the suspected crack region data inside the ice layer.

[0011] Furthermore, step S3 includes the following steps: Step S31: Obtain ice echo signal data corresponding to the adjacent flight path of the fixed-wing aircraft in the shallow region of the target by using the suspected area data of the gap inside the ice layer; Step S32: Based on the ice echo signal data corresponding to adjacent tracks, perform three-dimensional point cloud registration on the suspected area data of the ice gap inside the ice layer to generate three-dimensional point cloud data of the ice gap inside the ice layer. Step S33: Denoise and smooth the three-dimensional point cloud data of the ice layer fissures by using statistical outlier removal and bilateral filtering and Laplacian smoothing to obtain the preprocessed three-dimensional point cloud data of the ice layer fissures. Step S34: Perform 3D reconstruction of the ice fissure surface on the preprocessed 3D point cloud data of the ice fissure. Based on the Poisson surface reconstruction algorithm, construct the corresponding 3D surface model of the ice fissure through normal vector estimation and signed distance function. Calculate the 3D morphological geometric parameters of the ice fissure surface model to extract the length, width, depth, area, and volume parameters corresponding to the ice fissure, so as to obtain the 3D morphological geometric parameter data of the ice fissure. Step S35: Classify the crack types based on the three-dimensional morphological geometric parameter data of the ice layer cracks to obtain ice layer crack type classification data.

[0012] Furthermore, step S32 includes the following steps: Spatiotemporal alignment is performed on the ice echo signal data corresponding to adjacent tracks to calculate the time difference and spatial offset of adjacent tracks in the same detection area, and generate track spatiotemporal alignment parameters; based on the track spatiotemporal alignment parameters, coordinate mapping alignment processing is performed on the ice echo signal data corresponding to adjacent tracks to obtain spatiotemporally aligned ice echo signal data; Time-frequency domain cross-correlation analysis was performed on the spatiotemporally aligned ice layer echo signal data to calculate the similarity of signal amplitude, phase and frequency at the same spatial location from different perspectives, generate ice layer echo signal similarity matrix, and filter out corresponding matching signal point pairs by setting a similarity threshold to obtain initial signal matching point pair data; Geometric constraint verification was performed on the initial signal matching point pair data. An epipolar constraint model was constructed based on the radar detection angle and baseline distance of adjacent tracks. Mismatched point pairs that did not conform to the epipolar geometric relationship were eliminated. The random sampling consensus algorithm was used to iteratively optimize the remaining matching point pairs. The matching point set with the highest proportion of interior points was retained to obtain the optimized signal matching point pair data. Based on the optimized signal matching point pair data, three-dimensional point cloud registration is performed on the suspected area of ​​the ice layer internal fissure. The two-dimensional echo signal coordinates are converted into three-dimensional spatial coordinates based on the signal matching point pairs. The overlapping area corresponding to the adjacent track point clouds is selected as the reference benchmark to calculate the transformation matrix between the point clouds. At the same time, the corresponding ice layer internal fissure point cloud is iteratively updated by minimizing the distance error from the point to the plane, so as to obtain the three-dimensional point cloud data of the ice layer internal fissure.

[0013] Furthermore, step S4 includes the following steps: Step S41: Perform spatial distribution analysis on the ice crevice type classification data to calculate the spatial density, directional distribution and clustering characteristics of each type of ice crevice, and obtain the spatial distribution feature vector of each type of ice crevice. Step S42: Based on the spatial distribution feature vectors corresponding to each type of gap, predict the gap distribution evolution trend of the ice layer structure data. Use the spatial distribution feature vectors to train the spatiotemporal autoregressive moving average model and combine it with Monte Carlo simulation. At the same time, consider the evolution of ice flow velocity and temperature change to generate corresponding ice layer gap evolution trend data. Step S43: Based on the ice layer gap evolution trend data, perform three-dimensional reconstruction of the ice layer gap evolution of each sub-layer interface in the ice layer layer structure data to generate a three-dimensional fine model of the gap inside the ice layer. Step S44: Perform 3D visualization of the fine 3D model of the fissures inside the ice layer to obtain 3D visualization positioning data of the fissures inside the ice layer; Step S45: Generate an internal ice gap identification and positioning report based on the three-dimensional visualization positioning data of the internal ice gap, and send the internal ice gap identification and positioning report to the terminal.

[0014] The beneficial effects of this invention are: The ice layer internal gap identification and localization method based on airborne shallow ice-detecting radar proposed in this invention has the following advantages compared with the prior art: By using a fixed-wing aircraft equipped with ice-detecting radar to conduct comprehensive flight detection of the target shallow area, it can widely cover the target area and acquire multi-channel ice layer echo signal data. This method has high detection efficiency and is particularly suitable for detecting large-scale ice layers. The ice-detecting radar, through the acquisition of echo signals in different frequency bands, can perform detailed analysis of the internal structure of the ice layer from multiple angles and dimensions, thereby obtaining multi-channel ice layer echo data. Using time-frequency enhancement technology to process the echo signals can effectively improve the quality of the echo signals, making the details in the signals clearer, thus improving the accuracy of subsequent analysis. The ice layer layer identification process further analyzes the echo signals to distinguish different layers of the ice layer, obtaining layered structure data. This layered structure data not only reflects the basic characteristics of the ice layer but also provides strong data support for subsequent gap detection and structural analysis, laying the foundation for in-depth research on the internal physical properties of the ice layer. Secondly, by analyzing the interlayer interfaces of the ice layer's layered structure data, the interrelationships and transitional characteristics between different ice layers can be revealed, thereby extracting characteristic data of the internal interfaces of the ice layer. This analysis not only helps to further understand the layered structure of the ice layer but also enables the precise identification of differences in material characteristics between layers. For example, ice layers at different depths have different densities, ice particle structures, or bubble distributions. Through interface feature analysis, interlayer interface discontinuities or potential crack areas can be effectively identified. Based on these interface feature data, further detection of suspected crack areas is possible. The key to this step is that by detecting suspected crack areas, weak points within the ice layer can be identified in advance, providing target areas for subsequent crack type classification and three-dimensional geometric analysis. This detection method effectively utilizes the layered structure data of the ice layer, reduces human intervention, and can efficiently discover potential unstable areas of the ice layer, preventing risks caused by ice layer collapse or fracturing, thus enabling better and more effective identification and location of these tiny cracks. Then, by performing three-dimensional morphological geometric analysis on the suspected areas of fissures inside the ice layer, three-dimensional morphological parameter data of the fissures can be obtained. These geometric parameters include the depth, width, length, and shape of the fissures. By quantifying these parameters, the spatial distribution and morphological characteristics of fissures inside the ice layer can be more clearly understood. Three-dimensional morphological geometric analysis not only helps researchers accurately understand the spatial structure of fissures inside the ice layer, but also identifies the distribution trend of fissures, providing a scientific basis for subsequent fissure type classification. Based on this geometric data, further classification of fissure types can not only distinguish different types of micro-fissures (such as bubble fissures, cracks, pores, etc.), but also provide a reference for ice layer stability assessment.Finally, a refined reconstruction of the ice layer structure data based on ice fracture type classification data was performed to study the fracture distribution evolution. The core of this step lies in generating a detailed 3D model of the fractures within the ice layer by meticulously modeling the evolution of fracture distribution. This detailed 3D model accurately represents the spatial distribution and evolution of fractures within the ice layer, providing a dynamic reference for long-term ice layer monitoring. The generation of this model helps improve the accuracy of ice layer stability and safety assessments, especially under the influence of climate change or geological movements, where changes in internal ice fractures can significantly impact ice layer stability. Based on the visualized location representation of the fractures using the 3D fracture model, researchers can intuitively understand the fracture distribution within the ice layer, enabling them to make corresponding judgments and decisions. By promptly sending fracture identification and location reports to terminals, relevant decision-makers can quickly obtain information and take appropriate countermeasures, thereby effectively reducing potential safety risks. Attached Figure Description

[0015] Other features, objects, and advantages of the invention will become more apparent from the following detailed description of non-limiting embodiments with reference to the accompanying drawings: Figure 1 This is a schematic diagram of the steps of the method for identifying and locating gaps inside ice layers based on airborne shallow ice-detecting radar according to the present invention. Figure 2 for Figure 1 A detailed flowchart of step S1; Figure 3 for Figure 1 A detailed flowchart of step S2. Detailed Implementation

[0016] The technical system of the present invention will now be clearly and completely described with reference to the accompanying drawings. Obviously, the described embodiments are only some, not all, of the embodiments of the present invention. All other embodiments obtained by those skilled in the art based on the embodiments of the present invention without inventive effort are within the scope of protection of the present invention.

[0017] Furthermore, the accompanying drawings are merely illustrative of the invention and are not necessarily drawn to scale. The same reference numerals in the drawings denote the same or similar parts, and therefore repeated descriptions of them will be omitted. Some block diagrams shown in the drawings are functional entities and do not necessarily correspond to physically or logically independent entities. These functional entities can be implemented in software, in one or more hardware steps or integrated circuits, or in different network and / or processor systems and / or microcontroller systems.

[0018] It should be understood that although the terms "first," "second," etc., may be used herein to describe various units, these units should not be limited by these terms. These terms are used merely to distinguish one unit from another. For example, without departing from the scope of the exemplary embodiments, a first unit may be referred to as a second unit, and similarly, a second unit may be referred to as a first unit. The term "and / or" as used herein includes any and all combinations of one or more of the associated listed items.

[0019] To achieve the above objectives, please refer to Figures 1 to 3 This invention provides a method for identifying and locating internal ice gaps based on airborne shallow-layer ice-penetrating radar. The system includes the following steps: Step S1: Use a fixed-wing aircraft equipped with ice-detecting radar to conduct a comprehensive flight detection of the shallow area of ​​the target to obtain multi-channel ice layer echo signal data; perform time-frequency enhancement and ice layer layer identification on the multi-channel ice layer echo signal data to obtain ice layer layer structure data; Step S2: Perform interlayer interface analysis on the ice layer layer structure data to obtain the internal interface feature data of the ice layer; based on the internal interface feature data of the ice layer, perform crack detection on the ice layer layer structure data to obtain crack suspected area data inside the ice layer. Step S3: Perform three-dimensional morphological geometric analysis on the suspected area of ​​the ice layer's internal fissures to obtain the three-dimensional morphological geometric parameter data of the ice layer fissures; classify the fissure types based on the three-dimensional morphological geometric parameter data of the ice layer fissures to obtain the ice layer fissure type classification data; Step S4: Based on the classification data of ice layer gap types, perform fine reconstruction of the gap distribution evolution of ice layer layer structure data to generate a three-dimensional fine model of the gaps inside the ice layer; based on the three-dimensional fine model of the gaps inside the ice layer, perform visualization and localization expression and generate a report on the identification and localization of the gaps inside the ice layer, and send the report on the identification and localization of the gaps inside the ice layer to the terminal.

[0020] In the embodiments of this invention, please refer to Figure 1 The diagram shown is a flowchart illustrating the steps of the ice layer internal crevice identification and location method based on airborne shallow ice-penetrating radar according to the present invention. In this example, the ice layer internal crevice identification and location method based on airborne shallow ice-penetrating radar includes the following steps: Step S1: Use a fixed-wing aircraft equipped with ice-detecting radar to conduct a comprehensive flight detection of the shallow area of ​​the target to obtain multi-channel ice layer echo signal data; perform time-frequency enhancement and ice layer layer identification on the multi-channel ice layer echo signal data to obtain ice layer layer structure data; In this embodiment of the invention, a fixed-wing aircraft is equipped with an ice-detecting radar with a center frequency of 500MHz. The radar has an 8-channel phased array structure, with a channel spacing of 0.5 meters, a pulse width of 0.1 microseconds, and a pulse repetition frequency of 10kHz. The aircraft conducts a coverage flight over a 10 square kilometer shallow subsurface area (X=0-1000 meters, Y=0-100 meters) at an altitude of 300 meters and a speed of 120 km / h. The flight paths are parallel, spaced 50 meters apart, and each path is 1000 meters long. During detection, data is collected simultaneously from all 8 channels at a sampling rate of 200MHz, with each channel generating 2×102 data per second. 8 Each 16-bit binary sampling point records the echo amplitude and phase, continuously acquiring multi-channel ice layer echo signal data for 5 hours. Time-frequency enhancement is performed on the data: a short-time Fourier transform (256-point window, Hanning window, 50% overlap) is used to amplify components in the time-frequency matrix below 30% of the maximum amplitude by a factor of 2.5. Ice layer layer identification is performed using a deep belief network combined with a conditional random field model. The model contains three layers of restricted Boltzmann machines (nodes 32, 16, and 8). The linear chain conditional random field defines the transition probability between adjacent pixels (0.9 for the same layer, 0.1 for boundaries). After 500 training iterations, the loss function is reduced to 0.01, identifying the surface layer (0-5 meters), middle layer (5-15 meters), and bottom layer (15-20 meters), obtaining ice layer layer structure data including the three-dimensional coordinates and echo characteristics of each layer.

[0021] Step S2: Perform interlayer interface analysis on the ice layer layer structure data to obtain the internal interface feature data of the ice layer; based on the internal interface feature data of the ice layer, perform crack detection on the ice layer layer structure data to obtain crack suspected area data inside the ice layer. In this embodiment of the invention, interlayer interface analysis is performed on the ice layer layer structure data to calculate the characteristic parameters of the surface-middle layer (5 meters) and middle-bottom layer (15 meters) interfaces: the surface-middle layer interface has a thickness of 0.2 meters, ice crystal particle diameter of 0.1-0.3 mm (average 0.2 mm), density change of -0.07 g / cm³, and temperature gradient change of +0.2℃ / m; the middle-bottom layer interface has a thickness of 0.3 meters, particle diameter of 0.3-0.5 mm (average 0.4 mm), density change of +0.05 g / cm³, and temperature gradient change of +0.2℃ / m, forming the internal interface characteristic data of the ice layer. Based on these data, suspected gap areas are detected, and thresholds are set: the surface-middle layer interface has a particle density of <100 particles / mm² and a hardness of <2.3 MPa, and the middle-bottom layer interface has a particle density of <70 particles / mm² and a hardness of <3.0 MPa. Sampling was performed at 1m x 1m intervals along the interface. In the surface-middle layer (X=8-10m, Y=8-10m), the particle density was 95 particles / mm², and the hardness was 2.2 MPa. In the middle-bottom layer (X=15-17m, Y=15-17m), the particle density was 65 particles / mm², and the hardness was 2.9 MPa, both meeting the thresholds. The continuity of the verification area (area ≥ 9 square meters) was verified. Once passed, data on suspected fissure areas within the ice layer were obtained, and the three-dimensional coordinate range of the area was recorded, such as the surface-middle layer area (X=8-10m, Y=8-10m, Z=5m).

[0022] Step S3: Perform three-dimensional morphological geometric analysis on the suspected area of ​​the ice layer's internal fissures to obtain the three-dimensional morphological geometric parameter data of the ice layer fissures; classify the fissure types based on the three-dimensional morphological geometric parameter data of the ice layer fissures to obtain the ice layer fissure type classification data; In this embodiment of the invention, three-dimensional morphological geometric analysis is performed on the data of suspected areas with fissures inside the ice layer. Through spatiotemporal alignment of adjacent track echo signals (time difference 2 seconds, Y-axis offset 50 meters), 1200 matching point pairs are screened using time-frequency cross-correlation (similarity ≥ 0.8). After epipolar constraint (distance < 0.5 meters) and optimization using a random sampling consensus algorithm (interior point ratio 92%), triangulation (Z = (3 × 10⁻⁶) is used. 8(×Δt) / 2) generates a 3D point cloud (1200 points, accuracy 0.1m). After denoising (statistical outlier removal, 20-neighborhood, threshold 2.5σ) and smoothing (bilateral filtering, spatial 0.3m, value range 0.1m), a Poisson surface reconstruction (normal vector based on 8-neighborhood, voxel 0.4m) is used to generate a surface model. Calculated parameters: large horizontal fissure, length 115m (X=505-620m), width 50m (Y=5-55m), depth 2.5m (Z=1-3.5m), area 600m², volume 100m³, obtaining the 3D morphological geometric parameter data of the ice layer fissure. Based on parameter classification: length > 50m, width > 30m, depth > 1m, and volume > 50m³ are classified as large horizontal fissures. This fissure meets the standard, obtaining ice layer fissure type classification data, including type labels and parameters.

[0023] Step S4: Based on the classification data of ice layer gap types, perform fine reconstruction of the gap distribution evolution of ice layer layer structure data to generate a three-dimensional fine model of the gaps inside the ice layer; based on the three-dimensional fine model of the gaps inside the ice layer, perform visualization and localization expression and generate a report on the identification and localization of the gaps inside the ice layer, and send the report on the identification and localization of the gaps inside the ice layer to the terminal.

[0024] In this embodiment of the invention, the evolution of the fissure distribution is finely reconstructed based on the classification data of ice fissure types. A spatial distribution feature vector of large horizontal fissures [0.1 locations / km², 0°, 1] (density, direction, clustering) is calculated. The evolution over the next 5 years is predicted using a spatiotemporal autoregressive moving average model (ARIMA(2,1,2)) combined with Monte Carlo simulation (1000 samplings): length 115 meters, width 50 meters, depth 2.5 meters. An improved Poisson surface reconstruction (10-neighborhood normal vector, 0.4-meter voxel) is used to generate a fine 3D model of the internal ice fissures, containing 20,000 triangular facets with an edge error ≤0.2 meters, fitting the layered interface (deviation ≤0.3 meters). Visualization is achieved using volume rendering, with depth-based coloring (1-2 meters blue, 2-3 meters green, 3-3.5 meters red), transparency 0.3-0.6, 3D coordinate axis scale 10 meters, and multi-view rendering (top view, side view, cross-section, 1920×1080 pixels). Generate an identification and positioning report, including an introduction (region and method), a parameter table (10 parameters), an evolution diagram, a 3D diagram, and conclusions (key areas X=550-570 meters, Y=20-30 meters), in PDF format (Arial font, size 12, figures 300dpi). Transmit the report (10MB, checksum verification) to the terminal via TCP wired connection, store it automatically, and complete the process.

[0025] Furthermore, step S1 includes the following steps: Step S11: Use a fixed-wing aircraft equipped with ice-detecting radar to conduct a comprehensive flight detection of the shallow area of ​​the target to obtain multi-channel ice layer echo signal data; Step S12: Perform spatiotemporal synchronization and noise suppression processing on the multi-channel ice layer echo signal data to generate noise-suppressed ice layer echo signal data; Step S13: Perform time-frequency feature enhancement processing on the noise-suppressed ice echo signal data to obtain enhanced ice echo signal data; Step S14: Perform wavelet transform multi-scale decomposition on the enhanced ice layer echo signal data to obtain the ice layer echo coefficient data at different scales; perform multi-scale feature extraction on the ice layer echo coefficient data at different scales to extract time-domain features, frequency-domain features, and time-frequency joint features respectively, and perform feature dimensionality reduction and fusion through principal component analysis to obtain multi-scale feature data corresponding to the internal structure of the ice layer. Step S15: Based on multi-scale feature data, a deep belief network combined with a conditional random field model is used to construct an ice layer layer identification model, and the ice layer layer identification model is used to identify the ice layer layer structure data of the multi-scale feature data.

[0026] As an embodiment of the present invention, reference is made to... Figure 2 As shown, Figure 1 A detailed flowchart of step S1 is shown below. In this embodiment, step S1 includes the following steps: Step S11: Use a fixed-wing aircraft equipped with ice-detecting radar to conduct a comprehensive flight detection of the shallow area of ​​the target to obtain multi-channel ice layer echo signal data; In this embodiment of the invention, a fixed-wing aircraft is equipped with an ice-detecting radar with a center frequency of 500MHz. The radar antenna has an 8-channel phased array structure with a channel spacing of 0.5 meters, a transmitted pulse width of 0.1 microseconds, and a pulse repetition frequency of 10kHz. The aircraft conducts a comprehensive flight detection over a shallow target area (10 square kilometers) at an altitude of 300 meters and a speed of 120 km / h. The flight path is a parallel route with a spacing of 50 meters between each route, and each route is 10 kilometers long. The radar transmits electromagnetic waves that penetrate the ice layer and receives the reflected echoes. Data is collected simultaneously through all 8 channels at a sampling rate of 200MHz, with each channel generating 2×102 data per second. 8 Each sampling point uses 16-bit binary data to record the amplitude and phase information of the echo signal. The detection lasted for 5 hours, acquiring multi-channel ice layer echo signal data, including reflection characteristics of the ice surface, internal layers, and substrate. The echoes from the crevice region showed a sudden drop in amplitude and abrupt phase change.

[0027] Step S12: Perform spatiotemporal synchronization and noise suppression processing on the multi-channel ice layer echo signal data to generate noise-suppressed ice layer echo signal data; In this embodiment of the invention, spatiotemporal synchronization processing is performed on multi-channel ice echo signal data. Using the fourth channel out of eight as a reference, the data of each channel is aligned using timestamps to correct the 0.2 microsecond time offset caused by antenna position differences, ensuring that the echo signals from the same detection point coincide on the time axis. An adaptive noise cancellation algorithm is used to suppress noise. A reference noise sequence is constructed (1000 consecutive sampling points in the signal blank area), and the noise mean is calculated to be 0.02V with a variance of 0.001V². A signal mean subtraction operation is performed on each sampling point, and then a 5th-order Butterworth low-pass filter (cutoff frequency 50MHz) is used to filter high-frequency noise. After processing, the signal-to-noise ratio of the ice surface echo signal is increased from 15dB to 30dB. For example, the original signal amplitude at a sampling point is 0.82V (including 0.15V noise), and after processing, the amplitude is 0.67V, and the noise is reduced to 0.03V, generating noise-suppressed ice echo signal data while preserving the characteristic abrupt changes in the gap region.

[0028] Step S13: Perform time-frequency feature enhancement processing on the noise-suppressed ice echo signal data to obtain enhanced ice echo signal data; In this embodiment of the invention, time-frequency feature enhancement processing is performed on the noise-suppressed ice layer echo signal data. A short-time Fourier transform is used, with a sliding window length of 256 sampling points, a Hanning window function, and a 50% overlap rate. The time-frequency matrix for each window is calculated (time resolution 0.64 microseconds, frequency resolution 781.25kHz). Threshold enhancement is performed on the time-frequency matrix, setting the amplitude threshold to 30% of the maximum amplitude. Frequency components below the threshold are multiplied by 2.5 times, while components above the threshold remain unchanged, enhancing the low-frequency attenuation characteristics of the gap region. For example, in the time-frequency matrix corresponding to a certain gap location, the amplitude in the 50-100kHz frequency band is 0.12V (below the threshold of 0.15V), which is amplified to 0.3V, increasing the time-frequency difference between the gap and normal ice layer from 0.1V to 0.28V. The resulting enhanced ice layer echo signal data has more significant time-frequency characteristics of the internal gaps, providing clear input for subsequent decomposition.

[0029] Step S14: Perform wavelet transform multi-scale decomposition on the enhanced ice layer echo signal data to obtain the ice layer echo coefficient data at different scales; perform multi-scale feature extraction on the ice layer echo coefficient data at different scales to extract time-domain features, frequency-domain features, and time-frequency joint features respectively, and perform feature dimensionality reduction and fusion through principal component analysis to obtain multi-scale feature data corresponding to the internal structure of the ice layer. In this embodiment of the invention, wavelet transform multi-scale decomposition is performed on the enhanced ice echo signal data. The db4 wavelet basis is selected and the decomposition layer is 5 layers to obtain ice echo coefficient data at scale 1 (25-50MHz) to scale 5 (0.78-1.56MHz). Extracting time-domain features: Calculating the mean (e.g., mean coefficient of scale 3, 0.23V), variance (0.05V²), and peak factor (3.2) of the coefficients at each scale; extracting frequency-domain features: obtaining the center frequency (center frequency of scale 2, 25MHz), bandwidth (5MHz), and energy percentage (18%) of each scale through Fourier transform; extracting time-frequency joint features: calculating wavelet entropy (entropy value of scale 4, 0.82) and cross-correlation coefficient (correlation coefficient between scale 1 and scale 5, 0.15) of the data in 24 dimensions of the three types of features, performing principal component analysis, selecting the top 8 principal components with a cumulative contribution rate of 95%, and fusing them to obtain multi-scale feature data. The first principal component (weight 35%) contains time-domain peak and frequency-domain energy information, and the second principal component (weight 25%) contains time-frequency entropy and correlation coefficient information, clearly reflecting the differences in the internal structure of the ice layer, and finally obtaining the multi-scale feature data corresponding to the internal structure of the ice layer.

[0030] Step S15: Based on multi-scale feature data, a deep belief network combined with a conditional random field model is used to construct an ice layer layer identification model, and the ice layer layer identification model is used to identify the ice layer layer structure data of the multi-scale feature data.

[0031] In this embodiment of the invention, an ice layer layer recognition model is constructed based on multi-scale feature data. The deep belief network contains three Restricted Boltzmann Machines (RBMs). The first RBM input is 8-dimensional principal component features, with 32 hidden layer nodes; the second RBM input is 32-dimensional features, with 16 hidden layer nodes; and the third RBM outputs 8-dimensional abstract features. The Conditional Random Field model is a linear chain structure, using the RBM output features as input, and defining the state transition probabilities of adjacent pixels (0.9 within the same layer, 0.1 at the layer boundary). The model was trained using 1000 sets of labeled data (including 500 fissures and 1000 layer boundaries) and iterated 500 times until the loss function converged to 0.01. Multi-scale feature data was input into the model, and a deep belief network was used to extract low-energy features of the fissure region. Conditional random fields were used to correct the boundary ambiguity problem. The ice layer was identified as being divided into three layers (surface layer 0-5 meters, middle layer 5-15 meters, and bottom layer 15-20 meters). There are two fissures in the middle layer, located at coordinates (3000 meters, 2000 meters) at a depth of 8-9 meters and (7000 meters, 5000 meters) at a depth of 12-13 meters, respectively. The ice layer layer structure data was generated, and the fissure location error was less than 0.5 meters.

[0032] Furthermore, step S2 includes the following steps: Step S21: Identify ice layer interfaces in the ice layer layer structure data to obtain a set of ice layer layer structure interfaces; Step S22: Perform interlayer interface identification analysis on the different interlayer interfaces in the ice layer within the ice layer layer structure profile interface set to analyze the physical, chemical and structural characteristics between different interlayer interfaces, including interface thickness and particle distribution, and obtain the ice layer interlayer interface feature identification set. Step S23: Based on the ice layer interlayer interface feature identification set, perform spatiotemporal evolution analysis of ice layer interfaces between different interlayer interfaces in the ice layer layer structure profile interface set, measure the ice layer depth annual rings between different interlayer interfaces to infer the generation, development and degradation laws of different ice layer interfaces over time, and generate spatiotemporal evolution simulation data corresponding to different interlayer interfaces inside the ice layer. Step S24: Based on the spatiotemporal evolution simulation data corresponding to different interlayer profiles within the ice layer, perform hierarchical clustering of the interfaces between different interlayer profiles in the ice layer within the ice layer layer structure profile set to generate an internal interface hierarchical classification dataset of the ice layer; based on the internal interface hierarchical classification dataset of the ice layer, perform interlayer interface feature coupling analysis on the interface features between different interlayer profiles within the ice layer interlayer interface feature recognition set to obtain internal interface feature data of the ice layer. Step S25: Detect suspected gap regions in the ice layer layer structure data based on the internal interface feature data of the ice layer to obtain suspected gap region data inside the ice layer.

[0033] As an embodiment of the present invention, reference is made to... Figure 3 As shown, Figure 1 A detailed flowchart of step S2 is shown below. In this embodiment, step S2 includes the following steps: Step S21: Identify ice layer interfaces in the ice layer layer structure data to obtain a set of ice layer layer structure interfaces; In this embodiment of the invention, ice layer profile interfaces are identified by analyzing the layered structure data of the ice layer (including the three-dimensional coordinates and echo intensity of the surface layer (0-5 meters), the middle layer (5-15 meters), and the bottom layer (15-20 meters). A three-dimensional edge detection algorithm is used to calculate the echo intensity gradient at each depth point, with a gradient threshold set to 0.2V / meter (when the change in echo intensity between adjacent depth points exceeds this value, it is determined to be an interface). At a depth of 5 meters, the echo intensity drops sharply from 0.8V at the surface layer to 0.5V at the middle layer, with a gradient of 0.06V / meter (5-meter depth difference), exceeding the threshold, and is marked as a surface-middle layer profile interface. At a depth of 15 meters, the echo intensity rises from 0.5V at the middle layer to 0.7V at the bottom layer, with a gradient of 0.02V / meter, exceeding the threshold, and is marked as a middle-bottom layer profile interface. Planar fitting was performed on each cross-section. The fitted plane equation for the surface-middle layer interface was Z = 5.0 + 0.002X + 0.001Y, and for the middle-bottom layer interface it was Z = 15.0 - 0.001X + 0.003Y, with plane deviations ≤ 0.3 meters for both. A set of cross-sections representing the ice layer's layered structure was generated, containing the three-dimensional coordinate ranges (X = 0-100 meters, Y = 0-100 meters) and plane parameters of the two cross-sections, providing a benchmark for subsequent interlayer analysis.

[0034] Step S22: Perform interlayer interface identification analysis on the different interlayer interfaces in the ice layer within the ice layer layer structure profile interface set to analyze the physical, chemical and structural characteristics between different interlayer interfaces, including interface thickness and particle distribution, and obtain the ice layer interlayer interface feature identification set. In this embodiment of the invention, interlayer interface identification and analysis were performed on two interlayer interfaces within the ice layer layer structure profile set. An ultrasonic thickness gauge was used to measure the interface thickness; the surface-middle layer interface thickness was 0.2 meters, and the middle-bottom layer interface thickness was 0.3 meters. X-ray diffraction analysis of particle distribution showed that the ice crystal particle diameter at the surface-middle layer interface ranged from 0.1 to 0.3 mm, with an average of 0.2 mm and a particle density of 120 particles / mm²; the ice crystal particle diameter at the middle-bottom layer interface ranged from 0.3 to 0.5 mm, with an average of 0.4 mm and a particle density of 80 particles / mm². Regarding physical characteristics, the hardness values ​​at the interfaces were measured: the surface-middle layer interface had a hardness of 2.5 MPa, and the middle-bottom layer interface had a hardness of 3.2 MPa. Chemically, spectral analysis revealed an oxygen content of 32% at the surface-middle layer interface and 28% at the middle-bottom layer interface. These characteristic parameters are organized into a feature recognition set for interlayer interfaces of ice layers, recording the thickness, particle diameter distribution, density, hardness and oxygen content of each interface to form quantitative feature comparison data.

[0035] Step S23: Based on the ice layer interlayer interface feature identification set, perform spatiotemporal evolution analysis of ice layer interfaces between different interlayer interfaces in the ice layer layer structure profile interface set, measure the ice layer depth annual rings between different interlayer interfaces to infer the generation, development and degradation laws of different ice layer interfaces over time, and generate spatiotemporal evolution simulation data corresponding to different interlayer interfaces inside the ice layer. In this embodiment of the invention, spatiotemporal evolution analysis of ice layer interfaces is performed based on an ice layer interlayer interface feature identification set. Through ice layer depth annual ring measurements (sampling every 0.5 meters to detect annual ring spacing), the annual ring spacing of the surface-middle layer interface gradually increases from 0.3 meters to 0.5 meters from the bottom upwards, corresponding to a formation time of 10 years (an annual increase of 0.5 meters); the annual ring spacing of the middle-bottom layer interface increases from 0.4 meters at the bottom to 0.6 meters at the top, forming over 20 years. Combined with historical temperature data (average winter temperatures of -10℃ to -15℃ over the past 30 years), the interface formation pattern is simulated: the surface-middle layer interface grows at a rate of 0.05 meters / year in low-temperature years (-15℃) and 0.03 meters / year in high-temperature years (-10℃); the middle-bottom layer interface grows at rates of 0.07 meters / year (low-temperature) and 0.04 meters / year (high-temperature). The degradation pattern is characterized by a decrease in interface hardness of 0.2 MPa every 10 years. The generated spatiotemporal evolution simulation data includes the annual growth thickness, hardness changes, and particle density change curves of the two interfaces. For example, the growth thickness of the surface-middle layer interface in the 5th year is 0.25 meters, the hardness is 2.7 MPa, and the particle density is 115 particles / square millimeter.

[0036] Step S24: Based on the spatiotemporal evolution simulation data corresponding to different interlayer profiles within the ice layer, perform hierarchical clustering of the interfaces between different interlayer profiles in the ice layer within the ice layer layer structure profile set to generate an internal interface hierarchical classification dataset of the ice layer; based on the internal interface hierarchical classification dataset of the ice layer, perform interlayer interface feature coupling analysis on the interface features between different interlayer profiles within the ice layer interlayer interface feature recognition set to obtain internal interface feature data of the ice layer. In this embodiment of the invention, hierarchical clustering of two interlayer interfaces is performed based on spatiotemporal evolution simulation data, using the Euclidean distance clustering algorithm with a distance threshold of 0.1 (after feature parameter standardization). The Euclidean distance between the feature vectors of the surface-middle layer interface (thickness 0.2, particle density 120, hardness 2.5, oxygen content 32) and the middle-bottom layer interface (0.3, 80, 3.2, 28) is 0.15 > 0.1, resulting in two clusters. Based on the clustering results, interface feature coupling analysis is performed, and the correlation coefficients of the feature parameters are calculated: the thickness of the surface-middle layer interface is negatively correlated with particle density (-0.6), and hardness is positively correlated with oxygen content (0.7); the thickness of the middle-bottom layer interface is negatively correlated with particle density (-0.5), and hardness is positively correlated with oxygen content (0.6). The coupling analysis also includes the correlation of feature change rates, such as the ratio of the annual thickness growth rate of 0.01 meters and the annual hardness reduction rate of 0.02 MPa of the surface-middle layer interface being 0.5. These coupling relationships are organized into ice layer internal interface feature data, including the feature correlation and rate of change correlation parameters of each interface.

[0037] Step S25: Detect suspected gap regions in the ice layer layer structure data based on the internal interface feature data of the ice layer to obtain suspected gap region data inside the ice layer.

[0038] In this embodiment of the invention, suspected gap areas are detected based on the characteristic data of the internal interface of the ice layer. Detection thresholds are set as follows: particle density at the surface-middle layer interface < 100 particles / mm² and hardness < 2.3 MPa; particle density at the middle-bottom layer interface < 70 particles / mm² and hardness < 3.0 MPa. Both interfaces are scanned (sampled every 1 meter × 1 meter). At the surface-middle layer interface, in the region X = 8-10 meters and Y = 8-10 meters, the particle density is 95 particles / mm² and the hardness is 2.2 MPa, meeting the threshold conditions. At the middle-bottom layer interface, in the region X = 15-17 meters and Y = 15-17 meters, the particle density is 65 particles / mm² and the hardness is 2.9 MPa, also meeting the threshold conditions. Continuity verification is performed on these areas, requiring the area of ​​continuously meeting the conditions to be ≥ 9 square meters (3 × 3 meters). Both of the above areas have an area of ​​9 square meters, thus passing the verification. Data on suspected areas of internal ice fissures were generated, and the coordinate range and characteristic parameters of the areas were recorded, such as the average particle density of 92 particles / mm² and the average hardness of 2.1 MPa in the surface-middle layer, providing a basis for subsequent precise positioning.

[0039] Furthermore, step S21 includes the following steps: The layer interface features of the ice layer structure data are analyzed to obtain the density change, temperature gradient and pressure change values ​​corresponding to each layer interface in the ice layer. In this embodiment of the invention, layer interface feature analysis is performed on the ice layer layer structure data, which includes the three-dimensional coordinates and echo intensity information of the surface layer (0-5 meters), middle layer (5-15 meters), and bottom layer (15-20 meters). A density calculation model is used, and the density change value at each layer interface is calculated based on the linear relationship between echo intensity and density (density = echo intensity × 0.2 + 0.8 g / cm³): at the surface and middle layer interface (5-meter depth), the density abruptly changes from 0.92 g / cm³ to 0.85 g / cm³, a change of -0.07 g / cm³; at the middle and bottom layer interface (15-meter depth), the density abruptly changes from 0.85 g / cm³ to 0.90 g / cm³, a change of +0.05 g / cm³. Using temperature gradient sensor data (sampled once per meter), the temperature gradient values ​​are calculated: at the 5-meter interface, the upper temperature gradient is -0.5℃ / m (decreasing with depth), and the lower gradient is -0.3℃ / m, with a change of +0.2℃ / m; at the 15-meter interface, the upper temperature gradient is -0.3℃ / m, and the lower gradient is -0.1℃ / m, with a change of +0.2℃ / m. Based on the pressure calculation formula (pressure = density × gravitational acceleration × depth), the pressure changes are obtained: at the 5-meter interface, the upper pressure is 45120 Pa, and the lower pressure is 41650 Pa, with a change of -3470 Pa; at the 15-meter interface, the upper pressure is 124950 Pa, and the lower pressure is 132300 Pa, with a change of +7350 Pa, forming the characteristic parameter set corresponding to each sub-level interface.

[0040] Preferably, the layered structure data of the ice layer is used to identify the abrupt boundary points between layers based on the density change value, temperature gradient value and pressure change value corresponding to each layer interface in the ice layer, and the boundary boundary point set between each layer interface in the ice layer is identified based on the abrupt threshold set corresponding to density, temperature gradient and pressure. In this embodiment of the invention, abrupt change thresholds are set based on the characteristic parameters of each layer interface: the density change threshold is ±0.04 g / cm³ (an absolute value ≥0.04 g / cm³ is considered a sudden change), the temperature gradient change threshold is ±0.15℃ / m, and the pressure change threshold is ±3000 Pa. The ice layer layer structure data is sampled every 0.1 meters along the depth direction (0-20 meters), and the density, temperature gradient, and pressure values ​​at each sampling point are extracted. The parameter changes between adjacent sampling points are calculated. Within a depth range of 5 meters (4.9-5.1 meters), the density is 0.92 g / cm³ at 4.9 meters and 0.85 g / cm³ at 5.0 meters, with a change of -0.07 g / cm³ (≤-0.04); the temperature gradient is -0.5℃ / m at 4.9 meters and -0.3℃ / m at 5.0 meters, with a change of +0.2℃ / m (≥0.15); the pressure is 44668 Pa at 4.9 meters and 41650 Pa at 5.0 meters, with a change of -3018 Pa (≤-3000). All three parameters meet the abrupt change threshold, and are identified as interlayer abrupt change boundary points. Within a 15-meter depth range (14.9-15.1 meters), the density is 0.85 g / cm³ at 14.9 meters and 0.90 g / cm³ at 15.0 meters, a variation of +0.05 g / cm³ (≥0.04). The temperature gradient and pressure changes also meet the thresholds, identifying this as another boundary point. Traversing the entire depth range, a set of boundary points (X and Y are planar coordinates) including coordinates such as (5.0 meters, X1, Y1) and (15.0 meters, X2, Y2) is generated. Each boundary point is verified through abrupt changes in three parameters.

[0041] Preferably, the ice layer layer structure data is calibrated and identified based on the set of boundary demarcation points corresponding to the interfaces between the layers within the ice layer, so as to obtain the ice layer layer structure interface set.

[0042] In this embodiment of the invention, ice layer profile interfaces are identified and calibrated based on the boundary boundary point set. A three-dimensional surface fitting algorithm is used to spatially interpolate the boundary points of the same interface (e.g., all (5.0, X, Y) points at a depth of 5 meters). The interpolation method is Kriging interpolation, and a spherical model (range 50 meters, sill value 0.1) is selected as the variogram model. Depth values ​​are calculated for each 1-meter × 1-meter grid point in the planar coordinate system to ensure that the depth deviation between adjacent points is ≤0.1 meters. After fitting, the interface between the surface and middle layers is a continuous surface with an average depth of 5.0 meters and a maximum undulation of 0.3 meters; the interface between the middle and bottom layers has an average depth of 15.0 meters and a maximum undulation of 0.2 meters. The interface smoothness is verified (calculation of the radius of curvature, requiring ≥10 meters), and outliers with a radius of curvature <10 meters (a total of 3, which are then re-interpolated and corrected) are removed. The final result is a set of ice layer stratification profiles, which includes three-dimensional mesh data of two complete profiles (10,000 mesh points for each profile). The coordinates of each mesh point are accurate to 0.1 meters, which can clearly reflect the spatial distribution of the profiles and provide a stratification reference for subsequent crack location. Furthermore, step S25 includes the following steps: Step S251: Divide the ice layer layer structure data into small units inside the ice layer to obtain the small structural units inside each ice layer; In this embodiment of the invention, the ice layer structure data (including the three-dimensional coordinates and echo characteristics of the surface layer 0-5 meters, the middle layer 5-15 meters, and the bottom layer 15-20 meters) is divided into small units inside the ice layer. A regular grid division method is used to divide each ice layer into 1-meter × 1-meter intervals along the horizontal direction (X, Y axis) and 0.5-meter thickness along the vertical direction (Z axis), forming cubic small units. The spatial dimensions of each unit are 1 meter × 1 meter × 0.5 meters. The surface layer is divided into 10 vertical units (5m / 0.5m = 10 layers), and a horizontal 100m × 100m area is divided into 100 × 100 units, totaling 100 × 100 × 10 = 100,000 micro-units. The middle layer, 10m thick, is divided into 20 vertical units, with the same horizontal division as the surface layer, totaling 100 × 100 × 20 = 200,000 units. The bottom layer, 5m thick, is divided into 10 layers, totaling 100 × 100 × 10 = 100,000 units. Each unit is marked with unique three-dimensional coordinates (X = 1-100m, Y = 1-100m, Z = 0.5-20m, interval 0.5m). For example, the coordinates of the first unit in the surface layer are (1,1,0.5). By correlating the coordinates with echo signal data, the micro-structural units within each ice layer are obtained, with the unit boundary error controlled within ±0.05m.

[0043] Step S252: Based on the internal interface feature data of the ice layer, perform feature quantization between the micro-structural units of each ice layer to quantify and calculate the interface curvature, interface roughness and interface reflectivity between the micro-units of each ice layer. In this embodiment of the invention, the characteristics between micro-structural units within each ice layer are quantized based on the internal interface feature data of the ice layer (including the amplitude, phase, and three-dimensional morphological data of the echo signal). Interface curvature is calculated: the interface between adjacent units is fitted with a sphere using the least squares method, with curvature = 1 / sphere radius. The sphere radius fitted for the interface between the surface units (5,5,0.5) and (5,5,1.0) is 50 meters, with a curvature of 0.02 meters⁻¹; the sphere radius fitted for the interface between the middle units (10,10,6.0) and (10,10,6.5) is 100 meters, with a curvature of 0.01 meters⁻¹. Calculate interface roughness: Take the height values ​​of 100 sampling points on the interface and calculate the root mean square roughness. The interface height deviation between the surface unit (8,8,2.0) and (8,8,2.5) is 0.03-0.08 meters, and the root mean square roughness is 0.05 meters; the interface height deviation between the middle unit (12,12,8.0) and (12,12,8.5) is 0.01-0.04 meters, and the root mean square roughness is 0.02 meters. Calculate interface reflectivity: Calculate using the ratio of echo intensity to incident intensity. The interface reflectivity between the surface unit (6,6,3.0) and (6,6,3.5) is 45% (echo intensity 4500 / incident intensity 10000); the interface reflectivity between the middle unit (15,15,10.0) and (15,15,10.5) is 3200 / 10000 = 32%. For each adjacent unit pair, the specific values ​​of interface curvature, roughness, and reflectivity are calculated to form a quantitative feature set.

[0044] Step S253: Based on the interface curvature, interface roughness and interface reflectivity between the micro units inside each ice layer, evaluate the interface jump of the micro structural units inside each ice layer to obtain the order of the interface jump between the micro units inside each ice layer. In this embodiment of the invention, interface transition is evaluated based on the interface curvature, interface roughness, and interface reflectivity between the micro-units within each ice layer. The evaluation formula is: Interface transition order = (curvature × 20 + roughness × 10 + (1 - reflectivity / 100) × 5) × 2, where the unit of curvature is meters⁻¹ and the unit of roughness is meters. For example, the interface curvature of the surface unit (5,5,0.5) and (5,5,1.0) is 0.02 meters⁻¹, the roughness is 0.05 meters, and the reflectivity is 45%. Substituting these values ​​into the formula, we get: (0.02 × 20 + 0.05 × 10 + (1 - 45 / 100) × 5) × 2 = (0.4 + 0.5 + 2.75) × 2 = 3.65 × 2 = 7.3. The interface curvature of the middle layer unit (10,10,6.0) and (10,10,6.5) is 0.01 m⁻¹, roughness is 0.02 m, and reflectivity is 32%. The calculated values ​​are: (0.01×20+0.02×10+(1-32 / 100)×5)×2=(0.2+0.2+3.4)×2=3.8×2=7.6. Performing this calculation on all adjacent unit pairs yields an interface transition order ranging from 5.2 to 12.8, with an average order of 7.5 for the surface layer and 8.1 for the middle layer. A higher order indicates a more drastic change in interface characteristics, potentially suggesting the presence of gaps.

[0045] Step S254: Based on the transition order of the internal interface between the micro-units inside each ice layer, perform crack detection on the ice layer layer structure data to obtain crack suspected region data inside the ice layer.

[0046] In this embodiment of the invention, the ice layer structure data is used to detect suspected gap regions based on the order of interface transitions between micro-units within each ice layer. An order threshold is set: ≥8.0 for the surface layer and ≥8.5 for the middle layer (determined based on the order range corresponding to gaps from historical data statistics). The order of all adjacent unit pairs is filtered. In the surface layer, consecutive unit pairs from (8,8,2.0) to (8,10,3.0) have orders of 8.2, 8.5, and 8.3, respectively, all ≥8.0, and are marked as suspected regions. In the middle layer, consecutive unit pairs from (15,15,10.0) to (15,17,11.0) have orders of 8.6, 8.8, and 8.7, all ≥8.5, and are also marked as suspected regions. The marked suspected regions are then verified for continuity, requiring ≥5 consecutive adjacent unit pairs. Two regions pass verification in the surface layer (containing 6 and 8 consecutive unit pairs, respectively), and one region passes verification in the middle layer (containing 7 consecutive unit pairs). The generated data on suspected areas of internal ice cracks includes the three-dimensional coordinate range of each area. For example, the surface area is X=8-10 meters, Y=8-10 meters, and Z=2.0-3.0 meters, while the middle area is X=15-17 meters, Y=15-17 meters, and Z=10.0-11.0 meters. The area boundaries are accurate to 0.1 meters, providing a range basis for subsequent fine detection.

[0047] Furthermore, step S253 includes the following steps: Interface distribution coupling analysis is performed based on the interface curvature, interface roughness, and interface reflectivity between the micro-units inside each ice layer to generate the ice layer interface coupling distribution field between the micro-units inside each ice layer. In this embodiment of the invention, interface distribution coupling analysis is performed based on the interface curvature, interface roughness, and interface reflectivity between the corresponding micro-units (1m × 1m × 0.5m cubes) within each ice layer. Interface curvature is obtained through three-dimensional laser scanning; the interface curvature ranges from 0.02 to 0.05 m⁻¹ between the surface and middle layers, and from 0.01 to 0.03 m⁻¹ between the middle layer units. Interface roughness is calculated using root mean square height; the surface layer roughness is 0.03-0.08 m, and the middle layer roughness is 0.01-0.04 m. Interface reflectivity is calculated based on the echo intensity of the ice-detecting radar; the surface layer reflectivity is 30%-45%, and the middle layer reflectivity is 20%-35%. The coupling degree between units is calculated using a coupling coefficient matrix (curvature weight 0.4, roughness weight 0.3, reflectivity weight 0.3). For example, the coupling degree between surface unit A (curvature 0.03m⁻¹, roughness 0.05m, reflectivity 35%) and adjacent unit B (curvature 0.04m⁻¹, roughness 0.06m, reflectivity 38%) is 0.4×|0.03-0.04|+0.3×|0.05-0.06|+0.3×|35%-38%|=0.016. This calculation is performed on all adjacent units to generate the coupling distribution field of the ice layer interface, which is presented in the form of a three-dimensional mesh. The value of each mesh point represents the coupling degree between the corresponding units, ranging from 0.005 to 0.025. Regions with high coupling degree indicate significant differences in interface characteristics.

[0048] Preferably, the interface jump barrier is calculated for the ice layer interface coupling distribution field between the micro-units inside each ice layer, so as to calculate the height of the interface jump barrier at the interface of the micro-units inside different ice layers and obtain the interface jump barrier matrix inside the ice layer. In this embodiment of the invention, the interface transition barrier is calculated by performing interface abrupt change on the coupled distribution field of the ice layer interface. The barrier height formula is: Barrier height = Coupling degree × Density difference × Gravitational acceleration × Unit thickness, where the density difference is the density difference between adjacent units (in g / cm³), and the unit thickness is 0.5 meters. For example, the density difference between surface unit A (density 0.92 g / cm³) and unit B (density 0.89 g / cm³) is 0.03 g / cm³, the coupling degree is 0.016, and the barrier height = 0.016 × 0.03 × 9.8 × 0.5 = 0.016 × 0.147 = 0.002352 J / m². The density difference between middle-layer unit C (density 0.85 g / cm³) and unit D (density 0.87 g / cm³) is -0.02 g / cm³, the coupling degree is 0.012, and the barrier height is 0.012 × 0.02 × 9.8 × 0.5 = 0.012 × 0.098 = 0.001176 J / m². The barrier height is calculated for all adjacent unit pairs, generating an interface transition barrier matrix within the ice layer. This matrix is ​​a 100 × 100 × 20 three-dimensional matrix (corresponding to 100 × 100 planar units and 20 thickness units), with element values ​​ranging from 0.0008 to 0.0032 J / m². The barrier height at each unit interface is recorded.

[0049] Preferably, the abrupt change singularity of the interface jump barrier matrix inside the ice layer is identified to obtain the corresponding abrupt change singularity of the interface jump barrier between each small unit inside the ice layer. In this embodiment of the invention, singularities are identified by abruptly changing the potential barrier matrix of the ice layer's internal interface. The threshold for determining abrupt changes is set as a difference in barrier height between adjacent matrix elements ≥ 0.001 J / m². A sliding window method (3×3×1 window) is used to traverse the matrix, calculating the difference between the center element and its eight neighboring elements. For example, the barrier height at a position (50,50,10) in the matrix is ​​0.0025 J / m², and its neighboring element values ​​are 0.0012, 0.0013, and 0.0014 J / m², with differences ≥ 0.001 J / m², thus identifying it as a singularity. After traversing the entire matrix, 12 singularities are identified in the surface layer and 8 in the middle layer. Each singularity records its three-dimensional coordinates (X,Y,Z) and abrupt barrier height values; for example, the abrupt value at (50,50,10) is 0.0013 J / m². These singularities correspond to locations where interface features change abruptly, and may be associated with fissures within the ice layer.

[0050] Preferably, the intensity and distribution density of the ice layer interface jump singularity are obtained by identifying the abrupt change of the ice layer interface jump barrier between the micro-units inside each ice layer. Based on the intensity and distribution density of the ice layer interface jump singularity, the interface jump of each micro-structural unit inside the ice layer is evaluated to obtain the ice layer interface jump order between the micro-units inside each ice layer.

[0051] In this embodiment of the invention, the singularity intensity (i.e., the abrupt change in barrier height) and distribution density (number of singularities per 100 square meters) are obtained through the abrupt change of the potential barrier at the ice interface. The average intensity of the surface singularities is 0.0015 J / m², and the distribution density is 0.3 singularities / 100m²; the average intensity of the middle layer singularities is 0.0012 J / m², and the distribution density is 0.2 singularities / 100m². The interface abrupt change is evaluated based on the intensity and density, and the evaluation formula is: abrupt change order = (intensity / maximum intensity) × 10 + (density / maximum density) × 5, where the maximum intensity is 0.002 J / m² and the maximum density is 0.4 singularities / 100m². The transition order of a surface unit is calculated as (0.0015 / 0.002)×10 + (0.3 / 0.4)×5 = 7.5 + 3.75 = 11.25; the transition order of a middle unit is calculated as (0.0012 / 0.002)×10 + (0.2 / 0.4)×5 = 6 + 2.5 = 8.5. The transition order between the interface transitions of the small units within each ice layer is generated, ranging from 5 to 15. A higher order indicates a more intense interface transition, providing a quantitative basis for subsequent gap identification.

[0052] Furthermore, step S254 includes the following steps: The probability density of ice gaps is estimated based on the order of the interface transition between the micro-units within each ice layer, so as to obtain the probability density of ice gaps between the micro-units within each ice layer. In this embodiment of the invention, the probability density of ice gaps is estimated based on the transition order of the interface between the micro-units within each ice layer. The probability density function used is: gap probability density = transition order / maximum transition order × 0.8 + 0.2, where the maximum transition order is 15. The transition order of the interface between surface unit A and unit B is 11.25, and its gap probability density = 11.25 / 15 × 0.8 + 0.2 = 0.6 × 0.8 + 0.2 = 0.68; the transition order of the interface between middle unit C and unit D is 8.5, and the gap probability density = 8.5 / 15 × 0.8 + 0.2 ≈ 0.453 × 0.8 + 0.2 ≈ 0.362 + 0.2 = 0.562. The probability density was calculated for all adjacent cell pairs. The generated probability density of fissures within the ice layer ranged from 0.53 to 0.82, with an average probability density of 0.65 for the surface layer and 0.58 for the middle layer. A higher probability density value indicates a greater likelihood of fissures at the cell interface. For example, the cell interface with a transition order of 15 has a fissure probability density of 15 / 15 × 0.8 + 0.2 = 1.0, which is directly identified as a high-probability fissure region.

[0053] Preferably, the ice layer structure data is used to detect suspected gap regions based on the probability density of gaps between the micro-units inside each ice layer, so as to divide the suspected gap regions between the micro-units inside each ice layer by adaptive threshold segmentation, and obtain the initial suspected gap region data inside the ice layer. In this embodiment of the invention, the ice layer structure data is analyzed to detect suspected gap regions based on the probability density of gaps within the ice layer corresponding to the micro-units within each ice layer. An adaptive threshold segmentation algorithm is used, and the threshold is calculated as follows: Threshold = Average probability density of the ice layer × 1.2. For the surface layer, the average probability density is 0.65, so the threshold is 0.65 × 1.2 = 0.78; for the middle layer, the average probability density is 0.58, so the threshold is 0.58 × 1.2 = 0.696. The probability density of each unit interface is compared with the corresponding threshold, and areas exceeding the threshold are marked as suspected gap regions. In the surface layer, the probability density of units (50,50,5) and (50,51,5) is 0.82 > 0.78, and they are marked as suspected regions; in the middle layer, the probability density of units (30,30,10) and (30,31,10) is 0.72 > 0.696, and they are marked as suspected regions. After detection, three suspected areas were found in the surface layer and two suspected areas were found in the middle layer. Initial crack suspected area data inside the ice layer was generated, and each area contains continuous unit interface coordinates and probability density values.

[0054] Preferably, the initial suspected crack region data inside the ice layer is subjected to morphological processing to eliminate false crack regions using opening and closing operations, and region screening is performed in combination with geometric constraints to obtain suspected crack region data inside the ice layer.

[0055] In this embodiment of the invention, morphological processing is performed on the initial suspected crack region data inside the ice layer, and opening and closing operations are used to eliminate false crack regions: the opening operation uses a 3×3×1 cubic structuring element, first eroding and then expanding to remove isolated regions with an area < 3 units (e.g., a false region with an area of ​​2 units in the surface layer is eliminated); the closing operation uses the same structuring element, first expanding and then eroding to fill the holes in the region (e.g., a region with a hole of 1 unit in the middle layer is filled). Region selection is performed based on geometric constraints, requiring a region volume ≥ 5 units (1 m × 1 m × 0.5 m) and an aspect ratio ≤ 5. The remaining two regions in the surface layer have volumes of 8 units and 6 units respectively, with aspect ratios of 3 and 2; the remaining region in the middle layer has a volume of 7 units and an aspect ratio of 2.5, all of which meet the constraints. The final data obtained were the suspected areas of the ice layer's internal fissures, including the three-dimensional coordinate range of three suspected areas, such as the surface area coordinates (49-51, 49-51, 5) and the middle area coordinates (29-31, 29-31, 10). The boundary accuracy of each area was down to 0.1 meters, providing a reliable basis for subsequent precise positioning.

[0056] Furthermore, step S3 includes the following steps: Step S31: Obtain ice echo signal data corresponding to the adjacent flight path of the fixed-wing aircraft in the shallow region of the target by using the suspected area data of the gap inside the ice layer; In this embodiment of the invention, the suspected area of ​​the ice layer's internal fissure is marked as X=500-600 meters, Y=0-50 meters, and Z=1-3 meters. Based on this coordinate range, adjacent flight track data covering this area is extracted from the fixed-wing aircraft's flight track database. The adjacent flight tracks are designated A (flight path X=0-1000 meters, Y=0 meters) and B (flight path X=0-1000 meters, Y=50 meters), with a track spacing of 50 meters, a flight altitude of 300 meters, and a flight speed of 120 km / h. Echo signals corresponding to X=500-600 meters and Y=0 meters are extracted from track A, with sampling times from 100 to 105 seconds, containing 2000 sampling points; echo signals corresponding to X=500-600 meters and Y=50 meters are extracted from track B, with sampling times from 102 to 107 seconds, containing 2000 sampling points. Both types of echo signals record amplitude (range 0-5V), phase (0-2π radians) and timestamp information to ensure that the intercepted signal completely covers the suspected area of ​​the gap. Each sampling point corresponds to a spatial resolution of 0.5m × 0.5m × 0.1m, and the acquired ice layer echo signal data of adjacent tracks can completely reflect the echo characteristics of the suspected area.

[0057] Step S32: Based on the ice echo signal data corresponding to adjacent tracks, perform three-dimensional point cloud registration on the suspected area data of the ice gap inside the ice layer to generate three-dimensional point cloud data of the ice gap inside the ice layer. In this embodiment of the invention, three-dimensional point cloud registration is generated based on ice echo signal data corresponding to adjacent tracks. First, the echo signals of track A and track B are spatiotemporally aligned (time difference 2 seconds, Y-direction offset 50 meters). Coordinate mapping is then used to make the signals of the same spatial point overlap. Time-frequency domain cross-correlation analysis (short-time Fourier transform window length 1024 points) is used to calculate signal similarity, and 1200 matching point pairs with a comprehensive similarity ≥ 0.8 are selected. After epipolar constraint verification (epoch equation Y = 0.577X + 25, distance threshold 0.5 meters) and optimization using a random sampling consensus algorithm (1000 iterations, interior point ratio 92%), 994 reliable matching point pairs are obtained. Based on the principle of triangulation, the two-dimensional signal coordinates are converted to three-dimensional coordinates (X = (X_A + X_B) / 2, Y = (Y_A + Y_B) / 2, Z = (3 × 10^2 / 25)^2 / 25 ... 8(×Δt) / 2), such as Z=1.5 meters when Δt=0.01 microseconds for the matched point pair. Select the overlapping area of ​​X=500-600 meters and Y=0-50 meters as the benchmark, calculate the transformation matrix and iteratively optimize (error of 0.08 meters after 50 iterations), and generate 3D point cloud data of the ice layer internal fissures containing 1200 points, with coordinate accuracy of 0.1 meters and point spacing of 0.5 meters.

[0058] Step S33: Denoise and smooth the three-dimensional point cloud data of the ice layer fissures by using statistical outlier removal and bilateral filtering and Laplacian smoothing to obtain the preprocessed three-dimensional point cloud data of the ice layer fissures. In this embodiment of the invention, the 3D point cloud data of the ice layer's internal fissures is denoised and smoothed. A statistical outlier removal algorithm is used, with the number of neighboring points set to 20 and the standard deviation threshold set to 2.5. The average distance from each point to its 20 neighboring points is calculated, with an overall average distance of 0.3 meters and a standard deviation of 0.05 meters. Outliers with a distance > 0.3 + 2.5 × 0.05 = 0.425 meters (a total of 32, mostly signal interference points) are removed. Bilateral filtering is used, with a spatial domain standard deviation of 0.3 meters and a value domain standard deviation of 0.1 meters. A weighted average is applied to each point (the weight of neighboring points decreases with distance and intensity differences), keeping the point cloud edges clear while reducing noise. Laplacian smoothing is then applied, iterated 5 times, adjusting the point coordinates to the average position of the neighboring points in each iteration, controlling the smoothing amplitude to ≤ 0.05 meters / iteration. After processing, the average distance deviation of the point cloud decreased from 0.12 meters to 0.06 meters. For example, the coordinates of point (550,25,1.5) were corrected to (550.02,25.01,1.49) after processing. The preprocessed three-dimensional point cloud data of the ice layer fissures was obtained, retaining 1168 points. The point cloud density was uniform and the edge features were complete.

[0059] Step S34: Perform 3D reconstruction of the ice fissure surface on the preprocessed 3D point cloud data of the ice fissure. Based on the Poisson surface reconstruction algorithm, construct the corresponding 3D surface model of the ice fissure through normal vector estimation and signed distance function. Calculate the 3D morphological geometric parameters of the ice fissure surface model to extract the length, width, depth, area, and volume parameters corresponding to the ice fissure, so as to obtain the 3D morphological geometric parameter data of the ice fissure. In this embodiment of the invention, the surface of the ice fissure is reconstructed in three dimensions by preprocessing the three-dimensional point cloud data of the ice fissure. The Poisson surface reconstruction algorithm is used to first calculate the normal vector of each point (based on a plane fitted by 8 neighboring points, with the normal vector pointing outwards from the fissure). For example, the normal vector of point (550, 25, 1.5) is (0.02, 0.01, -0.9997). A signed distance function is constructed (negative inside the point cloud and positive outside), and the reconstruction resolution is set to 512×512×512 voxels (each voxel is 0.2m×0.2m×0.2m). The marchingcubes algorithm is used to extract isosurfaces, generating a three-dimensional surface model of the ice fissure (15,000 triangular facets). The three-dimensional geometric parameters of the model were calculated as follows: length is the distance along the longest axis (X-axis), ranging from 505 meters to 595 meters, with a length of 90 meters; width is the maximum distance perpendicular to the longest axis (Y-axis), ranging from 5 meters to 45 meters, with a width of 40 meters; depth is the maximum difference along the Z-axis, ranging from 1 meter to 3 meters, with a depth of 2 meters; area is the sum of the areas of the triangular facets of the surface model, with each facet having an area of ​​0.04 square meters, for a total area of ​​15000 × 0.04 = 600 square meters; volume is the number of voxels inside the model multiplied by the voxel volume, containing 12500 voxels, with a volume of 12500 × 0.008 = 100 cubic meters. The three-dimensional geometric parameter data of the ice layer fissures were obtained, and each parameter was recorded to an accuracy of 0.1 meters or 0.1 square meters.

[0060] Step S35: Classify the crack types based on the three-dimensional morphological geometric parameter data of the ice layer cracks to obtain ice layer crack type classification data.

[0061] In this embodiment of the invention, the crack types are classified based on the three-dimensional geometric parameters of the ice layer cracks. Classification criteria are set as follows: cracks with a length > 50 meters, width > 30 meters, depth > 1 meter, and volume > 50 cubic meters are classified as large horizontal cracks; cracks with a length of 30-50 meters, width 10-30 meters, depth 0.5-1 meter, and volume 10-50 cubic meters are classified as medium-sized oblique cracks; and cracks with a length < 30 meters, width < 10 meters, depth < 0.5 meters, and volume < 10 cubic meters are classified as small vertical cracks. The crack parameters to be classified are a length of 90 meters, a width of 40 meters, a depth of 2 meters, and a volume of 100 cubic meters, all of which meet the criteria for large horizontal cracks. Further verification is achieved by calculating the width-to-depth ratio (40 / 2 = 20) and the length-to-width ratio (90 / 40 = 2.25), showing a width-to-depth ratio > 10 and a length-to-width ratio > 2, matching the horizontal extension characteristics. The final classification data of ice layer fissure types was obtained, and the fissure was marked as a large horizontal fissure. The data included the parameter values ​​and ratios of the classification criteria, providing a type identifier for subsequent engineering assessments. Furthermore, step S32 includes the following steps: Spatiotemporal alignment is performed on the ice echo signal data corresponding to adjacent tracks to calculate the time difference and spatial offset of adjacent tracks in the same detection area, and generate track spatiotemporal alignment parameters; based on the track spatiotemporal alignment parameters, coordinate mapping alignment processing is performed on the ice echo signal data corresponding to adjacent tracks to obtain spatiotemporally aligned ice echo signal data; In this embodiment of the invention, the ice echo signal data corresponding to adjacent tracks are spatiotemporally aligned. The adjacent tracks are parallel routes with a 50-meter spacing, designated as track A (X=0-1000 meters, Y=0 meters) and track B (X=0-1000 meters, Y=50 meters), both detecting the same area (X=500-600 meters, Y=0-50 meters). By comparing timestamps, the detection time difference between track A and track B in this area is calculated to be 2 seconds, with a spatial offset of 0 meters in the X direction and 50 meters in the Y direction. Based on these spatiotemporal alignment parameters, the signal data of track B is shifted -50 meters along the Y-axis so that the Y coordinates of the two tracks coincide in the same detection area. The signal data of track B is truncated with a 2-second delay on the time axis to ensure that the detection data at the same time point correspond to the same spatial location. After processing, the signals at (550,0,5) in track A and (550,50,5) in track B are aligned in time and space to form spatiotemporally aligned ice layer echo signal data. The sampling time deviation of the signal at the same spatial point is ≤0.01 seconds, and the coordinate deviation is ≤0.1 meters.

[0062] Preferably, time-frequency domain cross-correlation analysis is performed on the spatiotemporally aligned ice layer echo signal data to calculate the similarity of signal amplitude, phase and frequency at the same spatial location from different perspectives, generate an ice layer echo signal similarity matrix, and filter out corresponding matching signal point pairs by setting a similarity threshold to obtain initial signal matching point pair data; In this embodiment of the invention, time-frequency domain cross-correlation analysis is performed on the spatiotemporally aligned ice echo signal data. A short-time Fourier transform (window length 1024 sampling points, Hanning window, overlap rate 50%) is used to convert the signal to the time-frequency domain. The similarity of signal amplitude, phase, and frequency between track A and track B at the same spatial location (X=550, Y=25, Z=5) is calculated. The amplitude similarity is calculated using the normalized cross-correlation coefficient, with a value of 0.85; the phase similarity is calculated using the average of the absolute values ​​of the phase differences, with a value of 0.12 radians; and the frequency similarity is calculated using the percentage of spectral overlap, with a value of 0.88. The overall similarity is (0.85 + (1 - 0.12 / π) + 0.88) / 3 = (0.85 + 0.96 + 0.88) / 3 = 0.897, generating an ice echo signal similarity matrix containing the overall similarity of all spatial locations. A similarity threshold of 0.8 was set, and signal point pairs with a comprehensive similarity of ≥0.8 were selected, such as continuous point pairs (550,25,5) and (551,25,5), to obtain the initial signal matching point pair data, which contained a total of 1200 matching point pairs. Each point pair recorded its three-dimensional coordinates and similarity value.

[0063] Preferably, the initial signal matching point pair data is geometrically constrained and verified. An epipolar constraint model is constructed based on the radar detection angle and baseline distance of adjacent tracks. Mismatched point pairs that do not conform to the epipolar geometric relationship are eliminated. The remaining matching point pairs are iteratively optimized using a random sampling consensus algorithm. The matching point set with the highest proportion of interior points is retained to obtain the optimized signal matching point pair data. In this embodiment of the invention, geometric constraints are verified on the initial signal matching point pair data. An epipolar constraint model is constructed based on the radar detection angle (30 degrees) and baseline distance (50 meters) of adjacent tracks. The epipolar equation is Y = 0.577X + 25 (calculated based on the relative positions of track A and track B). The distance from each matching point pair to the epipolar line is calculated, with a distance threshold of 0.5 meters. Mismatched point pairs with a distance > 0.5 meters (a total of 120) are removed. The random sampling consensus algorithm is used to iteratively optimize the remaining 1080 matching point pairs. In each iteration, 8 point pairs are randomly selected to calculate the transformation matrix, and the number of interior points that satisfy the distance error < 0.3 meters is counted. After 1000 iterations, a matching point set with an interior point ratio of 92% (994 point pairs) is retained. In the optimized signal matching point pair data, all point pairs satisfy the epipolar distance < 0.5 meters and the average distance error is 0.2 meters, ensuring that the matching point pairs conform to the geometric projection relationship and providing a reliable foundation for subsequent 3D point cloud generation.

[0064] Preferably, based on the optimized signal matching point pair data, three-dimensional point cloud registration is performed on the suspected area data of the internal cracks of the ice layer. The two-dimensional echo signal coordinates are converted into three-dimensional spatial coordinates based on the signal matching point pairs, and the overlapping area corresponding to the adjacent track point clouds is selected as the reference benchmark to calculate the transformation matrix between the point clouds. At the same time, the corresponding internal crack point cloud of the ice layer is iteratively updated by minimizing the distance error from the point to the plane, so as to obtain the three-dimensional point cloud data of the internal cracks of the ice layer.

[0065] In this embodiment of the invention, three-dimensional point cloud registration is performed on the suspected area of ​​the ice layer's internal fissures based on optimized signal matching point pair data. The two-dimensional echo signal coordinates are converted into three-dimensional spatial coordinates using the triangulation principle. The calculation formulas are: X = (X_A + X_B) / 2, Y = (Y_A + Y_B) / 2, Z = (c × Δt) / 2 (where c is the speed of electromagnetic wave propagation in ice, 3 × 10⁻⁶). 8 (m / s, Δt is the echo time difference). For example, the echo time of the matched point pair (550,25,5) on track A is 10 microseconds, and the echo time on track B is 10.02 microseconds, Δt = 0.01 microseconds, Z = (3 × 10⁻⁶ m / s, Δt is the echo time difference). 8 ×0.01×10⁻ 6 ) / 2=1.5 meters. The overlapping region (X=500-600 meters, Y=0-50 meters) corresponding to adjacent track point clouds was selected as the reference benchmark. The transformation matrix between point clouds (rotation angle 0 degrees, translation (0,0,0)) was calculated. Registration was iteratively updated by minimizing the distance error from the point to the plane (target error 0.1 meters). After 50 iterations, the error was reduced to 0.08 meters. The generated 3D point cloud data of the ice layer's internal fissures contained 1200 3D points in the suspected fissure region, with a coordinate accuracy of 0.1 meters, such as (550,25,1.5), (551,25,1.6), etc., clearly presenting the spatial distribution of the fissures.

[0066] Furthermore, step S4 includes the following steps: Step S41: Perform spatial distribution analysis on the ice crevice type classification data to calculate the spatial density, directional distribution and clustering characteristics of each type of ice crevice, and obtain the spatial distribution feature vector of each type of ice crevice. In this embodiment of the invention, spatial distribution analysis is performed on ice layer fissure type classification data. This data includes one large horizontal fissure (90 meters long, 40 meters wide, etc.) and three other types of fissures (two medium-sized oblique fissures and one small vertical fissure). Spatial density is calculated: the number of fissures of each type within a 10 square kilometer area is counted. The density of large horizontal fissures is 1 / 10 = 0.1 fissures / km², the density of medium oblique fissures is 2 / 10 = 0.2 fissures / km², and the density of small vertical fissures is 0.1 fissures / km². Directional distribution is analyzed: principal component analysis is used to determine the azimuth of the longest axis of each fissure. Large horizontal fissures extend along the X-axis (azimuth 0°), medium oblique fissures have azimuths of 30° and 45°, and a directional frequency histogram is calculated. The 0° direction accounts for 25%, 30° and 45° each account for 25%, and the remaining directions account for 25%. Clustering characteristics were determined using the DBSCAN algorithm (neighborhood radius 50 meters, minimum number of points 2). Large horizontal gaps and one medium-sized oblique gap were detected as clustering (40-meter spacing), while the rest were isolated. These features were quantized into spatial distribution feature vectors: the vector for large horizontal gaps was [0.1, 0°, 1] (density, main direction, cluster identifier); the vector for medium-sized oblique gaps was [0.2, 35°, 1]; and the vector for small vertical gaps was [0.1, 90°, 0]. The vector dimension was fixed at 3, and the values ​​were accurate to 0.1.

[0067] Step S42: Based on the spatial distribution feature vectors corresponding to each type of gap, predict the gap distribution evolution trend of the ice layer structure data. Use the spatial distribution feature vectors to train the spatiotemporal autoregressive moving average model and combine it with Monte Carlo simulation. At the same time, consider the evolution of ice flow velocity and temperature change to generate corresponding ice layer gap evolution trend data. In this embodiment of the invention, the evolution trend of fracture distribution is predicted based on the spatial distribution feature vector of ice layer structure data. A spatiotemporal autoregressive moving average model (ARIMA(2,1,2)) is used, with fracture feature vector data from the past 5 years (sampled once per year) as input. During model training, spatial density and orientation angle are used as independent variables, while ice flow velocity (0.5 m / year) and temperature change (0.2℃ increase per year) are used as covariates. The model is iterated 1000 times to reduce the mean square error to 0.01. Combined with Monte Carlo simulation (1000 samplings), with ice flow velocity fluctuation range set at ±0.1 m / year and temperature change at ±0.1℃ / year, the fracture evolution parameters for the next 10 years are calculated: large horizontal fractures increase in length by an average of 5 meters per year (reaching 140 meters in the 10th year), width by an average of 2 meters per year (reaching 60 meters), and depth by an average of 0.1 meters per year due to temperature increase (reaching 3 meters). Medium-sized oblique fractures shift along the ice flow direction, moving 0.5 meters per year, with a length increase of 3 meters per year. The generated ice layer fissure evolution trend data includes the three-dimensional coordinate range and morphological parameters for each year of the next 10 years. For example, the range of large horizontal fissures in the 5th year is X=505-595 meters, Y=5-45 meters, and Z=1-3.5 meters, ensuring that the parameter changes conform to physical laws.

[0068] Step S43: Based on the ice layer gap evolution trend data, perform three-dimensional reconstruction of the ice layer gap evolution of each sub-layer interface in the ice layer layer structure data to generate a three-dimensional fine model of the gap inside the ice layer. In this embodiment of the invention, the layered structure data of the ice layer is reconstructed in three dimensions based on the evolution trend data of ice layer fissures, and the evolution parameters for the next 5 years are selected (e.g., a large horizontal fissure with a length of 115 meters and a width of 50 meters). An improved Poisson surface reconstruction algorithm is used to first calculate the normal vector (based on a plane fitted with 10 neighboring points, pointing outwards) of the fissure region point cloud in the layered structure data. For example, the normal vector of point (550, 25, 2) is (0.03, 0.02, -0.999). A signed distance function is constructed, and the voxel resolution is increased to 256×256×256 (0.4 meters / voxel). A fine three-dimensional model of the internal fissures of the ice layer is generated by isosurface extraction, containing 20,000 triangular facets with an edge error ≤0.2 meters. The model needs to fit the layer interfaces (5 meters between the surface and middle layers, and 15 meters between the middle and bottom layers), with an intersection deviation ≤0.3 meters from the interface. After reconstruction, the model clearly shows the intersection of the gaps and the layered interfaces. For example, the cross-section at a depth of 5 meters is rectangular (115 meters × 50 meters). The generated three-dimensional fine model can accurately reflect the extension state of the gaps in each layer.

[0069] Step S44: Perform 3D visualization of the fine 3D model of the fissures inside the ice layer to obtain 3D visualization positioning data of the fissures inside the ice layer; In this embodiment of the invention, a three-dimensional visualization of a detailed three-dimensional model of the ice layer's internal fissures is achieved using a volume rendering algorithm. Color mapping rules are set: the fissure surface is layered and colored according to depth (blue for 1-2 meters, green for 2-3 meters, and red for 3-3.5 meters), with transparency decreasing from 0.6 to 0.3 as depth increases. The background ice layer is displayed in grayscale (value 0.8). Three-dimensional coordinate axes (red for X-axis, green for Y-axis, and blue for Z-axis) are added, with scale intervals of 10 meters and unit labels in meters. The model is rendered from multiple perspectives: a top view shows the horizontal distribution range (X=505-620 meters, Y=5-55 meters), a side view shows the depth changes (Z=1-3.5 meters), and a cross-sectional view (along X=560 meters) shows the intersection with the layered interfaces. The rendering resolution is 1920×1080 pixels, with a frame rate of 30 frames per second, ensuring no jagged edges during dynamic rotation. The generated 3D visualization positioning data of the internal cracks in the ice layer includes multi-view images and coordinate annotations. Each feature point (such as endpoints and center points) is labeled with an accuracy of 0.1 meters, which can intuitively distinguish the spatial relationship between the cracks and the ice layer structure.

[0070] Step S45: Generate an internal ice gap identification and positioning report based on the three-dimensional visualization positioning data of the internal ice gap, and send the internal ice gap identification and positioning report to the terminal.

[0071] In this embodiment of the invention, an ice layer internal fissure identification and location report is generated based on three-dimensional visualization positioning data. The report structure includes five parts: an introduction describing the detection area (10 square kilometers) and the method; a fissure parameter table listing 10 parameters of large horizontal fissures, such as length 115 meters and width 50 meters, with a three-dimensional coordinate range (X=505-620 meters, etc.); an evolution trend diagram showing the length and width growth curves over the next 10 years (data points for each year are connected by broken lines); a three-dimensional visualization map inserting rendered images from three perspectives, with key coordinates marked; and a conclusion indicating the key monitoring areas (X=550-570 meters, Y=20-30 meters). The report is in PDF format, with the text in Arial font, size 12, and the chart resolution at 300 dpi. The report is sent to the terminal via wired transmission using TCP protocol, with a data packet size of 10MB. Checksum verification ensures data integrity. After receiving the report, the terminal automatically stores it to a designated path in the same storage format as the report, completing the entire identification and location process.

[0072] Therefore, the embodiments should be considered as exemplary and non-limiting in all respects, and the scope of the invention is defined by the appended claims rather than the foregoing description. Thus, all variations falling within the meaning and scope of the equivalents of the application are intended to be included within the invention.

[0073] The above description is merely a specific embodiment of the present invention, enabling those skilled in the art to understand or implement the invention. Various modifications to these embodiments will be readily apparent to those skilled in the art, and the general principles defined herein may be implemented in other embodiments without departing from the spirit or scope of the invention. Therefore, the present invention is not to be limited to the embodiments shown herein, but is to be accorded the widest scope consistent with the principles and novel features of the invention herein.

Claims

1. An ice layer internal crack identification and positioning method based on an airborne shallow ice exploration radar, characterized in that, Includes the following steps: Step S1: Use a fixed-wing aircraft equipped with ice-detecting radar to conduct a comprehensive flight detection of the shallow area of ​​the target to obtain multi-channel ice layer echo signal data; perform time-frequency enhancement and ice layer layer identification on the multi-channel ice layer echo signal data to obtain ice layer layer structure data; Step S2: Perform interlayer interface analysis on the ice layer layer structure data to obtain the internal interface characteristic data of the ice layer; Based on the internal interface feature data of the ice layer, the suspected gap region is detected in the layered structure data of the ice layer, and the suspected gap region data of the ice layer is obtained. Step S3: Perform three-dimensional morphological geometric analysis on the suspected area of ​​the ice layer's internal fissures to obtain the three-dimensional morphological geometric parameter data of the ice layer fissures; Based on the three-dimensional morphological geometric parameter data of ice layer fissures, fissure types are classified to obtain ice layer fissure type classification data; Step S4: Based on the classification data of ice layer gap types, perform fine reconstruction of the gap distribution evolution of ice layer layer structure data to generate a three-dimensional fine model of the gaps inside the ice layer; based on the three-dimensional fine model of the gaps inside the ice layer, perform visualization and localization expression and generate a report on the identification and localization of the gaps inside the ice layer, and send the report on the identification and localization of the gaps inside the ice layer to the terminal.

2. The method according to claim 1, wherein, Step S1 includes the following steps: Step S11: Use a fixed-wing aircraft equipped with ice-detecting radar to conduct a comprehensive flight detection of the shallow area of ​​the target to obtain multi-channel ice layer echo signal data; Step S12: Perform spatiotemporal synchronization and noise suppression processing on the multi-channel ice layer echo signal data to generate noise-suppressed ice layer echo signal data; Step S13: Perform time-frequency feature enhancement processing on the noise-suppressed ice echo signal data to obtain enhanced ice echo signal data; Step S14: Perform wavelet transform multi-scale decomposition on the enhanced ice layer echo signal data to obtain the ice layer echo coefficient data at different scales; perform multi-scale feature extraction on the ice layer echo coefficient data at different scales to extract time-domain features, frequency-domain features, and time-frequency joint features respectively, and perform feature dimensionality reduction and fusion through principal component analysis to obtain multi-scale feature data corresponding to the internal structure of the ice layer. Step S15: Based on multi-scale feature data, a deep belief network combined with a conditional random field model is used to construct an ice layer layer identification model, and the ice layer layer identification model is used to identify the ice layer layer structure data of the multi-scale feature data.

3. The method for identifying and locating internal ice gaps based on airborne shallow-layer ice-penetrating radar according to claim 1, characterized in that, Step S2 includes the following steps: Step S21: Identify ice layer interfaces in the ice layer layer structure data to obtain a set of ice layer layer structure interfaces; Step S22: Perform interlayer interface identification analysis on the different interlayer interfaces in the ice layer within the ice layer layer structure profile interface set to analyze the physical, chemical and structural characteristics between different interlayer interfaces, including interface thickness and particle distribution, and obtain the ice layer interlayer interface feature identification set. Step S23: Based on the ice layer interlayer interface feature identification set, perform spatiotemporal evolution analysis of ice layer interfaces between different interlayer interfaces in the ice layer layer structure profile interface set, measure the ice layer depth annual rings between different interlayer interfaces to infer the generation, development and degradation laws of different ice layer interfaces over time, and generate spatiotemporal evolution simulation data corresponding to different interlayer interfaces inside the ice layer. Step S24: Based on the spatiotemporal evolution simulation data corresponding to different interlayer profiles within the ice layer, perform hierarchical clustering of the interfaces between different interlayer profiles in the ice layer within the ice layer layer structure profile set to generate an internal interface hierarchical classification dataset of the ice layer; based on the internal interface hierarchical classification dataset of the ice layer, perform interlayer interface feature coupling analysis on the interface features between different interlayer profiles within the ice layer interlayer interface feature recognition set to obtain internal interface feature data of the ice layer. Step S25: Detect suspected gap regions in the ice layer layer structure data based on the internal interface feature data of the ice layer to obtain suspected gap region data inside the ice layer.

4. The method for identifying and locating internal ice gaps based on airborne shallow-layer ice-penetrating radar according to claim 3, characterized in that, Step S21 includes the following steps: The layer interface features of the ice layer structure data are analyzed to obtain the density change, temperature gradient and pressure change values ​​corresponding to each layer interface in the ice layer. Based on the density change, temperature gradient and pressure change values ​​corresponding to each layer interface in the ice layer, the abrupt boundary points between the layers are identified. The boundary point set between each layer interface in the ice layer is identified by setting the abrupt threshold corresponding to density, temperature gradient and pressure. Based on the set of boundary demarcation points between the interfaces of each layer within the ice layer, the ice layer layer structure data is used to identify and calibrate the ice layer interface to obtain the ice layer layer structure interface set.

5. The method for identifying and locating internal ice gaps based on airborne shallow-layer ice-penetrating radar according to claim 3, characterized in that, Step S25 includes the following steps: Step S251: Divide the ice layer layer structure data into small units inside the ice layer to obtain the small structural units inside each ice layer; Step S252: Based on the internal interface feature data of the ice layer, perform feature quantization between the micro-structural units of each ice layer to quantify and calculate the interface curvature, interface roughness and interface reflectivity between the micro-units of each ice layer. Step S253: Based on the interface curvature, interface roughness and interface reflectivity between the micro units inside each ice layer, evaluate the interface jump of the micro structural units inside each ice layer to obtain the order of the interface jump between the micro units inside each ice layer. Step S254: Based on the transition order of the internal interface between the micro-units inside each ice layer, perform crack detection on the ice layer layer structure data to obtain crack suspected region data inside the ice layer.

6. The method for identifying and locating internal ice gaps based on airborne shallow-layer ice-penetrating radar according to claim 5, characterized in that, Step S253 includes the following steps: Interface distribution coupling analysis is performed based on the interface curvature, interface roughness, and interface reflectivity between the micro-units inside each ice layer to generate the ice layer interface coupling distribution field between the micro-units inside each ice layer. The interface jump barrier is calculated for the coupling distribution field of the ice layer interface between the micro units inside each ice layer, so as to calculate the height of the interface jump barrier at the interface of the micro units inside different ice layers and obtain the interface jump barrier matrix inside the ice layer. By identifying the abrupt change singularity of the interface jump barrier matrix inside the ice layer, the abrupt change singularity of the ice layer interface jump barrier between each small unit inside the ice layer is obtained. By obtaining the intensity and distribution density of the ice layer interface jump singularity corresponding to the abrupt change of the potential barrier between the small units inside each ice layer, and by evaluating the interface jump of the small structural units inside each ice layer based on the intensity and distribution density of the ice layer interface jump singularity, the order of the ice layer interface jump between the small units inside each ice layer can be obtained.

7. The method for identifying and locating internal ice gaps based on airborne shallow-layer ice-penetrating radar according to claim 5, characterized in that, Step S254 includes the following steps: The probability density of ice gaps is estimated based on the order of the interface transition between the micro-units within each ice layer, so as to obtain the probability density of ice gaps between the micro-units within each ice layer. Based on the probability density of gaps in the ice layer between the micro-units inside each ice layer, the ice layer layer structure data is used to detect suspected gap regions. Then, the suspected gap regions in the ice layer between the micro-units inside each ice layer are divided by adaptive threshold segmentation to obtain the initial suspected gap region data inside the ice layer. Morphological processing was performed on the initial suspected crack region data inside the ice layer to eliminate false crack regions using opening and closing operations, and region filtering was performed in combination with geometric constraints to obtain the suspected crack region data inside the ice layer.

8. The method for identifying and locating internal ice gaps based on airborne shallow-layer ice-penetrating radar according to claim 1, characterized in that, Step S3 includes the following steps: Step S31: Obtain ice echo signal data corresponding to the adjacent flight path of the fixed-wing aircraft in the shallow region of the target by using the suspected area data of the gap inside the ice layer; Step S32: Based on the ice echo signal data corresponding to adjacent tracks, perform three-dimensional point cloud registration on the suspected area data of the ice gap inside the ice layer to generate three-dimensional point cloud data of the ice gap inside the ice layer. Step S33: Denoise and smooth the three-dimensional point cloud data of the ice layer fissures by using statistical outlier removal and bilateral filtering and Laplacian smoothing to obtain the preprocessed three-dimensional point cloud data of the ice layer fissures. Step S34: Perform 3D reconstruction of the ice fissure surface on the preprocessed 3D point cloud data of the ice fissure. Based on the Poisson surface reconstruction algorithm, construct the corresponding 3D surface model of the ice fissure through normal vector estimation and signed distance function. Calculate the 3D morphological geometric parameters of the ice fissure surface model to extract the length, width, depth, area, and volume parameters corresponding to the ice fissure, so as to obtain the 3D morphological geometric parameter data of the ice fissure. Step S35: Classify the crack types based on the three-dimensional morphological geometric parameter data of the ice layer cracks to obtain ice layer crack type classification data.

9. The method for identifying and locating internal ice gaps based on airborne shallow-layer ice-penetrating radar according to claim 8, characterized in that, Step S32 includes the following steps: Spatiotemporal alignment is performed on the ice echo signal data corresponding to adjacent tracks to calculate the time difference and spatial offset of adjacent tracks in the same detection area, and generate track spatiotemporal alignment parameters; based on the track spatiotemporal alignment parameters, coordinate mapping alignment processing is performed on the ice echo signal data corresponding to adjacent tracks to obtain spatiotemporally aligned ice echo signal data; Time-frequency domain cross-correlation analysis was performed on the spatiotemporally aligned ice layer echo signal data to calculate the similarity of signal amplitude, phase and frequency at the same spatial location from different perspectives, generate ice layer echo signal similarity matrix, and filter out corresponding matching signal point pairs by setting a similarity threshold to obtain initial signal matching point pair data; Geometric constraint verification was performed on the initial signal matching point pair data. An epipolar constraint model was constructed based on the radar detection angle and baseline distance of adjacent tracks. Mismatched point pairs that did not conform to the epipolar geometric relationship were eliminated. The random sampling consensus algorithm was used to iteratively optimize the remaining matching point pairs. The matching point set with the highest proportion of interior points was retained to obtain the optimized signal matching point pair data. Based on the optimized signal matching point pair data, three-dimensional point cloud registration is performed on the suspected area of ​​the ice layer internal fissure. The two-dimensional echo signal coordinates are converted into three-dimensional spatial coordinates based on the signal matching point pairs. The overlapping area corresponding to the adjacent track point clouds is selected as the reference benchmark to calculate the transformation matrix between the point clouds. At the same time, the corresponding ice layer internal fissure point cloud is iteratively updated by minimizing the distance error from the point to the plane, so as to obtain the three-dimensional point cloud data of the ice layer internal fissure.

10. The method for identifying and locating internal ice gaps based on airborne shallow-layer ice-penetrating radar according to claim 1, characterized in that, Step S4 includes the following steps: Step S41: Perform spatial distribution analysis on the ice crevice type classification data to calculate the spatial density, directional distribution and clustering characteristics of each type of ice crevice, and obtain the spatial distribution feature vector of each type of ice crevice. Step S42: Based on the spatial distribution feature vectors corresponding to each type of gap, predict the gap distribution evolution trend of the ice layer structure data. Use the spatial distribution feature vectors to train the spatiotemporal autoregressive moving average model and combine it with Monte Carlo simulation. At the same time, consider the evolution of ice flow velocity and temperature change to generate corresponding ice layer gap evolution trend data. Step S43: Based on the ice layer gap evolution trend data, perform three-dimensional reconstruction of the ice layer gap evolution of each sub-layer interface in the ice layer layer structure data to generate a three-dimensional fine model of the gap inside the ice layer. Step S44: Perform 3D visualization of the fine 3D model of the fissures inside the ice layer to obtain 3D visualization positioning data of the fissures inside the ice layer; Step S45: Generate an internal ice gap identification and positioning report based on the three-dimensional visualization positioning data of the internal ice gap, and send the internal ice gap identification and positioning report to the terminal.

Images