Re: Re: Re: Re: Re: Re: FindMinimum and the minimum-radius circle

*To*: mathgroup at smc.vnet.net*Subject*: [mg50280] Re: [mg50266] Re: [mg50235] Re: [mg50213] Re: [mg50193] Re: Re: FindMinimum and the minimum-radius circle*From*: DrBob <drbob at bigfoot.com>*Date*: Tue, 24 Aug 2004 06:22:19 -0400 (EDT)*References*: <cfuq6g$653$1@smc.vnet.net> <200408180834.EAA08732@smc.vnet.net> <cg204q$msq$1@smc.vnet.net> <200408200857.EAA12448@smc.vnet.net> <200408210704.DAA24766@smc.vnet.net> <200408220419.AAA10264@smc.vnet.net> <200408231034.GAA14705@smc.vnet.net> <opsc6666wtiz9bcq@monster.cox-internet.com> <opsc68aguniz9bcq@monster.cox-internet.com> <opsc7f2og6iz9bcq@monster.cox-internet.com>*Reply-to*: drbob at bigfoot.com*Sender*: owner-wri-mathgroup at wolfram.com

I left this off that last post: Attributes[rationalize] = {Listable}; Bobby On Mon, 23 Aug 2004 16:45:50 -0500, DrBob <drbob at bigfoot.com> wrote: > Oops! Fred's routine returned a too-small circle for this data: > > data10={{0.00604914, -1.54375}, { > 0.201781, 0.7955}, {0.703343, -0.400924}, {-0.485319, -0.704386}, \ > {-0.710707, 0.905265}, {0.463783, 0.470773}, {-1.06043, 1.8205}, { > 0.104751, 0.235759}, {-1.89801, 1.03611}, {-0.330475, 0.931081}} > > The center returned is this: > > rationalize[x_] := Module[{n = Floor[Log[10, Abs[x]]]}, > Rationalize[x, 10^(n - 10)]] > N@smallestcircle@rationalize@data10 > > Circle[{-0.527188,0.138376},1.63862] > > smallestcircle@data10 > > Circle[{-0.527188,0.138376},1.63862] > > That center is midway between two points that probably should define a diameter of the correct circle--my plotter marked it as a diameter--but the radius is too small, so both points are outside the circle. > > Ten more samplings gave me another example, with the right center but a too-small radius: > > {{0.737908, -1.39895}, {-0.602274, -1.7662}, {-0.876949, 0.705707}, { > 0.697127, 0.393087}, {1.22223, 0.205725}, { > 0.483356, -0.22763}, {1.45802, 0.91891}, {1.73444, -0.814056}, {1.26108, \ > -1.69484}, {1.31919, 1.10769}} > > Bobby > > On Mon, 23 Aug 2004 13:57:42 -0500, DrBob <drbob at bigfoot.com> wrote: > >> Ah! A simple fix is to Rationalize the data before applying Fred's method. For instance: >> >> smallestcircle@Rationalize[data,10^-10] >> >> The second argument to Rationalize has to be large enough that all the data is Rational afterward, but small enough to lose nothing important in the approximation--which may be an art in itself, if we get really technical. >> >> Bobby >> >> On Mon, 23 Aug 2004 13:34:08 -0500, DrBob <drbob at bigfoot.com> wrote: >> >>> That loops infinitely on ten to twenty percent of my samples for n=3 or n=10, including the following triangle: >>> >>> sampler[n_Integer] := RandomArray[NormalDistribution[0, 1], {n, 2}] >>> data=sampler@3 >>> >>> {{1.48943,-0.0340306},{-0.743613,0.145276},{1.06966,1.15398}} >>> >>> The NMinimize method returns this center for the same data: >>> >>> circleFinder[data][[{3,2}]] >>> >>> {{0.388114,0.245001},1.13611} >>> >>> Bobby >>> >>> On Mon, 23 Aug 2004 06:34:59 -0400 (EDT), Fred Simons <f.h.simons at tue.nl> wrote: >>> >>>> The recent messages about the problem of finding the minimum radius circle >>>> are actually a discussion on the numerical aspects of the problem. The >>>> following remark is a little bit outside that scope. >>>> >>>> It is possible to construct that circle directly. I have to admit that I did >>>> not follow this thread very carefully, so if someone else already gave a >>>> similar construction, I strongly apologize. >>>> >>>> The idea is that the smallest circle is a circle that passes through three >>>> points or passes through two diametrical points (that has been remarked). A >>>> little bit more can be said: when the smallest circle passes through three >>>> points, then the center of the circle must be inside the triangle. So the >>>> construction runs as follows. First construct a circle that contains all >>>> points and passes through one point. Then move the center of the circle in >>>> the direction of that point until it passes through a second point. When >>>> these points are diametric, we are ready. Otherwise, move the center of the >>>> circle to the connecting line until it passes through a third point. Then >>>> test. If the center of the circle is outside the triangle, move it in the >>>> direction of the largest side of the triangle taking care that all points >>>> remain inside the circle, and so on. Here is an implementation: >>>> >>>> smallestcircle[points_] := >>>> Block[{m, p1, p2, p3, v, t, t0, pts = points}, >>>> m = Mean[points]; >>>> p1 = p2 = p3 = points[[Ordering[points, -1, >>>> (#1 - m) . (#1 - m) < (#2 - m) . (#2 - m) & ][[1]]]]; >>>> While[Length[Union[{p1, p2, p3}]] < 3 || >>>> (p1 + p2)/2 == m || (p1 - p2) . (p1 - p2) > >>>> (p1 - p3) . (p1 - p3) + (p2 - p3) . >>>> (p2 - p3), >>>> v = (p1 + p2)/2 - m; t = 1; p = p3; >>>> pts = >>>> Reap[Scan[If[v . (#1 - p1) < 0, >>>> Sow[#1]; t0 = ((m - #1) . (m - #1) - >>>> (m - p1) . (m - p1))/(2* >>>> v . (#1 - p1)); If[t0 < t, >>>> t = t0; p = #1]] & , pts]][[2]]; >>>> If[pts == {}, Break[], pts = pts[[1]]]; >>>> m = m + t*v; p3 = p; {p1, p2, p3} = >>>> If[(p3 - p1) . (p3 - p1) >= (p2 - p1) . >>>> (p2 - p1), {p1, p3, p2}, >>>> {p3, p2, p1}]]; >>>> Circle[m, Norm[m - p1]]] >>>> >>>> For exact coordinates we find the exact circle. This function is faster than >>>> the numerical ones: >>>> >>>> In[561]:= >>>> SeedRandom[1]; >>>> points = Partition[Table[Random[], {1000}], 2]; >>>> Timing[NMinimize[ >>>> Sqrt[Max[((x - #1[[1]])^2 + (y - #1[[2]])^ >>>> 2 & ) /@ points]], {x, y}]] >>>> Timing[smallestcircle[points]] >>>> >>>> Out[563]= >>>> {4.625*Second, {0.6628153793262672, >>>> {x -> 0.48683054240909296, >>>> y -> 0.46145830394871523}}} >>>> Out[564]= >>>> {0.14000000000000057*Second, >>>> Circle[{0.48684175891608683, >>>> 0.46146943938871915}, 0.6628150360165743]} >>>> >>>> >>>> Fred Simons >>>> Eindhoven University of Technology >>>> >>>> >>>> >>> >>> >>> >> >> >> > > > -- DrBob at bigfoot.com www.eclecticdreams.net

**References**:**Re: FindMinimum and the minimum-radius circle***From:*Thomas Burton <tburton@brahea.com>

**Re: Re: FindMinimum and the minimum-radius circle***From:*Thomas Burton <tburton@brahea.com>

**Re: Re: Re: FindMinimum and the minimum-radius circle***From:*"Janos D. Pinter" <jdpinter@hfx.eastlink.ca>

**Re: Re: Re: Re: FindMinimum and the minimum-radius circle***From:*Daniel Lichtblau <danl@wolfram.com>

**Re: Re: Re: Re: Re: FindMinimum and the minimum-radius circle***From:*"Fred Simons" <f.h.simons@tue.nl>