CN102779340B - Automatic corresponding method of feature point coordinates based on Delaunay triangulation - Google Patents

Abstract

The invention belongs to the field of machine vision and provides an automatic corresponding method of feature point coordinates based on Delaunay triangulation. The method comprises the following steps of: (1) binary segmenting and sorting an effective feature point image region; (2) calculating the pixel coordinate of the feature point through least square ellipse fitting; (3) performing Delaunay triangulation on the feature point set; (4) sorting effective subdivision triangles; (5) traversing the non-longest margins of the effective subdivision triangles and establishing a four-neighborhood information table of the feature point set; and (6) inquiring the four-neighborhood information table, searching and traversing the feature point set in four neighborhoods, calculating the world coordinate of each feature point during the traversing process to complete correspondence of the pixel coordinate of the feature point to the world coordinate. If an image has shooting inclined angle, distortion and feature point loss, the method provided by the invention can still effectively operate, so that the method is especially suitable for the online marking of an IC encapsulation vision positioning system.

Description

A kind of automatic corresponding method of unique point coordinate based on Delaunay triangulation

Technical field

The invention belongs to field of machine vision, relate to a kind of automatic corresponding method for video camera on-line proving process circular array target unique point coordinate.Unique point coordinate correspondence refers to the corresponding of unique point pixel coordinate and world coordinates.

Background technology

Vision location, because of advantages such as its noncontact, high precision, not damageds, is widely used on IC sealed in unit.Camera calibration is one of important step of vision location, for setting up contacting between image coordinate system and world coordinate system, proofreaies and correct various linearities, nonlinear images distortion.At present, conventional demarcation target has gridiron pattern type and round dot array type.Its orbicular spot array type target adopts circular image regional center point as unique point, has higher stated accuracy, on IC sealed in unit, is widely used.Target unique point array distribution and spacing (world coordinates) are known, its world coordinates easily obtains, yet unique point pixel extraction order is often inconsistent with the order of world coordinates, in order to realize calibration algorithm, must guarantee that the image coordinate of feature point for calibration is correctly corresponding with world coordinates.When taking uncalibrated image, existing method adopts the artificial mode of adjusting scaling board position just to make it to put with, and the ranks number of complicate statistics feature point mesh mostly.This has run counter to online, the full automatic service requirement of IC sealed in unit, has affected the work efficiency of equipment.In addition, in the situation that feature point mesh is densely distributed, also very easily there is mistake in the ranks number of complicate statistics feature point mesh.Therefore, exploitation has important theory significance and practical value without manual intervention and complicate statistics, the unique point coordinate corresponding method that adapts to the requirement of IC sealed in unit on-line operation.

For video camera on-line proving, reduce the demand of manual intervention, researcher has proposed the automatic corresponding method of various features point coordinate both at home and abroad: document " Fully automatic algorithm for region of interest location in camera calibration " (Optical Engineering, 2002, 41 (6): 1220-1226) propose the automatic correspondence that a kind of intelligent image of interest region method based on RADON conversion is realized unique point coordinate on array round dot target, but it can only be little at shooting inclination angle, effective under the little prerequisite of distortion (image non-linear distortion), document " Robust recognition or checkerboard pattern for camera calibration " (Optical Engineering, 2006,45 (9): 1-9) with article " characteristic points automatic extraction method based on round dot array target " (China Mechanical Engineering, 2010,21 (16): method 1906-1910) proposing can adapt to the situation that scaling board exists rotation, but need in target image region, add special warning triangle, Delaunay triangulation is the angle Optimized triangulations of a kind of in computational geometry theory " making the leg-of-mutton minimum angle of subdivision maximum ".If ρ is the arbitrary point set in plane,, in the inside of the leg-of-mutton circumscribed circle of each Delaunay subdivision of ρ, do not comprise any point in ρ.It is usually used to the division of planar point set, by the leg-of-mutton end points of traversal subdivision, can obtain " near point " information of any point in a set.Document " Automatic Grid Finding in Calibration Patterns Using Delaunay triangulation " (Technical Report, NRC-46487/ERB-1104, National Research Council, Canada, August, 2003) introduced first the automatic correspondence of Delaunay triangulation method realization character point coordinate, but the method is only demarcated target effectively and is depended on the azimuth information that three marked circle provide to gridiron pattern type.Document " a kind of automatic corresponding method of unique point in camera calibration " (photoelectron laser, 2011,22 (5): corresponding method 736-739) proposing is based on Delaunay triangulation, it allows feature point mesh to have incompleteness to a certain degree, but still needs to guarantee that it is complete that feature point mesh periphery has a limit at least.

Summary of the invention

Technical matters to be solved by this invention is: a kind of automatic corresponding method of unique point coordinate based on Delaunay triangulation is provided, the unique point coordinate corresponding method that the method proposes can automatically complete camera calibration without manual intervention in the situation that, simplify camera calibration process, met the demand of IC sealed in unit for video camera on-line proving.

A kind of automatic corresponding method of unique point coordinate based on Delaunay triangulation provided by the invention, is characterized in that, the method comprises the steps:

The 1st step: read uncalibrated image I, then it is carried out to binary segmentation, obtain bianry image I ';

The 2nd step: bianry image I ' is carried out to connected domain analysis, obtain set of image regions N, calculate the area value S of each image-region;

The 3rd step: extract the point coordinate of each image-region in set of image regions N, and calculate accordingly the circularities Δ R of each image-region, and radius value R; Bring Δ R and R into relative circularities Δ r that formula 2 calculates each image-region in N;

The 4th step: screening meets that " S is at [S in set of image regions N ₁, S ₂] in scope and Δ r < Δ r ₀" image-region, the image-region that obtains is formed to validity feature dot image regional ensemble N _e, wherein, S ₁, S ₂be respectively the high and low threshold value of predefined image-region area, Δ r ₀for the relative circularity threshold value of predefined image-region;

The 5th step: to validity feature dot image regional ensemble N _ein each validity feature dot image region n _ei, i representation feature point sequence number, the point coordinate extracting according to the 3rd step obtains an elliptic equation by the Fitting Calculation, and the elliptical center point being obtained by this equation calculation of parameter is image-region n _eicharacteristic of correspondence point P _i, its coordinate is (x _i, y _i);

The 6th step: establish by all unique point P _ithe point set forming is combined into P, and P is carried out to Delaunay triangulation, obtains subdivision triangle set τ;

The 7th step: for each subdivision triangle in subdivision triangle set τ, according to formula η=min (| ω ₁-90 ° |, | ω ₂-90 ° |, | ω ₃-90 ° |) calculate the leg-of-mutton form variations of subdivision angle η, wherein, ω ₁, ω ₂, ω ₃for the leg-of-mutton angle value of subdivision; Predefined form variations angle threshold value η ₀, in subdivision triangle set τ, screening meets " η≤η ₀" triangle, obtain effective subdivision triangle set τ _e;

The 8th step: travel through effective subdivision triangle set τ _ein each leg-of-mutton non-longest edge, the sequence number of two end points and deflection are joined respectively in neighbours' domain information structure of the other side, while finishing to traversal, can set up neighbours' domain information table of unique point set;

The 9th step: neighbours' domain information table of query characteristics point set, complete unique point P is put in an optional neighbours territory ₀as Seed Points, the set of neighbours territory recursive search traversal unique point; In ergodic process, determine current some P of every one deck recurrence _kwith respect to Seed Points P ₀feature point mesh coordinate (r _k, c _k), k represents the recurrence number of plies, obtains each unique point P in unique point set when traversal finishes _irelative mesh coordinate (r _i, c _i);

The 10th step: the row-coordinate minimum value r that finds the relative mesh coordinate of all unique points in unique point set P _min=min (r _i) and row coordinate minimum value c _min=min (c _i), calculate each unique point P _iworld coordinates (X _i, Y _i),

$\left\{\begin{array}{c}{X}_i=L\&times;(c_i-c_{\min })\\ Y_i=L\&times;(r_i-r_{\min })\end{array}\right.$

Wherein, L is the world coordinates value of unique point spacing.By each unique point P _ipixel coordinate (x _i, y _i) and world coordinates (X _i, Y _i) corresponding one by one, complete the correspondence of unique point coordinate.

It is automatically corresponding that the present invention can be applicable in IC sealed in unit on-line proving process circular array target unique point coordinate.In the on-line proving process of IC sealed in unit, scaling board is placed automatically by machine, due to the spacing of unique point often very little (being generally no more than 3mm), so be easy to cause Partial Feature dot image region to be blocked and form incomplete situation generation, when blocking when developing into feature point mesh periphery there is no a limit be complete, the method that background technology is mentioned is just no longer applicable.In addition, the method that background technology is mentioned also needs the ranks number of complicate statistics feature point mesh, and this has not only run counter to the requirement of equipment automatization operation but also has very easily made mistakes.ACDT provided by the invention (Automatic Correspondence based on Delaunay Triangulation) method has puted forth effort to solve above-mentioned problem under the prerequisite of inheriting previous method advantage, gets final product the correspondence of realization character point coordinate without any manual intervention.When there is the rotation of feature point mesh, incomplete problem because of scaling board positioning error, ACDT method still can work orderly, and in addition, the method is insensitive for the distortion that is no more than lens distortion peak demand for general measure (being no more than 1%).

The unique point coordinate corresponding method that the present invention proposes can automatically complete camera calibration without manual intervention and complicate statistics in the situation that, has simplified camera calibration process, has met the demand of IC sealed in unit for video camera on-line proving; Owing to having utilized array distribution feature and the local anti-distortion character of feature point set, the method for distortion, to take inclination angle insensitive, feature point mesh exist by target rotation, translation, caused block defect time still can work orderly; Than two kinds that mention in the background technology same unique point corresponding method based on Delaunay triangulation, the present invention proposes ACDT method and does not rely on the azimuth information that marked circle provides, without complicate statistics feature point mesh number, and in the situation that block very serious (it is complete that feature point mesh periphery does not have a limit), still can normally work, adaptability, robustness are better.

Accompanying drawing explanation

Fig. 1 represents the process flow diagram of ACDT method.

Fig. 2 (a) represents uncalibrated image, and Fig. 2 (b) represents the validity feature dot image region obtaining through screening.

Fig. 3 represents the process flow diagram of the validity feature dot image region filtering algorithm that example of the present invention proposes.

The Delaunay triangulation result of Fig. 4 (a) representation feature point set, Fig. 4 (b) represents effective subdivision triangle, Fig. 4 (c) represents the leg-of-mutton non-longest edge of effective subdivision.

Fig. 5 represents that feature point mesh neighbours domain information table that example of the present invention proposes sets up the process flow diagram of algorithm.

Fig. 6 represents to choose No. 10 points as Seed Points, the path of neighbours territory recurrence traversal unique point set.

Embodiment

The present invention utilizes the leg-of-mutton character of unique point set Delaunay subdivision to obtain each neighborhood of a point dot information, and the set of anti-distortion character neighbours territory, the part that utilizes feature point mesh recursive search traversal unique point, determine each point residing position on feature point mesh.

Two conditions of joint image region area of the present invention circularity relative to image-region are carried out the screening in validity feature dot image region.Wherein, image-region area can be approximately equal to the pixel number that forms this image-region, and the relative circularity of image-region can detect the minimum image region circle model about roundness evaluation in regulation > > according to GB1598-80 < < geometrical and toleranging and calculate.Existing similar approach is only screened by image-region area mostly, and other has added simple shape judgement, and it is accurate that the screening effect under the same terms is not so good as this method.

The present invention utilizes unique point and the common array distribution feature of image slices vegetarian refreshments, the neighbours territory recursive search ergodic algorithm that is usually used in image connectivity domain analysis is applied to the traversal of unique point set, and according to searching route, determines that each unique point is with respect to the position of Seed Points in ergodic process.

The present invention utilizes the anti-distortion in the part of feature point mesh, determines current some partial row's column direction angle and some position, current some neighbours territory.

The anti-distortion in part of feature point mesh: establish 1 P in feature point mesh ξ _i(i representation feature point sequence number) and its neighbours' territory point P _i1, P _i2, P _i3, P _i4the vector forming angle between middle adjacent vector is respectively α _i12, α _i23, α _i34, α _i41, according to the maximum distortion angle β of formula 1 defined feature dot grid ξ _max.Emulation experiment and application practice prove, is no more than measure with image capturing system peak demand (being not more than 1%) in the situation that β in pattern distortion _maxcan remain in a less scope (being no more than 10 °), this character is called the anti-distortion in part of feature point mesh.

β _max=max (| α _i12-90 ° |, | α _i23-90 ° |, | α _i34-90 ° |, | α _i41-90 ° |) formula 1

As shown in Figure 1, the concrete steps of the inventive method comprise:

The 1st step: read uncalibrated image I, cut apart I by binary segmentation method (as OSTU adaptive threshold dividing method etc.), obtain bianry image I ';

The 2nd step: to image I ' carry out connected domain (Blob) and analyze and to obtain set of image regions N, calculate the area value S of each image-region.

The 3rd step: use profile extraction algorithm (as eight neighborhood track algorithms etc.) to extract the point of each image-region in set of image regions N, and calculate accordingly the circularities Δ R of each image-region and radius value R thereof.Bring Δ R and R into relative circularities Δ r that formula 2 calculates each image-region in N.

$\mathrm{\&Delta;r}=\frac{\mathrm{\&Delta;R}}{R}$ Formula 2

Calculate circularities Δ R and radius value R and can stipulate in > > the minimum image region circule method about roundness evaluation according to GB1598-80 < < geometrical and toleranging-detections, also can adopt other circularity Modeling Calculation method, as minimum two-multiply law, minimum circumscribed circle method, minimum inscribe circule method etc.

The 4th step:, in set of image regions N, screening meets that " S is at [S ₁, S ₂] in scope and Δ r < Δ r ₀" image-region, the image-region that obtains is formed to validity feature dot image regional ensemble N _e.The high and low threshold value S of image-region area ₁, S ₂and the relative circularity threshold value of image-region Δ r ₀can set according to the imaging actual conditions of scaling board S ₁, S ₂generally get respectively pre-estimation unique point image-region area value obtain 1.2 times with 0.8 times, Δ r ₀span generally between 0.05 to 0.3.

The 5th step: to N _ein each validity feature dot image region n _ei(i representation feature point sequence number), obtains an elliptic equation according to its point coordinate (extract and obtain in the 3rd step) by the Fitting Calculation (as least square method), and the elliptical center point being obtained by this equation calculation of parameter is image-region n _eicharacteristic of correspondence point P _i, its image coordinate is (x _i, y _i).

The 6th step: establish by all unique point P _ithe point set forming is combined into P, and P is carried out to Delaunay triangulation, obtains subdivision triangle set τ.

The 7th step: set form variations angle threshold value η ₀(span is generally at 15 ° ~ 25 °), calculate the leg-of-mutton form variations of each subdivision angle η in τ according to formula 3, and in subdivision triangle set τ, screening meets " η≤η ₀" triangle, obtain effective subdivision triangle set τ _e, wherein, ω ₁, ω ₂, ω ₃for the leg-of-mutton angle value of subdivision;

η=min (| ω ₁-90 ° |, | ω ₂-90 ° |, | ω ₃-90 ° |) formula 3

The 8th step: travel through effective subdivision triangle set τ _ein each leg-of-mutton non-longest edge, to two of every limit end points Q ₁, Q ₂proceed as follows: by Q ₁some sequence number and deflection join Q ₂neighbours' domain information structure, and by Q ₂some sequence number and deflection join Q ₁neighbours' domain information structure.While finishing to traversal, can set up neighbours' domain information table of feature point mesh (being unique point set P).

Neighbours' domain information table is the structure array of a customization type, index sequence number in construction array has represented the unique point sequence number of its representative, wherein four long type variablees have recorded neighbours' territory point sequence number of point corresponding to current index sequence number, and other four double type variablees have recorded each neighbours territory and put corresponding deflection.

The 9th step: neighbours' domain information table of query characteristics dot grid, complete unique point P is put in an optional neighbours territory ₀as Seed Points, neighbours territory recursive search traversal feature point set.In ergodic process, determine current some P of every one deck recurrence _k(k represents the recurrence number of plies) is with respect to Seed Points P ₀feature point mesh coordinate (r _k, c _k) (being called for short " mesh coordinate relatively " below).While finishing to traversal, can obtain each unique point P in unique point set _irelative mesh coordinate (r _i, c _i).

The concrete grammar of neighbours territory recursive search traversal is as follows:

(9.1) calculate current some P of k layer recurrence _kpartial row's column direction angle: when k=0, the partial row of regulation Seed Points just, row just, negative, the row negative direction angle δ of row _0R+, δ _0C+, δ _0R-, δ _0C-be respectively 0 °, 90 °, 180 °, 270 °.When k > 0, establish from current some P of k-1 layer recurrence _kto P _kvector with level to the angle of right, be γ, according to P _krelative P _k-1residing position, calculates P according to formula 4 _kpartial row's column direction angle δ _kR+, δ _kC+, δ _kR-, δ _kC-, and they are normalized in the scope of 0 ~ 360 °.

formula 4

(9.2) determine current some P of k layer recurrence _ksome position, neighbours territory: according to formula 5, determine P _kneighbours' territory point P _kj(j represents neighbours' territory point sequence number) is with respect to P _kposition, ε wherein ₀for predefined deviation of directivity threshold value (generally in 20 °).

$\left\{\begin{array}{c}|{\mathrm{\&theta;}}_{\mathrm{kj}}-{\mathrm{\&delta;}}_{\mathrm{kR}+}|\&le;{\mathrm{\&epsiv;}}_0\&DoubleRightArrow;P_{\mathrm{kj}}=P_{\mathrm{kR}+}\\ |{\mathrm{\&theta;}}_{\mathrm{kj}}-{\mathrm{\&delta;}}_{\mathrm{kC}-}|\&le;{\mathrm{\&epsiv;}}_0\&DoubleRightArrow;P_{\mathrm{kj}}=P_{\mathrm{kC}-}\\ |{\mathrm{\&theta;}}_{\mathrm{kj}}-{\mathrm{\&delta;}}_{\mathrm{kR}-}|\&le;{\mathrm{\&epsiv;}}_0\&DoubleRightArrow;P_{\mathrm{kj}}=P_{\mathrm{kR}-}\\ |{\mathrm{\&theta;}}_{\mathrm{kj}}-{\mathrm{\&delta;}}_{\mathrm{kC}+}|\&le;{\mathrm{\&epsiv;}}_0\&DoubleRightArrow;P_{\mathrm{kj}}=P_{\mathrm{kC}+}\end{array}\right.$ Formula 5

θ _kjrepresent vector with the angle of level to right, P _kR+, P _kR+, P _kR+, P _kR+represent respectively P _kjust be expert at, the neighbours territory point of row just, in negative, the row negative direction of row.

(9.3) calculate current some P of k layer recurrence _krelative mesh coordinate: establish current some P of k-1 layer recurrence _k-1with respect to Seed Points P ₀mesh coordinate be (r _k-1, c _k-1), according to P _kwith respect to P _k-1residing position, calculates P according to formula 6 _kwith respect to Seed Points P ₀mesh coordinate (r _k, c _k).

$\left\{\begin{array}{c}{P}_k=P_{k-1R+}\&DoubleRightArrow;[r_k,c_k]=[r_{k-1},c_{k-1}]+\left[\mathrm{0,1}\right]\\ P_k=P_{k-1C+}\&DoubleRightArrow;[r_k,c_k]=[r_{k-1},c_{k-1}]+\left[\mathrm{1,0}\right]\\ P_k=P_{k-1R-}\&DoubleRightArrow;[r_k,c_k]=[r_{k-1},c_{k-1}]+[0,-1]\\ P_k=P_{k-1C-}\&DoubleRightArrow;[r_k,c_k]=[r_{k-1},c_{k-1}]+[-\mathrm{1,0}]\end{array}\right.$ Formula 6

The 10th step: the row-coordinate minimum value r that finds the relative mesh coordinate of all unique points in unique point set P _min=min (r _i) and row coordinate minimum value c _min=min (c _i), according to formula 7, calculate each unique point P _iworld coordinates (X _i, Y _i), wherein, L is the world coordinates value of unique point spacing.By each unique point P _ipixel coordinate (x _i, y _i) and world coordinates (X _i, Y _i) corresponding one by one, complete the correspondence of unique point coordinate.

$\left\{\begin{array}{c}{X}_i=L\&times;(c_i-c_{\min })\\ Y_i=L\&times;(r_i-r_{\min })\end{array}\right.$ Formula 7

Example:

The core concept of the ACDT method that example of the present invention proposes be by neighbours territory recursive search traversal method analogy conventional in image connectivity domain analysis be used for traversal and there is equally the unique point set of array distribution feature, and in ergodic process, determine the position of all unique points on feature point mesh, thereby complete the correspondence of unique point coordinate.

The line segment aggregate that the non-longest edge of effective subdivision triangle in feature point set Delaunay triangulation result forms equates the line segment aggregate forming with feature point mesh neighbours territory point, this character provides a kind of method for traveling through unique point: by traveling through the non-longest edge of effective subdivision triangle, set up neighbours' domain information table of feature point mesh, it is the self-defined structure body array that an element number equals unique point number, and the C language definition of structure is as follows:

In unique point set ergodic process, can utilize the anti-distortion in part of feature point mesh, inquiry neighbours domain information table is upgraded partial row corresponding to current unique point, column direction, and determine that based on this neighbours territory point of current unique point, with respect to its position, obtains unique point with respect to the mesh coordinate of Seed Points.

Below in conjunction with accompanying drawing, the specific embodiment of the present invention is further described, the ACDT method that application the present invention proposes is carried out the corresponding of unique point pixel coordinate and world coordinates to the uncalibrated image shown in Fig. 2 (a) (unique point spacing 0.8mm), and concrete steps are as follows:

The 1st step: use OSTU adaptive threshold method to carry out binary segmentation to Fig. 2 (a), carry out the screening of validity feature point in cutting apart the bianry image obtaining.Observe the diameter of discovery feature point imaging 85 pixel left and right, set the high and low threshold value of image-region area and be respectively 7000,5000, the relative circularity threshold value of image-region is 0.25, carry out the screening of validity feature dot image region, result is as shown in Fig. 2 (b), and the process flow diagram of screening process as shown in Figure 3.

The 2nd step: according to the point coordinate in validity feature dot image region, this image-region characteristic of correspondence point pixel coordinate of least square ellipse the Fitting Calculation, detailed process is as follows:

Suppose validity feature dot image regional ensemble N _ein certain image-region n _eithe point coordinate of (the corresponding unique point of central point of each image-region) is (x _il, y _il) (i representation feature point sequence number, l represents point sequence number), carry it in the general expression of elliptic equation, can form with elliptic equation general expression parameter A _i～F _ifor the system of equations of unknown number, the quadratic sum Z of this system of equations residual error _ibe shown below;

$Z_i=\underset{l}{\mathrm{\&Sigma;}}{({x_{\mathrm{il}}}^2{A}_i+{2x}_{\mathrm{il}}{y}_{\mathrm{il}}{B}_i+{y_{\mathrm{il}}}^2{C}_i+x_{\mathrm{il}}{D}_i+y_{\mathrm{il}}{E}_i+F_i)}^2$

Calculating makes Z _iwhile getting minimum duration, the least square solution a of system of equations parameter _i～f _i, bring following formula calculating into and can obtain n _eicorresponding regional center pixel coordinate (x _i, y _i), this coordinate is validity feature dot image region n _eicharacteristic of correspondence point coordinate;

$x_i=\frac{b_i{e}_i-c_i{d}_i}{a_i{c}_i-{b_i}^2},$ $y_i=\frac{b_i{d}_i-a_i{e}_i}{a_i{c}_i-{b_i}^2}$

To N _ein each region carry out aforesaid operations, can obtain the pixel coordinate of all unique points, the results are shown in Table 2.

The 3rd step: unique point set is carried out to Delaunay triangulation, and result is as shown in Fig. 4 (a).

The 4th step: setting form variations angle threshold value is 5 ° of effective subdivision triangles of screening, and result is as shown in Fig. 4 (b).

The 5th step: remove the leg-of-mutton longest edge of effective subdivision, result is as shown in Fig. 4 (c), travel through the leg-of-mutton non-longest edge of effective subdivision, the sequence number of two end points and deflection are joined respectively in neighbours' domain information structure of the other side, neighbours' domain information table of setting up feature point mesh, the results are shown in Table 1.The algorithm flow chart of this step as shown in Figure 5.

Table 1 feature point set neighbours domain information table

The 6th step: choose No. 10 points as the set of Seed Points neighbours territory recursive search traversal unique point, ergodic process as shown in Figure 6.

The 7th step: in the ergodic process of unique point set, determine that according to searching route each unique point, with respect to the mesh coordinate of Seed Points, the results are shown in Table 2.

The 8th step: according to unique point spacing, unique point is changed with respect to the mesh coordinate of Seed Points, obtained the world coordinates of unique point, complete unique point coordinate corresponding, the results are shown in Table 2

The intermediate result that table 2 unique point is corresponding and net result

The above is preferred embodiment of the present invention, but the present invention should not be confined to the disclosed content of this embodiment and accompanying drawing.So every, do not depart from the equivalence completing under spirit disclosed in this invention or revise, all falling into the scope of protection of the invention.

Claims (3)

1. the automatic corresponding method of unique point coordinate based on Delaunay triangulation, is characterized in that, the method comprises the steps:

The 1st step: read uncalibrated image I, then it is carried out to binary segmentation, obtain bianry image I ';

The 2nd step: bianry image I ' is carried out to connected domain analysis, obtain set of image regions N, calculate the area value S of each image-region;

The 3rd step: extract the point coordinate of each image-region in set of image regions N, and calculate accordingly the circularities △ R of each image-region, and radius value R; Bring △ R and R into relative circularities △ r that formula 2 calculates each image-region in N;

$\mathrm{\&Delta;r}=\frac{\mathrm{\&Delta;R}}{R}$ Formula 2

The 4th step: screening meets that " S is at [S in set of image regions N ₁, S ₂] in scope and △ r < △ r ₀" image-region, the image-region that obtains is formed to validity feature dot image regional ensemble N _e, wherein, S ₁, S ₂be respectively the high and low threshold value of predefined image-region area, △ r ₀for the relative circularity threshold value of predefined image-region;

The 5th step: to validity feature dot image regional ensemble N _ein each validity feature dot image region n _ei, i representation feature point sequence number, the point coordinate extracting according to the 3rd step obtains an elliptic equation by the Fitting Calculation, and the elliptical center point being obtained by this equation calculation of parameter is image-region n _eicharacteristic of correspondence point P _i, its coordinate is (x _i, y _i);

The 6th step: establish by all unique point P _ithe point set forming is combined into P, and P is carried out to Delaunay triangulation, obtains subdivision triangle set τ;

The 7th step: for each subdivision triangle in subdivision triangle set τ, according to formula

η＝min(|ω ₁-90°|,|ω ₂-90°|,|ω ₃-90°|)

Calculate the leg-of-mutton form variations of subdivision angle η, wherein, ω ₁, ω ₂, ω ₃for the leg-of-mutton angle value of subdivision; Predefined form variations angle threshold value η ₀, in subdivision triangle set τ, screening meets " η≤η ₀" triangle, obtain effective subdivision triangle set τ _e;

The 8th step: travel through effective subdivision triangle set τ _ein each leg-of-mutton non-longest edge, the sequence number of two end points and deflection are joined respectively in neighbours' domain information structure of the other side, while finishing to traversal, can set up neighbours' domain information table of unique point set;

The 9th step: neighbours' domain information table of query characteristics point set, complete unique point P is put in an optional neighbours territory ₀as Seed Points, the set of neighbours territory recursive search traversal unique point; In ergodic process, determine current some P of every one deck recurrence _kwith respect to Seed Points P ₀feature point mesh coordinate (r _k, c _k), k represents the recurrence number of plies, obtains each unique point P in unique point set when traversal finishes _irelative mesh coordinate (r _i, c _i);

The 10th step: the row-coordinate minimum value r that finds the relative mesh coordinate of all unique points in unique point set P _min=min (r _i) and row coordinate minimum value c _min=min (c _i), calculate each unique point P _iworld coordinates (X _i, Y _i),

$\left\{\begin{array}{c}{X}_i=L\&times;(c_i-c_{\min })\\ Y_i=L\&times;(r_i-r_{\min })\end{array}\right.$

Wherein, L is the world coordinates value of unique point spacing; By each unique point P _ipixel coordinate (x _i, y _i) and world coordinates (X _i, Y _i) corresponding one by one, complete the correspondence of unique point coordinate.

2. the automatic corresponding method of unique point coordinate described in a claim 1, it is characterized in that, in the 8th step, index sequence number in described neighbours' domain information construction array represents the unique point sequence number of representative, wherein four long type variablees record neighbours' territory point sequence number of point corresponding to current index sequence number, and other four double type variablees record each neighbours territory and put corresponding deflection.

3. the automatic corresponding method of unique point coordinate described in claim 1, is characterized in that, the 9th step, and the detailed process of described unique point set neighbours territory recursive search traversal method is:

(9.1) calculate current some P of k layer recurrence _kpartial row's column direction angle: when k=0, the partial row of regulation Seed Points just, row just, negative, the row negative direction angle δ of row _0R+, δ _0C+, δ _0R-, δ _0C-be respectively 0 °, 90 °, 180 °, 270 °; When k>0, establish from current some P of k-1 layer recurrence _k-1to P _kvector with level to the angle of right, be γ, according to P _krelative P _k-1residing position, calculates P according to formula I _kpartial row's column direction angle δ _kR+, δ _kC+, δ _kR-, δ _kC-, and they are normalized in the scope of 0～360 °;

formula I

(9.2) determine current some P of k layer recurrence _ksome position, neighbours territory: according to formula, II determines P _kneighbours' territory point P _kjwith respect to P _kposition, wherein, j represents neighbours' territory point sequence number, ε ₀for predefined deviation of directivity threshold value;

$\left\{\begin{array}{c}|{\mathrm{\&theta;}}_{\mathrm{kj}}-{\mathrm{\&delta;}}_{\mathrm{kR}+}|\&le;{\mathrm{\&epsiv;}}_0\&DoubleRightArrow;P_{\mathrm{kj}}=P_{\mathrm{kR}+}\\ |{\mathrm{\&theta;}}_{\mathrm{kj}}-{\mathrm{\&delta;}}_{\mathrm{kC}-}|\&le;{\mathrm{\&epsiv;}}_0\&DoubleRightArrow;P_{\mathrm{kj}}=P_{\mathrm{kC}-}\\ |{\mathrm{\&theta;}}_{\mathrm{kj}}-{\mathrm{\&delta;}}_{\mathrm{kR}-}|\&le;{\mathrm{\&epsiv;}}_0\&DoubleRightArrow;P_{\mathrm{kj}}=P_{\mathrm{kR}-}\\ |{\mathrm{\&theta;}}_{\mathrm{kj}}-{\mathrm{\&delta;}}_{\mathrm{kC}+}|\&le;{\mathrm{\&epsiv;}}_0\&DoubleRightArrow;P_{\mathrm{kj}}=P_{\mathrm{kC}+}\end{array}\right.$ Formula II

θ _kjrepresent vector with the angle of level to right, P _kR+, P _kC+, P _kR-, P _kC-represent respectively P _kjust be expert at, the neighbours territory point of row just, in negative, the row negative direction of row;

(9.3) calculate current some P of k layer recurrence _krelative mesh coordinate: establish current some P of k-1 layer recurrence _k-1with respect to Seed Points P ₀mesh coordinate be (r _k-1, c _k-1), according to P _kwith respect to P _k-1residing position, III calculates P according to the following formula _kwith respect to Seed Points P ₀mesh coordinate (r _k, c _k):

$\left\{\begin{array}{c}{P}_k=P_{k-1R+}\&DoubleRightArrow;[r_k,c_k]=[r_{k-1},c_{k-1}]+\left[\mathrm{0,1}\right]\\ P_k=P_{k-1C+}\&DoubleRightArrow;[r_k,c_k]=[r_{k-1},c_{k-1}]+\left[\mathrm{1,0}\right]\\ P_k=P_{k-1R-}\&DoubleRightArrow;[r_k,c_k]=[r_{k-1},c_{k-1}]+[0,-1]\\ P_k=P_{k-1C-}\&DoubleRightArrow;[r_k,c_k]=[r_{k-1},c_{k-1}]+[-\mathrm{1,0}]\end{array}\right.$ Formula III.