//4th order Linkwitz-Riley filters
//
//The filter is unstable for fast automation changes in the lower frequency range.
//Parameter interpolation and/or oversampling should fix this.
//
//The sum of the Linkwitz-Riley (Butterworth squared) HP and LP outputs, will result an all-pass filter
//at Fc and flat magnitude response - close to ideal for crossovers.

desc: (LO=1+2  HI=3+4)  (Mono)  Read info...

slider1:200<20,22000,1>Frequency (Hz)

@init
pi=22/7;

@slider
fc=slider1;

wc=2*pi*fc;
wc2=wc*wc;
wc3=wc2*wc;
wc4=wc2*wc2;
k=wc/tan(pi*fc/srate);
k2=k*k;
k3=k2*k;
k4=k2*k2;
sqrt2=sqrt(2);
sq_tmp1=sqrt2*wc3*k;
sq_tmp2=sqrt2*wc*k3;
a_tmp=4*wc2*k2+2*sq_tmp1+k4+2*sq_tmp2+wc4;

b1=(4*(wc4+sq_tmp1-k4-sq_tmp2))/a_tmp;
b2=(6*wc4-8*wc2*k2+6*k4)/a_tmp;
b3=(4*(wc4-sq_tmp1+sq_tmp2-k4))/a_tmp;
b4=(k4-2*sq_tmp1+wc4-2*sq_tmp2+4*wc2*k2)/a_tmp;

//================================================
// low-pass
//================================================
a0=wc4/a_tmp;
a1=4*wc4/a_tmp;
a2=6*wc4/a_tmp;
a3=a1;
a4=a0;

//=====================================================
// high-pass
//=====================================================
a0=k4/a_tmp;
a1=-4*k4/a_tmp;
a2=6*k4/a_tmp;
a3=a1;
a4=a0;

@sample
tempx=(spl0+spl1)*0.5;

tempy=a0*tempx+a1*xm1+a2*xm2+a3*xm3+a4*xm4-b1*ym1-b2*ym2-b3*ym3-b4*ym4;
xm4=xm3;
xm3=xm2;
xm2=xm1;
xm1=tempx;
ym4=ym3;
ym3=ym2;
ym2=ym1;
ym1=tempy;

spl0=spl0-tempy;
spl1=spl1-tempy;
spl2=tempy;
spl3=tempy;