Full-automatic calibration method for internal and external parameters of camera of unmanned aerial vehicle laser scanning system

Abstract

The invention relates to a full-automatic calibration method for internal and external parameters of a camera of an unmanned aerial vehicle laser scanning system, and the method comprises the following steps: extracting feature points and descriptors of each image obtained by the unmanned aerial vehicle laser scanning system, carrying out the flight strip division of the images, and matching homonymy points in the flight strips and between the flight strips; based on the matched homonymy points, using a bundle adjustment method to obtain internal parameters and distortion parameters of the camera; carrying out distortion correction on the image, and searching a three-dimensional point closest to the same-name light beam in the laser point cloud through a ray tracing method; and minimizing the distance between the three-dimensional points with the same name, and solving the external parameters of the camera. According to the method, the internal parameters and the external parameters of the camera can be calculated only through the image data, the point cloud data and the POS data collected by the unmanned aerial vehicle during aerial survey operation, control points and artificial markers do not need to be arranged, and the calibration process is fully automatic; dependence of calibration of internal parameters and external parameters of the camera of the unmanned aerial vehicle on a specific calibration object, a professional calibration site and man-machine interaction is eliminated.

Description

Technical Field

[0001] The present application relates to the technical field of UAV parameter calibration, and in particular to a method for fully automatic calibration of internal and external parameters of a camera of a UAV laser scanning system. Background Art

[0002] In traditional UAV photogrammetry, when image control points are included, the camera's intrinsic and extrinsic parameters can be obtained by self-calibration in the bundle adjustment with control points by manually marking the control points in the image. The calibration process requires manual marking of control points, which is a semi-automatic method and difficult to achieve full automation. In UAV laser scanning systems, control points are usually not measured, but 3D laser point cloud data is included. The camera's intrinsic and distortion parameters can be obtained through POS-assisted bundle adjustment, but the traditional bundle adjustment with control points cannot be used to fully automatically calibrate the camera's extrinsic parameters.

[0003] Indoor calibration of camera extrinsic parameters usually requires the use of a calibration object with known size. By establishing the correspondence between the known coordinates of the calibration object and its image points, a certain algorithm is used to obtain the intrinsic and extrinsic parameters of the camera model. For example, the "hand-eye" calibration method proposed by Radu Horaud and Fadi Dornaika in the paper "Hand-eye calibration" published in 1995, in which the "hand" (GPS / IMU) and "eye" (LiDAR / camera) are fixed on a machine. When the machine moves, the posture changes of the "hand" and "eye" meet certain constraints. Based on this, the coordinate transformation relationship between the "hand" and "eye" can be solved to solve the extrinsic parameters. The traditional camera calibration method needs to be calibrated in the laboratory, and the calibration process requires the use of a high-precision calibration plate. The calibration process is complicated and difficult for customers to complete the calibration independently after the laser scanning system is delivered to customers. Summary of the invention

[0004] The purpose of the present invention is to propose a fully automatic calibration method for the internal and external parameters of the camera of a UAV laser scanning system. Only the image data, point cloud data and POS data collected by the UAV during aerial survey operations can be used to solve the camera's internal and external parameters. No manual intervention or layout of control points is required, and the calibration of the UAV camera's internal and external parameters is freed from dependence on specific calibration objects and professional calibration sites.

[0005] The technical solution adopted by the present invention is: a method for fully automatic calibration of internal and external parameters of a camera of a UAV laser scanning system, comprising the following steps:

[0006] S1: Extract feature points and descriptors from each image acquired by the drone, divide the drone image into flight strips, and match the same-name points within and between flight strips;

[0007] S2: Based on the matched synonymous points in step S1, the intrinsic parameters and distortion parameters of the camera are obtained using bundle adjustment;

[0008] S3: After the image is corrected for distortion, the 3D point closest to the beam of the same name is searched in the 3D laser point cloud by ray tracing;

[0009] S4: Minimize the distance between 3D points with the same name and solve the camera extrinsic parameters.

[0010] Furthermore, the specific method of step S1 is to use the SIFT algorithm to extract SIFT feature points and descriptors of the image, and perform SIFT feature point matching between overlapping images within the flight strip and between flight strips to obtain image points with the same name.

[0011] Furthermore, the specific steps of step S2 are:

[0012] S201: Obtaining intrinsic parameters and distortion parameters of a camera by using bundle adjustment, wherein the intrinsic parameters include a principal point, a focal length, and distortion parameters of the camera, and the distortion parameters include radial distortion parameters and tangential distortion parameters;

[0013] S202: Perform distortion correction on the image using the obtained distortion parameters. The specific formula is:

[0014] x dt =x(1+k1r 2 +k2r 4 +k3r 6 )+2p1xy+p2(r 2 +2x 2 )

[0015] y dt =y(1+k1r 2 +k2r 4 +k3r 6 )+p1(r 2 +2y 2 )+2p2xy;

[0016] Among them, x dt Represents the horizontal coordinate of the image after distortion correction, y dt represents the ordinate of the image after distortion correction, x represents the abscissa of the image before distortion correction, y represents the ordinate of the image after distortion correction, k1 represents the first radial distortion parameter, k2 represents the second radial distortion parameter, k3 represents the third radial distortion parameter, p1 represents the first tangential distortion parameter, p2 represents the second tangential distortion parameter, and r represents the distance between the image point and the image principal point.

[0017] Furthermore, the specific steps of step S3 are:

[0018] S301: The camera photography centers in the two images are denoted as O1 and O2, and a pair of SIFT feature points f1 and f2 are selected in the two images, and a ray is emitted from the camera photography centers O1 and O2 respectively, passing through the SIFT feature points f1 and f2 respectively, and intersecting with the three-dimensional laser point cloud respectively;

[0019] S302: Using an octree structure, search from the three-dimensional laser point cloud for the three-dimensional point Q1 with the same name that is closest to the ray passing through the camera photography center O1 and the SIFT feature point f1, and the three-dimensional point Q2 with the same name that is closest to the ray passing through the camera photography center O2 and the SIFT feature point f2. The specific process is: take the smallest cube containing the three-dimensional laser point cloud as the root node, and then divide the cube into eight sub-cubes until each sub-cube contains at most a certain number of three-dimensional points; take the ray formed by the camera photography center and the SIFT feature point as the axis, and construct a cylinder with a radius of R; search for the point closest to the ray based on the octree structure, and determine whether the cylinder intersects with the sub-cube node of the octree. If the cylinder does not intersect with the node, skip the node and its sub-nodes. If the cylinder intersects with a sub-cube, traverse the sub-nodes of the node until the leaf node. For all points in the leaf node, calculate the distance d from each point to the ray, and select the point with the smallest distance as the final search result.

[0020] Furthermore, the specific steps of step S4 are:

[0021] S401: The images where the SIFT feature points f1 and f2 are located are recorded as M1 and M2 respectively, the time when the images M1 and M2 are acquired are recorded as t1 and t2, and the POS position of the camera at time t1 is recorded as Posture is recorded as The POS position of the camera at time t2 is recorded as Posture is recorded as The external parameters of the camera include the translation and rotation angle between the camera coordinate system and the inertial navigation coordinate system. The translation between the camera coordinate system and the inertial navigation coordinate system is (l x , l y , l z ), the rotation angle is recorded as in, Indicates the X coordinate of the camera's POS position at time t1, Indicates the Y coordinate of the camera's POS position at time t1, Indicates the Z coordinate of the camera's POS position at time t1, Indicates the roll angle of the camera's POS posture at time t1, Indicates the pitch angle of the camera's POS posture at time t1, Indicates the heading angle of the camera's POS posture at time t1, represents the X coordinate of the camera's POS position at time t2, Indicates the Y coordinate of the camera's POS position at time t2, Indicates the Z coordinate of the camera's POS position at time t2, Indicates the roll angle of the camera's POS posture at time t2, represents the pitch angle of the camera's POS posture at time t2, Indicates the heading angle of the camera's POS posture at time t2, l x Indicates the offset in the x direction in the external parameters, l y Indicates the offset in the y direction in the external parameters, l z Indicates the offset in the z direction in the external parameters, represents the rotation angle around the x-axis in the external parameters, ω represents the rotation angle around the y-axis in the external parameters, and k represents the rotation angle around the z-axis in the external parameters;

[0022] S402: According to the positioning principle of drone laser scanning, we can get:

[0023]

[0024] in, represents the coordinates of the three-dimensional point Q1 with the same name, represents the coordinates of the three-dimensional point Q2 with the same name, and R(·) represents the rotation matrix; represents the X coordinate of the 3D point Q1 with the same name, represents the Y coordinate of the 3D point Q1 with the same name, represents the Z coordinate of the 3D point Q1 with the same name, represents the X coordinate of the three-dimensional point Q2 with the same name, represents the Y coordinate of the 3D point Q2 with the same name, represents the Z coordinate of the three-dimensional point Q2 with the same name; represents the X coordinate of Q1 in the camera coordinate system, represents the Y coordinate of Q1 in the camera coordinate system, represents the Z coordinate of Q1 in the camera coordinate system, represents the X coordinate of Q2 in the camera coordinate system, represents the Y coordinate of Q2 in the camera coordinate system, Indicates the Z coordinate of Q2 in the camera coordinate system;

[0025] Substitute the initial external parameters and the coordinates of the same-name 3D point Q1 and the same-name 3D point Q2 into the above formula to calculate the coordinates of the same-name 3D point Q1 and the same-name 3D point Q2 in the camera coordinate system. and

[0026] S403: Obtain the distance between the three-dimensional point Q1 and the three-dimensional point Q2 with the same name by the formula in step S402

[0027] S404: By minimizing the distances between multiple pairs of 3D points with the same name and solving the external parameters, the specific formula is:

[0028]

[0029] Among them, x * represents the external parameters of the solution, represents the distance between the i-th pair of 3D points with the same name Q1Q2, N represents the number of 3D points with the same name, and argmin represents minimization; the formula is iteratively solved by the LM nonlinear optimization method, and the value of the external parameter at the j-th iteration is x j , the objective function f(x) is at x j The Jacobian matrix at is J(x j ), then the increment Δx of the parameter vector to be solved at this time is:

[0030] Δx=-[J(x j ) T J(x j )+μI] -1 J(x j ) T f(x j );

[0031] Among them, μ is the damping coefficient, which is adaptively and dynamically adjusted during the solution process, and I is the unit matrix;

[0032] After obtaining the increment of the parameter vector, update the parameters to be solved:

[0033] x j+1 =x j +Δx;

[0034] Repeat the above steps until the change in the parameter to be solved is less than the threshold or the maximum number of iterations is reached to obtain the optimal solution x * .

[0035] The beneficial effects of the present invention are:

[0036] The present invention calibrates the intrinsic parameters and distortion parameters of the camera through bundle adjustment, uses the camera distortion parameters to correct the distortion of the image, and uses laser point cloud as an aid to solve the extrinsic parameters of the camera by minimizing the distance between the three-dimensional point clouds corresponding to the same-name matching points of the image, thereby realizing simple and full-automatic calibration of the intrinsic and extrinsic parameters of the camera. Compared with the existing camera intrinsic and extrinsic parameter calibration methods in unmanned aerial vehicle laser scanning systems, which have the problems of low automation, many restrictions, and cumbersome calibration process, the present invention only needs the image data, point cloud data, and POS data collected by the unmanned aerial vehicle during aerial survey operations to fully automatically solve the camera intrinsic and extrinsic parameters through an algorithm without manual intervention, and at the same time gets rid of the dependence of the calibration of the intrinsic and extrinsic parameters of the unmanned aerial vehicle camera on specific calibration objects and professional calibration sites, and does not require manual layout of control points, saving manpower and time costs. BRIEF DESCRIPTION OF THE DRAWINGS

[0037] In order to more clearly illustrate the technical solutions in the embodiments of the present invention, the drawings required for use in the embodiments will be briefly introduced below. Obviously, the drawings described below are only some embodiments of the present application. For ordinary technicians in this field, other drawings can be obtained based on these drawings without paying creative work.

[0038] Figure 1 is a flow chart of a method according to an embodiment of the present invention;

[0039] Figure 2 This is a schematic diagram of the principle of calibration of internal and external parameters in an embodiment of the present invention;

[0040] Figure 3 Schematic diagram of the structure of an octree in an embodiment of the present invention;

[0041] Figure 4 A schematic diagram of the principle of ray tracing to search for three-dimensional points of the same name in an embodiment of the present invention;

[0042] Figure 5 Schematic diagram of the coloring effect of the point cloud using the parameters before and after calibration in an embodiment of the present invention Figure 1 ;

[0043] Figure 6 Schematic diagram of the coloring effect of the point cloud using the parameters before and after calibration in an embodiment of the present invention Figure 2 . DETAILED DESCRIPTION

[0044] In order to more clearly understand the above-mentioned purpose, features and advantages of the present invention, the present invention is further described in detail below in conjunction with the accompanying drawings and specific embodiments. In the following description, many specific details are set forth to facilitate a full understanding of the present invention, but the present invention can also be implemented in other ways different from those described herein, and therefore, the present invention is not limited to the limitations of the specific embodiments disclosed below.

[0045] like Figure 1 As shown, a method for fully automatic calibration of internal and external parameters of a camera of a UAV laser scanning system comprises the following steps:

[0046] S1: Extract feature points and descriptors from each image scanned by the drone, divide the drone images into flight strips, and match the same-name points within and between flight strips; specifically, the SIFT algorithm can be used to extract SIFT feature points and descriptors of the image, and SIFT feature points can be matched between overlapping images within and between flight strips to obtain the same-name image points of the image.

[0047] S2: Based on the matched synonymous points in step S1, use bundle adjustment to obtain the intrinsic parameters and distortion parameters of the camera; the specific steps are:

[0048] S201: Obtaining intrinsic parameters and distortion parameters of a camera by using bundle adjustment, wherein the intrinsic parameters include a principal point, a focal length, and distortion parameters of the camera, and the distortion parameters include radial distortion parameters and tangential distortion parameters.

[0049] S202: Perform distortion correction on the image using the obtained distortion parameters. The specific formula is:

[0050] x dt =x(1+k1r 2 +k2r 4 +k3r 6 )+2p1xy+p2(r 2 +2x 2 )

[0051] y dt =y(1+k1r 2 +k2r 4 +k3r 6 )+p1(r 2 +2y 2 )+2p2xy;

[0052] Among them, x dt Represents the horizontal coordinate of the image after distortion correction, y dtrepresents the ordinate of the image after distortion correction, x represents the abscissa of the image before distortion correction, y represents the ordinate of the image after distortion correction, k1 represents the first radial distortion parameter, k2 represents the second radial distortion parameter, k3 represents the third radial distortion parameter, p1 represents the first tangential distortion parameter, p2 represents the second tangential distortion parameter, and r represents the distance between the image point and the image principal point.

[0053] S3: After the image is distorted, the ray tracing method is used to search for the 3D point closest to the beam of the same name in the 3D laser point cloud; the specific steps are:

[0054] S301: Let the camera photography center in the two images be O1 and O2, and select a pair of SIFT feature points f1 and f2 in the two images, and emit a ray from the camera photography center O1 and O2 respectively, passing through the SIFT feature points f1 and f2 respectively, and intersecting with the 3D laser point cloud respectively. If there is no system, then theoretically the two rays should intersect with the 3D laser point cloud at the same point; if there is an error in the external parameters, then the two rays intersect with the 3D laser point cloud at Figure 2 The two different points shown are the same-name 3D points Q1 and Q2. The same-name 3D points Q1 and Q2 are functions of the image intrinsic parameters, POS pose and extrinsic parameters. Therefore, when the image intrinsic parameters and POS pose are known, the extrinsic parameters can be solved by minimizing the distance between the same-name 3D points.

[0055] S302: Using an octree structure to search through ray tracing from the 3D laser point cloud, the 3D point Q1 with the same name closest to the ray passing through the camera photography center O1 and the SIFT feature point f1 and the 3D point Q2 with the same name closest to the ray passing through the camera photography center O2 and the SIFT feature point f2 are searched. The specific process is as follows: Figure 3 As shown in , the smallest cube containing the 3D laser point cloud is taken as the root node, and then the cube is divided into eight sub-cubes until each sub-cube contains at most a certain number of 3D points. Figure 4 As shown, a cylinder with a radius of R is constructed with the ray formed by the camera photography center and the SIFT feature points as the axis; based on the octree structure, the point closest to the ray is searched to determine whether the cylinder intersects with the sub-cube node of the octree. If the cylinder does not intersect with the node, the node and its sub-nodes are skipped; if the cylinder intersects with a sub-cube, the sub-nodes of the node are traversed until the leaf node. For all points in the leaf node, the distance d from each point to the ray is calculated, and the point with the smallest distance is selected as the final search result.

[0056] S4: Minimize the distance between the three-dimensional points with the same name and solve the camera external parameters; the specific steps are:

[0057] S401: The images where the SIFT feature points f1 and f2 are located are recorded as M1 and M2 respectively, the time when the images M1 and M2 are acquired are recorded as t1 and t2, and the POS position of the camera at time t1 is recorded as Posture is recorded as The POS position of the camera at time t2 is recorded as Posture is recorded as The external parameters of the camera include the translation and rotation angle between the camera coordinate system and the inertial navigation coordinate system. The translation between the camera coordinate system and the inertial navigation coordinate system is (l x , l y , l z ), the rotation angle is recorded as in, Indicates the X coordinate of the camera's POS position at time t1, Indicates the Y coordinate of the camera's POS position at time t1, Indicates the Z coordinate of the camera's POS position at time t1, Indicates the roll angle of the camera's POS posture at time t1, Indicates the pitch angle of the camera's POS posture at time t1, Indicates the heading angle of the camera's POS posture at time t1, represents the X coordinate of the camera's POS position at time t2, Indicates the Y coordinate of the camera's POS position at time t2, Indicates the Z coordinate of the camera's POS position at time t2, Indicates the roll angle of the camera's POS posture at time t2, represents the pitch angle of the camera's POS posture at time t2, Indicates the heading angle of the camera's POS posture at time t2, l x Indicates the offset in the x direction in the external parameters, l y Indicates the offset in the y direction in the external parameters, l z Indicates the offset in the z direction in the external parameters, represents the rotation angle around the x-axis in the external parameters, ω represents the rotation angle around the y-axis in the external parameters, and k represents the rotation angle around the z-axis in the external parameters.

[0058] S402: According to the positioning principle of drone laser scanning, we can get:

[0059]

[0060] in, represents the coordinates of the three-dimensional point Q1 with the same name, represents the coordinates of the three-dimensional point Q2 with the same name, and R(·) represents the rotation matrix; represents the X coordinate of the 3D point Q1 with the same name, represents the Y coordinate of the 3D point Q1 with the same name, represents the Z coordinate of the 3D point Q1 with the same name, represents the X coordinate of the three-dimensional point Q2 with the same name, represents the Y coordinate of the 3D point Q2 with the same name, represents the Z coordinate of the three-dimensional point Q2 with the same name; represents the X coordinate of Q1 in the camera coordinate system, represents the Y coordinate of Q1 in the camera coordinate system, represents the Z coordinate of Q1 in the camera coordinate system, represents the X coordinate of Q2 in the camera coordinate system, represents the Y coordinate of Q2 in the camera coordinate system, Indicates the Z coordinate of Q2 in the camera coordinate system;

[0061] Substitute the initial external parameters and the coordinates of the same-name 3D point Q1 and the same-name 3D point Q2 into the above formula to calculate the coordinates of the same-name 3D point Q1 and the same-name 3D point Q2 in the camera coordinate system. and

[0062] S403: Obtain the distance between the three-dimensional point Q1 and the three-dimensional point Q2 with the same name by the formula in step S402

[0063] S404: By minimizing the distances between multiple pairs of 3D points with the same name and solving the external parameters, the specific formula is:

[0064]

[0065] Among them, x * represents the external parameters of the solution, represents the distance between the i-th pair of 3D points with the same name Q1Q2, N represents the number of 3D points with the same name, and argmin represents minimization; the formula is iteratively solved by the LM nonlinear optimization method, and the value of the external parameter at the j-th iteration is x j , the objective function f(x) is at x j The Jacobian matrix at is J(x j ), then the increment Δx of the parameter vector to be solved at this time is:

[0066] Δx=-[J(x j ) T J(x j )+μI] -1 J(x j ) T f(x j );

[0067] Among them, μ is the damping coefficient, which is adaptively and dynamically adjusted during the solution process, and I is the unit matrix;

[0068] After obtaining the increment of the parameter vector, update the parameters to be solved:

[0069] x j+1 =x j +Δx;

[0070] Repeat the above steps until the change in the parameter to be solved is less than the threshold or the maximum number of iterations is reached to obtain the optimal solution x * .

[0071] Figure 5 and Figure 6 Schematic diagram of the effect of coloring the point cloud using the parameters before and after calibration, where: Figure 5 and Figure 6 The left side of the figure is a schematic diagram of the effect of coloring the point cloud using the parameters before calibration, and the right side is a schematic diagram of the effect of coloring the point cloud using the parameters after calibration. By comparing the coloring effects before and after camera calibration, it can be seen that there are misalignments and discontinuities in the coloring results before calibration, while there are no misalignments and discontinuities in the coloring results using the parameters after calibration described in the embodiment of the present invention.

[0072] At present, the camera intrinsic parameters and extrinsic parameters of most unmanned aerial vehicle laser scanning systems are calibrated by the system provider. This approach has two major problems: (1) For the system provider, it is necessary to calibrate the camera extrinsic parameters of each unmanned aerial vehicle laser scanning system. The traditional method of calibration requires human-computer interaction, which is very labor-intensive; (2) For users, during the long-term use of the system, the camera's intrinsic parameters and extrinsic parameters may undergo slight changes, resulting in a decrease in data quality, while the user cannot perform calibration. The embodiment of the present invention only requires the image data, point cloud data and POS data collected by the unmanned aerial vehicle during aerial survey operations to calculate the camera's intrinsic parameters and extrinsic parameters, calibrates the camera's intrinsic parameters and distortion parameters through bundle adjustment, uses the camera's distortion parameters to correct the image's distortion, and uses the laser point cloud as an aid to calculate the camera's extrinsic parameters by minimizing the distance between the three-dimensional point clouds corresponding to the same-name matching points in the image, thereby achieving simple and fully automatic calibration of the camera's intrinsic parameters and extrinsic parameters. There is no need for manual intervention or deployment of control points, and the calibration of the internal and external parameters of the drone camera is no longer dependent on specific calibration objects or professional calibration sites, saving manpower and time costs.

[0073] By adopting the method described in the embodiment of the present invention, on the one hand, the system provider can efficiently calibrate the intrinsic and extrinsic parameters of the camera, or even provide the system directly to the user without calibrating the intrinsic and extrinsic parameters, and the user completes the calibration when collecting test data or training data, which greatly reduces the workload of the system provider and improves the work efficiency of the system provider; on the other hand, during the long-term use of the drone laser scanning system, the user can recalibrate at any time according to actual conditions, which can effectively improve the data quality of the drone laser scanning system.

[0074] The above description is only a preferred embodiment of the present invention and is not intended to limit the present invention. For those skilled in the art, the present invention may have various modifications and variations. Any modification, equivalent replacement, improvement, etc. made within the spirit and principle of the present invention shall be included in the protection scope of the present invention.

Claims

1. A method for fully automatic calibration of internal and external parameters of a drone laser scanning system camera, characterized in that: The steps include: S1: Extract feature points and descriptors from each image acquired by the drone, divide the drone image into flight strips, and match the same-name points within and between flight strips; S2: Based on the matched synonymous points in step S1, the intrinsic parameters and distortion parameters of the camera are obtained using bundle adjustment; S3: After the image is corrected for distortion, the 3D point closest to the beam of the same name is searched in the 3D laser point cloud by ray tracing; S4: Minimize the distance between 3D points with the same name and solve the camera extrinsic parameters.

2. The method for fully automatic calibration of internal and external parameters of a drone laser scanning system camera according to claim 1, characterized in that: The specific method of step S1 is to use the SIFT algorithm to extract SIFT feature points and descriptors of the image, and perform SIFT feature point matching between overlapping images within the flight strip and between flight strips to obtain image points with the same name.

3. The method for fully automatic calibration of internal and external parameters of a drone laser scanning system camera according to claim 2, characterized in that: The specific steps of step S2 are: S201: Obtaining intrinsic parameters and distortion parameters of a camera by using bundle adjustment, wherein the intrinsic parameters include a principal point, a focal length, and distortion parameters of the camera, and the distortion parameters include radial distortion parameters and tangential distortion parameters; S202: Perform distortion correction on the image using the obtained distortion parameters. The specific formula is: x dt =x(1+k1r 2 +k2r 4 +k3r 6 )+2p1xy+p2(r 2 +2x 2 ) <h2 style=";text-align:left;direction:ltr">y<h2 style=";text-align:left;direction:ltr"> dt <h2 style=";text-align:left;direction:ltr"> =y(1+k1r<h2 style=";text-align:left;direction:ltr"> 2 <h2 style=";text-align:left;direction:ltr"> +k2r<h2 style=";text-align:left;direction:ltr"> 4 <h2 style=";text-align:left;direction:ltr"> +k3r<h2 style=";text-align:left;direction:ltr"> 6 <h2 style=";text-align:left;direction:ltr"> )+p1(r<h2 style=";text-align:left;direction:ltr"> 2 <h2 style=";text-align:left;direction:ltr"> +2y<h2 style=";text-align:left;direction:ltr"> 2 <h2 style=";text-align:left;direction:ltr"> )+2p2xy; Among them, x dt Represents the horizontal coordinate of the image after distortion correction, y dt represents the ordinate of the image after distortion correction, x represents the abscissa of the image before distortion correction, y represents the ordinate of the image after distortion correction, k1 represents the first radial distortion parameter, k2 represents the second radial distortion parameter, k3 represents the third radial distortion parameter, p1 represents the first tangential distortion parameter, p2 represents the second tangential distortion parameter, and r represents the distance between the image point and the image principal point.

4. The method for fully automatic calibration of internal and external parameters of a drone laser scanning system camera according to claim 3, characterized in that: The specific steps of step S3 are: S301: The camera photography centers in the two images are denoted as O1 and O2, and a pair of SIFT feature points f1 and f2 are selected in the two images, and a ray is emitted from the camera photography centers O1 and O2 respectively, passing through the SIFT feature points f1 and f2 respectively, and intersecting with the three-dimensional laser point cloud respectively; S302: Using an octree structure, search from the three-dimensional laser point cloud for the three-dimensional point Q1 with the same name that is closest to the ray passing through the camera photography center O1 and the SIFT feature point f1, and the three-dimensional point Q2 with the same name that is closest to the ray passing through the camera photography center O2 and the SIFT feature point f2. The specific process is: take the smallest cube containing the three-dimensional laser point cloud as the root node, and then divide the cube into eight sub-cubes until each sub-cube contains at most a certain number of three-dimensional points; take the ray formed by the camera photography center and the SIFT feature point as the axis, and construct a cylinder with a radius of R; search for the point closest to the ray based on the octree structure, and determine whether the cylinder intersects with the sub-cube node of the octree. If the cylinder does not intersect with the node, skip the node and its sub-nodes. If the cylinder intersects with a sub-cube, traverse the sub-nodes of the node until the leaf node. For all points in the leaf node, calculate the distance d from each point to the ray, and select the point with the smallest distance as the final search result.

5. The method for fully automatic calibration of internal and external parameters of a drone laser scanning system camera according to claim 4, characterized in that: The specific steps of step S4 are: S401: The images where the SIFT feature points f1 and f2 are located are recorded as M1 and M2 respectively, the time when the images M1 and M2 are acquired are recorded as t1 and t2, and the POS position of the camera at time t1 is recorded as Posture is recorded as The POS position of the camera at time t2 is recorded as Posture is recorded as The external parameters of the camera include the translation and rotation angle between the camera coordinate system and the inertial navigation coordinate system. The translation between the camera coordinate system and the inertial navigation coordinate system is (l x , l y , l z ), the rotation angle is recorded as in, Indicates the X coordinate of the camera's POS position at time t1, Indicates the Y coordinate of the camera's POS position at time t1, Indicates the Z coordinate of the camera's POS position at time t1, Indicates the roll angle of the camera's POS posture at time t1, Indicates the pitch angle of the camera's POS posture at time t1, Indicates the heading angle of the camera's POS posture at time t1, represents the X coordinate of the camera's POS position at time t2, Indicates the Y coordinate of the camera's POS position at time t2, Indicates the Z coordinate of the camera's POS position at time t2, Indicates the roll angle of the camera's POS posture at time t2, represents the pitch angle of the camera's POS posture at time t2, Indicates the heading angle of the camera's POS posture at time t2, l x Indicates the offset in the x direction in the external parameters, l y Indicates the offset in the y direction in the external parameters, l z Indicates the offset in the z direction in the external parameters, represents the rotation angle around the x-axis in the external parameters, ω represents the rotation angle around the y-axis in the external parameters, and k represents the rotation angle around the z-axis in the external parameters; S402: According to the positioning principle of drone laser scanning, we can get: in, represents the coordinates of the three-dimensional point Q1 with the same name, represents the coordinates of the three-dimensional point Q2 with the same name, and R(·) represents the rotation matrix; represents the X coordinate of the 3D point Q1 with the same name, represents the Y coordinate of the 3D point Q1 with the same name, represents the Z coordinate of the 3D point Q1 with the same name, represents the X coordinate of the three-dimensional point Q2 with the same name, represents the Y coordinate of the 3D point Q2 with the same name, represents the Z coordinate of the three-dimensional point Q2 with the same name; represents the X coordinate of Q1 in the camera coordinate system, represents the Y coordinate of Q1 in the camera coordinate system, represents the Z coordinate of Q1 in the camera coordinate system, represents the X coordinate of Q2 in the camera coordinate system, represents the Y coordinate of Q2 in the camera coordinate system, Indicates the Z coordinate of Q2 in the camera coordinate system; Substitute the initial external parameters and the coordinates of the same-name 3D point Q1 and the same-name 3D point Q2 into the above formula to calculate the coordinates of the same-name 3D point Q1 and the same-name 3D point Q2 in the camera coordinate system. and S403: Obtain the distance between the three-dimensional point Q1 and the three-dimensional point Q2 with the same name by the formula in step S402 S404: By minimizing the distances between multiple pairs of 3D points with the same name and solving the external parameters, the specific formula is: Among them, x * represents the external parameters of the solution, represents the distance between the i-th pair of 3D points with the same name Q1Q2, N represents the number of 3D points with the same name, and argmin represents minimization; the formula is iteratively solved by the LM nonlinear optimization method, and the value of the external parameter at the j-th iteration is x j , the objective function f(x) is at x j The Jacobian matrix at is J(x j ), then the increment Δx of the parameter vector to be solved at this time is: Δx=-[J(x j ) T J(x j )+μI] -1 J(x j ) T f(x j ); Among them, μ is the damping coefficient, which is adaptively and dynamically adjusted during the solution process, and I is the unit matrix; After obtaining the increment of the parameter vector, update the parameters to be solved: x j+1 =x j +Δx; Repeat the above steps until the change in the parameter to be solved is less than the threshold or the maximum number of iterations is reached to obtain the optimal solution x * .

Images