Smallest enclosing circle problem: Difference between revisions

Line 358:

function circle_from2(point a, b)

-- return the smallest circle that intersects 2 points:

point ~~centre~~ = sq_div(sq_add(a,b),2~~) -- (ie midpoint~~)

point midpoint = sq_div(sq_add(a,b),2)

atom ~~radius~~ = distance(a,b)/2

atom halfdiameter = distance(a,b)/2

circle res = { ~~centre~~, ~~radius~~ }

circle res = { midpoint, halfdiameter }

return res

end function

Line 366:

function circle_from3(point a, b, c)

-- return a unique circle that intersects three points

~~atom~~ ~~{{aX,aY},{bX,bY},{cX,cY}}~~ = ~~{a,~~sq_sub(b,a)~~,sq_sub(c,a)}~~,

point bma = sq_sub(b,a),

B = ~~bX * bX + bY * bY~~,

cma = sq_sub(c,a)

~~C =~~ cX * cX + cY * cY,

atom {{aX,aY},{bX,bY},{cX,cY}} = {a,bma,cma},

D = ~~bX * cY - bY * cX~~

B = sum(sq_power(bma,2)),

~~point~~ ~~centre~~ = { (cY * B - bY * C) / (2 * D) ~~+ aX~~,

C = sum(sq_power(cma,2)),

(bX * C - cX * ~~B) / (~~2 * D) + aY }

D = (bX*cY - bY*cX)*2

~~circle~~ ~~res~~ = { ~~centre,~~ ~~distance(centre~~, ~~a) }~~

point centre = {(cY*B - bY*C)/D + aX,

(bX*C - cX*B)/D + aY }

atom radius = distance(centre,a) -- (=== b,c)

circle res = { centre, radius }

return res

end function

Line 396:

Line 399:

procedure welzl(sequence p)

~~printf(1,"Points~~ %v ==> ~~centre %v radius %.14g\n",prepend(~~welzlr(shuffle(p),{}~~),p)~~)

circle c = welzlr(shuffle(p),{})

string s = sprintf("centre %v radius %.14g",c)

printf(1,"Points %v ==> %s\n",{p,s})

end procedure

constant tests = {{ { 0, 0 }, { 0, 1 }, { 1, 0 } },

constant tests = {{{0, 0},{ 0, 1},{ 1,0}},

{ { 5, -2 }, { -3, -2 }, { -2, 5 }, { 1, 6 }, { 0, 2 } }}

{{5,-2},{-3,-2},{-2,5},{1,6},{0,2}}}

papply(tests,welzl)</lang>