ch03::RealWorldHaskell(2)

haskell

11/3の続き。正解を探したけど見つからなかったので本当に正しいかどうかは不明。


ex.3-9

data Direction = RTurn  -- 右回り
               | LTurn  -- 左回り
               | Line   -- 一直線
               deriving (Eq, Show)


ex.3-10
3点a,b,cとした時、ベクトルbaから見てベクトルbcが時計回り(RTurn)か、反時計回り(LTurn)か、同じか正反対の方向(Line)かを判定する。
ベクトルba,bcがなす角をそれぞれθ1,θ2として、α=θ1-θ2 を考えると、

  • 2nπ < α < π+n2π (n = 0,±1,±2,...)　の時、左回り(LTurn) ※ この時 sinα > 0
  • -π+2nπ < α < 2nπ (n = 0,±1,±2,...) の時、右回り(RTurn) ※ この時 sinα < 0
  • α = nπ (n = 0,±1,±2,...) の時、一直線(Line) ※ この時 sinα = 0

という様な判定ができそう。
α(=θ1-θ2)は加法定理から求めれば良い感じ。
sin(α)=sin(θ1-θ2)=sin(θ1)cos(θ2)-cos(θ1)sin(θ2)

-- ex.3-10
-- (補助データ型)
-- 2次元平面上の点
data Point = Point Double Double
             deriving (Eq, Show)

-- 2次元平面上のベクトル
data Vec = Vec Double Double
           deriving (Eq, Show)

-- (補助関数)
-- 2乗を返す
square :: Double -> Double
square x = x * x

-- 2次元平面上の2点から大きさ(長さ)1のベクトルを取得
mkVec :: Point -> Point -> Vec
mkVec (Point x1 y1) (Point x2 y2) = Vec ((x2 - x1)/r) ((y2 - y1)/r)
   where r = sqrt ((square (x2 - x1)) + (square (y2 - y1)))

-- ベクトルが成す角をθとしてsinθを取得
mySin :: Vec -> Double
mySin (Vec x y) = y / (sqrt (square x) + (square y))

-- ベクトルが成す角をθとしてcosθを取得
myCos :: Vec -> Double
myCos (Vec x y) = x / (sqrt (square x) + (square y))

-- 2つのベクトルが成す角をαとしてsinαを取得
diffSin :: Vec -> Vec -> Double
diffSin v1 v2 = (mySin v1)*(myCos v2) - (myCos v1)*(mySin v2)

-- (これが求めたい関数、やっと本命!)
-- 2次元平面上の3点から方向を取得
direction :: Point -> Point -> Point -> Direction
direction p1 p2 p3
   | ds > 0  = RTurn
   | ds < 0  = LTurn
   | ds == 0 = Line
   where ds = diffSin (mkVec p2 p1) (mkVec p2 p3)

Main> let p1 = Point 10 10       -- 第1象限の点
Main> let p2 = Point (-10) 10    -- 第2象限の点
Main> let p3 = Point (-10) (-10) -- 第3象限の点
Main> let p4 = Point 10 (-10)    -- 第4象限の点
Main> let p0 = Point 1 1         -- 原点付近の点
Main>
Main> direction p1 p0 p2         -- 左回りになって欲しい
LTurn                            -- Good!
Main> direction p2 p0 p3         -- 左回りになってほしい
LTurn                            -- Good!
Main> direction p4 p0 p3         -- 今度は右回り
RTurn                            -- Good!
Main> direction p3 p0 p1         -- 直線はどうか
Line                             -- Good!
Main> direction p2 p0 p2         -- ベクトルが同じなら?
Line                             -- Nice!
Main>

なんだかゴチャゴチャしていてイケてない。もっとスマートな解法がありそうだけどとりあえずこれでいいや。


ex.3-11
directionを使えば簡単。

-- ex.3-11
mkDirectionList :: [Point] -> [Direction]
mkDirection [] = []      -- リストが空の場合
mkDirection [_] = []     -- 要素数が1の場合
mkDirection [_,_] = []   -- 要素数が2の場合
mkDirection (p1:p2:p3:ps) = (direction p1 p2 p3) : mkDirectionList (p2:p3:ps)

Main> mkDirectionList []
[]
Main> mkDirectionList [p1]
[]
Main> mkDirectionList [p1,p2]
[]
Main> mkDirectionList [p1,p0,p3]
[LTurn]
Main> mkDirectionList [p1,p2,p3,p4,p0]
[RTurn,RTurn,RTurn]
Main> mkDirectionList [p0,p4,p3,p2,p1]
[LTurn,LTurn,LTurn]
Main>


ex.3-12
graham-scanアルゴリズムって昔C++の本で読んだような気がする。さっぱり忘れてた!

-- ex.3-12
import Data.List (sortBy)

-- (補助関数)
-- Point型のデータをx,yでソート
sortByXY :: [Point] -> [Point]
sortByXY xs = sortBy compByX (sortBy CompByY xs)
   where compByX :: Point -> Point -> Ordering
         compByX (Point x1 _) (Point x2 _)
            | x1 > x2   = GT
            | otherwise = LT
         compByY :: Point -> Point -> Ordering
         compByY (Point _ y1) (Point _ y2)
            | y1 > y2   = GT
            | otherwise = LT

-- Point型の点と成す角でソート(ex.3-10で定義したmySin,mkVecを使用)
sortByTheta :: Point -> [Point] -> [Point]
sortByTheta p xs = sortBy compByTheta xs
   where compByTheta :: Point -> Point -> Ordering
         compByTheta p1 p2
            | mySin (mkVec p p1) < mySin (mkVec p p2) = GT
            | otherwise                                = LT

-- sortByXY,sortByThetaを使ってgraham-scanで各点を評価する順にソート
sortForGraham :: [Point] -> [Point]
sortForGraham xs = h : (sortByTheta h t)
   where ss = sortByXY xs
         h  = head ss
         t  = tail ss

-- リスト内の連続する3点が全て左回り(LTurn)かどうか判定
leftTurnTest :: [Direction] -> Bool
leftTurnTest []      = False
leftTurnTest [LTurn] = True
leftTurnTest [RTurn] = False
leftTurnTest [Line]  = False
leftTurnTest (x:xs)
   | x == LTurn = leftTurnTest xs
   | otherwise  = False

-- (Graham-scan)
-- 連続する3点が左回り以外の時、真ん中の点を削除する
grahamCheck :: [Point] -> [Point]
grahamCheck [] = []
grahamCheck [_] = []
grahamCheck [p1,p2] = [p1,p2]
grahamCheck (p1:p2:p3:ps)
   | (direction p1 p2 p3) == RTurn = grahamCheck (p1:p3:ps)
   | (direction p1 p2 p3) == Line  = grahamCheck (p1:p3:ps)
   | otherwise                     = p1 : grahamCheck (p2:p3:ps)

-- 全ての連続する3点が左回りになるまでgrahamCheckを繰り返す(ex.3-11で定義したmkDirectionListを使用)
grahamScan :: [Point] -> [Point]
grahamScan xs
   | leftTurnTest (mkDirectionList xs) = xs
   | otherwise                         = grahamScan $ grahamCheck $ sortForGraham xs

-- sample data
a1  = Point (-3) 0
a2  = Point 0 3
a3  = Point 3 0
a4  = Point 0 (-3)
a5  = Point 2 2
a6  = ...          -- a6以降はa1〜a5の内部に収まるように設定

aList = [a2,a4,a1,a5,a3,...]

Main> grahamScan aList
[Point (-3.0) 0.0,Point 0.0 3.0,Point 2.0 2.0,Point 3.0 0.0,Point 0.0 (-3.0)]
Main>

ちゃんと動いてるっぽい。