A three-dimensional terrain modeling method and system

Abstract

The present application provides a three-dimensional terrain modeling method and system, the modeling method comprising: obtaining aerial images of a plurality of unmarked ground control points of an area to be modeled; generating a first sparse point cloud set and a rough three-dimensional terrain model of the area to be modeled based on the aerial images of each unmarked ground control point; automatically marking ground control points on each aerial image based on the rough three-dimensional terrain model; generating a second sparse point cloud set of the area to be modeled based on the aerial images of each marked ground control point; and generating a dense point cloud set and a refined three-dimensional terrain model of the area to be modeled based on the second sparse point cloud set. The modeling method provided by the present application effectively improves the marking speed while improving the marking accuracy of ground control points, significantly improving the efficiency of high-precision three-dimensional terrain modeling overall.

Description

Technical Field

[0001] The present application belongs to the field of aerial photogrammetry technology, relates to three-dimensional modeling technology based on aerial images, and specifically provides a three-dimensional terrain modeling method and system. Background Art

[0002] The 3D terrain modeling method based on aerial photogrammetry is a technology that uses photographic images and related data to reconstruct a 3D terrain model. This modeling method has the advantages of high precision, high efficiency and low cost. It is currently widely used in geographic information systems (GIS), urban planning, land resources and major infrastructure project surveys, environmental and natural disaster monitoring and other fields. With the development of drone technology, the application of 3D modeling based on photogrammetry in the data collection process is becoming more and more common, making real-time and large-scale terrain modeling possible.

[0003] Existing 3D terrain modeling methods based on aerial photogrammetry generally involve acquiring a collection of aerial images of the area to be modeled, matching and aerial triangulating each aerial image to create a sparse point cloud of the area to be modeled, and then using a densification algorithm to obtain a dense point cloud of the area to be modeled, thereby generating a 3D terrain digital model of the area to be modeled. To improve modeling accuracy, several ground control points (GCPs) can be deployed in the area to be modeled and marked on the aerial images. This provides geometric constraints for the iterative optimization of the point cloud coordinates during the subsequent aerial triangulation phase. Obviously, the ability to accurately calibrate the positions of the GCPs on the aerial images will significantly affect the measurement and modeling results.

[0004] At present, the main method for marking ground control points on aerial images is manual processing, that is, manually selecting each ground control point on each aerial image and establishing an association between the same ground control point in different aerial images. This operation not only consumes a lot of time, but also significantly reduces the modeling speed. In addition, there are also methods for automatically identifying ground control points through texture information on aerial images. However, the existing algorithms do not establish a matching relationship between different aerial images when automatically identifying ground feature points. As a result, even if the images of each ground feature point can be extracted from each aerial image through texture features, the same ground feature points on different aerial images cannot be associated. Generally, manual specification of the association between the ground control points on each aerial image is required. Otherwise, different ground control points on two or more aerial images may be mistakenly identified as the same one, resulting in an adverse effect on the subsequent aerial triangulation measurement results. Summary of the Invention

[0005] In order to solve the problems existing in the above-mentioned prior art, the present application provides a three-dimensional terrain modeling method through an embodiment, which includes the following steps:

[0006] S1, obtain aerial images of several unmarked ground control points in the area to be modeled;

[0007] S2, generating a first sparse point cloud set and a rough three-dimensional terrain model of the area to be modeled based on the aerial images of each unmarked ground control point;

[0008] S3, automatically marking and associating ground control points on each aerial image based on the rough three-dimensional terrain model;

[0009] S4, generating a second sparse point cloud set of the area to be modeled based on the aerial images of each marked ground control point;

[0010] S5: Generate a dense point cloud set and a fine three-dimensional terrain model of the area to be modeled based on the second sparse point cloud set.

[0011] Compared with the existing technology, the three-dimensional terrain modeling method provided by the present application adds an operation of generating a sparse point cloud using aerial photos without marked ground control points to the conventional steps of three-dimensional terrain modeling based on aerial photogrammetry, and uses this operation to establish a rough three-dimensional terrain model, and then uses the rough three-dimensional terrain model to automatically mark and associate ground control points. Through the above steps, the same ground control points can be automatically identified and accurately associated in different aerial photos without manual operation one by one, which effectively improves the marking speed while improving the marking accuracy of ground control points, thereby significantly improving the efficiency of high-precision three-dimensional terrain modeling as a whole.

[0012] Furthermore, step S2 includes the following steps:

[0013] S21, identifying multiple feature points in each aerial image;

[0014] S22, matching feature points in different aerial images;

[0015] S23, obtaining a first estimated value of a camera orientation element of each aerial image and a first estimated value of a three-dimensional coordinate of each feature point based on an unconstrained bundle adjustment;

[0016] S24, constructing a first sparse point cloud set based on the first estimated value of the three-dimensional coordinates of each feature point;

[0017] S25: Establish a rough three-dimensional terrain model of the area to be modeled based on the first sparse point cloud set.

[0018] Furthermore, step S3 includes the following steps:

[0019] S31, selecting a ground control point as a target ground control point;

[0020] S32, selecting a reference image and several evaluation images corresponding to the target ground control point;

[0021] S33, determining at least one annotation point from the image corresponding to the target ground control point on the reference image;

[0022] S34, projecting the marked point onto the surface of the rough three-dimensional terrain model and intercepting a first projection surface element;

[0023] S35, searching for the best position and best orientation of the first projection plane element;

[0024] S36, projecting an image within a reference window on the reference picture onto a first projection plane element having an optimal position and an optimal orientation to obtain a second projection plane element;

[0025] S37 , based on the expected image obtained by reprojecting the second projection bin onto each of the to-be-marked pictures, searching and determining the imaging position corresponding to the target ground control point on each of the to-be-marked pictures.

[0026] Preferably, the similarity between the image of the target ground control point on the reference picture and its orthographic image is greater than a preset first similarity threshold, or the deviation is less than a preset first deviation threshold.

[0027] Preferably, the similarity between the image of the target ground control point on the evaluation picture and its orthographic image or its image on the reference picture is greater than a preset second similarity threshold, or the deviation is less than a preset second deviation threshold, wherein the first similarity threshold is greater than or equal to the second similarity threshold, and the first deviation threshold is less than or equal to the second deviation threshold.

[0028] Furthermore, step S34 includes the following steps:

[0029] The marked point is projected onto the surface of the rough three-dimensional terrain model along a reference projection line corresponding to the reference image to obtain its projection point on the surface of the rough three-dimensional terrain, wherein the reference projection line of the marked point corresponding to the reference image is determined based on the camera orientation element of the reference image.

[0030] Furthermore, step S35 includes the following steps:

[0031] S351, reprojecting the first projection bin onto each evaluation image to obtain a reprojection window;

[0032] S352, determining the cross-correlation characteristics of the images in the reprojection windows of each evaluation picture;

[0033] S353, determining whether the mutual correlation feature satisfies the correlation criterion; if not, executing step S354 and returning to step S351; if satisfied, executing step S355;

[0034] S354, translating or rotating the first projection surface element along the reference projection line;

[0035] S355: Taking the current position and orientation of the first projection plane element as the optimal position and optimal orientation of the first projection plane element.

[0036] Furthermore, step S37 includes cyclically executing the following steps until all the aerial images that need to be marked are processed:

[0037] S371, selecting an aerial image of the area to be modeled as an image to be labeled;

[0038] S372, reprojecting the second projection bin onto the image to be labeled to obtain a desired image;

[0039] S373, extracting the edge of the desired image to construct a matching window;

[0040] S374 , searching for the best position of the matching window on the image to be marked, thereby marking the target ground control point on the image to be marked.

[0041] Preferably, in the process of searching for the best position of the matching window on the image to be marked, the matching window is only translated on the image to be marked.

[0042] Furthermore, step S4 includes the following steps:

[0043] S41, identifying multiple feature points in the aerial image of each marked ground control point and performing feature point matching;

[0044] S42, obtaining a second estimated value of a camera orientation element of each aerial image and a second estimated value of a three-dimensional coordinate of each feature point by a bundle adjustment with constraints, wherein the constraints include position constraints of the marked ground control points;

[0045] S43: Construct a second sparse point cloud set based on the second estimated value of the three-dimensional coordinate of each feature point.

[0046] The present application also provides a three-dimensional terrain modeling system, including a processor and a readable storage medium, wherein the readable storage medium stores an executable program, and when the executable program is executed by the processor, the aforementioned three-dimensional terrain modeling method can be implemented. BRIEF DESCRIPTION OF THE DRAWINGS

[0047] Figure 1A flowchart of an existing method for 3D terrain modeling based on aerial images;

[0048] Figure 2 A flowchart of a three-dimensional terrain modeling method provided according to an embodiment of the present application;

[0049] Figure 3A is a schematic diagram of an aerial image of an area to be modeled in some embodiments;

[0050] Figure 3B This is an enlarged aerial photo;

[0051] Figure 4 Flowchart for implementing step S2 according to an embodiment of the present application;

[0052] Figure 5 A schematic diagram of the principle of iteratively estimating camera orientation elements and three-dimensional coordinates of feature points using an unconstrained bundle adjustment method according to an embodiment of the present application;

[0053] Figure 6 This is a flowchart for implementing the iterative estimation of camera orientation elements and three-dimensional coordinates of feature points using the unconstrained bundle adjustment method according to an embodiment of the present application;

[0054] Figure 7 Flowchart for implementing step S3 according to an embodiment of the present application;

[0055] Figure 8A A schematic diagram of a reference image provided according to an embodiment of the present application;

[0056] Figure 8B A schematic diagram of an evaluation image provided according to an embodiment of the present application;

[0057] Figure 9 A schematic diagram of intercepting a first projection surface element on a rough three-dimensional model according to an embodiment of the present application;

[0058] Figure 10 A flowchart of searching for the optimal position and optimal orientation of a target ground control point according to an embodiment of the present application;

[0059] Figure 11 A schematic diagram of an operation for searching for the optimal position and optimal orientation of a target ground control point according to an embodiment of the present application;

[0060] Figure 12 Schematic diagram showing how the NCC coefficient of an image in a reprojection window varies with the size of the first projection bin according to an embodiment of the present application;

[0061] Figure 13 Flowchart for implementing step S37 according to an embodiment of the present application;

[0062] Figure 14 This is a schematic diagram of an operation of searching for the optimal position of a target ground control point matching window on a to-be-marked image according to an embodiment of the present application;

[0063] Figure 15 Schematic diagram of an operation for searching for an optimal position of a matching window on a to-be-marked image partially displaying target ground control points according to an embodiment of the present application;

[0064] Figure 16 A schematic diagram of a first sparse point cloud of an area to be modeled according to an embodiment of the present application;

[0065] Figure 17 A schematic diagram of a rough three-dimensional terrain model of an area to be modeled according to an embodiment of the present application;

[0066] Figure 18 A schematic diagram of a dense point cloud of an area to be modeled according to an embodiment of the present application;

[0067] Figure 19 A schematic diagram of a fine three-dimensional terrain model of an area to be modeled according to an embodiment of the present application;

[0068] Figure 20 A schematic diagram of a three-dimensional model of a region to be modeled according to an embodiment of the present application;

[0069] Figure 21 Schematic diagram of the architecture of a three-dimensional terrain modeling system provided according to an embodiment of the present application. DETAILED DESCRIPTION

[0070] Hereinafter, the present application will be further described based on preferred embodiments with reference to the accompanying drawings.

[0071] In addition, various components in the drawings are enlarged or reduced in size for ease of understanding, but this is not intended to limit the scope of protection of this application.

[0072] In the description of the embodiments of the present application, it should be noted that if the terms "upper", "lower", "inner", "outer" and the like indicate an orientation or positional relationship based on the orientation or positional relationship shown in the accompanying drawings, or are the orientation or positional relationship in which the products of the embodiments of the present application are usually placed when in use, they are only for the convenience of describing the present application and simplifying the description, and do not indicate or imply that the device or element referred to must have a specific orientation, be constructed and operated in a specific orientation, and therefore cannot be understood as a limitation on the present application. In addition, in the description of the present application, in order to distinguish different units, words such as first and second are used in this specification, but these are not limited by the order of manufacture, nor can they be understood as indicating or implying relative importance. Their names may be different in the detailed description and claims of the present application.

[0073] The vocabulary in this specification is used to illustrate the embodiments of the present application, but is not intended to limit the present application. It should also be noted that, unless otherwise clearly specified and limited, the terms "disposed", "connected", and "connected" should be understood in a broad sense. For example, it can be a fixed connection, a detachable connection, or an integral connection; it can be a mechanical connection, a direct connection, an indirect connection through an intermediate medium, or a communication between the two components. For those skilled in the art, the specific meanings of the above terms in this application can be specifically understood.

[0074] Figure 1 A conventional three-dimensional terrain modeling method based on aerial photogrammetry is presented. This method first obtains a set of aerial images of the area to be modeled, consisting of multiple aerial images that meet a certain overlap ratio, through aerial photography or other means. The captured ground control points (GCPs) are then marked on the aerial images. Furthermore, a sparse point cloud of the area to be modeled is obtained by matching and aerial triangulation of the individual aerial images. Finally, a dense point cloud of the area to be modeled is obtained using a densification algorithm, and a three-dimensional terrain digital model of the area to be modeled is generated from the dense point cloud. In this three-dimensional terrain modeling process, ground control points serve as basic data, providing geometric constraints for the iterative optimization of point cloud coordinates during the aerial triangulation phase. Clearly, the ability to accurately calibrate the positions of each ground control point in the aerial images and accurately associate the same ground control points in different images will significantly affect the measurement and modeling results.

[0075] The current mainstream method of manually marking ground control points requires manually selecting each ground control point on each aerial image and establishing associations between the same ground control point in different aerial images. This operation not only consumes a lot of time, but also significantly reduces the modeling speed.

[0076] In addition, there are also algorithms that automatically identify ground control points based on texture information in aerial images. However, the applicant found that the above-mentioned automatic recognition algorithm based on texture features of ground control point images has the following problems:

[0077] First, the same ground control point may have significant shape distortion and image area differences in aerial images taken from different perspectives. For example, a ground control point that appears square in one aerial image (taken from orthographic projection) may be compressed into a very thin rectangle in another aerial image (where the viewing angle is almost perpendicular to the normal direction of the ground control point). Due to the large difference in texture information between the two, neither correlation nor deep learning can resolve this large discrepancy.

[0078] Secondly, when a ground control point is imaged at the edge of an aerial image, the image may only contain a partial image of the ground control point. In this case, using conventional image recognition methods may result in incorrect recognition due to the lack of complete image information.

[0079] More importantly, the image of a ground control point in a certain aerial image may be wrongly associated with the image of another ground control point in other aerial images for some reason. The reason for this problem is that when using, for example Figure 1 In the existing modeling method shown, when marking the ground control points, the camera orientation information of each aerial image is not limited, resulting in matching based only on texture information. In this case, image distortion caused by perspective deformation or terrain structure factors of the shooting area may lead to matching errors. For example, there are some relatively similar buildings that are arranged periodically, causing the recognition algorithm to identify two different ground control points with different spatial positions but similar texture features in aerial images taken at two different locations as the same one.

[0080] It can be seen that if the matching relationship between different aerial images has not been established during the automatic identification of ground feature points, then even if the images of each ground feature point can be extracted from each aerial image through texture features, the same ground feature points in different aerial images cannot be associated. Generally, manual specification of the association between the ground control points in each aerial image is required. Otherwise, different ground control points in two or more aerial images may be mistakenly identified as the same one, which will have an adverse effect on the subsequent aerial triangulation measurement results.

[0081] The root cause of the above problems is that the existing automatic recognition algorithm lacks the necessary posture information when identifying ground control points. It is obviously impossible to fundamentally solve the problem simply by optimizing the feature extraction algorithm or deep learning network structure, or increasing the training level of the deep learning network.

[0082] In order to solve the problems existing in the above-mentioned prior art, this application improves the ground control point annotation and association method in the existing modeling process, thereby providing a new three-dimensional terrain modeling method. Figure 2 Flowchart of a three-dimensional terrain modeling method according to an embodiment of the present application, as shown in FIG. Figure 2 As shown, the method includes the following steps:

[0083] S1, obtain aerial images of several unmarked ground control points in the area to be modeled;

[0084] S2, generating a first sparse point cloud set and a rough three-dimensional terrain model of the area to be modeled based on the aerial images of each unmarked ground control point;

[0085] S3, automatically marking and associating ground control points on each aerial image based on the rough three-dimensional terrain model;

[0086] S4, generating a second sparse point cloud set of the area to be modeled based on the aerial images of each marked ground control point;

[0087] S5: Generate a dense point cloud set and a fine three-dimensional terrain model of the area to be modeled based on the second sparse point cloud set.

[0088] By comparison Figure 1 and Figure 2 It can be seen that compared with the existing technology, the three-dimensional terrain modeling method provided by the present application first uses aerial photos without marked ground control points (GCPs) to generate a sparse point cloud in step S2, and uses this operation to establish a rough three-dimensional terrain model. Then, in step S3, the rough three-dimensional terrain model is used to automatically mark and associate the ground control points. Through this step, there is no need for manual marking and association of each image, and the association relationship between the same ground control points in different aerial photos can be automatically identified and accurately established, thereby providing accurate constraints for subsequent aerial image registration and a second sparse point cloud generation operation through aerial triangulation.

[0089] The above steps are described in detail below with reference to the accompanying drawings.

[0090] <Obtain aerial images of the area to be modeled>

[0091] Step S1 is used to obtain an initial aerial picture of the area to be modeled, that is, an aerial picture without ground control point marking.

[0092] Figure 3A In a specific embodiment, several original aerial images obtained by aerial photography of an area to be modeled are shown. The aerial image collection comes from the aerial image collection (Dort-ZECHE dataset) of a building as a designated area published by the International Society for Photogrammetry and Remote Sensing (ISPRS).

[0093] Figure 3B One of the aerial photos is enlarged and displayed. Figure 3B As shown in the red circles, different aerial images obtained through aerial photography may contain images of 0, 1, 2 or more ground control points.

[0094] During aerial photogrammetry and 3D terrain modeling, ground control points can provide effective geometric constraints using their precise location information. Generally, they can be made of sheet materials or soft fabrics, with their shape and size determined according to the aerial photography altitude and shooting resolution (for example, a square of 50cm×50cm or 1m×1m). Patterns with obvious texture features and high contrast are printed or painted on the surface.

[0095] Furthermore, before conducting aerial photography of the area to be modeled, several ground control points can be placed at a certain density on the ground, buildings or natural objects in the area to be modeled according to the established aerial photography plan, and the spatial position of each ground control point can be measured with high precision using a total station, GPS, etc.

[0096] After completing the above settings, aerial photography of the area to be modeled can be carried out according to the pre-established aerial photography plan. The image data obtained from aerial photography can be multiple aerial pictures or aerial videos taken continuously. The above aerial images generally need to be converted according to the requirements of subsequent data processing, such as unifying the picture size, or converting continuous videos into multiple aerial pictures. The number of aerial pictures is generally determined by the scale of the area to be modeled, the area covered by each aerial picture, and the overlap of adjacent aerial pictures. The overlap between consecutive aerial pictures is usually set between 60% and 80% to ensure the reliability of subsequent feature matching.

[0097] <Establish a rough 3D terrain model of the area to be modeled>

[0098] Step S2 is used to establish a rough three-dimensional terrain model of the area to be modeled. This rough three-dimensional terrain model is used to provide "indexes" and constraints for the imaging position, shape, and texture features of each ground control point in different aerial images in the subsequent step S2, so that the same ground control points can be identified and associated in different aerial images without manual one-by-one association.

[0099] Figure 4 The specific implementation process of step S2 in some preferred embodiments is shown as follows: Figure 4 As shown, the rough three-dimensional terrain model can be established by the following steps:

[0100] Step S21, identifying multiple feature points in each aerial image;

[0101] Step S22, matching feature points in different aerial images;

[0102] Step S23, obtaining a first estimated value of a camera orientation element of each aerial image and a first estimated value of a three-dimensional coordinate of each feature point based on an unconstrained bundle adjustment;

[0103] Step S24, constructing a first sparse point cloud set based on the first estimated value of the three-dimensional coordinates of each feature point;

[0104] Step S25 : establishing a rough three-dimensional terrain model of the area to be modeled based on the first sparse point cloud set.

[0105] Specifically, first in step S21, feature points can be extracted from each aerial picture through a feature detection algorithm (such as SIFT, SURF or ORB, etc.); then in step S22, feature points of different aerial pictures can be matched. This step can use methods such as nearest neighbor search to match the same feature points in aerial pictures of different perspectives.

[0106] After completing the feature point matching of the aerial image, the bundle adjustment operation can be performed in step S23 to estimate the spatial position of each feature point. It should be emphasized that in the method provided in this application, since the ground control points are not marked and associated on the aerial image before executing step S23, this step adopts the bundle free network adjustment that does not use position coordinates for constraints.

[0107] Unconstrained bundle adjustment refers to a method for estimating the camera orientation elements (including internal orientation elements such as the camera's focal length, principal point, and distortion coefficient, as well as external orientation elements such as the camera's position and orientation) and the three-dimensional coordinates of each feature point of an aerial image without relying on external known control points.

[0108] Figure 5 and Figure 6 The principles and specific processes of iteratively estimating the three-dimensional coordinates of camera orientation elements and feature points using the unconstrained bundle adjustment method in some specific embodiments are respectively shown. Figure 6As shown, in step S231, the internal and external orientation elements of the camera can be initially estimated through the description document of the camera parameter specifications and the flight position and posture records of the aircraft during aerial photography, and then the three-dimensional coordinates of each feature point can be calculated based on the camera orientation elements and feature point matching results of each aerial picture through the triangulation method (step S232).

[0109] Then, in step S233, the first objective function value is calculated using the three-dimensional coordinates of each feature point and the camera orientation elements of each aerial image (generally, the objective function of optical adjustment is a robust function of the sum of squares of mapping errors, that is, the error is minimized after each feature point is mapped to the aerial image), and in step S234, it is determined whether the function value converges or the number of iterations reaches a preset upper limit. If the judgment result is no, the camera orientation elements are re-estimated and the three-dimensional coordinates of each feature point are recalculated through step S235, and then the first objective function value is recalculated. If the judgment result is yes, the first estimated value of the camera orientation element of each aerial image and the first estimated value of the three-dimensional coordinates of each feature point are determined based on the current iteration result (step S236).

[0110] After obtaining the first estimated value of the three-dimensional coordinates of each feature point, each feature point can be used as a point cloud point in step S24, and its position in the three-dimensional space can be set according to the first estimated value, thereby obtaining a first sparse point cloud set. Then, in step S25, the above-mentioned feature points in the sparse point cloud set are used as grid points, and algorithms such as Delaunay triangulation are used to construct a rough three-dimensional terrain model of the area to be modeled. Generally, the rough three-dimensional terrain model constructed in this step is composed of several small patches, and the vertex of each patch is a feature point.

[0111] In some preferred embodiments, the camera orientation elements of each aerial image obtained through the above steps, the sparse point cloud set of the area to be modeled, and the rough three-dimensional terrain model are all stored in a database for use in subsequent steps.

[0112] <Automatically mark and associate ground control points>

[0113] It is particularly important to point out that in conventional three-dimensional terrain modeling methods, since the spatial resolution of the sparse point cloud composed of feature points is too low, the three-dimensional model generated therefrom is relatively rough. Therefore, the three-dimensional grid model is generally not generated directly from the sparse point cloud. Instead, the sparse point cloud is densified and then the dense point cloud is used to reconstruct the three-dimensional terrain. However, in the embodiment of the present application, since the surface of the rough three-dimensional terrain model reflects the approximate three-dimensional shape of the terrain of the area to be modeled, that is, the position of each facet of its surface can roughly characterize the orientation and spatial coordinates of the area to be modeled at that location, it can be regarded as a reference for finding ground control points on each aerial image and establishing the correlation between the same ground control point on each aerial image, so that in step S3, the ground control points are automatically marked and associated with the help of the rough three-dimensional terrain model.

[0114] Figure 7 In some specific embodiments, the specific implementation steps of step S3 are shown, such as Figure 7 As shown, step S3 further includes the following steps:

[0115] Step S31, selecting a ground control point as a target ground control point;

[0116] Step S32, selecting a reference image and several evaluation images corresponding to the target ground control point;

[0117] Step S33, determining at least one annotation point from the image corresponding to the target ground control point on the reference image;

[0118] Step S34, projecting the marked points onto the surface of the rough three-dimensional terrain model and intercepting the first projection surface element;

[0119] Step S35, searching for the best position and best orientation of the first projection plane element;

[0120] Step S36, projecting the image within a reference window on the reference picture onto the first projection plane element having the optimal position and optimal orientation to obtain a second projection plane element;

[0121] Step S37 : Based on the expected image obtained by reprojecting the second projection bin onto each of the to-be-marked pictures, searching and determining the imaging position corresponding to the target ground control point on each of the to-be-marked pictures.

[0122] In some specific embodiments, when there are multiple ground control points that need to be automatically marked and associated, such as Figure 7 As shown, steps S31 to S37 may be executed cyclically until all operations on ground control points that need to be automatically identified and associated are completed.

[0123] The implementation of step S3 is described in detail below with reference to the accompanying drawings.

[0124] First, in step S31, a ground control point is selected as the target ground control point that needs to be automatically marked and associated on each picture to be marked. Then, based on the imaging conditions of the target ground control point on each aerial picture, a reference picture is selected, and several evaluation pictures are selected around it.

[0125] In order to ensure the accuracy of subsequent search results, the reference image should be selected so that the target ground control point is photographed at an orthographic angle as much as possible. That is, the aerial image whose image corresponding to the target ground control point is close to its orthographic shape should be selected as the reference image. For example, when the actual shape of the target ground control point is a square, the aerial image that captures the target ground control point and whose imaging shape is also basically a square should be selected from the aerial image collection.

[0126] In some preferred embodiments, the deviation or similarity can be used to evaluate the deviation or similarity between the shape of the image formed by a ground control point in the aerial image and the shape it should have when photographed according to the orthographic projection perspective. For example, the Hausdorff distance, average distance, or similarity ratio can be used to measure the geometric difference between the image and the ideal square. For another example, the difference between the length of each side (or perimeter) of the image formed by the GCP in the aerial image and the length of each side (or perimeter) of a standard square with the same area can be used as a deviation indicator:

[0127]

[0128] Among them, D is the deviation index, l i is the length of each side of the image formed by a GCP on the aerial picture, l square is the side length of a square of equal area to the image.

[0129] Obviously, a preset similarity threshold (first similarity threshold) or deviation threshold (first deviation threshold) can be set. When the similarity between the shape of the image formed by a ground control point on an aerial picture and the shape of its orthographic image is greater than the first similarity threshold, or the deviation is less than the first deviation threshold, the aerial picture is used as a reference picture with the ground control point as the target ground control point.

[0130] After determining the benchmark image of the target ground control point, several aerial images can be selected "around" the benchmark image as evaluation images. Here, "around" means that the evaluation images have a camera position and perspective similar to those of the benchmark image. Since aerial images are generally numbered consecutively according to the order in which they were taken, the evaluation images in the aerial image collection are generally several images that are adjacent to and before the benchmark image.

[0131] Figure 8A and Figure 8B Schematic diagrams of a reference image and its surrounding evaluation images selected from a plurality of aerial images after determining a target ground control point to be marked in a specific embodiment are respectively shown.

[0132] like Figure 8A As shown in , the evaluation picture may include a picture of the target ground control point taken at a nearly orthographic perspective, or a picture of the target ground control point taken at a more deviated perspective than the reference picture. Similarly, Figure 8B As shown, the evaluation pictures can also be screened according to the similarity or deviation between the shape of the image formed by the target ground control point on it and the shape of the orthographic image of the ground control point (or the shape of the image formed on the reference picture), that is, a second similarity threshold or a second deviation threshold can be preset, and only aerial pictures whose similarity with the shape of the orthographic image of the target ground control point (or the shape of the image formed on the reference picture) is greater than the second similarity threshold or the deviation is less than the second deviation threshold are selected as evaluation pictures.

[0133] Preferably, the image of the target ground control point on the reference image is closer to its orthographic image than the image of the surrounding evaluation image. Therefore, in some preferred embodiments, the first similarity threshold is greater than or equal to the second similarity threshold; in other preferred embodiments, the first deviation threshold is less than or equal to the second deviation threshold.

[0134] In some preferred embodiments, the number of evaluation pictures is not less than 6. For example, 3 aerial pictures can be selected in sequence on both sides of the reference picture as 6 evaluation pictures. In other preferred embodiments, the number of evaluation pictures is not more than 10. For example, 5 aerial pictures can be selected in sequence on both sides of the reference picture as 10 evaluation pictures.

[0135] Step S33 is used to mark the annotation point corresponding to the target ground control point on the reference image. The marking in this step can be performed manually, for example, manually selecting a pixel in the image of the target ground control point on the reference image (preferably the pixel corresponding to the center of the target ground control point) and using the pixel as the annotation point of the target ground control point on the reference image; or, it can also be performed semi-automatically in a human-computer interactive manner, for example, first manually defining a rough range, using an automatic recognition algorithm to identify the image of the ground control point therein and marking one of the pixels as the annotation point, and then manually adjusting the position of the annotation point.

[0136] It should be noted that although the operation of marking the target ground control point in step S33 occurs after the selection of the reference picture in step S32, the preferred principle for selecting the reference picture is whether a ground control point that is desired to be automatically marked can be imaged in an aerial picture in an orthographic or near-orthographic direction. The purpose is to facilitate accurate marking of the ground control point in the aerial picture. Therefore, for different target ground control points, different aerial pictures may need to be selected as their reference pictures.

[0137] After the target ground control points are accurately marked on the reference image, in step S34, Figure 9 As shown, first, the marked point of the target ground control point is projected onto the rough three-dimensional terrain model of the area to be modeled along its projection line corresponding to the reference image (in the embodiment of the present application, the projection line corresponding to the marked point on the reference image is called the reference projection line), wherein the reference projection line of the marked point can be determined by the camera orientation element of the reference image obtained in step S2 and stored in the database.

[0138] Camera orientation elements play an important role in fields such as photogrammetry, computer vision, and image processing. They include the camera's intrinsic orientation elements and extrinsic orientation elements. The intrinsic orientation elements describe the camera's internal parameters, which are related to the camera's geometric structure and optical properties. These parameters mainly include focal length, principal point, radial distortion coefficients, tangential distortion coefficients, etc., which are used to characterize how to project from a three-dimensional space point to a two-dimensional image plane; the extrinsic orientation elements describe the camera's spatial position and orientation during shooting, mainly including the camera's three-dimensional position coordinates and rotation angles (Roll, Pitch, Yaw or the corresponding rotation matrix / quaternion), which are used to characterize the camera's position and posture relative to the world coordinate system.

[0139] refer to Figure 5 Once the camera orientation elements of an aerial image are determined through aerial triangulation, for a three-dimensional space point, the straight line from the camera center to this space point is the projection line of the space point corresponding to the aerial image. The intersection of this projection line and the aerial image is the position of its image on the aerial image. Therefore, the projection line shows how the space point is mapped to the image through the internal and external orientation elements of the camera.

[0140] In order to clearly indicate the direction of the projection operation, in an embodiment of the present application, the operation of projecting one or several elements (such as pixels, geometric features, or texture features) from an aerial image along their corresponding projection lines onto a three-dimensional terrain model in the form of a mesh patch of the area to be modeled is called projection; correspondingly, the operation of projecting elements on the three-dimensional terrain model onto the aerial image along their projection lines is called reprojection.

[0141] As analyzed above, the rough three-dimensional model obtained by unconstrained bundle adjustment can roughly reflect the three-dimensional terrain of the area to be modeled, and the camera orientation element provides the direction of projection of each pixel on the reference image to its corresponding spatial point. Therefore, after projecting the marked point of the target ground control point onto the surface of the rough three-dimensional model according to its corresponding reference projection line on the reference image, a projection point will be obtained. The three-dimensional coordinates of the projection point and the orientation of the plane on which it is located on the rough three-dimensional model (i.e., the normal direction) represent the approximate spatial position and orientation of the target ground control point when it is arranged in the area to be modeled.

[0142] Then, if Figure 9 As shown, continue with step S34 to intercept a surface element around the projection point on the plane where the projection point is located on the rough groove three-dimensional model. In this application, this surface element is called the first projection surface element.

[0143] The first projection surface element can preferably be intercepted with the projection point as the center, or it can be adjusted accordingly according to the actual intersection of the various planes constituting the rough three-dimensional model (for example, when the surface element intercepted with the projection point as the center happens to span two facets constituting the rough three-dimensional model, the first projection surface element can be appropriately offset so that it is in the same plane). It is only necessary to ensure that the projection point is included in the intercepted surface element.

[0144] The shape and size of the first projection surface element can be selected according to the scale of the area to be modeled, the resolution parameters of the aerial image, etc. In some preferred embodiments, the first projection surface element can be in the shape of a square, with the projection point as the center of the square. In other embodiments, the shape of the first projection surface element can also be set accordingly according to the actual shape of the target ground control point, such as a rectangle, a circle, etc.

[0145] In an embodiment of the present application, the first projection surface element is equivalent to a Mask, which constructs a "window" with a limited shape and size (number of pixels) based on the approximate position and orientation information of the target ground control point reflected by the projection point on the rough three-dimensional model. Furthermore, this can be used as the starting point for the search, and an iterative search can be performed in step S35 to obtain the best estimate of the position and orientation of the target ground control point.

[0146] refer to Figure 10 as well as Figure 11 In some preferred embodiments of the present application, the search for the optimal position and orientation of the target ground control point is performed by the following steps:

[0147] Step S351 , reprojecting the first projection bin onto each evaluation image to obtain a reprojection window;

[0148] Step S352, determining the cross-correlation characteristics of the images in the reprojection windows of each evaluation picture;

[0149] Step S353, determining whether the mutual correlation feature satisfies the correlation criterion; if not, executing step S354 and returning to step S351; if satisfied, executing step S355;

[0150] Step S354, translating or rotating the first projection surface element along the reference projection line;

[0151] Step S355: taking the current position and orientation of the first projection plane element as the optimal position and optimal orientation of the first projection plane element.

[0152] Specifically, in step S351, the first projection plane is constructed according to its current position and orientation, and the projection lines of each point on it projecting onto each evaluation picture are constructed according to the camera orientation element information of each evaluation picture (determined in step S2). The projection lines corresponding to each evaluation picture are reprojected onto each evaluation picture. After reprojection onto the evaluation picture, the area covered on the evaluation picture is as follows: Figure 11 The edge of this region, shown in red, constitutes the reprojection window.

[0153] Obviously, during the first reprojection operation, the intersection of the first projection element and the reference projection line is the projection point on the rough three-dimensional model, and the position and orientation of the first projection element are the position and orientation of the element intercepted in step S34.

[0154] refer to Figure 11The reprojection window represents the expected position, size, and shape (i.e., the position of the reprojection window, the number of pixels covered, and the shape formed) of the target ground control point on the corresponding aerial image after it is photographed by the cameras at various positions when the target ground control point is in a certain alternative posture. Obviously, only the first projection surface element that is consistent with the actual posture of the target ground control point, the texture features of the image within the reprojection window obtained by reprojecting it to each evaluation image, have the greatest correlation with the target ground control point.

[0155] Since the rough 3D model has already obtained the approximate pose of the target ground control point, using it to determine the initial pose of the first projection bin can significantly reduce the computational complexity of the search. Only a slight translation or rotation around the initial pose is required to accurately estimate the true pose of the target ground control point. At the same time, for any alternative pose of the first projection bin, its position determines the position of the reprojection window on the evaluation image, and its orientation determines the shape and size of the reprojection window on the evaluation image. It can be seen that the texture features of the image within the reprojection window are simultaneously restricted and constrained by factors such as position, pose, shape, and area. Any reprojection window generated according to any pose that deviates from the actual pose by more than a certain degree will cause the image it covers to deviate significantly from the true imaging result of the target ground control point on the evaluation image in one or more of the following: position, shape, and area, and further cause a significant decrease in the correlation of the image texture within the reprojection window. Therefore, the above method can greatly improve the accuracy of estimating the true pose of the target ground control point.

[0156] When evaluating the correlation of images in each reprojection window in step S352, an image correlation evaluation index known to those skilled in the art can be selected to calculate the correlation of images in the reprojection windows on any two evaluation pictures. For example, an optional correlation evaluation index is the normalized cross-correlation coefficient (NCC), or other indicators for image correlation can also be selected.

[0157] In some preferred embodiments, the reference picture is also used as the evaluation picture, that is, the first projection bin is also reprojected to the reference picture, and the image in its reprojection window is used to participate in the correlation evaluation.

[0158] After completing the correlation calculation of the images in each reprojection window on each evaluation picture, it is judged in step S353 whether the calculation results meet the correlation criterion, and step S354 or step S355 is selected to be executed based on the judgment result.

[0159] In some optional embodiments, the correlation of the images in each reprojection window may be statistically analyzed first to obtain any one or more of the statistical characteristics such as the mean, median, maximum, or minimum, and then determine whether the statistical characteristics are greater than or equal to a preset correlation threshold. If the determination result is negative, the first projection plane is translated or rotated in step S354 (see Figure 11 , the translation or rotation of the first projection surface element is always performed with the intersection of the first projection surface element and the reference projection line as the movement base point, that is, the projection result of the annotation point of the target ground control point to its spatial point is always "bound" to the reference projection line, and its translation and rotation are always performed along the reference projection line), and then return to step S351 to re-project and recalculate the correlation. When the judgment result is yes, execute step S355, and use the position and orientation (normal direction) of the first projection surface element at this time as its optimal position and optimal orientation. The position and orientation will be regarded as the most likely estimated results of the actual position and orientation of the target control point in the area to be modeled, and provide effective constraints when subsequently searching and associating the target control point on the image to be marked.

[0160] Figure 12 The figure shows how the NCC coefficient of the image in each reprojection window changes with the size of the first projection element in a specific embodiment. In the figure, the first projection elements are all square elements, and the horizontal axis represents the number of pixels contained in each side of the element.

[0161] like Figure 12 As shown in the figure, in the area where the first projection bin size is small, the NCC value generally shows a gradual upward trend as the bin size increases, indicating that the matching accuracy improves as the bin size increases. However, as the size further increases, the NCC value no longer increases significantly and even fluctuates (such as the fluctuation between 11×11 pixels and 17×17 pixels).

[0162] Analyzing the texture features of each evaluation image reveals that this correlation metric no longer significantly improves, and may even decrease, after the bin size increases to a certain extent. This indicates that as the number of pixels used for correlation calculation increases, additional image noise may be introduced, resulting in a decrease in the contribution of the increased bin size to pose optimization, or even a negative impact. Therefore, the size of the first projection bin cannot be too small, otherwise there will be insufficient information to evaluate the cross-correlation during the subsequent search process. Preferably, the lower limit of the size of the first projection bin is 3×3 pixels; at the same time, the size of the first projection bin does not need to be too large, as further increasing the size will no longer improve the cross-correlation. In some preferred embodiments, the upper limit of the size of the first projection bin can be determined by calculating the maximum value of the correlation statistics of the images within each reprojection window.

[0163] After determining the actual position of the target ground control point, its actual position can be used as a geometric constraint. Obviously, using this constraint, the search for the target ground control point on the image to be marked can be effectively limited to a specific area, thereby significantly improving the efficiency of automatic search. More importantly, using this geometric constraint, the image of the ground control point automatically searched and marked on the image to be marked must be the image of the target ground control point. That is to say, while automatically searching, the association of the same ground feature point on different aerial images is also completed, thereby fundamentally avoiding the problem of identifying the images of different ground control points on two aerial images as the same ground control point.

[0164] like Figure 7 As shown, searching for target ground control points on the image to be marked is achieved through steps S36 and S37.

[0165] Step S36 is used to generate a second projection surface element: the image edge corresponding to the target ground control point can be generated on the reference image through manual operation or manual plus automatic algorithm interaction, thereby obtaining a reference window; the image in the reference window is projected along the reference projection line onto the first projection surface element with the optimal position and optimal orientation, thereby obtaining a second projection surface element.

[0166] Step S37 is used to automatically search for target ground control points on each image to be marked. As analyzed above, since this step searches for the specified target ground control points on each image to be marked, it also establishes the association between the same ground control points on different aerial images.

[0167] In some specific embodiments, such as Figure 13 As shown, step S37 includes cyclically executing the following steps until all the aerial pictures that need to be marked are processed:

[0168] S371, selecting an aerial image of the area to be modeled as an image to be labeled;

[0169] S372, reprojecting the second projection bin onto the image to be labeled to obtain a desired image;

[0170] S373, extracting the edge of the desired image to construct a matching window;

[0171] S374 , searching for the best position of the matching window on the image to be marked, thereby marking the target ground control point on the image to be marked.

[0172] Specifically, in step S372, the second projection bin can be reprojected onto the image to be labeled along the projection line corresponding to the image to be labeled to obtain its expected image on the image to be labeled. Obviously, the position and orientation of the second projection bin are consistent with the first projection bin having the optimal position and orientation. The difference between the two is that the second projection bin adds the expected shape and texture feature information of the corresponding image of the target ground control point in the reference image when projected onto the target ground control point's true position. Therefore, when the multiple information with this specific position, orientation, shape, and texture features is further reprojected onto the image to be labeled, the resulting expected image can roughly define the expected imaging position of the target ground control point in the image to be labeled, as well as the expected shape, size, and texture of the image. Utilizing this strongly constrained information, the search for the target ground control point eliminates the possibility of identifying different ground control points in two aerial images as the same. Furthermore, the matching window only needs to be translated and slid within a small area around the expected image, eliminating the need for the large-scale search and rotation operations required in conventional image matching retrieval, which undoubtedly greatly improves the search speed.

[0173] Preferably, the search for the optimal position of the matching window in step S374 can be implemented using a least squares image matching algorithm, such as Figure 14 As shown in the figure, the matching window can be slid translationally on the image to be marked according to the set sliding step size, and the image within the matching window is extracted, and the minimum error square calculation is performed on its texture features and the texture features of the expected image. The above process is repeated until the position with the minimum error is found. The image surrounded by the matching window at this position represents the position where the target ground control point is imaged, thereby completing the automatic marking of the target ground control point.

[0174] like Figure 14 As shown, when the desired image is completely inside the picture to be marked, the matching window can be constructed only by the edge of the extracted desired image; Figure 15 The illustrated embodiment further demonstrates the search for the optimal position of a matching window on an image to be marked that only partially displays the target ground control point. When the target ground control point is located at the edge of an aerial image, only a portion of it may be imaged on the aerial image. For conventional automatic labeling methods based solely on texture information, the texture features used for search and comparison are obviously missing, causing difficulty in recognition. However, for the method of the present application, since the desired image is imaged at the edge of the aerial image according to its expected position, and only a portion of it is located on the aerial image, the edge of the desired image and the edge of the aerial image itself can be used to jointly determine the matching window, and its optimal position can be searched on the image to be marked. That is, through the method of the present application, it is possible to complete the identification and marking of the target ground control point using only partial information on the image to be marked.

[0175] Back to Figure 7 After completing steps S31 to S37 once, the automatic identification and association of a target ground control point is complete. You can then select another ground control point as the target ground control point and repeat steps S31 to S37. By repeatedly performing the above operations, you can complete the automatic labeling and association of all or a specific number of ground control points.

[0176] <High-precision modeling of 3D terrain>

[0177] Back to Figure 2 After completing the operation of the above step S3, step S4 can be performed to generate a sparse point cloud set again using the aerial photos that have been marked and associated with the ground control points (i.e., a constrained bundle adjustment operation is performed under the constraints provided by the ground control points to obtain a second sparse point cloud).

[0178] In some specific embodiments, step S4 includes the following steps:

[0179] S41, identifying multiple feature points in the aerial image of each marked ground control point and performing feature point matching;

[0180] S42, obtaining a second estimated value of a camera orientation element of each aerial image and a second estimated value of a three-dimensional coordinate of each feature point by a bundle adjustment with constraints, wherein the constraints include position constraints of the marked ground control points;

[0181] S43: Construct a second sparse point cloud set based on the second estimated value of the three-dimensional coordinate of each feature point.

[0182] The main difference between step S4 and step S2 is that, the bundle adjustment operation at this moment uses the positional information of the ground control points that have been marked and associated as a constraint condition, specifically, because the true three-dimensional coordinates of each ground control point have been obtained by accurate measurement, its position of imaging on each aerial picture has also been marked and associated by step S3, therefore, the positional information of the ground control point can be used as a constraint condition, in the iterative estimation process of bundle adjustment, the estimated result of the three-dimensional coordinate of each feature point is constantly optimized, thereby finally obtaining the second estimated value of the three-dimensional coordinate of each feature point and generating the second sparse point cloud.

[0183] After obtaining the second sparse point cloud, a dense point cloud can be generated in step S5 to obtain a fine three-dimensional terrain model of the area to be modeled. The specific implementations of the above steps S4 and S5 are well known to those skilled in the art and will not be repeated here.

[0184] <Specific Example>

[0185] This embodiment adopts Figure 3A and Figure 3B The aerial image collection shown is used Figure 2 The method shown is used to model the buildings with high precision.

[0186] Figure 16 、 Figure 17 Schematic diagrams of the three-dimensional sparse point cloud and coarse mesh model of the building generated during the implementation process are respectively shown. In the process of automatically marking and associating each ground control point, the size of the first projection element is 11×11 pixels. Table 1 below lists the statistics of the marking results of 10 ground control points.

[0187] Table 1 Statistics of automatic marking results of ground control points

[0188]

[0189] Table 1 shows that when the reference window size is 15×15, a total of 482 points are successfully labeled, with a success rate of approximately 67.3%. The root mean square error between the coordinates of these 482 labeled points and their true coordinates is 0.564 pixels. When the window size is 25×25, a total of 683 points are successfully labeled, with a success rate of approximately 95.5%. The root mean square error between the coordinates of these 683 labeled points and their true coordinates is 0.549 pixels. Statistical results show that the above automatic labeling and association steps can automatically label more than 700 ground control points in less than one minute. As the reference window size increases to 25×25 pixels, the labeling success rate exceeds 95%. At the same time, the overall accuracy of the labeling of each ground control point reaches the sub-pixel level, demonstrating a comprehensive improvement in labeling speed and accuracy compared to manual labeling.

[0190] Figure 18 、 Figure 19 The dense point cloud and fine three-dimensional terrain model of the building are shown respectively. Figure 20 To obtain the final 3D model after mapping the fine model, Figures 18 to 20 It can be seen that the three-dimensional terrain modeling method provided in this application can efficiently and accurately complete three-dimensional terrain modeling based on aerial images.

[0191] Some embodiments of the present application also provide a three-dimensional terrain modeling system, such as Figure 21 As shown, the system includes a processor and a readable storage medium, wherein the readable storage medium stores an executable program, and when the executable program is executed by the processor, the aforementioned three-dimensional terrain modeling method can be implemented.

[0192] The above is a detailed introduction to the specific implementation methods of the present application. For those skilled in the art, several improvements and modifications can be made to the present application without departing from the principles of the present application. These improvements and modifications also fall within the scope of protection of the claims of the present application.

Claims

1. A three-dimensional terrain modeling method, characterized in that: The following steps are involved: S1, obtain aerial images of several unmarked ground control points in the area to be modeled; S2, generating a first sparse point cloud set and a rough three-dimensional terrain model of the area to be modeled based on the aerial images of each unmarked ground control point; S3, automatically marking and associating ground control points on each aerial image based on the rough three-dimensional terrain model; S4, generating a second sparse point cloud set of the area to be modeled based on the aerial images of each marked ground control point; S5, generating a dense point cloud set and a fine three-dimensional terrain model of the area to be modeled based on the second sparse point cloud set; Step S3 includes the following steps: S31, selecting a ground control point as a target ground control point; S32, selecting a reference image and several evaluation images corresponding to the target ground control point; S33, determining at least one annotation point from the image corresponding to the target ground control point on the reference image; S34, projecting the marked point onto the surface of the rough three-dimensional terrain model and intercepting a first projection surface element; S35, searching for the best position and best orientation of the first projection plane element; S36, projecting an image within a reference window on the reference picture onto a first projection plane element having an optimal position and an optimal orientation to obtain a second projection plane element; S37 , based on the expected image obtained by reprojecting the second projection bin onto each of the to-be-marked pictures, searching and determining the imaging position corresponding to the target ground control point on each of the to-be-marked pictures.

2. The three-dimensional terrain modeling method according to claim 1, characterized in that: Step S2 includes the following steps: S21, identifying multiple feature points in each aerial image; S22, matching feature points in different aerial images; S23, obtaining a first estimated value of a camera orientation element of each aerial image and a first estimated value of a three-dimensional coordinate of each feature point based on an unconstrained bundle adjustment; S24, constructing a first sparse point cloud set based on the first estimated value of the three-dimensional coordinates of each feature point; S25: Establish a rough three-dimensional terrain model of the area to be modeled based on the first sparse point cloud set.

3. The three-dimensional terrain modeling method according to claim 1, characterized in that: The similarity between the image of the target ground control point on the reference image and its orthographic image is greater than a preset first similarity threshold, or the deviation is less than a preset first deviation threshold; The similarity between the image of the target ground control point on the evaluation picture and its orthographic image or its image on the reference picture is greater than a preset second similarity threshold, or the deviation is less than a preset second deviation threshold, wherein the first similarity threshold is greater than or equal to the second similarity threshold, and the first deviation threshold is less than or equal to the second deviation threshold.

4. The three-dimensional terrain modeling method according to claim 1, characterized in that: Step S34 includes the following steps: The marked point is projected onto the surface of the rough three-dimensional terrain model along a reference projection line corresponding to the reference image to obtain its projection point on the surface of the rough three-dimensional terrain, wherein the reference projection line of the marked point corresponding to the reference image is determined based on the camera orientation element of the reference image.

5. The three-dimensional terrain modeling method according to claim 1, characterized in that: Step S35 includes the following steps: S351, reprojecting the first projection bin onto each evaluation image to obtain a reprojection window; S352, determining the cross-correlation characteristics of the images in the reprojection windows of each evaluation picture; S353, determining whether the mutual correlation feature satisfies the correlation criterion; if not, executing step S354 and returning to step S351; if satisfied, executing step S355; S354, translating or rotating the first projection surface element along the reference projection line; S355: Taking the current position and orientation of the first projection plane element as the optimal position and optimal orientation of the first projection plane element.

6. The three-dimensional terrain modeling method according to claim 1, characterized in that: Step S37 includes cyclically executing the following steps until all the aerial images that need to be marked are processed: S371, selecting an aerial image of the area to be modeled as an image to be labeled; S372, reprojecting the second projection bin onto the image to be labeled to obtain a desired image; S373, extracting the edge of the desired image to construct a matching window; S374 , searching for the best position of the matching window on the image to be marked, thereby marking the target ground control point on the image to be marked.

7. The three-dimensional terrain modeling method according to claim 6, characterized in that: In the process of searching for the best position of the matching window on the image to be marked, the matching window is only translated on the image to be marked.

8. The three-dimensional terrain modeling method according to claim 1, characterized in that: Step S4 includes the following steps: S41, identifying multiple feature points in the aerial image of each marked ground control point and performing feature point matching; S42, obtaining a second estimated value of a camera orientation element of each aerial image and a second estimated value of a three-dimensional coordinate of each feature point by a bundle adjustment with constraints, wherein the constraints include position constraints of the marked ground control points; S43: Construct a second sparse point cloud set based on the second estimated value of the three-dimensional coordinate of each feature point.

9. A three-dimensional terrain modeling system, comprising a processor and a readable storage medium, characterized in that: The readable storage medium stores an executable program, and when the executable program is executed by the processor, the three-dimensional terrain modeling method according to claim 1 can be implemented.

Images