qm-dsp
1.8
|
00001 /* -*- c-basic-offset: 4 indent-tabs-mode: nil -*- vi:set ts=8 sts=4 sw=4: */ 00002 00003 /* 00004 QM DSP library 00005 Centre for Digital Music, Queen Mary, University of London. 00006 00007 This program is free software; you can redistribute it and/or 00008 modify it under the terms of the GNU General Public License as 00009 published by the Free Software Foundation; either version 2 of the 00010 License, or (at your option) any later version. See the file 00011 COPYING included with this distribution for more information. 00012 */ 00013 00014 #include "KaiserWindow.h" 00015 00016 #include "maths/MathUtilities.h" 00017 00018 KaiserWindow::Parameters 00019 KaiserWindow::parametersForTransitionWidth(double attenuation, 00020 double transition) 00021 { 00022 Parameters p; 00023 p.length = 1 + (attenuation > 21.0 ? 00024 ceil((attenuation - 7.95) / (2.285 * transition)) : 00025 ceil(5.79 / transition)); 00026 p.beta = (attenuation > 50.0 ? 00027 0.1102 * (attenuation - 8.7) : 00028 attenuation > 21.0 ? 00029 0.5842 * pow(attenuation - 21.0, 0.4) + 0.07886 * (attenuation - 21.0) : 00030 0); 00031 return p; 00032 } 00033 00034 static double besselTerm(double x, int i) 00035 { 00036 if (i == 0) { 00037 return 1; 00038 } else { 00039 double f = MathUtilities::factorial(i); 00040 return pow(x/2, i*2) / (f*f); 00041 } 00042 } 00043 00044 static double bessel0(double x) 00045 { 00046 double b = 0.0; 00047 for (int i = 0; i < 20; ++i) { 00048 b += besselTerm(x, i); 00049 } 00050 return b; 00051 } 00052 00053 void 00054 KaiserWindow::init() 00055 { 00056 double denominator = bessel0(m_beta); 00057 bool even = (m_length % 2 == 0); 00058 for (int i = 0; i < (even ? m_length/2 : (m_length+1)/2); ++i) { 00059 double k = double(2*i) / double(m_length-1) - 1.0; 00060 m_window.push_back(bessel0(m_beta * sqrt(1.0 - k*k)) / denominator); 00061 } 00062 for (int i = 0; i < (even ? m_length/2 : (m_length-1)/2); ++i) { 00063 m_window.push_back(m_window[int(m_length/2) - i - 1]); 00064 } 00065 }