A multi-modal remote sensing image registration method based on multi-scale template matching

Abstract

This invention belongs to the field of remote sensing image registration technology, specifically relating to a multimodal remote sensing image registration method based on multi-scale template matching. First, multimodal remote sensing image pairs are acquired, and a feature description map of each image is constructed. Multiple layer-by-layer downsampling operations are then performed on these maps to construct Gaussian pyramids for the two images. Next, feature points of the image pairs are extracted, and the Gaussian pyramids of the images to be registered are rotated in different directions. Based on the extracted feature points, a template matching algorithm is used to match the remote sensing image pairs in different rotation directions, selecting the direction with the most matching points as the principal direction. Finally, the spatial transformation relationship of the remote sensing image pairs is determined based on the feature matching results corresponding to the principal direction. This invention can reliably handle general multi-source remote sensing image registration problems, solving the problem that existing methods cannot effectively handle images with severe nonlinear radiation distortion, and exhibits strong robustness and reliability.

Description

Technical Field

[0001] This invention belongs to the field of remote sensing image registration technology, specifically relating to a multimodal remote sensing image registration method based on multi-scale template matching. Background Technology

[0002] With the rapid development of aerospace and sensor technologies, remote sensing methods have become increasingly diversified, enabling the acquisition of diverse remote sensing data, such as optical, depth information, LiDAR, infrared, synthetic aperture radar (SAR), and topographic data. Combining multiple data sources to form multi-source remote sensing data is beneficial for better applications in change detection, target localization, and environmental and disaster monitoring. Multimodal image matching is key to achieving complementary advantages and leveraging the strengths of multi-source data.

[0003] As shown in Figures 1(a) and 1(b), there may be differences in displacement, rotation, scale, radiometrics, and noise between remote sensing image data from different sources. Therefore, it is necessary to achieve multi-distance remote sensing image registration. Registration methods mainly include three categories: deep learning-based methods, feature-based methods, and template-based methods.

[0004] Deep learning-based methods are currently a highly advanced technology, and many related works have achieved good results on corresponding datasets, but they are difficult to apply in real-world scenarios. On the one hand, there is currently no universal multimodal image dataset for model training, making it impossible to train reliable and stable parameters; on the other hand, well-trained models have poor adaptability to different real-world data and insufficient generalization ability. However, in recent years, with the rapid development of artificial intelligence, deep learning-based convolutional neural networks (CNNs) have developed rapidly in the field of computer vision, proposing algorithms such as LIFT, SuperPoint, DELF, and D2-Net, which also have some adaptability to nonlinear radiation distortion.

[0005] The specific process of feature-based methods involves first extracting image features, and then matching based on the similarity between feature descriptors of image feature points. The extracted features should be robust, stable, and repeatable; these features can be points, lines, or surfaces. Specifically, some researchers have proposed the PSO-SIFT algorithm, which offers an improved method for calculating image gradients, enhancing robustness to nonlinear radiometric distortion. Others have proposed the RIFT algorithm, which employs phase consistency and innovatively introduces the maximum index map descriptor (MIM). It first extracts feature points from the phase-consistent map and then uses the MIM for feature matching, achieving better results. However, due to the unstable localization accuracy of the extracted feature points, the matching accuracy of feature-based methods is generally lower than that of template-based methods.

[0006] Template-based methods match by comparing the similarity metrics of two selected regions. This is a commonly used image processing method in the field of computer vision. The main process is to determine the matching region in the target image that is most similar to the template image. The goal is to find corresponding image points in two or more images containing the same scene, and then find the transformation relationship between the images to achieve the alignment of the reference frames of the two images. Specifically, it mainly includes the following three steps: (1) First, design a search strategy to determine the candidate matching window regions corresponding to all template images; (2) Select the feature vector of the corresponding window region on the designed feature description map; (3) Calculate the similarity between the feature vector of the template window and the feature vector of the candidate window using the designed similarity metric method. Finally, the candidate window with the highest score is regarded as the correct matching region. At present, template matching plays an important role in face recognition, visual positioning and other fields. However, no matter how good the template-based method is, it has extremely high requirements for the initial position. For example, CFOG uses the points calculated by the RPC parameters of the satellite as the initial matching reference control points. However, at present, many different sources of data lack spatial mapping relationships. If the position deviation of the initial reference control point is too large or there is no initial reference control point, the matching performance will drop sharply. Its versatility in practice is reduced because it is sensitive to geometric differences and its speed is limited by the size of the selected template.

[0007] In summary, because the same ground feature may be located at different locations, and because different sensors can cause significant differences in the appearance of the same ground feature, there is severe nonlinear radiometric distortion (NRD). These significant differences pose serious challenges to the existing image matching algorithms described above, leading to a substantial decrease in the matching performance of most algorithms and making it difficult to meet the requirements of multi-source data fusion. Therefore, it is crucial to research a reliable multimodal image matching algorithm to provide accurate multimodal observation data for the same scene. Summary of the Invention

[0008] The purpose of this invention is to provide a multimodal remote sensing image registration method based on multi-scale template matching, in order to solve the problem that existing multimodal image registration methods cannot effectively handle multimodal images with severe nonlinear radiation distortion.

[0009] To address the aforementioned technical problems, this invention provides a multimodal remote sensing image registration method based on multi-scale template matching, comprising the following steps:

[0010] 1) Obtain multimodal remote sensing image pairs, namely a reference image and an image to be registered; construct a feature description map for each image, and perform multiple layer-by-layer downsampling operations on it to construct a Gaussian pyramid for the two images;

[0011] 2) Extract feature points from the image pairs and rotate the Gaussian pyramid of the images to be registered in different directions; based on the extracted feature points, use a template matching algorithm to match the remote sensing image pairs in different rotation directions, and select the direction with the most matching points as the main direction.

[0012] 3) Determine the spatial transformation relationship of the remote sensing image pair based on the feature matching results corresponding to the main direction.

[0013] The beneficial effects of the above technical solution are as follows: Addressing the problem that traditional image registration algorithms cannot effectively handle multimodal images with severe nonlinear radiometric distortion, this invention first constructs a Gaussian pyramid for the multimodal remote sensing image pair. Then, it rotates the Gaussian pyramid of the image to be registered in different directions, matching the remote sensing image pair in different rotation directions to find the optimal direction as the principal direction. Based on the feature matching results corresponding to the principal direction, it determines the spatial transformation relationship of the remote sensing image pair, achieving robustness to rotated images. This enables automatic registration even with poor multimodal image quality and severe geometrical differences and nonlinear radiometric distortion. Thus, while maintaining the high accuracy of the template matching method, it overcomes the drawback of poor adaptability to geometrical differences, automatically completing multimodal remote sensing image registration under most conditions. Furthermore, in the specific processing, a feature description map is first generated for the original image, and then it is downsampled layer by layer to obtain a series of feature description maps at different resolutions. Compared to generating feature description maps separately for images at different scales, this significantly improves computational efficiency. Overall, the method of the present invention can reliably handle general multi-source remote sensing image registration problems, including images from different sensors, at different times, from different perspectives, and at different resolutions, and has strong robustness and reliability.

[0014] Preferably, the constructed feature description map is a frequency domain convolutional map. The process of constructing the frequency domain convolutional map includes: ① performing wavelet decomposition on the input image to obtain low-frequency wavelet components and high-frequency wavelet components on the horizontal, vertical, and diagonal lines; ② using the Log-Gabor function to obtain the amplitude components of the low-frequency wavelet components at different scales, which are used as features of the low-frequency wavelet components at different scales; ③ normalizing the features at different scales and reconstructing the low-frequency wavelet components using the normalization results; ④ performing wavelet reconstruction based on the reconstructed low-frequency wavelet components and the high-frequency wavelet components on the horizontal, vertical, and diagonal lines from step ① to obtain the frequency domain convolutional map.

[0015] The beneficial effects of the above technical solution are: using frequency domain convolutional maps can more accurately and reliably describe the feature information of multimodal remote sensing images, and significantly reduce nonlinear radiation distortion between multimodal images.

[0016] Preferably, the feature points extracted in step 2) are Shi-Tomasi corner points.

[0017] The beneficial effects of the above technical solution are as follows: Shi-Tomasi corner point extraction is a relatively stable and computationally simple corner point extraction method. Based on this method, the extracted feature points are concentrated in the areas of the image with rich structural information, thereby improving the registration accuracy.

[0018] Preferably, the extracted Shi-Tomasi corner points are Shi-Tomasi corner points after distribution homogenization processing, and the distribution homogenization processing is performed using an adaptive non-maximum suppression algorithm. The formula for calculating the threshold for the variable length of the suppression region is:

[0019]

[0020]

[0021] In the formula, a j The threshold for suppressing variable region length; m is the number of feature points to be extracted; H I and W I These represent the height and width of the image, respectively.

[0022] The beneficial effects of the above technical solution are as follows: considering that Shi-Tomasi corner points are concentrated in areas with rich structural information in the image, the distribution uniformization of the corner point extraction results helps to reduce the local error of the image transformation model and improve the accuracy of image matching.

[0023] Preferably, in step 2), the template matching algorithm requires multi-scale matching using the following method: performing template matching between a certain level of the Gaussian pyramid of the reference image and each level of the Gaussian pyramid of the image to be registered, in order to find the hierarchical correspondence between the Gaussian pyramid of the reference image and the Gaussian pyramid of the image to be registered; and then finding the matching search area in the template matching algorithm process based on the hierarchical correspondence.

[0024] The beneficial effects of the above technical solution are: the above processing can ensure that the correct matching search area is found while ensuring processing efficiency.

[0025] Preferably, in step 3), outlier feature points in the feature matching results need to be eliminated, and the spatial transformation relationship of the remote sensing image pairs needs to be determined based on the feature matching results after the outlier feature points are eliminated.

[0026] The beneficial effect of the above technical solution is that removing outlier feature points can improve registration accuracy.

[0027] Preferably, the MAGSAC++ method is used for outlier feature point elimination.

[0028] The beneficial effect of the above technical solution is that the MAGSAC++ method can improve the elimination accuracy.

[0029] Preferably, in step 3), when determining the spatial transformation relationship of the remote sensing image pair, it is also necessary to combine the feature matching results of the adjacent directions of the main direction.

[0030] The beneficial effects of the above technical solution are as follows: considering that the rotation difference of the matched image is not strictly in the main direction in most cases, and the correct match may occur in the main direction and the front and back directions, the feature matching results in the main direction and the feature matching results in the direction adjacent to the main direction are merged to form the best matching result.

[0031] Preferably, the similarity measurement method used in the template matching algorithm is as follows: the SSD between the template window corresponding to a certain feature description map included in the Gaussian pyramid of the reference image and the image to be registered is determined using the following formula:

[0032]

[0033] In the formula, p represents the position of a pixel in the feature description map; T i (s) is the SSD similarity metric function between the template window on the feature descriptor map M1 of the reference image and the candidate matching window after displacement s at the corresponding position on the feature descriptor map M2 of the image to be registered; the matching function is defined as:

[0034]

[0035] In the formula, s i F represents the defined matching function; 2D and These represent forward and backward 2D FFTs, respectively. Describe the complex conjugate of the 2DFFT; determine the minimum T. i (s) to determine the best match between windows on M1 and M2.

[0036] The beneficial effects of the above technical solution are: the similarity measurement method described above can improve the matching speed.

[0037] Preferably, the multimodal remote sensing image pair in step 1) is a denoised multimodal remote sensing image pair.

[0038] The beneficial effects of the above technical solution are: noise reduction processing can improve registration accuracy. Attached Figure Description

[0039] Figure 1(a) and Figure 1(b) are two remote sensing images from different sources;

[0040] Figure 2 This is a flowchart of the multimodal remote sensing image registration method based on multi-scale template matching of the present invention.

[0041] Figure 3(a) and Figure 3(b) show the effects before and after homogenization when the feature point is 800, respectively.

[0042] Figure 4 This is a diagram illustrating the generation framework of the frequency domain convolutional map in this invention;

[0043] Figure 5(1a) is the original optical image;

[0044] Figures 5(1b), 5(1c), 5(1d), and 5(1e) are the gradient plot, phase coherence plot, MIM plot, and FDCM plot of Figure 5(1a), respectively.

[0045] Figure 5(2a) is the original map image;

[0046] Figures 5(2b), 5(2c), 5(2d), and 5(2e) are the gradient plot, phase coherence plot, MIM plot, and FDCM plot of Figure 5(2a), respectively.

[0047] Figures 6(1a) and 6(1b) are optical and SAR images from the first set of comparative experiments, respectively.

[0048] Figures 6(1c), 6(1d), 6(1e), and 6(1f) are comparison results of the NCC, HOPC, CFOG, and FDCM similarity measurement experiments in Figures 6(1a) and 6(1b), respectively.

[0049] Figures 6(2a) and 6(2b) are optical and infrared images from the second set of comparative experiments, respectively.

[0050] Figures 6(2c), 6(2d), 6(2e), and 6(2f) are comparison results of the NCC, HOPC, CFOG, and FDCM similarity measurement experiments in Figures 6(2a) and 6(2b), respectively.

[0051] Figure 7 This is a scale-space pyramid diagram designed in the MSTM of this invention;

[0052] Figure 8 This is a diagram of the multi-scale template matching strategy in the MSTM of this invention;

[0053] Figures 9(a), 9(b), 9(c), 9(d), 9(e), and 9(f) represent the optical-depth dataset, optical-infrared dataset, optical-map dataset, optical-optical dataset, optical-SAR dataset, and day-night dataset, respectively.

[0054] Figures 10(a) and 10(b) show the matching success rate and matching time for different template sizes at each level, respectively.

[0055] Figure 11 It is the number of correct matches of the rotation-invariant MSTM algorithm on multiple types of rotated image pairs between -90° and 90°.

[0056] Figure 12 These are NCMs plots of MSTM under different types of multi-source remote sensing images (rotation angle from 0° to 359°);

[0057] Figures 13(a), 13(b), 13(c), 13(d), 13(e), 13(f), 13(g), and 13(h) show the MSTM matching performance results at rotation angles of 0°, 45°, 90°, 135°, 180°, 225°, 270°, and 315°, respectively.

[0058] Figures 14(a) and 14(b) show the performance results of MSTM at scales of 1:1.5 and 1:2, respectively.

[0059] Figures 15(1a), 15(1b), 15(1c), 15(1d), 15(1e), and 15(1f) are the matching results of image pairs selected in the first multimodal dataset using CFOG, Ms-HLMO, HAPCG, RIFT, SuperGlue, RIFT, and the MSTM of the present invention, respectively.

[0060] Figures 15(2a), 15(2b), 15(2c), 15(2d), 15(2e), and 15(2f) are the matching results of image pairs selected in the second multimodal dataset using CFOG, Ms-HLMO, HAPCG, RIFT, SuperGlue, RIFT, and the MSTM of the present invention, respectively.

[0061] Figures 15(3a), 15(3b), 15(3c), 15(3d), 15(3e), and 15(3f) are the matching results of image pairs selected in the third multimodal dataset using CFOG, Ms-HLMO, HAPCG, RIFT, SuperGlue, RIFT, and the MSTM of this invention, respectively.

[0062] Figures 15(4a), 15(4b), 15(4c), 15(4d), 15(4e), and 15(4f) are the matching results of image pairs selected in the fourth multimodal dataset using CFOG, Ms-HLMO, HAPCG, RIFT, SuperGlue, RIFT, and the MSTM of the present invention, respectively.

[0063] Figures 15(5a), 15(5b), 15(5c), 15(5d), 15(5e), ​​and 15(5f) are the matching results of image pairs selected in the fifth multimodal dataset using CFOG, Ms-HLMO, HAPCG, RIFT, SuperGlue, RIFT, and the MSTM of the present invention, respectively.

[0064] Figures 15(6a), 15(6b), 15(6c), 15(6d), 15(6e), and 15(6f) are the matching results of image pairs selected in the sixth multimodal dataset using CFOG, Ms-HLMO, HAPCG, RIFT, SuperGlue, RIFT, and the MSTM of the present invention, respectively.

[0065] Figures 16(a), 16(b), 16(c), 16(d), 16(e), and 16(f) are NCM diagrams for different datasets.

[0066] Figures 17(a), 17(b), 17(c), 17(d), 17(e), and 17(f) are RMSE plots for different datasets.

[0067] Figures 18(1a) and 18(1b) are respectively a checkerboard pattern display of the registration accuracy of optical-depth image pair 1 and a display of the image fusion effect;

[0068] Figures 18(1c) and 18(1d) are respectively a checkerboard pattern display of the registration accuracy of optical-depth image pair 2 and a display of the image fusion effect;

[0069] Figures 18(1e) and 18(1f) are respectively a checkerboard pattern display of the registration accuracy of optical-depth image pair 3 and an image fusion effect display;

[0070] Figures 18(1g) and 18(1h) are respectively a checkerboard pattern display of the registration accuracy of optical-depth image pair 4 and a display of the image fusion effect;

[0071] Figures 18(2a) and 18(2b) are respectively a checkerboard diagram showing the registration accuracy of optical-infrared image pair 1 and a diagram showing the image fusion effect;

[0072] Figures 18(2c) and 18(2d) are respectively a checkerboard pattern display of the registration accuracy of the optical-infrared image pair 2 and a display of the image fusion effect;

[0073] Figures 18(2e) and 18(2f) are respectively a checkerboard pattern display of the registration accuracy of optical-infrared image pair 3 and a display of the image fusion effect;

[0074] Figures 18(2g) and 18(2h) are respectively a checkerboard pattern display of the registration accuracy of optical-infrared image pair 4 and a display of the image fusion effect;

[0075] Figures 18(3a) and 18(3b) are respectively a checkerboard pattern display of the registration accuracy of optical-map image pair 1 and a display of the image fusion effect;

[0076] Figures 18(3c) and 18(3d) are respectively a checkerboard pattern display of the registration accuracy of optical-map image pair 2 and a display of the image fusion effect;

[0077] Figures 18(3e) and 18(3f) are respectively a checkerboard pattern display of the registration accuracy of optical-map image pair 3 and a display of the image fusion effect;

[0078] Figures 18(3g) and 18(3h) are respectively a checkerboard pattern display of the registration accuracy of optical-map image pair 4 and a display of the image fusion effect;

[0079] Figures 18(4a) and 18(4b) are respectively a checkerboard diagram showing the registration accuracy of optical-optical image pair 1 and a diagram showing the image fusion effect;

[0080] Figures 18(4c) and 18(4d) are respectively a checkerboard pattern display of the registration accuracy of optical-optical image pairs 2 and a display of the image fusion effect;

[0081] Figures 18(4e) and 18(4f) are respectively a checkerboard diagram showing the registration accuracy of optical-optical image pairs 3 and a diagram showing the image fusion effect;

[0082] Figures 18(4g) and 18(4h) are respectively a checkerboard diagram showing the registration accuracy of optical-optical image pairs 4 and a diagram showing the image fusion effect;

[0083] Figures 18(5a) and 18(5b) are respectively a checkerboard diagram showing the registration accuracy of optical-SAR image pair 1 and a diagram showing the image fusion effect;

[0084] Figures 18(5c) and 18(5d) are respectively a checkerboard diagram showing the registration accuracy of optical-SAR image pair 2 and a diagram showing the image fusion effect;

[0085] Figures 18(5e) and 18(5f) are respectively a checkerboard diagram showing the registration accuracy of optical-SAR image pair 3 and a diagram showing the image fusion effect;

[0086] Figures 18(5g) and 18(5h) are respectively a checkerboard display of the registration accuracy of optical-SAR image pair 4 and a display of the image fusion effect;

[0087] Figures 18(6a) and 18(6b) are respectively a checkerboard pattern display of the registration accuracy of day-night image pair 1 and a display of the image fusion effect;

[0088] Figures 18(6c) and 18(6d) are respectively a checkerboard pattern display of the registration accuracy of day-night image pairs 2 and an image fusion effect display;

[0089] Figures 18(6e) and 18(6f) are respectively a checkerboard pattern display of the registration accuracy of day-night image pairs 3 and an image fusion effect display;

[0090] Figures 18(6g) and 18(6h) are respectively a checkerboard pattern display of the registration accuracy of day-night image pairs 4 and an image fusion effect display. Detailed Implementation

[0091] To make the objectives, technical solutions, and advantages of the present invention clearer, the present invention will be further described in detail below with reference to the accompanying drawings and embodiments.

[0092] The framework of the multimodal remote sensing image registration method (referred to as MSTM) based on multi-scale template matching of the present invention is as follows: Figure 2 As shown, the overall approach is as follows: First, the input images are preprocessed, primarily for image denoising. Then, feature point detection is performed on the denoised images. Specifically, the Shi-Tomasi corner detection algorithm can be used to extract feature points from the multimodal image pairs, providing initial control points for subsequent template window determination. Next, the MSTM proposed in this invention is used to construct feature description maps and perform template matching under a multi-scale strategy, ultimately determining the accurate spatial transformation relationship between image pairs. Specifically, under the multi-scale strategy, a Gaussian pyramid is constructed to create the multi-scale space of the feature description map. The specific processing steps are as follows: First, construct the lowest-level frequency domain convolutional map (FDCM). Then, perform multiple downsampling operations on the FDCM, with the scale space from top to bottom. Next, extract the feature templates corresponding to the feature points of each FDCM layer, rotate the feature templates multiple times, and aggregate all rotated feature templates to generate an all-directional aggregated feature vector. Then, match the window with the highest similarity at the same position in the target image as the correct match, and perform outlier removal to remove incorrect matches. Thus, through a multi-scale strategy, find the correct corresponding scale level and accurately transform the matrix parameters level by level. Finally, determine the accurate spatial transformation relationship between image pairs based on the final calculated accurate transformation matrix parameters.

[0093] Specifically, the specific contributions of this invention include: (1) A multi-scale template matching method (MSTM) for multimodal image registration is designed to solve the multimodal image registration problem with severe nonlinear radiation distortion. MSTM can reliably handle general multi-source remote sensing image registration problems, including images from different sensors, different time phases, different viewpoints, and different resolutions, and has strong robustness and reliability. (2) A wavelet transform-based frequency domain convolutional graph (FDCM) is proposed, which can more accurately and reliably describe the feature information of multimodal images and significantly reduce the nonlinear radiation distortion between multimodal images; and a rotation-invariant omnidirectional aggregation feature vector is designed, which can complete the automatic registration work when the multimodal image quality is poor and there is severe nonlinear radiation distortion due to geometric differences. (3) To address the geometric difference problem of multimodal image data, a scale-space-based multi-scale template matching strategy is applied to MSTM, which solves the problem of finding the mapping relationship between the template window and the candidate matching window, and achieves good results in solving geometric differences and optimizing template matching speed.

[0094] The implementation process of the entire method will be described in more detail below.

[0095] Step 1: Obtain a multimodal image pair, which is referred to as the reference image and the image to be registered, respectively; and preprocess the obtained image pair, which mainly involves image denoising.

[0096] Step two involves constructing frequency domain convolutional maps (FDCM) of the reference image and the image to be registered as feature description maps, and then performing multiple layer-by-layer downsampling processes on these maps to construct Gaussian pyramids for the two images. Specifically:

[0097] Due to severe nonlinear radiative distortion among multimodal images, direct matching using image information is difficult to achieve good results. Therefore, most excellent multimodal image matching algorithms, such as HOPC, CFOG, and RIFT, describe the features of the images using certain methods to generate a feature map, increasing the similarity between multimodal images and greatly improving the matching performance of the algorithm. As mentioned above, compared with other description methods such as gradient information, phase-consistent maps (PCs) are more resistant to nonlinear radiative distortion among multimodal images, and most of the aforementioned excellent algorithms are based on PC models to construct feature maps. However, when significant nonlinear radiative distortion exists, the PC model may be affected, and the generated feature map may not be clear enough in terms of feature description or may lack some geometric structure. This invention constructs a new frequency domain convolutional map (FDCM) based on the PC model using wavelet transform. Wavelet transform can decompose an image into low-frequency and high-frequency components. The low-frequency components contain the overall geometric structure information of the image, while the high-frequency components contain the local detail information of the image. Most of the nonlinear radiative distortion in the image is contained in the local details of the image. The low-frequency image after wavelet transform is significantly less affected by nonlinear radiative distortion of local image details, thus greatly enhancing its robustness to nonlinear radiative distortion. The FDCM generation process is as follows: Figure 4 As shown.

[0098] First, the input image is decomposed using wavelet decomposition, calculated as follows:

[0099] [C N H s V s D s ] = WT N (I)

[0100] In the formula, I is the input image; N represents the scale in the wavelet decomposition process; s is the number of layers; C N This represents the low-frequency wavelet components after decomposition, corresponding to Figure 4 LL;H s V s D s These represent the high-frequency wavelet components on the horizontal, vertical, and diagonal lines, respectively, corresponding to... Figure 4 HL, LH, and HH in the text.

[0101] Good noise suppression and edge extraction are key to the superior performance of the feature descriptor map. Log-Gabor wavelets are used in the construction of the FDCM. They can describe the geometric structure information of an image at multiple scales and from multiple perspectives. Typically, the Log-Gabor wavelet is represented as follows:

[0102]

[0103] In the formula, (ρ,θ) represents logarithmic polar coordinates; s and o represent the scale and direction of the Log-Gabor wavelet, respectively; (ρ s ,θ so ) represents the center frequency of the Log-Gabor wavelet; σ ρ ,σ θ Let ρ and θ represent the bandwidths, respectively. The Log-Gabor wavelet can be decomposed into even-symmetric and odd-symmetric filters in the spatial domain, defined as follows:

[0104] L(x,y,s,o)=L even (x,y,s,o)+iL odd (x,y,s,o)

[0105] In the formula, the real part L even (x,y,s,o) denotes an even-symmetric wavelet; the imaginary part L odd (x,y,s,o) represents an odd-symmetric wavelet.

[0106] Then, for the low-frequency component C generated by wavelet decomposition... N , will C N With even-symmetric wavelet L even (x,y,s,o) and odd-symmetric wavelet L odd Convolution is performed on (x,y,s,o) to generate E respectively. so (x,y) and O so (x,y) is as follows:

[0107] E so (x,y)=C N (x,y)*L even (x,y,s,o)

[0108] O so (x,y)=C N (x,y)*L odd (x,y,s,o)

[0109] The amplitude of the image in the frequency domain is:

[0110]

[0111] The multi-scale Log-Gabor features in a specified direction are:

[0112]

[0113] Next, the feature map after normalization and dimensionality reduction of the multi-directional Log-Gabor features is used as the reconstructed low-frequency wavelet component:

[0114]

[0115]

[0116] In the formula, ε represents a very small constant.

[0117] Then, wavelet reconstruction is performed on the four sub-images (including the high-frequency wavelet components and the reconstructed low-frequency wavelet components) to generate the final frequency domain convolutional image:

[0118] f(x,y)=(IWT N (F,H s V s D s ))

[0119] Where f(x,y) is the FDCM; IWT N (·) represents the reconstruction function.

[0120] Visualization comparison of FDCM with other typical feature description maps, for example Figures 5(1a) to 5(2e) As shown in Figures 5(1a) and 5(2a), the original data consists of visible light-map image pairs, exhibiting significant nonlinear radiometric distortion. For comparison purposes, the two images have been manually registered, essentially eliminating geometric differences. The image gradient feature description maps generated using the following formula are shown in Figures 5(1b) and 5(2b). This is the method used by most algorithms. However, due to the influence of nonlinear radiometric distortion, the feature description is discontinuous and unstable, resulting in poor algorithm robustness.

[0121]

[0122]

[0123] The PC maps shown in Figures 5(1c) and 5(2c) can describe the geometric information of an image, but they are easily affected by noise, and edges may be blurred and ghosted. The MIM maps shown in Figures 5(1d) and 5(2d) are feature descriptor maps of the RIFT algorithm. Based on the PC map, they select the principal direction of local image features for description. Although they are more robust to nonlinear radiometric distortion, they lose the geometric structure information of the image and are less robust to some low-quality images. The FDCM maps shown in Figures 5(1e) and 5(2e) are not only very robust to nonlinear radiometric distortion, but also clearly describe the geometric structure information of the image, showing good performance on complex and low-quality images. Multi-scale local feature information is used in MSTM to describe the feature point matching window for template matching.

[0124] Additionally, it should be noted that generating FDCM feature descriptors for images at different scales is time-consuming, and the FDCM feature descriptors generated from the original images are not high-dimensional features. Therefore, if... Figure 7 As shown, the present invention directly downsamples the data layer by layer to obtain a series of feature description maps with different resolutions.

[0125] Step 3: Extract feature points from the feature pyramid of the image pair and perform distribution homogenization on the extracted feature points.

[0126] Specifically:

[0127] 1) Shi-Tomasi corner extraction is performed as feature points. Harris corner detection is one of the most stable corner extraction algorithms, and the Shi-Tomasi corner extraction algorithm is an improved version of the Harris corner detector. The corner response value calculation of this algorithm is simpler, and it can achieve better results than the Harris operator in most cases. The Shi-Tomasi corner response value for each image is calculated as follows:

[0128] corner = min(λ1, λ2)

[0129]

[0130] In the formula, λ1 and λ2 are two eigenvalues ​​of matrix M; G x G y σ represents the gradient values ​​of the image along the x and y directions, respectively; W is a Gaussian sliding window with variance σ. Pixels with corner response values ​​greater than a threshold are considered reliable feature points of the multimodal image.

[0131] 2) Perform feature point distribution homogenization. Previous research suggests that uniformly distributed feature points help reduce local errors in image transformation models and improve image matching accuracy. However, the Shi-Tomasi corner detector extracts feature points based on image intensity, resulting in most extracted feature points clustering in regions rich in structural information. To overcome the problem of uneven feature point distribution, an efficient adaptive non-maximum suppression algorithm (ANMS) is employed. This algorithm quickly and uniformly divides the extracted feature points, avoiding block-based feature point extraction and significantly improving algorithm performance. The suppression region side length threshold 'a' is used. h The calculation is as follows:

[0132]

[0133]

[0134] In the formula, H I and W IHere, represents the height and width of the image, respectively; m represents the number of feature points to be extracted. Therefore, when the input image size changes, the feature points can be re-divided according to the image size and the number of feature points to be extracted, making the feature point distribution more uniform, as shown in Figures 3(a) and 3(b).

[0135] Step four: Based on the extracted feature points, template matching is used to match image pairs and determine their spatial transformation relationship. The specific implementation steps are as follows:

[0136] 1) Extract multi-directional feature description maps from the images to be registered in each image pair, and aggregate these feature description maps to construct an omnidirectional aggregated feature description vector for matching windows adapted to different directions. Specifically:

[0137] To extract feature description vectors for candidate matching windows in different directions, it is necessary to rotate the series of feature description maps (FDCM) generated from the original image n times (actually rotating the Gaussian pyramid composed of a series of feature description maps) to achieve full-directional coverage and construct directional feature description maps. For the r-th direction, the rotation transformation formula is:

[0138]

[0139] In the formula, Rot represents the rotation operation on the image; FDCM 0 is the original, unrotated feature map; m is the rotation angle of the r-th feature map.

[0140] Directional feature description map FDCM obtained by rotation r Keypoints can be obtained using the Shi-Tomasi algorithm. Feature Description Map (FDCM) r The key point at position (i,j) can be represented as The directional feature description vector of the corresponding candidate matching window can be represented as: Image FDCM r The set of all feature point orientation feature description vectors in the set can be described as This is called the r-direction feature description vector of the image. In the aggregation model, the feature description map FDCM aggregates the feature keypoints and the omnidirectional feature description vectors of the matching window.

[0141] 2) Based on the extracted feature points, a template matching algorithm is used to match remote sensing image pairs in different rotation directions. The direction with the best matching effect is selected as the main direction, and the feature matching results of the main direction and the feature matching results of the directions adjacent to the main direction are merged to form the best matching result. Specifically, this includes Steps 1 to 4.

[0142] Step 1: Perform multi-scale matching on the reference image and the image to be registered to find the correct matching search region corresponding to the feature points on the reference image in the image to be registered. Specifically, first, perform template matching between the highest level of the reference image and all levels of the image to be registered to find the correct corresponding level relationship. Based on this, complete the entire multi-scale matching by following the conventional matching process. The complete matching process is as follows: Figure 8 As shown.

[0143] The biggest challenge of template matching-based methods lies in finding the correct template matching search region corresponding to a feature point. For a feature point in the reference image, a full-image search can be performed in the image to be matched. While this approach might find the correct search region, it is extremely time-consuming. This is especially true when dealing with a large number of feature points, where the search time increases exponentially, significantly degrading matching performance. Furthermore, when dealing with two images with scale differences, even if the template matching windows correspond correctly, inconsistencies in the local feature descriptions within the template windows can lead to unreliable matching results or even matching failures. To address these two key issues, eliminate dependence on initial control points, and achieve scale robustness, the multi-scale template matching strategy mentioned in this step was designed.

[0144] Furthermore, after establishing the feature description scale space, template matching on low-level feature maps using templates with large windows is time-consuming, while using templates with small windows makes it difficult to find the correct corresponding region. However, using template matching with larger windows on high-level feature maps can shorten the matching time and improve the matching success rate. Although the location of the correctly matched point may not be very precise, after outlier removal using MAGSAC++, a relatively reliable transformation model can still be estimated using the remaining points. This model is then passed to the next level for template matching. Because of the initial transformation model, the template matching search window can be appropriately reduced, further shortening the search time and improving the matching accuracy. This process continues until the lowest level is reached to complete the entire multi-scale matching. When scale scaling exists, the levels of the reference image and the image to be matched may be inconsistent; therefore, the multi-scale matching method described in this step is used.

[0145] The feature description scale space of this invention is similar to the Gaussian scale space, the most classic and widely used in SIFT. The biggest difference is that SIFT provides scale attribute information for feature points based on the scale space and uses it for feature extraction and multi-scale feature description, with all scale features participating in matching simultaneously. The scale space of this invention does not involve scale information in feature matching, but rather plays a role in searching the corresponding level, thus effectively solving the problem of scale differences.

[0146] Step 2: After finding the correct matching search region, template matching is performed on the search region of the image to be matched using the feature point location. Specifically, a similarity metric based on SSD is used to transform to the frequency domain, and Fast Fourier Transform (FFT) is used for template matching to find corresponding points with similarity greater than the threshold as the correct matching points.

[0147] Let M1 and M2 represent the 2D feature description maps (FDCM) of the reference image and the image to be registered, respectively. The SSD between corresponding template windows on the feature description map is defined as:

[0148]

[0149] In the formula, p represents the position of a pixel in the feature description map; T i (s) represents the candidate matching window on the feature descriptor map M1 of the reference image, shifted by s to the corresponding position on the feature descriptor map M2 of the image to be registered. The SSD similarity metric function between the two windows is defined as follows: When T reaches its minimum value, optimal matching between the template window on M1 and the candidate matching window on M2 is achieved.

[0150]

[0151] In the formula, s i This represents the displacement vector between M1 and M2 when the SSD similarity metric function of the corresponding template windows reaches its minimum value. Expanding and simplifying the above equation, the displacement vector s depends only on the parameter in the following equation:

[0152]

[0153] Convolution in the spatial domain is equivalent to multiplication in the frequency domain. Using 2D FFT to accelerate spatial domain convolution calculations, the above equation can be defined as:

[0154]

[0155] In the formula, F 2D and These represent forward and backward 2D FFTs, respectively. This represents the complex conjugate of a 2D FFT.

[0156] Figures 6(1a) to 6(2f) Similarity graphs for various similarity measurement methods are shown. The test images, Figures 6(1a) and 6(1b), are optical-SAR image pairs, and Figures 6(2a) and 6(2b) are optical-infrared image pairs. The NCC, HOPC, CFOG, and FDCM similarity measurement methods were used respectively, and the results are shown below. Figures 6(1d) to 6(1f) , Figures 6(2d) to 6(2f)As shown, the NCC method failed to find the correct matching relationship; HOPC could detect the correspondence, but it suffered from severe noise interference; both CFOG and FDCM could detect the correct corresponding region, but CFOG was affected by interference from other peaks, which easily led to mismatches; FDCM had smoother peaks and better robustness.

[0157] Step 3: After obtaining the matching point pairs in each direction, score the number of matching points obtained in the r-th direction. r This is recorded as the matching score. r The highest direction is determined as the principal characteristic direction r. main The formula is:

[0158]

[0159] Since the rotation difference of the matched images is not strictly in the principal direction r in most cases... main In this case, correct matches may occur in both the main direction and the forward / backward directions. Therefore, this invention merges the feature matching results of the direction adjacent to the main direction to form the optimal matching result.

[0160] Step 4: The feature point positions in the best matching result are the rotated Feature Description Map (FDCM). r The coordinates on the map are not the positions on the original feature map, so they need to be rotated inversely to map them back to the coordinates on the original feature map:

[0161]

[0162] In the formula, Key points for the orientation features after rotation; Rot -1 For the feature description map FDCM r The reverse rotation transformation outputs the result. These are the coordinates on the original feature description map.

[0163] 3) Since there are still a large number of outliers in the best matching results that need to be removed, otherwise it will affect the parameter estimation of the transformation model, this step uses the MAGSAC++ method with high accuracy to eliminate outliers.

[0164] 4) Based on the correct corresponding scale level, transform the matrix parameters step by step with precision, and determine the accurate spatial transformation relationship between image pairs based on the final calculated precise transformation matrix parameters.

[0165] This completes the multimodal remote sensing image registration method based on multi-scale template matching of the present invention.

[0166] To verify the matching performance of the MSTM proposed in this invention, several typical multi-source remote sensing image data were selected, and the following data were compared with CFOG (Y. Ye, J. Shan, L. Bruzzone, and L. Shen, “Robust registration of multimodal remote sensing images based on structural similarity,” IEEE Transactions on Geoscience and Remote Sensing, vol. 55, no. 5, pp. 2941–2958, 2017.), SuperGlue (P.-E. Sarlin, D. DeTone, T. Malisiewicz, and A. Rabinovich, “Superglue: Learning feature matching with graph neural networks,” in Proceedings of the IEEE / CVF conference on computer vision and pattern recognition, 2020, pp. 4938–4947.), and RIFT (J. Li, Q. Hu, and M. Ai, “Rift: Multi-modal image matching based on radiation-variation insensitive feature transform,” IEEE Transactions on Image A comparative experiment was conducted on several effective matching algorithms, including Processing, vol.29, pp.3296–3310, 2019., HAPCG (Yao Yongxiang, Zhang Yongjun, Wan Yi, Liu Xinyi, Guo Haoyu. Heterogeneous image matching considering anisotropic weighted moment and absolute phase direction [J / OL]. Journal of Wuhan University (Information Science Edition): 1-13 [2021-04-02].), Ms-HLMO (C.Gao, W.Li, R.Tao, and Q.Du, “Ms-hlmo: Multiscale histogram of local main orientation for remotesensing image registration,” IEEE Transactions on Geoscience and RemoteSensing, vol.60, pp.1–14, 2022.).The computer used in the experiment was a Lenovo Y9000K laptop with an i7-10875H CPU, a GeForce RTX 2080 graphics card (8GB VRAM), 64GB of RAM, and Windows 10 64-bit operating system.

[0167] To ensure fairness in the comparative experiments, the default optimal parameters were used for all algorithms. The parameters for SuperGlue, RIFT, HAPCG, and Ms-HLMO were consistent with the versions published by their authors. CFOG was removed because the authors' version required known initial control point information, making it unsuitable for multimodal data lacking known information; the remaining parameters remained unchanged.

[0168] (1) Data set and evaluation metrics.

[0169] Six real-world multimodal image datasets were collected for experiments, including optical-depth, optical-infrared, optical-electronic map, optical-optical, optical-SAR, and day-night datasets. The image sizes ranged from 400×400 pixels to 1000×1000 pixels. The optical-optical dataset included different seasons and viewing angles. These multimodal remote sensing images cover almost all scenarios for multimodal image matching. Each dataset contains 10 image pairs, for a total of 60 multimodal image pairs. Most images exhibit significant nonlinear radiometric distortion, making them highly representative and valuable for experimental testing, effectively validating the performance of the multimodal matching algorithm. Figures 9(a) to 9(f) Several sample pairs from the dataset are shown.

[0170] First, feature matching was performed on the images of the six datasets. Matching pairs with an error of less than 3 pixels were considered correct matches. The number of correct matches (NCM), success rate (SR), and root mean square error (RMSE) were used as evaluation metrics, and their calculation formulas are shown below:

[0171]

[0172] In the formula, and H is the coordinates of a pair of matching points; H is the ground reality transformation between image pairs (calculated by manually selecting points); N is the total number of all matching points.

[0173]

[0174] In the formula, NCM represents the number of correct matches; From the reference image via the transformation matrix H It is calculated; These are the actual coordinates of the feature points on the image to be matched.

[0175]

[0176] (2) Parameter settings.

[0177] The main parameters of the MSTM algorithm proposed in this invention include the template size of the highest and lowest levels, and the number of rotation directions in the omnidirectional aggregated feature vector. To determine the optimal parameters for each level of template, we conducted comparative experiments, the results of which are shown in Figures 10(a) and 10(b). Figure 10(a) shows that when the highest level template size is 200 pixels, the matching success rate is higher when the lowest level template size is greater than 100 pixels. Figure 10(b) shows that the matching time increases exponentially with the template size. Based on the above analysis, we determined the highest level template size J1 to be 200 pixels and the lowest level template size J2 to be 100 pixels. When encountering images with matching difficulties, the lowest level template size can be increased to 200 pixels.

[0178] To determine the number of directions in the omnidirectional aggregated feature vector of the MSTM algorithm, we randomly selected one image pair from each of the six datasets for a test experiment. The results are as follows: Figure 11 As shown. Keep one image unchanged, rotate the other image from -90° to 90°, and calculate the number of correctly matched points after rotation. From... Figure 11 As can be seen, the MSTM algorithm without rotation invariance has an adaptability of around ±15° for rotated image pairs, with the number of directions n set to 24. Theoretically, this parameter can provide good adaptability for rotations in all directions.

[0179] (3) Invariance test.

[0180] The key to template-based multimodal image matching lies in whether the template matching based on the constructed feature description map has good robustness against nonlinear radiation distortion. Based on analysis, rotational and scale differences are relatively common types of distortion. Therefore, we designed experiments to test whether the algorithm has good invariance in these two aspects and whether it can handle various multimodal image data.

[0181] ① Rotation Invariance. To test the good adaptability of this method to rotation, for image pairs in the above dataset that have no rotational difference, one image is fixed while the other is rotated. The rotation angle is 0° to 359°, with an interval of 1°, generating 360 image pairs per pair. Then, a matching algorithm is performed on each image pair, and... Figure 12 The corresponding NCMs for the dataset are plotted in the image. Matching results at different rotation angles are shown below. Figures 13(a) to 13(h) As shown. From Figure 12It can be observed that the NCMs for all six scenarios change with the rotation angle, but remain stable within a certain range. No matching failures occurred at any angle, demonstrating that the proposed algorithm has good adaptability to rotation. Figures 13(a) to 13(h) It can be seen that, for image pairs at various rotation angles, the matching results of the method of the present invention are sufficient to meet the needs of image registration and other functions.

[0182] ② Scale Invariance. To test the adaptability of the method of this invention to images with certain scale differences, the same strategy was used to fix one image in an image pair and scale the other image. The scaling ratio was 1-2 with an interval of 0.2, generating 5 image pairs for each image pair. Some matching results are shown in Figure 14(a) and Figure 14(b). When images have different scale differences, the method can still obtain accurate matching results. This is because we adopted a multi-scale matching strategy in the algorithm, which enhances the matching performance for images with scale differences and has a certain adaptability to multimodal image data with certain scale differences.

[0183] (4) Qualitative analysis.

[0184] Select one image pair from each of the six multimodal datasets for testing, such as... Figures 15(1a) to 15(6f) As shown in Figures 15(1a) and 15(1e), these images contain translation and minor rotation variations; Figures 15(1b) and 15(1d) contain only translation variations; and Figures 15(1c) and 15(1f) contain translation, minor rotation, and scaling variations. Because these image pairs are multimodal image data from different imaging devices, they exhibit severe nonlinear distortion. Therefore, testing on these image pairs is very challenging. Figures 15(1a) to 15(6f) Matching results for CFOG, Ms-HLMO, HAPCG, RIFT, SuperGlue, and the proposed MSTM were plotted respectively.

[0185] The results show that the CFOG algorithm failed to match any images, achieving a SR accuracy of 0. Although CFOG performs pixel-by-pixel feature description, possessing strong feature description capabilities, its requirement for initial reference control points causes a sharp performance drop on images with displacement, rotation, and scale differences. The Ms-HLMO algorithm failed to match the second and third image pairs, achieving an SR accuracy of 50%. The HAPCG algorithm failed to match the first, fifth, and sixth image pairs, also achieving an SR accuracy of 50%. While two image pairs failed to match, HAPCG performed well on successfully matched pairs. This is because HAPCG constructs an anisotropic scale space and uses logarithmic polar coordinates for feature description, resulting in strong robustness. The SuperGlue algorithm failed to match the third and fourth image pairs, achieving an SR accuracy of 66.6%. Although SuperGlue employs deep learning, combining attention mechanisms and graph neural networks, it exhibits strong robustness to multimodal data. However, due to the lack of multimodal datasets in its training dataset, its matching performance still degrades on some multimodal datasets. The RIFT algorithm failed to match on the third image pair, achieving an SR accuracy of 83.3%. Because the RIFT algorithm incorporates a MIM map based on the phase-consistent map, it is more robust to nonlinear radiometric distortions. Furthermore, it analyzes the effect of rotation on the MIM, achieving rotation invariance in feature description, thus improving the algorithm's performance. In contrast, the MSTM algorithm proposed in this invention successfully matched on all six image pairs, achieving an SR accuracy of 100%.

[0186] MSTM outperforms the aforementioned popular algorithms in matching multimodal images, primarily due to: ① MSTM employs an FDCM feature description map, which exhibits better feature robustness compared to other algorithms, laying a crucial foundation for subsequent matching. ② The adoption of a multi-scale matching strategy enhances the algorithm's robustness to multimodal images, ensuring matching accuracy.

[0187] (5) Quantitative analysis.

[0188] Figures 16(a) to 16(f)The quantitative results of the number of NCMs are shown, along with the results of six methods on six multimodal datasets. The results indicate that the CFOG algorithm outperforms the other datasets on the optical-depth dataset. This is likely because the displacement and scale differences between depth and optical images are small, allowing for successful matching without initial reference control points. The Ms-HLMO algorithm shows relatively balanced performance across all datasets, thanks to its construction of PMOM, use of the GGLOH descriptor, and multi-scale matching strategy, which significantly enhances stability. The SuperGlue algorithm outperforms the other three algorithms on the optical-depth dataset, possibly because it is trained using a combination of RGB images and depth data. Its performance on the optical-SAR dataset is similar to the other three algorithms. This is partly due to the severe nonlinear radiometric distortion between SAR and optical images, making matching inherently difficult; and partly because this dataset has the largest gap with SuperGlue's training dataset, causing SuperGlue's training parameters to be unsuitable for this dataset. The HAPCG algorithm outperforms all other algorithms on some datasets, primarily because it utilizes anisotropic filtering to perform nonlinear diffusion on images, constructing an anisotropic weighted moment scale space. Then, it extends the phase consistency model, establishing an absolute phase consistency directional gradient and combining it with a log-polar coordinate descriptor template to create an absolute phase consistency directional gradient histogram, significantly enhancing the robustness of the descriptor. The RIFT algorithm performs relatively consistently across all datasets because it incorporates a MIM map based on the phase consistency map, making it more robust to nonlinear radiometric distortion. Furthermore, it achieves rotation invariance in feature description, contributing to its superior performance. In contrast, the proposed MSTM algorithm successfully matched all images across six datasets, with almost all image pairs having NCMs far greater than 250. This algorithm exhibits stable and robust matching performance, demonstrates good adaptability to nonlinear radiometric distortion, and outperforms other algorithms.

[0189] Table 1 Comparison of Matching Success Rates

[0190]

[0191] Table 1 summarizes the matching success rates (SRs) of the six methods for each dataset. As shown in the table, the CFOG algorithm performs best on the optical-depth dataset with an SR of 90%; the Ms-HLMO algorithm performs worst on the optical-SAR dataset with an SR of only 50%; the HAPCG algorithm achieves 100% SR on multiple datasets; the SuperGlue algorithm performs worst on the optical-SAR dataset with an SR of only 70%, which is inferior to traditional algorithms; RIFT performs consistently across all datasets and is generally successful in matching; the MSTM algorithm achieves 100% SR on all datasets. The average SRs of CFOG, Ms-HLMO, HAPCG, SuperGlue, and RIFT across all six datasets are 60%, 80%, 86.6%, 88.3%, and 98.3%, respectively. Compared to the CFOG algorithm, the MSTM algorithm shows a 40% performance improvement. Table 2 and... Figures 17(a) to 17(f) The RMSE of six methods across all datasets is presented. Since some algorithms have lower SR on certain datasets, only the RMSE of successfully matched images is calculated. In terms of RMSE, feature-based matching methods have slightly lower accuracy; HAPCG and RIFT algorithms are both lower than CFOG and MSTM algorithms. CFOG is currently the most accurate multimodal matching algorithm. Compared to CFOG, MSTM's RMSE accuracy is basically consistent with CFOG, and even slightly higher on some datasets, demonstrating superior performance.

[0192]

[0193] Based on the qualitative and quantitative test results above, the following conclusions can be drawn: Each part of the MSTM algorithm is designed specifically for nonlinear radiation distortion, including feature point extraction, feature descriptor mapping, template matching, omnidirectional feature vector aggregation, and multi-scale strategies. Therefore, this algorithm has good adaptability to nonlinear radiation distortion, achieving excellent NCM and RMSE accuracies on all six datasets, surpassing other algorithms. Thus, the method described in this invention is a superior multimodal matching method.

[0194] (6) Registration effect.

[0195] This section applies MSTM to image registration and image fusion. MSTM is used for image matching, and after calculating the transformation matrix, the two images are mapped onto the same reference frame to generate the corrected image. Tests were conducted on six datasets, and the visualization results are shown below. Figures 18(1a) to 18(6h) As shown, all images were well registered and fused, without displacement, ghosting, or blurring, further demonstrating that the MSTM algorithm achieves high accuracy and good distribution in matching.

[0196] (7) Time analysis.

[0197] Table 3 shows the average runtime of each algorithm across all datasets. CFOG, Ms-HLMO, HAPCG, and RIFT were computed using Matlab 2021a, while SuperGlue and MSTM were computed using Python. The table shows that MSTM is slightly slower than CFOG and SuperGlue, taking about one-fifth of the time of HAPCG and RIFT. Analyzing the algorithm principles, the significant time consumption of HAPCG is mainly due to the substantial time spent constructing the nonlinear diffusion scale space and generating high-dimensional feature vectors, leading to significant time requirements for subsequent matching. The main time consumption of the RIFT algorithm is in the construction of the MIM graph and the implementation of rotation invariance. Ms-HLMO is the most time-consuming because it constructs a PMOM based on the MIM graph and then employs the GGLOH descriptor and a multi-scale matching strategy. The MSTM algorithm employs a multi-scale matching strategy, which should result in a significantly longer runtime than other algorithms. However, because MSTM reduces the dimensionality of multi-dimensional feature vectors when constructing the feature map, it greatly reduces the search time for the template matching window. Furthermore, the matching results from the upper pyramid can be shared with the lower pyramid, allowing for a rough estimation of the mapping relationship and saving considerable time. While MSTM is rotation-invariant, it takes some time to generate omnidirectional aggregated feature vectors, thus its runtime is slightly longer than CFOG and SuperGlue.

[0198] Table 3 Comparison of Time for Various Algorithms

[0199]

[0200] In summary, addressing the inability of traditional image matching algorithms to effectively handle multimodal images with severe nonlinear radiation distortion, this invention proposes a novel multiscale template matching algorithm (MSTM) for multimodal image matching. To address this issue, a novel frequency domain convolutional map (FDCM) based on wavelet transform and phase consistency is first proposed to construct a feature description map, significantly reducing nonlinear radiation distortion between multimodal images. Simultaneously, robustness to rotated images is achieved by designing rotation-invariant omnidirectional aggregated feature vectors. Finally, a multiscale template matching strategy improves matching performance on multimodal images with displacement and scale changes. To improve the algorithm's time efficiency, most complex calculations are performed in the frequency domain. Furthermore, experimental results on six multimodal image datasets demonstrate that this method can achieve accurate and robust matching results between multimodal images. Qualitative and quantitative evaluations confirm that this method outperforms several mainstream multimodal image matching algorithms in terms of matching accuracy, matching success rate, and matching time.

Claims

1. A multimodal remote sensing image registration method based on multi-scale template matching, characterized in that, Includes the following steps: 1) Acquire multimodal remote sensing image pairs, namely a reference image and an image to be registered; Each image is decomposed into low-frequency wavelet components, as well as high-frequency wavelet components in the horizontal, vertical, and diagonal directions. The Log-Gabor function is used to obtain the amplitude components of the low-frequency wavelet components at different scales, which serve as the features of the low-frequency wavelet components at different scales. Features at different scales are normalized, and low-frequency wavelet components are reconstructed using the normalization results. Wavelet reconstruction is performed based on the reconstructed low-frequency wavelet components and the high-frequency wavelet components on the horizontal, vertical, and diagonal lines to obtain the frequency domain convolutional map of each image. Multiple layer-by-layer downsampling operations are performed on the frequency domain convolutional maps of the two images to construct the Gaussian pyramids of the two images. 2) Extract feature points from the image pairs and rotate the Gaussian pyramid of the images to be registered in different directions; based on the extracted feature points, use a template matching algorithm to match the remote sensing image pairs in different rotation directions, and select the direction with the most matching points as the main direction. 3) Determine the spatial transformation relationship of the remote sensing image pair based on the feature matching results corresponding to the main direction.

2. The multimodal remote sensing image registration method based on multi-scale template matching according to claim 1, characterized in that, In step 2), the rotation transformation formula for rotating the Gaussian pyramid of the image to be registered in different directions is as follows: ； in This is the feature description map for the r-th direction. This indicates a rotation operation on the image; This is the original, non-rotated feature description map; For the first The rotation angle of the feature description map, where n is the total number of rotations.

3. The multimodal remote sensing image registration method based on multi-scale template matching according to claim 1, characterized in that, The feature points extracted in step 2) are Shi-Tomasi corner points.

4. The multimodal remote sensing image registration method based on multi-scale template matching according to claim 3, characterized in that, The extracted Shi-Tomasi corner points are those that have undergone uniform distribution processing, and the uniform distribution processing is performed using an adaptive non-maximum suppression algorithm. The formula for calculating the threshold for the variable length of the suppression region is as follows: ； ； In the formula, To suppress the threshold of region length variation; The number of feature points to be extracted; These represent the height and width of the image, respectively.

5. The multimodal remote sensing image registration method based on multi-scale template matching according to claim 1, characterized in that, In step 2), the template matching algorithm requires multi-scale matching using the following method: perform template matching between a certain level of the Gaussian pyramid of the reference image and each level of the Gaussian pyramid of the image to be registered, in order to find the hierarchical correspondence between the Gaussian pyramid of the reference image and the Gaussian pyramid of the image to be registered; and then find the matching search region in the template matching algorithm process based on the hierarchical correspondence.

6. The multimodal remote sensing image registration method based on multi-scale template matching according to claim 1, characterized in that, In step 3), outlier feature points in the feature matching results need to be eliminated, and the spatial transformation relationship of the remote sensing image pairs needs to be determined based on the feature matching results after the outlier feature points are eliminated.

7. The multimodal remote sensing image registration method based on multi-scale template matching according to claim 6, characterized in that, The MAGSAC++ method is used to eliminate outlier features.

8. The multimodal remote sensing image registration method based on multi-scale template matching according to claim 1, characterized in that, In step 3), when determining the spatial transformation relationship of remote sensing image pairs, it is also necessary to combine the feature matching results of adjacent directions of the main direction.

9. The multimodal remote sensing image registration method based on multi-scale template matching according to claim 1, characterized in that, The similarity measurement method used in the template matching algorithm is as follows: The SSD between the template window corresponding to a certain feature description map included in the Gaussian pyramid of the reference image and the image to be registered is determined using the following formula: ； In the formula, It represents the position of a pixel in the feature description map; It is a feature description map of the reference image. The template window and the feature description map of the image to be registered. displacement of the corresponding position The SSD similarity metric function between subsequent candidate matching windows; the matching function is defined as: ； In the formula, This represents the defined matching function; and These represent forward and backward 2D FFTs, respectively. Represents the complex conjugate of a 2D FFT; Determine the minimum To determine and The best match between the upper windows.

10. The multimodal remote sensing image registration method based on multi-scale template matching according to any one of claims 1 to 9, characterized in that, The multimodal remote sensing image pair in step 1) is the denoised multimodal remote sensing image pair.

Images