LPF
Overview
LPF
implements a fairly common 2-pole resonant
low-pass IIR filter.
Tangled Files
lpf.c
and lpf.h
.
#ifndef SK_LPF_H
#define SK_LPF_H
#ifndef SKFLT
#define SKFLT float
#endif
<<typedefs>>
<<funcdefs>>
#ifdef SK_LPF_PRIV
<<structs>>
#endif
#endif
#include <math.h>
#define SK_LPF_PRIV
#include "lpf.h"
#ifndef M_PI
#define M_PI 3.14159265358979323846
#endif
<<funcs>>
Struct and Initialization
The struct for LPF is called sk_lpf
.
typedef struct sk_lpf sk_lpf;
struct sk_lpf {
<<sk_lpf>>
};
LPF is initialized with sk_lpf_init
.
void sk_lpf_init(sk_lpf *lpf, int sr);
void sk_lpf_init(sk_lpf *lpf, int sr)
{
<<init>>
}
Being a 2-pole filter requires 2 samples of memory, and 3 filter coefficients: 2 beta coefficients (the recursive part of the filter) and one alpha coefficient (applied to the input).
SKFLT y[2];
SKFLT a0, b1, b2;
lpf->y[0] = 0.0;
lpf->y[1] = 0.0;
lpf->a0 = 0.0;
lpf->b1 = 0.0;
lpf->b2 = 0.0;
A update
flag is used to implement parameter caching. If either
cutoff or frequency is set, it will cause the filter
coefficients to update themselves.
int update;
lpf->update = 1;
The constant tpidsr
is a classic variable shorthand for
"two pi divided by the sampling rate", and is used to shave
off some arithmetic operations when computing filter
coefficients.
SKFLT tpidsr;
lpf->tpidsr = 2.0 * M_PI / (SKFLT)sr;
Parameters
Cutoff Frequency
The cutoff frequency is set with sk_lpf_cutoff
, in
units of Hz.
void sk_lpf_cutoff(sk_lpf *lpf, SKFLT cutoff);
The update flag is set to indicate that the next
call to sk_lpf_tick
needs to update coefficients
before computing.
void sk_lpf_cutoff(sk_lpf *lpf, SKFLT cutoff)
{
lpf->cutoff = cutoff;
lpf->update = 1;
}
Cutoff the variable is defined and set at init time below:
SKFLT cutoff;
sk_lpf_cutoff(lpf, 1000); /* arbitrary default */
Q Factor
The Q factor is set with sk_lpf_q
.
The Q factor controls how much resonance there is. Should be a positive non-zero value. A value of 1.0 is passive, greater than 1.0 will cause resonance. Less than 1? Not sure.
void sk_lpf_q(sk_lpf *lpf, SKFLT q);
Similar to sk_lpf_cutoff
, sk_lpf_q
will set the
update
flag, which will update filter coefficients
next time the tick function sk_lpf_tick
is computed.
Note that Q cannot be exactly 0 since it is used in a division operation, so it is set to be a reasonably small value.
void sk_lpf_q(sk_lpf *lpf, SKFLT q)
{
if (q < 0.001) q = 0.001;
lpf->q = q;
lpf->update = 1;
}
The Q variable is set and initialized below:
SKFLT q;
sk_lpf_q(lpf, 1.0); /* arbitrary default */
Computation
A single sample of audio is computed with
sk_lpf_tick
.
SKFLT sk_lpf_tick(sk_lpf *lpf, SKFLT in);
SKFLT sk_lpf_tick(sk_lpf *lpf, SKFLT in)
{
SKFLT out;
SKFLT y0;
out = 0.0;
<<check_and_update_coefficients>>
<<compute_difference_equation>>
return out;
}
Before a filter is computed, the update
flag is
checked to see if it is set. If so, the filter
coefficients must be updated.
To compute new coefficients, frequency must be converted
from cycles-per-second to radians-per-sample. Multiplying
by
gives radians, then
dividing by the sampling rate (or multiplying by the
inverse, big T
, converts
from seconds to radians. In our function, this gets
smooshed together in a constant called tpidsr
.
The intermediate variables C
and D
are computed next,
followed by the alpha and beta filter coefficients. These
look like the remains from some bilinear transform from
an S-plane filter, but I'm not sure.
The update flag is then reset back to zero.
if (lpf->update) {
SKFLT C, D;
SKFLT freq;
SKFLT qres;
qres = (1.0 / lpf->q);
if (qres < 0.001) qres = 0.001;
/* convert to radians/sample */
freq = lpf->cutoff * lpf->tpidsr;
/* intermediates */
D = tan(freq * qres * 0.5);
C = (1.0 - D) / (1.0 + D);
lpf->b1 = (1.0 + C) * cos(freq);
lpf->b2 = -C;
lpf->a0 = (1.0 + C - lpf->b1) * 0.25;
lpf->update = 0;
}
With all the up-to-date coefficients, computing filter is a matter of computing the difference equation and updating the filter memory.
The output itself seems to be effectively boosting the filter memory. I am not sure why it is doing that.
y0 = lpf->a0*in + lpf->b1*lpf->y[0] + lpf->b2*lpf->y[1];
out = y0 + 2.0*lpf->y[0] + lpf->y[1];
lpf->y[1] = lpf->y[0];
lpf->y[0] = y0;