DSPRelated.com
Forums

Scilab IIR filter, Bode plot trouble

Started by msvilans.dsp March 17, 2011
Dear Comp.DSP,

I am using Scilab 5.3.0 (64-bit) to experiment with digital IIR filter
design.

When I create a Bode plot of the filter response, the magnitudes don't
appear to make sense. I created a low pass filter, and plotted the
magnitudes on a linear scale, using the frmag() function. The response for
the passband frequencies was at 1.0, as expected. In a Bode plot of the
same function, the passband was sitting at 20dB, not 0dB as I would have
expected.

Also, the filter passband has a dip due to ripple, but this is not visible
in the Bode plot at all.

This image illustrates what I mean:
http://www.aeonyx.ca/public/2011.03.16-iir-bode-plot-01.png

Code to generate the image is appended to this message.

Any help explaining what is happening would be greatly appreciated.

Best regards,
Markus.

-- Begin Scilab code --



On 03/17/2011 06:19 AM, msvilans.dsp wrote:
> Dear Comp.DSP, > > I am using Scilab 5.3.0 (64-bit) to experiment with digital IIR filter > design. > > When I create a Bode plot of the filter response, the magnitudes don't > appear to make sense. I created a low pass filter, and plotted the > magnitudes on a linear scale, using the frmag() function. The response for > the passband frequencies was at 1.0, as expected. In a Bode plot of the > same function, the passband was sitting at 20dB, not 0dB as I would have > expected. > > Also, the filter passband has a dip due to ripple, but this is not visible > in the Bode plot at all. > > This image illustrates what I mean: > http://www.aeonyx.ca/public/2011.03.16-iir-bode-plot-01.png > > Code to generate the image is appended to this message. > > Any help explaining what is happening would be greatly appreciated. > > Best regards, > Markus. > > -- Begin Scilab code --
?? I don't see the appended code... Certainly the Bode plot doesn't appear to be correct. If you plot the output from frmag() as a log log (i.e., instead of using 'plot2d(f, x)', use 'plot2d(f, x, logflags='ll')'), then it'll be the magnitude portion of your bode. -- Tim Wescott Wescott Design Services http://www.wescottdesign.com Do you need to implement control loops in software? "Applied Control Theory for Embedded Systems" was written for you. See details at http://www.wescottdesign.com/actfes/actfes.html
>> Any help explaining what is happening would be greatly appreciated. >> >> Best regards, >> Markus. >> >> -- Begin Scilab code -- > ?? > >I don't see the appended code... > >Certainly the Bode plot doesn't appear to be correct. If you plot the >output from frmag() as a log log (i.e., instead of using 'plot2d(f, x)', >use 'plot2d(f, x, logflags='ll')'), then it'll be the magnitude portion >of your bode. >
Hi Tim, Sorry, I must have neglected to append the code. I've included it below. The plot2d() function will not allow me to plot in log-log. It complains that all values must be non-negative. I assume that there must be some zero values, which cause it to fail. (Any suggestions to a Scilab newbie on how to get around this?) The Bode plot is handy because it just makes the plot, although it looks wrong. Regards, Markus -- Begin Scilab code -- filterType = 'ellip'; cutoff1 = 0.3; cutoff2 = 0.4; passbandRipple = 0.05; stopbandRipple = 0.025; [cells, fact, zzeros, zpoles] = eqiir('lp', filterType, ... 2 * %pi * [ cutoff1, cutoff2 ], ... passbandRipple, stopbandRipple); // Generate the equivalent linear system of the filter num = fact * real(poly(zzeros,'s')); den = real(poly(zpoles,'s')); elatf = syslin('c', num, den); // Get filter frequency response [hf, fr] = frmag(num, den, 1024); clf subplot(1, 2, 1); plot2d(fr, abs(hf)); xgrid(2); xtitle("IIR Filter Response"); subplot(1, 2, 2); bode(elatf, 0.0000001, 10000); xtitle("IIR Filter Response Bode Plot"); -- End Scilab code --
On 03/17/2011 06:24 PM, msvilans.dsp wrote:
>>> Any help explaining what is happening would be greatly appreciated. >>> >>> Best regards, >>> Markus. >>> >>> -- Begin Scilab code -- >> ?? >> >> I don't see the appended code... >> >> Certainly the Bode plot doesn't appear to be correct. If you plot the >> output from frmag() as a log log (i.e., instead of using 'plot2d(f, x)', >> use 'plot2d(f, x, logflags='ll')'), then it'll be the magnitude portion >> of your bode. >> > > Hi Tim, > > Sorry, I must have neglected to append the code. I've included it below. > > The plot2d() function will not allow me to plot in log-log. It complains > that all values must be non-negative. I assume that there must be some zero > values, which cause it to fail. (Any suggestions to a Scilab newbie on how > to get around this?) The Bode plot is handy because it just makes the plot, > although it looks wrong. > > Regards, > Markus > > -- Begin Scilab code -- > filterType = 'ellip'; > cutoff1 = 0.3; > cutoff2 = 0.4; > passbandRipple = 0.05; > stopbandRipple = 0.025; > > [cells, fact, zzeros, zpoles] = eqiir('lp', filterType, ... > 2 * %pi * [ cutoff1, cutoff2 ], ... > passbandRipple, stopbandRipple); > > // Generate the equivalent linear system of the filter > num = fact * real(poly(zzeros,'s')); > den = real(poly(zpoles,'s')); > elatf = syslin('c', num, den); > > // Get filter frequency response > [hf, fr] = frmag(num, den, 1024); > > clf > subplot(1, 2, 1); > plot2d(fr, abs(hf)); > xgrid(2); > xtitle("IIR Filter Response"); > > subplot(1, 2, 2); > bode(elatf, 0.0000001, 10000); > xtitle("IIR Filter Response Bode Plot"); > -- End Scilab code -- >
Ah ha! I thought this might be the case, but I wasn't sure enough to say anything -- the 'eqiir' filter function is giving you a filter design in the z domain, but you're making a filter out of them in the Laplace domain. Try making your polynomials in 'z', and use 'd' for the time domain in your syslin call. See if that helps. Let us know how it works out. -- Tim Wescott Wescott Design Services http://www.wescottdesign.com Do you need to implement control loops in software? "Applied Control Theory for Embedded Systems" was written for you. See details at http://www.wescottdesign.com/actfes/actfes.html
>> -- Begin Scilab code -- >> filterType = 'ellip'; >> cutoff1 = 0.3; >> cutoff2 = 0.4; >> passbandRipple = 0.05; >> stopbandRipple = 0.025; >> >> [cells, fact, zzeros, zpoles] = eqiir('lp', filterType, ... >> 2 * %pi * [ cutoff1, cutoff2 ], ... >> passbandRipple, stopbandRipple); >> >> // Generate the equivalent linear system of the filter >> num = fact * real(poly(zzeros,'s')); >> den = real(poly(zpoles,'s')); >> elatf = syslin('c', num, den); >> >> // Get filter frequency response >> [hf, fr] = frmag(num, den, 1024); >> >> clf >> subplot(1, 2, 1); >> plot2d(fr, abs(hf)); >> xgrid(2); >> xtitle("IIR Filter Response"); >> >> subplot(1, 2, 2); >> bode(elatf, 0.0000001, 10000); >> xtitle("IIR Filter Response Bode Plot"); >> -- End Scilab code -- >> > >Ah ha! I thought this might be the case, but I wasn't sure enough to >say anything -- the 'eqiir' filter function is giving you a filter >design in the z domain, but you're making a filter out of them in the >Laplace domain. > >Try making your polynomials in 'z', and use 'd' for the time domain in >your syslin call. See if that helps. > >Let us know how it works out.
Hi, Tim, Your suggestion worked. The resulting plots look like this: http://www.aeonyx.ca/public/2011.03.19-iir-bode-plot-02.png I've included the modified Scilab code below. Thanks very much for your help! It seems it's time for me to brush up on my z- and Laplace transforms, and how to properly use them in Scilab. Best regards, Markus -- Begin Scilab code -- filterType = 'ellip'; cutoff1 = 0.3; cutoff2 = 0.4; passbandRipple = 0.05; stopbandRipple = 0.025; [cells, fact, zzeros, zpoles] = eqiir('lp', filterType, ... 2 * %pi * [ cutoff1, cutoff2 ], ... passbandRipple, stopbandRipple); // Generate the equivalent linear system of the filter num = fact * real(poly(zzeros, 'z')); den = real(poly(zpoles, 'z')); elatf = syslin('d', num, den); // Get filter frequency response [hf, fr] = frmag(num, den, 1024); clf subplot(1, 2, 1); plot2d(fr, abs(hf)); xgrid(2); xtitle("IIR Filter Response"); subplot(1, 2, 2); bode(elatf, 0.01, 10000); xtitle("IIR Filter Response Bode Plot"); -- End Scilab code --
On Sat, 19 Mar 2011 17:34:52 -0500, msvilans.dsp wrote:

>>> -- Begin Scilab code -- >>> filterType = 'ellip'; >>> cutoff1 = 0.3; >>> cutoff2 = 0.4; >>> passbandRipple = 0.05; >>> stopbandRipple = 0.025; >>> >>> [cells, fact, zzeros, zpoles] = eqiir('lp', filterType, ... >>> 2 * %pi * [ cutoff1, cutoff2 ], ... >>> passbandRipple, stopbandRipple); >>> >>> // Generate the equivalent linear system of the filter num = fact * >>> real(poly(zzeros,'s')); >>> den = real(poly(zpoles,'s')); >>> elatf = syslin('c', num, den); >>> >>> // Get filter frequency response >>> [hf, fr] = frmag(num, den, 1024); >>> >>> clf >>> subplot(1, 2, 1); >>> plot2d(fr, abs(hf)); >>> xgrid(2); >>> xtitle("IIR Filter Response"); >>> >>> subplot(1, 2, 2); >>> bode(elatf, 0.0000001, 10000); >>> xtitle("IIR Filter Response Bode Plot"); -- End Scilab code -- >>> >>> >>Ah ha! I thought this might be the case, but I wasn't sure enough to >>say anything -- the 'eqiir' filter function is giving you a filter >>design in the z domain, but you're making a filter out of them in the >>Laplace domain. >> >>Try making your polynomials in 'z', and use 'd' for the time domain in >>your syslin call. See if that helps. >> >>Let us know how it works out. > > > Hi, Tim, > > Your suggestion worked. > > The resulting plots look like this: > http://www.aeonyx.ca/public/2011.03.19-iir-bode-plot-02.png > > I've included the modified Scilab code below. > > Thanks very much for your help! > It seems it's time for me to brush up on my z- and Laplace transforms, > and how to properly use them in Scilab. > > Best regards, > Markus > > -- Begin Scilab code -- > filterType = 'ellip'; > cutoff1 = 0.3; > cutoff2 = 0.4; > passbandRipple = 0.05; > stopbandRipple = 0.025; > > [cells, fact, zzeros, zpoles] = eqiir('lp', filterType, ... > 2 * %pi * [ cutoff1, cutoff2 ], ... > passbandRipple, stopbandRipple); > > // Generate the equivalent linear system of the filter num = fact * > real(poly(zzeros, 'z')); den = real(poly(zpoles, 'z')); > elatf = syslin('d', num, den); > > // Get filter frequency response > [hf, fr] = frmag(num, den, 1024); > > clf > subplot(1, 2, 1); > plot2d(fr, abs(hf)); > xgrid(2); > xtitle("IIR Filter Response"); > > subplot(1, 2, 2); > bode(elatf, 0.01, 10000); > xtitle("IIR Filter Response Bode Plot"); -- End Scilab code --
You might find that your Bode plot looks better if you use the call bode(elatf, 0.01, 0.499) This "cheats" the thing a little bit, so that it doesn't investigate the depth of the null at f = fs/2. The alternate way is a pain in the behind, but you can get a handle to the figure (g = gcf()), then spelunk through the various child object handles until you find the magnitude Bode plot and change it's data range so that it only shows down to -30 or -50dB -- then you'll have a chance of seeing the passband detail. -- http://www.wescottdesign.com