Remote sensing image region registration method based on complex gradient description

Abstract

The invention discloses a remote sensing image region registration method based on complex gradient description. The remote sensing image region registration method comprises the following steps: S1, geocoding an acquired remote sensing image; s2, determining positions of homonymy points on the to-be-matched image through rough matching transformation, and establishing a search area by taking the positions as centers; s3, gradient information in the SAR image is extracted by adopting an ROEWA operator, and gradient information in the optical image is extracted by adopting a Sobel operator; s4, a template area of the image is taken out, zero padding is carried out on the periphery of the template area, the template area and the search area are made to be consistent in size, and frequency domain cross-correlation operation is carried out on the template area and the search area; s5, eliminating abnormal values by using a random sample consistency algorithm; and S6, describing a conversion relation between the images by adopting an affine transformation model. The method not only can be used for registration of SAR / optical remote sensing images, but also can be used for registration of other multi-mode remote sensing images.

Description

Technical Field

[0001] The present invention relates to the technical field of image region registration, and more specifically, to a remote sensing image region registration method based on complex gradient description. Background Art

[0002] Image registration, as a key task, aims to align images captured at different times and by different sensors, and is crucial for multimodal remote sensing applications. There have been many studies on remote sensing image registration. Traditional registration methods can be roughly divided into two categories: feature-based registration algorithms and region-based registration algorithms.

[0003] Feature-based methods usually do not rely on the geographic reference information of the images as the initial constraint condition. These methods usually first detect significant features (such as point features, line features, and surface features) between the images, and then identify corresponding points by describing the detected features. However, although this method overcomes the radiation differences between optical images and SAR images to a certain extent, it is difficult to obtain high-precision matching performance for images with obvious geometric distortions.

[0004] Region-based remote sensing registration methods are divided into methods based on the spatial domain and the frequency domain according to the different image representation domains. In the spatial domain, the sum of squared differences (SSD), normalized cross correlation (NCC), and mutual information (MI) are the most representative similarity metrics. However, the above traditional region-based registration methods usually evaluate the similarity of the image intensity information to find the optimal geometric transformation parameters to maximize or minimize the selected similarity metric, thereby driving the registration process. However, these similarity metrics are more sensitive to the significant non-linear radiation differences commonly existing in multimodal remote sensing images and cannot provide robust registration results for multimodal remote sensing images.

[0005] Currently, the more mainstream method is the structure-based regional registration method. According to whether it is a pixel-by-pixel structural feature, the structure-based methods can be divided into sparse structure feature-based methods and dense structure feature-based methods. For structure-based methods, sparse structure feature-based methods only perform calculations in the designed sparse sampling grid, which means that they can only use similarity metrics in the spatial domain to evaluate similarity and identify corresponding points, while similarity metrics in the frequency domain are not applicable. Dense structure feature-based methods can use similarity metrics in both the spatial domain and the frequency domain for matching, showing obvious effectiveness and advantages in dealing with significant non-linear radiation differences between multi-modal remote sensing images, and can meet the requirements of many current applications. Although dense structure feature-based methods can effectively overcome significant non-linear radiation differences between images, this method of mapping the intensity information of images to high-dimensional spatial pixel-level feature representations through pixel-by-pixel descriptors to drive the registration process greatly consumes computing resources, not only has high requirements for hardware but also is difficult to meet the needs of real-time operation. Summary of the Invention

[0006] The purpose of the present invention is to provide a remote sensing image regional registration method based on complex gradient description to overcome the defects of the prior art.

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

[0008] A remote sensing image regional registration method based on complex gradient description, comprising the following steps:

[0009] S1. Perform geocoding on the acquired remote sensing images;

[0010] S2. Use the block Harris operator to extract feature points in the geocoded remote sensing optical image, calculate the Harris feature values at each pixel in the reference image, divide the reference image into n×n non-overlapping grid regions, extract the k points with larger feature values in each grid as feature points, and a total of k×n×n feature points will be obtained for the entire reference image. Determine the positions of corresponding points on the image to be matched through rough matching transformation and establish a search area centered on this position;

[0011] S3. Extract the gradient information in the geocoded remote sensing optical image using the ROEWA operator and the Sobel operator. Calculate the gradient gx in the x-axis direction and the gradient gy in the y-axis direction obtained by the calculation using complex number representation. Open a statistical window for the 3×3 neighborhood around the pixel of the geocoded remote sensing optical image, add the complex number gradients of the pixels, obtain the initial feature description at the pixel (x, y), perform 2D Gaussian filtering on the two-dimensional feature description corresponding to the complex number feature description. This regional feature description aggregates all gradient information into one dimension to obtain a refined description of the feature;

[0012] S4. Take out the template area of the image, pad zeros around the template area to make the size of the template area the same as that of the search area, and perform a frequency domain cross-correlation operation on the template area and the search area;

[0013] S5. Use the random sample consensus algorithm to eliminate outliers;

[0014] S6. Use an affine transformation model to describe the conversion relationship between images.

[0015] Furthermore, in step S1, an initial RPC model and ground elevation data are used during the geocoding process to perform rough positioning and terrain correction on the SAR image.

[0016] Furthermore, in step S2, k = 1 to 5, and the number of grids n is calculated according to the total number of points expected to be obtained.

[0017] Furthermore, in step S3, the gradient g x in the x-axis direction and the gradient g y in the y-axis direction calculated using complex number representation are calculated by the formula: g complex = g x + ig y .

[0018] Furthermore, in step S4, the formula for performing the frequency domain CC operation on the zero-padded template area I 1 and the search area S is:

[0019] CFG_cc(x,y) = FFT -1 (FFT(I 1 (x,y))·[FFT(S(x,y))] * )

[0020] where FFT -1 represents the inverse Fourier transform operation, FFT represents the Fourier transform, and [] * represents taking the conjugate. The position with the highest similarity between the two at the maximum value position of the image is the coordinate of the center of the template area in the search area

[0021] Compared with the prior art, the advantages of the present invention are as follows: The present invention can not only be used for the registration of SAR / optical remote sensing images, but also for the registration of other multi-modal remote sensing images. BRIEF DESCRIPTION OF THE DRAWINGS

[0022] In order to more clearly illustrate the technical solutions in the embodiments of the present invention or the prior art, the following will briefly introduce the drawings required for the description of the embodiments or the prior art. Obviously, the drawings in the following description are only some embodiments of the present invention. For those of ordinary skill in the art, without creative efforts, other drawings can be obtained based on these drawings.

[0023] Figure 1 It is a schematic flow chart of the remote sensing image region registration method based on complex gradient description of the present invention.

[0024] Figure 2 It is a schematic diagram of the principle of gradient extraction and region feature description in the present invention.

[0025] Figure 3 It is a schematic diagram of the principle of zero-padding for the template region in the present invention.

[0026] Figure 4 In the present invention, a is the geocorrected optical image, and b is the geocorrected SAR image.

[0027] Figure 5 It is the corrected checkerboard fusion image in the present invention. DETAILED DESCRIPTION OF THE EMBODIMENTS

[0028] The following will elaborate on the preferred embodiments of the present invention in conjunction with the drawings, so that the advantages and features of the present invention can be more easily understood by those skilled in the art, thereby making a clearer and more definite definition of the protection scope of the present invention.

[0029] Refer to Figure 1 As shown, this embodiment discloses a remote sensing image region registration method based on complex gradient description, including the following steps:

[0030] Step S1: Considering the imaging mechanism of the SAR sensor, perform geocoding on the acquired remote sensing image.

[0031] In this embodiment, during the geocoding process, the initial RPC model and ground elevation data are adopted to roughly locate and correct the terrain of the SAR image, eliminating significant differences in scale and rotation.

[0032] Step S2: Use the block Harris operator to extract feature points from the geocoded remote sensing optical image, calculate the Harris eigenvalue at each pixel in the reference image, divide the reference image into n×n non-overlapping grid regions, extract the k points with larger eigenvalues in each grid as feature points, and the entire reference image will obtain k×n×n feature points. Regarding the selection of the number of divided grids n and the number of feature points k in the grid, generally empirically set k to 1 to 5 points, and then calculate the number of grids n according to the total number of points expected. After completing the extraction of feature points in the reference image, determine the position of the corresponding points on the image to be matched through a rough matching transformation (usually both SAR images and optical images contain geographical information, and a rough transformation model can be established through geographical information) and establish a search area centered on this position.

[0033] Step S3: Use the ROEWA operator to extract the gradient information in the geocoded remote sensing optical image and the Sobel operator to extract the gradient information in the geocoded remote sensing optical image, and calculate the x-axis direction gradient g x and the y-axis direction gradient g y using complex number representation. The formula is: g complex = g x + ig y .

[0034] After obtaining the complex gradient description, open a statistical window for the 3×3 neighborhood around the pixel of the geocoded remote sensing optical image, add the complex gradients of these 9 pixels, obtain the initial feature description at the pixel (x, y), perform 2D Gaussian filtering on the two-dimensional feature description corresponding to this complex feature description. This regional feature description aggregates all gradient information into one dimension, obtains a fine description of the feature, and compared with HOG-like descriptions, the constructed feature dimension is lower and the computational resources consumed are lower. The specific principle is as Figure 2 shown.

[0035] Step S4: In this embodiment, normalized cross-correlation is a traditional similarity measurement method. It is necessary to slide the template region in the search region and calculate the similarity function at each position to obtain the similarity surface and then obtain the matching position. However, this method requires multiple iterative runs and has a high computational complexity. In order to achieve fast similarity measurement of feature descriptions, this embodiment proposes a frequency-domain cross-correlation similarity measurement method based on CFG, named CFG_cc similarity measurement.

[0036] First, take out the template region of the image, fill zeros around the template region to make the size of the template region the same as that of the search region, as Figure 3 shown, where I is the template region and I 1 is the template region after filling zeros.

[0037] For the padded template region I 1 and the search region S, the formula for performing frequency-domain CC operation is as follows:

[0038] CFG_cc(x,y) = FFT -1 (FFT(I 1 (x,y)) · [FFT(S(x,y))] * )

[0039] Where FFT -1 represents the inverse Fourier transform operation, FFT represents the Fourier transform, and * represents taking the conjugate. The position with the highest similarity between the two at the maximum value position of the image is the coordinate of the center of the template region in the search region

[0040] Since frequency-domain cross-correlation transfers the correlation operation in the spatial domain to the frequency domain and replaces spatial convolution with frequency-domain multiplication, the operation efficiency can be greatly improved. At the same time, since the proposed method uses complex numbers to describe gradient information, the image gradient information is completely retained, and fine similarity discrimination between the template region and the search region can be achieved.

[0041] Step S5: After completing the matching of all feature point template regions, several matching point pairs can be obtained. However, there may still be some incorrect matching points in the matching point pairs at this time, which need to be eliminated. In this embodiment, the Random Sample Consensus (RANSAC) algorithm is used to eliminate outliers.

[0042] Step S6: Since the image was geocoded in step S1, an affine transformation model is used to describe the transformation relationship between images.

[0043] In this embodiment, Table 1 shows the comparison between this embodiment and other methods. The quantitative results are evaluated by two criteria: Root Mean Square Error (RMSE) and Correct Matching Rate (CMR). Among them, the threshold of CMR is set to 2.5 pixels. The proposed method shows the lowest RMSE of 1.86. The RMSE of AWOG and HOPC is slightly higher than that of the proposed method, and CFOG performs the worst with an RMSE of 4.69. At the same time, the CMR of this embodiment is much higher than that of other methods. In addition, the running speed of the method in this embodiment is the fastest.

[0044] Table 1

[0045] This embodiment The proposed method HOPC CFOG AWOG RMSE(m) 1.996090 2.422635 4.694291 2.260054 CMR 0.843750 0.672897 0.572727 0.696078 Time(s) 0.374128 16.078963 2.038008 1.843254

[0046] In this embodiment, complex numbers are used to describe the gradient information of an image, and a planar gradient feature description, named complex gradients feature (CGF), is designed. This feature description can not only represent the gradient intuitively and comprehensively, but also has significant computational advantages. This is because the real part and the imaginary part of the complex number can be used to represent the gradient information in the x-axis and y-axis directions of the image respectively, and a single complex number can represent both the direction and the magnitude of the gradient at the same time, thus achieving a full and compact representation of the image gradient information. In addition, due to the two-dimensional planar structure represented by complex numbers, compared with other dense feature description methods, it not only reduces the amount of calculation, saves a large amount of hardware storage space, but also is more convenient for frequency domain calculation. In terms of similarity measurement, this embodiment designs a CGF_cc similarity measurement method. Based on the characteristics of the CGF feature descriptor, this method uses frequency domain cross-correlation as a similarity measurement means to achieve fast, accurate and robust matching of features.

[0047] The present invention can be used not only for the registration of SAR / optical remote sensing images, but also for the registration of other multi-modal remote sensing images.

[0048] Although the embodiments of the present invention are described in conjunction with the accompanying drawings, the patent owner can make various deformations or modifications within the scope of the appended claims. As long as it does not exceed the protection scope described in the claims of the present invention, it should be within the protection scope of the present invention.

Claims

1. A remote sensing image region registration method based on complex gradient description, characterized in that: The following steps are involved: S1. geocoding the acquired remote sensing images; S2. Use the block Harris operator to extract feature points in the geocoded remote sensing optical image, calculate the Harris eigenvalue at each pixel in the reference image, divide the reference image into n×n non-overlapping grid areas, extract k points with larger eigenvalues ​​in each grid as feature points, and the entire reference image will obtain k×n×n feature points. Determine the position of the same-name point on the image to be matched through coarse matching transformation and establish a search area centered on the position; S3, using ROEWA operator to extract the gradient information in the geocoded remote sensing optical image and Sobel operator to extract the gradient information in the geocoded remote sensing optical image, and using complex number representation to calculate the x-axis direction gradient g x and the y-axis gradient g y , a statistical window is opened in the 3×3 neighborhood around the pixel of the geocoded remote sensing optical image, the complex gradient of the pixel is added to obtain the initial feature description at the pixel (x, y), and a 2D Gaussian filter is performed on the two-dimensional feature description corresponding to the complex feature description. The regional feature description aggregates all gradient information into one dimension to obtain a detailed description of the feature; S4, taking out the template area of ​​the image, padding the template area with zeros around it to make the template area consistent with the search area in size, and performing frequency domain cross-correlation operation on the template area and the search area; S5, using random sample consistency algorithm to eliminate outliers; S6. Use affine transformation model to describe the transformation relationship between images.

2. The remote sensing image region registration method based on complex gradient description according to claim 1, characterized in that: In the step S1, the initial RPC model and ground elevation data are used in the geocoding process to perform rough positioning and terrain correction on the SAR image.

3. The remote sensing image region registration method based on complex gradient description according to claim 1, characterized in that: In step S2, k=1-5, and the number of grids n is calculated according to the total number of points expected to be obtained.

4. The remote sensing image region registration method based on complex gradient description according to claim 1, characterized in that: In step S3, the x-axis direction gradient g is calculated by using a complex number. x and the y-axis gradient g y The formula is: g complex =g x +ig y .

5. The remote sensing image region registration method based on complex gradient description according to claim 1, characterized in that: The formula for performing frequency domain CC operation on the zero-filled template area I1 and the search area S in step S4 is: CFG_cc(x,y)=FFT -1 (FFT(I1(x,y))·[FFT(S(x,y))] * ) Among them, FFT -1 stands for inverse Fourier transform operation, FFT stands for Fourier transform, [] * The position with the highest similarity between the maximum image position and the maximum image position is the coordinate of the center of the template area in the search area.

Images