[Ruby] 三角形の外心・重心・内心・垂心を計算するプログラム ※ブラウザで動く版 http://hhiro.net/triangle/

triangle.rb

#!/usr/bin/env ruby

# 三角形の外心・重心・内心・垂心を計算するプログラム
# (C) 2012 H.Hiro(Maraigue) <main[at]hhiro.net>
# MIT Licenseの元で自由に利用可

require "matrix"

# 頂点を表すクラス
class Point
  # コンストラクタ
  def initialize(x = 0.0, y = 0.0)
    @x = (x.kind_of?(Integer) ? x.to_f : x)
    @y = (y.kind_of?(Integer) ? y.to_f : y)
  end
  
  def self.[](x, y)
    Point.new(x, y)
  end
  
  attr_accessor :x, :y
  
  # 比較
  def ==(other)
    @x == other.x && @y == other.y
  end
  
  # selfからdir方向に長さlenだけ進んだ点を得る
  # dirはラジアンで指定し、x軸の正の向きを0、
  # そこから反時計回りに増えていくとする
  def proceed(len, dir)
    Point[@x + len * Math.cos(dir), @y + len * Math.sin(dir)]
  end
  
  # 中点
  def midpoint(other)
    Point.new((@x + other.x) * 0.5, (@y + other.y) * 0.5)
  end
  
  # 角の二等分線
  # （2つの直線「self - target1」と「self - target2」の）
  def angle_bisection(target1, target2)
    dir = (Line[self, target1].direction + Line[self, target2].direction) / 2.0
    Line[self, self.proceed(1.0, dir)]
  end
  
  # selfを通るlineへの垂線
  def perpendicular(line)
    Line[self, self.proceed(1.0, line.direction + Math::PI / 2.0)]
  end
  
  # デバッグ用表示
  def to_s
    sprintf("(%f, %f)", @x, @y)
  end
  
  def inspect
    "<Point: #{to_s}>"
  end
end

# 線分（あるいは直線）を表すクラス
class Line
  # コンストラクタ
  def initialize(p1, p2)
    @p1 = p1
    @p2 = p2
  end
  
  def self.[](p1, p2)
    Line.new(p1, p2)
  end
  
  attr_accessor :p1, :p2
  
  # 長さ
  def length
    Math.sqrt((p1.x - p2.x) ** 2 + (p1.y - p2.y) ** 2)
  end
  
  # p1からp2に向かう向き（ラジアン）
  def direction
    if p1 == p2
      0.0
    else
      Math.atan2(p2.y - p1.y, p2.x - p1.x)
    end
  end
  
  # 中点
  def midpoint
    p1.midpoint(p2)
  end
  
  # 垂直二等分線
  def segment_bisection
    mp = midpoint
    half_len = length / 2
    dir = direction
    
    new_p1 = mp.proceed(half_len, dir + Math::PI / 2)
    new_p2 = mp.proceed(half_len, dir - Math::PI / 2)
    Line.new(new_p1, new_p2)
  end
  
  # 交点
  # ただし、selfもotherも「線分」ではなく「直線」とみなす
  # （すなわち、平行しているのでなければ必ず交わる）。
  def crossing(other)
    # l1 = (x1, y1)-(x2, y2), l2 = (x3, y3)-(x4, y4) のとき
    # (y2 - y1)X + (x1 - x2)Y = x1y2 - x2y1
    # (y4 - y3)X + (x3 - x4)Y = x3y4 - x4y3 を解けばよい。
    
    a = Matrix[[@p2.y - @p1.y, @p1.x - @p2.x],
               [other.p2.y - other.p1.y, other.p1.x - other.p2.x]]
    b = Matrix[[@p1.x * @p2.y - @p2.x * @p1.y],
               [other.p1.x * other.p2.y - other.p2.x * other.p1.y]]
    r = a.inv * b
    Point[r[0, 0], r[1, 0]]
  end
  
  # デバッグ用表示
  def inspect
    "<Line: #{@p1}-#{@p2}>"
  end
end

# 三角形を表すクラス
class Triangle
  # コンストラクタ
  def initialize(p1, p2, p3)
    @p1 = p1
    @p2 = p2
    @p3 = p3
  end
  
  def self.[](p1, p2, p3)
    Triangle.new(p1, p2, p3)
  end
  
  # 外心
  def circumcenter
    Line[@p1, @p2].segment_bisection.crossing(Line[@p1, @p3].segment_bisection)
  end
  
  # 重心
  def centroid
    Line[Line[@p1, @p2].midpoint, @p3].crossing(Line[Line[@p2, @p3].midpoint, @p1])
  end
  
  # 内心
  def incircle
    @p1.angle_bisection(@p2, @p3).crossing(@p2.angle_bisection(@p1, @p3))
  end
  
  # 垂心
  def orthocenter
    @p1.perpendicular(Line[@p2, @p3]).crossing(@p2.perpendicular(Line[@p1, @p3]))
  end
end

if $0 == __FILE__
  puts "---- Line(1) ----"
  l1 = Line.new(Point[1,0], Point[0,1])
  l2 = l1.segment_bisection
  p l1
  p l2
  p l1.crossing(l2) # 理論値は(0.5, 0.5)
  
  puts "---- Line(2) ----"
  l1 = Line.new(Point[1,0], Point[2,0])
  l2 = l1.segment_bisection
  p l1
  p l2
  p l1.crossing(l2) # 理論値は(1.5, 0.0)
  
  puts "---- Triangle ----"
  p1 = Point[3, 0]
  p2 = Point[-3, 0]
  p3 = Point[0, 9]
  t = Triangle[p1, p2, p3]
  p t.circumcenter # 理論値は(0.0, 4.0)
  p t.centroid     # 理論値は(0.0, 3.0)
  p t.incircle     # 理論値は(0.0, -1+√10)
  p t.orthocenter  # 理論値は(0.0, 1.0)
end