c@119: /* -*- c-basic-offset: 4 indent-tabs-mode: nil -*-  vi:set ts=8 sts=4 sw=4: */
c@119: /*
c@119:     Constant-Q library
c@119:     Copyright (c) 2013-2014 Queen Mary, University of London
c@119: 
c@119:     Permission is hereby granted, free of charge, to any person
c@119:     obtaining a copy of this software and associated documentation
c@119:     files (the "Software"), to deal in the Software without
c@119:     restriction, including without limitation the rights to use, copy,
c@119:     modify, merge, publish, distribute, sublicense, and/or sell copies
c@119:     of the Software, and to permit persons to whom the Software is
c@119:     furnished to do so, subject to the following conditions:
c@119: 
c@119:     The above copyright notice and this permission notice shall be
c@119:     included in all copies or substantial portions of the Software.
c@119: 
c@119:     THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
c@119:     EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF
c@119:     MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND
c@119:     NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHOR BE LIABLE FOR ANY
c@119:     CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF
c@119:     CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION
c@119:     WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
c@119: 
c@119:     Except as contained in this notice, the names of the Centre for
c@119:     Digital Music; Queen Mary, University of London; and Chris Cannam
c@119:     shall not be used in advertising or otherwise to promote the sale,
c@119:     use or other dealings in this Software without prior written
c@119:     authorization.
c@119: */
c@119: 
c@119: #include "KaiserWindow.h"
c@119: 
c@122: #include "MathUtilities.h"
c@119: 
c@119: KaiserWindow::Parameters
c@119: KaiserWindow::parametersForTransitionWidth(double attenuation,
c@119: 					   double transition)
c@119: {
c@119:     Parameters p;
c@119:     p.length = 1 + (attenuation > 21.0 ?
c@119: 		    ceil((attenuation - 7.95) / (2.285 * transition)) :
c@119: 		    ceil(5.79 / transition));
c@119:     p.beta = (attenuation > 50.0 ? 
c@119: 	      0.1102 * (attenuation - 8.7) :
c@119: 	      attenuation > 21.0 ? 
c@119: 	      0.5842 * pow(attenuation - 21.0, 0.4) + 0.07886 * (attenuation - 21.0) :
c@119: 	      0);
c@119:     return p;
c@119: }
c@119: 
c@119: static double besselTerm(double x, int i)
c@119: {
c@119:     if (i == 0) {
c@119: 	return 1;
c@119:     } else {
c@119: 	double f = MathUtilities::factorial(i);
c@119: 	return pow(x/2, i*2) / (f*f);
c@119:     }
c@119: }
c@119: 
c@119: static double bessel0(double x)
c@119: {
c@119:     double b = 0.0;
c@119:     for (int i = 0; i < 20; ++i) {
c@119: 	b += besselTerm(x, i);
c@119:     }
c@119:     return b;
c@119: }
c@119: 
c@119: void
c@119: KaiserWindow::init()
c@119: {
c@119:     double denominator = bessel0(m_beta);
c@119:     bool even = (m_length % 2 == 0);
c@119:     for (int i = 0; i < (even ? m_length/2 : (m_length+1)/2); ++i) {
c@119: 	double k = double(2*i) / double(m_length-1) - 1.0;
c@119: 	m_window.push_back(bessel0(m_beta * sqrt(1.0 - k*k)) / denominator);
c@119:     }
c@119:     for (int i = 0; i < (even ? m_length/2 : (m_length-1)/2); ++i) {
c@119:         m_window.push_back(m_window[int(m_length/2) - i - 1]);
c@119:     }
c@119: }