-- | Filter coefficient calculations.
module Sound.Sc3.Common.Math.Filter where

-- | SOS filter coefficients, (a0,a1,a2,b1,b2)
type SosCoef t = (t,t,t,t,t)

-- | Butterworth low pass or high pass SOS filter coefficients, (a0,a1,a2,b1,b2).
bw_lpf_or_hpf_coef :: Floating n => Bool -> n -> n -> SosCoef n
bw_lpf_or_hpf_coef :: forall n. Floating n => Bool -> n -> n -> SosCoef n
bw_lpf_or_hpf_coef Bool
is_hpf n
sample_rate n
f =
    let f' :: n
f' = n
f forall a. Num a => a -> a -> a
* forall a. Floating a => a
pi forall a. Fractional a => a -> a -> a
/ n
sample_rate
        c :: n
c = if Bool
is_hpf then forall a. Floating a => a -> a
tan n
f' else n
1.0 forall a. Fractional a => a -> a -> a
/ forall a. Floating a => a -> a
tan n
f'
        c2 :: n
c2 = n
c forall a. Num a => a -> a -> a
* n
c
        s2c :: n
s2c = forall a. Floating a => a -> a
sqrt n
2.0 forall a. Num a => a -> a -> a
* n
c
        a0 :: n
a0 = n
1.0 forall a. Fractional a => a -> a -> a
/ (n
1.0 forall a. Num a => a -> a -> a
+ n
s2c forall a. Num a => a -> a -> a
+ n
c2)
        a1 :: n
a1 = if Bool
is_hpf then -n
2.0 forall a. Num a => a -> a -> a
* n
a0 else n
2.0 forall a. Num a => a -> a -> a
* n
a0
        a2 :: n
a2 = n
a0
        b1 :: n
b1 = if Bool
is_hpf then n
2.0 forall a. Num a => a -> a -> a
* (n
c2 forall a. Num a => a -> a -> a
- n
1.0) forall a. Num a => a -> a -> a
* n
a0 else n
2.0 forall a. Num a => a -> a -> a
* (n
1.0 forall a. Num a => a -> a -> a
- n
c2) forall a. Num a => a -> a -> a
* n
a0
        b2 :: n
b2 = (n
1.0 forall a. Num a => a -> a -> a
- n
s2c forall a. Num a => a -> a -> a
+ n
c2) forall a. Num a => a -> a -> a
* n
a0
    in (n
a0,n
a1,n
a2,n
b1,n
b2)

-- | Two place infinite impulse response filter coefficients, (a0,b1,b2)
type Iir2Coef t = (t,t,t)

-- | rlpf coefficients, (a0,b1,b2).
rlpf_coef :: Floating n => (n -> n -> n) -> (n,n,n) -> Iir2Coef n
rlpf_coef :: forall n. Floating n => (n -> n -> n) -> (n, n, n) -> (n, n, n)
rlpf_coef n -> n -> n
max_f (n
radians_per_sample,n
f,n
rq) =
    let qr :: n
qr = n -> n -> n
max_f n
0.001 n
rq
        pf :: n
pf = n
f forall a. Num a => a -> a -> a
* n
radians_per_sample
        d :: n
d = forall a. Floating a => a -> a
tan (n
pf forall a. Num a => a -> a -> a
* n
qr forall a. Num a => a -> a -> a
* n
0.5)
        c :: n
c = (n
1.0 forall a. Num a => a -> a -> a
- n
d) forall a. Fractional a => a -> a -> a
/ (n
1.0 forall a. Num a => a -> a -> a
+ n
d)
        b1 :: n
b1 = (n
1.0 forall a. Num a => a -> a -> a
+ n
c) forall a. Num a => a -> a -> a
* forall a. Floating a => a -> a
cos n
pf
        b2 :: n
b2 = forall a. Num a => a -> a
negate n
c
        a0 :: n
a0 = (n
1.0 forall a. Num a => a -> a -> a
+ n
c forall a. Num a => a -> a -> a
- n
b1) forall a. Num a => a -> a -> a
* n
0.25
    in (n
a0,n
b1,n
b2)

-- | resonz coefficients, (a0,b1,b2).
resonz_coef :: Floating n => (n,n,n) -> Iir2Coef n
resonz_coef :: forall n. Floating n => (n, n, n) -> (n, n, n)
resonz_coef (n
radians_per_sample,n
f,n
rq) =
    let ff :: n
ff = n
f forall a. Num a => a -> a -> a
* n
radians_per_sample
        b :: n
b = n
ff forall a. Num a => a -> a -> a
* n
rq
        r :: n
r = n
1.0 forall a. Num a => a -> a -> a
- n
b forall a. Num a => a -> a -> a
* n
0.5
        two_r :: n
two_r = n
2.0 forall a. Num a => a -> a -> a
* n
r
        r2 :: n
r2 = n
r forall a. Num a => a -> a -> a
* n
r
        ct :: n
ct = (n
two_r forall a. Num a => a -> a -> a
* forall a. Floating a => a -> a
cos n
ff) forall a. Fractional a => a -> a -> a
/ (n
1.0 forall a. Num a => a -> a -> a
+ n
r2)
        b1 :: n
b1 = n
two_r forall a. Num a => a -> a -> a
* n
ct
        b2 :: n
b2 = forall a. Num a => a -> a
negate n
r2
        a0 :: n
a0 = (n
1.0 forall a. Num a => a -> a -> a
- n
r2) forall a. Num a => a -> a -> a
* n
0.5
    in (n
a0,n
b1,n
b2)

-- * Pinking

type Pinking_Param t = (t, t, t, t, t, t, t)

{- | Sample rate variable pinking filter.

https://www.musicdsp.org/en/latest/Filters/76-pink-noise-filter.html
-}
pinking_filter_freq_48000 :: Fractional t => Pinking_Param t
pinking_filter_freq_48000 :: forall t. Fractional t => Pinking_Param t
pinking_filter_freq_48000 = (t
4752.456, t
4030.961, t
2784.711, t
1538.461, t
357.681, t
70, t
30)

pinking_filter_freq_96000 :: Fractional t => Pinking_Param t
pinking_filter_freq_96000 :: forall t. Fractional t => Pinking_Param t
pinking_filter_freq_96000 = (t
8227.219, t
8227.219, t
6388.570, t
3302.754, t
479.412, t
151.070, t
54.264)

pinking_filter_freq_192000 :: Fractional t => Pinking_Param t
pinking_filter_freq_192000 :: forall t. Fractional t => Pinking_Param t
pinking_filter_freq_192000 = (t
9211.912, t
8621.096, t
8555.228, t
8292.754, t
518.334, t
163.712, t
240.241)

-- > pinking_filter_coef 48000 pinking_filter_freq_48000
pinking_filter_coef :: Floating t => t -> Pinking_Param t -> Pinking_Param t
pinking_filter_coef :: forall t. Floating t => t -> Pinking_Param t -> Pinking_Param t
pinking_filter_coef t
sr (t
f0, t
f1, t
f2, t
f3, t
f4, t
f5, t
f6) =
  let f :: t -> t
f t
n = forall a. Floating a => a -> a
exp ((- t
2) forall a. Num a => a -> a -> a
* forall a. Floating a => a
pi forall a. Num a => a -> a -> a
* t
n forall a. Fractional a => a -> a -> a
/ t
sr)
  in (t -> t
f t
f0, t -> t
f t
f1, t -> t
f t
f2, t -> t
f t
f3, t -> t
f t
f4, t -> t
f t
f5, t -> t
f t
f6)

pinking_filter_next :: Floating t => Pinking_Param t -> Pinking_Param t -> t -> (t, Pinking_Param t)
pinking_filter_next :: forall t.
Floating t =>
Pinking_Param t -> Pinking_Param t -> t -> (t, Pinking_Param t)
pinking_filter_next (t
k0, t
k1, t
k2, t
k3, t
k4, t
k5, t
k6) (t
m0, t
m1, t
m2, t
m3, t
m4, t
m5, t
m6) t
white =
  let (t
b0, t
b1, t
b2, t
b3, t
b4, t
b5, t
b6) = ((t
k0 forall a. Num a => a -> a -> a
* t
white) forall a. Num a => a -> a -> a
+ (t
k0 forall a. Num a => a -> a -> a
* t
m0)
                                     ,(t
k1 forall a. Num a => a -> a -> a
* t
white) forall a. Num a => a -> a -> a
+ (t
k1 forall a. Num a => a -> a -> a
* t
m1)
                                     ,(t
k2 forall a. Num a => a -> a -> a
* t
white) forall a. Num a => a -> a -> a
+ (t
k2 forall a. Num a => a -> a -> a
* t
m2)
                                     ,(t
k3 forall a. Num a => a -> a -> a
* t
white) forall a. Num a => a -> a -> a
+ (t
k3 forall a. Num a => a -> a -> a
* t
m3)
                                     ,(t
k4 forall a. Num a => a -> a -> a
* t
white) forall a. Num a => a -> a -> a
+ (t
k4 forall a. Num a => a -> a -> a
* t
m4)
                                     ,(t
k5 forall a. Num a => a -> a -> a
* t
white) forall a. Num a => a -> a -> a
+ (t
k5 forall a. Num a => a -> a -> a
* t
m5)
                                     ,(t
k6 forall a. Num a => a -> a -> a
* t
white) forall a. Num a => a -> a -> a
+ (t
k6 forall a. Num a => a -> a -> a
* t
m6))
      pink :: t
pink = t
b0 forall a. Num a => a -> a -> a
+ t
b1 forall a. Num a => a -> a -> a
+ t
b2 forall a. Num a => a -> a -> a
+ t
b3 forall a. Num a => a -> a -> a
+ t
b4 forall a. Num a => a -> a -> a
+ t
b5 forall a. Num a => a -> a -> a
+ t
white forall a. Num a => a -> a -> a
- t
b6
  in (t
pink, (t
b0, t
b1, t
b2, t
b3, t
b4, t
b5, t
b6))