From: Rick DeNatale <rick.denatale@...>
Date: 2008-02-21T04:10:50+09:00
Subject: Re: The Smallest Circle (#157)

Here's my solution.

I also eschewed doing any google research and started from scratch.  I
think I've taken a slightly different approach than any already
described.

The basic approach is to start with a circle containing only the first
point, and then add the points one by one changing the circle as
necessary.  So

First point
   Set center to the point and radius to zero.

Subsequent points.
   If the point is already in the circle then the point is simply
added to the collection of points and the circle is unchanged.  The
point is in the circle if the distance from the center of the circle
is less than the radius (within some epsilon to deal with floating
point comparisons).

  If the point is further from the center than the radius, the we know
that it must lie on the circumference of the new circle which will
contain all of the points examined thus far.  We now find the set of
points which are furthest from the new point (again within some
epsilon).  There will of course be one or more of these points.

  Now the center of the new circle is the average of those furthest
points plus the new points, and the radius is the distance between the
center and the new point.

The code is in several files:

==== lib/circle.rb =====
class Circle < Struct.new(:center, :radius)
  def to_s
    "{center:#{center}, radius:#{radius}}"
  end
end

==== lib/point.rb =====
class Point < Struct.new(:x, :y)

  def self.random
    new(rand, rand)
  end

  def to_s
    "Point(#{x}, #{y})"
  end

  def inspect
    to_s
  end

  def -(aPoint)
    Point.new(x - aPoint.x, y - aPoint.y)
  end

  def dist_squared(other_point)
    diff = self - other_point
    diff.x * diff.x + diff.y * diff.y
  end

  def dist(other_point)
    Math.sqrt(dist_squared(other_point))
  end

  def ==(other_point)
    x = other_point.x && y == other_point.y
  end

  def center(other_point)
    Point.new((x + other_point.x)/2, (y + other_point.y)/2)
  end

end

==== lib/point_collection.rb ====
require File.dirname(__FILE__) + '/point'
require File.dirname(__FILE__) + '/circle'

class PointCollection
  attr_reader :radius, :center

  def initialize
    @points = []
    @radius = 0
  end

  def circle
    Circle.new(center, radius)
  end

  def add_point(a_point)
    @points << a_point
    recalc_after_adding(a_point)
  end

  def empty?
    @points.empty?
  end

  private
  # Recalculate the center and radius if a point has been added.
  # If the new point is within the current circle then no change
  # is needed.
  # If it isn't then, it must be on the circumference of the new circle
  # In this case find the point or points furthest from the new point,
  # The new center will be the center of the new point plus these far points
  # The new radius will be the distance of the new point from the new center.
  def recalc_after_adding(a_point)
    if @center
      if a_point.dist_squared(center) > @radius * @radius
        circ_points = (farthest_from(a_point) + [a_point]).uniq
        @center = center_of(circ_points)
        @radius = @center.dist(a_point)
      end
    else
      @center = a_point
    end
  end

  # Given an array of points find the center.
  def center_of(an_array_of_points)
    x_sum = y_sum = 0.0
    an_array_of_points.each do |point|
        x_sum += point.x
        y_sum += point.y
    end
    Point.new(x_sum / an_array_of_points.length, y_sum /
an_array_of_points.length)
  end

  def farthest_from(a_point)
    max_dist_squared = 0.0
    delta = 0.000001
    farthest = []
    p @points
    @points.each do |point|
      dist_squared = point.dist_squared(a_point)
      diff = dist_squared - max_dist_squared
      if diff.abs <= delta
        max_dist_squared = dist_squared if diff > 0
        farthest << point
      elsif diff > delta
        max_dist_squared = dist_squared
        farthest = [ point ]
      end
    end
    farthest
  end
end

==== spec/point_spec.rb
require File.dirname(__FILE__) + '/../lib/point'

def reflect(point)
  Point.new(-point.x, -point.y)
end

describe Point, "(4, 5)" do

  before(:each) do
    @point = Point.new(4, 5)
  end

  it "should have a dist_squared of 25.0 from point(1, 1)" do
    @point.dist_squared(Point.new(1,1)).should == 25.0
  end

  it "should have a dist of 5 from point(1,1)" do
    @point.dist(Point.new(1,1)).should be_close(5.0, 0.0000001)

  end

  it "should equal a new point with the same coordinates" do
    @point.should == Point.new(4, 5)
  end

  it "should be its own center when paired with self" do
    @point.center(@point).should == @point
  end

  it "should have a center of zero when paired with its reflection" do
    @point.center(reflect(@point)).should == Point.new(0,0)
  end
end

==== spec/point_collection_spec.rb =====
require File.dirname(__FILE__) + '/../lib/point'
require File.dirname(__FILE__) + '/../lib/point_collection'

describe PointCollection do
  before(:each) do
    @pc = PointCollection.new
  end

  it "should be empty initially" do
    @pc.should be_empty
  end

  it "should should have a radius of zero" do
    @pc.radius.should == 0
  end

  it "should not have a center" do
    @pc.center.should == nil
  end

  describe PointCollection, "with one point" do
    before(:each) do
      @first_point = Point.new(1,1)
      @pc.add_point(@first_point)
    end

    it "should not be empty" do
      @pc.should_not be_empty
    end

    it "should have the first point as the center" do
      @pc.center.should == @first_point
    end
  end

  describe PointCollection, "with two points" do
    before(:each) do
      @first_point = Point.new(1,1)
      @pc.add_point(@first_point)
      @second_point = Point.new(5, 5)
      @pc.add_point(@second_point)
    end

    it "should have the center of the two points as its center" do
      @pc.center.should == @first_point.center(@second_point)
    end

    it "should have a radius half the distance between the two points" do
      @pc.radius.should be_close(@first_point.dist(@second_point) / 2.0, 0.001)
    end

    describe PointCollection, "adding a point inside" do

      before(:each) do
        @old_center, @old_radius = @pc.center, @pc.radius
        @pc.add_point(Point.new(4,4))
      end

      it "should have the same radius as before" do
        @pc.radius.should == @old_radius
      end

      it "should have the same center as before" do
        @pc.center.should == @old_center
      end
    end

    describe PointCollection, "adding a point ouside" do
      before(:each) do
        @old_center, @old_radius = @pc.center, @pc.radius
        @third_point = Point.new(6,6)
        @pc.add_point(@third_point)
      end

      it "should calculate the radius correctly" do
        @pc.radius.should
be_close(@first_point.dist(@third_point)/2.0, 0.0001)
      end
    end
  end
end
==============

This much was done about 5 days ago,  I've spent time sporadically
since then getting rmagick working with X11 on OS X Leopard so that I
could visually convince myself that this was working.

Here's the final piece, with the visualization code adapted from
someone elses submission (sorry but I can't seem to find it).

==== lib/main.rb =====
#! /usr/bin/env ruby

require 'rubygems'
require File.dirname(__FILE__) + '/point'
require File.dirname(__FILE__) + '/circle'
require File.dirname(__FILE__) + '/point_collection'

def generate_samples(n)
  (1..n).map { Point.random }
end

def encircle(points)
  point_collection = PointCollection.new
  points.each do |point|
    point_collection.add_point(point)
  end
  point_collection.circle
end

def scale(value)
  $size / 6 + value * $size / 1.5
end

def draw_circle(circle, color=false)
  center = circle.center
  radius = circle.radius
  $gc.stroke_color(color.to_s) if color
  $gc.circle(scale(center.x), scale(center.y), scale(center.x),
scale(center.y+radius))
end

def draw_point(center, color=false)
  $gc.stroke_color(color.to_s) if color
  $gc.circle(scale(center.x), scale(center.y), scale(center.x),
scale(center.y)+2)
end

def draw_line(p0, p1, color=false)
  $gc.stroke_color(color.to_s) if color
  $gc.line(scale(p0.x), scale(p0.y), scale(p1.x), scale(p1.y))
end


require 'RMagick'
$size = 600
$gc = Magick::Draw.new
$gc.fill_opacity(0)
$gc.stroke_width(2)
$gc.stroke_opacity(1)

samples = generate_samples(ARGV[0].to_i)
circle = encircle(samples)
puts circle

samples.each { |p| draw_point(p, :black) }
draw_circle(circle)
canvas = Magick::ImageList.new
canvas.new_image($size, $size)
$gc.draw(canvas)
canvas.display
=======


-- 
Rick DeNatale

My blog on Ruby
http://talklikeaduck.denhaven2.com/