DSP.Estimation.Frequency.PerMax

Description

This module implements an algorithm to maximize the peak value of a DFT/FFT. It is based off an aticle by Mark Sullivan from Personal Engineering Magazine.

Maximizes

S(w) = 1/N * sum(k=0,N-1) |x[k] * e^(-jwk)|^2

which is equivalent to solving

S'(w) = Im{X(w) * ~Y(w)} = 0

where

X(w) = sum(k=0,N-1) (x[k] * e^(-jwk)) Y(w) = X'(w) = sum(k=0,N-1) (k * x[k] * e^(-jwk))

This algorithm used the bisection method for finding the zero of a function. The search area is +- half a bin width.

Regula falsi requires an additional (x,f(x)) pair which is expensive in this case. Newton's method could be used but requires S''(w), which takes twice as long to caculate as S'(w). Brent's method may be best here, but it also requires three (x,f(x)) pairs

Synopsis

# Documentation

Arguments

 :: (Ix a, Integral a, RealFloat b) => Array a (Complex b) X[k] -> a k -> b w

Discrete frequency periodigram maximizer