CN110910495A - Three-dimensional geometric form recovery method for dome structure - Google Patents

Three-dimensional geometric form recovery method for dome structure Download PDF

Info

Publication number
CN110910495A
CN110910495A CN201911042921.3A CN201911042921A CN110910495A CN 110910495 A CN110910495 A CN 110910495A CN 201911042921 A CN201911042921 A CN 201911042921A CN 110910495 A CN110910495 A CN 110910495A
Authority
CN
China
Prior art keywords
coordinates
intersection
section
geological boundary
point
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
CN201911042921.3A
Other languages
Chinese (zh)
Other versions
CN110910495B (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 Panzhi Geographic Information Industry Research Institute Co Ltd
Nanjing Normal University
Original Assignee
Nanjing Panzhi Geographic Information Industry Research Institute Co Ltd
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 Panzhi Geographic Information Industry Research Institute Co Ltd, Nanjing Normal University filed Critical Nanjing Panzhi Geographic Information Industry Research Institute Co Ltd
Priority to CN201911042921.3A priority Critical patent/CN110910495B/en
Publication of CN110910495A publication Critical patent/CN110910495A/en
Application granted granted Critical
Publication of CN110910495B publication Critical patent/CN110910495B/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)
  • Image Generation (AREA)

Abstract

The invention discloses a method for recovering a three-dimensional geometric form of a dome structure, which specifically comprises the following steps: (1) establishing a three-dimensional coordinate system by taking the centroid of any layer of geological boundary of the dome structure as an origin O, calculating the coordinates of all break points of the layer of geological boundary, and storing the coordinates into a point set G; (2) generating n sections which are perpendicular to the xOy plane and pass through the O point based on the equal included angle, and extracting intersection points of the sections and the geological boundary to obtain an intersection point set JP; (3) generating n second-order Bezier curves according to two intersection points of each section and the geological boundary and corresponding occurrence data at the intersection points, and calculating the coordinates of the top points of the domes; (4) for each intersection in the set JP, generating 2n Bezier curves; (5) interpolating between adjacent Bezier curves to generate a contour line set based on a Morphing method, and constructing an irregular triangular network between the adjacent contour lines; (6) and circularly executing the steps until all layers are recovered. The automatic restoration method effectively realizes the automatic restoration of the dome structure landform, and has the advantages of good restoration effect and high automation degree.

Description

Three-dimensional geometric form recovery method for dome structure
Technical Field
The invention relates to the field of geographic information technology application, in particular to a method for recovering a three-dimensional geometric form of a dome structure.
Background
The dome structure is a structural deformation pattern causing surface elevation, and the plane projection form of the dome structure usually shows an approximate circle or an ellipse, the diameter of which is from several kilometers to hundreds of kilometers, such as sheep mountains in the united states, Huang Ling anticline in the three gorges in China, and the like. The formation of the vault is a lot of factors, but is usually caused by the invasion of rock slurry or the superposition of two sets of anticline in the direction vertical to the hinge direction, and the result is that the originally horizontal sedimentary rock stratum on the upper part of the vault is subjected to bending deformation. These formations are susceptible to erosion and often form annular single-sided hills around their periphery, with steep slopes towards the center of the vault. On the geological map, the fornix core part exposes the older stratum, and the older stratum is gradually renewed outwards and distributed in a shape of nearly concentric circles. The vault structure deformation process is difficult to dissect through the plane geometric shape, so that the three-dimensional geometric shape before the vault is eroded is recovered, the understanding of the structure deformation process can be facilitated, the landform evolution information such as the denudation amount, the swelling rate and the like can be acquired, and the method has important academic value.
Disclosure of Invention
The purpose of the invention is as follows: the invention provides a method for recovering a three-dimensional geometric form of a dome structure, aiming at the problems in the prior art.
The technical scheme is as follows: the method for recovering the three-dimensional geometrical morphology of the dome structure comprises the following steps:
(1) establishing an O-xyz coordinate system in a three-dimensional space by taking the northward direction as an x axis and the eastern direction as a y axis, taking the centroid of any layer of geological boundary of a vault structure as an origin O and the upward direction perpendicular to the xOy plane as a z axis, calculating the coordinates of all break points of the layer of geological boundary, and storing the coordinates into a point set G;
(2) generating n sections which are perpendicular to the xOy plane and pass through the O points based on equal included angles according to the point set G, extracting 2n intersection points of the n sections and the geological boundary, and storing an intersection point set JP, wherein n is a positive integer;
(3) generating n second-order Bezier curves according to two intersection points of each section and the geological boundary and corresponding occurrence data at the intersection points, and calculating to obtain a dome vertex V coordinate;
(4) for each intersection point in the set JP, generating 2n Bezier curves based on the corresponding position dip angle, the dome vertex V and the tangent line at the vertex V;
(5) interpolation is carried out between adjacent Bezier curves to generate a contour line set based on a Morphing method, an irregular triangular net is constructed between adjacent contour lines, and the recovery of the current layer of the dome structure is completed;
(6) and (5) circularly executing the steps (1) to (5) until all the layers are recovered.
Further, the step (1) specifically comprises:
(1-1) taking the north direction as an x axis and the east direction as a y axis, selecting an arbitrary point with the same height as a horizontal plane as an origin O ', and taking the upward direction perpendicular to the xO ' y plane as a z axis, and establishing an O ' -xyz coordinate system in a three-dimensional space;
(1-2) extracting the coordinates of all break points of any layer of geological boundary of the dome structure to be recovered in an O' -xyz coordinate system according to the geological map and DEM data, and storing a point set A ═ p in sequencei(pxi,pyi,pzi) I ═ 1, …, m }, where m is the number of points in the geological boundary element, (pxi,pyi,pzi) Represents a point piX, y, z coordinate values in the O' -xyz coordinate system;
(1-3) establishing an O-xyz coordinate system in a three-dimensional space by taking the centroid of the geological boundary as an origin O, and calculating the coordinates g of all break points of the geological boundary of the layer in the O-xyz coordinate system according to the following formulai(xi,yi,zi) And storing the point set G ═ Gi(xi,yi,zi)|i＝1,…,m}：
Figure BDA0002253339720000021
Further, the step (2) specifically comprises:
(2-1) generating n sections which are perpendicular to the xOy plane and pass through the O point in the clockwise direction according to the principle that the included angle between any two adjacent sections is equal by taking y as a first section and storing a section set SP as { SP ═ SPt|t＝1,…,n}；
(2-2) dividing each of the generated sections along the z-axis to generateForming 2n half-sections, and adding a half-section set P { P } in a clockwise order from the northj1, …,2n, where the half-section PjThe included angle theta between the projection line on the xOy plane and the positive direction of the x-axisjComprises the following steps:
Figure BDA0002253339720000022
(2-3) calculating the origin O and the element G in the point set G based on the following formulaiFormed ray OgiAngle with positive direction of x-axisi；
Figure BDA0002253339720000023
In the formula, giRepresents the ith element in the point set G, m represents the number of elements in the point set G, xi,yiEach represents giCoordinate values;
(2-4) calculating the range of the included angle between the line segment formed by any two adjacent points in the point set G and the positive direction of the x axis based on the following formula:
Figure BDA0002253339720000031
in the formula, rangeiRepresents element G in point set Gi、gi+1Line segment L of compositioniRange of angle with positive direction of x-axis, rmini＝min(anglei,anglei+1),rmaxi＝max(anglei,anglei+1)，anglei+1Represents Ogi+1The included angle with the positive direction of the x axis;
(2-5) calculating coordinates of intersection points of each element in the half-section set P and the geological boundary based on the following formula:
Figure BDA0002253339720000032
in the formula, xj,yj,zjRepresenting an element P in a half-section set PjIntersection Jp with geological boundaryjIs equal to 1, 2n, i is e to [1, m ]]And satisfy rangeiIncluding thetaj,Ai＝yi+1-yi,Bi＝xi+1-xi,Ci＝xi+1yi-xiyi+1，xi,yi,ziDenotes giThree-dimensional spatial coordinate value of (2), xi+1,yi+1,zi+1Denotes gi+1Three-dimensional space coordinate values of (a);
(2-6) storing the intersection point coordinates of each element in the half-section set P and the geological boundary into an intersection point set JP ═ Jpj(xj,yj,zj) 1, …,2n, where the midpoint Jp is in the settAnd Jpt+nIs a profile SPtTwo intersection points with the geological boundary, t 1.
Further, the step (3) specifically comprises:
(3-1) transforming the coordinates of the intersection of each section plane and the geological boundary generated in the step (2) according to the intersection set JP based on the following formula to obtain the coordinates when the coordinates are mapped to the coordinate system x 'Oy' corresponding to the section plane:
Figure BDA0002253339720000033
in the formula, xt,yt,ztIndicating the t-th point Jp in the intersection set JPtIs also a section plane SPtThree-dimensional spatial coordinates, x, of the first intersection with the geological boundaryt+n,yt+n,zt+nRepresents the t + n point Jp in the intersection set JPt+nIs also a section plane SPtThree-dimensional spatial coordinates, x ', of a second intersection with the geological boundary't1,y′t1Denotes JptMapping to a profile SPtTwo-dimensional coordinates of (c), x't2,y′t2Denotes Jpt+nMapping to a profile SPtTwo-dimensional coordinates of (a);
(3-2) drawing an inclination angle line at two intersection points of each section and the geological boundary respectively, and calculating the intersection point coordinates of the two inclination angle lines according to the following formula:
Figure BDA0002253339720000041
in formula (II), x't0,y′t0Shown in cross section SPtTwo points of intersection Jp with the geological boundaryt、Jpt+nAt the intersection Jp of two lines of inclination drawn respectivelyt0Two dimensional coordinates, α for JptAngle of inclination of (d), β denotes Jpt+nThe dip angle is n represents the number of the sections;
(3-3) according to Jpt、Jpt0、Jpt+nDetermining a second-order Bezier curve, and calculating the highest point T of the second-order Bezier curve according to the following formulatThree-dimensional space coordinates in an O-xyz coordinate system:
Figure BDA0002253339720000042
in the formula (I), the compound is shown in the specification,
Figure BDA0002253339720000043
(3-4) according to the highest point T of all second-order Bezier curvestThe coordinates V (0,0, topz) of the vertex of the dome before denudation are calculated, wherein the topz is all the coordinates
Figure BDA0002253339720000045
Average value of (a).
Further, the step (4) specifically comprises:
(4-1) transforming the coordinates of the intersection of each section plane and the geological boundary generated in the step (2) and the coordinates of the dome apex according to the intersection set JP based on the following formula to obtain the coordinates when the coordinates are mapped to the coordinate system x 'Oy' corresponding to the section plane:
Figure BDA0002253339720000044
x′V＝0,y′V＝topz
in the formula, xt,yt,ztIndicating the t-th point Jp in the intersection set JPtIs also a section plane SPtThree-dimensional spatial coordinates, x, of the first intersection with the geological boundaryt+n,yt+n,zt+nRepresents the t + n point Jp in the intersection set JPt+nIs also a section plane SPtThree-dimensional spatial coordinates, x ', of a second intersection with the geological boundary't1,y′t1Denotes JptMapping to a profile SPtTwo-dimensional coordinates of (c), x't2,y′t2Denotes Jpt+nMapping to a profile SPtTwo-dimensional coordinates of (c), x'V,y′VRepresenting mapping of dome apex to profile SPtTopz represents the z-axis value in the three-dimensional space coordinate of the dome apex;
(4-2) drawing an inclination angle line at two intersection points of each section and the geological boundary respectively, drawing a tangent line at the dome vertex V, and calculating the coordinates of the intersection points of the tangent line and the two inclination angle lines according to the following formula as shown in figure 3:
Figure BDA0002253339720000051
in formula (II), x'tA,y′tARepresenting the tangent at the dome apex V and the point JptTwo-dimensional coordinate, x 'at intersection A of the lines of inclination'tB,y′tBRepresenting the tangent at the dome apex V and the point Jpt+nTwo-dimensional coordinates of the intersection B of the lines of inclination, α denotes JptAngle of inclination of (d), β denotes Jpt+nThe dip angle is n represents the number of the sections;
(4-3) respectively generating two second-order Bezier curves BZ at two sides of the dome vertext1、BZt2：
BZt1:
Figure BDA0002253339720000052
BZt2:
Figure BDA0002253339720000053
Wherein t is 1, n, x ', y' respectively represent x 'axis and y' axis coordinate variables on a coordinate system x 'Oy', and BZt1From points V, A, JptDetermination of BZt2From points V, B, Jpt+nDetermining;
(4-4) dividing each section along the z-axis to obtain 2n half sections, and converting coordinates of each point on the Bezier curves corresponding to all the half sections into an O-xyz coordinate system according to the following formula, thereby completing the generation of 2n Bezier curves:
Figure BDA0002253339720000054
in formula (II), x'j,y′jDenotes the j-th half section P obtained by divisionjCorresponding Bezier curve BZjCoordinate variable on the coordinate system x 'Oy', thetajShowing a half section PjThe included angle between the projection line on the xOy plane and the positive direction of the x-axis, xB,yB,zBRepresents the curve BZjCurve BZTD obtained by converting to O-xyz coordinate systemjWherein the section SPjThe half section obtained after the division is PjAnd Pj+nA half section PjThe corresponding Bezier curve is BZj1I.e. BZjAfter coordinate transformation, BZTDjHalf section Pj+nThe corresponding Bezier curve is BZj2I.e. BZj+nAfter coordinate transformation, BZTDj+n。
Further, the step (5) specifically comprises:
(5-1) Bezier curve BZTD generated according to two adjacent sections of fornixjAnd BZTDj+1And geological boundary GL between two sectionsjAs a condition, Morphing interpolation is performed, j is 1, …, n-1; the specific method comprises the following steps: using dome vertex V as initial geological boundary SRC, GLjDEST, BZTD as a target geological boundaryjAnd BZTDj+1Performing Morping interpolation between every two sections as a first constraint FC and a last constraint LC to form a transitional geological boundary with the number of the sections specified by a user;
and (5-2) after connecting the corresponding points, determining an irregular triangular net formed among contour lines according to the shortest diagonal principle, completing the restoration of the current layer of the dome structure, and then adjusting the x and y coordinates of the current layer of the dome to the horizontal plane, namely performing translation opposite to that in the step (1-3).
Has the advantages that: compared with the prior art, the invention has the following remarkable advantages: the method can effectively recover the dome-shaped landform, and solves the problem that no good dome-shaped landform recovery algorithm exists at present.
Drawings
FIG. 1 is a schematic flow diagram of one embodiment of the present invention;
FIG. 2 is a geological map of the experimental area in this example;
FIG. 3 is the DEM, geological boundary and its break points input in this embodiment;
FIG. 4 is a schematic diagram of a profile generating Bezier curve;
FIG. 5 is a schematic diagram of the top line, vertex, intersection, etc. elements of the cross section;
FIG. 6 is a three-dimensional and two-dimensional representation of Bezier curves and geological boundaries generated in the present example ((a) three-dimensional representation, (b) xOy plane representation, (c) yOz plane representation, and (d) xOz plane representation);
FIG. 7 is a schematic diagram of the basic elements of the Morphing process;
FIG. 8 is a schematic diagram of constructing an irregular triangulation network;
FIG. 9 shows the results of stratum recovery in examples (a top view and (b) perspective view).
Detailed Description
The embodiment discloses a method for restoring the three-dimensional geometric form of a dome structure, which comprises the following steps as shown in fig. 1:
step 1: and establishing an O-xyz coordinate system in a three-dimensional space by taking the northward direction as an x axis, the eastern direction as a y axis, the centroid of any layer of geological boundary of the dome structure as an origin O and the upward direction perpendicular to the xOy plane as a z axis, calculating the coordinates of all break points of the layer of geological boundary, and storing the coordinates into a point set G.
The experimental subject of this example was a vault structure developed in the ancient strata in Hubei province Ziguo county. The data source is geological map and DEM data (30 m resolution) of the region with a scale of 1:20 ten thousand, which are respectively shown in FIG. 2 and FIG. 3. The inner and outer boundary geological boundary data of a certain stratum of the dome are extracted according to a geological map (as shown in fig. 7), and each stratum inclination angle is 45 degrees. In this embodiment, 6 cross sections are provided.
The method specifically comprises the following steps:
(1-1) taking the north direction as an x axis and the east direction as a y axis, selecting an arbitrary point with the same height as a horizontal plane as an origin O ', and taking the upward direction perpendicular to the xO ' y plane as a z axis, and establishing an O ' -xyz coordinate system in a three-dimensional space;
(1-2) extracting the coordinates of all break points of any layer of geological boundary of the dome structure to be recovered in an O' -xyz coordinate system according to the geological map and DEM data, and storing a point set A ═ p in sequencei(pxi,pyi,pzi) I ═ 1, …, m }, where m is the number of points in the geological boundary element, (pxi,pyi,pzi) Represents a point piX, y, z coordinate values in the O' -xyz coordinate system;
(1-3) establishing an O-xyz coordinate system in a three-dimensional space by taking the centroid of the geological boundary as an origin O, and calculating the coordinates g of all break points of the geological boundary of the layer in the O-xyz coordinate system according to the following formulai(xi,yi,zi) And storing the point set G ═ Gi(xi,yi,zi)|i＝1,…,m}：
Figure BDA0002253339720000071
Step 2: and according to the point set G, generating n sections which are perpendicular to the xOy plane and pass through the O points based on the equal included angles, extracting 2n intersection points of the n sections and the geological boundary, and storing an intersection point set JP, wherein n is a positive integer.
The method specifically comprises the following steps:
(2-1) generating n sections which are perpendicular to the xOy plane and pass through the O point in the clockwise direction according to the principle that the included angle between any two adjacent sections is equal, taking y as a first section, andstoring the profile set SP ═ SPt|t＝1,…,n}；
(2-2) dividing each generated cross section along the z-axis to generate 2n half cross sections, and adding a half cross section set P in a clockwise order from the north direction to { P ═ P { (P) }j1, …,2n, where the half-section PjThe included angle theta between the projection line on the xOy plane and the positive direction of the x-axisjComprises the following steps:
Figure BDA0002253339720000081
θj∈[0,2π)，θjis the angle of clockwise rotation along the x-axis;
(2-3) calculating the origin O and the element G in the point set G based on the following formulaiFormed ray OgiAngle with positive direction of x-axisi；
Figure BDA0002253339720000082
In the formula, giRepresents the ith element in the point set G, m represents the number of elements in the point set G, xi,yiEach represents giCoordinate values;
(2-4) calculating the range of the included angle between the line segment formed by any two adjacent points in the point set G and the positive direction of the x axis based on the following formula:
Figure BDA0002253339720000083
in the formula, rangeiRepresents element G in point set Gi、gi+1Line segment L of compositioniRange of angle with positive direction of x-axis, rmini＝min(anglei,anglei+1),rmaxi＝max(anglei,anglei+1)，anglei+1Represents Ogi+1The included angle with the positive direction of the x axis;
(2-5) calculating coordinates of intersection points of each element in the half-section set P and the geological boundary based on the following formula:
Figure BDA0002253339720000084
in the formula, xj,yj,zjRepresenting an element P in a half-section set PjIntersection Jp with geological boundaryjIs equal to 1, 2n, i is e to [1, m ]]And satisfy rangeiIncluding thetaj,Ai＝yi+1-yi,Bi＝xi+1-xi,Ci＝xi+1yi-xiyi+1，xi,yi,ziDenotes giThree-dimensional spatial coordinate value of (2), xi+1,yi+1,zi+1Denotes gi+1Three-dimensional space coordinate values of (a);
(2-6) storing the intersection point coordinates of each element in the half-section set P and the geological boundary into an intersection point set JP ═ Jpj(xj,yj,zj) 1, …,2n, where the midpoint Jp is in the settAnd Jpt+nIs a profile SPtTwo intersection points with the geological boundary, t 1.
In this embodiment, all intersection coordinates are obtained by taking n as 6. All intersection coordinates in the clockwise direction from the x-axis are: (3393.15362372841,0,605.094404264791), (3390.78069807738,1957.66814879796,760.933831183495), (1949.90353157613,3377.33198654784,906.783044368775), (0,3257.85887014261,1042.50000636182), (-1299.09890681403,2250.10531065908,990.384219627612), (-2071.95295401048,1196.24259574619,979.442328482689), (-3287.76889180167,4.02635564968115E-13,1083.13663798333), (-3147.75429689418, -1817.35679065466,1081.71255498761), (-1921.07303496226, -3327.39610160518,1232.57227005481), (0, -2043.62213296309,709.400003467141), (864.585146566349, -1497.5054013223,817.870484683787), (1809.40841259259, -1044.6624340843,835.809199717015).
And step 3: and generating n second-order Bezier curves according to two intersection points of each section and the geological boundary and corresponding occurrence data at the intersection points, and calculating to obtain the V coordinate of the dome vertex.
The method specifically comprises the following steps:
(3-1) transforming the coordinates of the intersection of each section plane and the geological boundary generated in the step (2) according to the intersection set JP based on the following formula to obtain the coordinates when the coordinates are mapped to the coordinate system x 'Oy' corresponding to the section plane:
Figure BDA0002253339720000091
in the formula, xt,yt,ztIndicating the t-th point Jp in the intersection set JPtIs also a section plane SPtThree-dimensional spatial coordinates, x, of the first intersection with the geological boundaryt+n,yt+n,zt+nRepresents the t + n point Jp in the intersection set JPt+nIs also a section plane SPtThree-dimensional spatial coordinates, x ', of a second intersection with the geological boundary't1,y′t1Denotes JptMapping to a profile SPtTwo-dimensional coordinates of (c), x't2,y′t2Denotes Jpt+nMapping to a profile SPtTwo-dimensional coordinates of (a);
(3-2) drawing an inclination angle line at two intersection points of each section and the geological boundary respectively, and calculating the intersection point coordinates of the two inclination angle lines according to the following formula:
Figure BDA0002253339720000101
in formula (II), x't0,y′t0Shown in cross section SPtTwo points of intersection Jp with the geological boundaryt、Jpt+nAt the intersection Jp of two lines of inclination drawn respectivelyt0Two dimensional coordinates, as shown in FIG. 4, α denotes JptAngle of inclination of (d), β denotes Jpt+nThe dip angle is n represents the number of the sections;
(3-3) according to Jpt、Jpt0、Jpt+nDetermining a second-order Bezier curve, and calculating the highest point T of the second-order Bezier curve according to the following formulatThree-dimensional space coordinates in an O-xyz coordinate system:
Figure BDA0002253339720000102
in the formula (I), the compound is shown in the specification,
Figure BDA0002253339720000103
(3-4) according to the highest point T of all second-order Bezier curvestThe coordinates V (0,0, topz) of the vertex of the dome before denudation are calculated, wherein the topz is all the coordinates
Figure BDA0002253339720000105
Average value of (a). In the present example, topz is determined as 3358.1907.
And 4, step 4: for each intersection in set JP, 2n Bezier curves are generated based on the corresponding position dip, dome apex V, and tangent at apex V.
The method specifically comprises the following steps:
(4-1) transforming the coordinates of the intersection of each section plane and the geological boundary generated in the step (2) and the coordinates of the dome apex according to the intersection set JP based on the following formula to obtain the coordinates when the coordinates are mapped to the coordinate system x 'Oy' corresponding to the section plane:
Figure BDA0002253339720000104
x′V＝0,y′V＝topz
in the formula, xt,yt,ztIndicating the t-th point Jp in the intersection set JPtIs also a section plane SPtThree-dimensional spatial coordinates, x, of the first intersection with the geological boundaryt+n,yt+n,zt+nRepresents the t + n point Jp in the intersection set JPt+nIs also a section plane SPtThree-dimensional spatial coordinates, x ', of a second intersection with the geological boundary't1,y′t1Denotes JptMapping to a profile SPtTwo-dimensional coordinates of (c), x't2,y′t2Denotes Jpt+nMapping to a profile SPtTwo ofDimension coordinate, x'V,y′VRepresenting mapping of dome apex to profile SPtTopz represents the z-axis value in the three-dimensional space coordinate of the dome apex;
(4-2) drawing an inclination angle line at each section and two intersection points of the geological boundary, drawing a tangent line at the dome vertex V, and calculating the coordinates of the intersection points of the tangent line and the two inclination angle lines according to the following formula as shown in FIG. 5:
Figure BDA0002253339720000111
in formula (II), x'tA,y′tARepresenting the tangent at the dome apex V and the point JptTwo-dimensional coordinate, x 'at intersection A of the lines of inclination'tB,y′tBRepresenting the tangent at the dome apex V and the point Jpt+nTwo-dimensional coordinates of the intersection B of the lines of inclination, α denotes JptAngle of inclination of (d), β denotes Jpt+nThe dip angle is n represents the number of the sections;
(4-3) respectively generating two second-order Bezier curves BZ at two sides of the dome vertext1、BZt2：
BZt1:
Figure BDA0002253339720000112
BZt2:
Figure BDA0002253339720000113
Wherein t is 1, n, x ', y' respectively represent x 'axis and y' axis coordinate variables on a coordinate system x 'Oy', and BZt1From points V, A, JptDetermination of BZt2From points V, B, Jpt+nDetermining;
(4-4) dividing each section along the z-axis to obtain 2n half sections, and converting coordinates of each point on the Bezier curves corresponding to all the half sections into an O-xyz coordinate system according to the following formula, thereby completing the generation of 2n Bezier curves:
Figure BDA0002253339720000114
in formula (II), x'j,y′jDenotes the j-th half section P obtained by divisionjCorresponding Bezier curve BZjCoordinate variable on the coordinate system x 'Oy', thetajShowing a half section PjThe included angle between the projection line on the xOy plane and the positive direction of the x-axis, xB,yB,zBRepresents the curve BZjCurve BZTD obtained by converting to O-xyz coordinate systemjWherein the section SPjThe half section obtained after the division is PjAnd Pj+nA half section PjThe corresponding Bezier curve is BZj1I.e. BZjAfter coordinate transformation, BZTDjHalf section Pj+nThe corresponding Bezier curve is BZj2I.e. BZj+nAfter coordinate transformation, BZTDj+n。
In this embodiment, the Bezier curves generated by the inner and outer layer boundary lines in each three-dimensional space are shown in fig. 6.
And 5: and (3) interpolating between adjacent Bezier curves to generate a contour line set based on a Morphing method, and constructing an irregular triangular network between the adjacent contour lines to finish the recovery of the current layer of the dome structure.
The method specifically comprises the following steps:
(5-1) Bezier curve BZTD generated according to two adjacent sections of fornixjAnd BZTDj+1And geological boundary GL between two sectionsjAs a condition, Morphing interpolation is performed, j is 1, …, n-1; the specific method comprises the following steps: using dome vertex V as initial geological boundary SRC, GLjDEST, BZTD as a target geological boundaryjAnd BZTDj+1Performing Morping interpolation between every two sections as a first constraint FC and a last constraint LC to form a transitional geological boundary with the number of the sections specified by a user, as shown in FIG. 7;
(5-2) after connecting the corresponding points (as shown in fig. 8, A and C between adjacent contour lines are corresponding points, B and D are corresponding points), determining an irregular triangular net formed between the contour lines according to the shortest diagonal principle, and completing the current layer recovery of the dome structure. The x, y coordinates of the current level of the dome are then adjusted to the horizontal plane, i.e., translated in the opposite direction as in (1-3).
Step 6: and circularly executing the steps 1-5 until all layers are recovered. The final restoration result is shown in fig. 9, and fig. 9 shows the effect of overlapping the dome restoration result and the original DEM.
According to the embodiment, the method can effectively recover the dome-shaped landform, and solves the problem that no good dome-shaped landform recovery algorithm exists at present.
While the invention has been described in connection with what is presently considered to be the most practical and preferred embodiment, it is to be understood that the invention is not to be limited to the disclosed embodiment, but on the contrary, is intended to cover various modifications and equivalent arrangements included within the spirit and scope of the appended claims.

Claims (6)

1. A method for recovering the three-dimensional geometry of a dome, the method comprising:
(1) establishing an O-xyz coordinate system in a three-dimensional space by taking the northward direction as an x axis and the eastern direction as a y axis, taking the centroid of any layer of geological boundary of a vault structure as an origin O and the upward direction perpendicular to the xOy plane as a z axis, calculating the coordinates of all break points of the layer of geological boundary, and storing the coordinates into a point set G;
(2) generating n sections which are perpendicular to the xOy plane and pass through the O points based on equal included angles according to the point set G, extracting 2n intersection points of the n sections and the geological boundary, and storing an intersection point set JP, wherein n is a positive integer;
(3) generating n second-order Bezier curves according to two intersection points of each section and the geological boundary and corresponding occurrence data at the intersection points, and calculating to obtain a dome vertex V coordinate;
(4) for each intersection point in the set JP, generating 2n Bezier curves based on the corresponding position dip angle, the dome vertex V and the tangent line at the vertex V;
(5) interpolation is carried out between adjacent Bezier curves to generate a contour line set based on a Morphing method, an irregular triangular net is constructed between adjacent contour lines, and the recovery of the current layer of the dome structure is completed;
(6) and (5) circularly executing the steps (1) to (5) until all the layers are recovered.
2. The method for restoring a three-dimensional geometry of a dome structure according to claim 1, wherein: the step (1) specifically comprises the following steps:
(1-1) taking the north direction as an x axis and the east direction as a y axis, selecting an arbitrary point with the same height as a horizontal plane as an origin O ', and taking the upward direction perpendicular to the xO ' y plane as a z axis, and establishing an O ' -xyz coordinate system in a three-dimensional space;
(1-2) extracting the coordinates of all break points of any layer of geological boundary of the dome structure to be recovered in an O' -xyz coordinate system according to the geological map and DEM data, and storing a point set A ═ p in sequencei(pxi,pyi,pzi) I ═ 1, …, m }, where m is the number of points in the geological boundary element, (pxi,pyi,pzi) Represents a point piX, y, z coordinate values in the O' -xyz coordinate system;
(1-3) establishing an O-xyz coordinate system in a three-dimensional space by taking the centroid of the geological boundary as an origin O, and calculating the coordinates g of all break points of the geological boundary of the layer in the O-xyz coordinate system according to the following formulai(xi,yi,zi) And storing the point set G ═ Gi(xi,yi,zi)|i＝1,…,m}：
Figure FDA0002253339710000011
3. The method for restoring a three-dimensional geometry of a dome structure according to claim 1, wherein: the step (2) specifically comprises the following steps:
(2-1) generating n sections which are perpendicular to the xOy plane and pass through the O point in the clockwise direction according to the principle that the included angle between any two adjacent sections is equal by taking y as a first section and storing a section set SP as { SP ═ SPt|t＝1,…,n}；
(2-2) dividing each generated cross section along the z-axis to generate 2n half cross sections, and adding a half cross section set P in a clockwise order from the north direction to { P ═ P { (P) }j1, …,2n, where the half-section PjThe included angle theta between the projection line on the xOy plane and the positive direction of the x-axisjComprises the following steps:
Figure FDA0002253339710000021
(2-3) calculating the origin O and the element G in the point set G based on the following formulaiFormed ray OgiAngle with positive direction of x-axisi；
Figure FDA0002253339710000022
In the formula, giRepresents the ith element in the point set G, m represents the number of elements in the point set G, xi,yiEach represents giCoordinate values;
(2-4) calculating the range of the included angle between the line segment formed by any two adjacent points in the point set G and the positive direction of the x axis based on the following formula:
Figure FDA0002253339710000023
in the formula, rangeiRepresents element G in point set Gi、gi+1Line segment L of compositioniRange of angle from the positive direction of the x-axis, rmini＝min(anglei,anglei+1),rmaxi＝max(anglei,anglei+1)，anglei+1Represents Ogi+1The included angle with the positive direction of the x axis;
(2-5) calculating coordinates of intersection points of each element in the half-section set P and the geological boundary based on the following formula:
Figure FDA0002253339710000024
in the formula, xj,yj,zjRepresenting an element P in a half-section set PjIntersection Jp with geological boundaryjIs equal to 1, 2n, i is e to [1, m ]]And satisfy rangeiIncluding thetaj,Ai＝yi+1-yi,Bi＝xi+1-xi,Ci＝xi+1yi-xiyi+1，xi,yi,ziDenotes giThree-dimensional spatial coordinate value of (2), xi+1,yi+1,zi+1Denotes gi+1Three-dimensional space coordinate values of (a);
(2-6) storing the intersection point coordinates of each element in the half-section set P and the geological boundary into an intersection point set JP ═ Jpj(xj,yj,zj) 1, …,2n, where the midpoint Jp is in the settAnd Jpt+nIs a profile SPtTwo intersection points with the geological boundary, t 1.
4. The method for restoring a three-dimensional geometry of a dome structure according to claim 1, wherein: the step (3) specifically comprises the following steps:
(3-1) transforming the coordinates of the intersection of each section plane and the geological boundary generated in the step (2) according to the intersection set JP based on the following formula to obtain the coordinates when the coordinates are mapped to the coordinate system x 'Oy' corresponding to the section plane:
Figure FDA0002253339710000031
in the formula, xt,yt,ztIndicating the t-th point Jp in the intersection set JPtIs also a section plane SPtThree-dimensional spatial coordinates, x, of the first intersection with the geological boundaryt+n,yt+n,zt+nRepresents the t + n point Jp in the intersection set JPt+nIs also a section plane SPtThree-dimensional spatial coordinates, x ', of a second intersection with the geological boundary't1,y′t1Denotes JptMapping to a profile SPtTwo-dimensional coordinates of (c), x't2,y′t2Denotes Jpt+nMapping to a profile SPtTwo-dimensional coordinates of (a);
(3-2) drawing an inclination angle line at two intersection points of each section and the geological boundary respectively, and calculating the intersection point coordinates of the two inclination angle lines according to the following formula:
Figure FDA0002253339710000032
in formula (II), x't0,y′t0Shown in cross section SPtTwo points of intersection Jp with the geological boundaryt、Jpt+nAt the intersection Jp of two lines of inclination drawn respectivelyt0Two dimensional coordinates, α for JptAngle of inclination of (d), β denotes Jpt+nThe dip angle is n represents the number of the sections;
(3-3) according to Jpt、Jpt0、Jpt+nDetermining a second-order Bezier curve, and calculating the highest point T of the second-order Bezier curve according to the following formulatThree-dimensional space coordinates in an O-xyz coordinate system:
Figure FDA0002253339710000041
in the formula (I), the compound is shown in the specification,
Figure FDA0002253339710000042
(3-4) according to the highest point T of all second-order Bezier curvestThe coordinates V (0,0, topz) of the vertex of the dome before denudation are calculated, wherein the topz is all the coordinates
Figure FDA0002253339710000045
Average value of (a).
5. The method for restoring a three-dimensional geometry of a dome structure according to claim 1, wherein: the step (4) specifically comprises the following steps:
(4-1) transforming the coordinates of the intersection of each section plane and the geological boundary generated in the step (2) and the coordinates of the dome apex according to the intersection set JP based on the following formula to obtain the coordinates when the coordinates are mapped to the coordinate system x 'Oy' corresponding to the section plane:
Figure FDA0002253339710000043
x′V＝0,y′V＝topz
in the formula, xt,yt,ztIndicating the t-th point Jp in the intersection set JPtIs also a section plane SPtThree-dimensional spatial coordinates, x, of the first intersection with the geological boundaryt+n,yt+n,zt+nRepresents the t + n point Jp in the intersection set JPt+nIs also a section plane SPtThree-dimensional spatial coordinates, x ', of a second intersection with the geological boundary't1,y′t1Denotes JptMapping to a profile SPtTwo-dimensional coordinates of (c), x't2,y′t2Denotes Jpt+nMapping to a profile SPtTwo-dimensional coordinates of (c), x'V,y′VRepresenting mapping of dome apex to profile SPtTopz represents the z-axis value in the three-dimensional space coordinate of the dome apex;
(4-2) drawing an inclination angle line at two intersection points of each section and the geological boundary respectively, drawing a tangent line at the vertex V of the dome, and calculating the intersection point coordinates of the tangent line and the two inclination angle lines according to the following formula:
Figure FDA0002253339710000044
in formula (II), x'tA,y′tARepresenting the tangent at the dome apex V and the point JptTwo-dimensional coordinate, x 'at intersection A of the lines of inclination'tB,y′tBRepresenting the tangent at the dome apex V and the point Jpt+nTwo-dimensional coordinates of the intersection B of the lines of inclination, α denotes JptAngle of inclination of (d), β denotes Jpt+nThe dip angle is n represents the number of the sections;
(4-3) respectively generating two second-order Bezier curves BZ at two sides of the dome vertext1、BZt2：
BZt1:
Figure FDA0002253339710000051
BZt2:
Figure FDA0002253339710000052
Wherein t is 1, n, x ', y' respectively represent x 'axis and y' axis coordinate variables on a coordinate system x 'Oy', and BZt1From points V, A, JptDetermination of BZt2From points V, B, Jpt+nDetermining;
(4-4) dividing each section along the z-axis to obtain 2n half sections, and converting coordinates of each point on the Bezier curves corresponding to all the half sections into an O-xyz coordinate system according to the following formula, thereby completing the generation of 2n Bezier curves:
Figure FDA0002253339710000053
in formula (II), x'j,y′jDenotes the j-th half section P obtained by divisionjCorresponding Bezier curve BZjCoordinate variable on the coordinate system x 'Oy', thetajShowing a half section PjThe included angle between the projection line on the xOy plane and the positive direction of the x-axis, xB,yB,zBRepresents the curve BZjCurve BZTD obtained by converting to O-xyz coordinate systemjWherein the section SPjThe half section obtained after the division is PjAnd Pj+nA half section PjThe corresponding Bezier curve is BZj1I.e. BZjAfter coordinate transformation, BZTDjHalf section Pj+nThe corresponding Bezier curve is BZj2I.e. BZj+nAfter coordinate transformation, BZTDj+n。
6. The method for restoring a three-dimensional geometry of a dome structure according to claim 1, wherein: the step (5) specifically comprises the following steps:
(5-1) Bezier curve BZTD generated according to two adjacent sections of fornixjAnd BZTDj+1And geological boundary GL between two sectionsjAs a conditionCarrying out Morping interpolation, wherein j is 1, …, n-1; the specific method comprises the following steps: using dome vertex V as initial geological boundary SRC, GLjDEST, BZTD as a target geological boundaryjAnd BZTDj+1Performing Morping interpolation between every two sections as a first constraint FC and a last constraint LC to form a transitional geological boundary with the number of the sections specified by a user;
and (5-2) after connecting the corresponding points, determining an irregular triangular net formed among contour lines according to the shortest diagonal principle, completing the restoration of the current layer of the dome structure, and then adjusting the x and y coordinates of the current layer of the dome to the horizontal plane.
CN201911042921.3A 2019-10-30 2019-10-30 Three-dimensional geometric form restoration method of dome structure Active CN110910495B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201911042921.3A CN110910495B (en) 2019-10-30 2019-10-30 Three-dimensional geometric form restoration method of dome structure

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201911042921.3A CN110910495B (en) 2019-10-30 2019-10-30 Three-dimensional geometric form restoration method of dome structure

Publications (2)

Publication Number Publication Date
CN110910495A true CN110910495A (en) 2020-03-24
CN110910495B CN110910495B (en) 2023-05-26

Family

ID=69814734

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201911042921.3A Active CN110910495B (en) 2019-10-30 2019-10-30 Three-dimensional geometric form restoration method of dome structure

Country Status (1)

Country Link
CN (1) CN110910495B (en)

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN112581558A (en) * 2020-10-26 2021-03-30 南京师范大学 Model construction method and system for intrusion structure in map-cut geological profile

Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101750626A (en) * 2008-12-16 2010-06-23 中国石油天然气集团公司 Data acquisition designing method in three-dimensional seismic physical simulation
WO2010148625A1 (en) * 2009-06-24 2010-12-29 中国石油集团川庆钻探工程有限公司 Block model constructing method for complex geological structures
CN104077810A (en) * 2014-06-26 2014-10-01 中南大学 Offset fault restoration modeling method based on TIN
CN105093311A (en) * 2015-06-29 2015-11-25 成都理工大学 Stratum denudation thickness measurement method for superposition basin multi-phase difference upheaval region
CN105388524A (en) * 2014-09-04 2016-03-09 中国石油化工股份有限公司 Complex geologic body precise description taking triangular surface grid as limiting condition

Patent Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101750626A (en) * 2008-12-16 2010-06-23 中国石油天然气集团公司 Data acquisition designing method in three-dimensional seismic physical simulation
WO2010148625A1 (en) * 2009-06-24 2010-12-29 中国石油集团川庆钻探工程有限公司 Block model constructing method for complex geological structures
CN104077810A (en) * 2014-06-26 2014-10-01 中南大学 Offset fault restoration modeling method based on TIN
CN105388524A (en) * 2014-09-04 2016-03-09 中国石油化工股份有限公司 Complex geologic body precise description taking triangular surface grid as limiting condition
CN105093311A (en) * 2015-06-29 2015-11-25 成都理工大学 Stratum denudation thickness measurement method for superposition basin multi-phase difference upheaval region

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN112581558A (en) * 2020-10-26 2021-03-30 南京师范大学 Model construction method and system for intrusion structure in map-cut geological profile
CN112581558B (en) * 2020-10-26 2024-04-19 南京师范大学 Model construction method and system for intrusion structure in cut geological section

Also Published As

Publication number Publication date
CN110910495B (en) 2023-05-26

Similar Documents

Publication Publication Date Title
Owen A survey of unstructured mesh generation technology.
De Smith Determination of gradient and curvature constrained optimal paths
CN108053477B (en) Numerical processing method for deformation in pipeline
CN116152461B (en) Geological modeling method, device, computer equipment and computer readable storage medium
CN102881048B (en) Point cloud clipping-based generation method for spatial curved surface
CN110910495B (en) Three-dimensional geometric form restoration method of dome structure
CN112489209B (en) Collapse rock disaster scene reappearing method
CN109215124A (en) The construction method of underground engineering 3D grid model under a kind of complex geological condition
CN115630478A (en) Discrete fracture network generation method based on iterative function system
CN101587597A (en) Construction method based on the geologic rule constraint complex-structure blocky geologic model
CN108986212B (en) Three-dimensional virtual terrain LOD model generation method based on crack elimination
CN116152512A (en) Height measuring and calculating method based on building shadow restoration
CN114777745A (en) Inclined evidence obtaining modeling method based on unscented Kalman filtering
JP6025615B2 (en) Image drawing device
CN110910499B (en) Method and device for constructing geological environment carrier fault model based on Revit software
CN110084886A (en) A kind of geological space restored method for taking the constraint of geology volume morphing-occurrence-toughness into account
CN104834005B (en) A kind of method for determining geological interface
Fua Cartographic applications of model-based optimization
CN113269886B (en) Slope three-dimensional digital twin model building method based on multi-source data fusion
CN113240812B (en) Ultra-thin manganese ore body three-dimensional modeling method based on incremental simulation
CN113821850B (en) Geological boundary optimization method for steep slope oblique photography model by utilizing point offset technology
CN110992488B (en) Inclined crack grid dividing method based on embedded discrete crack model
CN113110429A (en) Minimum lasting formation generation and control method of multi-robot system under visual field constraint
CN105760601A (en) Method for automatically matching ore body boundary lines in sequence section
CN114298949A (en) Multi-region fusion track generation method

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