Coalfield aquifer spatial distribution analysis method based on surveying and mapping fusion

Abstract

The present application relates to the technical field of geological surveying and mapping, in particular to a coalfield aquifer spatial distribution analysis method based on surveying and mapping fusion, which comprises: performing spatial registration on the original surveying and mapping data set of a target region to obtain spatialized surveying and mapping data; obtaining spatial difference partitions through three-dimensional reconstruction of aquifer lithology and thickness and correlation mapping of water conductivity performance, determining boundary information by combining aquifer outcrop and structure point overlay analysis, and constructing a digital spatial structure; updating the structure through multi-source real-time monitoring data parameter inversion, interpreting water level anomaly spatial characteristics, positioning risk areas in combination with water conductivity performance partition, inferring trends by fusing water inflow records, and generating a risk distribution information map; the present application can improve the efficiency of coalfield aquifer spatial distribution analysis based on surveying and mapping fusion.

Description

Technical Field

[0001] This invention relates to the field of geological surveying and mapping technology, and in particular to a method for analyzing the spatial distribution of coalfield aquifers based on surveying and mapping fusion. Background Technology

[0002] Existing analyses of the spatial distribution of coalfield aquifers mostly employ a single mapping data processing mode. Multi-source mapping information cannot achieve unified spatial registration and efficient fusion. Core attribute data such as aquifer lithology and thickness are difficult to reconstruct in three dimensions. The correlation mapping between lithology, thickness, and water conductivity lacks a systematic method. The spatial difference zoning results have low accuracy and cannot truly reflect the spatial distribution characteristics of aquifers.

[0003] Traditional analysis methods fail to accurately overlay spatial zoning results with aquifer outcrop and structural point data, resulting in biased aquifer boundary definitions. The construction of digital spatial structures lacks a dynamic correction mechanism, and multi-source monitoring data cannot be effectively combined with spatial models for parameter inversion. The identification of water level anomalies, location of risk areas, and projection of risk distribution are time-consuming. Overall, the analysis results and timeliness cannot meet the needs of practical applications. Therefore, improving the efficiency of spatial distribution analysis of coalfield aquifers has become an urgent problem to be solved. Summary of the Invention

[0004] This invention provides a method for analyzing the spatial distribution of coalfield aquifers based on mapping fusion, in order to solve the problems mentioned in the background art.

[0005] To achieve the above objectives, the present invention provides a method for analyzing the spatial distribution of coalfield aquifers based on surveying and mapping fusion, comprising:

[0006] A1. Spatial registration is performed on the original mapping dataset of the target area to obtain the spatialized mapping dataset of the target area;

[0007] A2. Reconstruct the lithology and thickness attribute data of the aquifer in the spatialized mapping dataset in three dimensions, and associate and map the reconstructed three-dimensional spatial distribution trend information with the water-conducting performance attribute data of the aquifer in the spatialized mapping dataset to obtain the spatial difference zoning information of the target area.

[0008] A3. Perform spatial overlay analysis on the spatial difference zoning information and the aquifer outcrop and structural point data in the original mapping dataset to obtain the boundary information set of the target area;

[0009] A4. Construct a digital spatial structure of the target area based on the boundary information set and the spatialized mapping dataset;

[0010] A5. Perform parameter inversion between the multi-source real-time monitoring data of the target area and the digital spatial structure to obtain the updated spatial structure of the target area;

[0011] A6. Interpret the comparison results between water level observation data and theoretical water level distribution information from the updated spatial structure to generate the spatial characteristics of water level anomalies in the target area, and extract the risk area location information of the target area by combining the water-conducting performance zoning information in the updated spatial structure.

[0012] A7. Perform trend extrapolation between the risk area location information and the water inflow record data in the updated spatial structure to construct a risk distribution information map of the target area.

[0013] In a preferred embodiment, the step of spatially registering the original mapping dataset of the target area to obtain a spatialized mapping dataset of the target area includes:

[0014] Obtain the original mapping dataset of the target area and parse the original coordinate reference information of each data source in the original mapping dataset. The original coordinate reference information includes geodetic coordinate system parameters and map projection parameters.

[0015] Based on the geographical extent of the target area, the projection parameters of the standard spatial reference frame are extracted to generate a unified target coordinate reference system for the target area.

[0016] Based on the unified target coordinate reference system, the original coordinate reference information is spatially unified to obtain a dataset in the transition coordinate system of the target region.

[0017] Spatial feature points in the dataset under the transition coordinate system are spatially matched with the standard geographic base map of the target area to obtain a spatialized mapping dataset of the target area.

[0018] In a preferred embodiment, the step of performing three-dimensional reconstruction of the aquifer lithology and thickness attribute data in the spatialized mapping dataset, and then mapping the reconstructed three-dimensional spatial distribution trend information with the aquifer hydraulic conductivity attribute data in the spatialized mapping dataset to obtain the spatial difference zoning information of the target area includes:

[0019] Extract the aquifer lithology description data and aquifer thickness numerical data corresponding to the sampling points from the spatialized mapping dataset to generate a discrete point attribute set for the target area;

[0020] Spatial interpolation is performed on the sampling points in the discrete point attribute set to construct a continuous voxel model of the aquifer thickness in the target region.

[0021] Using the aquifer lithology description data as labels, values ​​are assigned to the aquifer thickness continuum voxel model on a voxel-by-voxel basis to obtain the three-dimensional aquifer structure of the target region.

[0022] The lithological variation gradient and thickness variation gradient at each voxel location are analyzed from the three-dimensional aquifer structure, and the lithological variation gradient and thickness variation gradient are tensor synthesized to obtain the three-dimensional spatial distribution trend information of the target area.

[0023] The three-dimensional spatial distribution trend information is correlated point by point with the aquifer water conductivity attribute data extracted from the spatialized mapping dataset to obtain the correspondence table of the target area;

[0024] Based on the correspondence table, the three-dimensional aquifer structure is divided into regions to obtain the spatial difference partitioning information of the target region.

[0025] In a preferred embodiment, the step of resolving the lithological variation gradient and thickness variation gradient at each voxel location from the three-dimensional aquifer structure, and then performing tensor synthesis of the lithological variation gradient and thickness variation gradient to obtain the three-dimensional spatial distribution trend information of the target region includes:

[0026] Traverse each voxel in the three-dimensional aquifer structure, and within a preset neighborhood range, obtain the lithological label sequence and thickness value sequence of each adjacent voxel of the current voxel.

[0027] The lithological variation gradient vector of the current voxel is obtained by measuring the difference between the lithological label of the current voxel and the lithological label sequence of the adjacent voxels.

[0028] The thickness value of the current voxel is compared with the thickness value sequence of the adjacent voxels by performing neighborhood difference to obtain the thickness change gradient vector of the current voxel.

[0029] The lithological variation gradient vector and the thickness variation gradient vector are bilinearly aggregated to obtain the local tensor synthesis result of the current voxel;

[0030] The local tensor synthesis results of all voxels in the three-dimensional aquifer structure are arranged and combined to obtain the three-dimensional spatial distribution trend information of the target region.

[0031] In a preferred embodiment, the step of spatially overlaying the spatially differentiated zoning information with the aquifer outcrop and structural point data in the original mapping dataset to obtain the boundary information set of the target area includes:

[0032] The spatial extent outline of each partition unit in the spatial difference partitioning information is used as the partition boundary outline set of the target region.

[0033] Locate the spatial coordinates of aquifer outcrops and the spatial coordinates and structural type attributes of structural points in the original mapping dataset, and compile them into a discrete geological feature point set for the target area;

[0034] The spatial overlay determination is performed between the partition boundary contour set and the discrete geological feature point set to obtain the overlay associated point set of the target region.

[0035] Based on the structural type attributes in the set of superimposed associated points, identify the partition boundary segments adjacent to the aquifer outcrop and specific structural points, and merge the partition boundary segments into the boundary fragment set of the target region;

[0036] The boundary fragment set is topologically stitched together to obtain the boundary information set of the target region.

[0037] In a preferred embodiment, constructing the digital spatial structure of the target area based on the boundary information set and the spatial mapping dataset includes:

[0038] The spatial coordinate sequence of the boundary lines in the boundary information set is spatially enclosed to obtain the closed spatial constraint domain of the target region.

[0039] The aquifer attribute data in the spatialized mapping dataset is mapped to discrete grid nodes within the closed spatial constraint domain to construct the initial field for attribute assignment of the target region.

[0040] Based on the boundary line, the grid cells in the initial field of the attribute assignment are subjected to boundary adaptation and clipping to obtain the boundary-fitting grid of the target region;

[0041] The boundary-fitting mesh is assembled into a digital spatial structure for the target region.

[0042] In a preferred embodiment, the step of performing parameter inversion between the multi-source real-time monitoring data of the target area and the digital spatial structure to obtain the updated spatial structure of the target area includes:

[0043] The dynamic water level change sequence and water quality response characteristics in the multi-source real-time monitoring data of the target area are compared with the corresponding initial hydrogeological parameters in the digital spatial structure to obtain the parameter sensitivity index set of the target area.

[0044] Based on the parameter sensitivity index set, the weights of the non-homogeneous parameter fields in the digital spatial structure are adjusted.

[0045] Based on the adjusted heterogeneous parameter field, the boundary constraints of the digital spatial structure are dynamically modified to obtain the intermediate transitional spatial structure of the target region;

[0046] The intermediate transition space structure is logically consistent with the historical water inrush event records in the multi-source real-time monitoring data. If the verification is successful, the intermediate transition space structure is determined as the updated space structure of the target area.

[0047] In a preferred embodiment, the step of interpreting the comparison results between water level observation data and theoretical water level distribution information from the updated spatial structure to generate spatial characteristics of water level anomalies in the target area, and combining this with the water-conducting performance zoning information in the updated spatial structure to extract risk area location information of the target area, includes:

[0048] Read water level observation data and corresponding theoretical water level distribution information from the updated spatial structure, and perform grid-by-grid difference calculation on the water level observation data and corresponding theoretical water level distribution information to obtain the water level residual dataset of the target area.

[0049] Based on the sign and magnitude of the residual values ​​at each grid point in the water level residual dataset, the water level residual dataset is divided into a positive residual subset and a negative residual subset.

[0050] Spatial clustering analysis is performed on the positive residual subset and the negative residual subset respectively, and the analyzed positive residual spatial blocks and negative residual spatial blocks are jointly labeled as the spatial features of water level anomalies in the target area;

[0051] Extract water-conducting performance zoning information from the updated spatial structure. The water-conducting performance zoning information includes a high water-conducting performance region, a medium water-conducting performance region, and a low water-conducting performance region.

[0052] The positive residual spatial blocks in the spatial features of the water level anomaly are spatially overlaid with the high water conductivity regions in the water conductivity zoning information, and the overlapping regions are identified as the risk area location information of the target area.

[0053] In a preferred embodiment, the formula for calculating the water level residual value of each grid point in the water level residual dataset is as follows:

[0054] ;

[0055] In the formula, For the first Water level residual value at each grid point For the first Water level observation data at each grid point For the first Theoretical water level distribution information for each grid point The preset spatial smoothing coefficient, For the first The set of neighboring grid points of a grid point The number of grid points in the neighboring grid point set. For the first Water level observation data from neighboring grid points, For the first Theoretical water level distribution information for each neighboring grid point.

[0056] In a preferred embodiment, the step of performing trend extrapolation between the risk area location information and the updated water inflow record data in the spatial structure to construct a risk distribution information map of the target area includes:

[0057] Based on the spatial location of each risk block in the risk area location information, extract the water inflow record data at the corresponding location from the updated spatial structure to generate a risk-related water inflow time series sequence for the target area;

[0058] Trend fitting is performed on the time series of risk-related water inflow to extract the direction and intensity of the water inflow change trend for each risk block;

[0059] The direction and intensity of the water inflow change trend are spatially correlated with the risk area location information to construct a risk distribution information map of the target area.

[0060] Compared with the prior art, the present invention has the following beneficial effects:

[0061] 1. This invention, through spatial registration processing of original surveying data, can quickly establish a surveying dataset with a unified spatial benchmark, achieving efficient integration and standardized conversion of multi-source surveying information. Relying on 3D reconstruction technology of aquifer lithology and thickness data, it can accurately construct the 3D spatial structure of aquifers, clearly analyze the spatial distribution patterns of aquifers, and then, through correlation mapping with hydraulic performance data, quickly generate accurate spatial difference zoning results. Combined with spatial overlay analysis of aquifer outcrop and structural point data, it can accurately define aquifer boundary information, construct a digital spatial structure that conforms to actual geological conditions, and through parameter inversion of multi-source real-time monitoring data, it can achieve dynamic updates of the spatial structure, allowing model parameters to match the actual hydrogeological conditions in real time, comprehensively improving the accuracy of spatial modeling and the smoothness of data processing.

[0062] 2. This invention can quickly interpret the discrepancies between observed water level data and theoretical water level distribution, automatically generate spatial characteristics of water level anomalies, and accurately pinpoint risk areas by combining water-conducting performance zoning information. Furthermore, through trend extrapolation from recorded water inflow data, it efficiently constructs a risk distribution information map of coalfield aquifers. The entire technology achieves fully automated operation from data fusion and model construction to risk identification and map generation, significantly improving the efficiency of spatial distribution analysis of coalfield aquifers based on surveying and mapping fusion, while ensuring the accuracy and reliability of the analysis results. This method can provide complete technical support and data basis for refined management of coalfield aquifers and early prevention and control of water hazards, effectively meeting the practical application needs of coalfield geological exploration and safe production. Attached Figure Description

[0063] Figure 1 This is a flowchart illustrating a method for analyzing the spatial distribution of coalfield aquifers based on mapping fusion, provided in an embodiment of the present invention.

[0064] The realization of the objective, functional features and advantages of the present invention will be further explained in conjunction with the embodiments and with reference to the accompanying drawings. Detailed Implementation

[0065] It should be understood that the specific embodiments described herein are merely illustrative of the invention and are not intended to limit the invention.

[0066] This application provides a method for analyzing the spatial distribution of coalfield aquifers based on surveying and mapping fusion. The executing entity of this method includes, but is not limited to, at least one of the following electronic devices that can be configured to execute the method provided in this application: a server, a terminal, etc. In other words, the method can be executed by software or hardware installed on a terminal device or a server device. The server includes, but is not limited to, a single server, a server cluster, a cloud server, or a cluster of cloud servers. The server can be an independent server or a cloud server providing basic cloud computing services such as cloud services, cloud databases, cloud computing, cloud functions, cloud storage, network services, cloud communication, middleware services, domain name services, security services, content delivery networks (CDNs), and big data and artificial intelligence platforms.

[0067] Reference Figure 1 The diagram shown is a flowchart illustrating a method for analyzing the spatial distribution of coalfield aquifers based on surveying and mapping fusion, according to an embodiment of the present invention. In this embodiment, the method for analyzing the spatial distribution of coalfield aquifers based on surveying and mapping fusion includes:

[0068] A1. Spatial registration is performed on the original mapping dataset of the target area to obtain the spatialized mapping dataset of the target area;

[0069] In this embodiment of the invention, the step of spatially registering the original mapping dataset of the target area to obtain the spatialized mapping dataset of the target area includes:

[0070] Obtain the original mapping dataset of the target area and parse the original coordinate reference information of each data source in the original mapping dataset. The original coordinate reference information includes geodetic coordinate system parameters and map projection parameters.

[0071] Based on the geographical extent of the target area, the projection parameters of the standard spatial reference frame are extracted to generate a unified target coordinate reference system for the target area.

[0072] Based on the unified target coordinate reference system, the original coordinate reference information is spatially unified to obtain a dataset in the transition coordinate system of the target region.

[0073] Spatial feature points in the dataset under the transition coordinate system are spatially matched with the standard geographic base map of the target area to obtain a spatialized mapping dataset of the target area.

[0074] The data collection range is locked by defining the boundary of the target area. Within this range, ground-based and aerial surveying equipment is used to collect corresponding spatial location and topographic information, forming a complete original surveying dataset. The data sources, such as satellite imagery, laser point cloud, and terrain vector data, contained in the original surveying dataset are decomposed one by one. The geodetic coordinate system parameters and map projection parameters corresponding to each type of data source are identified and extracted one by one. The above parameters are then summarized to form the original coordinate reference information of each data source.

[0075] Based on the pre-defined latitude and longitude range of the target area, the standard spatial reference frame type of the area is determined. The projection parameters that match the target area are extracted from the nationally unified standard spatial reference frame document. The extracted projection parameters are combined with the preset central meridian scale coefficient coordinate origin value to generate a unified target coordinate reference system applicable to the target area.

[0076] Using a unified target coordinate reference system as the benchmark for spatial transformation, the original coordinate reference information obtained in the early stage is transformed parameter by parameter. The geodetic coordinate system parameters of different data sources are uniformly transformed to the geodetic benchmark corresponding to the unified target coordinate reference system. The map projection parameters of different data sources are uniformly adjusted to the projection rules corresponding to the unified target coordinate reference system. The spatial benchmark of all data sources is uniformly processed to form a dataset under the transition coordinate system corresponding to the target area.

[0077] Extract all feature point information used to identify spatial location from the dataset under the transition coordinate system. Position each spatial feature point according to the position rules corresponding to the unified target coordinate reference system. Compare the position of each positioned spatial feature point with the corresponding ground feature point on the standard geographic base map of the target area. When the position deviation between the spatial feature point and the corresponding ground feature point on the standard geographic base map is less than the preset position matching threshold, it is determined to be a successful match. All successfully matched spatial feature points are integrated with the corresponding surveying and mapping data to form a spatialized surveying and mapping dataset of the target area.

[0078] The beneficial effects are that by parsing the coordinate reference information and unifying the spatial benchmark of the original surveying and mapping dataset, the differences in coordinate systems between different data sources can be eliminated, ensuring the consistency of surveying and mapping data under the same spatial benchmark. By accurately matching the spatial feature points with the standard geographic base map, the spatial positioning accuracy of the surveying and mapping data can be improved, ensuring that the final spatialized surveying and mapping dataset can completely and accurately reflect the real geographic spatial form of the target area. At the same time, the entire data processing flow is clear and controllable, can be directly reproduced, and will not have data deviations or information loss, thus meeting the application requirements of high-precision spatialized surveying and mapping data.

[0079] A2. Reconstruct the lithology and thickness attribute data of the aquifer in the spatialized mapping dataset in three dimensions, and associate and map the reconstructed three-dimensional spatial distribution trend information with the water-conducting performance attribute data of the aquifer in the spatialized mapping dataset to obtain the spatial difference zoning information of the target area.

[0080] In this embodiment of the invention, the step of performing three-dimensional reconstruction of the aquifer lithology and thickness attribute data in the spatialized mapping dataset, and then associating and mapping the reconstructed three-dimensional spatial distribution trend information with the aquifer hydraulic conductivity attribute data in the spatialized mapping dataset to obtain the spatial difference zoning information of the target area includes:

[0081] Extract the aquifer lithology description data and aquifer thickness numerical data corresponding to the sampling points from the spatialized mapping dataset to generate a discrete point attribute set for the target area;

[0082] Spatial interpolation is performed on the sampling points in the discrete point attribute set to construct a continuous voxel model of the aquifer thickness in the target region.

[0083] Using the aquifer lithology description data as labels, values ​​are assigned to the aquifer thickness continuum voxel model on a voxel-by-voxel basis to obtain the three-dimensional aquifer structure of the target region.

[0084] The lithological variation gradient and thickness variation gradient at each voxel location are analyzed from the three-dimensional aquifer structure, and the lithological variation gradient and thickness variation gradient are tensor synthesized to obtain the three-dimensional spatial distribution trend information of the target area.

[0085] The three-dimensional spatial distribution trend information is correlated point by point with the aquifer water conductivity attribute data extracted from the spatialized mapping dataset to obtain the correspondence table of the target area;

[0086] Based on the correspondence table, the three-dimensional aquifer structure is divided into regions to obtain the spatial difference partitioning information of the target region.

[0087] The process of resolving the lithological variation gradient and thickness variation gradient at each voxel location from the three-dimensional aquifer structure, and then performing tensor synthesis of the lithological variation gradient and thickness variation gradient to obtain the three-dimensional spatial distribution trend information of the target region includes:

[0088] Traverse each voxel in the three-dimensional aquifer structure, and within a preset neighborhood range, obtain the lithological label sequence and thickness value sequence of each adjacent voxel of the current voxel.

[0089] The lithological variation gradient vector of the current voxel is obtained by measuring the difference between the lithological label of the current voxel and the lithological label sequence of the adjacent voxels.

[0090] The thickness value of the current voxel is compared with the thickness value sequence of the adjacent voxels by performing neighborhood difference to obtain the thickness change gradient vector of the current voxel.

[0091] The lithological variation gradient vector and the thickness variation gradient vector are bilinearly aggregated to obtain the local tensor synthesis result of the current voxel;

[0092] The local tensor synthesis results of all voxels in the three-dimensional aquifer structure are arranged and combined to obtain the three-dimensional spatial distribution trend information of the target region.

[0093] The spatial coordinates of all sampling points within the target area are determined. Based on these coordinates, the lithological description information and specific thickness values ​​of the aquifer corresponding to each sampling point are filtered one by one from the spatialized mapping dataset. The filtered information is checked one by one to ensure that each sampling point corresponds to unique lithological description data and thickness data. All sampling points and their corresponding two types of data are associated and integrated to form a discrete point attribute set of the target area.

[0094] The voxel side length of the target area is set to 1 meter. The spatial range of the interpolation filling is determined to be completely consistent with the target area. Based on the aquifer thickness data of all sampling points in the discrete point attribute set, spatial interpolation filling is performed on the blank areas between sampling points. During interpolation, the thickness value of each blank position is estimated by referring to the thickness values ​​of the three nearest sampling points to ensure that the thickness value after filling is continuous and conforms to the geological distribution law of the area. After filling, a continuous voxel model of aquifer thickness in the target area is constructed. Each voxel in this model contains a unique aquifer thickness value.

[0095] Using the aquifer lithology description data corresponding to each sampling point in the discrete point attribute set as lithology labels, the lithology label of each sampling point is assigned to the voxel where the sampling point is located in the aquifer thickness continuous voxel model. For voxels not directly covered by sampling points, the lithology label is determined based on the lithology labels of the three nearest sampling points around the voxel. If the lithology labels of the three sampling points are consistent, the label is directly assigned; otherwise, the lithology label of the nearest sampling point is used. After assigning the corresponding lithology labels to all voxels, the three-dimensional aquifer structure of the target area is obtained. This structure contains both the thickness value and the lithology label of each voxel.

[0096] The preset neighborhood range is set to a 3×3×3 voxel region around the current voxel. Each voxel in the three-dimensional aquifer structure is traversed in the order from shallow to deep, from left to right, and from top to bottom. For each current voxel, the lithology tags of all adjacent voxels in its neighborhood range are extracted one by one and arranged in the order of the neighboring voxels to form a lithology tag sequence. At the same time, the thickness values ​​of all adjacent voxels in the neighborhood range are extracted and arranged in the same order to form a thickness value sequence, ensuring that each current voxel can obtain a complete lithology tag sequence and thickness value sequence.

[0097] The lithological label of the current voxel is compared with each label in the lithological label sequence of the adjacent voxels one by one. If the descriptions are completely consistent, it is recorded as 0; if the descriptions are inconsistent, it is recorded as 1. All comparison results are arranged in the original order of the lithological label sequence to form an ordered set of numerical combinations. This set of numerical combinations is the lithological change gradient vector of the current voxel, which can intuitively reflect the changes in the lithology around the current voxel.

[0098] The thickness value of the current voxel is calculated by subtracting the thickness value of the corresponding adjacent voxel from the thickness value of the current voxel. All the calculated differences are arranged in the original order of the thickness value sequence to form an ordered combination of differences. This combination of differences is the thickness change gradient vector of the current voxel, which can intuitively reflect the change of thickness around the current voxel.

[0099] First, adjust the lithological variation gradient vector and the thickness variation gradient vector to be sequences of the same length. Take the values ​​at corresponding positions in the two vectors, and multiply the lithological difference value and the thickness difference value at each position to obtain the aggregated value at each position. Then, arrange the aggregated values ​​at all positions in the original order to form a new set of ordered numerical combinations. This numerical combination is the local tensor synthesis result of the current voxel, which can simultaneously reflect the comprehensive variation characteristics of the lithology and thickness around the current voxel.

[0100] Based on the spatial coordinates of each voxel in the three-dimensional aquifer structure, the local tensor synthesis results corresponding to each voxel are arranged and combined in the order of spatial position from shallow to deep, from left to right, and from top to bottom. The local tensor synthesis results of all voxels are integrated into a complete three-dimensional data set. This data set can clearly present the variation trend of aquifer lithology and thickness in three-dimensional space in the target area, which is the three-dimensional spatial distribution trend information of the target area.

[0101] Extract the aquifer conductivity attribute data corresponding to each sampling point from the spatial mapping dataset, clarify the correspondence between the spatial coordinates of each sampling point and the conductivity attribute data, then compare the spatial coordinates of each voxel in the three-dimensional spatial distribution trend information with the spatial coordinates of all sampling points, find the sampling point that is closest to the spatial position of each voxel, bind the conductivity attribute data of the sampling point with the local tensor synthesis result of the corresponding voxel, and record the spatial position, local tensor synthesis result and conductivity attribute data of all voxels one by one to obtain the correspondence table of the target area.

[0102] A region segmentation criterion is set: when the difference in the hydraulic conductivity attribute data of two adjacent voxels is less than 0.1 m / d, and the difference in the corresponding position values ​​in the local tensor synthesis results of the two voxels is less than 1, the two voxels are determined to belong to the same region; if this condition is not met, they are determined to be different regions. Based on this criterion, combined with the hydraulic conductivity attribute data of each voxel in the corresponding relationship table and the local tensor synthesis results, all voxels in the three-dimensional aquifer structure are divided into regions. Voxels that meet the same criterion are divided into one region, and each region is assigned a unique partition identifier. After summarizing the division results of all regions, the spatial difference partitioning information of the target region is obtained.

[0103] The beneficial effects are that by performing layered processing on spatialized mapping datasets, a three-dimensional aquifer structure is gradually constructed and three-dimensional spatial distribution trend information is extracted. Combined with aquifer water conductivity attribute data, regional segmentation is completed. This can accurately present the spatial distribution characteristics of lithology, thickness, and water conductivity of aquifers in the target area, eliminate correlation deviations between different data, ensure the accuracy and reproducibility of spatial difference zoning information, and clearly distinguish aquifer areas with different characteristics. This provides accurate and reliable technical support for subsequent aquifer development and utilization, water resource assessment, and geological disaster prevention and control. The entire process is clear and the steps are well-defined, which can effectively avoid data loss or bias problems and improve data utilization efficiency.

[0104] A3. Perform spatial overlay analysis on the spatial difference zoning information and the aquifer outcrop and structural point data in the original mapping dataset to obtain the boundary information set of the target area;

[0105] In this embodiment of the invention, the step of performing spatial overlay analysis on the spatial difference zoning information and the aquifer outcrop and structural point data in the original mapping dataset to obtain the boundary information set of the target area includes:

[0106] The spatial extent outline of each partition unit in the spatial difference partitioning information is used as the partition boundary outline set of the target region.

[0107] Locate the spatial coordinates of aquifer outcrops and the spatial coordinates and structural type attributes of structural points in the original mapping dataset, and compile them into a discrete geological feature point set for the target area;

[0108] The spatial overlay determination is performed between the partition boundary contour set and the discrete geological feature point set to obtain the overlay associated point set of the target region.

[0109] Based on the structural type attributes in the set of superimposed associated points, identify the partition boundary segments adjacent to the aquifer outcrop and specific structural points, and merge the partition boundary segments into the boundary fragment set of the target region;

[0110] The boundary fragment set is topologically stitched together to obtain the boundary information set of the target region.

[0111] The unique partition identifier of each partition unit in the spatial difference partition information is identified. The spatial range boundary of each partition unit is extracted one by one. During the extraction, the outline is drawn along the outermost voxel edge of each partition unit. During the outlining process, the intersection of the edges of each two adjacent outer voxels is used as the node of the outline. All nodes are connected in a clockwise order to form a complete closed outline of each partition unit. This ensures that each partition unit corresponds to a unique closed outline and that the outline can accurately cover the entire spatial range of the partition unit. The closed outlines of all partition units are summarized and integrated to form the partition boundary outline set of the target area. This set contains the complete spatial range outline of each partition unit and the corresponding partition identifier.

[0112] The criteria for determining aquifer outcrops are defined as points marked as aquifers directly exposed on the surface in the original mapping data. The criteria for determining structural points are points marked as faults, folds, and joints in the original mapping data. Based on these two criteria, all aquifer outcrops and structural points are selected one by one from the original mapping dataset. The spatial coordinates of each aquifer outcrop and the spatial coordinates and corresponding structural type attributes of each structural point are extracted. All extracted data are checked one by one to ensure that the coordinates of each aquifer outcrop are unique and that the coordinates and structural type attributes of each structural point correspond one-to-one, with no missing or mismatched data. The coordinates of all aquifer outcrops, the coordinates of all structural points, and the structural type attributes are all collected and organized to form a discrete geological feature point set for the target area.

[0113] A spatial overlay determination threshold of 0.5 meters is set. This threshold is set in conjunction with the mapping accuracy requirements of the target area and is used to determine the spatial relationship between contour lines and geological feature points. The coordinates of all contour lines in the partition boundary contour set are compared with the spatial coordinates of all points in the discrete geological feature point set. The straight-line distance between each contour line node and each geological feature point is calculated. When the straight-line distance between a certain contour line node and a certain geological feature point is less than 0.5 meters, it is determined that the node and the geological feature point have achieved spatial overlay. The coordinates of the overlaid contour line node, the corresponding partition identifier, the coordinates of the corresponding geological feature point, and the structural type attribute (if it is a structural point) are bound and recorded. After all the binding records are summarized, an overlay association point set of the target area is formed, ensuring that each overlay association record contains complete information on contour nodes, partitions, and geological feature points.

[0114] The scope of specific structural points is defined as structural points with structural type attributes identified as faults and folds. First, the associated records superimposed on aquifer outcrops are selected from the set of superimposed associated points. Based on the coordinates of the contour line nodes in the records, the two adjacent nodes on the boundary contour line of the zone where the node is located are found. Connecting the node with the two adjacent nodes forms two zone boundary line segments, which are the zone boundary line segments adjacent to the aquifer outcrops. Then, the associated records superimposed on specific structural points are selected from the set of superimposed associated points. The adjacent nodes on the contour line of the corresponding node are found using the same method, forming the zone boundary line segments adjacent to the specific structural points. All the selected zone boundary line segments are sorted, and duplicate line segments are removed. The criterion for determining duplicate line segments is that the coordinates of the two endpoints of the two line segments are completely consistent. All the remaining zone boundary line segments are classified and merged according to the zone identifier to form the boundary segment set of the target area. This set contains all zone boundary line segments adjacent to aquifer outcrops and specific structural points, as well as the corresponding zone identifiers.

[0115] The topology splicing criteria are set as follows: if the distance between the endpoints of two boundary line segments is less than 0.5 meters and they belong to the same partition identifier, the two line segments are considered to be splicable. The line segments in the boundary segment set are grouped according to the partition identifier. Within each group, the endpoint coordinates of the line segments are compared one by one to find the line segments that meet the splicing criteria. The splicable line segments are then connected into a complete line segment after the endpoints are aligned. During splicing, it is ensured that the line segments have the same direction, no intersections or overlaps, and no breaks. The spliced ​​line segments are checked. If there are any isolated line segments that are not spliced, their endpoints are checked again to see if they meet the splicing criteria. This process continues until all splicable line segments in each group are spliced. All spliced ​​complete line segments, unspliced ​​isolated line segments, corresponding partition identifiers, and associated geological feature point information are integrated to form the boundary information set of the target area. This set contains all complete line segments of the partition boundaries related to geological feature points within the target area and their related attribute information.

[0116] The beneficial effects are that by extracting the boundary contours of the zones and aggregating discrete geological feature points, combined with clear spatial overlay judgment criteria and topological splicing rules, the zone boundaries are accurately associated with aquifer outcrops and specific structural points, effectively eliminating boundary splicing deviations and data redundancy, ensuring the accuracy and integrity of boundary fragments. The resulting boundary information set can clearly present the spatial relationship between the zone boundaries of the target area and geological structures and aquifer outcrops, improving the accuracy and reproducibility of boundary information. This provides accurate and reliable boundary data support for subsequent work such as aquifer boundary delineation, geological structure analysis, rational development and utilization of water resources, and geological disaster prevention and control in the target area. The entire process is clear and the steps are well-defined, effectively avoiding boundary omissions or incorrect splicing, and ensuring the scientific nature and reliability of subsequent related work.

[0117] A4. Construct a digital spatial structure of the target area based on the boundary information set and the spatialized mapping dataset;

[0118] In this embodiment of the invention, constructing a digital spatial structure of the target area based on the boundary information set and the spatial mapping dataset includes:

[0119] The spatial coordinate sequence of the boundary lines in the boundary information set is spatially enclosed to obtain the closed spatial constraint domain of the target region.

[0120] The aquifer attribute data in the spatialized mapping dataset is mapped to discrete grid nodes within the closed spatial constraint domain to construct the initial field for attribute assignment of the target region.

[0121] Based on the boundary line, the grid cells in the initial field of the attribute assignment are subjected to boundary adaptation and clipping to obtain the boundary-fitting grid of the target region;

[0122] The boundary-fitting mesh is assembled into a digital spatial structure for the target region.

[0123] Extract the spatial coordinate sequence of each boundary line from the boundary information set. Use this coordinate sequence as the core basis for the closed loop. Connect all coordinate nodes one by one in a clockwise direction. During the connection process, ensure that the line segments between nodes do not intersect and that the first and last nodes are precisely closed. Set a complete spatial range in the internal region of the closed loop. This range serves as the constraint boundary and can completely surround all relevant aquifer voxels and geological feature points in the target area. Integrate all spatial voxels that meet the constraint conditions to obtain the closed spatial constraint domain of the target area. The boundary of this constraint domain completely coincides with the boundary line in the boundary information set.

[0124] The layout of discrete grid nodes within the closed spatial constraint domain is determined, with a preset spacing of 2 meters. Grid nodes are uniformly distributed within the closed spatial constraint domain based on this spacing, ensuring that the spatial coordinates of each grid node are unique and ordered. Various aquifer attribute data, such as aquifer thickness values, lithological labels, and hydraulic conductivity data, contained in the spatial mapping dataset are matched with the corresponding grid nodes. The matching criterion is that the spatial coordinates of the grid node are completely consistent with the spatial coordinates of the sampling point. For grid nodes with sampling points, the corresponding aquifer attribute data is directly assigned. For grid nodes without sampling points, spatial interpolation is performed based on the aquifer attribute data of the three nearest surrounding grid nodes, and the calculation result is used as the aquifer attribute data for that node. Attribute assignment is completed for all grid nodes, constructing an initial attribute assignment field for the target area. This initial field completely covers the aquifer attribute information of all nodes within the closed spatial constraint domain.

[0125] Using the boundary line of the boundary information set as the clipping reference, the spatial positional relationship between each grid cell and the boundary line in the initial field of attribute assignment is detected one by one. When all or part of the area of ​​the grid cell is located outside the boundary line, it is determined that the grid cell needs to be adapted and clipped. During clipping, the intersection of the boundary line and the grid cell is used as the new boundary node. The grid cell is divided according to the direction of the boundary line. The grid part located inside the boundary line is retained as the effective grid cell. For the grid cell that is crossed by the boundary line, the boundary fitting shape is reconstructed according to its relative position with the boundary line. The same clipping operation is performed on all grid cells that need to be clipped. Finally, the boundary fitting grid body of the target area is obtained with a shape that completely fits the boundary of the target area and has no grid cells that exceed the boundary range.

[0126] All valid grid cells within the boundary-fitted grid are assembled in an orderly manner according to their respective spatial coordinates. During assembly, it is ensured that the edges of adjacent grid cells completely overlap and the nodes are fully connected, without any gaps or overlaps. During the assembly process, each grid cell is assigned a unique structural identifier, and its corresponding aquifer attribute data and spatial morphology information are bound together. After completing the assembly operation of all grid cells, a digital entity with complete spatial morphology and detailed attribute information is formed. This entity can accurately reflect the spatial distribution and geological characteristics of the aquifer in the target area, which is the digital spatial structure of the target area.

[0127] The beneficial effects include: forming a closed spatial constraint domain by spatial enclosure, which can accurately define the research scope of the target area and avoid interference from irrelevant data; constructing an initial field for attribute assignment by mapping aquifer attribute data to discrete grid nodes, which can achieve the ordering and continuity of attribute data; obtaining a boundary-fitting grid body through boundary adaptation and trimming, which can ensure that the grid structure completely matches the actual boundary and eliminate errors caused by the grid exceeding the boundary; and forming a digital spatial structure by assembly, which can realize the digital representation and accurate modeling of the aquifer in the target area. The entire process is logically clear, the steps are well-defined, and it is highly reproducible. It can provide high-precision and complete digital basic data support for subsequent aquifer numerical simulation, water resource assessment, and geological disaster early warning, effectively improving the accuracy and reliability of related analysis and applications.

[0128] A5. Perform parameter inversion between the multi-source real-time monitoring data of the target area and the digital spatial structure to obtain the updated spatial structure of the target area;

[0129] In this embodiment of the invention, the step of performing parameter inversion between the multi-source real-time monitoring data of the target area and the digital spatial structure to obtain the updated spatial structure of the target area includes:

[0130] The dynamic water level change sequence and water quality response characteristics in the multi-source real-time monitoring data of the target area are compared with the corresponding initial hydrogeological parameters in the digital spatial structure to obtain the parameter sensitivity index set of the target area.

[0131] Based on the parameter sensitivity index set, the weights of the non-homogeneous parameter fields in the digital spatial structure are adjusted.

[0132] Based on the adjusted heterogeneous parameter field, the boundary constraints of the digital spatial structure are dynamically modified to obtain the intermediate transitional spatial structure of the target region;

[0133] The intermediate transition space structure is logically consistent with the historical water inrush event records in the multi-source real-time monitoring data. If the verification is successful, the intermediate transition space structure is determined as the updated space structure of the target area.

[0134] The specific scope of multi-source real-time monitoring data is clearly defined, including the water level dynamic change sequence recorded hourly for 30 consecutive days from all monitoring points within the target area, as well as water quality response characteristic data such as pH, turbidity, and dissolved oxygen at corresponding monitoring points within the same monitoring period. Initial hydrogeological parameters, including initial water level benchmark values ​​and initial water quality benchmark values, are extracted from the digital spatial structure for each monitoring point. The water level value for each time period in the water level dynamic change sequence is compared with the corresponding initial water level benchmark value, and the water level difference for each time period is calculated. The values ​​of each water quality index in the water quality response characteristics are then analyzed. Each water quality indicator is compared with the corresponding initial water quality benchmark value. The numerical difference of each indicator is calculated, and a difference judgment standard is set. When the water level difference is greater than 0.1 meters or the water quality indicator difference exceeds the preset benchmark range (i.e., pH value exceeds 6.5-8.5, turbidity exceeds 0-5 NTU, and dissolved oxygen exceeds 5-10 mg / L), the difference, the corresponding monitoring point location, and the parameter type are recorded as a sensitivity indicator. After summarizing all records that meet the conditions, a parameter sensitivity indicator set for the target area is obtained. It is ensured that each indicator in this set includes the difference, monitoring point location, parameter type, and corresponding comparison standard.

[0135] The specific types of heterogeneous parameter fields in the digital spatial structure are identified, including aquifer hydraulic conductivity parameter fields, lithological distribution parameter fields, and thickness variation parameter fields. The initial weight coefficient for each parameter field is set to 1. Based on the index types in the parameter sensitivity index set, the number of sensitivity indexes corresponding to each parameter field is counted one by one. The proportion of the number of sensitivity indexes corresponding to each parameter field to the total number of sensitivity indexes is calculated, and a weight adjustment rule is set. If the proportion of the number of sensitivity indexes corresponding to a certain parameter field is greater than 30%, the weight coefficient of that parameter field is adjusted to 1.5. If the proportion is less than or equal to 30%, the weight coefficient remains unchanged at 1. The weight adjustment of all heterogeneous parameter fields is completed according to this rule. After adjustment, the weight coefficient of each parameter field is checked one by one to ensure that the weight adjustment result matches the sensitivity of the parameter sensitivity index set, thus completing the weight adjustment of the heterogeneous parameter fields in the digital spatial structure.

[0136] Based on the spatial coordinate sequence of the boundary lines in the boundary information set, and combined with the adjusted heterogeneous parameter field, the boundary constraints of the digital spatial structure are dynamically corrected. The boundary constraints include the spatial coordinate accuracy constraints of the boundary lines and the attribute constraints of the grid cells inside the boundary. During the correction, the adaptability of the current boundary constraints of the digital spatial structure with the adjusted heterogeneous parameter field is checked one by one. When the parameter field weight corresponding to a certain boundary line is adjusted, if the spatial coordinates of the boundary line deviate from the grid node coordinates corresponding to the adjusted parameter field by more than 0.3 meters, the coordinate nodes of the boundary line are fine-tuned to a position where the deviation from the grid node coordinates is less than 0.3 meters. At the same time, the attribute constraint rules of the grid cells inside the boundary are adjusted to ensure that the attribute assignments of the grid cells are consistent with the weights of the adjusted heterogeneous parameter field. All boundary constraints are checked and corrected one by one. After the correction is completed, the intermediate transitional spatial structure of the target area is obtained. The boundary constraints of this structure are fully adapted to the adjusted heterogeneous parameter field, and the core attributes and spatial morphology of the digital spatial structure are preserved.

[0137] The specific content of historical water inrush event records in multi-source real-time monitoring data is clarified, including the spatial coordinates of the water inrush event, the water inrush period, the water inrush flow rate, and the corresponding changes in hydrogeological parameters. The grid cell attributes and hydrogeological parameters of the corresponding water inrush event location in the intermediate transition space structure are logically verified to be consistent with the relevant data in the historical water inrush event records. The verification criteria are that the initial water level at this location in the intermediate transition space structure deviates from the water level before the historical water inrush event by less than 0.1 meters, the water conductivity parameters deviate from the historical records by less than 0.05 m / d, and the parameter change trend at this location in the structure is completely consistent with the parameter change trend in the historical water inrush event. The adaptability of each historical water inrush event record to the intermediate transition space structure is verified one by one. If all historical water inrush event records meet the verification criteria, the verification is deemed successful, and the intermediate transition space structure is determined as the updated spatial structure of the target area. If there are records that do not meet the verification criteria, the non-homogeneous parameter field and boundary constraints are readjusted until all historical water inrush event records meet the verification criteria.

[0138] The beneficial effects are as follows: by accurately comparing multi-source real-time monitoring data with the initial hydrogeological parameters of the digital spatial structure, a set of parameter sensitivity indicators is generated, enabling targeted weight adjustments for heterogeneous parameter fields. Combined with dynamic boundary line correction, an intermediate transitional spatial structure is formed. Logical consistency is then verified through historical water inrush event records. This effectively improves the fit between the digital spatial structure and actual hydrogeological conditions, eliminates deviations between initial parameters and real-time monitoring data, and ensures the accuracy and reliability of the updated spatial structure. It can accurately reflect the dynamic changes of aquifers in the target area, providing high-precision and dynamic digital support for subsequent aquifer numerical simulation, water inrush risk prediction, rational water resource regulation, and geological disaster prevention. The entire process is clearly defined and reproducible, effectively avoiding parameter deviations and insufficient structural adaptability, and ensuring the scientific validity and effectiveness of subsequent related work.

[0139] A6. Interpret the comparison results between water level observation data and theoretical water level distribution information from the updated spatial structure to generate the spatial characteristics of water level anomalies in the target area, and extract the risk area location information of the target area by combining the water-conducting performance zoning information in the updated spatial structure.

[0140] In this embodiment of the invention, the step of interpreting the comparison results between water level observation data and theoretical water level distribution information from the updated spatial structure to generate spatial characteristics of water level anomalies in the target area, and combining this with the water-conducting performance zoning information in the updated spatial structure to extract risk area location information of the target area, includes:

[0141] Read water level observation data and corresponding theoretical water level distribution information from the updated spatial structure, and perform grid-by-grid difference calculation on the water level observation data and corresponding theoretical water level distribution information to obtain the water level residual dataset of the target area.

[0142] Based on the sign and magnitude of the residual values ​​at each grid point in the water level residual dataset, the water level residual dataset is divided into a positive residual subset and a negative residual subset.

[0143] Spatial clustering analysis is performed on the positive residual subset and the negative residual subset respectively, and the analyzed positive residual spatial blocks and negative residual spatial blocks are jointly labeled as the spatial features of water level anomalies in the target area;

[0144] Extract water-conducting performance zoning information from the updated spatial structure. The water-conducting performance zoning information includes a high water-conducting performance region, a medium water-conducting performance region, and a low water-conducting performance region.

[0145] The positive residual spatial blocks in the spatial features of the water level anomaly are spatially overlaid with the high water conductivity regions in the water conductivity zoning information, and the overlapping regions are identified as the risk area location information of the target area.

[0146] The formula for calculating the water level residual value of each grid point in the water level residual dataset is as follows:

[0147] ;

[0148] In the formula, For the first Water level residual value at each grid point For the first Water level observation data at each grid point For the first Theoretical water level distribution information for each grid point The preset spatial smoothing coefficient, For the first The set of neighboring grid points of a grid point The number of grid points in the neighboring grid point set. For the first Water level observation data from neighboring grid points, For the first Theoretical water level distribution information for each neighboring grid point.

[0149] For the first The water level observation data for each grid point comes from the actual water level values ​​collected by water level monitoring equipment deployed in the target area at preset monitoring times. The deployment locations of the monitoring equipment correspond one-to-one with the grid points in the target area, with each grid point corresponding to an independent water level monitoring device. The water level data collected by the equipment undergoes data validity verification, and abnormal values ​​caused by equipment failure or signal interference are removed before being directly used as the water level observation data for the corresponding grid point. This data is the basic input data for calculating the water level residual value. The collection process follows the standardized specifications for hydrogeological monitoring in the target area, and the data source is traceable and reproducible.

[0150] For the first The theoretical water level distribution information for each grid point is generated based on hydrogeological survey data and hydrological model calculations for the target area. The hydrogeological survey data includes fundamental parameters such as aquifer thickness, permeability coefficient, and recharge boundary conditions. These parameters are input into a pre-defined hydrological numerical model, which calculates the theoretical water level value for each grid point according to the grid division rules for the target area. This theoretical water level value is the theoretical water level distribution for the [number]th grid point. The theoretical water level distribution information for each grid point is generated based on the validated hydrogeological model of the target area, and the data results are highly consistent with the actual hydrogeological conditions of the target area.

[0151] The preset spatial smoothing coefficient is determined based on the hydrogeological spatial variation characteristics of the target area and the statistical analysis of historical water level monitoring data. First, grid-by-grid water level observation data and corresponding theoretical water level data of the target area over the past five years are collected. The spatial variation function of water level residuals at different spatial distances is calculated, and the spatial correlation range of water level residuals is obtained. Combined with the empirical value of the smoothing characteristics of aquifer water level spatial distribution in the field of coalfield hydrogeology, the value range of the spatial smoothing coefficient is determined to be 0.3 to 0.7. Finally, a fixed spatial smoothing coefficient value is selected according to the uniformity of aquifer distribution in the target area. This coefficient is used to adjust the influence of neighboring grid points on the residual of the current grid point. The value selection process is based on historical data statistics and field experience, and has clear basis and reproducibility.

[0152] For the first The set of neighboring grid points of the nth grid point is determined based on the grid spatial arrangement rules of the target region and the preset neighboring grid point judgment criteria. The target region adopts a regular rectangular grid point division method, with the nth grid point as the neighboring grid point. Taking the first grid point as the center, select the grid points in its four directions (up, down, left, and right) as neighboring grid points, and number and summarize these adjacent grid points to form the first grid point. The set of neighboring grid points of the i-th grid point, the size of which is related to the grid resolution of the target region, and all grid points in the set are neighboring grid points of the i-th grid point. Each grid point has a direct spatial adjacency relationship, and the process of determining the set follows a unified spatial location determination rule with no ambiguity in the determination conditions.

[0153] This represents the number of grid points in the neighboring grid point set. This value is directly determined by the number of grid point numbers in the neighboring grid point set; that is, by counting the total number of grid points in the neighboring grid point set. The counting process involves counting each grid point number in the neighboring grid point set one by one, and using the counted number as the total number. The specific value, The value is directly related to the rules for determining the set of neighboring grid points. The statistical process is an automated count, and the result is a clear integer with no fuzzy values.

[0154] For the first Water level observation data from the nearest grid point, the source of which is the first... The water level observation data for each grid point comes from the same source: the actual water level values ​​collected at the same preset monitoring time by water level monitoring equipment deployed at the corresponding grid points within the target area. This data is then compared with the data from the first grid point. The data validity verification process is the same for all grid point water level observation data. After removing outliers, it is used as the first... Water level observation data from neighboring grid points, the first The value of the water level observation data of the nearest grid point is related to the first The water level observation data at each grid point belong to the same monitoring batch and are consistent in both time and space.

[0155] For the first The theoretical water level distribution information of the nth neighboring grid points, the method of which generates this information is the same as that of the nth grid point. The theoretical water level distribution information for each grid point is generated in a completely consistent manner. Based on the hydrogeological survey data of the target area and the same hydrological numerical model, the theoretical water level distribution information for each grid point is calculated. The theoretical water level value of the nth neighboring grid point, the nth The theoretical water level distribution information of the nearest grid points and the first The theoretical water level distribution information of each grid point is theoretical data generated by the same hydrological model under the grid point division rules of the target area, and has a unified generation basis and data consistency.

[0156] For the first The water level residual value at the nth grid point, this value is derived from the nth grid point. The water level observation data at each grid point and the first The difference in theoretical water level distribution information at each grid point is calculated by adding the product of the spatial smoothing coefficient and the average residual of neighboring grid points. The calculation process is as follows: first calculate the... The observed and theoretical water level difference for each grid point is calculated, then the average of the observed and theoretical water level differences for all grid points within the neighboring grid point set is calculated. This average is multiplied by the spatial smoothing coefficient and then compared with the observed and theoretical water level difference for the first grid point. The differences between each of the 1 and 2 grid points are added together to obtain the final result. The water level residual value at the nth grid point, the th The calculation process of the water level residual value of each grid point integrates the influence of the water level deviation of a single grid point and the water level deviation of spatially adjacent grid points, and is the core indicator reflecting the comprehensive deviation of the grid point water level in the target area.

[0157] The water level observation data corresponding to each grid point is read sequentially from the updated spatial structure according to the preset grid point numbering order. At the same time, the theoretical water level distribution information that completely corresponds to the grid point position is read. The water level observation data and theoretical water level distribution information corresponding to each position are matched one-to-one. The water level observation data is the minuend and the theoretical water level distribution information is the subtrahend. The numerical subtraction operation is performed at the same grid point position. The result of the subtraction is used as the water level residual value corresponding to the grid point. The water level residual values ​​calculated from all grid points are integrated according to the original grid point spatial arrangement to form the water level residual dataset of the target area.

[0158] The residual value corresponding to each grid point in the water level residual dataset is compared with the value of zero. Grid points with residual values ​​greater than zero are assigned to the positive residual subset, grid points with residual values ​​less than zero are assigned to the negative residual subset, and grid points with residual values ​​equal to zero are not included in any residual subset. This completes the partitioning operation of the water level residual dataset.

[0159] The planar coordinates and corresponding residual values ​​of all grid points in the positive residual subset are selected as the basis for analysis. The grid points are classified according to the spatial distance between them and the similarity of their residual values. Grid points that are spatially adjacent and whose residual values ​​are in the same range are aggregated into the same spatial unit to form a positive residual spatial block. The planar coordinates and corresponding residual values ​​of all grid points in the negative residual subset are selected as the basis for analysis. They are aggregated according to the same classification rules of spatial distance and numerical similarity to form a negative residual spatial block. The positive and negative residual spatial blocks are uniformly labeled. The labeling result is the spatial feature of water level anomaly in the target area.

[0160] The water-conducting performance-related partition data is extracted from the updated spatial structure according to the preset regional division boundaries. Based on the preset water-conducting performance judgment threshold, the partition data is divided into three types of regions: regions with water-conducting performance values ​​higher than the first preset threshold are marked as high water-conducting performance regions; regions with water-conducting performance values ​​between the first and second preset thresholds are marked as medium water-conducting performance regions; and regions with water-conducting performance values ​​lower than the second preset threshold are marked as low water-conducting performance regions. The above three types of regions together constitute the water-conducting performance partition information.

[0161] The spatial boundary range of the positive residual spatial block within the spatial characteristics of water level anomalies is compared with the spatial boundary range of the high water conductivity area within the water conductivity zoning information. Each grid point within the positive residual spatial block is determined to be within the spatial boundary of the high water conductivity area. If a grid point is simultaneously within both the positive residual spatial block and the high water conductivity area, it is determined to be an overlapping area. All the determined overlapping areas are integrated, and the integrated area information is the risk area location information of the target area.

[0162] This formula is directly related to the generation of water level residual datasets and the identification of spatial features of water level anomalies. The formula calculates the first... The grid-point water level residual values ​​serve as the core components of the water level residual dataset. Each grid-point water level residual value accurately reflects the comprehensive deviation between the observed and theoretical water levels at that grid point. By integrating all grid-point water level residual values, a complete water level residual dataset is formed. This dataset provides a quantitative basis for subsequently dividing the dataset into positive and negative residual subsets. The division into positive and negative residual subsets is based solely on the comparison between the water level residual value and zero. Spatial clustering analysis, on the other hand, divides the dataset into blocks based on the spatial distribution characteristics of the water level residual values. The formula provides accurate water level residual values, ensuring a quantitative basis for identifying spatial features of water level anomalies. Meanwhile, the identification process of risk area location information is based on the superposition of positive residual spatial blocks and high water conductivity areas. The division of positive residual spatial blocks depends on the positive or negative attributes of the water level residual values. Therefore, the formula provides core quantitative data support for the entire process of water level anomaly identification and risk area location. The calculation results of the formula directly determine the analysis objects and judgment criteria of all subsequent steps, ensuring the logical coherence and data accuracy of the entire technical solution.

[0163] The beneficial effects are that by calculating the water level residuals point by point and dividing them into positive and negative residual subsets, the distribution characteristics of water level deviations can be accurately located. Spatial clustering analysis can clearly identify the spatial distribution pattern of water level anomalies. Combined with the water-conducting performance zoning information, spatial overlay analysis can achieve accurate location of risk areas. The entire process is based on clear numerical comparison and spatial location determination rules, which can stably reproduce the analysis results, effectively improve the accuracy and reliability of hydrogeological risk identification, and avoid identification errors caused by fuzzy judgments. This provides clear and quantifiable risk basis for hydrological management in the target area.

[0164] A7. Perform trend extrapolation between the risk area location information and the water inflow record data in the updated spatial structure to construct a risk distribution information map of the target area.

[0165] In this embodiment of the invention, the step of performing trend extrapolation between the risk area location information and the water inflow record data in the updated spatial structure to construct a risk distribution information map of the target area includes:

[0166] Based on the spatial location of each risk block in the risk area location information, extract the water inflow record data at the corresponding location from the updated spatial structure to generate a risk-related water inflow time series sequence for the target area;

[0167] Trend fitting is performed on the time series of risk-related water inflow to extract the direction and intensity of the water inflow change trend for each risk block;

[0168] The direction and intensity of the water inflow change trend are spatially correlated with the risk area location information to construct a risk distribution information map of the target area.

[0169] Based on the spatial boundary coordinates of each risk block in the risk area location information, the spatial range of each risk block is determined. This spatial range is based on grid coordinates and accurately corresponds to the grid positions in the updated spatial structure. According to the spatial range of each risk block, all water inflow records at the corresponding grid positions are extracted one by one from the updated spatial structure. The water inflow records are collected at a preset monitoring time interval, which is fixed at once a day. During extraction, it is ensured that the water inflow records at the corresponding positions of each risk block are complete, without missing or abnormal data. All water inflow records corresponding to each risk block are arranged in chronological order of monitoring time, from the earliest monitoring time to the latest monitoring time. Each risk block corresponds to a set of water inflow data arranged by time. The time series data of water inflow corresponding to all risk blocks together form the risk-related water inflow time series sequence of the target area.

[0170] Trend fitting was performed on the time series data corresponding to each risk block in the risk-related water inflow time series. The fitting process involved sequentially matching the water inflow data of each risk block according to the monitoring time, with the monitoring time as the horizontal axis and the water inflow value as the vertical axis, and connecting all data points sequentially to form a continuous change curve. The direction of the water inflow change trend was determined by observing the overall trend of the curve. If the curve shows an overall upward trend from the earliest monitoring time to the latest monitoring time, it is considered an upward trend; if it shows an overall downward trend, it is considered a downward trend. If the absolute value of the cumulative difference in the amplitude of the curve's fluctuations is less than 20 cubic meters per day, it is considered a stable trend. The trend strength was determined by... By calculating the difference in water inflow between all two adjacent monitoring times, the cumulative difference is obtained by summing all the differences. A positive cumulative difference with an absolute value greater than or equal to 50 cubic meters per day indicates a strong upward trend; a positive cumulative difference with an absolute value between 20 and 50 cubic meters per day indicates a weak upward trend; a negative cumulative difference with an absolute value greater than or equal to 50 cubic meters per day indicates a strong downward trend; a negative cumulative difference with an absolute value between 20 and 50 cubic meters per day indicates a weak downward trend; and a cumulative difference with an absolute value less than 20 cubic meters per day indicates a stable trend. Based on this, the direction and intensity of the water inflow change trend in each risk area can be extracted.

[0171] The spatial coordinates of each risk block are correlated one-to-one with the corresponding trend direction and intensity of water inflow. First, the spatial boundary and center coordinates of each risk block in the risk area location information are determined. Then, the trend direction and intensity information of each risk block are bound to the center coordinates of that risk block. A unique identifier is set for different trend directions: red for an upward trend, blue for a downward trend, and green for a stable trend. Different shades of identifier are set for different trend intensities: dark for a strong trend, light for a weak trend, and standard for a stable trend. All risk blocks bound with trend information are accurately superimposed onto the base spatial map of the target area according to their spatial location, ensuring that the location of each risk block is completely consistent with the actual spatial coordinates. The resulting map, which includes the location of risk blocks, the trend direction of water inflow, and the trend intensity, is the risk distribution information map of the target area.

[0172] The beneficial effects are as follows: by extracting the time series data of water inflow at the corresponding location of the risk area, a precise correlation between the risk area and the water inflow change is achieved. The trend fitting process adopts clear judgment criteria, which can accurately extract the trend direction and intensity of water inflow change. The spatial correlation operation deeply integrates trend information with spatial location. The constructed risk distribution information map can intuitively present the water inflow change characteristics of each risk block in the target area. The entire process is standardized and the judgment criteria are clear. The results can be stably reproduced, which effectively makes up for the limitations of relying solely on water level anomalies to identify risks. It improves the comprehensiveness and accuracy of hydrogeological risk assessment and provides an intuitive and quantifiable reference for the formulation of risk prevention and control measures in the target area.

[0173] It will be apparent to those skilled in the art that the present invention is not limited to the details of the exemplary embodiments described above, and that the present invention can be implemented in other specific forms without departing from the spirit or essential characteristics of the present invention.

[0174] This application embodiment can acquire and process relevant data based on artificial intelligence technology. Artificial intelligence is the theory, method, technology, and application system that uses digital computers or machines controlled by digital computers to simulate, extend, and expand human intelligence, perceive the environment, acquire knowledge, and use that knowledge to obtain optimal results.

[0175] Finally, it should be noted that 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 preferred embodiments, those skilled in the art should understand that modifications or equivalent substitutions can be made to the technical solutions of the present invention without departing from the spirit and scope of the technical solutions of the present invention.

Claims

1. A method for analyzing the spatial distribution of coalfield aquifers based on surveying and mapping fusion, characterized in that, The method includes: A1. Spatial registration is performed on the original mapping dataset of the target area to obtain the spatialized mapping dataset of the target area; A2. Reconstruct the lithology and thickness attribute data of the aquifer in the spatialized mapping dataset in three dimensions, and associate and map the reconstructed three-dimensional spatial distribution trend information with the water-conducting performance attribute data of the aquifer in the spatialized mapping dataset to obtain the spatial difference zoning information of the target area. A3. Perform spatial overlay analysis on the spatial difference zoning information and the aquifer outcrop and structural point data in the original mapping dataset to obtain the boundary information set of the target area; A4. Construct a digital spatial structure of the target area based on the boundary information set and the spatialized mapping dataset; A5. Perform parameter inversion between the multi-source real-time monitoring data of the target area and the digital spatial structure to obtain the updated spatial structure of the target area; A6. Interpret the comparison results between water level observation data and theoretical water level distribution information from the updated spatial structure to generate the spatial characteristics of water level anomalies in the target area, and extract the risk area location information of the target area by combining the water-conducting performance zoning information in the updated spatial structure. A7. Perform trend extrapolation between the risk area location information and the water inflow record data in the updated spatial structure to construct a risk distribution information map of the target area.

2. The method for spatial distribution analysis of coalfield aquifers based on surveying and mapping fusion as described in claim 1, characterized in that, The spatial registration of the original mapping dataset of the target area to obtain the spatialized mapping dataset of the target area includes: Obtain the original mapping dataset of the target area and parse the original coordinate reference information of each data source in the original mapping dataset. The original coordinate reference information includes geodetic coordinate system parameters and map projection parameters. Based on the geographical extent of the target area, the projection parameters of the standard spatial reference frame are extracted to generate a unified target coordinate reference system for the target area. Based on the unified target coordinate reference system, the original coordinate reference information is spatially unified to obtain a dataset in the transition coordinate system of the target region. Spatial feature points in the dataset under the transition coordinate system are spatially matched with the standard geographic base map of the target area to obtain a spatialized mapping dataset of the target area.

3. The method for spatial distribution analysis of coalfield aquifers based on surveying and mapping fusion as described in claim 1, characterized in that, The step involves performing three-dimensional reconstruction of the aquifer lithology and thickness attribute data in the spatialized mapping dataset, and then mapping the reconstructed three-dimensional spatial distribution trend information with the aquifer hydraulic conductivity attribute data in the spatialized mapping dataset to obtain the spatial difference zoning information of the target area, including: Extract the aquifer lithology description data and aquifer thickness numerical data corresponding to the sampling points from the spatialized mapping dataset to generate a discrete point attribute set for the target area; Spatial interpolation is performed on the sampling points in the discrete point attribute set to construct a continuous voxel model of the aquifer thickness in the target region. Using the aquifer lithology description data as labels, values ​​are assigned to the aquifer thickness continuum voxel model on a voxel-by-voxel basis to obtain the three-dimensional aquifer structure of the target region. The lithological variation gradient and thickness variation gradient at each voxel location are analyzed from the three-dimensional aquifer structure, and the lithological variation gradient and thickness variation gradient are tensor synthesized to obtain the three-dimensional spatial distribution trend information of the target area. The three-dimensional spatial distribution trend information is correlated point by point with the aquifer water conductivity attribute data extracted from the spatialized mapping dataset to obtain the correspondence table of the target area; Based on the correspondence table, the three-dimensional aquifer structure is divided into regions to obtain the spatial difference partitioning information of the target region.

4. The method for spatial distribution analysis of coalfield aquifers based on surveying and mapping fusion as described in claim 3, characterized in that, The process of resolving the lithological variation gradient and thickness variation gradient at each voxel location from the three-dimensional aquifer structure, and then performing tensor synthesis of the lithological variation gradient and thickness variation gradient to obtain the three-dimensional spatial distribution trend information of the target region includes: Traverse each voxel in the three-dimensional aquifer structure, and within a preset neighborhood range, obtain the lithological label sequence and thickness value sequence of each adjacent voxel of the current voxel. The lithological variation gradient vector of the current voxel is obtained by measuring the difference between the lithological label of the current voxel and the lithological label sequence of the adjacent voxels. The thickness value of the current voxel is compared with the thickness value sequence of the adjacent voxels by performing neighborhood difference to obtain the thickness change gradient vector of the current voxel. The lithological variation gradient vector and the thickness variation gradient vector are bilinearly aggregated to obtain the local tensor synthesis result of the current voxel; The local tensor synthesis results of all voxels in the three-dimensional aquifer structure are arranged and combined to obtain the three-dimensional spatial distribution trend information of the target region.

5. The method for spatial distribution analysis of coalfield aquifers based on surveying and mapping fusion as described in claim 1, characterized in that, The step of spatially overlaying the spatially differentiated zoning information with the aquifer outcrop and structural point data in the original mapping dataset to obtain the boundary information set of the target area includes: The spatial extent outline of each partition unit in the spatial difference partitioning information is used as the partition boundary outline set of the target region. Locate the spatial coordinates of aquifer outcrops and the spatial coordinates and structural type attributes of structural points in the original mapping dataset, and compile them into a discrete geological feature point set for the target area; The spatial overlay determination is performed between the partition boundary contour set and the discrete geological feature point set to obtain the overlay associated point set of the target region. Based on the structural type attributes in the set of superimposed associated points, identify the partition boundary segments adjacent to the aquifer outcrop and specific structural points, and merge the partition boundary segments into the boundary fragment set of the target region; The boundary fragment set is topologically stitched together to obtain the boundary information set of the target region.

6. The method for spatial distribution analysis of coalfield aquifers based on surveying and mapping fusion as described in claim 1, characterized in that, The step of constructing a digital spatial structure of the target area based on the boundary information set and the spatial mapping dataset includes: The spatial coordinate sequence of the boundary lines in the boundary information set is spatially enclosed to obtain the closed spatial constraint domain of the target region. The aquifer attribute data in the spatialized mapping dataset is mapped to discrete grid nodes within the closed spatial constraint domain to construct the initial field for attribute assignment of the target region. Based on the boundary line, the grid cells in the initial field of the attribute assignment are subjected to boundary adaptation and clipping to obtain the boundary-fitting grid of the target region; The boundary-fitting mesh is assembled into a digital spatial structure for the target region.

7. The method for spatial distribution analysis of coalfield aquifers based on surveying and mapping fusion as described in claim 1, characterized in that, The step of performing parameter inversion between the multi-source real-time monitoring data of the target area and the digital spatial structure to obtain the updated spatial structure of the target area includes: The dynamic water level change sequence and water quality response characteristics in the multi-source real-time monitoring data of the target area are compared with the corresponding initial hydrogeological parameters in the digital spatial structure to obtain the parameter sensitivity index set of the target area. Based on the parameter sensitivity index set, the weights of the non-homogeneous parameter fields in the digital spatial structure are adjusted. Based on the adjusted heterogeneous parameter field, the boundary constraints of the digital spatial structure are dynamically modified to obtain the intermediate transitional spatial structure of the target region; The intermediate transition space structure is logically consistent with the historical water inrush event records in the multi-source real-time monitoring data. If the verification is successful, the intermediate transition space structure is determined as the updated space structure of the target area.

8. The method for spatial distribution analysis of coalfield aquifers based on surveying and mapping fusion as described in claim 1, characterized in that, The comparison results between the water level observation data and the theoretical water level distribution information interpreted from the updated spatial structure are used to generate the spatial characteristics of water level anomalies in the target area. Combined with the water-conducting performance zoning information in the updated spatial structure, the risk area location information of the target area is extracted, including: Read water level observation data and corresponding theoretical water level distribution information from the updated spatial structure, and perform grid-by-grid difference calculation on the water level observation data and corresponding theoretical water level distribution information to obtain the water level residual dataset of the target area. Based on the sign and magnitude of the residual values ​​at each grid point in the water level residual dataset, the water level residual dataset is divided into a positive residual subset and a negative residual subset. Spatial clustering analysis is performed on the positive residual subset and the negative residual subset respectively, and the analyzed positive residual spatial blocks and negative residual spatial blocks are jointly labeled as the spatial features of water level anomalies in the target area; Extract water-conducting performance zoning information from the updated spatial structure. The water-conducting performance zoning information includes a high water-conducting performance region, a medium water-conducting performance region, and a low water-conducting performance region. The positive residual spatial blocks in the spatial features of the water level anomaly are spatially overlaid with the high water conductivity regions in the water conductivity zoning information, and the overlapping regions are identified as the risk area location information of the target area.

9. The method for spatial distribution analysis of coalfield aquifers based on surveying and mapping fusion as described in claim 8, characterized in that, The formula for calculating the water level residual value of each grid point in the water level residual dataset is as follows: ； In the formula, For the first Water level residual value at each grid point For the first Water level observation data at each grid point For the first Theoretical water level distribution information for each grid point The preset spatial smoothing coefficient, For the first The set of neighboring grid points of a grid point The number of grid points in the neighboring grid point set. For the first Water level observation data from neighboring grid points, For the first Theoretical water level distribution information for each neighboring grid point.

10. The method for spatial distribution analysis of coalfield aquifers based on surveying and mapping fusion as described in claim 1, characterized in that, The step of performing trend extrapolation between the risk area location information and the updated water inflow record data in the spatial structure to construct a risk distribution information map of the target area includes: Based on the spatial location of each risk block in the risk area location information, extract the water inflow record data at the corresponding location from the updated spatial structure to generate a risk-related water inflow time series sequence for the target area; Trend fitting is performed on the time series of risk-related water inflow to extract the direction and intensity of the water inflow change trend for each risk block; The direction and intensity of the water inflow change trend are spatially correlated with the risk area location information to construct a risk distribution information map of the target area.

Images