A coal mine subsidence surveying and mapping system

Abstract

The present invention relates to the technical field of subsidence surveying and mapping, and specifically to a coal mine subsidence surveying and mapping system, including: a data acquisition module, a subsidence boundary recognition module, a subsidence evolution module, a subsidence coupling module, and a subsidence evaluation module; associating the coal mine subsidence images pair with the subsidence area according to the horizontal displacement of the coal mine subsidence images at multiple time windows; extracting the subsidence markers of the coal mine subsidence images pair, classifying the subsidence area with the subsidence markers, and outputting the subsidence vector boundary; dividing grid cells according to the subsidence vector boundary, identifying the deformation points and the deformation components in multiple directions within each grid cell, and forming a continuous deformation surface; outputting the subsidence evolution trend through the integral probability and the associated position of each deformation point in the continuous deformation surface; coupling the subsidence evolution trend to identify the regional deformation characteristics of the subsidence evolution trend; setting the subsidence surveying and mapping levels of each subsidence area according to the regional deformation characteristics; improving the accuracy and efficiency of subsidence surveying and mapping.

Description

Technical Field

[0001] The present invention relates to the technical field of subsidence surveying and mapping, and specifically to a coal mine subsidence surveying and mapping system. Background Art

[0002] Quickly, accurately, and comprehensively monitoring the subsidence of the mining area is of great significance for reasonably guiding coal resource mining activities, preventing mining subsidence disasters, and protecting the ecological environment from being damaged. Traditional surface subsidence detection technologies are prone to ignoring the deformation boundaries and the relative trend boundaries of each deformation rate in the subsiding areas, resulting in conflicts between the positions of each boundary line and the data of its settlement changing with time when dividing the subsidence boundary line, reducing the regional accuracy of subsidence surveying and mapping.

[0003] For example, Chinese Patent Publication No. CN112577470A discloses a method and system for fusing UAV and InSAR to monitor the dynamic subsidence basin of a mining area, belonging to the technical field of mining area subsidence monitoring. The method includes the following steps: respectively obtaining the UAV subsidence basin image and the InSAR subsidence basin image of the monitored area by using UAV and InSAR; obtaining the fusion boundary line of the UAV subsidence basin image and the InSAR subsidence basin image; wherein comparing the subsidence contour lines of the UAV subsidence basin image and the InSAR subsidence basin image pairwise, and taking the two closest subsidence contour lines as the fusion boundary line; based on the fusion boundary line, fusing the UAV subsidence basin image and the InSAR subsidence basin image to obtain a complete high-precision subsidence basin image.

[0004] For example, Chinese Patent Publication No. CN112629485A discloses a method for monitoring the surface subsidence of a mine, including the following steps: arranging a surface movement observation station, before the underground project excavation, according to the excavation size and excavation depth, arranging a cross-shaped monitoring line on the surface, and the intersection points of the cross-shaped monitoring line are grid monitoring points; setting permanent monitoring points according to the position of the working face cutting eye; preliminarily determining the key monitoring range on the surface according to the underground project excavation parameters; installing wireless signal receiving devices at each monitoring point; determining the monitoring frequency scheme.

[0005] The prior art illustrates that the boundary lines of two forms of images can be combined to identify the subsiding area, and the monitoring points are set by the excavation size and excavation depth. However, when the prior art combines and describes multiple images, it is difficult to effectively associate the deformations and display contents at multiple time points only through the form of boundary lines, resulting in the inability to track the change trend of the subsiding area. At the same time, when identifying its subsidence boundary, it is necessary to set the boundary according to its changing form to quickly identify the risk level of the subsiding area and make decisions on the subsidence in a timely manner, improving the ability to warn of subsidence. Summary of the Invention

[0006] To solve the above technical problems, the technical solution adopted by the present invention is as follows: A coal mine subsidence mapping system, comprising: a data acquisition module, which is used to obtain pairs of subsidence images of the mining area, and associate the pairs of subsidence images of the mining area with the subsidence area according to the horizontal displacements of the subsidence images of the mining area in multiple time windows.

[0007] A subsidence boundary recognition module, which is used to extract subsidence marks in the pairs of subsidence images of the mining area, classify the subsidence area according to the temporal changes of the subsidence marks, and output a subsidence vector boundary.

[0008] A subsidence evolution module, which is used to divide grid cells according to the subsidence vector boundary, identify deformation points and deformation components in multiple directions within each grid cell to form a continuous deformation surface; obtain a subsidence evolution path through the integral probability and associated positions of the deformation points within the continuous deformation surface, and output a subsidence evolution trend based on the deformation rate on the subsidence evolution path.

[0009] A subsidence coupling module, which is used to couple the subsidence evolution trend according to the positions of the grid cells, and identify the regional deformation characteristics of the subsidence evolution trend in the corresponding time window based on the deformation gradient of adjacent grid cells.

[0010] A subsidence evaluation module, which is used to integrate the regional deformation characteristics, and set the subsidence mapping levels of each subsidence area according to the subsidence expansion rate and direction of each grid cell.

[0011] The beneficial effects of the present invention are as follows: First, by obtaining pairs of subsidence images of the mining area, marking the ground objects in the mining area, identifying the horizontal displacement according to the deformation rate of the existing deformation in multiple time windows, and using the center points of multiple points as the horizontal displacement, the number of processing points for the subsidence area is reduced, and the situation of data distortion easily occurring during the gradual matching of multiple points is avoided; after marking multiple data in the pairs of subsidence images of the mining area, classifying the subsidence area according to the data content corresponding to the subsidence marks and horizontal displacements, and describing different forms of subsidence vector boundaries, which is convenient for subsequent staff to perform data backtracking and subsidence prediction, and prevent the occurrence of construction risks in the mining area. And the accuracy of identifying different subsidence boundaries is improved.

[0012] Second, by dividing grid cells according to the maximum subsidence direction and the minimum subsidence direction, and extracting the subsidence evolution path according to the integral probability and associated positions in the grid points, the obtained subsidence evolution path can include the subtle trends in the horizontal and vertical directions, improving the accuracy of the changes of each deformation point during the subsidence evolution process, providing a data basis for subsequent risk assessment, and the subsidence evolution path can show the position of the subsidence evolution in the corresponding area, which is convenient for observation and monitoring at multiple positions.

[0013] III. After identifying the distribution of the subsidence evolution trend in each grid cell, the present invention identifies the corresponding nodes of the deformation rate and the subsidence evolution path for its adjacent grid cells to obtain the subsidence patterns corresponding to each region. Then, the differences in multiple time windows are obtained in the form of a symmetric difference set, so as to quickly evaluate the risk level of the subsidence area, help managers adjust the mining strategy in a timely manner, and reduce potential risks. BRIEF DESCRIPTION OF THE DRAWINGS

[0014] The present invention will be further described below with reference to the drawings and embodiments.

[0015] Figure 1 is a system framework diagram of a coal mine subsidence surveying and mapping system.

[0016] Figure 2 is a schematic flow diagram of a data acquisition module of a coal mine subsidence surveying and mapping system.

[0017] Figure 3 is a schematic flow diagram of a subsidence boundary identification module of a coal mine subsidence surveying and mapping system.

[0018] Figure 4 is a schematic flow diagram of a subsidence evolution module of a coal mine subsidence surveying and mapping system.

[0019] Figure 5 is a schematic flow diagram of a subsidence coupling module of a coal mine subsidence surveying and mapping system. DETAILED DESCRIPTION OF THE EMBODIMENTS

[0020] The embodiments of the present invention will be described in detail below. The described embodiments are exemplary and are only used to explain the present invention and should not be construed as a limitation to the present invention. For those not specified in the embodiments regarding specific technologies or conditions, they shall be carried out according to the technologies or conditions described in the literature in this field or according to the product specifications.

[0021] Refer to Figure 1 , a coal mine subsidence surveying and mapping system, comprising: a data acquisition module, a subsidence boundary identification module, a subsidence evolution module, a subsidence coupling module, and a subsidence evaluation module; wherein, the output end of the data acquisition module is connected to the subsidence boundary identification module, the output end of the subsidence boundary identification module is connected to the subsidence evolution module, the output end of the subsidence evolution module is connected to the subsidence coupling module, and the output end of the subsidence coupling module is connected to the subsidence evaluation module.

[0022] The data acquisition module is used to obtain pairs of coal mine subsidence images and associate the pairs of coal mine subsidence images with the subsidence area according to the horizontal displacement of the coal mine subsidence images in multiple time windows.

[0023] The subsidence boundary recognition module is used to extract subsidence markers in the subsidence image pairs of the mining area, classify the subsidence areas based on the temporal changes of the subsidence markers, and output the subsidence vector boundary.

[0024] The subsidence evolution module is used to divide grid cells according to the subsidence vector boundary, identify the deformation points and deformation components in multiple directions within each grid cell to form a continuous deformation surface; obtain the subsidence evolution path through the integral probability and associated positions of the deformation points within the continuous deformation surface, and output the subsidence evolution trend based on the deformation rate on the subsidence evolution path.

[0025] The subsidence coupling module is used to couple the subsidence evolution trend according to the positions of the grid cells, and identify the regional deformation characteristics of the subsidence evolution trend in the corresponding time window with the deformation gradient of adjacent grid cells.

[0026] The subsidence evaluation module is used to integrate the regional deformation characteristics, and set the subsidence surveying and mapping levels of each subsidence area according to the rate and direction of subsidence expansion of each grid cell.

[0027] Preferably, the subsidence image pairs of the mining area are represented as multiple images collected for the area where ground subsidence needs to be identified, including the collected InSAR data and spectral satellite data. When these two types of data are used, multiple images are taken based on one area, and the multiple images taken are combined into an image pair. The spectral satellite data is used to identify the composition of surface substances, soil vegetation, etc. in the area image, while the InSAR data identifies the deformation at the corresponding positions to illustrate the subsidence area and relative boundary of the area.

[0028] When obtaining the subsidence image pairs of the mining area, it is also necessary to preprocess the images, including but not limited to data registration. Use geocoding tools (such as GDAL, ENVI) to perform spatial registration on multi-source data to ensure that all data is in the same geographic coordinate system; for spectral satellite data, radiometric normalization processing is required to eliminate the influence of factors such as the atmosphere and solar altitude angle; for InSAR, a phase unwrapping algorithm is required to obtain the deformation amount that needs to be identified in the current subsidence area.

[0029] As Figure 2 shown, the implementation method of the data acquisition module includes obtaining the historical subsidence image pairs of the current subsidence area, comparing the image boundaries of the historical subsidence image pairs with the current subsidence image pairs to obtain the aligned subsidence image pairs.

[0030] At this time, based on the historical mining area subsidence images of the captured content, the mining area subsidence image pairs collected within multiple time periods are positionally calibrated to identify the subsided position and boundary of the current subsidence area; this step belongs to the form of image preprocessing and is only used here to combine multiple images as the content that can calculate and identify mining area subsidence; as for the deformation corresponding to InSAR data and the ground object shape corresponding to spectral satellite data contained in the mining area subsidence image pairs, the data in the two images are combined according to the geographical coordinates of the image acquisition.

[0031] According to the aligned mining area subsidence image pairs, calculate the horizontal displacement of the monitoring points within the subsidence area, and use the horizontal displacement of the monitoring points to associate the mining area subsidence image pairs with the subsidence area.

[0032] Preferably, the horizontal displacement of the monitoring points within the subsidence area represents the horizontal change in the position of the monitoring points marked in the mining area subsidence image pairs each time they are marked. For example, in the historical mining area subsidence image pairs, point A is marked on the image, and point A has a settlement of 5 mm. When the mining area subsidence image pairs are collected next time, the ground position corresponding to point A has a horizontal displacement. Then, the horizontal displacement that appears at this position is associated with the mining area subsidence image pairs collected next time to identify its horizontal displacement, and combined with the deformation field identified in the image, it explains the displacement conditions of the ground at the positions of the monitoring points in multiple directions.

[0033] When setting monitoring points here to identify horizontal displacement, the implementation methods for calculating the horizontal displacement of the monitoring points within the subsidence area include: using the orthophoto back-projection technology to obtain the pixel coordinates of each monitoring point in the mining area subsidence pair.

[0034] Match the monitoring points of adjacent time series, use the center point of the monitoring points as the target feature point, calculate the matching degree of the target feature point, and use the distance of the target feature point with the maximum matching degree in the adjacent time series as the horizontal displacement of the monitoring points within the subsidence area.

[0035] Among them, the calculation method of the matching degree is obtained by calculating the Pearson correlation coefficient according to the pixel distances of the multiple monitoring points corresponding to the target feature point in the adjacent time series.

[0036] The reason for choosing the center point here is to reduce the number of points used in the calculation of the horizontal displacement. At the same time, choosing the center point of the monitoring points is to select relatively stable points to calculate the horizontal displacement. Because the monitoring points are generally set with the ground objects existing on the ground as markers, when the ground objects change, the positions of the monitoring points are likely to change greatly, resulting in data distortion when using a single monitoring point for step-by-step matching. Using the center point is to prevent data distortion, and choosing the center point can avoid the calculation of a large number of redundant monitoring points, ensuring that the horizontal displacement on the ground can be identified in the presence of noise.

[0037] Preferably, the center point of the monitoring points is extracted according to the area where multiple monitoring points are arranged. For example, the monitoring points are arranged in a certain small area in the form of triangles, squares, or polygons, and the extracted center point is the center point of multiple monitoring points in this area.

[0038] Preferably, the above-obtained horizontal displacement is based on the joint solution of the ascending and descending orbit InSAR for the eastward and northward components, and is obtained by taking the square root of the sum of the squares of the eastward and northward deformation rates. It is a comprehensive index representing the magnitude of the surface horizontal movement speed calculated using the Pythagorean theorem in a two-dimensional plane.

[0039] In an embodiment of the present invention, the above-mentioned subsidence vector boundary is used to describe the gradient of the change of subsidence marks at the boundary. The extracted subsidence marks include bare soil, vegetation, cracks, and deformation marks, and the subsidence area is classified according to the proportion of the corresponding area change of the subsidence marks.

[0040] For example, after extracting the settlement amount from the mining area subsidence image pair, the cumulative settlement amount, the change ratio of bare soil, the vegetation index, and the crack length growth rate in a fixed time period are used as its subsidence marks, and it is divided into various types of areas such as active subsidence areas, ecological damage areas, fault zone expansion areas, and stable areas according to the subsidence marks.

[0041] The settlement amount here can be extracted according to the vertical deformation component of the InSAR data in the mining area subsidence image pair, and the corresponding contents of bare soil, vegetation, and cracks are extracted using spectral satellite data. At the same time, in order to ensure that the divided areas can reflect the deformation in the horizontal and vertical directions, the monitoring points corresponding to the horizontal displacement are introduced as subsidence marks, and the label combination of the subsidence marks is regarded as the basis for the subsidence area division vector boundary; for example, deformation + increase in bare soil → active subsidence area, deformation + vegetation degradation → ecological damage area, crack + mineral exposure → fault zone expansion area, low deformation + material stability → stable area, high settlement + high horizontal displacement → composite subsidence area, low settlement + high horizontal displacement → horizontal slip area; these contents can be represented as shown in Table 1. It should be noted that these combinations only represent some parts that can be combined when subsidence is recognized in the image, and the data for combination listed below are only for illustration. Other values can be used in actual judgment to select the subsidence marks for combination.

[0042] Table 1. Schematic diagram of mark combinations

[0043]

[0044] In Table 1, the ways of various subsidence identification combinations are illustrated. For example, the deformation rate described in the first row can be the deformation rate in any direction. When the deformation rate in any direction exceeds a certain value, severe deformation is likely to occur. "mm / yr" represents millimeters per year, and the proportion of bare soil change indicates the loss of surface soil or the exposure of rock strata, which is used to illustrate the intensity of surface change. When it exceeds a certain value later, it can be considered an active subsidence area. For the second row, it shows that when the vegetation index decreases and the cumulative settlement in this area reaches a certain level, it indicates the damage to the ecological environment caused by the subsidence in this area, which belongs to the area that needs ecological restoration. The vegetation index generally ranges from -1 to 1, and the normal value is generally greater than 0.2. At this time, choosing 0.3 is to judge the form of its index change to illustrate the corresponding impact of mining in the mining area on the environment. For the third row, it shows whether there are rapidly expanding cracks during ground subsidence. The calcite index is the index calculated when identifying ground markers from spectral satellite data, which is used to describe that the proportion of mineral spectral characteristics exceeds the threshold, indicating the activity of the fault zone; corresponding disasters may occur in this area, and measurement needs to be strengthened. For the fourth row, it shows that when the deformation rate in any direction is less than a certain value and the proportion of bare soil change is small, it indicates that there is no significant risk in this area in the near future, and the monitoring frequency can be reduced, etc. For the fifth row, it further differentiates the existing subsidence and illustrates it with the deformation rates in the horizontal and vertical directions. This label combination can be used as a secondary label for the larger scenarios described in the first four rows to illustrate the form of drastic changes in the deformation in this area. Similarly, for the sixth and seventh rows, it should be noted that the standard deviation of the displacement direction is calculated based on the eastward component and northward component of the identified deformation rate in the corresponding area, and the arctangent function is calculated with the ratio of the northward component to the eastward component to obtain the displacement direction angle. Then, the standard deviation is obtained by checking multiple deformation positions within the area of the main label of the label combination. According to these data, the subsidence area can be further classified, and the characteristics of each area can be analyzed to identify the expected evolution trend of the current subsidence area, which is convenient for subsequent staff to perform data backtracking and subsidence prediction to prevent mining area construction risks.

[0045] As Figure 3 shown, the implementation method of the subsidence boundary identification module includes: using the deformation amount, bare soil, and vegetation time-series change values of the subsidence area as subsidence markers, performing feature encoding on the subsidence markers, and obtaining the regional subsidence feature set.

[0046] According to the preset determination logic of the regional subsidence feature set, performing feature clustering on the subsidence markers to obtain the marker combination set.

[0047] According to the positional inclusion relationship of each element in the marker combination set, connecting the elements of the same category in the marker combination set to obtain the classified label area; outputting the boundary line connecting the classified label areas as the subsidence vector boundary.

[0048] The preset determination logic here can be set according to the triggering conditions described in Table 1 and the priorities represented by the classification labels. For example, the priority of the active subsidence area is the highest. Then, the data extracted from the set of regional subsidence features is clustered to obtain multiple groups of data with different label combinations such as high deformation amounts and low vegetation indices. Then, the clustering results are used as the set of flag combinations to illustrate the data set of various subsidence flag combinations. After that, the positions corresponding to these data sets in the mining area subsidence image pair are marked to illustrate whether there is an inclusion relationship between the data of each cluster. Then, the areas of the same category are connected to form continuous subsidence areas. After that, the boundary lines of these continuous subsidence areas are used as the output subsidence vector boundaries to identify the change trends of subsidence areas of different categories over time. This step takes into account the spatial continuity of the subsidence areas, avoids misjudgment of isolated points, and makes the boundary extraction more in line with the actual geological situation.

[0049] Preferably, the implementation method for feature clustering of subsidence flags includes calculating the within-class sum of squares under different numbers of clusters. This sum of squares is calculated based on the sum of squares of the corresponding data in the current cluster. For example, when the label combination includes the deformation amount and the vegetation index, the within-class sum of squares of the deformation amount and the vegetation index during clustering are calculated respectively. Then, a relationship graph between the number of clusters and the within-class sum of squares is plotted, with the number of clusters on the horizontal axis and the within-class sum of squares on the vertical axis. The point with the largest curvature of the curve in the relationship graph is selected, and the corresponding number of clusters and clustering data are used as the set of flag combinations output at this time. At this time, since the label combination includes the deformation amount and the vegetation index, the within-class sum of squares after clustering of the two are normalized and then added together, and the relationship graph is plotted with the added within-class sum of squares, and finally the set of flag combinations is output.

[0050] At the same time, when outputting the boundary lines connecting the classified label areas, since there will be partial overlaps and gaps in the line segments connected by the elements in the set of flag combinations, it is necessary to combine the boundary lines of the classified label areas they connect. For example, according to the trend of the boundary lines of the classified label areas, the boundary lines of each classified label area are combined, and the nearest boundary lines are judged. When the distance between the nearest boundary lines is less than the preset interval, the nearest boundary lines are combined, and the combined boundary lines are used as the boundary lines connecting the classified label areas. At this time, the preset interval can be set to 1.5 times the data spatial resolution. For example, if the resolution of the currently collected mining area subsidence image pair is 1 meter, the preset interval is 1.5 meters, which is used to identify the areas covered by different combined labels on these images.

[0051] In one embodiment of the present invention, deformation components in multiple directions are used to describe the deformation variables of the deformation point in the three-dimensional XYZ directions; that is, the XY axis on the image is taken as its horizontal direction, and the Z axis is taken as its vertical direction. The data in the corresponding vertical direction is used to illustrate the relevant deformation variables in the InSAR data to explain the size of the subsidence generated in the current area.

[0052] When dividing grid cells, the subsidence evolution module sets them based on the maximum and minimum subsidence directions within the subsidence vector boundary. The maximum subsidence direction indicates the direction of most significant subsidence in the area corresponding to the subsidence vector boundary, and grid cells are set along this direction. The minimum subsidence direction then identifies the trend in the direction with the least subsidence. This approach ensures that changes in key deformation directions are fully expressed, reducing data distortion caused by improper grid division. This allows for more accurate reflection of actual surface deformation during subsidence monitoring.

[0053] Preferably, the maximum subsidence direction and the minimum subsidence direction are set based on the deformation rate value at the vertical direction of the subsidence vector boundary. If there are two adjacent points on the subsidence vector boundary with the largest difference in deformation rate, the area covered by the subsidence vector boundary is divided into multiple grid units according to the vertical direction corresponding to the line connecting these two points. The minimum subsidence direction is set in the same way, and is also set according to the difference in deformation rate in the subsidence vector boundary. At this time, if the grid units divided by the two directions overlap, the corresponding grid units are merged and combined to serve as the divided grid units.

[0054] like Figure 4 As shown in FIG, the implementation method of the subsidence evolution module includes: setting grid cells according to the maximum subsidence direction and the minimum subsidence direction of the subsidence vector boundary, and adaptively dividing each grid cell. When dividing, a 500m×500m grid can be generated based on the area of the area included in the maximum subsidence direction and the minimum subsidence direction, or the grid cells can be obtained by dividing the area by equal size.

[0055] The deformation points and deformation components stored in each grid cell are detected, and the deformation components are combined with the position index of the deformation points to set a continuous deformation surface. The continuous deformation surface combines the deformation components of the deformation points through Kriging interpolation to form data in the form of a surface, thereby converting sparse points into a continuous map, which is convenient for describing the deformation position, range and corresponding intensity in the area corresponding to the subsidence vector boundary. Then, the deformation rate and other values of these deformation points in the grid cell are used to identify the relationship between adjacent grid cells, so as to obtain the associated positions with subsidence association and the integral probability of deformation at each point, and finally describe its deformation trend in multiple directions.

[0056] It should be noted that when describing the subsidence evolution trend in multiple directions, the deformation component in the vertical direction is mainly used to illustrate the correlation between its associated positions, while the two components in the horizontal direction mainly illustrate the relative content of the integral probability in the horizontal direction when the horizontal displacement expands, indicating the total probability of the deformation amount at each point under horizontal pushing.

[0057] Perform probability statistics on each deformation point in the continuous deformation surface, calculate the integral probability that the deformation component in the horizontal direction of each deformation point exceeds a given threshold; connect numerical intervals based on the integral probability of each deformation point to set the horizontal deformation path.

[0058] Preferably, the integral probability is used to quantify the cumulative probability that the deformation amount in a certain area exceeds a given threshold, and is used to reflect the distribution of its subsidence risk. Generally, the given threshold for a single grid cell will use the average value of the deformation amount under normal conditions as its given threshold. Assuming that the deformation amount follows a Gaussian distribution process, at this time, the form of Gaussian distribution can be combined with the form of Bayesian probability to describe the deformation amount at each deformation point.

[0059] Then the integral probability of the deformation point is expressed as: ; where represents the deformation component in the horizontal direction of deformation point i greater than the given threshold of the integral probability value, represents pi, represents the standard deviation of the deformation component in the horizontal direction, and this standard deviation represents the standard deviation of the deformation points in the continuous deformation surface where the current deformation point is located. represents the average value of the deformation component in the horizontal direction, and the average value is obtained in the same way as the standard deviation; Denotes the exponential constant. The integral probability of the deformation points calculated at this time is based on the deformation values within the grid cells where the deformation points are located and the values of the given thresholds to determine the trends generated by each deformation point during deformation. The integral probability calculated in this part represents the key part of the probability estimation. When the value of this integral probability is large, such as when the value of a certain deformation point is close to 1, it means that the deformation value at this position may exceed the set threshold, and there may be relatively serious subsidence problems in this area, such as goaf collapse, surface subsidence caused by groundwater level decline, etc. The surface stability of this area is poor, and this area needs to be marked. If this integral probability approaches a medium probability, such as between 0.3 and 0.7, it indicates that the deformation value in this area may exceed the given threshold and may be in the transition zone between the rapid subsidence area and the stable area. The surface deformation is developing but at a slow speed. This part of the area needs to be observed. If this integral probability is close to 0, for example, if it is described as less than 0.3, it means that there is no obvious subsidence phenomenon in these areas, and there is no need to deliberately identify the subsidence evolution trend. At the same time, when calculating the integral probability of each deformation point, its deformation components can also be ignored, and the integral probability can be directly calculated using the corresponding deformation in the InSAR data to illustrate the correlation in all directions.

[0060] Preferably, when connecting numerical intervals according to the integral probability of each deformation point, the integral probability is divided into three intervals according to the magnitude of the value. The deformation points corresponding to the high probability, low probability, and medium probability in the three intervals are used to sequentially select the nearest deformation points according to the position index for connection to form a horizontal deformation path to represent the path under different horizontal deformation types. The thresholds for dividing these intervals can use 0.7 and 0.3 to divide its integral probability, or use the average value of the thresholds set as medium probability, low probability, and high probability in historical data as the form for dividing its integral probability.

[0061] Afterwards, when identifying its associated positions, hotspot analysis is performed using the deformation components in the vertical direction of each deformation point to identify the existing subsidence centers or diffusion paths to illustrate the existing paths and trends. The deformation classification in the vertical direction generally refers to its settlement amount. After combining multiple points identified by the horizontal displacement, the subsidence trend covering the corresponding details in both directions can be combined, which can improve the sensitivity of the system to different deformation types, thereby achieving more accurate mine area management.

[0062] That is, the implementation method of the subsidence evolution module also includes: calculating the Moran index corresponding to each deformation point based on the distance between each deformation point in the continuous deformation surface and the deformation components in the vertical direction; screening each deformation point based on the Moran index to obtain the associated positions corresponding to each deformation point; at this time, the associated positions are used to represent multiple deformation points that have spatial correlation with the current deformation point.

[0063] After connecting the deformation points contained in the associated positions, an associated position path is formed. The intersection of the associated position path and the horizontal deformation path is used as the subsidence evolution path, and the deformation points on the subsidence evolution path are converted into subsidence evolution trend output.

[0064] Preferably, the Moran index can be expressed as ;in, The Moran index represents the deformation point. Here, it is calculated based on the deformation component in the vertical direction. For deformation point j, other deformation points in the continuous deformation surface corresponding to deformation point i are selected. The value range of i and j is 1 to N. Represents the number of deformation points; Represents the weight of deformation points i and j, which is set as the inverse of the distance between deformation points i and j; represents the value of the deformation classification in the vertical direction of the deformation point i, represents the value of the deformation classification in the vertical direction of the deformation point j, Represents the average value of the vertical component of the deformation point. The Moran index is used to cluster deformation points in the vertical direction. The Moran index ranges from -1 to 1. When it is close to 1, it shows a strong positive correlation, which is prone to concentrated subsidence and clustering in which the core of the concentrated subsidence expands to the periphery. When it is close to 0, the subsidence distribution is random. When it is close to -1, there is a strong negative correlation, and subsidence and stable parts are staggered, which is prone to geological boundaries and local instability. When using the Moran index to screen deformation points, any deformation point is used as the current deformation point to be processed. Based on the Moran index, multiple deformation points corresponding to its positive correlation and multiple deformation points corresponding to its negative correlation are obtained respectively. The calculation is iterated until the calculated Moran index values are close to the range of -1 or 1. These multiple deformation points with positive and negative correlation are regarded as associated positions. For example, the range of these two values is close to ±1, divided by ±0.6, or the average of the positive and negative correlation values determined in the sampling historical data is set.

[0065] When obtaining the associated position of each deformation point, the implementation method includes: based on the Moran index, obtaining the positive correlation points and negative correlation points of the current deformation point, and outputting the positive correlation points and negative correlation points as the associated positions. When judging the positive correlation points and negative correlation points, the average value of the positive correlation and negative correlation of each deformation point can be used to select the positive correlation and negative correlation areas corresponding to the current deformation point based on historical data. These areas are used as the main areas for subsidence judgment to find the locations where unstable boundaries exist and the locations where subsidence is obviously concentrated. After combining these locations with the deformation generated in the horizontal direction, a relatively complete subsidence evolution trend can be obtained, thereby indicating the areas where subsidence changes exist in the current area, which is convenient for predicting the subsidence areas that will exist under subsequent mining.

[0066] In the subsidence evolution trend diagram, it is displayed according to the magnitude of the deformation amount that appears on the map along the subsidence evolution path. Its horizontal axis can be the geographical location coordinates or the distance along the subsidence evolution path, and then the vertical axis shows the deformation rate. It is also possible to use two vertical axes respectively, such as the YZ axis, to display the horizontal deformation rate and the vertical deformation rate respectively.

[0067] In one embodiment of the present invention, the regional deformation characteristics are used to describe the subsidence patterns presented under the deformation of multiple grid cells in multiple time windows, such as subsidence forms like bowl-shaped, trough-shaped, step-shaped, etc. At the same time, the subsidence stage will be described according to the time window, such as the initial stage, the active stage, and the decline stage.

[0068] When the subsidence coupling module couples the subsidence evolution trend, it aggregates adjacent grid cells based on the deformation rate and the subsidence evolution direction of the subsidence evolution trend, describes the subsidence pattern in the time domain corresponding to the time window of the subsidence evolution trend, and then outputs the aggregated partial features as its regional deformation characteristics.

[0069] Such as Figure 5 As shown, the implementation method of the subsidence coupling module includes: extracting the deformation rate and the subsidence evolution path in each grid cell of the subsidence evolution trend, and taking the direction consistency and the deformation rate of the subsidence evolution path in adjacent grid cells as the deformation gradient.

[0070] Use the deformation gradient to query the subsidence pattern and output it as the regional deformation characteristics according to the subsidence stage corresponding to the time window.

[0071] Preferably, after setting the angles of each point in the subsidence evolution path in adjacent grid cells in the form of its map coordinate system, the average value of the cosine values obtained by subtracting the angles of each point in the subsidence evolution path in adjacent grid cells is used as the recognized direction consistency.

[0072] Preferably, when obtaining the subsidence evolution trend, the deformation points on each grid cell are displayed as a continuous deformation surface to illustrate the subsidence evolution process in the continuous deformation surface. Then, the subsidence evolution trend at this time is mapped onto each grid cell to identify the positions with obvious subsidence on the grid cells divided according to the maximum subsidence direction and the minimum subsidence direction. After that, the subsidence patterns represented by these contents are output. For example, if the direction consistency between two adjacent regions is relatively small, and the deformation rate is the highest at the center of these two regions, and then it is verified that the subsidence evolution path decreases outward in a radial form, it indicates that there is a bowl-shaped subsidence in these adjacent regions, and the direction consistency at this time may be less than 0.2. If the deformation rates of two adjacent grid cells increase along a specific direction and the value of their direction consistency is greater than 0.8, it can be considered as a trough-shaped subsidence. For a stepped subsidence, it is when the difference in deformation rates between adjacent grid cells is greater than 5 mm / km and the value of the direction consistency is less than 0.5, then it is considered a stepped subsidence. All these contents can be obtained by querying and retrieving data according to the currently identified deformation gradient to determine what type of subsidence situation the adjacent grid cells belong to in the subsidence evolution trend. After that, when outputting the regional deformation characteristics under multiple time windows, the deformation rate and other data calculated within the current time window are used to calculate the overlap ratio with the data of the next time window. When the overlap ratio exceeds 70%, it is considered that the subsidence stages of adjacent time windows are consistent, and the corresponding subsidence pattern is output as the regional deformation characteristic. Otherwise, the deformation gradient of the subsidence pattern is used as its regional deformation characteristic for output.

[0073] In an embodiment of the present invention, the integrated regional deformation characteristics will combine the change trends of the subsidence evolution trend in adjacent grid cells to obtain, for each grid cell after being divided according to the subsidence vector boundary, in which direction the boundary expands under the influence of time and the expansion rate, so as to obtain the subsidence amount, deformation rate, affected area, etc. corresponding to the subsidence mapping level of the regional subsidence under the current subsidence mapping, which is convenient for subsequent adjustment of the location or intensity of the mining construction in the mining area based on the mapping data.

[0074] When integrating the regional deformation characteristics, the symmetric difference set of the regional deformation characteristics is used to identify the subsidence expansion rate and direction. The symmetric difference set is expressed as the part of the union of the regional deformation characteristics minus the intersection, which is used to identify the different parts in multiple time windows. The subsidence mapping level is defined by the subsidence amount, deformation rate, and affected area at the positions of each part, indicating the existence of high risk, low risk, and medium risk, which helps the mining area managers adjust the location or intensity of the mining construction to reduce the risks brought by subsidence.

[0075] The implementation method of the subsidence evaluation module further includes: obtaining the symmetric difference set of the regional deformation features, and setting the subsidence mapping level according to the path direction of the symmetric difference set, using the subsidence amount, deformation rate, and influence area at the corresponding position as retrieval conditions.

[0076] As shown in Table 2, the corresponding subsidence mapping levels can be expressed as follows.

[0077] Table 2. Schematic Diagram of Subsidence Mapping Levels

[0078]

[0079] In Table 2, the data that needs to be identified in the corresponding area of the symmetric difference set is described. These data are only used for adaptive display. When the type of the area where subsidence is collected changes, the required values such as subsidence amount, deformation rate, and influence area will be adaptively adjusted; at the same time, the direction of the symmetric difference set is actually extracted from multiple subsidence evolution paths existing in the subsidence evolution trend after mapping the data of the difference set, or the direction of multiple paths formed by the integral probability and associated positions of the deformation points in the subsidence evolution module is selected for identification. After calculating the subsidence amount, deformation rate, and influence area of the small areas included in these directions, the situations of other combinations of regional deformation features are continuously traversed, so as to complete the entire process of identifying the subsidence in the current mining area, and all the identified multiple data indicators are sent to the external data terminal. According to the data of the external data terminal, the deformation of the corresponding area can be viewed, which helps the staff to timely check in which dimensions there are problems in the current area, and improves the accuracy and efficiency of the subsidence warning in the mining area.

[0080] Preferably, the influence area refers to the area where the symmetric difference set exists in the path direction. This area is the calculation content directly obtained when the symmetric difference set is subtracted. The deformation rate is obtained based on the deformation rate measured at the deformation points existing in the symmetric difference set, and the subsidence amount is also inferred based on the value of the deformation rate in the vertical direction.

[0081] Although the embodiments of the present invention have been shown and described above, it can be understood that the above embodiments are exemplary and should not be construed as limiting the present invention. Those of ordinary skill in the art can make changes, modifications, substitutions, and variations to the above embodiments within the scope of the present invention, and still be covered by the protection scope of the present invention.

Claims

1. A coal mine subsidence surveying and mapping system, characterized in that, Including: A data acquisition module, which is used to obtain pairs of subsidence images of a mining area, and associate the pairs of subsidence images of the mining area with the subsidence area according to the horizontal displacements of the subsidence images of the mining area in multiple time windows; A subsidence boundary recognition module, which is used to extract subsidence markers in the pairs of subsidence images of the mining area, classify the subsidence area according to the temporal changes of the subsidence markers, and output a subsidence vector boundary; A subsidence evolution module, which is used to divide grid cells according to the subsidence vector boundary, identify deformation points and deformation components in multiple directions within each grid cell, and form a continuous deformation surface; [[ID=G4]]Obtain a subsidence evolution path through the integral probability and associated positions of each deformation point in the continuous deformation surface, and output a subsidence evolution trend based on the deformation rate on the subsidence evolution path; A subsidence coupling module, which is used to couple the subsidence evolution trend according to the positions of the grid cells, and identify the regional deformation characteristics of the subsidence evolution trend in the corresponding time window based on the deformation gradients of adjacent grid cells; A subsidence evaluation module, which is used to integrate the regional deformation characteristics, and set the subsidence surveying and mapping levels of each subsidence area according to the rates and directions of subsidence expansion of each grid cell.

2. The coal mine subsidence surveying and mapping system according to claim 1, wherein, The implementation method of the data acquisition module includes: Obtain pairs of historical subsidence images of the current subsidence area, compare the image boundaries of the pairs of historical subsidence images of the mining area with the pairs of current subsidence images of the mining area, and obtain the aligned pairs of subsidence images of the mining area; According to the aligned pairs of subsidence images of the mining area, calculate the horizontal displacements of the monitoring points in the subsidence area, and associate the pairs of subsidence images of the mining area with the subsidence area by using the horizontal displacements of the monitoring points.

3. A coal mine subsidence surveying and mapping system according to claim 2, wherein, The implementation method of calculating the horizontal displacements of the monitoring points in the subsidence area includes: Using the orthophoto inverse projection technology, obtain the pixel coordinates of each monitoring point in the subsidence pair of the mining area; Match the monitoring points in adjacent time series, use the center points of the monitoring points as target feature points, calculate the matching degrees of the target feature points, and use the distances of the target feature points with the maximum matching degrees in the adjacent time series as the horizontal displacements of the monitoring points in the subsidence area.

4. A coal mine subsidence surveying and mapping system according to claim 1, characterized in that, The implementation method of the subsidence boundary recognition module includes: Use the deformation amount, bare soil, and vegetation temporal change values of the subsidence area as subsidence markers, perform feature encoding on the subsidence markers, and obtain a set of regional subsidence characteristics; According to the preset determination logic of the set of regional subsidence characteristics, perform feature clustering on the subsidence markers, and obtain a set of marker combinations; According to the position inclusion relationship of each element in the set of marker combinations, connect the elements of the same category in the set of marker combinations to obtain a classified label area; output the boundary line connected by the classified label area as a subsidence vector boundary.

5. A coal mine subsidence surveying and mapping system according to claim 4, characterized in that, The implementation method of performing feature clustering on the subsidence markers includes; By calculating the within-class sum of squares under different numbers of clusters, draw a relationship diagram of the number of clusters and the within-class sum of squares, with the number of clusters on the horizontal axis and the within-class sum of squares on the vertical axis, select the point with the largest curvature of the curve in the relationship diagram, and use the corresponding number of clusters and clustering data of this point as the output set of marker combinations.

6. The coal mine subsidence mapping system according to claim 4, characterized in that, When outputting the boundary line of the classified label area, it also includes: Combine the boundary lines of each classification label area according to the trend of the boundary lines of the classification label area, and make a judgment based on the closest boundary line. When the distance to the closest boundary line is less than the preset interval, combine the closest boundary lines, and use the combined boundary line as the boundary line connecting the classification label areas.

7. A coal mine subsidence surveying and mapping system according to claim 1, characterized in that, The implementation method of the subsidence evolution module includes: Set grid cells according to the maximum subsidence direction and the minimum subsidence direction of the subsidence vector boundary; Detect the deformation points and deformation components stored in each grid cell, combine the deformation components according to the position index of the deformation points, and set a continuous deformation surface; Perform probability statistics on each deformation point in the continuous deformation surface, calculate the integral probability that the deformation component in the horizontal direction of each deformation point exceeds a given threshold; connect the numerical intervals according to the integral probability of each deformation point, and set a horizontal deformation path; Calculate the Moran index corresponding to each deformation point based on the distance between each deformation point in the continuous deformation surface and the deformation component in the vertical direction; screen each deformation point based on the Moran index to obtain the associated position corresponding to each deformation point; After connecting the deformation points included in the associated position, form an associated position path, use the intersection route of the associated position path and the horizontal deformation path as the subsidence evolution path, and convert the deformation points on the subsidence evolution path into subsidence evolution trends for output.

8. A coal mine subsidence surveying and mapping system according to claim 7, characterized in that, The implementation method of obtaining the associated position of each deformation point includes: Based on the Moran index, obtain the positive correlation points and negative correlation points of the current deformation point, and output the positive correlation points and negative correlation points as the associated positions.

9. The coal mine subsidence surveying and mapping system according to claim 1, characterized in that, The implementation method of the subsidence coupling module includes: Extract the deformation rate and subsidence evolution path of the subsidence evolution trend in each grid cell, and use the direction consistency and deformation rate of the subsidence evolution path in adjacent grid cells as the deformation gradient; Use the deformation gradient to query the subsidence mode, and output it as the regional deformation characteristics according to the subsidence stage corresponding to the time window.

10. A coal mine subsidence surveying and mapping system according to claim 1, characterized in that, The implementation method of the subsidence evaluation module also includes: Obtain the symmetric difference set of the regional deformation characteristics, and set the subsidence mapping level based on the path trend of the symmetric difference set, using the settlement amount, deformation rate, and influence area at the corresponding position as the retrieval conditions.

Images