made window functions more flexible

This commit is contained in:
Ahmet Inan 2018-03-05 17:44:07 +01:00
commit 375d8dea41

118
window.hh
View file

@ -12,101 +12,95 @@ Copyright 2018 Ahmet Inan <inan@aicodix.de>
namespace DSP {
template <typename TYPE>
struct Window
{
virtual TYPE operator () (int, int) = 0;
};
template <int TAPS, typename TYPE>
class Rect
class Taps
{
TYPE w[TAPS];
public:
Rect()
Taps(Window<TYPE> *func)
{
for (int n = 0; n < TAPS; ++n)
w[n] = TYPE(1);
w[n] = (*func)(n, TAPS);
}
inline TYPE operator () (int n) { return n >= 0 && n < TAPS ? w[n] : 0; }
inline operator const TYPE * () const { return w; }
};
template <int TAPS, typename TYPE>
class Hann
template <typename TYPE>
struct Rect : public Window<TYPE>
{
TYPE w[TAPS];
public:
Hann()
{
for (int n = 0; n < TAPS; ++n)
w[n] = TYPE(0.5) * (TYPE(1) - std::cos(Const<TYPE>::TwoPi() * TYPE(n) / TYPE(TAPS - 1)));
}
inline TYPE operator () (int n) { return n >= 0 && n < TAPS ? w[n] : 0; }
inline operator const TYPE * () const { return w; }
TYPE operator () (int n, int N) { return n >= 0 && n < N ? 1 : 0; }
};
template <int TAPS, typename TYPE>
class Hamming
template <typename TYPE>
struct Hann : public Window<TYPE>
{
TYPE w[TAPS];
public:
Hamming()
TYPE operator () (int n, int N)
{
for (int n = 0; n < TAPS; ++n)
w[n] = TYPE(0.54) - TYPE(0.46) * std::cos(Const<TYPE>::TwoPi() * TYPE(n) / TYPE(TAPS - 1));
return TYPE(0.5) * (TYPE(1) - std::cos(Const<TYPE>::TwoPi() * TYPE(n) / TYPE(N - 1)));
}
inline TYPE operator () (int n) { return n >= 0 && n < TAPS ? w[n] : 0; }
inline operator const TYPE * () const { return w; }
};
template <int TAPS, typename TYPE>
class Lanczos
template <typename TYPE>
struct Hamming : public Window<TYPE>
{
TYPE w[TAPS];
TYPE sinc(TYPE x)
TYPE operator () (int n, int N)
{
return TYPE(0.54) - TYPE(0.46) * std::cos(Const<TYPE>::TwoPi() * TYPE(n) / TYPE(N - 1));
}
};
template <typename TYPE>
class Lanczos : public Window<TYPE>
{
static TYPE sinc(TYPE x)
{
return TYPE(0) == x ? TYPE(1) : std::sin(Const<TYPE>::Pi() * x) / (Const<TYPE>::Pi() * x);
}
public:
Lanczos()
TYPE operator () (int n, int N)
{
for (int n = 0; n < TAPS; ++n)
w[n] = sinc(TYPE(2 * n) / TYPE(TAPS - 1) - TYPE(1));
return sinc(TYPE(2 * n) / TYPE(N - 1) - TYPE(1));
}
inline TYPE operator () (int n) { return n >= 0 && n < TAPS ? w[n] : 0; }
inline operator const TYPE * () const { return w; }
};
template <int TAPS, typename TYPE>
class Blackman
template <typename TYPE>
class Blackman : public Window<TYPE>
{
TYPE w[TAPS];
TYPE a0, a1, a2;
public:
Blackman(TYPE a0, TYPE a1, TYPE a2)
{
for (int n = 0; n < TAPS; ++n)
w[n] = a0 - a1 * std::cos(Const<TYPE>::TwoPi() * TYPE(n) / TYPE(TAPS - 1)) + a2 * std::cos(Const<TYPE>::FourPi() * TYPE(n) / TYPE(TAPS - 1));
}
Blackman(TYPE a0, TYPE a1, TYPE a2) : a0(a0), a1(a1), a2(a2) {}
Blackman(TYPE a) : Blackman((TYPE(1) - a) / TYPE(2), TYPE(0.5), a / TYPE(2)) {}
// "exact Blackman"
Blackman() : Blackman(TYPE(7938) / TYPE(18608), TYPE(9240) / TYPE(18608), TYPE(1430) / TYPE(18608)) {}
inline TYPE operator () (int n) { return n >= 0 && n < TAPS ? w[n] : 0; }
inline operator const TYPE * () const { return w; }
};
template <int TAPS, typename TYPE>
class Gauss
{
TYPE w[TAPS];
public:
Gauss(TYPE o)
TYPE operator () (int n, int N)
{
for (int n = 0; n < TAPS; ++n)
w[n] = std::exp(- TYPE(0.5) * std::pow((TYPE(n) - TYPE(TAPS - 1) / TYPE(2)) / (o * TYPE(TAPS - 1) / TYPE(2)), TYPE(2)));
return a0 - a1 * std::cos(Const<TYPE>::TwoPi() * TYPE(n) / TYPE(N - 1)) + a2 * std::cos(Const<TYPE>::FourPi() * TYPE(n) / TYPE(N - 1));
}
inline TYPE operator () (int n) { return n >= 0 && n < TAPS ? w[n] : 0; }
inline operator const TYPE * () const { return w; }
};
template <int TAPS, typename TYPE>
class Kaiser
template <typename TYPE>
class Gauss : public Window<TYPE>
{
TYPE w[TAPS];
TYPE o;
public:
Gauss(TYPE o) : o(o) {}
TYPE operator () (int n, int N)
{
return std::exp(- TYPE(0.5) * std::pow((TYPE(n) - TYPE(N - 1) / TYPE(2)) / (o * TYPE(N - 1) / TYPE(2)), TYPE(2)));
}
};
template <typename TYPE>
class Kaiser : public Window<TYPE>
{
TYPE a;
/*
i0() implements the zero-th order modified Bessel function of the first kind:
https://en.wikipedia.org/wiki/Bessel_function#Modified_Bessel_functions:_I%CE%B1,_K%CE%B1
@ -115,7 +109,7 @@ class Kaiser
We obviously can't use the factorial here, so let's get rid of it:
$= 1 + \left(\frac{x}{2 \cdot 1}\right)^2 + \left(\frac{x}{2 \cdot 1}\cdot \frac{x}{2 \cdot 2}\right)^2 + \left(\frac{x}{2 \cdot 1}\cdot \frac{x}{2 \cdot 2}\cdot \frac{x}{2 \cdot 3}\right)^2 + .. = 1 + \sum_{m=1}^\infty \left(\prod_{n=1}^m \frac{x}{2n}\right)^2$
*/
TYPE i0(TYPE x)
static TYPE i0(TYPE x)
{
Kahan<TYPE> sum(1.0);
TYPE val = 1.0;
@ -130,13 +124,11 @@ class Kaiser
return sum();
}
public:
Kaiser(TYPE a)
Kaiser(TYPE a) : a(a) {}
TYPE operator () (int n, int N)
{
for (int n = 0; n < TAPS; ++n)
w[n] = i0(Const<TYPE>::Pi() * a * std::sqrt(TYPE(1) - std::pow(TYPE(2 * n) / TYPE(TAPS - 1) - TYPE(1), TYPE(2)))) / i0(Const<TYPE>::Pi() * a);
return i0(Const<TYPE>::Pi() * a * std::sqrt(TYPE(1) - std::pow(TYPE(2 * n) / TYPE(N - 1) - TYPE(1), TYPE(2)))) / i0(Const<TYPE>::Pi() * a);
}
inline TYPE operator () (int n) { return n >= 0 && n < TAPS ? w[n] : 0; }
inline operator const TYPE * () const { return w; }
};
}