From: Rob Gaddi on
On Mon, 15 Mar 2010 13:11:42 -0500
Vladimir Vassilevsky <nospam(a)nowhere.com> wrote:

>
>
> Rob Gaddi wrote:
>
> > Hey y'all --
> >
> > So I've been working with some filter topologies trying to do
> > programmable low-pass filtering over a 100,000:1 Fc range in a
> > single stage, no pre-decimation. (That's unfortunately
> > non-negotiable for reasons above my control).
> >
> > The second order sections that I'm using are basically derived as
> > digital equivalents of analog state variable filters. Flow diagram
> > is:
> >
> > x[n] -->--+-->--k1-->--{A}-->--+-->--k2-->--{A}-->---- y[n]
> > | | |
> > | | -z^-1
> > | | |
> > \----------<---------\------------<-------/
> >
> > Where each {A} is an accumulator, structured as:
> >
> > -->--+------>---------->--
> > | |
> > \--<--z^-1--<--/
> >
> > There are no registers in the forward path, so y[n] contains (among
> > other things) k1*k2*x[n]. The most important aspect, however, is
> > that the negative feedback around the overall loop is -1, which
> > means that the section has unity gain at DC. As compared to, for
> > instance, direct form sections, where the gains become ludicrous at
> > low cut-off frequencies.
> >
> > I'm trying to find some literature on this sucker in hopes that
> > someone smarter than I has some more insight into it, just to make
> > sure I know what I'm getting myself into. But I can't figure out
> > for the life of me what it's called. It's sorta-kinda-similairish
> > to a Chamberlin topology, but does all of the work in only two
> > multiplies rather than 3.
> >
> > In short, to name a thing is to have Google power over it. Anyone
> > know where to find this one?
>
> As mentioned by Tim, this structure is a lot like what is used in
> control systems.
>
> This won't save precision though, you will have to do 64-bit
> accumulation. What it could help with is the precision of the
> coefficients. You can probably get by 16 bits.
>

48 bits was the plan. I'm in an FPGA so there's no benefit to sticking
with powers of 2. But wider than is easy to work with, certainly.

> However, you can accomplish exactly the same result with the direct
> form feedback path.
>
> y[n] = x[n] + B1*y[n-1] + B2*y[n-2]
>
> B1 ~ 2
> B2 ~ 1
>
> Let B1 = 2 - alpha
> Let B2 = 1 - beta
>
> Then alpha and beta could be 16-bit numbers.
>
>
> Vladimir Vassilevsky
> DSP and Mixed Signal Design Consultant
> http://www.abvolt.com
>

Still not wide enough. I asked Octave to crank me out two 2nd order
Butterworths, one for a cutoff of 1/250000, the other for 100000/250000
(my extremes, 1 Hz and 100 kHz with a 500kHz sample rate). Taking a
look at what would be B1 there, for the low frequency case it's
1.7772e-05, and for the high it's 1.6305. That's a better than 90000:1
ratio. I think that on some level, if you really need to accurately
express math with a horrible wide dynamic range, you've got to have
numbers with a horrible wide dynamic range.

The other thing that was killing me on the direct form implementations,
however, was the gain. Once again relying on Octave, I get poles for
that 1 Hz filter at (1 - 8.8858e-06 +/- 8.8857e-06*j). That makes the
DC gain ~1e+10, which is better than 33 bits of gain to have to pull
back down each stage, and radically dependent on the coefficient
quantization.

--
Rob Gaddi, Highland Technology
Email address is currently out of order
From: Tim Wescott on
Rob Gaddi wrote:
> On Mon, 15 Mar 2010 13:11:42 -0500
> Vladimir Vassilevsky <nospam(a)nowhere.com> wrote:
>
>>
>> Rob Gaddi wrote:
>>
>>> Hey y'all --
>>>
>>> So I've been working with some filter topologies trying to do
>>> programmable low-pass filtering over a 100,000:1 Fc range in a
>>> single stage, no pre-decimation. (That's unfortunately
>>> non-negotiable for reasons above my control).
>>>
>>> The second order sections that I'm using are basically derived as
>>> digital equivalents of analog state variable filters. Flow diagram
>>> is:
>>>
>>> x[n] -->--+-->--k1-->--{A}-->--+-->--k2-->--{A}-->---- y[n]
>>> | | |
>>> | | -z^-1
>>> | | |
>>> \----------<---------\------------<-------/
>>>
>>> Where each {A} is an accumulator, structured as:
>>>
>>> -->--+------>---------->--
>>> | |
>>> \--<--z^-1--<--/
>>>
>>> There are no registers in the forward path, so y[n] contains (among
>>> other things) k1*k2*x[n]. The most important aspect, however, is
>>> that the negative feedback around the overall loop is -1, which
>>> means that the section has unity gain at DC. As compared to, for
>>> instance, direct form sections, where the gains become ludicrous at
>>> low cut-off frequencies.
>>>
>>> I'm trying to find some literature on this sucker in hopes that
>>> someone smarter than I has some more insight into it, just to make
>>> sure I know what I'm getting myself into. But I can't figure out
>>> for the life of me what it's called. It's sorta-kinda-similairish
>>> to a Chamberlin topology, but does all of the work in only two
>>> multiplies rather than 3.
>>>
>>> In short, to name a thing is to have Google power over it. Anyone
>>> know where to find this one?
>> As mentioned by Tim, this structure is a lot like what is used in
>> control systems.
>>
>> This won't save precision though, you will have to do 64-bit
>> accumulation. What it could help with is the precision of the
>> coefficients. You can probably get by 16 bits.
>>
>
> 48 bits was the plan. I'm in an FPGA so there's no benefit to sticking
> with powers of 2. But wider than is easy to work with, certainly.
>
>> However, you can accomplish exactly the same result with the direct
>> form feedback path.
>>
>> y[n] = x[n] + B1*y[n-1] + B2*y[n-2]
>>
>> B1 ~ 2
>> B2 ~ 1
>>
>> Let B1 = 2 - alpha
>> Let B2 = 1 - beta
>>
>> Then alpha and beta could be 16-bit numbers.
>>
>>
>> Vladimir Vassilevsky
>> DSP and Mixed Signal Design Consultant
>> http://www.abvolt.com
>>
>
> Still not wide enough. I asked Octave to crank me out two 2nd order
> Butterworths, one for a cutoff of 1/250000, the other for 100000/250000
> (my extremes, 1 Hz and 100 kHz with a 500kHz sample rate). Taking a
> look at what would be B1 there, for the low frequency case it's
> 1.7772e-05, and for the high it's 1.6305. That's a better than 90000:1
> ratio. I think that on some level, if you really need to accurately
> express math with a horrible wide dynamic range, you've got to have
> numbers with a horrible wide dynamic range.
>
> The other thing that was killing me on the direct form implementations,
> however, was the gain. Once again relying on Octave, I get poles for
> that 1 Hz filter at (1 - 8.8858e-06 +/- 8.8857e-06*j). That makes the
> DC gain ~1e+10, which is better than 33 bits of gain to have to pull
> back down each stage, and radically dependent on the coefficient
> quantization.
>
As long as your coefficients come directly from your denominator
polynomial, you're going to be radically dependent on their
quantization. The pole positions are just plain highly sensitive to the
last term in the denominator.

Unless you can come up with a topology (like the state-space one that I
suggested) that frees you from the (pole value)^2 nature of the trailing
coefficient, you're stuck with that.

--
Tim Wescott
Control system and signal processing consulting
www.wescottdesign.com
From: Tim Wescott on
Tim Wescott wrote:
> Rob Gaddi wrote:
>> Hey y'all --
>>
>> So I've been working with some filter topologies trying to do
>> programmable low-pass filtering over a 100,000:1 Fc range in a single
>> stage, no pre-decimation. (That's unfortunately non-negotiable for
>> reasons above my control).
>>
>> The second order sections that I'm using are basically derived as
>> digital equivalents of analog state variable filters. Flow diagram is:
>>
>> x[n] -->--+-->--k1-->--{A}-->--+-->--k2-->--{A}-->---- y[n]
>> | | |
>> | | -z^-1
>> | | |
>> \----------<---------\------------<-------/
>>
>> Where each {A} is an accumulator, structured as:
>>
>> -->--+------>---------->--
>> | |
>> \--<--z^-1--<--/
>>
>> There are no registers in the forward path, so y[n] contains (among
>> other things) k1*k2*x[n]. The most important aspect, however, is that
>> the negative feedback around the overall loop is -1, which means that
>> the section has unity gain at DC. As compared to, for instance, direct
>> form sections, where the gains become ludicrous at low cut-off
>> frequencies.
>>
>> I'm trying to find some literature on this sucker in hopes that someone
>> smarter than I has some more insight into it, just to make sure I know
>> what I'm getting myself into. But I can't figure out for the life of
>> me what it's called. It's sorta-kinda-similairish to a Chamberlin
>> topology, but does all of the work in only two multiplies rather than 3.
>>
>> In short, to name a thing is to have Google power over it. Anyone know
>> where to find this one?
>
> In the control systems world that'd be a controllability canonical form,
> only missing the forward terms from the input and the first state that
> are necessary to make the numerator come out right.
>
> Note: This is exactly the answer you asked for, yet it is useless for
> your purposes. I must be an engineer*.
>
> Dunno what -- if anything -- it's called in the DSP world.
>
> * http://www.anvari.org/fun/Job/Engineer_and_Manager.html
>
Oops -- it's _not_ any canonical form. My bad.

Further, if I've done my math right (check me), it's transfer function is

k1 k2
H(z) = ----------------------------
z^2 - (1 + k2 - k1 k2)z + k2

As your poles approach 1, this will become very sensitive to the value
of k2. Basically, the trailing coefficient needs to have about twice
the precision needed to distinguish your desired pole position from 1,
so if you're talking a pole pair at 1-10^-6, then you'd need 12
_decimal_ digits to do the job.

--
Tim Wescott
Control system and signal processing consulting
www.wescottdesign.com
From: Rob Gaddi on
On Mon, 15 Mar 2010 14:44:52 -0700
Tim Wescott <tim(a)seemywebsite.now> wrote:

> Tim Wescott wrote:
> > Rob Gaddi wrote:
> >> Hey y'all --
> >>
> >> So I've been working with some filter topologies trying to do
> >> programmable low-pass filtering over a 100,000:1 Fc range in a
> >> single stage, no pre-decimation. (That's unfortunately
> >> non-negotiable for reasons above my control).
> >>
> >> The second order sections that I'm using are basically derived as
> >> digital equivalents of analog state variable filters. Flow
> >> diagram is:
> >>
> >> x[n] -->--+-->--k1-->--{A}-->--+-->--k2-->--{A}-->---- y[n]
> >> | | |
> >> | | -z^-1
> >> | | |
> >> \----------<---------\------------<-------/
> >>
> >> Where each {A} is an accumulator, structured as:
> >>
> >> -->--+------>---------->--
> >> | |
> >> \--<--z^-1--<--/
> >>
> >> There are no registers in the forward path, so y[n] contains (among
> >> other things) k1*k2*x[n]. The most important aspect, however, is
> >> that the negative feedback around the overall loop is -1, which
> >> means that the section has unity gain at DC. As compared to, for
> >> instance, direct form sections, where the gains become ludicrous
> >> at low cut-off frequencies.
> >>
> >> I'm trying to find some literature on this sucker in hopes that
> >> someone smarter than I has some more insight into it, just to make
> >> sure I know what I'm getting myself into. But I can't figure out
> >> for the life of me what it's called. It's sorta-kinda-similairish
> >> to a Chamberlin topology, but does all of the work in only two
> >> multiplies rather than 3.
> >>
> >> In short, to name a thing is to have Google power over it. Anyone
> >> know where to find this one?
> >
> > In the control systems world that'd be a controllability canonical
> > form, only missing the forward terms from the input and the first
> > state that are necessary to make the numerator come out right.
> >
> > Note: This is exactly the answer you asked for, yet it is useless
> > for your purposes. I must be an engineer*.
> >
> > Dunno what -- if anything -- it's called in the DSP world.
> >
> > * http://www.anvari.org/fun/Job/Engineer_and_Manager.html
> >
> Oops -- it's _not_ any canonical form. My bad.
>
> Further, if I've done my math right (check me), it's transfer
> function is
>
> k1 k2
> H(z) = ----------------------------
> z^2 - (1 + k2 - k1 k2)z + k2
>
> As your poles approach 1, this will become very sensitive to the
> value of k2. Basically, the trailing coefficient needs to have about
> twice the precision needed to distinguish your desired pole position
> from 1, so if you're talking a pole pair at 1-10^-6, then you'd need
> 12 _decimal_ digits to do the job.
>

I've been getting
k1 k2
H(z) = ----------------------------------
z^2 - (k1 k2 + k2 - 2)z + 1 - k2

Which for z=1 collapses to (k1 k2)/(k1 k2) = 1. That makes sense
because if the filter is stable, then at DC the values going into both
accumulators must be 0, meaning that x[n] = y[n-1].

--
Rob Gaddi, Highland Technology
Email address is currently out of order
From: Tim Wescott on
Rob Gaddi wrote:
> On Mon, 15 Mar 2010 14:44:52 -0700
> Tim Wescott <tim(a)seemywebsite.now> wrote:
>
>> Tim Wescott wrote:
>>> Rob Gaddi wrote:
>>>> Hey y'all --
>>>>
>>>> So I've been working with some filter topologies trying to do
>>>> programmable low-pass filtering over a 100,000:1 Fc range in a
>>>> single stage, no pre-decimation. (That's unfortunately
>>>> non-negotiable for reasons above my control).
>>>>
>>>> The second order sections that I'm using are basically derived as
>>>> digital equivalents of analog state variable filters. Flow
>>>> diagram is:
>>>>
>>>> x[n] -->--+-->--k1-->--{A}-->--+-->--k2-->--{A}-->---- y[n]
>>>> | | |
>>>> | | -z^-1
>>>> | | |
>>>> \----------<---------\------------<-------/
>>>>
>>>> Where each {A} is an accumulator, structured as:
>>>>
>>>> -->--+------>---------->--
>>>> | |
>>>> \--<--z^-1--<--/
>>>>
>>>> There are no registers in the forward path, so y[n] contains (among
>>>> other things) k1*k2*x[n]. The most important aspect, however, is
>>>> that the negative feedback around the overall loop is -1, which
>>>> means that the section has unity gain at DC. As compared to, for
>>>> instance, direct form sections, where the gains become ludicrous
>>>> at low cut-off frequencies.
>>>>
>>>> I'm trying to find some literature on this sucker in hopes that
>>>> someone smarter than I has some more insight into it, just to make
>>>> sure I know what I'm getting myself into. But I can't figure out
>>>> for the life of me what it's called. It's sorta-kinda-similairish
>>>> to a Chamberlin topology, but does all of the work in only two
>>>> multiplies rather than 3.
>>>>
>>>> In short, to name a thing is to have Google power over it. Anyone
>>>> know where to find this one?
>>> In the control systems world that'd be a controllability canonical
>>> form, only missing the forward terms from the input and the first
>>> state that are necessary to make the numerator come out right.
>>>
>>> Note: This is exactly the answer you asked for, yet it is useless
>>> for your purposes. I must be an engineer*.
>>>
>>> Dunno what -- if anything -- it's called in the DSP world.
>>>
>>> * http://www.anvari.org/fun/Job/Engineer_and_Manager.html
>>>
>> Oops -- it's _not_ any canonical form. My bad.
>>
>> Further, if I've done my math right (check me), it's transfer
>> function is
>>
>> k1 k2
>> H(z) = ----------------------------
>> z^2 - (1 + k2 - k1 k2)z + k2
>>
>> As your poles approach 1, this will become very sensitive to the
>> value of k2. Basically, the trailing coefficient needs to have about
>> twice the precision needed to distinguish your desired pole position
>> from 1, so if you're talking a pole pair at 1-10^-6, then you'd need
>> 12 _decimal_ digits to do the job.
>>
>
> I've been getting
> k1 k2
> H(z) = ----------------------------------
> z^2 - (k1 k2 + k2 - 2)z + 1 - k2
>
> Which for z=1 collapses to (k1 k2)/(k1 k2) = 1. That makes sense
> because if the filter is stable, then at DC the values going into both
> accumulators must be 0, meaning that x[n] = y[n-1].
>
Oops again.

Well, if I haven't completely destroyed my credibility:

You're getting some gain in precision by the topology then: if you
consider that your zero-order coefficient needs to be b0 = d^2, where d
is the real part of your pole pair. Then if you take d = 1-ed, b0 = 1 -
2ed + ed^2. You _must_ capture that ed^2 part, but if you use some
floating point-ish thing to represent the number, then you just need
enough bits to encompass the difference between 2ed and 2ed - ed^2.

You still need oodles of precision in your states, though.

--
Tim Wescott
Control system and signal processing consulting
www.wescottdesign.com