module Graphics.Triangulation.Triangulation where import Graphics.Formats.Collada.ColladaTypes import Graphics.Formats.Collada.Transformations (cycleNeighbours,cycleN) import Data.Array (Array(..), listArray, (!)) import Debug.Trace import List type Points = Array Int (Float,Float) type TriangulationFunction = Points -> [Int] -> [(Int,Int,Int)] data Tree = Node Int Int [Tree] type F2 = (Float,Float) instance Show Tree where show (Node c p tree) = "Node " ++ (show c) ++ " " ++ (show p) ++ "[" ++ (concat(map show tree)) ++ "]" -- | since there are a lot of triangulation algorithms -- a triangulation function can be passed triangulate :: TriangulationFunction -> Geometry -> Geometry triangulate f (Geometry name prims (Vertices vname ps ns)) = Geometry name (map triPoly prims) (Vertices vname ps ns) where triPoly (LP (LinePrimitive pIndices nIndices tex col)) = PL (LinePrimitive (tri pIndices) (normals pIndices nIndices) tex col) -- TO DO: other patterns tri pIndices = map (\(x,y,z) -> [x,y,z]) (concat (map (f arr) pIndices) ) normals pIndices nIndices = replicate (length (concat pIndices)) (head nIndices) -- TO DO: Why not (tri pIndices) arr = listArray (0,l-1) $ map (\(x,y,z) -> (x,z)) ps l = length ps -- | some triangulation algorithms on't support polygons with holes -- These polygons with (nested) holes have to be cut so that they consist of only one outline -- I.e. the chars a,b,d,e,g,o,p,q contain holes tat have to be deleted. deleteHoles :: Geometry -> Geometry deleteHoles (Geometry name prims (Vertices vname ps ns)) = Geometry name newPrims (Vertices vname ps ns) where newPrims = zipWith3 (\pInd tex col -> LP (LinePrimitive pInd pInd tex col)) flattenedTrees (map t prims) (map c prims) flattenedTrees = zipWith flatten trees indices arr = listArray (0,l-1) $ map (\(x,y,z) -> (x,z)) ps l = length ps trees = map (generateTrees arr insidePoly) indices pI (LP (LinePrimitive pIndices nIndices tex col)) = pIndices t (LP (LinePrimitive pIndices nIndices tex col)) = tex c (LP (LinePrimitive pIndices nIndices tex col)) = col indices = map pI prims flatten :: [Tree] -> [[Int]] -> [[Int]] flatten [] is = [] flatten ((Node c poly []):ts) is = (alternate c (pdir (is!!poly)) (is!!poly)) : (flatten ts is) flatten ((Node c poly ps):ts) is = (embed arr (flatten ps is) (alternate c (pdir (is!!poly)) (is!!poly))): (flatten ts is) pdir poly = polygonDirection (map (arr!) poly) -- |cut a polygon at a good position and insert the contained hole-polygon with opposite direction embed :: Points -> [[Int]] -> [Int] -> [Int] embed _ [] poly = poly embed points (s:sub_polys) poly = embed points sub_polys ((take (n+1) poly) ++ s ++ [head s] ++ (drop n poly)) where n = fst (rotatePoly (head s) points poly) -- |make sure that direction (clockwise or ccw) of polygons alternates depending on the nesting number c of poly alternate :: Int -> Bool -> [Int] -> [Int] alternate c b poly | (b && (even c)) || (not b && (odd c)) = poly | otherwise = reverse poly -- |f should be the funtion to test "contains" -- the trees then are the hierarchy of containedness of outlines generateTrees :: Points -> (Points -> [Int] -> [Int] -> Bool) -> [[Int]] -> [Tree] generateTrees points f [] = [] generateTrees points f ps = (treesList points containedPolys []) ++ (map (\x -> Node 0 x []) separateOutlines) where containedPolys = filter (\[p0,p1] -> f points (ps!!p0) (ps!!p1)) (combi ++ (map reverse combi)) combi = combinationsOf 2 [0..((length ps)-1)] -- all 2-subsets i.e. [[0,1],[0,2],[1,2]] separateOutlines = ([0..((length ps)-1)]) \\ (nub $ concat containedPolys) -- separate outlines don't contain other outlines treesList :: Points -> [[Int]] -> [Tree] -> [Tree] treesList points [] trees = trees treesList points ([x,y]:cs) trees = treesList points cs (insertTrees points [x,y] trees) insertTrees :: Points -> [Int] -> [Tree] -> [Tree] insertTrees points [x,y] trees | or (map fst ins) = map snd ins | otherwise = (map snd ins) ++ [ Node 0 y [Node 1 x []] ] where ins = map (insertTree points [x,y]) trees insertTree :: Points -> [Int] -> Tree -> (Bool, Tree) insertTree points [x,y] (Node c i []) | y == i = (True, Node c i [Node (c+1) x []] ) | otherwise = (False, Node c i []) insertTree points [x,y] (Node c i trees) | y == i = (True, Node c i ((Node (c+1) x []):trees) ) | otherwise = (b, Node c i (map snd subtrees)) where subtrees = map (insertTree points [x,y]) trees b = or (map fst subtrees) -- subsets of size k -- borrowed from David Amos' library: HaskellForMaths combinationsOf 0 _ = [[]] combinationsOf _ [] = [] combinationsOf k (x:xs) = map (x:) (combinationsOf (k-1) xs) ++ combinationsOf k xs -- |how many positions to rotate a polygon until the start point is nearest to some other point -- call i.e. with nearest (3,4) [(0,0),(1,2), ... ] 0 0 rotatePoly :: Int -> Points -> [Int] -> (Int,Float) rotatePoly p points poly = (fst tup, snd tup) where tup = nearest (points!p) (map (points!) poly) (-1) 0 0 nearest :: F2 -> [F2] -> Float -> Int -> Int -> (Int,Float) nearest _ [] dist l ml = (ml,dist) nearest (x0,y0) ((x1,y1):ps) dist l ml | (newDist < dist) || (dist < 0) = nearest (x0,y0) ps newDist (l+1) l | otherwise = nearest (x0,y0) ps dist (l+1) ml where newDist = (x0-x1)*(x0-x1)+(y0-y1)*(y0-y1) -- | returns True iff the first point of the first polygon is inside the second poylgon insidePoly :: Points -> [Int] -> [Int] -> Bool insidePoly _ [] _ = False insidePoly _ _ [] = False insidePoly points poly1 poly2 = pointInside (points!(head poly1)) (map (points!) poly2) -- | A point is inside a polygon if it has an odd number of intersections with the boundary (Jordan Curve theorem) pointInside :: F2 -> [F2] -> Bool pointInside (x,y) poly = (length intersectPairs) `mod` 2 == 1 where intersectPairs = [ p | p <- allPairs, positiveXAxis p, aboveBelow p] --, specialCases p] allPairs = cycleNeighbours poly positiveXAxis p = (x0 p) > x || (x1 p) > x -- intersect with positive x-axis -- only lines with one point above + one point below can intersect aboveBelow p = (((y0 p)> y && (y1 p)< y) || ((y0 p) < y && (y1 p) > y)) specialCases p = (((dir1 p) > 0 && (dir2 p) <= 0) || ((dir1 p) <= 0 && (dir2 p) > 0))-- cross product for special cases dir1 p = cross ((x1 p)-(x0 p),(y1 p)-(y0 p)) (1,0) dir2 p = cross ((x1 p)-(x0 p),(y1 p)-(y0 p)) (x-(x0 p),y-(y0 p)) cross (a0,b0) (a1,b1) = a0*b1 - a1*b0 x0 p = fst (head p) x1 p = fst (last p) y0 p = snd (head p) y1 p = snd (last p) -- | the direction of a polygon can be obtained by looking at a maximal point -- returns True if counterclockwise -- False if clockwise polygonDirection :: [F2] -> Bool polygonDirection poly | dir > 0 = True | dir < 0 = False | dir ==0 = (fst (p!!lminus) > fst (p!!lplus)) || (snd (p!!lminus) < snd (p!!lplus)) where p = nub poly dir = area2 (p!!lminus) (p!!l) (p!!lplus) l = maxim p 0 0 (0,0) lminus = (l-1) `mod` (length p) lplus = (l+1) `mod` (length p) -- the index of the right-/upmost point maxim [] count ml (mx,my) = ml maxim ((x,y):xs) count ml (mx,my) | (x > mx) || (x >= mx && (y > my)) = maxim xs (count+1) count (x,y) | otherwise = maxim xs (count+1) ml (mx,my) isRightTurnOrOn m x p = (area2 m x p) <= 0 isLeftTurn :: F2 -> F2 -> F2 -> Bool isLeftTurn m x p = (area2 m x p) > 0 area2 (x2,y2) (x0,y0) (x1,y1) = (x1-x0)*(y2-y0)-(x2-x0)*(y1-y0)