CN108876861B

CN108876861B - Stereo matching method for extraterrestrial celestial body patrolling device

Info

Publication number: CN108876861B
Application number: CN201810513673.5A
Authority: CN
Inventors: 李海超; 李莹; 陈亮; 顾征; 李峰; 辛蕾
Original assignee: China Academy of Space Technology CAST
Current assignee: China Academy of Space Technology CAST
Priority date: 2018-05-25
Filing date: 2018-05-25
Publication date: 2020-10-20
Anticipated expiration: 2038-05-25
Also published as: CN108876861A

Abstract

The invention relates to a stereo matching method of an extraterrestrial celestial body patrolling device, belonging to the stereo vision matching technology. The method comprises the following steps: acquiring left and right images by left and right cameras of the extraterrestrial celestial body patrolling device; performing epipolar line correction on the left and right images, and respectively constructing left and right pyramids of L layers for the left and right images after the epipolar line correction; performing bidirectional stereo matching on the two images of the L-th layer of the left pyramid and the right pyramid to obtain left and right disparity maps; determining the parallax search range of each pixel of the L-1 layer by using the parallax map of the L-1 layer, and performing bidirectional stereo matching on the two images of the L-1 layer to obtain left and right parallax maps of the L-1 layer; and then, circulating until the parallax search range of each pixel of the 1 st layer is determined by using the parallax map of the 2 nd layer, and performing bidirectional stereo matching on the two images of the 1 st layer to obtain the left and right parallax maps of the 1 st layer.

Description

Stereo matching method for extraterrestrial celestial body patrolling device

Technical Field

The invention relates to a stereo matching method of an extraterrestrial celestial object inspection device, in particular to a multi-resolution fast stereo matching method of the extraterrestrial celestial object inspection device, and belongs to the technical field of computer vision and deep space exploration.

Background

Extraterrestrial celestial body patroller exploration has become an important element of deep space exploration. Since extraterrestrial celestial bodies have the characteristics of unknown working environment, unstructured working environment, complex terrain, large working range, remote distance from the earth and the like, a series of problems of autonomous environment perception, path planning and the like of the patrol device need to be mainly solved. The rapid three-dimensional reconstruction of the terrain on the surface of the extraterrestrial celestial body is a key problem for ensuring the patrol device to smoothly advance in the extraterrestrial body.

At present, the extraterrestrial celestial body patrol device mainly realizes three-dimensional reconstruction of a terrain by a laser radar method and a stereoscopic vision method. The laser radar method can generate high-precision depth information, and has the advantages of small calculation amount compared with stereoscopic vision, but has the defects of large power consumption, large volume and weight and the like of equipment. For example, laser striping machines were used for obstacle detection on the united states "songera" train landed in 7 months 1997 due to the processing speed limitations of the processor for radiation hardening at that time. However, it has been estimated that tens of thousands of dollar costs are incurred for every 1 kg of weight added to space. The stereoscopic vision system has the characteristics of large visual field, low power consumption, small volume and mass, low cost and the like, and meets the emission and use requirements of the carried equipment with the mass, the energy consumption and the volume as small as possible. Particularly, with the improvement of the processing speed of the processor and the progress of the stereoscopic vision technology, the stereoscopic vision has wide application prospect in the field of extraterrestrial celestial body detection. For example, U.S. Mars, courage and opportunity landings in 2004, utilize stereo vision to render three-dimensional maps; the curio Mars train landing Mars in 8 months of 2012 also adopts stereoscopic vision as the main technical means of obstacle avoidance, path planning and navigation positioning; the three-dimensional reconstruction of the unknown environment on the lunar surface is realized by adopting a stereoscopic vision technology, and the autonomous obstacle avoidance based on the stereoscopic vision is realized.

The method based on stereoscopic vision is to perform stereoscopic matching on two or more images obtained by a binocular stereoscopic camera to obtain three-dimensional information of ground terrain. Scharstein et al mainly divide the stereo matching method into four steps (see D.Scharstein, R.Szeliski. A. taxomony and evaluation of dense two-frame stereo response algorithm. IJCV,47(1-3):7-42,2002): (1) calculating pixel point matching cost; (2) aggregation of matching costs; (3) disparity calculation and optimization (sub-pixel); (4) and (4) performing parallax refinement, and dividing the matching method into a local matching method and a global optimization method according to an optimization theory. The local matching method calculates the matching cost by using a matching window around the pixel, and finds a minimum matching cost value as a parallax value through the optimization idea of 'Winner Take All-WTA'. The method mainly considers problems from local parts, so that mismatching is easy to generate, but the method has the advantage in calculation speed and can adapt to occasions with high real-time requirements. The global optimization method is mainly based on the energy minimization idea, and the globally optimal parallax image is obtained. Generally, the global optimization problem can be expressed by a Markov Random Field (MRF), but the MRF model Field is a two-dimensional image, so that the MRF-based energy function global optimization is an NP problem, and although some approximate optimization methods such as a graph cut algorithm and a confidence propagation algorithm exist, the computation efficiency is still low, and the method is not suitable for binocular stereo matching of the extraterrestrial celestial body rover. The dynamic programming is to optimize a single scanning line, and the dynamic programming algorithm is an algorithm with higher real-time performance in the global algorithm. But due to the one-dimensional nature of the dynamic programming algorithm, there are tomographic phenomena between different scan lines. Therefore, to overcome the above problems, the SGM (Semi-global) algorithm generalizes the dynamic planning of a single scan line to multiple directions, thereby expressing approximately global relationships and introducing disparity continuity constraints in the dynamic planning (see Hirschmuller h., Stereo processing by Semi-global planning and physical information, ieee transformations and Pattern analysis and Machine Intelligence,30(2): 328-. To overcome the possible banding in the weak texture region by the SGM algorithm, facciio et al improved the SGM algorithm and proposed the MGM (More global matching) algorithm (see facciio g., Franchis c., MGM: a signalizing More global matching for Stereo Vision. british Machine Vision Conference, September 2015.).

Although the SGM and MGM algorithms achieve better effects, on one hand, in the extraterrestrial celestial environment of unknown environment, since the parallax search range of the scene cannot be predicted, if the conventional SGM or MGM algorithm is adopted, a relatively large parallax search range is used for all pixels of the whole image. However, if the parallax search range is set to be too large, more calculation time is consumed, more memory space is occupied, and meanwhile, false matching may be caused; if the parallax search range is set too small, the correct parallax is not obtained. Meanwhile, the imaging resolution of the extraterrestrial celestial body patrolling device is higher and higher, and huge calculation amount is brought. On the other hand, both the SGM and MGM processing algorithms use the same set of penalty parameters P in the entire image₁And P₂Resulting in poor matching accuracy in the parallax discontinuous area and the occlusion area.

In order to solve the problems of the two aspects, on one hand, the invention establishes a multi-resolution fast stereo matching strategy of the extraterrestrial celestial body patroller, and the generated low-resolution matching result can be propagated step by step to restrict the parallax searching range of a high-resolution image, so that the acceleration essence is to reduce the parallax searching range and simultaneously reduce the occurrence probability of mismatching; on the other hand, in order to overcome the defect that the traditional SGM or MGM algorithm has poor matching precision in a parallax discontinuous area and an occlusion area, a strategy for adaptively adjusting penalty parameters according to edge information is provided.

Disclosure of Invention

The technical problem to be solved by the invention is as follows: the method is small in parallax search range and high in calculation speed, and parallax discontinuous areas and occlusion areas have good matching accuracy.

The technical solution of the invention is as follows:

a stereo matching method for an extraterrestrial celestial body patroller comprises the following steps:

s1, acquiring a stereoscopic image pair of the extraterrestrial celestial body scene by using a left camera and a right camera of the extraterrestrial celestial body patrolling device, wherein the image acquired by the left camera is defined as a left image, and the image acquired by the right camera is defined as a right image;

s2, performing epipolar line correction on the left image and the right image in the step S1 respectively, constructing an L-layer pyramid, namely a left pyramid, on the left image subjected to epipolar line correction, and constructing an L-layer pyramid, namely a right pyramid, on the right image subjected to epipolar line correction;

l is a natural number, and is preferably 2-5;

the pyramid is preferably a Gaussian pyramid;

s3, performing bidirectional stereo matching on the image of the left pyramid highest layer, namely the L-th layer, and the image of the right pyramid highest layer, namely the L-th layer, which are constructed in the step S2 to obtain a disparity map of the L-th layer of the left pyramid and a disparity map of the L-th layer of the right pyramid;

the method for performing bidirectional stereo matching on the image at the highest layer, i.e., the L-th layer, of the left pyramid and the image at the highest layer, i.e., the L-th layer, of the right pyramid, which are constructed in the step S2, specifically includes the following steps:

s31: extracting SIFT feature points in the highest-level image of the left pyramid, and the extraction is called as leftSIFT feature points, namely extracting SIFT feature points in the highest-layer image of the right pyramid, namely right SIFT feature points, and matching the extracted left SIFT feature points with the right SIFT feature points to obtain n pairs of matching points; the coordinate of the characteristic point corresponding to the ith matching point in the highest layer image of the left pyramid in the n pairs of matching points is recorded as (x)_i,y_i) And the coordinates of the feature point corresponding to the ith pair of matching points in the n pairs of matching points in the image of the highest layer of the right pyramid are recorded as (x)_i',y_i'), i is 0,1, …, n-1, n is more than or equal to 3; the SIFT (Scale Invariant feature transform, i.e., Scale Invariant feature transform);

s32: d is calculated from the i-th pair of matching points among the n pairs of matching points obtained in step S31_i＝x_i'-x_iObtaining a spatial point corresponding to the ith characteristic point in the highest-layer image of the left pyramid, and recording the coordinate of the spatial point as (x)_i,y_i,d_i) Obtaining n space points by the n feature points in the left pyramid highest-layer image, performing space plane fitting on the n space points obtained by the n feature points in the left pyramid highest-layer image, and calculating to obtain the parallax search range of all pixel points in the left pyramid highest-layer image;

d is calculated from the i-th pair of matching points among the n pairs of matching points obtained in step S31_i'＝x_i-x_i', obtaining a spatial point corresponding to the ith characteristic point in the highest-layer image of the right pyramid, and recording the coordinate of the spatial point as (x)_i',y_i',d_i'), obtaining n space points by n feature points in the image of the highest layer of the right pyramid, performing space plane fitting on the n space points obtained by the n feature points in the image of the highest layer of the right pyramid, and calculating to obtain the parallax search range of all pixel points in the image of the highest layer of the right pyramid;

the method for performing spatial plane fitting on n spatial points obtained from n feature points in the highest-layer image of the left pyramid and calculating to obtain the parallax search range with pixel points in the highest-layer image of the left pyramid comprises the following steps:

s321, performing space plane fitting on n space points obtained from the n feature points in the highest-level image of the left pyramid, and expressing an equation of an obtained space plane Π as D ═ aX + bY + c;

wherein a, b and c represent coefficients in a space plane equation; (X, Y, D) is the coordinate of any point on the space plane pi, and D is parallax;

fitting the n spatial points using a least squares method to determine the coefficients in the spatial plane equation, i.e.

According to

Obtaining space plane coefficients a, b and c;

s322, regarding the coordinate (x) of the ith feature point in the n feature points in the highest-level image of the left pyramid_i,y_i) Calculating the parallax D corresponding to the ith characteristic point according to the equation D ═ aX + bY + c of the space plane pi_i ^ΠCalculating the parallax difference d of the ith feature point_i ^Π-d_i(ii) a For n parallax difference values corresponding to n feature points in the highest-layer image of the left pyramid, counting the maximum value D in the n parallax difference values_maxAnd the minimum value D_min；

S323, for any pixel point p (x, y) in the highest-level image of the left pyramid, calculating a parallax D corresponding to the pixel point p (x, y) according to an equation D ═ aX + bY + c of the spatial plane Π^Π(x, y), calculating d_min(x,y)＝d^Π(x,y)+D_min-Δd，d_max(x,y)＝d^Π(x,y)+D_max+ Δ d, where Δ d is a reserved parallax margin, preferably 2-5, and the parallax search range of the pixel point p (x, y) is determined as [ d [_min(x,y),d_max(x,y)]Obtaining a parallax searching range with pixel points in the highest layer image of the left pyramid by using the integer in the left pyramid;

and (3) performing space plane fitting on n space points obtained from the n feature points in the highest-layer image of the right pyramid, and repeating the steps S321-S323 according to the same principle to obtain the parallax search range of all pixel points in the highest-layer image of the right pyramid.

S33: performing Census transformation on all pixel points in the left pyramid highest-layer image to obtain Census sequences of all pixel points in the left pyramid highest-layer image, and performing Census transformation on all pixel points in the right pyramid highest-layer image to obtain Census sequences of all pixel points in the right pyramid highest-layer image; calculating to obtain the matching cost of all the pixel points in the left pyramid highest-layer image within the parallax search range of all the pixel points in the left pyramid highest-layer image according to the obtained Census sequences of all the pixel points in the left pyramid highest-layer image and the Census sequences of all the pixel points in the right pyramid highest-layer image; calculating the matching cost of all pixel points in the right pyramid highest-layer image within the parallax search range of all pixel points in the right pyramid highest-layer image according to the obtained Census sequences of all pixel points in the left pyramid highest-layer image and Census sequences of all pixel points in the right pyramid highest-layer image;

the Census conversion method for all the pixel points in the left pyramid highest-level image to obtain the Census sequences of all the pixel points in the left pyramid highest-level image comprises the following steps:

s331, establishing a Census window by taking any pixel point in the highest-layer image of the left pyramid as a center, comparing the gray value magnitude relation between pixel points except the center point and the gray value of the center point in the established Census window, marking the pixel point of which the gray value is less than or equal to the gray value of the center point as 0, otherwise marking the pixel point as 1, and finally connecting the comparison result into 0/1 binary code streams according to bits, wherein the Census sequences are called as current pixel points;

s332, repeating the step S331 until all pixel points in the highest-layer image of the left pyramid are processed, and obtaining Census sequences of all pixel points in the highest-layer image of the left pyramid;

and repeating the steps S331-S332 according to the same principle by the method for performing Census conversion on all the pixel points in the right pyramid highest-level image to obtain the Census sequences of all the pixel points in the right pyramid highest-level image.

The method for calculating the matching cost of all the pixel points in the highest-layer image of the left pyramid in the parallax search range of all the pixel points in the highest-layer image of the left pyramid according to the obtained Census sequences of all the pixel points in the highest-layer image of the left pyramid and the Census sequences of all the pixel points in the highest-layer image of the right pyramid comprises the following steps:

s333, calculating the Hamming distance of the Census sequence of the pixel point p (x, y) in the left pyramid highest-level image and the pixel point q (x-d, y) in the right pyramid highest-level image, wherein d is parallax, and d belongs to [ d [ d ]_min(x,y),d_max(x,y)]And d is an integer, and obtaining a matching cost C (p, d) of a pixel point p (x, y) in the highest-layer image of the left pyramid when the parallax is d;

s334, repeating the step S333 to obtain the matching cost of the pixel point p (x, y) in the highest-layer image of the left pyramid in different parallaxes d;

s335, repeating the steps S333-S334 until all pixel points in the highest-layer image of the left pyramid are processed, and obtaining the matching cost of all pixel points in the highest-layer image of the left pyramid;

and repeating the steps S333-S335 according to the same principle by the method for calculating the matching cost of all the pixel points in the left pyramid highest-level image according to the parallax search range of all the pixel points in the left pyramid highest-level image, so as to obtain the matching cost of all the pixel points in the right pyramid highest-level image.

S34: performing matching cost aggregation on all pixel points in the left pyramid highest-layer image according to the matching cost of all pixel points in the left pyramid highest-layer image obtained in the step S33 to obtain a matching cost aggregation value of the left pyramid highest-layer image, and obtaining an initial disparity map of the left pyramid highest-layer image according to the obtained matching cost aggregation value of the left pyramid highest-layer image; performing matching cost aggregation on all pixel points in the right pyramid highest-layer image according to the matching costs of all pixel points in the right pyramid highest-layer image obtained in the step S33 to obtain a matching cost aggregation value of the right pyramid highest-layer image, and obtaining an initial disparity map of the right pyramid highest-layer image according to the obtained matching cost aggregation value of the right pyramid highest-layer image;

the method for aggregating the matching costs of all the pixel points in the left pyramid highest-layer image according to the matching costs of all the pixel points in the left pyramid highest-layer image obtained in step S33 to obtain a matching cost aggregate value, and obtaining the initial disparity map of the left pyramid highest-layer image according to the obtained matching cost aggregate value includes:

s341, Canny edge feature detection is carried out on the left pyramid highest-layer image to obtain a plurality of edges; the points on the edge are called edge points, and for any edge point, the range [ d ] is searched according to the parallax of the pixel point corresponding to the edge point_min(x,y),d_max(x,y)]Calculating a difference d_max(x,y)-d_min(x, y), called the parallax span of the edge point, calculating the average value of the parallax spans of all the edge points on each edge, if the average value is less than the set threshold, removing the current edge, if the average value is not less than the set threshold, keeping the current edge; processing a plurality of edges obtained by detecting the left pyramid highest-layer image to obtain qualified edges in the left pyramid highest-layer image; the set threshold value is preferably a real number in a range of 2-8;

s342, performing matching cost aggregation on all pixel points in the left pyramid highest-layer image according to the matching costs of all pixel points in the left pyramid highest-layer image to obtain a matching cost aggregation value;

calculating a matching cost aggregation value S (p, d) of the pixel point p (x, y) when the parallax is d:

in the formula, r represents the propagation direction of a pixel point p (x, y); l is_r(p, d) represents the matching cost of the pixel point p (x, y) along the direction r;

wherein R represents m directions derived from the propagation direction R, and R ∈ { R ∈ {₁,…,r_m}，m≥1，r₁In the same direction as r, i.e. r₁R; p-R represents the previous pixel point of the pixel point p along the direction R; l is_r(p-R, d) represents the matching cost of the pixel p-R when the disparity is d, L_r(p-R, d +1) represents the matching cost of the pixel p-R when the disparity is d +1, L_r(p-R, d-1) represents the matching cost of the pixel point p-R when the parallax is d-1; d 'is the disparity of pixel p-R, whose disparity search range d' ∈ [ d ∈ ]/_min(p-R),d_max(p-R)]，

Indicating that pixel p-R is in the disparity search range d_min(p-R),d_max(p-R)]The inner minimum path matching cost; p₁ ^R、P₂ ^RDenotes a penalty parameter, P₂ ^R>P₁ ^R；

The propagation direction r is 8 directions,

r∈{(-1,0),(1,0),(0,1),(0,-1),(-1,-1),(1,-1),(1,1),(-1,1)}，

or in the 4 directions, the direction of the light beam,

r∈{(-1,0),(1,0),(0,1),(0,-1)}；

m is 4, 2 or 1;

when m is 4 and the propagation directions R are 8, R derived from the 8 propagation directions R is sequentially represented as:

R∈{(-1,0),(0,-1),(-1,-1),(1,-1)}、R∈{(1,0),(0,1),(1,1),(-1,1)}、

R∈{(0,1),(-1,0),(-1,1),(-1,-1)}、R∈{(0,-1),(1,0),(1,-1),(1,1)}

R∈{(-1,-1),(1,-1),(0,-1),(1,0)}、R∈{(1,-1),(1,1),(1,0),(0,1)}、

R∈{(1,1),(-1,1),(0,1),(-1,0)}、R∈{(-1,1),(-1,-1),(-1,0),(0,-1)}

when m is 4 and the propagation directions R are 4, R derived from the 4 propagation directions R is sequentially represented as:

R∈{(-1,0),(0,-1),(-1,-1),(1,-1)}、R∈{(1,0),(0,1),(1,1),(-1,1)}、

R∈{(0,1),(-1,0),(-1,1),(-1,-1)}、R∈{(0,-1),(1,0),(1,-1),(1,1)}

when m is 2 and the propagation directions R are 8, R derived from the 8 propagation directions is sequentially represented as:

R∈{(-1,0),(0,-1)}、R∈{(1,0),(0,1)}、R∈{(0,1),(-1,0)}、R∈{(0,-1),(1,0)}、

R∈{(-1,-1),(1,-1)}、R∈{(1,-1),(1,1)}、R∈{(1,1),(-1,1)}、R∈{(-1,1),(-1,-1)}

when m is 2 and the propagation directions R are 4, R derived from the 4 propagation directions is sequentially represented as: r belongs to { (-1,0), (0, -1) }, R belongs to { (1,0), (0,1) }, R belongs to { (0,1), (-1,0) }, R belongs to { (0, -1), (1,0) };

when m is 1 and the propagation directions R are 8, R derived from the 8 propagation directions is sequentially represented as:

R∈{(-1,0)}、R∈{(1,0)}、R∈{(0,1)}、R∈{(0,-1)}、

R∈{(-1,-1)}、R∈{(1,-1)}、R∈{(1,1)}、R∈{(1,1)}

when m is 1 and the propagation directions R are 4, R derived from the 4 propagation directions is sequentially represented as: r belongs to { (-1,0) }, R belongs to { (0,1) }, R belongs to { (0, -1) };

if the current pixel point p (x, y) is an edge point on the edge meeting the conditions in the highest image of the left pyramid, the punishment parameter corresponding to p (x, y) is as follows:

if the current pixel point p (x, y) is not the edge point on the edge meeting the conditions in the highest-layer image of the left pyramid, the punishment parameter corresponding to p (x, y) is as follows:

wherein, P₁Preferably 8, P₂Preferably 32;

processing other pixel points in the left pyramid highest-layer image according to the step S342 to obtain matching cost aggregate values of all the pixel points in the left pyramid highest-layer image;

s343, calculating an initial disparity map of the highest-layer image of the left pyramid: calculating the parallax value at the point p (x, y) according to the winner-king strategy:

processing other pixel points in the left pyramid highest-layer image according to the step S343 to obtain parallax values of all the pixel points in the left pyramid highest-layer image, and forming an initial parallax map of the left pyramid highest-layer image;

and repeating the steps S341 to S343 according to the same principle by the method for obtaining the right initial parallax map according to the obtained matching cost aggregation value to obtain the initial parallax map of the right pyramid highest layer image.

S35: performing sub-pixel parallax calculation on the initial parallax map of the left pyramid highest-layer image obtained in the step S34 to obtain a sub-pixel parallax map of the left pyramid highest-layer image; performing sub-pixel parallax calculation on the initial parallax map of the right pyramid highest-layer image obtained in the step S34 to obtain a sub-pixel parallax map of the right pyramid highest-layer image; and then carrying out left-right consistency detection processing on the obtained sub-pixel disparity map of the left pyramid highest-layer image and the obtained sub-pixel disparity map of the right pyramid highest-layer image to obtain a disparity map of the left pyramid highest-layer image and a disparity map of the right pyramid highest-layer image.

The sub-pixel parallax calculation adopts a quadratic curve fitting method;

the left-right consistency detection processing method comprises the following steps:

setting upThe sub-pixel parallax map of the highest-level image of the left pyramid is marked as D_lr(x, y), the sub-pixel parallax map of the highest layer image of the right pyramid is marked as D_rl(x, y), if any pixel point in the sub-pixel disparity map of the highest-layer image of the left pyramid satisfies | D_lr(x,y)-D_rl(x+D_lr(x,y),y)|<And T, if T is 1, the pixel point is a left and right consistency point, otherwise, the pixel point is a left and right inconsistency point.

S4, determining the parallax search range of all pixel points in the L-1 layer image of the left pyramid by using the parallax map of the L layer of the left pyramid obtained in the step S3; determining the parallax search range of all pixel points in the L-1 layer image of the right pyramid by using the parallax map of the L layer of the right pyramid obtained in the step S3; performing left-to-right stereo matching on the L-1 layer image of the left pyramid and the L-1 layer image of the right pyramid according to the determined parallax search ranges of all pixel points in the L-1 layer image of the left pyramid, and performing right-to-left stereo matching on the L-1 layer image of the left pyramid and the L-1 layer image of the right pyramid according to the determined parallax search ranges of all pixel points in the L-1 layer image of the right pyramid to obtain a parallax image of the L-1 layer of the left pyramid and a parallax image of the L-1 layer of the right pyramid;

the method for determining the parallax search range of all pixel points in the L-1 layer image of the left pyramid by using the parallax map of the L layer of the left pyramid obtained in the step S3 includes:

s41: to any pixel point p in the L-1 layer image of the left pyramid_L-1(x_L-1,y_L-1) Calculating

According to the parallax map of the L-th layer of the left pyramid obtained in the step S3, pixel points p in the L-th layer image of the left pyramid are processed_L(x_L,y_L) The corresponding parallax is denoted as d_L(x_L,y_L) With p_L(x_L,y_L) Establishing a window of 3 x 3 pixels size for the center;

s42: if all the pixels in the window are left-right consistency points, calculating the corresponding parallax of all the pixels in the windowMaximum value d_{L max}And the minimum value d of parallax_{L min}：

d_{L max}＝max{d_L(x_L+i,y_L+j)|i＝-1,0,1,j＝-1,0,1}

d_{L min}＝min{d_L(x_L+i,y_L+j)|i＝-1,0,1,j＝-1,0,1}

Pixel point p in L-1 layer image of left pyramid_L-1(x_L-1,y_L-1) Is set to [2d ]_{L min}-,2d_{L max}+]The integer in (1) is the rest, and is preferably 1-3;

s43: if there is a left-right inconsistency point for the pixels in the window, then at y_LLine, with p_L(x_L,y_L) As a starting point, searching leftwards until a left-right consistency point p is obtained_l(x_L-x_l,y_L) Searching right until left and right consistency points p are obtained_r(x_L+x_r,y_L) (ii) a At the y-th_L-1 line with p_L(x_L,y_L-1) as a starting point, searching to the left until a left and right consistency point p is obtained_ul(x_L-x_ul,y_L-1), search right until left and right consistency points p are obtained_ur(x_L+x_ur,y_L-1); at the y-th_L+1 lines, with p_d(x_L,y_L+1) as starting point, searching left until obtaining left and right consistency point p_dl(x_L-x_dl,y_L+1), searching right until left and right consistency points p are obtained_dr(x_L+x_dr,y_L+1)；

Obtaining the point p according to the disparity map of the L-th layer of the left pyramid obtained in the step S3_ul,p_ur,p_l,p_r,p_dl,p_drMaximum value d 'of corresponding parallax'_LmaxAnd minimum value of d'_LminThe pixel point p in the L-1 layer image of the left pyramid is processed_L-1(x_L-1,y_L-1) Is set to

The whole number of (1);

s44: repeating the steps of S41-S43 until all the pixel points in the L-1 layer image of the left pyramid are processed, and obtaining the parallax search range of all the pixel points in the L-1 layer image of the left pyramid;

and (5) repeating the steps S41-S44 according to the same principle by using the method for determining the parallax search range of all the pixel points in the L-1 layer image of the right pyramid, which is obtained in the step S3, by using the parallax map of the L-1 layer of the right pyramid, so as to obtain the parallax search range of all the pixel points in the L-1 layer image of the right pyramid.

The method comprises the following steps of performing left-to-right stereo matching on the L-1 layer image of the left pyramid and the L-1 layer image of the right pyramid according to the determined parallax search ranges of all pixel points in the L-1 layer image of the left pyramid, and performing right-to-left stereo matching on the L-1 layer image of the left pyramid and the L-1 layer image of the right pyramid according to the determined parallax search ranges of all pixel points in the L-1 layer image of the right pyramid to obtain a left parallax image of the L-1 layer of the left pyramid and a right parallax image of the L-1 layer of the right pyramid, wherein the specific steps are as follows:

s45: performing Census conversion on all pixel points in the L-1 layer image of the left pyramid to obtain Census sequences of all pixel points in the L-1 layer image of the left pyramid, and performing Census conversion on all pixel points in the L-1 layer image of the right pyramid to obtain Census sequences of all pixel points in the L-1 layer image of the right pyramid; calculating the matching cost of all pixel points in the L-1 layer image of the left pyramid in the parallax search range of all pixel points in the L-1 layer image of the left pyramid according to the obtained Census sequences of all pixel points in the L-1 layer image of the left pyramid and the Census sequences of all pixel points in the L-1 layer image of the right pyramid; calculating the matching cost of all pixel points in the L-1 layer image of the right pyramid according to the obtained Census sequences of all pixel points in the L-1 layer image of the left pyramid and the Census sequences of all pixel points in the L-1 layer image of the right pyramid and in the parallax searching range of all pixel points in the L-1 layer image of the right pyramid;

s46: performing matching cost aggregation on all pixel points in the L-1 layer image of the left pyramid according to the matching cost of all the pixel points in the L-1 layer image of the left pyramid obtained in the step S45 to obtain a matching cost aggregation value of the L-1 layer image of the left pyramid, and obtaining an initial parallax map of the L-1 layer image of the left pyramid according to the obtained matching cost aggregation value of the L-1 layer image of the left pyramid; performing matching cost aggregation on all pixel points in the L-1 layer image of the right pyramid according to the matching cost of all pixel points in the L-1 layer image of the right pyramid obtained in the step S45 to obtain a matching cost aggregation value of the L-1 layer image of the right pyramid, and obtaining an initial parallax map of the L-1 layer image of the right pyramid according to the obtained matching cost aggregation value of the L-1 layer image of the right pyramid;

s47: performing sub-pixel parallax calculation on the initial parallax map of the L-1 layer image of the left pyramid obtained in the step S46 to obtain a sub-pixel parallax map of the L-1 layer image of the left pyramid; performing sub-pixel parallax calculation on the initial parallax map of the L-1 layer image of the right pyramid obtained in the step S46 to obtain a sub-pixel parallax map of the L-1 layer image of the right pyramid; and then carrying out left-right consistency detection processing on the obtained sub-pixel disparity map of the L-1 layer image of the left pyramid and the sub-pixel disparity map of the L-1 layer image of the right pyramid to obtain a disparity map of the L-1 layer image of the left pyramid and a disparity map of the L-1 layer image of the right pyramid.

S5, determining the parallax search range of all pixel points in the L-2 layer image of the left pyramid by using the parallax map of the L-1 layer of the left pyramid obtained in the step S4; determining the parallax search range of all pixel points in the L-2 layer image of the right pyramid by using the parallax map of the L-1 layer of the right pyramid obtained in the step S4; performing left-to-right stereo matching on the L-2 layer image of the left pyramid and the L-2 layer image of the right pyramid according to the determined parallax search ranges of all pixel points in the L-2 layer image of the left pyramid, and performing right-to-left stereo matching on the L-2 layer image of the left pyramid and the L-2 layer image of the right pyramid according to the determined parallax search ranges of all pixel points in the L-2 layer image of the right pyramid to obtain a parallax image of the L-2 layer of the left pyramid and a parallax image of the L-2 layer of the right pyramid;

……

s6, determining the parallax search range of all pixel points in the left pyramid layer 1 image (namely the left image after epipolar line correction) by using the parallax map of the left pyramid layer 2; determining the parallax search range of all pixel points in the image of the 1 st layer of the right pyramid (namely the right image after epipolar line correction) by using the parallax image of the 2 nd layer of the right pyramid; and performing stereo matching from left to right on the left pyramid layer-1 image and the right pyramid layer-1 image according to the determined parallax search ranges of all pixel points in the left pyramid layer-1 image, and performing stereo matching from right to left on the left pyramid layer-1 image and the right pyramid layer-1 image according to the determined parallax search ranges of all pixel points in the right pyramid layer-1 image to obtain a parallax image of the left pyramid layer-1 and a parallax image of the right pyramid layer-1.

Compared with the prior art, the invention has the beneficial effects that:

(1) according to the characteristic that the parallax is in inverse proportion to the distance of an observer, namely, the parallax value of a point close to the inspection device in the left image and the right image obtained by the inspection device is large, and the parallax value of a point far away from the inspection device is relatively small, the parallax searching range is estimated by using a space plane fitting method, the parallax searching range is greatly reduced, the parallax searching range does not need to be manually specified, and the method is suitable for being applied to extraterrestrial bodies with unknown environments.

(2) The invention provides a pyramid stereo matching strategy from coarse to fine, and determines the parallax searching range of each pixel in the lower layer of image according to the parallax value of the higher layer of image of the pyramid, thereby greatly reducing the parallax searching range and the matching calculation time, and obtaining higher processing speed for the high-resolution inspector image.

(3) In order to overcome the problem that the conventional method has poor matching accuracy in a parallax discontinuous area and an occlusion area, the invention provides the method for extracting the edge characteristics of the image and adaptively selecting a penalty coefficient aiming at the area containing the edge, and the result shows that the high matching accuracy is obtained in the parallax discontinuous area and the occlusion area.

Drawings

FIG. 1 is a flow chart of a stereo matching method of the extraterrestrial celestial body rover of the present invention;

FIG. 2A is a left image obtained by the left camera of the lunar rover after epipolar line correction, the image size being 947 × 696 pixels;

FIG. 2B is a right image of the lunar vehicle, obtained by the right camera of the present invention, after epipolar line correction, the image size being 947 × 696 pixels;

FIG. 3 is a flowchart of obtaining a left disparity map and a right disparity map from a left and right pyramid top-level image according to the present invention;

FIG. 4A is a SIFT matching result of the top-level image of the left and right pyramids according to the embodiment of the present invention;

fig. 4B is a spatial plane obtained by fitting SIFT points in the left pyramid highest-level image and spatial points formed by parallax according to the embodiment of the present invention;

FIG. 4C is a Canny edge detection diagram of the left pyramid top level image;

fig. 4D is a left disparity map corresponding to a highest-level image of the left pyramid in the embodiment of the present invention;

FIG. 4E is the disparity map shown in FIG. 4D after it has been superimposed with an edge;

FIG. 5A is a diagram illustrating a disparity search range of a lower layer image according to a disparity map of a higher layer image in accordance with the present invention;

FIG. 5B is a schematic diagram of the pixel disparity search range when left and right inconsistent points exist in a 3 × 3 region according to the present invention;

FIG. 6A is an edge map of left pyramid layer 2 image extraction with an image size of 473X 348 pixels in accordance with an embodiment of the present invention;

fig. 6B is a disparity map corresponding to a layer 2 image of the left pyramid in accordance with the embodiment of the present invention;

FIG. 6C is the disparity map shown in FIG. 6B after it has been superimposed with an edge;

FIG. 7A is an edge diagram of the left pyramid layer 1 image extraction according to an embodiment of the present invention;

fig. 7B is a disparity map corresponding to a layer 1 image of the left pyramid in accordance with the embodiment of the present invention;

FIG. 7C is the disparity map shown in FIG. 7B after being superimposed with an edge;

fig. 7D is a disparity map obtained by the conventional MGM method after edge superimposition.

Detailed Description

In order to provide an embodiment of image matching of an extraterrestrial celestial object rover, the embodiment of the invention is described below with reference to the accompanying drawings.

The invention provides a stereo matching method of an extraterrestrial celestial body patrolling device, which comprises the following steps of:

and S1, acquiring a stereo image pair of the extraterrestrial celestial body scene by using the left camera and the right camera of the extraterrestrial celestial body rover, wherein the image acquired by the left camera is defined as a left image, and the image acquired by the right camera is defined as a right image.

And S2, performing epipolar line correction on the left image and the right image in the step S1 respectively, constructing an L-layer pyramid, namely a left pyramid, on the left image subjected to epipolar line correction, and constructing an L-layer pyramid, namely a right pyramid, on the right image subjected to epipolar line correction, wherein L is preferably 2-5, and the pyramid in the invention is preferably a Gaussian pyramid.

Fig. 2A and 2B show left and right images obtained by left and right cameras of a lunar rover after epipolar line correction, the image sizes of which are 947 × 696 pixels, and in the embodiment of the present invention, a left gaussian pyramid and a right gaussian pyramid of 3 layers are respectively constructed for the pair of images, so that the image size of the highest layer (i.e., layer 3) is 236 × 174 pixels.

S3, performing bidirectional stereo matching on the image of the left pyramid highest layer, i.e., the L-th layer, and the image of the right pyramid highest layer, i.e., the L-th layer, to obtain a disparity map of the left pyramid L-th layer and a disparity map of the right pyramid L-th layer, where a main flow is as shown in fig. 3, and the specific steps are as follows in this embodiment:

s31: SIFT (Scale Invariant Feature Transform) the invention applies SIFT algorithm to realize the Feature point extraction and matching of two images at the highest layer, and mainly comprises the steps of Scale space extreme value detection, key point positioning, key point local direction calculation, descriptor generation, SIFT matching, RANSAC mismatching elimination and the like.

According to the embodiment of the invention, SIFT feature points (left SIFT feature points) in the highest-level image of the left pyramid and SIFT feature points (right SIFT feature points) in the highest-level image of the right pyramid are extracted, the extracted left SIFT feature points and the extracted right SIFT feature points are matched to obtain n pairs of matching points, and the coordinates of the feature points of the ith pair of matching points in the highest-level image of the left pyramid are marked as (x)_i,y_i) And the coordinates of the feature point of the ith pair of matching points in the right pyramid highest-level image in the n pairs of matching points are recorded as (x)_i',y_i'), i is 0,1, …, n-1, n is not less than 3. Fig. 4A shows SIFT feature point extraction results and matching results of the top-level image of the left and right pyramids, which includes 260 pairs of matching points.

S32: after the SIFT feature points are matched, the conventional method is to determine the parallax search range of the whole image according to the SIFT matching result. Under a binocular stereo camera imaging device, the parallax is inversely proportional to the distance of an observer. The extraterrestrial celestial body patrolling device images the front, and the two obtained images have the following characteristics: the parallax value of the scene point close to the patrol device in the image is larger, and the parallax value of the scene point far away from the patrol device in the image is relatively small. Therefore, if each pixel in the whole image is assigned with the same parallax search range according to the conventional method, the calculation amount is enormous, but most of the calculation is redundant, and mismatching may be caused.

D is calculated from the i-th pair of matching points among the n pairs of matching points obtained in step S31_i＝x_i'-x_iObtaining a spatial point corresponding to the ith characteristic point in the highest-layer image of the left pyramid, and recording the coordinate of the spatial point as (x)_i,y_i,d_i) And n space points are obtained by n feature points in the highest layer image of the left pyramid. Then, performing space plane fitting on n space points obtained by n feature points in the left pyramid highest-layer image, and calculating to obtain a parallax search range of all pixel points of the left pyramid highest-layer image, wherein the specific steps comprise:

s321, performing space plane fitting on n space points obtained from the n feature points in the highest-level image of the left pyramid, and obtaining an equation of a space plane Π, which is expressed as D ═ aX + bY + c.

Wherein a, b and c represent coefficients in a space plane equation; (X, Y, D) is the coordinate of any point on the space plane Π, and D is the parallax corresponding to the coordinate (X, Y).

To minimize S, one should satisfy:

the spatial plane coefficients a, b, c can be obtained.

Fig. 4B is a spatial plane obtained by fitting 260 spatial points corresponding to 260 SIFT feature points in the left pyramid top-level image with 260 matching points, where a is-0.0052, B is-0.226, and c is 30.5.

S322, regarding the coordinate (x) of the ith feature point in the n feature points in the highest-level image of the left pyramid_i,y_i) Calculating the parallax D corresponding to the ith characteristic point according to the equation D ═ aX + bY + c of the space plane pi_i ^ΠCalculating the parallax difference value of the ith characteristic point as d_i ^Π-d_i(ii) a For n parallax difference values corresponding to n feature points in the highest-layer image of the left pyramid, counting the maximum value D in the n parallax difference values_maxAnd minimum D_min。

S323, for any pixel p (x, y) in the highest-level image of the left pyramid, calculating the parallax D of the pixel p (x, y) according to the equation D ═ aX + bY + c of the spatial plane Π^Π(x, y), calculating d_min(x,y)＝d^Π(x,y)+D_min-Δd，d_max(x,y)＝d^Π(x,y)+D_max+ Δ d, where Δ d is a reserved parallax margin, preferably 2-5, and the parallax search range of the pixel point p (x, y) is determined as [ d [_min(x,y),d_max(x,y)]The integral number in the left pyramid is obtained to obtain the parallax search range with pixel points in the highest-layer image of the left pyramidAnd (5) enclosing.

D is calculated from the i-th pair of matching points among the n pairs of matching points obtained in step S31_i'＝x_i-x_i', obtaining a spatial point corresponding to the ith characteristic point in the highest-layer image of the right pyramid, and recording the coordinate of the spatial point as (x)_i',y_i',d_i') obtaining n space points by n feature points in the image of the highest layer of the right pyramid, carrying out space plane fitting on the n space points obtained by the n points in the image of the highest layer of the right pyramid, and calculating to obtain the parallax search range of all pixel points of the image of the highest layer of the right pyramid.

S33: performing Census transformation on all pixel points in the left pyramid highest-layer image to obtain Census sequences of all pixel points in the left pyramid highest-layer image, and performing Census transformation on all pixel points in the right pyramid highest-layer image to obtain Census sequences of all pixel points in the right pyramid highest-layer image; calculating to obtain the matching cost of all the pixel points in the left pyramid highest-layer image within the parallax search range of all the pixel points in the left pyramid highest-layer image according to the obtained Census sequences of all the pixel points in the left pyramid highest-layer image and the Census sequences of all the pixel points in the right pyramid highest-layer image; and calculating the matching cost of all the pixel points in the image of the highest layer of the right pyramid in the parallax search range of all the pixel points in the image of the highest layer of the right pyramid according to the obtained Census sequences of all the pixel points in the image of the highest layer of the left pyramid and the Census sequences of all the pixel points in the image of the highest layer of the right pyramid.

The Census conversion method for all pixel points in the left pyramid highest-level image to obtain Census sequences of all pixel points in the left pyramid highest-level image comprises the following steps:

s331, establishing a Census window by taking any pixel point in the highest-layer image of the left pyramid as a center, preferably 3 × 3 pixels or 5 × 5 pixels in the size of the window, comparing the size relationship between pixel points except the center point and the gray value of the center point in the established Census window, marking the pixel point with the gray value less than or equal to the gray value of the center point as 0, otherwise marking the pixel point as 1, and finally bitwise connecting the comparison result into 0/1 binary code streams which are called Census sequences of the current pixel points.

And S332, repeating the step S331 until all the pixel points in the highest-layer image of the left pyramid are processed, and obtaining Census sequences of all the pixel points in the highest-layer image of the left pyramid.

According to the Census sequences of all the pixel points in the left pyramid highest-layer image and the Census sequences of all the pixel points in the right pyramid highest-layer image, in the parallax search range of all the pixel points in the left pyramid highest-layer image, the method for calculating the matching cost of all the pixel points in the left pyramid highest-layer image comprises the following steps:

s333, calculating the Hamming distance of the Census sequence of the pixel point p (x, y) in the left pyramid highest-level image and the pixel point q (x-d, y) in the right pyramid highest-level image, wherein d is parallax, and d belongs to [ d [ d ]_min(x,y),d_max(x,y)]And d is an integer, and the matching cost C (p, d) of the pixel point p (x, y) in the highest-layer image of the left pyramid when the parallax is d is obtained.

And S334, repeating the step S333 to obtain the matching cost of the pixel point p (x, y) in the highest-layer image of the left pyramid in different parallaxes d.

S34: and performing matching cost aggregation on all pixel points in the left pyramid highest-layer image according to the matching cost of all pixel points in the left pyramid highest-layer image obtained in the step S33 to obtain a matching cost aggregation value of the left pyramid highest-layer image, and obtaining an initial disparity map of the left pyramid highest-layer image according to the obtained matching cost aggregation value of the left pyramid highest-layer image.

Based on the calculated matching cost, an energy function is constructed as follows:

in the above energy function, p is the pixel in the image, N_pIs a neighborhood of pixel P, P₁、P₂Penalty factor for parallax change in neighborhood of pixel P and P₂>P₁>0，P₁Penalizing minor parallax variations, P₂Punishment of larger parallax variation, D is a dense parallax map of the image to be matched, if T [ ·]If the parameter in the function is true, 1 is returned, otherwise 0 is returned.

To solve the problem of global minimization, the SGM algorithm and the MGM algorithm convert the global problem into a multi-directional dynamic programming optimization calculation strategy, and both of them show excellent effects. However, the two methods do not consider the sudden parallax change, and a penalty parameter with a constant value is used in all directions of all pixel points, so that the parallax calculation of the two algorithms in the area near the sudden parallax change is incorrect. Generally, since extraterrestrial celestial bodies tend to have a large number of rocks, there is usually a sharp sudden change in parallax, and the position of the sudden change in parallax in the parallax map occurs at the edge. Therefore, the embodiment of the invention firstly extracts the edge information of the left image and the right image; then, in the calculation process of matching cost aggregation, different penalty parameters are adopted at the edge position with large parallax change.

S341, because the Canny edge detection operator has the advantages of low error rate edge detection, optimal positioning, no false edge caused by image noise and the like, the Canny edge detection operator is firstly used for detecting the edges of the left and right images, and the method mainly comprises the steps of smoothing the images by a Gaussian filter, calculating the amplitude and the direction of the gradient by finite difference of first-order partial derivatives, applying non-maximum suppression to the gradient amplitude, detecting and connecting the edges by a dual-threshold algorithm and the like.

Canny edge feature detection is carried out on the highest layer image of the left pyramid, and a plurality of edges are obtained; the points on the edge are called edge points, and for any edge point, the range [ d ] is searched according to the parallax of the pixel point corresponding to the edge point_min(x,y),d_max(x,y)]Calculating a difference d_max(x,y)-d_min(x, y), called the parallax span of the edge point, calculating the average value of the parallax spans of all the edge points on the current strip edge, if the average value is less than the set threshold value, removing the current strip edge, if the average value is not less than the set threshold value, keeping the current strip edge; and processing a plurality of edges obtained from the left pyramid highest-layer image to obtain the eligible edges in the left pyramid highest-layer image. The set threshold value is preferably a real number in the range of 2-8.

Fig. 4C is a Canny edge detection diagram of the image of the highest layer of the left pyramid, and it can be seen that the place where the parallax mutation is severe at the boundary between the rock and the lunar soil on the moon corresponds to the edge position in the image. When calculating the average value of the parallax span, different thresholds can be set according to the pyramid layer where the image is located, in the embodiment of the present invention, the threshold value set in the 3 rd layer image is 2.5.

And S342, performing matching cost aggregation on all pixel points in the left pyramid highest-layer image according to the matching costs of all pixel points in the left pyramid highest-layer image to obtain a matching cost aggregation value.

Because the dynamic programming has strong correlation constraint to one direction in the calculation process of the direction and stripes are easy to generate in the implementation, the matching cost aggregation is calculated by the sum of the matching costs on the optimization paths in a plurality of directions. Calculating a matching cost aggregation value S (p, d) of the pixel point p (x, y) when the parallax is d:

in the formula, r represents the propagation direction of the pixel point p (x, y), and r belongs to { (-1,0), (1,0), (0,1), (0, -1), (-1, -1), (1, -1), (1,1), (-1,1) } when the propagation direction r takes 8 directions, or r belongs to { (-1,0), (1,0), (0,1), (0, -1) } when the propagation direction r takes 4 directions.

L_r(p, d) represents the matching cost of the pixel point p (x, y) along the direction r;

wherein R represents m directions derived from the propagation direction R, and R ∈ { R ∈ {₁,…,r_m}，m≥1，r₁In the same direction as r, i.e. r₁Where m according to the invention can be chosen to be 4, 2 or 1.

R∈{(-1,0),(0,-1),(-1,-1),(1,-1)}、R∈{(1,0),(0,1),(1,1),(-1,1)}、

R∈{(0,1),(-1,0),(-1,1),(-1,-1)}、R∈{(0,-1),(1,0),(1,-1),(1,1)}

R∈{(-1,-1),(1,-1),(0,-1),(1,0)}、R∈{(1,-1),(1,1),(1,0),(0,1)}、

R∈{(1,1),(-1,1),(0,1),(-1,0)}、R∈{(-1,1),(-1,-1),(-1,0),(0,-1)}

R∈{(-1,0),(0,-1),(-1,-1),(1,-1)}、R∈{(1,0),(0,1),(1,1),(-1,1)}、

R∈{(0,1),(-1,0),(-1,1),(-1,-1)}、R∈{(0,-1),(1,0),(1,-1),(1,1)}

when m is 2 and the propagation directions R are 4, R derived from the 4 propagation directions is sequentially represented as:

R∈{(-1,0),(0,-1)}、R∈{(1,0),(0,1)}、R∈{(0,1),(-1,0)}、R∈{(0,-1),(1,0)}；

R∈{(-1,0)}、R∈{(1,0)}、R∈{(0,1)}、R∈{(0,-1)}、

R∈{(-1,-1)}、R∈{(1,-1)}、R∈{(1,1)}、R∈{(1,1)}

when m is 1 and the propagation directions R are 4, R derived from the 4 propagation directions is sequentially represented as: r belongs to { (-1,0) }, R belongs to { (0,1) }, R belongs to { (0, -1) }.

p-R represents the previous pixel point of the pixel point p along the direction R; l is_r(p-R, d) represents the matching cost of the pixel p-R when the disparity is d, L_r(p-R, d +1) represents the matching cost of the pixel p-R when the disparity is d +1, L_r(p-R, d-1) represents the matching cost of the pixel point p-R when the parallax is d-1; d 'is the disparity of pixel p-R, whose disparity search range d' ∈ [ d ∈ ]/_min(p-R),d_max(p-R)]，

Indicating that pixel p-R is in the disparity search range d_min(p-R),d_max(p-R)]The inner minimum path matching cost; p₁ ^R、P₂ ^RDenotes a penalty parameter, P₂ ^R>P₁ ^R。

Parameter P₁ ^RAnd P₂ ^RThe choice of (2) has a large influence on the final matching result, and the magnitude of (2) determines the smoothness of the parallax. If the penalty parameter is large, the parallax is smooth, but edge details cannot be reserved in the parallax abrupt change area, and the edge positioning accuracy corresponding to the parallax is low. The invention ensures the edge details of the barrierAnd selecting a penalty parameter according to the edge information when the edge is accurate. For a flat area, namely a parallax non-abrupt area, a larger penalty parameter is adopted to strengthen the parallax smoothness constraint; and for the non-flat area at the parallax abrupt change position, a smaller penalty parameter is adopted to reserve the parallax abrupt change at the edge, and the positioning accuracy of the abrupt change edge is improved.

wherein, P₁And P₂Is referred to the value of the MGM method, P₁Preferably 8, P₂Preferably 32.

And processing other pixel points in the left pyramid highest-layer image according to the step S342 to obtain the matching cost aggregate value of all the pixel points in the left pyramid highest-layer image.

and processing other pixel points in the left pyramid highest-layer image according to the step S343 to obtain the parallax values of all the pixel points in the left pyramid highest-layer image, and forming an initial parallax image of the left pyramid highest-layer image.

And (3) sub-pixel parallax calculation: the initial disparity maps of the left and right pyramids are all at a pixel level, which causes severe aliasing phenomenon on some continuous planes, especially in scenes such as large-inclination planes and curved surfaces which are not directly opposite to a camera. Therefore, there is a need to improve the matching accuracy from the pixel level to the sub-pixel level so that the parallax on the target surface presents a natural smooth transition. The invention adopts a quadratic curve fitting method to obtain a parallax map with sub-pixel precision.

Left and right consistency detects the initial left and right disparity maps to realize the detection of the occlusion areas and the mismatching points, which requires that the disparities of the left and right disparity maps are consistent.

Setting the sub-pixel parallax map of the highest layer image of the left pyramid as D_lr(x, y), the sub-pixel parallax map of the highest layer image of the right pyramid is marked as D_rl(x, y) for the left pyramidAny pixel point in the sub-pixel disparity map of the highest layer image if | D is satisfied_lr(x,y)-D_rl(x+D_lr(x,y),y)|<And T, if T is 1, the pixel point is a left and right consistency point, otherwise, the pixel point is a left and right inconsistency point.

The present embodiment is obtained when the propagation direction R is 8 directions and m in the derivative direction R is 4. Fig. 4D is a disparity map corresponding to the highest-level image of the left pyramid in the method of the present invention, and fig. 4E is a disparity map after the disparity map shown in fig. 4D is superimposed on a plurality of edges satisfying a disparity span condition.

S4, determining the parallax search range of all pixel points in the L-1 layer image of the left pyramid by using the parallax map of the L layer of the left pyramid obtained in the step S3; determining the parallax search range of all pixel points in the L-1 layer image of the right pyramid by using the parallax map of the L layer of the right pyramid obtained in the step S3; and performing stereo matching from left to right on the L-1 layer image of the left pyramid and the L-1 layer image of the right pyramid according to the determined parallax search ranges of all pixel points in the L-1 layer image of the left pyramid, and performing stereo matching from right to left on the L-1 layer image of the left pyramid and the L-1 layer image of the right pyramid according to the determined parallax search ranges of all pixel points in the L-1 layer image of the right pyramid to obtain a parallax image of the L-1 layer of the left pyramid and a parallax image of the L-1 layer of the right pyramid.

The method for determining the parallax search range of all pixel points in the L-1 layer image of the left pyramid by using the parallax map of the L layer obtained in the step S3 includes:

s42: if all the pixels in the window are left-right consistency points, as shown in fig. 5A, the disparity maximum value d corresponding to all the pixels in the window is calculated_{L max}And the minimum value d of parallax_{L min}：

d_{L max}＝max{d_L(x_L+i,y_L+j)|i＝-1,0,1,j＝-1,0,1}

d_{L min}＝min{d_L(x_L+i,y_L+j)|i＝-1,0,1,j＝-1,0,1}

Pixel point p in L-1 layer image of left pyramid_L-1(x_L-1,y_L-1) Is set to [2d ]_{L min}-,2d_{L max}+]The integer of (1) to (3) is preferable as the remainder.

S43: if there are left and right inconsistency points for the pixels in the window, as shown in FIG. 5B, then at the y-th_LLine, with p_L(x_L,y_L) As a starting point, searching leftwards until a left-right consistency point p is obtained_l(x_L-x_l,y_L) Searching right until left and right consistency points p are obtained_r(x_L+x_r,y_L) (ii) a At the y-th_L-1 line with p_L(x_L,y_L-1) as a starting point, searching to the left until a left and right consistency point p is obtained_ul(x_L-x_ul,y_L-1), search right until left and right consistency points p are obtained_ur(x_L+x_ur,y_L-1); at the y-th_L+1 lines, with p_d(x_L,y_L+1) as starting point, searching left until obtaining left and right consistency point p_dl(x_L-x_dl,y_L+1), searching right until left and right consistency points p are obtained_dr(x_L+x_dr,y_L+1)。

Is an integer of (1).

S44: and repeating the steps of S41-S43 until all the pixel points in the L-1 layer image of the left pyramid are processed, and obtaining the parallax search range of all the pixel points in the L-1 layer image of the left pyramid.

The method comprises the following steps of performing left-to-right stereo matching on an L-1 layer image of a left pyramid and an L-1 layer image of a right pyramid according to the determined parallax search ranges of all pixel points in the L-1 layer image of the left pyramid, and performing right-to-left stereo matching on the L-1 layer image of the left pyramid and the L-1 layer image of the right pyramid according to the determined parallax search ranges of all pixel points in the L-1 layer image of the right pyramid to obtain a left parallax image of the L-1 layer of the left pyramid and a right parallax image of the L-1 layer of the right pyramid, wherein the specific steps are as follows:

Fig. 6A is an edge map extracted from a left pyramid layer 2 image, the image size is 473 × 348 pixels, fig. 6B is a left disparity map corresponding to the left pyramid layer 2 image obtained by the method of the present invention, fig. 6C is a disparity map obtained by superimposing the disparity map shown in fig. 6B and an edge satisfying a disparity span condition, and the threshold value set in the averaging process of the disparity spans of the layer is 5.

S5, determining the parallax search range of all pixel points in the L-2 layer image of the left pyramid by using the parallax map of the L-1 layer of the left pyramid obtained in the step S4; determining the parallax search range of all pixel points in the L-2 layer image of the right pyramid by using the parallax map of the L-1 layer of the right pyramid obtained in the step S4; and performing stereo matching from left to right on the L-2 layer image of the left pyramid and the L-2 layer image of the right pyramid according to the determined parallax search ranges of all pixel points in the L-2 layer image of the left pyramid, and performing stereo matching from right to left on the L-2 layer image of the left pyramid and the L-2 layer image of the right pyramid according to the determined parallax search ranges of all pixel points in the L-2 layer image of the right pyramid to obtain a parallax image of the L-2 layer of the left pyramid and a parallax image of the L-2 layer of the right pyramid.

Then, the process is circularly carried out until the parallax image of the left pyramid layer 2 is used for determining the parallax searching range of all pixel points in the left pyramid layer 1 image (namely the left image after epipolar line correction); determining the parallax search range of all pixel points in the image of the 1 st layer of the right pyramid (namely the right image after epipolar line correction) by using the parallax image of the 2 nd layer of the right pyramid; and performing stereo matching from left to right on the left pyramid layer-1 image and the right pyramid layer-1 image according to the determined parallax search ranges of all pixel points in the left pyramid layer-1 image, and performing stereo matching from right to left on the left pyramid layer-1 image and the right pyramid layer-1 image according to the determined parallax search ranges of all pixel points in the right pyramid layer-1 image to obtain a parallax image of the left pyramid layer-1 and a parallax image of the right pyramid layer-1.

Fig. 7A is an edge map extracted from the left pyramid layer 1 image (left image after epipolar line correction, i.e., image shown in fig. 2A), fig. 7B is a disparity map corresponding to the left pyramid layer 1 image in the method of the present invention, fig. 7C is a disparity map obtained by superimposing the disparity map shown in fig. 7B and an edge satisfying a condition, where when an edge satisfying a condition is obtained, the threshold set in the layer 1 image is 7.5, and fig. 7D is a disparity map obtained by superimposing the disparity map obtained by the conventional MGM method and an edge satisfying a disparity condition. As can be seen by comparing fig. 7C and 7D, the conventional MGM method is inaccurate in calculating the parallax of the area near the edge and inaccurate in positioning the sudden change of parallax, but the method of the present invention can accurately position the sudden change of parallax, so that the accurate parallax can be obtained in the area near the edge.

And (5) analyzing the computational complexity. The conventional MGM method requires D times of calculation for each pixel of the original image with resolution W × H, where D is the parallax span, i.e. the calculation complexity is O (W × H × D), where the parallax search range of the whole original image is [ -40,150], and the parallax span D is 190 pixels, so the calculation complexity is represented as O (947 × 696 × 190). For the computational complexity analysis of the method of the present invention, it is necessary to count the disparity search ranges of all pixels of all layers in the pyramid, where the matching computational complexity of the two images of the highest layer (layer 3) is O (734903), the matching computational complexity of the two images of layer 2 is O (611766), and the matching computational complexity of the two images of layer 1 (original) is O (2457204). The ratio of the computational complexity of the two methods is:

Claims

1. a stereo matching method for an extraterrestrial celestial body patrolling device is characterized by comprising the following steps:

s1, acquiring two images of the extraterrestrial celestial body scene by using a left camera and a right camera of the extraterrestrial celestial body patrolling device, wherein the image acquired by the left camera is defined as a left image, and the image acquired by the right camera is defined as a right image;

circularly carrying out;

s6, determining the parallax search range of all pixel points in the left pyramid layer 1 image, namely the epipolar line corrected left image, by using the parallax map of the left pyramid layer 2; determining the parallax search range of all pixel points in the right pyramid layer 1 image, namely the right image after epipolar line correction, by using the parallax map of the right pyramid layer 2; performing left-to-right stereo matching on the left pyramid layer-1 image and the right pyramid layer-1 image according to the determined parallax search ranges of all pixel points in the left pyramid layer-1 image, and performing right-to-left stereo matching on the left pyramid layer-1 image and the right pyramid layer-1 image according to the determined parallax search ranges of all pixel points in the right pyramid layer-1 image to obtain a left pyramid layer-1 parallax image and a right pyramid layer-1 parallax image;

in the step S2, L is a natural number and takes a value of 2-5; the pyramid is a Gaussian pyramid;

in the step S3, the method for performing bidirectional stereo matching on the image at the highest layer of the left pyramid, i.e., the L-th layer, and the image at the highest layer of the right pyramid, i.e., the L-th layer, constructed in the step S2 includes the following specific steps:

s31: extracting SIFT feature points in the highest-layer image of the left pyramid, wherein the SIFT feature points are called left SIFT feature points, extracting SIFT feature points in the highest-layer image of the right pyramid, wherein the SIFT feature points are called right SIFT feature points, and matching the extracted left SIFT feature points and the right SIFT feature points to obtain n pairs of matching points; the coordinate of the characteristic point corresponding to the ith matching point in the highest layer image of the left pyramid in the n pairs of matching points is recorded as (x)_i,y_i) And the coordinates of the feature point corresponding to the ith pair of matching points in the n pairs of matching points in the image of the highest layer of the right pyramid are recorded as (x)_i',y_i')，i＝0,1,…,n-1，n≥3；

2. The stereo matching method for the extraterrestrial celestial body rover according to claim 1, wherein: in step S32, the method for obtaining the parallax search range with pixels in the left pyramid highest-level image by performing spatial plane fitting on n spatial points obtained from n feature points in the left pyramid highest-level image includes:

According to

Obtaining space plane coefficients a, b and c;

S323, for any pixel point p (x, y) in the highest-level image of the left pyramid, calculating a parallax D corresponding to the pixel point p (x, y) according to an equation D ═ aX + bY + c of the spatial plane Π^Π(x, y), calculating d_min(x,y)＝d^Π(x,y)+D_min-Δd，d_max(x,y)＝d^Π(x,y)+D_max+ Δ d, where Δ d is the reserved disparity margin, and the disparity search range for pixel point p (x, y) is determined to be [ d [ d ] ]_min(x,y),d_max(x,y)]Obtaining a parallax searching range with pixel points in the highest layer image of the left pyramid by using the integer in the left pyramid;

3. The stereo matching method for the extraterrestrial celestial body rover according to claim 1, wherein: in step S33, Census conversion is performed on all the pixels in the left pyramid highest-level image, and a Census sequence of all the pixels in the left pyramid highest-level image is obtained by the method including:

and repeating the steps S331-S332 according to the same principle to obtain the Census sequences of all the pixel points in the right pyramid highest-level image.

4. The stereo matching method of the extraterrestrial celestial body rover as claimed in claim 1 or 2, wherein: in step S33, the method for calculating the matching cost of all the pixels in the left pyramid highest-level image in the parallax search range of all the pixels in the left pyramid highest-level image according to the Census sequences of all the pixels in the left pyramid highest-level image and the Census sequences of all the pixels in the right pyramid highest-level image includes:

and (4) calculating the matching cost of all the pixel points in the highest-layer image of the left pyramid according to the parallax searching range of all the pixel points in the highest-layer image of the left pyramid, and repeating the steps S333-S335 according to the same principle to obtain the matching cost of all the pixel points in the highest-layer image of the right pyramid.

5. The stereo matching method for the extraterrestrial celestial body rover according to claim 1, wherein: in step S34, the method for aggregating the matching costs of all the pixel points in the left pyramid highest-level image according to the matching costs of all the pixel points in the left pyramid highest-level image obtained in step S33 to obtain a matching cost aggregate value includes:

s341, Canny edge feature detection is carried out on the left pyramid highest-layer image to obtain a plurality of edges; the points on the edge are called edge points, and for any edge point, the range [ d ] is searched according to the parallax of the pixel point corresponding to the edge point_min(x,y),d_max(x,y)]Calculating a difference d_max(x,y)-d_min(x, y), called the parallax span of the edge point, calculating the average value of the parallax spans of all the edge points on each edge, if the average value is less than the set threshold, removing the current edge, if the average value is not less than the set threshold, keeping the current edge; processing a plurality of edges obtained by detecting the left pyramid highest-layer image to obtain qualified edges in the left pyramid highest-layer image;

wherein R represents m directions derived from the propagation direction R, and R ∈ { R ∈ {₁,…,r_m}，m≥1，r₁In the same direction as r, i.e. r₁R; p-R represents the previous pixel point of the pixel point p along the direction R; l is_r(p-R, d) represents the matching cost of the pixel p-R when the disparity is d, L_r(p-R, d +1) represents the matching cost of the pixel p-R when the disparity is d +1, L_r(p-R, d-1) represents the matching cost of the pixel point p-R when the parallax is d-1; d 'is the disparity of pixel p-R, whose disparity search range d' ∈ [ d ∈ ]/_min(p-R),d_max(p-R)]And d' is an integer, and is,

indicating that pixel p-R is in the disparity search range d_min(p-R),d_max(p-R)]The inner minimum path matching cost; p₁ ^R、P₂ ^RDenotes a penalty parameter, P₂ ^R＞P₁ ^R；

wherein, P₁Is 8, P₂Is 32;

6. The stereo matching method of the extraterrestrial celestial body rover as claimed in claim 5, wherein: the propagation direction r is 8 directions,

r∈{(-1,0),(1,0),(0,1),(0,-1),(-1,-1),(1,-1),(1,1),(-1,1)}，

or in the 4 directions, the direction of the light beam,

r∈{(-1,0),(1,0),(0,1),(0,-1)}；

m is 4, 2 or 1;

R∈{(-1,0),(0,-1),(-1,-1),(1,-1)}、R∈{(1,0),(0,1),(1,1),(-1,1)}、

R∈{(0,1),(-1,0),(-1,1),(-1,-1)}、R∈{(0,-1),(1,0),(1,-1),(1,1)}

R∈{(-1,-1),(1,-1),(0,-1),(1,0)}、R∈{(1,-1),(1,1),(1,0),(0,1)}、

R∈{(1,1),(-1,1),(0,1),(-1,0)}、R∈{(-1,1),(-1,-1),(-1,0),(0,-1)}

R∈{(-1,0),(0,-1),(-1,-1),(1,-1)}、R∈{(1,0),(0,1),(1,1),(-1,1)}、

R∈{(0,1),(-1,0),(-1,1),(-1,-1)}、R∈{(0,-1),(1,0),(1,-1),(1,1)}

R∈{(-1,0)}、R∈{(1,0)}、R∈{(0,1)}、R∈{(0,-1)}、

R∈{(-1,-1)}、R∈{(1,-1)}、R∈{(1,1)}、R∈{(1,1)}

7. The stereo matching method for the extraterrestrial celestial body rover according to claim 1, wherein: in the step S35, a quadratic curve fitting method is used for the sub-pixel parallax calculation;

setting the sub-pixel parallax map of the highest layer image of the left pyramid as D_lr(x, y), the sub-pixel parallax map of the highest layer image of the right pyramid is marked as D_rl(x, y), if any pixel point in the sub-pixel disparity map of the highest-layer image of the left pyramid satisfies | D_lr(x,y)-D_rl(x+D_lr(x,y),y)|<And T, if T is 1, the pixel point is a left and right consistency point, otherwise, the pixel point is a left and right inconsistency point.

8. The stereo matching method for the extraterrestrial celestial body rover according to claim 1, wherein: in step S4, the method for determining the parallax search range of all the pixels in the L-1 th layer image of the left pyramid by using the L-th layer parallax map of the left pyramid obtained in step S3 includes:

s41: to any pixel point p in the L-1 layer image of the left pyramid_L-1(x_L-1,y_L-1)，Computing

s42: if all the pixels in the window are left-right consistency points, calculating the parallax maximum value d corresponding to all the pixels in the window_LmaxAnd the minimum value d of parallax_Lmin：

d_Lmax＝max{d_L(x_L+i,y_L+j)|i＝-1,0,1,j＝-1,0,1}

d_Lmin＝min{d_L(x_L+i,y_L+j)|i＝-1,0,1,j＝-1,0,1}

Pixel point p in L-1 layer image of left pyramid_L-1(x_L-1,y_L-1) Is set to [2d ]_Lmin-,2d_Lmax+]The integer of (1) is the rest;

The whole number of (1);

9. The stereo matching method for the extraterrestrial celestial body rover according to claim 1, wherein: in step S4, the left pyramid L-1 layer image and the right pyramid L-1 layer image are stereoscopically matched from left to right according to the determined parallax search ranges of all the pixels in the left pyramid L-1 layer image, and the left pyramid L-1 layer image and the right pyramid L-1 layer image are stereoscopically matched from right to left according to the determined parallax search ranges of all the pixels in the right pyramid L-1 layer image, so as to obtain a left parallax map of the left pyramid L-1 layer and a right parallax map of the right pyramid L-1 layer, which specifically include the following steps: