lut.hpp

Go to the documentation of this file.

// sbl/dsp/lut.hpp — Lookup table interpolation functions
//
// These functions read LUTE-generated constexpr tables. The table pointer
// must point to an array of (Size + guard_points) entries.
//
// Phase convention: uint32_t phase where the full 32-bit range maps to one
// cycle of the table. Upper bits select the table index, lower bits provide
// sub-sample fractional precision.
//
// Inspired by MI stmlib/utils/dsp.h (MIT). Generalized for arbitrary
// power-of-2 table sizes via constexpr template parameters.
 
#ifndef SBL_DSP_MATH_LUT_HPP_
#define SBL_DSP_MATH_LUT_HPP_
 
#include <cstdint>
 
namespace sbl::dsp::lut {
 
// ─── Compile-time helpers ──────────────────────────────────────────────
 
namespace detail {
 
constexpr uint8_t log2_of(uint16_t n) {
    uint8_t r = 0;
    while (n > 1) { n >>= 1; ++r; }
    return r;
}
 
} // namespace detail
 
// ═══════════════════════════════════════════════════════════════════════
// Float lookup functions
//
// For float tables in [-1.0, 1.0]. Same phase convention (uint32_t full
// range = one cycle), dramatically simpler math — no Q16 fractions, no
// int64 intermediates.
// ═══════════════════════════════════════════════════════════════════════
 
// ─── Float linear interpolation ──────────────────────────────────────
 
/// @note All functions in sbl::dsp::lut are ISR-safe — bounded computation, no I/O.
 
template<uint16_t Size>
inline float lookup_linear(const float* table, uint32_t phase) {
    constexpr uint8_t bits = detail::log2_of(Size);
    static_assert(bits >= 1 && bits <= 16, "Size must be power of 2 in [2, 65536]");
 
    uint32_t idx = phase >> (32 - bits);
    float frac = static_cast<float>(phase & ((1u << (32 - bits)) - 1))
               * (1.0f / static_cast<float>(1u << (32 - bits)));
    return table[idx] + (table[idx + 1] - table[idx]) * frac;
}
 
// ─── Float cubic Hermite interpolation ───────────────────────────────
 
template<uint16_t Size>
inline float lookup_cubic(const float* table, uint32_t phase) {
    constexpr uint8_t bits = detail::log2_of(Size);
    static_assert(bits >= 2 && bits <= 16, "Size must be power of 2 in [4, 65536]");
 
    uint32_t idx = phase >> (32 - bits);
    float f = static_cast<float>(phase & ((1u << (32 - bits)) - 1))
            * (1.0f / static_cast<float>(1u << (32 - bits)));
 
    float xm1 = table[idx > 0 ? idx - 1 : Size - 1];
    float x0  = table[idx];
    float x1  = table[idx + 1];
    float x2  = table[idx + 2];
 
    // Catmull-Rom coefficients
    float c = (x1 - xm1) * 0.5f;
    float v = x0 - x1;
    float w = c + v;
    float a = w + v + (x2 - x0) * 0.5f;
    float b_neg = w + a;
 
    // Horner's method: ((a*f - b_neg)*f + c)*f + x0
    return ((a * f - b_neg) * f + c) * f + x0;
}
 
// ─── Float crossfade between two tables ──────────────────────────────
//
// balance: 0.0 = 100% table_a, 1.0 = 100% table_b.
 
template<uint16_t Size>
inline float crossfade(const float* table_a, const float* table_b,
                       uint32_t phase, float balance) {
    float a = lookup_linear<Size>(table_a, phase);
    float b = lookup_linear<Size>(table_b, phase);
    return a + (b - a) * balance;
}
 
} // namespace sbl::dsp::lut
 
#endif // SBL_DSP_MATH_LUT_HPP_