module Synthesizer.Plain.Filter.LinearPredictive where

import Synthesizer.Plain.Analysis (scalarProduct)

import qualified Data.List.Match as ListMatch
import qualified Data.List as List

import qualified Algebra.Field    as Field

import NumericPrelude.Numeric
import NumericPrelude.Base
import Prelude ()


{- |
Determine optimal filter coefficients and residue by adaptive approximation.
The number of initial filter coefficients is used as filter order.
-}
approxCoefficients :: Field.C a =>
   a -> [a] -> [a] -> [(a,[a])]
approxCoefficients :: forall a. C a => a -> [a] -> [a] -> [(a, [a])]
approxCoefficients a
k [a]
mask0 [a]
xs =
   let infixes :: [[a]]
infixes = ([a] -> [a]) -> [[a]] -> [[a]]
forall a b. (a -> b) -> [a] -> [b]
map ([a] -> [a] -> [a]
forall b a. [b] -> [a] -> [a]
ListMatch.take [a]
mask0) ([a] -> [[a]]
forall a. [a] -> [[a]]
List.tails [a]
xs)
       targets :: [a]
targets = [a] -> [a] -> [a]
forall b a. [b] -> [a] -> [a]
ListMatch.drop [a]
mask0 [a]
xs
   in  ((a, [a]) -> ([a], a) -> (a, [a]))
-> (a, [a]) -> [([a], a)] -> [(a, [a])]
forall b a. (b -> a -> b) -> b -> [a] -> [b]
scanl
          (\(a
_,[a]
mask) ([a]
infx,a
target) ->
              let residue :: a
residue = a
target a -> a -> a
forall a. C a => a -> a -> a
- [a] -> [a] -> a
forall y. C y => T y -> T y -> y
scalarProduct [a]
mask [a]
infx
                  norm2 :: a
norm2 = [a] -> [a] -> a
forall y. C y => T y -> T y -> y
scalarProduct [a]
infx [a]
infx
              in  (a
residue,
                   [a]
mask [a] -> [a] -> [a]
forall a. C a => a -> a -> a
+ (a -> a) -> [a] -> [a]
forall a b. (a -> b) -> [a] -> [b]
map ((a
ka -> a -> a
forall a. C a => a -> a -> a
*a
residuea -> a -> a
forall a. C a => a -> a -> a
/a
norm2)a -> a -> a
forall a. C a => a -> a -> a
*) [a]
infx))
          (a
forall a. C a => a
zero,[a]
mask0) ([[a]] -> [a] -> [([a], a)]
forall a b. [a] -> [b] -> [(a, b)]
zip [[a]]
infixes [a]
targets)
{-
mapM print $ take 20 $ drop 2000 $ approxCoefficients (1::Double) [0,0,0,0.1] (iterate (1+) 100)


mapM print $ take 20 $ drop 10000 $ approxCoefficients (0.2::Double) [0.1,0] (map sin (iterate (0.03+) 0))

must yield coefficients [-1, 2*cos(0.03::Double)]
-}