{-# LANGUAGE
    MultiParamTypeClasses,
    FlexibleInstances, FlexibleContexts,
    UndecidableInstances
  #-}

module Data.Random.Distribution.Triangular where

import Data.Random.RVar
import Data.Random.Distribution
import Data.Random.Distribution.Uniform

-- |A description of a triangular distribution - a distribution whose PDF
-- is a triangle ramping up from a lower bound to a specified midpoint
-- and back down to the upper bound.  This is a very simple distribution
-- that does not generally occur naturally but is used sometimes as an
-- estimate of a true distribution when only the range of the values and
-- an approximate mode of the true distribution are known.
data Triangular a = Triangular {
    -- |The lower bound of the triangle in the PDF (the smallest number the distribution can generate)
    forall a. Triangular a -> a
triLower  :: a,
    -- |The midpoint of the triangle (also the mode of the distribution)
    forall a. Triangular a -> a
triMid    :: a,
    -- |The upper bound of the triangle (and the largest number the distribution can generate)
    forall a. Triangular a -> a
triUpper  :: a}
    deriving (Triangular a -> Triangular a -> Bool
forall a. Eq a => Triangular a -> Triangular a -> Bool
forall a. (a -> a -> Bool) -> (a -> a -> Bool) -> Eq a
/= :: Triangular a -> Triangular a -> Bool
$c/= :: forall a. Eq a => Triangular a -> Triangular a -> Bool
== :: Triangular a -> Triangular a -> Bool
$c== :: forall a. Eq a => Triangular a -> Triangular a -> Bool
Eq, Int -> Triangular a -> ShowS
forall a. Show a => Int -> Triangular a -> ShowS
forall a. Show a => [Triangular a] -> ShowS
forall a. Show a => Triangular a -> String
forall a.
(Int -> a -> ShowS) -> (a -> String) -> ([a] -> ShowS) -> Show a
showList :: [Triangular a] -> ShowS
$cshowList :: forall a. Show a => [Triangular a] -> ShowS
show :: Triangular a -> String
$cshow :: forall a. Show a => Triangular a -> String
showsPrec :: Int -> Triangular a -> ShowS
$cshowsPrec :: forall a. Show a => Int -> Triangular a -> ShowS
Show)

-- |Compute a triangular distribution for a 'Floating' type.
floatingTriangular :: (Floating a, Ord a, Distribution StdUniform a) => a -> a -> a -> RVarT m a
floatingTriangular :: forall a (m :: * -> *).
(Floating a, Ord a, Distribution StdUniform a) =>
a -> a -> a -> RVarT m a
floatingTriangular a
a a
b a
c
    | a
a forall a. Ord a => a -> a -> Bool
> a
b     = forall a (m :: * -> *).
(Floating a, Ord a, Distribution StdUniform a) =>
a -> a -> a -> RVarT m a
floatingTriangular a
b a
a a
c
    | a
b forall a. Ord a => a -> a -> Bool
> a
c     = forall a (m :: * -> *).
(Floating a, Ord a, Distribution StdUniform a) =>
a -> a -> a -> RVarT m a
floatingTriangular a
a a
c a
b
    | Bool
otherwise = do
        let p :: a
p = (a
cforall a. Num a => a -> a -> a
-a
b)forall a. Fractional a => a -> a -> a
/(a
cforall a. Num a => a -> a -> a
-a
a)
        a
u <- forall a (m :: * -> *). Distribution StdUniform a => RVarT m a
stdUniformT
        let d :: a
d   | a
u forall a. Ord a => a -> a -> Bool
>= a
p    = a
a
                | Bool
otherwise = a
c
            x :: a
x   | a
u forall a. Ord a => a -> a -> Bool
>= a
p    = (a
u forall a. Num a => a -> a -> a
- a
p) forall a. Fractional a => a -> a -> a
/ (a
1 forall a. Num a => a -> a -> a
- a
p)
                | Bool
otherwise = a
u forall a. Fractional a => a -> a -> a
/ a
p
-- may prefer this: reusing u costs resolution, especially if p or 1-p is small and c-a is large.
--        x <- stdUniform
        forall (m :: * -> *) a. Monad m => a -> m a
return (a
b forall a. Num a => a -> a -> a
- ((a
1 forall a. Num a => a -> a -> a
- forall a. Floating a => a -> a
sqrt a
x) forall a. Num a => a -> a -> a
* (a
bforall a. Num a => a -> a -> a
-a
d)))

-- |@triangularCDF a b c@ is the CDF of @realFloatTriangular a b c@.
triangularCDF :: RealFrac a => a -> a -> a -> a -> Double
triangularCDF :: forall a. RealFrac a => a -> a -> a -> a -> Double
triangularCDF a
a a
b a
c a
x
    | a
x forall a. Ord a => a -> a -> Bool
< a
a
    = Double
0
    | a
x forall a. Ord a => a -> a -> Bool
<= a
b
    = forall a b. (Real a, Fractional b) => a -> b
realToFrac ((a
x forall a. Num a => a -> a -> a
- a
a)forall a b. (Num a, Integral b) => a -> b -> a
^(Int
2 :: Int) forall a. Fractional a => a -> a -> a
/ ((a
c forall a. Num a => a -> a -> a
- a
a) forall a. Num a => a -> a -> a
* (a
b forall a. Num a => a -> a -> a
- a
a)))
    | a
x forall a. Ord a => a -> a -> Bool
<= a
c
    = forall a b. (Real a, Fractional b) => a -> b
realToFrac (a
1 forall a. Num a => a -> a -> a
- (a
c forall a. Num a => a -> a -> a
- a
x)forall a b. (Num a, Integral b) => a -> b -> a
^(Int
2 :: Int) forall a. Fractional a => a -> a -> a
/ ((a
c forall a. Num a => a -> a -> a
- a
a) forall a. Num a => a -> a -> a
* (a
c forall a. Num a => a -> a -> a
- a
b)))
    | Bool
otherwise
    = Double
1

instance (RealFloat a, Ord a, Distribution StdUniform a) => Distribution Triangular a where
    rvarT :: forall (n :: * -> *). Triangular a -> RVarT n a
rvarT (Triangular a
a a
b a
c) = forall a (m :: * -> *).
(Floating a, Ord a, Distribution StdUniform a) =>
a -> a -> a -> RVarT m a
floatingTriangular a
a a
b a
c
instance (RealFrac a, Distribution Triangular a) => CDF Triangular a where
    cdf :: Triangular a -> a -> Double
cdf  (Triangular a
a a
b a
c) = forall a. RealFrac a => a -> a -> a -> a -> Double
triangularCDF a
a a
b a
c