A coastal mineral resource multi-source heterogeneous data fusion evaluation method and system

Abstract

The application belongs to the technical field of computer data processing, and specifically discloses a kind of coastal mineral resources multi-source heterogeneous data fusion evaluation method and system, comprising: first, through water depth compensation and quadtree grid, the absolute space topology of land and water multi-source data is aligned;Subsequently, based on the water attenuation model, the optical frequency domain component is separated, and the acoustic backscattering texture gradient field is used to guide enhancement, to eliminate muddy water interference and generate clear seabed physical feature map;Then, the drilling grade data is extracted, combined with the seabed wave impedance to construct anisotropic conduction field, to drive the spatial diffusion of discrete data to generate continuous mineral probability field;Finally, the isosurface is extracted and the three-dimensional closed connectivity is detected, and the three-dimensional ore body boundary model without topological conflict is output.The application improves the level of coastal mineral automation prediction and the accuracy of three-dimensional modeling.

Description

Technical Field

[0001] This invention belongs to the field of computer data processing technology and relates to a method and system for the fusion and evaluation of multi-source heterogeneous data on coastal mineral resources. Background Technology

[0002] In existing technologies, coastal mineral resource assessment is a comprehensive task that heavily relies on multi-source survey data. This data encompasses land-based elevation measurements, marine multibeam bathymetry and backscattering data, geological borehole core analysis reports, and satellite or aerial remote sensing optical imagery. These data come from different sources, resulting in significant differences in their spatial reference frames, elevation datums, data formats, and semantic representations; they are typically heterogeneous data. To conduct joint analysis to delineate seabed ore bodies, it is necessary to utilize computer systems to effectively integrate, mine, and model this massive amount of multimodal data.

[0003] Currently, the common practice in the industry for processing this type of data is to manually register land elevation and sea depth data separately, and then use interpolation tools in GIS software to directly interpolate discrete borehole grade information into a two-dimensional or three-dimensional grade field using distance-inverse weighting or ordinary kriging. Simultaneously, optical remote sensing images are used for visual interpretation of the seabed sediment, or general image dehazing algorithms are used to improve underwater visibility. These operations are often independent, and the inherent physical connections between the data are not fully utilized. For example, the guiding role of seabed acoustic texture on optical image sharpness, or the constraint relationship between seabed topography and geological structure on mineralization distribution, are usually ignored in existing processes.

[0004] However, the aforementioned existing technical solutions have the following limitations: First, the lack of unified land and sea elevation benchmarks and the dramatic topographic changes near the coastline mean that direct stitching can lead to elevation jumps and topological conflicts, making manual registration inefficient and unable to guarantee accuracy over a wide area. Second, due to the scattering and absorption of suspended matter in seawater, high-frequency details in shallow sea optical remote sensing images are severely lost, and conventional image enhancement methods cannot distinguish between geological edges and water noise, resulting in enhanced images that cannot reliably reflect the continuous physical characteristics of the seabed. Third, spatial interpolation of discrete borehole data often uses isotropic mathematical methods, which do not consider the differences in physical properties of actual geological bodies in different directions, often resulting in mineralization distribution fields that exhibit smooth diffusion or singular isolated high-value areas that do not match the actual geological structure. Finally, the orebody boundaries extracted through isosurfaces often contain unclosed patches and topological errors, which can lead to incorrect volumetric calculations when directly used for geological modeling and reserve calculations, failing to meet the accuracy requirements of mineral exploration. Summary of the Invention

[0005] To overcome the aforementioned deficiencies of the existing technology and to achieve the above objectives, the present invention proposes the following technical solution: a method for evaluating multi-source heterogeneous data fusion of coastal mineral resources, comprising: S1, acquiring land-based elevation vector point data, multibeam acoustic water depth raster data and geological borehole text data, constructing a water depth compensation surface based on historical tidal extreme values, inverting the absolute elevation of the multibeam acoustic water depth raster data, and establishing a quadtree variable resolution grid for topological alignment to obtain an absolutely spatially aligned data field.

[0006] S2. Obtain shallow sea area optical remote sensing images that are spatially matched with the data field aligned with the absolute space. Transform the optical remote sensing images to the frequency domain using Fourier transform and separate the high-frequency components affected by water turbidity to obtain the optical high-frequency matrix to be compensated and the low-frequency background components.

[0007] S3. Extract the seabed backscattering texture gradient field from the multibeam acoustic grid data in the absolute spatial alignment data field. Use the seabed backscattering texture gradient field as the guiding matrix to perform directional edge enhancement on the optical high-frequency matrix to be compensated. Merge the enhanced high-frequency matrix with the low-frequency background component and perform inverse Fourier transform to obtain a clear seabed continuous physical feature map.

[0008] S4. Extract mineral composition keywords and grade values ​​from geological borehole text data. Combine the borehole opening coordinates and sampling depth, map the grade values ​​onto a quadtree variable resolution grid to obtain a three-dimensional discrete point source vector.

[0009] S5. Calculate the continuous rate of change of spatial wave impedance of a clear continuous physical feature map of the seabed, construct an anisotropic transmission coefficient field, take the three-dimensional discrete point source vector as the initial source term, and perform spatial iterative diffusion along the anisotropic transmission coefficient field to generate a continuous mineral distribution probability field.

[0010] S6. Extract isosurfaces from the probability field of continuous mineral distribution, and check the three-dimensional closed connectivity of the extracted isosurfaces in a quadtree variable resolution grid to generate a three-dimensional ore body boundary model without topological conflicts as the evaluation result.

[0011] The second aspect of the present invention provides a multi-source heterogeneous data fusion evaluation system for coastal mineral resources, comprising: a spatial alignment module for acquiring land-based elevation vector point data, multibeam acoustic water depth raster data and geological borehole text data, constructing a water depth compensation surface based on historical tidal extreme values, inverting the absolute elevation of the multibeam acoustic water depth raster data, and establishing a quadtree variable resolution grid for topological alignment to obtain an absolute spatially aligned data field.

[0012] The frequency domain decomposition module is used to acquire shallow sea optical remote sensing images that are spatially matched with the data field aligned with the absolute space. The optical remote sensing images are converted to the frequency domain through Fourier transform, and the high-frequency components affected by water turbidity are separated to obtain the optical high-frequency matrix to be compensated and the low-frequency background components.

[0013] The guided enhancement module is used to extract the seabed backscattering texture gradient field of the multibeam acoustic grid data in the absolute spatial alignment data field. The seabed backscattering texture gradient field is used as the guiding matrix to perform directional edge enhancement on the optical high-frequency matrix to be compensated. The enhanced high-frequency matrix is ​​merged with the low-frequency background component and subjected to inverse Fourier transform to obtain a clear seabed continuous physical feature map.

[0014] The point source mapping module is used to extract mineral composition keywords and grade values ​​from geological borehole text data. Combined with borehole opening coordinates and sampling depth, the grade values ​​are mapped onto a quadtree variable resolution grid to obtain a three-dimensional discrete point source vector.

[0015] The iterative diffusion module is used to calculate the continuous rate of change of spatial wave impedance in a clear continuous physical feature map of the seabed, construct an anisotropic transmission coefficient field, use the three-dimensional discrete point source vector as the initial source term, and perform spatial iterative diffusion along the anisotropic transmission coefficient field to generate a continuous mineral distribution probability field.

[0016] The orebody modeling module is used to extract isosurfaces from the probability field of continuous mineral distribution, and to check the three-dimensional closed connectivity of the extracted isosurfaces in a quadtree variable resolution grid, generating a three-dimensional orebody boundary model without topological conflicts as the evaluation result.

[0017] Compared with the prior art, the beneficial effects of the present invention are as follows: (1) The present invention realizes the benchmark unification and collaborative analysis of multi-source heterogeneous data such as land-based elevation, multi-beam acoustics, geological boreholes and optical remote sensing. It restores the details of the seabed optical image under turbid water by means of the water body suspended matter attenuation model and acoustic texture guidance. It generates a continuous mineral distribution probability field through anisotropic diffusion constrained by geological features and outputs a three-dimensional ore body boundary model without topological conflicts, which effectively improves the computer automation processing level of data integration and prediction evaluation in coastal mineral exploration.

[0018] (2) This invention constructs a water depth compensation surface by combining historical tidal extreme values ​​and uses quadtree variable resolution for topological alignment, which eliminates the inconsistency of spatial position and replaces the inefficient and difficult-to-guarantee manual registration, thus realizing high-precision fusion of land-sea heterogeneous data in vertical reference and absolute spatial position.

[0019] (3) This invention utilizes the water suspension attenuation function derived from actual measurements to accurately separate the optical frequency domain components, and employs seabed backscattering texture gradient field guidance enhancement. This method can effectively suppress water scattering noise and avoid introducing artificial artifacts while maintaining the true geological and physical edges, thus improving the clarity and reliable interpretability of underwater optical images.

[0020] (4) This invention constructs an anisotropic transmission field by utilizing clear changes in seabed wave impedance, and diffuses grade information along the direction of physical geological structure to generate a continuous mineral distribution probability field with clear physical meaning; further, by detecting the three-dimensional closed connectivity and topological correction of the isosurface of the probability field, it directly outputs a ore body boundary model without holes or free boundaries, ensuring the scientific nature of the three-dimensional ore body prediction and providing a reliable data body that can be directly used by computers for subsequent accurate resource estimation. Attached Figure Description

[0021] To more clearly illustrate the technical solutions of the embodiments of the present invention, the accompanying drawings used in the description of the embodiments will be briefly introduced below. Obviously, the drawings described below are only some embodiments of the present invention. For those skilled in the art, other drawings can be obtained based on these drawings without creative effort.

[0022] Figure 1 This is a schematic diagram of the implementation steps of the method of the present invention.

[0023] Figure 2 This is a schematic diagram of the attenuation trend of the water suspended matter attenuation function model in the frequency domain according to the present invention.

[0024] Figure 3 This is a schematic diagram of the attenuation function of the anisotropic conduction coefficient as a function of the wave impedance gradient in this invention.

[0025] Figure 4 This is a schematic diagram comparing the number of triangular facets in the three-dimensional closed connectivity detection results of the isosurface of the present invention.

[0026] Figure 5 This is a schematic diagram of the system module connections of the present invention. Detailed Implementation

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

[0028] Please see Figure 1As shown, the present invention proposes a multi-source heterogeneous data fusion evaluation method for coastal mineral resources, including: S1, acquiring land-based elevation vector point data, multibeam acoustic water depth raster data and geological borehole text data, constructing a water depth compensation surface based on historical tidal extreme values, inverting the absolute elevation of the multibeam acoustic water depth raster data, and establishing a quadtree variable resolution grid for topological alignment to obtain an absolutely spatially aligned data field.

[0029] In a preferred embodiment, the step of constructing a depth compensation surface based on historical tidal extremes, inverting the absolute elevation of the multibeam acoustic depth raster data, and establishing a quadtree variable resolution grid for topological alignment includes: acquiring historical tidal extreme value records of the area covered by the land-based elevation vector point data, generating a depth compensation surface reflecting tidal changes; subtracting the relative depth of each grid cell in the multibeam acoustic depth raster data from the compensation value at the corresponding position of the depth compensation surface to invert the absolute elevation; constructing a quadtree variable resolution grid covering the coastal transition area, and within the same grid level, performing spatial topological snapping between the inverted acoustic raster pixels and the land-based elevation vector points to eliminate elevation jumps and output an absolutely spatially aligned data field.

[0030] In a further preferred embodiment, the historical tidal extreme value record includes long-term tidal level observation data obtained from coastal tide gauge stations, and the water depth compensation surface is obtained by interpolating the highest and lowest tidal levels in the long-term tidal level observation data.

[0031] Specifically, land-based elevation vector point data, multibeam acoustic depth grid data, and geological borehole text data covering the coastal transition zone are obtained from coastal surveying departments. To ensure spatial matching of multi-source data, all data with spatial location attributes are uniformly projected onto the same plane projection coordinate system, such as the CGCS2000 coordinate system. Land-based elevation vector point data refers to the longitude, latitude, and absolute elevation values ​​of the land surface recorded in discrete point form. The absolute elevation value is referenced to the 1985 National Height Datum and obtained from field leveling or GNSS fitting. Multibeam acoustic depth grid data is generated by transmitting and receiving acoustic pulses through a shipborne multibeam echo sounder, followed by beamforming and repositioning calculations to produce regular gridded depth information. Each grid cell stores the relative depth value with the local mean sea level as the reference surface, measured in meters and with positive values ​​facing the seabed. Geological borehole text data refers to unstructured or semi-structured text reports generated from borehole logging during geological exploration. These reports contain information such as borehole opening coordinates, sampling depth ranges, lithological descriptions, and mineral grade analysis.

[0032] After obtaining the above data, a depth compensation surface is constructed based on historical tidal extremes, and then the absolute elevation of the multibeam acoustic depth raster data is retrieved. Historical tidal extreme records are obtained from coastal tide gauge stations. These stations continuously record the changes in tide level over time, forming a long-term tide level observation sequence measured using a unified absolute elevation benchmark. For each tide gauge station, the highest and lowest tide levels for all recorded periods are extracted from the long-term tide level observation data. The highest tide level is the elevation of the highest point of the tide observed at that station, and the lowest tide level is the elevation of the lowest point of the tide observed. The arithmetic mean of the sum of the highest and lowest tide levels is taken as the approximate elevation of the mean sea level at that station's location. This processing is based on the assumption of coastal tidal wave symmetry, with the center of the long-term extreme value envelope reflecting the location with the highest probability of mean sea level. (Tide gauge station) The approximate average sea level elevation at that location is denoted as Its value is obtained by summing the historical highest tide level and the historical lowest tide level given by the elevation datum after harmonic analysis and correction of the long-term tide level observation records of the station, and taking the arithmetic mean, in meters.

[0033] Discrete data obtained from all tide gauge stations The value is used to generate a continuous surface covering the range of land-based elevation vector point data and multibeam acoustic depth raster data through inverse distance-weighted spatial interpolation. This surface is defined as the depth compensation surface. At any position on the depth compensation surface... Compensation amount at the location Represented as: in, This represents the total number of tide gauge stations. It is a location To the The horizontal distance between the measuring points of the tide gauge stations is calculated using Euclidean distance in a unified plane projection coordinate system, with the unit being meters. When the interpolated position coincides with a certain tide gauge station, the distance term of that tide gauge station dominates the weight, and the surface naturally passes through the known point.

[0034] The absolute elevation of a point is obtained by subtracting the relative depth from the corresponding location in the multibeam acoustic depth raster data from the compensation amount at the same position on the depth compensation surface. The relative depth, provided by the multibeam system after sound velocity, attitude, and tidal level corrections, represents the vertical depth from the local mean sea level to the seabed, in meters. The compensation amount represents the vertical offset from the local mean sea level to a preset absolute elevation reference surface, also in meters. The resulting absolute elevation is the absolute elevation of the seabed point relative to this absolute elevation reference, with consistent dimensions and no dimensional conflicts. This calculation is performed pixel-by-pixel through the entire multibeam acoustic depth raster data to form the inverted acoustic absolute elevation raster.

[0035] Next, a quadtree variable-resolution mesh is constructed for topological alignment to integrate complex local elevation differences between land and underwater terrain. When constructing the quadtree variable-resolution mesh, the coverage area is determined by the bounding rectangle of the land-based elevation vector point data and the inverted acoustic absolute elevation raster, buffered outwards by 200 meters from the coastline. The initial mesh side length is set to 256 meters, i.e., level 0. For each mesh cell, the set of land-based vector point elevations and acoustic raster cell absolute elevations falling within that cell is statistically analyzed, and the elevation range is calculated, i.e., the difference between the maximum and minimum values. If the range exceeds a threshold of 1.5 meters and the current mesh side length is greater than the minimum mesh side length of 8 meters, the cell is divided into four equal parts to generate a more refined level. Hierarchical sub-units, with side lengths halved. The threshold of 1.5 meters is derived from the analysis of 200 sets of measured coastal topographic profiles, which can distinguish obvious intertidal micro-topographic changes. Continuous subdivision continues until the range does not exceed the limit or the side length reaches 8 meters, thus forming a quadtree variable resolution grid from coarse to fine, with the grid cell side length halved at each level as the number of levels increases.

[0036] Spatial topological snapping is performed within the same grid level to eliminate elevation jumps. For any pixel on the inverted acoustic absolute elevation grid... The cell is located within a quadtree variable resolution grid based on its center coordinates. All land-based elevation vector points are searched within this cell, and the pixel is calculated. To each vector point horizontal distance If the minimum horizontal distance is less than the adsorption threshold of 3 meters, then the pixel... Adsorption elevation By satisfying The elevation values ​​of all vector points in meters are obtained by weighting the reciprocal of the square of the distance: in, The number of vector points within the adsorption threshold; It is a vector point The absolute elevation value, in meters; The unit is meters, ensuring the weighted result is in meters. If there are no land-based elevation vector points with a horizontal distance of less than 3 meters within the cell, then the pixel... Preserve its original absolute elevation The adsorption threshold of 3 meters is set based on the density of conventional feature survey points from the tide gauge station to the shoreline, and can distinguish between real topographical abrupt changes and false jumps caused by inconsistencies in the measurement benchmark.

[0037] After pixel-by-pixel adsorption, the acoustic pixel elevations in shallow water areas automatically align with adjacent land control points, replacing the previously rigid land-sea elevation transition steps with a continuous transition. The adsorbed elevation values ​​are then written back to the raster, replacing the original values, and stored within the same spatial reference frame along with the land-based elevation vector points and borehole coordinates from the geological borehole text data, resulting in an absolutely spatially aligned data field. Within this data field, elevation information at any location remains consistent in both vertical reference and spatial position, without any topological conflicts.

[0038] For example, a total of 1562 land-based elevation vector point data points were acquired for a certain bedrock coastal area. The spacing between the measurement points was approximately 20 meters in flat areas and approximately 5 meters in the intertidal reef area. Multibeam acoustic depth raster data had a raster resolution of 4 meters, 1500 rows and 1800 columns, with depth values ​​ranging from 0.2 meters at the shore to 23.5 meters in the underwater shallows. Geological borehole text data included data from 12 boreholes. Three tide gauge stations, numbered T1, T2, and T3, were located along the coast of this area. Station T1 recorded a maximum tide level of 3.82 meters and a minimum tide level of -2.46 meters; station T2 recorded a maximum tide level of 4.15 meters and a minimum tide level of -2.71 meters; and station T3 recorded a maximum tide level of 3.94 meters and a minimum tide level of -2.59 meters. The arithmetic mean sea level of each station was approximately 0.68 meters, 0.72 meters, and 0.675 meters, respectively. Pixel Q in the existing multibeam grid is located at 121.8123°E, 30.7231°N, with a relative water depth of 8.4 meters. The horizontal distances from Q to three tide gauge stations are measured to be 1250 meters, 980 meters, and 2150 meters, respectively. Using the inverse distance weighted interpolation algorithm, combined with the horizontal distances from Q to the three tide gauge stations, the compensation at this location is calculated to be approximately 0.701 meters. Subtracting this from the relative water depth, the absolute elevation of pixel Q is approximately -7.699 meters.

[0039] When constructing the quadtree variable resolution mesh, the initial cell side length was 256 meters, covering the area. Within a certain mesh cell containing cell Q, the elevation range was statistically determined to be 5.2 meters, exceeding the threshold of 1.5 meters. Further subdivision was performed to a sub-cell with a side length of 64 meters. Within this sub-cell, there were two land-based elevation vector points: V1 with an absolute elevation of 1.34 meters and V2 with an absolute elevation of 1.45 meters. Before absorption, the elevation of the acoustically inverted cell Q was -7.699 meters, forming a steep jump of approximately 9 meters with the vector points. The horizontal distance from Q to V1 was measured to be 2.4 meters, and the horizontal distance to V2 was 2.9 meters, both less than the absorption threshold of 3 meters. Performing an inverse distance squared weighted calculation, the absorbed elevation of cell Q was found to be approximately 1.384 meters. This value falls precisely between the elevations of adjacent land-based points, removing the original jump and resulting in a smooth and natural terrain appearance at the land-water interface. After this adsorption process, the acoustic grid and the land-based vector points are fused into a continuous elevation surface. Combined with the borehole coordinate information in the geological borehole text data, this forms an absolutely spatially aligned data field, providing a unified axis and vertical reference for subsequent shallow-sea optical remote sensing image fusion and ore body modeling.

[0040] S2. Obtain shallow sea area optical remote sensing images that are spatially matched with the data field aligned with the absolute space. Transform the optical remote sensing images to the frequency domain using Fourier transform and separate the high-frequency components affected by water turbidity to obtain the optical high-frequency matrix to be compensated and the low-frequency background components.

[0041] In a preferred embodiment, separating the high-frequency components affected by water turbidity to obtain the optical high-frequency matrix to be compensated and the low-frequency background components includes: constructing a water suspended matter attenuation function model corresponding to the optical remote sensing image of the shallow sea area, which expresses the attenuation characteristics of the optical signal in the spatial frequency domain under different water turbidity; comparing the frequency domain image after Fourier transform with the water suspended matter attenuation function model to separate the optical high-frequency components attenuated due to water scattering and absorption as the optical high-frequency matrix to be compensated, while retaining the low-frequency background components that have not been attenuated.

[0042] In a further preferred embodiment, the water body suspended matter decay function model is an exponential decay function, the independent variables of which include a spatial frequency component and a spatial cutoff frequency; wherein, the spatial cutoff frequency is inversely proportional to the product of the water depth and turbidity at the corresponding location.

[0043] Specifically, acquire shallow-sea optical remote sensing images that spatially match the absolutely spatially aligned data field. These shallow-sea optical remote sensing images refer to digital images obtained by optical sensors mounted on airborne or spacecraft platforms, imaging underwater coastal areas in the visible to near-infrared bands. Image pixels store quantized values ​​of the radiance reflected by ground features. If the original data are dimensionless digital quantized values, they must be pre-calibrated to radiance values ​​to ensure dimensional uniformity in the subsequent physical model. Their spatial extent completely overlaps with the absolutely spatially aligned data field constructed in step S1 in both the projected coordinate system and grid position. Spatial matching is achieved as follows: Using the highest resolution level of the quadtree variable-resolution grid in the absolutely spatially aligned data field as a reference, coastline, island / reef edges, and obvious underwater topographic feature points are extracted as corresponding control points. Polynomial geometric fine correction is performed on the optical remote sensing image, with at least 12 pairs of correction control points evenly distributed around the image's perimeter and center. Bilinear interpolation is used for resampling to ensure that the coordinates and size of the resampled optical pixels are consistent with the acoustic raster pixels in the absolutely spatially aligned data field. For multispectral images with 8-meter or 4-meter resolution, the pixel alignment error is controlled within 0.5 pixels. The spatially matched shallow sea area optical remote sensing image is denoted as... ,in and This represents the row and column number of the pixel in the image coordinate system. The radiance at the corresponding location is expressed in watts per square meter per steradian micrometer, but it is treated as a pure numerical sequence in subsequent processing, and its dimensions remain unchanged after Fourier transform.

[0044] The spatially matched optical remote sensing image of the shallow sea area was transformed to the frequency domain using a standard two-dimensional discrete Fourier transform, resulting in a frequency domain matrix containing complex elements. After shifting the zero-frequency component to the center of the matrix, a spectral image symmetrically distributed around the zero frequency is obtained. Based on the total number of rows M, the total number of columns N, and the pixel ground side length Δs, the dimensionless frequency domain index (u,v) is converted into a physical space frequency in units of "per period per meter" to ensure the dimensional consistency of the subsequent physical model.

[0045] Following frequency domain transformation, a water suspension attenuation function model corresponding to the optical remote sensing images of the shallow sea area was constructed. This model expresses the attenuation characteristics of optical signals in the spatial frequency domain under different water turbidity levels. It is a filter function that decreases with increasing spatial frequency, used to characterize the attenuation of high-frequency details in underwater imaging caused by scattering and absorption effects from suspended sediment, plankton, etc. To construct this model, the optical attenuation coefficients of water bodies with different turbidity were measured at 25 stations along the coast using an on-site optical profiling system. Turbidity was expressed in scattering turbidity units (NTU). Simultaneously, modulation transfer function (MTF) curves were obtained at different underwater depths using a standard whiteboard and a resolution target. Each MTF curve was fitted to an exponential decay form to determine the spatial cutoff frequency when the MTF decreases to 1 / e. Statistical analysis revealed that in turbid water bodies dominated by sediment, the spatial cutoff frequency... With water depth (unit: meter) and turbidity There is an inverse proportional relationship between the units (NTU). An empirical formula was derived by fitting 200 sets of measured data: in, These are empirical constants determined by regression analysis and calibrated based on these 200 sets of measured data. The dimension is NTU. Therefore, for any pixel in the spatially matched optical image corresponding to the ground location, the absolute elevation is extracted from the absolute spatial alignment data field through interpolation and converted into the actual water depth at that location by combining it with the tide data acquired during the optical image acquisition. Simultaneously, turbidity data from discrete water quality sampling points acquired synchronously with optical remote sensing images are obtained. A continuous turbidity surface is generated through spatial interpolation to extract the turbidity at the corresponding locations. Then the spatial cutoff frequency at that location can be calculated. Then, establish an attenuation function that varies with frequency: in and This represents the spatial frequency component. The function's range is between 0 and 1, and it is dimensionless. When the radial spatial frequency... When the frequency approaches 0, it approaches 1; when the frequency is much greater than the cutoff frequency, it approaches 0, reflecting the low-pass filtering effect of the water body.

[0046] Please see Figure 2 As shown, Figure 2 This diagram illustrates the attenuation trend of the suspended solids attenuation function model in the water body according to the present invention in the frequency domain. As shown in the figure, the model uses the normalized radial spatial frequency as the independent variable. With increasing frequency, the normalized attenuation function value decreases exponentially. When the spatial frequency reaches or exceeds the cutoff frequency... At that time, the high-frequency details affected by scattering and absorption by turbid water are strongly attenuated. The attenuation curve intuitively reflects the low-pass filtering effect of the water body, providing physical model support for the subsequent accurate separation of high and low frequency components.

[0047] The frequency domain image after Fourier transform is compared with the water suspended matter attenuation function model to separate the high-frequency components affected by water turbidity, thus obtaining the optical high-frequency matrix to be compensated and the low-frequency background components. The comparison operation is implemented by frequency-wise multiplication: for frequency coordinates... frequency domain values Multiply by the decay function at the corresponding position of that pixel As a low-frequency background component Because the attenuation function of the low-frequency part is close to 1, it is hardly weakened; at the same time, Multiply by the complement The optical high-frequency matrix to be compensated is obtained. This section includes high-frequency details and noise that are attenuated by scattering and absorption by water. The specific formula is as follows: The above operation is performed point-by-point on the spectrum of each pixel in the image. Since the water depth and turbidity may be different in each pixel, The spatial domain changes are gradual, so a block-based strategy is adopted in actual processing. The image is divided into small blocks of 32×32 pixels. Within each block, a uniform attenuation function is calculated using the water depth and turbidity of the central pixel to balance computational efficiency and transition continuity. In this way, the original frequency domain image is decomposed into two complex matrices of the same size: It retains low-frequency information such as large areas of continuous seabed background and uneven lighting. This includes high-frequency components such as seabed texture and target boundaries that are suppressed by turbid water. The elements in both matrices retain the dimensions of the original radiance frequency domain values. Combining them yields the complete matrix. .

[0048] For example, continuing with the bedrock coastal zone area in step S1 of the embodiment, underwater digital topography can be obtained from the absolutely spatially aligned data field. A spatially matched optical remote sensing image of the shallow sea area is acquired in the same region. The image is a single-band blue-green band image with a pixel ground side length of 4 meters and a row and column number of... , After radiometric correction, the radiance value of a certain pixel A in the image. The unit is watts per square meter per steradian micrometer. Pixel A corresponds to a position in the absolutely spatial aligned data field, at which the water depth is calculated from the aligned data field. The turbidity was measured in meters through simultaneous on-site sampling and interpolation. NTU. Substituting the above parameters into the empirical formula for spatial cutoff frequency, the spatial cutoff frequency at this location is calculated to be approximately 0.02345 cycles / meter. This cutoff frequency means that approximately 0.023 cycles can be resolved per meter, corresponding to a wavelength of approximately 42.6 meters, indicating that water turbidity severely attenuates details smaller than tens of meters. Now, a Fourier transform is performed on the small block containing pixel A, with a block size of 32×32 pixels and a spatial range of 128 meters × 128 meters. For a certain frequency location within the block, its physical spatial frequency is... The period / meter, when substituted, yields an attenuation function value of approximately 0.033 at that frequency, indicating that only 3.3% of this high-frequency component is retained after passing through the water, meaning it is almost completely lost. Correspondingly, the complement... If the complex value in the frequency domain at that frequency position, after being subjected to a discrete Fourier transform, is... (Modulus approximately 5.07) If the unit is imaginary, then the value of that frequency retained in the low-frequency background component can be calculated using the high-low frequency separation formula. The corresponding value stored in the optical high-frequency matrix to be compensated is... This high-frequency component is the part that needs to be enhanced using multi-beam acoustic backscattering texture in subsequent steps. This separation is performed sequentially on all frequencies within a block, and the results from all blocks are combined to obtain the full-image optical high-frequency matrix to be compensated. and low-frequency background components These two matrices maintain spatial matching with the absolute spatial alignment data field output from step S1, providing accurate high-frequency missing components and a stable background for texture-guided fusion in step S3.

[0049] S3. Extract the seabed backscattering texture gradient field from the multibeam acoustic grid data in the absolute spatial alignment data field. Use the seabed backscattering texture gradient field as the guiding matrix to perform directional edge enhancement on the optical high-frequency matrix to be compensated. Merge the enhanced high-frequency matrix with the low-frequency background component and perform inverse Fourier transform to obtain a clear seabed continuous physical feature map.

[0050] In a preferred embodiment, the step of extracting the seabed backscattering texture gradient field from the multibeam acoustic raster data in the absolute spatial alignment data field, and using the seabed backscattering texture gradient field as a guiding matrix to perform directional edge enhancement on the optical high-frequency matrix to be compensated, includes: obtaining backscattering intensity information from the multibeam acoustic raster data in the absolute spatial alignment data field; calculating the spatial rate of change of the backscattering intensity information to generate a gradient field characterizing the direction and intensity of the seabed micro-topographic texture, as the seabed backscattering texture gradient field; using the seabed backscattering texture gradient field as the guiding matrix for guided filtering, retaining edge components in the optical high-frequency matrix to be compensated that are consistent with the gradient field direction, suppressing water scattering noise that is inconsistent with the gradient field direction, and obtaining the enhanced high-frequency matrix.

[0051] Specifically, the seabed backscattering texture gradient field of the multibeam acoustic grid data is extracted from the absolute spatial alignment data field output in step S1, and this is used as a guiding matrix to perform directional edge enhancement on the optical high-frequency matrix to be compensated obtained in step S2. The enhanced high-frequency matrix is ​​merged with the low-frequency background component, and then a clear seabed continuous physical feature map is generated through inverse Fourier transform.

[0052] The absolute spatially aligned data field contains multibeam acoustic grid data. When this grid data is acquired by the multibeam bathymetry system in step S1, each pixel synchronously stores two inherent attributes: absolute elevation value and backscatter intensity information. The backscatter intensity information is obtained by the multibeam bathymetry system recording the acoustic energy of the seabed reflected echoes simultaneously with the transmission of acoustic pulses. After acoustic propagation loss compensation, seabed incident angle correction, and Lambert's law normalization, it is converted into a backscatter intensity related to the physical characteristics of the seabed, expressed in decibels. The backscatter intensity grid of the multibeam acoustic grid data has the same grid resolution and spatial reference as the inverted absolute elevation grid, with a one-to-one pixel correspondence. The backscatter intensity grid is extracted separately from the absolute spatially aligned data field and denoted as... ,in and These are row and column coordinates. This represents the backscattering intensity, measured in dB. To perform pixel-level matching with the optical image, spatial relationships are utilized... Perform bilinear interpolation resampling to ensure that its spatial resolution and pixel alignment are completely consistent with the shallow sea area optical remote sensing image after spatial matching in step S2, i.e., the pixel ground edge length is the same. meters, image size is OK List.

[0053] Based on this, the spatial rate of change of backscattering intensity information is calculated, generating a gradient field characterizing the direction and intensity of seabed micro-topographic texture. For the resampled backscattering intensity image... The convolution operation is performed using the commonly used 3-pixel × 3-pixel Sobel horizontal and vertical operators. At the pixel position... At each location, the horizontal gradient components were calculated. and vertical gradient components Specifically, the gradient magnitude is calculated as the sum of the products of each pixel value in the neighborhood and its corresponding operator element, with the result in dB per meter. Based on these two components, the square root of the sum of the squares of the horizontal and vertical gradient components is taken to calculate the gradient magnitude matrix pixel by pixel. Simultaneously, the arctangent of the ratio of the vertical gradient component to the horizontal gradient component is calculated, and the gradient direction matrix is ​​calculated pixel by pixel. The arctangent function returns the value in radians, and the gradient direction ranges from -π to π. and Together they form the gradient field of the seabed backscattering texture. This reflects the drastic changes in acoustic response caused by seabed boundaries and micro-topographic undulations. This indicates the main direction of the change. For example, in a sand wave development area, the backscattering intensity changes gently along the sand wave ridge line but changes drastically perpendicular to the ridge line. The corresponding gradient direction will be perpendicular to the ridge line and point to the direction of the steepest undulating gradient.

[0054] Next, extract the optical high-frequency matrix to be compensated from the output of step S2. and low-frequency background components Both matrices are complex matrices in the frequency domain, with dimensions consistent with the optical remote sensing image of the shallow sea area. To utilize the seabed backscattering texture gradient field for directional edge enhancement of high-frequency components, the optical high-frequency matrix to be compensated in the frequency domain is first transformed to the spatial domain using a two-dimensional discrete inverse Fourier transform, yielding a spatial domain high-frequency image. This spatial domain high-frequency image is a real-valued matrix, characterizing the details, edges, and noise in the optical image attenuated due to water scattering and absorption. Its values ​​represent radiance residuals, with units consistent with the original shallow-sea optical remote sensing image: watts per square meter per steradian micrometer. Similarly, its spatial gradient information is also calculated to obtain the gradient magnitude of the optical high-frequency image. and direction The calculation formula is exactly the same as that for backscattering gradient calculation, only the input is replaced with... This yields the gradient components in units of radiance per meter.

[0055] We now use the gradient field of the seabed backscattering texture as the guiding matrix for guided filtering of high-frequency images in the spatial domain. Perform directional enhancement guided filtering. The core idea of ​​guided filtering is to enhance the directionality of filtering within a local window. Inside, assuming the target high-frequency image It is a guide image The linear transformation, i.e., the locally linear model, is expressed as: , where p represents a local window Spatial location index of each pixel within the cell; This represents the feature magnitude of the guiding image G at pixel p. Here, the backscattered gradient magnitude is used. As a guide image Furthermore, adaptive weights defined by directional consistency are introduced into the filtering cost function. This ensures that the filtering process preserves optical high-frequency edges aligned with the backscattering texture direction, while smoothing out water scattering noise with inconsistent directions. For the pixel to be enhanced... A square window centered The window side length is set to 9 pixels, and the local linear coefficients are defined according to the standard guided filtering. and This is obtained by solving the following weighted least squares problem: In the formula, Display window Each pixel within; It is a high-frequency image in the spatial domain at the pixel level The value; It is the backscattered texture gradient field at the pixel The amplitude; and These represent local windows. The linear scaling factor and translation intercept within the range are obtained by solving the above weighted least squares cost function analytically. To prevent Excessively large regularization parameters should be set within the window. 0.01 times the variance. Weight The directional consistency factor is defined as the square of the cosine of the angle between the optical high-frequency gradient direction and the backscattering texture gradient direction. in This represents the optical high-frequency gradient direction angle at pixel p; This represents the angle of the backscattered texture gradient direction at pixel p. When the two gradient directions are aligned, the angle difference is close to 0 or... When the cosine square is close to 1, the fitting error of that pixel is fully accounted for, forcing the local linear model to retain the high-frequency changes at that point; when the directions are orthogonal, the weight is close to 0, the pixel error is ignored, and the output value tends to the window mean, i.e., smoothing the noise. The smoothing intensity is controlled in conjunction with the window size. The dimensions are consistent with the variance of the guiding image, ensuring homogeneity of dimensions. Solving the above problem yields the center of each window. and Then, for all covered pixels The average coefficient of the window is used to calculate the enhanced pixel value. : in, , , For all containing cells Local window center The set, Represents a set The number of elements in the image. The enhanced image yields a high-frequency image in the spatial domain. It retains the edge features in the original optical high frequency that are consistent with the acoustic texture gradient structure, and effectively suppresses random noise caused by water scattering.

[0056] Enhanced spatial domain high-frequency image The enhanced high-frequency matrix is ​​obtained by converting it back to the frequency domain using a standard two-dimensional discrete Fourier transform. .

[0057] The enhanced high-frequency matrix The low-frequency background components retained in step S2 The corresponding positions in the frequency domain are added together to complete the merging, resulting in the merged frequency domain matrix. The merged frequency domain matrix It fully includes texture-guided enhanced high-frequency information and low-frequency background information unaffected by turbid water. Finally, for Performing a standard two-dimensional discrete inverse Fourier transform yields a clear, continuous physical feature map of the seabed in the spatial domain. Its pixel values ​​have been restored to radiance measurements, maintaining the same physical units as the original shallow sea optical remote sensing image. The image clearly shows the continuous landforms and seabed boundaries, and has eliminated the blurring and loss of details caused by water turbidity.

[0058] For example, continuing with the bedrock coastal zone example in steps S1 and S2, the absolutely spatially aligned data field has been constructed. From this data field, a backscatter intensity grid aligned with the optical image is extracted, where the backscatter intensity value corresponds to pixel A in the S2 example. dB. Using the Sobel operator, the horizontal gradient component of pixel A and its neighborhood is approximately 1.6 dB / m, and the vertical gradient component is approximately -0.9 dB / m. The gradient magnitude can then be calculated. It is 1.836 dB / m, gradient direction The value is -0.512 radians. Simultaneously, the optical high-frequency matrix to be compensated obtained from step S2... After inverse Fourier transform, the high-frequency value at pixel A in the spatial domain (Radiance unit), its own gradient amplitude ,direction Radius. Visible. and The difference is about 0.025 radians, which is almost identical. This indicates that the optical high-frequency edge direction at this location is highly consistent with the acoustic texture gradient direction and should be preserved.

[0059] In a 9x9 window containing pixel A Inside, assuming there exists another pixel B with an optical high-frequency gradient direction of 1.2 radians, while the backscattering texture gradient direction is still close to -0.5 radians, the directional difference is 1.7 radians. Substituting these values ​​into the directional consistency factor formula, the weights are calculated... Approximately 0.017, it hardly contributes to determining the local linearity coefficients. The value is obtained through weighted least squares calculation. and Enhanced output at pixel A It is close to the original value of 1.52, while at pixel B... Then, regressing to the window mean results in significant smoothing. After processing the entire image, [the following steps are taken]. Perform a Fourier transform, then combine with... By merging and inversely transforming, a clear and continuous physical feature map of the seabed can be obtained. .exist In the image, the bedrock trenches and gravel boundaries, which were previously blurred in the optical image of turbid water, become clearly discernible. The boundary orientation is consistent with the hard substrate edge revealed by the multibeam backscattering texture, verifying the effectiveness of directional edge enhancement. This image will serve as the basis for calculating the continuous rate of change of spatial wave impedance in step S5.

[0060] S4. Extract mineral composition keywords and grade values ​​from geological borehole text data. Combine the borehole opening coordinates and sampling depth, map the grade values ​​onto a quadtree variable resolution grid to obtain a three-dimensional discrete point source vector.

[0061] In a preferred embodiment, the step of extracting mineral composition keywords and grade values ​​from geological borehole text data, combining borehole opening coordinates and sampling depth, and mapping the grade values ​​onto a quadtree variable resolution grid to obtain a three-dimensional discrete point source vector includes: using natural language processing named entity recognition technology to extract the composition keywords corresponding to the target mineral and the grade values ​​associated with the composition keywords from the geological borehole text data; obtaining the borehole opening coordinates and sampling depth, using the borehole opening coordinates and sampling depth as spatial locations, setting the extracted grade values ​​as attributes at that spatial location, forming a three-dimensional discrete point source vector, wherein the spatial reference of the three-dimensional discrete point source vector is consistent with the absolute space alignment data field.

[0062] Specifically, from the geological borehole text data obtained in step S1, natural language processing (NLP) named entity recognition technology is used to extract the component keywords corresponding to the target minerals and the grade values ​​associated with the component keywords. Geological borehole text data consists of unstructured or semi-structured text reports generated from borehole logging during geological exploration. These reports describe stratigraphic lithology, tectonic phenomena, and mineral analysis results in natural language paragraphs. Mineral analysis results are often presented in the form of "grade value of a certain mineral %" or "average grade value of a certain element g / t". NLP named entity recognition technology refers to an information extraction method that uses a pre-trained sequence labeling model to locate entity boundaries and classify entity types from a text sequence. In this step, entity types are divided into two categories: "mineral component keywords" and "grade values". Mineral component keywords are the standard or petrological names of the target mineral in the text, such as "magnetite," "chalcopyrite," and "gold." Grade values ​​are the grade values ​​obtained from borehole analysis that are immediately adjacent to or semantically related to the mineral in the text.

[0063] The named entity recognition model employs an architecture combining a bidirectional long short-term memory network (LSTM) and a conditional random field (CRF). This architecture excels at capturing long-distance lexical dependencies commonly found in geological reports and effectively handles complex sentence structures. The input layer consists of Chinese character-level embedding vectors. After encoding contextual features through the LSTM layer, the CRF layer outputs the globally optimal entity label sequence. The training dataset comprises 300 full-text coastal mineral geological borehole reports extracted from the geological archives of coastal provinces. These reports were manually annotated using the character-level "BIO" labeling method. The specific annotation rules are: the first character of the mineral name is labeled "B-MIN", the remaining characters are labeled "I-MIN", the first character of the grade number is labeled "B-GRADE", subsequent numbers or decimal points are labeled "I-GRADE", and other characters are labeled "O". The trained model achieves an F1 score of 0.937 for mineral name entity recognition and 0.961 for grade value recognition, meeting the requirements for production applications.

[0064] For each geological borehole text data set, the model inputs a complete text paragraph and outputs a series of "mineral composition keywords - grade values" pairs. Each pair also includes the sampling depth range corresponding to the paragraph description. When the text contains records like "magnetite grade 42.7%, depth 12.30m to 13.80m," the model identifies the mineral composition keyword "magnetite" and the grade value "42.7," and infers the sampling depth range to be 12.30m to 13.80m based on preset depth number recognition rules that include keywords such as "depth," "from," and "to," along with numerical analysis. The midpoint of this depth range is taken as the sampling depth for that grade sample, in meters. If the same depth range describes multiple minerals, they are extracted separately, each forming an independent sample. The grade value is recorded in the raw dimensions given in the laboratory report, commonly as a mass percentage or grams per ton.

[0065] To unify the dimensions of attributes in subsequent diffusion calculations, if the grade value The original unit is grams per ton, then it can be converted to a dimensionless mass percentage using the following formula. : In the formula, This refers to the grade value in grams per ton, with 10,000 as the conversion factor. This is based on the fact that 1% equals 10,000 grams per ton. After conversion... This is a pure numerical percentage. If the original grade value is already a percentage, retain the original value. This process ensures that all grade values ​​are uniformly expressed as dimensionless mass percentages.

[0066] Next, the borehole coordinates and sampling depth are obtained. During the design and construction of each borehole, the geodetic coordinates of the borehole opening are measured using GNSS static positioning. These coordinates are then converted to Cartesian coordinates using the same projection transformation formula as the data field aligned with absolute space. The unit is meters. The orifice coordinates and the sampling depth obtained in the previous stage are then compared. The vertical component, representing the three-dimensional spatial location, is used to extract and convert the grade value. As an attribute at this spatial location, grid mapping is performed. The mapping operation uses the quadtree variable resolution grid constructed in step S1 as the planar index structure. This grid shares the same spatial reference and grid system with all raster data processed in steps S1, S2, and S3. For each grade sample, based on the orifice plane coordinates... Starting from the finest level of the quadtree variable resolution mesh, the leaf node mesh cells containing the given point are searched up level by level until a unique planar adsorption cell is determined. The coordinates of the center point of this cell are marked as... This adsorption operation normalizes the precise measurement location of the borehole to the center of the grid cell, ensuring the alignment of raster cells in the absolute spatial alignment data field.

[0067] Each grade sample forms a three-dimensional spatial point after adsorption. Its spatial location is defined as The attribute value carried by this point is the grade value. All samples from all boreholes are processed according to this rule to obtain a series of three-dimensional spatial sampling points with grade attributes. When multiple sampling points fall into the same quadtree variable resolution planar grid cell, and the difference in their sampling depths is less than the vertical adsorption threshold of 0.5 meters, these points are merged into one point in the vertical direction. The sampling depth of the merged point is the arithmetic mean of the depths of each point, and the grade value is also the arithmetic mean of the grade values ​​of each point. The vertical adsorption threshold of 0.5 meters is set based on the 1-meter sampling interval commonly used in coastal placer and rock mineral geological exploration, and takes into account the error of laboratory sampling, which can distinguish different mineralization layers. All the retained points after merging constitute a three-dimensional discrete point source vector, which consists of a series of quadruples recording spatial location and attributes. The structure records the discrete distribution of grade data in three-dimensional space. The spatial reference plane coordinate system and elevation datum of the three-dimensional discrete point source vectors are completely consistent with the absolute spatial alignment data field and can be directly used as the internal source term of the Poisson diffusion equation in step S5.

[0068] For example, based on the aforementioned bedrock coastal zone example, step S1 obtains text data from 12 geological boreholes. The logging text for borehole ZK05 includes: "Second sampling segment, from 8.20m to 9.80m depth, lithology copper-bearing magnetite quartzite, magnetite grade 36.5%, chalcopyrite grade 1.22%". After processing this text, the named entity recognition model outputs the paired terms "magnetite-36.5" and "chalcopyrite-1.22", with the sampling depth range being 8.20m to 9.80m, and the midpoint depth... Meters. Both grade values ​​are percentages and do not require conversion. The coordinates of borehole ZK05 are... rice, The coordinates fall within an 8-meter-side cell of the finest level of the quadtree variable resolution grid, and the center coordinates of this cell are... rice, Meters. This generates two points in the three-dimensional discrete point source vector: and .

[0069] Another borehole, ZK08, recorded a gold-bearing layer: "Depth 21.50m to 22.20m, gold grade 8.3g / t". Gold grade is measured in grams per ton (g / t), converted to a percentage. Sampling midpoint depth meters, the center coordinates of the adsorption grid cell are rice, Rice, forming point .

[0070] A total of 87 grade samples were extracted from 12 boreholes in the area. After planar grid adsorption and vertical merging, 65 three-dimensional discrete point sources were finally formed. In the four-tuple records, the spatial coordinate unit is meters, and the grade is a dimensionless mass percentage. The maximum iron grade within the coverage area is 42.7 g / t, the maximum copper grade is 2.15 g / t, and the maximum gold grade is 12.6 g / t. These three-dimensional discrete point source vectors will serve as the basic control data for subsequent construction of a continuous mineral distribution probability field, and their spatial distribution accurately reflects the mineralization locations revealed by existing exploration projects.

[0071] S5. Calculate the continuous rate of change of spatial wave impedance of a clear continuous physical feature map of the seabed, construct an anisotropic transmission coefficient field, take the three-dimensional discrete point source vector as the initial source term, and perform spatial iterative diffusion along the anisotropic transmission coefficient field to generate a continuous mineral distribution probability field.

[0072] In a preferred embodiment, the step of calculating the continuous rate of change of spatial wave impedance of a clear continuous physical feature map of the seabed, constructing an anisotropic conduction coefficient field, and using a three-dimensional discrete point source vector as an initial source term to perform spatial iterative diffusion along the anisotropic conduction coefficient field includes: calculating the wave impedance value of each point in space based on the reflection intensity distribution of the clear continuous physical feature map of the seabed; obtaining the gradient of the change of wave impedance value in each direction in space to generate an anisotropic conduction coefficient field related to the spatial location, which characterizes the conduction difference of the geological medium in different directions; and using a three-dimensional discrete point source vector as an initial source term to perform spatial iterative diffusion along the anisotropic conduction coefficient field, including: establishing a Poisson diffusion equation, using the three-dimensional discrete point source vector as an internal source term of the equation, and using the anisotropic conduction coefficient field as the diffusion coefficient, performing iterative calculations until the diffusion energy distribution reaches a preset convergence condition to obtain a continuous mineral distribution probability field.

[0073] Specifically, based on the reflection intensity distribution of the clear, continuous physical feature map of the seabed output in step S3, the wave impedance value at each point in space is calculated. (Clear, continuous physical feature map of the seabed) The comprehensive optical reflectance radiance of the seabed surface was recorded, and the spatial variation of its pixel values ​​simultaneously carries information about the differences in sediment type and shallow strata properties. Wave impedance is a physical quantity describing the acoustic properties of a medium in seismic exploration and acoustic detection, defined as the product of the medium density per unit area and the P-wave velocity, with units of kilograms per square meter per second. However, in the context of optical remote sensing, the relative variation of wave impedance values ​​reflects the differences in medium properties better than their absolute values. Therefore, a linear mapping is used for the conversion from optical reflectance radiance to wave impedance. An initial wave impedance benchmark value for soft seabed surface sediments is set. Kilograms per square second (kg / m² / s), this value is obtained by averaging the statistical results of indoor acoustic measurements of shallow marine silty sediments. Mapping coefficient between reflected radiance and wave impedance. The coefficient, kilograms per square meter per second, was calibrated using least-squares fitting from 12 sets of simultaneously acquired underwater hyperspectral reflectance data and shallow seismic profile wave impedance inversion results. To maintain consistency with the discretized grid, any horizontal grid point in a clear, continuous seabed physical feature map was used. ,in For the column index corresponding to the east-west direction, For the row index corresponding to the north-south direction, the corresponding wave impedance value Calculate using the following formula: In the formula, To clearly visualize the continuous physical features of the seabed at grid points The radiance value, in units of watts per square meter per steradian micrometer; and These are the maximum and minimum global radiance values ​​of the image, respectively, with units of and . Same; the ratio term is the normalized radiance, dimensionless; The dimensions of the product term are such that they have the same dimension as the wave impedance, kilograms per square meter per second, resulting in a consistent overall dimension.

[0074] The above calculations only yield the wave impedance distribution at the seabed surface. Considering that mineral assessment requires three-dimensional spatial information, the two-dimensional wave impedance is extended downwards along the depth direction to construct a three-dimensional wave impedance volume. The extension rule is: for any horizontal grid point... Starting from the absolute elevation of the seabed, proceed vertically downwards at depth intervals. A calculation layer is selected each time, until the maximum depth is reached that traverses all known borehole sampling depths; the layer is numbered. , No. Layer depth , This is the absolute seabed elevation at that location in the absolutely spatially aligned data field obtained in step S1. The vertical profile of the three-dimensional wave impedance volume at this location is set as a constant extension, i.e. This means that optical remote sensing images only reflect the properties of the seabed surface, and the deep wave impedance remains consistent with the seabed for the time being without the constraint of drilling. Subsequent diffusion calculations will readjust the deep probability distribution.

[0075] Based on a three-dimensional wave impedance volume, the gradient of the wave impedance value in each direction in space is calculated to generate an anisotropic conduction coefficient field. For a three-dimensional lattice... wave impedance The gradient components along the horizontal two directions (corresponding to the east and north directions respectively) and the vertical direction are approximated using the central difference: in, The horizontal pixel side length ; For vertical interlayer spacing, the dimensions of the three gradient components are all kilograms per square meter per second per meter. rice.

[0076] The anisotropic conduction coefficient field consists of principal conduction coefficients in three orthogonal directions, denoted as […]. , and There is an inverse relationship between the conduction coefficient and the wave impedance gradient: directions with drastic changes in wave impedance indicate geological interfaces or faults, where mineralization diffusion and connectivity are poor, and the conduction coefficient should be small; directions with gentle changes in wave impedance indicate homogeneous rock layers or bedding planes, and the conduction coefficient is large. Based on this physical principle, the principal conduction coefficients in the three directions are determined by the following formula: In the formula, The substrate diffusion coefficient is . For each iteration step, the maximum convergence value is selected based on a two-dimensional coastal sand mineral diffusion model numerical experiment, while ensuring numerical stability. It is an exponential function with the natural constant e as its base. It is the gradient normalization constant. The units are kilograms per square meter per second per meter, with dimensions consistent with the gradient components. They are determined by the upper quartile of the histogram of the wave impedance volume gradient distribution in the study area, ensuring that the conduction coefficient rapidly decays to less than one-third of the base value when the gradient value is in a high partition. The anisotropic conduction coefficient field is stored in the form of three scalar fields with uniform three-dimensional volume dimensions.

[0077] Please see Figure 3 As shown, Figure 3 This diagram illustrates the attenuation function of the anisotropic conduction coefficient as a function of the acoustic impedance gradient. As shown, the normalized principal conduction coefficient is inversely proportional to the absolute value of the normalized acoustic impedance gradient. In low-gradient regions where acoustic impedance changes are gradual, the system assigns a larger conduction coefficient to achieve strong diffusion smoothing, promoting the extension of grade information. However, in high-gradient regions where acoustic impedance changes drastically, such as at geological interfaces or faults, the conduction coefficient rapidly decays to the baseline value, achieving weak diffusion edge preservation. This mechanism ensures that the spatial diffusion of discrete data is strictly controlled within the actual physical and geological structure.

[0078] Next, the three-dimensional discrete point source vector obtained in step S4 is used as the internal source term of the equation, and spatial iterative diffusion is performed along the anisotropic transmission coefficient field. The three-dimensional discrete point source vector consists of 65 spatial points with quality attributes, which need to be transformed into an initial source term distribution on a three-dimensional grid first. Traverse all quadtree variable resolution mesh cells, for each source point It is assigned to the grid point with the nearest spatial location. , where horizontal grid index From plane coordinates Correspondence determined, Closest in the depth direction The stratum number. The initial source term value for this grid point is set to the grade value. If multiple source points are assigned to the same grid point, the arithmetic mean of the number of items is taken. Grid points without assigned source points are initialized to zero. Initial source items. It is a three-dimensional matrix with its horizontal dimension clearly showing the continuous physical features of the seabed. Figure 1 The depth dimension covers the range from the seabed to the maximum drilling depth.

[0079] A Poisson diffusion equation is established and iteratively evolved. The difference form of the diffusion equation on a discrete three-dimensional volume is:

[0080] In the formula, and The first Second and third The next iteration step grid point The diffusion field value, with dimensions consistent with grade, is a dimensionless mass percentage; The iteration time step is set to 1.0, used to control the magnitude of each update; the interface conduction coefficient is the harmonic average of the principal conduction coefficients of two adjacent grid points, i.e., the terms with a 1 / 2 fractional subscript in the formula, such as... For example, To ensure continuous flux; the principal conduction coefficient is as follows: , Equals represents the intrinsic conductivity of the corresponding lattice point in the corresponding direction. This value is based on the medium properties of the three-dimensional volume itself or a priori parameters.

[0081] The above difference equations are applied only to blank grid points, with the source point being forcibly assigned a fixed boundary condition. For the source grid point mapped in step S4, it is treated as a first-type fixed boundary condition throughout the entire iterative evolution, i.e., forcibly assigned... This allows it to continuously transmit quality information to its surroundings without accumulating or diverging in its own value.

[0082] Iteration from Beginning, initial field The settings have already been configured. Each update step calculates the entire field using the difference equations described above. The spatial boundary treatment of the diffusion equation is as follows: a zero-flux boundary is used in the horizontal direction, meaning the conduction coefficient at the outer boundary of the computational domain is set to zero, preventing mass from escaping the computational area; a zero-flux boundary is used on the top surface of the seabed, excluding air and water above the seabed from calculation; a zero-flux boundary is used on the bottom interface, placed 20 meters below the maximum borehole depth to prevent information loss from the bottom. The convergence criterion is set as the maximum absolute change in the field value between two consecutive iterations being less than a threshold. This threshold, determined through prior numerical experiments, ensures that diffusion is fully complete after the observed millimeter change in grade number. Iteration stops upon reaching convergence, and the total number of iterations is recorded. Output the final field As a continuous mineral distribution probability field, this field has a dimensionless mass percentage and records the mineralization probability intensity obtained by anisotropic diffusion extrapolation guided by physical properties from known borehole mineralization information.

[0083] For example, continuing with the aforementioned bedrock coastal zone example, from a pixel C in a clear map of continuous physical features of the seabed... Watts per square meter per micrometer, maximum value of the entire image area , Normalized radiance is wave impedance kilograms per square meter per second.

[0084] Seabed Elevation of Pixel C in Absolute Space Aligned Data Field Meters, extending downwards to the strata ,depth Meters. This grid point belongs to a hard substrate region near a fault zone, with a horizontal wave impedance gradient. kilograms per square meter per second per meter, vertical gradient Kilograms per square meter per second per meter. The anisotropic conductivity coefficient is calculated using the exponential formula: The square of meters per iteration step; The square meter is measured per iteration step. This shows that the change in wave impedance in the horizontal direction is greater than that in the vertical direction, where the conduction coefficient is slightly smaller, and diffusion preferentially proceeds along the vertical and nearly parallel bedding directions.

[0085] Step S4 gives a point in the three-dimensional discrete point source vector. Corresponding grid points Initial source term value The diffusion equation is passed through... After one iteration, the maximum change in the entire field decreased to Below the convergence threshold Diffusion ceases. The final continuous mineral distribution probability field maintains its highest value at the source point and extends outward along the stratigraphic direction where the wave impedance changes gently. The probability value decreases rapidly at the location where the fault zone passes, forming a mineralization probability distribution cloud map consistent with the geological structure constraints, providing a continuous input field for isosurface extraction in step S6.

[0086] S6. Extract isosurfaces from the probability field of continuous mineral distribution, and check the three-dimensional closed connectivity of the extracted isosurfaces in a quadtree variable resolution grid to generate a three-dimensional ore body boundary model without topological conflicts as the evaluation result.

[0087] In a preferred embodiment, generating a topologically conflict-free 3D ore body boundary model includes: determining the industrial mining probability threshold of the target mineral type; extracting isosurfaces with probability values ​​equal to the industrial mining probability threshold from a continuous mineral distribution probability field; performing 3D closed connectivity detection on the extracted isosurfaces under a spatial index constructed by a quadtree variable resolution mesh; deleting isolated patches that do not form closed volumes; and generating a topologically conflict-free 3D ore body boundary model.

[0088] Specifically, the industrial mining probability threshold for the target mineral is determined from the continuous mineral distribution probability field output in step S5, and isosurfaces with probability values ​​equal to this threshold are extracted. The continuous mineral distribution probability field is a three-dimensional regular grid covering the coastal transition zone. The value stored at each grid point represents the probability intensity that the target mineral exists at that spatial location and its grade meets the economically mineable standard. This field is generated by anisotropic diffusion from the three-dimensional discrete point source vector in step S4, and its dimension is a dimensionless mass percentage. The industrial mining probability threshold is a boundary grade value defined for the target mineral. Ore reserves above this value have economic value for mining. It is set based on a comprehensive consideration of concentrate market prices, comprehensive mining and beneficiation costs, and beneficiation recovery rates, using the net present value break-even method to determine the lower limit of the grade, expressed as a dimensionless mass percentage. For seafloor iron ore, the industrial mining probability threshold is set to 20; for copper ore, it is set to 0.5; and for gold ore, it is set to... This is an approximate value after converting 1.0 gram per ton to a mass percentage. In this step, the target mineral type is specified by the user when the evaluation task is initiated, and the system automatically retrieves the corresponding industrial mining probability threshold, denoted as... .

[0089] From the probability field of continuous mineral distribution Extracting isosurfaces The process is achieved by traversing the three-dimensional volume data cell by cell. The continuous mineral distribution probability field is considered as a cubic array of eight adjacent grid points, the three-dimensional dimensions of which are identical to the three-dimensional impedance volume constructed in step S5. The probability value at each cell vertex is known. For probability values ​​crossing a threshold... An edge, i.e., the probability value of one vertex is greater than The other end is smaller An isosurface must intersect an edge. The location of the intersection point is determined using linear interpolation, where the two vertices of the edge are denoted as . and Spatial coordinates are denoted as follows: and The corresponding probability value is and Intersection coordinates Calculate using the following formula:

[0090] In the formula, , , All values ​​are dimensionless, and the result of subtraction or division is also dimensionless. Multiplying by the coordinate difference yields a length unit consistent with the spatial coordinates, making the formula dimensionally homogeneous. The intersection points are calculated by traversing all edges crossing the threshold within the cube cell. These intersection points are then connected to form triangular faces according to the pre-established lookup table rules for 256 vertex states in the classic Marching Cubes algorithm. The set of all triangular faces constitutes the threshold. The corresponding isosurface triangular mesh. This mesh approximates the ore body boundary in the form of three-dimensional facets, with each triangular facet vertex carrying its three-dimensional spatial coordinates.

[0091] Based on the extracted isosurface triangular mesh, three-dimensional closed connectivity detection needs to be performed under the spatial index of the quadtree variable resolution mesh. Using the quadtree variable resolution mesh established in step S1, the three-dimensional volume space of the continuous mineral distribution probability field is mapped onto the mesh plane, and attached to the corresponding mesh cells along the vertical direction at depth levels. In this way, each triangular facet of the isosurface can be located to a unique quadtree leaf node cell based on the horizontal coordinates of its centroid. Through the parent-child adjacency relationship between cells, potential adjacent faces of the facet can be found in a constant amount of time, accelerating the subsequent construction of the connectivity graph.

[0092] 3D closed connectivity detection is achieved by constructing a triangular patch adjacency graph. Each triangle has three edges. If two triangles share the same edge (i.e., the spatial coordinates of the two endpoints of the edge completely coincide within a preset floating-point precision tolerance), then they are adjacent patches. All triangles are traversed, and a mapping from edge identifier to triangle index is established for each edge. The edge identifier is represented by an ascending or descending ordered pair of the indices of the two vertices of the edge, ensuring that edges sharing the same vertex receive the same identifier. After establishing the mapping, each triangle can quickly query its three adjacent triangle lists, forming an undirected adjacency graph.

[0093] Next, a breadth-first search is used to identify all connected components on the adjacency graph. Starting from any unvisited triangle, all reachable triangles are traversed along adjacency relationships. All triangles within the same connected component share the same connected component label. After obtaining all connected components, each connected component is checked to see if it is a closed volume. The criterion for a closed volume is: within the connected component, each edge must be referenced by exactly two triangles, i.e., there are no boundary edges referenced by only one triangle. All triangles contained in the connected component are scanned, and the number of times each edge is referenced is counted. If an edge has a reference count of 1, it indicates that the boundary is open, and the connected component cannot form a finite closed volume. It is marked as an isolated face and deleted.

[0094] After deleting isolated facets that do not form closed volumes, all remaining closed connected components are merged, outputting a topologically conflict-free 3D orebody boundary model. This topologically conflict-free 3D orebody boundary model is a 2D manifold surface composed of triangular facets. It is completely closed, without holes, free boundaries, or self-intersections, and is strictly extracted according to the isosurface of the continuous mineral distribution probability field. Its spatial reference elevation datum is consistent with the absolute spatial alignment data field. This model can be directly used for reserve estimation, 3D visualization, and mining design.

[0095] For example, continuing with the aforementioned bedrock coastal zone assessment example, the target mineral is magnetite, and the industrial mining probability threshold is... The value is set to 20. In the continuous mineral distribution probability field output by step S5, the probability value of the core ore section controlled by the borehole reaches the highest level of 42.7, decreasing outwards and dropping below 20 at the edge of the ore body. Isosurface extraction uses linear interpolation to traverse all approximately... Each cube unit generates approximately [number] triangular facets. After constructing the adjacency graph under the quadtree variable resolution grid index, connected component detection found a total of 4 connected components, of which component 1 contains A series of triangles, with no free edges after boundary checks, constitute a closed ore body; components 2, 3, and 4 each contain approximately , and These are triangular patches, each with a boundary edge, and are considered isolated facets. These isolated facets are mostly located in seafloor depressions or localized enrichment areas of shallow weathered layers, and are not connected to the main ore body. Because they do not form a closed volume, they are automatically deleted. Please refer to [link to relevant documentation]. Figure 4 As shown, Figure 4This diagram illustrates the comparison of the number of triangular facets in the three-dimensional closed connectivity detection results of the isosurface of this invention. As shown, due to the significant differences in the volume of each component, the vertical axis is displayed on a logarithmic coordinate. Among the four extracted connected components, connected component 1 has a large number of triangular facets and, after detection, has no free boundaries, so it is determined to be a closed ore body and is retained; while connected components 2, 3, and 4 have fewer facets and, due to the presence of boundary edges that prevent closure, are determined to be isolated facets and are discarded. This processing logic effectively filters out local interfering facets, ensuring that the final generated three-dimensional ore body boundary model has no topological conflicts. The final generated three-dimensional ore body boundary model without topological conflicts contains only component 1, depicting a layered magnetite ore body boundary approximately 860 meters long, 420 meters wide, and with a maximum vertical thickness of 38 meters. Its surface is topologically closed and can be directly input into geological modeling software for resource quantity calculation.

[0096] Please see Figure 5 As shown, the second aspect of the present invention provides a multi-source heterogeneous data fusion evaluation system for coastal mineral resources, comprising: a spatial alignment module, a frequency domain decomposition module, a guided enhancement module, a point source mapping module, an iterative diffusion module, and an ore body modeling module.

[0097] The spatial alignment module is used to acquire land-based elevation vector point data, multibeam acoustic water depth raster data, and geological borehole text data. It constructs a water depth compensation surface based on historical tidal extreme values, inverts the absolute elevation of the multibeam acoustic water depth raster data, and establishes a quadtree variable resolution grid for topological alignment to obtain an absolute spatial alignment data field.

[0098] The frequency domain decomposition module is used to acquire shallow sea area optical remote sensing images that are spatially matched with the data field aligned with the absolute space. The optical remote sensing images are converted to the frequency domain through Fourier transform, and the high-frequency components affected by water turbidity are separated to obtain the optical high-frequency matrix to be compensated and the low-frequency background components.

[0099] The guided enhancement module is used to extract the seabed backscattering texture gradient field of the multibeam acoustic grid data in the absolute spatial alignment data field. The seabed backscattering texture gradient field is used as the guiding matrix to perform directional edge enhancement on the optical high-frequency matrix to be compensated. The enhanced high-frequency matrix is ​​merged with the low-frequency background component and subjected to inverse Fourier transform to obtain a clear seabed continuous physical feature map.

[0100] The point source mapping module is used to extract mineral composition keywords and grade values ​​from geological borehole text data, and combine the borehole opening coordinates and sampling depth to map the grade values ​​onto a quadtree variable resolution grid to obtain a three-dimensional discrete point source vector.

[0101] The iterative diffusion module is used to calculate the continuous rate of change of spatial wave impedance of a clear continuous physical feature map of the seabed, construct an anisotropic transmission coefficient field, take the three-dimensional discrete point source vector as the initial source term, and perform spatial iterative diffusion along the anisotropic transmission coefficient field to generate a continuous mineral distribution probability field.

[0102] The ore body modeling module is used to extract isosurfaces from the probability field of continuous mineral distribution, and to check the three-dimensional closed connectivity of the extracted isosurfaces in a quadtree variable resolution grid, generating a three-dimensional ore body boundary model without topological conflicts as the evaluation result.

[0103] Each of the modules can be implemented in whole or in part through software, hardware, or a combination thereof. It supports hardware embedded in or independent of the processor in the computer device, and also supports software stored in the memory of the computer device, so that the processor can call and execute the operations corresponding to each of the above modules.

[0104] The above embodiments are only used to illustrate the technical solutions of the present invention, and are not intended to limit it. Although the present invention has been described in detail with reference to the foregoing embodiments, those skilled in the art should understand that modifications can still be made to the technical solutions described in the foregoing embodiments, or equivalent substitutions can be made to some of the technical features. Such modifications or substitutions do not cause the essence of the corresponding technical solutions to deviate from the spirit and scope of the technical solutions of the embodiments of the present invention, and should all be included within the protection scope of the present invention.

Claims

1. A method for evaluating coastal mineral resources through multi-source heterogeneous data fusion, characterized in that, include: S1. Acquire land-based elevation vector point data, multibeam acoustic water depth raster data, and geological borehole text data. Construct a water depth compensation surface based on historical tidal extreme values, invert the absolute elevation of the multibeam acoustic water depth raster data, and establish a quadtree variable resolution grid for topological alignment to obtain an absolute spatially aligned data field. S2. Obtain the optical remote sensing image of the shallow sea area that is spatially matched with the data field aligned with the absolute space. Transform the optical remote sensing image to the frequency domain through Fourier transform and separate the high-frequency component affected by water turbidity to obtain the optical high-frequency matrix to be compensated and the low-frequency background component. S3. Extract the seabed backscattering texture gradient field from the multibeam acoustic grid data in the absolute spatial alignment data field. Use the seabed backscattering texture gradient field as the guiding matrix to perform directional edge enhancement on the optical high-frequency matrix to be compensated. Merge the enhanced high-frequency matrix with the low-frequency background component and perform inverse Fourier transform to obtain a clear seabed continuous physical feature map. S4. Extract mineral composition keywords and grade values ​​from geological borehole text data. Combine the borehole opening coordinates and sampling depth, map the grade values ​​to a quadtree variable resolution grid to obtain a three-dimensional discrete point source vector. S5. Calculate the continuous rate of change of spatial wave impedance of a clear continuous physical feature map of the seabed, construct an anisotropic transmission coefficient field, take the three-dimensional discrete point source vector as the initial source term, and perform spatial iterative diffusion along the anisotropic transmission coefficient field to generate a continuous mineral distribution probability field. S6. Extract isosurfaces from the probability field of continuous mineral distribution, and check the three-dimensional closed connectivity of the extracted isosurfaces in a quadtree variable resolution grid to generate a three-dimensional ore body boundary model without topological conflicts as the evaluation result.

2. The method for evaluating coastal mineral resources through multi-source heterogeneous data fusion according to claim 1, characterized in that, The process of constructing a water depth compensation surface based on historical tidal extremes, inverting the absolute elevation of multibeam acoustic water depth raster data, and establishing a quadtree variable resolution grid for topological alignment includes: acquiring historical tidal extreme value records of the area covered by land-based elevation vector point data, and generating a water depth compensation surface that reflects tidal changes. The absolute elevation is obtained by subtracting the relative water depth of each grid cell in the multibeam acoustic water depth grid data from the compensation value of the corresponding position on the water depth compensation surface. A quadtree variable resolution grid covering the coastal transition zone is constructed. Within the same grid level, the inverted acoustic raster pixels and the land-based elevation vector points are spatially topologically snapped together to eliminate elevation jumps and output an absolutely spatially aligned data field.

3. The method for evaluating coastal mineral resources through multi-source heterogeneous data fusion according to claim 2, characterized in that, The historical tidal extreme value record includes long-term tidal level observation data obtained from coastal tide gauge stations, and the water depth compensation surface is obtained by interpolating the highest and lowest tidal levels in the long-term tidal level observation data.

4. The method for evaluating coastal mineral resources through multi-source heterogeneous data fusion according to claim 1, characterized in that, The separation of high-frequency components affected by water turbidity yields the optical high-frequency matrix to be compensated and the low-frequency background component, including: A water suspension attenuation function model corresponding to optical remote sensing images of shallow sea areas was constructed. This model expresses the attenuation characteristics of optical signals in the spatial frequency domain under different water turbidity. The frequency domain image after Fourier transform is compared with the water suspended matter attenuation function model to separate the optical high-frequency components that are attenuated due to water scattering and absorption, which are used as the optical high-frequency matrix to be compensated, while the low-frequency background components that are not attenuated are retained.

5. The method for evaluating coastal mineral resources by fusing multi-source heterogeneous data according to claim 4, characterized in that, The water suspended matter decay function model is an exponential decay function, and the independent variables of the exponential decay function include the spatial frequency component and the spatial cutoff frequency. The spatial cutoff frequency is inversely proportional to the product of the water depth and turbidity at the corresponding location.

6. The method for evaluating coastal mineral resources through multi-source heterogeneous data fusion according to claim 1, characterized in that, The extraction of the seabed backscattering texture gradient field from the multibeam acoustic raster data in the absolute spatial alignment data field, and the use of the seabed backscattering texture gradient field as a guiding matrix for directional edge enhancement of the optical high-frequency matrix to be compensated, includes: Extract backscattering intensity information from multibeam acoustic raster data from an absolutely spatially aligned data field; The spatial variation rate of backscattering intensity information is calculated to generate a gradient field characterizing the direction and intensity of seabed micro-topographic texture, which serves as the seabed backscattering texture gradient field. The gradient field of the seabed backscattering texture is used as the guiding matrix for the guided filter. The edge components in the high-frequency optical matrix to be compensated that are consistent with the gradient field direction are retained, while water scattering noise that is inconsistent with the gradient field direction is suppressed, resulting in the enhanced high-frequency matrix.

7. The method for evaluating coastal mineral resources through multi-source heterogeneous data fusion according to claim 1, characterized in that, The process involves extracting mineral composition keywords and grade values ​​from geological borehole text data, combining the borehole opening coordinates and sampling depth, and mapping the grade values ​​onto a quadtree variable resolution grid to obtain a three-dimensional discrete point source vector, including: Using natural language processing named entity recognition technology, we extract the component keywords corresponding to the target minerals and the grade values ​​associated with the component keywords from geological borehole text data; Obtain the borehole coordinates and sampling depth of the geological borehole. Using the borehole coordinates and sampling depth as spatial locations, set the extracted grade values ​​as attributes at these spatial locations to form a three-dimensional discrete point source vector. The spatial reference of the three-dimensional discrete point source vector is consistent with the absolute spatial alignment data field.

8. The method for evaluating coastal mineral resources through multi-source heterogeneous data fusion according to claim 1, characterized in that, The calculation of the continuous rate of change of spatial wave impedance in a clear, continuous physical feature map of the seabed, the construction of an anisotropic transmission coefficient field, and the use of a three-dimensional discrete point source vector as the initial source term, along with spatial iterative diffusion along the anisotropic transmission coefficient field, include: Based on the reflection intensity distribution of a clear, continuous physical feature map of the seabed, the wave impedance value at each point in space is calculated. The gradient of wave impedance value in various spatial directions is calculated to generate an anisotropic conduction coefficient field related to spatial location. This field characterizes the conduction differences of the geological medium in different directions. A Poisson diffusion equation is established, with the three-dimensional discrete point source vector as the internal source term of the equation and the anisotropic conduction coefficient field as the diffusion coefficient. Iterative calculations are performed until the diffusion energy distribution reaches the preset convergence condition, thus obtaining the continuous mineral distribution probability field.

9. The method for evaluating coastal mineral resources through multi-source heterogeneous data fusion according to claim 1, characterized in that, The generation of a topologically conflict-free 3D orebody boundary model includes: Determine the industrial mining probability threshold for the target mineral type, and extract isosurfaces with probability values ​​equal to the industrial mining probability threshold from the continuous mineral distribution probability field; Under the spatial index of the quadtree variable resolution mesh, the extracted isosurfaces are subjected to three-dimensional closed connectivity detection, isolated patches that do not form closed volumes are deleted, and a three-dimensional ore body boundary model without topological conflicts is generated.

10. A multi-source heterogeneous data fusion and evaluation system for coastal mineral resources, characterized in that, include: The spatial alignment module is used to acquire land-based elevation vector point data, multibeam acoustic water depth raster data, and geological borehole text data. It constructs a water depth compensation surface based on historical tidal extremes, inverts the absolute elevation of the multibeam acoustic water depth raster data, and establishes a quadtree variable resolution grid for topological alignment to obtain an absolute spatial alignment data field. The frequency domain decomposition module is used to acquire shallow sea area optical remote sensing images that are spatially matched with the data field aligned with the absolute space. The optical remote sensing images are converted to the frequency domain through Fourier transform, and the high-frequency components affected by water turbidity are separated to obtain the optical high-frequency matrix to be compensated and the low-frequency background components. The guided enhancement module is used to extract the seabed backscattering texture gradient field of the multibeam acoustic grid data in the absolute spatial alignment data field. The seabed backscattering texture gradient field is used as the guiding matrix to perform directional edge enhancement on the optical high frequency matrix to be compensated. The enhanced high frequency matrix is ​​merged with the low frequency background component and subjected to inverse Fourier transform to obtain a clear seabed continuous physical feature map. The point source mapping module is used to extract mineral composition keywords and grade values ​​from geological borehole text data. Combined with borehole opening coordinates and sampling depth, the grade values ​​are mapped to a quadtree variable resolution grid to obtain a three-dimensional discrete point source vector. The iterative diffusion module is used to calculate the continuous rate of change of spatial wave impedance of a clear continuous physical feature map of the seabed, construct an anisotropic transmission coefficient field, take the three-dimensional discrete point source vector as the initial source term, and perform spatial iterative diffusion along the anisotropic transmission coefficient field to generate a continuous mineral distribution probability field. The orebody modeling module is used to extract isosurfaces from the probability field of continuous mineral distribution, and to check the three-dimensional closed connectivity of the extracted isosurfaces in a quadtree variable resolution grid, generating a three-dimensional orebody boundary model without topological conflicts as the evaluation result.

Images