Are there any "good" low-computational-cost approximations for parameterized Kaiser window generation suitable for small systems/languages that do not include Bessel functions in their standard library?

(maybe only to 8 to 12-bit accuracy, and something that can be run in Basic on an Apple II, et.al.)

  • 1
    I don't know how good this technique is, but CORDIC has been adapted to compute the Kaiser window: link.springer.com/article/10.1007%2Fs11265-013-0781-z – MBaz Feb 18 '17 at 21:59
  • Apple II? Computational wise LUT will probably be the lowest costing approach. Only the storage is the concern. Or may be you can truncate the polynomial Bessel series with enough terms for 12 bit accuracy. – Fat32 Feb 18 '17 at 23:46
  • 1
    Kaiser is an approximation of the DPSS (discrete prolate spheroidal sequence) a.k.a. Slepian window. So do you need Kaiser in itself or would you accept something with similar performance? – Olli Niemitalo Mar 4 '17 at 8:49

You could try the exponential window:

$$w_n=\frac{ \exp \left[\alpha \sqrt{1-\left(\frac{n-M}{M}\right)^2}\right]}{\exp(\alpha)}$$ $$\alpha=-427.5*10^{-6} A_s^2+0.1808*A_s-3.516$$

or the hyperbolic cosine window:

$$w_n=\frac{\cosh \left[ \alpha \sqrt{1-\left(\frac{n-M}{M}\right)^2} \right]}{\cosh(\alpha)}$$ $$\alpha=-325.1*10^{-6}*A_s^2+0.1677*A_s-3.149$$

Both have similar results with themselves and the Kaiser window, while being less computational intensive.

If you want to go with Robert's answer, you could simplify the Bessel function by using Horner's method, or by applying a Gauss-Kronrod (or whatever other) integration. For example, I got this formula (after some reductions&co):

$$f(x)=\frac{\left[ ~~ \cosh(x) ~+~ 2*(\cosh(0.970941817426052*x) ~+~ \cosh(0.8854560256532099*x) ~+~ \cosh(0.7485107481711011*x) ~+~ \cosh(0.5680647467311558*x) ~+~ \cosh(0.3546048870425356*x) ~+~ \cosh(0.120536680255323*x)) \right]}{13}$$

For $\beta=12.26526$ ($A_s=120dB$), results in a these values:

$I_0(\beta)=24430.40185694905$ $f(\beta)=24430.4018627702$

which, IMHO, is pretty good and it's faster than using giga-km long numbers you'd deal with in the case of the classic formula, either as a sum or as Horner -- either way, I think even long double would not be enough for 20 terms (which is how I've truncated it).

just implement the Bessel function. it's easy.

i always like my window definitions centered about zero, since pretty much all of them are even symmetry.

i'll do this in discrete-time, but it's essentially the same thing in continuous-time.

Kaiser window:

$$ w[n] \triangleq \begin{cases} \frac{1}{I_0(\beta)} I_0\left(\beta \sqrt{1 - \left(\frac{n}{M/2}\right)^2} \right) \quad \quad & |n| \le M/2 \\ 0 & |n|>M/2 \\ \end{cases} $$

$I_0(x)$ is the 0th-order modified Bessel function of the 1st kind. $M+1$ is the number of non-zero samples or FIR taps (the FIR filter order is $M$ and, in my centered and symmetrical case, must be even). $\beta$ is a "shape parameter" and O&S recommend this heuristic:

$$ \beta = \begin{cases} 0.1102 \cdot (A-8.7) & A>50 \\ 0.5842 \cdot (A-21)^{2/5} + 0.07886 \cdot (A-21) \quad & 21 \le A \le 50 \\ 0.0 & A<21 \\ \end{cases}$$

$$ M = 2 \left\lceil \frac{A-8}{4.57 \cdot \Delta\omega} \right\rceil $$

$A$ is the desired stopband attention in dB and $\Delta\omega$ is the desired width of the transition band in normalized angular frequency.

finally, the Bessel is evaluated as:

$$ I_0(x) = 1 \ + \ \lim_{K \to \infty} \ \sum\limits_{k=1}^{K} \left(-\frac{x^2}{4}\right)^{k} (k!)^{-2} $$

when you evaluate this with a computer, pick a $K$ decently large (my guess is that $K=8$ may be good enough) and evaluate the summation starting with $k=K$ and work it backwards to $k=1$ to keep numerical accuracy. you might want to use Horner's method.

$$ I_0(x) \approx 1 + x^2\left( \tfrac{-1}{(1!)^2 \, 4^1} + x^2\left(\tfrac{+1}{(2!)^2 \, 4^2} + x^2\left(... + \, x^2\left(\tfrac{(-1)^{K-1}}{((K-1)!)^2 \, 4^{K-1}} + x^2 \tfrac{(-1)^K}{(K!)^2 \, 4^K} \right) \right) \right) \right) $$

you can evaluate all of the $(k!)^{-2}$ in advance with a short table.

Your Answer

 
discard

By clicking "Post Your Answer", you acknowledge that you have read our updated terms of service, privacy policy and cookie policy, and that your continued use of the website is subject to these policies.

Not the answer you're looking for? Browse other questions tagged or ask your own question.