Chamberlin state variable filter transfer function

Hi,
I am comparing coefficient quantization errors for second order IIR filter structures (for fixed point implementation). I have found the Chamberlin state variable structure to perform well at low frequencies and have seen in the article "The modified Chamberlin and Zölzer filter structures" that the quantized pole distribution is dense in the low frequency area of the unit circle in the z-plane.

I have implemented the filter (in software) based on the diagram in:
http://www.earlevel.com/Digital%20Audio/StateVar.html

Now I want to derive the transfer function from the diagram and find out how f and q (see the link) are related to the b-coefficients in:

H(z) a0 + a1*z^-1 + a2*z^-2
= ------------------------
1 - b1*z^-1 - b2*z^-2

Then I will find the relation to the poles and make my own quantized pole distribution plot.
From the diagram I can see that the following simple implementation can be used (working matlab code, only HP filteret samples are saved):
low = 0;
high = 0;
band = 0;
y = zeros(length(x),1);

for i = 1:length(x)
low = low + f * band;
high = x(i) - low - q * band;
band = f * high + band;
y(i) = high;
end

I have then tried rewriting it to something like this:

ylp(n) = ylp(n-1)+f*ybp(n-1)
yhp(n) = x(n) - ylp(n)-q*ybp(n-1)
ybp(n) = f*yhp(n)+ybp(n-1)

but the problem is I now have three equations with more than three unknowns and also I need to get some x[n-1] etc. into the equation. I have tried doing some substitution and for instance turned up with:

ybp(n) = -f*ylp(n-1)+f*x(n)-(f^2+f*q-1)*ybp(n-1)

Then by looking at the diagram I tried to express ylp(n-1) by some by some x(n-...), f and q's but as it is recursive I can't (don't know how to) do that.
There must be some mathematical trick or sequency knowledge I should know?

Best regards,
Kim Start a New Thread

Kim-

Chamberlin's book is from 1985, and addresses 16-bit microprocessors. So I have to
ask what is your purpose... Are you limited to 16-bit precision? How low are your
required filter frequencies? If your sampling rate is 48 kHz then are your filter
frequencies less than 100 Hz?

If so, then trying to "fine tune" quantized coefficients will make little difference
and you will still face instability. Only higher precision yields stability. On a
16-bit uP or DSP this typically means implementing double or even quad precision
(usually in optimized asm lang). As one example, on a TI C55x DSP, we use 64-bit
precision for high accuracy IIR filters with fc less than 100 Hz (48 kHz sampling
rate).

I might suggest that unless you're incorporating into your MATLAB calculations very
precise simulation of your final device quantization characteristics, then the MATLAB
calcs won't do much good.

-Jeff

> I am comparing coefficient quantization errors for second order IIR filter
> structures (for fixed point implementation). I have found the Chamberlin
> state variable structure to perform well at low frequencies and have seen
> in the article "The modified Chamberlin and Zölzer filter structures"
> that the quantized pole distribution is dense in the low frequency area
> of the unit circle in the z-plane.
>
> I have implemented the filter (in software) based on the diagram in:
> http://www.earlevel.com/Digital%20Audio/StateVar.html
>
> Now I want to derive the transfer function from the diagram and find out
> how f and q (see the link) are related to the b-coefficients in:
>
> H(z) a0 + a1*z^-1 + a2*z^-2
> = ------------------------
> 1 - b1*z^-1 - b2*z^-2
>
> Then I will find the relation to the poles and make my own quantized pole
> distribution plot.
> From the diagram I can see that the following simple implementation can
> be used (working matlab code, only HP filteret samples are saved):
> low = 0;
> high = 0;
> band = 0;
> y = zeros(length(x),1);
>
> for i = 1:length(x)
> low = low + f * band;
> high = x(i) - low - q * band;
> band = f * high + band;
> y(i) = high;
> end
>
> I have then tried rewriting it to something like this:
>
> ylp(n) = ylp(n-1)+f*ybp(n-1)
> yhp(n) = x(n) - ylp(n)-q*ybp(n-1)
> ybp(n) = f*yhp(n)+ybp(n-1)
>
> but the problem is I now have three equations with more than three unknowns and also I need to get some x[n-1] etc. into the equation. I have tried doing some substitution and for instance turned up with:
>
> ybp(n) = -f*ylp(n-1)+f*x(n)-(f^2+f*q-1)*ybp(n-1)
>
> Then by looking at the diagram I tried to express ylp(n-1) by some by some x(n-...), f and q's but as it is recursive I can't (don't know how to) do that.
> There must be some mathematical trick or sequency knowledge I should know?
>
> Best regards,
> Kim Start a New Thread

Hi,

I have sent Jeff an email describing what I need it for and my problem has now been solved. It turned out that there was an easy solution to my problem - I should have done the calculations in the z domain then you get three equations and three unknowns and it is easy to rearrange the equations to find the A and B coefficients.
The answer:

http://groups.google.com/group/comp.dsp/browse_thread/thread/b0753153b2506a06/e268fe8f0b4b1deb?hle268fe8f0b4b1deb

/Kim Start a New Thread

Kim-

> I have sent Jeff an email describing what I need it for and my problem
> has now been solved. It turned out that there was an easy solution to
> my problem - I should have done the calculations in the z domain then
> you get three equations and three unknowns and it is easy to rearrange
> the equations to find the A and B coefficients.
> The answer:
>
> http://groups.google.com/group/comp.dsp/browse_thread/thread/b0753153b2506a06/e268fe8f0b4b1deb?hle268fe8f0b4b1deb

The approach sounds good for your application (20 Hz highpass, 44.1 kHz Fs). Our
applications are usually along the lines of speech compression and equalizers, and we
typically have to implement IIR filters with coefficients generated by a variety of
methods. So we need DSP libraries that "just work", no matter where the coefficients
come from and what fc/Fs.

In our low frequency cases (10 to 50 Hz), passband has to be flat immediately on both
sides of the cutoffs -- this seems to rule out Chamberlin. Was this not a problem in
your application?

-Jeff Start a New Thread

You might want to look into the work of Mitra, Dana Massie, and others
who have developed a 2nd order structure wherein an all pass filter is
embedded in an adder-subtractive setup. It turns out to be very
insensitive to coefficient errors, is good at low frequencies, and has
cofficients which are nearly uncoupled from each other. I've used it
for parametric EQ and it's good. If those references don't allow you to
find it, let me know and I'll dig up real citations,

Chris Moore
k...@am3d.com wrote:
> Hi,
>
> I have sent Jeff an email describing what I need it for and my problem
> has now been solved. It turned out that there was an easy solution to my
> problem - I should have done the calculations in the z domain then you
> get three equations and three unknowns and it is easy to rearrange the
> equations to find the A and B coefficients.
> The answer:
>
> http://groups. google.com/ group/comp. dsp/browse_ thread/thread/
> b0753153b2506a06 /e268fe8f0b4b1de b?hle268fe8f 0b4b1deb
> /Kim
> Start a New Thread