diff --git a/README.md b/README.md index c88a39f..ddcc0a8 100644 --- a/README.md +++ b/README.md @@ -122,7 +122,10 @@ Some everyday helpers: * [probability density function](https://en.wikipedia.org/wiki/Probability_density_function) of the [normal distribution](https://en.wikipedia.org/wiki/Normal_distribution) * [sinc function](https://en.wikipedia.org/wiki/Sinc_function) * [delta function](https://en.wikipedia.org/wiki/Dirac_delta_function) -* [decibel function](https://en.wikipedia.org/wiki/Decibel) + +### [decibel.hh](decibel.hh) + +[Decibel](https://en.wikipedia.org/wiki/Decibel) calculation helpers. ### [resampler.hh](resampler.hh) diff --git a/decibel.hh b/decibel.hh new file mode 100644 index 0000000..256f041 --- /dev/null +++ b/decibel.hh @@ -0,0 +1,44 @@ +/* +Decibel calculation helpers + +Copyright 2018 Ahmet Inan +*/ + +#pragma once + +namespace DSP { + +template +TYPE decibel(TYPE v) +{ +#if 0 + return TYPE(10) * log10(v); +#else + static constexpr TYPE + // 2*1024*10/l(10) + scale = 8894.350989378597430295120259412072085389250678859091274024013479236L, + inv1 = TYPE(1) / TYPE(1), + inv3 = TYPE(1) / TYPE(3), + inv5 = TYPE(1) / TYPE(5), + inv7 = TYPE(1) / TYPE(7), + inv9 = TYPE(1) / TYPE(9), + inv11 = TYPE(1) / TYPE(11), + inv13 = TYPE(1) / TYPE(13); + TYPE x = v; + for (int i = 0; i < 10; ++i) + x = sqrt(x); + x = (x - TYPE(1)) / (x + TYPE(1)); + TYPE xx = x * x; + return scale * x * (xx * (xx * (xx * (xx * (xx * (xx * + inv13 + inv11) + inv9) + inv7) + inv5) + inv3) + inv1); +#endif +} + +template +TYPE idecibel(TYPE dB) +{ + return pow(TYPE(10), dB / TYPE(10)); +} + +} + diff --git a/utils.hh b/utils.hh index f5165bb..8a004e5 100644 --- a/utils.hh +++ b/utils.hh @@ -40,31 +40,5 @@ TYPE delta(TYPE x) return TYPE(0) == x ? TYPE(1) : TYPE(0); } -template -TYPE decibel(TYPE v) -{ -#if 0 - return TYPE(10) * log10(v); -#else - static constexpr TYPE - // 2*1024*10/l(10) - scale = 8894.350989378597430295120259412072085389250678859091274024013479236L, - inv1 = TYPE(1) / TYPE(1), - inv3 = TYPE(1) / TYPE(3), - inv5 = TYPE(1) / TYPE(5), - inv7 = TYPE(1) / TYPE(7), - inv9 = TYPE(1) / TYPE(9), - inv11 = TYPE(1) / TYPE(11), - inv13 = TYPE(1) / TYPE(13); - TYPE x = v; - for (int i = 0; i < 10; ++i) - x = sqrt(x); - x = (x - TYPE(1)) / (x + TYPE(1)); - TYPE xx = x * x; - return scale * x * (xx * (xx * (xx * (xx * (xx * (xx * - inv13 + inv11) + inv9) + inv7) + inv5) + inv3) + inv1); -#endif -} - }