Overground and underground integrated three-dimensional modeling anti-mold-penetration processing method
By processing data through a Geographic Information System (GIS), the fusion of high-precision surface models and low-precision underground models is automated, solving the problem of model penetration and ensuring the model's authenticity and morphological characteristics, thus achieving an automated process without model penetration.
Patent Information
- Authority / Receiving Office
- CN · China
- Patent Type
- Applications(China)
- Current Assignee / Owner
- Filing Date
- 2025-12-26
- Publication Date
- 2026-03-31
AI Technical Summary
In digital twin and smart geology applications, when high-precision surface models are fused with low-precision underground models, clipping can easily occur, affecting visualization and spatial analysis. Existing technologies are inefficient and the processing quality is difficult to guarantee.
The data is processed using a Geographic Information System (GIS). Through steps such as neighborhood minimum value calculation, feature area detection, and intelligent smoothing, the problem of clipping through the model is automatically handled to ensure that the underground model does not penetrate the surface model while maintaining the true shape of the model.
It achieves automated processing without model penetration, effectively handles complex features, maintains model morphology, ensures model authenticity, and realizes a complete automated process from data input to result verification.
Smart Images
Figure CN121767590A_ABST
Abstract
Description
Technical Field
[0001] This invention relates to the interdisciplinary fields of computer 3D modeling, geographic information systems, and geological engineering, and particularly to a method for preventing clipping during integrated above-ground and underground 3D modeling. More specifically, it is an intelligent processing method for solving the clipping phenomenon that occurs when 3D models of different accuracies are merged. Background Technology
[0002] In applications such as digital twins and smart geology, constructing an integrated 3D model from the surface to the subsurface is of significant value. Integrated 3D modeling of the surface and subsurface requires merging a high-precision surface model with a low-precision subsurface 3D geological model. However, the data sources and accuracy of the two types of models differ significantly: surface models are typically acquired using technologies such as lidar and oblique photogrammetry, generating detailed models with centimeter-level accuracy that accurately reflect topographic relief details; in contrast, subsurface 3D geological models are constructed based on limited geological map data, borehole data, and geophysical exploration data, with an accuracy of only meters or even lower.
[0003] When these two types of models with vastly different levels of accuracy are fused in the same three-dimensional space, the lower accuracy and inference errors of the underground three-dimensional geological model often result in unavoidable protrusions on its surface that penetrate the high-precision surface three-dimensional model, creating a clipping phenomenon. This geometric conflict not only affects the visualization effect but also renders the model unusable for accurate spatial analysis and engineering calculations.
[0004] Existing solutions suffer from two main drawbacks: first, they rely on manual editing, adjusting penetration vertices one by one, which is inefficient and makes it difficult to guarantee processing quality; second, strategies that involve lowering the overall underground 3D geological model or raising the overall surface 3D model, while avoiding penetration, significantly increase the gap between the above-ground and underground, thus compromising the model's realism and application value. Therefore, there is an urgent need for a technical solution that can automatically and accurately handle the penetration problem while preserving the model's true shape to the greatest extent possible. Summary of the Invention
[0005] The technical problem to be solved by the present invention is to provide an integrated above-ground and underground modeling method that can automatically process massive amounts of three-dimensional data, strictly prevent model penetration, and maintain the authenticity of underground geological morphology, thereby overcoming the defects of the above-mentioned prior art.
[0006] To solve the above-mentioned technical problems, the present invention adopts the following technical solution:
[0007] A method for preventing mold penetration in integrated above-ground and underground 3D modeling, using a Geographic Information System (GIS) to process data, includes the following steps: S01, Data Conversion: Convert the surface DSM to ASC format data; S02, Neighborhood Minimum Calculation: Determine the thinning ratio and neighborhood radius, and preliminarily calculate the elevation value of the new terrain; S03, Feature Region Detection: Identifying high-risk feature regions by fusing curvature and gradient features; S04, Subsurface Elevation Generation: Generating subsurface elevation based on neighborhood minimum and characteristic regions; S05, Perform intelligent smoothing: Intelligent smoothing is performed using feature-weighted anisotropic diffusion; S06, Strictly enforced constraints: For detected penetration points, the elevation value is forcibly corrected to the minimum value; S07, Output Data: Output the thinned ASCII format data.
[0008] Based on the above technical solutions, the present invention may be supplemented with the following technical means to better or more specifically solve the technical problems to be solved by the present invention:
[0009] When executing S02, for each thinned subsurface model grid cell, calculate the minimum elevation value within its corresponding original high-precision surface area and surrounding neighborhood: Let the horizontal grid size of the original 3D surface model be... The horizontal grid size of the target underground three-dimensional geological model is ,in, , ; For each target grid The corresponding original region boundary is: ; ; ; ; In the formula: h and w are the sampling intervals for the underground horizontal grid, where h = m / p, indicating that in the east-west direction, each underground horizontal grid corresponds to h consecutive original surface horizontal grids; w = n / q, indicating that in the north-south direction, each underground horizontal grid corresponds to w consecutive original surface horizontal grids; h × w is the thinning factor of the underground 3D geological model's horizontal grid. r is the neighborhood radius, a dimensionless expansion factor. When determining the original search area, the original projection area of the current underground horizontal grid is used as the center, and r×h original surface horizontal grids are expanded outward in the east-west direction, and r×w original surface horizontal grids are expanded outward in the north-south direction. This constitutes the neighborhood range, and the minimum original surface elevation value within this neighborhood range is... ; In the formula, Hsurface H is the surface elevation value. surface (x,y) represents the elevation value of any point (x,y) on the Earth's surface, within the neighborhood.
[0010] Furthermore, the value of r ranges from 1 to 3. When r=1, it means that the search range covers the current grid and its immediate surrounding grids. As the value of r increases, the search range expands accordingly, which can further avoid the potential clipping risk from terrain undulations further away.
[0011] Furthermore, when performing S03, the detected feature regions include local geometric feature regions formed by abrupt changes in terrain or land features.
[0012] Furthermore, the local geometric feature region is a curvature feature. Or the region where the gradient feature G exceeds the set threshold.
[0013] Furthermore, when executing S03, the calculation formula for the curvature feature is as follows: ; In the formula, Indicates curvature characteristics; It is the second partial derivative, which represents the rate of change of the east-west slope. The larger the value, the more curved the east-west slope. It is the second partial derivative, which represents the rate of change of the north-south slope. The larger the value, the more curved the north-south slope. The formula for calculating the gradient feature is as follows: ; In the formula, G represents the gradient feature; It is the first-order partial derivative, representing the rate of change of elevation H in the east-west direction. The larger the absolute value, the steeper the slope in the east-west direction; It is the first-order partial derivative, representing the rate of change of elevation H in the north-south direction. The larger the absolute value, the steeper the slope in the north-south direction; The calculation formula for feature region fusion is as follows: ; In the formula, and As the normalization factor, The feature region recognition coefficient, The range is 0-1, with a value close to 1 indicating a high-risk area and a value close to 0 indicating a low-risk, flat area.
[0014] Furthermore, when executing S04, the formula for calculating the initial underground elevation is as follows: ; The formula for calculating feature region adjustment is as follows: ; In the formula, H underground This is the initial underground elevation. H represents the adjusted underground elevation of the characteristic region. min This represents the lowest surface elevation within the adjacent area. This is the safety threshold.
[0015] Furthermore, when executing S05, the smoothing weight calculation formula is as follows: ; The anisotropic diffusion equation is as follows: ; In the formula, It is the first-order partial derivative, representing the rate of change of anisotropic diffusion; s represents the smoothing weight; and F is the feature region recognition coefficient. This represents a vector differential operator used to calculate the gradient of a scalar field or the divergence of a vector field; when When applied to the elevation scalar field H, it constitutes a gradient calculation; Let H represent the gradient of the elevation field H. This is a vector whose direction points in the direction of the most drastic change in surface elevation at that point, and its magnitude represents the rate of elevation change at that point. The modulus is the value of G, which is the geometric basis for identifying the steepness of the terrain and the characteristic regions.
[0016] Furthermore, in the feature region, that is, in the region where the feature region recognition coefficient F > 0.6, the smoothing weight is smaller to preserve terrain details; in the flat region, that is, in the region where F < 0.3, the smoothing weight is larger to enhance the smoothness of the model.
[0017] Furthermore, a smaller smoothing weight refers to the case where s < 0.58, while a larger smoothing weight refers to the case where s > 0.79.
[0018] Compared with the prior art, the present invention has the following beneficial technical effects:
[0019] First, it absolutely ensures no penetration: This invention uses neighborhood minimum value calculation to mathematically ensure that the underground model does not penetrate the surface model.
[0020] Second, it can effectively handle complex features: This invention is specifically designed with detection and optimization algorithms for complex geometric features such as building corners, terrain grooves, and ditches.
[0021] Third, good shape preservation: This invention uses intelligent smoothing technology with feature region perception to maintain the shape characteristics of the model while avoiding clipping.
[0022] Fourth, fully automated processing: This invention uses a Geographic Information System (GIS) to process data, realizing a complete automated process from data input to result verification. Attached Figure Description
[0023] Figure 1 This is a flowchart of the overall processing method of the present invention;
[0024] Figure 2 This is a schematic diagram of the neighborhood minimum value calculation in this invention;
[0025] Figure 3 This is a schematic diagram of the distribution of the feature regions in this invention;
[0026] Figure 4 Is with Figure 3 A corresponding diagram of the smoothed weight distribution;
[0027] Figure 5 This is a comparison image of feature-aware smoothing processing in this invention.
[0028] Figure 6 This is a schematic diagram of penetration detection and correction in this invention;
[0029] Figure 7 These are application example renderings (clipping) obtained using traditional methods.
[0030] Figure 8 It is obtained by using the method of the present invention, and Figure 7 Corresponding application example screenshots (without clipping);
[0031] Figure 9 yes Figure 7 A magnified view of a section (clipping);
[0032] Figure 10 yes Figure 8 A magnified view of a portion (without clipping), this portion is related to... Figure 9 The corresponding local area is shown. Detailed Implementation
[0033] To facilitate a thorough understanding of this technical solution by those skilled in the art, an embodiment of the present invention is described below in conjunction with the accompanying drawings.
[0034] Figure 1 This is a flowchart of the overall processing method of the present invention, which shows the complete processing flow of the present invention from data input to result output.
[0035] Figure 2This is a schematic diagram of the neighborhood minimum value calculation in this invention; the diagram shows the correspondence between the original high-precision grid and the thinned grid, reflecting the process of neighborhood range expansion and minimum value extraction.
[0036] Figure 3 This is a schematic diagram of the distribution of feature regions in this invention. The diagram shows the identification results of feature regions such as building corners, terrain depressions, and ditches. The gray boxes in the diagram represent areas with large surface undulations, and the white boxes represent areas with small surface undulations.
[0037] Figure 4 Is with Figure 3 The corresponding smoothing weight distribution diagram shows that the value 0.3 in the diagram represents the smoothing effect (weak smoothing effect) in areas with large surface undulations, and the values 0.6 and 1.0 represent the smoothing effect in areas with small surface undulations. Among them, the value 0.6 represents the medium smoothing effect, and the value 1.0 represents the strong smoothing effect.
[0038] Figure 5 This is a comparison chart of feature-aware smoothing processing in this invention. The left side of the chart, from top to bottom, represents the model curves of the building corner area (F≈1, s≈0.3), ditch area (F≈0.6, s≈0.58), and flat area (F≈0.1, s≈0.93) before processing (before correction). The right side of the chart, from top to bottom, represents the model curves of the building corner area, ditch area, and flat area after processing (after correction). F is the feature region recognition coefficient, and s represents the smoothing weight.
[0039] Figure 6 This is a schematic diagram of penetration detection and correction in this invention; the diagram illustrates the detection method and correction principle of penetration points: for the detected penetration points, the elevation value is forcibly corrected to the minimum value.
[0040] Figures 7 to 10 These are two sets of comparison images, specifically... Figure 7 and Figure 8 It is a set of comparison pictures. Figure 9 and Figure 10 Is with Figure 7 and Figure 8 Another set of related comparison charts. Figure 7 These are application example images obtained using traditional methods (there is clipping, see the pink part in the image). Figure 8 It is obtained by using the method of the present invention, and Figure 7 The corresponding application example image (without clipping) shows a landslide in Gaolan Mountain, Lanzhou City, and the surrounding town area. Figure 9 yes Figure 7 A magnified view of a section (exhibiting clipping). Figure 10 yes Figure 8A magnified view of a portion of the image (no clipping), this portion is related to... Figure 9 The corresponding local area is shown.
[0041] like Figure 1 and combined Figures 2-6 As shown, a method for preventing mold penetration in integrated above-ground and underground 3D modeling, using a Geographic Information System (GIS) to process data, is characterized by including the following steps:
[0042] S01, Data Conversion: Convert the Digital Surface Model (DSM) data to ASC format, that is, convert the Digital Surface Model (DSM) data from its original format to ASC (American Standard Code for Information Interchange) text format.
[0043] S02, Neighborhood Minimum Calculation: Determine the thinning ratio and neighborhood radius, and preliminarily calculate the elevation value of the new terrain.
[0044] S03, Feature Region Detection: Identify feature regions with high risk of mold penetration by fusing curvature features and gradient features.
[0045] S04, Underground Elevation Generation: Generates underground elevation based on neighborhood minimum and feature region.
[0046] S05, Perform intelligent smoothing: Intelligent smoothing is performed using feature-weighted anisotropic diffusion.
[0047] S06, Strict Constraint Execution: For the detected penetration point, the elevation value is forcibly corrected to the minimum value.
[0048] S07, Output Data: Output the thinned ASCII format data.
[0049] In this embodiment, when executing S02, for each thinned underground model grid cell, the minimum elevation value of its corresponding original high-precision surface area and surrounding neighborhood is calculated. The specific calculation method is as follows:
[0050] Let the horizontal grid size of the original 3D surface model be... The horizontal grid size of the target underground three-dimensional geological model is ,in, , ; For each target grid The corresponding original region boundary is: ; ; ; ; In the formula: h and w are the sampling intervals for the underground horizontal grid, where h = m / p, indicating that in the east-west direction, each underground horizontal grid corresponds to h consecutive original surface horizontal grids; w = n / q, indicating that in the north-south direction, each underground horizontal grid corresponds to w consecutive original surface horizontal grids; h × w is the thinning factor of the underground 3D geological model's horizontal grid; regarding the thinning factor, as... Figure 2 As shown, the surface has a 10×10 grid (top of the image), and the underground becomes a 1×1 grid (bottom middle of the image). The thinning factor is (10 / 1)×(10 / 1)=100 times. Let h×w represent the area thinning factor, and h and w alone represent the thinning factors in the east-west and north-south directions. The neighborhood radius is a dimensionless expansion factor. When determining the original search area, the original projection area of the current underground horizontal grid is used as the center, and the area is expanded outward by r×h original surface horizontal grids in the east-west direction and by r×w original surface horizontal grids in the north-south direction, thus forming the neighborhood range. The minimum original surface elevation value within this neighborhood range is... The calculation formula is as follows: ; In the formula, H surface H is the surface elevation value. surface (x,y) represents the elevation value of any point (x,y) on the Earth's surface, within the neighborhood.
[0051] In this embodiment, the typical range of the neighborhood radius r is 1 to 3. When r=1, it means that the search range covers the current grid and its immediate surrounding grids. As the value of r increases, the search range expands accordingly, which can further avoid the potential clipping risk from terrain undulations further away.
[0052] In this embodiment, when executing S03, the detected feature regions include local geometric feature regions formed by abrupt changes in terrain or landforms. For example, building corners, recesses, ditches, and walls are typical local geometric feature regions. Local geometric feature regions are also curvature features. Or the region where the gradient feature G exceeds a set threshold. The curvature feature... The calculation formula is as follows: ; In the formula, It is the second partial derivative, which represents the rate of change of the east-west slope. The larger the value, the more curved the east-west slope. It is the second partial derivative, representing the rate of change of the north-south slope. The larger the value, the more curved the north-south slope.
[0053] The formula for calculating the gradient feature G is as follows: ; In the formula, It is the first-order partial derivative, representing the rate of change of elevation H in the east-west direction. The larger the absolute value, the steeper the slope in the east-west direction; It is the first-order partial derivative, representing the rate of change of elevation H in the north-south direction. The larger the absolute value, the steeper the slope in the north-south direction; The calculation formula for feature region fusion is as follows: ; In the formula, and As the normalization factor, The feature region recognition coefficient, The range is 0-1, with values close to 1 indicating a high-risk region and values close to 0 indicating a low-risk, flat region. In the actual implementation of this invention, two empirical thresholds, 0.3 and 0.6, can be set: F < 0.3 represents a flat region, 0.3-0.6 represents an intermediate region, and F ≥ 0.6 represents a high-risk region.
[0054] In this embodiment, when executing S04, the formula for calculating the initial underground elevation is as follows: ; The formula for calculating feature region adjustment is as follows: ; In the formula, H underground This is the initial underground elevation. H represents the adjusted underground elevation of the characteristic region. min This represents the lowest surface elevation within the adjacent area. This is the safety threshold. The numerical value is related to the model accuracy; the lower the model accuracy, the better. The larger the value, the higher the model accuracy. The smaller the value, the better. The value generally does not exceed the size of the underground grid; for example, if the underground grid is 10 meters, then... The value can be selected between 0 and 10 meters based on the model's accuracy.
[0055] When executing S05 in this embodiment, the smoothing weight calculation formula is as follows: ; The anisotropic diffusion equation is as follows: ; In the formula, It is the first-order partial derivative, representing the rate of change of anisotropic diffusion; s represents the smoothing weight; and F is the feature region recognition coefficient. This represents a vector differential operator used to calculate the gradient of a scalar field or the divergence of a vector field; when When applied to the elevation scalar field H, it constitutes a gradient calculation; Let H represent the gradient of the elevation field H. This is a vector whose direction points in the direction of the most drastic change in surface elevation at that point, and its magnitude represents the rate of elevation change at that point. The modulus is the value of G, which is the geometric basis for identifying the steepness of the terrain and the characteristic regions.
[0056] In feature regions, i.e., regions where the feature region recognition coefficient F > 0.6, the smoothing weight is smaller to preserve terrain details; in flat regions, i.e., regions where F < 0.3, the smoothing weight is larger to enhance model smoothness. The smaller smoothing weight refers to the case where s < 0.58, and the larger smoothing weight refers to the case where s > 0.79.
[0057] In this embodiment, when executing S06, the formula for forced correction of the detected penetration point is: H corrected =H min — ; In the formula, This indicates a forced correction of the underground elevation value.
[0058] The above, in conjunction with the attached diagram ( Figures 1-6 This paper provides a detailed description of one embodiment of the present invention. It should also be noted that the technical solution of the present invention has been verified in a landslide area and surrounding towns in Gaolan Mountain, Lanzhou City. Figure 8 and Figure 10 As shown, the three-dimensional geological model obtained by applying the technical solution of this invention does not exhibit clipping, while the three-dimensional geological model obtained by applying traditional methods exhibits more severe clipping (see [link]). Figure 7 and Figure 9 In summary, in the process of integrated 3D modeling of above-ground and underground environments, this invention can not only effectively handle complex features, but also maintain good morphological features of the model while avoiding clipping, and realize a complete automated process from data input to result verification.
Claims
1. A method for preventing penetration of a mold in a three-dimensional modeling of an above-ground and underground integration, using a geographic information system to process data, characterized by, The method comprises the following steps: S01, data conversion: converting the ground DSM into ASC format data; S02, neighborhood minimum value calculation: determining the thinning ratio and neighborhood radius, and preliminarily calculating the elevation value of the new terrain; S03, feature region detection: identifying the feature region of high penetration risk by fusing the curvature feature and the gradient feature; S04, underground elevation generation: generating the underground elevation based on the neighborhood minimum value and the feature region; S05, intelligent smoothing: performing intelligent smoothing by using feature-weighted anisotropic diffusion; S06, strict constraint execution: for the detected penetration point, the elevation value is forced to be corrected to the minimum value; S07, output data: outputting the thinned ASC format data.
2. The above-and-underground integrated three-dimensional modeling penetration prevention processing method according to claim 1, characterized by: When S02 is performed, for each underground model grid cell after thinning, the minimum elevation value in the corresponding original high-precision ground area and the surrounding neighborhood range is calculated: The horizontal grid size of the original surface three-dimensional model is The horizontal grid size of the target underground three-dimensional geological model is wherein, , ; For each target grid whose corresponding original region boundary is: ； ； ； ； In the formula, h and w are the underground horizontal grid sampling intervals, wherein h=m / p represents that in the east-west direction, each underground horizontal grid corresponds to h consecutive original ground horizontal grids; w=n / q represents that in the north-south direction, each underground horizontal grid corresponds to w consecutive original ground horizontal grids; h×w is the thinning multiple of the underground three-dimensional geological model horizontal grid; The value range of r is 1 to 3, when r=1, it represents that the search range covers the current grid and the adjacent surrounding grid, and as the value of r increases, the search range correspondingly expands. r is a neighborhood radius, which is a dimensionless expansion multiple, and in determining the original search area, the original projection area of the current underground horizontal grid is taken as the center, r x h original ground surface horizontal grids are expanded outward in the east-west direction, and r x w original ground surface horizontal grids are expanded outward in the south-north direction, thereby forming a neighborhood range, and the original ground surface minimum elevation value in the neighborhood range ; In the formula, H surface is a ground elevation value, H surface (x, y) represents an elevation value of any point (x, y) on the ground, (x, y) within the neighborhood range.
3. The above-and-underground integrated three-dimensional modeling penetration prevention processing method according to claim 2, characterized by: When S03 is performed, the detected feature region includes a local geometric feature region formed due to the mutation of the terrain or the ground object.
4. The above-and-underground integrated 3D modeling penetration prevention processing method according to claim 1, characterized in that: When S03 is performed, the calculation formula of the curvature feature is as follows:
5. The above-and-underground integrated 3D modeling penetration prevention processing method according to claim 4, characterized in that: The local geometric feature region is a curvature feature or a region where a gradient feature G exceeds a set threshold.
6. The above-and-underground integrated 3D modeling penetration prevention processing method according to claim 1, characterized in that: The calculation formula of the gradient feature is as follows: ； In the formula, represents the curvature characteristic; is a second-order partial derivative, representing the change speed of the east-west direction slope, The greater the value of, the greater the degree of curvature in the east-west direction; is a second-order partial derivative, representing the change speed of the north-south direction slope, The greater the value of, the greater the degree of curvature in the north-south direction; The calculation formula of the feature region fusion is as follows: ； In the formula, G represents a gradient feature; is a first-order partial derivative, representing a change rate of the height H in the east-west direction, The greater the absolute value of G is, the steeper the terrain is in the east-west direction. is the first-order partial derivative, representing the rate of change of the elevation H in the north-south direction, The greater the absolute value of the first-order partial derivative, the steeper the slope in the north-south direction. When S04 is performed, the calculation formula of the initial underground elevation is as follows: ； wherein and is a normalization factor, F is a feature region identification coefficient, F ranges from 0-1, with closer to 1 indicating a high risk feature region and closer to 0 indicating a low risk flat region.
7. The above-and-underground integrated 3D modeling penetration prevention processing method according to claim 6, characterized in that: The calculation formula of the feature region adjustment is as follows: ； When S05 is performed, the smoothing weight calculation formula is as follows: ； wherein, H underground is the initial subsurface elevation, is the adjusted subsurface elevation for the feature area, H min is the minimum ground surface elevation in the neighborhood, is a safety threshold.
8. The above-and-underground integrated 3D modeling penetration prevention processing method according to claim 6, characterized in that: The anisotropic diffusion equation is as follows: ； In the feature region, that is, in the region where the feature region identification coefficient F>0.6, the smoothing weight is small, so as to maintain the terrain details; in the flat region, that is, in the region where F<0.3, the smoothing weight is large, so as to enhance the model smoothness. ； In the formula, is the first-order partial derivative, indicating the rate of change of anisotropic diffusion, s represents the smoothing weight, and F is a feature region identification coefficient, represents a vector differential operator, which is used to calculate the gradient of a scalar field or the divergence of a vector field; when acts on the elevation scalar field H, that is, a gradient operation is formed; represents the gradient of the elevation field H, which is a vector, the direction of which points to the direction of the most severe change of the surface elevation at the point, and the modulus length thereof represents the rate of change of the elevation at the point, the modulus of is the numerical value of the G, which is the geometric basis for identifying the steepness of the terrain and the feature region.
9. The above-and-underground integrated three-dimensional modeling penetration prevention processing method according to claim 8, characterized by: The small smoothing weight refers to the case that s<0.58, and the large smoothing weight refers to the case that s>0.
79.
10. The above-and-underground integrated three-dimensional modeling penetration prevention processing method of claim 9, wherein: