special.h

Design notes

Types 1

Functions 11

• fninline double besselI0(double x)
  Modified Bessel function of the first kind, order zero.
  More

  The Kaiser window is defined in terms of it and almost nothing else here is, which is why it is easy to write twice.

  A series rather than a table: it converges in a couple of dozen terms for the arguments a window uses -- beta is rarely above ten -- and a table would be more code than it removes, with a resolution to argue about. Computed once per window, never per sample.

  special::besselI0 · special.h:35
• fninline double kaiser(double r, double beta)
  The Kaiser window's shape at a position, r running -1 to 1 across it.
  More

  The window BOTH callers actually wanted, rather than the Bessel function they each derived it from separately.

  special::kaiser · special.h:60
• fninline double ellipticK(double k)
  K(k), the complete elliptic integral of the first kind -- the quarter period -- by the arithmetic-geometric mean: K = pi / (2 AGM(1, k')). The AGM converges quadratically, so this is exact to the last digit in a handful of iterations. Infinite at k = 1.
  special::ellipticK · special.h:127
• fninline void ellipticIntegralK(double k, double& K, double& Kprime)
  K(k) and K'(k) = K(k'), together: the ratio K'/K is the imaginary period that the elliptic design's degree equation is written in.
  special::ellipticIntegralK · special.h:141
• fninline std::complex<double> cde(std::complex<double> u, double k)
  cd(u K, k): the Jacobi cd function, argument in quarter periods.
  special::cde · special.h:147
• fninline std::complex<double> sne(std::complex<double> u, double k)
  sn(u K, k): the Jacobi sn function, argument in quarter periods.
  special::sne · special.h:156
• fninline std::complex<double> asne(std::complex<double> w, double k)
  The inverse of sne: a u with sn(u K, k) = w.
  More

  Each Gauss step is inverted on the way DOWN the sequence -- the quadratic k(n) w^2 s^2 - (1 + k(n)) s + w = 0 solved for s, in the rationalised form that does not divide by w -- and at the bottom u = (2/pi) asin. sn is doubly periodic, so there are infinitely many answers; this is the one the principal square roots and arcsine lead to, which for real w in [-1, 1] is the real u in [-1, 1], and for real w in (1, 1/k] -- the stopband the design maps -- is 1 + i v with v in [0, K'/K].

  special::asne · special.h:175
• fninline float tanhApprox(float x)
  tanh, as a [9/8] Padé approximant.
  More

  The rational form keeps the S-curve exactly ODD, which a polynomial fit does not, and odd symmetry is what stops a waveshaper generating even harmonics it was never asked for. Beyond ±5 it clamps, where real tanh is already within 1e-4 of ±1.

  The famous short version, x(27+x²)/(27+9x²), is only good to 2e-2: fine for a soft clipper, nowhere near transparent. The measured error of THIS one is asserted in the test suite rather than claimed here.

  fastmath::tanhApprox · special.h:228
• fninline float exp2Approx(float x)
  2^x for x in [-1, 1], the range a per-sample pitch or gain modulation actually spans. Outside it, falls back to the exact function rather than returning nonsense.
  More

  Degree 7, because 2^x = e^(x ln2) and the series coefficients are (ln2)^k/k!: truncating at degree 4 (a first version here) left 1e-3 of error at the ends, which is audible as a stepped gain.

  fastmath::exp2Approx · special.h:242
• fninline float sinTurns(float x)
  sin(2*pi*x) with x in TURNS, for a phase that already runs 0..1.
  More

  An oscillator holds its phase as a fraction of a cycle, so the library call costs a multiply by 2*pi before it even starts. Taking turns removes that and lets the polynomial be written directly in the variable the caller has.

  Why it is worth having at all: a polyphonic synth runs this once per voice per sample. Eight voices at 48 kHz is 384,000 calls a second for one note held down, and std::sin is correctly rounded to the last bit -- an accuracy nobody can hear and everybody pays for.

  Degree 13, which is the Taylor series of sin(2*pi*x) truncated where its error falls below the noise floor of 24-bit audio. The measured error is asserted in the test suite rather than claimed here, because a stated accuracy nobody checks is a stated accuracy that drifts.

  Accurate over the WHOLE turn, not a quadrant: the wrap below folds any phase into [-0.5, 0.5) first, so a caller does not have to know the domain. An approximation used outside its range is a bug, not a speed-up.

  fastmath::sinTurns · special.h:271
• fninline float dbToGainApprox(float db)
  Decibels to linear gain, for control-rate use.
  fastmath::dbToGainApprox · special.h:292