module Triangulation where

import Array

------------------------------------------------------------------------------
-- Polygon
------------------------------------------------------------------------------

type Point = (Double, Double)
type Polygon = Array Int Point

polygon :: [Point] -> Polygon
polygon ps = array (1, length ps) (zip [1..] ps)

dist :: Polygon -> Int -> Int -> Double
dist p i j | adjacent p i j = 0.0
dist p i j | otherwise      = dist' (p!i) (p!j)

dist' :: Point -> Point -> Double
dist' (x,y) (x',y') = sqrt ((x - x')^2 + (y - y')^2)

adjacent :: Polygon -> Int -> Int -> Bool
adjacent p i j = d == 1 || d == n - 1
    where
        d = abs (i - j)
        n = snd (bounds p)


------------------------------------------------------------------------------
-- Triangulation
------------------------------------------------------------------------------

triCosts :: Polygon -> Array (Int, Int) Double
triCosts p = cost
    where
        cost = array ((1,2), (n,n))
                     ([ ((i, s),
                         minimum [ cost ! (i, k + 1)
                                   + cost ! (i +% k, s - k)
                                   + dist p i (i +% k)
                                   + dist p (i +% k) (i +% (s - 1))
                                 | k <- [1..s-2]
                                 ]
                        )
                      | i <- [1..n], s <- [4..n]
                      ] 
                      ++ [ ((i, s), 0.0) | i <- [1..n], s <- [2..3] ]
                     )

        n = snd (bounds p)

        i +% k = (i + k - 1) `mod` n + 1


triCost :: Polygon -> Double
triCost p = triCosts p ! (1, snd (bounds p))


triCostsTrace :: Polygon -> Array (Int, Int) (Int, Double)
triCostsTrace p = cost
    where
        cost = array ((1,2), (n,n))
                     ([ ((i, s),
                         minimumTrace [ (k,
					 snd (cost ! (i, k + 1))
                                         + snd (cost ! (i +% k, s - k))
                                         + dist p i (i +% k)
                                         + dist p (i +% k) (i +% (s - 1))
                                        )
                                      |  k <- [1..s-2]
                                      ]
                        )
                      | i <- [1..n], s <- [4..n]
                      ] 
                      ++ [ ((i, s), (1, 0.0)) | i <- [1..n], s <- [2..3] ]
                     )

        n = snd (bounds p)

        i +% k = (i + k - 1) `mod` n + 1

        minimumTrace = foldr1 (\ix@(_,x) jy@(_,y) -> if x < y then ix else jy)


ppCosts cost = putStrLn (unlines (map show (assocs cost)))

------------------------------------------------------------------------------
-- Test data
------------------------------------------------------------------------------

-- Isosceles triangle
poly1 = polygon [(0.0, 0.0), (5.0, 10.0), (10.0, 10.0)]

-- Square
poly2 = polygon [(0.0, 0.0), (0, 10.0), (10.0, 10.0), (10.0, 0.0)]

poly3 = polygon [(2.0, 0.0), (0.0, 5.0), (4.0, 10.0), (6.0, 5.0)]

poly4 = polygon [(4.0, 0.0), (0.0, 5.0), (2.0, 10.0),
                 (6.0, 10.0), (10.0, 5.0), (8.0, 0.0)]

poly5 = polygon [(4.0, 0.0), (0.0, 5.0), (2.0, 10.0),
                 (8.0, 8.0), (10.0, 2.0), (8.0, 0.0)]