Forums

Please check my filter!

Started by NightHawk July 18, 2005
Im trying to implement an 8th order Butterworth low pass filter in Matlab.
I know you can design filters in matlab, but because of the specific
implementation of my project, im doing with no matlab functions in an m
file. So, the problem im having is that the input before filtering is
EXACTLY the same as the output... what am I doing wrong?

The sampling frecuency is 44100 hz.
The filter cutoff is 60 hz(I know its low!!)

t=0:1/Fs:1;
x=sin(2*pi*40*t)+sin(2*pi*1000*t);

y=[1:1:length(x)];

for n=0:1:length(x)
    y=0;
end

for n=9:1:length(x)
    n

y(n) = (  1 * x(n- 8));...
     + (  8 * x(n- 7));...
     + ( 28 * x(n- 6));...
     + ( 56 * x(n- 5));...
     + ( 70 * x(n- 4));...
     + ( 56 * x(n- 3));...
     + ( 28 * x(n- 2));...
     + (  8 * x(n- 1));...
     + (  1 * x(n- 0));...

     + ( -0.9571275625 * y(n- 8));...
     + (  7.6989603129 * y(n- 7));...
     + (-27.0940698453 * y(n- 6));...
     + ( 54.4854090892 * y(n- 5));...
     + (-68.4806794719 * y(n- 4));...
     + ( 55.0855564533 * y(n- 3));...
     + (-27.6942305609 * y(n- 2));...
     + (  7.9561815852 * y(n- 1));

end

plot(0:1/Fs:1,x)
plot(0:1/Fs:1,y)


		
This message was sent using the Comp.DSP web interface on
www.DSPRelated.com
>Im trying to implement an 8th order Butterworth low pass filter in
Matlab.
>I know you can design filters in matlab, but because of the specific >implementation of my project, im doing with no matlab functions in an m >file. So, the problem im having is that the input before filtering is >EXACTLY the same as the output... what am I doing wrong? > >The sampling frecuency is 44100 hz. >The filter cutoff is 60 hz(I know its low!!) > >t=0:1/Fs:1; >x=sin(2*pi*40*t)+sin(2*pi*1000*t); > >y=[1:1:length(x)]; > >for n=0:1:length(x) > y=0; >end > >for n=9:1:length(x) > n > >y(n) = ( 1 * x(n- 8));... > + ( 8 * x(n- 7));... > + ( 28 * x(n- 6));... > + ( 56 * x(n- 5));... > + ( 70 * x(n- 4));... > + ( 56 * x(n- 3));... > + ( 28 * x(n- 2));... > + ( 8 * x(n- 1));... > + ( 1 * x(n- 0));... > > + ( -0.9571275625 * y(n- 8));... > + ( 7.6989603129 * y(n- 7));... > + (-27.0940698453 * y(n- 6));... > + ( 54.4854090892 * y(n- 5));... > + (-68.4806794719 * y(n- 4));... > + ( 55.0855564533 * y(n- 3));... > + (-27.6942305609 * y(n- 2));... > + ( 7.9561815852 * y(n- 1)); > >end > >plot(0:1/Fs:1,x) >plot(0:1/Fs:1,y) > > > >This message was sent using the Comp.DSP web interface on >www.DSPRelated.com >
I been playing around with the pocedures in matalb... and have restructured the equations inside the FOR loop, and now their are all in the same line. Still.. no good, however, I now have another problem. The output of the filter is now a horizontal line intersected with a verticl line... and when I examine the output vector, from sample 22288 onwards is full of NaNs, how can this be if all im doing is adding numbers? This message was sent using the Comp.DSP web interface on www.DSPRelated.com
NightHawk wrote:
>>Im trying to implement an 8th order Butterworth low pass filter in > > Matlab. > >>I know you can design filters in matlab, but because of the specific >>implementation of my project, im doing with no matlab functions in an m >>file. So, the problem im having is that the input before filtering is >>EXACTLY the same as the output... what am I doing wrong? >> >>The sampling frecuency is 44100 hz. >>The filter cutoff is 60 hz(I know its low!!) >> >>t=0:1/Fs:1; >>x=sin(2*pi*40*t)+sin(2*pi*1000*t); >> >>y=[1:1:length(x)]; >> >>for n=0:1:length(x) >> y=0; >>end >> >>for n=9:1:length(x) >> n >> >>y(n) = ( 1 * x(n- 8));... >> + ( 8 * x(n- 7));... >> + ( 28 * x(n- 6));... >> + ( 56 * x(n- 5));... >> + ( 70 * x(n- 4));... >> + ( 56 * x(n- 3));... >> + ( 28 * x(n- 2));... >> + ( 8 * x(n- 1));... >> + ( 1 * x(n- 0));... >> >> + ( -0.9571275625 * y(n- 8));... >> + ( 7.6989603129 * y(n- 7));... >> + (-27.0940698453 * y(n- 6));... >> + ( 54.4854090892 * y(n- 5));... >> + (-68.4806794719 * y(n- 4));... >> + ( 55.0855564533 * y(n- 3));... >> + (-27.6942305609 * y(n- 2));... >> + ( 7.9561815852 * y(n- 1)); >> >>end >> >>plot(0:1/Fs:1,x) >>plot(0:1/Fs:1,y) >> >> >> >>This message was sent using the Comp.DSP web interface on >>www.DSPRelated.com >> > > > > I been playing around with the pocedures in matalb... and have > restructured the equations inside the FOR loop, and now their are all in > the same line. Still.. no good, however, I now have another problem. The > output of the filter is now a horizontal line intersected with a verticl > line... and when I examine the output vector, from sample 22288 onwards is > full of NaNs, how can this be if all im doing is adding numbers?
I see some multiplications, too. Jerry -- Engineering is the art of making what you want from things you can get. �����������������������������������������������������������������������
I could not get the filter to work, so what I did was to design a second
order filter and apply it 4 times to simulate a 8th order filter, it seems
to work, but I wonder if there are any problems or disavantages for using
the a filter many times. Is there a problem or is it good?
		
This message was sent using the Comp.DSP web interface on
www.DSPRelated.com
"NightHawk" <nightstalker_099@yahoo.com> wrote in message 
news:rpWdnYiLyOyt-UHfRVn-qg@giganews.com...
>I could not get the filter to work, so what I did was to design a second > order filter and apply it 4 times to simulate a 8th order filter, it seems > to work, but I wonder if there are any problems or disavantages for using > the a filter many times. Is there a problem or is it good?
That is a very good way to implement a high-order IIR filter. Breaking it into sections lessens problems of round-off error which can lead to instability. Keep in mind to get a true Butterworth response, the 4 individual filters each need to be different. Cascading 4 identical 2nd order Butterworth filters doesn't give you an 8th order Butterworth filter. This was discussed in detail here: http://tinyurl.com/bp22k Also, to answer your previous question, you can get NaNs in a filter because of the feedback involved. If the filter is unstable, the output will grow continually until it overflows the number system. With fixed point, this usually means clipping at full scale. With floating point like Matlab uses, this means NaNs.
NightHawk wrote:
> I could not get the filter to work, so what I did was to design a second > order filter and apply it 4 times to simulate a 8th order filter, it seems > to work, but I wonder if there are any problems or disavantages for using > the a filter many times. Is there a problem or is it good? > > This message was sent using the Comp.DSP web interface on > www.DSPRelated.com
That's an eighth-order filter, but not one that takes advantage of the degrees of freedom available to it. You have a cascade of four second-order filters; that's good. Identical second-order sections are not very good. A direct implementation of an eighth-order filter is likely to be numerically unstable; that's consistent with NANs. Decompose the filter into four second-order sections. (Partial fraction decomposition is likely to help with the math.) Then cascade those second-order sections and you will have the filter you want. Jerry -- Engineering is the art of making what you want from things you can get. &#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;
"Jerry Avins" <jya@ieee.org> wrote in message 
news:3eqdna7-g54760HfRVn-3w@rcn.net...
> NightHawk wrote: >> I could not get the filter to work, so what I did was to design a second >> order filter and apply it 4 times to simulate a 8th order filter, it seems >> to work, but I wonder if there are any problems or disavantages for using >> the a filter many times. Is there a problem or is it good? >> >> This message was sent using the Comp.DSP web interface on >> www.DSPRelated.com > > That's an eighth-order filter, but not one that takes advantage of the degrees > of freedom available to it. You have a cascade of four second-order filters; > that's good. Identical second-order sections are not very good. > > A direct implementation of an eighth-order filter is likely to be numerically > unstable; that's consistent with NANs. Decompose the filter into four > second-order sections. (Partial fraction decomposition is likely to help with > the math.) Then cascade those second-order sections and you will have the > filter you want.
Sense the OP is working in Matlab, he can use the function TF2SOS to easily convert to second order sections. This assumes the Signal Processing Toolbox. If not, Roots can be used to at least factor the numerator and denominator polynomials.
How did you get the filter coefficients? I am unable to get the frequency
response of the filter as you mentioned (60Hz Cut Off). This is what i
tried..

Fs = 44100;
b = [1 8 28 56 70 56 28 8 1];
a = [7.9561815852 -27.6942305609 55.0855564533 -68.4806794719
54.4854090892 -27.0940698453 7.6989603129 -0.9571275625];

figure;freqz(b,a,8192,Fs);

where 'a' and 'b' are from your matlab code.

Regards,
Ajith


>NightHawk wrote: >> I could not get the filter to work, so what I did was to design a
second
>> order filter and apply it 4 times to simulate a 8th order filter, it
seems
>> to work, but I wonder if there are any problems or disavantages for
using
>> the a filter many times. Is there a problem or is it good? >> >> This message was sent using the Comp.DSP web interface on >> www.DSPRelated.com > >That's an eighth-order filter, but not one that takes advantage of the >degrees of freedom available to it. You have a cascade of four >second-order filters; that's good. Identical second-order sections are >not very good. > >A direct implementation of an eighth-order filter is likely to be >numerically unstable; that's consistent with NANs. Decompose the filter >into four second-order sections. (Partial fraction decomposition is >likely to help with the math.) Then cascade those second-order sections >and you will have the filter you want. > >Jerry >-- >Engineering is the art of making what you want from things you can get. >&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095; >
This message was sent using the Comp.DSP web interface on www.DSPRelated.com
Ajith wrote:
> How did you get the filter coefficients? I am unable to get the frequency > response of the filter as you mentioned (60Hz Cut Off). This is what i > tried..
... Did I mention a particular filter? I wrote:
>>That's an eighth-order filter, but not one that takes advantage of the >>degrees of freedom available to it. You have a cascade of four >>second-order filters; that's good. Identical second-order sections are >>not very good. >> >>A direct implementation of an eighth-order filter is likely to be >>numerically unstable; that's consistent with NANs. Decompose the filter >>into four second-order sections. (Partial fraction decomposition is >>likely to help with the math.) Then cascade those second-order sections >>and you will have the filter you want.
Jerry -- Engineering is the art of making what you want from things you can get. &#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;&#2013266095;
Sorry Jerry. My message is the reply for the first message
("NightHawk's"). I just replied to the last message in the thread.

Regards,
Ajith

Jerry Avins wrote:
> Ajith wrote: > > How did you get the filter coefficients? I am unable to get the frequen=
cy
> > response of the filter as you mentioned (60Hz Cut Off). This is what i > > tried.. > > ... > > Did I mention a particular filter? I wrote: > > >>That's an eighth-order filter, but not one that takes advantage of the > >>degrees of freedom available to it. You have a cascade of four > >>second-order filters; that's good. Identical second-order sections are > >>not very good. > >> > >>A direct implementation of an eighth-order filter is likely to be > >>numerically unstable; that's consistent with NANs. Decompose the filter > >>into four second-order sections. (Partial fraction decomposition is > >>likely to help with the math.) Then cascade those second-order sections > >>and you will have the filter you want. > > Jerry > -- > Engineering is the art of making what you want from things you can get. > =AF=AF=AF=AF=AF=AF=AF=AF=AF=AF=AF=AF=AF=AF=AF=AF=AF=AF=AF=AF=AF=AF=AF=AF=
=AF=AF=AF=AF=AF=AF=AF=AF=AF=AF=AF=AF=AF=AF=AF=AF=AF=AF=AF=AF=AF=AF=AF=AF=AF= =AF=AF=AF=AF=AF=AF=AF=AF=AF=AF=AF=AF=AF=AF=AF=AF=AF=AF=AF=AF=AF=AF