{-# LANGUAGE TypeFamilies #-}
{-# LANGUAGE TypeOperators #-}
{-# LANGUAGE Rank2Types #-}
{- |
<http://arxiv.org/abs/0911.5171>
-}
module Synthesizer.LLVM.CausalParameterized.Helix (
   -- * time and phase control based on the helix model
   static,
   staticPacked,
   dynamic,
   dynamicLimited,

   -- * useful control curves
   zigZag,
   zigZagPacked,
   zigZagLong,
   zigZagLongPacked,
   ) where

import qualified Synthesizer.LLVM.CausalParameterized.ProcessValue as CausalPV
import qualified Synthesizer.LLVM.CausalParameterized.ProcessPacked as CausalPS
import qualified Synthesizer.LLVM.CausalParameterized.ProcessPrivate
                                                              as CausalPrivP
import qualified Synthesizer.LLVM.CausalParameterized.Process as CausalP
import qualified Synthesizer.LLVM.CausalParameterized.Functional as Func
import qualified Synthesizer.LLVM.Parameterized.SignalPacked as SigPS
import qualified Synthesizer.LLVM.Parameterized.SignalPrivate as SigP
import qualified Synthesizer.LLVM.CausalParameterized.RingBufferForward as RingBuffer
import qualified Synthesizer.LLVM.Frame.SerialVector as Serial
import qualified Synthesizer.LLVM.Simple.Value as Value
import qualified Synthesizer.LLVM.Interpolation as Ip
import qualified Synthesizer.LLVM.Parameter as Param
import Synthesizer.LLVM.CausalParameterized.Functional (($&), (&|&), )
import Synthesizer.LLVM.CausalParameterized.Process (($*), ($<), )
import Synthesizer.LLVM.Simple.Value ((%>), (%>=), (?), (??), )

import qualified Synthesizer.LLVM.Storable.Vector as SVU
import qualified Data.StorableVector as SV

import qualified LLVM.Extra.ScalarOrVector as SoV
import qualified LLVM.Extra.Vector as Vector
import qualified LLVM.Extra.Arithmetic as A
import qualified LLVM.Extra.Memory as Memory
import qualified LLVM.Extra.MaybeContinuation as MaybeCont
import LLVM.Extra.Class (MakeValueTuple, ValueTuple, )

import qualified LLVM.Core as LLVM
import LLVM.Core (CodeGenFunction, Value, IsSized, IsFloating, )

import qualified Type.Data.Num.Decimal as TypeNum
import Type.Data.Num.Decimal (D1, )

import Data.Word (Word32, )

import Foreign.ForeignPtr (touchForeignPtr, )
import Foreign.Storable (Storable, )

import Control.Arrow (first, (<<<), (^<<), (<<^), )
import Control.Category (id, )
import Control.Applicative (liftA2, )
import Control.Functor.HT (unzip, )
import Data.Traversable (mapM, )
import Data.Tuple.HT (mapFst, )

import qualified Algebra.Field as Field
import qualified Algebra.Ring as Ring

import NumericPrelude.Numeric hiding (splitFraction, )
import NumericPrelude.Base hiding (unzip, zip, mapM, id, )

import Prelude ()


{- |
Inputs are @(shape, phase)@.

The shape parameter is limited at the beginning and at the end
such that only available data is used for interpolation.
Actually, we allow almost one step less than possible,
since the right boundary of the interval of admissible @shape@ values is open.
-}
static ::
   (Storable vh, MakeValueTuple vh, ValueTuple vh ~ v, Memory.C v,
    Ip.C nodesStep, Ip.C nodesLeap,
    SoV.RationalConstant a, SoV.Fraction a,
    Storable a, MakeValueTuple a, ValueTuple a ~ Value a,
    Memory.FirstClass a, Memory.Stored a ~ am, IsSized am,
    LLVM.NumberOfElements a ~ D1) =>
   (forall r. Ip.T r nodesLeap (Value a) v) ->
   (forall r. Ip.T r nodesStep (Value a) v) ->
   Param.T p Int ->
   Param.T p a ->
   Param.T p (SV.Vector vh) ->
   CausalP.T p (Value a, Value a) v
static ipLeap ipStep periodInt period vec =
   let period32 = Param.word32 periodInt
       cellMargin = combineMarginParams ipLeap ipStep periodInt
   in  interpolateCell ipLeap ipStep
       <<<
       first (peekCell cellMargin period32 vec)
       <<<
       flattenShapePhaseProc period32 period
       <<<
       first
          (limitShape cellMargin period32
              (Param.word32 $ fmap SV.length vec))


staticPacked ::
   (Storable vh, MakeValueTuple vh, ValueTuple vh ~ ve,
    Serial.Element v ~ ve, Memory.C ve,
    Ip.C nodesStep, Ip.C nodesLeap,
    Serial.Size (nodesLeap (nodesStep v)) ~ n,
    Serial.C (nodesLeap (nodesStep v)),
    Serial.Element (nodesLeap (nodesStep v)) ~
       nodesLeap (nodesStep (Serial.Element v)),
    TypeNum.Positive n,
    SoV.RationalConstant a, SoV.Fraction a, Vector.Real a,
    Storable a, MakeValueTuple a, ValueTuple a ~ Value a,
    Memory.FirstClass a, Memory.Stored a ~ am, IsSized am,
    LLVM.IsPrimitive a, LLVM.IsPrimitive am) =>
   (forall r. Ip.T r nodesLeap (Serial.Value n a) v) ->
   (forall r. Ip.T r nodesStep (Serial.Value n a) v) ->
   Param.T p Int ->
   Param.T p a ->
   Param.T p (SV.Vector vh) ->
   CausalP.T p (Serial.Value n a, Serial.Value n a) v
staticPacked ipLeap ipStep periodInt period vec =
   let period32 = Param.word32 periodInt
       cellMargin = combineMarginParams ipLeap ipStep periodInt
   in  interpolateCell ipLeap ipStep
       <<<
       first (CausalPS.pack
          (peekCell (fmap elementMargin cellMargin) period32 vec))
       <<<
       flattenShapePhaseProcPacked period32 period
       <<<
       first
          (limitShapePacked cellMargin period32
              (Param.word32 $ fmap SV.length vec))


{- |
In contrast to 'dynamic' this one ends
when the end of the manipulated signal is reached.
-}
dynamicLimited ::
   (Ip.C nodesStep, Ip.C nodesLeap,
    A.Additive v, Memory.C v,
    SoV.RationalConstant a, SoV.Fraction a,
    LLVM.CmpRet a, LLVM.CmpResult a ~ Bool,
    Storable a, MakeValueTuple a, ValueTuple a ~ Value a,
    Memory.FirstClass a, Memory.Stored a ~ am, IsSized am,
    LLVM.NumberOfElements a ~ D1) =>
   (forall r. Ip.T r nodesLeap (Value a) v) ->
   (forall r. Ip.T r nodesStep (Value a) v) ->
   Param.T p Int ->
   Param.T p a ->
   SigP.T p v ->
   CausalP.T p (Value a, Value a) v
dynamicLimited ipLeap ipStep periodInt period sig =
   dynamicGen
      (\cellMargin (skips, fracs) ->
         let windows =
                RingBuffer.trackSkip (fmap Ip.marginNumber cellMargin) sig $& skips
         in  (windows,
              CausalP.delay1Zero $& skips,
              CausalP.delay1Zero $& fracs))
      ipLeap ipStep periodInt period

{- |
If the time control exceeds the end of the input signal,
then the last waveform is locked.
This is analogous to 'static'.
-}
dynamic ::
   (Ip.C nodesStep, Ip.C nodesLeap,
    A.Additive v, Memory.C v,
    SoV.RationalConstant a, SoV.Fraction a,
    LLVM.CmpRet a, LLVM.CmpResult a ~ Bool,
    Storable a, MakeValueTuple a, ValueTuple a ~ Value a,
    Memory.FirstClass a, Memory.Stored a ~ am, IsSized am,
    LLVM.NumberOfElements a ~ D1) =>
   (forall r. Ip.T r nodesLeap (Value a) v) ->
   (forall r. Ip.T r nodesStep (Value a) v) ->
   Param.T p Int ->
   Param.T p a ->
   SigP.T p v ->
   CausalP.T p (Value a, Value a) v
dynamic ipLeap ipStep periodInt period sig =
   dynamicGen
      (\cellMargin (skips, fracs) ->
         let {-
             For conformance with 'static'
             we stop one step before the definite end.
             We achieve this by using a buffer
             that is one step longer than necessary.
             -}
             ((running, actualSkips), windows) =
                mapFst unzip $ unzip $
                RingBuffer.trackSkipHold
                   (fmap (succ . Ip.marginNumber) cellMargin) sig $& skips
             holdFracs =
                CausalPV.zipWithSimple (\r fr -> r ? (fr, 1))
                $&
                running &|& (CausalP.delay1Zero $& fracs)
         in  (windows, actualSkips, holdFracs))
      ipLeap ipStep periodInt period

dynamicGen ::
   (Ip.C nodesStep, Ip.C nodesLeap,
    A.Additive v, Memory.C v,
    SoV.RationalConstant a, SoV.Fraction a,
    LLVM.CmpRet a, LLVM.CmpResult a ~ Bool,
    Storable a, MakeValueTuple a, ValueTuple a ~ Value a,
    Memory.FirstClass a, Memory.Stored a ~ am, IsSized am,
    LLVM.NumberOfElements a ~ D1) =>
   (Param.T p (Ip.Margin (nodesLeap (nodesStep v))) ->
    (Func.T p (Value a, Value a) (Value Word32),
     Func.T p (Value a, Value a) (Value a)) ->
    (Func.T p (Value a, Value a) (RingBuffer.T v),
     Func.T p (Value a, Value a) (Value Word32),
     Func.T p (Value a, Value a) (Value a))) ->
   (forall r. Ip.T r nodesLeap (Value a) v) ->
   (forall r. Ip.T r nodesStep (Value a) v) ->
   Param.T p Int ->
   Param.T p a ->
   CausalP.T p (Value a, Value a) v
dynamicGen limitMaxShape ipLeap ipStep periodInt period =
   let period32 = Param.word32 periodInt
       cellMargin = combineMarginParams ipLeap ipStep periodInt
       minShape =
          Param.word32 $ fmap fst $
          liftA2 shapeMargin cellMargin periodInt

   in  Func.withArgs $ \(shape, phase) ->
          let (windows, skips, fracs) =
                 limitMaxShape cellMargin $
                 unzip (integrateFrac $& (limitMinShape minShape $& shape))
              (offsets, shapePhases) =
                 unzip
                    (flattenShapePhaseProc period32 period $&
                       (constantFromWord32 minShape + fracs)
                       &|&
                       (CausalP.osciCoreSync $&
                          phase
                          &|&
                          negate
                             (CausalPV.map (flip (/)) period $&
                                (CausalP.mapSimple LLVM.inttofp $& skips))))
          in  interpolateCell ipLeap ipStep $&
                 (CausalP.map (uncurry . cellFromBuffer) period32
                  $&
                  windows
                  &|&
                  offsets)
                 &|&
                 shapePhases

constantFromWord32 ::
   (IsFloating a, LLVM.NumberOfElements a ~ TypeNum.D1) =>
   Param.T p Word32 -> Func.T p inp (Value a)
constantFromWord32 x =
   Func.fromSignal
      (CausalP.mapSimple LLVM.inttofp $* SigP.constant x)

limitMinShape ::
   (Storable a, Memory.FirstClass a, Memory.Stored a ~ am, IsSized am,
    LLVM.NumberOfElements a ~ TypeNum.D1,
    LLVM.IsFloating a, LLVM.CmpRet a, LLVM.CmpResult a ~ Bool) =>
   Param.T p Word32 ->
   CausalP.T p (Value a) (Value a)
limitMinShape xLim =
   CausalPV.mapAccum
      (\_ x lim -> (x%>=lim) ? ((x-lim,zero), (zero,lim-x)))
      (Value.lift1 LLVM.inttofp) (return ()) xLim

integrateFrac ::
   (IsFloating a,
    Memory.FirstClass a, Memory.Stored a ~ am, IsSized am,
    LLVM.NumberOfElements a ~ TypeNum.D1) =>
   CausalP.T p (Value a) (Value Word32, Value a)
integrateFrac =
   CausalP.mapAccumSimple
      (\a (_n,frac) -> do
         s <- splitFraction =<< A.add a frac
         return (s, s))
      (return (A.zero, A.zero))


interpolateCell ::
   (Ip.C nodesStep, Ip.C nodesLeap) =>
   (forall r. Ip.T r nodesLeap a v) ->
   (forall r. Ip.T r nodesStep a v) ->
   CausalP.T p (nodesLeap (nodesStep v), (a, a)) v
interpolateCell ipLeap ipStep =
   CausalP.mapSimple
      (\(nodes, (leap,step)) ->
         ipLeap leap =<< mapM (ipStep step) nodes)

cellFromBuffer ::
   (Memory.C a, Ip.C nodesLeap, Ip.C nodesStep) =>
   Value Word32 ->
   RingBuffer.T a ->
   Value Word32 ->
   CodeGenFunction r (nodesLeap (nodesStep a))
cellFromBuffer periodInt buffer offset =
   Ip.indexNodes
      (Ip.indexNodes (flip RingBuffer.index buffer) A.one)
      periodInt offset

elementMargin ::
   Ip.Margin (nodesLeap (nodesStep v)) ->
   Ip.Margin (nodesLeap (nodesStep (Serial.Element v)))
elementMargin (Ip.Margin x y) = Ip.Margin x y

peekCell ::
   (Storable a, MakeValueTuple a, ValueTuple a ~ value, Memory.C value,
    Ip.C nodesLeap, Ip.C nodesStep) =>
   Param.T p (Ip.Margin (nodesLeap (nodesStep value))) ->
   Param.T p Word32 ->
   Param.T p (SV.Vector a) ->
   CausalP.T p (Value Word32) (nodesLeap (nodesStep value))
peekCell margin period32 vec =
   Param.with (Param.word32 $ fmap Ip.marginOffset margin) $ \getOffset valueOffset ->
   Param.with period32 $ \getPeriod valuePeriod -> CausalPrivP.Cons
      (\(p,off,per) () n () -> MaybeCont.lift $ do
         offset <- A.sub n (valueOffset off)
         nodes <-
            Ip.loadNodes (Ip.loadNodes Memory.load A.one) (valuePeriod per)
               =<< LLVM.getElementPtr p (offset, ())
         return (nodes, ()))
      (return ())
      (return . flip (,) ())
      (const $ const $ return ())
      (\p ->
         let (fp,ptr,_l) = SVU.unsafeToPointers $ Param.get vec p
         in  return (fp, (ptr, getOffset p, getPeriod p)))
      touchForeignPtr


flattenShapePhaseProc ::
   (IsFloating a, SoV.Fraction a, SoV.RationalConstant a,
    Storable ah, MakeValueTuple ah, ValueTuple ah ~ Value a,
    Memory.FirstClass a, Memory.Stored a ~ am, IsSized am,
    LLVM.NumberOfElements a ~ D1) =>
   Param.T p Word32 ->
   Param.T p ah ->
   CausalP.T p
      (Value a, Value a)
      (Value Word32, (Value a, Value a))
flattenShapePhaseProc period32 period =
   CausalP.map
      (\(perInt, per) (shape, phase) ->
         flattenShapePhase perInt per shape phase)
      (liftA2 (,) period32 period)

flattenShapePhaseProcPacked ::
   (IsFloating a, Vector.Real a, SoV.RationalConstant a,
    Storable ah, MakeValueTuple ah, ValueTuple ah ~ Value a,
    Memory.FirstClass a, Memory.Stored a ~ am, IsSized am,
    TypeNum.Positive n) =>
   Param.T p Word32 ->
   Param.T p ah ->
   CausalP.T p
      (Serial.Value n a, Serial.Value n a)
      (Serial.Value n Word32,
       (Serial.Value n a, Serial.Value n a))
flattenShapePhaseProcPacked period32 period =
   CausalP.map
      (\(perInt, per) (Serial.Cons shape, Serial.Cons phase) -> do
         perIntVec <- SoV.replicate perInt
         perVec <- SoV.replicate per
         (i, (leap, step)) <-
            flattenShapePhase perIntVec perVec shape phase
         return (Serial.Cons i, (Serial.Cons leap, Serial.Cons step)))
      (liftA2 (,) period32 period)

flattenShapePhase ::
   (IsFloating a, SoV.Fraction a, SoV.RationalConstant a,
    LLVM.NumberOfElements a ~ LLVM.NumberOfElements i,
    LLVM.IsInteger i) =>
   Value i ->
   Value a ->
   Value a -> Value a ->
   CodeGenFunction r (Value i, (Value a, Value a))
flattenShapePhase = Value.unlift4 $ \periodInt period shape phase ->
   let qLeap = Value.lift1 A.fraction $ shape/period - phase
       (n,qStep) =
          unzip $ Value.lift1 splitFraction $
          {-
          If 'shape' is correctly limited,
          the value is always non-negative algebraically,
          but maybe not numerically.
          -}
          Value.max zero $
          shape - qLeap * Value.lift1 LLVM.inttofp periodInt
   in  (n,(qLeap,qStep))

{- |
You must make sure, that the argument is non-negative.
-}
splitFraction ::
   (IsFloating a, LLVM.IsInteger i,
    LLVM.NumberOfElements a ~ LLVM.NumberOfElements i) =>
   Value a -> CodeGenFunction r (Value i, Value a)
splitFraction x = do
   n <- LLVM.fptoint x
   frac <- A.sub x =<< LLVM.inttofp n
   return (n, frac)


limitShape ::
   (Memory.FirstClass t, Memory.Stored t ~ tm, IsSized tm,
    IsFloating t, SoV.Real t,
    LLVM.NumberOfElements t ~ LLVM.NumberOfElements i,
    MakeValueTuple i, ValueTuple i ~ Value i,
    Ring.C i, LLVM.IsInteger i, SoV.IntegerConstant i,
    Storable i, Memory.FirstClass i, Memory.Stored i ~ im, IsSized im,
    Ip.C nodesStep, Ip.C nodesLeap) =>
   Param.T p (Ip.Margin (nodesLeap (nodesStep value))) ->
   Param.T p i ->
   Param.T p i ->
   CausalP.T p (Value t) (Value t)
limitShape margin periodInt len =
   CausalPV.zipWithSimple (Value.limit . unzip)
   $<
   limitShapeSignal margin periodInt len

limitShapePacked ::
   (Memory.FirstClass t, Memory.Stored t ~ tm,
    IsSized tm, LLVM.IsPrimitive tm,
    LLVM.SizeOf tm ~ tmsize,
    IsFloating t, LLVM.IsPrimitive t, Vector.Real t,
    TypeNum.Positive n,
    Ip.C nodesStep, Ip.C nodesLeap) =>
   Param.T p (Ip.Margin (nodesLeap (nodesStep value))) ->
   Param.T p Word32 ->
   Param.T p Word32 ->
   CausalP.T p (Serial.Value n t) (Serial.Value n t)
limitShapePacked margin periodInt len =
   CausalPV.zipWithSimple
      (\minmax shape ->
         let (minShape,maxShape) = unzip minmax
         in  Value.limit
                (Value.lift1 Serial.upsample minShape,
                 Value.lift1 Serial.upsample maxShape)
                shape)
   $<
   limitShapeSignal margin periodInt len

limitShapeSignal ::
   (Memory.FirstClass t, Memory.Stored t ~ tm, IsSized tm,
    IsFloating t,
    LLVM.NumberOfElements t ~ LLVM.NumberOfElements i,
    MakeValueTuple i, ValueTuple i ~ Value i,
    Ring.C i, LLVM.IsInteger i, SoV.IntegerConstant i,
    Storable i, Memory.FirstClass i, Memory.Stored i ~ im, IsSized im,
    Ip.C nodesStep, Ip.C nodesLeap) =>
   Param.T p (Ip.Margin (nodesLeap (nodesStep value))) ->
   Param.T p i ->
   Param.T p i ->
   SigP.T p (Value t, Value t)
limitShapeSignal margin periodInt len =
   SigP.Cons
      (\minMax () () -> return (minMax, ()))
      (return ())
      (\(minShapeInt, maxShapeInt) -> do
         minShape <- LLVM.inttofp minShapeInt
         maxShape <- LLVM.inttofp maxShapeInt
         return ((minShape, maxShape), ()))
      (const $ const $ return ())
      (\p -> return ((),
         shapeLimits
            (Param.get margin p)
            (Param.get periodInt p)
            (Param.get len p)))
      (const $ return ())


_limitShape ::
   (Ring.C th, Storable th, MakeValueTuple th, ValueTuple th ~ t,
    Memory.C t, A.Real t,
    Ip.C nodesStep, Ip.C nodesLeap) =>
   Ip.Margin (nodesLeap (nodesStep value)) ->
   Param.T p th ->
   Param.T p th ->
   CausalP.T p t t
_limitShape margin periodInt len =
   CausalPrivP.Cons
      (\(minShape,maxShape) () shape () -> MaybeCont.lift $ do
         limited <- A.min maxShape =<< A.max minShape shape
         return (limited, ()))
      (return ())
      (\minmax -> return (minmax, ()))
      (const $ const $ return ())
      (\p ->
         return
            ((),
             shapeLimits margin
                (Param.get periodInt p)
                (Param.get len p)))
      (const $ return ())

shapeLimits ::
   (Ip.C nodesLeap, Ip.C nodesStep, Ring.C t) =>
   Ip.Margin (nodesLeap (nodesStep value)) ->
   t ->
   t ->
   (t, t)
shapeLimits margin periodInt len =
   case shapeMargin margin periodInt of
      (leftMargin, rightMargin) ->
         (leftMargin, len - rightMargin)

_shapeLimits ::
   (Ip.C nodesLeap, Ip.C nodesStep,
    IsFloating t, LLVM.NumberOfElements t ~ D1) =>
   Ip.Margin (nodesLeap (nodesStep value)) ->
   Value.T (Value Word32) ->
   Value.T (Value t) ->
   (Value.T (Value t), Value.T (Value t))
_shapeLimits margin periodInt len =
   let (leftMargin, rightMargin) = shapeMargin margin periodInt
   in  (Value.lift1 LLVM.inttofp leftMargin,
        len - Value.lift1 LLVM.inttofp rightMargin)

shapeMargin ::
   (Ip.C nodesLeap, Ip.C nodesStep, Ring.C i) =>
   Ip.Margin (nodesLeap (nodesStep value)) ->
   i -> (i, i)
shapeMargin margin periodInt =
   let leftMargin = fromIntegral (Ip.marginOffset margin) + periodInt
       rightMargin = fromIntegral (Ip.marginNumber margin) - leftMargin
   in  (leftMargin, rightMargin)

combineMarginParams ::
   (Ip.C nodesStep, Ip.C nodesLeap) =>
   (forall r. Ip.T r nodesLeap a v) ->
   (forall r. Ip.T r nodesStep a v) ->
   Param.T p Int ->
   Param.T p (Ip.Margin (nodesLeap (nodesStep v)))
combineMarginParams ipLeap ipStep periodInt =
   fmap
      (combineMargins (Ip.toMargin ipLeap) (Ip.toMargin ipStep))
      periodInt

combineMargins ::
   Ip.Margin (nodesLeap value) ->
   Ip.Margin (nodesStep value) ->
   Int ->
   Ip.Margin (nodesLeap (nodesStep value))
combineMargins marginLeap marginStep periodInt =
   Ip.Margin {
      Ip.marginNumber =
         Ip.marginNumber marginStep +
         Ip.marginNumber marginLeap * periodInt,
      Ip.marginOffset =
         Ip.marginOffset marginStep +
         Ip.marginOffset marginLeap * periodInt
   }


{- |
@zigZagLong loopStart loopLength@
creates a curve that starts at 0
and is linear until it reaches @loopStart+loopLength@.
Then it begins looping in a ping-pong manner
between @loopStart+loopLength@ and @loopStart@.
It is useful as @shape@ control for looping a sound.
Input of the causal process is the slope (or frequency) control.
Slope values must not be negative.

*Main> Sig.renderChunky SVL.defaultChunkSize (Causal.take 25 <<< Helix.zigZagLong 6 10 $* 2) () :: SVL.Vector Float
VectorLazy.fromChunks [Vector.pack [0.0,1.999999,3.9999995,6.0,8.0,10.0,12.0,14.0,15.999999,14.000001,12.0,10.0,7.999999,6.0,8.0,10.0,12.0,14.0,16.0,14.0,11.999999,9.999998,7.999998,6.0000024,8.000002]]
-}
zigZagLong ::
   (Storable a, MakeValueTuple a, ValueTuple a ~ Value a,
    Memory.FirstClass a, Memory.Stored a ~ am, IsSized am,
    SoV.Fraction a, IsFloating a, SoV.RationalConstant a, LLVM.CmpRet a,
    Field.C a) =>
   Param.T p a ->
   Param.T p a ->
   CausalP.T p (Value a) (Value a)
zigZagLong =
   zigZagLongGen (CausalP.fromSignal . SigP.constant) zigZag

zigZagLongPacked ::
   (Storable a, MakeValueTuple a, ValueTuple a ~ Value a,
    SoV.Fraction a, SoV.RationalConstant a, Vector.Real a,
    Memory.FirstClass a, Memory.Stored a ~ am, IsSized am,
    LLVM.IsPrimitive am,
    Field.C a,
    (n TypeNum.:*: LLVM.SizeOf am) ~ amsize,
    TypeNum.Positive amsize,
    TypeNum.Positive n) =>
   Param.T p a ->
   Param.T p a ->
   CausalP.T p (Serial.Value n a) (Serial.Value n a)
zigZagLongPacked =
   zigZagLongGen (CausalP.fromSignal . SigPS.constant) zigZagPacked

zigZagLongGen ::
   (A.RationalConstant al, A.Field al, Field.C a) =>
   (Param.T p a -> CausalP.T p al al) ->
   (Param.T p a -> CausalP.T p al al) ->
   Param.T p a ->
   Param.T p a ->
   CausalP.T p al al
zigZagLongGen constant zz prefix loop =
   zz (negate $ prefix/loop) * constant loop + constant prefix
   <<<
   id / constant loop

{- |
@zigZag start@ creates a zig-zag curve with values between 0 and 1, inclusively,
that is useful as @shape@ control for looping a sound.
Input of the causal process is the slope (or frequency) control.
Slope values must not be negative.
The start value must be at most 2 and may be negative.
-}
zigZag ::
   (Storable a, MakeValueTuple a, ValueTuple a ~ Value a,
    Memory.FirstClass a, Memory.Stored a ~ am, IsSized am,
    SoV.Fraction a, IsFloating a, SoV.RationalConstant a,
    LLVM.CmpRet a) =>
   Param.T p a ->
   CausalP.T p (Value a) (Value a)
zigZag start =
   CausalPV.mapSimple (\x -> 1-abs (1-x))
   <<<
   CausalPV.mapAccum
      (\_ d t0 ->
         let t1 = t0+d
         in  (t0, wrap (curry . (??)) t1))
      id (return ()) start

zigZagPacked ::
   (Storable a, MakeValueTuple a, ValueTuple a ~ Value a,
    Memory.FirstClass a, Memory.Stored a ~ am, IsSized am,
    SoV.Fraction a, Vector.Real a, IsFloating a, SoV.RationalConstant a,
    LLVM.CmpRet a,
    TypeNum.Positive n) =>
   Param.T p a ->
   CausalP.T p (Serial.Value n a) (Serial.Value n a)
zigZagPacked start =
   Serial.Cons
   ^<<
   CausalPV.mapSimple (\x -> 1 - abs (1-x))
   <<<
   CausalPV.mapAccum
      (\_ d t0 ->
         let (t1, cum) = unzip $ Value.lift2 Vector.cumulate t0 d
         {-
         Vector.select can be replaced by (??)
         once vector select is implemented by LLVM.
         -}
         in  (wrap (Value.lift3 Vector.select) cum, t1))
      id (return ()) start
   <<^
   (\(Serial.Cons v) -> v)

wrap ::
   (SoV.RationalConstant a, IsFloating a, SoV.Fraction a, LLVM.CmpRet a) =>
   (Value.T (Value (LLVM.CmpResult a)) ->
    Value.T (Value a) ->
    Value.T (Value a) ->
    Value.T (Value a)) ->
   Value.T (Value a) -> Value.T (Value a)
wrap select a = select (a%>0) (2 * Value.fraction (a/2)) a