Circuits Syst Signal Process (2010) 29: 51–64 DOI 10.1007/s00034-009-9119-2 L O W P O W E R D I G I TA L F I LT E R S
Digital Filter Design Using Computer Algebra Systems ´ Miroslav D. Lutovac · Jelena D. Certi´ c· Ljiljana D. Mili´c
Received: 1 May 2008 / Revised: 17 October 2008 / Published online: 20 June 2009 © Birkhäuser Boston 2009
Abstract Digital filter design can be performed very efficiently using modern computer tools. The drawback of the numeric-based tools is that they usually generate a tremendous amount of numeric data, and the user might easily lose insight into the phenomenon being investigated. The computer algebra systems successfully overcome some problems encountered in the traditional numeric-only approach. In this paper, we introduce an original approach to algorithm development and digital filter design using a computer algebra system. The main result of the paper is the development of an algorithm for an infinite impulse response (IIR) filter design that, theoretically, is impossible to be implemented using the traditional approach. We present a step-by-step procedure which includes derivations of closed-form expressions for (1) the transfer functions of the implemented digital filter which contains the algebraic loop; (2) the closed-form expression for computing the number of requested iteration steps; and (3) the error function representing the difference of the output sample values of the new filter and that of the conventional filter. We demonstrate how one can use some already-known multiplierless digital filter as a start-up filter to design a new digital filter whose passband edge frequency can be simply moved by using a single parameter. As a result, we obtain a multiplierless IIR filter, which belongs to the family of low-power digital filters where multipliers are replaced with a small number of adders and shifters. This work was partially supported by the Ministry of Science of Serbia under Grant 110002 TR. M.D. Lutovac State University of Novi Pazar, Vuka Karadži´ca bb, 36300 Novi Pazar, Serbia e-mail:
[email protected] ´ c · L.D. Mili´c () J.D. Certi´ School of Electrical Engineering, University of Belgrade, Bulevar Kralja Aleksandra 73, 11000 Belgrade, Serbia e-mail:
[email protected] ´ c J.D. Certi´ e-mail:
[email protected]
52
Circuits Syst Signal Process (2010) 29: 51–64
Keywords Symbolic signal processing · Algebraic loop · Algorithm development · Modeling · Implementation · Multiplierless IIR filter 1 Introduction Traditional digital signal processing is based on numeric computations that might be inefficient to satisfy the needs of new efficient systems [7]. Symbolic techniques can be used to complement traditional numeric algorithms. They can be used to effectively search for new solutions, see e.g. [3], where a new low-sensitivity structure has been derived. Nowadays, computer programs enable a knowledge about the design of algorithms, and can be used to employ that knowledge in symbolic manipulations. The design process is the activity of taking algorithms expressed in an algorithmic form and turning them into working systems, say a computer source code. The numbers are only one of the possible symbols that the computer can handle, and a computer program can be viewed as a set of instructions for manipulating symbols. The use of computer algebra systems (CASs) and symbolic processing can help us gain insight into how some system works, although this system is preferred to experiment with numeric simulations. The current leaders in implementing CASs (Mathematica by Wolfram Research [14] and Symbolic Math Toolbox of M ATLAB and Maple [9]) are extremely powerful in doing symbolic mathematics for technical computing. A new theory called algebraic signal processing (ASP) theory is proposed based on the linear signal model defined as a triple algebra-module-mapping in which filtering, spectrum analysis, the Fourier transform, and other notions follow from the corresponding ASP concepts [12, 13], where authors explain how to design signal models from the specification of a special filter. Using a CAS [6] the transfer function matrix of a complex multiple-input multiple-output (MIMO) discrete system can be derived in terms of the system parameters kept as symbols. For a class of nonlinear systems, it was demonstrated how to use a CAS to obtain a closed-form expression of the output signal for a known symbolic stimulus [6]. The purpose of this paper is threefold: first, to demonstrate how a CAS can be used in examining the structure and performances of an existing digital filter; second, to use the CAS in finding solutions for transforming an existing infinite impulse response (IIR) half-band filter to another digital filter that retains the passband/stopband performances of the start-up half-band filter; and third, to use the CAS in order to develop a multiplierless IIR filter with an adjustable cutoff frequency. This paper is organized as follows. Following this introduction, in Sect. 2, we describe how CAS can be used for verification of a previously published digital filter, i.e., find a transfer function in terms of symbolic multiplier coefficients, plot a magnitude response for numeric values of filter parameters, and generate a computer code for processing. In Sect. 3 we introduce the digital filter with an algebraic loop and show that the passband edge can be simply controlled using a single parameter. Since the digital filter with algebraic loop cannot be realized in a traditional way, Sect. 4 is devoted to the design of a digital filter without an algebraic loop using a CAS, and to the derivation of the closed-form expressions of filter coefficients. In Sect. 5 we describe that digital filters with algebraic loops can be implemented by inserting a delay
Circuits Syst Signal Process (2010) 29: 51–64
53
element into a loop, and we present how the CAS can be used to derive the transfer function and the error function as a closed-form expression in terms of the filter parameters and of the number of iterations. Section 6 is devoted to application of the digital filter with an algebraic loop. The concluding remarks are given in Sect. 7.
2 Multiplierless Half-Band IIR Filters and Symbolic Analysis The multiplierless half-band IIR filters are useful for implementation using programmable logic devices with low-power consumption [4]. The straightforward design of a multiplierless half-band IIR filter is based on a first-order sensitivity analysis combined with the choice of the filter transfer function and the appropriate realization structure [10]. One of the unavoidable steps in deriving new solutions is the analysis of previously published filters such as that given in [10]. The verification process of published results can be time consuming, especially due to possible printing errors in papers, or when using results from a number of references. As an example, let us derive the transfer function and generate a computer source code of the half-band IIR digital filter from Example 2 reported in paper [10] using the first-order all-pass filter structure presented in Fig. 1 of paper [5] to realize four second-order filter sections. Instead of using the time-consuming manual derivations, we apply CAS [6] and [14] to draw the schematic of the filter, derive the transfer function, and generate a computer source code. Let us consider an all-pass section first. The schematic of the all-pass section can be drawn using an easy point-and-click interface (graphical user interface, GUI). When drawing a new element, the corresponding textual description is automatically generated. The textual schematic description contains all details for drawing, solving, simulating, and implementing the system [6]. As an example, Fig. 1 displays the resulting second-order section produced by the GUI when starting from the first-order section from [5], and replacing a unit delay with the two delays.
Fig. 1 Schematic of the second-order section drawn using GUI, and the automatically generated textual description
54
Circuits Syst Signal Process (2010) 29: 51–64
Fig. 2 Schematic of the ninth-order half-band filter generated using schematic of the second-order section, and the automatically derived transfer function Table 1 Coefficients of the example multiplierless half-band filter Symbol
b1
b2
b3
b4
Value
1/22 + 1/23 + 1/26
1/23 − 1/28
1 − 1/23 + 1/26
1/2 + 1/23 + 1/25 + 1/27
The transfer function is obtained using a command DiscreteSystemTransferFunction[schHB0]. It finds the output Y in terms of the input X using the textual description of the schematic, schHB0, and defined functions in the z domain of the two-input adder, the multiplier (−b), and the delay. The result (the filter transfer function) can be viewed using the traditional form b + z−2 . 1 + bz−2
(1)
The larger schematic can be generated by joining: (1) the existing schematics, say schHB0, (2) their copies with changed coordinates of elements and changed values of
the system parameters (for example multiplier b becomes multiplier b1 ), and (3) new schematics that are drawn using the GUI. Figure 2 shows the resulting schematic for the ninth-order half-band filter. In Fig. 2, the derived transfer function is expressed in terms of the symbolic values of the multipliers {b1 , b2 , b3 , b4 }, and the reserved symbol for the complex frequency z. The process of verification is completed when the magnitude response is computed for the given values. Figure 3 displays the magnitude response of the filter from [10] whose multiplier values are shown in Table 1. The magnitude response exhibits an extremely low passband ripple, and achieves the requirement of 46 dB stopband attenuation for f > 0.28. Notice that due to the coefficient quantization, the magnitude response is not strictly elliptic.
Circuits Syst Signal Process (2010) 29: 51–64
55
Fig. 3 Magnitude response of the ninth-order multiplierless half-band filter of Fig. 2 and Table 1
Fig. 4 Automatically generated computer code with the usage
The last step in verification is to generate a computer code, similar to the function filter in M ATLAB, which can be used for processing. The code is automatically evaluated when using the command that generates the code and the usage explaining the format of function inputs and outputs, see Fig. 4. The computer code has a sequence of equations that are evaluated for each input sample as shown in Fig. 4. Therefore, the schematic is not seen as a static picture. It is a symbolic object that contains all details for drawing, symbolic solving, analyzing, and implementing the filter. The knowledge embedded in the schematic object is used according to the task that is requested. In the example considered above, a larger digital filter is generated
56
Circuits Syst Signal Process (2010) 29: 51–64
when using four copies of the basic second-order filter section whose multiplier value b is set to the four different values {b1 , b2 , b3 , b4 }. 3 Digital Filter with Algebraic Loop A new lowpass filter can be derived from the start-up half-band filter when applying the lowpass-to-lowpass transformation [1, 2, 11] z−1 →
a + z−1 , 1 + az−1
a = 0.
(2)
Actually, the above lowpass-to-lowpass transformation is the first-order all-pass transfer function. According to [1, 2, 11], the parameter a is expressible in terms of the 3 dB cutoff frequency, a=−
cos(ω3 dB ) . 1 + sin(ω3 dB )
(3)
In that case a new filter is obtained, which retains the passband/stopband magnitude responses of the start-up half-band filter; only the 3 dB cutoff is moved from the position π/2 to the desired position denoted by ω3 dB . The above transformation can be implemented on the schematic, as shown in Fig. 5. Here, the second-order section has been drawn using the GUI by replacing each delay with the first-order section, whose transfer function is given in (2). The transfer function of the second-order sections is automatically derived from the schematic of Fig. 5(b), resulting in a 2 + b + (−2ba − 2a)z−1 + (ba 2 + 1)z−2 . ba 2 + (−2ba − 2a)z−1 + (a 2 + b)z−2 + 1
(4)
The advantage of the new digital filter is that all coefficients are multiplierless, the b-coefficients are the same as in [10], see Table 1, and the a-coefficient can be simply selected from the set of values in such a way that the 3 dB-cutoff frequency is from the required set of frequencies, which is determined from (3) for a being a power of two (i.e. a = ±2−p , where p is an integer), or a simple combination of powers of two [5–8]. For example, for 0.36 ≤ f 3 dB ≤ 0.37, we can choose f 3 dB ≈ 0.364 and a = 2−2 + 2−3 . However, the all-pass sections of Fig. 5 contain the algebraic loops. In numerical simulations, algebraic loops may be unsolvable. An algebraic loop generally occurs when an input port with direct feedthrough is driven by the output of the same block, either directly, or by a feedback path through other blocks with direct feedthrough [9]. When a model contains an algebraic loop, a numeric system, such as Simulink [9], calls a loop solving routine at each time step. The loop solver performs iterations to determine the solution to the problem (if it can) [9]. To solve the problem, for example, to solve an algebraic equation or constraint of the form F (z) = 0, the Simulink loop solver uses Newton’s method with weak line search and rank-one updates to a Jacobian matrix of partial derivatives. Although the method is robust, it is possible to create loops for which the loop solver does not converge
Circuits Syst Signal Process (2010) 29: 51–64
57
Fig. 5 Schematics of filters with algebraic loops (thick lines): (a) the ninth-order filter generated using schematic of ninth-order half-band filter and lowpass-to-lowpass transformation and (b) the basic second-order section
without a good initial guess for the algebraic states z. The basic idea of this paper is to show that the algebraic loop can be overcome in such a way that it becomes a part of the computer code that implements the filtering function—not to be solved by a numeric system such as Simulink. Symbolic analysis of systems is inherently immune to the problem imposed by algebraic loops. Notice that the transfer function of the filter with algebraic loops exists, but this filter cannot be implemented because of the algebraic loops. If we try to automatically generate a code using a CAS [6], a warning message is shown: DiscreteSystemImplementationEquations::algloop: The system cannot be implemented because it has algebraic loops. In the next sections, we show that using CASs the transfer function derived in (4) can be used to develop: (a) a new digital filter that has no algebraic loop; and (b) a numerical algorithm that implements a digital filter with algebraic loop using several iteration steps.
4 Digital Filter Without Algebraic Loop Let us rearrange automatically derived equation (4) in a more suitable form, H (z) =
(a 2 + b) − 2a(1 + b)z−1 + (1 + a 2 b)z−2 . (1 + a 2 b) − 2a(1 + b)z−1 + (a 2 + b)z−2
(5)
58
Circuits Syst Signal Process (2010) 29: 51–64
It is to be observed that the transfer function of a second-order all-pass section derived from the filter with the algebraic loop has the term (1 + a 2 b) instead of 1 in the denominator. Therefore, the transfer function from (5) should be rewritten in the new form: H (z) =
a 2 +b − 2a(1+b) z−1 + z−2 1+a 2 b 1+a 2 b . −1 + a 2 +b z−2 1 − 2a(1+b) z 2 2 1+a b 1+a b
(6)
The above function can be rearranged to the form of the second-order all-pass filter transfer function as used in [5, 7] and [8]: H (z) =
bnew − anew (1 + bnew )z−1 + z−2 , 1 − anew (1 + bnew )z−1 + bnew z−2
(7a)
where anew and bnew are given by bnew =
a2 + b , 1 + a2b
anew =
2a . 1 + a2
(7b)
Knowing the basic filter structure and using the GUI, a new digital filter is designed as illustrated in Fig. 6. The filter of Fig. 6 has no algebraic loop, and therefore it can be implemented. The drawback of this filter is that all b-coefficients are changed, and poles are moved towards the unit circle, i.e. bnew > b. Obviously, bnew coefficients are no more multiplierless, and coefficient quantization can produce undesired effects in the fixed point implementation.
Fig. 6 Schematics of (a) the ninth-order filter generated with the adjustable passband edge using schematic of the ninth-order half-band filter and (b) the basic second-order section
Circuits Syst Signal Process (2010) 29: 51–64
59
Fig. 7 Magnitude response of the ninth-order filter with adjusted 3 dB-cutoff frequency to f 3 dB = 0.364
The verification of the structure can be done by invoking the function that generates the transfer function directly from the textual description of the schematic, as shown in Fig. 7. This example shows that, starting from some known solution, the CAS can be used to derive a new digital filter and required filter sections including the values of the filter coefficients in terms of the known parameters. Here, the start-up solution is the multiplierless half-band filter, the desired sections are the all-pass sections of Fig. 6, and the existing known parameters are the half-band filter coefficients and the desired 3 dB-cutoff frequency. The resulting digital filter without the algebraic loop is not multiplierless, and it may fail to satisfy the filter specification after coefficient quantization.
5 Implementation of Digital Filter with Algebraic Loop In this section, we consider the possibility of implementing the all-pass second-order section with the algebraic loop, and also the overall digital filter based on this type of all-pass section. The reason is the attractive property that the filter sections share the common constant value a, which is determined by the 3 dB-cutoff frequency only, see Fig. 5 and equation (3). The remaining constants in the structure of Fig. 5 are the unchanged multiplierless constants b of the start-up half-band filter. In order to implement a digital filter with the algebraic loop, we propose a new computation procedure: We cut the loop and insert the unit delay element. Consequently, during the time between two input samples, we use a higher sampling rate and compute an approximate value at the filter output, similarly to the multirate filtering operation. Even more, there are some filter implementations in which the output of the multiplier is computed in two or more steps in order to implement multiplier coefficients with higher accuracy [4]. The processing speed can be much higher than the data rate, permitting us to perform some operations between two input samples. Let us draw a modified schematic of the basic second-order section used in filters with algebraic loops, see Fig. 8. Since we need to compute the output sample before the arrival of a new input sample, the input X is equal to the last input to the section and the inputs X d1 and X d2 are equal to the existing states of the delay elements. Therefore, this structure is a three-input three-output system. This filter section no longer has the algebraic loop,
60
Circuits Syst Signal Process (2010) 29: 51–64
Fig. 8 Schematics of the basic second-order section used in filters with algebraic loop, and automatically derived computer code
and therefore, it can be implemented. The question is: How many iterations must we process in order to achieve an acceptable error? The closed-form expression can be derived if we use the symbolic processing. That is, we generate the implementation code from the schematic and process with symbolic inputs X, X d1 , and X d2 under the assumption that the unknown state of the delay element is d. In two processing steps, we obtain two recurrence equations. We can use the CAS command Reduce to eliminate unknown state d, thus arriving at one recurrence equation. The CAS can solve recurrence relations for a known initial state, say d = 0. In this way, we derive a closed-form expression of the value at the output in terms of the symbolic input values, and of the number of iterations. The exact output value is a limiting value when the number of iterations n becomes very large. The limiting value for n → ∞ should be equal to that obtained with an alternative filter that has no algebraic loop. The difference between the exact output value and the value at the nth iteration is the error function. In this example, the error function is expressible in the form e(n) =
(b + a 2 )(−ba 2 )n−1 . 1 + ba 2
(8)
Finally, for specified acceptable error emax , and for known filter parameters a and b, the minimal number of iterations nmin can be computed using the following expression: log( emax +emax ba 2 ) b+a 2 nmin = Ceiling + 1, (9) log(a 2 b) where the function Ceiling[x] gives the smallest integer greater than or equal to x. Can we derive the transfer function of this type of filter as a closed-form expression, or can we compute it from the numeric values only? In this case we are an-
Circuits Syst Signal Process (2010) 29: 51–64
61
alyzing a second-order IIR filter section for obtaining a closed-form expression of the transfer function. Hence, we should derive the time response in terms of the current input and output signals and of the two previous values of the input and output signals. Symbolic signal processing can provide the time response. Starting from any initial time, processing three successive symbolic input samples, for unknown states of two delays, say d1 and d2 , we set three recurrence equations. Eliminating two unknown states d1 and d2 , we derive a recurrence equation with coefficients in terms of parameters a, and b, and of the number of iterations n. In this way, the closed-form expression of the transfer function is computed, n H (z) = 1 − −a 2 b ×
a 2 +b − 2a(1+b) (−1 + (−a 2 b)n )z−1 + (−1 + (−a 2 b)n )2 z−2 1+a 2 b 1+a 2 b . 2 b)n )z−1 + a 2 +b (−1 + (−a 2 b)n )2 z−2 1 − 2a(1+b) (−1 + (−a 2 2 1+a b 1+a b
(10)
Notice that the coefficients of the derived second-order section are not symmetric as in the case of the conventional all-pass filter. Using the second-order sections from (10), an overall ninth-order filter is generated. The three examples shown in Fig. 9 illustrate the effects of the number of iterations on the overall filter magnitude response. The b coefficients are those from Table 1, and the number of iterations performed in the corresponding second-order sections, for the maximal error of emax = 1/218 , is given Table 2. In spite of the fact that the second-order sections are not exact all-passes, the magnitude response after several iterations becomes almost the same as the magnitude response of the implementation without algebraic loop. It is not necessary to perform the same number of iterations in all second-order sections. The smaller number of iterations is required in sections with smaller b, see Fig. 9 and Table 2, Example (c). The final step is to convert the CAS code into the code of some other program, say M ATLAB. This is a trivial step, because the automatically generated code of CAS
Fig. 9 Magnitude response of the ninth-order filter with algebraic loops using different number of iterations per input sample as given in Table 2
62
Circuits Syst Signal Process (2010) 29: 51–64
Table 2 Number of iterations in the second-order sections b1 , b2 , b3 , and b4 from Table 1, and a = 2−2 + 2−3 Example
Section with b1
Section with b2
Section with b3
Section with b4
(a)
n1 = 8
n2 = 8
n3 = 8
n4 = 8
(b)
n1 = 7
n2 = 7
n3 = 7
n4 = 7
(c)
n1 = 5
n2 = 4
n3 = 8
n4 = 6
is displayed using traditional input form and then imported as a sequence of commands in the new M ATLAB function. The numeric values in M ATLAB agree with those computed in Mathematica. The main benefit of the implementation of filters with algebraic loops is that we can use the multiplierless start-up filter design and adjust the edge frequency by simultaneously changing only one parameter in all filter sections. There is no need for computing new filter coefficients, and problems with coefficient quantization are avoided. The drawback of the proposed approach is the increased number of additions between two input samples. From the power consumption point of view, the main benefit is the replacement of general multipliers with a small number of adders, for example b1 = 1/22 + 1/23 + 1/26 , b2 = 1/23 − 1/28 , b3 = 1 − 1/23 + 1/26 , and b4 = 1/2 + 1/23 + 1/25 + 1/27 . In all first- and second-order sections, the common coefficient a can be simply selected from the set of values in such a way that the 3 dB-cutoff frequency is from the required set of frequencies, which is determined from (3) for a being a power of two (i.e. a = ±2−p , where p is an integer), or a simple combination of powers of two [5, 8, 10]. The number of iterations depends only on the maximal error and the coefficient values (a and b); for example, n1 = 5, n2 = 4, n3 = 8, and n4 = 6, for emax = 1/218 as in Example (c) of Table 2 and Fig. 9. The shorter word lengths can be used, and thus the power consumption can be lower. Examples of Fig. 9 demonstrate that starting from a multiplierless half-band filter we obtain the multiplierless filter with an adjustable cutoff frequency. The way to deal with the algebraic loop is not quite as new. Simulink [6] uses a similar method (for example, loop solver using Newton’s method, or Algebraic Constraint blocks that specify initial guesses). But, in hardware implementations we cannot use Simulink. The number of iterations is not known in advance, and the transfer function is not known because it is a function of the number of iterations, see Fig. 9. Therefore, the iterations should be included in the implementation code. Although the iterated operations can be viewed as processing with higher sampling rate, the overall filter generates an output sample for each input sample. Alternatively, the iterated operations can be viewed as another way of computing multiplication of a sample with a constant using several operations of addition. In this way, the actual value of a constant is a function of other filter parameters and of the number of iterations, and the overall transfer function can be simply derived as in the case of a single-rate system. This paper presents the direct procedure of developing implementation code, which is how to design a filter, with an additional parameter—the appropriate number of iterations in each second-order section. The CAS is used to automate the design
Circuits Syst Signal Process (2010) 29: 51–64
63
procedure, the software Mathematica for basic manipulations with symbols, and the application package SchematicSolver [6] for GUI and some DSP functions. The design procedure is not given in Mathematica and SchematicSolver.
6 Fixed Point Implementation The proposed digital filter with the algebraic loop is compared with the filter implemented according to the approximate procedure presented in [2]. We use the start-up multiplierless ninth-order half-band filter of Figs. 2 and 5 to generate three complementary lowpass/highpass filter pairs with the 3 dB-crossover frequency placed at 0.05 (ω3 dB = 0.1π). Therefore, we consider: (1) a filter pair with the constants computed according to the exact formulae given in [11]; (2) a filter pair without algebraic loops whose constants are computed with the use of the approximate formulae derived in [2] and the processing evaluated in the fixed point arithmetic; and (3) the filter pair with the all-pass sections with the algebraic loops proposed in this paper, and the fixed point signal processing. Figure 10 displays the magnitude responses of three designs. The filter with the algebraic loop is practically insensitive to the coefficient quantization because b-coefficients are not changed after frequency transformation. The above example shows the benefits of CAS. The new digital filter is designed using a straightforward procedure and all necessary closed-form relations are derived. In this way, by specifying the number of iterations, the filter performances can achieve the desired specifications. The particular benefits are the filter structure and the proposed computation procedure that provides the change of the filter cutoff frequency (the crossover frequency in the case of a filter pair) by adjusting the value of only one parameter. This property may be efficiently used in digital filters (filter pairs) with a variable cutoff (crossover) frequency.
Fig. 10 Magnitude responses of filters in fixed point implementation
64
Circuits Syst Signal Process (2010) 29: 51–64
7 Conclusion This paper presents a vision of the role of symbolic computations in signal processing. The paper provides illustrative application examples appropriate to algorithm development, modeling, and implementation of digital filters. Digital filters were highlighted as visualized algorithms by means of blockdiagram representations. The schematic has been established as a symbolic object containing all details for drawing, symbolic solving, analyzing, and implementing the system. It has been shown how a computer algebra system (CAS) analyzes the schematic as the symbolic object and how it identifies symbolic system parameters to be used, for example, as coefficients of a new filter structure. The knowledge embedded in the schematic object has been used to generate the implementation code and to derive the transfer function. For the special type of digital filters with algebraic loops, it was demonstrated how to use the CAS to obtain the closed-form expression for the transfer function and for the output signal, and how to derive the error function in terms of the filter parameters and of the number of iteration steps. The presented examples demonstrate the possibility of the multiplierless implementation of IIR digital filters (filter pairs) with variable cutoff (crossover) frequency. Acknowledgement tions.
The authors are grateful to the anonymous reviewers for their constructive sugges-
References 1. A.C. Constantinides, Spectral transformations for digital filters. Proc. IEEE 117, 1585–1590 (1970) ´ c, L. Mili´c, Signal processor implementation of a low-pass/high-pass IIR digital filter with 2. J. Certi´ variable cut-off frequency. EUROCON 2005, 1618–1621 (2005) 3. G. Jovanovic-Dolecek, S.K. Mitra, Symbolic sensitivity analysis of the new second-order IIR structure. ING. Invest. Tecnol. IX/1, 59–65 (2008) 4. B. Lutovac, M.D. Lutovac, Design and VHDL description of multiplierless half-band IIR filter. Int. J. ˝ 56, 348–350 (2002) Electron. Commun. AEU 5. M.D. Lutovac, L.D. Mili´c, Design of computationally efficient elliptic IIR filters with a reduced number of shift-and-add operations in multipliers. IEEE Trans. Signal Process. 45, 2422–2430 (1997) 6. M.D. Lutovac, D.V. Toši´c, SchematicSolver Version 2 7. M.D. Lutovac, D.V. Toši´c, Symbolic analysis and design of control systems using Mathematica. Int. J. Control 79(11), 1368–1381 (2006). Special Issue on Symbolic Computing in Control 8. M.D. Lutovac, D.V. Toši´c, B.L. Evans, Filter Design for Signal Processing Using MATLAB and Mathematica. Prentice-Hall, New York (2001) 9. MATLAB Version 7, The MathWorks, Inc., Natick, MA, 2005 10. L. Mili´c, M. Lutovac, Design of multiplierless elliptic IIR filters with a small quantization error. IEEE Trans. Signal Process. 47, 469–479 (1999) 11. L.D. Milic, T. Saramäki, Three classes of IIR complementary filter pairs with an adjustable crossover frequency, in Proc. IEEE Int. Symp. Circuits and Systems, 2003, ISCAS’03, vol. IV, pp. 145–148 (2003) 12. M. Püschel, J. Moura, Algebraic signal processing theory: 1-D space. IEEE Trans. Signal Process. 56, 3586–3599 (2008) 13. M. Püschel, J. Moura, Algebraic signal processing theory: Foundation and 1-D time. IEEE Trans. Signal Process. 56, 3572–3585 (2008) 14. S. Wolfram, The Mathematica Book (Cambridge University Press, Cambridge, 2003)
