module DSP.Estimation.Spectral.ARMA (arma_mywe) where
import Data.Array
import Data.Complex
import DSP.Correlation
arma_mywe ::
(RealFloat b, Integral i, Ix i) =>
Array i (Complex b) -> i -> i -> Array i (Complex b)
arma_mywe :: forall b i.
(RealFloat b, Integral i, Ix i) =>
Array i (Complex b) -> i -> i -> Array i (Complex b)
arma_mywe Array i (Complex b)
x i
p i
q = Array i (Complex b)
a'
where r :: Array i (Complex b)
r = forall i e. Ix i => (i, i) -> [(i, e)] -> Array i e
array (i
qforall a. Num a => a -> a -> a
-i
2forall a. Num a => a -> a -> a
*i
pforall a. Num a => a -> a -> a
+i
1,i
qforall a. Num a => a -> a -> a
+i
p) [ (i
k, forall a b.
(Ix a, Integral a, RealFloat b) =>
Array a (Complex b) -> a -> Complex b
rxx_u Array i (Complex b)
x i
k) | i
k <- [(i
qforall a. Num a => a -> a -> a
-i
2forall a. Num a => a -> a -> a
*i
pforall a. Num a => a -> a -> a
+i
1)..(i
qforall a. Num a => a -> a -> a
+i
p)] ]
a' :: Array i (Complex b)
a' = forall i e. Ix i => (i, i) -> [(i, e)] -> Array i e
array (i
1,i
p) [ (i
k, Array (i, i) (Complex b)
aforall i e. Ix i => Array i e -> i -> e
!(i
p,i
k)) | i
k <- [i
1..i
p] ]
a :: Array (i, i) (Complex b)
a = forall i e. Ix i => (i, i) -> [(i, e)] -> Array i e
array ((i
1,i
1),(i
p,i
p)) [ ((i
k,i
i), i -> i -> Complex b
ak i
k i
i) | i
k <- [i
1..i
p], i
i <- [i
1..i
k] ]
b :: Array (i, i) (Complex b)
b = forall i e. Ix i => (i, i) -> [(i, e)] -> Array i e
array ((i
1,i
1),(i
pforall a. Num a => a -> a -> a
-i
1,i
pforall a. Num a => a -> a -> a
-i
1)) [ ((i
k,i
i), i -> i -> Complex b
bk i
k i
i) | i
k <- [i
1..(i
pforall a. Num a => a -> a -> a
-i
1)], i
i <- [i
1..i
k] ]
rho :: Array i (Complex b)
rho = forall i e. Ix i => (i, i) -> [(i, e)] -> Array i e
array (i
1,i
pforall a. Num a => a -> a -> a
-i
1) [ (i
k, i -> Complex b
rhok i
k) | i
k <- [i
1..(i
pforall a. Num a => a -> a -> a
-i
1)] ]
ak :: i -> i -> Complex b
ak i
1 i
1 = -Array i (Complex b)
rforall i e. Ix i => Array i e -> i -> e
!(i
qforall a. Num a => a -> a -> a
+i
1) forall a. Fractional a => a -> a -> a
/ Array i (Complex b)
rforall i e. Ix i => Array i e -> i -> e
!i
q
ak i
k i
i | i
iforall a. Eq a => a -> a -> Bool
==i
k = -(Array i (Complex b)
rforall i e. Ix i => Array i e -> i -> e
!(i
qforall a. Num a => a -> a -> a
+i
k) forall a. Num a => a -> a -> a
+ forall (t :: * -> *) a. (Foldable t, Num a) => t a -> a
sum [ Array (i, i) (Complex b)
aforall i e. Ix i => Array i e -> i -> e
!(i
kforall a. Num a => a -> a -> a
-i
1,i
l) forall a. Num a => a -> a -> a
* Array i (Complex b)
rforall i e. Ix i => Array i e -> i -> e
!(i
qforall a. Num a => a -> a -> a
+i
kforall a. Num a => a -> a -> a
-i
l) | i
l <- [i
1..(i
kforall a. Num a => a -> a -> a
-i
1)] ] ) forall a. Fractional a => a -> a -> a
/ Array i (Complex b)
rhoforall i e. Ix i => Array i e -> i -> e
!(i
kforall a. Num a => a -> a -> a
-i
1)
| Bool
otherwise = Array (i, i) (Complex b)
aforall i e. Ix i => Array i e -> i -> e
!(i
kforall a. Num a => a -> a -> a
-i
1,i
i) forall a. Num a => a -> a -> a
+ Array (i, i) (Complex b)
aforall i e. Ix i => Array i e -> i -> e
!(i
k,i
k) forall a. Num a => a -> a -> a
* Array (i, i) (Complex b)
bforall i e. Ix i => Array i e -> i -> e
!(i
kforall a. Num a => a -> a -> a
-i
1,i
kforall a. Num a => a -> a -> a
-i
i)
bk :: i -> i -> Complex b
bk i
1 i
1 = -Array i (Complex b)
rforall i e. Ix i => Array i e -> i -> e
!(i
qforall a. Num a => a -> a -> a
-i
1) forall a. Fractional a => a -> a -> a
/ Array i (Complex b)
rforall i e. Ix i => Array i e -> i -> e
!i
q
bk i
k i
i | i
iforall a. Eq a => a -> a -> Bool
==i
k = -(Array i (Complex b)
rforall i e. Ix i => Array i e -> i -> e
!(i
qforall a. Num a => a -> a -> a
-i
k) forall a. Num a => a -> a -> a
+ forall (t :: * -> *) a. (Foldable t, Num a) => t a -> a
sum [ Array (i, i) (Complex b)
bforall i e. Ix i => Array i e -> i -> e
!(i
kforall a. Num a => a -> a -> a
-i
1,i
l) forall a. Num a => a -> a -> a
* Array i (Complex b)
rforall i e. Ix i => Array i e -> i -> e
!(i
qforall a. Num a => a -> a -> a
-i
kforall a. Num a => a -> a -> a
-i
l) | i
l <- [i
1..(i
kforall a. Num a => a -> a -> a
-i
1)] ] ) forall a. Fractional a => a -> a -> a
/ Array i (Complex b)
rhoforall i e. Ix i => Array i e -> i -> e
!(i
kforall a. Num a => a -> a -> a
-i
1)
| Bool
otherwise = Array (i, i) (Complex b)
bforall i e. Ix i => Array i e -> i -> e
!(i
kforall a. Num a => a -> a -> a
-i
1,i
i) forall a. Num a => a -> a -> a
+ Array (i, i) (Complex b)
bforall i e. Ix i => Array i e -> i -> e
!(i
k,i
k) forall a. Num a => a -> a -> a
* Array (i, i) (Complex b)
aforall i e. Ix i => Array i e -> i -> e
!(i
kforall a. Num a => a -> a -> a
-i
1,i
kforall a. Num a => a -> a -> a
-i
i)
rhok :: i -> Complex b
rhok i
1 = (Complex b
1 forall a. Num a => a -> a -> a
- Array (i, i) (Complex b)
aforall i e. Ix i => Array i e -> i -> e
!(i
1,i
1) forall a. Num a => a -> a -> a
* Array (i, i) (Complex b)
bforall i e. Ix i => Array i e -> i -> e
!(i
1,i
1)) forall a. Num a => a -> a -> a
* Array i (Complex b)
rforall i e. Ix i => Array i e -> i -> e
!i
q
rhok i
k = (Complex b
1 forall a. Num a => a -> a -> a
- Array (i, i) (Complex b)
aforall i e. Ix i => Array i e -> i -> e
!(i
k,i
k) forall a. Num a => a -> a -> a
* Array (i, i) (Complex b)
bforall i e. Ix i => Array i e -> i -> e
!(i
k,i
k)) forall a. Num a => a -> a -> a
* Array i (Complex b)
rhoforall i e. Ix i => Array i e -> i -> e
!(i
kforall a. Num a => a -> a -> a
-i
1)