{-# LANGUAGE UndecidableInstances #-}
--------------------------------------------------------------------------------
-- |
-- Module      :  HGeometry.Direction
-- Copyright   :  (C) Frank Staals
-- License     :  see the LICENSE file
-- Maintainer  :  Frank Staals
--
-- Compuations that have to do with directions
--
--------------------------------------------------------------------------------
module HGeometry.Direction
  ( CanComputeNormalVector(..)
  , uniformDirection
  , uniformUpwardDirectionWrt
  , module HGeometry.Direction.Cardinal
  ) where

import Control.Lens
import HGeometry.Ball
import HGeometry.Direction.Cardinal
import HGeometry.Number.Radical
import HGeometry.Point
import HGeometry.Properties
import HGeometry.Triangle
import HGeometry.Vector
import System.Random.Stateful

--------------------------------------------------------------------------------

class (r ~ NumType geom
      ) => CanComputeNormalVector geom r | geom -> r where
  {-# mINIMAL normalVectorAt #-}
  -- | Given a point on the object, and the object. Compute the outward normal vector at
  -- the given point. The normal is returned as a unit vector.
  --
  -- pre: the query point lies on (the surface of) the object. (This is not checked)
  normalUnitVectorAt     :: ( Point_ point d r
                            , Has_ Metric_ d r
                            , d ~ Dimension geom
                            , Radical r, Fractional r)
                         => point -> geom -> Vector d r
  normalUnitVectorAt point
q geom
g = Vector d r -> Vector d r
forall vector (d :: Nat) r.
(Metric_ vector d r, Radical r, Fractional r) =>
vector -> vector
signorm (Vector d r -> Vector d r) -> Vector d r -> Vector d r
forall a b. (a -> b) -> a -> b
$ point -> geom -> Vector d r
forall point (d :: Nat).
(Point_ point d r, d ~ Dimension geom, Num r) =>
point -> geom -> Vector d r
forall geom r point (d :: Nat).
(CanComputeNormalVector geom r, Point_ point d r,
 d ~ Dimension geom, Num r) =>
point -> geom -> Vector d r
normalVectorAt point
q geom
g

  -- | Given a point on the object, and the object. Compute the outward normal vector at
  -- the given point. No Guarantees are given about the length of the resulting vector.
  --
  -- pre: the query point lies on (the surface of) the object. (This is not checked)
  normalVectorAt :: ( Point_ point d r, d ~ Dimension geom, Num r)
                 => point -> geom -> Vector d r

instance Point_ center d r => CanComputeNormalVector (Ball center) r where
  normalVectorAt :: forall point (d :: Nat).
(Point_ point d r, d ~ Dimension (Ball center), Num r) =>
point -> Ball center -> Vector d r
normalVectorAt point
q Ball center
b = (point
qpoint -> Getting (Point d r) point (Point d r) -> Point d r
forall s a. s -> Getting a s a -> a
^.Getting (Point d r) point (Point d r)
forall point (d :: Nat) r.
Point_ point d r =>
Lens' point (Point d r)
Lens' point (Point d r)
asPoint) Point d r -> Point d r -> Vector d r
forall point (d :: Nat) r.
(Affine_ point d r, Num r) =>
point -> point -> Vector d r
.-. (Ball center
bBall center
-> Getting (Point d r) (Ball center) (Point d r) -> Point d r
forall s a. s -> Getting a s a -> a
^.(center -> Const (Point d r) center)
-> Ball center -> Const (Point d r) (Ball center)
forall geom point. HasCenter geom point => Lens' geom point
Lens' (Ball center) center
center((center -> Const (Point d r) center)
 -> Ball center -> Const (Point d r) (Ball center))
-> ((Point d r -> Const (Point d r) (Point d r))
    -> center -> Const (Point d r) center)
-> Getting (Point d r) (Ball center) (Point d r)
forall b c a. (b -> c) -> (a -> b) -> a -> c
.(Point d r -> Const (Point d r) (Point d r))
-> center -> Const (Point d r) center
forall point (d :: Nat) r.
Point_ point d r =>
Lens' point (Point d r)
Lens' center (Point d r)
asPoint)

instance Point_ vertex 3 r => CanComputeNormalVector (Triangle vertex) r where
  normalVectorAt :: forall point (d :: Nat).
(Point_ point d r, d ~ Dimension (Triangle vertex), Num r) =>
point -> Triangle vertex -> Vector d r
normalVectorAt point
_ (Triangle vertex
u vertex
v vertex
w) = Vector 3 r -> Vector 3 r -> Vector 3 r
forall r. Num r => Vector 3 r -> Vector 3 r -> Vector 3 r
cross (vertex
v vertex -> vertex -> Vector 3 r
forall point (d :: Nat) r.
(Affine_ point d r, Num r) =>
point -> point -> Vector d r
.-. vertex
u) (vertex
w vertex -> vertex -> Vector 3 r
forall point (d :: Nat) r.
(Affine_ point d r, Num r) =>
point -> point -> Vector d r
.-. vertex
u)
  -- FIXME: shoudn't we check the orientation of the triangle: I think it may point the
  -- wrong way now


--------------------------------------------------------------------------------
-- * Generating directions


-- | Generates a unit vector corresponding to a direction.
--
-- The general idea is to generate a random vector in the [-1,1]^d, if it lies inside the
-- unit ball, scale it to be length 1, otherwise, retry.
uniformDirection     :: ( StatefulGen gen m
                        , Fractional r, Ord r, Radical r, UniformRange r
                        , Has_ Metric_ d r, Ord (Vector d r)
                        ) => gen -> m (Vector d r)
uniformDirection :: forall gen (m :: * -> *) r (d :: Nat).
(StatefulGen gen m, Fractional r, Ord r, Radical r, UniformRange r,
 Has_ Metric_ d r, Ord (Vector d r)) =>
gen -> m (Vector d r)
uniformDirection gen
gen = m (Vector d r)
go
  where
    generate' :: IxValue vector -> vector
generate' IxValue vector
x = (Int -> IxValue vector) -> vector
forall vector (d :: Nat) r.
Vector_ vector d r =>
(Int -> r) -> vector
generate (IxValue vector -> Int -> IxValue vector
forall a b. a -> b -> a
const IxValue vector
x)
    go :: m (Vector d r)
go = do v <- (Vector d r, Vector d r) -> gen -> m (Vector d r)
forall a g (m :: * -> *).
(UniformRange a, StatefulGen g m) =>
(a, a) -> g -> m a
forall g (m :: * -> *).
StatefulGen g m =>
(Vector d r, Vector d r) -> g -> m (Vector d r)
uniformRM (IxValue (Vector d r) -> Vector d r
forall {vector}.
Vector_ vector (Dimension vector) (IxValue vector) =>
IxValue vector -> vector
generate' (-r
1), IxValue (Vector d r) -> Vector d r
forall {vector}.
Vector_ vector (Dimension vector) (IxValue vector) =>
IxValue vector -> vector
generate' r
IxValue (Vector d r)
1) gen
gen
            if quadrance v <= 1 then pure $ signorm v
                                else go

-- | Generates an upward direction (with respect to the given normal/up vector) uniformly
-- at random.
uniformUpwardDirectionWrt            :: ( StatefulGen gen m
                                        , Fractional r, Ord r, Radical r, UniformRange r
                                        , Has_ Metric_ d r, Ord (Vector d r)
                                        ) => Vector d r -> gen -> m (Vector d r)
uniformUpwardDirectionWrt :: forall gen (m :: * -> *) r (d :: Nat).
(StatefulGen gen m, Fractional r, Ord r, Radical r, UniformRange r,
 Has_ Metric_ d r, Ord (Vector d r)) =>
Vector d r -> gen -> m (Vector d r)
uniformUpwardDirectionWrt Vector d r
normal gen
gen = gen -> m (Vector d r)
forall gen (m :: * -> *) r (d :: Nat).
(StatefulGen gen m, Fractional r, Ord r, Radical r, UniformRange r,
 Has_ Metric_ d r, Ord (Vector d r)) =>
gen -> m (Vector d r)
uniformDirection gen
gen m (Vector d r) -> (Vector d r -> Vector d r) -> m (Vector d r)
forall (f :: * -> *) a b. Functor f => f a -> (a -> b) -> f b
<&> \Vector d r
dir ->
  if (Vector d r
normal Vector d r -> Vector d r -> r
forall vector (d :: Nat) r.
(Metric_ vector d r, Num r) =>
vector -> vector -> r
`dot` Vector d r
dir) r -> r -> Bool
forall a. Ord a => a -> a -> Bool
>= r
0 then Vector d r
dir else Vector d r -> Vector d r
forall r vector (d :: Nat).
(Num r, Vector_ vector d r) =>
vector -> vector
negated Vector d r
dir