Multispectral image registration method based on improved SIFT (Scale Invariant Feature Transform) features

Abstract

The invention relates to a multispectral image registration method, in particular to a multispectral image registration method based on improved SIFT (Scale Invariant Feature Transform) features. The method comprises the following steps: step 1, preprocessing multispectral remote sensing image data, wherein the multispectral remote sensing image data comprises a reference multispectral image and a to-be-registered multispectral image; step 2, extracting image gradient features based on the preprocessed multispectral remote sensing image data; extracting edge corner point features based on image gradient features; step 3, extracting an improved SIFT feature based on the image gradient feature and the edge corner feature; 4, performing homonymy point matching calculation on the reference multispectral image and the multispectral image to be registered based on the improved SIFT features to obtain a transformation matrix; and step 5, converting the multispectral image to be registered into the reference multispectral image based on the transformation matrix to obtain a final registration result. According to the method, the number of effective feature points in the image registration process can be increased, and rapid and accurate image registration is realized.

Description

Technical Field

[0001] This invention relates to a multispectral image registration method, specifically a multispectral image registration method based on improved SIFT features. Background Technology

[0002] With the continuous development of sensor technology, multispectral imaging technology has also made significant progress. Multispectral imaging technologies based on filters, gratings, and prisms have all achieved substantial advancements, enabling the acquisition of increasingly more multispectral data of the same scene from different sensors and angles. Image registration technology, on the other hand, is the process of aligning and superimposing images of the same scene under different conditions (such as different scenes, lighting, angles, times, and devices), ensuring that the two images are geometrically aligned, i.e., maintaining consistency in their spatial coordinates. This technology provides a data foundation for subsequent image stitching, fusion, and many other applications, and can be widely used in image processing, computer vision, and many other fields.

[0003] In the image registration process, one image is called the reference image and the other is called the image to be registered. The data of the two come from different detectors, but they are imaging the same scene. Therefore, there is complementary and redundant information. Registration can integrate image information from different sensors to obtain richer details of the target area and better interpret the scene area.

[0004] The registration process comprises four basic steps: feature detection, feature matching, transform model estimation, and image resampling and transformation. Feature detection: Detects salient or unique regions in both the original and target images, including closed boundary regions, edges, contours, line intersections, corners, etc. For ease of processing, these features can be represented by "representative points," such as centers, line endpoints, and feature points; these points are collectively referred to as control points. Feature matching: Establishes the correspondence between features in the reference and target images, using feature descriptors, similarity measures, and spatial relationships between features to determine the accuracy of registration. Transform model estimation: Estimates the type and parameters of the mapping model to align the reference and target images. The parameters of the mapping model are calculated by establishing the correspondence between features. Image resampling and transformation: Transforms the target image using the mapping model.

[0005] Image registration methods can be broadly categorized into two types: grayscale-based registration methods and feature-based registration methods. Grayscale-based registration methods primarily calculate the similarity of pixel values ​​between the reference image and the image to be registered to determine the transformation matrix from the image to the reference image. Common methods include those based on cross-correlation and mutual information, grayscale correlation, and gradient mutual information. These methods directly utilize image grayscale values, effectively improving the accuracy and robustness of image registration, but the overall process is computationally intensive. Feature-based registration methods achieve final registration by extracting texture and grayscale feature information between images. These include Hough transform-based registration methods, Harris detection-based registration methods, and LSD-based registration methods. These methods utilize extracted image features to establish matching relationships, significantly reducing the computational load and effectively shortening the computation time. However, they suffer from insufficient expressiveness in multispectral image feature extraction and difficulties in matching corresponding points in the images. Summary of the Invention

[0006] The purpose of this invention is to address the technical problems of insufficient expressiveness of multispectral image feature extraction and difficulty in matching corresponding points in existing registration methods, and to provide a multispectral image registration method based on improved SIFT features.

[0007] To achieve the above objectives, the technical solution adopted by the present invention is as follows:

[0008] A multispectral image registration method based on improved SIFT features, characterized by the following steps:

[0009] Step 1: Preprocess the multispectral remote sensing image data, which includes a reference multispectral image and a multispectral image to be registered.

[0010] Step 2: Based on the preprocessed multispectral remote sensing image data, extract image gradient features; then, based on the image gradient features, extract edge corner features.

[0011] Step 3: Extract improved SIFT features based on image gradient features and edge corner features;

[0012] Step 4: Based on the improved SIFT features, perform corresponding point matching calculations between the reference multispectral image and the multispectral image to be registered to obtain the transformation matrix;

[0013] Step 5: Transform the multispectral image to be registered into the reference multispectral image based on the transformation matrix to obtain the final registration result and complete the multispectral image registration.

[0014] Furthermore, in step 1, the preprocessing of the multispectral remote sensing image data specifically involves:

[0015] Y1, Multispectral Band Selection

[0016] Principal component analysis algorithm is used to extract the first N principal components from multispectral remote sensing image data as representative spectral bands; N≥3;

[0017] Y2, Data Normalization Operation

[0018] Data normalization is performed on all extracted representative spectral segments to normalize them to a set range, so as to facilitate subsequent feature extraction.

[0019] Furthermore, step 2 specifically involves:

[0020] A. Extracting image gradient features

[0021] A1. Using preprocessed multispectral remote sensing image data as the bottom layer of the Gaussian gradient pyramid.

[0022] A2. Extract the gradient magnitude map in Gaussian scale space.

[0023]

[0024] In the formula: σ is the scale parameter of the Gaussian scale space; and These are the Gaussian scale spatial horizontal and vertical components of scale σ, respectively;

[0025] A3. In the Gaussian gradient magnitude diagram Continue extracting gradient direction With amplitude Complete the construction of the Gaussian gradient pyramid for the current layer:

[0026]

[0027] In the formula: and Gaussian gradient magnitude plots Gradient components in the horizontal and vertical directions;

[0028] A4. Change the Gaussian spatial scale parameter σ according to the set value, and complete the construction of the remaining layers of Gaussian gradient pyramids in the same way as A1 to A3, until the number of iterations meets the condition: log2(min(width,height))+1; where width and height are the length and width of the original reference multispectral image and the multispectral image to be registered, respectively.

[0029] B. Extract the edge corner features of each layer in the Gaussian gradient pyramid using the following method.

[0030] B1. Extract the edge image of the multispectral remote sensing image data in the current layer Gaussian gradient pyramid;

[0031] B2. Based on the edge image, corner points are extracted to obtain the contour set C:

[0032]

[0033] In the formula: c j Let j represent the j-th contour in the contour set C; n is the contour c. j The number of pixels in N; S P is the number of contours in the contour set C; i j For contour c j Each pixel in the array, i = 1, 2, ..., n;

[0034] B3. First, extract each contour c from the contour set C. j All pixels Pj:

[0035]

[0036] In the formula: and These are the coordinates of the pixels;

[0037] Then, around each pixel Pj, extract two redundant pixels spaced 2 pixels apart on the left and right. and

[0038]

[0039] Where: G σ For Gaussian weighted kernels; λ L and λ R These are two redundant pixels, one on the left and one on the right.

[0040] Finally, each pixel Pj and its corresponding redundant pixels are... and The median vector of the triangle is used as the feature of the edge corners.

[0041] Furthermore, step 3 specifically involves:

[0042] Step 3.1: For each pixel of all edge corner features in the second to penultimate layers of the Gaussian gradient pyramid, calculate the maximum value in its surrounding 3×3 neighborhood and the 3×3 neighborhood in the two layers above and below, in the direction from bottom to top. These maximum values ​​are taken as the extreme points of each pixel. The extreme points contain magnitude, direction and scale information.

[0043] Step 3.2: Divide each layer of Gaussian gradient pyramid into M fan-shaped regions along the circumference. For all pixels in the 4×4 neighborhood around each extreme point, count their orientation information and count the orientation of all edge corner features corresponding to each pixel, i.e. the number of pixels falling into the M fan-shaped regions. This yields the final improved SIFT feature with a feature scale of 4×4×M.

[0044] Furthermore, step 4 specifically involves:

[0045] Step 4.1: Calculate the positional difference e between the improved SIFT features corresponding to the reference multispectral image and the multispectral image to be registered using the following formula. p (i) Scale difference e s (i) and angular difference e o (i):

[0046]

[0047] In the formula: (x i ,y i ) and (x' i ,y' i ) represent the coordinates of the i-th feature point in the reference image and the coordinates of the i-th feature point in the image to be registered, respectively; σ i and σ' i These are the scale information of the reference image and the image to be registered, respectively; θ i and θ i ′ represent the feature angles of the reference image and the image to be registered, respectively; μ represents the parameters of the transformation matrix T; abs represents the absolute value; T is the transformation matrix;

[0048] Step 4.2: The transformation matrix T can be obtained by calculating the following minimized objective function:

[0049] E(i)=λ1e p (i)+λ2e s (i)+λ3e o (i)

[0050] In the formula: λ k Let λ1 be the penalty coefficient, k = 1, 2, 3, and λ1 + λ2 + λ3 = 1.

[0051] Furthermore, in step 5, the multispectral image to be registered is transformed into the reference multispectral image using the following formula based on the transformation matrix:

[0052] DATA = T·DATA′;

[0053] In the formula: DATA′ is the original multispectral data, and DATA is the transformed multispectral data.

[0054] The beneficial effects of this invention are:

[0055] 1. This invention reduces the computational and storage requirements of the entire registration process by performing band selection and amplitude normalization on multispectral data. Furthermore, it normalizes spectral data captured by different detectors and multispectral cameras, significantly improving the speed and accuracy of the registration process. This solves the problem that existing image registration algorithms cannot address the high dimensionality and varying spectral bands in multispectral data, leading to high computational costs and poor registration results due to significant differences in the characteristics of different spectral bands.

[0056] 2. This invention first effectively enhances the texture, edge and other details of the image in the Gaussian gradient space, and then extracts the edge corner features of the image in this space. Finally, it realizes an improved SIFT feature extraction based on edge corner features and image gradient features, which increases the number of effective feature points in the image registration process and solves the problem of insufficient expressiveness of multispectral image feature extraction in existing registration methods.

[0057] 3. The position, scale, and orientation objective function model proposed in this invention effectively estimates the transformation matrix parameters of the image to be registered to the reference image, achieving fast and accurate image registration and solving the problem of difficult matching of corresponding points in existing registration methods. Attached Figure Description

[0058] Figure 1 This is a flowchart of an embodiment of the present invention;

[0059] Figure 2 In step 1 of this embodiment of the invention, aerial multispectral remote sensing image data of a certain area taken from different angles are selected. a is a reference multispectral image, and b is the multispectral image to be registered.

[0060] Figure 3 This is a schematic diagram of the edge corner feature extraction step in step 2 of this embodiment of the invention;

[0061] Figure 4 This is a schematic diagram of the edge corner features in step 2 of an embodiment of the present invention;

[0062] Figure 5 Here are the edge corner feature maps in step 2 of this embodiment of the invention: a is the edge corner feature map extracted based on the reference multispectral image, and b is the edge corner feature map extracted based on the multispectral image to be registered.

[0063] Figure 6 Here are the gradient feature maps in step 2 of this embodiment of the invention: a is the gradient feature map extracted based on the reference multispectral image, and b is the gradient feature map extracted based on the multispectral image to be registered.

[0064] Figure 7 This is a diagram showing the matching effect of corresponding points in step 4 of this embodiment of the invention;

[0065] Figure 8 This is a diagram showing the effect of transforming the reference multispectral image and the multispectral image to be registered after estimating the parameters of the transformation matrix for matching corresponding points in step 5 of this embodiment of the invention; a is the reference multispectral image, and b is the effect of transforming the multispectral image to be registered.

[0066] Figure 9 This is a registration result diagram of an embodiment of the present invention. Detailed Implementation

[0067] To make the objectives, advantages, and features of the present invention clearer, the following detailed description of a multispectral image registration method based on improved SIFT features, in conjunction with the accompanying drawings and specific embodiments, will further clarify the advantages and features of the present invention. The advantages and features of the present invention will become clearer according to the following specific embodiments.

[0068] See Figure 1 The specific implementation steps of the multispectral image registration method based on improved SIFT features in this embodiment are as follows:

[0069] Step 1, for example Figure 2 The multispectral remote sensing image data shown is preprocessed. The multispectral remote sensing image data includes a reference multispectral image and a multispectral image to be registered.

[0070] Specifically, Y1, multispectral band selection

[0071] Principal Component Analysis (PCA) is used to extract the first three principal components from the multispectral remote sensing image data as representative spectral segments. These representative segments effectively reduce the computational load of subsequent processing while maintaining registration accuracy. In other embodiments of this invention, the first N principal components can be extracted as representative spectral segments, where N ≥ 3, as needed.

[0072] Y2, Data Normalization Operation

[0073] Data normalization was performed on all extracted representative spectral segments to normalize them to the range of [0, 255] to facilitate subsequent feature extraction.

[0074] Step 2, see Figure 6 Based on preprocessed multispectral remote sensing image data, image gradient features are extracted; see [link / reference]. Figures 3 to 5 Then, based on the image gradient features, the edge corner features are extracted;

[0075] A. Extracting image gradient features

[0076] A1. Using preprocessed multispectral remote sensing image data as the bottom layer of the Gaussian gradient pyramid.

[0077] A2. Extracting the gradient magnitude map in Gaussian scale space using a Sobel filter.

[0078]

[0079] In the formula: σ is the scale parameter of the Gaussian scale space; and These are the Gaussian scale spatial horizontal and vertical components of scale σ, respectively;

[0080] A3. In the Gaussian gradient magnitude diagram Continue extracting gradient direction With amplitude Complete the construction of the Gaussian gradient pyramid for the current layer:

[0081]

[0082] In the formula: and Gaussian gradient magnitude plots Gradient components in the horizontal and vertical directions;

[0083] A4. Change the Gaussian spatial scale parameter σ according to the set value. In this embodiment, the Gaussian spatial scale parameter σ of the next layer is changed according to the rule that it is twice that of the previous layer. That is, the Gaussian spatial scale parameter σ of each layer of Gaussian gradient pyramid is σ, 2σ, 4σ... from the bottom layer to the top layer. Complete the construction of the remaining layers of Gaussian gradient pyramid in the same way as A1 to A3 until the number of iterations meets the condition: log2(min(width,height))+1; where width and height are the length and width of the original reference multispectral image and the multispectral image to be registered, respectively.

[0084] B. Extract the edge corner features of each layer in each Gaussian gradient pyramid using the following method: B1. Use the Canny corner detection method to extract the edge image of the multispectral remote sensing image data in the current layer of the Gaussian gradient pyramid.

[0085] B2. Using an ergodic search method, corner points are extracted from the edge image to obtain the contour set C:

[0086]

[0087] In the formula: c j Let j represent the j-th contour in the contour set C; n is the contour c. jThe number of pixels in N; S P is the number of contours in the contour set C; i j For contour c j Each pixel in the array, i = 1, 2, ..., n;

[0088] B3. First, extract each contour c from the contour set C. j All pixels Pj:

[0089]

[0090] In the formula: and These are the coordinates of the pixels;

[0091] Then, around each pixel Pj, extract two redundant pixels spaced 2 pixels apart on the left and right. and

[0092]

[0093] Where: G σ The purpose of using a Gaussian weighted kernel is to smooth the sampling at corner points; λ L and λ R These are two redundant pixels, one on the left and one on the right.

[0094] Finally, each pixel Pj and its corresponding redundant pixels are... and The median vector of the triangle is used as the feature of the edge corner points.

[0095] Step 3: Extract improved SIFT features based on image gradient features and edge corner features.

[0096] Step 3.1: Following the bottom-up direction, for each pixel of all edge corner features in the second to second-to-last layers of the Gaussian gradient pyramid, calculate the maximum value in its surrounding 3×3 neighborhood and the 3×3 neighborhood in the two layers above and below it. These values ​​are taken as the extreme points of each pixel, and the magnitude, direction, and scale information of the extreme points can be obtained.

[0097] Step 3.2: Divide each layer of Gaussian gradient pyramid into 18 fan-shaped regions along the circumference. For all pixels in the 4×4 neighborhood around each extreme point, count their orientation information and count the orientation of all edge corner features corresponding to each pixel, i.e. the number of pixels falling into the 18 fan-shaped regions, to obtain the final improved SIFT feature with a feature scale of 4×4×18.

[0098] Step 4, see Figure 7Based on the improved SIFT features, corresponding point matching calculations are performed between the reference multispectral image and the multispectral image to be registered to obtain the transformation matrix.

[0099] Each improved SIFT feature contains positional, scale, and angular information. Therefore, the positional difference e between the improved SIFT features corresponding to the reference multispectral image and the multispectral image to be registered is first calculated using the following formula. p (i) Scale difference e s (i) and angular difference e o (i):

[0100]

[0101] In the formula: (x i ,y i ) and (x' i ,y' i ) represent the coordinates of the i-th feature point in the reference image and the coordinates of the i-th feature point in the image to be registered, respectively; σ i and σ' i These are the scale information of the reference image and the image to be registered, respectively; θ i and θ i ′ represent the feature angles of the reference image and the image to be registered, respectively; μ represents the parameters of the transformation matrix T; abs represents the absolute value; T is the transformation matrix;

[0102] Secondly, the transformation matrix T can be obtained by calculating the following minimized objective function:

[0103] E(i)=λ1e p (i)+λ2e s (i)+λ3e o (i)

[0104] In the formula: λ k Let λ1 be the penalty coefficient, k = 1, 2, 3, and λ1 + λ2 + λ3 = 1.

[0105] Step 5, see Figure 8 The following formula is used to transform the multispectral image to be registered into the reference multispectral image based on the transformation matrix T, so as to obtain the final registration result. Figure 9 As shown, multispectral image registration is completed:

[0106] DATA=T·DATA′

[0107] In the formula: DATA′ is the original multispectral data, and DATA is the transformed multispectral data.

Claims

1. A multispectral image registration method based on improved SIFT features, characterized in that, Includes the following steps: Step 1: Preprocess the multispectral remote sensing image data, which includes a reference multispectral image and a multispectral image to be registered. Step 2: Based on the preprocessed multispectral remote sensing image data, extract image gradient features; then, based on the image gradient features, extract edge corner features. Step 3: Extract improved SIFT features based on image gradient features and edge corner features; Step 4: Based on the improved SIFT features, perform corresponding point matching calculations between the reference multispectral image and the multispectral image to be registered to obtain the transformation matrix; Step 5: Transform the multispectral image to be registered into the reference multispectral image based on the transformation matrix to obtain the final registration result and complete the multispectral image registration.

2. The multispectral image registration method based on improved SIFT features according to claim 1, characterized in that, In step 1, the preprocessing of the multispectral remote sensing image data specifically includes: Y1, Multispectral Band Selection Principal component analysis algorithm is used to extract the first N principal components from multispectral remote sensing image data as representative spectral bands; N≥3; Y2, Data Normalization Operation Data normalization is performed on all extracted representative spectral segments to normalize them to a set range, so as to facilitate subsequent feature extraction.

3. The multispectral image registration method based on improved SIFT features according to claim 2, characterized in that: The N=3 is set to a range of [0,255].

4. A multispectral image registration method based on improved SIFT features according to any one of claims 1-3, characterized in that, Step 2 is as follows: A. Extracting image gradient features A1. Using preprocessed multispectral remote sensing image data as the bottom layer of the Gaussian gradient pyramid. A2, at the bottom layer Based on this, extract the gradient magnitude map in the Gaussian scale space. In the formula: σ is the scale parameter of the Gaussian scale space; and These are the Gaussian scale spatial horizontal and vertical components of scale σ, respectively; A3. In the Gaussian gradient magnitude diagram Continue extracting gradient direction With amplitude Complete the construction of the Gaussian gradient pyramid for the current layer: In the formula: and Gaussian gradient magnitude plots Gradient components in the horizontal and vertical directions; A4. Change the Gaussian spatial scale parameter σ according to the set value, and complete the construction of the remaining layers of Gaussian gradient pyramids in the same way as A1 to A3, until the number of iterations meets the condition: log2(min(width,height))+1; where width and height are the length and width of the original reference multispectral image and the multispectral image to be registered, respectively. B. Extract the corresponding edge corner features in each layer of the Gaussian gradient pyramid using the following method. B1. Extract the edge image of the multispectral remote sensing image data in the current layer Gaussian gradient pyramid; B2. Based on the edge image, corner points are extracted to obtain the contour set C: In the formula: c j Let j represent the j-th contour in the contour set C; n is the contour c. j The number of pixels in N; S The number of contours in the contour set C; For contour c j Each pixel in the array, i = 1, 2, ..., n; B3. First, extract each contour c from the contour set C. j All pixels P j : In the formula: and These are the coordinates of the pixels; Then at each pixel P j Extract two redundant pixels from the surrounding area, spaced 2 pixels apart. and Where: G σ For Gaussian weighted kernels; λ L and λ R These are two redundant pixels, one on the left and one on the right. Finally, each pixel P j and corresponding redundant pixels and The median vector of the triangle is used as the feature of the edge corners.

5. The multispectral image registration method based on improved SIFT features according to claim 4, characterized in that: In A2, the gradient magnitude map in the Gaussian scale space is extracted using a Sobel filter. In B1, the edge image of multispectral remote sensing image data is extracted using the Canny corner detection method. In B2, the corner points of the edge image are extracted using the ergodic search method to obtain the contour set C.

6. The multispectral image registration method based on improved SIFT features according to claim 4, characterized in that, Step 3 specifically involves: Step 3.1: For each pixel of all edge corner features in the second to penultimate layers of the Gaussian gradient pyramid, calculate the maximum value in its surrounding 3×3 neighborhood and the 3×3 neighborhood in the two layers above and below, in the direction from bottom to top. These maximum values ​​are taken as the extreme points of each pixel. The extreme points contain magnitude, direction and scale information. Step 3.2: Divide each layer of Gaussian gradient pyramid into M fan-shaped regions along the circumference. For all pixels in the 4×4 neighborhood around each extreme point, count their orientation information and count the orientation of all edge corner features corresponding to each pixel, i.e. the number of pixels falling into the M fan-shaped regions. This yields the final improved SIFT feature with a feature scale of 4×4×M.

7. The multispectral image registration method based on improved SIFT features according to claim 6, characterized in that: The value of M is 18.

8. The multispectral image registration method based on improved SIFT features according to claim 1, characterized in that, Step 4 specifically involves: Step 4.1: Calculate the positional difference e between the improved SIFT features corresponding to the reference multispectral image and the multispectral image to be registered using the following formula. p (i) Scale difference e s (i) and angular difference e o (i): In the formula: (x i ,y i ) and (x' i ,y' i ) represent the coordinates of the i-th feature point in the reference image and the coordinates of the i-th feature point in the image to be registered, respectively; σ i and σ' i These are the scale information for the reference image and the image to be registered, respectively; θ i and θ i ′ represent the feature angles of the reference image and the image to be registered, respectively; μ represents the parameters of the transformation matrix T; abs represents the absolute value; T is the transformation matrix; Step 4.2: The transformation matrix T can be obtained by calculating the following minimized objective function: E(i)=λ1e p (i)+λ2e s (i)+λ3e o (i) In the formula: λ k Let λ1 be the penalty coefficient, k = 1, 2, 3, and λ1 + λ2 + λ3 = 1.

9. A multispectral image registration method based on improved SIFT features according to claim 8, characterized in that: In step 5, the multispectral image to be registered is transformed into the reference multispectral image using the following formula based on the transformation matrix: DATA = T·DATA′; In the formula: DATA′ is the original multispectral data, and DATA is the transformed multispectral data.

Images