グラハム走査

昨日は，包装アルゴリズムで凸包を求めるプログラムに取り組んでみましたが，今日は，グラハム(Graham)走査で作ってみました。
my_canvas.rbについては，昨日と同じです。

graham.rb

グラハム走査については，以下を参考にしました。

• 『アルゴリズムC〈第2巻〉探索・文字列・計算幾何』

# -*- coding: utf-8 -*-
require 'my_canvas'

# http://d.hatena.ne.jp/rubyco/20060609/zipwith
module Enumerable
  def zip_with(*others, &block)
    zip(*others).map &block
  end
end

# 線分p1p2と水平な直線がなす角度を求める(0〜360)
def theta(p1, p2)
  dx = p2[0] - p1[0]; ax = dx.abs
  dy = p2[1] - p1[1]; ay = dy.abs
  t = (ax + ay == 0) ? 0 : dy.to_f / (ax + ay)
  if (dx < 0)
    t = 2 - t
  elsif (dy < 0)
    t = 4 + t
  end
  return t * 90.0
end

def ccw(p0, p1, p2)
  dx1 = p1[0] - p0[0]; dy1 = p1[1] - p0[1]
  dx2 = p2[0] - p0[0]; dy2 = p2[1] - p0[1]
  return +1 if (dx1 * dy2 > dy1 * dx2)
  return -1 if (dx1 * dy2 < dy1 * dx2)
  return -1 if ((dx1 * dx2 < 0) || (dy1 * dy2 < 0))
  return +1 if ((dx1 * dx1 + dy1 * dy1) < (dx2 * dx2 + dy2 * dy2))
  return 0
end

# グラハム(Graham)走査
def grahamscan(points)
  path = []
  min = 0
  pts = Marshal.load(Marshal.dump(points)) # deep copy しておく。
  pts.each_with_index do |elem, i|         # y座標の最小の点を探す
    min = i if elem[1] < pts[min][1]
  end
  pts.each_with_index do |elem, i| # 一番右下の点を求める
    if elem[1] == pts[min][1]
      min = i if elem[0] > pts[min][0]
    end
  end
  start = pts.delete_at(min)
  path << start

  # start地点と各点となす角を求め，その角が小さい順にソートする。
  ths = pts.map {|ele| theta(path.last, ele)} 
  pts = pts.zip_with(ths) {|a, b| a << b}.sort_by {|ele| ele[2]}

  path << pts.shift
  pts.each do |ele|
    idx = path.length - 1
    while idx >= 1 && ccw(path[idx-1], path[idx], ele) >= 0
      idx = idx - 1
    end
    if idx == 0
      path << ele               # 凸包に追加
    else
      path.pop                  # 凸包から点を消去
      redo
    end
  end
  return path << start
end

points = [[ 3,  9], [11,  1], [ 6,  8], [ 4,  3],
          [ 5, 15], [ 8, 11], [ 1,  6], [ 7,  4],
          [ 9,  7], [14,  5], [10, 13], [16, 14],
          [15,  2], [13, 16], [ 3, 12], [12, 10]]

solved = grahamscan(points)


canvas = MyCanvas.new(20, 20)

points.each do |p|
  canvas.draw_point(p, 0.1, "#FF3030")
end

pre = solved.shift
canvas.draw_point(pre, 0.1)
solved.each do |p|
  canvas.draw_line(pre, p)
  canvas.draw_point(p, 0.2, "#FF3030")
  pre = p
end

canvas.write_to_png("graham.png")

実行結果