{-# LANGUAGE NoImplicitPrelude #-} {- | Copyright : (c) Henning Thielemann 2008 License : GPL Maintainer : synthesizer@henning-thielemann.de Stability : provisional Portability : requires multi-parameter type classes -} module Synthesizer.Plain.Filter.Recursive.MovingAverage (sumsStaticInt, modulatedFrac, ) where import qualified Synthesizer.Plain.Signal as Sig -- import qualified Synthesizer.Plain.Modifier as Modifier import qualified Synthesizer.Plain.Filter.Recursive.Integration as Integration import Synthesizer.Plain.Filter.NonRecursive (delay, ) import qualified Algebra.Module as Module import qualified Algebra.RealField as RealField -- import qualified Algebra.Field as Field -- import qualified Algebra.Ring as Ring import qualified Algebra.Additive as Additive import Control.Monad.Fix (fix) import Data.List (tails) import NumericPrelude.Base import NumericPrelude.Numeric {- | Like 'Synthesizer.Plain.Filter.NonRecursive.sums' but in a recursive form. This needs only linear time (independent of the window size) but may accumulate rounding errors. @ ys = xs * (1,0,0,0,-1) \/ (1,-1) ys * (1,-1) = xs * (1,0,0,0,-1) ys = xs * (1,0,0,0,-1) + ys * (0,1) @ -} sumsStaticInt :: (Additive.C v) => Int -> Sig.T v -> Sig.T v sumsStaticInt n xs = fix (\ys -> let (xs0,xs1) = splitAt n xs in (xs0 ++ (xs1-xs)) + (zero:ys)) {- staticInt :: (Module.C a v, Additive.C v) => Int -> Sig.T v -> Sig.T v staticInt n xs = -} {- Sum of a part of a vector with negative sign for reverse order. It adds from @from@ (inclusively) to @to@ (exclusively), that is, it sums up @abs (to-from)@ values -} _sumFromTo :: (Additive.C v) => Int -> Int -> Sig.T v -> v _sumFromTo from to = if from <= to then sum . take (to-from) . drop from else negate . sum . take (from-to) . drop to {- It would be a nice approach to interpolate not just linearly at the borders but in a way that the cut-off frequency is perfectly suppressed. However suppression depends on the phase shift of the wave. Actually, we could use a complex factor, but does this help? -} sumFromToFrac :: (RealField.C a, Module.C a v) => a -> a -> Sig.T v -> v sumFromToFrac from to xs = let (fromInt, fromFrac) = splitFraction from (toInt, toFrac) = splitFraction to in case compare fromInt toInt of EQ -> (to-from) *> (xs !! fromInt) LT -> sum $ zipWith id (((1-fromFrac) *>) : replicate (toInt-fromInt-1) id ++ (toFrac *>) : []) $ drop fromInt xs GT -> negate $ sum $ zipWith id (((1-toFrac) *>) : replicate (fromInt-toInt-1) id ++ (fromFrac *>) : []) $ drop toInt xs {- | Sig.T a must contain only non-negative elements. -} sumDiffsModulated :: (RealField.C a, Module.C a v) => a -> Sig.T a -> Sig.T v -> Sig.T v sumDiffsModulated d ds = -- prevent negative d's since 'drop' cannot restore past values zipWith3 sumFromToFrac ((d+1) : ds) (map (1+) ds) . init . init . tails . (zero:) {- zipWith3 sumFromToFrac (d : map (subtract 1) ds) ds . init . tails -} _sumsModulated :: (RealField.C a, Module.C a v) => Int -> Sig.T a -> Sig.T v -> Sig.T v _sumsModulated maxDInt ds xs = let maxD = fromIntegral maxDInt posXs = sumDiffsModulated 0 ds xs negXs = sumDiffsModulated maxD (map (maxD-) ds) (delay maxDInt xs) in Integration.run (posXs - negXs) {- | Shift sampling points by a half sample period in order to preserve signals for window widths below 1. -} sumsModulatedHalf :: (RealField.C a, Module.C a v) => Int -> Sig.T a -> Sig.T v -> Sig.T v sumsModulatedHalf maxDInt ds xs = let maxD = fromIntegral maxDInt d0 = maxD+0.5 delXs = delay maxDInt xs posXs = sumDiffsModulated d0 (map (d0+) ds) delXs negXs = sumDiffsModulated d0 (map (d0-) ds) delXs in Integration.run (posXs - negXs) {- *Synthesizer.Plain.Filter.NonRecursive> movingAverageModulated 10 (replicate 10 (3::Double) ++ [1.1,2.2,2.6,0.7,0.1,0.1]) (repeat (1::Double)) [0.5,0.6666666666666666,0.8333333333333333,1.0,1.0,1.0,1.0,1.0,1.0,1.0,0.9999999999999999,1.0,0.9999999999999998,0.999999999999999,0.9999999999999942,0.9999999999999942] -} modulatedFrac :: (RealField.C a, Module.C a v) => Int -> Sig.T a -> Sig.T v -> Sig.T v modulatedFrac maxDInt ds xs = zipWith (\d y -> recip (2*d) *> y) ds $ sumsModulatedHalf maxDInt ds xs