{-# LANGUAGE TypeFamilies #-}
module Stochastic.Distributions(
  UniformBase(rDouble)
  ,stdBase
  ,Empirical(..)
  ,mkEmpirical
  ) where

import System.Random
import Control.Monad.State.Lazy
import Stochastic.Tools

data UniformBase = UniformBase {
  rDouble :: (Double, UniformBase)
}


stdGen2Uni gen = UniformBase {
  rDouble = mapTuple
            (id)
            (stdGen2Uni)
            (randomR (0,1) gen)
  }

stdBase s = stdGen2Uni (mkStdGen s)


data Empirical = Empirical {
  degreesOfFreedom :: Int,
  cdf  :: Double -> Double,
  cdf' :: Double -> Double
  }

empiricalCDF' hist x = u - ((c - x)/s)
  where
    interval = head $ filter (\y -> cum_frequence y > x) $ hist
    u = upper_bound interval
    s = slope interval
    c = cum_frequence interval

  
                 
empiricalCDF hist x
  | x <  (lower_bound $ head hist) = 0
  | x >= (upper_bound $ head $ reverse hist) = 1
  | otherwise =
    f $ filter (\y -> lower_bound y <= x  && upper_bound y >= x) hist
  where
    f [] =
      error $ "x "++(show x)++" not in histogram range\n" ++
      (foldr (\h str -> (show (lower_bound h, upper_bound h, frequency h)) ++ "\n" ++ str) "" hist)
    f (y:ys) =
      let part = (x - (lower_bound y)) in
      let step = part * slope y        in
      (cum_frequence y) + (step * (rel_frequence y))
    


mkEmpirical :: [Double] -> Empirical
mkEmpirical samples = Empirical {
  degreesOfFreedom = length h,
  cdf  = empiricalCDF h,
  cdf' = empiricalCDF' h
  }
  where
    h = fIHistogram samples
    count :: Double
    count = fromInteger . toInteger $ length samples