changeset 399:ac828725aeae resample

Update for Kaiser parameters changes; ensure Kaiser window is generated odd-length so with exactly 1 at the "centre"

--- a/src/may/stream/filter.yeti	Thu Sep 26 14:47:27 2013 +0100
+++ b/src/may/stream/filter.yeti	Thu Sep 26 15:04:45 2013 +0100
@@ -441,12 +441,20 @@
  * and transition bandwidth Hz at the given sampleRate. The filter
  * will contain an odd number of samples.
  */
-kaiserSincFilterFor { n, alpha, bandwidth, sampleRate } =
-   (kw = win.kaiserForBandwidth alpha bandwidth sampleRate;
-    kw = if (vec.length kw) % 2 == 0 then vec.concat [vec.zeros 1, kw] else kw fi;
-    sw = win.sinc (n*2) (vec.length kw);
+kaiserSincFilterFor { n, attenuation, bandwidth, samplerate } =
+   (pp = win.kaiserParameters (ByFrequency { attenuation, bandwidth, samplerate });
+    length = if pp.length % 2 == 0 then pp.length + 1 else pp.length fi;
+    kw = win.kaiser (pp with { length });
+    sw = win.sinc (n*2) length;
     println "kaiserSincFilterFor: filter length is \(vec.length sw)";
-    bf.multiply sw kw);
+
+    f = bf.multiply sw kw;
+    
+    println "sinc peak is \(bf.max sw) at index \(bf.maxindex sw)";
+    println "kaiser peak is \(bf.max kw) at index \(bf.maxindex kw)";
+    println "filter peak is \(bf.max f) at index \(bf.maxindex f)";
+
+    f);
 
 bandpassed f0 f1 attenuation bandwidth s =
    (rate = s.sampleRate;
@@ -517,6 +525,7 @@
         nzc = 13;
         attenuation = 80;
         filter = kaiserSincWindow nzc factor attenuation;
+//        println "interpolated: filter length \(vec.length filter)";
         out = delayedBy (- (nzc * factor))
            (convolvedWith [Framesize 1024]
                (mat.newMatrix (RowMajor ()) (map \filter [1..s.channels]))
@@ -595,7 +604,7 @@
 
     //!!! ponder
 //    filt = kaiserSincFilterFor { n, alpha = 80, bandwidth = 10, sampleRate = lower };
-    filt = kaiserSincFilterFor { n, alpha = 80, bandwidth = 1, sampleRate = n };
+    filt = kaiserSincFilterFor { n, attenuation = 80, bandwidth = 1, samplerate = higher };
 
 //\() (    plot.plot [Vector filt] );
 
@@ -646,7 +655,7 @@
 
     reconstructOne ch =
        (var out = 0;
-        var series = mat.getRow ch buffer;
+        var series = mat.getRow ch buffer; //!!! check this has length flen
         for [0..(vec.length series)-1] do i:
             x = vec.at series i;
             if x != 0 then
@@ -658,6 +667,8 @@
     readOne () =
        (result = mat.newColumnVector
            (vec.fromList (map reconstructOne [0 .. s.channels-1]));
+//        println "result = \(result)";
+//print ".";
         m = sourceRate / g;
         buffer := mat.concat (Horizontal ())
            [mat.columnSlice buffer m (flen - m),