{-# LANGUAGE NoImplicitPrelude #-}
{-# LANGUAGE MultiParamTypeClasses #-}
{-# LANGUAGE FlexibleInstances #-}
module Synthesizer.Plain.Filter.Recursive.AllpassPoly where
import qualified Algebra.Module as Module
import qualified Algebra.RealTranscendental as RealTrans
import qualified Algebra.Transcendental as Trans
import qualified Algebra.Field as Field
import qualified Algebra.ZeroTestable as ZeroTestable
import Number.Complex (cis,(+:),real,imag)
import qualified Number.Complex as Complex
import Orthogonals(Scalar,one_ket_solution)
import qualified Prelude as P
import NumericPrelude.Numeric
import NumericPrelude.Base
newtype Parameter a = Parameter [a]
deriving Int -> Parameter a -> ShowS
[Parameter a] -> ShowS
Parameter a -> String
(Int -> Parameter a -> ShowS)
-> (Parameter a -> String)
-> ([Parameter a] -> ShowS)
-> Show (Parameter a)
forall a. Show a => Int -> Parameter a -> ShowS
forall a. Show a => [Parameter a] -> ShowS
forall a. Show a => Parameter a -> String
forall a.
(Int -> a -> ShowS) -> (a -> String) -> ([a] -> ShowS) -> Show a
$cshowsPrec :: forall a. Show a => Int -> Parameter a -> ShowS
showsPrec :: Int -> Parameter a -> ShowS
$cshow :: forall a. Show a => Parameter a -> String
show :: Parameter a -> String
$cshowList :: forall a. Show a => [Parameter a] -> ShowS
showList :: [Parameter a] -> ShowS
Show
shiftParam :: (Scalar a, P.Fractional a, Trans.C a) =>
Int -> a -> a -> Parameter a
shiftParam :: forall a.
(Scalar a, Fractional a, C a) =>
Int -> a -> a -> Parameter a
shiftParam Int
order a
weight a
phase =
let
normalVector :: [a]
normalVector = (a -> a) -> [a] -> [a]
forall a b. (a -> b) -> [a] -> [b]
map a -> a
forall a. C a => a -> a
negate
((Int -> a) -> [Int] -> [a]
forall a b. (a -> b) -> [a] -> [b]
map (a -> Int -> a -> Int -> Int -> a
forall a. C a => a -> Int -> a -> Int -> Int -> a
scalarProdScrewExp a
weight Int
order a
phase Int
0) [Int
1..Int
order])
normalMatrix :: [[a]]
normalMatrix = (Int -> [a]) -> [Int] -> [[a]]
forall a b. (a -> b) -> [a] -> [b]
map (\Int
j ->
(Int -> a) -> [Int] -> [a]
forall a b. (a -> b) -> [a] -> [b]
map (a -> Int -> a -> Int -> Int -> a
forall a. C a => a -> Int -> a -> Int -> Int -> a
scalarProdScrewExp a
weight Int
order a
phase Int
j) [Int
1..Int
order]) [Int
1..Int
order]
in [a] -> Parameter a
forall a. [a] -> Parameter a
Parameter ([[a]] -> [a] -> [a]
forall a. (Scalar a, Fractional a) => [[a]] -> [a] -> [a]
one_ket_solution [[a]]
normalMatrix [a]
normalVector)
makePhase :: (RealTrans.C a, ZeroTestable.C a) => Parameter a -> a -> a
makePhase :: forall a. (C a, C a) => Parameter a -> a -> a
makePhase (Parameter [a]
ks) a
frequency =
let omega :: a
omega = a
2a -> a -> a
forall a. C a => a -> a -> a
*a
forall a. C a => a
pi a -> a -> a
forall a. C a => a -> a -> a
* a
frequency
omegas :: [a]
omegas = (a -> a) -> a -> [a]
forall a. (a -> a) -> a -> [a]
iterate (a
omegaa -> a -> a
forall a. C a => a -> a -> a
+) a
omega
denom :: T a
denom = T a
1T a -> T a -> T a
forall a. C a => a -> a -> a
+[T a] -> T a
forall a. C a => [a] -> a
sum ((a -> a -> T a) -> [a] -> [a] -> [T a]
forall a b c. (a -> b -> c) -> [a] -> [b] -> [c]
zipWith (\a
k a
w -> a
ka -> a -> a
forall a. C a => a -> a -> a
*a -> a
forall a. C a => a -> a
cos a
w a -> a -> T a
forall a. a -> a -> T a
+: a
ka -> a -> a
forall a. C a => a -> a -> a
*a -> a
forall a. C a => a -> a
sin a
w) [a]
ks [a]
omegas)
in a
2 a -> a -> a
forall a. C a => a -> a -> a
* T a -> a
forall a. (C a, C a) => T a -> a
Complex.phase T a
denom a -> a -> a
forall a. C a => a -> a -> a
- a
omegaa -> a -> a
forall a. C a => a -> a -> a
*(Int -> a
forall a b. (C a, C b) => a -> b
fromIntegral ([a] -> Int
forall a. [a] -> Int
forall (t :: * -> *) a. Foldable t => t a -> Int
length [a]
ks))
scalarProdScrewExp :: Trans.C a => a -> Int -> a -> Int -> Int -> a
scalarProdScrewExp :: forall a. C a => a -> Int -> a -> Int -> Int -> a
scalarProdScrewExp a
r Int
order a
phase Int
k Int
j =
let (a
intCos,a
intSin) = a -> Int -> (a, a)
forall a. C a => a -> Int -> (a, a)
integrateScrewExp a
r (Int
kInt -> Int -> Int
forall a. C a => a -> a -> a
+Int
jInt -> Int -> Int
forall a. C a => a -> a -> a
-Int
order)
in a
2 a -> a -> a
forall a. C a => a -> a -> a
* ((a, a) -> a
forall a b. (a, b) -> a
fst (a -> Int -> (a, a)
forall a. C a => a -> Int -> (a, a)
integrateScrewExp a
r (Int
kInt -> Int -> Int
forall a. C a => a -> a -> a
-Int
j)) a -> a -> a
forall a. C a => a -> a -> a
-
(a -> a
forall a. C a => a -> a
cos a
phase a -> a -> a
forall a. C a => a -> a -> a
* a
intCos a -> a -> a
forall a. C a => a -> a -> a
+ a -> a
forall a. C a => a -> a
sin a
phase a -> a -> a
forall a. C a => a -> a -> a
* a
intSin))
screwProd :: Trans.C a => Int -> a -> Int -> Int -> a -> a
screwProd :: forall a. C a => Int -> a -> Int -> Int -> a -> a
screwProd Int
order a
phase Int
k Int
j a
omega =
let z0 :: T a
z0 = a -> T a
forall a. C a => a -> T a
cis (Int -> a
forall a b. (C a, C b) => a -> b
fromIntegral Int
k a -> a -> a
forall a. C a => a -> a -> a
* a
omega) T a -> T a -> T a
forall a. C a => a -> a -> a
-
a -> T a
forall a. C a => a -> T a
cis a
phase T a -> T a -> T a
forall a. C a => a -> a -> a
* a -> T a
forall a. C a => a -> T a
cis (Int -> a
forall a b. (C a, C b) => a -> b
fromIntegral (Int
orderInt -> Int -> Int
forall a. C a => a -> a -> a
-Int
k) a -> a -> a
forall a. C a => a -> a -> a
* a
omega)
z1 :: T a
z1 = a -> T a
forall a. C a => a -> T a
cis (Int -> a
forall a b. (C a, C b) => a -> b
fromIntegral Int
j a -> a -> a
forall a. C a => a -> a -> a
* a
omega) T a -> T a -> T a
forall a. C a => a -> a -> a
-
a -> T a
forall a. C a => a -> T a
cis a
phase T a -> T a -> T a
forall a. C a => a -> a -> a
* a -> T a
forall a. C a => a -> T a
cis (Int -> a
forall a b. (C a, C b) => a -> b
fromIntegral (Int
orderInt -> Int -> Int
forall a. C a => a -> a -> a
-Int
j) a -> a -> a
forall a. C a => a -> a -> a
* a
omega)
in T a -> a
forall a. T a -> a
real T a
z0 a -> a -> a
forall a. C a => a -> a -> a
* T a -> a
forall a. T a -> a
real T a
z1 a -> a -> a
forall a. C a => a -> a -> a
+ T a -> a
forall a. T a -> a
imag T a
z0 a -> a -> a
forall a. C a => a -> a -> a
* T a -> a
forall a. T a -> a
imag T a
z1
integrateScrewExp :: Trans.C a => a -> Int -> (a,a)
integrateScrewExp :: forall a. C a => a -> Int -> (a, a)
integrateScrewExp a
r Int
kInt =
let k :: a
k = Int -> a
forall a b. (C a, C b) => a -> b
fromIntegral Int
kInt
q :: a
q = (a -> a
forall a. C a => a -> a
exp (a
2a -> a -> a
forall a. C a => a -> a -> a
*a
forall a. C a => a
pia -> a -> a
forall a. C a => a -> a -> a
*a
r) a -> a -> a
forall a. C a => a -> a -> a
- a
1) a -> a -> a
forall a. C a => a -> a -> a
/ (a
ra -> Integer -> a
forall a. C a => a -> Integer -> a
^Integer
2 a -> a -> a
forall a. C a => a -> a -> a
+ a
ka -> Integer -> a
forall a. C a => a -> Integer -> a
^Integer
2)
in (a
ra -> a -> a
forall a. C a => a -> a -> a
*a
q, -a
ka -> a -> a
forall a. C a => a -> a -> a
*a
q)
integrateNum :: (Field.C a, Module.C a v) => Int -> (a,a) -> (a->v) -> v
integrateNum :: forall a v. (C a, C a v) => Int -> (a, a) -> (a -> v) -> v
integrateNum Int
n (a
lo,a
hi) a -> v
f =
let xs :: [a]
xs = (Int -> a) -> [Int] -> [a]
forall a b. (a -> b) -> [a] -> [b]
map (\Int
k -> a
lo a -> a -> a
forall a. C a => a -> a -> a
+ (a
hia -> a -> a
forall a. C a => a -> a -> a
-a
lo) 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)
[Int
1..(Int
nInt -> Int -> Int
forall a. C a => a -> a -> a
-Int
1)]
in ((a
hia -> a -> a
forall a. C a => a -> a -> a
-a
lo) a -> a -> a
forall a. C a => a -> a -> a
/ Int -> a
forall a b. (C a, C b) => a -> b
fromIntegral Int
n) a -> v -> v
forall a v. C a v => a -> v -> v
*>
((v -> v -> v) -> v -> [v] -> v
forall b a. (b -> a -> b) -> b -> [a] -> b
forall (t :: * -> *) b a.
Foldable t =>
(b -> a -> b) -> b -> t a -> b
foldl v -> v -> v
forall a. C a => a -> a -> a
(+) ((a
1a -> a -> a
forall a. C a => a -> a -> a
/a
2 a -> a -> a
forall a. a -> a -> a
`asTypeOf` a
lo) a -> v -> v
forall a v. C a v => a -> v -> v
*> (a -> v
f a
lo v -> v -> v
forall a. C a => a -> a -> a
+ a -> v
f a
hi))
((a -> v) -> [a] -> [v]
forall a b. (a -> b) -> [a] -> [b]
map a -> v
f [a]
xs))