57409dbc4b9677388c289bab48f72172054a5409
1 <?xml version="1.0" encoding="UTF-8"?>
2 <refentry version="5.0-subset Scilab" xml:id="DesignEllipticFilter"
3           xmlns="http://docbook.org/ns/docbook"
5           xmlns:xi="http://www.w3.org/2001/XInclude"
6           xmlns:svg="http://www.w3.org/2000/svg"
7           xmlns:mml="http://www.w3.org/1998/Math/MathML"
8           xmlns:html="http://www.w3.org/1999/xhtml"
9           xmlns:db="http://docbook.org/ns/docbook">
10   <refnamediv>
11     <refname>How to design an elliptic filter</refname>
13     <refpurpose>How to design an elliptic filter (analog and
14     digital)</refpurpose>
15   </refnamediv>
17   <refsection>
18     <title>Description</title>
20     <para>The goal is to design a simple analog and digital elliptic
21     filter.</para>
22   </refsection>
24   <refsection>
25     <title>Designing an analog elliptic filter</title>
27     <para>There are several possibilities to design an elliptic lowpass
30     zeros of the filter. Once we have got these poles and zeros, we will have
31     to translate this representation into a <link
34     <para>And then, the filter can be represented in bode plot.</para>
36     <programlisting>
37 // analog elliptic (Bessel), order 2, cutoff 1 Hz
38 Epsilon = 3;  // ripple of filter in pass band (0&lt;epsilon&lt;1)
39 A       = 60; // attenuation of filter in stop band (A&gt;1)
40 OmegaC  = 10; // pass band cut-off frequency in Hertz
41 OmegaR  = 50; // stop band cut-off frequency in Hertz
43 // Generate the filter
44 [_zeros,pols,gain] = zpell(3,60,10,50);
46 // Generate the equivalent linear system of the filter
47 num   = gain * real(poly(_zeros,'s'));;
48 den   = real(poly(pols,'s'));
49 elatf = syslin('c',num,den);
51 // Plot the resulting filter
52 bode(elatf,0.01,100);
53 title('Analog Elliptic filter');
54  </programlisting>
56     <para>Bode plot is only suited for analog filters.</para>
58     <mediaobject>
59       <imageobject>
60         <imagedata align="center"
61                    fileref="../../images/analog_elliptic_filter.png" />
62       </imageobject>
63     </mediaobject>
65     <para>If you want to design a highpass, bandpass or bandstop filter, you
66     can first design a lowpass and then transfrom this lowpass filter using
68   </refsection>
70   <refsection>
71     <title>Designing a digital elliptic filter</title>
73     <para>Now, let's focus on how to produce a digital lowpass elliptic
74     filter.</para>
76     <para>We can produce two kinds of digital filters:</para>
78     <itemizedlist>
79       <listitem>
80         <para>an IIR (Infinite Impulse Response).</para>
82         <para>To compute such a filter, we can use the following
83         functions:</para>
85         <itemizedlist>
86           <listitem>
88           </listitem>
90           <listitem>
92           </listitem>
93         </itemizedlist>
94       </listitem>
96       <listitem>
97         <para>a FIR (Finite Impulse Response).</para>
99         <para>To compute such a filter, we can use the following
100         functions:</para>
102         <itemizedlist>
103           <listitem>
105           </listitem>
107           <listitem>
109           </listitem>
111           <listitem>
113           </listitem>
115           <listitem>
117           </listitem>
118         </itemizedlist>
119       </listitem>
120     </itemizedlist>
122     <para>For our demonstration, we will use the <link
125     <programlisting role="example">
126 Order   = 2; // The order of the filter
127 Fs      = 1000; // The sampling frequency
128 Fcutoff = 40;   // The cutoff frequency
130 // We design a low pass elliptic filter
131 hz = iir(Order,'lp','ellip',[Fcutoff/Fs/2 0],[0.1 0.1]);
133 // We compute the frequency response of the filter
134 [frq,repf]=repfreq(hz,0:0.001:0.5);
135 [db_repf, phi_repf] = dbphi(repf);
137 // And plot the bode like representation of the digital filter
138 subplot(2,1,1);
139 plot2d(Fs*frq,db_repf);
140 xtitle('Obtained Frequency Response (Magnitude)');
141 subplot(2,1,2);
142 plot2d(Fs*frq,phi_repf);
143 xtitle('Obtained Frequency Response (Phase in degree)');
144  </programlisting>
146     <para>Here is the representation of the digital elliptic filter.</para>
148     <mediaobject>
149       <imageobject>
150         <imagedata align="center"
151                    fileref="../../images/digital_elliptic_filter.png" />
152       </imageobject>
153     </mediaobject>
155     <para>To represent the filter in phase and magnitude, we need first to
156     convert the discrete impulse response into magnitude and phase using the
158     a set of normalized frequencies.</para>
159   </refsection>
161   <refsection>
162     <title>Filtering a signal using the digital filter</title>
164     <para>Designing a filter is a first step. Once done, this filter will be
165     used to transform a signal. To get rid of some noise for example.</para>
167     <para>In the following examples, we will filter a gaussian noise.</para>
169     <programlisting>
170 rand('normal');
171 Input = rand(1,1000); // Produce a random gaussian noise
172 t     = 1:1000;
174 sl= tf2ss(hz); // From transfert function to syslin representation
175 y = flts(Input,sl); // Filter the signal
177 subplot(2,1,1);
178 plot(t,Input);
179 xtitle('The gaussian noise','t','y');
180 subplot(2,1,2);
181 plot(t,y);
182 xtitle('The filtered gaussian noise','t','y');
183  </programlisting>
185     <para>Here is the representation of the signal before and after
186     filtering.</para>
188     <mediaobject>
189       <imageobject>
190         <imagedata fileref="../../images/digital_filtered_noise.png" />
191       </imageobject>
192     </mediaobject>
194     <para>As we can see in the result, the high frequencies of the noise have
195     been removed and it remains only the low frequencies. The signal is still
196     noisy, but it contains mainly low frequencies.</para>
197   </refsection>
199   <refsection>
200     <title>Filtering a signal using the analog filter</title>
202     <para>To filter a signal using an analog filter, we have two
203     strategies:</para>
205     <itemizedlist>
206       <listitem>
207         <para>transform the analog filter into a discrete one using the <link
209       </listitem>
211       <listitem>
213         the signal</para>
214       </listitem>
215     </itemizedlist>
220     <programlisting>
221 rand('normal');
222 Input = rand(1,1000); // Produce a random gaussian noise
223 n     = 1:1000; // The sample index
225 eldtf = dscr(elatf,1/100); // Discretization of the linear filter
226 y = flts(Input,eldtf); // Filter the signal
228 subplot(2,1,1);
229 plot(n,Input);
230 xtitle('The gaussian noise','n','y');
231 subplot(2,1,2);
232 plot(n,y);
233 xtitle('The filtered gaussian noise','n','y');
234  </programlisting>
236     <para>Here is the representation of the signal before and after filtering
240     <mediaobject>
241       <imageobject>
242         <imagedata fileref="../../images/analog_filtered_noise.png" />
243       </imageobject>
244     </mediaobject>
248     <programlisting>
249 rand('normal');
250 Input = rand(1,1000); // Produce a random gaussian noise
251 t     = 1:1000;
252 t     = t*0.01; // Convert sample index into time steps
254 y = csim(Input,t,elatf); // Filter the signal
256 subplot(2,1,1);
257 plot(t,Input);
258 xtitle('The gaussian noise','t','y');
259 subplot(2,1,2);
260 plot(t,y);
261 xtitle('The filtered gaussian noise','t','y');
262  </programlisting>
264     <para>Here is the representation of the signal before and after filtering
267     <mediaobject>
268       <imageobject>
269         <imagedata fileref="../../images/analog_filtered_noise_csim.png" />
270       </imageobject>
271     </mediaobject>
278   </refsection>
280   <refsection>
281     <title>See Also</title>
283     <simplelist type="inline">