CN111415415A - Automatic construction method of three-dimensional map cutting geological profile - Google Patents

Automatic construction method of three-dimensional map cutting geological profile Download PDF

Info

Publication number
CN111415415A
CN111415415A CN202010242546.3A CN202010242546A CN111415415A CN 111415415 A CN111415415 A CN 111415415A CN 202010242546 A CN202010242546 A CN 202010242546A CN 111415415 A CN111415415 A CN 111415415A
Authority
CN
China
Prior art keywords
dimensional
point
stratum
line
geological
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN202010242546.3A
Other languages
Chinese (zh)
Other versions
CN111415415B (en
Inventor
沈言根
李安波
黄键初
闾国年
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Nanjing Normal University
Original Assignee
Nanjing Normal University
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Nanjing Normal University filed Critical Nanjing Normal University
Priority to CN202010242546.3A priority Critical patent/CN111415415B/en
Publication of CN111415415A publication Critical patent/CN111415415A/en
Application granted granted Critical
Publication of CN111415415B publication Critical patent/CN111415415B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T17/00—Three dimensional [3D] modelling, e.g. data description of 3D objects
    • G06T17/05—Geographic models

Landscapes

  • Physics & Mathematics (AREA)
  • Engineering & Computer Science (AREA)
  • Geometry (AREA)
  • Software Systems (AREA)
  • Remote Sensing (AREA)
  • Computer Graphics (AREA)
  • General Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Geophysics And Detection Of Objects (AREA)
  • Processing Or Creating Images (AREA)

Abstract

The invention discloses an automatic construction method of a three-dimensional map cutting geological section, which comprises the steps of (1) reading geological section line data, a geological map and a DEM to form a section line set P L, a geological boundary line set Geo L ine, a ground plane set GeoPolygon and a grid data set GeoDEM, (2) obtaining any section line from the set P L, obtaining intersection points of the section line and all geological boundary lines based on the set Geo L ine to form a three-dimensional ground intersection point set PP, (3) obtaining a stratum attitude set AT and a stratum attribute set MK of a current section line based on the GeoPolygon, (4) obtaining bottom intersection points of a section bottom side line corresponding to the current section line and all ground lines based on the set GeoDEM and the set PP to form a three-dimensional bottom intersection point set EP, (5) constructing the three-dimensional map cutting geological section based on the PP and the set EP, (6) circulating the steps (2) - (5) to obtain all three-dimensional map cutting sections.

Description

Automatic construction method of three-dimensional map cutting geological profile
Technical Field
The invention relates to the field of three-dimensional modeling and geology, in particular to an automatic construction method of a three-dimensional map cutting geological profile.
Background
The map-cutting geological profile is a geological profile which is compiled by selecting a certain direction on a geological map and using a projection method according to various geological and geographic elements and a certain scale. The two-dimensional map cutting section map is matched with the geological plane map, and the working experience and the space imagination of geological experts are combined, so that the spreading rule and the penetration and cutting relation of the geological body and the geological structure on the three-dimensional space can be recognized by human beings to a certain extent.
However, the expression mode of geological information is too abstract, which is not beneficial to information acquisition and rule cognition of geological space, and is difficult to meet the requirements of three-dimensional geological modeling and expression. The three-dimensional map cutting geological profile can support three-dimensional geological information expression to a certain degree while keeping the characteristics of the traditional two-dimensional map cutting profile, and can provide a data base for profile-based three-dimensional geological modeling application. Therefore, the method for automatically constructing the three-dimensional map cutting geological profile is developed, and has important practical value and research significance.
Disclosure of Invention
The purpose of the invention is as follows: aiming at the problems in the prior art, the invention provides an automatic construction method of a three-dimensional map cutting geological profile, which has high automation degree and good expression effect.
The technical scheme is as follows: the automatic construction method of the three-dimensional map cutting geological profile comprises the following steps:
(1) reading the map-cut geological section line data, the geological map and the DEM to form a section line set P L, a geological boundary set Geo L ine, a stratum set GeoPolygon and a grid data set GeoDEM;
(2) acquiring any section line from the section line set P L, and acquiring the intersection points of the section line and all geological boundary lines based on the geological boundary line set Geo L ine to form a three-dimensional earth surface intersection point set PP;
(3) acquiring a stratum attitude set AT and a stratum attribute set MK of the current section line based on the stratum surface set GeoPolygon;
(4) based on the raster data set GeoDEM and the three-dimensional earth surface intersection point set PP, acquiring bottom intersection points of the profile bottom side line corresponding to the current section line and all the stratum lines to form a three-dimensional bottom intersection point set EP;
(5) constructing a three-dimensional map cutting geological profile based on the three-dimensional surface intersection point set PP and the three-dimensional bottom intersection point set EP;
(6) and (5) circulating the steps (2) to (5) until all section lines are processed, and obtaining all three-dimensional map cutting geological profiles.
Further, the step (1) specifically comprises:
(1-1) reading and cutting geological section line data, and storing the section line data in a section line set P L ═ pli1,2, …, L I } where pliI represents the ith section line, L I represents the number of the section lines;
(1-2) extracting geological boundary data and ground level data from the geological map, and respectively storing the geological boundary data and the ground level data into a geological boundary set Geo L ine and a ground level set geoPolygon;
and (1-3) extracting raster data from the DEM and storing the raster data into a raster data set GeoDEM.
Further, the step (2) specifically comprises:
(2-1) obtaining any one of the section lines pliReading all points on the section line generates a sectional segmented line end point set TP ═ TPf(xf,yf) 1,2, …, PF, where tpfPoint f on the cross-sectional line, PF represents the number of endpoints;
(2-2) for all the end points in the set TP, the azimuth angle of the segment line formed by any two adjacent end points is calculated according to the following formula, and the set AG of azimuth angles is stored as { β {s|s＝1,2,…,PF-1}：
Figure BDA0002433032480000021
in the formula ,(xa,ya)、(xb,yb) Coordinates representing two adjacent end points, respectively, βsRepresents the azimuth of the segment line s;
(2-3) extraction of the hatching pIiIntersection points with all geological boundary lines in the geological boundary line set Geo L ine, and storing the stratum attitude of the corresponding position of each intersection point as the attribute into an intersection point set InPoint;
(2-4) the intersection point set InPoint is defined as a three-dimensional surface intersection point set PP ═ PPj1,2, …, PJ }, where ppjRepresents the jth intersectionThe point number and PJ represent the number of the intersection points in the three-dimensional surface intersection point set PP.
Further, the step (3) specifically comprises:
(3-1) acquiring stratum attitude of all points in the three-dimensional surface intersection point set PP to obtain an attitude set AT ═ αn(ρ,θ,)|n＝1,2,…,AN}，αn(rho, theta) represents the n-th stratum attitude, rho is the stratum inclination, theta is the stratum inclination and is the stratum trend, and AN represents the number of attitude;
(3-2) obtaining the Current section line pliAll segmented lines of (a), forming a line set Se L ine;
(3-3) storing the stratum unique number GID of the stratum surface set GeoPolygon as the attribute of the corresponding line in the line set Se L ine, and establishing the association between GeoPolygon and Se L ine;
and (3-4) acquiring the GIDs of all the lines from Se L ine and storing the GIDs in the stratum attribute set MK.
Further, the step (4) specifically comprises:
(4-1) at the current section line pliSampling at intervals according to a preset distance d to obtain a sampling point set IP (IP) { IP }m1,2, …, PM, where ipmThe m point of sampling, PM is the number of sampling points;
(4-2) according to the section line pliThe run of (c) is sequenced, the section line pl isiThe starting point, the end point, the three-dimensional earth surface intersection point set PP and the sampling point set IP are merged to form an earth surface point set
Figure BDA0002433032480000031
wherein ,
Figure BDA0002433032480000032
representing gp at point kkPJ is the number of intersection points in the three-dimensional earth surface intersection point set PP;
(4-3) endowing all points in the earth surface point set GP with elevation values according to the raster data set GeoDEM;
(4-4) calculating the section line pl in a two-dimensional coordinate systemiForming a set ZV by the vertical coordinate of the intersection point of the bottom edge line of the corresponding section and all the stratum lines;
(4-5) acquiring three-dimensional bottom intersection points of all the stratum lines and the bottom side line of the section based on the three-dimensional ground surface intersection point set PP, the azimuth angle set AG of the section line and the set ZV to form a three-dimensional bottom intersection point set
Figure BDA0002433032480000033
Indicating the coordinates of the jth point in the EP.
Further, the step (4-4) specifically comprises:
(4-4-1) converting the surface point set GP from three-dimensional coordinates to two-dimensional planar coordinates using the following formula:
Figure BDA0002433032480000034
in the formula ,
Figure BDA0002433032480000035
representing GP of k point in the surface point set GPkIs determined by the three-dimensional coordinates of (a),
Figure BDA0002433032480000036
denotes gpkTwo-dimensional coordinates of (a);
(4-4-2) adopting the following formula to downwards move all the points in the surface point set GP by a preset distance h to form a point set of a bottom side line
Figure BDA0002433032480000037
Figure BDA0002433032480000038
(4-4-3) calculating the apparent inclination angle from the attitude set AT by using the following formula
Figure BDA0002433032480000039
Figure BDA00024330324800000310
Figure BDA00024330324800000311
Figure BDA00024330324800000312
Where ρ is the stratigraphic inclination, θ is the stratigraphic dip,
Figure BDA0002433032480000041
is the azimuth of the corresponding segment line;
(4-4-5) two-dimensional plane coordinates based on the intersection of the Earth's surface and the corresponding apparent inclination
Figure BDA0002433032480000042
Calculating to obtain the coordinates of the bottom points of the stratum line by adopting the following formula to form a set of the bottom points of the stratum line
Figure BDA0002433032480000043
Figure BDA0002433032480000044
in the formula ,
Figure BDA0002433032480000045
represents the jth bottom point of the horizon
Figure BDA0002433032480000046
The coordinates of the position of the object to be imaged,
Figure BDA0002433032480000047
representing the jth point pp in the set of surface intersectionsjPJ represents the number of points in the surface intersection set, le represents the point dpjAnd ppjThe preset distance of (a);
(4-4-6) based on Point dpjAnd ppjConstructing a formation line sljForming a set of layer lines S L ═ { sl ═j|j＝1,2,…,PJ}；
(4-4-7) creating a profile bottom edge line cl based on the two-dimensional coordinates of the first point of the surface point set GP, the point set BP of the bottom edge line and the two-dimensional coordinates of the last point in the surface point set GP;
(4-4-8) the set of formation lines S L is traversed to obtain the ordinate value of the intersection of each formation line and the section bottom edge cl, and the ordinate value is stored in the set ZV ═ ZVj1,2, …, PJ.
Further, the step (4-5) specifically comprises:
(4-5-1) calculating the relative elevations of each point in the surface intersection point set and the corresponding point in the three-dimensional bottom intersection point set according to the following formula:
Figure BDA0002433032480000048
in the formula ,hejRepresenting the jth point pp in the set of surface intersectionsjAnd the corresponding point epjThe relative elevation of the light source (c),
Figure BDA0002433032480000049
denotes ppjZ-axis coordinate of (1), zvjRepresents the jth value in the set ZV;
(4-5-2) calculating three-dimensional coordinates of the three-dimensional bottom intersection point according to the following equation
Figure BDA00024330324800000410
Figure BDA0002433032480000051
Further, the step (5) specifically comprises:
(5-1) moving down the surface point set GP by a preset distance h according to the following formula to obtain a three-dimensional bottom point set
Figure BDA0002433032480000052
Figure BDA0002433032480000053
(5-2) constructing a MultiPatch polyhedral element based on the three-dimensional surface intersection point set PP and the three-dimensional bottom intersection point set EP, and generating a three-dimensional cutting geological profile set SectionMap;
(5-3) creating an ID field in the SectionMap attribute table, and storing stratum attribute values corresponding to the three-dimensional stratum surface elements in the ID field;
and (5-4) transferring other attribute information of the GeoPolygon into an attribute table of the SectionMap according to the peer relationship between the ID field in the SectionMap and the GID field in the GeoPolygon in the stratum level set.
Further, the step (5-2) specifically comprises:
(5-2-1) respectively inserting the first point and the last point in the surface point set GP to the first position and the last position in the three-dimensional surface intersection point set PP;
(5-2-2) respectively inserting the first point and the last point in the three-dimensional bottom point set WP into the first position and the last position in the three-dimensional bottom intersection point set EP;
(5-2-3) acquiring any two adjacent points PP in the three-dimensional earth surface intersection point set PP in the earth surface point set GPjAnd ppj+1All upper points in between, constitute an upper point set UP ═ { UP ═ UPl|l＝1,2,…,PL}；
(5-2-4) acquiring any two adjacent points EP in the three-dimensional bottom intersection point set EP in the three-dimensional bottom point set WPjAnd epj+1All the lower points in between constitute the lower point set XP ═ { XP ═l|l＝1,2,…,PU}；
(5-2-5) arranging UP and XP in a clockwise order to form a stratum vertex set V ═ Vr|r＝1，2，…，PL+PU}；
(5-2-6) sequentially constructing MultiPatch polyhedral elements in a clockwise order based on the stratum-level vertex set V to form a three-dimensional stratum plane few；
(5-2-7) performing (5-2-3) - (5-2-6) circularly to obtain a three-dimensional stratum plane set FE ═ FEw|w＝1，2，…，PJ+1}；
(5-2-8) tessellating each three-dimensional horizon FE in three-dimensional horizon set FEwAnd constructing a single three-dimensional map cutting geological profile to obtain a three-dimensional map cutting geological profile set SectionMap.
Has the advantages that: the invention improves the three-dimensional expression effect of the map-cutting geological profile and has higher automation degree through links such as data loading, intersection point calculation, stratum attitude and attribute acquisition, three-dimensional bottom intersection point calculation, three-dimensional map-cutting profile construction and the like.
Drawings
FIG. 1 is section line and geological map data used in the present example;
fig. 2 is DEM data employed in the present embodiment;
FIG. 3 is a flow chart of an embodiment of the present invention;
FIG. 4 is a schematic two-dimensional cross-sectional view provided by the present invention;
FIG. 5 is a three-dimensional slice geological profile collection constructed in accordance with the present invention;
FIG. 6 is a diagram of the relative position of a collection of geological profiles and a DEM provided in three-dimensional view in accordance with the present invention.
Detailed Description
To explain the technical solution of the present invention in further detail, the experimental data of this embodiment uses the geological data of the tangram of Nanjing city (fig. 1) and the DEM data of the tangram of Nanjing city (fig. 2), and the projection coordinate system used in the experimental data is WGS _1984_ UTM _ Zone _ 50N. The following further description is provided by describing a specific embodiment in conjunction with the accompanying drawings.
As shown in fig. 3, the embodiment provides an automatic construction method of a three-dimensional map-cut geological section, which specifically includes the following steps:
(1) and reading the map cutting geological section line data, the geological map (including the stratum boundary lines and the stratum) and the DEM to form a section line set P L, a geological boundary line set Geo L ine, a stratum layer set GeoPolygon and a grid data set GeoDEM.
The method specifically comprises the following steps:
(1-1) reading and cutting geological section line data, and storing the section line data in a section line set P L ═ pli1,2, …, L I } where pliThe number of section lines is represented by the ith section line L I, and in the embodiment, L I is 3;
(1-2) extracting geological boundary data and ground level data from the geological map, and respectively storing the geological boundary data and the ground level data into a geological boundary set Geo L ine and a ground level set geoPolygon;
and (1-3) extracting raster data from the DEM and storing the raster data into a raster data set GeoDEM.
(2) Any cross line is obtained from the cross line set P L, and the intersection points of the cross line and all geological boundaries are obtained based on the geological boundary set Geo L ine, so that a three-dimensional surface intersection point set PP is formed.
The method specifically comprises the following steps:
(2-1) obtaining any one of the section lines pliReading all points on the section line generates a sectional segmented line end point set TP ═ TPf(xf，yf) 1,2, …, PF, where tpfPoint f on the cross-sectional line, PF represents the number of endpoints;
(2-2) for all the end points in the set TP, the azimuth angle of the segment line formed by any two adjacent end points is calculated according to the following formula, and the set AG of azimuth angles is stored as { β {s|s＝1,2,…,PF-1}：
Figure BDA0002433032480000071
in the formula ,(xa,ya)、(xb,yb) Coordinates representing two adjacent end points, respectively, βsRepresents the azimuth of the segment line s; in the present embodiment, when i is 1, PF is 2, xa＝153224.48,ya＝150759.99,xb＝153193.79,,yb＝150579.27,β＝3.31。
(2-3) extraction of the hatching line pl based on the Arcgis Engine APIiIntersection points with all geological boundary lines in the geological boundary line set Geo L ine, and storing the stratum attitude of the corresponding position of each intersection point as the attribute into an intersection point set InPoint;
(2-4) the intersection point set InPoint is defined as a three-dimensional surface intersection point set PP ═ PPj1,2, …, PJ }, where ppjThe j-th intersection number is shown, PJ represents the number of intersections in the three-dimensional surface intersection set PP, and when i is 1, PJ is 6.
(3) And acquiring a stratum attitude set AT and a stratum attribute set MK of the current section line based on the stratum surface set GeoPolygon.
The method specifically comprises the following steps:
(3-1) acquiring stratum attitude of all points in the three-dimensional surface intersection point set PP to obtain an attitude set AT ═ αn(ρ,θ，)|n＝1，2，…，AN}，αn(rho, theta) represents the n-th stratum attitude, rho is the stratum inclination, theta is the stratum inclination and is the stratum trend, and AN represents the number of attitude; when i is 1, AN is 6;
(3-2) calling Arcgis Engine API to obtain the current section line pliAll segmented lines of (a), forming a line set Se L ine;
(3-3) storing the stratum unique codes GID of the stratum surface set GeoPolygon as the attributes of corresponding lines in a line set Se L ine, and establishing the association between GeoPolygon and Se L ine;
(3-4) the GID of all lines is obtained from Se L ine and stored in the stratum attribute set MK.
(4) And acquiring bottom intersections of the profile bottom side lines corresponding to the current section lines and all the stratum lines based on the raster data set GeoDEM and the three-dimensional ground surface intersection point set PP to form a three-dimensional bottom intersection point set EP.
The method specifically comprises the following steps:
(4-1) at the current section line pliSampling at intervals according to a preset distance d to obtain a sampling point set IP (IP) { IP }m1,2, …, PM, where ipmThe m point of sampling, PM is the number of sampling points; in the embodiment, d is 4, and when i is 1, PM is 158;
(4-2) according to the section line pliThe run of (c) is sequenced, the section line pl isiThe starting point, the end point, the three-dimensional earth surface intersection point set PP and the sampling point set IP are merged to form an earth surface point set
Figure BDA0002433032480000081
wherein ,
Figure BDA0002433032480000082
representing gp at point kkIs PJ is three-dimensionalThe number of intersection points in the table intersection point set PP; when i is 1 in the present embodiment, PM + PJ +2 is 165;
(4-3) loading GeoDEM, calling an Arcgis Engine API, and endowing all points in the surface point set GP with elevation values according to the raster data set GeoDEM;
(4-4) calculating the section line pl in a two-dimensional coordinate systemiForming a set ZV by the vertical coordinate of the intersection point of the bottom edge line of the corresponding section and all the stratum lines;
(4-5) acquiring three-dimensional bottom intersection points of all the stratum lines and the bottom side line of the section based on the three-dimensional ground surface intersection point set PP, the azimuth angle set AG of the section line and the set ZV to form a three-dimensional bottom intersection point set
Figure BDA0002433032480000083
Indicating the coordinates of the jth point in the EP.
Wherein, the step (4-4) specifically comprises the following steps:
(4-4-1) converting the surface point set GP from three-dimensional coordinates to two-dimensional planar coordinates using the following formula:
Figure BDA0002433032480000084
in the formula ,
Figure BDA0002433032480000085
representing GP of k point in the surface point set GPkIs determined by the three-dimensional coordinates of (a),
Figure BDA0002433032480000086
denotes gpkTwo-dimensional coordinates of (a);
(4-4-2) adopting the following formula to downwards move all the points in the surface point set GP by a preset distance h to form a point set of a bottom side line
Figure BDA0002433032480000091
Wherein h is 20;
Figure BDA0002433032480000092
(4-4-3) calculating the apparent inclination angle from the attitude set AT by using the following formula
Figure BDA0002433032480000093
Figure BDA0002433032480000094
Figure BDA0002433032480000095
Figure BDA0002433032480000096
Where ρ is the stratigraphic inclination, θ is the stratigraphic dip,
Figure BDA0002433032480000097
is the azimuth of the corresponding segment line;
(4-4-5) two-dimensional plane coordinates based on the intersection of the Earth's surface and the corresponding apparent inclination
Figure BDA0002433032480000098
Calculating to obtain the coordinates of the bottom points of the stratum line by adopting the following formula to form a set of the bottom points of the stratum line
Figure BDA0002433032480000099
Figure BDA00024330324800000910
in the formula ,
Figure BDA00024330324800000911
represents the jth bottom point dp of the horizonjThe coordinates of the position of the object to be imaged,
Figure BDA00024330324800000912
representing the jth point pp in the set of surface intersectionsjTwo-dimensional plane coordinates of (1), PJ representing a set of surface intersectionsThe number of midpoints, le, represents a point dpjAnd ppjThe preset distance of (a) is 20;
(4-4-6) based on Point dpjAnd ppjConstructing a formation line sljForming a set of layer lines S L ═ { sl ═j|j＝1,2,…,PJ}；
(4-4-7) calling an Arcgis Engine API, and creating a profile bottom edge line cl based on the two-dimensional coordinates of the first point of the surface point set GP, the point set BP of the bottom edge line and the two-dimensional coordinates of the last point in the surface point set GP;
(4-4-8) the set of formation lines S L is traversed to obtain the ordinate value of the intersection of each formation line and the section bottom edge cl, and the ordinate value is stored in the set ZV ═ ZVj1,2, …, PJ. In the present embodiment, a schematic diagram of the ordinate of the bottom intersection point is shown in fig. 4.
Wherein, the step (4-5) specifically comprises:
(4-5-1) calculating the relative elevations of each point in the surface intersection point set and the corresponding point in the three-dimensional bottom intersection point set according to the following formula:
Figure BDA0002433032480000101
in the formula ,hejRepresenting the jth point pp in the set of surface intersectionsjAnd the corresponding point epjThe relative elevation of the light source (c),
Figure BDA0002433032480000102
denotes ppjZ-axis coordinate of (1), zvjRepresents the jth value in the set ZV;
(4-5-2) calculating three-dimensional coordinates of the three-dimensional bottom intersection point according to the following equation
Figure BDA0002433032480000103
Figure BDA0002433032480000104
(5) And constructing a three-dimensional map cutting geological profile based on the three-dimensional surface intersection point set PP and the three-dimensional bottom intersection point set EP.
The method specifically comprises the following steps:
(5-1) moving down the surface point set GP by a preset distance h according to the following formula to obtain a three-dimensional bottom point set
Figure BDA0002433032480000105
Figure BDA0002433032480000106
(5-2) constructing a MultiPatch polyhedral element according to the three-dimensional surface intersection point set PP and the three-dimensional bottom intersection point set EP based on the Arcgis Engine API, and generating a three-dimensional cutting geological section set section map;
(5-3) creating an ID field in the SectionMap attribute table, and storing stratum attribute values corresponding to the three-dimensional stratum surface elements in the ID field;
and (5-4) transferring other attribute information of the GeoPolygon into an attribute table of the SectionMap according to the peer relationship between the ID field in the SectionMap and the GID field in the GeoPolygon in the stratum level set.
Wherein, the step (5-2) specifically comprises the following steps:
(5-2-1) respectively inserting the first point and the last point in the surface point set GP to the first position and the last position in the three-dimensional surface intersection point set PP;
(5-2-2) respectively inserting the first point and the last point in the three-dimensional bottom point set WP into the first position and the last position in the three-dimensional bottom intersection point set EP;
(5-2-3) acquiring any two adjacent points PP in the three-dimensional earth surface intersection point set PP in the earth surface point set GPjAnd ppj+1All upper points in between, constitute an upper point set UP ═ { UP ═ UPl1,2, …, P L }, in this embodiment, when j is 1, P L is 21;
(5-2-4) acquiring any two adjacent points EP in the three-dimensional bottom intersection point set EP in the three-dimensional bottom point set WPjAnd epj+1All the lower points in between constitute the lower point set XP ═ { XP ═l1,2, …, PU }; in the present embodiment, when j is 1, PU is 4;
(5-2-5) arranging UP and XP in a clockwise order to form a stratum vertex set V ═ Vr|r＝1,2，…，PL+PU}；
(5-2-6) sequentially constructing MultiPatch polyhedral elements in a clockwise order based on the stratum-level vertex set V to form a three-dimensional stratum plane few；
(5-2-7) performing (5-2-3) - (5-2-6) circularly to obtain a three-dimensional stratum plane set FE ═ FEw|w＝1,2,…,PJ+1}；
(5-2-8) tessellating each three-dimensional horizon FE in three-dimensional horizon set FEwAnd constructing a single three-dimensional map cutting geological profile to obtain a three-dimensional map cutting geological profile set SectionMap.
(6) And (5) circulating the steps (2) to (5) until all section lines are processed, and obtaining all three-dimensional map cutting geological profiles. In this embodiment, a set of three-dimensional slice geological profiles is constructed as shown in FIG. 5.
In this embodiment, the relative positions of the three-dimensional sliced geological profile collection and the DEM are as shown in fig. 6. After the construction is finished, the three-dimensional map cutting geological profile model in the ESRI TIN format can be exported to be a three-dimensional model file in the obj, FBX and other formats based on the Arcgis Engine API.
In the embodiment of the invention, partial GIS operation is provided based on the Arcgis Engine API, and related steps can also use the APIs of software such as SuperMap, Arcgis Object and the like to carry out corresponding GIS operation.

Claims (9)

1. An automatic construction method of a three-dimensional cutting geological section is characterized by comprising the following steps:
(1) reading the map-cut geological section line data, the geological map and the DEM to form a section line set P L, a geological boundary set Geo L ine, a stratum set GeoPolygon and a grid data set GeoDEM;
(2) acquiring any section line from the section line set P L, and acquiring the intersection points of the section line and all geological boundary lines based on the geological boundary line set Geo L ine to form a three-dimensional earth surface intersection point set PP;
(3) acquiring a stratum attitude set AT and a stratum attribute set MK of the current section line based on the stratum surface set GeoPolygon;
(4) based on the raster data set GeoDEM and the three-dimensional earth surface intersection point set PP, acquiring bottom intersection points of the profile bottom side line corresponding to the current section line and all the stratum lines to form a three-dimensional bottom intersection point set EP;
(5) constructing a three-dimensional map cutting geological profile based on the three-dimensional surface intersection point set PP and the three-dimensional bottom intersection point set EP;
(6) and (5) circulating the steps (2) to (5) until all section lines are processed, and obtaining all three-dimensional map cutting geological profiles.
2. The method of automatically constructing a three-dimensional sliced geological profile of claim 1, wherein: the step (1) specifically comprises the following steps:
(1-1) reading and cutting geological section line data, and storing the section line data in a section line set P L ═ pli1,2, …, L I } where pliI represents the ith section line, L I represents the number of the section lines;
(1-2) extracting geological boundary data and ground level data from the geological map, and respectively storing the geological boundary data and the ground level data into a geological boundary set Geo L ine and a ground level set geoPolygon;
and (1-3) extracting raster data from the DEM and storing the raster data into a raster data set GeoDEM.
3. The method of automatically constructing a three-dimensional sliced geological profile of claim 1, wherein: the step (2) specifically comprises the following steps:
(2-1) obtaining any one of the section lines pliReading all points on the section line generates a sectional segmented line end point set TP ═ TPf(xf,yf) 1,2, …, PF, where tpfPoint f on the cross-sectional line, PF represents the number of endpoints;
(2-2) for all the end points in the set TP, the azimuth angle of the segment line formed by any two adjacent end points is calculated according to the following formula, and the set AG of azimuth angles is stored as { β {s|s＝1,2,…,PF-1}：
Figure FDA0002433032470000011
in the formula ,(xa,ya)、(xb,yb) Coordinates representing two adjacent end points, respectively, βsRepresents the azimuth of the segment line s;
(2-3) extraction of the hatching pIiIntersection points with all geological boundary lines in the geological boundary line set Geo L ine, and storing the stratum attitude of the corresponding position of each intersection point as the attribute into an intersection point set InPoint;
(2-4) the intersection point set InPoint is defined as a three-dimensional surface intersection point set PP ═ PPj1,2, …, PJ }, where ppjThe j-th intersection point sequence number is shown, and PJ shows the number of intersection points in the three-dimensional surface intersection point set PP.
4. The method of automatically constructing a three-dimensional sliced geological profile of claim 1, wherein: the step (3) specifically comprises the following steps:
(3-1) acquiring stratum attitude of all points in the three-dimensional surface intersection point set PP to obtain an attitude set AT ═ αn(ρ，θ，)|n＝1，2，…，AN}，αn(rho, theta) represents the n-th stratum attitude, rho is the stratum inclination, theta is the stratum inclination and is the stratum trend, and AN represents the number of attitude;
(3-2) obtaining the Current section line pliAll segmented lines of (a), forming a line set Se L ine;
(3-3) storing the stratum unique number GID of the stratum surface set GeoPolygon as the attribute of the corresponding line in the line set Se L ine, and establishing the association between GeoPolygon and Se L ine;
and (3-4) acquiring the GIDs of all the lines from Se L ine and storing the GIDs in the stratum attribute set MK.
5. The method of automatically constructing a three-dimensional sliced geological profile of claim 1, wherein: the step (4) specifically comprises the following steps:
(4-1) at the current section line pliSampling at intervals according to a preset distance d to obtain a sampling point set IP (IP) { IP }m1,2, …, PM, where ipmThe m point of sampling, PM is the number of sampling points;
(4-2) according to the section line pliThe run of (c) is sequenced, the section line pl isiThe starting point, the end point, the three-dimensional earth surface intersection point set PP and the sampling point set IP are merged to form an earth surface point set
Figure FDA0002433032470000021
wherein ,
Figure FDA0002433032470000022
representing gp at point kkPJ is the number of intersection points in the three-dimensional earth surface intersection point set PP;
(4-3) endowing all points in the earth surface point set GP with elevation values according to the raster data set GeoDEM;
(4-4) calculating the section line pl in a two-dimensional coordinate systemiForming a set ZV by the vertical coordinate of the intersection point of the bottom edge line of the corresponding section and all the stratum lines;
(4-5) acquiring three-dimensional bottom intersection points of all the stratum lines and the bottom side line of the section based on the three-dimensional ground surface intersection point set PP, the azimuth angle set AG of the section line and the set ZV to form a three-dimensional bottom intersection point set
Figure FDA0002433032470000023
Figure FDA0002433032470000024
Indicating the coordinates of the jth point in the EP.
6. The method of automatically constructing a three-dimensional sliced geological profile of claim 5, wherein: the step (4-4) specifically comprises the following steps:
(4-4-1) converting the surface point set GP from three-dimensional coordinates to two-dimensional planar coordinates using the following formula:
Figure FDA0002433032470000031
in the formula ,
Figure FDA0002433032470000032
representing GP of k point in the surface point set GPkIs determined by the three-dimensional coordinates of (a),
Figure FDA0002433032470000033
denotes gpkTwo-dimensional coordinates of (a);
(4-4-2) adopting the following formula to downwards move all the points in the surface point set GP by a preset distance h to form a point set of a bottom side line
Figure FDA0002433032470000034
Figure FDA0002433032470000035
(4-4-3) calculating the apparent inclination angle from the attitude set AT by using the following formula
Figure FDA0002433032470000036
Figure FDA0002433032470000037
Figure FDA0002433032470000038
Figure FDA0002433032470000039
Wherein ρ is a stratigraphic inclination, θ is a stratigraphic dip, β'sIs the azimuth of the corresponding segment line;
(4-4-5) two-dimensional plane coordinates based on the intersection of the Earth's surface and the corresponding apparent inclination
Figure FDA00024330324700000310
Calculating to obtain the coordinates of the bottom points of the stratum line by adopting the following formula to form a set of the bottom points of the stratum line
Figure FDA00024330324700000311
Figure FDA00024330324700000312
in the formula ,
Figure FDA00024330324700000313
represents the jth bottom point dp of the horizonjThe coordinates of the position of the object to be imaged,
Figure FDA00024330324700000314
representing the jth point pp in the set of surface intersectionsjPJ represents the number of points in the surface intersection set, le represents the point dpjAnd ppjThe preset distance of (a);
(4-4-6) based on Point dpjAnd ppjConstructing a formation line sljForming a set of layer lines S L ═ { sl ═j|j＝1,2，…,PJ}；
(4-4-7) creating a profile bottom edge line cl based on the two-dimensional coordinates of the first point of the surface point set GP, the point set BP of the bottom edge line and the two-dimensional coordinates of the last point in the surface point set GP;
(4-4-8) the set of formation lines S L is traversed to obtain the ordinate value of the intersection of each formation line and the section bottom edge cl, and the ordinate value is stored in the set ZV ═ ZVj1,2, …, PJ.
7. The method of automatically constructing a three-dimensional sliced geological profile of claim 5, wherein: the step (4-5) specifically comprises the following steps:
(4-5-1) calculating the relative elevations of each point in the surface intersection point set and the corresponding point in the three-dimensional bottom intersection point set according to the following formula:
Figure FDA0002433032470000041
in the formula ,hejRepresenting sets of surface intersectionsThe j point pp in the mergerjAnd the corresponding point epjThe relative elevation of the light source (c),
Figure FDA0002433032470000042
denotes ppjZ-axis coordinate of (1), zvjRepresents the jth value in the set ZV;
(4-5-2) calculating three-dimensional coordinates of the three-dimensional bottom intersection point according to the following equation
Figure FDA0002433032470000043
Figure FDA0002433032470000044
8. The method of automatically constructing a three-dimensional sliced geological profile of claim 1, wherein: the step (5) specifically comprises the following steps:
(5-1) moving down the surface point set GP by a preset distance h according to the following formula to obtain a three-dimensional bottom point set
Figure FDA0002433032470000045
Figure FDA0002433032470000046
(5-2) constructing a MultiPatch polyhedral element based on the three-dimensional surface intersection point set PP and the three-dimensional bottom intersection point set EP, and generating a three-dimensional cutting geological profile set SectionMap;
(5-3) creating an ID field in the SectionMap attribute table, and storing stratum attribute values corresponding to the three-dimensional stratum surface elements in the ID field;
and (5-4) transferring other attribute information of the GeoPolygon into an attribute table of the SectionMap according to the peer relationship between the ID field in the SectionMap and the GID field in the GeoPolygon in the stratum level set.
9. The method of automatically constructing a three-dimensional sliced geological profile of claim 1, wherein: the step (5-2) specifically comprises the following steps:
(5-2-1) respectively inserting the first point and the last point in the surface point set GP to the first position and the last position in the three-dimensional surface intersection point set PP;
(5-2-2) respectively inserting the first point and the last point in the three-dimensional bottom point set WP into the first position and the last position in the three-dimensional bottom intersection point set EP;
(5-2-3) acquiring any two adjacent points PP in the three-dimensional earth surface intersection point set PP in the earth surface point set GPjAnd ppj+1All upper points in between, constitute an upper point set UP ═ { UP ═ UPl|l＝1，2，…，PL}；
(5-2-4) acquiring any two adjacent points EP in the three-dimensional bottom intersection point set EP in the three-dimensional bottom point set WPjAnd epj+1All the lower points in between constitute the lower point set XP ═ { XP ═l|l＝1，2，…，PU}；
(5-2-5) arranging UP and XP in a clockwise order to form a stratum vertex set V ═ Vr|r＝1，2，…，PL+PU}；
(5-2-6) sequentially constructing MultiPatch polyhedral elements in a clockwise order based on the stratum-level vertex set V to form a three-dimensional stratum plane few；
(5-2-7) performing (5-2-3) - (5-2-6) circularly to obtain a three-dimensional stratum plane set FE ═ FEw|w＝1，2，…，PJ+1}；
(5-2-8) tessellating each three-dimensional horizon FE in three-dimensional horizon set FEwAnd constructing a single three-dimensional map cutting geological profile to obtain a three-dimensional map cutting geological profile set SectionMap.
CN202010242546.3A 2020-03-31 2020-03-31 Automatic construction method of three-dimensional map cutting geological section Active CN111415415B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202010242546.3A CN111415415B (en) 2020-03-31 2020-03-31 Automatic construction method of three-dimensional map cutting geological section

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202010242546.3A CN111415415B (en) 2020-03-31 2020-03-31 Automatic construction method of three-dimensional map cutting geological section

Publications (2)

Publication Number Publication Date
CN111415415A true CN111415415A (en) 2020-07-14
CN111415415B CN111415415B (en) 2023-05-09

Family

ID=71493482

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202010242546.3A Active CN111415415B (en) 2020-03-31 2020-03-31 Automatic construction method of three-dimensional map cutting geological section

Country Status (1)

Country Link
CN (1) CN111415415B (en)

Cited By (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111968231A (en) * 2020-08-14 2020-11-20 成都理工大学 Three-dimensional stratum modeling method based on geological map section
CN112233238A (en) * 2020-10-28 2021-01-15 南京师范大学 Bedrock geological body three-dimensional model construction method and device based on map-cut geological parallel section
CN112233230A (en) * 2020-09-08 2021-01-15 南京师范大学 Three-dimensional model construction method and device for fault structure in map-cut geological profile
CN112734926A (en) * 2021-01-29 2021-04-30 南京师范大学 Automatic generation method of map-cutting geological profile facing unconsolidated formation coverage area
CN113393578A (en) * 2021-06-02 2021-09-14 南京师范大学 Construction method and device of geological profile three-dimensional model

Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN105469443A (en) * 2014-09-30 2016-04-06 中国地质调查局发展研究中心 Method for generating three-dimensional geological map based on geological route (PRB) process double modeling
CN106934858A (en) * 2017-03-14 2017-07-07 中国地质科学院矿产资源研究所 Three-dimensional geological modeling method and system for scale region of mining area
CN109753707A (en) * 2018-12-25 2019-05-14 核工业北京地质研究院 A method of stratigraphic boundary, which is extracted, using section of exploration line carries out three-dimensional modeling

Patent Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN105469443A (en) * 2014-09-30 2016-04-06 中国地质调查局发展研究中心 Method for generating three-dimensional geological map based on geological route (PRB) process double modeling
CN106934858A (en) * 2017-03-14 2017-07-07 中国地质科学院矿产资源研究所 Three-dimensional geological modeling method and system for scale region of mining area
CN109753707A (en) * 2018-12-25 2019-05-14 核工业北京地质研究院 A method of stratigraphic boundary, which is extracted, using section of exploration line carries out three-dimensional modeling

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
周良辰 等: "平面地质图的三维地质体建模方法研究", 《地球信息科学学报》 *
高士娟 等: "基于平面地质图的地质体三维建模", 《地质找矿论丛》 *

Cited By (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111968231A (en) * 2020-08-14 2020-11-20 成都理工大学 Three-dimensional stratum modeling method based on geological map section
CN112233230A (en) * 2020-09-08 2021-01-15 南京师范大学 Three-dimensional model construction method and device for fault structure in map-cut geological profile
CN112233230B (en) * 2020-09-08 2024-02-27 南京师范大学 Three-dimensional model construction method and device for fault structure in cut geological section
CN112233238A (en) * 2020-10-28 2021-01-15 南京师范大学 Bedrock geological body three-dimensional model construction method and device based on map-cut geological parallel section
CN112233238B (en) * 2020-10-28 2024-02-27 南京师范大学 Bedrock geologic body three-dimensional model construction method and device based on graph cut geological parallel section
CN112734926A (en) * 2021-01-29 2021-04-30 南京师范大学 Automatic generation method of map-cutting geological profile facing unconsolidated formation coverage area
CN112734926B (en) * 2021-01-29 2024-04-19 南京师范大学 Automatic generation method of cut geological section for loose layer coverage area
CN113393578A (en) * 2021-06-02 2021-09-14 南京师范大学 Construction method and device of geological profile three-dimensional model

Also Published As

Publication number Publication date
CN111415415B (en) 2023-05-09

Similar Documents

Publication Publication Date Title
CN111415415A (en) Automatic construction method of three-dimensional map cutting geological profile
CN104574511B (en) A kind of quick progressive three-dimensional geological modeling method
CN106446910B (en) Complex geological curved surface feature extraction and reconstruction method
CN102385067B (en) Drawing method for isoline containing reverse fault
CN102521884A (en) Three-dimensional roof reconstruction method based on LiDAR data and ortho images
CN110032771B (en) DEM accurate cutting method considering local detail characteristics for opencast coal mine
CN103970837B (en) Discontinuous DEM classified manufacturing method based on urban land and vertical planning
CN109003330A (en) A kind of three dimensional contour line method based on basement rock boundary constraint
CN102013114B (en) Microstation v8i-based city rapid-modeling method
CN116152461A (en) Geological modeling method, device, computer equipment and computer readable storage medium
CN112434360B (en) Three-dimensional modeling and profile mapping method for complex rock-soil body
CN106875479A (en) A kind of automatic horizontally-placed method of digital elevation model
Englert et al. One-eye stereo system for the acquisition of complex 3D building descriptions
CN106875330B (en) Method for rotating plane model into spherical model
Taillandier Automatic building reconstruction from cadastral maps and aerial images
CN105205311B (en) A kind of complex number algorithm that stratigraphic section is drawn
CN110363844B (en) Coal mine tunnel three-dimensional modeling method and system
CN112233230A (en) Three-dimensional model construction method and device for fault structure in map-cut geological profile
CN110942510A (en) Three-dimensional rapid combination method for planar geological map
CN109858160A (en) A kind of modeling method of the rail traffic geological information model based on BIM technology
CN114463494A (en) Automatic topographic feature line extracting algorithm
Vassallo et al. Digging in excavation diaries: digital re-assessment of stratigraphy in 3D GIS. The sanctuary of Ayia Irini, Cyprus
CN113393578A (en) Construction method and device of geological profile three-dimensional model
CN114036608A (en) Method for converting geological profile into ultimate balance calculation model
CN108875755B (en) Automatic horizontal rock stratum extraction method based on curve parallel characteristics

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant