Aerial photography data building height calculation method based on 3D Gaussian technology

Abstract

The invention provides an aerial photography data building height calculation method based on a 3D Gaussian technology. The method comprises the steps that an aerial photography picture is acquired; performing data preprocessing; geographic coordinate conversion is carried out; all aerial photos are taken, the positions of the aerial photos in a local coordinate system before conversion are calculated, and a conversion function fbest (x, y, z) is applied to enable the local coordinate system to be aligned with a target coordinate system; inputting the internal parameters, the external parameters and the sparse point cloud of the converted aerial photograph into a 3DGS algorithm, and performing scene reconstruction; acquiring a click position of a user; and calculating the building height of the user click position. According to the aerial photography data building height calculation method based on the 3D Gaussian technology, the 3DGS technology is introduced into urban building three-dimensional modeling, and interactive browsing and automatic height measurement tasks are achieved while the advantages of the 3DGS technology are considered for high-definition image data collected by an unmanned aerial vehicle.

Description

TECHNICAL FIELD

[0001] The application relates to the technical field of aerial building height calculation, and in particular to a method for calculating the height of a building based on 3D Gaussian technology. BACKGROUND

[0002] With the rise of smart city construction, three-dimensional modeling of city-level large scenes has wide application value in many fields, including city roaming, autonomous driving and navigation and tour guidance. Traditional three-dimensional modeling methods, such as oblique photography modeling, express scenes by displaying, and describe scene information through elements such as geometry, texture and material. This method has a long modeling process, produces a large amount of data, and in order to achieve high fidelity, often requires manual participation in a large amount of editing work, and the data production efficiency is low.

[0003] 3D Gaussian is a new breakthrough in the field of multi-view synthesis, which represents the geometry and scene in 3D space as a large number of Gaussian distributions, and then projects and renders pixel colors on these Gaussian distributions through the parallel computing power of GPU to generate realistic 2D images. First, the colmap software is used to calculate the camera intrinsic and extrinsic parameters of each photo to generate a sparse point cloud in the local coordinate system. Then, the 3DGS technology starts from sparse points and represents the scene with a three-dimensional Gaussian distribution, which retains the ideal properties of continuous body radiation field for scene optimization, while avoiding unnecessary calculations in blank space; according to the camera intrinsic and extrinsic parameters, the 3D Gaussian distribution is projected onto a two-dimensional plane, and then the 3DGS uses a fast visibility-aware rendering algorithm to calculate the color of each pixel, which allows real-time rendering, supports anisotropic splatting, and accelerates training. According to the loss calculated from the inference image and the real image, the parameters of the 3DGS are optimized through back propagation; at the same time, the three-dimensional Gaussian function is optimized by interleaving or density control, especially the anisotropic covariance, to achieve accurate representation of the scene.

[0004] However, using 3D Gaussian technology for building height measurement based on aerial data has the following problems: the previous height measurement method has the problems of low precision, low efficiency, high cost and great difficulty, and is difficult to promote; the model established by the 3DGS algorithm has high fidelity and simplifies the modeling steps, but the model established by the 3DGS algorithm is only a visual model similar to the real world, but the model cannot directly recover the real spatial distance. The 3DGS modeling algorithm based on depth map has the problems of unclear boundary under the current view and inaccurate depth acquisition at the boundary, which also directly affects the accuracy of spatial distance calculation.

[0005] Patent application CN119762983A discloses a method, apparatus, device, and storage medium for identifying building height changes based on remote sensing spatiotemporal information. The method includes: obtaining the building's bottom outline using a vector map; obtaining the building's top outline at different times using multi-temporal high-resolution optical remote sensing images; obtaining the projected displacement of the building's top and bottom outlines based on the Euclidean distance between the geometric center points of the building's bottom and top outlines; calculating the relative projected displacement ratio of multiple buildings and the absolute projected displacement ratio of a single building based on the projected displacements of the building's top and bottom outlines at different times; calculating the satellite inclination tangent ratio based on the satellite observation angles at different times; and identifying the building's height change based on the relative projected displacement ratio of multiple buildings, the absolute projected displacement ratio of a single building, and the satellite inclination tangent ratio. A drawback of this method is that for buildings with irregular shapes and complex top structures, it is difficult to accurately determine the geometric center point of the top outline, which may lead to errors in the projected displacement calculation and affect the accuracy of height change identification.

[0006] Patent application CN 119131119A discloses a method, system, storage medium, and computer program product for determining building height. The method includes: determining a top predicted bounding box and multiple bottom candidate bounding boxes of a building to be identified in a remote sensing image; determining the bottom predicted bounding box of the building to be identified from among the multiple bottom candidate bounding boxes based on the geometric features of the top and bottom predicted bounding boxes; and analyzing the pixel offset between the top and bottom predicted bounding boxes to determine the height of the building to be identified. The pixel offset represents the pixel dimension distance from the center point of the top predicted bounding box to the center point of the bottom predicted bounding box, and the pixel offset is obtained from the remote sensing image. The disadvantages of this method are that it only analyzes remote sensing images from a single viewpoint, making it difficult to handle occluded or overlapping buildings; and weather conditions affect the stability of top / bottom detection. Summary of the Invention

[0007] To address the aforementioned technical issues, this invention proposes a method for calculating building height based on aerial photography data using 3D Gaussian technology. This method introduces 3DGS technology into urban building 3D modeling and, while leveraging the advantages of 3DGS technology, enables interactive browsing and automatic height measurement for high-definition image data collected by drones.

[0008] The purpose of this invention is to provide a method for calculating building height based on aerial photography data using 3D Gaussian technology, including acquiring aerial photographs and the following steps:

[0009] Step 1: Perform data preprocessing;

[0010] Step 2: Perform geographic coordinate transformation;

[0011] Step 3: Take all aerial photographs, calculate their positions in the local coordinate system before transformation, and apply the transformation function f_best (x,y,z) to align the local coordinate system with the target coordinate system;

[0012] Step 4: Input the intrinsic parameters, extrinsic parameters, and sparse point cloud data of the converted aerial photographs into the 3DGS algorithm for scene reconstruction;

[0013] Step 5: Obtain the user's click location;

[0014] Step 6: Calculate the building height at the user's clicked location.

[0015] Preferably, step 1 includes downsampling the aerial photograph and estimating the camera's intrinsic and extrinsic parameters.

[0016] In any of the above solutions, step 2 preferably includes the following sub-steps:

[0017] Step 21: Extract the latitude, longitude, and altitude from the EXIF ​​information of the aerial photograph and convert them to the target coordinate system;

[0018] Step 22: Randomly select several aerial photographs, calculate their positions in the local coordinate system before transformation, then extract the corresponding transformed coordinates, and calculate the similarity transformation function f(x,y,z);

[0019] Step 23: Take all aerial photographs, calculate their positions in the local coordinate system before the transformation, and input them into the similarity transformation function f(x,y,z) to obtain the transformed coordinates X', Y', Z';

[0020] Step 24: For the same aerial photograph, calculate the distance between the converted coordinates X', Y', Z' and X, Y, Z. If the distance is less than the threshold, it is recorded as an interior point; if the distance is greater than or equal to the threshold, it is recorded as an exterior point. Count the number of interior points.

[0021] Step 25: When the number of interior points no longer increases, record the transformation function f_best(x,y,z); otherwise, go to step 22.

[0022] In any of the above solutions, the preferred method includes obtaining the transformed coordinates, as shown in the formula:

[0023] X=(N+alt)*cos_lat*cos_lon

[0024] Y=(N+alt)*cos_lat*sin_lon

[0025] Z=(N*(1-e2_)+alt)*sin_lat

[0026] Where X, Y, and Z are the transformed coordinates, alt is the altitude, cos_lat is the cosine of the latitude of the aerial photograph, cos_lon is the cosine of the longitude of the aerial photograph, sin_lat is the sine of the latitude of the aerial photograph, sin_lon is the sine of the longitude of the aerial photograph, e2_ is the square of the first eccentricity of the Earth ellipsoid, and N is the radius of curvature of the prime mover and trochanter. The formula is:

[0027]

[0028] Where α is the Earth's semi-major axis, approximately 6,378,137 meters.

[0029] In any of the above schemes, the preferred method is to continue the transformation to obtain the result coordinates, as shown in the formula:

[0030] R*(point_i-point_0)

[0031] Where point_i are the coordinates Xi, Yi, Zi of each point after transformation, point_0 are the coordinates X0, Y0, Z0 of the point corresponding to the first aerial photograph after transformation, and R is the rotation matrix, with the formula:

[0032]

[0033] Here, sin_lat_point_0 applies the sine function to the latitude of the first aerial photograph, sin_lon_point_0 applies the sine function to the latitude of the first aerial photograph, cos_lat_point_0 applies the cosine function to the latitude of the first aerial photograph, and cos_lon_point_0 applies the cosine function to the latitude of the first aerial photograph.

[0034] In any of the above schemes, step 22 preferably includes estimating the transformation matrix between the two point sets using the Umeyama algorithm, i.e., solving for R and t to optimize the following function:

[0035]

[0036] The transformation matrix is ​​obtained by solving:

[0037]

[0038] Where n is the total number of points, i is the i-th point, and p i Let P be the i-th point in the source point set. i Let be the i-th point in the target point set, c be the scaling factor, R be the rotation matrix, t be the translation matrix, and T be the transformation matrix between the two point sets.

[0039] In any of the above solutions, step 6 preferably includes the following sub-steps:

[0040] Step 61: Using the user's click location as the tracking point, rotate the viewpoint around the tracking point to synthesize multiple depth maps from the rooftop viewpoint;

[0041] Step 62: For each depth map, use a confidence function to filter the depths and discard erroneous depths;

[0042] Step 63: Convert all pixels with the correct depth into point clouds;

[0043] Step 64: Take the maximum z-coordinate in the point cloud, subtract the building's bottom height, and obtain the building height.

[0044] In any of the above solutions, step 61 preferably includes the following sub-steps:

[0045] Step 611: Obtain the synthesized viewpoint information;

[0046] Step 612: Construct the rotation matrix;

[0047] Step 613: Calculate the camera position, i.e., the translation vector;

[0048] Step 614: Construct the complete camera_to_world matrix, which includes the rotation matrix and translation vector.

[0049] In any of the above schemes, it is preferred that the viewpoint information includes a forward vector, a right vector, and a down vector. The forward vector is a unit vector pointing from the camera to the target point. The right vector is obtained by the cross product of the world_up vector and the forward vector. The correct down vector is obtained by the cross product of the forward vector and the right vector.

[0050] In any of the above schemes, it is preferred to arrange the right vector, down vector, and forward vector in columns to obtain the rotation matrix.

[0051] In any of the above schemes, the preferred option is that the formula for the translation vector is:

[0052] T=[camera_x, camera_y, camera_z]

[0053] camera_x = target_x + distance * sin(angle) * cos(pitch_angle)

[0054] camera_y = target_y + distance * sin(pitch_angle)

[0055] camera_z = target_z + distance * cos(angle) * cos(pitch_angle)

[0056] Where camera_x, camera_y, and camera_z are the x, y, and z coordinates of the camera position, target_x, target_y, and target_z are the x, y, and z coordinates of the tracking point, distance is the distance from the camera position to the tracking point set by the user, angle is the horizontal angle, and pitch_angle is the downward angle.

[0057] In any of the above schemes, the preferred option is that the confidence level is calculated using the following formula:

[0058]

[0059] Where N_confidence is the ordered sequence of all Gaussian points sorted by depth, and α m_confidence Let α be the opacity at the m-th Gaussian point. n_confidence Let be the opacity of the nth Gaussian point, m_confidence be the mth Gaussian point in the ordered sequence, and n_confidence be the nth Gaussian point in the ordered sequence.

[0060] In any of the above schemes, the preferred option is that the confidence function formula is:

[0061]

[0062] Where threshold=0.999.

[0063] In any of the above schemes, it is preferred that when the confidence function of a certain pixel is calculated to be 1, the depth synthesis is accurate and the depth corresponding to that pixel is retained; when the confidence function of a certain pixel is calculated to be 0, the depth synthesis is inaccurate and the depth corresponding to that pixel is discarded.

[0064] In any of the above solutions, step 63 preferably includes the following sub-steps:

[0065] Step 631: For each pixel (u, v) in the depth map, its corresponding depth value is d. Calculate the 3D coordinates (x, y) of this point in the camera coordinate system. c ,y c ,zc );

[0066] Step 632: Transfer the point (x) in the camera coordinate system to... c ,y c ,z c Convert ) to a point in world coordinates (x) w ,y w ,z w ).

[0067] Preferably, in any of the above schemes, the three-dimensional coordinates (x) c ,y c ,z c The formula for calculating ) is

[0068]

[0069] Where K is the camera intrinsic parameter and d is the depth value corresponding to pixel (u, v).

[0070] In any of the above schemes, the preferred option is that point (x) w ,y w ,z w The conversion formula is:

[0071]

[0072] Where c2w is the transformation matrix from the camera coordinate system to the world coordinate system.

[0073] This invention proposes a method for calculating building height based on aerial photography data using 3D Gaussian technology. The multi-view height finding (highest point) process traverses every view near the rooftop, eliminating the need for manual searching, and is highly universal, enabling automatic finding of the highest point. Attached Figure Description

[0074] Figure 1 This is a flowchart of a preferred embodiment of the method for calculating building height based on aerial photography data according to the present invention using 3D Gaussian technology.

[0075] Figure 2 This is a flowchart of another preferred embodiment of the method for calculating building height based on aerial photography data according to the present invention using 3D Gaussian technology.

[0076] Figure 3 This is a schematic diagram of an embodiment of the user-selected building according to the method for calculating building height based on aerial photography data using 3D Gaussian technology of the present invention.

[0077] Figure 4 This is a schematic diagram of an embodiment of the viewpoint rotation method for calculating building height based on aerial photography data according to the present invention using 3D Gaussian technology. Detailed Implementation

[0078] The present invention will be further described below with reference to the accompanying drawings and specific embodiments.

[0079] Example 1

[0080] like Figure 1 As shown, step 1000 is executed to obtain aerial photographs.

[0081] Step 1100 is executed to perform data preprocessing, including downsampling the aerial photographs and estimating camera intrinsic and extrinsic parameters.

[0082] Execute step 1200 to perform geographic coordinate transformation, which includes the following sub-steps:

[0083] Execute step 1210 to extract the latitude, longitude, and altitude from the EXIF ​​information of the aerial photograph and convert them to the target coordinate system. The conversion method includes obtaining the converted coordinates using the following formula:

[0084] X=(N+alt)*cos_lat*cos_lon

[0085] Y=(N+alt)*cos_lat*sin_lon

[0086] Z=(N*(1-e2_)+alt)*sin_lat

[0087] Where X, Y, and Z are the transformed coordinates, alt is the altitude, cos_lat is the cosine of the latitude of the aerial photograph, cos_lon is the cosine of the longitude of the aerial photograph, sin_lat is the sine of the latitude of the aerial photograph, sin_lon is the sine of the longitude of the aerial photograph, e2_ is the square of the first eccentricity of the Earth ellipsoid, and N is the radius of curvature of the prime mover and trochanter. The formula is:

[0088]

[0089] Where α is the Earth's semi-major axis, approximately 6,378,137 meters.

[0090] The transformation method includes further transformation to obtain the resulting coordinates, as shown in the formula:

[0091] R*(point_i-point_0)

[0092] Where point_i are the coordinates Xi, Yi, Zi of each point after transformation, point_0 are the coordinates X0, Y0, Z0 of the point corresponding to the first aerial photograph after transformation, and R is the rotation matrix, with the formula:

[0093]

[0094] Here, sin_lat_point_0 applies the sine function to the latitude of the first aerial photograph, sin_lon_point_0 applies the sine function to the latitude of the first aerial photograph, cos_lat_point_0 applies the cosine function to the latitude of the first aerial photograph, and cos_lon_point_0 applies the cosine function to the latitude of the first aerial photograph.

[0095] Execute step 1220, randomly select several aerial photographs, calculate their positions in the local coordinate system before transformation, then extract the corresponding transformed coordinates, and calculate the similarity transformation function f(x,y,z).

[0096] The transformation matrix between two sets is estimated using the Umeyama algorithm, i.e., solving for R and t to optimize the following function:

[0097]

[0098] The transformation matrix is ​​obtained by solving:

[0099]

[0100] Where n is the total number of points, i is the i-th point, and p i Let P be the i-th point in the source point set. i Let be the i-th point in the target point set, c be the scaling factor, R be the rotation matrix, t be the translation matrix, and T be the transformation matrix between the two point sets.

[0101] Execute step 1230, take all aerial photographs, calculate their positions in the local coordinate system before transformation, and input them into the similarity transformation function f(x,y,z) to obtain the transformed coordinates X', Y', Z'.

[0102] Execute step 1240: For the same aerial photograph, calculate the distance between the converted coordinates X', Y', Z' and X, Y, Z. If the distance is less than the threshold, it is recorded as an interior point; if the distance is greater than or equal to the threshold, it is recorded as an exterior point. Count the number of interior points.

[0103] Execute step 1250. When the number of interior points no longer increases, record the transformation function f_best(x,y,z); otherwise, proceed to step 1220.

[0104] Execute step 1300, take all aerial photographs, calculate their positions in the local coordinate system before transformation, and apply the transformation function f_best (x,y,z) to align the local coordinate system with the target coordinate system.

[0105] Step 1400 involves inputting the intrinsic and extrinsic parameters of the converted aerial photograph, along with the sparse point cloud data, into the 3DGS algorithm for scene reconstruction.

[0106] Perform step 1500 to obtain the user's click location.

[0107] Execute step 1600 to calculate the building height at the user-clicked location, including the following sub-steps:

[0108] Step 1610 involves using the user's click location as a tracking point, rotating the viewpoint around the tracking point, and synthesizing multiple depth maps from the rooftop perspective. This includes the following sub-steps:

[0109] Execute step 1611 to obtain the synthesized viewpoint information, which includes a forward vector, a right vector, and a down vector. The forward vector is a unit vector pointing from the camera to the target point. The right vector is obtained by the cross product of the world_up vector and the forward vector. The correct down vector is obtained by the cross product of the forward vector and the right vector.

[0110] Perform step 1612 to construct the rotation matrix by arranging the right vector, down vector, and forward vector column-wise.

[0111] Execute step 1613 to calculate the camera position, i.e., the translation vector. The formula for the translation vector is:

[0112] T=[camera_x, camera_y, camera_z]

[0113] camera_x = target_x + distance * sin(angle) * cos(pitch_angle)

[0114] camera_y = target_y + distance * sin(pitch_angle)

[0115] camera_z = target_z + distance * cos(angle) * cos(pitch_angle)

[0116] Where camera_x, camera_y, and camera_z are the x, y, and z coordinates of the camera position, target_x, target_y, and target_z are the x, y, and z coordinates of the tracking point, distance is the distance from the camera position to the tracking point set by the user, angle is the horizontal angle, and its value is equal to 360 / the number of intervals specified by the user; pitch_angle is the pitch angle, which is 45 degrees.

[0117] Execute step 1614 to construct the complete camera_to_world matrix, which includes the rotation matrix and translation vector.

[0118] Execute step 1620: For each depth map, use a confidence function to filter the depths, discarding erroneous depths. The confidence calculation formula is:

[0119]

[0120] Where N_confidence is the ordered sequence of all Gaussian points sorted by depth, α m_confidence Let α be the opacity at the m-th Gaussian point. n_confidence Let be the opacity of the nth Gaussian point, m_confidence be the mth Gaussian point in the ordered sequence, and n_confidence be the nth Gaussian point in the ordered sequence.

[0121] The confidence function formula is:

[0122]

[0123] Where threshold=0.999.

[0124] When the confidence function of a pixel is calculated to be 1, the depth synthesis is accurate, and the depth corresponding to that pixel is retained; when the confidence function of a pixel is calculated to be 0, the depth synthesis is inaccurate, and the depth corresponding to that pixel is discarded.

[0125] Execute step 1630 to convert all pixels with the correct depth into a point cloud, including the following sub-steps:

[0126] Execute step 1631: For each pixel (u, v) in the depth map, its corresponding depth value is d. Calculate the 3D coordinates (x, y) of that point in the camera coordinate system. c ,y c ,z c The three-dimensional coordinates (x) c ,y c ,z c The formula for calculating ) is

[0127]

[0128] Where K is the camera intrinsic parameter and d is the depth value corresponding to pixel (u, v).

[0129] Execute step 1632 to set the point (x) in the camera coordinate system.c ,y c ,z c Convert ) to a point in world coordinates (x) w ,y w ,z w ), point (x) w ,y w ,z w The conversion formula is:

[0130]

[0131] Where c2w is the transformation matrix from the camera coordinate system to the world coordinate system.

[0132] Execute step 1640, take the maximum z-coordinate in the point cloud, subtract the building's bottom height, and obtain the building height.

[0133] Example 2

[0134] This invention discloses a method for calculating building height from aerial photography data based on 3D Gaussian technology. It utilizes the 3DGS algorithm for scene reconstruction, requiring only drone aerial images. A multi-view height-finding (highest point) process is then developed, traversing the vicinity of the building rooftop and using information from multiple perspectives to find the highest point, achieving high accuracy, low cost, and minimal complexity. The 3DGS model directly reconstructs the true spatial distance. A coordinate system transformation function is used to convert the aerial image shooting location (local coordinate system) calculated by the 3DGS algorithm to a globally unified coordinate system (unit: meters). An outlier removal algorithm is used to improve the robustness of the transformation function, achieving accurate distance recovery. During depth map synthesis, outlier depth removal is performed: based on the cumulative opacity of pixels, the confidence level of each pixel's depth value is obtained, and a confidence filtering function is designed to filter depth values, removing those with large errors.

[0135] like Figure 2 As shown, a method for calculating building height based on aerial photography data using 3D Gaussian technology includes the following steps:

[0136] Step 1: Data Preprocessing

[0137] Data preprocessing includes downsampling aerial photographs (based on actual hardware and rendering requirements), estimating camera intrinsic and extrinsic parameters, etc.

[0138] Step 2: Geographic coordinate system transformation

[0139] Step 2.1: Extract the latitude, longitude, and altitude from the EXIF ​​information of the aerial photograph, and directly convert them to the target coordinate system using the following formula:

[0140] X=(N+alt)*cos_lat*cos_lon

[0141] Y=(N+alt)*cos_lat*sin_lon

[0142] Z=(N*(1-e2_)+alt)*sin_lat

[0143] Where X, Y, and Z are the transformed coordinates; the suffixes lat and lon represent latitude; the prefixes sin and cos represent the sine and cosine functions; alt represents the altitude in meters; e²_ is the square of the first eccentricity of the Earth ellipsoid, approximately 0.00669438002290; and N is the radius of curvature of the prime mover, calculated using the formula:

[0144]

[0145] Where α is the Earth's semi-major axis, approximately 6,378,137 meters.

[0146] Then, continue the transformation to obtain the resulting coordinates:

[0147] R*(point_i-point_0)

[0148] Where, point_i are the coordinates Xi, Yi, Zi of each point after the above transformation, point_0 are the coordinates X0, Y0, Z0 of the point corresponding to the first aerial photograph after the above transformation, and R is the rotation matrix, with the formula as follows:

[0149] [-sin_lon, cos_lon, 0.,-sin_lat * cos_lon, -sin_lat * sin_lon, cos_lat,cos_lat * cos_lon, cos_lat * sin_lon, sin_lat;]

[0150] In this context, the suffix lat represents the latitude of the first aerial photograph, the suffix lon represents the latitude of the first aerial photograph, the prefix sin represents the sine function, and the prefix cos represents the cosine function.

[0151] Step 2.2: Randomly select several aerial photographs, calculate their positions in the local coordinate system before transformation, then extract the corresponding transformed coordinates (results of Step 2.1), and calculate the similarity transformation function f(x,y,z).

[0152] Step 2.3: Take all aerial photographs, calculate their positions in the local coordinate system before transformation, and input them into the similarity transformation function f(x,y,z) to obtain the transformed coordinates X',Y',Z'.

[0153] Step 2.4: For the same aerial photograph, calculate the distance between the converted coordinates X', Y', Z' and X, Y, Z. If the distance is less than the threshold, it is recorded as an "interior point"; if the distance is greater than or equal to the threshold, it is recorded as an "outer point". Count the number of "interior points".

[0154] Step 2.5: If the number of "interior points" no longer increases, record the transformation function f_best(x,y,z) and go to step 3; otherwise, go to step 2.2.

[0155] Step 3: Take all aerial photographs, calculate their positions in the local coordinate system before transformation, and apply the transformation function f_best (x,y,z) to align the local coordinate system with the target coordinate system.

[0156] Step 4: Input the intrinsic parameters, extrinsic parameters, and sparse point cloud data of the converted aerial photographs into the 3DGS algorithm for scene reconstruction.

[0157] Step 5: The system reads the user's click location, such as... Figure 3 As shown.

[0158] Step 6: Altitude Measurement Process

[0159] Step 6.1: As Figure 4 As shown, the location clicked by the user is used as the tracking point, and the viewpoint rotates around the tracking point to synthesize multiple depth maps from the rooftop perspective.

[0160] The synthesized viewpoint information requires the following parameters:

[0161] forward vector: the unit vector pointing from the camera to the target point.

[0162] The right vector is obtained by the cross product of the world_up vector and the forward vector.

[0163] The down vector is obtained by the cross product of the forward and right vectors.

[0164] Construct the rotation matrix: Arrange the right, down, and forward vectors column-wise.

[0165] Calculate the camera position, i.e., the translation vector.

[0166] Construct the complete camera_to_world matrix: including rotation matrix + translation vector

[0167] Step 6.2: [Original Technology] For each depth map, use a confidence function to filter the depths and discard erroneous depths.

[0168] The confidence level is calculated using the following formula:

[0169]

[0170] Where N_confidence is the ordered sequence of all Gaussian points sorted by depth, and α m_confidence Let α be the opacity at the m-th Gaussian point. n_confidence Let be the opacity of the nth Gaussian point, where m_confidence is the opacity of the mth Gaussian point in the ordered sequence, and n_confidence is the opacity of the nth Gaussian point in the ordered sequence. When the number of Gaussian points is sufficiently large, mathematically speaking, the upper limit of the sum is 1. Therefore, if the cumulative opacity reaches 0.999, the confidence level of the depth value reaches 0.999.

[0171] The confidence function formula is:

[0172]

[0173] The threshold is typically set to 0.999. When the confidence level of a pixel is input into the formula, if the result is 1, the depth synthesis is accurate, and the depth corresponding to that pixel is retained; if the result is 0, the depth synthesis is inaccurate, and the depth corresponding to that pixel is discarded.

[0174] Step 6.3: Convert all pixels with the correct depth into point clouds.

[0175] For each pixel (u, v) in the depth map, its corresponding depth value is d. Then, the 3D coordinates (x, y) of that point in the camera coordinate system are... c ,y c ,z c It can be calculated using the following formula:

[0176]

[0177] Where K is the camera intrinsic parameter, (x c ,y c ,z c () represents the three-dimensional coordinates in the camera coordinate system.

[0178] Then, the point (x) in the camera coordinates c ,y c ,z c Convert ) to a point in world coordinates (x) w ,y w ,z w ).

[0179]

[0180] Where c2w is the transformation matrix from the camera coordinate system to the world coordinate system.

[0181] Step 6.4: Take the maximum z-coordinate in the point cloud, subtract the building's bottom height, and obtain the building height.

[0182] Step 6.5: The user chooses whether to continue. If they choose to continue, proceed to step 6.1; otherwise, the process ends.

[0183] To better understand this invention, specific embodiments have been described in detail above, but these are not intended to limit the invention. Any simple modifications made to the above embodiments based on the technical essence of this invention still fall within the scope of this invention. Each embodiment in this specification focuses on its differences from other embodiments; similar or identical parts between embodiments can be referred to mutually. For system embodiments, since they basically correspond to method embodiments, the description is relatively simple; relevant parts can be referred to the descriptions of the method embodiments.

Claims

1. A method for calculating building height based on aerial photographic data using 3D Gaussian technology, comprising acquiring aerial photographs, characterized in that, It also includes the following steps: Step 1: Perform data preprocessing; Step 2: Perform geographic coordinate transformation; Step 3: Take all aerial photographs, calculate their positions in the local coordinate system before transformation, and apply the transformation function f_best(x,y,z) to align the local coordinate system with the target coordinate system; Step 4: Input the intrinsic parameters, extrinsic parameters, and sparse point cloud data of the converted aerial photographs into the 3DGS algorithm for scene reconstruction; Step 5: Obtain the user's click location; Step 6: Calculate the building height at the user-clicked location, including the following sub-steps: Step 61: Using the user's click location as the tracking point, rotate the viewpoint around the tracking point to synthesize multiple depth maps from the rooftop viewpoint; Step 62: For each depth map, use a confidence function to filter the depths and discard erroneous depths; Step 63: Convert all pixels with the correct depth into point clouds; Step 64: Take the maximum z-coordinate in the point cloud, subtract the building's bottom height, and obtain the building height.

2. The method for calculating building height based on aerial photography data using 3D Gaussian technology as described in claim 1, characterized in that, Step 1 includes downsampling the aerial photographs and estimating the camera's intrinsic and extrinsic parameters.

3. The method for calculating building height based on aerial photography data using 3D Gaussian technology as described in claim 2, characterized in that, Step 2 includes the following sub-steps: Step 21: Extract the latitude, longitude, and altitude from the EXIF ​​information of the aerial photograph and convert them to the target coordinate system; Step 22: Randomly select several aerial photographs, calculate their positions in the local coordinate system before transformation, then extract the corresponding transformed coordinates, and calculate the similarity transformation function f(x,y,z); Step 23: Take all aerial photographs, calculate their positions in the local coordinate system before the transformation, and input them into the similarity transformation function f(x,y,z) to obtain the transformed coordinates X', Y', Z'; Step 24: For the same aerial photograph, calculate the distance between the converted coordinates X', Y', Z' and X, Y, Z. If the distance is less than the threshold, it is recorded as an interior point; if the distance is greater than or equal to the threshold, it is recorded as an exterior point. Count the number of interior points. Step 25: When the number of interior points no longer increases, record the transformation function f_best(x,y,z); otherwise, go to step 22.

4. The method for calculating building height based on aerial photography data using 3D Gaussian technology as described in claim 3, characterized in that, The transformation method includes obtaining the transformed coordinates, using the following formula: X=(N+alt)*cos_lat*cos_lon Y=(N+alt)*cos_lat*sin_lon Z=(N*(1-e2_)+alt)*sin_lat Where X, Y, and Z are the transformed coordinates, alt is the altitude, cos_lat is the cosine of the latitude of the aerial photograph, cos_lon is the cosine of the longitude of the aerial photograph, sin_lat is the sine of the latitude of the aerial photograph, sin_lon is the sine of the longitude of the aerial photograph, e2_ is the square of the first eccentricity of the Earth ellipsoid, and N is the radius of curvature of the prime mover and trochanter. The formula is: ， Where α is the Earth's semi-major axis.

5. The method for calculating building height based on aerial photography data using 3D Gaussian technology as described in claim 4, characterized in that, The transformation method includes further transformation to obtain the resulting coordinates, as shown in the formula: R*(point_i-point_0) Where point_i are the coordinates Xi, Yi, Zi of each point after transformation, point_0 are the coordinates X0, Y0, Z0 of the point corresponding to the first aerial photograph after transformation, and R is the rotation matrix, with the formula: ， Here, sin_lat_point_0 applies the sine function to the latitude of the first aerial photograph, sin_lon_point_0 applies the sine function to the latitude of the first aerial photograph, cos_lat_point_0 applies the cosine function to the latitude of the first aerial photograph, and cos_lon_point_0 applies the cosine function to the latitude of the first aerial photograph.

6. The method for calculating building height based on aerial photography data using 3D Gaussian technology as described in claim 5, characterized in that, Step 61 includes the following sub-steps: Step 611: Obtain the synthesized viewpoint information; Step 612: Construct the rotation matrix; Step 613: Calculate the camera position, i.e., the translation vector; Step 614: Construct the complete camera_to_world matrix, which includes the rotation matrix and translation vector.

7. The method for calculating building height based on aerial photography data using 3D Gaussian technology as described in claim 6, characterized in that, The viewpoint information includes a forward vector, a right vector, and a down vector. The forward vector is a unit vector pointing from the camera to the target point. The right vector is obtained by the cross product of the world_up vector and the forward vector. The correct down vector is obtained by the cross product of the forward vector and the right vector.

8. The method for calculating building height based on aerial photography data using 3D Gaussian technology as described in claim 7, characterized in that, The rotation matrix is ​​obtained by arranging the right vector, down vector, and forward vector column by column.

9. The method for calculating building height based on aerial photography data using 3D Gaussian technology as described in claim 8, characterized in that, The confidence level is calculated using the following formula: ， Where N_confidence is the ordered sequence of all Gaussian points sorted by depth, and α m_confidence Let α be the opacity at the m-th Gaussian point. n_confidence Let be the opacity of the nth Gaussian point, m_confidence be the mth Gaussian point in the ordered sequence, and n_confidence be the nth Gaussian point in the ordered sequence.

10. The method for calculating building height based on aerial photography data using 3D Gaussian technology as described in claim 9, characterized in that, The confidence function formula is: ， Where threshold=0.999.

Images