{-# LANGUAGE NoImplicitPrelude #-}
{-# LANGUAGE MultiParamTypeClasses #-}
{-# LANGUAGE FlexibleInstances #-}
{- |
Copyright   :  (c) Henning Thielemann 2008
License     :  GPL

Maintainer  :  synthesizer@henning-thielemann.de
Stability   :  provisional
Portability :  requires multi-parameter type classes

Butterworth lowpass and highpass
-}
module Synthesizer.Plain.Filter.Recursive.Butterworth (
   Parameter,
   causal,
   causalPole,
   highpassCausalPole, highpassPole,
   lowpassCausalPole,  lowpassPole,
   modifier,
   parameter,
   partialParameter,
   runPole,

   -- used in Dimensional.Causal.FilterParameter
   checkedHalf,
   -- used in LLVM.Filter.Butterworth
   partialRatio,
   makeSines,
   ) where

import Synthesizer.Plain.Filter.Recursive (Passband(Lowpass,Highpass), Pole(Pole))
import qualified Synthesizer.Plain.Filter.Recursive.SecondOrderCascade as Cascade
import qualified Synthesizer.Plain.Filter.Recursive.SecondOrder as Filt2
import qualified Synthesizer.Plain.Signal   as Sig
import qualified Synthesizer.Plain.Modifier as Modifier
import qualified Synthesizer.Causal.Process as Causal
import Control.Arrow ((>>>), )

import qualified Algebra.Module                as Module
import qualified Algebra.Transcendental        as Trans
import qualified Algebra.Ring                  as Ring

import qualified Data.StorableVector as SV
import Foreign.Storable (Storable)

import NumericPrelude.Numeric
import NumericPrelude.Base



sineList, _sineListSlow, sineListFast :: (Trans.C a) => a -> [a]
sineList :: forall a. C a => a -> [a]
sineList = forall a. C a => a -> [a]
sineListFast

_sineListSlow :: forall a. C a => a -> [a]
_sineListSlow a
x =
   forall a b. (a -> b) -> [a] -> [b]
map forall a. C a => a -> a
sin forall a b. (a -> b) -> a -> b
$ forall a b. (a -> b) -> [a] -> [b]
map (a
xforall a. C a => a -> a -> a
*) forall a b. (a -> b) -> a -> b
$ forall a. (a -> a) -> a -> [a]
iterate (a
2forall a. C a => a -> a -> a
+) a
1

sineListFast :: forall a. C a => a -> [a]
sineListFast a
x =
   let sinx :: a
sinx  = forall a. C a => a -> a
sin a
x
       cos2x :: a
cos2x = a
2 forall a. C a => a -> a -> a
- a
4forall a. C a => a -> a -> a
*a
sinxforall a. C a => a -> a -> a
*a
sinx
       --cos2x = 2 * cos (2*x)
       sines :: [a]
sines = (-a
sinx) forall a. a -> [a] -> [a]
: a
sinx forall a. a -> [a] -> [a]
:
                  forall a b c. (a -> b -> c) -> [a] -> [b] -> [c]
zipWith (\a
y1 a
y0 -> a
cos2x forall a. C a => a -> a -> a
* a
y0 forall a. C a => a -> a -> a
- a
y1) [a]
sines (forall a. [a] -> [a]
tail [a]
sines)
   in  forall a. [a] -> [a]
tail [a]
sines

makeSines :: (Trans.C a) =>
   Int -> [a]
makeSines :: forall a. C a => Int -> [a]
makeSines Int
order =
   forall a. Int -> [a] -> [a]
take (String -> Int -> Int
checkedHalf String
"makeSines" Int
order) (forall a. C a => a -> [a]
sineList (forall a. C a => a
pi forall a. C a => a -> a -> a
/ forall a b. (C a, C b) => a -> b
fromIntegral (Int
2forall a. C a => a -> a -> a
*Int
order)))

partialRatio :: (Trans.C a) =>
   Int -> a -> a
partialRatio :: forall a. C a => Int -> a -> a
partialRatio Int
order a
ratio =
   (a
1forall a. C a => a -> a -> a
/a
ratioforall a. C a => a -> Integer -> a
^Integer
2forall a. C a => a -> a -> a
-a
1) forall a. C a => a -> a -> a
** (-a
1 forall a. C a => a -> a -> a
/ forall a b. (C a, C b) => a -> b
fromIntegral (Int
2forall a. C a => a -> a -> a
*Int
order))



_partialLowpassParameterInstable, partialLowpassParameter :: (Trans.C a) =>
   a -> a -> a -> Filt2.Parameter a

{- must handle infinite values when 'freq' approaches 0.5 -}
_partialLowpassParameterInstable :: forall a. C a => a -> a -> a -> Parameter a
_partialLowpassParameterInstable a
ratio a
freq a
sinw =
   let wc :: a
wc    = a
ratio forall a. C a => a -> a -> a
* forall a. C a => a -> a
tan (forall a. C a => a
piforall a. C a => a -> a -> a
*a
freq)
       sinw2 :: a
sinw2 = a
2 forall a. C a => a -> a -> a
* a
wc forall a. C a => a -> a -> a
* a
sinw
       wc2 :: a
wc2   = a
wc forall a. C a => a -> a -> a
* a
wc
       denom :: a
denom = a
wc2 forall a. C a => a -> a -> a
+ a
sinw2 forall a. C a => a -> a -> a
+ a
1
       c0 :: a
c0    = a
wc2 forall a. C a => a -> a -> a
/ a
denom
   in  forall a. a -> a -> a -> a -> a -> Parameter a
Filt2.Parameter a
c0 (a
2forall a. C a => a -> a -> a
*a
c0) a
c0
          (a
2forall a. C a => a -> a -> a
*(a
1forall a. C a => a -> a -> a
-a
wc2)forall a. C a => a -> a -> a
/a
denom) ((-a
wc2forall a. C a => a -> a -> a
+a
sinw2forall a. C a => a -> a -> a
-a
1)forall a. C a => a -> a -> a
/a
denom)

-- using ratio disallows simplification by trigonometric Pythagoras' theorem
partialLowpassParameter :: forall a. C a => a -> a -> a -> Parameter a
partialLowpassParameter a
ratio a
freq =
   let phi :: a
phi      = forall a. C a => a
piforall a. C a => a -> a -> a
*a
freq
       rsin2phi :: a
rsin2phi = a
ratio forall a. C a => a -> a -> a
* forall a. C a => a -> a
sin (a
2forall a. C a => a -> a -> a
*a
phi)
       cosphi2 :: a
cosphi2  = forall a. C a => a -> a
cos a
phi forall a. C a => a -> Integer -> a
^ Integer
2
       c0d :: a
c0d      = (a
ratio forall a. C a => a -> a -> a
* forall a. C a => a -> a
sin a
phi) forall a. C a => a -> Integer -> a
^ Integer
2
       d1d :: a
d1d      = (a
cosphi2 forall a. C a => a -> a -> a
- a
c0d) forall a. C a => a -> a -> a
* a
2
       ratsin :: a
ratsin   = a
cosphi2 forall a. C a => a -> a -> a
+ a
c0d
   in  \a
sinw ->
          let rsinsin :: a
rsinsin = a
rsin2phi forall a. C a => a -> a -> a
* a
sinw
              denom :: a
denom   = a
rsinsin forall a. C a => a -> a -> a
+ a
ratsin
              d2d :: a
d2d     = a
rsinsin forall a. C a => a -> a -> a
- a
ratsin
              c0 :: a
c0      = a
c0d forall a. C a => a -> a -> a
/ a
denom
              d1 :: a
d1      = a
d1d forall a. C a => a -> a -> a
/ a
denom
              d2 :: a
d2      = a
d2d forall a. C a => a -> a -> a
/ a
denom
          in  forall a. a -> a -> a -> a -> a -> Parameter a
Filt2.Parameter a
c0 (a
2forall a. C a => a -> a -> a
*a
c0) a
c0 a
d1 a
d2


-- * use second order filter parameters for control

type Parameter a = Cascade.Parameter a

{-# INLINE parameter #-}
parameter ::
   (Trans.C a, Storable a) =>
   Passband -> Int -> Pole a -> Parameter a
parameter :: forall a.
(C a, Storable a) =>
Passband -> Int -> Pole a -> Parameter a
parameter Passband
kind Int
order =
   -- I hope that the 'let' is floated out of a 'map'
   let sinesVec :: Vector a
sinesVec = forall a. Storable a => [a] -> Vector a
SV.pack (forall a. C a => Int -> [a]
makeSines Int
order)
   in  \ (Pole a
ratio a
freq) ->
           forall a. Vector (Parameter a) -> Parameter a
Cascade.Parameter forall a b. (a -> b) -> a -> b
$
           forall a b.
(Storable a, Storable b) =>
(a -> b) -> Vector a -> Vector b
SV.map
              (\a
sinw ->
                 forall a. C a => Passband -> a -> a -> a -> Parameter a
partialParameter Passband
kind (forall a. C a => Int -> a -> a
partialRatio Int
order a
ratio) a
sinw a
freq) forall a b. (a -> b) -> a -> b
$
           Vector a
sinesVec

partialParameter ::
   Trans.C a =>
   Passband -> a -> a -> a -> Filt2.Parameter a
partialParameter :: forall a. C a => Passband -> a -> a -> a -> Parameter a
partialParameter Passband
kind a
partRatio a
sinw a
freq =
   forall a. C a => Passband -> (a -> Parameter a) -> a -> Parameter a
Filt2.adjustPassband Passband
kind
      (forall a b c. (a -> b -> c) -> b -> a -> c
flip (forall a. C a => a -> a -> a -> Parameter a
partialLowpassParameter a
partRatio) a
sinw)
      a
freq

{-# INLINE modifier #-}
modifier ::
   (Ring.C a, Module.C a v, Storable a, Storable v) =>
   Int ->
   Modifier.Simple (Cascade.State v) (Parameter a) v v
modifier :: forall a v.
(C a, C a v, Storable a, Storable v) =>
Int -> Simple (State v) (Parameter a) v v
modifier =
   forall a v.
(C a, C a v, Storable a, Storable v) =>
Int -> Simple (State v) (Parameter a) v v
Cascade.modifier

{-# INLINE causal #-}
causal :: (Ring.C a, Module.C a v, Storable a, Storable v) =>
   Int ->
   Causal.T (Parameter a, v) v
causal :: forall a v.
(C a, C a v, Storable a, Storable v) =>
Int -> T (Parameter a, v) v
causal Int
order =
   forall a v.
(C a, C a v, Storable a, Storable v) =>
Int -> T (Parameter a, v) v
Cascade.causal (String -> Int -> Int
checkedHalf String
"causal" Int
order)


{-# INLINE checkedHalf #-}
checkedHalf :: String -> Int -> Int
checkedHalf :: String -> Int -> Int
checkedHalf String
funcName Int
order =
   let (Int
order2,Int
r) = forall a. C a => a -> a -> (a, a)
divMod Int
order Int
2
   in  if Int
rforall a. Eq a => a -> a -> Bool
==Int
0
         then Int
order2
         else forall a. HasCallStack => String -> a
error forall a b. (a -> b) -> a -> b
$ String
"Butterworth."forall a. [a] -> [a] -> [a]
++String
funcNameforall a. [a] -> [a] -> [a]
++String
": order must be even"

{-
lowpassCausal, highpassCausal :: (Trans.C a, Module.C a v) =>
   Int -> Causal.T (Parameter a, v) v
lowpassCausal  = causal Lowpass
highpassCausal = causal Highpass

lowpass, highpass :: (Trans.C a, Module.C a v) =>
   Int -> Sig.T (Parameter a) -> Sig.T v -> Sig.T v
lowpass  = run Lowpass
highpass = run Highpass
-}


-- * directly use frequency as control parameter

{- |
When called as @runPole kind order ratio freqs@,
the filter amplifies frequency 0 with factor 1
and frequency @freq@ with factor @ratio@.

It uses the frequency and ratio information directly
and thus cannot benefit from efficient parameter interpolation
(asynchronous run of a ControlledProcess).
-}
runPole :: (Trans.C a, Module.C a v) =>
   Passband -> Int -> Sig.T a -> Sig.T a -> Sig.T v -> Sig.T v
runPole :: forall a v.
(C a, C a v) =>
Passband -> Int -> T a -> T a -> T v -> T v
runPole Passband
kind Int
order T a
ratios T a
freqs =
   let makePartialFilter :: a -> T v -> T v
makePartialFilter a
s =
          forall a v. (C a, C a v) => T (Parameter a) -> T v -> T v
Filt2.run forall a b. (a -> b) -> a -> b
$
          forall a b c. (a -> b -> c) -> [a] -> [b] -> [c]
zipWith
             (\a
ratio a
freq ->
                forall a. C a => Passband -> a -> a -> a -> Parameter a
partialParameter Passband
kind (forall a. C a => Int -> a -> a
partialRatio Int
order a
ratio) a
s a
freq)
             T a
ratios T a
freqs
   in  forall (t :: * -> *) b a.
Foldable t =>
(b -> a -> b) -> b -> t a -> b
foldl forall b c a. (b -> c) -> (a -> b) -> a -> c
(.) forall a. a -> a
id (forall a b. (a -> b) -> [a] -> [b]
map forall {v}. C a v => a -> T v -> T v
makePartialFilter (forall a. C a => Int -> [a]
makeSines Int
order))

causalPole :: (Trans.C a, Module.C a v) =>
   Passband -> Int -> Causal.T (Pole a, v) v
causalPole :: forall a v. (C a, C a v) => Passband -> Int -> T (Pole a, v) v
causalPole Passband
kind Int
order =
   let {-# INLINE makePartialFilter #-}
       makePartialFilter :: a -> T (Pole a, c) c
makePartialFilter a
s =
          forall (a :: * -> * -> *) b c d.
Arrow a =>
a b c -> a (b, d) (c, d)
Causal.first
             (forall a b. (a -> b) -> T a b
Causal.map (\(Pole a
ratio a
freq) ->
                forall a. C a => Passband -> a -> a -> a -> Parameter a
partialParameter Passband
kind (forall a. C a => Int -> a -> a
partialRatio Int
order a
ratio) a
s a
freq))
          forall {k} (cat :: k -> k -> *) (a :: k) (b :: k) (c :: k).
Category cat =>
cat a b -> cat b c -> cat a c
>>>
          forall a v. (C a, C a v) => T (Parameter a, v) v
Filt2.causal
   in  forall c x. [T (c, x) x] -> T (c, x) x
Causal.chainControlled forall a b. (a -> b) -> a -> b
$ forall a b. (a -> b) -> [a] -> [b]
map forall {a} {c}. (C a, C a c) => a -> T (Pole a, c) c
makePartialFilter forall a b. (a -> b) -> a -> b
$ forall a. C a => Int -> [a]
makeSines Int
order


lowpassCausalPole, highpassCausalPole :: (Trans.C a, Module.C a v) =>
   Int -> Causal.T (Pole a, v) v
lowpassCausalPole :: forall a v. (C a, C a v) => Int -> T (Pole a, v) v
lowpassCausalPole  = forall a v. (C a, C a v) => Passband -> Int -> T (Pole a, v) v
causalPole Passband
Lowpass
highpassCausalPole :: forall a v. (C a, C a v) => Int -> T (Pole a, v) v
highpassCausalPole = forall a v. (C a, C a v) => Passband -> Int -> T (Pole a, v) v
causalPole Passband
Highpass

lowpassPole, highpassPole :: (Trans.C a, Module.C a v) =>
   Int -> Sig.T a -> Sig.T a -> Sig.T v -> Sig.T v
lowpassPole :: forall a v. (C a, C a v) => Int -> T a -> T a -> T v -> T v
lowpassPole  = forall a v.
(C a, C a v) =>
Passband -> Int -> T a -> T a -> T v -> T v
runPole Passband
Lowpass
highpassPole :: forall a v. (C a, C a v) => Int -> T a -> T a -> T v -> T v
highpassPole = forall a v.
(C a, C a v) =>
Passband -> Int -> T a -> T a -> T v -> T v
runPole Passband
Highpass