aboutsummaryrefslogtreecommitdiffstats
diff options
context:
space:
mode:
authorChris Robinson <[email protected]>2023-11-07 15:04:51 -0800
committerChris Robinson <[email protected]>2023-11-07 15:04:51 -0800
commitd44112514d84614af6c16ec48d9450706da3b335 (patch)
tree9b4fa89c1918b6523b3c85aacb5fec92d943498d
parent34cde4b55f87d775e79c1f4cbd6b14d129b1a949 (diff)
Precompute a value used multiple times
-rw-r--r--common/polyphase_resampler.cpp28
1 files changed, 9 insertions, 19 deletions
diff --git a/common/polyphase_resampler.cpp b/common/polyphase_resampler.cpp
index 14f7e40d..23102b9d 100644
--- a/common/polyphase_resampler.cpp
+++ b/common/polyphase_resampler.cpp
@@ -3,6 +3,7 @@
#include <algorithm>
#include <cmath>
+#include <numeric>
#include "alnumbers.h"
#include "opthelpers.h"
@@ -67,23 +68,11 @@ constexpr double BesselI_0(const double x)
*
* k = 2 i / M - 1, where 0 <= i <= M.
*/
-double Kaiser(const double b, const double k)
+double Kaiser(const double beta, const double k, const double besseli_0_beta)
{
if(!(k >= -1.0 && k <= 1.0))
return 0.0;
- return BesselI_0(b * std::sqrt(1.0 - k*k)) / BesselI_0(b);
-}
-
-// Calculates the greatest common divisor of a and b.
-constexpr uint Gcd(uint x, uint y)
-{
- while(y > 0)
- {
- const uint z{y};
- y = x % y;
- x = z;
- }
- return x;
+ return BesselI_0(beta * std::sqrt(1.0 - k*k)) / besseli_0_beta;
}
/* Calculates the size (order) of the Kaiser window. Rejection is in dB and
@@ -124,11 +113,11 @@ constexpr double CalcKaiserBeta(const double rejection)
* p -- gain compensation factor when sampling
* f_t -- normalized center frequency (or cutoff; 0.5 is nyquist)
*/
-double SincFilter(const uint l, const double b, const double gain, const double cutoff,
- const uint i)
+double SincFilter(const uint l, const double beta, const double besseli_0_beta, const double gain,
+ const double cutoff, const uint i)
{
const double x{static_cast<double>(i) - l};
- return Kaiser(b, x / l) * 2.0 * gain * cutoff * Sinc(2.0 * cutoff * x);
+ return Kaiser(beta, x/l, besseli_0_beta) * 2.0 * gain * cutoff * Sinc(2.0 * cutoff * x);
}
} // namespace
@@ -137,7 +126,7 @@ double SincFilter(const uint l, const double b, const double gain, const double
// that's used to cut frequencies above the destination nyquist.
void PPhaseResampler::init(const uint srcRate, const uint dstRate)
{
- const uint gcd{Gcd(srcRate, dstRate)};
+ const uint gcd{std::gcd(srcRate, dstRate)};
mP = dstRate / gcd;
mQ = srcRate / gcd;
@@ -160,11 +149,12 @@ void PPhaseResampler::init(const uint srcRate, const uint dstRate)
// calculating the left offset to avoid increasing the transition width.
const uint l{(CalcKaiserOrder(180.0, width)+1) / 2};
const double beta{CalcKaiserBeta(180.0)};
+ const double besseli_0_beta{BesselI_0(beta)};
mM = l*2 + 1;
mL = l;
mF.resize(mM);
for(uint i{0};i < mM;i++)
- mF[i] = SincFilter(l, beta, mP, cutoff, i);
+ mF[i] = SincFilter(l, beta, besseli_0_beta, mP, cutoff, i);
}
// Perform the upsample-filter-downsample resampling operation using a