{-# LANGUAGE NoImplicitPrelude #-}
module Synthesizer.Storable.Repair (removeClicks, ) where

import qualified Data.StorableVector.Lazy.Pattern as SVP
import qualified Data.StorableVector.Lazy as SVL
import qualified Data.StorableVector as SV

import Control.Monad (guard, )

import qualified Data.List as List
import Data.Tuple.HT (mapFst, mapSnd, )
import Data.Maybe.HT (toMaybe, )
import Data.Word (Word8, )

import Foreign.Storable (Storable, )

import qualified Algebra.RealField as RealField
import qualified Algebra.Field     as Field
import qualified Algebra.Additive  as Additive

import NumericPrelude.Numeric
import NumericPrelude.Base


{-
could use Generic.Control.line
-}
ramp ::
   (Storable a, Field.C a) =>
   Int -> (a, a) -> SVL.Vector a
ramp :: forall a. (Storable a, C a) => Int -> (a, a) -> Vector a
ramp Int
n (a
y0, a
y1) =
   ChunkSize -> Int -> (Int -> a) -> Vector a
forall a. Storable a => ChunkSize -> Int -> (Int -> a) -> Vector a
SVL.sampleN ChunkSize
SVL.defaultChunkSize Int
n ((Int -> a) -> Vector a) -> (Int -> a) -> Vector a
forall a b. (a -> b) -> a -> b
$ \Int
k ->
      (a
y0 a -> a -> a
forall a. C a => a -> a -> a
* Int -> a
forall a b. (C a, C b) => a -> b
fromIntegral (Int
nInt -> Int -> Int
forall a. C a => a -> a -> a
-Int
k) a -> a -> a
forall a. C a => a -> a -> a
+ a
y1 a -> a -> a
forall a. C a => a -> a -> a
* Int -> a
forall a b. (C a, C b) => a -> b
fromIntegral Int
k)
         a -> a -> a
forall a. C a => a -> a -> a
/ Int -> a
forall a b. (C a, C b) => a -> b
fromIntegral Int
n

svHead :: (Storable a) => SVL.Vector a -> a
svHead :: forall a. Storable a => Vector a -> a
svHead =
   a -> (a -> Vector a -> a) -> Vector a -> a
forall a b.
Storable a =>
b -> (a -> Vector a -> b) -> Vector a -> b
SVL.switchL ([Char] -> a
forall a. HasCallStack => [Char] -> a
error [Char]
"head: empty storable vector") a -> Vector a -> a
forall a b. a -> b -> a
const

{- |
@less-than@ comparison where 'Nothing' means \"unbounded\".
-}
maybeLT :: Ord a => Maybe a -> a -> Bool
maybeLT :: forall a. Ord a => Maybe a -> a -> Bool
maybeLT Maybe a
mx a
y =
   Bool -> (a -> Bool) -> Maybe a -> Bool
forall b a. b -> (a -> b) -> Maybe a -> b
maybe Bool
False (a -> a -> Bool
forall a. Ord a => a -> a -> Bool
<a
y) Maybe a
mx


type Jump = Word8

up, down, smooth :: Jump
up :: Jump
up = Jump
1
down :: Jump
down = -Jump
1
smooth :: Jump
smooth = Jump
0

{- |
both @up@ and @down@ threshold must be positive numbers.
-}
splitAtJump ::
   (Storable a, Additive.C a, Ord a) =>
   Int ->
   Maybe a -> Maybe a ->
   SVL.Vector a -> (SVL.Vector a, (Jump, SVL.Vector a))
splitAtJump :: forall a.
(Storable a, C a, Ord a) =>
Int
-> Maybe a -> Maybe a -> Vector a -> (Vector a, (Jump, Vector a))
splitAtJump Int
slopeDistance Maybe a
thresholdUp Maybe a
thresholdDown Vector a
xs =
   (\(Vector Jump
beforeJump, Vector Jump
afterJump) ->
      (Vector a -> (Jump, Vector a))
-> (Vector a, Vector a) -> (Vector a, (Jump, Vector a))
forall b c a. (b -> c) -> (a, b) -> (a, c)
mapSnd ((,) (Jump -> (Jump -> Vector Jump -> Jump) -> Vector Jump -> Jump
forall a b.
Storable a =>
b -> (a -> Vector a -> b) -> Vector a -> b
SVL.switchL Jump
up (\Jump
d Vector Jump
_ -> Jump
d) Vector Jump
afterJump)) ((Vector a, Vector a) -> (Vector a, (Jump, Vector a)))
-> (Vector a, Vector a) -> (Vector a, (Jump, Vector a))
forall a b. (a -> b) -> a -> b
$
      LazySize -> Vector a -> (Vector a, Vector a)
forall a.
Storable a =>
LazySize -> Vector a -> (Vector a, Vector a)
SVP.splitAt (Vector Jump -> LazySize
forall a. Vector a -> LazySize
SVP.length Vector Jump
beforeJump) Vector a
xs) ((Vector Jump, Vector Jump) -> (Vector a, (Jump, Vector a)))
-> (Vector Jump, Vector Jump) -> (Vector a, (Jump, Vector a))
forall a b. (a -> b) -> a -> b
$
   (Jump -> Bool) -> Vector Jump -> (Vector Jump, Vector Jump)
forall a.
Storable a =>
(a -> Bool) -> Vector a -> (Vector a, Vector a)
SVL.span (Jump
smoothJump -> Jump -> Bool
forall a. Eq a => a -> a -> Bool
==) (Vector Jump -> (Vector Jump, Vector Jump))
-> Vector Jump -> (Vector Jump, Vector Jump)
forall a b. (a -> b) -> a -> b
$
   (a -> a -> Jump) -> Vector a -> Vector a -> Vector Jump
forall a b c.
(Storable a, Storable b, Storable c) =>
(a -> b -> c) -> Vector a -> Vector b -> Vector c
SVL.zipWith
      (\a
x0 a
x1 ->
         case (Maybe a -> a -> Bool
forall a. Ord a => Maybe a -> a -> Bool
maybeLT Maybe a
thresholdUp (a
x1a -> a -> a
forall a. C a => a -> a -> a
-a
x0),
               Maybe a -> a -> Bool
forall a. Ord a => Maybe a -> a -> Bool
maybeLT Maybe a
thresholdDown (a
x0a -> a -> a
forall a. C a => a -> a -> a
-a
x1)) of
            (Bool
True, Bool
_) -> Jump
up
            (Bool
_, Bool
True) -> Jump
down
            (Bool, Bool)
_ -> Jump
smooth) Vector a
xs (Vector a -> Vector Jump) -> Vector a -> Vector Jump
forall a b. (a -> b) -> a -> b
$
   Vector a -> (Vector a -> a -> Vector a) -> Vector a -> Vector a
forall a b.
Storable a =>
b -> (Vector a -> a -> b) -> Vector a -> b
SVL.switchR Vector a
forall a. Storable a => Vector a
SVL.empty
      (\Vector a
ys a
y ->
         Vector a -> Vector a -> Vector a
forall a. Storable a => Vector a -> Vector a -> Vector a
SVL.append Vector a
ys (Vector a -> Vector a) -> Vector a -> Vector a
forall a b. (a -> b) -> a -> b
$
         [Vector a] -> Vector a
forall a. Storable a => [Vector a] -> Vector a
SVL.fromChunks [Int -> a -> Vector a
forall a. Storable a => Int -> a -> Vector a
SV.replicate Int
slopeDistance a
y]) (Vector a -> Vector a) -> Vector a -> Vector a
forall a b. (a -> b) -> a -> b
$
   Int -> Vector a -> Vector a
forall a. Storable a => Int -> Vector a -> Vector a
SVL.drop Int
slopeDistance Vector a
xs

chopAtJumps ::
   (Storable a, Additive.C a, Ord a) =>
   Int -> Int ->
   Maybe a -> Maybe a ->
   SVL.Vector a -> (SVL.Vector a, [(Jump, SVP.Vector a)])
chopAtJumps :: forall a.
(Storable a, C a, Ord a) =>
Int
-> Int
-> Maybe a
-> Maybe a
-> Vector a
-> (Vector a, [(Jump, Vector a)])
chopAtJumps
      Int
slopeDistance Int
minSpikeDistance
      Maybe a
thresholdUp Maybe a
thresholdDown =
   ((Jump, Vector a) -> [(Jump, Vector a)])
-> (Vector a, (Jump, Vector a)) -> (Vector a, [(Jump, Vector a)])
forall b c a. (b -> c) -> (a, b) -> (a, c)
mapSnd
      (((Jump, Vector a) -> Maybe ((Jump, Vector a), (Jump, Vector a)))
-> (Jump, Vector a) -> [(Jump, Vector a)]
forall b a. (b -> Maybe (a, b)) -> b -> [a]
List.unfoldr
          (\(Jump
dir,Vector a
rest) ->
             Bool
-> ((Jump, Vector a), (Jump, Vector a))
-> Maybe ((Jump, Vector a), (Jump, Vector a))
forall a. Bool -> a -> Maybe a
toMaybe (Bool -> Bool
not (Bool -> Bool) -> Bool -> Bool
forall a b. (a -> b) -> a -> b
$ Vector a -> Bool
forall a. Storable a => Vector a -> Bool
SVL.null Vector a
rest)
                (let (Vector a
ys,Vector a
zs) = Int -> Vector a -> (Vector a, Vector a)
forall a. Storable a => Int -> Vector a -> (Vector a, Vector a)
SVL.splitAt Int
minSpikeDistance Vector a
rest
                 in  (Vector a -> (Jump, Vector a))
-> (Vector a, (Jump, Vector a))
-> ((Jump, Vector a), (Jump, Vector a))
forall a c b. (a -> c) -> (a, b) -> (c, b)
mapFst ((,) Jump
dir (Vector a -> (Jump, Vector a))
-> (Vector a -> Vector a) -> Vector a -> (Jump, Vector a)
forall b c a. (b -> c) -> (a -> b) -> a -> c
. Vector a -> Vector a -> Vector a
forall a. Storable a => Vector a -> Vector a -> Vector a
SVL.append Vector a
ys) ((Vector a, (Jump, Vector a))
 -> ((Jump, Vector a), (Jump, Vector a)))
-> (Vector a, (Jump, Vector a))
-> ((Jump, Vector a), (Jump, Vector a))
forall a b. (a -> b) -> a -> b
$
                     Int
-> Maybe a -> Maybe a -> Vector a -> (Vector a, (Jump, Vector a))
forall a.
(Storable a, C a, Ord a) =>
Int
-> Maybe a -> Maybe a -> Vector a -> (Vector a, (Jump, Vector a))
splitAtJump Int
slopeDistance Maybe a
thresholdUp Maybe a
thresholdDown Vector a
zs)))
   ((Vector a, (Jump, Vector a)) -> (Vector a, [(Jump, Vector a)]))
-> (Vector a -> (Vector a, (Jump, Vector a)))
-> Vector a
-> (Vector a, [(Jump, Vector a)])
forall b c a. (b -> c) -> (a -> b) -> a -> c
.
   Int
-> Maybe a -> Maybe a -> Vector a -> (Vector a, (Jump, Vector a))
forall a.
(Storable a, C a, Ord a) =>
Int
-> Maybe a -> Maybe a -> Vector a -> (Vector a, (Jump, Vector a))
splitAtJump Int
slopeDistance Maybe a
thresholdUp Maybe a
thresholdDown


{- |
Usage:
   @removeClicks
      slopeDistance maxSpikeWidth minSpikeDistance
      thresholdUp thresholdDown@

@slopeDistance@ is the distance of samples in which we analyse differences.
The smoother the spike slope the larger @slopeDistance@ must be.

@slopeDistance@ should be smaller than the _minimal_ spike width.
@maxSpikeWidth@ should be smaller than the minimal spike distance.
Spike distance is measured from one spike beginning to the next one.

@thresholdUp@ is the minimal difference of two samples at @slopeDistance@
that are to be considered an upward jump.
@thresholdDown@ is for downward jumps.
If a threshold is 'Nothing' then jumps in this direction are ignored.
You should only use this if you are very sure
that spikes with the according sign do not occur.
Otherwise the algorithm will be confused
by the jump in reverse direction at the end of the spike.

Example: @removeClicks 1 5 20 (Just 0.1) (Just 0.1)@.

The algorithm works as follows:
Chop the signal at jumps.
Then begin at a certain distance behind the jump
and search backwards for the matching jump at the end of the spike.
If the borders of a spike are found this way,
then they are connected by a linear ramp.
-}
removeClicks ::
   (Storable a, RealField.C a) =>
   Int -> Int -> Int ->
   Maybe a -> Maybe a ->
   SVL.Vector a -> SVL.Vector a
removeClicks :: forall a.
(Storable a, C a) =>
Int -> Int -> Int -> Maybe a -> Maybe a -> Vector a -> Vector a
removeClicks
      Int
slopeDistance Int
maxSpikeWidth Int
minSpikeDistance
      Maybe a
thresholdUp Maybe a
thresholdDown =
   [Vector a] -> Vector a
forall a. Storable a => [Vector a] -> Vector a
SVL.concat
   ([Vector a] -> Vector a)
-> (Vector a -> [Vector a]) -> Vector a -> Vector a
forall b c a. (b -> c) -> (a -> b) -> a -> c
.
   (Vector a -> [Vector a] -> [Vector a])
-> (Vector a, [Vector a]) -> [Vector a]
forall a b c. (a -> b -> c) -> (a, b) -> c
uncurry (:)
   ((Vector a, [Vector a]) -> [Vector a])
-> (Vector a -> (Vector a, [Vector a])) -> Vector a -> [Vector a]
forall b c a. (b -> c) -> (a -> b) -> a -> c
.
   ([(Jump, Vector a)] -> [Vector a])
-> (Vector a, [(Jump, Vector a)]) -> (Vector a, [Vector a])
forall b c a. (b -> c) -> (a, b) -> (a, c)
mapSnd
      (((Jump, Vector a) -> Vector a) -> [(Jump, Vector a)] -> [Vector a]
forall a b. (a -> b) -> [a] -> [b]
map
          (\(Jump
dir, Vector a
chunk) ->
             (Vector a -> Vector a -> Vector a)
-> (Vector a, Vector a) -> Vector a
forall a b c. (a -> b -> c) -> (a, b) -> c
uncurry Vector a -> Vector a -> Vector a
forall a. Storable a => Vector a -> Vector a -> Vector a
SVL.append
             ((Vector a, Vector a) -> Vector a)
-> (Vector a -> (Vector a, Vector a)) -> Vector a -> Vector a
forall b c a. (b -> c) -> (a -> b) -> a -> c
.
             (Vector a -> Vector a)
-> (Vector a, Vector a) -> (Vector a, Vector a)
forall a c b. (a -> c) -> (a, b) -> (c, b)
mapFst
                ((\(Vector a
ys, ~(Jump
_dir,Vector a
click)) ->
                    Vector a -> Vector a -> Vector a
forall a. Storable a => Vector a -> Vector a -> Vector a
SVL.append
                       (Int -> (a, a) -> Vector a
forall a. (Storable a, C a) => Int -> (a, a) -> Vector a
ramp (Vector a -> Int
forall a. Vector a -> Int
SVL.length Vector a
click)
                           (Vector a -> a
forall a. Storable a => Vector a -> a
svHead Vector a
chunk, Vector a -> a
forall a. Storable a => Vector a -> a
svHead Vector a
click))
                       (Vector a -> Vector a
forall a. Storable a => Vector a -> Vector a
SVL.reverse Vector a
ys))
                 ((Vector a, (Jump, Vector a)) -> Vector a)
-> (Vector a -> (Vector a, (Jump, Vector a)))
-> Vector a
-> Vector a
forall b c a. (b -> c) -> (a -> b) -> a -> c
.
                 Int
-> Maybe a -> Maybe a -> Vector a -> (Vector a, (Jump, Vector a))
forall a.
(Storable a, C a, Ord a) =>
Int
-> Maybe a -> Maybe a -> Vector a -> (Vector a, (Jump, Vector a))
splitAtJump Int
slopeDistance
                    (Bool -> Maybe ()
forall (f :: * -> *). Alternative f => Bool -> f ()
guard (Jump
dirJump -> Jump -> Bool
forall a. Eq a => a -> a -> Bool
==Jump
up)   Maybe () -> Maybe a -> Maybe a
forall a b. Maybe a -> Maybe b -> Maybe b
forall (m :: * -> *) a b. Monad m => m a -> m b -> m b
>> Maybe a
thresholdUp)
                    (Bool -> Maybe ()
forall (f :: * -> *). Alternative f => Bool -> f ()
guard (Jump
dirJump -> Jump -> Bool
forall a. Eq a => a -> a -> Bool
==Jump
down) Maybe () -> Maybe a -> Maybe a
forall a b. Maybe a -> Maybe b -> Maybe b
forall (m :: * -> *) a b. Monad m => m a -> m b -> m b
>> Maybe a
thresholdDown)
                 (Vector a -> (Vector a, (Jump, Vector a)))
-> (Vector a -> Vector a)
-> Vector a
-> (Vector a, (Jump, Vector a))
forall b c a. (b -> c) -> (a -> b) -> a -> c
.
                 Vector a -> Vector a
forall a. Storable a => Vector a -> Vector a
SVL.reverse)
             ((Vector a, Vector a) -> (Vector a, Vector a))
-> (Vector a -> (Vector a, Vector a))
-> Vector a
-> (Vector a, Vector a)
forall b c a. (b -> c) -> (a -> b) -> a -> c
.
             Int -> Vector a -> (Vector a, Vector a)
forall a. Storable a => Int -> Vector a -> (Vector a, Vector a)
SVL.splitAt Int
maxSpikeWidth
             (Vector a -> Vector a) -> Vector a -> Vector a
forall a b. (a -> b) -> a -> b
$
             Vector a
chunk))
   ((Vector a, [(Jump, Vector a)]) -> (Vector a, [Vector a]))
-> (Vector a -> (Vector a, [(Jump, Vector a)]))
-> Vector a
-> (Vector a, [Vector a])
forall b c a. (b -> c) -> (a -> b) -> a -> c
.
   Int
-> Int
-> Maybe a
-> Maybe a
-> Vector a
-> (Vector a, [(Jump, Vector a)])
forall a.
(Storable a, C a, Ord a) =>
Int
-> Int
-> Maybe a
-> Maybe a
-> Vector a
-> (Vector a, [(Jump, Vector a)])
chopAtJumps
      Int
slopeDistance Int
minSpikeDistance
      Maybe a
thresholdUp Maybe a
thresholdDown