CN104809692A - Interpolation method for DEM (Digital Elevation Model) generated by contour lines - Google Patents

Interpolation method for DEM (Digital Elevation Model) generated by contour lines Download PDF

Info

Publication number
CN104809692A
CN104809692A CN201510237207.5A CN201510237207A CN104809692A CN 104809692 A CN104809692 A CN 104809692A CN 201510237207 A CN201510237207 A CN 201510237207A CN 104809692 A CN104809692 A CN 104809692A
Authority
CN
China
Prior art keywords
level line
dem
contour
image
contour lines
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
CN201510237207.5A
Other languages
Chinese (zh)
Other versions
CN104809692B (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.)
Shandong Linyi Tobacco Co Ltd
Original Assignee
Shandong Linyi Tobacco Co Ltd
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 Shandong Linyi Tobacco Co Ltd filed Critical Shandong Linyi Tobacco Co Ltd
Priority to CN201510237207.5A priority Critical patent/CN104809692B/en
Publication of CN104809692A publication Critical patent/CN104809692A/en
Application granted granted Critical
Publication of CN104809692B publication Critical patent/CN104809692B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Landscapes

  • Image Generation (AREA)
  • Image Analysis (AREA)

Abstract

The invention discloses an interpolation method for a DEM (Digital Elevation Model) generated by contour lines. The interpolation method comprises the following steps: 1, acquiring and loading vector contour line data meeting specific conditions, wherein the specific conditions are: each contour line is spatially and continuously distributed and cannot be interrupted, and odd contour lines and even contour lines exist at intervals; 2, performing morphological corrosion reconstruction and morphological expansion reconstruction to generate a lower plateau image and an upper plateau image respectively; 3, performing obstacle Euclidean distance conversion to acquire a geodesic path and a geodesic distance, and performing pipelined linear interpolation to generate the DEM. The generated DEM has the characteristics of high accuracy and high fidelity and is generated by using the linear interpolation method, the DEM accuracy obtained under the condition that only contour line data are used is equal to or higher than the DEM accuracy obtained by using a map algebra method, manual intervention is not needed, the DEM is not limited by iterations, an adaptive algorithm is truly achieved, and the interpolation method can be applied to various industries.

Description

The interpolating method of a kind of contour lines creation DEM
Technical field
The present invention relates to the interpolating method of a kind of digital elevation model (DEM), specifically the interpolating method of a kind of contour lines creation DEM, belong to cartography and Geographical Information Sciences technical field.
Background technology
Digital elevation model (Digital Elevation Model, DEM) is a kind of discrete numeral expression to earth topographical surface landforms.Since the 1950's, the later stage was suggested, DEM receives great concern, and is widely used in fields such as mapping, civil engineering work, geology, mine engineering, Landscape Architecture, highway layout, flood control, agricultural, planning, military engineering, aircraft and battlefield emulation.Meanwhile, in the ascendant all the time to the research of DEM, the interpolating method wherein generating DEM is one of emphasis of research.Generally speaking, different pieces of information source needs different interpolating methods to generate DEM.At present, the data generating DEM are mainly derived from topomap, remotely-sensed data (both comprised space flight and aviation image, and comprised again synthetic aperture radar interferometry data and laser radar data), ground survey, existing DEM etc.; Wherein, obtaining DEM from topomap is a kind of method be most widely used at present.The Mapping departments of China just utilize the digital adventure appearance of 1:1 ten thousand, 1:5 ten thousand and 1:25 ten thousand engineer's scale (a kind of carrier of topomap) to generate the DEM of multiple resolution respectively.
Usually, when obtaining DEM by topomap, based on isocontour distribution characteristics, three kinds of modes are had to generate DEM: level line discretize, the interpolation of contours and level line build Di Luoni TIN (Triangulated Irregular Network, TIN).Level line discretization method essence level line is regarded as the data of irregular distribution, do not consider the characteristic of level line itself, and this causes the DEM generated may occur some exceptions; Steepest gradient (streamline) interpolation algorithm generating DEM based on contour line data that Lebel and Olson proposes, interpolation theory is fairly simple, but because digitized level line is far from papery topomap level line intuitively, therefore, the method implements also exists many problems; Because the direct TIN built by level line exists " flat triangular " (namely, level triangle) problem, therefore, current engineering generally adopts the method building TIN based on level line and additional " characteristic " (as topographic structure line and characteristic number strong point such as mountain top point, depression points, saddle point etc.) in producing.
Proposed the interpolating method of much new contour lines creation DEM in recent years, wherein the achievement in research of the people such as Hu Peng, Hu Hai is more representative.Their argument foundation: above-mentioned " characteristic " is isocontour dual form in essence, is not necessary; And in engineering is produced, be difficult to the density of controlling feature data to balance precision and the workload of DEM.Based on this, Hu Peng, Hu Hai etc. successively discuss the method utilizing map algebra directly to be generated DEM by the interpolation of contours, that is: MADEM.Map algebra essentially, is a kind of image manipulation be based upon on range conversion operating basis; It is used for interpolation when generating DEM, not only do not need extra supplemental characteristic data, and the DEM generated has higher precision, meets the DEM precision evaluation standard of " elevation orderisomorphism ".At present, the DEM based on map algebra exploitation generates software and has been applied to beginning in the production run of China 1:5 ten thousand DEM, and its efficiency and precision obtain the accreditation of research and production department gradually.
But also there are urgently improvements based on the interpolating method of map algebra.Specifically, because the method asks for half-interval contour (namely to two adjacent contour intervals from equal line) C by iteration l/ 2, C l/ 4, C l/ 8, C l/ 16, C l/ 32 ... (C lfor basic contour interval isocontour on topomap) generate DEM, namely ask for the border of two adjacent isocontour Voronoi figure iteratively and give this border by two isocontour mean values; To when dividing unnecessarily, be boundary's (wherein, n is maximum iteration time) with 1/2n+1Voronoi figure, corresponding elevation is composed in layering again.This is also a kind of linear interpolation method in essence.But the method needs to pre-enter maximum iteration time, when having exceeded maximum iteration time, the point be not interpolated into all gives the increment of last half contour interval.Visible, maximum iteration time is that the method interpolation generates the crucial index of in DEM one.But because this index have employed the mode of manual intervention, do not accomplish self-adaptation, this limits the engineering application of the method to a certain extent.
Summary of the invention
For the deficiencies in the prior art, the invention provides the interpolating method of a kind of contour lines creation DEM, it need not any manual intervention, generates DEM by linear interpolation method, and the DEM precision obtained is equal to or higher than the DEM precision obtained by map algebraic method.
The present invention solves the technical scheme that its technical matters takes: the interpolating method of a kind of contour lines creation DEM, is characterized in that, comprise the following steps:
Step one, obtain and load the vector contour line data meeting specified conditions, described specified conditions are: each level line spatially continuous distribution, can not interrupt, and odd number level line and even number contour interval occur;
Step 2, carries out morphological erosion and rebuilds and the morphological dilations reconstruction lower high original image of generation and upper high original image respectively;
Step 3, carries out obstacle Euclidean Distance Transform and obtains geodetic path and geodesic distance, and the linear interpolation carried out along streamline generates DEM.
Preferably, in step, level line thematic data is separated to obtain the contour line data of vector form by the level line in digitizing papery topomap or from digital line is drawn, and give corresponding height value to each level line, require each level line spatially continuous distribution and must not interrupting without reason; In the level line obtained, level line is odd number level line in any two adjacent contours, another level line is even number level line, and odd number level line and even number contour interval occur.
Preferably, in step 2, rasterizing operation generation level line image is carried out to vector level line, to level line image carry out morphological erosion reconstruction and morphological dilations rebuild generate one respectively under high original image and upper high original image.
Preferably, in step 3, obtain the length to level line on it and lower isocontour geodetic path of any one pixel on level line image by the isocontour range conversion of the odd number being barrier with even number level line with the isocontour range conversion of even number two obstacle Euclidean Distance Transforms that odd number level line is barrier respectively, the linear interpolation carried out along streamline generates DEM.
Further, described step one specifically comprises the following steps: (1) obtains vector contour line data; (2) level line of process fracture, and give its corresponding height value to level line; (3) the odd even feature of adjacent contour height value is checked.
Further, described step 2 specifically comprises the following steps: (1) vector level line rasterizing; (2) assisted image is generated; (3) carry out morphological erosion to rebuild and morphological dilations reconstruction.
Further, described step 3 specifically comprises the following steps: (1) odd, even isocontour obstacle Euclidean Distance Transform; (2) linear interpolation along streamline is carried out.
Beneficial effect of the present invention is as follows:
(1) according to the concept of the geodesic distance of mathematical morphology, convert two kinds of image operations based on morphological reconstruction and distance of obstacle, the present invention proposes a kind of based on morphological reconstruction and the integrated DEM interpolating method of obstacle Euclidean Distance Transform.Interpolating method is the key factor affecting digital elevation model precision, obtains the process of DEM from level line, and the selection of interpolating method is a very important sport technique segment, and it has important impact to the precision generating DEM.In the present invention, morphological reconstruction is for obtaining higher level line corresponding to arbitrary pixel and lower isocontour height value, obstacle Euclidean Distance Transform is for obtaining this point to upper and lower two isocontour geodesic distances, linear interpolation is used to obtain the height value of this point, suitable tight of whole flow process and have science.
(2) compared with " MADEM " interpolating method proposed with Hu Peng, Hu Hai, when only using contour line data, the result that the DEM precision that the present invention generates equals or obtains higher than " MADEM ".But " MADEM " interpolating method that Hu Peng, Hu Hai propose needs manual intervention, is subject to the restriction of iterations, limit the engineering application of the method; And the present invention adopts DEM interpolating method not only without the need to manual intervention, and be not subject to the restriction of iterations, really achieve adaptive algorithm, all trades and professions can be applied.
Accompanying drawing explanation
Fig. 1 is method flow diagram of the present invention;
Fig. 2 (a) to Fig. 2 (g) for isocontour elevation is distributed between 630-920m, long and wide when being respectively 1400m × 1100m based on the schematic diagram of morphological reconstruction and the integrated interpolation of contours of obstacle Euclidean Distance Transform;
Fig. 2 (a) is the schematic diagram of level line image;
Fig. 2 (b) is the schematic diagram of lower high original image;
Fig. 2 (c) is the schematic diagram of upper high original image;
Fig. 2 (d) is for strange level line is to the schematic diagram of even isocontour geodesic distance;
Fig. 2 (e) is for even level line is to the schematic diagram of strange isocontour geodesic distance;
The schematic diagram of the DEM that Fig. 2 (f) generates for interpolation;
The schematic diagram of the level line that Fig. 2 (g) generates for DEM and former level line fit figure;
Fig. 3 is the upper and lower isocontour schematic diagram of certain pixel p on level line image;
Fig. 4 is the schematic diagram being obtained geodesic distance by obstacle Euclidean Distance Transform;
Fig. 4 (a) is not for carrying out the former figure (white is entity, and black is barrier) of obstacle Euclidean Distance Transform;
Fig. 4 (b) is for carrying out the result schematic diagram after obstacle Euclidean Distance Transform;
Fig. 5 (a) and (b) are the schematic diagram of the DEM generated based on " MADEM " interpolating method of Hu Peng, the schematic diagram of DEM when Fig. 5 (a) iterations is 10 times, the schematic diagram of DEM when Fig. 5 (b) iterations is 30 times;
The difference schematic diagram of the DEM contrast of " MADEM " interpolating method acquisition of Fig. 6 interpolating method of the present invention and Hu Peng;
Fig. 7 (a) to Fig. 7 (c) for isocontour elevation is distributed between 240-940m, long and wide when being respectively 4400m × 3100m based on the schematic diagram of morphological reconstruction and the integrated interpolation of contours of obstacle Euclidean Distance Transform;
The strange level line of Fig. 7 (a) is to the schematic diagram of even isocontour geodesic distance (eigenwert is between 0 to 15.11m);
The schematic diagram of the DEM (brightness represents the size of elevation) that Fig. 7 (b) interpolation generates;
The level line (black) that Fig. 7 (c) interpolation DEM generates and former level line fit figure and partial enlarged drawing thereof.
Embodiment
For clearly demonstrating the technical characterstic of this programme, below by embodiment, and in conjunction with its accompanying drawing, the present invention will be described in detail.Disclosing hereafter provides many different embodiments or example is used for realizing different structure of the present invention.Of the present invention open in order to simplify, hereinafter the parts of specific examples and setting are described.In addition, the present invention can in different example repeat reference numerals and/or letter.This repetition is to simplify and clearly object, itself does not indicate the relation between discussed various embodiment and/or setting.It should be noted that parts illustrated in the accompanying drawings are not necessarily drawn in proportion.Present invention omits the description of known assemblies and treatment technology and process to avoid unnecessarily limiting the present invention.
For the needs generating high precision, high-fidelity DEM based on the interpolation of contours, the invention provides the interpolating method of a kind of contour lines creation DEM based on morphological reconstruction and range conversion, especially a kind of upper plateau based on morphological reconstruction, lower plateau image generating method, a kind of geodetic path based on obstacle Euclidean Distance Transform and geodesic distance generation method is proposed, with the DEM utilizing the interpolation of contours to generate high precision, high-fidelity.
As shown in Figure 1, the interpolating method of a kind of contour lines creation DEM based on morphological reconstruction and range conversion of the present invention, it is characterized in that, arbitrary level line spatially continuous distribution, can not to interrupt and odd number level line and even number contour interval occur, said method comprising the steps of:
Step one, obtains and loads the vector contour line data meeting specified conditions;
Step 2, carries out morphological erosion reconstruction and morphological dilations reconstruction generates the two panel height original images such as lower high original image, upper high original image respectively;
Step 3, carries out obstacle Euclidean Distance Transform and obtains geodetic path and geodesic distance, and the linear interpolation carried out along streamline generates DEM.
Further, described step one specifically comprises the following steps:
(1) vector contour line data is obtained
Original topomap, can by the data of digital acquisition vector form if papery topomap; If raw data itself provides with the form of digital adventure appearance, then can extract its level line.
(2) level line of process fracture, and give its corresponding height value to level line
To the requirement of vector level line demand fulfillment space continuity.That is, in a figure amplitude range, a level line must not interrupt without reason.Otherwise, for the level line of fracture, need the level line being connected fracture by artificial attended operation, to obtain complete level line.Its corresponding height value is given to every bar level line.
(3) the odd even feature of adjacent contour height value is checked
According to the size of height value, all level lines are sorted, and from minimum value, carry out odd number or the isocontour Type division of even number.In principle, in two adjacent height value, an odd number, an even number must be had, and odd number level line and even number contour interval occur.
Further, described step 2 specifically comprises the following steps:
(1) vector level line rasterizing
Rasterizing operation is carried out to vector level line, generates a level line image, as shown in Fig. 2 (a).On level line image, the corresponding height value of this level line of the pixel on each level line is filled, and is an invalid value (as shown in Figure 3), is designated as contour map as CL for the pixel filling on what level line not in office.In Fig. 3, the pixel p between white portion surround by two level lines, the level line that height value is larger is expressed as C up (), the level line that height value is less is expressed as C l(p); Meanwhile, p be mountain top point or concave point time can only by C u(p) or C lp () surrounded, now C u(p)=C l(p).
For each pixel, if from two panel height original image P below contour map can be determined as CL l(being designated as: lower high original image, as Suo Shi Fig. 2 (c)) and P u(being designated as: upper high original image, as Suo Shi Fig. 2 (b)), then two isocontour elevations that arbitrary pixel is corresponding are known, that is:
P l(p)=CL[C l(p)] (1)
P u(p)=CL[C u(p)]
Wherein CL [] is corresponding isocontour height value.
(2) assisted image is generated
With above-mentioned level line image for reference, generate an assisted image M, and M is identical with codomain with the level line image size of rasterizing, the rule of generation is:
Wherein h maxfor being greater than a number of this level line image maximal value.
(3) carry out morphological erosion to rebuild and morphological dilations reconstruction
Being marking image with M, take CL as mask image, carries out morphological erosion and rebuilds and obtain lower high original image P l; In like manner, being marking image with CL, take M as mask image, carries out morphological dilations and rebuilds the upper high original image P of acquisition u.
Morphological erosion reconstruction and morphological dilations are rebuild two kinds of geodetics operated from mathematical morphology and are converted.In mathematical morphology, the fundamental operation of geodetic conversion comprises geodesic dilation and geodetic corrodes two operators.Geodesic dilation relates to two width images: marking image and mask image, and the size of this two width image is identical with field of definition, but the value of each pixel of mask image must be greater than or equal to the value of pixel of the same name corresponding to marking image.The implementation procedure of geodesic dilation is: utilize basic isotropic structure element to make dilation operation to marking image, require in calculating process that the result images obtained must remain under mask image, namely mask image plays a part limit markers image expansion and spreads.
In like manner, geodetic corrosion requires that mask image must be less than or equal to marking image, its implementation procedure is: utilize basic isotropic structure element to make erosion operation to marking image, require in calculating process that the result images obtained must remain on mask image, namely mask image plays a part the contraction of limit markers image.
The geodesic dilation of bounded image or geodetic corrode conversion, through the circulation always convergence of certain number of times, namely until marking image expansion or contraction are stoped by mask image completely; Once, the value of arbitrary pixel of marking image no longer changes in now recycle, from marking image, carry out morphological reconstruction just based on this principle to mask image.From marking image f, carry out expansion reconstruction to mask image g (f≤g) be expressed as it is defined as f relative to the geodesic dilation of g until stable:
R g δ ( f ) = δ g ( i ) ( f ) - - - ( 3 )
Wherein i is the number of times of Shi Xunhuan.In like manner, from marking image f, carry out corrosion reconstruction to mask image g (f>=g) to be expressed as it is defined as the geodetic corrosion of f relative to g until stable:
R g ϵ ( f ) = ϵ g ( i ) ( f ) - - - ( 4 )
Wherein i is the number of times of Shi Xunhuan.
The high original image of certain contour line data using said method to obtain is as shown in Fig. 2 (b), (c).Further, described step 3 specifically comprises the following steps:
(1) odd, even isocontour obstacle Euclidean Distance Transform
On contour map picture, connect each pixel to be obtained to the length in level line on it and lower isocontour geodetic path by two obstacle Euclidean Distance Transforms: the isocontour range conversion of odd number that (1) is barrier with even number level line, as shown in Fig. 2 (d); (2) with the isocontour range conversion of even number that odd number level line is barrier, as shown in Fig. 2 (e).
Range conversion calculates also identifier space point to concentrate each point to arrive conversion or the process of reference body distance, is divided into Euclidean Distance Transform and taxi distance conversion.What the former identified is Euclidean distance, is applicable to the figure of natural form; What the latter identified is manhatton distance, and be applicable to the figure of conformation of rules, Euclidean Distance Transform is refered in particular in the range conversion in the present invention.
Euclidean space with obstacle G (Ω) refers in Euclidean space to there is barrier set (as shown in Figure 4), that is: Ω={ w 1, w 2..., w k, wherein, w ketc. being barrier figure, and mutually disjoint.If P1, P2 are two points in plane.We are from one of them point, at the most only with Ω={ w 1, w 2..., w kborder adjoin, and arrive another and put the bee-line that obtains, be called the penetrate thing set omega={ w between P1, P2 under Euclidean space with obstacle 1, w 2..., w kdistance, referred to as the distance of obstacle between P1, P2, be denoted as D b(P1, P2).
Visible, when the image space outside geodetic mask is designated as barrier, obstacle Euclidean Distance Transform is carried out to pending a certain entity, just can obtain the geodesic distance of arbitrfary point in this entity to mask, as shown in Figure 4.In Fig. 4 (a), black represents barrier, and white represents current entity, and grey represents geodetic mask; In Fig. 4 (b), other colors represent the geodesic distance size of distance entity.
In the present invention, the techniqueflow of Euclidean Distance Transform is as follows:
1. the distance value of all spatial point being composed is an enough large positive number M 1, in addition the distance value of all entity points being composed is 0.
2. each grid of sequential access, also namely line number from 0,1,2 ... access successively, row number also from 0,1,2 ... access successively, and rewrite the square value of each point according to the following formula:
de 2(i,j)=min(de 2 1(i,j),de 2 2(i,j),de 2 3(i,j),de 2 4(i,j),de 2(i,j)) (5)
3. backward accesses each grid, and rewrites the square value of each grid point:
de 2(i,j)=min(de 2 5(i,j),de 2 6(i,j),de 2 7(i,j),de 2 8(i,j),de 2(i,j)) (6)
4. each point square distance value is transferred to distance value:
c(i,j)=[de 2(i,j)] 1/2(7)
In the present invention, the techniqueflow of obstacle Euclidean Distance Transform can with reference to the techniqueflow of above-mentioned Euclidean Distance Transform, but need to do following change:
1. at former algorithm, " 1. the distance value of all spatial point being composed is an enough large positive number M 1, in addition the distance value of all entity points being composed is 0." basis on increase: the distance value of all obstacle object points is composed as another enough big integer M 2, and establish M 2>>M 1.
2. at original algorithm 2., 3., 4. in step to have a square distance value de 2(i, j) is M 2point all do not process, namely do not rewrite their square distance value yet.
3. original algorithm 3., 4. all de in step 2the square distance value of (i, j) neighborhood middle distance initial point is M 2point all do not participate in comparing.
(2) linear interpolation along streamline is carried out
Along from upper level line C up () is to lower level line C l(p) and cross p geodetic path, utilize linear interpolation to calculate the interpolation of each pixel p.As shown in Fig. 2 (d), (e), the pixel p be positioned in figure between upper and lower level line corresponds to connection two isocontour geodetic path and can be determined by following method: from upper level line C up () is to lower level line C l(p) and cross the shortest path of p by from p to C l(p) and p to C up the geodetic path composition of (), this shortest path is also " streamline " or be called " line of steepest descent ".
Usually upper level line C is supposed in earth science research u(p) and lower level line C lp the difference in height between () is equally distributed along geodetic path, the interpolation H (p) of p can by C uthe height value P of (p) u(p) and C lthe height value P of (p) lp () weighting obtains:
H ( p ) = P u ( p ) * d M [ ( p , C l ( p ) ] + P l ( p ) * d M [ ( p , C u ( p ) ] d M [ ( p , C l ( p ) ] + d M [ ( p , C u ( p ) ] - - - ( 8 )
Wherein, d m[(p, C u(p)] be from p to C uthe geodesic distance of (p), d m[(p, C l(p)] be from p to C lthe geodesic distance of (p).Fig. 2 (f) illustrates the DEM generated through process of the present invention by the level line of Fig. 2 (a).
Needs particularly point out, for contour line data interpolation, MADEM side's ratio juris that interpolating method and Hu Peng etc. that the present invention generates DEM propose is on all four: along the linear interpolation of line of steepest descent, and the interpolating method that the present invention generates DEM is only realize identical object from another approach, and only have part point (refer to " and to again minute unnecessarily time, with 1/2 n+1voronoi figure is boundary, and corresponding elevation is composed in layering " those particular points) the vertical accuracy precision of elevation that may obtain higher than MADEM method.
The present invention creatively proposes the interpolating method of a kind of contour lines creation DEM based on morphological reconstruction and range conversion, especially a kind of upper plateau based on morphological reconstruction, lower plateau image generating method, a kind of geodetic path based on obstacle Euclidean Distance Transform and geodesic distance generation method is proposed, with the DEM utilizing the interpolation of contours to generate high precision, high-fidelity.
Interpolating method is the key factor affecting digital elevation model precision, obtains the process of DEM from level line, and the selection of interpolating method is a very important sport technique segment, and it has important impact to the precision generating DEM.In the present invention, morphological reconstruction is for obtaining higher level line corresponding to arbitrary pixel and lower isocontour height value, obstacle Euclidean Distance Transform is for obtaining this point to upper and lower two isocontour geodesic distances, linear interpolation is used to obtain the height value of this point, suitable tight of whole flow process and have science.
In order to verify the validity of the interpolating method of the contour lines creation DEM based on morphological reconstruction and range conversion, two scape experimental datas are used to carry out the interpolation experiment of DEM.Experimental data have employed the contour line data of domestic Mount Taishan, Tai'an, Shandong Province 2 the regional area 1:1 ten thousand in mountain range, and contour interval is 10 meters.Meanwhile, in order to verify both differences, using " MADEM " interpolating method of Hu Peng and interpolating method of the present invention to generate resolution respectively for first experiment data is the DEM of 1m, and second experimental data only have employed the interpolating method of the present invention's proposition.
The isocontour elevation of first experiment data is distributed between 630-920m, long and is widely respectively 1400m × 1100m.Use the DEM of the inventive method acquisition as shown in Fig. 2 (f), its resolution is 1m.In addition, DEM generates the level line that contour interval is 10m thus, sees Fig. 2 (g), and this level line can accurately approach original level line, this demonstrate that the correctness of this method and generates the high precision of DEM.Meanwhile, employ " MADEM " interpolating method of Hu Peng, respectively through 10 times and 30 iteration, generate two DEM, as shown in Fig. 5 (a), (b).Fig. 5 illustrates, use the DEM that the interpolating method of Hu Peng obtains, the number of times of iteration is more, more approaches the result that the present invention generates the interpolating method acquisition of DEM.
In order to compare the quality of the DEM that two kinds of interpolating methods obtain, reducing computing is carried out to the DEM obtained, as shown in Figure 6.As seen in Figure 6, the DEM of two kinds of interpolating method acquisitions is in most of pixel, and height value is equal, but on partial pixel, there is certain elevation difference, and maximum difference about 0.3 meter, larger discrepancy is mainly distributed in away from isocontour position.Visible, the interpolating method of Hu Peng obtains the precision of DEM in the restriction of part point by iterations; And the interpolating method that the present invention generates DEM is by the impact of this constraint condition, and the quality of the DEM generated higher than or equal the DEM that interpolating method that Hu Peng proposes obtains.The DEM quality of two kinds of method acquisitions tends to consistent reason and is, although the implementation procedure of two kinds of interpolating methods is different, two kinds of interpolating methods are inherently to adjacent two isocontour linear interpolations along streamline.
Second isocontour elevation of experimental data is distributed between 240-940m, long and is widely respectively 4400m × 3100m.The DEM of data genaration is shown in Fig. 7 (b) thus, and its resolution is 0.5m; Meanwhile, Fig. 7 (a) shows the result of strange isocontour distance of obstacle conversion.The elevation that DEM generates thus is spaced apart the level line of 10m and former level line fit situation is shown in Fig. 7 (c).Find by observing our, both goodnesses of fit are fairly good, and the situation of this and first scene is similar.The description of test of above-mentioned two scenes, the present invention generates DEM interpolating method and has very high fidelity.
Visible, the interpolation of contours method of patent of the present invention achieves the DEM of high precision, high-fidelity, compared with " MADEM " interpolating method proposed with Hu roc, when only using contour line data, the result that the DEM precision that patented method of the present invention generates equals or obtains higher than Hu roc method.But " MADEM " interpolating method that Hu Peng, Hu Hai propose needs manual intervention, is subject to the restriction of iterations, limit the engineering application of the method; And the present invention adopts DEM interpolating method not only without the need to manual intervention, and be not subject to the restriction of iterations, really achieve adaptive algorithm, all trades and professions can be applied.
The above is the preferred embodiment of the present invention, and for those skilled in the art, under the premise without departing from the principles of the invention, can also make some improvements and modifications, these improvements and modifications are also regarded as protection scope of the present invention.

Claims (7)

1. an interpolating method of contour lines creation DEM, is characterized in that, comprises the following steps:
Step one, obtain and load the vector contour line data meeting specified conditions, described specified conditions are: each level line spatially continuous distribution, can not interrupt, and odd number level line and even number contour interval occur;
Step 2, carries out morphological erosion and rebuilds and the morphological dilations reconstruction lower high original image of generation and upper high original image respectively;
Step 3, carries out obstacle Euclidean Distance Transform and obtains geodetic path and geodesic distance, and the linear interpolation carried out along streamline generates DEM.
2. the interpolating method of a kind of contour lines creation DEM according to claim 1, it is characterized in that, in step, level line thematic data is separated to obtain the contour line data of vector form by the level line in digitizing papery topomap or from digital line is drawn, and give corresponding height value to each level line, require each level line spatially continuous distribution and must not interrupting without reason; In the level line obtained, level line is odd number level line in any two adjacent contours, another level line is even number level line, and odd number level line and even number contour interval occur.
3. the interpolating method of a kind of contour lines creation DEM according to claim 1, it is characterized in that, in step 2, rasterizing operation generation level line image is carried out to vector level line, to level line image carry out morphological erosion reconstruction and morphological dilations rebuild generate one respectively under high original image and upper high original image.
4. the interpolating method of a kind of contour lines creation DEM based on morphological reconstruction and range conversion according to claim 1, it is characterized in that, in step 3, obtain the length to level line on it and lower isocontour geodetic path of any one pixel on level line image by the isocontour range conversion of the odd number being barrier with even number level line with the isocontour range conversion of even number two obstacle Euclidean Distance Transforms that odd number level line is barrier respectively, the linear interpolation carried out along streamline generates DEM.
5. the interpolating method of a kind of contour lines creation DEM according to any one of Claims 1-4, it is characterized in that, described step one specifically comprises the following steps: (1) obtains vector contour line data; (2) level line of process fracture, and give its corresponding height value to level line; (3) the odd even feature of adjacent contour height value is checked.
6. the interpolating method of a kind of contour lines creation DEM according to any one of Claims 1-4, it is characterized in that, described step 2 specifically comprises the following steps: (1) vector level line rasterizing; (2) assisted image is generated; (3) carry out morphological erosion to rebuild and morphological dilations reconstruction.
7. the interpolating method of a kind of contour lines creation DEM according to any one of Claims 1-4, it is characterized in that, described step 3 specifically comprises the following steps: (1) odd, even isocontour obstacle Euclidean Distance Transform; (2) linear interpolation along streamline is carried out.
CN201510237207.5A 2015-05-12 2015-05-12 A kind of interpolating method of contour lines creation DEM Active CN104809692B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201510237207.5A CN104809692B (en) 2015-05-12 2015-05-12 A kind of interpolating method of contour lines creation DEM

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201510237207.5A CN104809692B (en) 2015-05-12 2015-05-12 A kind of interpolating method of contour lines creation DEM

Publications (2)

Publication Number Publication Date
CN104809692A true CN104809692A (en) 2015-07-29
CN104809692B CN104809692B (en) 2019-03-12

Family

ID=53694497

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201510237207.5A Active CN104809692B (en) 2015-05-12 2015-05-12 A kind of interpolating method of contour lines creation DEM

Country Status (1)

Country Link
CN (1) CN104809692B (en)

Cited By (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN106056671A (en) * 2016-05-20 2016-10-26 西安电子科技大学 Contour line broken-line connection method
CN107449404A (en) * 2017-09-12 2017-12-08 中煤航测遥感集团有限公司 DLG collecting methods and device
CN107464272A (en) * 2017-07-19 2017-12-12 南京理工大学 The interpolation method of central diffusion type meteorological causes isopleth based on key point
CN110325820A (en) * 2017-02-24 2019-10-11 赫尔环球有限公司 Height map for indoor positioning service
CN110412617A (en) * 2019-09-06 2019-11-05 李娜 It is a kind of based on the unmanned plane rescue mode of self feed back laser radar scanning and application
CN110458954A (en) * 2019-07-31 2019-11-15 阿里巴巴集团控股有限公司 A kind of contour lines creation method, device and equipment

Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20050104884A1 (en) * 2002-04-22 2005-05-19 Kansei Iwata Digital altimetric map drawing method and device
CN102682479A (en) * 2012-04-13 2012-09-19 国家基础地理信息中心 Method for generating three-dimensional terrain feature points on irregular triangulation network
CN103868503A (en) * 2014-03-26 2014-06-18 西安煤航信息产业有限公司 Water area elevation processing method for making DEM (Dynamic Effect Model) based on INSAR (Interferometric Synthetic Aperture Radar)
CN104331870A (en) * 2014-11-25 2015-02-04 武汉大学 Structured multi-scale reconstruction method suitable for regular grid DEM

Patent Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20050104884A1 (en) * 2002-04-22 2005-05-19 Kansei Iwata Digital altimetric map drawing method and device
CN102682479A (en) * 2012-04-13 2012-09-19 国家基础地理信息中心 Method for generating three-dimensional terrain feature points on irregular triangulation network
CN103868503A (en) * 2014-03-26 2014-06-18 西安煤航信息产业有限公司 Water area elevation processing method for making DEM (Dynamic Effect Model) based on INSAR (Interferometric Synthetic Aperture Radar)
CN104331870A (en) * 2014-11-25 2015-02-04 武汉大学 Structured multi-scale reconstruction method suitable for regular grid DEM

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
吴艳兰: ""地貌三维综合的地图代数模型和方法研究"", 《中国优秀博硕士学位论文全文数据库 (硕士) 基础科学辑》 *
沈晶等: ""用形态学重建方法进行机载LiDAR数据滤波"", 《武汉大学学报·信息科学版》 *
沈晶等: ""集成距离变换和区域邻接图生成Delaunay三角网的方法研究"", 《武汉大学学报·信息科学版》 *

Cited By (10)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN106056671A (en) * 2016-05-20 2016-10-26 西安电子科技大学 Contour line broken-line connection method
CN106056671B (en) * 2016-05-20 2018-11-16 西安电子科技大学 A kind of contour connection of broken lines method
CN110325820A (en) * 2017-02-24 2019-10-11 赫尔环球有限公司 Height map for indoor positioning service
CN107464272A (en) * 2017-07-19 2017-12-12 南京理工大学 The interpolation method of central diffusion type meteorological causes isopleth based on key point
CN107449404A (en) * 2017-09-12 2017-12-08 中煤航测遥感集团有限公司 DLG collecting methods and device
CN107449404B (en) * 2017-09-12 2019-12-06 中煤航测遥感集团有限公司 DLG data acquisition method and device
CN110458954A (en) * 2019-07-31 2019-11-15 阿里巴巴集团控股有限公司 A kind of contour lines creation method, device and equipment
CN110458954B (en) * 2019-07-31 2023-06-06 创新先进技术有限公司 Contour line generation method, device and equipment
CN110412617A (en) * 2019-09-06 2019-11-05 李娜 It is a kind of based on the unmanned plane rescue mode of self feed back laser radar scanning and application
CN110412617B (en) * 2019-09-06 2023-06-30 李娜 Unmanned aerial vehicle rescue method based on self-feedback laser radar scanning and application

Also Published As

Publication number Publication date
CN104809692B (en) 2019-03-12

Similar Documents

Publication Publication Date Title
CN104809692A (en) Interpolation method for DEM (Digital Elevation Model) generated by contour lines
CN102122395B (en) Adaptive scale DEM (digital elevation model) modeling method capable of keeping terrain features
CN110544305B (en) Terrain steeple ridge line information fusion method constructed by DEM (digital elevation model) oriented to regular grids
CN111028335B (en) Point cloud data block surface patch reconstruction method based on deep learning
CN108733711A (en) Distribution line space length acquisition methods based on three-dimension GIS technology
CN104851133A (en) Image self-adaptive grid generation variational method
CN113963050A (en) Method and system for calculating earth volume based on point cloud
CN107492130B (en) Non-continuous equivalent cloud picture generation method based on layer rendering
CN115795402B (en) Variational method-based multi-source precipitation data fusion method and system
CN117171855A (en) Hilly area flow field model modeling method based on Delaunay triangulation
Xu et al. Methods for the construction of DEMs of artificial slopes considering morphological features and semantic information
CN106156376A (en) A kind of UHV transmission line electromagnetic interference Three-dimensional Simulation System
CN113012286A (en) Method for constructing road DEM (digital elevation model) based on dense point cloud data
CN114969944A (en) High-precision road DEM construction method
Liao et al. Topological relationship-based flow direction modeling: stream burning and depression filling
CN107045488A (en) A kind of surface boundary member computational methods based on Bezier curved surface quadrangle parametric equations
Miao et al. Numerical calculation of channel dredging volume using 3D digital stratum model
Abdel-Bary Ebrahim et al. Planning road networks in new cities using GIS: The case of New Sohag, Egypt
CN112950778B (en) Construction method and system of rasterized curved surface for measuring landform wrinkles
Xu et al. 3D GIS based long-distance project construction progress visualization
CN107992702B (en) Method for converting various BIM model data into DEM data
Zhang et al. Research on the calculation method of engineering earth volume under complex terrain conditions of Yellow River dam bank
Hovad et al. Slope based Grid Creation using Interpolation of LIDAR Data Sets.
Chen et al. Design and implementation of flight simulation visual database based on worldwide database of Mantis
Zhou et al. New high-fidelity terrain modeling method constrained by terrain semanteme

Legal Events

Date Code Title Description
C06 Publication
PB01 Publication
EXSB Decision made by sipo to initiate substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant