Hello all,

I've been experimenting with filters (again) and have learned how to use
Matlab.  The following two Matlab functions work fairly well for evaluating
the frequency response a tuned resonant low pass filter.

The variable "band" controls the band-pass nature of the filter to some
extent.  A band of 1 has a pass band gain of 1.  Values of band less than 1
(but greater than 0) have reduced gain in the pass band.  The Moog VCF
behaves somewhat like this.

x is the angle, fco is the frequency cut off, rez controls the amount of
resonance (rez should be >~.5) and sr is the sample rate.

The 1+exp(fco/11000) term came about because of a problem with resonance.
The peak resonance stays pretty flat at low frequencies with the value of
resonance approximating the gain of the peak.  Above ~3kHz the resonance
value starts to take off and attenuation of the high frequencies lessens.
The 1+exp(fco/11000) fixes the climb in resonance pretty well.  The reduced
attenuation of high frequencies was not fixed, but who can hear up there
anyway?  Over sampling would help to solve these problems.

As far as I can tell it does not self oscillate.  I think there needs to be
a zero in the denominator to do that.  The tuning is only approximate but I
think it is close.  (In other words fco is in Hz).

I don't have time to post a Csound version right now but it should be
pretty easy to implement using John Fitch's nlfilter.  I haven't been able
to figure out how to generate an equivalent high pass version yet.

function [omega] = rezfilter(x, fco, rez, band, sr)
fqc=66/441*sr;
c=fqc/fco;
rez2=rez/(1+exp(fco/11000));
a0=c/rez2+c^2;
a1=-(1+rez2*(1-band))+c/rez2+2*c^2;
a2=c^2;
omega=abs(1./(a0-a1*exp(-i*x)+a2*exp(-2*i*x)));

function [] = rzgraph(rez, band)

x0=-4:.05:.56;
x=10.^x0;
y0=rezfilter(x, 30, rez, band, 44100);
y1=rezfilter(x, 100, rez, band, 44100);
y2=rezfilter(x, 300, rez, band, 44100);
y3=rezfilter(x, 1000, rez, band, 44100);
y4=rezfilter(x, 3000, rez, band, 44100);
y5=rezfilter(x, 10000, rez, band, 44100);
y6=rezfilter(x, 10, rez, band, 44100);
semilogx(x, y0,'b',x, y1,'b',x, y2,'b',x, y3,'b',x, y4,'b',x, y5,'b',x,
y6,'b');