From: Kostya Shishkov <kostya.shishkov@gmail.com>
Date: Sat, 27 Oct 2018 16:57:11 +0000 (+0200)
Subject: dsp: implement Kaiser-Bessel derived window generation
X-Git-Url: https://git.nihav.org/?p=nihav.git;a=commitdiff_plain;h=9669f269bf8c09420b64dfe1347581bbfdf12e4a

dsp: implement Kaiser-Bessel derived window generation
---

diff --git a/src/dsp/window.rs b/src/dsp/window.rs
index d59a005..42b6b6d 100644
--- a/src/dsp/window.rs
+++ b/src/dsp/window.rs
@@ -4,7 +4,7 @@ use std::f32::consts;
 pub enum WindowType {
     Square,
     Sine,
-    KaiserBessel,
+    KaiserBessel(f32),
 }
 
 pub fn generate_window(mode: WindowType, scale: f32, size: usize, half: bool, dst: &mut [f32]) {
@@ -23,8 +23,30 @@ pub fn generate_window(mode: WindowType, scale: f32, size: usize, half: bool, ds
                     dst[n] = (((n as f32) + 0.5) * param).sin() * scale;
                 }
             },
-        WindowType::KaiserBessel => {
-unimplemented!();
+        WindowType::KaiserBessel(alpha) => {
+                let dlen = if half { size as f32 } else { (size as f32) * 0.5 };
+                let alpha2 = ((alpha * consts::PI / dlen) * (alpha * consts::PI / dlen)) as f64;
+
+                let mut kb: Vec<f64> = Vec::with_capacity(size);
+                let mut sum = 0.0;
+                for n in 0..size {
+                    let b = bessel_i0(((n * (size - n)) as f64) * alpha2);
+                    sum += b;
+                    kb.push(sum);
+                }
+                sum += 1.0;
+                for n in 0..size {
+                    dst[n] = (kb[n] / sum).sqrt() as f32;
+                }
             },
     };
 }
+
+fn bessel_i0(inval: f64) -> f64 {
+    let mut val: f64 = 1.0;
+    for n in (1..64).rev() {
+        val *= inval / ((n * n) as f64);
+        val += 1.0;
+    }
+    val
+}