Subpixel determination of imperfect circles characteristics

Pattern Recognition 41 (2007) 250–271www.elsevier.com/locate/pr

Subpixel determination of imperfect circles characteristics

Fabrice Mairessea,∗, Tadeusz Sliwaa, Stéphane Binczakb, Yvon Voisina

aUniversité de Bourgogne, Le2i – CNRS UMR 5158, Route des plaines de l’Yonne, BP 16, 89010 Auxerre, FrancebUniversité de Bourgogne, Le2i – CNRS UMR 5158, Aile de l’ingénieur, BP 47870, Dijon, France

Received 5 April 2006; received in revised form 4 January 2007; accepted 10 March 2007

Abstract

This article deals with the problem of the determination of characteristics of imperfect circular objects in discrete images, namely the radiusand center coordinates. To limit distortion, a multi-level method based on active contours was developed. Its originality is to furnish a set ofgeometric envelopes in one pass, with a correspondence between grayscale and a regularity scale. The adequacy of this approach was testedwith several methods, among them is the Radon-based method. More particularly, this study indicates the relevance of the use of active contourscombined with a Radon transform-based method which was improved using a fitting considering the discrete implementation of the Radontransform.� 2007 Pattern Recognition Society. Published by Elsevier Ltd. All rights reserved.

Keywords: Discrete circles; Geometric noise; Radon transform; Model fitting; Active contours

1. Simulated defaults and discretization problem

Picture data given by an image sensor is a discrete dataset taken from continuous parts of natural objects. Thus, theproblem of discretization is intrinsic in any physical image.The study of simple geometric forms such as lines or circlesbecomes difficult, especially if they are noised or deformed.Many studies are concerned with the several variants of discretelines. The point of interest of this article is the problem ofsubpixel measurement of imperfect discrete circles.

A known continuous circle can resemble several different ob-jects once discretized. The easiest example to understand is thepossible choice between the two discretized circles surround-ing an Euclidian circle, i.e. the “inner” or the “outer” (Fig. 1)shape for circles. Another simple example consists of takingpixels containing a part of the circle.

A unique precise circle has several discretized shape repre-sentations. It is not surprising to find some common definitionswith which to create discretized circles [1–3]. A natural

∗ Corresponding author. Tel.: +33 3 86 49 28 53; fax: +33 3 86 49 28 50.E-mail addresses: [email protected] (F. Mairesse),

[email protected] (T. Sliwa), [email protected](S. Binczak), [email protected] (Y. Voisin).

0031-3203/$30.00 � 2007 Pattern Recognition Society. Published by Elsevier Ltd. All rights reserved.doi:10.1016/j.patcog.2007.03.010

question is what definition is closest to the discrete circle pro-duced by the image acquisition chain? It is obvious that it de-pends on many parameters such as the scene, optics and illu-mination.

General 3D modeling by a direct or a reverse approach usingvirtual environments is a very difficult problem and can remainill defined. For this article observations will be restrained toin ideal conditions. The circles will consist of holes presentin homogenous or low textured plane parts placed in the focalplane under gaussian conditions in a dark box with illuminationequipment producing the simplest possible illumination model:the intensity of the light corresponding to a pixel is proportionalto the length of the contour observed by the pixel.

So, sample images will be obtained by discretization of acontinuous circle by the principle of intersection, i.e. binarycircles defined by the presence or absence of a circle part in thepixel’s field of view (FOV). It is agreed that discretization isdone as on a CCD sensor with square pixels without separations.

However, our goal is not to study the way to discretize a circlebut instead the way to return from a discretized circle to theEuclidian circle characteristics, mainly the center coordinatesand the radius, and the reachable precision of their subpixelmeasure [4], notably in the case of distortions and noise onthe contour. Indeed, it is not uncommon in image processing

http://www.elsevier.com/locate/pr

mailto:[email protected]




F. Mairesse et al. / Pattern Recognition 41 (2007) 250–271 251

to be confronted with a problem of distorted or noisy circularshapes. Noise is an intrinsic problem in every signal acquisitionchain, especially when the signal corresponds to images wherelight noise can vary greatly. For the distortion problem, in thisarticle, the point of interest is only on interior deformation oftenproduced by some mechanical production processes like crossdrilling which lightly touches the opposite surface of entranceof a loose material. The result of the drilling operation leavessome materials that occult a part of the drilling.

It is important to note that the study is not focused in contourextraction but in characteristic measurements. In the same way,the study considers distorted but complete circles. Indeed, evenif least-squares, circular Hough transform or other methodslike curve fitting presented [5–7] can be applied on incompletecircles, but it is not the case of the presented Radon-basedalgorithms. We also do not consider the presence of noise inthe background as it requires a pre-processing step that is notin the scope of this article. The studied case of noised circlescan be considered as the result of a noise reduction step.

The authors considered a test set produced by modeling acertain variety of perfect and imperfect discrete circles. Thecharacterization of an imperfect aspect of a circle can, for ex-ample, be given by the (4� area)/perimeter circular coefficientthat must be close to 1. The test set is composed of circles ofvarious sizes and positions, and can be divided into four groups:

Perfect discrete circles described by the following equations(Fig. 2):

x = E[x0 + r cos(�)],y = E[y0 + r sin(�)] (1)

E[x] being the nearest integer of x.Noised discrete circles which are derived from Eq. (1) in

order to produce imprecision at a hole’s edge [8,9] and are

Fig. 1. From left to right: Euclidian circle on a discretized frame, “inner”discretized circle and “outer” discretized circle.

Fig. 2. From left to right: perfect discretized circle, noised perfect discretized circle, distorted discretized circle and noised and distorted discretized circle.

described by the following equations (Fig. 2):

x = E[x0 + (r − Anoise) cos(�)],y = E[y0 + (r − Anoise) sin(�)]. (2)

Anoise is a random variable of centered Gaussian density prob-ability of 1/�

√2� e−(1/2)x2/�2

,As can be seen, the authors opted for a simple noise following

a normal centered law.Distorted discrete circles which simulate interior distortion

are described by the following equations (Fig. 2):

x = E[x0 + (r − adistortion) cos(�)],y = E[y0 + (r − adistortion) sin(�)] (3)

adistortion = a| sin(�) sin(��) e−�(�+�)4 |, where a, � and � rep-resent parameters of this deformation.

The authors choose a non-symmetric geometric model withinterior convex parts, but other models can be used.

Noised and distorted circles combine noise (Eq. (2)) anddistortion (Eq. (3)) (see Fig. 2).

The relative scale unit is the pixel. Subpixel sizes and shiftswill have an influence on the measurement error. This can beexplained by the modification of the value of the nearest integerwith a non-integer shift (Fig. 3). While varying the parametersdefining the different groups, we obtain a statistically significantset of simulated circles in order to compare the efficiency ofour measurement tools.

Fig. 3. Distortion of a perfect discretized circle due to subpixel horizontalshifts of 0.125 pixels.

252 F. Mairesse et al. / Pattern Recognition 41 (2007) 250–271

2. Circle parameters determination methods

This part essentially deals with the parameter estimationtools. Several approaches were studied in order to obtain thebest precision. Among them barycentric, circular Hough trans-form and least-mean squares based tools can be considered asthe classical methods. Another more innovative approach, us-ing Radon transform, furnishes globally good results.

The Radon transform is well known in image processingfor finding lines in a picture. The main idea is to analyze acircle through its tangents. Indeed, a circle can be completelydescribed by two tangents.

Consider a circle whose two describing tangents are paral-lel. The center of the circle is located at the intersection of theperpendicular line going through the contact point of each tan-gent of the circle and a fourth line parallel to the tangents andlocated mid-distance (“radius”) (Fig. 4). However, it is obviousthat this situation is only possible in the case of perfect cir-cles where this intersection really corresponds with the centerof the circle. That is why the study of other cases requires sev-eral pairs of tangents in order to find the most probable center.So, to a certain extent, while the characteristics of a lightly de-formed circle can be determined as we will see in Section 2.3,this kind of technique is not enough for heavy deformations. Inorder to avoid this limitation, we explored reducing the impactof deformations on circles. Our goal is to make the shape morecircular without modifying the center position. An approach ofthe shape at a coarse scale leads to a more regular contour. Thenatural question that comes into one’s mind is how to choosethe best scale for our purpose. In order to answer to this ques-tion, the use of a multi-level method is an adapted solution thatleads us to combine it with an active contours method in orderto obtain a multi-level contour image in a single pass as we willsee in Section 3. The application of algorithms on this imageor level-fixed extracted image will give results that justify theexplored method.

Moreover, the implementation of the Radon transform fixesthe center of an image as its reference point and needs to dis-cretize the angle parameter. This is a problem in a measurementprocess because the error on measures will not be the same ifthe shape is near to or distant from the reference point. In or-der to avoid this problem, a preliminary centering can be done.

Fig. 4. Tangents surrounding a circle.

For each shape, the barycenter is computed allowing us to centerthe form in a new smaller picture. Thus, we can reduce problemswith contour positioning and resources.

In order to validate our way of thinking, calculations withand without centering have been done. It is important to notethat the reference method will be the barycenter. As a reminder,we recall that this study takes into consideration holes withinner defects.

2.1. Least-mean squares methods for circles

In the circle estimation domain, least-mean squares basedalgorithms are numerous and form one of the two main cate-gories of circle estimators (with Hough being the other one).Our least-square methods are well-known methods [10]. Thesecond algorithm implies normalization [10,11] due to the in-troduction of a fourth parameter.

2.2. Circular Hough transform

The 3D circular Hough transform is also a very commontool for circle estimation. As its name denotes, the 3D circularHough transform needs an accumulator array of three dimen-sions; two for coordinates and one for radius space domain.Our approach leads us to consider an n−2D array where n isthe number of radii explored. This consumes less memory re-sources.

Let us focus on the Hough transform for circle forms. Cur-rently, there are several algorithms of this type, grouped ac-cording to their simplicity of programming or their resourcesrequirements:

The standard version uses tangents and normals. More pre-cisely, for each pixel of the contour, the tangent at that point isfirst computed in order to determine the normal. Then, the cur-vature is also computed for the side of the center in relation tothe contour. In an accumulator array, the point of the normal ata distance of one radius from the contour and where the curva-ture is positive is increased by 1. It is important to note that it isnot the local curvature which is considered, but a more globalcurvature that avoid discretization problems. In an ideal case,we obtain a maximum value of P where P is the number ofpixel contours [12]. This algorithm is not easy to implement,and furthermore, it requires a large amount of resources.

A simplification of the previous algorithm can be done. In-stead of computing the curvature to determine the side of thecenter, a double incrementing is made on each side of the con-tour. It leads to, in an ideal case, an accumulator containing amaximum point equal to P surrounded by a circle whose ra-dius is equal to twice the radius of the real circle. Admittedly,local maxima appear but they do not interfere with finding theglobal maximum, i.e. the center of the circle, even in cases ofdistorted or noised circular forms [13].

A difficulty remains concerning the determination of tangentson a discrete grid. For a given pixel, there can be several validtangents which differ only in a small shift of orientation. Con-sidering all valid tangents of a pixel, we can describe an arc


of possible circle centers. The maximum intersection of arcs isthen considered as the circle center. The orientation of tangentsdoes not matter as long as the distance of the projection doesnot change. Indeed, in a perfect case, the intersection of arcsoccurs at a single point. Since the distance of the projection ismore important, it is thus possible to consider a complete circlein order to find this intersection. In this way, the computationof tangents is no longer necessary, thus reducing computations.This is the algorithm that the authors have chosen.

All the above described methods still need the exact valueof the radius in order to give a correct answer. However, thisvalue is the unknown parameter that we are looking for. Twosolutions are then possible. One can explore all possible radiiin a large domain and then take the radius that gives the bestmaximum value. Otherwise, one can use another method toobtain a first approximation of the radius and improve thisvalue by finding the best Hough 3D maximum accumulator.There was no indication about which one would lead to thebest results, therefore, both were investigated.

In the third method, we study a discrete modeling of an idealform (i.e. circle). So, the principle that the accumulation is per-formed exactly at the center of the model is incorrect. Indeed,the discretization of circles can lead to an accumulation overa neighborhood around the real center position. The real max-imum accumulation is then not easily found as the accumula-tor is discrete. So, it is easy to find the accumulator maximumbut not necessarily the exact accumulation position. It is thesame for radii which are discrete. Our approach is to followthe evolution of the maximum discrete accumulation at eachradius and interpolate the real radius value. By nonlinear least-squares method, this step is done in order to furnish a subpixelapproximation of the circle radius (Fig. 5). Then, a 3D circular

Fig. 5. Fitting for radius estimation.

Fig. 6. From left to right: initial circle, 3D circular Hough accumulator (normalized by the accumulator maximum) at, respectively, r = 5, 11.5, 13.5, and 14.5.

Hough transformation is computed with this non-integer valueof the radius. Eventually, an integer value of the center’s coor-dinates is found.

The problem of the previous method is that it is very time-consuming because of the computation of the 3D circularHough transform over a large domain. Another way to proceedis to approach the radius value by increasing or decreasing thevalue of the radius parameter of the Hough transform in orderto find the adequate radius. This is possible thanks to the dis-cretized framework of an image. Nevertheless, geometricallyspeaking, there are not more than two intersections betweencircles except if they are identical. The radius (and the center)will be found when projected circles have a common inter-section point that is the center of the initial circle. It is thischaracteristic that will be used to optimize the 3D circularHough transform. As ideally there is only one common inter-section point, it is obvious that on a discretized framework thesame global effect can be seen with an advantage: globally,each time the estimation of the radius gets better, the maximumvalue of the accumulator increases. Following the maximumaccumulation value of projected circles will lead to the highestvalue that is associated with the best radius estimation. So, thenearer the intersections are to the center, the bigger the max-imum is in the accumulator (Fig. 6). It remains to follow theevolution of the maximum value over some radius in order toknow if the real radius value is inferior or superior. It is thenenough to follow the same reasoning until finding the greatestvalue in the accumulator. This described method can reduceresource usage.

2.3. Envelope by Radon and diameters accumulation

The Radon transform [14,15] is a transformation that allowsone to compute the projections of a data matrix, like an imagefor example, according to a certain angle �. To that aim, itcomputes the integral along the line orthogonal to � and returnsa 1D signal. It is a common tool in tomography as the inversetransform reconstructs the object from projections. The Radontransform is given by

R�(x′)=

∫ +∞

−∞f (x′ cos �−y′ sin �, x′ sin �+y′ cos �) dy′ (4)

where

[x′, y′] =[

cos � sin �− sin � cos �

] [x

y

]

and f an image.


This transformation is also well known for finding lines ina picture. Actually, an important value of the integral along aprojection is significant of an important amount of non-zerovalue pixels aligned which can be considered as a segment orportions of a line. The Radon transform, also known as the1D Hough transform, is thus a commonly known tool for lineidentification, so it is not, a priori, a method adapted to circularforms [16,17].

It can however be used for circle detection [18,19], to findchords of a circle, and to determine center position. As previ-ously said, we are interested in the search for whole tangentsof a circle and, because of the discretization, the tangents canbe confused with a small arc of a circle. The description of acircle by its tangents has already been done [20,21]. The Radontransform will give us a set of possible lines which are presentin a picture for each discrete angle (Fig. 7).

More precisely, the Radon transform converts an image(x, y) into a new domain (�, �). In this new domain, it is pos-sible to pair parallel tangents in order to, after further steps,obtain the circle center. This reasoning can be applied in thecase of a perfect circle since a discrete line passing through themaximum of non-zero value pixels will be a discrete tangent.This is obviously not the case for any discrete contour as onecan see in Fig. 8. A distortion can consequently lead to miss-ing the oriented lines we wish to use, i.e. the tangents. Suchan error clearly produces inaccuracies but our results indicatethat these are not as large as the errors given by the least-meansquares methods.

2.3.1. Initial algorithmAs previously explained, the first step of the algorithm is to

describe a circle by its tangents. To this end, the authors use

Fig. 7. Circle envelope reconstruction by numerous tangents.

Fig. 8. Oriented lines passing by the maximum of active pixels.

the Radon transform. This way we obtain an image where thex-axis corresponds to the orientation of the integral line, they-axis to the � parameter and the intensity level to the lineintegral value. The larger this value is, the more pixels the linecontains. Once this step is done, it is necessary to identify thelines that are tangents.

Considering the discretization phenomenon, the contact be-tween a circle and a tangent need not necessarily reduce to onepixel and is most of the time composed of several pixels. Thatis why the authors consider the line which contains the mostpixels as the most probable tangent. It is thus a question offinding maxima in the computed Radon image.

In order to select each tangent pair of the same orientation, itis necessary to find the two maxima which correspond to thesetangents. However, it is not possible to directly consider max-ima as several maxima could be found that do not necessarilycorrespond to these tangents. To that end, the authors split theRadon signal into two parts in the middle according to the �parameter. Indeed, the integral line is equal to zero where itdoes not intersect the circle and non-null when it does. So fora given orientation, the Radon transform has a height of twicethe radius of the circle. It is then possible to isolate the contri-bution of each half circle for a given orientation. This splittingis done by a barycentric calculation. Now, the Radon image issplit into two parts as can be seen in Fig. 9. It is now possibleto identify, for each orientation, maxima in these two parts inorder to obtain all the circle’s tangents (Fig. 11). The pairing oftangents of the same orientation is simple to do as the Radonrepresentation is given in (�, �) (Figs. 10, 11).

It is obvious that, for noised or distorted circles, the tangentsfound can perhaps not correspond to tangents of the ideal circle.This is due to a perturbation (induced by the noise or distor-tion) in the Radon signal that is translated by a global maximapositions displacement (Fig. 12). However, from now on eachpair of lines found is considered as tangents of the evaluatedcontour.

Once each pair of tangents is computed, the next calcula-tion step is to determine the line having the same orienta-tion and passing mid-distance between each pair of tangents.It is, in fact, a question of lines that contains radius for anideal circle. The intersection of these lines leads to the circle


Fig. 9. Top: initial Radon signal. Bottom: upper and lower part of the initial signal separated by a barycentric method.

ρ(θ)>0

Tangent contact point

Tangent contact point

θ oriented line passing

through the center of the

image

ρ(θ)<0

Fig. 10. � positive and negative domains considering the circle center as the center of the image.

center position. This can be done in the Radon image by find-ing the mean of � parameters and value for each tangent pair(Figs. 11 and 12). This gives us a new Radon signal corre-sponding to lines passing through the center circle in the initialimage.

The final step is to compute the accumulator array. This canbe done by inversing the Radon transform [22] (Fig. 13). Themaximum of this array corresponds to the most probable centerposition since it corresponds to the position where there are themost number of line intersections (Fig. 19).


Fig. 11. Top: initial Radon signal of a perfect circle, Middle: upper and lower parts of the signal and the envelope obtained by maximum search. Bottom:complete envelope of the signal and average curve describing the circle center.


θ (degrees)

0 20 40 60 80 100 120 140 160

-100

-80

-60

-40

-20

0

20

40

60

80

100

5

10

15

20

25

Rθ (x′)

x′

θ (degrees)

Rθ (x′)

0 20 40 60 80 100 120 140 160

-100

-80

-60

-40

-20

0

20

40

60

80

102

4

6

8

10

12

14

16

18

20

22

θ (degrees)

Rθ (x′)

0 20 40 60 80 100 120 140 160

-100

-80

-60

-40

-20

0

20

40

60

80

1002

4

6

8

10

12

14

16

18

20

22

θ (degrees)

0 20 40 60 80 100 120 140 160

-100

-80

-60

-40

-20

0

20

40

60

80

100

5

10

15

20

25

θ (degrees)

Rθ (x′) Rθ (x′)

0 20 40 60 80 100 120 140 160

-100

-80

-60

-40

-20

0

20

40

60

80

100

5

10

15

20

25

θ (degrees)

Rθ (x′)

0 20 40 60 80 100 120 140 160

-100

-80

-60

-40

-20

0

20

40

60

80

100

5

10

15

20

25

θ (degrees)

Rθ (x′)

0 20 40 60 80 100 120 140 160

-100

-80

-60

-40

-20

0

20

40

60

80

1002

4

6

8

10

12

14

16

18

20

22

x′x′

x′ x′x′

x′

Fig. 12. Top: initial Radon signal of a noised circle, Middle: upper and lower parts of the signal and the envelope obtained by maximum search. Bottom:complete envelope of the signal and average curve describing the circle center.


Fig. 13. Center determination by accumulation of lines.

The radius value calculation is based on the variation of the� parameter between each tangent pair.

So, the center position is found at single pixel-precision.This is not quite as precise as other methods. Our goal is toobtain a subpixel estimation of each parameter with the pur-pose of a better comparison between our approach and theclassical tools of circle estimation. To that end, the approxi-mation on the circle center position obtained during the searchof the maximum accumulation can be avoided by determiningcircle parameters in the (�, �) domain as we can see in 2.3.3.Moreover, in the case of a noisy circle, there is a possible wayof improvement to be explored. A noisy circle shows its con-tour in the form of a cloud of points which are more or lessdistant from the real shape. The deletion of the farthest pointsallows us to obtain the closest of the points which well definethe circle. This might increase precision in circle parameterdetermination.

2.3.2. Algorithm including an unlikely point elimination (UPE)Tangent determination in the case of noisy or distorted

circles is not optimal as the hypothesis that lines passing bythe high number of active pixels are tangents of the idealcircle is not valid anymore. The dispersion of contour pointsat the ideal circle because of noise or distortion implies that“pseudo-tangent” positions can change compared to theirideal positions, i.e. for a given position, the orientation of thepseudo-tangent is not the same, causing problems in the tan-gents pairing. In order to reduce the induced errors, the dataset is pared of the points whose position seems to be far froman estimated circle position. The goal of this computing is tosuppress the points which are too far from the estimated circlecenter.

The authors opted for a model of �� after three successiveiterations, where � is the standard deviation. The parameter �is chosen in the triplet [1; 2; 3] that means 27 possibilities.

The aim of this algorithm is to discover if a set [�1 �2 �3] isparticularly relevant for the suppression of unlikely contourpoints. This might be useful for noised circles for which thecontour is not well defined. The algorithm is based on the factthat a point in the image domain will be described in the (�, �)

domain by the equation

� = x cos(�) − y sin(�). (5)

Obviously, it is also the case for the center of the circle. FromEq. (5), it is possible to determine the distance of computedpoints of the discretized Radon transform at the ideal curve,and then the standard deviation. The iterative process is thencarried out (Fig. 20).

2.3.3. Algorithm considering fitting aspect (FA)As we have seen, a convenient way to avoid the pixel ap-

proximation of the accumulation image of the basic algorithmis to directly manipulate data in the (�, �) domain, i.e. withoutreturning to the image domain. The aim is to obtain the bestapproximation of the center of a circle. The computed set ofpoints is expected to be the description of the center of thecircle in the Radon domain and, as a point in the image do-main, it is described by a curve following equation (1.5). Thesearch for the two optimum parameters set, which best fitsthe data, provides a pair of real numbers that correspond tothe coordinates of the circle center. More precisely, these noninteger numbers are a subpixel estimation of the circle center.The method used for the fitting is a trust region algorithm [23]for nonlinear least mean squares. Fig. 14 is an illustration ofthat fitting on both distorted and noisy circles and on a perfectone.

A complementary process is used to increase measurementprecision. It is similar to the UPE (2.3.2) as it is based onthe suppression of fitting unlikely points. The distance ofpoints from the ideal estimated circle is computed thanks to


0 20 40 60 80 100 120 140 160 180-40

-30

-20

-10

0

10

0 20 40 60 80 100 120 140 160 180

-40

-30

-20

-10

0

10

20

0 20 40 60 80 100 120 140 160 180-40

-30

-20

-10

0

10

0 20 40 60 80 100 120 140 160 180

-40

-30

-20

-10

0

10

20

Center position (-15,35)

Center estimation (-15.07,34.42)

Center position(-15,35)






Fig. 14. Left: fitting for perfect circle by FA (top) and FADE (bottom) methods. Right: fitting for distorted and noised circle by FA (top) and FADE (bottom)method.

Fig. 15. Top-left: initial noised circle, from top-middle to bottom-right: suppressed pixels according to the number of iterations.

the standard deviation reported at the calculated circle po-sition. The method used is an unlikely method of type ��.In order to study effects of this approach, no constraint ismade on the number of suppressed points except to keep atleast two points in order to complete the fitting but in prac-tice it is considered complete after nine iterations. This is

an iterative process which deletes several points in a cycleto slowly converge to a solution. An illustration is presentedin Fig. 15 which represents this point deletion process. Letus specify that the calculation is directly done in the Houghdomain by considering the distance of a point on the curve(Fig. 21).


2.3.4. Algorithm considering fitting aspect with discretizationeffect (FADE)

In order to obtain better results, the authors decided to takeinto consideration the discretization of the circle in the fittingmodel. Currently, Eq. (5) does not exactly correspond to thedata set even for a perfect circle. This is due to the discretizedaspect of the representation in the (�, �) domain. In order toconsider this aspect in the fitting, one needs to change Eq. (5).The empirical choice the authors carried out is to consider thefirst terms of the Fourier series expansion of Eq. (5), whilekeeping the previous fitting method and considering a triangleform of the signal thanks to the two additional parameters (c& d). The development is done on the fourth order. Tests withhigher orders give worse results. It seems that at higher ordersthe addition of terms with strong variation hinders the follow-ing of the discretized steps. This gives rise to the followingequation:

y = X + c + sin(�(2X + 1))

�+ sin(2�(2X + 1))

2�

+ sin(3�(2X + 1))

3�+ sin(4�(2X + 1))

4�, (6)

where X = a cos x − b sin x + d; (a, b) the center coordinatesof the estimated circle; and c and d are real numbers.

The same process of points’ suppression is used except thatit needs four points in order to complete the determination of(a, b, c, d).

In order to be convinced of the difference between FA andFADE methods, the reader is invited to compare results pro-vided in Fig. 15 where the discretized step effect is clearlyvisible on perfect circle, and noised and distorted circles.

3. Active contours multi-level approximations

As can be noticed, several approaches are possible for com-pensating imperfections in circles. Among them, one can quotefitting methods and suppression of problematic points whichhave been used successfully. Considering this last method, thedifficulty is to determine which pixel of the contour belongs tothe form or to the noise and distortion i.e. if the pixel must bepreserved or eliminated. However, in order to limit pixel elim-ination errors that can lead from a better to a worse estimationof the center position, a solution can be arrived at to have amore global approach of the contour instead of a contour pixelapproach. Indeed, the authors opted for a pre-processing on thewhole contour in order to adapt data to the philosophy of theRadon-based methods. This pre-processing provides a familyof new contours obtained by an active contours approach inthe aim of rounding the initial form especially for noised ordistorted circles.

3.1. Physical approach

In order to have a better approximation of circle parameters,a pre-processing of the initial image, before applying measure-ment tools can be a source of precision improvement. The aim

is to modify initial data in order to create a new data set thatwill be more efficiently processed than the initial one. So, themeasurement methods will be applied on that modified set ofdata, instead of raw data, in order to improve global precisionof measures.

The basic idea is to obtain a new family of envelopes fromthe initial circular shape. This family is composed of severalcontours which represent the circular form viewed at differentscales. These new contours follow more or less the originalform according to the scale (Fig. 16).

The concept of granularity of a contour expresses the no-tion that the authors wish to highlight. The granularity, in theauthors’ minds, translates the impression of roughness of acontour. The more irregular is the contour, the larger the gran-ularity. The study of the form at different scales amounts tostudying the contour at several granularities. This way, a setof several new contours of the initial form is furnished as theenvelope of a contour is not the same for all granularity levels.It narrows down to having a family of envelopes following rea-sonably well the contour of the object. The advantage of thatmethod, compared to a convex envelope, is that the contouris rounded without impinging too much on the position of thecenter of the object i.e. a distortion is not worsened. That way,the authors hope that the deformation will be partly compen-sated thanks to the rounding effect, improving the final results(Fig. 17).

Regarding the family of envelopes, a granularity parameteris empirically fixed to an intended value of 2 which means thatthe approximation will not mould the circle but, on the contrary,will form a global envelope of the form. This value is chosenexperimentally but can be related to the radius of the circle,for example. The multi-level contour is then calculated aroundthis granularity in order to give a precise grayscale (256 levels)picture that translates the evolution of the form following thegranularity. Each gray level corresponds to a scale, i.e. to agranularity. A simple thresholding allows getting the contourapproach of the object at a fixed scale, i.e. granularity. Finally,the edge of the threshold is determined and used as a newcontour (Fig. 18).

The contour is determined thanks to the progress of an ini-tial marker from the border of the image to the center of thepicture. This marker spreads itself out depending on a front.This propagation is allowed or not according to the granular-ity which varies according to the magnitude of the granularityparameter. If a propagation path is too reduced compared tothe granularity parameter, the propagation path is stopped forthis scale, although the propagation can still continue at smallerscales. So an advantage of this method, in comparison withclassic active contours, is that the calculation can be made inonly one pass, to have the complete family of envelopes arounda granularity/scale. Moreover, the algorithm is well adapted foran electronic implementation. This implies that the computa-tion can be delegated to specific hardware. For the moment, thecomputation of the active contours consumes a large amount ofresources which is why a hardware implementation is neces-sary for industrial applications. Also note that this method canbe used for multi-scale classification.


Fig. 16. Top: initial image (left) and final active contours image (right). Bottom: several thresholds of the final active contours image (i.e. several views atdifferent scale) with inner and outer markers.

Fig. 17. From left to right: initial image, final active contours image, and two thresholds of the final active contours image (i.e. several views at differentscale) with outer markers only.

Fig. 18. From left to right: initial image, edge of final active contours image, and edges of two thresholds of the final active contours image (i.e. several viewsat different scale) with outer markers only.

3.2. Theory

Let us consider a bi-dimensional regular discrete (N, M)

length grid � on which the following bistable diffusive systemis defined:

dvn,m

dt= Dn,m[vn−1,m + vn+1,m + vn,m−1 + vn,m+1 − 4vn,m]

− vn,m(a − vn,m)(1 − vn,m), (7)

where Dn,m is a local diffusion parameter and, a, a thresholdparameter.

The system is completed by the Neumann conditions (zero-flux conditions) on the border �� of the definition domain �,so that

�vn,m

�= 0 if (n, m) ∈ ��, (8)

where �/� denotes the outer derivative boundary.Note that it is a discrete version of the FitzHugh–Nagumo

partial differential equation (PDE). The presented active con-tours method is based on a propagation of topological wavesthat are based on this equation, so that an image evolves through


BINARY IMAGE

BINARY IMAGE IN RADON DOMAIN

ESTIMATION OF RADIUS POSITION

ESTIMATION OF CENTER POSITION BY

MAXIMUM ACCUMULATION

UPPER PART LOWER PART

UPPER TANGENTS LOWER TANGENTS

Radon transform

Barycenter

Maximum intensity

Mean of positions

Inverse Radon transform

Fig. 19. Scheme of the initial Radon algorithm.

the dynamic system

E(A, B) :

⎧⎪⎪⎪⎪⎨⎪⎪⎪⎪⎩

dvn,m

dt= D(B)

(vn−1,m + vn+1,m + vn,m−1 + vn,m+1 − 4.vn,m) − fa(vn,m)

with fa(v) = v(a − v)(1 − v) for a <1

2,

D(B) = 1 + tanh(20.B − 12)

2and v|t=0 = A

⎫⎪⎪⎪⎪⎬⎪⎪⎪⎪⎭

, (9)

where vn,m corresponds to the intensity of the image at pixel(n, m), and A is an initial marker from which propagationarises. f (v) is a cubic function which sets two attracting val-ues v = 0 and 1. Due to the fact that � < 1

2 , only the stablestate v=1 can propagate. The diffusive parameter controls boththe velocity of the traveling wave and the path of propagation,so that when D < D∗ >1, the path is forbidden, (with pa-rameter D∗, only function of �). When D > D∗, a minimumsize of allowed path D > Dcritical is required to let the wave

propagate, as D also defines the wavelength of the travelingwave. Applying relevant information from the image to be pro-cessed (i.e. defining B in the system), one can obtain controlledtraveling waves in the medium leading to a stable convergingstate where the wavefronts are stopped in the vicinity of forbid-den pixels or in-between them if they are close enough. Thisresult leads to obtaining active curves connecting different pix-els depending on the local diffusive parameter.


BINARY IMAGE



ESTIMATION OF CENTER POSITION



Radon transform

Barycenter

Maximum intensity

Mean of positions

Unlikely points elimination

Fig. 20. Scheme of the UPE algorithm.

4. Study results

4.1. Information on study conditions

We focus on cases with noise amplitude of 5 and distortionamplitude of 10 (for a radius ranging between 25 and 40 pixels).We consider that it furnishes sufficiently noised and distortedcircles for our purpose. So, in the rest of this article, we analyzethese cases and specify when noise or distortion varies. Withoutexplicit changes, the following conclusions (cf. 3.3) are madeon a model whose characteristics are a noise amplitude of 5and distortion amplitude of 10. The granularity parameter is setto 2 for the following results.

Algorithms were developed and tested with MatLab 7.0.4.on a Windows XP Pro platform (Pentium M 1.7 GHz, 640 MoRAM). Computation times are however not representative oftrue performance as the algorithms have not been already opti-mized, compiled or targeted. In its actual form, active contours

are very time-consuming because of the propagation of thewave front, but combined with the centering (smaller picture)it only takes few seconds. It is the same case for Radon-basedmethods. So, algorithms in their present state are not readyfor real-time processing, but further optimized implementationmay be applied in real production conditions.

4.2. Individual results analysis

4.2.1. Least-mean squares methodsAs previously stated, the least-mean squares method is a very

common and efficient tool for circle recognition. It was ob-vious that its results will be very precise for a perfect circle.One can see that mean distance error is lower than 0.015 pix-els with a standard deviation of the same order. In the sameway, the addition of a symmetric noise clearly modifies but notdramatically results as the mean error is lower than 0.15 andhas a standard deviation of half of that value. Difficulties in


BINARY IMAGE



ESTIMATION OF CENTER POSITION



Radon transform

Barycenter

Maximum intensity

Mean of positions

Fitting

Fig. 21. Scheme of the FA/FADE algorithm.

getting a well-positioned center appear with the processing ofdistorted circles as, for a distortion of 10, the mean error withleast-squares methods reaches more than 2.1 pixels (� ∼ 0.09).The noised and distorted case simply adds previous errors oncircle center position to increase a little the mean error to 2.25pixels (� ∼ 0.14).

The use of active contours associated with the least-squaresfor measurements does not seem to be profitable for perfectcircles as the contour is less precise. In the other cases, theimprovement can be interesting as the precision can only beincreased ∼ 0.15 pixels (counting on standard deviation vari-ation). It is clear that the use of active contours, for the least-squares methods, is helpful only for a very precise measure-ment. The improvement is visible but not at the level of whatit can be as can be seen in the rest of this article. The center-ing of the form does not bring anything here as the position ofthe contour in the picture does not intervene in the calculation

of least-squares approximations. The difference between thenormalized and non-normalized versions is inconsequential, atleast for a noise amplitude of 5.

4.2.2. Barycentric methodThe barycentric method, which is used to center the form in

the picture in order to improve precision with Radon-based al-gorithms, is classically less precise than least-squares approach.However, the combination with active contours can greatly im-prove results in perfect or noised circles. For a perfect dis-cretized contour, the precision order is about 0.36 pixel butreaches about 0.05 pixel because the active contours combinedwith a depletion of the standard deviation from 0.21 to 0.7. Ina noised context, results are not as impressive as the precisionof a simple barycenter varies a lot in the [0, 1] pixel domain.There, active contours do not always improve results but rarelyworsen precision to more than 0.1 pixels. The interesting fact is


that the precision is generally lower than 0.5 pixels (noise=5).It can be very helpful for industrial inspection to be sure thatthe precision can be delimited. Globally, the mean error fallsfrom 0.4 to 0.24 pixels with active contours (threshold of 27)

and the standard deviation from 0.2 to 0.12. Thus we can hopefor a little improvement with the use of active contours in anoised contour. The case of distortion is problematic as the pre-cision with active contours gives practically the same result asthat without, respectively, ∼ 1.8 pixels with a standard devia-tion of 0.11 and ∼ 1.9 pixels error with a standard deviation of0.28. One can see that the biggest improvement occurs in thestandard deviation and not in the precision. The best explana-tion for this constancy is that the deformation is the same onthe test set so even if positions are different, envelopes mustbe practically the same so the error position tends toward thesame value. Oddly, the mean error with noise and deformationis lower than the mean error with only distortion ( ∼ 1.7 vs 1.9)for a standard deviation of the same order. This can be due tothe noise but the sample set is too small to be sure of the expla-nation. The active contours do not bring anything in this caseas the best results (threshold of 27) are equal to those withoutactive contours.

4.2.3. 3D circular Hough transformThe 3D circular Hough transform is not a subpixel process.

So it is obvious that the precision offered is often less thana subpixel version. We applied the active contours before the3D circular Hough transform for study purposes. As we sus-pected, it is not helpful but quite the opposite. This behavioris predictable. The 3D circular Hough transform uses an ac-cumulator where a circle is summed for each point of initialpicture. The active contours give us several envelopes but theinterior part is full as a disk. So there is a loss of precision onthe accumulator maximum due to the addition of pixels. Theuse of an internal marker in addition to the external marker toform a ring envelope could be a better approach, but we donot think it will increase precision because of reasons previ-ously discussed, and that will stay true for a thick envelope.The precision obtained is much worse than those obtained withleast-squares or barycentric methods. The mean error is of 2.45for perfect case, 1.8 for noised contour, 4.38 for distorted cir-cles, and 4.14 for noised and distorted circles. One can see thatthe addition of noise increases the precision. This phenomenoncould be due to the compensation of the discretization of theform but this is not the purpose of this article. This method isclearly not precise enough for metrological measurements.

4.2.4. Initial Radon transform-based methodThe initial implementation of the Radon transform estimator

is not subpixel so the error position is the square root of theinteger. The precision for a perfect circle is then 0, 1, or

√2.

The precision obtained is quite good (0.61 pixel error, � ∼ 0.53)but it is a perfect case. The use of active contours does not addanything as they have little impact on the results (non-subpixelapproach effect). The addition of noise (amplitude 5) impliesa degradation of precision which reaches a boundary of

√5

even if mean error little varies (0.9 pixel error, � ∼ 0.65). Here,the active contours are helpful as the position error varies togive a mean error and a standard deviation of 0.52 pixel fora threshold of 26. For distorted circles, the precision is worsethan with least-square methods with an error of 2.89 pixels(� ∼ 2.12). Again, the active contours are very helpful as theposition error falls to one pixel error (1.08 pixel error, � ∼ 0.90,threshold of 26). This is a great improvement as the error isless than that of least-square methods in distorted case. Thecase of distorted and noised circles, as one could expect, is thecombination of the precision on noised and distorted contours.Here, the position error reaches nearly 3.5 pixels. That is a verybad result compared to those previously obtained. However,thanks to the active contours (threshold of 26), it is possible toimprove results up to those of least-squares methods, i.e. nearly2 pixels error. The problem of the precision of the measurementremains as the standard deviation is around 1.5 pixel.

4.2.5. UPE methodAll the following Radon-transform based algorithms give

subpixel results. The first subpixel method that we will studyis a combination of the Radon transform and an unlikely pointsuppression (UPE), with and without active contours. The ad-vantage is that the precision will not depend on a discrete accu-mulator for center position. For perfect circles, the best resultsgive a 0.45 pixels error (� ∼ 0.16). The use of active contoursincreases by little the precision, with a mean error of 0.31 pix-els and a standard deviation of 0.12. Globally, it seems to bejudicious to make a choice of � = [3, 3, 3] (among the tripletsthat give the same precision) combined with active contoursthreshold at 26 that gives the more reliable results. But theseresults are not at the level of least-squares based results. Thecase of the noised contour is problematic. It is not possible toclearly identify a good triplet �. Indeed, without active con-tours, a precision of 0.54 pixel for � ∼ 0.27 can be obtainedwith � = [2, 1, 3]. The use of active contours is also a prob-lem as � is not necessarily the same depending on the level ofapproximation. An � of [1, 1, 3] or [2, 1, 3] gives better resultsfor low thresholds. On the contrary, � = [3, 1, 3] gives betterresults for high threshold as well as the best precision for athreshold of 27 with a mean error of 0.34 pixels and a stan-dard deviation of 0.12. One can see that the best precision isobtained with a threshold of 25 and an � of [1, 1, 3] (� ∼ 0.14).Distorted circles take advantage of the new way of calculationas the best precision is at 1.73 pixels with a � of 0.55. This is agreat improvement since the first method gives a mean error of2.89 (� ∼ 2.12). Active contours are also useful for improvingresults. Indeed, the precision can fall to nearly 1.68 pixel error(� ∼ 0.16) with, for example, � = [3, 1, 3] and � = [1, 1, 3] ora 1.45 pixel error (� ∼ 0.4) with � = [1, 1, 2]. It is the samething for noised and distorted circles except that the problemof finding an optimum � is even more difficult. Indeed, severaltriplets give very close mean errors and standard deviations.This way, it is possible to choose each triplet whose last integeris 3 as they give very similar results at the end (one can see thecumulative error among the different level of approximations).This is due to the unlikely point suppression process. A dras-


tic suppression at the beginning of the process is not really aproblem as the value of � is big, that is to say that only pixelswhich are far from the computed center are deleted. But, � de-creases with the iteration process and, this way, the eliminationis done on closer and closer points. That is why it is necessaryto have an � that compensates the decreasing of � in order tonot suppress too many pixels.

4.2.6. FA methodThe second vision of subpixel Radon transform combines

the fitting with an unlikely point suppression (FA). The studyof a perfect circle shows an improvement of results as the pre-cision is around 0.45 pixel error (� ∼ 0.15). The process ofpoint suppression is, here, not really interesting as it inducesvery little improvement ( ∼ 0.002 pixel), which is logical as thecontour is “perfectly” discretized. Active contours improve re-sults and the best result is obtained for a threshold of 26 with0.31 pixel error (� ∼ 0.16). This time, the suppression of pointsclearly decreases the precision. But the least-squares methodgives more precisely the position of the center. The addition ofnoise (amplitude 5) on the contour globally leads to an approx-imately half pixel error on the centers position (� ∼ 0.27). It ismore than least-squares error but a great improvement can beseen compared to the previous version. The use of an envelopesfamily mostly improves the precision of center position in or-der to furnish a precision at 0.33–0.34 pixel for a standard de-viation of 0.12–0.13 (threshold of 25 and 27). Oddly, the pointsuppression process does not improve the precision although itis a question of a favorable case. Seeing these results, it seemsunnecessary for noised contours. An explanation would be thatthe noise is symmetric and can, this way, be compensated forin the computing. Adding distortion to perfect circles revealsthe propensity of this algorithm to make abstraction of defor-mation. Indeed, this version is already better than the previousversion, but combined with active contours, it becomes morereliable than least-squares algorithms. Here, the elimination ofpoints is completely justified as the precision obtained is lessthan 1 pixel error for a distortion of 10. Although this is withoutenvelope computing, the mean error is of 1.5 pixel for �=0.83(8th iteration), and a precision of 0.69 pixels (� ∼ 0.15) can beobtained with active contours (threshold of 27 and 7th iteration)or a precision of 0.64 pixels (� ∼ 0.21) for a threshold of 27

and the 8th iteration. In the case of noised and distorted circles,one can see that, most of the time, a precision of 1.17 pixels(� ∼ 0.63) is available thanks to active contours and point sup-pression process (threshold of 27 and 9th iteration). Comparedto least-squares that furnish results with an error of 1.8–2.2 pix-els (considering the standard deviation), it is a great enhance-ment. But the problem is to define the level of the thresholdingfor envelope determination. Point suppression processes canlower algorithm efficiency, but, globally, it clearly improves re-sults even if the standard deviation increases.

4.2.7. FADE methodThe last Radon method, compared to the FA version, takes

into account the discretized phenomenon of the computed

Radon transform (FADE). One can see that the obtained pre-cision is now better than that of the before. The mean error isof 0.14 pixels (� ∼ 0.2). As previously explained, the unlikelypoint elimination does not improve results as there are nounlikely points to suppress. Active contours do not bring anyimprovement on perfect circles as they lead to an approxima-tion that is necessarily more perturbed than the perfect form.So a lack of precision can be seen when active contours areapplied. Taking into account the discretization problem provesits interest even if the least-squares methods give better preci-sion results. However, as one can see, the precision obtained(1.27, � ∼ 0.24) for noised circles is worse than those withoutthe aspect of discretization. Even if the use of active contourscan partially compensate this problem (0.36, � ∼ 0.24), it ispreferable not to consider the discretization. It is the same fora distorted circle, but this is more sensitive on the standarddeviation whose value is more than the double compared to theprevious version (2.42, � ∼ 2.18 compared to 1.50, � ∼ 0.83before). Indeed, the point suppression process and the activecontours are not sufficiently helpful to compensate the loss ofprecision (1.99, � ∼ 1.6 compared to 0.69, � ∼ 0.15 before).The same problem is visible for noised and distorted circles.The precision is close to the previous version (3.32 pixel com-pared to 3.03) but the standard deviation is twice the previous( ∼ 1.09) if point elimination is not considered. With activecontours, the precision reaches 2.44 pixel error (� ∼ 1.36) atbest. This is far from the previously obtained precision i.e.1.17 pixels (� ∼ 0.63). As can be seen, to the contrary ofwhat it seems in 2.3.2.b, taking into account the discretizationphenomenon is only interesting for perfect circles, and least-squares methods are clearly better in this case. This algorithmis not the best for measurement of distorted and/or noisedcircles. One may prefer the third Radon implementation.

4.3. Global results analysis

Seeing these results, it is obvious that the active contoursare very helpful when they are combined with advanced Radontransform estimator. Surely, the improvement depends on thechoice of a level approximation, but it can be noticed that thethresholding in the middle of the domain gives the best resultfor distorted cases (noised or not). If the position is well definedwithout active contours, the “blurring effect” due to an envelopecauses a little imprecision. In return, if the center position isnot precise, the use of active contours improves precision. Inthe same way, the unlikely point elimination gives noticeablegains of precision for the third Radon-based method, even if itis to the detriment of the standard deviation.

The use of a pre-centering of the form thanks to the barycen-ter is without effect for least-squares methods and, obviously,for the barycenter method. Its aim is to improve the precisionof Radon-based calculus. Indeed, the reference point for theRadon transform is the center of the picture and the calculationof a centered form would give better results (cf. Fig. 8). Forthe first implementation of the Radon transform, results are notso much improved as they are quite the same, with or with-


Table 1Results of barycentric method

Mean error of the barycentric method given in pixels

Threshold Perfect circle Noised circle Distorted circle Noised and distorted circle

Original 0.35986 0.39971 1.9195 1.6698Active contours 2 0.049311 0.97638 1.8133 2.4061

4 0.049311 0.85625 1.8133 2.29128 0.049311 0.72115 1.8133 2.2213

16 0.049311 0.57251 1.8133 2.087432 0.049311 0.38491 1.8133 2.046364 0.050374 0.27618 1.8223 1.8944

128 0.049589 0.24229 1.8241 1.8636

Table 2Results of least-squares method

Mean error of the least-squares method given in pixels



4 0.022057 0.12014 2.1448 2.02868 0.022057 0.11924 2.1448 2.0472

16 0.022057 0.11253 2.1448 2.072732 0.022057 0.11003 2.1448 2.093864 0.025468 0.10829 2.1367 2.1188

128 0.025539 0.11306 2.1377 2.1219

Table 3Results of normalized least-squares method

Mean error of the normalized least-squares method given in pixels



4 0.022054 0.12076 2.1475 2.02788 0.022054 0.11967 2.1475 2.0468

16 0.022054 0.1127 2.1475 2.073432 0.022054 0.11005 2.1475 2.094664 0.025464 0.10831 2.1392 2.1209

128 0.025538 0.11306 2.1401 2.124

out active contours, except for some values which vary a lot.But, for the other Radon-based methods, one can see a greatchange. This change does not necessarily concern the preci-sion obtained but rather the precision in general. Indeed, with-out pre-centering, the precision can vary a lot from one itera-tion cycle to the following iteration (e.g. from over 7 to lowerthan 0.5 pixels). The advantage of the pre-centering is to fur-nish more homogeneous precision along the iteration process(mostly lower than 0.5 pixel from one iteration to the follow-ing one). Nevertheless, the best precision is not always foundwith this method but the variation from the best result withpre-centering is often marginal (Tables 1–8 ).

Regarding the radius (Table 9), we obtain five values whichare given by Radon transform methods (the same way of cal-culus), least-squares methods (normalized or not, that is to saytwo values), and 3D circular Hough transform (subpixel andclassic version, that is to say two values). We propose studyingeach kind of circle to determine which radius determination isthe best for each type. First, we take an interest in the perfectdiscretized circle. It is obvious seeing the results that the preci-sion is maximal with least-square approximation, this withoutpre-processing by active contours. Indeed, active contours, inthis case, add imprecision as the contour is “perfectly” repre-sented. We can note that the non-normalized version is a little


Table 4Results of the initial Radon-based method

Mean error of the initial Radon-based method given in pixels



4 0.64247 0.586 1.2151 2.62778 0.64247 0.60038 1.2151 2.5572

16 0.64247 0.602 1.2151 2.415532 0.64247 0.54876 1.2151 2.595764 0.63675 0.51657 1.0754 2.3413

128 0.60914 0.54705 1.3182 2.0354

Table 5Results of the UPE Radon-based method

Mean error of the UPE Radon-based method given in pixels


Original 1st best precision 0.44843 {1, 4, 7, 10, 16, 19, 25} 0.53691 {16} 1.7334 {26} 2.5157 {20}2nd best precision 0.44989 {13, 22} 0.53986 {4} 2.0813 {25} 2.5185 {2}

Active contours 2 1st best precision 0.40309 {1, 4, 7, 10, 13, 16, 19, 22, 25} 0.38202 {25} 1.7597 {7, 16, 25} 1.9546 {23}2nd best precision 0.4885 {2, 11, 20} 0.38525 {7} 1.778 {26} 1.9633 {5}

4 1st best precision 0.40309 {1, 4, 7, 10, 13, 16, 19, 22, 25} 0.36178 {25} 1.7597 {7, 16, 25} 1.8678{23}2nd best precision 0.4885 {2, 11, 20} 0.37136 {16} 1.778 {26} 1.9788 {8}

8 1st best precision 0.40309 {1, 4, 7, 10, 13, 16, 19, 22, 25} 0.34517 {16} 1.7597 {7, 16, 25} 1.971 {23}2nd best precision 0.4885 {2, 11, 20} 0.34666 {25} 1.778 {26} 1.9822 {8}



64 1st best precision 0.31065 {1, 4, 7, 10, 13, 16, 19, 22, 25} 0.35546 {7} 1.4578 {26} 1.8298 {26}2nd best precision 0.46558 {2, 11, 20} 0.35927 {25} 1.6803 {7, 16, 25} 1.8651 {4}

128 1st best precision 0.38347 {1, 4, 7, 10, 13, 16, 19, 22, 25} 0.34059 {7} 1.6022 {26} 1.8021 {7}2nd best precision 0.48285 {2, 11, 20} 0.34207 {25} 1.6632 {17} 1.8059 {16}

Reference [�1, �2, �3] Reference [�1, �2, �3] Reference [�1, �2, �3]{1} [3,3,3] {10} [2,3,3] {19} [1,3,3]{2} [3,3,2] {11} [2,3,2] {20} [1,3,2]{3} [3,3,1] {12} [2,3,1] {21} [1,3,1]{4} [3,2,3] {13} [2,2,3] {22} [1,2,3]{5} [3,2,2] {14} [2,2,2] {23} [1,2,2]{6} [3,2,1] {15} [2,2,1] {24} [1,2,1]{7} [3,1,3] {16} [2,1,3] {25} [1,1,3]{8} [3,1,2] {17} [2,1,2] {26} [1,1,2]{9} [3,1,1] {18} [2,1,1] {27} [1,1,1]

more precise than the normalized one. Then, for noised circles,the active contours are useful as it rounds out the noised con-tour to form a well defined new family contour. One can seethat the precision is better for a certain level than others. Thislevel corresponds to a threshold of 25. The normalized versionof least-squares proves that it is better than the classic versionfor radius estimation when noise is added. The circular Houghmethod is also improved by active contours. For distorted andfor noised and distorted circles, active contours significantlyimproves the precision for a threshold of 27 to furnish subpixelradius estimation. The estimation of radius by Radon-based

method is interesting as it is also under a pixel precision but ata lower level than least-squares and circular Hough transform.Indeed, the computation of the radius has not already been fur-ther studied and needs some improvements, notably the takinginto account of deleted points in the calculations. It remains thatthe active contours approach implies a noticeable improvementin the estimation of the radius.

Regarding the influence of distortion, a better precision is ob-tained for a distortion of amplitude 5 that is obviously normal asthe distortion is less than before. But, depending on the methodused, the variation of the precision is not the same. Based on


Table 6Results of the FA Radon-based method

Mean error of the FA Radon-based method given in pixels


Original Precision without iteration 0.44946 0.52831 2.0436 3.0316Best precision (number of iteration) 0.44724 (2) 0.52831 (0) 1.5025 (8) 2.8024 (9)

Active contours 2 Precision without iteration 0.40495 0.37664 1.6016 2.0918Best precision (number of iteration) 0.40452 (2) 0.37137 (5) 0.77922 (8) 1.7336 (9)

4 Precision without iteration 0.40495 0.36582 1.6016 2.0457Best precision (number of iteration) 0.40452 (2) 0.36452 (3) 0.77922 (8) 1.643 (9)






Table 7Results of the FADE Radon-based method

Mean error of the FADE Radon-based method given in pixels


Original Precision without iteration 0.14009 1.2697 2.8697 3.3246Best precision (number of iteration) 0.14009 (0) 1.2697 (0) 2.4252 (5) 3.2846 (2)

Active contours 2 Precision without iteration 0.1997 0.65548 2.4793 2.7986Best precision (number of iteration) 0.1997 (0) 0.65548 (0) 2.1916 (5) 2.7773 (3)







Table 8Results of circular Hough method

Mean error of the circular Hough method given in pixels



4 2.9639 3.8698 5.4979 6.13038 2.9639 3.9188 5.4979 6.1806

16 2.9639 3.9119 5.4979 6.380532 2.9639 3.7139 5.4979 6.067764 3.1736 3.735 5.4388 6.1115

128 2.7204 3.6151 5.5259 5.961


Table 9Results on the radius: [LS]: least-squares method, [NLS]: normalized least-squares method, [CH]: circular Hough method and [R]: Radon-based method

Mean error of the radius given in pixels


Original Without centring 0.03992 [LS] 2.4143 [CH] 1.4214 [CH] 5.3286 [CH]With centring 0.03992 [LS] 2.4143 [CH] 1.4214 [CH] 5.3286 [CH]

Active contours 2 Without centring 0.54683 [R] 0.21429 [CH] 0.42143 [CH] 2.1558 [NLS]With centring 0.66825 [R] 0.21429 [CH] 0.42143 [CH] 2.1558 [NLS]

4 Without centring 0.54683 [R] 0.21429 [CH] 0.42143 [CH] 2.109 [NLS]With centring 0.66825 [R] 0.21429 [CH] 0.42143 [CH] 2.109 [NLS]

8 Without centring 0.54683 [R] 0.21429 [CH] 0.42143 [CH] 2.0493 [NLS]With centring 0.66825 [R] 0.21429 [CH] 0.42143 [CH] 2.0493 [NLS]

16 Without centring 0.54683 [R] 0.17404 [NLS] 0.42143 [CH] 1.9533 [NLS]With centring 0.66825 [R] 0.17404 [NLS] 0.42143 [CH] 1.9533 [NLS]

32 Without centring 0.54683 [R] 0.060943 [NLS] 0.42143 [CH] 1.6571 [CH]With centring 0.66825 [R] 0.060943 [NLS] 0.42143 [CH] 1.6571 [CH]

64 Without centring 0.80516 [R] 0.2524 [LS] 0.51074 [NLS] 1.5286 [CH]With centring 0.8869 [R] 0.2524 [LS] 0.51074 [NLS] 1.5286 [CH]

128 Without centring 1.0584 [LS] 0.49837 [R] 0.15228 [LS] 0.70374 [NLS]With centring 1.0584 [LS] 0.50052 [R] 0.15228 [LS] 0.70374 [NLS]

Table 10Overview of results method

Method Precision Computation velocity

Perfect Noised Distorted Noised and distorted

Least-mean squares + + + + + + − − + + +Normalized least-mean squares + + + + + + − − + + +barycenter + + − − + + +3D circular Hough transform − − − − − − − − − − −Initial Radon-based method+ active contours − − + + +UPC Radon-based method+ active contours + + + + −FA Radon-based method+ active contours + + + + + + + + −FADE Radon-based method+ active contours + + ++ ++ −

some results, it seems that subpixel Radon-based methods areless sensitive than other methods such as least-squares ones atan increase of distortion. On the contrary, least-squares meth-ods are less sensitive to noise variation (amplitude = 10) thanRadon-based methods but active contours partially succeeds incorrecting for the loss of precision. A comparative overview onprecision and computation velocity is available in Table 10.

5. Conclusions

In this study, we developed or improved several methods fordetermining circle characteristics, and compared them to clas-sic algorithms in circle recognition. Results are rather goodand very promising in a noised or/and distorted context, par-ticularly, for Radon transform-based methods whose results, inits more reliable version (FA), trump those of other methods.Taking into account the discretization effect appears to greatlyimprove global efficiency in a perfect circle but if some per-

turbations due to noise or distortion are present, the precisionis lower if discretization is considered. The use of active con-tours, in order to produce a family of rounded envelopes, ap-pears judicious in order to avoid or limit some noise/distortioneffects. This is a promising path of research as it can be elec-tronically combined with a sensor for a real-time implementa-tion. For the time being, the granularity parameter is set to 2so it still remains to automatically determine the level of ap-proximation according to the noise/distortion of the contour inorder to obtain the best precision.

References

[1] E. Andres, Discrete circles, rings and spheres, Comput. Graph. 18 (5)(1994) 695–706.

[2] J.E. Bresenham, A linear algorithm for incremental digital display ofcircular arcs, Commun. ACM 20 (2) (1977) 100–106.

[3] M.L.V. Pitteway, Integer circles, etc—some further thoughts, Comput.Graph. Image Process. 3 (1974) 262–265.


[4] F.L. Chen, S.W. Lin, Subpixel estimation of circle parameters usingorthogonal circular detector, Comput. Vision and Image Understanding78 (2) (2000) 206–221.

[5] R. Hanek, M. Beetz, The contracting curve density algorithm: fittingparametric curve models to images using local self-adapting separationcriteria, Int. J. Comput. Vis. 59 (3) (2004) 233–258.

[6] R. Hanek, The contracting curve density algorithm and its application tomodel-based image segmentation, IEEE Computer Society Conferenceon Computer Vision and Pattern Recognition (CVPR 2001), 2001, pp.797–804.

[7] L.R. Williams, K.K. Thornber, A comparison of measures for detectingnatural shapes in clustered backgrounds, Int. J. Comput. Vision 34 (2–3)(2000) 81–96.

[8] E.E. Zelniker, I.V.L. Clarkson, Approximate maximum-likelihoodestimation of circle parameters using a phase-coded kernel, XII EuropeanSignal Processing Conference, 2004, pp. 617–620.

[9] E.E. Zelniker, I.V.L. Clarkson, A statistical analysis of least-squares circlecentre estimation, IEEE International Symposium on Signal Processingand Information Technology (ISSPIT 2003), 2003, pp. 114–117.

[10] W. Gander, G.H. Golub, R. Strebel, Least-squares fitting of circles andellipses, BIT 34 (1994) 558–578.

[11] N. Chernov, C. Lesort, Least-squares fitting of circles, J math. Imaging.Vis. 23 (3) (2005) 239–252.

[12] D.S. Luo, P. Smart, J.E.S. Macleod, Circular hough transform forroundness measurement of objects, Pattern Recogn. 28 (11) (1995)1745–1749.

[13] S. Chiu, J. Liaw, An effective voting method for circle detection, PatternRecogn. Lett. 26 (2) (2005) 121–133.

[14] R.N. Bracewell, Two-Dimensional Imaging, Prentice-Hall, EnglewoodCliffs, NJ, 1995, pp. 505–537.

[15] J.S. Lim, Two-Dimensional Signal and Image Processing, Prentice-Hall,Englewood Cliffs, NJ, 1990, pp. 42–45.

[16] R.O. Duda, P.E. Hart, Use of the Hough transformation to detect linesand curves in pictures, ACM 15 (1972) 11–15.

[17] Y.T. Ching, Detecting line segments in an image: a new implementationfor Hough transform, Pattern Recogn. Lett. 22 (2001) 421–429.

[18] D. Ioannou, W. Huda, A.F. Lane, Circle recognition through a 2D HoughTransform and radius histogramming, Image Vis. Comput. 17 (1999)15–26.

[19] H. Kim, J. Kim, A two-step circle detection algorithm from theintersecting chords, Pattern Recogn. Lett. 22 (2001) 787–798.

[20] E. Kim, M. Haseyama, H. Kitajima, A new fast and robust circleextraction algorithm, 15th International Conference on Vision Interface(VI 2002), 2002, pp. 439–446.

[21] E. Kim, M. Haseyama, H. Kitajima, Extraction of circles from arcssegmented using short straight lines, International Conference onInformation Technology & Applications (ICITA 2002), 2002.

[22] A.C. Kak, M. Slaney, Principles of Computerized Tomographic Imaging,IEEE Press, New York, 1998.

[23] M.A. Branch, T.F. Coleman, Y. Li, A subspace, interior, and conjugategradient method for large-scale bound-constrained minimizationproblems, SIAM J. Sci. Comput. 21 (1) (1999) 1–23.

About the Author—FABRICE MAIRESSE is a Ph.D. student at the University of Burgundy (France) and is also a member of the Image Processing group ofthe Laboratory LE2I (“Laboratoire d’Electronique, Informatique et Image”). In 2004, he received a master engineering degree from ISTASE (“Institut supérieurdes techniques avancées de Saint-Etienne”) and a master research degree (of the University of Saint-Etienne) in image processing. His research interests areapplications of industrial control by image processing. His interests are metrology and mosaicking.

About the Author—TADEUSZ SLIWA received his Ph.D. degree in 2003 from the University of Burgundy, France. Since 2004, he has been an associateprofessor at the University of Burgundy, France. He is also a member of the Image Processing group of the Laboratory LE2I (“Laboratoire d’Electronique,Informatique et Image”). His research interests are in the applications of image processing.

About the Author—STÉPHANE BINCZAK received his Ph.D. degree in 1999 from the University of Burgundy, France. Since 2000, he has been an associateprofessor at the University of Burgundy, France. He is also a member of the Image Processing group of the Laboratory Le2i (“Laboratoire d’Electronique,Informatique et Image”). His main interest is nonlinear signal processing.

About the Author—YVON VOISIN received his Ph.D. degree in 1993 from the University of Franche-Comté, France. Since 2005, he has been a full professorat the University of Burgundy, France. He is also a member of the Image Processing group of the Laboratory Le2i (“Laboratoire d’Electronique, Informatiqueet Image”). His interests are the application of artificial vision and 3-D reconstruction and motion analysis.

Documents

Subpixel determination of imperfect circles characteristics