CN112529945B - Multi-view three-dimensional ISAR scattering point set registration method - Google Patents
Multi-view three-dimensional ISAR scattering point set registration method Download PDFInfo
- Publication number
- CN112529945B CN112529945B CN202011290620.5A CN202011290620A CN112529945B CN 112529945 B CN112529945 B CN 112529945B CN 202011290620 A CN202011290620 A CN 202011290620A CN 112529945 B CN112529945 B CN 112529945B
- Authority
- CN
- China
- Prior art keywords
- point
- points
- point set
- pairs
- target
- Prior art date
- Legal status (The legal status 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 status listed.)
- Active
Links
Images
Classifications
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T7/00—Image analysis
- G06T7/30—Determination of transform parameters for the alignment of images, i.e. image registration
- G06T7/33—Determination of transform parameters for the alignment of images, i.e. image registration using feature-based methods
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/10—Image acquisition modality
- G06T2207/10032—Satellite or aerial image; Remote sensing
- G06T2207/10044—Radar image
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
- G06T2207/00—Indexing scheme for image analysis or image enhancement
- G06T2207/20—Special algorithmic details
- G06T2207/20048—Transform domain processing
- G06T2207/20056—Discrete and fast Fourier transform, [DFT, FFT]
Landscapes
- Engineering & Computer Science (AREA)
- Computer Vision & Pattern Recognition (AREA)
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- Image Analysis (AREA)
Abstract
The invention discloses a registration method of a multi-view three-dimensional ISAR scattering point set, which comprises the following steps: obtaining a source point set and a target point set of an object based on three-dimensional ISAR imaging; performing surface fitting on each point in the source point set and each point in the target point set respectively to extract a curvature value, and selecting points with curvature values meeting the sorting requirement from a plurality of neighborhood scales as characteristic points of the source point set and the target point set respectively; performing initial registration on the source point set and the target point set: combining the coordinate value root mean square error and the distance root mean square error evaluation function to obtain the best matching four point pairs, and substituting the coordinates of the matching point pairs into a singular value decomposition method to calculate the transformation relation among the point sets; and obtaining a global optimal solution by adaptively changing the iteration step length by utilizing the best matching four point pairs and an iterative nearest neighbor algorithm based on an adaptive threshold value so as to make the point set converge to the global optimal solution. The method can effectively find the matching point pair with higher matching degree and can improve the registration precision.
Description
Technical Field
The invention belongs to the technical field of three-dimensional ISAR imaging, and particularly relates to a registration method of a multi-view three-dimensional ISAR scattering point set.
Background
Three-dimensional ISAR (Inverse Synthetic Aperture Radar) imaging can obtain a three-dimensional ISAR scattering point set of a target, and compared with a two-dimensional image, the three-dimensional ISAR scattering point set can reflect the real size and shape of the target, so that necessary information support is provided for target parameter extraction and target identification.
A scattering point set obtained by three-dimensional ISAR imaging is composed of a series of scattered points containing coordinate information, the structure of the scattering point set is similar to that of a three-dimensional laser point cloud, the three-dimensional ISAR scattering point set is sparser, overlapping areas among different visual angles are fewer, and a registration method aiming at the multi-visual angle three-dimensional ISAR scattering point set is rarely appeared at present. The method develops the three-dimensional ISAR scattering point set registration method research on the basis of the existing research results of three-dimensional laser point cloud registration.
When the laser point cloud data are registered, firstly, the registration is realized by a human-computer interaction method. Before point cloud data are generated, special marks are pasted on the surface of a target, so that point cloud data obtained from different angles comprise mark points, corresponding points are conveniently searched to obtain a transformation matrix of the point cloud, and registration is achieved. However, this method is time consuming and requires object matching, and is not suitable for registration of non-cooperative objects in three-dimensional ISAR imaging. One type of method for realizing automatic registration of three-dimensional point cloud data is principal component analysis, and the method obtains three principal axis vectors by calculating a covariance matrix of the point cloud data so as to obtain a change matrix to realize registration. However, the method requires more overlapping areas among the point cloud data to be registered, and the registration effect is poor when the target is shielded.
The three-dimensional point cloud descriptor is one of the widely used three-dimensional point cloud registration. In the method, partial overlapping areas exist among all groups of point clouds, the proper descriptor can calculate the characteristics of each point and the adjacent areas, and the matching point pair set among all groups of point clouds can be obtained by comparing the characteristics, so that the transformation matrix among the point clouds is calculated, and the point cloud registration is realized. The three-dimensional point cloud descriptor needs to be robust and distinguishable, and the core problem is to find the three-dimensional point cloud descriptor which can effectively capture the point cloud local set characteristics and has invariant translation and rotation and how to find and match the corresponding points by using the selected descriptor.
Wherein, the 'point signature descriptor' represents the minimum distance from a point on the three-dimensional point cloud to a certain reference plane. The point signature descriptor has translation and rotation invariance and is a simple method for representing the local structure of the curved surface. The method is complex in calculation process, large in calculation amount when point signatures are constructed, sensitive to noise especially when reference vectors are extracted, and possibly invalid when a local curved surface is a plane or a spherical surface, so that corresponding points among point clouds in the areas cannot be found.
Johnson et al propose a point cloud registration method based on a curved local histogram (called a spin image). Each spinning image is a local surface descriptor, and can be obtained by arbitrary two-dimensional calculation in the three-dimensional cylindrical coordinates of a certain neighborhood point. The spin image is a two-dimensional histogram with rotational and translational invariance. When the matching point pairs of the two point clouds are searched, the matching point pairs can be judged and determined according to the correlation degree between the spinning images. Then, an initial transformation matrix can be calculated from the corresponding point pairs for point cloud registration. Carmichael et al applied this method to three-dimensional datasets that handle large non-uniform sampling. Correlators also propose optimization methods such as multi-resolution and spherical spin images, which can improve the time efficiency of the algorithm and adaptively select algorithm parameters (such as two-dimensional histogram width). The main disadvantage of this type of algorithm is the low resolving power, which can generate a large number of mis-matched pairs.
Feldmar et al propose a registration algorithm based on principal curvatures of points of a point cloud. The algorithm firstly calculates two principal curvatures and principal directions of each point in the point cloud, and then selects points with similar principal curvatures and principal directions between the point clouds as matching point pairs by setting a threshold value to obtain a three-dimensional transformation matrix. The method is also suitable under the condition that any prior information and topological structure about the point cloud data are unknown, but because the method carries out matching point pair screening according to the characteristics of the single points and the screening process is simple, more wrong matching point pairs can be obtained.
Most local surface descriptors are defined by a single point. Mian proposes a registration algorithm based on the descriptors of the surface tensor. The descriptors are defined at a pair of directed points, and the change of the surface is measured by a high-dimensional surface histogram (tensor), and each descriptor contains hundreds to thousands of elements. Tensor matching between point clouds is performed by using a modified set hashing algorithm. However, this tensor-based approach requires resampling the surface before constructing the local descriptors, and may erroneously change the surface topology, resulting in a decrease in registration accuracy.
The Iterative Closest Point (ICP) algorithm is a fine registration algorithm proposed by Besl and McKay. The objective of the algorithm is to improve the initial registration result iteration to improve the overall quality of the registration. The quality of registration is usually determined by the sum of the distances of the corresponding points in the two sets of point clouds. Since there may be many local minima in the parameter range searched by the ICP algorithm, which makes the algorithm extremely easy to converge to a local optimal solution, the initial transformation of the point cloud must be as close as possible to the true value, so that an initial registration (e.g. point signature method or spin image method) is required before using the ICP to optimize the configuration result.
In conclusion, the algorithm for realizing dense laser point cloud registration is mature, but the three-dimensional ISAR scattering point set is sparse, the point sets with different visual angles have few overlapped areas, the method for directly applying point cloud registration is easier to generate mismatching point pairs, and the precise registration is difficult to converge to the global optimal solution.
Disclosure of Invention
In order to solve the problems in the prior art, the invention provides a multi-view three-dimensional ISAR scattering point set registration method. The technical problem to be solved by the invention is realized by the following technical scheme:
the invention provides a registration method of a multi-view three-dimensional ISAR scattering point set, which comprises the following steps:
s1: obtaining a source point set and a target point set of an object based on three-dimensional ISAR imaging;
s2: performing surface fitting on each point in the source point set and each point in the target point set respectively to extract a curvature value, and selecting points with the curvature values meeting the sorting requirement from a plurality of neighborhood scales as feature points of the source point set and the target point set respectively;
s3: initially registering the source point set and the target point set: combining coordinate values with a distance root mean square error evaluation function to obtain optimal matching four point pairs, and substituting coordinates of the matching point pairs into a singular value decomposition method to calculate a transformation relation among point sets;
s4: and obtaining a global optimal solution by adaptively changing the iteration step length by utilizing the optimal matching four-point pairs and an iterative nearest neighbor algorithm based on an adaptive threshold value so as to converge the point set to the global optimal solution.
In one embodiment of the present invention, the S2 includes:
s2.1: for the source set S, according to a point p in the set i And k thereof a Total k of neighborhoods a +1 point, obtaining each item coefficient of the curved surface by adopting a least square quadric surface fitting method, and further calculating to obtain any point p i Main curvature K of 1 ，K 2 And a curvature value f pi ；
S2.2: filtering plane points: for a point p in the source point set S data i When the curvature value f of the point pi When the value is larger than the plane point screening threshold value xi, reserving, otherwise, filtering;
s2.3: extracting candidate feature points by neighborhood curvature sorting: for a point P in the source point set S i The point and k thereof a In the neighborhood, co-k a And +1 points are sorted from large to small according to the curvature value, and the sorting result is recorded as sr Q When point p is i Has a curvature value of sr Q When there are first h, p may be i If not, directly eliminating to obtain a feature point set as P a ；
S2.4: solving intersection of multiple neighborhood scale candidate points, and extracting feature points: selecting another domain scale, here let neighborhood scale k b =6, repeating the steps S2.1 to S2.3 to obtain a candidate feature point set P b And finally can be obtainedFeature point set P of source point set S s ＝P a ∩P b ；
S2.5: and for the target point set D, repeating the steps from S2.1 to S2.4 to obtain a characteristic point set P of the target point set D d 。
In one embodiment of the present invention, in the S2.3, the sorting result sr Q Comprises the following steps:
sr Q ＝sort(f Q )
wherein f is Q And (3) representing the curvature value of k +1 points in the Q set, and sort (·) representing an inverse ordering operator.
In one embodiment of the invention, in S2.3, the candidate point feature p i Satisfies the following conditions:
f pi ∈Top h (sr Q )h∈[1,k+1]
wherein, top h (sr Q ) Is sr Q The first h curvature values, f pi For candidate point feature p i The curvature value of (a).
In one embodiment of the present invention, the S3 includes:
s3.1: two point pairs are generated: set of feature points P in source set S s Two points are arbitrarily selected to form two point pairs { s i ,s j Is traversed by s i ,s j The point pairs are compared with a characteristic point set P in a target point set D d Selecting { d } from the set of potential matching points in (1) i ,d j Makes the distance between the two pairs of matching points root-mean-square s i -s j |-|d i -d j Minimum, | traverse the feature point set P s All the point pairs of (1) are combined, in O (n) 2 ) Get two-point pair set E under the complexity of 2 And E is 2 Arranging according to the ascending order of the DRMS of each point pair;
s3.2: two point pairs are combined into four point pairs: after sorting E 2 Sequentially selecting a group of two point pairs e i Then E is again 2 Middle traverse e i Two-point pair e without coincident points j Calculate e i 、e j Selecting the DRMS with the smallest E in the DRMS of the corresponding points of the middle four groups j Form a set of four point pairs { e i ,e j Find a set of four points at a timeTo time, at E 2 Removing all two point pairs containing the four points, and traversing E 2 Find all four-point pair sets E 4 A 1 is mixing E 4 The DRMS of the four point pairs of each group are arranged in an ascending order, and the four point pairs of the top 10 percent are taken to form a group E 4 ' A set;
s3.3: calculation of E 4 ' transformation relationship of each set of four point pairs { R, t }: will E 4 Substituting the four-point pairs of the four groups of coordinates into a singular value decomposition method to obtain the transformation relation { R, t } of the four-point pairs of the four groups;
s3.4: calculating CRMS errors among the transformed point sets according to the transformation relation of the four point pairs of each group: at E 4 Taking a set of four point pairs and corresponding transformation relations { R, t }, according to p d (i)＝Rp s (i) + t, transforming each point in the original point set S into the coordinate system of the target point set D, marking as S ', regarding one point S ' (i) in S ', using KD-tree to search the nearest neighbor point D (i) of S ' (i) in the target point set D as the corresponding point, traversing S ' to find the corresponding points of all the points, combining CRMS evaluation function, calculating the coordinate value root mean square Error Error transformed according to the set of four point pairs 2 (S', D), traverse E 4 ' calculating the root mean square error of the coordinate values corresponding to each four point pairs.
S3.5: selection of E 4 ' four point pairs with the smallest root mean square error of the coordinates are taken as the best matching four point pairs.
In one embodiment of the present invention, in S3.4, the coordinate value root mean square Error is 2 (S', D) is:
where n is the number of S' point sets.
In one embodiment of the present invention, the S4 includes:
s4.1: searching the nearest neighbor point of each point S' in the target point set D to form corresponding points, calculating the Euclidean distance of each group of corresponding points, and counting the distance to be less than a set threshold value lambda s The number of corresponding points of (1), is noted as c 1 ；
S4.2: computing correspondencesMean value mu and variance sigma of Euclidean distances of points, the distances are in the interval [ lambda ] 1 ,λ 2 ]Put the corresponding points in the set S i ,D i In which S is i Storing a source point, D, of the ith iteration corresponding points i Storage target point, λ 1 ＝μ-ασ，λ 2 = μ + β σ is a set adaptive distance threshold, where μ and σ are a mean value and a variance of distances of corresponding points of the source point set in the target point set, and α and β are both adjustment coefficients;
s4.3: calculation of S using singular value decomposition i ,D i The corresponding transformation relation { R } i ,t i And will P 1 'transforming to the coordinate system of the target point set D again, and recording as S';
s4.4: searching the nearest neighbor point of each point S' in the target point set D, and counting the distance to be less than lambda s The number of corresponding points of (1), is noted as c 2 ；
S4.5: judgment c 1 And c 2 Size of (c) if 1 ＜c 2 Let S ' = ' S ', c 1 ＝c 2 Repeating S4.2 to S4.5 if c 1 ≥c 2 Then the best registration position between the point sets has been reached and the iteration is ended.
Compared with the prior art, the invention has the beneficial effects that:
1. the registration method of the multi-view three-dimensional ISAR scattering point set solves the problem that three-dimensional ISAR scattering point sets with different views are difficult to register. The method provides a multi-scale neighborhood feature point extraction algorithm based on curvature sorting, extracts uniformly distributed feature points, solves the problem that the extracted feature points are distributed in a centralized manner and are difficult to reflect the overall features of a point set, and effectively improves the feature point extraction effect.
2. Aiming at the problem of point pair mismatching between two point sets due to less coincident areas, a matching point pair screening algorithm based on CRMS and DRMS is provided, the algorithm can effectively find matching point pairs with higher matching degree, and the initial registration precision is improved.
3. In order to solve the problem that the iterative nearest neighbor algorithm is easy to fall into local optimum in the fine registration stage, the iterative nearest neighbor algorithm based on the self-adaptive threshold is provided, and the global optimum solution is obtained by adaptively changing the iteration step length.
The present invention will be described in further detail with reference to the accompanying drawings and examples.
Drawings
Fig. 1 is a flowchart of a multi-view three-dimensional ISAR scattering point set registration method according to an embodiment of the present invention;
fig. 2a to 2F are graphs of feature point extraction effects of source point sets and target point sets of boeing 737, apache and F14 obtained by using a multi-scale neighborhood feature point extraction algorithm based on curvature sorting;
fig. 3a to 3c are graphs of the effect of the best matching four-point pair in the boeing 737, apache and F14 feature point sets obtained by the screening algorithm using the matching point pair based on CRMS and DRMS according to the embodiment of the present invention;
fig. 4a to 4c are initial registration result graphs of point sets of boeing 737, apache and F14 after a screening algorithm;
fig. 5a to 5c are graphs of fine registration results of the initial registration of the point sets of boeing 737, apache and F14 in sequence.
Detailed Description
To further illustrate the technical means and effects of the present invention adopted to achieve the predetermined object, the following detailed description will be made on a multi-view three-dimensional ISAR scattering point set registration method according to the present invention with reference to the accompanying drawings and the detailed description thereof.
The foregoing and other technical matters, features and effects of the present invention will be apparent from the following detailed description of the embodiments, which is to be read in connection with the accompanying drawings. The technical means and effects of the present invention adopted to achieve the predetermined purpose can be more deeply and specifically understood through the description of the specific embodiments, however, the attached drawings are provided for reference and description only and are not used for limiting the technical scheme of the present invention.
It is noted that, herein, relational terms such as first and second, and the like may be used solely to distinguish one entity or action from another entity or action without necessarily requiring or implying any actual such relationship or order between such entities or actions. Also, the terms "comprises," "comprising," or any other variation thereof, are intended to cover a non-exclusive inclusion, such that an article or device that comprises a list of elements does not include only those elements but may include other elements not expressly listed. Without further limitation, an element defined by the phrase "comprising one of 8230, and" comprising 8230does not exclude the presence of additional like elements in an article or device comprising the element.
Referring to fig. 1, fig. 1 is a flowchart of a multi-view three-dimensional ISAR scattering point set registration method according to an embodiment of the present invention. The method comprises the following steps:
s1: a source point set S and a target point set D of the object are obtained based on three-dimensional ISAR imaging.
Specifically, the echo signals are subjected to translational compensation, FRFT (fractional Fourier transform) is used for performing parameter estimation on each component linear frequency modulation signal, coordinates of scattering points are determined according to parameters, and a three-dimensional ISAR imaging result consisting of a series of scattering points is obtained and serves as a three-dimensional ISAR scattering source point set S. Similarly, three-dimensional ISAR scattering point sets of the same object from different perspectives are obtained as a target point set D.
S2: and performing surface fitting on each point in the source point set S and the target point set D to extract a curvature value, and selecting points with the curvature values meeting the sorting requirement from a plurality of neighborhood scales as the characteristic points of the source point set S and the target point set D respectively.
The specific substeps of step S2 include:
s2.1: for the source set S, according to a point p in the set i And k thereof a Total k of neighborhoods a +1 point, obtaining each item coefficient of the curved surface by adopting a least square quadratic surface fitting method, and further calculating to obtain any point p i Main curvature K of 1 ，K 2 And a curvature value f pi . Where k is selected a ＝4。
S2.2: and filtering out plane points.
Because a large number of points in a plane area exist in three-dimensional ISAR scattering point set data, the points are not characterizedObviously, the curvature calculation result is small and close to 0, the assistance to the characteristic point extraction process is not large, a large amount of time is wasted, and the plane points can be filtered by setting a sufficiently small curvature threshold. For a point p in the S data of the source point set i When the curvature value f of the point pi And keeping when the value is larger than the plane point screening threshold value xi, and filtering otherwise. In the present embodiment, the plane point filtering threshold ξ takes 0.1.
S2.3: and (5) sorting and extracting candidate characteristic points by neighborhood curvature.
After calculating the curvature value of each point, for one point P in the source point set S i Let this point and the points in its neighborhood (denoted as Q set, total k) a +1 point) are sorted from large to small according to the curvature value, and the sorting result is recorded as sr Q When point p is i Has a curvature value of sr Q When there are first h, p may be i If not, directly eliminating to obtain a feature point set as P a 。
Further, the ranking results are:
sr Q ＝sort(f Q )
wherein, f Q And (3) the curvature value of k +1 points in the Q set, and sort (·) is an inverse ordering operator.
In the present embodiment, the candidate point feature p i Satisfies the following conditions:
f pi ∈Top h (sr Q )h∈[1,k+1]
wherein, top h (sr Q ) Is sr Q The first h curvature values of (i.e. p) i Should be the point of Q where h is large before curvature, f pi As a candidate point feature p i The curvature value of (a). The parameter h may control the number of candidate feature points. Generally, the size of h is positively correlated with the number of candidate feature points obtained by extraction, and when the number of candidate feature points is large, a better registration effect is easily obtained, but the calculation amount is greatly increased. Therefore, the actual use needs to be flexibly adjusted according to the size of the point set and the number k of the neighborhoods. In this embodiment, h =2 is taken, that is, two values are selected as candidate feature values.
S2.4: and solving intersection of the multiple neighborhood scale candidate points, and extracting characteristic points.
When local curved surface summation is carried out, the number k of the neighborhoods influences the fitting effect, and further the calculation results of the curvature and the curvature of each point are changed. The change in curvature at each point in turn affects whether the point can become a characteristic point. But representative feature points should have higher curvature values in a plurality of neighborhood ranges. Therefore, by changing the value of the neighborhood number k, refitting the local curved surface and calculating the curvature of each point, obtaining a plurality of groups of candidate feature point sets through the ordering of the plane point green area and the neighborhood curvature, and obtaining the intersection to obtain the final feature point.
In this step, another domain scale is selected, where the neighborhood scale k is assigned b =6, repeating the steps S2.1 to S2.3 to obtain a candidate feature point set P b Finally, a feature point set P of the source point set S can be obtained s ＝P a ∩P b 。
S2.5: and for the target point set D, repeating the steps from S2.1 to S2.4 to obtain a characteristic point set P of the target point set D d 。
Referring to fig. 2a to 2F, fig. 2a to 2F are feature point extraction effect diagrams of source point sets and target point sets of a boeing 737, apache and F14 obtained by using a multi-scale neighborhood feature point extraction algorithm based on curvature sorting, wherein fig. 2a and 2b are feature point extraction effect diagrams of the source point set and the target point set of the boeing 737; FIG. 2c and FIG. 2d are graphs of the effect of feature point extraction for the Apache source point set and the target point set; fig. 2e and 2F are graphs of the feature point extraction effect of the F14 source point set and the target point set.
S3: performing an initial registration on the set of source points S and the set of target points D.
Combining Coordinate value Root Mean square Error (CRMS) and Distance Root Mean square Error evaluation function (DRMS) to obtain the best matching four point pairs; substituting the matched four-point pairs into a singular value decomposition method to calculate a transformation relation { R, t } between the source point set S and the target point set D, wherein R is a rotation matrix, t is a translation vector, and p is a substitution formula d (i)＝Rp s (i) + t transforming the coordinates of the source set S into the set D of destination points, where p s (i) Is the ith point three-dimensional coordinate column in the source point set SVector, p d (i) And (4) three-dimensional coordinate column vectors of the ith point in the target point set D.
The specific substeps of step S3 include:
s3.1: two point pairs are generated.
Set of feature points P in source set S s Two points are arbitrarily selected to form two point pairs { s i ,s j Is traversed by s i ,s j Point pair feature point set P in target point set D d The set of potential matching points in (1) is selected as { d } i ,d j Makes the distance between the two pairs of matching points root-mean-square s i -s j |-|d i -d j And | l is minimal. Traverse feature point set P s All the point pairs of (1) are combined, and can be in O (n) 2 ) Obtaining a two-point pair set E under the complexity of 2 And E is 2 The DRMS errors for each pair are sorted in ascending order.
S3.2: the two-point pairs combine four-point pairs.
After sorting E 2 Sequentially selecting a group of two point pairs e i Then E is again 2 In traverse e i Two-point pair e without coincident points j Calculating e i 、e j Selecting the DRMS with the smallest E j Form a set of four point pairs { e i ,e j }: each time a set of four-point pairs is found, it needs to be at E 2 Removing all two point pairs containing the four points, and traversing E 2 Find all four-point pair sets E 4 (at this time E 2 Empty). Will E 4 The distances of the four point pairs are arranged in ascending order of root mean square error. In this embodiment, four-point pairs in each group are arranged in order of DRMS values, and the four-point pair composition E of the top 10% is taken 4 ' set.
S3.3: calculation of E 4 ' the transformation relationship of each set of four point pairs { R, t }. Will E 4 And substituting the four-point pairs of the four groups of coordinates into a singular value decomposition method to obtain the transformation relation { R, t } of the four-point pairs of the four groups.
S3.4: and calculating CRMS errors among the point sets after transformation according to the transformation relation of the four point pairs of each group.
Specifically, at E 4 ' taking a group of four point pairs and corresponding transformation relation { RT, transforming each point in the source point set S into the coordinate system of the target point set D, and marking as a transformation point set S' of the source point set S. For one point S '(i) (three-dimensional coordinate vector) in S', the nearest neighbor point D (i) of S '(i) is searched for as a corresponding point in the target point set D using the KD-tree, and corresponding points of all the points are found by traversing S'. And calculating the root mean square error of the coordinate values after transformation according to the set of four-point pairs by combining the CRMS evaluation function.
Further, the coordinate value root mean square Error 2 (S', D) is:
where n is the number of S' point sets.
S3.5: selection of E 4 ' four point pairs with the smallest root mean square error of the coordinates are taken as the best matching four point pairs.
It should be noted that repeating S3.1-S3.5 can generate pairs E between point sets 8 Or even sixteen to E 16 However, the number of the matching point pairs is not doubled, the algorithm complexity is increased, and the transformation relation among the point sets can be determined by four pairs of matching points, so that the four-point pair set E is obtained by the method 4 The continued combining is stopped.
The results of this embodiment are shown in fig. 3a to 3c, and fig. 3a to 3c are graphs of the effect of the best matching four-point pair in the boeing 737, apache and F14 feature point sets obtained by the matching point pair screening algorithm based on CRMS and DRMS according to this embodiment of the present invention, and the source point is slightly shifted to better show the matching effect of the point pair.
S3.6: substituting the optimal matching four point pairs into a singular value decomposition method again to calculate a transformation relation { R, t } among the point sets, transforming each point in the source point set S into a coordinate system of the target point set D, wherein the initial registration result refers to fig. 4a to 4c, and fig. 4a to 4c are initial registration result graphs of the point sets of the Boeing 737, the Apache and the F14 after a screening algorithm.
And 4, step 4: performing fine registration on the source point set S and the target point set D: and obtaining a global optimal solution by adaptively changing the iteration step length by utilizing the optimal matching four point pairs and an iterative nearest neighbor algorithm based on an adaptive threshold value so as to make the point set converge to the global optimal solution.
Setting an adaptive distance threshold λ 1 ＝μ-ασ、λ 2 = μ + β σ, where λ 1 Is used for removing error corresponding points with larger distance difference in the nearest neighbor corresponding point set during each iterative registration, namely lambda 2 And the method is used for screening out corresponding points with small distances and increasing the iteration step length so that the point set can be converged to the global optimal solution more quickly. Mu and sigma are respectively the mean value and variance of the distance of the corresponding points of the source point set in the target point set, alpha and beta are both regulating coefficients, and the experience value is taken. In this example, α =0.4 and β =0.2 are taken.
The specific sub-steps of step 4 include:
s4.1: and searching the nearest neighbor point of each point S' in the target point set D to form a corresponding point. Calculating Euclidean distance of corresponding points of each group, and counting the distance to be smaller than a threshold lambda s The number of corresponding points of (1), is noted as c 1 In this embodiment, λ s ＝0.05。
S4.2: calculating the mean value mu and the variance sigma of the Euclidean distance of the corresponding points, and enabling the distance to be in the interval [ lambda ] 1 ,λ 2 ]Put the corresponding points in the set S i ,D i In which S is i Storing a source point, D, of the ith iteration corresponding points i The target point is stored.
S4.3: calculation of S using singular value decomposition i ,D i Corresponding transformation relation R i ,t i And transforming S 'into the coordinate system of D again, and recording as S'.
S4.4: counting that the distance between the S' and the D is less than lambda s Number of corresponding points of (2), noted as c 2 In this embodiment, the threshold λ s ＝0.05。
S4.5: if c is 1 ＜c 2 That is, the number of corresponding points increases after iteration, indicating that the point set is proceeding in a direction of convergence, and making S '= S', c 1 ＝c 2 And repeating S4.2 to S4.5. If c is 1 ≥c 2 And reducing the number of corresponding points, namely the iteration reaches the optimal registration position among the point sets, and ending the iteration. Please refer to fig. 5a to 5c, which5a to 5c are graphs of fine registration results of point sets of boeing 737, apache and F14 after initial registration.
The registration method of the multi-view three-dimensional ISAR scattering point set solves the problem that three-dimensional ISAR scattering point sets with different view angles are difficult to register, the method provides a multi-scale neighborhood feature point extraction algorithm based on curvature sorting, extracts uniformly distributed feature points, solves the problem that the extracted feature points are concentrated in distribution and difficult to reflect the overall features of the point set, and effectively improves the feature point extraction effect. Aiming at the problem of point pair mismatching between two point sets due to less coincident areas, a matching point pair screening algorithm based on CRMS and DRMS is provided, the algorithm can effectively find matching point pairs with higher matching degree, and the initial registration precision is improved. In order to solve the problem that the iterative nearest neighbor algorithm is easy to fall into local optimum in the fine registration stage, the iterative nearest neighbor algorithm based on the self-adaptive threshold is provided, and the global optimum solution is obtained by adaptively changing the iteration step length.
The foregoing is a further detailed description of the invention in connection with specific preferred embodiments and it is not intended to limit the invention to the specific embodiments described. For those skilled in the art to which the invention pertains, several simple deductions or substitutions can be made without departing from the spirit of the invention, and all shall be considered as belonging to the protection scope of the invention.
Claims (5)
1. A multi-view three-dimensional ISAR scattering point set registration method is characterized by comprising the following steps:
s1: obtaining a source point set and a target point set of an object based on three-dimensional ISAR imaging;
s2: performing surface fitting on each point in the source point set and each point in the target point set respectively to extract a curvature value, and selecting points with curvature values meeting the sorting requirement from a plurality of neighborhood scales as feature points of the source point set and the target point set respectively;
s3: performing initial registration on the source point set and the target point set: combining coordinate values with a distance root mean square error evaluation function to obtain optimal matching four point pairs, and substituting coordinates of the matching point pairs into a singular value decomposition method to calculate a transformation relation among point sets;
s4: obtaining a global optimal solution by adaptively changing the iteration step length by utilizing the optimal matching four point pairs and an iterative nearest neighbor algorithm based on an adaptive threshold value so as to make the point set converge to the global optimal solution;
the S3 comprises the following steps:
s3.1: two point pairs are generated: set of feature points P in set of source points S s Two points are arbitrarily selected to form two point pairs { s i ,s j Is traversed by s i ,s j Point pair feature point set P in target point set D d Selecting { d } from the set of potential matching points in (1) i ,d j Makes the distance between the two pairs of matching points root-mean-square s i -s j |-|d i -d j Minimum, | traverse the feature point set P s All pairs of points of (2) are combined, at O (n) 2 ) Obtaining a two-point pair set E under the complexity of 2 And E is 2 Arranging according to the DRMS ascending order of each point pair;
s3.2: two point pairs are combined into four point pairs: after sorting E 2 Sequentially selecting a group of two point pairs e i Then, at E 2 Middle traverse e i Two-point pair e without coincident points j Calculate e i 、e j Selecting the DRMS with the smallest E j Form a set of four point pairs { e i ,e j Every time a set of four-point pairs is found, at E 2 Removing all two point pairs containing the four points, and traversing E 2 Find all four-point pair sets E 4 A 1 is mixing E 4 The DRMS of the four point pairs of each group are arranged in an ascending order, and the four point pairs of the top 10 percent are taken to form a group E 4 ' set;
s3.3: calculation of E 4 ' the transformation relation of each set of four point pairs { R, t }: will E 4 Substituting the four-point pairs of coordinates into a singular value decomposition method to obtain the transformation relation { R, t } of the four-point pairs;
s3.4: calculating CRMS errors among the transformed point sets according to the transformation relation of the four point pairs of each group: at E 4 Taking a set of four point pairs and corresponding transformation relations { R, t }, according to p d (i)＝Rp s (i) + t, transforming each point in the original point set S into the coordinate system of the target point set D, marking as S ', regarding one point S ' (i) in S ', using KD-tree to search the nearest neighbor point D (i) of S ' (i) in the target point set D as the corresponding point, traversing S ' to find the corresponding points of all the points, combining CRMS evaluation function, calculating the coordinate value root mean square Error Error transformed according to the set of four point pairs 2 (S', D), traverse E 4 ' calculating coordinate value root mean square error corresponding to each four point pairs;
s3.5: selection of E 4 Taking the four point pairs with the minimum mean square error of the central coordinates as the best matching four point pairs;
the S4 comprises the following steps:
s4.1: searching nearest neighbor points of S' points in the target point set D to form corresponding points, calculating Euclidean distances of the corresponding points of each group, and counting the distance to be smaller than a set threshold lambda s The number of corresponding points of (1), is noted as c 1 ；
S4.2: calculating the mean value mu and the variance sigma of the Euclidean distance of the corresponding points, and enabling the distance to be in the interval [ lambda ] 1 ,λ 2 ]Put the corresponding points in the set S i ,D i In which S is i Storing a source point, D, of the ith iteration corresponding points i Storage target point, λ 1 ＝μ-ασ，λ 2 = μ + β σ is a set adaptive distance threshold, where μ and σ are a mean value and a variance of distances of corresponding points of the source point set in the target point set, and α and β are both adjustment coefficients;
s4.3: calculation of S using singular value decomposition i ,D i Corresponding transformation relation R i ,t i And will P 1 'transforming to the coordinate system of the target point set D again, and recording as S';
s4.4: searching the nearest neighbor point of each point S' in the target point set D, wherein the statistical distance is less than lambda s Number of corresponding points of (2), noted as c 2 ；
S4.5: judgment c 1 And c 2 Size of (c) if 1 <c 2 Let S' = S ", c 1 ＝c 2 Repeating S4.2 to S4.5 if c 1 ≥c 2 Then the best registration position between the point sets has been reached, and the stack is endedAnd (4) generation.
2. The multi-view three-dimensional ISAR scattering point set registration method of claim 1, wherein the S2 comprises:
s2.1: for the source set S, according to a point p in the set i And k thereof a Total k of neighborhoods a +1 point, obtaining each item coefficient of the curved surface by adopting a least square quadric surface fitting method, and further calculating to obtain any point p i Main curvature K of 1 ，K 2 And a curvature value f pi ；
S2.2: filtering out plane points: for a point p in the source point set S data i When the curvature value f of the point pi When the value is larger than the plane point screening threshold value xi, reserving, otherwise, filtering;
s2.3: and (3) extracting candidate feature points by neighborhood curvature sorting: for a point P in the source point set S i The point and k thereof a In-neighborhood co-k a And +1 points are sorted from large to small according to the curvature value, and the sorting result is recorded as sr Q When point p is i Has a curvature value of sr Q When there are first h, p can be i If not, directly eliminating to obtain a feature point set as P a ；
S2.4: solving intersection of the multiple neighborhood scale candidate points, and extracting characteristic points: selecting another neighborhood dimension, here let neighborhood dimension k b =6, repeating the steps S2.1 to S2.3 to obtain a candidate feature point set P b Finally, a feature point set P of the source point set S can be obtained s ＝P a ∩P b ；
S2.5: and for the target point set D, repeating the steps from S2.1 to S2.4 to obtain a characteristic point set P of the target point set D d 。
3. The multi-view three-dimensional ISAR scattering point set registration method of claim 2, wherein in S2.3, the sorting result sr is Q Comprises the following steps:
sr Q ＝sort(f Q )
wherein, f Q Curvature for k +1 points in the Q setThe value sort (·) is the reverse order sort operator.
4. The multi-view three-dimensional ISAR scattering point set registration method of claim 2, wherein in S2.3, the candidate point feature p is i Satisfies the following conditions:
f pi ∈Top h (sr Q )h∈[1,k+1]
wherein, top h (sr Q ) Is sr Q The first h curvature values, f pi As a candidate point feature p i The curvature value of (a).
5. The multi-view three-dimensional ISAR scattering point set registration method of claim 4, wherein in S3.4, the coordinate values are root mean square Error (ERror) 2 (S', D) is:
where n is the number of S' point sets.
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202011290620.5A CN112529945B (en) | 2020-11-17 | 2020-11-17 | Multi-view three-dimensional ISAR scattering point set registration method |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202011290620.5A CN112529945B (en) | 2020-11-17 | 2020-11-17 | Multi-view three-dimensional ISAR scattering point set registration method |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN112529945A CN112529945A (en) | 2021-03-19 |
| CN112529945B true CN112529945B (en) | 2023-02-21 |
Family
ID=74980820
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202011290620.5A Active CN112529945B (en) | 2020-11-17 | 2020-11-17 | Multi-view three-dimensional ISAR scattering point set registration method |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN112529945B (en) |
Families Citing this family (3)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN113256692B (en) * | 2021-05-14 | 2024-04-02 | 西安交通大学 | Rigid body registration method and device based on self-adaptive neighborhood weight learning |
| CN113835082B (en) * | 2021-11-26 | 2022-03-11 | 中国电子科技集团公司信息科学研究院 | Method and device for distributed radar track association |
| CN114937122A (en) * | 2022-06-16 | 2022-08-23 | 黄冈强源电力设计有限公司 | Rapid three-dimensional model reconstruction method for cement fiberboard house |
Citations (18)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN102663753A (en) * | 2012-04-12 | 2012-09-12 | 浙江大学 | Body surface rebuilding method and body surface rebuilding device based on multiple visual angles of binocular stereoscopic vision |
| CN102880737A (en) * | 2012-07-30 | 2013-01-16 | 北京信息科技大学 | Workpiece registering method and system based on flexible assembly |
| CN103616688A (en) * | 2013-12-16 | 2014-03-05 | 中国科学院电子学研究所 | Method for improving quality of three-dimensional interferometric inverse synthetic aperture radar image |
| CN104134216A (en) * | 2014-07-29 | 2014-11-05 | 武汉大学 | Laser point cloud auto-registration method and system based on 16-dimension feature description |
| CN104535978A (en) * | 2014-12-19 | 2015-04-22 | 西安工程大学 | Three-dimensional InISAR image registration and fusion method based on mutual information |
| CN105279762A (en) * | 2015-11-20 | 2016-01-27 | 北京航空航天大学 | An oral cavity soft and hard tissue CT sequence and three-dimensional grid model registration method |
| CN105302874A (en) * | 2015-10-09 | 2016-02-03 | 苏州盛景信息科技股份有限公司 | Geographic cloud data based spatial matching method |
| CN106447708A (en) * | 2016-10-10 | 2017-02-22 | 吉林大学 | OCT eye fundus image data registration method |
| CN107886529A (en) * | 2017-12-06 | 2018-04-06 | 重庆理工大学 | A kind of point cloud registration method for three-dimensional reconstruction |
| CN108279643A (en) * | 2018-01-29 | 2018-07-13 | 西南交通大学 | A kind of workpiece posture method of adjustment based on measurement point Yu adaptive differential evolution algorithm |
| CN109146001A (en) * | 2018-09-14 | 2019-01-04 | 西安电子科技大学 | Multi-angle of view ISAR image interfusion method |
| CN109345620A (en) * | 2018-08-13 | 2019-02-15 | 浙江大学 | Merge the improvement ICP object under test point cloud method of quick point feature histogram |
| CN109472816A (en) * | 2018-09-17 | 2019-03-15 | 西北大学 | A kind of point cloud registration method |
| CN109541597A (en) * | 2018-12-12 | 2019-03-29 | 中国人民解放军国防科技大学 | Multi-station radar ISAR image registration method |
| CN109919984A (en) * | 2019-04-15 | 2019-06-21 | 武汉惟景三维科技有限公司 | A kind of point cloud autoegistration method based on local feature description's |
| CN110148165A (en) * | 2019-05-17 | 2019-08-20 | 西安电子科技大学 | A kind of three-dimensional interference ISAR method for registering images based on particle group optimizing |
| CN111080684A (en) * | 2019-12-12 | 2020-04-28 | 哈尔滨工程大学 | Point cloud registration method for point neighborhood scale difference description |
| CN111896954A (en) * | 2020-08-06 | 2020-11-06 | 华能澜沧江水电股份有限公司 | Corner reflector coordinate positioning method for shipborne SAR image |
-
2020
- 2020-11-17 CN CN202011290620.5A patent/CN112529945B/en active Active
Patent Citations (18)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN102663753A (en) * | 2012-04-12 | 2012-09-12 | 浙江大学 | Body surface rebuilding method and body surface rebuilding device based on multiple visual angles of binocular stereoscopic vision |
| CN102880737A (en) * | 2012-07-30 | 2013-01-16 | 北京信息科技大学 | Workpiece registering method and system based on flexible assembly |
| CN103616688A (en) * | 2013-12-16 | 2014-03-05 | 中国科学院电子学研究所 | Method for improving quality of three-dimensional interferometric inverse synthetic aperture radar image |
| CN104134216A (en) * | 2014-07-29 | 2014-11-05 | 武汉大学 | Laser point cloud auto-registration method and system based on 16-dimension feature description |
| CN104535978A (en) * | 2014-12-19 | 2015-04-22 | 西安工程大学 | Three-dimensional InISAR image registration and fusion method based on mutual information |
| CN105302874A (en) * | 2015-10-09 | 2016-02-03 | 苏州盛景信息科技股份有限公司 | Geographic cloud data based spatial matching method |
| CN105279762A (en) * | 2015-11-20 | 2016-01-27 | 北京航空航天大学 | An oral cavity soft and hard tissue CT sequence and three-dimensional grid model registration method |
| CN106447708A (en) * | 2016-10-10 | 2017-02-22 | 吉林大学 | OCT eye fundus image data registration method |
| CN107886529A (en) * | 2017-12-06 | 2018-04-06 | 重庆理工大学 | A kind of point cloud registration method for three-dimensional reconstruction |
| CN108279643A (en) * | 2018-01-29 | 2018-07-13 | 西南交通大学 | A kind of workpiece posture method of adjustment based on measurement point Yu adaptive differential evolution algorithm |
| CN109345620A (en) * | 2018-08-13 | 2019-02-15 | 浙江大学 | Merge the improvement ICP object under test point cloud method of quick point feature histogram |
| CN109146001A (en) * | 2018-09-14 | 2019-01-04 | 西安电子科技大学 | Multi-angle of view ISAR image interfusion method |
| CN109472816A (en) * | 2018-09-17 | 2019-03-15 | 西北大学 | A kind of point cloud registration method |
| CN109541597A (en) * | 2018-12-12 | 2019-03-29 | 中国人民解放军国防科技大学 | Multi-station radar ISAR image registration method |
| CN109919984A (en) * | 2019-04-15 | 2019-06-21 | 武汉惟景三维科技有限公司 | A kind of point cloud autoegistration method based on local feature description's |
| CN110148165A (en) * | 2019-05-17 | 2019-08-20 | 西安电子科技大学 | A kind of three-dimensional interference ISAR method for registering images based on particle group optimizing |
| CN111080684A (en) * | 2019-12-12 | 2020-04-28 | 哈尔滨工程大学 | Point cloud registration method for point neighborhood scale difference description |
| CN111896954A (en) * | 2020-08-06 | 2020-11-06 | 华能澜沧江水电股份有限公司 | Corner reflector coordinate positioning method for shipborne SAR image |
Non-Patent Citations (4)
| Title |
|---|
| Automatic registration of partial overlap three-dimensional surfaces;Yong Zhuo等;《2010 International Conference on Mechanic Automation and Control Engineering》;20100803;1-4 * |
| 三维点云数据配准方法研究;陈阳;《中国优秀硕士学位论文全文数据库 信息科技辑》;20150515;第2015年卷(第5期);I138-857 * |
| 基于激光雷达的车前地形复现关键技术研究;郭强;《中国优秀硕士学位论文全文数据库 工程科技Ⅱ辑》;20191015;第2019年卷(第10期);C035-12 * |
| 逆向工程中点云数据自动配准方法研究;金龙;《中国优秀硕士学位论文全文数据库 信息科技辑》;20170115;第2017年卷(第1期);I138-785 * |
Also Published As
| Publication number | Publication date |
|---|---|
| CN112529945A (en) | 2021-03-19 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| CN112529945B (en) | Multi-view three-dimensional ISAR scattering point set registration method | |
| CN109887015B (en) | Point cloud automatic registration method based on local curved surface feature histogram | |
| Kumar et al. | Image based leaf segmentation and counting in rosette plants | |
| CN110287873B (en) | Non-cooperative target pose measurement method and system based on deep neural network and terminal equipment | |
| JP5705147B2 (en) | Representing 3D objects or objects using descriptors | |
| JP6216508B2 (en) | Method for recognition and pose determination of 3D objects in 3D scenes | |
| Scharr et al. | Leaf segmentation in plant phenotyping: a collation study | |
| JP5677798B2 (en) | 3D object recognition and position and orientation determination method in 3D scene | |
| CN106778823B (en) | Automatic identification method for reading of pointer instrument | |
| CN111795704A (en) | Method and device for constructing visual point cloud map | |
| CN111553409B (en) | Point cloud identification method based on voxel shape descriptor | |
| CN113128610B (en) | Industrial part pose estimation method and system | |
| CN109697692B (en) | Feature matching method based on local structure similarity | |
| CN108376408A (en) | A kind of three dimensional point cloud based on curvature feature quickly weights method for registering | |
| CN110930456A (en) | Three-dimensional identification and positioning method of sheet metal part based on PCL point cloud library | |
| CN107862708A (en) | A kind of SAR and visible light image registration method | |
| CN109242774A (en) | A kind of flat-type part point cloud method based on hyperspace invariant features | |
| CN114358166B (en) | Multi-target positioning method based on self-adaptive k-means clustering | |
| CN116612165A (en) | Registration method for large-view-angle difference SAR image | |
| CN114332172A (en) | Improved laser point cloud registration method based on covariance matrix | |
| Srivastava et al. | Drought stress classification using 3D plant models | |
| CN106951873B (en) | Remote sensing image target identification method | |
| CN116704497B (en) | Rape phenotype parameter extraction method and system based on three-dimensional point cloud | |
| CN112884057B (en) | Point cloud data-based three-dimensional curved surface quality classification method and system and storage medium | |
| CN113554704A (en) | Electronic component positioning method based on improved SURF algorithm |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| PB01 | Publication | ||
| PB01 | Publication | ||
| SE01 | Entry into force of request for substantive examination | ||
| SE01 | Entry into force of request for substantive examination | ||
| GR01 | Patent grant | ||
| GR01 | Patent grant |