CN110211043B - A Registration Method Based on Grid Optimization for Panoramic Image Stitching - Google Patents

Abstract

The invention belongs to the technical field of digital images, and particularly relates to a grid optimization-based registration method suitable for free view image stitching. In recent years, panoramic images have received great attention due to their application prospects in the fields of virtual reality, medical imaging, and the like. Image registration mainly acquires the position parameters of each image in the panorama. The traditional image registration needs to execute a series of operations such as feature point extraction, matching, homography matrix solving, camera parameter correction and the like. The invention uses a registration method based on grid optimization to replace the classical registration technology, not only has obvious improvement in speed, but also is suitable for splicing images with large parallax photographed by free view angles. In the method, ORB rapid feature extraction is used for acquiring feature points of images, a strategy of step-by-step coarse-fine matching is introduced, and finally three constraint terms based on grid optimization are introduced to obtain optimal registration parameters between the images by a method of minimizing an error function.

Description

A Registration Method Based on Grid Optimization for Panoramic Image Stitching

technical field

The invention belongs to the technical field of digital image processing, and in particular relates to an image registration method based on grid optimization for panorama image mosaic.

Background technique

With the rapid development of information technology and the improvement of people's living standards, people's demand for high-quality panoramic images is also increasing. Whether it is virtual reality, panoramic live broadcast, medical imaging, or assisted driving, panoramic images are very good. application prospects. Functionally, people can get more information from panoramic images and get a better visual experience.

In order to obtain wide-angle images, the traditional method is to use a wide-angle lens, such as a fisheye lens, which can capture nearly 180 degrees of the entire scene on the horizontal plane, but this method has at least the following three shortcomings: First, it will introduce distortion visible to the naked eye and distortion, the second is that the resolution is reduced due to the large shooting angle, and the third is that the wide-angle lens is expensive.

Under such a background and demand, image mosaic technology came into being. This technology aims to combine a group of image sequences with overlapping areas and small viewing angles, through a series of processing and operations, and finally obtain a picture that contains all scene information. Panoramic images with large viewing angle characteristics and ultra-high resolution. Since each lens is only responsible for shooting a part of the scene, if a high-definition lens is used, the resolution of the details will be higher, and the normal lens will also avoid the distortion introduced by the wide-angle.

The traditional image stitching technology can be divided into the following steps: image acquisition, SIFT feature extraction, feature matching, homography matrix solution, camera parameter correction and image fusion. A homography matrix can only align one plane, if there is parallax and stretching in the initial image, very serious ghosting and skewing will appear in the mosaic obtained by using traditional methods, and even the mosaic will fail (the resulting image is completely distorted) .

Contents of the invention

The purpose of the present invention is to overcome the shortcomings of the current technology and provide an image registration method based on grid optimization for panoramic image mosaic.

The invention is suitable for image registration of free viewing angles, and exhibits excellent stability for image scale transformation, brightness transformation and rotation transformation. Therefore, this invention can not only be applied to panoramic monitoring with a fixed lens, but also to drone surveys with moving and rotating lenses.

For the acquisition of the image sequence, the present invention uses three methods: one, the fixed camera rotates and shoots to obtain the surrounding panorama, and the image sequence obtained like this is the most regular; Yes, the image sequence obtained in this way will have transformations such as rotation, shaking, and translation; third, the UAV acquires images at high altitude, the camera position is fixed but the shooting direction will be adjusted, and the obtained images have scale, rotation transformation and high resolution.

The grid optimization-based image registration method for panoramic image mosaic provided by the present invention, the specific steps are as follows:

(1) Apply ORB [1] for fast feature extraction. Use the FAST [2] algorithm for feature detection, and then generate the relevant descriptor of the feature point based on the improved BRIEF [3] descriptor, which contains the scale information, position and direction information of the feature point.

(2) Use the high-dimensional tree (KD tree) and the optimal node first algorithm (BBF) for rough matching of feature points, and then use the formula (1) to perform a ratio test on the obtained matching points, p is the current feature point, which is defined as p The ratio of the Hamming distance between the nearest feature point p _best-closed and the second adjacent feature point p _{second-closed} is less than a threshold (ratio), the threshold is generally taken as 0.65; apply the formula (2) for cross-testing; traverse the image I The feature points in the image J, find the matching point in the image J, recorded as M _I->J , then traverse the feature points in the image J, find the corresponding point in the image I, recorded as M _J->I , cross The test considers a pair of correct matches only when the two correspond to each other;

M _I&J ＝M _I->J ∩M _J->I (2)

(3) Then, perform multi-layer RANSAC [4] screening on the matching point pairs after fine matching, and filter the interior point sets of feature point pairs in multiple planes of the image, so that the number of final interior point sets accounts for the total number of matching pairs More than 80% of the matching information is retained to the greatest extent.

表示网格中某一匹配点的初始坐标(x,y,1)，/>

表示变换后的坐标(x′,y′,1)，均为三维，因此矩阵H的维度为3×3。(4) Through MDLT (Moving Direct Linear Transformation), the feature matching points are mapped to vertex matching points with more regular and uniform distribution. Divide the image into dense grids, and each grid corresponds to a homography matrix projection transformation, as shown in formula (3), where

Indicates the initial coordinates (x,y,1) of a matching point in the grid, />

Indicates the transformed coordinates (x′, y′, 1), all of which are three-dimensional, so the dimension of the matrix H is 3×3.

May wish to set H as:

Then formula (3) can be decomposed to get:

Organized into the format of matrix multiplication, it is:

a _i h = 0 (6)

Among them, a _i represents the 2×9 matrix formed by the i-th pair of matching points, and h is the column vector format of H (the dimension is 9×1):

h＝[h ₁ h ₂ h ₃ h ₄ h ₅ h ₆ h ₇ h ₈ h ₉ ] ^T

通过公式(4)计算；/>

为第i对匹配点对第k个网络的影响权值，由该匹配点与网格中心的距离决定。公式(4)中x_k表示第k个网格的中心坐标，x_i表示第i对匹配点的坐标，μ为调整参数。Taking all M pairs of matching points into consideration, as shown in formula (8), the transformation matrix h is solved by minimizing the square error e _k , where A is the transformation matrix corresponding to M pairs of matching points, and its dimension is 2M×9; W _k is The weight diagonal matrix of each pair of matching points, the dimension is 2M×2M, each element

Calculated by formula (4); />

is the influence weight of the i-th pair of matching points on the k-th network, which is determined by the distance between the matching point and the grid center. In formula (4), x _k represents the center coordinate of the k-th grid, x _i represents the coordinates of the i-th pair of matching points, and μ is an adjustment parameter.

e _k ＝argmin||W _k Ah|| ² st||h||＝1 (8)

Then use the calculated homography matrix to find the matching points of the grid vertices in another image. The distribution of the matching points of the vertices is even and regular, and they are used as the matching points in the grid optimization, which can effectively reduce the amount of calculation.

(5) Three constraint terms based on grid optimization are introduced: the coordinate alignment constraint term for overlapping regions, the local similarity constraint term for non-overlapping regions, and the global similarity constraint term for global structural consistency. First set several flags, store the matching relationship between images into the set T, and set the set of matching point pairs obtained after image i and image j are mapped in step (4) as M _ij , since the image has been divided into grids, Let V _i and E _i be the set of vertices and edges of the mesh in image i.

The coordinate alignment constraint item is shown in formula (10), and its function is to ensure that the coordinates of the corresponding matching points are as consistent as possible after grid optimization, and reduce the alignment error in the overlapping area of adjacent images. Among them, m(p) returns the matching point of feature point p on another image, and Ψ(p) uses the linear combination of four grid vertex coordinates to represent the position of feature point p.

表示边/>

的相似变换，具体计算如公式(12)，其中/>

和/>

是关于顶点变量的线性组合，主要模拟旋转和大小变换【5】；/>

和/>

表示原图像上同一条边的两个顶点位置，v_j′ⁱ和v′_k ^k表示网格优化后的顶点位置，E_i表示网格的边集合。The local similarity constraint item is shown in formula (11), and its main purpose is to ensure that the length and direction of the same vector edge will not change much before and after mesh optimization. Since the projection matrix is mainly suitable for overlapping areas, similar transformation matrices are introduced in very overlapping areas

means side />

The similarity transformation of , the specific calculation is as formula (12), where />

and />

It is a linear combination of vertex variables, mainly simulating rotation and size transformation [5]; />

and />

Indicates the two vertex positions of the same edge on the original image, v _j ′ ⁱ and v′ _k ^k represent the vertex positions after grid optimization, and E _i represents the edge set of the grid.

为每条边的权值，从重叠区域到非重叠区域渐变，距离重叠区域越远，权重越大，定义为公式(14)。s_i定义为图像i的尺度量，可通过光束平差法估计图像i的相机参数获得；θ_i定义为图像i对于基准图像的旋转量，本发明定义LSD【6】特征线检测得到的线特征之间的夹角的平均值为两幅图像之间的旋转角；参数/>

和/>

在公式(12)中已经有所提及。The global similarity constraint term is shown in formula (13), which aims to improve the overall structural consistency in the image sequence. in

is the weight of each edge, gradually changing from the overlapping area to the non-overlapping area, the farther away from the overlapping area, the greater the weight, defined as formula (14). s _i is defined as the scale quantity of image i, which can be obtained by estimating the camera parameters of image i through bundle adjustment method; θ _i is defined as the rotation amount of image i relative to the reference image, and the present invention defines the line obtained by LSD [6] characteristic line detection The average of the angles between the features is the rotation angle between the two images; parameter />

and />

It has been mentioned in formula (12).

是共享边/>

的网格集合(1个或2个，取决于是否为边界边)，Mⁱ表示图像i的重叠区域的全部网格的结合；d(q_k,Mⁱ)是一个函数，用以计算集合/>

内的网格q_k到重叠区域的距离；R_i和C_i表示图像i网格的行数和列数。Wherein, η and λ are adjustment parameters, which are generated by experiments. In the embodiment, η is 6, and λ is 20.

is the shared side/>

The set of grids (1 or 2, depending on whether it is a boundary edge), M ⁱ represents the combination of all grids in the overlapping area of image i; d(q _k , M ⁱ ) is a function to calculate the set />

The distance from the grid q _k within to the overlapping area; R _i and C _i represent the number of rows and columns of the grid of image i.

Equation (15) combines the three constraints, and the coordinate value of each pixel in the panorama after grid optimization can be obtained by minimizing this equation, and the image registration is completed.

In the formula, γ is the adjustment coefficient of the local similarity constraint item, which is 0.54 in the embodiment.

The above is the basic steps of the present invention, as shown in Figure 1.

The invention can effectively eliminate the ghost and tilt in the result image, improve the registration accuracy, and reduce the time consumption of registration at the same time.

Description of drawings

FIG. 1 is a flowchart of the architecture of the present invention.

Figure 2 shows the registration results and mosaic panorama of the image sequence captured by the hand-held camera.

Figure 3 shows the registration results and mosaic panorama of the high-resolution image sequence captured by the drone at high altitude.

Detailed ways

Below by example, further describe the method of the present invention in conjunction with accompanying drawing.

For a group of test image sequences A, use the method of the present invention to carry out image registration, the specific process is as follows:

1. Use the ORB algorithm to perform fast feature extraction on each image in A, and obtain the feature points of each image;

2. Introduce the K-D tree in the rough matching of feature points. Since the descriptor of each feature point is a 256-bit binary string, it can be used as 32-dimensional data, with 8 bits per dimension. It is applied to the K-D tree construction, and then found based on BBF search The points closest to this feature point are used as a pair of matching points;

3. Apply ratio test and cross test to the obtained matching point pairs, eliminate the wrong matches, and keep the remaining ones as fine matching point pairs;

4. Perform multi-layer RANAC screening on fine matching point pairs, and filter the interior point sets of feature point pairs in multiple planes of the image, so that the number of final interior point sets accounts for more than 80% of the total number of matching pairs, and the maximum retention matching information;

5. Through MDLT (Moving Direct Linear Transformation), the feature belt matching points are mapped to vertex matching points with more regular and uniform distribution. Divide the image into dense grids, each grid corresponds to a homography matrix projection transformation, and then apply the obtained homography matrix to find the matching points of the grid vertices in another image, the vertex matching points are evenly distributed, Using it as a matching point in grid optimization can effectively reduce the amount of calculation;

6. Introduce three constraints based on grid optimization: coordinate alignment constraints for overlapping regions, local similarity constraints for non-overlapping regions, and global similarity constraints for global structural consistency. By minimizing the error function To obtain the coordinate position of the pixel point on the panorama.

In summary, the image registration of the test image sequence A is completed, and the obtained registration results and the final mosaic panorama are shown in Figures 2 and 3 .

references:

【1】Oriented concise feature extraction algorithm, reference: E.Rublee, V.Rabaud, K.Konolige, etal.ORB: An efficient alternative to SIFTor SURF[C].International Conference on Computer Vision,2011:2564-2571.

【2】Fast corner detection algorithm, reference: E.Rosten, T.Drummond. Machine learning for high-speed corner detection [C]. European Conference on computer Vision, 2006: 430-443.

【3】Concise descriptors, references: M.Calonder, V.Lepetit, C.Strecha, et al.BRIEF: Binary Robust Independent Elementary Features[C].European Conference on Computer Vision.2010:778-792.

【4】Random sampling consensus algorithm, reference: M.A.Fischler, Bolles R.C..Random sample consensus: a paradigm for model fitting with applications to image analysis and automatic cartography[J].Communications of the ACM,1981,24(6):381 -395.

【5】For details of similar transformation, please refer to: T.Igarashi, Y.Igarashi. Implementing as-rigid-as-possible shape manipulation and surface flattening[J].Journal ofGraphics Gpu&Game Tools,2009.

【6】Abbreviation of Line Segment Detector, a line feature detection algorithm, for details, please refer to: R.G.V.Gioi, J.Jakubowicz, J.M.Morel, er al.LSD: a line segment detector[J].ImageProcessing On Line,2012,2 :35-55.

Claims (1)

1. The grid optimization-based image registration method for panoramic image stitching is characterized by comprising the following specific steps of:

(1) Performing rapid feature extraction by using ORB; performing feature detection by using a FAST algorithm, and then generating related descriptors of feature points based on the improved BRIEF descriptor, wherein the related descriptors comprise scale information, position and direction information of the feature points;

(2) Coarse matching of feature points is carried out by adopting a K-D tree and an optimal node priority algorithm, then the obtained matching points are subjected to ratio test by using a formula (1), wherein p is the current feature point and is defined as the nearest neighbor feature point p of p _best-closed And next adjacent feature point p _{second-closed} The ratio of hamming distances of (a) is less than a threshold ratio; applying the formula (2) to perform cross test; traversing the characteristic points in the image I, searching the matched points in the image J, and marking the points as M _I-＞J Then traversing the characteristic points in the image J to find the corresponding points in the image I, which are marked as M _J-＞I Cross testing considers that only when the two correspond to each other, the parties are a pair of correct matches;

M _I&J ＝M _I-＞J ∩M _J-＞I (2)

M _I&J indicating that image I and image J are the correct matches;

(3) Then, carrying out multilayer RANSAC screening on the matched point pairs after the fine matching, and screening the inner point sets of the characteristic point pairs in a plurality of planes of the image, so that the number of the final inner point sets accounts for more than 80% of the total matched pairs, and reserving the matching information to the maximum extent;

(4) Mapping the feature matching points into vertex matching points with more regular and uniform distribution through mobile direct linear conversion; dividing the image into dense grids, each grid corresponding to a homography matrix projective transformation, as shown in equation (3), wherein

Initial coordinates (x, y, 1) representing a certain matching point in the grid, are->

The transformed coordinates (x ', y', 1) are all three-dimensional, so the dimension of the matrix H is 3 x 3;

let H be:

the formula (3) is formed as follows:

the data are arranged into a matrix multiplication format, namely:

a _i h＝0 (6)

wherein a is _i Represents a 2X 9 matrix formed by the ith pair of matching points, H is the column direction of HVolume format, dimension 9×1:

taking all M pairs of matching points into account, as shown in equation (8), by minimizing the squared error e _k Solving a transformation matrix h, wherein A is a transformation matrix corresponding to M pairs of matching points, and the dimension is 2M multiplied by 9; w (W) _k For each pair of matching points, the dimension is 2M×2M, each element

Calculated by equation (4); />

The influence weight of the ith pair of matching points on the kth network is determined by the distance between the matching points and the grid center; x in formula (4) _k Represents the center coordinates, x, of the kth grid _i Representing the coordinates of the ith pair of matching points, wherein mu is an adjustment parameter;

e _k ＝argmin||W _k Ah|| ² s.t.||h||＝1 (8)

then, finding matching points of grid vertexes in another image by using the homography matrix obtained by calculation, and uniformly distributing the matching points of the vertexes to serve as the matching points in grid optimization;

(5) Three constraint terms based on grid optimization are introduced: aligning constraint terms for coordinates of overlapping regions, local similarity constraint terms for non-overlapping regions, and global similarity constraint terms for global structural consistency; firstly, setting a plurality of marks, storing the matching relation between images into a set T, and setting a set of matching point pairs obtained by mapping an image i and an image j in the step (4) as M _ij Since the image is divided into grids, set V _i And E is _i Vertex sets and edge sets for the mesh in image i;

the coordinate alignment constraint term is shown in a formula (10), so that the coordinates of the grid optimized matching points are ensured to be as consistent as possible, and the alignment error of the overlapping area of the adjacent images is reduced; wherein m (p) returns a matching point of the feature point p on the other image, and ψ (p) represents the position of the feature point p with a linear combination of 4 mesh vertex coordinates;

the local similarity constraint term is shown in a formula (11), so that the length and the direction of the same vector edge before and after grid optimization are not greatly changed; since the projection matrix is mainly suitable for the overlapping area, the similar transformation matrix is introduced in the very overlapping area

Representing edge->

Is calculated as equation (12), wherein +.>

And->

Linear combination of vertex variables, mainly modeling rotation and size transformation; />

And->

Representing the positions of two vertexes of the same edge on the original image, v _j ^′i And v _k ′ ^k Representing the vertex position after grid optimization, E _i Representing a set of edges of a mesh;

the global similarity constraint term is shown in a formula (13) and aims to improve the consistency of the whole structure in the image sequence; wherein the method comprises the steps of

For the weight of each edge, gradually changing from an overlapping area to a non-overlapping area, wherein the farther the distance from the overlapping area is, the larger the weight is, and the formula (14) is defined; s is(s) _i The method comprises the steps of defining a ruler measure of an image i, and estimating camera parameters of the image i through a beam adjustment method to obtain the ruler measure; θ _i Defining the rotation quantity of the image i relative to a reference image, wherein the average value of included angles between line features detected according to LSD feature lines is the rotation angle between the two images;

where eta and lambda are the adjustment parameters,

is a shared edge->

Grid set, M of (2) ⁱ A combination of all meshes representing overlapping areas of image i; d (q) _k ,M ⁱ ) Is a function ofComputing set->

Inner grid q _k Distance to the overlap region; r is R _i And C _i Representing the number of rows and columns of the image i grid;

combining the three constraint terms, as shown in a formula (15), minimizing the formula to obtain coordinate values of each pixel point in the panorama after grid optimization, and completing image registration;

wherein, gamma is the adjustment coefficient of the local similarity constraint term.

Images