--------------------------------------------------------------------------------
-- |
-- Module      :  Algorithms.Geometry.ConvexHull.Naive
  Copyright    :   ( C )  
-- License     :  see the LICENSE file
  Maintainer   :    
--------------------------------------------------------------------------------
module Algorithms.Geometry.ConvexHull.Naive( ConvexHull
                                           , lowerHull', lowerHullAll

                                           , isValidTriangle, upperHalfSpaceOf
                                           ) where

import           Control.Lens
import           Data.Ext
import           Data.Foldable (toList)
import           Geometry.HalfSpace
import           Geometry.HyperPlane
import           Geometry.Line
import           Geometry.Point
import           Geometry.Triangle
import           Geometry.Vector
import           Data.Intersection(intersects)
import           Data.List.NonEmpty (NonEmpty(..))
import           Data.List (find)
import           Data.Maybe (isNothing)
import           Data.Util
--------------------------------------------------------------------------------

type ConvexHull d p r = [Triangle 3 p r]

-- | Computes the lower hull without its vertical triangles.
--
  pre : The points are in general position . In particular , no four
-- points should be coplanar.
--
  running time :  )
lowerHull' :: forall r p. (Ord r, Fractional r, Show r)
           => NonEmpty (Point 3 r :+ p) -> ConvexHull 3 p r
lowerHull' = filter (not . isVertical) . lowerHullAll
  where
    toPt2   :: Point 3 r :+ p ->  Point 2 r :+ p
    toPt2 p = p&core %~ projectPoint

    isVertical (Triangle p q r) = ccw (toPt2 p) (toPt2 q) (toPt2 r) == CoLinear

-- | Generates a set of triangles to be used to construct a complete
-- convex hull. In particular, it may contain vertical triangles.
--
  pre : The points are in general position . In particular , no four
-- points should be coplanar.
--
  running time :  )
lowerHullAll                 :: forall r p. (Ord r, Fractional r, Show r)
                             => NonEmpty (Point 3 r :+ p) -> ConvexHull 3 p r
lowerHullAll (toList -> pts) = let mkT (Three p q r) = Triangle p q r in
    [ t | t <- mkT <$> uniqueTriplets pts, isNothing (isValidTriangle t pts) ]



_killOverlapping :: (Ord r, Fractional r) => [Triangle 3 p r] -> [Triangle 3 p r]
_killOverlapping = foldr keepIfNotOverlaps []
  where
    keepIfNotOverlaps t ts | any (t `overlaps`) ts = ts
                           | otherwise             = t:ts

overlaps :: (Fractional r, Ord r) => Triangle 3 p1 r -> Triangle 3 p2 r -> Bool
t1 `overlaps` t2 = upperHalfSpaceOf t1 == upperHalfSpaceOf t2 && False



-- | Tests if this is a valid triangle for the lower envelope. That
-- is, if all point lie above the plane through these points. Returns
-- a Maybe; if the result is a Nothing the triangle is valid, if not
-- it returns a counter example.
--
-- >>> let t = (Triangle (ext origin) (ext $ Point3 1 0 0) (ext $ Point3 0 1 0))
-- >>> isValidTriangle t [ext $ Point3 5 5 0]
-- Nothing
-- >>> let t = (Triangle (ext origin) (ext $ Point3 1 0 0) (ext $ Point3 0 1 0))
-- >>> isValidTriangle t [ext $ Point3 5 5 (-10)]
-- Just (Point3 5 5 (-10) :+ ())
isValidTriangle   :: (Num r, Ord r)
                  => Triangle 3 p r -> [Point 3 r :+ q] -> Maybe (Point 3 r :+ q)
isValidTriangle t = find (\a -> not $ (a^.core) `intersects` h)
  where
    h = upperHalfSpaceOf t


-- | Computes the halfspace above the triangle.
--
-- >>> upperHalfSpaceOf (Triangle (ext $ origin) (ext $ Point3 10 0 0) (ext $ Point3 0 10 0))
  HalfSpace { _ boundingPlane = HyperPlane { _ inPlane = Point3 0 0 0 , _ normalVec = Vector3 0 0 100 } }
upperHalfSpaceOf                  :: (Ord r, Num r) => Triangle 3 p r -> HalfSpace 3 r
upperHalfSpaceOf (Triangle p q r) = HalfSpace h
  where
    h' = from3Points (p^.core) (q^.core) (r^.core)
    c  = p&core.zCoord -~ 1
    h  = if (c^.core) `liesBelow` h' then h' else h'&normalVec %~ ((-1) *^)
    a `liesBelow` plane = (a `onSideUpDown` plane) == Below

------------------------------------------------------------------------------
 |
 Module      :  Algorithms.Geometry.ConvexHull.Naive
 License     :  see the LICENSE file
------------------------------------------------------------------------------
------------------------------------------------------------------------------
 | Computes the lower hull without its vertical triangles.

 points should be coplanar.

 | Generates a set of triangles to be used to construct a complete
 convex hull. In particular, it may contain vertical triangles.

 points should be coplanar.

 | Tests if this is a valid triangle for the lower envelope. That
 is, if all point lie above the plane through these points. Returns
 a Maybe; if the result is a Nothing the triangle is valid, if not
 it returns a counter example.

 >>> let t = (Triangle (ext origin) (ext $ Point3 1 0 0) (ext $ Point3 0 1 0))
 >>> isValidTriangle t [ext $ Point3 5 5 0]
 Nothing
 >>> let t = (Triangle (ext origin) (ext $ Point3 1 0 0) (ext $ Point3 0 1 0))
 >>> isValidTriangle t [ext $ Point3 5 5 (-10)]
 Just (Point3 5 5 (-10) :+ ())
 | Computes the halfspace above the triangle.

 >>> upperHalfSpaceOf (Triangle (ext $ origin) (ext $ Point3 10 0 0) (ext $ Point3 0 10 0))
  Copyright    :   ( C )  
  Maintainer   :    
module Algorithms.Geometry.ConvexHull.Naive( ConvexHull
                                           , lowerHull', lowerHullAll

                                           , isValidTriangle, upperHalfSpaceOf
                                           ) where

import           Control.Lens
import           Data.Ext
import           Data.Foldable (toList)
import           Geometry.HalfSpace
import           Geometry.HyperPlane
import           Geometry.Line
import           Geometry.Point
import           Geometry.Triangle
import           Geometry.Vector
import           Data.Intersection(intersects)
import           Data.List.NonEmpty (NonEmpty(..))
import           Data.List (find)
import           Data.Maybe (isNothing)
import           Data.Util

type ConvexHull d p r = [Triangle 3 p r]

  pre : The points are in general position . In particular , no four
  running time :  )
lowerHull' :: forall r p. (Ord r, Fractional r, Show r)
           => NonEmpty (Point 3 r :+ p) -> ConvexHull 3 p r
lowerHull' = filter (not . isVertical) . lowerHullAll
  where
    toPt2   :: Point 3 r :+ p ->  Point 2 r :+ p
    toPt2 p = p&core %~ projectPoint

    isVertical (Triangle p q r) = ccw (toPt2 p) (toPt2 q) (toPt2 r) == CoLinear

  pre : The points are in general position . In particular , no four
  running time :  )
lowerHullAll                 :: forall r p. (Ord r, Fractional r, Show r)
                             => NonEmpty (Point 3 r :+ p) -> ConvexHull 3 p r
lowerHullAll (toList -> pts) = let mkT (Three p q r) = Triangle p q r in
    [ t | t <- mkT <$> uniqueTriplets pts, isNothing (isValidTriangle t pts) ]



_killOverlapping :: (Ord r, Fractional r) => [Triangle 3 p r] -> [Triangle 3 p r]
_killOverlapping = foldr keepIfNotOverlaps []
  where
    keepIfNotOverlaps t ts | any (t `overlaps`) ts = ts
                           | otherwise             = t:ts

overlaps :: (Fractional r, Ord r) => Triangle 3 p1 r -> Triangle 3 p2 r -> Bool
t1 `overlaps` t2 = upperHalfSpaceOf t1 == upperHalfSpaceOf t2 && False



isValidTriangle   :: (Num r, Ord r)
                  => Triangle 3 p r -> [Point 3 r :+ q] -> Maybe (Point 3 r :+ q)
isValidTriangle t = find (\a -> not $ (a^.core) `intersects` h)
  where
    h = upperHalfSpaceOf t


  HalfSpace { _ boundingPlane = HyperPlane { _ inPlane = Point3 0 0 0 , _ normalVec = Vector3 0 0 100 } }
upperHalfSpaceOf                  :: (Ord r, Num r) => Triangle 3 p r -> HalfSpace 3 r
upperHalfSpaceOf (Triangle p q r) = HalfSpace h
  where
    h' = from3Points (p^.core) (q^.core) (r^.core)
    c  = p&core.zCoord -~ 1
    h  = if (c^.core) `liesBelow` h' then h' else h'&normalVec %~ ((-1) *^)
    a `liesBelow` plane = (a `onSideUpDown` plane) == Below
