Fixed Point IIR implementation

  • Permalink
  • submit to reddit
  • Email
  • Follow


Hi pals,

I would like to implement a IIR Biquad filter using the fixed point 
arithmetics...
Hence to reduced the intermediate states I plan to use the following trick:
s(k) = x(k) -a1*s(k-1) -a2*s(k-2)
y(n) = b0*s(k) + b1*s(k-1) + b2*s(k-2)

The I can compute each y, saving only 2 states (s(k-1) a�nd s(k-2))...

BUT ... in order to scale my coefficients or input, I need to know what 
are the boundaries of the s(k) serie... in order to find a proper Fixed 
Point format for my coefficients..
So my questions are :

* for which conditions, s(k) is bounded
* If the s(k) is bounded, what are the boundaries?

Thanks in advance for your help (and sorry for my poor englsih)

Fred.
0
Reply fred_nach (6) 11/30/2005 11:35:55 AM

See related articles to this posting

Fred Nach wrote:
> Hi pals,
> 
> I would like to implement a IIR Biquad filter using the fixed point 
> arithmetics...
> Hence to reduced the intermediate states I plan to use the following trick:
> s(k) = x(k) -a1*s(k-1) -a2*s(k-2)
> y(n) = b0*s(k) + b1*s(k-1) + b2*s(k-2)
> 
> The I can compute each y, saving only 2 states (s(k-1) a�nd s(k-2))...
> 
> BUT ... in order to scale my coefficients or input, I need to know what 
> are the boundaries of the s(k) serie... in order to find a proper Fixed 
> Point format for my coefficients..
> So my questions are :
> 
> * for which conditions, s(k) is bounded
> * If the s(k) is bounded, what are the boundaries?
> 
> Thanks in advance for your help (and sorry for my poor englsih)

Your English is good enough, I think. The s[k] are data you supply. They 
can't be larger or smaller than the representation allows, but they may 
have smaller limits.

I see no feedback in the formulas you give. Is that correct?

Jerry
-- 
Engineering is the art of making what you want from things you can get.
�����������������������������������������������������������������������
0
Reply jya (12870) 11/30/2005 4:17:06 PM

Fred Nach wrote:
>
> I would like to implement a IIR Biquad filter using the fixed point
> arithmetics...
> Hence to reduced the intermediate states I plan to use the following tric=
k:
> s(k) =3D x(k) -a1*s(k-1) -a2*s(k-2)
> y(n) =3D b0*s(k) + b1*s(k-1) + b2*s(k-2)
>
> The I can compute each y, saving only 2 states (s(k-1) a=E9nd s(k-2))...

this is the Direct 2 Form (sometimes called the "Direct II Canonical
Form").  unrecommeded for floating-point (if the Q is quite high, you
end up subtracting numbers very close to each other and  losing
precision) and *highly* unrecommended for fixed-point (saturation
clipping will be your lot).

try the Direct 1 Form:

    y[n] =3D b0*x[n] + b1*x[n-1] + b2*x[n-2] - a1*y[n-1] - a2*y[n-2]

every product on the right side of the =3D sign is a double precision
fixed point and all of those products should be added together in
double precision.  this is trivial in the Mot 56K and the ADI SHArC and
maybe in the TI fixed-pointers.

> BUT ... in order to scale my coefficients or input, I need to know what
> are the boundaries of the s(k) serie... in order to find a proper Fixed
> Point format for my coefficients..

your coefs are defined strictly in terms of the frequency response you
want.  the gross range of a1 is -2 to +2 and a2 is -1 to +1 for any
stable biquad filter.  the b0, b1, b2 coefs can be nearly anything but
you will have to choose a range and possibly do some scaling using the
arithmetic shift operation on the result before saving to y[n].

> So my questions are :
>
> * for which conditions, s(k) is bounded
> * If the s(k) is bounded, what are the boundaries?

that is a purely artificial construct.  are your fixed-point signal
values considered 16 bit signed integers?  then it's  -32768 <=3D x[n] <
+32768 .  are they normalized?  then it's -1 <=3D x[n]< +1.  their range
is whatever you want them to be.

> Thanks in advance for your help (and sorry for my poor englsih)

as Jerry said, your English is fine.

r b-j

0
Reply rbj (4087) 11/30/2005 5:57:01 PM

> I see no feedback in the formulas you give. Is that correct?

The aX coeffecients should be the feedback terms and the bX
coefficients should be the feed forward terms (from drawing a direct
form II simulation diagram).  The s(k) terms are the inputs to the
delay buffers, which will be a function of the input, feed forward and
feedback terms.

I see the question he is asking, how to determine if these terms
overflow the fixed point number format in use, and how to determine a
proper scaling value to keep it from happening.  Unfortunately, I don't
know how to answer this question, but in the project I am presently
working on I have encountered the same problem, so I am hoping someone
has a good answer.  I got around the problem (? maybe ?) by using a
direct form I structure.

I have seen Matlab scripts that run an impulse response till it
stabilizes and then sum up the impulse responses and use the maximum
value as a scale factor, though I don't know if this is a correct
solution.

0
Reply no_spam_me2 (215) 11/30/2005 6:01:31 PM

Noway2 wrote:
>>I see no feedback in the formulas you give. Is that correct?
> 
> 
> The aX coeffecients should be the feedback terms and the bX
> coefficients should be the feed forward terms (from drawing a direct
> form II simulation diagram).  The s(k) terms are the inputs to the
> delay buffers, which will be a function of the input, feed forward and
> feedback terms.
> 
> I see the question he is asking, how to determine if these terms
> overflow the fixed point number format in use, and how to determine a
> proper scaling value to keep it from happening.  Unfortunately, I don't
> know how to answer this question, but in the project I am presently
> working on I have encountered the same problem, so I am hoping someone
> has a good answer.  I got around the problem (? maybe ?) by using a
> direct form I structure.
> 
> I have seen Matlab scripts that run an impulse response till it
> stabilizes and then sum up the impulse responses and use the maximum
> value as a scale factor, though I don't know if this is a correct
> solution.

I expect that the b's are the denominator terms. Input samples are s[n], 
output and feedback samples are y[n]. I see no coefficients operating on 
any y[n].

I see what he's asking too. His problem lies with internal states, and I 
don't know how to deal with that in the abstract. If there's a general 
case, I haven't seen it.

Jerry
-- 
Engineering is the art of making what you want from things you can get.
�����������������������������������������������������������������������
0
Reply jya (12870) 11/30/2005 6:36:47 PM

Jerry Avins wrote:
>
> I expect that the b's are the denominator terms.

no, Jerry, the a's are in the denominator.  they reversed that
convention a while back.  in the 70's it was a_k on top and b_k on
bottom, but, i think because they think that poles are more important,
most texts have it switched around now.

r b-j

0
Reply rbj (4087) 11/30/2005 7:32:10 PM

robert bristow-johnson wrote:
> Jerry Avins wrote:
> 
>>I expect that the b's are the denominator terms.
> 
> 
> no, Jerry, the a's are in the denominator.  they reversed that
> convention a while back.  in the 70's it was a_k on top and b_k on
> bottom, but, i think because they think that poles are more important,
> most texts have it switched around now.
> 
> r b-j

I knew that too, but old habits sometimes bubble up. The important point 
was that Fred N. has no y terms to the right of an equal sign, so I 
don't see any recursion. Somewhere along the line I wrote that. Twice, I 
think.

Jerry
-- 
Engineering is the art of making what you want from things you can get.
�����������������������������������������������������������������������
0
Reply jya (12870) 11/30/2005 9:21:45 PM

Jerry Avins wrote:
> robert bristow-johnson wrote:
> > Jerry Avins wrote:
> >
> >>I expect that the b's are the denominator terms.
> >
> > no, Jerry, the a's are in the denominator.
>
> I knew that too, but old habits sometimes bubble up. The important point
> was that Fred N. has no y terms to the right of an equal sign, so I
> don't see any recursion.

there is recursion with the intermediate s[k] signal:

Fred Nach wrote:

> I would like to implement a IIR Biquad filter using the fixed point
> arithmetics...
> Hence to reduced the intermediate states I plan to use the following trick:
> s(k) = x(k) -a1*s(k-1) -a2*s(k-2)
> y(n) = b0*s(k) + b1*s(k-1) + b2*s(k-2)

this is the canonical Direct Form 2.  poles before zeros.  problem is,
if there is a high Q filter, those poles will amplify the signal so
much that it will saturate s[k] before the zeros get to beat the signal
back down.  DF2 ain't particularly good, especially for fixed-point
arithmetic.

r b-j

0
Reply rbj (4087) 11/30/2005 10:27:56 PM

On Wed, 30 Nov 2005 11:01:31 -0800, Noway2 wrote:
> I have seen Matlab scripts that run an impulse response till it
> stabilizes and then sum up the impulse responses and use the maximum
> value as a scale factor, though I don't know if this is a correct
> solution.

It's pretty close.  This is the main reason why these "direct form 2"
implmentations are not preferred.  A direct form 1 implmentation will use
more memory but (a) memory is relatively cheap and (b) the values stored
will be as constrained as the input and output values.

-- 
Andrew

0
Reply andrew-newspost (542) 12/1/2005 12:10:25 AM

robert bristow-johnson wrote:
> Jerry Avins wrote:
> 
>>robert bristow-johnson wrote:
>>
>>>Jerry Avins wrote:
>>>
>>>
>>>>I expect that the b's are the denominator terms.
>>>
>>>no, Jerry, the a's are in the denominator.
>>
>>I knew that too, but old habits sometimes bubble up. The important point
>>was that Fred N. has no y terms to the right of an equal sign, so I
>>don't see any recursion.
> 
> 
> there is recursion with the intermediate s[k] signal:

OK. Back to the books with my thinking cap on.

Jerry
-- 
Engineering is the art of making what you want from things you can get.
�����������������������������������������������������������������������
0
Reply jya (12870) 12/1/2005 4:19:12 AM

Jerry Avins wrote:
> robert bristow-johnson wrote:
> > Jerry Avins wrote:
....
> >>I knew that too, but old habits sometimes bubble up. The important point
> >>was that Fred N. has no y terms to the right of an equal sign, so I
> >>don't see any recursion.
> >
> >
> > there is recursion with the intermediate s[k] signal:
>
> OK. Back to the books with my thinking cap on.

no need.  here's a simple way to look at it:

   s[k] =    x[k] - a1*s[k-1] - a2*s[k-2]
   y[k] = b0*s[k] + b1*s[k-1] + b2*s[k-2]

think of this as two filters in cascade.  the first with input x[k] and
output s[k].  the second with input s[k] and output y[k].

the first is a all-pole filter (actually there are 2 zeros at the
origin, but don't worry about them) and the second is a all zero
filter, an FIR (like an FIR, there are 2 poles at the origin, but the
other 2 zeros kill 'em).  the transfer function of the first:

    S(z)/X(z) = 1/(1 + a1*z^-1 + a2*z-2)

the transfer function of the second

    Y(z)/S(z) = b0 + b1*z^-1 + b2*z^-2

and the transfer function of the whole sheee-bang is

   H(z) = Y(z)/X(z) = ( Y(z)/S(z) ) * ( S(z)/X(z) )

and the reason it has only two states instead of four is that the
states for the first are identical to the states of the second (both
delays on s[k]) so they can be combined into a single delay line.

 r b-j

0
Reply rbj (4087) 12/1/2005 7:09:00 AM

Hi all !

Firstly thanks a lot for your fast and usefull help ;-)

It seems that my previous post was not clear enough so here are some 
extra details :

* The filter structure is indeed the From II that robert B-J has 
described very clearly.

* a1 and a2 are hence the DENOMINATORS coeff (I was born in the 80's !!)
and ... b1 and b2 are NUMERATOR (no b0 cause it is supposed scaled)

* The input signal is x(n) and it is a bounded signal : PCM 16 bits so 
it varies from -32768 to 32767

* The s(k) are the internal states that I save between each sample 
processing...

* The filter is indeed stable hence  |a1| <= 2 and |a2| <= 1

I think I understood that this implementation is not recommended because 
it compute first the feedback terms and leads to big variations...
But is there a way to retrieve the boudary of s(k) mathematicaly ?

I think that as long as a1 and a2 are given for a stable filter s(k) 
should be bounded...

For my program, this implementation is very interesting cause i code in 
Assembly and It save a lot of MPIS to only read 2 states variable (s1 
and s2) then 4 (x1 x2 y1 y2)...

Fred.


robert bristow-johnson a �crit :
> Fred Nach wrote:
> 
>>I would like to implement a IIR Biquad filter using the fixed point
>>arithmetics...
>>Hence to reduced the intermediate states I plan to use the following trick:
>>s(k) = x(k) -a1*s(k-1) -a2*s(k-2)
>>y(n) = b0*s(k) + b1*s(k-1) + b2*s(k-2)
>>
>>The I can compute each y, saving only 2 states (s(k-1) a�nd s(k-2))...
> 
> 
> this is the Direct 2 Form (sometimes called the "Direct II Canonical
> Form").  unrecommeded for floating-point (if the Q is quite high, you
> end up subtracting numbers very close to each other and  losing
> precision) and *highly* unrecommended for fixed-point (saturation
> clipping will be your lot).
> 
> try the Direct 1 Form:
> 
>     y[n] = b0*x[n] + b1*x[n-1] + b2*x[n-2] - a1*y[n-1] - a2*y[n-2]
> 
> every product on the right side of the = sign is a double precision
> fixed point and all of those products should be added together in
> double precision.  this is trivial in the Mot 56K and the ADI SHArC and
> maybe in the TI fixed-pointers.
> 
> 
>>BUT ... in order to scale my coefficients or input, I need to know what
>>are the boundaries of the s(k) serie... in order to find a proper Fixed
>>Point format for my coefficients..
> 
> 
> your coefs are defined strictly in terms of the frequency response you
> want.  the gross range of a1 is -2 to +2 and a2 is -1 to +1 for any
> stable biquad filter.  the b0, b1, b2 coefs can be nearly anything but
> you will have to choose a range and possibly do some scaling using the
> arithmetic shift operation on the result before saving to y[n].
> 
> 
>>So my questions are :
>>
>>* for which conditions, s(k) is bounded
>>* If the s(k) is bounded, what are the boundaries?
> 
> 
> that is a purely artificial construct.  are your fixed-point signal
> values considered 16 bit signed integers?  then it's  -32768 <= x[n] <
> +32768 .  are they normalized?  then it's -1 <= x[n]< +1.  their range
> is whatever you want them to be.
> 
> 
>>Thanks in advance for your help (and sorry for my poor englsih)
> 
> 
> as Jerry said, your English is fine.
> 
> r b-j
> 
0
Reply fred_nach (6) 12/1/2005 3:18:16 PM

Me again ...

Ok, so i've just launch a Matlab plot of the filter derived from the 
s(k) sequence :

s(k) = x(k) -a1*s(k-1) -a2*s(k-2)
hence considering x(k) as input and s(k) as output
H(Z) = 1/(1 + a1*z^-1 + a2*z^-2)

the frequency response give a curve that decrease gradualy, *BUT* start 
witha gain of about 80 dB ... hence it means (as long as i 
understand...) that for almost constant x(k) (=low freq) the gain is 
very important so s(k) will take big values ... (gain of 80dB is like 
multiplying by 10 000 !) hence it for a Fixed Point architecture it 
requires ceil(log2(10000)) = 14 guard bits !!!

So yes... it's definitively not a good solution to implement a filter 
like that ... and I will do it using the direct form 1.

Salut !

Fred
0
Reply fred_nach (6) 12/1/2005 4:06:39 PM

Fred Nach wrote:
>
> It seems that my previous post was not clear enough

it was clear enough for me.

....

> * The filter is indeed stable hence  |a1| <= 2 and |a2| <= 1
>
> I think I understood that this implementation is not recommended because
> it compute first the feedback terms and leads to big variations...

it means that s[n] is too large to fit in you 16 bit words.  if you
represent them as 32 bit words, then you have to do double precision
multiplies when you multiply be b0, b1, b2.

> But is there a way to retrieve the boudary of s[n] mathematicaly ?
>
> I think that as long as a1 and a2 are given for a stable filter s[n]
> should be bounded...

it is.  the worst case is that if a1 and a2 are known, you can obtain
the impulse response for the all-pole filter

      H_s(z) = 1/(1 + a1*z^-1 + a2*z^-2)

call that impulse response h_s[n].

the worst case for the maximum value for |s[n]| is what happens if x[n]
shares the same sign as s[n] (for each value of n) and is at it at its
maximum magnitude.  from the convolution summation

      s[n] = SUM{ h_s[k] * x[n-k] }
              n

      max |s[n]| = max |x[n]| * SUM{ |h_s[n]| }
       n            n            n           (max's and SUM over all n)

you might find that max gain, SUM{ |h_s[n]| }, to be a very large
value.  much larger than 1. if you were to reduce you input to such a
degree that you guarantee that s[n] never saturates, you will likely
reduce

> For my program, this implementation is very interesting cause i code in
> Assembly and It save a lot of MIPS to only read 2 states variable (s1
> and s2) then 4 (x1 x2 y1 y2)...

there is no reason for the DF1 to require more MIPS in the inner loop
than the DF2.  at least for cascaded biquads.  if you're doing only one
biquad, that is your overall filter order is 2, then it will cost a
little more (but the cost is worth it).

the DF1 requires 2 additional states but, if you cascade 2nd order
biquads, the two output states can be shared with the 2 input states of
the next section.  so if you have N sections(an order 2N filter), the
number of states you need for DF1 is 2N+2 only 2 more than DF2.  not a
very high price to pay.

the DF2 is useless in your case (assuming that your filter will have
some decent Q making the poles very close to the unit circle).  but, if
you like learning this yourself (we call this "learning the hard way"),
feel free to implement it.  but you will find that either you
pre-scaling of the input will drive your signal into the noise floor
or,if you do not pre-scale the input, your intermediate state signal,
s[n], will saturate.  only with double precision (that will cost more
MIPS than DF1) will you avoid that.

r b-j

0
Reply rbj (4087) 12/1/2005 4:28:55 PM

Fred Nach wrote:
> Me again ...
>
> Ok, so i've just launch a Matlab plot of the filter derived from the
> s(k) sequence :
>
> s(k) = x(k) -a1*s(k-1) -a2*s(k-2)
> hence considering x(k) as input and s(k) as output
> H(Z) = 1/(1 + a1*z^-1 + a2*z^-2)
>
> the frequency response give a curve that decrease gradualy, *BUT* start
> witha gain of about 80 dB ... hence it means (as long as i
> understand...) that for almost constant x(k) (=low freq) the gain is
> very important so s(k) will take big values ... (gain of 80dB is like
> multiplying by 10 000 !) hence it for a Fixed Point architecture it
> requires ceil(log2(10000)) = 14 guard bits !!!
>
> So yes... it's definitively not a good solution to implement a filter
> like that ... and I will do it using the direct form 1.

our posts "crossed in the mail".  i'm glad that you "get it" now.  :-)

when using the DF1, make sure that all five terms (b0*x[n] b1*x[n-1],
etc) are added together in a double precision accumulator before
quantizing that back to a single precision output.  also, consider
"noise shaping" or "fraction saving" (Google Groups the latter term or
look up the DC Blocking filter in fixed-point at dspguru.com).  with
only 16 bits, you'll probably need it.

r b-j

0
Reply rbj (4087) 12/1/2005 4:35:02 PM

On 1 Dec 2005 08:35:02 -0800, "robert bristow-johnson"
<rbj@audioimagination.com> wrote:
>
>Fred Nach wrote:
>> Me again ...
>>
>> Ok, so i've just launch a Matlab plot of the filter derived from the
>> s(k) sequence :
>>
>> s(k) = x(k) -a1*s(k-1) -a2*s(k-2)
>> hence considering x(k) as input and s(k) as output
>> H(Z) = 1/(1 + a1*z^-1 + a2*z^-2)
>>
>> the frequency response give a curve that decrease gradualy, *BUT* start
>> witha gain of about 80 dB ... hence it means (as long as i
>> understand...) that for almost constant x(k) (=low freq) the gain is
>> very important so s(k) will take big values ... (gain of 80dB is like
>> multiplying by 10 000 !) hence it for a Fixed Point architecture it
>> requires ceil(log2(10000)) = 14 guard bits !!!
>>
>> So yes... it's definitively not a good solution to implement a filter
>> like that ... and I will do it using the direct form 1.
>
>our posts "crossed in the mail".  i'm glad that you "get it" now.  :-)
>
>when using the DF1, make sure that all five terms (b0*x[n] b1*x[n-1],
>etc) are added together in a double precision accumulator before
>quantizing that back to a single precision output.  also, consider
>"noise shaping" or "fraction saving" (Google Groups the latter term or
>look up the DC Blocking filter in fixed-point at dspguru.com).  with
>only 16 bits, you'll probably need it.
>

Isn't this a case where a lattice filter implementaion helps?
0
Reply just5217 (63) 12/1/2005 4:46:28 PM

Just Cocky wrote:

> Isn't this a case where a lattice filter implementaion helps?

the "normalized ladder" architecture is supposed to address this
problem.  it preserves the mean power in those internal states but does
not guarantee against saturation. in the DF1, you can scale down b0,
b1, b2 to such a point that the summer (and output y[n]) does not
saturate (and then scale it up later, if you can get away with it).

r b-j

0
Reply rbj (4087) 12/1/2005 5:07:45 PM

Just Cocky wrote:

> Isn't this a case where a lattice filter implementaion helps?

the "normalized ladder" architecture is supposed to address this
problem.  it preserves the mean power in those internal states but does
not guarantee against saturation. in the DF1, you can scale down b0,
b1, b2 to such a point that the summer (and output y[n]) does not
saturate (and then scale it up later, if you can get away with it).

r b-j

0
Reply rbj (4087) 12/1/2005 6:02:21 PM

robert bristow-johnson wrote:
> Jerry Avins wrote:
> 
>>robert bristow-johnson wrote:
>>
>>>Jerry Avins wrote:

   ...

>>OK. Back to the books with my thinking cap on.
> 
> 
> no need.  here's a simple way to look at it:

   ...

Thanks much. Part of my problem was misreading the original equations, 
the ones with () instead of [].

Jerry
-- 
Engineering is the art of making what you want from things you can get.
�����������������������������������������������������������������������
0
Reply jya (12870) 12/1/2005 6:29:17 PM

"Andrew Reilly" <andrew-newspost@areilly.bpc-users.org> wrote in message 
news:pan.2005.12.01.00.10.23.587775@areilly.bpc-users.org...
> On Wed, 30 Nov 2005 11:01:31 -0800, Noway2 wrote:
>> I have seen Matlab scripts that run an impulse response till it
>> stabilizes and then sum up the impulse responses and use the maximum
>> value as a scale factor, though I don't know if this is a correct
>> solution.
>
> It's pretty close.  This is the main reason why these "direct form 2"
> implmentations are not preferred.  A direct form 1 implmentation will use
> more memory but (a) memory is relatively cheap and (b) the values stored
> will be as constrained as the input and output values.

Even though memory is cheap, sometimes I find the overhead of the extra memory 
stores/fetches slows the filter down a bit (depending on the architecture). 


0
Reply jon99_harris7 (307) 12/3/2005 6:16:47 AM
comp.dsp 19739 articles. 22 followers. Post

19 Replies
138 Views

Similar Articles

[PageSpeed] 22

  • Permalink
  • submit to reddit
  • Email
  • Follow


Reply:

Similar Artilces:

IIR fixed-point implementation
HI everyone, I was wondering if anyone can recommend me good books or papers on fixed-point IIR implementation strategies on FPGAs or CPUs. Especially regarding fc versus Fs and coefficient quantization. I have issues with an order-1 DC-removal filter with Fs = 100 MHz and Fcut = 250 kHz. Basically, I need lots of precision in the intermediate calculations which impact the number of bits for the multiplier and thus the maximum frequency of my FPGA design. Best regards On 11/14/2011 1:15 PM, Benjamin Couillard wrote: > HI everyone, > > I was wondering if anyone can recommend me good...

IIR biquad implementation in fixed point / integer
Hello, I tried implementing a LPF biquad (that works perfectly in floating point) with: a = [1.00000, -1.90887, 0.91129], b = [6.1150e-04, 1.2230e-03, 6.1150e-04] in fixed point / integer, with 16 bit coefficients, 10 bit input signal and a 32 bit accumulator. I suspect the coefficient resolution is not enough. I tried: - the direct form 1 and direct form 2 transposed structures, but the output is zero if the input signal is not very high -- because of the small b coefficients, the signal is attenuated very much in the FIR portion of the filter and is lost; int biquad_ab(d, a0...

Fixed point implementation of 4'th order IIR filters
Hi Does anyone have some guidelines on how to implement a 4'th order low-pass Butterworth IIR filter in fixed point. My cut-off frequency is relatively close to the DC frequency so high precision is needed for the coefficients. What about realization structure and so on! I have implemented the bit-flipping algorithm in http://www.cmsa.wmin.ac.uk/~artur/pdf/Paper16.pdf for quantization of coefficients and it indeed works, but does some other techniques allow for further reductions of number of bits used to represent the coefficients. Thomas I've found the normalized lattice-ladder...

couple of questions on practical IIR filter implementation on fixed point systems
Hello. I am a professional programmer, but I am doing DSP just as a hobby. Three years ago I designed a morse / PSK31 decoder on a 16bit no- multiplier controller. After rouhgly 3 years of absence in the DSP field I am back to design a HAM shortwave software defined radio using a beautiful new part from Texas Instruments - TLV320AIC3254 codec with miniDSP core. The part could do all the signal processing I need with just a 5x5mm board space including earphone driver. That's an engineer's dream come true. The goal is to downconvert a SSB voice channel down to audio with a quadrature m...

floating point to fixed point conversion-in implementation of navigation system algorithm
Hi, I have a floating point c code to implement a navigation system algorithm. In order to implement it in the 64x+ fixed point DSP processor, i need to have the code in fixed point. I learnt that there are few processor specific IQmath libraries aiding this purpose. Is it practical to write the code using these libraries or should i manually write codes for each operation by scaling the inputs for each line of the code, in which case i have to create look up table for the trigonometric functions etc which is a harder process to begin with.. Thanks, Divya On 2/13/12 11:38 AM, divya_di...

fixed point implementation
Hi, i want to implement the equation p3 = p3*0.5/cos(3*PI/8); into 16bit processor, where p3 is a 16 bit value; any pointers as to how to go about it -vijay n_vijay_anand wrote: > Hi, > > i want to implement the equation > > p3 = p3*0.5/cos(3*PI/8); > > into 16bit processor, where p3 is a 16 bit value; > > any pointers as to how to go about it It looks like '0.5/cos(3*PI/8)' is a constant, so simplify it. cos(3*PI/8) = 0.99978861 0.5/0.99978861 = 0.50010571 Depending on the algorithm you might get away with simply dividing by 2: p3 >...

LMS fixed point implementation
Hi everyone! I need a fixed point implementation for the LMS algorithm;I have to implement it on a TI dsp board (c6711 or a c5000); can anyone suggest me where I can find this code? thank you guys! f Use this: It is a single step for fixed point LMS. c_ptr points to coefficient, r_ptr points to history. muErr = e * mu; mu = 2048(trial), using fixed mu; accA, accB are at least 32 bit accumulators. /* Macro: Performs single step of filter and adapt */ #define ADAPT_AND_FILTER(accA, accB, muErr, c_ptr, r_ptr) \ accB = 16384; \ accB += ((i3...

fixed point IIR coeffs
Hi i am using built-in Matlab functions to design IIR fixed point filter. For example: [b a] = cheby1(2,6,0.25) hd = dfilt.df1(b,a) set(hd,'arithmetic','fixed') Is there a simple method to check the coefficient values for fixed arithmetic? Precisely the ones which will appear in HDL code ( generatehdl(hd) ). Because i want to make a parametrized filter with a coefficient lookup table for different cutoff frequencies. ...

IIR with fixed-point numbers?
What�s the best way to determine the required number of bits necessary in a fixed-point IIR-filter? To avoid overflow during calculations. "Elektro" <blabla@bredband.net> wrote in news:4293a2ac$0$79465$14726298@news.sunsite.dk: > What�s the best way to determine the required number of bits necessary > in a fixed-point IIR-filter? To avoid overflow during calculations. > > > The number of bits needed in an IIR filter depends on many factors including the coefficients of the filter, the desired performance of the filter and the implementation of the filt...

fixed FFT point implementation woes
Hi, I am trying to create a mex file that does fixed point FFT (512 point?). I am generally ok with creating a C code in single precision to compute the the result and have ensured that it works fine as i get the correct result. (i read it back and plotted it in matlab.) the problem comes in when i try to convert it to fixed point. My input data is 24 bits integer. In addition, I take the sine and cosine twiddles to also be 24 bits... I do the cos(2*pi*n/N) and sin(2*pi/N) for n from 0 to 255 the multply the result by ((2^23)-1) and store it in a look up table. What i really dun qu...

IIR Filter fixed point DSP
Hi I am trying to design a IIR filter using a 16-bit fixed point DSP. I nee to calculate the coeff in the form of second order sections. I would lik to do this with MATLAB. Does anyone have some MATLAB code that will d this. I found a routine on Texas Instrument site but as well as the coef it creates some scaling factors for the inputs. I dont want this as I hav seen other programs such as Filter Express that just give the Coeff. Thanks J This message was sent using the Comp.DSP web interface o www.DSPRelated.com OK. You have several kinds of IIR filter, for instance elliptic (the be...

bit flipping for fixed point implementation!
Hi :) Has anyone tried to implement the algorithms presented in the following paper? http://www.cmsa.wmin.ac.uk/~artur/pdf/Paper16.pdf I'd be happy to see the code then if possible Best Regards Thomas ...

Goertzel fixed point implementation question
Hello, I have found fixed point implementation of Goertzel algorithm in boo “Real – Time Digital Signal Processing – Implementations an Applications”. I will paste here a part of code responsible bo computing recursive path. void gFilter (short *d, short in, short coef) { long AC0; d[0] = in >> 7; // Get input with scale down by 7 bit AC0 = (long) d[1] * coef; AC0 >>= 14; AC0 -= d[2]; d[0] += (short)AC0; d[2] = d[1]; // Update delay-line buffer d[1] = d[0]; } Book explains that this operation (AC0 >>= 14) aligns multiplicatio product to be...

About range estimation for fixed point implementation
I am reading a paper about fixed point implementation. The author has given an example for input range estimation for multi user detection. It is said "From the characteristic of the Gold code, we know that the maximum value of cross-correlation coefficients is the auto correlation of any particular spreading sequence,i.e., rang R is R=2x(2^r-1), where the spreading gain is 2^r-1". I can understand the maximum value part but I don't know how to deduce the range estimation is 2 x (spreading gain). Is there anybody explain to me? Many thanks. Zhi <threeinchnail@gmail.com> ...

Fixed point number hardware implementation
Hai, I want to know which is the right way of implementing and usage of fixed point number data types in hardware(industry standard)..I have referred various FIR implementations where they are mostly handling filter coefficients as integer(truncating from fixed or floating point using MATLAB) or binary.Is it difficult to handle and implement real(fraction) numbers i.e.,filter coefficients values directly in the hardware? for example: sample Filter coefficients generated by FDA tool: fixed point=0.211944580078125 or 16-bit signed integer= 13890 or fixed point binary =0011011001000010 all t...

FFT resolution in fixed-point implementation
What is the effect of lowered resolution of intermediate results on an FFT? I am considering a number of implementation choices for a fixed-point FFT. The time-series data is inherently 16-bit integers (audio data). I want to do a 8192-point FFT on this data. Because of the need to explicitly handle overflows (in the butterfly operations) in a fixed-point implementation, it is convenient to limit the resolution of the intermediate results. There is a definite tradeoff between speed and how I check for overflow. I am also limited by the C language. Ideally I would like to multiply two...

sin/cos implementation in fixed-point
Hello, Here's my problem. I converted some code from floating-point to fixed point. Currently, for sin/cos/tan I'm using the floating point math.h functions. I have two base formats q.22 and q.28 that I use and play with. I was initially using 1024 values sin/cos tables (in q.28 format). However, I'm finding that this is not giving me good-enough precision. I'm getting abosolute differences of about 0.004-0.008 between the floating point sin/cos and the lookup tables. The alogrithm (I'm using) seems very sensitive and I'd like to achieve an absolute error of 0.000001 ...

Fixed Point implementation of C exponent function
Hi, I tried to google a lot searching for fixed point implementation of exponent function y = exp(x) (e to the power x) but could not find it. Also I think there is no post on this group also. I found some CORDIC implementations but those were in C++. I need lightweight implementation in C. I want to use this function in an image processing algorithm to be implemented on a VLIW-DSP which doesn't have a floating point unit. I will appreciate in anyone here can suggest a plain C - code which handles both positive and negative arguments (x<0 and x<=0) ,with reasonable accu...

IIR filter on a fixed set of data points
I have a fixed 3600 data points(samples) and I want to apply a notch 2nd order Butterworth IIR Notch filter to remove the 1/4Fs frequency, the results start to converge(or filter start to work) after about 250 or so samples. The issue is the it is a continuous stream of data and I need to remove the harmonics on the initial 250 samples. What can I do to filter the whole data set? thank you, Andrew ahgu <honestgold@n_o_s_p_a_m.gmail.com> wrote: > I have a fixed 3600 data points(samples) and I want to apply a notch 2nd > order Butterworth IIR Notch filter to re...

Fixed-point implementation of levinson durbin algorithm
Hi, Any links to a fixed point implementation of the Levinson Durbin algorithm ...something like this: void levdur(pLpcCoeff, pAutoCorr, pReflec, nOrder) where pLpcCoeff is a pointer to the lpc coefficients where pAutoCorr is a pointer to the autocorrelation coefficients (range -40 to 40) where pReflec is a pointer to the reflection coefficients where nOrder is 10 The autocorrelation sequence is a 32bit array where values are in Q15 format. Thanks >Hi, > >Any links to a fixed point implementation of the Levinson Durbin >algorithm ...something like this: > >void levdur(...

iLBC codec: Fixed-point implementation of iCBSearch ??
Hello, Anybody out there who made a fixed point implementation of iCBSearch ?? John ...

multirate implementation on fix point DSP and efficiency
Hye, I have a question to try to find a way to improve my processing. I use a multirate filter with decimation. There is a decimation by 40. One FIr is used then decimation, and the IIR with decimation, then a second IIR with decimation, then again one IIR and decimation and at last one FIr . I use a fix point DSP. The filter are designed in 16 bits inputs / output (32 bits for computation). It appears that the results is that i have in output a 12bits resoutio and 4 lsb bits are lost during the stages. How can i improve this? Should i change the IIR design? or thedecimation? thank yo...

Implementing an 8 bit fixed point register
Hello to all I'm about to write a simulator for a microcontroller in python (why python? because I love it!!!) but I have a problem. The registry of this processor are all 8 bit long (and 10 bit for some other strange register) and I need to simulate the fixed point behaviour of the register, and to access the single bit. f.x. (this is a pseudo python session, only for understanding) >>> reg1 = fixed_int(8) >>> reg2 = fixed_int(10) >>> reg1[0].set() or >>> reg1[0] = 1 # or True? how to rapresent a binary bit >>> reg1[0] 1 >>> reg1[1...

Maximizing dynamic range of fixed point IIR filter
Let's say we need to implement H(z) = P(z)/Q(z) in the fixed point. The typical implementation would be a cascade of biquads. So we factor P(z) and Q(z) and distribute poles and zeroes between the stages. The dynamic range of a filter is limited by overflow at the top, and by quantization artifacts at the bottom. We can try all variants of assignment of poles and zeroes to different stages to maximize the dynamic range from rms quantization noise to full scale sine wave at the "worst" frequency. So far so good. However, this doesn't tell if some stage of f...