WO2019042232A1 - 一种快速鲁棒的多模态遥感影像匹配方法和系统 - Google Patents

一种快速鲁棒的多模态遥感影像匹配方法和系统 Download PDF

Info

Publication number
WO2019042232A1
WO2019042232A1 PCT/CN2018/102271 CN2018102271W WO2019042232A1 WO 2019042232 A1 WO2019042232 A1 WO 2019042232A1 CN 2018102271 W CN2018102271 W CN 2018102271W WO 2019042232 A1 WO2019042232 A1 WO 2019042232A1
Authority
WO
WIPO (PCT)
Prior art keywords
pixel
matching
feature
map
point
Prior art date
Application number
PCT/CN2018/102271
Other languages
English (en)
French (fr)
Inventor
叶沅鑫
Original Assignee
西南交通大学
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by 西南交通大学 filed Critical 西南交通大学
Priority to US16/639,046 priority Critical patent/US11244197B2/en
Publication of WO2019042232A1 publication Critical patent/WO2019042232A1/zh

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06V—IMAGE OR VIDEO RECOGNITION OR UNDERSTANDING
    • G06V20/00—Scenes; Scene-specific elements
    • G06V20/10—Terrestrial scenes
    • G06V20/13—Satellite images
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F18/00—Pattern recognition
    • G06F18/20—Analysing
    • G06F18/22—Matching criteria, e.g. proximity measures
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00—Image analysis
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06V—IMAGE OR VIDEO RECOGNITION OR UNDERSTANDING
    • G06V10/00—Arrangements for image or video recognition or understanding
    • G06V10/40—Extraction of image or video features
    • G06V10/42—Global feature extraction by analysis of the whole pattern, e.g. using frequency domain transformations or autocorrelation
    • G06V10/431—Frequency domain transformation; Autocorrelation
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06V—IMAGE OR VIDEO RECOGNITION OR UNDERSTANDING
    • G06V10/00—Arrangements for image or video recognition or understanding
    • G06V10/40—Extraction of image or video features
    • G06V10/44—Local feature extraction by analysis of parts of the pattern, e.g. by detecting edges, contours, loops, corners, strokes or intersections; Connectivity analysis, e.g. of connected components
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06V—IMAGE OR VIDEO RECOGNITION OR UNDERSTANDING
    • G06V10/00—Arrangements for image or video recognition or understanding
    • G06V10/40—Extraction of image or video features
    • G06V10/46—Descriptors for shape, contour or point-related descriptors, e.g. scale invariant feature transform [SIFT] or bags of words [BoW]; Salient regional features
    • G06V10/462—Salient features, e.g. scale invariant feature transforms [SIFT]
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06V—IMAGE OR VIDEO RECOGNITION OR UNDERSTANDING
    • G06V10/00—Arrangements for image or video recognition or understanding
    • G06V10/70—Arrangements for image or video recognition or understanding using pattern recognition or machine learning
    • G06V10/74—Image or video pattern matching; Proximity measures in feature spaces
    • G06V10/761—Proximity, similarity or dissimilarity measures

Definitions

  • the invention relates to the technical field of satellite image processing, in particular to a multi-modal remote sensing image matching method and system, which is used for automatic matching of multi-modal images such as visible light, infrared, laser radar, contract aperture radar and map.
  • Image matching is the process of identifying the same-named point between two or more images. It is the basic pre-processing step of many remote sensing image analysis such as image fusion, change detection, image mosaic, etc. Its matching accuracy is important for subsequent analysis work. influences.
  • remote sensing satellite sensors are equipped with global positioning system and inertial navigation system, which can directly locate and coarsely match remote sensing images, eliminating obvious rotation and scale changes between images, so that there is only a certain amount between images (such as within tens of pixels). The difference in translation.
  • there are significant nonlinear radiation differences between multi-modal remote sensing images visible light, infrared, lidar and contract aperture radar, etc.
  • the matching methods of multi-modal remote sensing images can be mainly divided into: feature matching and region matching.
  • Feature matching is to achieve image matching by the similarity of salient features between images.
  • Commonly used features include point feature methods, line feature methods, and surface features.
  • local invariant features such as Scale Invariant Feature Transform (SIFT) and shape context have also been applied in remote sensing image matching.
  • SIFT Scale Invariant Feature Transform
  • shape context have also been applied in remote sensing image matching.
  • these methods usually need to extract features with high repetition rate between images, and for multi-modal remote sensing with significant radiation differences, the repetition rate of feature extraction is often low, so they are automatic for multi-modal remote sensing images.
  • SIFT Scale Invariant Feature Transform
  • the region-based method mainly adopts the template matching strategy, and uses a certain matching similarity measure as a criterion to perform the same-name point identification between images.
  • the selection of the similarity measure is crucial, which directly affects the subsequent matching accuracy.
  • Commonly used similarity measures include the sum of squared gray-scale differences, normalized correlation coefficients, and mutual information.
  • these similarity measures use the similarity of gray scales between images to identify the same name.
  • due to the large difference in gray scale information between multi-modal remote sensing images they are not suitable for multi-modal remote sensing images. Automatic matching.
  • HOG Histogram of Orientated Gradient
  • LSS Local Self-Similarity
  • the present invention constructs a fast and robust multi-modal remote sensing image matching framework, which can integrate various local feature descriptors for automatic matching of multi-modal remote sensing images.
  • the framework first extracts local feature descriptors such as HOG, LSS or Speeded-Up Robust Features (SURF) for each pixel of the image to form a dense pixel-by-pixel feature representation map to reflect the common structures, shapes and textures between the images. .
  • SURF Speeded-Up Robust Features
  • a fast matching similarity measure is established in the frequency domain by using the three-dimensional Fourier transform, and the template matching strategy is used to identify the same name.
  • CFOG Channel Feature of Orientated Gradient
  • the object of the present invention is to overcome the deficiencies of the conventional matching method and provide a fast and robust multi-modal remote sensing image matching framework.
  • the framework utilizes pixel-by-pixel feature representation techniques to extract features such as structure, shape and texture shared between images, and establishes a fast matching similarity measure based on it, which can quickly and accurately between multi-modal remote sensing images. Get a large number of evenly distributed points of the same name.
  • CFOG Channel Feature of Orientated Gradient
  • the present invention provides a fast and robust multi-modal remote sensing image matching method, comprising the following steps:
  • the step D includes the following steps: calculating local feature descriptors such as HOG, LSS, or SURF for each pixel for the image data in the matching area, and arranging the feature vectors corresponding to each pixel in the Z direction to form a three-dimensional image. Pixel-by-pixel feature representation.
  • step D may further be to construct a directional gradient channel feature (CFOG) in the matching area, which specifically includes the following steps:
  • a Gaussian filter is used to convolute the three-dimensional directional gradient map to obtain a feature map.
  • Reusing the one-dimensional filter [1, 2, 1] to map the features in the Z direction Perform a convolution operation to obtain a feature map
  • the step D is to construct a directional gradient channel feature (CFOG) in the matching area, and the specific calculation process includes the following steps:
  • the gradients in the horizontal direction (X direction) and the vertical direction (Y direction) are calculated using the one-dimensional filters [-1, 0, 1] and [-1, 0, 1] T , respectively.
  • ⁇ represents the direction of the gradient of the quantization, and abs represents the absolute value. The symbol indicates that the value is positive when it is taken, otherwise it is 0;
  • the g ⁇ in each direction is superimposed to form a three-dimensional gradient g o , and then a convolution operation is performed on the g o in the X and Y directions using a two-dimensional Gaussian filter of standard ⁇ to obtain a feature map.
  • a convolution operation is performed on the g o in the X and Y directions using a two-dimensional Gaussian filter of standard ⁇ to obtain a feature map.
  • Reuse one-dimensional filter [1, 2, 1] in the Z direction Perform a convolution operation to obtain a feature map
  • Feature map Each pixel in the Z corresponds to a feature vector f i in the Z direction, traverses each pixel, normalizes its feature vector v i , and obtains the final directional gradient channel feature map.
  • the normalized calculation formula as follows:
  • ⁇ is a number that avoids dividing by zero.
  • the gradient direction is evenly divided into 18 equal parts in a range of 360°, the angular interval of each aliquot is 20°, and the values of ⁇ are ⁇ 0°, 20°, ..., 340° ⁇ .
  • the pixel-by-pixel feature expression map is converted into the frequency domain by using the three-dimensional fast Fourier transform, and the correlation operation is performed to obtain the similarity map, and the position of the maximum value of the similarity map is the matching position of the image.
  • the step E includes the following calculation steps:
  • step D the pixel-by-pixel feature representation graphs D 1 and D 2 of the regions AreaW 1i and AreaW 2i are respectively obtained, and D 1 is used as a template to slide at D 2 , and the sum of squares of the feature vector differences between D 1 and D 2 is utilized. Matching as a similarity measure;
  • c represents the pixel coordinates in the feature representation map
  • v represents the offset between D 1 and D 2
  • S i represents the sum of the squares of the eigenvector differences between D 1 and D 2 , when S i takes the minimum
  • the offset v i between D 1 and D 2 that is, the matching position, is obtained, and the calculation formula is:
  • the fast Fourier transform is used in the frequency domain to improve its computational efficiency, and the similarity function based on the Fourier transform is obtained as follows:
  • F and F -1 represent the fast Fourier forward transform and the inverse transform, respectively
  • F * represents the complex conjugate of F.
  • D 1 and D 2 are three-dimensional feature representation maps, according to the principle of convolution, it is necessary to use the three-dimensional fast Fourier transform to calculate the formula (7), and the final similarity function is:
  • 3DF and 3DF -1 represent a three-dimensional fast Fourier forward transform and an inverse transform, respectively, and 3DF * represents a complex conjugate of 3DF.
  • the present invention also provides a fast and robust multi-modal remote sensing image matching system, the system comprising the following units:
  • a pre-processing unit for determining resolution information of the reference image and the input image. If the two images have the same resolution, the subsequent units are entered, and if the resolutions are different, the two images are sampled to the same resolution;
  • a feature extraction unit configured to construct a pixel-by-pixel feature representation map in the matching area
  • the initial matching unit is used to establish a fast similarity measure for the same-name point matching by using the three-dimensional Fourier transform on the basis of the pixel-by-pixel feature representation map; for the same-named point, the local extremum of the similarity map is obtained.
  • the feature extraction unit is configured to calculate a local feature descriptor of each pixel for the image data in the matching region, and arrange the feature vectors corresponding to each pixel in the Z direction to form a three-dimensional pixel-by-pixel feature representation.
  • the initial matching unit converts the pixel-by-pixel feature representation map into the frequency domain by using a three-dimensional fast Fourier transform, and performs a correlation operation to obtain a similarity map, where the maximum value of the similarity map is the matching position of the image.
  • the present invention constructs a fast and robust multi-modal remote sensing image matching framework, and extracts a local feature descriptor of each pixel (such as HOG, LSS or SURF, etc.) to form a pixel-by-pixel feature representation map. It can better reflect the common structure, shape and texture of multi-modal remote sensing images. Based on the pixel-by-pixel feature representation, a fast matching similarity measure is established by using three-dimensional Fourier transform.
  • the framework can quickly and accurately acquire a large number of uniformly distributed points of the same name among multi-modal images, which can effectively improve the actual production efficiency of matching and meet the needs of operational operation, and it is a general framework that can integrate each Local feature description for image matching (not limited to descriptors such as HOG, LSS or SURF).
  • the pixel-by-pixel feature representation is formed by calculating local features such as HOG, LSS or SURF for each pixel, and is a dense feature representation technique. This is different from traditional HOG, LSS or SURF descriptors, which only perform feature construction in a sparse sampling grid (not for each pixel), or use the neighborhood information of the extracted feature points to calculate features. A relatively sparse feature representation.
  • the pixel-by-pixel feature representation technique of the present invention can better and more accurately reflect common structures, shapes, and textures among multi-modal remote sensing images, and the matching performance is more robust, and the three-dimensional based on the method is proposed.
  • the fast similarity measure of the Fourier transform enables fast matching between multimodal images.
  • the present invention constructs a novel pixel-by-pixel feature descriptor-direction gradient channel feature (CFOG) based on the direction gradient information of the image, which is superior to pixel-by-pixel in terms of matching efficiency and precision.
  • CFOG pixel-by-pixel feature descriptor-direction gradient channel feature
  • the present invention establishes a fast matching similarity measure in the frequency domain by using the three-dimensional Fourier transform. Compared with the similarity measures commonly used in the spatial domain, such as the squared difference of gray scales, the normalized correlation coefficient and mutual information, the calculation efficiency is higher and the calculation effect is better.
  • FIG. 2 is a schematic diagram of pixel-by-pixel feature representation of the present invention.
  • Figure 3 is a process of constructing a directional gradient channel feature (CFOG) in the present invention.
  • Figure 1 shows a fast and robust multi-modal remote sensing image matching method, including the following steps:
  • step A according to the resolution information of the reference image and the input image, it is determined whether the resolutions of the two images are consistent, and if the images are consistent, the subsequent processing is performed, and if the images are inconsistent, the two images are sampled and processed according to the same resolution.
  • Step B using a block strategy, using a Harris or Forstner operator to extract a large number of uniformly distributed feature points on the reference image, including:
  • the reference image is divided into n ⁇ n square grids that do not overlap each other.
  • the Harris or Forstner eigenvalues of each pixel are calculated, and the eigenvalues are sorted from large to small, and the eigenvalues are selected.
  • a larger k pixels are used as feature points.
  • k feature points can be detected in each grid, and n ⁇ n ⁇ k feature points are provided on the entire image.
  • the values of n and k can be set according to actual needs.
  • the feature extraction of the image may also be performed by other operators, which is not limited in the present invention.
  • Step D constructing a pixel-by-pixel feature representation map in the areas AreaW 1i and AreaW 2i , see FIG.
  • step D includes the following steps:
  • the local feature descriptors such as HOG, LSS, or SURF for each pixel in the region are calculated, and the feature vectors corresponding to each pixel are arranged in the Z direction to form a three-dimensional pixel-by-pixel feature representation map.
  • this embodiment does not limit this.
  • step D is to construct a Channel Feature of Orientated Gradient (CFOG) in the areas AreaW 1i and AreaW 2i .
  • CFOG Channel Feature of Orientated Gradient
  • ⁇ denotes the direction of the quantized gradient, where the gradient direction is evenly divided into 18 equal parts in the range of 360°, and the angular interval of each aliquot is 20°, so the value of ⁇ is ⁇ 0°, 20 °, whil,340° ⁇ .
  • Abs denotes an absolute value, the purpose of which is to convert gradient information with a gradient direction of [180°, 360°) between [0°, 180°), which can reduce the gradient inversion caused by multimodal images. Match the impact. The symbol indicates that the value is positive when it is taken, otherwise it is taken as 0.
  • Each of the pixels corresponds to a feature vector f i in the Z direction.
  • Each pixel is traversed, and its feature vector v i is normalized to further eliminate the influence of illumination changes, and the final CFOG feature map is obtained.
  • the normalized formula is as follows:
  • ⁇ is a number that avoids dividing by zero.
  • Step E based on the pixel-by-pixel feature representation graph, uses a three-dimensional Fourier transform to establish a fast similarity measure for the same-name point matching, including the following steps:
  • step D the pixel-by-pixel feature representation graphs D 1 and D 2 corresponding to the regions AreaW 1i and AreaW 2i respectively can be obtained, and D 1 is used as a template to slide on D 2 , and the feature vector difference between them can be utilized.
  • the sum of squares is matched as a measure of similarity.
  • the formula for summing the difference sum is simplified below, and the matching process is accelerated by using a three-dimensional fast Fourier transform.
  • c represents the pixel coordinates in the feature representation map
  • v represents the offset between D 1 and D 2
  • S i represents the sum of the squares of the feature vector differences between D 1 and D 2 .
  • the fast Fourier transform is used in the frequency domain to improve its computational efficiency.
  • the similarity function based on the Fourier transform is:
  • F and F -1 represent the fast Fourier forward transform and the inverse transform, respectively, and F * represents the complex conjugate of F.
  • D 1 and D 2 are three-dimensional feature expression maps, according to the principle of convolution, it is necessary to calculate equation (7) using a three-dimensional fast Fourier transform.
  • the final similarity function is:
  • 3DF and 3DF -1 represent a three-dimensional fast Fourier forward transform and an inverse transform, respectively, and 3DF * represents a complex conjugate of 3DF;
  • Step F for the above-mentioned point pair ⁇ P 1i (x, y), P * 2i (x, y) ⁇ , using a binary quadratic polynomial for local interpolation
  • the present invention also provides a fast and robust multi-modal remote sensing image matching system, the system comprising the following units:
  • a pre-processing unit for determining resolution information of the reference image and the input image. If the two images have the same resolution, the subsequent units are entered, and if the resolutions are different, the two images are sampled to the same resolution;
  • a feature extraction unit configured to construct a pixel-by-pixel feature representation map in the matching area
  • the initial matching unit is used to establish a fast similarity measure for the same-name point matching by using the three-dimensional Fourier transform on the basis of the pixel-by-pixel feature representation map; for the same-named point, the local extremum of the similarity map is obtained.
  • the feature extraction unit is configured to calculate a local feature descriptor of each pixel for the image data in the matching region, and arrange the feature vectors corresponding to each pixel in the Z direction to form a three-dimensional pixel-by-pixel feature representation.
  • the initial matching unit converts the pixel-by-pixel feature representation map into the frequency domain by using a three-dimensional fast Fourier transform, and performs a correlation operation to obtain a similarity map, where the maximum value of the similarity map is the matching position of the image.
  • the directional gradient channel feature (CFOG) constructed by the present invention is a novel pixel-by-pixel feature representation technique, which is superior to the pixel-by-pixel feature representation based on descriptors such as HOG, LSS and SURF in matching efficiency and precision. technology.
  • the technical scheme of the invention can make up for the deficiencies of the traditional matching method for the nonlinear radiation difference between multi-modal images, and can effectively solve the matching of multi-modal remote sensing data such as visible light, infrared, laser radar, synthetic aperture radar and map. problem.
  • the technical solution proposed by the present invention is a general technical framework, and can integrate various local feature descriptors (not limited to CFOG, HOG, LSS, SURF, etc.) for image matching.
  • the invention is not limited to the specific embodiments described above.
  • the invention extends to any new feature or any new combination disclosed in this specification, as well as any novel method or process steps or any new combination disclosed.

Abstract

本发明公开了一种快速、鲁棒的多模态遥感影像匹配方法和系统，可以整合各种局部特征描述符进行多模态遥感影像自动匹配。首先对影像的每个像素提取梯度方向直方图（Histogram of Oriented Gradient，HOG）,局部自相似（local self-similarity，LSS）或Speeded-Up Robust Features（SURF）等局部特征描述符，形成逐像素的特征表达图。然后基于该特征表达图，利用三维傅里叶变换在频率域建立一种快速的匹配相似性测度。最后采用模板匹配的策略进行同名点识别。另外，针对所发明的匹配方法和系统，还提出了一种新的逐像素特征表达技术，名为梯度方向特征通道（channel features of orientated gradients，CFOG），它在匹配性能和计算效率方面要优于基于HOG、LSS和SURF等描述符的逐像素特征表达方式。本发明能有效克服可见光、红外、激光雷达、合成孔径雷达以及地图等多模态影像间的非线性辐射差异，在影像间快速、精确地识别出同名点，实现影像的自动匹配。

Description

一种快速鲁棒的多模态遥感影像匹配方法和系统 技术领域
本发明涉及卫星影像处理技术领域,尤其是一种多模态遥感影像匹配方法和系统，用于可见光、红外、激光雷达、合同孔径雷达以及地图等多模态影像的自动匹配。
背景技术
影像匹配是在两幅或多福影像间识别同名点的过程，它是诸多遥感影像分析如影像融合，变化检测，影像镶嵌等的基本预处理步骤，其匹配精度对后续的分析工作产生重要的影响。目前遥感卫星传感器装载有全球定位系统和惯性导航系统，可对遥感影像进行直接定位和粗匹配，消除影像间明显的旋转和尺度变化，使影像间仅存在一定量(如几十个像素内)的平移差异。尽管如此，由于影像的成像机理的不同，多模态遥感影像间(可见光，红外，激光雷达和合同孔径雷达等)存在显著的非线性辐射差异，导致同名点的自动匹配仍然非常具有挑战性。
目前多模态遥感影像的匹配方法主要可分为：特征匹配和区域匹配。特征匹配是通过影像间显著特征的相似性来实现影像的匹配。常用的特征包括了点特征的方法，线特征的方法和面特征。最近局部不变性特征如Scale Invariant Feature Transform(SIFT)，shape context等在遥感影像匹配中也得到了一定的应用。但这些方法通常需要在影像间提取出具有高重复率的特征，而对于具有显著辐射差异的多模态遥感而言，特征提取的重复率往往较低，因此它们对于多模态遥感影像的自动匹配还存在一定的局限性。
基于区域的方法主要是采用模板匹配的策略，以某种匹配相似性测度为准则，在影像间进行同名点识别。在此过程中，相似性测度的选择至关重要，直接影响到后续的匹配精度。常用的相似性测度包括了灰度差平方和、归一化相关系数和互信息等。但这些相似性测度都是利用影像间灰度的相似性进行同名点识别，而由于多模态遥感影像间的灰度信息存在较大的差异，所以它们无法较好适用于多模态遥感影像的自动匹配。相比于灰度信息，影像的结构和形状属性具有较高的相似性，而且相关研究利用梯度方向直方图(Histogram of Orientated Gradient,HOG)、局部自相似(Local Self-Similarity,LSS)等局部描述符提取影像的结构和形状特征，并在此基础上建立相似性测度进行影像匹配，提高了匹配性能。尽管如此，HOG和LSS只是在一个稀疏的采样格网(不是针对每个像素)内进行特征构建，或者利用所提取特征点的邻域信息计算特征，是一种相对稀疏的特征表达方式，难以很好地反映多模态影像间的共有属性，而且其计算效率较低。
鉴于此，本发明构建了一种快速、鲁棒的多模态遥感影像匹配框架，该框架可以整合 各种局部特征描述符进行多模态遥感影像自动匹配。该框架首先对影像的每个像素提取HOG，LSS或Speeded-Up Robust Features(SURF)等局部特征描述符，形成稠密的逐像素特征表达图，来反映影像间共有的结构、形状和纹理等属性。然后基于该特征表达图，利用三维傅里叶变换在频率域建立一种快速的匹配相似性测度，并采用模板匹配的策略进行同名点识别，另外针对所发明的框架，构建了一种基于方向梯度特征的逐像素特征描述符，名为方向梯度通道特征(Channel Feature of Orientated Gradient,CFOG)。它在匹配性能和计算效率方面要优于逐像素的HOG，LSS和SURF等特征表达技术。
发明内容
本发明的发明目的在于：克服传统匹配方法的不足，提供了一种快速鲁棒的多模态遥感影像匹配框架。该框架利用通过逐像素的特征表达技术来提取影像间共有的结构、形状和纹理等特征，并在其基础上建立了快速的匹配相似性测度，可在多模态遥感影像间快速、精确地获取大量分布均匀的同名点。另外针对所发明的框架，构建了一种新颖的逐像素特征表达技术，名为方向梯度通道特征(Channel Feature of Orientated Gradient,CFOG)。
一方面，本发明提供了一种快速鲁棒的多模态遥感影像匹配方法，包括下列步骤：
A.判断参考影像和输入影像的分辨率信息，如果两幅影像具有相同的分辨率，则进行后续处理，如果分辨率不同，则将两幅影像采样为同样的分辨率；
B.采用分块的策略，在参考影像上检测出一系列分布均匀的特征点，记为P 1i(i＝1,2,3,……,N)，以点P 1i为中心选取模板区域AreaW 1i；
C.根据遥感影像自身提供的地理坐标信息，预测点集P 1i(i＝1,2,3,….,N)在输入影像上所对应的匹配区域AreaW 2i；
D.在匹配区域内构建逐像素的特征表达图；
E.在逐像素的特征表达图的基础上，利用三维傅里叶变换建立一种快速的相似性测度进行同名点匹配；
F.对于获得的同名点，对其相似性图进行局部极值拟合，求解出匹配点的亚像素位置；
G.重复步骤(C)—(F),遍历P 1i(i＝1,2,3,…,N)的每一个点，得到具有亚像素精度的同名点对{PD 1i(x,y)，PD 2i(x,y)}(i＝1,2,3,…,N)；
H.剔除{PD 1i(x,y)，PD 2i(x,y)}(i＝1,2,3,…,N)中误差较大的同名点对，获取最终的同名点对{PID 1i(x,y)，PID 2i(x,y)}(i＝1,2,3,…,S)。
其中，所述步骤D包括以下步骤：对于匹配区域内的影像数据，计算每个像素的HOG、LSS或者SURF等局部特征描述符，将每个像素对应的特征向量在Z方向进行排列，形成三 维的逐像素特征表达图。
进一步的，所述步骤D还可以为在匹配区域内构建方向梯度通道特征(CFOG)，具体包括以下步骤：
D1.对于匹配区域内的影像数据，计算每个像素在各个方向的梯度信息，形成三维的方向梯度图；
D2.在水平X方向和垂直Y方向，利用高斯滤波器对三维的方向梯度图做卷积运算，得到获得特征图
Figure PCTCN2018102271-appb-000001
再利用一维滤波器[1,2,1]在Z方向对特征图
Figure PCTCN2018102271-appb-000002
进行卷积运算，得到特征图
Figure PCTCN2018102271-appb-000003
D3.对特征图
Figure PCTCN2018102271-appb-000004
进行归一化操作，获得最终的方向梯度通道特征图。
其中所述步骤D为在匹配区域内构建方向梯度通道特征(CFOG)，具体计算过程包括如下步骤：
对于区域内的所有像素，分别利用一维滤波器[-1,0,1]和[-1,0,1] T计算它们在水平方向(X方向)和垂直方向(Y方向)的梯度g x和g y；
利用g x和g y计算它们在各个方向的梯度值g θ，计算公式如下：
Figure PCTCN2018102271-appb-000005
式中，θ表示量化的梯度方向，abs表示取绝对值，
Figure PCTCN2018102271-appb-000006
符号表示值为正时取本身，否则取0；
将各个方向的g θ叠置在一起，形成三维方向梯度图g o，然后在X和Y方向利用标准为σ的二维高斯滤波器对g o进行卷积运算获得特征图
Figure PCTCN2018102271-appb-000007
再利用一维滤波器[1,2,1]在Z方向对
Figure PCTCN2018102271-appb-000008
进行卷积运算得到特征图
Figure PCTCN2018102271-appb-000009
特征图
Figure PCTCN2018102271-appb-000010
中的每一个像素在Z方向上都对应了一个特征向量f i，遍历每个像素，对其特征向量v i进行归一化操作，得到最终的方向梯度通道特征图，归一化的计算公式如下：
Figure PCTCN2018102271-appb-000011
式中，ε是一个避免除零的数。
其中，所述梯度方向在360°的范围内均匀划分为18个等份，每个等份的角度间隔为20°，θ的取值为{0°，20°，……，340°}。
进一步的，所述步骤E为利用三维快速傅里叶变换将方逐像素特征表达图转换到频率域，并进行相关运算获得相似性图，取相似性图最大值的位置为影像的匹配位置。
其中所述步骤E具体包括以下计算步骤：
经过步骤D后分别得到区域AreaW 1i和AreaW 2i的逐像素的特征表达图D 1和D 2，将D 1作为模板在D 2进行滑动，利用D 1和D 2之间的特征向量差平方和作为相似性测度进行匹配；
D 1和D 2之间的差平方和计算公式为：
Figure PCTCN2018102271-appb-000012
式中，c表示特征表达图中的像素坐标，v表示D 1和D 2之间的偏移量，S i表示D 1和D 2之间的特征向量差平方和，当S i取得最小值时，将获得D 1和D 2之间的偏移量v i，即匹配位置，计算公式为：
Figure PCTCN2018102271-appb-000013
对公式(4)进行展开得:
Figure PCTCN2018102271-appb-000014
在公式(5)中，由于第一项和第二项的值接近于常数，所以当第三项的值最大时，公式(5)将获得最小值，因此，相似性函数可重新定义为：
Figure PCTCN2018102271-appb-000015
式中，
Figure PCTCN2018102271-appb-000016
是一个卷积运算；
考虑到频率域下的点乘运算等同于空间域下的卷积运算，因此在频率域利用快速傅里叶变换来提高其计算效率，得到基于傅里叶变换的相似性函数为：
Figure PCTCN2018102271-appb-000017
式中，F和F -1分别表示快速傅里叶正向变换和逆向变换，F *表示F的复数共轭。由于D 1和D 2是三维特征表达图，根据卷积的原理，需要利用三维快速傅里叶变换来计算公式(7)，得到最终的相似性函数为：
Figure PCTCN2018102271-appb-000018
式中，3DF和3DF -1分别表示三维快速傅里叶正向变换和逆向变换，3DF *表示3DF的复数共轭。
另一方面，本发明还提供了一种快速鲁棒的多模态遥感影像匹配系统，该系统包括下列单元：
预处理单元，用于判断参考影像和输入影像的分辨率信息，如果两幅影像具有相同的分辨率，则进入后续单元，如果分辨率不同，则将两幅影像采样为同样的分辨率；
模板区域选取单元，用于采用分块的策略，在参考影像上检测出一系列分布均匀的特征点，记为P 1i(i＝1,2,3,…,N)，以点P 1i为中心选取模板区域AreaW 1i；
匹配区域选取单元，用于根据遥感影像自身提供的地理坐标信息，预测点集P 1i(i＝1,2,3,…,N)在输入影像上所对应的匹配区域AreaW 2i；
特征提取单元，用于在匹配区域内构建逐像素特征表达图；
初匹配单元，用于在逐像素特征表达图的基础上，利用三维傅里叶变换建立一种快速相似性测度进行同名点匹配；对于获得的同名点，对其相似性图进行局部极值拟合，求解出匹配点的亚像素位置；重复上述单元的操作,遍历P 1i(i＝1,2,3,…,N)的每一个点，得到具有亚像素精度的同名点对{PD 1i(x,y)，PD 2i(x,y)}(i＝1,2,3,…,N)；
匹配点筛选单元，用于剔除{PD 1i(x,y)，PD 2i(x,y)}(i＝1,2,3,…,N)中误差较大的同名点对，获取最终的同名点对{PID 1i(x,y)，PID 2i(x,y)}(i＝1,2,3,…,S)。
进一步的，所述特征提取单元用于对于匹配区域内的影像数据，计算每个像素的局部特征描述符，将每个像素所对应的特征向量在Z方向进行排列，形成三维的逐像素特征表达图。
进一步的，所述初匹配单元是利用三维快速傅里叶变换将逐像素特征表达图转换到频率域，并进行相关运算获得相似性图，取相似性图最大值的位置为影像的匹配位置。
综上所述，由于采用了上述技术方案，本发明的有益效果是：
(1)本发明构建了一种快速、鲁棒的多模态遥感影像匹配框架，通过提取每个像素的局部特征描述符(如HOG，LSS或者SURF等)，形成逐像素的特征表达图，能较好地反映多模态遥感影像间的共有结构、形状和纹理等属性，并在逐像素的特征表达图的基础上，利用三维傅里叶变换建立了一种快速的匹配相似性测度。该框架可以快速精确、自动的在多模态影像间获取大量分布均匀的同名点，可有效提高匹配的实际生产效率，满足业务化运行的需求，而且它是一种通用的框架，可以整合各种局部特征描述进行影像匹配(不局限 于HOG，LSS或SURF等描述符)。
(2)在构建的框架中，逐像素的特征表达是通过计算每个像素的HOG，LSS或SURF等局部特征描述形成的，是一种稠密的特征表达技术。这不同于传统的HOG，LSS或SURF等描述符，它们只是在一个稀疏的采样格网(不是针对每个像素)内进行特征构建，或者利用所提取特征点的邻域信息计算特征，是一种相对稀疏的特征表达方式。相比而言，本发明的逐像素特征表达技术能更好，更精确地反映多模态遥感影像间的共有结构、形状和纹理等属性，匹配性能更稳健，而结合本方法提出的基于三维傅里叶变换的快速相似性测度，可实现多模态影像间的快速匹配。
(3)本发明针对所构建的框架，利用影像的方向梯度信息构建了一种新颖的逐像素特征描述符-方向梯度通道特征(CFOG)，它在匹配效率和精度方面都优于逐像素的HOG、LSS和SURF等特征表达方式。
(4)本发明在逐像素的特征表达基础上，利用三维傅里叶变换，在频率域建立了一种快速的匹配相似性测度。相比于空间域常用的相似性测度如灰度差平方和、归一化相关系数和互信息等，计算效率更高，计算效果更优。
(5)大量的实验结果表明，对于平坦地区的影像，匹配的总体精度可达到1个像素以内，而对于山区和城区的影像，匹配的总体精度可达到2.5个像素以内。对于大尺寸的遥感影像(超过20000×20000像素)，可在30秒内完成影像匹配。而且与目前流行的商业遥感影像软件(ENVI和ERDAS)相比，所发明的方法在匹配精度和计算效率都具有优势。
附图说明
图1是本发明的整体流程图。
图2是本发明的逐像素特征表达示意图。
图3是本发明中方向梯度通道特征(CFOG)的构建过程。
具体实施方式
为了使本领域的人员更好地理解本发明的技术方案，下面结合本发明的附图，对本发明的技术方案进行清楚、完整的描述，基于本申请中的实施例，本领域普通技术人员在没有做出创造性劳动的前提下所获得的其它类同实施例，都应当属于本申请保护的范围。
如图1，为一种快速鲁棒的多模态遥感影像匹配方法，包括以下步骤：
步骤A，根据参考影像和输入影像的分辨率信息，判断两幅影像的分辨率是否一致，一致则进行后续处理，不一致则将这两幅影像按照同样分辨率的进行采样处理。
步骤B，采用分块的策略，利用Harris或Forstner算子，在参考影像上提取大量的分 布均匀的特征点，具体包括：
将参考影像划分为n×n个互不重叠的方形格网，在每个格网内，计算每个像素点的Harris或Forstner特征值，并将特征值进行从大到小排序，选取特征值较大的k个像素作为特征点。这样可以在每个格网内检测出k个特征点，整幅影像上则拥有n×n×k个特征点。n和k的值可根据实际需求进行设定。记参考影像上检测的特征点集为P 1i(i＝1,2,3,….,N)。
在其他实施例中也可以采用其他算子进行图像的特征提取，本发明对此不进行限定。
步骤C，根据遥感影像自身提供的地理坐标信息，预测点集P 1i(i＝1,2,3,….,N)在输入影像上所对应的匹配搜索区域，具体包括以下步骤：
(1)提取P 1i(i＝1,2,3,….,N)的一个点P 1i(x,y)，x和y表示点P 1i的图像坐标，以点P 1i(x,y)为中心选取大小为M×M的模板区域AreaW 1i，并获取该点所对应的地理坐标Geo i；
(2)根据地理坐标Geo i并结合输入影像的地理坐标信息，计算其在输入影像上所对应的图像坐标P 2i(x,y)，并以P 2i(x,y)为中心确定一个大小也为M×M的方形窗口作为匹配搜索区域AreaW 2i。
步骤D，在区域AreaW 1i和AreaW 2i内构建逐像素的特征表达图，参见图2。
在一个实施例中，步骤D包括以下步骤：
计算对于区域内的每个像素的HOG、LSS或者SURF等局部特征描述符，将每个像素所对应的特征向量在Z方向进行排列，形成三维的逐像素特征表达图。对于使用何种具体的局部特征描述符，本实施例中对此不做限定。
在另一个实施例中，步骤D为在区域AreaW 1i和AreaW 2i内构建方向梯度通道特征(Channel Feature of Orientated Gradient,CFOG)，参见图3，具体包括以下步骤：
(1)对于区域内的所有像素，分别利用一维滤波器[-1,0,1]和[-1,0,1] T计算它们在水平方向(X方向)和垂直方向(Y方向)的梯度g x和g y.
(2)利用g x和g y计算它们在各个方向的梯度值g θ，计算公式如下：
Figure PCTCN2018102271-appb-000019
式中，θ表示量化的梯度方向，这里将梯度方向在360°的范围内均匀划分为18个等份，每个等份的角度间隔为20°，因此θ的取值为{0°，20°，……，340°}。abs表示取绝对值，其目的是将梯度方向为[180°，360°)的梯度信息转换到[0°，180°)之间，这样可以减少由于多模态影像间梯度反向所造成的匹配影响。
Figure PCTCN2018102271-appb-000020
符号表示值为正时取本身，否则取0。
(3)将各个方向的g θ叠置在一起，形成是三维的方向梯度图g o。然后在X和Y方向利 用标准为σ的二维高斯滤波器对g o进行卷积运算获得特征图
Figure PCTCN2018102271-appb-000021
再利用一维滤波器[1,2,1]在Z方向对
Figure PCTCN2018102271-appb-000022
进行卷积运算得到特征图
Figure PCTCN2018102271-appb-000023
(4)特征图
Figure PCTCN2018102271-appb-000024
中的每一个像素在Z方向上都对应了一个特征向量f i。遍历每个像素，对其特征向量v i进行归一化操作，进一步消除光照变化的影响，得到最终的CFOG特征图。归一化的计算公式如下：
Figure PCTCN2018102271-appb-000025
式中，ε是一个避免除零的数。
步骤E，在逐像素特征表达图的基础上，利用三维傅里叶变换建立一种快速的相似性测度进行同名点匹配，具体包括以下步骤：
(1)经过步骤D后，可以分别得到区域AreaW 1i和AreaW 2i对应的逐像素特征表达图D 1和D 2，将D 1作为模板在D 2进行滑动，可利用它们之间的特征向量差平方和作为相似性测度进行匹配。下面将差平方和的公式进行简化，并利用三维快速傅里叶变换加速其匹配过程。
D 1和D 2之间的差平方和计算公式为：
Figure PCTCN2018102271-appb-000026
式中，c表示特征表达图中的像素坐标，v表示D 1和D 2之间的偏移量，S i表示D 1和D 2之间的特征向量差平方和。当S i取得最小值时，将获得D 1和D 2之间的偏移量v i，即匹配位置，计算公式为：
Figure PCTCN2018102271-appb-000027
对公式(4)进行展开得:
Figure PCTCN2018102271-appb-000028
在式(5)中，由于第一项和第二项的值接近于常数，所以当第三项的值最大时，式(5)将获得最小值。因此，相似性函数可重新定义为：
Figure PCTCN2018102271-appb-000029
式中，
Figure PCTCN2018102271-appb-000030
是一个卷积运算。
考虑到频率域下的点乘运算等同于空间域下的卷积运算，这里在频率域利用快速傅里叶变换来提高其计算效率。基于傅里叶变换的相似性函数为：
Figure PCTCN2018102271-appb-000031
式中，F和F -1分别表示快速傅里叶正向变换和逆向变换，F *表示F的复数共轭。
由于D 1和D 2是三维特征表达图，根据卷积的原理，需要利用三维快速傅里叶变换来计算公式(7)。最终的相似性函数为：
Figure PCTCN2018102271-appb-000032
式中，3DF和3DF -1分别表示三维快速傅里叶正向变换和逆向变换，3DF *表示3DF的复数共轭；
(2)利用公式(8)进行匹配的过程为，首先对D 1进行三维傅里叶变换得到3DF(D 1)，并对D 2进行三维傅里叶变换并取复数共轭得到3DF *(D 2)，然后将3DF(D 1)和3DF *(D 2)进行点乘运算，并对运算结果进行三维傅里叶逆变换将获得D 1和D 2之间的相似性图，相似性图最大值的位置则对应了D 1和D 2之间的偏移量v i，也就是点P 1i(x,y)和点P 2i(x,y)之间的偏移量。将v i在X和Y方向上的偏移量记为(Δx,Δy)，则与点P 1i(x,y)对应的同名点为P 2i(x-Δx,y-Δy)，记为P * 2i(x,y)，获得的同名点对则为{P 1i(x,y)，P * 2i(x,y)}。
步骤F，对以上的同名点对{P 1i(x,y)，P* 2i(x,y)}，利用二元二次多项式进行局部插值获
取亚像素精度，具体包括以下步骤：
(1)以点P* 2i(x,y)为中心选取3×3像素的局部小窗口，并统计窗口内所有像素的相似性值，然后根据最小二乘的原理，采用二元二次多项式建立相似性值与像素位置的对应关系；
(2)对二元二次多项式求偏导，求解出偏导为0的位置，即获得亚像素精度的同名点对{PD 1i(x,y)，PD 2i(x,y)}；
步骤G，重复步骤C-F，遍历P 1i(i＝1,2,3,….,N)中的每一个点，获得具有亚像素精度的同名点对{PD 1i(x,y)，PD 2i(x,y)}(i＝1,2,3,….,N)。
步骤H，剔除{PD 1i(x,y)，PD 2i(x,y)}(i＝1,2,3,….,N)中误差较大的同名点对，获取最终的同名点对{PID 1i(x,y)，PID 2i(x,y)}(i＝1,2,3,…,S)，并用于进行匹配，具体包括以下步骤：
(1)基于最小二乘的原理，利用{PD 1i(x,y)，PD 2i(x,y)}(i＝1,2,3,….,N)的坐标建立投影变换模型。
(2)计算同名点对的均方根误差和残差，并剔除残差最大的同名点对。
(3)重复以上两个步骤，直到均方根误差小于1.5个像素，得到最终的具有亚像素精度的同名点对{PID 1i(x,y)，PID 2i(x,y)}(i＝1,2,3,…,S)。
在另一个实施例中，本发明还提供了一种快速鲁棒的多模态遥感影像匹配系统，该系统包括下列单元：
预处理单元，用于判断参考影像和输入影像的分辨率信息，如果两幅影像具有相同的分辨率，则进入后续单元，如果分辨率不同，则将两幅影像采样为同样的分辨率；
模板区域选取单元，用于采用分块的策略，在参考影像上检测出一系列分布均匀的特征点，记为P 1i(i＝1,2,3,…,N)，以点P 1i为中心选取模板区域AreaW 1i；
匹配区域选取单元，用于根据遥感影像自身提供的地理坐标信息，预测点集P 1i(i＝1,2,3,…,N)在输入影像上所对应的匹配区域AreaW 2i；
特征提取单元，用于在匹配区域内构建逐像素特征表达图；
初匹配单元，用于在逐像素特征表达图的基础上，利用三维傅里叶变换建立一种快速相似性测度进行同名点匹配；对于获得的同名点，对其相似性图进行局部极值拟合，求解出匹配点的亚像素位置；重复上述单元的操作,遍历P 1i(i＝1,2,3,…,N)的每一个点，得到具有亚像素精度的同名点对{PD 1i(x,y)，PD 2i(x,y)}(i＝1,2,3,…,N)；
匹配点筛选单元，用于剔除{PD 1i(x,y)，PD 2i(x,y)}(i＝1,2,3,…,N)中误差较大的同名点对，获取最终的同名点对{PID 1i(x,y)，PID 2i(x,y)}(i＝1,2,3,…,S)。
进一步的，所述特征提取单元用于对于匹配区域内的影像数据，计算每个像素的局部特征描述符，将每个像素所对应的特征向量在Z方向进行排列，形成三维的逐像素特征表达图。
进一步的，所述初匹配单元是利用三维快速傅里叶变换将逐像素特征表达图转换到频率域，并进行相关运算获得相似性图，取相似性图最大值的位置为影像的匹配位置。
以上为本发明具体实施方式的说明，通过本发明所构建的匹配框架，可利用各种局部特征描述符(如HOG、LSS或SURF等)进行有效的逐像素特征表达，从而有效地反映影像间共有结构、形状以及纹理等属性，并在此基础上利用三维傅里叶变换，在频率域建立了一种快速的匹配相似性测度，相比于空间域常用的相似性测度如归一化相关系数和互信息等，计算效率更高。另外，本发明所构建的方向梯度通道特征(CFOG)，是一种新颖的逐像素特征表达技术，它在匹配效率和精度上都优于基于HOG、LSS和SURF等描述符的逐像素特征 表达技术。本发明的技术方案能弥补传统匹配方法对于多模态影像间非线性辐射差异较为敏感的不足，可有效地解决了可见光、红外、激光雷达、合成孔径雷达以及地图等多模态遥感数据的匹配难题。
本发明所提出的技术方案是一种通用的技术框架，可以整合各种局部特征描述符(不限于CFOG、HOG，LSS和SURF等)进行影像匹配。
本发明并不局限于前述的具体实施方式。本发明扩展到任何在本说明书中披露的新特征或任何新的组合，以及披露的任一新的方法或过程的步骤或任何新的组合。

Claims (10)

  1. 一种快速鲁棒的多模态遥感影像匹配方法，其特征在于包括下列步骤：
    A.判断参考影像和输入影像的分辨率信息，如果两幅影像具有相同的分辨率，则进行步骤B，如果分辨率不同，则将两幅影像采样为同样的分辨率；
    B.采用分块的策略，在参考影像上检测出一系列分布均匀的特征点，记为P 1i(i＝1,2,3,…,N)，以点P 1i为中心选取模板区域AreaW 1i；
    C.根据遥感影像自身提供的地理坐标信息，预测点集P 1i(i＝1,2,3,…,N)在输入影像上所对应的匹配区域AreaW 2i；
    D.在匹配区域内构建逐像素特征表达图；
    E.在逐像素特征表达图的基础上，利用三维傅里叶变换建立一种快速相似性测度进行同名点匹配；
    F.对于获得的同名点，对其相似性图进行局部极值拟合，求解出匹配点的亚像素位置；
    G..重复步骤C-F,遍历P 1i(i＝1,2,3,…,N)的每一个点，得到具有亚像素精度的同名点对{PD 1i(x,y)，PD 2i(x,y)}(i＝1,2,3,…,N)；
    H.剔除{PD 1i(x,y)，PD 2i(x,y)}(i＝1,2,3,…,N)中误差较大的同名点对，获取最终的同名点对{PID 1i(x,y)，PID 2i(x,y)}(i＝1,2,3,…,S)。
  2. 如权利要求1所述的多模态遥感影像匹配方法，其特征在于，所述步骤D包括如下步骤：对于匹配区域内的影像数据，计算每个像素的局部特征描述符，将每个像素对应的特征向量在Z方向进行排列，形成三维的逐像素特征表达图。
  3. 如权利要求2所述的多模态遥感影像匹配方法，其特征在于，所述局部特征描述符为HOG、LSS或SURF的一种。
  4. 如权利要求1所述的多模态遥感影像匹配方法，其特征在于，所述步骤D为在匹配区域内构建方向梯度通道特征，具体包括以下步骤：
    D1.对于匹配区域内的影像数据，计算每个像素在各方向的梯度信息，形成三维的方向梯度图；
    D2.在水平X方向和垂直Y方向，利用高斯滤波器对三维的方向梯度图做卷积运算，得到获得特征图
    Figure PCTCN2018102271-appb-100001
    再利用一维滤波器[1,2,1]在Z方向对特征图
    Figure PCTCN2018102271-appb-100002
    进行卷积运算，得到特征图
    Figure PCTCN2018102271-appb-100003
    D3.对特征图
    Figure PCTCN2018102271-appb-100004
    进行归一化操作，获得最终的方向梯度通道特征图。
  5. 如权利要求4所述的多模态遥感影像匹配方法，其特征在于，所述步骤D中的构建方向梯度通道特征，具体包括以下计算步骤：
    对于区域内的所有像素，分别利用一维滤波器[-1,0,1]和[-1,0,1] T计算它们在水平方向(X方向)和垂直方向(Y方向)的梯度g x和g y；
    利用g x和g y计算它们在各个方向的梯度值g θ，计算公式如下：
    Figure PCTCN2018102271-appb-100005
    式中，θ表示量化的梯度方向，abs表示取绝对值，
    Figure PCTCN2018102271-appb-100006
    符号表示值为正时取本身，否则取0；
    将各个方向的g θ叠置在一起，形成三维方向梯度图g o，然后在X和Y方向利用标准为σ的二维高斯滤波器对g o进行卷积运算获得特征图
    Figure PCTCN2018102271-appb-100007
    再利用一维滤波器[1,2,1]在Z方向对
    Figure PCTCN2018102271-appb-100008
    进行卷积运算得到特征图
    Figure PCTCN2018102271-appb-100009
    特征图
    Figure PCTCN2018102271-appb-100010
    中的每一个像素在Z方向上都对应了一个特征向量f i，遍历每个像素，对其特征向量v i进行归一化操作，得到最终的方向梯度通道特征图，归一化的计算公式如下：
    Figure PCTCN2018102271-appb-100011
    式中，ε是一个避免除零的数。
  6. 如权利要求1所述的多模态遥感影像匹配方法，其特征在于，所述步骤E是利用三维快速傅里叶变换将逐像素特征表达图转换到频率域，并进行相关运算获得相似性图，取相似性图最大值的位置为影像的匹配位置。
  7. 如权利要求6所述的多模态遥感影像匹配方法，其特征在于，所述步骤E具体包括以下计算步骤：
    经过步骤D后分别得到区域AreaW 1i和AreaW 2i的逐像素特征表达图D 1和D 2，将D 1作为模板在D 2进行滑动，利用D 1和D 2之间的特征向量差平方和作为相似性测度进行匹配；
    D 1和D 2之间的差平方和计算公式为：
    Figure PCTCN2018102271-appb-100012
    式中，c表示特征表达图中的像素坐标，v表示D 1和D 2之间的偏移量，S i表示D 1和 D 2之间的特征向量差平方和，当S i取得最小值时，将获得D 1和D 2之间的偏移量v i，即匹配位置，计算公式为：
    Figure PCTCN2018102271-appb-100013
    对公式(4)进行展开得:
    Figure PCTCN2018102271-appb-100014
    在公式(5)中，由于第一项和第二项的值接近于常数，所以当第三项的值最大时，式(5)将获得最小值，因此，相似性函数可重新定义为：
    Figure PCTCN2018102271-appb-100015
    式中，
    Figure PCTCN2018102271-appb-100016
    是一个卷积运算；
    考虑到频率域下的点乘运算等同于空间域下的卷积运算，因此在频率域利用快速傅里叶变换来提高其计算效率，得到基于傅里叶变换的相似性函数为：
    Figure PCTCN2018102271-appb-100017
    式中，F和F -1分别表示快速傅里叶正向变换和逆向变换，F *表示F的复数共轭。由于D 1和D 2是三维特征图，根据卷积的原理，需要利用三维快速傅里叶变换来计算公式(7)，得到最终的相似性函数为：
    Figure PCTCN2018102271-appb-100018
    式中，3DF和3DF -1分别表示三维快速傅里叶正向变换和逆向变换，3DF *表示3DF的复数共轭。
  8. 一种快速鲁棒的多模态遥感影像匹配系统，其特征在于包括下列单元：
    预处理单元，用于判断参考影像和输入影像的分辨率信息，如果两幅影像具有相同的分辨率，则进入后续单元，如果分辨率不同，则将两幅影像采样为同样的分辨率；
    模板区域选取单元，用于采用分块的策略，在参考影像上检测出一系列分布均匀的特征点，记为P 1i(i＝1,2,3,…,N)，以点P 1i为中心选取模板区域AreaW 1i；
    匹配区域选取单元，用于根据遥感影像自身提供的地理坐标信息，预测点集P 1i(i＝1,2,3,…,N)在输入影像上所对应的匹配区域AreaW 2i；
    特征提取单元，用于在匹配区域内构建逐像素特征表达图；
    初匹配单元，用于在逐像素特征表达图的基础上，利用三维傅里叶变换建立一种快速相似性测度进行同名点匹配；对于获得的同名点，对其相似性图进行局部极值拟合，求解出匹配点的亚像素位置；重复上述单元的操作,遍历P 1i(i＝1,2,3,…,N)的每一个点，得到具有亚像素精度的同名点对{PD 1i(x,y)，PD 2i(x,y)}(i＝1,2,3,…,N)；
    匹配点筛选单元，用于剔除{PD 1i(x,y)，PD 2i(x,y)}(i＝1,2,3,…,N)中误差较大的同名点对，获取最终的同名点对{PID 1i(x,y)，PID 2i(x,y)}(i＝1,2,3,…,S)。
  9. 如权利要求8所述的多模态遥感影像匹配系统，其特征在于，所述特征提取单元用于对于匹配区域内的影像数据，计算每个像素的局部特征描述符，将每个像素对应的特征向量在Z方向进行排列，形成三维的逐像素特征表达图。
  10. 如权利要求8所述的多模态遥感影像匹配方法，其特征在于，所述初匹配单元是利用三维快速傅里叶变换将逐像素特征表达图转换到频率域，并进行相关运算获得相似性图，取相似性图最大值的位置为影像的匹配位置。
PCT/CN2018/102271 2017-08-31 2018-08-24 一种快速鲁棒的多模态遥感影像匹配方法和系统 WO2019042232A1 (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
US16/639,046 US11244197B2 (en) 2017-08-31 2018-08-24 Fast and robust multimodal remote sensing image matching method and system

Applications Claiming Priority (2)

Application Number Priority Date Filing Date Title
CN201710773285.6A CN107563438B (zh) 2017-08-31 2017-08-31 一种快速鲁棒的多模态遥感影像匹配方法和系统
CN201710773285.6 2017-08-31

Publications (1)

Publication Number Publication Date
WO2019042232A1 true WO2019042232A1 (zh) 2019-03-07

Family

ID=60978464

Family Applications (1)

Application Number Title Priority Date Filing Date
PCT/CN2018/102271 WO2019042232A1 (zh) 2017-08-31 2018-08-24 一种快速鲁棒的多模态遥感影像匹配方法和系统

Country Status (3)

Country Link
US (1) US11244197B2 (zh)
CN (1) CN107563438B (zh)
WO (1) WO2019042232A1 (zh)

Cited By (9)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111091133A (zh) * 2019-11-05 2020-05-01 西安建筑科技大学 一种基于sift算法的青铜器金文图像识别方法
CN111797744A (zh) * 2020-06-28 2020-10-20 武汉大学 一种基于共现滤波算法的多模态遥感图像匹配方法
CN112241957A (zh) * 2020-11-03 2021-01-19 南京多基观测技术研究院有限公司 一种仿生高分遥感影像变化检测方法
CN112633304A (zh) * 2019-09-23 2021-04-09 中国科学院沈阳自动化研究所 一种鲁棒的模糊图像匹配方法
CN112883850A (zh) * 2021-02-03 2021-06-01 湖北工业大学 一种基于卷积神经网络的多视角空天遥感图像匹配方法
CN113344030A (zh) * 2021-05-11 2021-09-03 淮阴工学院 一种基于判决相关分析的遥感图像特征融合方法及系统
CN113591933A (zh) * 2021-07-07 2021-11-02 中国人民解放军海军航空大学 一种基于相关性度量的遥感图像变化检测方法及系统
CN113763274A (zh) * 2021-09-08 2021-12-07 湖北工业大学 一种联合局部相位锐度定向描述的多源图像匹配方法
CN116343051A (zh) * 2023-05-29 2023-06-27 山东景闰工程研究设计有限公司 一种基于遥感影像的地质环境监测方法及系统

Families Citing this family (32)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
GB2555675B (en) * 2016-08-05 2019-05-08 Secr Defence Method and apparatus for generating an enhanced digital image of a physical object or environment
CN107563438B (zh) 2017-08-31 2019-08-30 西南交通大学 一种快速鲁棒的多模态遥感影像匹配方法和系统
CN108537833B (zh) * 2018-04-18 2022-06-21 昆明物理研究所 一种红外图像快速拼接方法
US11164015B2 (en) * 2018-05-08 2021-11-02 Ford Global Technologies, Llc Simultaneous diagnosis and shape estimation from a perceptual system derived from range sensors
US11361511B2 (en) * 2019-01-24 2022-06-14 Htc Corporation Method, mixed reality system and recording medium for detecting real-world light source in mixed reality
CN110163896B (zh) * 2019-03-29 2023-02-03 西南电子技术研究所（中国电子科技集团公司第十研究所） Sar图像匹配方法
CN112258561B (zh) * 2019-07-22 2023-08-25 复旦大学 针对图像拼接的匹配点获取方法
CN111190151B (zh) * 2020-01-14 2021-10-08 中国石油大学(华东) 扫描模式下多模态小卫星sar的系统参数设计及发射功率优化方法
CN111666858A (zh) * 2020-05-29 2020-09-15 中国科学院地理科学与资源研究所 基于单木识别的森林遥感图像配准方法及系统
CN112070745B (zh) * 2020-09-08 2022-11-22 中国科学院空天信息创新研究院 一种快速有效的遥感影像分割非监督评价方法
CN112183585A (zh) * 2020-09-08 2021-01-05 西安建筑科技大学 一种基于多特征测度的青铜器铭文相似性度量方法
CN112149728B (zh) * 2020-09-22 2023-11-17 成都智遥云图信息技术有限公司 一种快速的多模态图像模板匹配方法
CN112381748A (zh) * 2020-11-16 2021-02-19 上海亨临光电科技有限公司 基于纹理特征点的太赫兹与可见光图像配准方法及装置
CN112381864A (zh) * 2020-12-08 2021-02-19 兰州交通大学 一种基于对极几何的多源多尺度高分辨率遥感影像自动配准技术
CN112396643A (zh) * 2020-12-08 2021-02-23 兰州交通大学 一种尺度不变特征与几何特征融合的多模态高分影像配准方法
CN112559534B (zh) * 2020-12-24 2021-11-26 山东省国土测绘院 一种遥感影像数据归档管理系统及方法
CN112634335A (zh) * 2020-12-25 2021-04-09 清华大学 面向非线性辐射畸变的稳健遥感影像特征点对提取方法
CN113326724B (zh) * 2021-02-07 2024-02-02 海南长光卫星信息技术有限公司 一种遥感影像变化检测方法、装置、设备及可读存储介质
CN112950680B (zh) * 2021-02-20 2022-07-05 哈尔滨学院 一种卫星遥感图像配准方法
CN113034556B (zh) * 2021-03-19 2024-04-16 南京天巡遥感技术研究院有限公司 一种频域相关半稠密遥感影像匹配方法
CN113343747B (zh) * 2021-03-30 2023-05-09 西南电子技术研究所（中国电子科技集团公司第十研究所） 多模态图像鲁棒匹配vns方法
CN113313002A (zh) * 2021-05-24 2021-08-27 清华大学 一种基于神经网络的多模态遥感影像特征提取方法
CN113256653B (zh) * 2021-05-25 2023-05-09 南京信息工程大学 一种面向高层地物的异源高分遥感影像配准方法
CN113408615B (zh) * 2021-06-16 2022-04-12 中国石油大学(华东) 一种基于光学卫星遥感影像的船只自动匹配方法
CN114119686A (zh) * 2021-11-24 2022-03-01 刘文平 空间布局相似计算的多源遥感影像配准方法
CN115953600B (zh) * 2023-03-08 2023-07-21 中国测绘科学研究院 基于多方向滤波通道特征的多模态影像匹配方法及系统
CN116049482B (zh) * 2023-03-31 2023-10-20 北京数慧时空信息技术有限公司 基于时空域的遥感影像统筹方法
CN116091706B (zh) * 2023-04-07 2023-06-20 山东建筑大学 多模态遥感影像深度学习匹配的三维重建方法
CN116521927B (zh) * 2023-06-30 2024-02-13 成都智遥云图信息技术有限公司 基于网络地图瓦片的遥感影像匹配方法与系统
CN117076704B (zh) * 2023-10-16 2024-01-23 北京中遥数创科技有限公司 一种图像质量分析评价方法及系统
CN117237684B (zh) * 2023-11-13 2024-01-16 人工智能与数字经济广东省实验室(深圳) 一种像素级别的区域匹配方法及装置
CN117761695B (zh) * 2024-02-22 2024-04-30 中国科学院空天信息创新研究院 基于自适应分区sift的多角度sar三维成像方法

Citations (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103020945A (zh) * 2011-09-21 2013-04-03 中国科学院电子学研究所 一种多源传感器的遥感图像配准方法
CN104021556A (zh) * 2014-06-13 2014-09-03 西南交通大学 一种基于几何结构相似性的异源遥感影像配准方法
CN104574347A (zh) * 2013-10-24 2015-04-29 南京理工大学 基于多源遥感数据的在轨卫星图像几何定位精度评价方法
CN104867126A (zh) * 2014-02-25 2015-08-26 西安电子科技大学 基于点对约束和三角形网的有变化区域的合成孔径雷达图像配准方法
CN105427298A (zh) * 2015-11-12 2016-03-23 西安电子科技大学 基于各向异性梯度尺度空间的遥感图像配准方法
CN106971404A (zh) * 2017-03-20 2017-07-21 西北工业大学 一种鲁棒surf无人机彩色遥感图像配准方法
CN107563438A (zh) * 2017-08-31 2018-01-09 西南交通大学 一种快速鲁棒的多模态遥感影像匹配方法和系统
CN107909018A (zh) * 2017-11-06 2018-04-13 西南交通大学 一种稳健的多模态遥感影像匹配方法和系统

Family Cites Families (11)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
GB2303513B (en) * 1995-07-21 1999-09-29 Innovision Res Ltd Improved correlation processing for motion estimation
US7471086B2 (en) * 2007-04-20 2008-12-30 General Electric Company Magnetic resonance imaging visualization method and system
US8537337B2 (en) * 2009-12-22 2013-09-17 Weyerhaeuser Nr Company Method and apparatus for analyzing tree canopies with LiDAR data
US8712156B2 (en) * 2010-01-10 2014-04-29 Bronstein Bronstein Kimmel Technologies Ltd. Comparison of visual information
US9390327B2 (en) * 2013-09-16 2016-07-12 Eyeverify, Llc Feature extraction and matching for biometric authentication
US20170243084A1 (en) * 2015-11-06 2017-08-24 The Regents Of The University Of California Dsp-sift: domain-size pooling for image descriptors for image matching and other applications
JP6245590B1 (ja) * 2016-06-20 2017-12-13 公立大学法人大阪市立大学 皮膚診断装置、皮膚状態出力方法、プログラムおよび記録媒体
GB2555675B (en) * 2016-08-05 2019-05-08 Secr Defence Method and apparatus for generating an enhanced digital image of a physical object or environment
US20180053040A1 (en) * 2016-08-19 2018-02-22 Umm Al-Qura University System and method for 3d local surface matching
US10607364B2 (en) * 2017-10-31 2020-03-31 Honeywell International Inc. Position and attitude determination method and system using edge images
WO2019199579A1 (en) * 2018-04-09 2019-10-17 The Board Of Trustees Of The Leland Stanford Junior University Method of in situ gene sequencing

Patent Citations (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN103020945A (zh) * 2011-09-21 2013-04-03 中国科学院电子学研究所 一种多源传感器的遥感图像配准方法
CN104574347A (zh) * 2013-10-24 2015-04-29 南京理工大学 基于多源遥感数据的在轨卫星图像几何定位精度评价方法
CN104867126A (zh) * 2014-02-25 2015-08-26 西安电子科技大学 基于点对约束和三角形网的有变化区域的合成孔径雷达图像配准方法
CN104021556A (zh) * 2014-06-13 2014-09-03 西南交通大学 一种基于几何结构相似性的异源遥感影像配准方法
CN105427298A (zh) * 2015-11-12 2016-03-23 西安电子科技大学 基于各向异性梯度尺度空间的遥感图像配准方法
CN106971404A (zh) * 2017-03-20 2017-07-21 西北工业大学 一种鲁棒surf无人机彩色遥感图像配准方法
CN107563438A (zh) * 2017-08-31 2018-01-09 西南交通大学 一种快速鲁棒的多模态遥感影像匹配方法和系统
CN107909018A (zh) * 2017-11-06 2018-04-13 西南交通大学 一种稳健的多模态遥感影像匹配方法和系统

Cited By (15)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN112633304A (zh) * 2019-09-23 2021-04-09 中国科学院沈阳自动化研究所 一种鲁棒的模糊图像匹配方法
CN112633304B (zh) * 2019-09-23 2023-07-25 中国科学院沈阳自动化研究所 一种鲁棒的模糊图像匹配方法
CN111091133A (zh) * 2019-11-05 2020-05-01 西安建筑科技大学 一种基于sift算法的青铜器金文图像识别方法
CN111091133B (zh) * 2019-11-05 2023-05-30 西安建筑科技大学 一种基于sift算法的青铜器金文图像识别方法
CN111797744A (zh) * 2020-06-28 2020-10-20 武汉大学 一种基于共现滤波算法的多模态遥感图像匹配方法
CN112241957A (zh) * 2020-11-03 2021-01-19 南京多基观测技术研究院有限公司 一种仿生高分遥感影像变化检测方法
CN112883850B (zh) * 2021-02-03 2023-06-09 湖北工业大学 一种基于卷积神经网络的多视角空天遥感图像匹配方法
CN112883850A (zh) * 2021-02-03 2021-06-01 湖北工业大学 一种基于卷积神经网络的多视角空天遥感图像匹配方法
CN113344030A (zh) * 2021-05-11 2021-09-03 淮阴工学院 一种基于判决相关分析的遥感图像特征融合方法及系统
CN113344030B (zh) * 2021-05-11 2023-11-03 淮阴工学院 一种基于判决相关分析的遥感图像特征融合方法及系统
CN113591933A (zh) * 2021-07-07 2021-11-02 中国人民解放军海军航空大学 一种基于相关性度量的遥感图像变化检测方法及系统
CN113591933B (zh) * 2021-07-07 2024-04-09 中国人民解放军海军航空大学 一种基于相关性度量的遥感图像变化检测方法及系统
CN113763274B (zh) * 2021-09-08 2023-06-06 湖北工业大学 一种联合局部相位锐度定向描述的多源图像匹配方法
CN113763274A (zh) * 2021-09-08 2021-12-07 湖北工业大学 一种联合局部相位锐度定向描述的多源图像匹配方法
CN116343051A (zh) * 2023-05-29 2023-06-27 山东景闰工程研究设计有限公司 一种基于遥感影像的地质环境监测方法及系统

Also Published As

Publication number Publication date
CN107563438A (zh) 2018-01-09
US11244197B2 (en) 2022-02-08
CN107563438B (zh) 2019-08-30
US20200226413A1 (en) 2020-07-16

Similar Documents

Publication Publication Date Title
WO2019042232A1 (zh) 一种快速鲁棒的多模态遥感影像匹配方法和系统
CN111028277B (zh) 基于伪孪生卷积神经网络的sar和光学遥感图像配准方法
Chen et al. Building change detection with RGB-D map generated from UAV images
Yang et al. Automated registration of dense terrestrial laser-scanning point clouds using curves
Li et al. LNIFT: Locally normalized image for rotation invariant multimodal feature matching
US9396583B2 (en) Method of modelling buildings on the basis of a georeferenced image
CN107909018B (zh) 一种稳健的多模态遥感影像匹配方法和系统
CN109903313B (zh) 一种基于目标三维模型的实时位姿跟踪方法
WO2017049994A1 (zh) 一种高光谱图像角点检测方法与系统
US10043097B2 (en) Image abstraction system
Palenichka et al. Automatic extraction of control points for the registration of optical satellite and LiDAR images
CN107767456A (zh) 一种基于rgb‑d相机的物体三维重建方法
Troglio et al. Automatic extraction of ellipsoidal features for planetary image registration
Lee et al. Accurate registration using adaptive block processing for multispectral images
Han et al. LiDAR point cloud registration by image detection technique
JP2023530449A (ja) 空中と地上の位置合わせのためのシステムおよび方法
Parmehr et al. Automatic registration of optical imagery with 3d lidar data using local combined mutual information
Jin et al. Registration of UAV images using improved structural shape similarity based on mathematical morphology and phase congruency
CN112767459A (zh) 基于2d-3d转换的无人机激光点云与序列影像配准方法
CN112991388B (zh) 基于光流跟踪预测和凸几何距离的线段特征跟踪方法
Zhang et al. Integrating smartphone images and airborne lidar data for complete urban building modelling
Barazzetti et al. Automatic registration of multi-source medium resolution satellite data
Tahoun et al. Satellite image matching and registration: A comparative study using invariant local features
Barazzetti et al. Change detection from very high resolution satellite time series with variable off-nadir angle
KR101775124B1 (ko) 위성영상의 위치 정확도 개선을 위한 자동영상처리시스템 및 방법

Legal Events

Date Code Title Description
121 Ep: the epo has been informed by wipo that ep was designated in this application

Ref document number: 18850986

Country of ref document: EP

Kind code of ref document: A1

NENP Non-entry into the national phase

Ref country code: DE

122 Ep: pct application non-entry in european phase

Ref document number: 18850986

Country of ref document: EP

Kind code of ref document: A1