content: {-# LANGUAGE Arrows #-}
  {-# LANGUAGE FlexibleContexts #-}
  {-# LANGUAGE RecordWildCards #-}
  {-# LANGUAGE TypeFamilies #-}
  
  {- |
  Interpolation buffers.
  -}
  module FRP.Rhine.ResamplingBuffer.Interpolation where
  
  -- containers
  import Data.Sequence
  
  -- simple-affine-space
  import Data.VectorSpace
  
  -- time-domain
  import Data.TimeDomain (Diff)
  
  -- rhine
  import FRP.Rhine.ClSF
  import FRP.Rhine.ResamplingBuffer
  import FRP.Rhine.ResamplingBuffer.KeepLast
  import FRP.Rhine.ResamplingBuffer.Util
  
  -- | A simple linear interpolation based on the last calculated position and velocity.
  linear ::
    ( Monad m
    , Clock m cl1
    , Clock m cl2
    , VectorSpace v s
    , Num s
    , s ~ Diff (Time cl1)
    , s ~ Diff (Time cl2)
    ) =>
    -- | The initial velocity (derivative of the signal)
    v ->
    -- | The initial position
    v ->
    ResamplingBuffer m cl1 cl2 v v
  linear initVelocity initPosition =
    (derivativeFrom initPosition &&& clId) &&& timeInfoOf sinceInit
      ^->> keepLast ((initVelocity, initPosition), 0)
        >>-^ proc ((velocity, lastPosition), sinceInit1) -> do
          sinceInit2 <- timeInfoOf sinceInit -< ()
          let diff = sinceInit2 - sinceInit1
          returnA -< lastPosition ^+^ diff *^ velocity
  
  {- |
  sinc-Interpolation, or Whittaker-Shannon-Interpolation.
  
  The incoming signal is strictly bandlimited
  by the frequency at which @cl1@ ticks.
  Each incoming value is hulled in a sinc function,
  these are added and sampled at @cl2@'s ticks.
  In order not to produce a space leak,
  the buffer only remembers the past values within a given window,
  which should be chosen much larger than the average time between @cl1@'s ticks.
  -}
  sinc ::
    ( Monad m
    , Clock m cl1
    , Clock m cl2
    , VectorSpace v s
    , Ord s
    , Floating s
    , s ~ Diff (Time cl1)
    , s ~ Diff (Time cl2)
    ) =>
    -- | The size of the interpolation window
    --   (for how long in the past to remember incoming values)
    s ->
    ResamplingBuffer m cl1 cl2 v v
  sinc windowSize =
    historySince windowSize
      ^->> keepLast empty >>-^ proc as -> do
        sinceInit2 <- sinceInitS -< ()
        returnA -< vectorSum $ mkSinc sinceInit2 <$> as
    where
      mkSinc sinceInit2 (TimeInfo {..}, as) =
        let t = pi * (sinceInit2 - sinceInit) / sinceLast
         in (sin t / t) *^ as
      vectorSum = foldr (^+^) zeroVector
  
  -- TODO Do we want to give initial values?
  
  {- | Interpolates the signal with Hermite splines,
     using 'threePointDerivative'.
  
     Caution: In order to calculate the derivatives of the incoming signal,
     it has to be delayed by two ticks of @cl1@.
     In a non-realtime situation, a higher quality is achieved
     if the ticks of @cl2@ are delayed by two ticks of @cl1@.
  -}
  cubic ::
    ( Monad m
    , VectorSpace v s
    , Floating v
    , Eq v
    , Fractional s
    , s ~ Diff (Time cl1)
    , s ~ Diff (Time cl2)
    ) =>
    ResamplingBuffer m cl1 cl2 v v
  {- FOURMOLU_DISABLE -}
  cubic =
    ((delay zeroVector &&& threePointDerivative) &&& (sinceInitS >-> delay 0))
      >-> (clId &&& delay (zeroVector, 0))
     ^->> keepLast ((zeroVector, 0), (zeroVector, 0))
     >>-^ proc (((dv, v), t1), ((dv', v'), t1')) -> do
       t2 <- sinceInitS -< ()
       let
         t        = (t1 - t1') / (t2 - t1')
         tsquared = t ^ 2
         tcubed   = t ^ 3
         vInter   = ( 2 * tcubed - 3 * tsquared     + 1) *^  v'
                ^+^ (     tcubed - 2 * tsquared + t    ) *^ dv'
                ^+^ (-2 * tcubed + 3 * tsquared        ) *^  v
                ^+^ (     tcubed -     tsquared        ) *^ dv
       returnA -< vInter
  {- FOURMOLU_ENABLE -}
  