Drawing circles connecting four points

dotcircle.rb

#!/usr/bin/env ruby

# dotcircle.rb : Drawing circles connecting four points
#   by takehikom

# This program invokes ImageMagick's "convert" command.

class DotcircleDrawer
  def initialize
    @bgcolor = "#f0fff0"
    @dot_int = 100
    @margin = @dot_int / 2
    @image_side = @margin * 2 + @dot_int * 2
    @prefix = "dotcircle"
    @circle_color = "blue"
  end

  def start
    draw_base
    draw_circle1
    draw_circle2
    draw_circle3
    draw_circle4
    draw_circle5
  end

  def draw_base
    @filename_base = "#{@prefix}-base.png"

    command = "convert -size #{@image_side}x#{@image_side} \"xc:#{@bgcolor}\""
    command += " -fill black -stroke none -draw \""
    3.times do |i|
      y = @margin + @dot_int * i
      3.times do |j|
        x = @margin + @dot_int * j
        x2 = x + 6
        command += " circle #{x},#{y},#{x2},#{y}"
      end
    end
    command += "\" -quality 92 #{@filename_base}"

    puts command
    system command
  end

  def draw_circle(filename_out, x, y, r, draw_command = nil)
    command = "convert #{@filename_base}"
    command += " " + draw_command if draw_command
    vx1, vy1 = image_coordinate(x, y)
    vx2, vy2 = image_coordinate(x + r, y)
    command += " -fill none -stroke #{@circle_color} -strokewidth 3"
    command += " -draw \"circle #{vx1},#{vy1},#{vx2},#{vy2}\""
    command += " -quality 92 #{filename_out}"
    puts command
    system command
  end

  def image_coordinate(*a)
    a.map {|v| @margin + v * @dot_int}
  end

  def draw_circle1
    draw_circle("#{@prefix}1.png", 0.5, 0.5, Math::sqrt(2) * 0.5)

    vx1, vy1 = image_coordinate(0, 0)
    vx2, vy2 = image_coordinate(1, 0)
    vx3, vy3 = image_coordinate(0, 1)
    vx4, vy4 = image_coordinate(1, 1)
    com = " -fill none -stroke gray70 -strokewidth 2 -draw \""
    com += "line #{[vx1, vy1, vx4, vy4].join(',')}"
    com += " line #{[vx2, vy2, vx3, vy3].join(',')}"
    com += "\""
    draw_circle("#{@prefix}1a.png", 0.5, 0.5, Math::sqrt(2) * 0.5, com)
  end

  def draw_circle2
    draw_circle("#{@prefix}2.png", 1, 1, 1)

    vx1, vy1 = image_coordinate(1, 0)
    vx2, vy2 = image_coordinate(0, 1)
    vx3, vy3 = image_coordinate(2, 1)
    vx4, vy4 = image_coordinate(1, 2)
    com = " -fill none -stroke gray70 -strokewidth 2 -draw \""
    com += "line #{[vx1, vy1, vx4, vy4].join(',')}"
    com += " line #{[vx2, vy2, vx3, vy3].join(',')}"
    com += "\""
    draw_circle("#{@prefix}2a.png", 1, 1, 1, com)
  end

  def draw_circle3
    draw_circle("#{@prefix}3.png", 1, 1, Math::sqrt(2))

    vx1, vy1 = image_coordinate(0, 0)
    vx2, vy2 = image_coordinate(2, 0)
    vx3, vy3 = image_coordinate(0, 2)
    vx4, vy4 = image_coordinate(2, 2)
    com = " -fill none -stroke gray70 -strokewidth 2 -draw \""
    com += "line #{[vx1, vy1, vx4, vy4].join(',')}"
    com += " line #{[vx2, vy2, vx3, vy3].join(',')}"
    com += "\""
    draw_circle("#{@prefix}3a.png", 1, 1, Math::sqrt(2), com)
  end

  def draw_circle4
    draw_circle("#{@prefix}4.png", 0.5, 1, Math::sqrt(5) * 0.5)

    vx1, vy1 = image_coordinate(0, 0)
    vx2, vy2 = image_coordinate(1, 0)
    vx3, vy3 = image_coordinate(0, 2)
    vx4, vy4 = image_coordinate(1, 2)
    com = " -fill none -stroke gray70 -strokewidth 2 -draw \""
    com += "line #{[vx1, vy1, vx4, vy4].join(',')}"
    com += " line #{[vx2, vy2, vx3, vy3].join(',')}"
    com += "\""
    draw_circle("#{@prefix}4a.png", 0.5, 1, Math::sqrt(5) * 0.5, com)
  end

  def draw_circle5
    draw_circle("#{@prefix}5.png", 0.5, 1.5, Math::sqrt(10) * 0.5)

    vxo, vyo = image_coordinate(0.5, 1.5)
    vx1, vy1 = image_coordinate(0, 0)
    vx2, vy2 = image_coordinate(1, 0)
    vx3, vy3 = image_coordinate(2, 1)
    vx4, vy4 = image_coordinate(2, 2)
    vx5, vy5 = image_coordinate(0, 1)
    vx6, vy6 = image_coordinate(1, 1)
    vx7, vy7 = image_coordinate(0, 2)
    vx8, vy8 = image_coordinate(1, 2)
    com = " -fill none -stroke gray70 -strokewidth 2 -draw \""
    com += "line #{[vx5, vy5, vx8, vy8].join(',')}"
    com += " line #{[vx6, vy6, vx7, vy7].join(',')}"
    com += " line #{[vxo, vyo, vx1, vy1].join(',')}"
    com += " line #{[vxo, vyo, vx2, vy2].join(',')}"
    com += " line #{[vxo, vyo, vx3, vy3].join(',')}"
    com += " line #{[vxo, vyo, vx4, vy4].join(',')}"
    com += "\""
    draw_circle("#{@prefix}5a.png", 0.5, 1.5, Math::sqrt(10) * 0.5, com)
  end
end

class DotcircleCoordinate
  def initialize(size = nil)
    @size = size
  end

  def start
    init_dots(@size || 3)
    puts @p.map {|v| "%s(%d,%d)" % [v.label, v.x, v.y]}.join(" ")
    @c = {}

    @p.combination(4).each do |w, x, y, z|
      coord_a = [w, x, y, z].map {|v| [v.x, v.y]}
      circle_result = circle(*coord_a)
      if circle_result
        label4 = [w, x, y, z].map {|v| v.label}.join
        puts "%s => x=%1.3f, y=%1.3f, r=%1.3f" % ([label4] + circle_result)
        @c[label4] = circle_result
      end
    end
  end

  def init_dots(s)
    @@pos ||= Struct.new("Point", :label, :x, :y)
    case s
    when 3
    @p = [
      @@pos.new("A", 0, 0),
      @@pos.new("B", 1, 0),
      @@pos.new("C", 2, 0),
      @@pos.new("D", 0, 1),
      @@pos.new("E", 1, 1),
      @@pos.new("F", 2, 1),
      @@pos.new("G", 0, 2),
      @@pos.new("H", 1, 2),
      @@pos.new("I", 2, 2),
      ]
    when 4
    @p = [
      @@pos.new("A", 0, 0),
      @@pos.new("B", 1, 0),
      @@pos.new("C", 2, 0),
      @@pos.new("D", 3, 0),
      @@pos.new("E", 0, 1),
      @@pos.new("F", 1, 1),
      @@pos.new("G", 2, 1),
      @@pos.new("H", 3, 1),
      @@pos.new("I", 0, 2),
      @@pos.new("J", 1, 2),
      @@pos.new("K", 2, 2),
      @@pos.new("L", 3, 2),
      @@pos.new("M", 0, 3),
      @@pos.new("N", 1, 3),
      @@pos.new("O", 2, 3),
      @@pos.new("P", 3, 3),
      ]
    else
      raise
    end
  end

  # 2つの座標が同一のとき真を返す
  def samepoint(*a)
    return false if a.length < 2
    a.combination(2).each do |x, y|
      x_y = [y[0].to_f - x[0].to_f, y[1].to_f - x[1].to_f]
      return true if zero?(x_y[0]) && zero?(x_y[1])
    end
    false
  end

  # 3つの座標が同一直線上にあるとき真を返す
  def colinear(*a)
    return false if a.length < 3
    a.combination(3).each do |x, y, z|
      x_y = [y[0].to_f - x[0].to_f, y[1].to_f - x[1].to_f]
      y_z = [z[0].to_f - y[0].to_f, z[1].to_f - y[1].to_f]
      return true if zero?(x_y[1]) && zero?(y_z[1])
      return true if zero?(x_y[0] / x_y[1] - y_z[0] / y_z[1])
    end
    false
  end

  # 座標を全て通る円の[中心x,中心y,半径]を返す
  # 円が一意に定まらない場合にはnilを返す
  def circle(*a)
    return nil if a.length < 3
    return nil if samepoint(*a)
    return nil if colinear(*a)

    # 最初の3点から，[中心x,中心y,半径]を求める
    # 垂直二等分線の方程式 https://mathwords.net/suityokunitobun
    # (x2-x1)x + (y2-y1)y = (1/2)(x2^2-x1^2+y2^2-y1^2)
    c1 = (a[1][0] - a[0][0]).to_f
    puts "DEBUG: c1 = (#{a[1][0]} - #{a[0][0]}).to_f = #{c1}" if $DEBUG
    c2 = (a[1][1] - a[0][1]).to_f
    puts "DEBUG: c2 = (#{a[1][1]} - #{a[0][1]}).to_f = #{c2}" if $DEBUG
    c3 = (a[1][0]**2 - a[0][0]**2 + a[1][1]**2 - a[0][1]**2) * 0.5
    puts "DEBUG: c3 = (#{a[1][0]}**2 - #{a[0][0]}**2 + #{a[1][1]}**2 - #{a[0][1]}**2) * 0.5 = #{c3}" if $DEBUG
    c4 = (a[2][0] - a[1][0]).to_f
    puts "DEBUG: c4 = (#{a[2][0]} - #{a[1][0]}).to_f= #{c4}" if $DEBUG
    c5 = (a[2][1] - a[1][1]).to_f
    puts "DEBUG: c5 = (#{a[2][1]} - #{a[1][1]}).to_f = #{c5}" if $DEBUG
    c6 = (a[2][0]**2 - a[1][0]**2 + a[2][1]**2 - a[1][1]**2) * 0.5
    puts "DEBUG: c6 = (#{a[2][0]}**2 - #{a[1][0]}**2 + #{a[2][1]}**2 - #{a[1][1]}**2) * 0.5 = #{c6}" if $DEBUG
    puts "DEBUG: #{c1}x + #{c2}y = #{c3}, #{c4}x + #{c5}y = #{c6}" if $DEBUG
    # c1*x+c2*y=c3, c4*x+c5*5=c6を連立させてx, yを求める
    # 2元連立1次方程式の簡便解法 http://nicomar.style.coocan.jp/OPL/004.pdf
    d =  c1 * c5 - c2 * c4
    raise if zero?(d)
    x = (c3 * c5 - c2 * c6) / d
    y = (c1 * c6 - c3 * c4) / d
    puts "DEBUG: x=#{x}, y=#{y}" if $DEBUG
    # 1番目の座標を用いて半径を求める
    r = Math::sqrt((a[0][0] - x)**2 + (a[0][1] - y)**2)
    puts "DEBUG: r = Math::sqrt((#{a[0][0]} - #{x})**2 + (#{a[0][1]} - #{y})**2) = #{r}" if $DEBUG
    # 2番目以降の全ての点について，円周上にあるか確かめる
    a.each_with_index do |xy, i|
      next if i == 0
      r2 = Math::sqrt((xy[0] - x)**2 + (xy[1] - y)**2)
      puts "DEBUG: r2 = Math::sqrt((#{xy[0]} - #{x})**2 + (#{xy[1]} - #{y})**2) = #{r2}" if $DEBUG
      return false if !zero?(r2 - r)
    end

    [x, y, r]
  end

  def zero?(v)
    v.abs <= 1e-6
  end

  def colinear_print(*a)
    puts "colinear(#{a.inspect})... #{colinear(*a)}"
  end

  def samepoint_print(*a)
    puts "samepoint(#{a.inspect})... #{samepoint(*a)}"
  end

  def circle_print(*a)
    puts "circle(#{a.inspect})... #{circle(*a)}"
  end

  def self.test
    nc = self.new

    nc.colinear_print([0, 0], [1, 1], [2, 2])
    nc.colinear_print([0, 0], [1, 1], [2, 1])
    nc.colinear_print([0, 0], [1, 1], [2, 1], [2, 2])
    nc.colinear_print([0, 0], [1, 1], [2, 2], [2, 1])
    nc.colinear_print([1, 0], [0, 0], [2, 0])
    nc.colinear_print([2, 2], [0, 2], [1, 2])
    nc.colinear_print([0, 2], [0, 2], [0, 2])
    nc.colinear_print([0, 2], [0, 2], [0, -2])
    nc.colinear_print([0, 2], [0, 2], [0, -2], [0, 2])

    nc.samepoint_print([0, 0], [1, 1], [2, 2])
    nc.samepoint_print([0, 0], [1, 1], [0, 0], [2, 2])

    nc.circle_print([1, 0], [0, 1], [2, 1], [1, 2])
    nc.circle_print([1, 0], [0, 1], [2, 1], [1, 3])
    nc.circle_print([0, 0], [1, 0], [1, 1], [0, 1])
    nc.circle_print([0, 0], [2, 0], [2, 2], [0, 2])
    nc.circle_print([0, 0], [1, 0], [1, 2], [0, 2])
    nc.circle_print([0, 0], [1, 0], [2, 1], [2, 2])
  end
end

if __FILE__ == $0
  if ARGV.empty?
    DotcircleCoordinate.new(3).start
    exit
  end

  ARGV.each do |param|
    case param
    when /^-*d/i
      DotcircleDrawer.new.start
    when /^-*t/i
      DotcircleCoordinate.test
    when /^-*3/
      DotcircleCoordinate.new(3).start
    when /^-*4/
      DotcircleCoordinate.new(4).start
    when /^-*a/i
      DotcircleDrawer.new.start
      DotcircleCoordinate.test
      DotcircleCoordinate.new(3).start
      DotcircleCoordinate.new(4).start
    when /^-*h/i
      puts "usage: ruby #{$0} [3|4|draw|test|all|help]"
    end
  end
end