Unmanned aerial vehicle blind landing area autonomous detection method based on monocular vision

Abstract

The invention relates to the technical field of information technology and environmental perception, and discloses an unmanned aerial vehicle blind landing area autonomous detection method based on monocular vision, which can perform data processing only by using optical video data acquired by monocular vision without relying on landmark identification and ground priori knowledge. An elevation model of a flight scene is established, a three-dimensional terrain is analyzed autonomously, and various terrain factors are comprehensively considered to determine a safe landing area.

Description

Technical Field

[0001] The present invention relates to the fields of information technology and environmental perception technology, and in particular to an autonomous detection method for a blind landing area of ​​a UAV based on monocular vision. Background Art

[0002] The autonomous landing capability of drones is crucial for many application scenarios such as field search and rescue, reconnaissance, environmental monitoring, geological exploration, and disaster relief. They need to accurately identify, approach and land in a safe area without human intervention.

[0003] The Global Navigation Satellite System can provide three-dimensional spatial position parameters for drones. However, relying solely on GNSS navigation cannot obtain terrain and obstacle information in the flight area. In addition, in environments such as forests, valleys, and high-density high-rise buildings, GNSS positioning and navigation signals will be interfered with, and erroneous positioning and navigation data will affect the reliability of the system.

[0004] LiDAR can obtain 3D point cloud data of the scene. By designing a candidate plane extraction algorithm for point cloud data, a flat area that meets the landing requirements can be selected. However, such sensors are expensive, heavy, and not suitable for small drones. In contrast, optical cameras have the advantages of light weight and low cost, making them suitable as environmental perception sensors on drones.

[0005] Vision-based UAV autonomous landing technology has been developed for many years. Currently, it can be divided into two categories according to technical means: methods based on visual identification and methods based on ground feature extraction and matching.

[0006] (1) Methods based on visual landmark recognition. This type of method uses computer vision algorithms to identify specific landmarks placed on the landing area, and solves the flight posture parameters of the drone through operations such as feature extraction and matching. However, the landmark-based drone landing area detection method requires the design of different indicator landmarks according to different types of application scenarios, and the detection algorithm also needs to be designed according to the geometric characteristics of the indicator landmarks. Therefore, the drone can only achieve safe landing in a given scenario. This type of algorithm has poor stability and robustness, and is not suitable for use in scenarios where many types of indicator landmarks cannot be placed in advance, such as post-disaster rescue and emergency landing of drones due to malfunctions.

[0007] (2) Methods based on ground feature extraction and matching. This type of method uses the drone's visual system to obtain images of the landing area, extract ground features (such as edges, corners, etc.), match them with ground templates or previously established maps, and then find the ideal landing area based on the matching results. This type of method does not require pre-setting of landmarks. However, this type of method has high requirements for the accuracy of the ground template. If there is a large difference between the ground template and the actual scene, the matching result may be wrong, resulting in the inability to accurately find the ideal landing area.

[0008] When dealing with the challenge of safe landing of drones, in addition to considering the position and attitude information of the drone, the three-dimensional terrain of the landing site is also a crucial aspect of environmental perception. During the landing process, the drone needs to fully understand the terrain characteristics of the landing site, including terrain height, slope, obstacle distribution and other information to ensure the accuracy and safety of the landing.

[0009] The problem to be solved by the invention is that when a drone detects a landing area, it only relies on GNSS navigation signals, which limits the acquisition of terrain information, and the stability of the signal will be seriously affected in some complex environments. The robustness of the drone landing area detection method based on landmarks is poor, and specific landmarks need to be pre-designed and placed according to different terrains. It is not suitable for post-disaster rescue, emergency landing of drones due to malfunctions, and other scenarios where multiple indicating landmarks cannot be placed in advance. The method based on ground feature extraction and matching has high requirements on the accuracy of the ground template. Summary of the invention

[0010] The purpose of the present invention is to design a real-time autonomous landing area detection method for UAV based on monocular vision. The method described can be independent of landmark identification and ground prior knowledge, and only uses optical video data collected by monocular vision for data processing; establishes an elevation model of the flight scene, autonomously analyzes the three-dimensional terrain, and comprehensively considers various terrain factors to determine a safe landing area.

[0011] The technical solution of the present invention is achieved in this way:

[0012] A method for autonomously detecting blind landing areas of unmanned aerial vehicles based on monocular vision comprises the following steps:

[0013] Use the calibration board to calibrate the relative position and posture of the camera and the parameters within the image to obtain the calibration parameters;

[0014] The camera continuously acquires sequence images of the flight scene, covering the entire pre-landing area;

[0015] Using the synchronous positioning and mapping algorithm, the relative position and real-time flight trajectory of the drone are restored through the epipolar geometric relationship between the sequence images, and the initial three-dimensional point cloud of the flight scene is reconstructed;

[0016] The initial 3D point cloud is resampled using an octree-based hierarchical data structure to generate a voxel grid map with connectivity;

[0017] The voxel grid map is color-coded according to the height information in each node cell to build an elevation model of the current flight scene;

[0018] According to the elevation model, the landing area selection criteria including the flatness measurement and area area of ​​the landing area are formulated to obtain the nearest landing area location of the UAV.

[0019] As a further solution of the present invention: the parameter calibration method is specifically as follows:

[0020] Calibration is performed using more than 12 calibration plate images taken by the camera from different angles and positions;

[0021] The parameters obtained after the calibration operation include the camera's principal distance, principal point, distortion, relative position and attitude angle parameters.

[0022] As a further solution of the present invention: the calibration plate can be a checkerboard calibration plate and a solid origin calibration plate.

[0023] As a further solution of the present invention: the initial three-dimensional point cloud reconstruction method is specifically:

[0024] 1) Calculate the epipolar geometric relationship between the sequence images;

[0025] Let I 1 and I 2 are two adjacent images, p 1 and p 2 Represents the projection pixel of point P in the physical space on the two images, p 1 and p 2 The coordinates of are obtained by transforming the coordinates of the physical space point P in the camera coordinate system according to the camera intrinsic parameter matrix k:

[0026]

[0027] Among them, z 1 and z 2 They are p 1 and p 2 The depth value of P′ is the depth of the spatial point P at I 1 The coordinates in the camera coordinate system, R, t is I 1 and I 2 The geometric transformation relationship of the spatial coordinate system between two adjacent images;

[0028] Multiply the left side of the above equation by k -1 get:

[0029]

[0031] Substituting (3) into (4), we get:

[0032] z 2 k -1 p 2 =Rz 1 k -1 p 1 +t (5)

[0033] Multiply both sides of the above equation by [t] × have to:

[0034] z 2 [t] × k -1 p 2 =Rz 1 [t] × k -1 p 1 (6)

[0035] Where [t] × is the antisymmetric matrix of t;

[0036] Multiply (6) on the left have to:

[0037]

[0038] in, So the left side of the equation is 0, and we can get the following equation:

[0039]

[0040] Therefore, the epipolar geometry constraints are obtained as follows:

[0041]

[0042] Where F = k -T [t] × R -1 is the basic matrix;

[0043] The essential matrix E is determined by the camera's external parameters and has nothing to do with the camera's internal parameters. Let the unit-scale normalized vectors of the spatial point P in the image space coordinate system of the two images be expressed as P 0 , P 1 ; According to the spatial geometric transformation relationship between the two images:

[0044] P 1 =RP 0 +t (10)

[0045] Multiply both sides of the equation by the cross product of t and then by P. 1 The dot product of , we can get:

[0046] 0=P 1 t×RP 0 (11)

[0047] Right now:

[0048] (P 1 ) T [t] × RP 0 =0 (12)

[0049] (P 1 ) T EP 0 =0 (13)

[0050] Where E = [t] × R represents the essential matrix, which contains the relative position and posture transformation relationship of the two images; the relationship between the essential matrix E and the basic matrix F is:

[0051] E=k T Fk (14)

[0052] 2) Further calculate the relative parameters of the image position and attitude of the UAV during flight;

[0053] First, initialize the image position and scene structure:

[0054] E is decomposed into an antisymmetric matrix S (equivalent to [t] × ) and the product of a rotation matrix R:

[0055] E=SR (15)

[0056] If S is a real antisymmetric matrix, then S = UBU T , where B is a block diagonal matrix; two temporary matrices are introduced:

[0057] and

[0058] In the case of indefinite sign, Z = diag(1, 1, 0)W, and in the case of indefinite scale, S = Udiag(1, 1, 0)WU T ;

[0059] Perform singular value decomposition on E:

[0060] E=SR=Udiag(1,1,0)(WU T R) (16)

[0061] Then, the camera matrix is ​​extracted from the essential matrix to determine the relative parameters of the image position and posture:

[0062] E=UΣV T ,∑=diag(σ 1 , σ 2 , σ 3 ) (17)

[0063] Among them, σ 1 =σ 2 ≠0 and σ 3 = 0, U, V∈SO(3);

[0064] In the case of a constant factor difference, S-Udiag(1,1,0)U T , assuming that the rotation matrix R can be decomposed into R = UXV T , where X is a rotation matrix, then:

[0065] Udiag(1,1,0)V T =E=SR=(UZU T )(UXV T )=U(ZX)V T (18)

[0066] Then we have ZX = diag(1, 1, 0), and because X is a rotation matrix, we have X = W or X = W T ;

[0067] If the SVD decomposition of E is Udiag(1,1,0)V T :

[0068]

[0069] At this point there are 4 possible camera matrices, namely:

[0070]

[0071] Only one solution can satisfy the requirement that the reconstructed point is in front of the camera, that is, there is only one reliable solution. The relative parameters of the image position and posture are determined by testing whether a point is in front of the camera.

[0072] 3) Using the position and attitude parameters of the image and the matching relationship between pixels, the three-dimensional point cloud data of the flight scene is reconstructed using the triangulation method;

[0073] According to epipolar geometry, we have:

[0074] s 2 p 2 =s 1 R 1 +t (21)

[0075] Among them, p 1 , p 2 are the normalized coordinates of the feature points on the two images, s 1 ,s 2 is the depth of the two feature points; multiply both sides of the above formula by

[0076]

[0077] According to the above formula, we can get s 1 ,s 2 ; Due to the influence of noise, the obtained R,t may not accurately make the above formula equal to 0, so it is actually solved by least squares;

[0078] The result of triangulation is the real-time calculation of the coordinates of the three-dimensional points in the physical world corresponding to the image matching pixels, thereby reconstructing the three-dimensional point cloud data of the flight scene.

[0079] As a further solution of the present invention: the octree-based hierarchical data structure resamples the three-dimensional point cloud data, comprising the steps of:

[0080] 1) Performing octree partitioning on the three-dimensional space containing the three-dimensional point cloud data to obtain a hierarchical octree data structure, each node of which corresponds to a three-dimensional cube cell;

[0081] First, determine the bounding box in three-dimensional space, which is the smallest cube that contains all the point cloud data;

[0082] Secondly, divide the root node into eight child nodes, each of which corresponds to one eighth of the volume of the root node; for each child node, check whether it contains enough points; if the number of points contained exceeds the threshold, further divide the child node into eight smaller child nodes; repeat this process until the number of points in the node is less than or equal to the threshold;

[0083] Finally, the point cloud data is assigned to the corresponding octree nodes. After recursive partitioning, the points will eventually be distributed in the leaf nodes of the octree.

[0084] 2) According to the three-dimensional cube cell corresponding to each node of the octree data structure, the average value of the z coordinates of the n three-dimensional points falling into it is used as the center height value of the node cell;

[0085]

[0086] in, is the center height value of the j-th node cell;

[0087] 3) Using the Gaussian probability model statistical method, the center height value of each node of the octree data structure is statistically calculated with the height values ​​of the six nodes above, below, left, right, and front and back, to obtain the statistical standard deviation value of each center node; the formula for calculating the statistical standard deviation value of the z coordinates of the i-node cells is:

[0088]

[0089] in, is the center height value of the central node, and σ is the statistical standard deviation value of the node.

[0090] As a further solution of the present invention: the elevation model is specifically established as follows: each node cell contains the altitude information of the location, the voxel grid map is color-coded according to the altitude information in each node cell, and the elevation model of the current flight scene is established;

[0091] The node cell is a volume pixel in three-dimensional space, each node cell represents a small area in space, and each node cell contains the height information of the location, which refers to the position or height of the ground or terrain in the vertical direction; the color coding is to encode the height information according to a gradient color band; the generated elevation map is a type of map that uses color to represent the height changes of the terrain; with the incremental construction of the point cloud, the elevation map can be dynamically updated and maintained.

[0092] As a further solution of the present invention: it is necessary to first screen out the horizontal plane area according to the elevation model, specifically:

[0093] In the generated scene elevation model, the RANSAC algorithm is used to extract the plane structure, and the horizontal plane area is screened out based on the normal vector direction constraints of each plane;

[0094] The input of the RANSAC algorithm is the center of each node cell in the elevation map; the constraint condition is the direction of the normal vector of each plane. The relationship between each plane and the horizontal plane is measured by calculating the angle between the normal vector of each plane and the unit vector perpendicular to the horizontal plane. When the angle is close to or less than 20 degrees, the plane is identified as a horizontal plane.

[0095] As a further solution of the present invention: the landing area selection criteria including the landing area flatness measurement and the area of ​​the area are formulated as follows:

[0096] The flatness measurement of the horizontal plane area refers to the calculation based on the statistical standard deviation of the cells contained in the horizontal plane area, and the sum of the standard deviations of all cells is calculated. If the sum result is less than the threshold, it is considered to meet the flatness measurement requirements; the calculation formula of the flatness measurement F is:

[0097]

[0098] The area of ​​the horizontal plane region refers to the projection area of ​​the ground occupied by the cells calculated based on the projection area of ​​the cells contained in the horizontal plane region. Only the regions with an area greater than the area threshold are considered to meet the area requirement.

[0099] Among all horizontal plane areas that meet the flatness measurement requirements and area requirements, the area with the largest plane area is selected as the final UAV landing area detection result.

[0100] As a further solution of the present invention: the acquisition of the landing area position is specifically as follows:

[0101] Update the landing area in real time as the drone flies, ensuring that the pre-selected landing area at the current moment is closest to the drone to improve landing efficiency and safety;

[0102] The landing area selection location is updated with the real-time flight of the drone. This means that when an area with the same plane area appears in the updated elevation map, the landing area selection result will be dynamically adjusted to the landing area location closest to the drone.

[0103] Beneficial effects of the present invention:

[0104] 1. The present invention can complete the detection of the landing area under the condition of satellite signal denial; it does not need to rely on landmarks and prior knowledge of the ground; it can establish an elevation model of the flight scene, independently analyze the three-dimensional terrain, and comprehensively consider various terrain factors to determine a safe landing area.

[0105] 2. The present invention designs a three-dimensional terrain elevation modeling step to model the flight scene of the UAV, and uses the elevation model to analyze the geometric changes of the three-dimensional terrain.

[0106] 3. The present invention designs a landing area selection strategy, which comprehensively considers multiple types of terrain factors to complete the detection of a safe landing area.

[0107] 4. The present invention breaks through the limitations of satellite navigation signals and ground prior knowledge, and can autonomously analyze changes in the surrounding environment and detect landing areas in real time. Its results can provide support for environmental perception-related application fields such as blind landing of drones in the wild, and has good market application prospects.

[0108] The present invention will be further described below in conjunction with the accompanying drawings and embodiments. BRIEF DESCRIPTION OF THE DRAWINGS

[0109] Figure 1 It is a schematic diagram of the process of the present invention;

[0110] Figure 2 It is a schematic diagram of plane structure extraction and horizontal ground screening of the present invention;

[0111] Figure 3 A schematic diagram of a three-dimensional terrain analysis performed by a full-area traversal of the present invention;

[0112] Figure 4 This is a schematic diagram of autonomous detection of the landing area of ​​the UAV of the present invention;

[0113] Figure 5 It is a schematic diagram of updating the landing area detection results at different times of the present invention. DETAILED DESCRIPTION

[0114] The technical solutions in the embodiments of the present invention will be described clearly and completely below in conjunction with the accompanying drawings in the embodiments of the present invention.

[0115] See attached Figure 1-5 The present invention provides a method for autonomously detecting a blind landing zone of a UAV based on monocular vision, which comprises the following steps:

[0116] Use the calibration board to calibrate the relative position and posture of the camera and the parameters within the image to obtain the calibration parameters;

[0117] The camera continuously acquires sequence images of the flight scene, covering the entire pre-landing area;

[0118] Using the synchronous positioning and mapping algorithm, the relative position and real-time flight trajectory of the drone are restored through the epipolar geometric relationship between the sequence images, and the initial three-dimensional point cloud of the flight scene is reconstructed;

[0119] The initial 3D point cloud is resampled using an octree-based hierarchical data structure to generate a voxel grid map with connectivity;

[0120] The voxel grid map is color-coded according to the height information in each node cell to build an elevation model of the current flight scene;

[0121] According to the elevation model, the landing area selection criteria including the flatness measurement and area area of ​​the landing area are formulated to obtain the nearest landing area location of the UAV.

[0122] The parameter calibration method is specifically as follows:

[0123] Calibration is performed using more than 12 calibration plate images taken by the camera from different angles and positions;

[0124] The parameters obtained after the calibration operation include the camera's principal distance, principal point, distortion, relative position and attitude angle parameters.

[0125] The calibration plate can be a checkerboard calibration plate or a solid origin calibration plate.

[0126] The initial three-dimensional point cloud reconstruction method is specifically as follows:

[0127] 1) Calculate the epipolar geometric relationship between the sequence images;

[0128] Let I 1 and I 2 are two adjacent images, p 1 and p 2 Represents the projection pixel of point P in the physical space on the two images, p 1 and p 2 The coordinates of are obtained by transforming the coordinates of the physical space point P in the camera coordinate system according to the camera intrinsic parameter matrix k:

[0129]

[0130] Among them, z 1 and z 2 They are p 1 and p 2 The depth value of P′ is the depth of the spatial point P at I 1 The coordinates in the camera coordinate system, R, t is I 1 and I 2 The geometric transformation relationship of the spatial coordinate system between two adjacent images;

[0131] Multiply the left side of the above equation by k -1 get:

[0132]

[0133] Substituting (3) into (4), we get:

[0134] z 2 k -1 p 2 =Rz 1 k -1 p 1 +i (5)

[0135] Multiply both sides of the above equation by [t] × have to:

[0136] z 2 [t] × p -2 p 2 =Rz 1 [t] × k -1 p 1 (6)

[0137] Where [t] × is the antisymmetric matrix of t;

[0138] Multiply (6) on the left have to:

[0139]

[0140] in, So the left side of the equation is 0, and we can get the following equation:

[0141]

[0142] Therefore, the epipolar geometry constraints are obtained as follows:

[0143]

[0144] Where F = k -T [t] × R -1 is the basic matrix;

[0145] The essential matrix E is determined by the camera's external parameters and has nothing to do with the camera's internal parameters. Let the unit-scale normalized vectors of the spatial point P in the image space coordinate system of the two images be expressed as P 0 , P 1 ; According to the spatial geometric transformation relationship between the two images:

[0146] P 1 =RP 0 +t (10)

[0147] Multiply both sides of the equation by the cross product of t and then by P. 1 The dot product of , we can get:

[0148] 0=P 1 t×RP 0 (11)

[0149] Right now:

[0150] (P 1 )T[t] × RP 0 =0 (12)

[0151] (P 1 ) T EP 0 =0 (13)

[0152] Where E = [T] × R represents the essential matrix, which contains the relative position and posture transformation relationship of the two images; the relationship between the essential matrix E and the basic matrix F is:

[0153] E=k T Fk (14)

[0154] 2) Further calculate the relative parameters of the image position and attitude of the UAV during flight;

[0155] First, initialize the image position and scene structure:

[0156] E is decomposed into an antisymmetric matrix S (equivalent to [t] × ) and the product of a rotation matrix R:

[0157] E=SR (15)

[0158] If S is a real antisymmetric matrix, then S = UBU T , where B is a block diagonal matrix; two temporary matrices are introduced:

[0159] and

[0160] In the case of indefinite sign, Z = diag(1, 1, 0)W, and in the case of indefinite scale, S = Udiag(1, 1, 0)WU T ;

[0161] Perform singular value decomposition on E:

[0162] E=SR=Udiag(1,1,0)(WU T R) (16)

[0163] Then, the camera matrix is ​​extracted from the essential matrix to determine the relative parameters of the image position and posture:

[0164] E=U∑V T ,∑=diag(σ 1 , σ 2 , σ 3 ) (17)

[0165] Among them, σ 1 =σ 2 ≠0 and σ s = 0, U, V∈SO(3);

[0166] In the case of a constant factor difference, S = Udiag (1, 1, 0) U T , assuming that the rotation matrix R can be decomposed into R = UXV T , where X is a rotation matrix, then:

[0167] Udiag(1,1,0)V T =E=SR=(UZU T )(UXV T )=U(ZX)V T (18)

[0168] Then we have ZX = diag(1, 1, 0), and because X is a rotation matrix, we have X = W or X = W T ;

[0169] If the SVD decomposition of E is Udiag(1,1,0)V T :

[0170]

[0171] At this point there are 4 possible camera matrices, namely:

[0172]

[0173] Only one solution can satisfy the requirement that the reconstructed point is in front of the camera, that is, there is only one reliable solution. The relative parameters of the image position and posture are determined by testing whether a point is in front of the camera.

[0174] 3) Using the position and attitude parameters of the image and the matching relationship between pixels, the three-dimensional point cloud data of the flight scene is reconstructed using the triangulation method;

[0175] According to epipolar geometry, we have:

[0176] s 2 p 2 =s 1 R 1 +t (21)

[0177] Among them, p 1 , p 2 are the normalized coordinates of the feature points on the two images, s 1 ,s 2 is the depth of the two feature points; multiply both sides of the above formula by

[0178]

[0179] According to the above formula, we can get s 1 s 2 ; Due to the influence of noise, the obtained R,t may not accurately make the above formula equal to 0, so it is actually solved by least squares;

[0180] The result of triangulation is the real-time calculation of the coordinates of the three-dimensional points in the physical world corresponding to the image matching pixels, thereby reconstructing the three-dimensional point cloud data of the flight scene.

[0181] The octree-based hierarchical data structure resamples the three-dimensional point cloud data, comprising the steps of:

[0182] 1) Performing octree partitioning on the three-dimensional space containing the three-dimensional point cloud data to obtain a hierarchical octree data structure, each node of which corresponds to a three-dimensional cube cell;

[0183] First, determine the bounding box in three-dimensional space, which is the smallest cube that contains all the point cloud data;

[0184] Secondly, divide the root node into eight child nodes, each of which corresponds to one eighth of the volume of the root node; for each child node, check whether it contains enough points; if the number of points contained exceeds the threshold, further divide the child node into eight smaller child nodes; repeat this process until the number of points in the node is less than or equal to the threshold;

[0185] Finally, the point cloud data is assigned to the corresponding octree nodes. After recursive partitioning, the points will eventually be distributed in the leaf nodes of the octree.

[0186] 2) According to the three-dimensional cube cell corresponding to each node of the octree data structure, the average value of the 2D coordinates of the n three-dimensional points falling into it is used as the center height value of the node cell;

[0187]

[0188] in, is the center height value of the j-th node cell;

[0189] 3) Using the Gaussian probability model statistical method, the center height value of each node of the octree data structure is statistically calculated with the height values ​​of the six nodes above, below, left, right, and front and back, to obtain the statistical standard deviation value of each center node; the formula for calculating the statistical standard deviation value of the z coordinates of the i-node cells is:

[0190]

[0191] in, is the center height value of the central node, and σ is the statistical standard deviation value of the node.

[0192] The elevation model is established specifically as follows: each node cell contains the altitude information of the location, the voxel grid map is color-coded according to the altitude information in each node cell, and the elevation model of the current flight scene is established;

[0193] The node cell is a volume pixel in three-dimensional space, each node cell represents a small area in space, and each node cell contains the height information of the location, which refers to the position or height of the ground or terrain in the vertical direction; the color coding is to encode the height information according to a gradient color band; the generated elevation map is a type of map that uses color to represent the height changes of the terrain; with the incremental construction of the point cloud, the elevation map can be dynamically updated and maintained.

[0194] It is necessary to first screen out the horizontal plane area according to the elevation model, specifically:

[0195] In the generated scene elevation model, the RANSAC algorithm is used to extract the plane structure, and the horizontal plane area is screened out based on the normal vector direction constraints of each plane;

[0196] The input of the RANSAC algorithm is the center of each node cell in the elevation map; the constraint condition is the direction of the normal vector of each plane. The relationship between each plane and the horizontal plane is measured by calculating the angle between the normal vector of each plane and the unit vector perpendicular to the horizontal plane. When the angle is close to or less than 20 degrees, the plane is identified as a horizontal plane.

[0197] The landing area selection criteria, including the landing area flatness measurement and area area, are specifically:

[0198] The flatness measurement of the horizontal plane area refers to the calculation based on the statistical standard deviation of the cells contained in the horizontal plane area, and the sum of the standard deviations of all cells is calculated. If the sum result is less than the threshold, it is considered to meet the flatness measurement requirements; the calculation formula of the flatness measurement F is:

[0199]

[0200] The area of ​​the horizontal plane region refers to the projection area of ​​the ground occupied by the cells calculated based on the projection area of ​​the cells contained in the horizontal plane region. Only the regions with an area greater than the area threshold are considered to meet the area requirement.

[0201] Among all horizontal plane areas that meet the flatness measurement requirements and area requirements, the area with the largest plane area is selected as the final UAV landing area detection result.

[0202] The acquisition of the landing area position is specifically as follows:

[0203] Update the landing area in real time as the drone flies, ensuring that the pre-selected landing area at the current moment is closest to the drone to improve landing efficiency and safety;

[0204] The landing area selection location is updated with the real-time flight of the drone. This means that when an area with the same plane area appears in the updated elevation map, the landing area selection result will be dynamically adjusted to the landing area location closest to the drone.

[0205] This scheme is explained with examples. Figure 1 The present invention discloses a real-time autonomous landing area detection method for a UAV based on monocular vision, comprising:

[0206] 1. Before flying, perform a camera function test on the drone to ensure that the camera is working properly and the image quality is good. Adjust necessary settings such as exposure, focus, etc.

[0207] Perform a low-altitude test flight before flying to observe whether the camera works normally and whether the image is clear and stable.

[0208] During the flight, the camera's operating status and image quality are continuously monitored. Adjustments are made as needed to ensure the camera is always in optimal condition.

[0209] 2. Use the calibration equipment to calibrate the camera's internal parameters. You can choose to use a checkerboard calibration plate or a solid origin calibration plate. The parameters obtained after the calibration operation include the camera's principal distance, principal point, and distortion. The camera takes images of the calibration plate at different positions and angles, detects the checkerboard corner points each time, collects 10 to 20 images of the calibration plate, and calculates the calibration results. After a calibration is completed, the calibration parameters can be used repeatedly without recalibrating each time without changing the camera structure.

[0210] 3. Obtain in real time the sequence images of the ground scene collected by the drone during flight.

[0211] 4. Tightly couple the synchronous positioning and mapping algorithm with the subsequent elevation modeling module. Using computer vision technology, feature points are extracted from the image sequence, and the correspondence between images is determined by matching feature points. This step is implemented using feature descriptors (such as SIFT, SURF, etc.) and matching algorithms. Using the matched feature points, the basic matrix or essential matrix between images is calculated to restore the epipolar geometric relationship between images and determine the relative position and attitude of the drone. Based on the position and attitude parameters of the image and the matching relationship between pixels, the three-dimensional point cloud data of the flight scene is reconstructed using the triangulation method.

[0212] 5. When the drone receives the landing command, the initial 3D point cloud is resampled using an octree-based hierarchical data structure to generate a connected voxel grid map.

[0213] 1) Performing octree partitioning on the three-dimensional space containing the three-dimensional point cloud data to obtain a hierarchical octree data structure, each node of which corresponds to a three-dimensional cube cell;

[0214] 2) According to the three-dimensional cube cell corresponding to each node of the octree data structure, the coordinate values ​​of all three-dimensional points falling into it are averaged, and the z coordinate value of the average value is used as the average height value of the node cell to obtain the resampled value of the node;

[0215] 3) Using the Gaussian probability model statistical method, the center height value of each node of the octree data structure is statistically calculated with the height values ​​of the six nodes above, below, left, right, and front and back, to obtain the statistical standard deviation value of each center node.

[0216] 6. Each node cell contains the altitude information of that location. The voxel grid map is color-coded according to the altitude information in each node cell to establish the elevation model of the current flight scene.

[0217] 7. Use the center of each node cell in the constructed elevation model as the input of the RANSAC algorithm. According to the density of the point cloud, the threshold and parameters of the RANSAC algorithm are reasonably set, and the plane structure is extracted using the RANSAC algorithm; for each fitted plane, its parameters (normal vector, intercept, etc.) are recorded; among all the extracted planes, the horizontal plane area is screened out according to the normal vector direction constraint of each plane: for each extracted plane, its angle with the horizontal plane is calculated, and the plane with an angle close to or less than 20 degrees is considered to be a horizontal plane area. See Figure 2 .

[0218] 8. Considering the safety of drone landing, formulate landing area selection criteria including landing area flatness measurement and area size.

[0219] The flatness measure is calculated based on the statistical standard deviation of the cells contained in the horizontal plane area: the standard deviation of all cells is summed up, and if the sum is less than the threshold, it is considered to meet the flatness measure requirements;

[0220] Calculate the area of ​​the region based on the projected area of ​​the cells contained in the horizontal plane area: Calculate the projected area of ​​the ground occupied by the cells, and only the regions with an area greater than the area threshold are considered to meet the region area requirements;

[0221] Among all horizontal plane areas that meet the flatness measurement requirements and area requirements, the area with the largest plane area is selected as the final UAV landing area detection result.

[0222] 9. Traverse the entire horizontal plane area and analyze various terrain factors in each area.

[0223] The above-mentioned full domain traversal of the horizontal plane area first defines a node cell, forms a set of cells within a certain distance threshold around it, creates a rectangular box with a specific edge length, and moves the rectangular box on the ground with a certain iteration step length to perform full domain traversal of the horizontal plane area. Figure 3 .

[0224] The calculation and analysis of the terrain factors of each area is to calculate the flatness measurement and area of ​​all areas in the elevation map according to the safe landing standard;

[0225] The flatness measure is obtained by calculating the statistical standard deviation of all cells contained in the rectangular frame after each movement; the area of ​​the region is obtained by calculating the projected area of ​​all cells contained in the rectangular frame.

[0226] 10. Among all the horizontal plane areas that meet the flatness measurement requirements and area requirements, select the area with the largest plane area as the final UAV landing area detection result.

[0227] 11. If Figure 4 As shown in the figure, the real-time flight trajectory of the UAV, the ground scene elevation model, and the landing area detection results are integrated and unified, and real-time monitoring and analysis are carried out through visualization software.

[0228] 12. By constantly monitoring the real-time flight status of the drone and updating the elevation map, ensure that the selected landing area is always the area that meets the safe landing standards and is closest to the drone, thereby improving the efficiency and safety of landing. Figure 5 .

[0229] So far, the purpose of the present invention has been achieved.

[0230] The above are only preferred embodiments of the present invention and are not intended to limit the present invention. Any modifications, equivalent substitutions, improvements, etc. made within the spirit and principles of the present invention should be included in the protection scope of the present invention.

Claims

1. A method for autonomous detection of blind landing areas of unmanned aerial vehicles based on monocular vision, characterized in that: The steps include: Use the calibration board to calibrate the relative position and posture of the camera and the parameters within the image to obtain the calibration parameters; The camera continuously acquires sequence images of the flight scene, covering the entire pre-landing area; Using the synchronous positioning and mapping algorithm, the relative position and real-time flight trajectory of the drone are restored through the epipolar geometric relationship between the sequence images, and the initial three-dimensional point cloud of the flight scene is reconstructed; The initial 3D point cloud is resampled using an octree-based hierarchical data structure to generate a voxel grid map with connectivity; The voxel grid map is color-coded according to the height information in each node cell to build an elevation model of the current flight scene; According to the elevation model, the landing area selection criteria including the flatness measurement and area area of ​​the landing area are formulated to obtain the nearest landing area location of the UAV.

2. The method for autonomously detecting blind landing areas of unmanned aerial vehicles based on monocular vision according to claim 1 is characterized in that: The parameter calibration method is specifically as follows: Calibration is performed using more than 12 calibration plate images taken by the camera from different angles and positions; The parameters obtained after the calibration operation include the camera's principal distance, principal point, distortion, relative position and attitude angle parameters.

3. The method for autonomously detecting blind landing areas of unmanned aerial vehicles based on monocular vision according to claim 2 is characterized in that: The calibration plate can be a checkerboard calibration plate or a solid origin calibration plate.

4. The method for autonomously detecting blind landing areas of unmanned aerial vehicles based on monocular vision according to claim 2 is characterized in that: The initial three-dimensional point cloud reconstruction method is specifically as follows: 1) Calculate the epipolar geometric relationship between the sequence images; Assume I1 and I2 are two adjacent images, p1 and p2 represent the projected pixel points of point P in the physical space on the two images, and the coordinates of p1 and p2 are obtained by transforming the coordinates of point P in the physical space in the camera coordinate system according to the camera intrinsic parameter matrix k: Among them, z1 and z2 are the depth values ​​of p1 and p2 respectively, P′ is the coordinate of the spatial point P in the I1 camera coordinate system, and R,t is the geometric transformation relationship of the spatial coordinate system between the two adjacent images of I1 and I2; Multiply the left side of the above equation by k -1 get: Substituting (3) into (4), we get: z2k -1 p2=Rz1k -1 p1+i (5) Multiply both sides of the above equation by [t] x have to: z2[t] x k-1p2=Rz1[t] x k-1 p1 (6) Where [t] x is the antisymmetric matrix of t; Multiply (6) on the left have to: in, So the left side of the equation is 0, and we can get the following equation: Therefore, the epipolar geometry constraints are obtained as follows: Where F = k -T [t] x R -1 is the basic matrix; The essential matrix E is determined by the camera's external parameters and has nothing to do with the camera's internal parameters. Let the unit-scale normalized vectors of the spatial point P in the image space coordinate system of the two images be expressed as P 0 , P 1 ; According to the spatial geometric transformation relationship between the two images: P 1 =RP 0 +t (10) Multiply both sides of the equation by the cross product of t and then by P. 1 The dot product of , we can get: 0=P 1 t×RP 0 (11) Right now: (P 1 ) T [t] x RP 0 =0 (12) (P 1 ) T EP 0 =0 (13) Where E = [t] x R represents the essential matrix, which contains the relative position and posture transformation relationship of the two images; the relationship between the essential matrix E and the basic matrix F is: E=k T F k (14) 2) Further calculate the relative parameters of the image position and attitude of the UAV during flight; First, initialize the image position and scene structure: E is decomposed into an antisymmetric matrix S (equivalent to [t] x ) and the product of a rotation matrix R: E=SR (15) If S is a real antisymmetric matrix, then S = UBU T , where B is a block diagonal matrix; two temporary matrices are introduced: and In the case of indefinite sign, Z = diag(1, 1, 0)W, and in the case of indefinite scale, S = Udiag(1, 1, , 0)WU T ; Perform singular value decomposition on E: E=SR=Udiag(1,1,0)(WU T R) (16) Then, the camera matrix is ​​extracted from the essential matrix to determine the relative parameters of the image position and posture: E=U∑V T ,∑=diag(σ1,σ2,σ3) (17) Where o1=o2≠C and σ3=0, U, V∈SO(3); In the case of a constant factor difference, S = Udiag (1, 1, 0) U T , assuming that the rotation matrix R can be decomposed into R = UXV T , where X is a rotation matrix, then: Udiag(1,1,0)V T =B=SR=(UZU T )(UXV T )=U(ZX)V T (18) Then we have ZX = diag(1, 1, 0), and because X is a rotation matrix, we have X = W or X = W T ; If the SVD decomposition of E is Udiag(1,1,0)V T : At this point there are 4 possible camera matrices, namely: Only one solution can satisfy the requirement that the reconstructed point is in front of the camera, that is, there is only one reliable solution. The relative parameters of the image position and posture are determined by testing whether a point is in front of the camera. 3) Using the position and attitude parameters of the image and the matching relationship between pixels, the three-dimensional point cloud data of the flight scene is reconstructed using the triangulation method; According to epipolar geometry, we have: s2p2=s1Rp1+t (21) Among them, p1 and p2 are the normalized coordinates of the feature points on the two images, and s1 and s2 are the depths of the two feature points. Multiply both sides of the above formula by According to the above formula, s1s2 is obtained; due to the influence of noise, the obtained R, t may not accurately make the above formula equal to 0, so it is actually solved by least squares; The result of triangulation is the real-time calculation of the coordinates of the three-dimensional points in the physical world corresponding to the image matching pixels, thereby reconstructing the three-dimensional point cloud data of the flight scene.

5. The method for autonomously detecting blind landing areas of unmanned aerial vehicles based on monocular vision according to claim 4 is characterized in that: The octree-based hierarchical data structure resamples the three-dimensional point cloud data, comprising the steps of: 1) Performing octree partitioning on the three-dimensional space containing the three-dimensional point cloud data to obtain a hierarchical octree data structure, each node of which corresponds to a three-dimensional cube cell; First, determine the bounding box in three-dimensional space, which is the smallest cube that contains all the point cloud data; Secondly, divide the root node into eight child nodes, each of which corresponds to one eighth of the volume of the root node; for each child node, check whether it contains enough points; if the number of points contained exceeds the threshold, further divide the child node into eight smaller child nodes; repeat this process until the number of points in the node is less than or equal to the threshold; Finally, the point cloud data is assigned to the corresponding octree nodes. After recursive partitioning, the points will eventually be distributed in the leaf nodes of the octree. 2) According to the three-dimensional cube cell corresponding to each node of the octree data structure, the average value of the z coordinates of the n three-dimensional points falling into it is used as the center height value of the node cell; in, is the center height value of the j-th node cell; 3) Using the Gaussian probability model statistical method, the center height value of each node of the octree data structure is statistically calculated with the height values ​​of the six nodes above, below, left, right, and front and back, to obtain the statistical standard deviation value of each center node; the formula for calculating the statistical standard deviation value of the z coordinates of the i-node cells is: in, is the center height value of the central node, and σ is the statistical standard deviation value of the node.

6. The method for autonomously detecting blind landing areas of unmanned aerial vehicles based on monocular vision according to claim 5 is characterized in that: The elevation model is established specifically as follows: each node cell contains the altitude information of the location, the voxel grid map is color-coded according to the altitude information in each node cell, and the elevation model of the current flight scene is established; The node cell is a volume pixel in three-dimensional space, each node cell represents a small area in space, and each node cell contains the height information of the location, which refers to the position or height of the ground or terrain in the vertical direction; the color coding is to encode the height information according to a gradient color band; the generated elevation map is a type of map that uses color to represent the height changes of the terrain; with the incremental construction of the point cloud, the elevation map can be dynamically updated and maintained.

7. The method for autonomously detecting blind landing zones of unmanned aerial vehicles based on monocular vision according to claim 6 is characterized in that: It is necessary to first screen out the horizontal plane area according to the elevation model, specifically: In the generated scene elevation model, the RANSAC algorithm is used to extract the plane structure, and the horizontal plane area is screened out based on the normal vector direction constraints of each plane; The input of the RANSAC algorithm is the center of each node cell in the elevation map; the constraint condition is the direction of the normal vector of each plane. The relationship between each plane and the horizontal plane is measured by calculating the angle between the normal vector of each plane and the unit vector perpendicular to the horizontal plane. When the angle is close to or less than 20 degrees, the plane is identified as a horizontal plane.

8. The method for autonomously detecting blind landing areas of unmanned aerial vehicles based on monocular vision according to claim 7 is characterized in that: The landing area selection criteria, including the landing area flatness measurement and area area, are specifically: The flatness measurement of the horizontal plane area refers to the calculation based on the statistical standard deviation of the cells contained in the horizontal plane area, and the sum of the standard deviations of all cells is calculated. If the sum result is less than the threshold, it is considered to meet the flatness measurement requirements; the calculation formula of the flatness measurement F is: The area of ​​the horizontal plane region refers to the projection area of ​​the ground occupied by the cells calculated based on the projection area of ​​the cells contained in the horizontal plane region. Only the regions with an area greater than the area threshold are considered to meet the area requirement. Among all horizontal plane areas that meet the flatness measurement requirements and area requirements, the area with the largest plane area is selected as the final UAV landing area detection result.

9. The method for autonomous detection of blind landing areas of unmanned aerial vehicles based on monocular vision according to claim 8, characterized in that: The acquisition of the landing area position is specifically as follows: Update the landing area in real time as the drone flies, ensuring that the pre-selected landing area is closest to the drone at the current moment, so as to improve the efficiency and safety of landing; The landing area selection location is updated with the real-time flight of the drone. This means that when an area with the same plane area appears in the updated elevation map, the landing area selection result will be dynamically adjusted to the landing area location closest to the drone.

Images