- Email: [email protected]
Using Methodology for MATLAB Designing the First-Order Chebyshev Analogue and IIR Digital Filters
Using Methodology for MATLAB Designing the First-Order Chebyshev Analogue and IIR Digital Filters Petr Vojcinak ∗ Martin Pies ∗∗ Jiri Koziorek ∗∗∗ ∗
Department of Cybernetics and Biomedical Engineering, VSB-TUO, Czech Republic (e-mail: [email protected]) ∗∗ Department of Cybernetics and Biomedical Engineering, VSB-TUO, Czech Republic (e-mail: [email protected]) ∗∗∗ Department of Cybernetics and Biomedical Engineering, VSB-TUO, Czech Republic (e-mail: [email protected]). Abstract: This paper deals with using some MATLAB functions and tools for designing the first-order analogue Chebyshev filters and IIR Chebyshev filters. The first part of this paper is focused on a design of an analogue filter via Chebyshev approximation approach, including features and mathematical background of this iso-extremal approximation, approximation of a normalized low-pass (NLP) filter, and some mathematical formulas for calculating its fundamental parameters, such as constructing a tolerance scheme, an order of the Chebyshev approximation, poles of the NLP’s transfer function, characteristic equation, and group delay. Due to some frequency transformation formulas implemented in MATLAB, un-normalized forms of frequency functions for low-pass (LP), high-pass (HP), band-pass (BP), and band-stop (BS) analogue filters are also available. The second part of this paper describes a design of the firstorder IIR Chebyshev filter via Filter Visualization Tool (FVT), and functions, implemented in Signal Processing Toolbox, whereas a conversion from the analogue form into the digital form is done and discussed for bilinear transformation only. Keywords: Chebyshev analogue filter, Chebyshev approximation, Chebyshev polynomials, IIR digital filter, MATLAB, Normalized low-pass filter, NLP’s characteristic equation, NLP’s group delay, NLP’s transfer function, Un-normalized filters. 1. INTRODUCTION 1.1 Chebyshev Approximation Let us consider C as a space of continuous real functions, defined in a closed interval ha, bi and having a usual norm value, thus; Rektorys K. et al. (2000a): kf k = max |f (x)| . (1) x∈ha,bi
Because nth degree polynomials evidently constitute some finite-dimensional subspace of the space C , there is this (0) nth degree polynomial Qn , related to each function f ∈ C, to satisfy this equality, thus; Rektorys K. et al. (2000b):
En (f ) ≡ f − Q(0) (2) n = inf kf − Qn k .
The infimum is taken into account for all the polynomials of nth or lower degree (Qn ). This polynomial is called as the best approximation polynomial or the Chebyshev approximation polynomial. We often call it as “minimax” approximation, because maximal error is minimal in the interval ha, bi; Rektorys K. et al. (2000b).
If the polynomial Qn is the best regular approximation polynomial, then there are (n + 2) nodes to satisfy this equality; Rektorys K. et al. (2000b):
f (xi ) − Qn (xi ) = α(−1)i kf − Qn k . Where
i α
(3)
indexer for i = 0, . . . , m for all i, it is equalled to a value of +1 or −1.
Some feature, mentioned in (3), means, that maximal and minimal values (i. e. extrema) of the function f − Qn alternately occur in nodes x0 , x1 , · · · , xm . This characteristic feature is called as Chebyshev alternate feature and corresponding group of nodes is called as Chebyshev alternate group; Rektorys K. et al. (2000b). In practical use, we have to find some nth degree polynomial (the highest-power coefficient equals to 1), which has the smallest (it is related to the space C norm) difference from zero in the interval ha, bi. If we reformulate the previous sentence, then we found the (n − 1)th degree polynomial, which is a polynomial of the best regular approximation of function xn in the interval ha, bi. In case of desired polynomial Tn (x), in the interval ha, bi a number of consecutive nodes, in which Tn acquires the value of kTn k at alternate sign, is not smaller than (n + 1). For x ∈ h−1, +1i, the desired polynomial could be written in this form, thus: Tn (x) = kTn k cos [n arccos (x)]
= 2−n+1 cos [n arccos (x)] .
(4)
Finally, for −1 = x0 < x1 < · · · < xn = +1 and k = 0, · · · , n, a form of the desired (Chebyshev) polynomial is given by (5); Rektorys K. et al. (2000b): (n − k) π Tn (xk ) = kTn k cos = (−1)n−k kTn k . (5) n 1.2 Chebyshev Polynomials Chebyshev polynomials are a sequence of classic orthogonal polynomials, which are related to de Moivre’s formula and which can be defined recursively. We usually distinguish between two kinds of the Chebyshev polynomials, thus: • the first kind Tn (x) • the second kind Un (x)
These Chebyshev polynomials are the nth degree polynomials and their sequence composes some polynomial sequence. In practical use, the Chebyshev polynomials are important in approximation theory, whereas the roots of the first kind of the Chebyshev polynomials are used in polynomial interpolation, see (1), (2), and (3). This paper is focused on describing and using mathematical background for the first kind only. Generally, both kinds of the Chebyshev polynomials are orthogonal; in case of the first kind, we assume: • interval x ∈ h−1, +1i,
• weight function w (x) = 1 − x2
− 21
.
Based on the Sturm-Liouville theory and a form of an orthonormal basis in the Hilbert space, we assume this ˇas J. et al. (1977), Rektorys K. et al. formula; Nec (2000a): Z+1 Tn (x) Tm (x) w (x) dx = Fmn (x) .
−1
m 6= n 0 2 Fmn (x) = π/2 = kTn (x)k m = n, n 6= 0. π = kT0 (x)k2 m=n=0
(6)
The Chebyshev polynomials are also the solutions of these second-order ordinary differential equations (ODE): d2 d {f (x)} − x {f (x)} + n2 f (x) = 0. (7) 1 − x2 2 dx dx Where f (x) ≡ y = Tn (x) this ODE’s solution. d2 d 1 − x2 {f (x)} − 3x {f (x)} + dx2 dx +n (n + 2) f (x) = 0.
(8)
Where f (x) ≡ y = Un (x) this ODE’s solution.
These ODEs are special cases of the Sturm-Liouville differential equations. In case of the Chebyshev polynomials’ features, we assume definition, representations (polynomial, trigonometric, and integral), orthogonality, generating functions, recurrence relation, roots, and extrema. Because the Chebyshev polynomials are special cases of ultraspherical (also Gegenbauer) polynomials, which themselves are special cases of hypergeometric (also Jacobi)
polynomials, Gegenbauer differential equation has this form, see forms of (7) and (8): d2 d 1 − x2 {f (x)} − Ax {f (x)} + Bnf (x) = 0 (9) dx2 dx Where A = 2α + 1, B = 2α + n. In case of the trigonometric representation of the first kind Chebyshev polynomials, let us consider these formulas, based on Euler’s formula (and de Moivre’s formula). For |x| ≤ 1 (cosine function), see the form of (4): Tn [x = cos (Φ)] = cos (nΦ) = cos [n arccos (x)] .
(10)
For |x| ≥ 1 (hyperbolic cosine function): Tn [x = cosh (Φ)] = cos (nΦ) = cosh [n arg cosh (x)] . (11) Equations (10) and (11) show, that different approaches to defining the Chebyshev polynomials lead to different explicit formulas – it is very useful for designing a Chebyshev filter not only in pass band, see (10), but also in stop band, see (11). 2. DESIGNING THE FIRST-ORDER ANALOGUE CHEBYSHEV FILTER The Chebyshev polynomials are also the solutions of Pell equation, which is any Diophantine equation having this form, thus: x2 − ny 2 = 1. (12) If we consider some features of the Euler’s formula, the de Moivre’s formula, and both kinds of the Chebyshev polynomials, then we can write, thus: 2 Tn2 (x) − f (x) Un−1 (x) = Tn2 (x) − 2 (13) x2 − 1 d − {T (x)} = 1. n 2 n dx After simple adjustments, we get a solution, related to (10), in this form: √ √ n n x + x2 − 1 + x − x2 − 1 Tn (x) = . (14) 2 If we consider the features of cosine function and complex exponential function, when we get this desired result, compare with (10): n n eiΦ + e−iΦ Tn [x = cos (Φ)] = = cos (nΦ) . (15) 2 Designing the first-order analogue Chebyshev filter is based on the Chebyshev approximation, which uses socalled first Chebyshev approximation method. In this case, we find some polynomial solution in open interval Ω = (−1, +1), where Ω is a normalized radial frequency. This solution has to approximate zero as best (as possible) and with regular difference; Dav´ıdek V. et al. (2006). We assume some approximate differential equation; in this case it is the Pell equation, converted into this form, compare with (13): 2 d 2 1−x {Tn (x)} = n2 1 − Tn2 (x) . (16) dx
case of this filter kind, the poles are always located on an ellipse, see Fig. 2. 2.2 Practical Example for the NLP Filter Let us consider these values of normalized radial frequencies and primary parameters, thus: Fig. 1. Normalized low-pass filter: Tolerance scheme (left; including the primary parameters) and standard tolerance scheme (right; including the secondary parameters). 2.1 Approximating the Normalized Low-pass (NLP) Filter This approximation methodology is based on the first kind Chebyshev polynomials and focused on the normalized low-pass (NLP) filter only. Due to frequency normalization, it is possible to convert ideal LP’s requirements into a NLP prototype. Naturally, it is also possible to transform requirements for other types of analogue filters (i. e. low-pass, high-pass, band-pass, and band-stop) into the NLP filter, namely via other frequency formulas and/or using MATLAB functions, implemented in Signal Processing Toolbox or Filter Visualization Tool. The standard results of this approximation approach are as follows; Dav´ıdek V. et al. (2006): ˆ NLP (p), • transfer function H • characteristic function ϕNLP (p), • group delay τg (Ω).
There are defined and shown some ideal requirements for frequency-magnitude characteristic of the ideal low-pass (LP) filter. Radial frequency ωp is an un-normalized cutoff radial frequency of designed filter’s pass-band, where magnitude values at positive radial frequencies are given by this formula, thus: 1 : ω ∈ h0, ω i ˆ (17) H (jω) LP = 0 : ω ∈ hω , p+∞) s
Analogously to the ideal low-pass filter, there are defined and shown some ideal requirements for normalized frequency-magnitude characteristic of the ideal normalized low-pass (NLP) filter. Normalized radial frequency Ωp , naturally equalled to 1, is a normalized cut-off radial frequency of designed filter’s pass-band, where magnitude values at positive normalized radial frequencies are given by this formula, thus: 1 : Ω ∈ h0, Ω i = h0, 1i ˆ (18) HNLP (jΩ) = 0 : Ω ∈ hΩ , p+∞) s
Because frequency-magnitude characteristics of the ideal NLP filter need not meet all the feasibility requirements, it is required to consider so-called tolerance scheme or standard tolerance scheme, see Fig. 1; Dav´ıdek V. et al. (2006).
At approximating the Chebyshev NLP filter, there are exact analytic formulas to determine not only values of the standard results, but also values of other essential parameters, such as an approximation order, a number of extrema, and poles of the NLP’s transfer function. In
• For the normalized radial frequencies: ap = 1 [dB], as = 20 [dB] .
• For the primary parameters:
Ωp = 1 [−], Ωs = 2.15 [−] .
(19) (20)
From these input parameters, we are able to extract and calculate some output parameters. For the approximation order; this value is also round to the nearest higher integer value: arg cosh (Ωp /k1 ) ∼ n≥ (21) = 2.1141 ⇒ n = 3. arg cosh (Ωp /k) Result of (21) shows the approximation order is odd. Now, we must consider all the n-dependent formulas for odd order only, and we must also recalculate the other parameters for new integer value of n. For the number of pass band extrema and their values at odd order, we get: √ n−1 2·1−1 3 µ= = 1 : Ω01 = cos = . (22) 2 2·3 2 Through the calculating ellipse parameters (i. e. the primary half-axis a = 0.4942 and the secondary half-axis b = 1.1154), the poles (and their placement on the ellipse, see Fig. 2) of the NLP’s transfer function are given by this formula: 2µ − 1 2µ − 1 π + jb cos π . (23) pµ = −a sin 2n 2n Where µ is indexer for µ = 1, · · · , n. In MATLAB, we can create some algorithm, which calculates these poles via (23). Some MATLAB code can be written like this, whereas the poles are located in both half-planes of the normalized complex area: a = 0.5*(Clen1 - Clen2); b = 0.5*(Clen1 + Clen2); for mu = 1:n PolyLPolorov(mu) = -a*sin((2*mu - 1)*pi/(2*n))... + 1i*b*cos((2*mu - 1)*pi/(2*n)); PolyPPolorov(mu) = a*sin((2*mu - 1)*pi/(2*n))... + 1i*b*cos((2*mu - 1)*pi/(2*n)); end PolyLPolorov(1:n) PolyPPolorov(1:n)
In case of the transfer function (at odd order), its calculated form is as follows: 0.4913 ˆ NLP (p) ∼ H . (24) = 3 2 p + 0.9883p + 1.2384p + 0.4913 Where p is normalized Laplace operator.
Fig. 2. Normalized low-pass filter: Elliptic placement of three poles of the NLP’s transfer function ˆ NLP (p = jΩ) in the left half-plane of the normalized H complex area (Σ + jΩ).
Fig. 4. Normalized low-pass filter: Normalized radial frequency-group delay characteristics of τg1 (Ω) (obˆ 1 (p = jΩ)) and τg2 (Ω) (obtained by tained by H ˆ 2 (p = jΩ)); the Cartesian coordinates. H In case of the group delay (at odd order and considering the generalized approach, related to the transfer function and its first derivation), its calculated form is as follows (in radians): τg (Ω) ∼ =
0.9883Ω4 − 0.2500Ω2 + 0.6084 . Ω6 − 1.5001Ω4 + 0.5625Ω2 + 0.2414
(28)
Due to MATLAB possibilities and features, it is possible to plot not only normalized radial frequency-magnitude characteristic, but also other characteristics, such as real part and imaginary part of the transfer function, phase characteristic, and group delay in 3D graphs. The MATLAB code can be written like this:
Fig. 3. Normalized low-pass filter: Normalized radial ˆ frequency-magnitude characteristics of H1 (p = jΩ) ˆ (obtained analytically) and H 2 (p = jΩ) (obtained by implemented MATLAB function of cheb(n,Rp)); the Cartesian coordinates and tolerance scheme limits. In case of the normalized radial frequency-magnitude characteristic (at odd order and p = jΩ), its calculated form is as follows: 0.4913 ˆ .(25) HNLP (jΩ) ∼ =√ 6 4 Ω − 1.5001Ω + 0.5625Ω2 + 0.2414
In case of the normalized radial frequency-phase characteristic (at odd order and considering the first quadrant), its calculated form is as follows: Ω3 − 1.2384Ω ϕNLP (Ω) = arctg . (26) 0.4913 − 0.9883Ω2
In case of the characteristic function (based on the characteristic equation), its calculated form is as follows (at odd order n = 3, p = jΩ and j 3 = −j): φNLP (jΩ) ⇒ j 3 T3 (Ω) = 4(jΩ)3 − 3 j 3 Ω . (27)
B0 = 0.4913; A0 = B0; A1 = 1.2384; A2 = 0.9883; A3 = 1; dSigma = 0.01; dOmega = dSigma; [Delta, Sigma] = meshgrid(-1:dSigma:0, -2:dOmega:2); s = Delta + 1i*Sigma; FceZPrenos = B0./(A3*s.^3 + A2*s.^2 + A1*s + A0); FceZSkupinove = - real((-3*s.^2 - 1.9766*s - 1.2384)./... (s.^3 + 0.9883*s.^2 + 1.2384*s + 0.4913)); Imag = imag(FceZPrenos); Real = real(FceZPrenos); mesh(Delta, Sigma, unwrap(atan2(Imag, Real)));
This MATLAB code plots the 3D graph of the normalized radial frequency-phase characteristic at unwrapped mode. For this transfer function form, mentioned in (24), it is also possible to plot it in 3D graph. 2.3 Practical Example for Un-normalized Analogue Filters Due to MATLAB functions, implemented in Signal Processing Toolbox, it is possible to convert the normalized low-pass filter into the un-normalized analogue filters, see summarization in table 1.
Fig. 5. Normalized low-pass filter: 3D graph of the normalized radial frequency-magnitude characteristic, where three poles of the transfer function are shown in the left half-plane of the normalized complex area (Σ + jΩ); the Cartesian coordinates. Where
n Rp b a Wo Bw
approx. order pass-band ripple TF nominator TF denominator cut-off rad. freq. freq. bandwidth
n, Rp = 10−0.05ap , ˆNLP (p), B AˆNLP (p), ω0 , B ≡ Bm .
In case of Wo, this is not only the cut-off radial frequency (in case of the LP and HP filters), but also a central radial frequency (in case of the BP and BS filters), given by this formula (geometric average): p √ ω0 = ωmd ωmh = 2π fmd fmh . (29) For example, let us consider these values of un-normalized radial frequencies (in radians per second):
Fig. 6. Un-normalized low-pass filter: Radial frequencyˆ LP (s = jω); the Cartesian magnitude characteristic H coordinates. ω0 = 2000π rad · s−1 ⇒ f0 = 1000 [Hz] .
Where
HP
Filter Type NLP LP HP
BP
BS
Conversion Formulas Ωp = 1 Ωs ω p = ω 0 Ωp = ω 0 ω s = ω 0 Ωs ωp = ω0 /Ωp = ω0 ωs = ω 0 /Ωs 2 / ω B Ωs = ω02 − ωpd pd m Bm = |ωmh − ωmd | / (2π) Bp = ωph − ωpd / (2π) ωmh < ωph ωmd > ωpd Bm < Bp 2 Ωs = ωpd Bm / ω02 − ωpd Bm = |ωmh − ωmd | / (2π) Bp = ωph − ωpd / (2π) ωmh > ωph ωmd < ωpd Bm > Bp
MATLAB Function che1abp(n,Rp) lp2lp(b,a,Wo) lp2hp(b,a,Wo)
ωp = ω0 Ωp = ω0 = 2000π, ωs = ω0 Ωs = 2.15ω0 = 4300π, ωp = ω0 /Ωp = ω0 = 2000π, ωs = ω0 /Ωs = ω0 /2.15 ∼ = 930.23π.
For designed LP filter, the radial frequency-magnitude characteristic and corresponding MATLAB code are as follows, thus: [NumHDPs, DenHDPs] = lp2lp(NumHs, DenHs, omega0); HDPs = tf(NumHDPs, DenHDPs) [HDP, omegaDP] = freqs(NumHDPs, DenHDPs, omega); plot(omegaDP, abs(HDP))
For the BP filter at ω0 = 2000π rad · s−1 :
For the LP and HP filters: Table 1. Summarization of Conversion Formulas and Corresponding Callings of MATLAB ˇek J. (2009),THE MATHFunctions; Dolec WORKS Inc. (2011)
LP
(30)
ωmd = 1500π, ωmh = 2500π . ωpd = 1200π, ωph = 2800π Where
(31)
Bm = |ωmh − ωmd | / (2π) = 500,
Bp = |ωph − ωpd | / (2π) = 800, 2 Ωp = ω02 − ωmd / (ωmd Bm ) ∼ = 1.17, 2 2 Ωs = ω0 − ωpd / (ωpd Bm ) ∼ = 2.13.
For designed BP filter, the radial frequency-magnitude characteristic and corresponding MATLAB code are as follows, thus: [NumHPPs, DenHPPs] = lp2bp(NumHs, DenHs, omega0, Bomega); HPPs = tf(NumHPPs, DenHPPs)
lp2bp(b,a,Wo,Bw)
lp2bs(b,a,Wo,Bw)
[HPP, omegaPP] = freqs(NumHPPs, DenHPPs, omega); plot(omegaPP, abs(HPP))
Finally, for the BS filter at ω0 = 2000π rad · s−1 : ωpd = 1500π, ωph = 2500π . ωmd = 1200π, ωmh = 2800π
(32)
ˆ z H
−1
−1 2 · π 1 − z ˆ s → . =H −1 f tg fps π 1 + z
(35)
For example, let us consider an un-normalized IIR low−1 pass filter with ω = ω = 2000π rad · s and ωs = 0 p 8000π rad · s−1 . Then radial frequency ratio equals to: pω =
ω0 1 2000π = = 0.25. = ωs 8000π 4
(36)
For this designed filter, the “normalized” radial frequencymagnitude characteristic and corresponding MATLAB code are as follows, thus: For the analogue NLP filter (s-plane): Fig. 7. Un-normalized band-pass filter: Radial frequencyˆ BP (s = jω); the Cartesian magnitude characteristic H coordinates. 3. DESIGNING THE FIRST-ORDER DIGITAL (IIR) CHEBYSHEV FILTER Transfer function of IIR filter (IIR - Infinite Impulse Response) is given by this formula; Dav´ıdek V. et al. (2006):
ˆ z H
−1
=
M P
b0
bm z −m
m=0 N P
1+
n=0
= an z −n
M Q
µ=1 N Q
ν=1
1 − z0µ z −1
. (33)
(1 − z∞ν z −1 )
Designing the IIR filters is based on finding the coefficients ˆ z −1 and denominator Aˆ z −1 of this of nominator B ˆ z −1 transfer function, see (33), or finding zeros of B and poles of Aˆ z −1 to satisfy requirements of some tolerance scheme, whereas this design process is related to the quadrature of magnitude characteristic, and the group delay. Traditional approach is based on acceptance of analogue filters’ approximations, because some methodology of the NLP approximation has been developed yet. Well, this paper tries to demonstrate this idea. In MATLAB, there are two ways, how to approach bilinear transformation, thus: • bilinear(num, den, fs) - prewarped mode is not used to indicate “match” frequency, where num is a nominator of transfer function, den is a denominator of transfer function, and fs is a sampling frequency; THE MATHWORKS Inc. (2011): −1 2 1 − z −1 ˆ z −1 = H ˆ s → 2fs 1 − z = . (34) H 1 + z −1 Ts 1 + z −1 • bilinear(num, den, fs, fp) - fp parameter specifies prewarping, which indicates “match” frequency, for which frequency responses (before and after mapping) match exactly; this parameter is identical with edge of analogue filter’s pass-band; THE MATHWORKS Inc. (2011):
[NumHNDPs, DenHNDPs] = cheby1(n, ap, Omegap, ’low’, ’s’); [HNDPs, omegaNDPs] = freqs(NumHNDPs, DenHNDPs);
For the digital NLP filter (z-plane; IIR; prewarped mode is not considered): [NumHNDPz, DenHNDPz] = cheby1(n, ap,... 2*atan(Omegap/(2*fvz))/pi, ’low’, ’z’); [HNDPz, omegaNDPz] = freqz(NumHNDPz, DenHNDPz);
For the digital LP filter (z-plane; IIR; prewarped mode is not considered): [NumHDP1z, DenHDP1z] = cheby1(n, ap,... 2*atan(omega0*Omegap/(2*fvz))/pi, ’low’, ’z’); [HDP1z, omegaDP1z] = freqz(NumHDP1z, DenHDP1z);
Or (it is an equivalent code; using bilinear transformation, prewarped mode is not considered): [NumHDP2z, DenHDP2z] = bilinear(NumHDPs, DenHDPs, fvz); [HDP2z, omegaDP2z] = freqz(NumHDP2z, DenHDP2z);
For the digital LP filter (z-plane; IIR; prewarped mode is considered): [NumHDP3z, DenHDP3z] = cheby1(n, ap,... 2*atan(omega0*Omegap*tan(pi*f0/fvz)/(2*pi*f0))/pi,... ’low’, ’z’); [HDP3z, omegaDP3z] = freqz(NumHDP3z, DenHDP3z);
Or (it is an equivalent code; using bilinear transformation, prewarped mode is not considered): [NumHDP4z, DenHDP4z] = bilinear(NumHDPs, DenHDPs, fvz, f0); [HDP4z, omegaDP4z] = freqz(NumHDP4z, DenHDP4z);
Fig. 8 shows the “normalized” radial frequency-magnitude characteristics of the designed IIR low-pass filter. For its detailed analysis (e. g. phase delay, group delay, TF nominator, TF denominator, TF poles, and structure), we can use a MATLAB tool called Filter Visualization Tool (fvtool command). In this case, “normalized” means dividing the interval h0, πi by π to get the interval h0, 1i. If we consider the prewarped mode and sampling radial frequency ωs = 2π, then we get (after normalization): ωs = ω0 /pω = 2π ⇒ ω0 /π = pω ωs /π = 0.5.
(37)
ter), s-plane (un-normalized analogue filters), and z-plane (IIR filters; normalized digital low-pass filter and unnormalized digital filters). Designing the NLP filter has five important results (approximation order, characteristic equation, characteristic function, transfer function, and group delay), and other results in form of tolerance scheme, magnitude characteristics, phase characteristics, TF poles, roots, and extrema. All these formulas are used to solve practical examples, which include calculations, 2D graphs, 3D graphs, and some fragments of corresponding MATLAB codes. Due to MATLAB&Simulink possibilities, it is also possible to implement this approach, based on Dav´ıdek V. et al. (2006), not only in Simulink, but also in the field of FIR filters or other kinds of analogue filters. ACKNOWLEDGEMENTS Fig. 8. Un-normalized IIR low-pass filter: “Normalized” radial frequency-magnitude characteristics with the prewarped mode (light curve) and without this mode (dark curve); the Cartesian coordinates. In comparison with the non-prewarped mode, value of ω0 equals to 0.42. So, result of (36) is not considered. 4. CONCLUSIONS This paper deals with using some mathematical background - the Chebyshev approximation, Chebyshev polynomials, selected MATLAB functions from Signal Processing Toolbox, and other tool (e. g. Filter Visualization Tool) to design the first-order (also kind I) analogue Chebyshev filters and IIR Chebyshev filters. Because this paper is mainly centred on the mathematic viewpoint, there are lots of fundamental (maybe essential) formulas (de Moivre, Euler, also Jacobi, Gegenbauer, or Rodrigues). Features of the Chebyshev polynomials are also mentioned and discussed in forms of definition, representations, roots, and extrema. The second part of this paper is focused on designing the first-order Chebyshev filters in p-plane (the NLP fil-
This work is supported by project SP2012/111, named “Data Acquisition and Processing from Large Distributed Systems II” of Student Grant Agency (VSB - Technical University of Ostrava). REFERENCES ˇek M. (2006). AnaDav´ıdek V., Laipert M., and Vlc ˇ logov´e a ˇc´ıslicov´e filtry. CVUT. ˇek J. (2009). Modern´ı uˇcebnice elektroniky 6. d´ıl Dolec Kmitoˇctov´e filtry, gener´ atory sign´ al˚ u a pˇrevodn´ıky dat. BEN - technick´a literatura. ˇas J. et al. (1977). Aplikovan´ Nec a matematika I (A aˇz L). SNTL, Prague. Rektorys K. et al. (2000a). Pˇrehled uˇzit´e matematiky I. Prometheus. Rektorys K. et al. (2000b). Pˇrehled uˇzit´e matematiky II. Prometheus. THE MATHWORKS Inc. (2011). Signal processing toolboxTM : User’s guide. Version 6.16 (Release 2011b), offline.
Department of Cybernetics and Biomedical Engineering, VSB-TUO, Czech Republic (e-mail: [email protected]) ∗∗ Department of Cybernetics and Biomedical Engineering, VSB-TUO, Czech Republic (e-mail: [email protected]) ∗∗∗ Department of Cybernetics and Biomedical Engineering, VSB-TUO, Czech Republic (e-mail: [email protected]). Abstract: This paper deals with using some MATLAB functions and tools for designing the first-order analogue Chebyshev filters and IIR Chebyshev filters. The first part of this paper is focused on a design of an analogue filter via Chebyshev approximation approach, including features and mathematical background of this iso-extremal approximation, approximation of a normalized low-pass (NLP) filter, and some mathematical formulas for calculating its fundamental parameters, such as constructing a tolerance scheme, an order of the Chebyshev approximation, poles of the NLP’s transfer function, characteristic equation, and group delay. Due to some frequency transformation formulas implemented in MATLAB, un-normalized forms of frequency functions for low-pass (LP), high-pass (HP), band-pass (BP), and band-stop (BS) analogue filters are also available. The second part of this paper describes a design of the firstorder IIR Chebyshev filter via Filter Visualization Tool (FVT), and functions, implemented in Signal Processing Toolbox, whereas a conversion from the analogue form into the digital form is done and discussed for bilinear transformation only. Keywords: Chebyshev analogue filter, Chebyshev approximation, Chebyshev polynomials, IIR digital filter, MATLAB, Normalized low-pass filter, NLP’s characteristic equation, NLP’s group delay, NLP’s transfer function, Un-normalized filters. 1. INTRODUCTION 1.1 Chebyshev Approximation Let us consider C as a space of continuous real functions, defined in a closed interval ha, bi and having a usual norm value, thus; Rektorys K. et al. (2000a): kf k = max |f (x)| . (1) x∈ha,bi
Because nth degree polynomials evidently constitute some finite-dimensional subspace of the space C , there is this (0) nth degree polynomial Qn , related to each function f ∈ C, to satisfy this equality, thus; Rektorys K. et al. (2000b):
En (f ) ≡ f − Q(0) (2) n = inf kf − Qn k .
The infimum is taken into account for all the polynomials of nth or lower degree (Qn ). This polynomial is called as the best approximation polynomial or the Chebyshev approximation polynomial. We often call it as “minimax” approximation, because maximal error is minimal in the interval ha, bi; Rektorys K. et al. (2000b).
If the polynomial Qn is the best regular approximation polynomial, then there are (n + 2) nodes to satisfy this equality; Rektorys K. et al. (2000b):
f (xi ) − Qn (xi ) = α(−1)i kf − Qn k . Where
i α
(3)
indexer for i = 0, . . . , m for all i, it is equalled to a value of +1 or −1.
Some feature, mentioned in (3), means, that maximal and minimal values (i. e. extrema) of the function f − Qn alternately occur in nodes x0 , x1 , · · · , xm . This characteristic feature is called as Chebyshev alternate feature and corresponding group of nodes is called as Chebyshev alternate group; Rektorys K. et al. (2000b). In practical use, we have to find some nth degree polynomial (the highest-power coefficient equals to 1), which has the smallest (it is related to the space C norm) difference from zero in the interval ha, bi. If we reformulate the previous sentence, then we found the (n − 1)th degree polynomial, which is a polynomial of the best regular approximation of function xn in the interval ha, bi. In case of desired polynomial Tn (x), in the interval ha, bi a number of consecutive nodes, in which Tn acquires the value of kTn k at alternate sign, is not smaller than (n + 1). For x ∈ h−1, +1i, the desired polynomial could be written in this form, thus: Tn (x) = kTn k cos [n arccos (x)]
= 2−n+1 cos [n arccos (x)] .
(4)
Finally, for −1 = x0 < x1 < · · · < xn = +1 and k = 0, · · · , n, a form of the desired (Chebyshev) polynomial is given by (5); Rektorys K. et al. (2000b): (n − k) π Tn (xk ) = kTn k cos = (−1)n−k kTn k . (5) n 1.2 Chebyshev Polynomials Chebyshev polynomials are a sequence of classic orthogonal polynomials, which are related to de Moivre’s formula and which can be defined recursively. We usually distinguish between two kinds of the Chebyshev polynomials, thus: • the first kind Tn (x) • the second kind Un (x)
These Chebyshev polynomials are the nth degree polynomials and their sequence composes some polynomial sequence. In practical use, the Chebyshev polynomials are important in approximation theory, whereas the roots of the first kind of the Chebyshev polynomials are used in polynomial interpolation, see (1), (2), and (3). This paper is focused on describing and using mathematical background for the first kind only. Generally, both kinds of the Chebyshev polynomials are orthogonal; in case of the first kind, we assume: • interval x ∈ h−1, +1i,
• weight function w (x) = 1 − x2
− 21
.
Based on the Sturm-Liouville theory and a form of an orthonormal basis in the Hilbert space, we assume this ˇas J. et al. (1977), Rektorys K. et al. formula; Nec (2000a): Z+1 Tn (x) Tm (x) w (x) dx = Fmn (x) .
−1
m 6= n 0 2 Fmn (x) = π/2 = kTn (x)k m = n, n 6= 0. π = kT0 (x)k2 m=n=0
(6)
The Chebyshev polynomials are also the solutions of these second-order ordinary differential equations (ODE): d2 d {f (x)} − x {f (x)} + n2 f (x) = 0. (7) 1 − x2 2 dx dx Where f (x) ≡ y = Tn (x) this ODE’s solution. d2 d 1 − x2 {f (x)} − 3x {f (x)} + dx2 dx +n (n + 2) f (x) = 0.
(8)
Where f (x) ≡ y = Un (x) this ODE’s solution.
These ODEs are special cases of the Sturm-Liouville differential equations. In case of the Chebyshev polynomials’ features, we assume definition, representations (polynomial, trigonometric, and integral), orthogonality, generating functions, recurrence relation, roots, and extrema. Because the Chebyshev polynomials are special cases of ultraspherical (also Gegenbauer) polynomials, which themselves are special cases of hypergeometric (also Jacobi)
polynomials, Gegenbauer differential equation has this form, see forms of (7) and (8): d2 d 1 − x2 {f (x)} − Ax {f (x)} + Bnf (x) = 0 (9) dx2 dx Where A = 2α + 1, B = 2α + n. In case of the trigonometric representation of the first kind Chebyshev polynomials, let us consider these formulas, based on Euler’s formula (and de Moivre’s formula). For |x| ≤ 1 (cosine function), see the form of (4): Tn [x = cos (Φ)] = cos (nΦ) = cos [n arccos (x)] .
(10)
For |x| ≥ 1 (hyperbolic cosine function): Tn [x = cosh (Φ)] = cos (nΦ) = cosh [n arg cosh (x)] . (11) Equations (10) and (11) show, that different approaches to defining the Chebyshev polynomials lead to different explicit formulas – it is very useful for designing a Chebyshev filter not only in pass band, see (10), but also in stop band, see (11). 2. DESIGNING THE FIRST-ORDER ANALOGUE CHEBYSHEV FILTER The Chebyshev polynomials are also the solutions of Pell equation, which is any Diophantine equation having this form, thus: x2 − ny 2 = 1. (12) If we consider some features of the Euler’s formula, the de Moivre’s formula, and both kinds of the Chebyshev polynomials, then we can write, thus: 2 Tn2 (x) − f (x) Un−1 (x) = Tn2 (x) − 2 (13) x2 − 1 d − {T (x)} = 1. n 2 n dx After simple adjustments, we get a solution, related to (10), in this form: √ √ n n x + x2 − 1 + x − x2 − 1 Tn (x) = . (14) 2 If we consider the features of cosine function and complex exponential function, when we get this desired result, compare with (10): n n eiΦ + e−iΦ Tn [x = cos (Φ)] = = cos (nΦ) . (15) 2 Designing the first-order analogue Chebyshev filter is based on the Chebyshev approximation, which uses socalled first Chebyshev approximation method. In this case, we find some polynomial solution in open interval Ω = (−1, +1), where Ω is a normalized radial frequency. This solution has to approximate zero as best (as possible) and with regular difference; Dav´ıdek V. et al. (2006). We assume some approximate differential equation; in this case it is the Pell equation, converted into this form, compare with (13): 2 d 2 1−x {Tn (x)} = n2 1 − Tn2 (x) . (16) dx
case of this filter kind, the poles are always located on an ellipse, see Fig. 2. 2.2 Practical Example for the NLP Filter Let us consider these values of normalized radial frequencies and primary parameters, thus: Fig. 1. Normalized low-pass filter: Tolerance scheme (left; including the primary parameters) and standard tolerance scheme (right; including the secondary parameters). 2.1 Approximating the Normalized Low-pass (NLP) Filter This approximation methodology is based on the first kind Chebyshev polynomials and focused on the normalized low-pass (NLP) filter only. Due to frequency normalization, it is possible to convert ideal LP’s requirements into a NLP prototype. Naturally, it is also possible to transform requirements for other types of analogue filters (i. e. low-pass, high-pass, band-pass, and band-stop) into the NLP filter, namely via other frequency formulas and/or using MATLAB functions, implemented in Signal Processing Toolbox or Filter Visualization Tool. The standard results of this approximation approach are as follows; Dav´ıdek V. et al. (2006): ˆ NLP (p), • transfer function H • characteristic function ϕNLP (p), • group delay τg (Ω).
There are defined and shown some ideal requirements for frequency-magnitude characteristic of the ideal low-pass (LP) filter. Radial frequency ωp is an un-normalized cutoff radial frequency of designed filter’s pass-band, where magnitude values at positive radial frequencies are given by this formula, thus: 1 : ω ∈ h0, ω i ˆ (17) H (jω) LP = 0 : ω ∈ hω , p+∞) s
Analogously to the ideal low-pass filter, there are defined and shown some ideal requirements for normalized frequency-magnitude characteristic of the ideal normalized low-pass (NLP) filter. Normalized radial frequency Ωp , naturally equalled to 1, is a normalized cut-off radial frequency of designed filter’s pass-band, where magnitude values at positive normalized radial frequencies are given by this formula, thus: 1 : Ω ∈ h0, Ω i = h0, 1i ˆ (18) HNLP (jΩ) = 0 : Ω ∈ hΩ , p+∞) s
Because frequency-magnitude characteristics of the ideal NLP filter need not meet all the feasibility requirements, it is required to consider so-called tolerance scheme or standard tolerance scheme, see Fig. 1; Dav´ıdek V. et al. (2006).
At approximating the Chebyshev NLP filter, there are exact analytic formulas to determine not only values of the standard results, but also values of other essential parameters, such as an approximation order, a number of extrema, and poles of the NLP’s transfer function. In
• For the normalized radial frequencies: ap = 1 [dB], as = 20 [dB] .
• For the primary parameters:
Ωp = 1 [−], Ωs = 2.15 [−] .
(19) (20)
From these input parameters, we are able to extract and calculate some output parameters. For the approximation order; this value is also round to the nearest higher integer value: arg cosh (Ωp /k1 ) ∼ n≥ (21) = 2.1141 ⇒ n = 3. arg cosh (Ωp /k) Result of (21) shows the approximation order is odd. Now, we must consider all the n-dependent formulas for odd order only, and we must also recalculate the other parameters for new integer value of n. For the number of pass band extrema and their values at odd order, we get: √ n−1 2·1−1 3 µ= = 1 : Ω01 = cos = . (22) 2 2·3 2 Through the calculating ellipse parameters (i. e. the primary half-axis a = 0.4942 and the secondary half-axis b = 1.1154), the poles (and their placement on the ellipse, see Fig. 2) of the NLP’s transfer function are given by this formula: 2µ − 1 2µ − 1 π + jb cos π . (23) pµ = −a sin 2n 2n Where µ is indexer for µ = 1, · · · , n. In MATLAB, we can create some algorithm, which calculates these poles via (23). Some MATLAB code can be written like this, whereas the poles are located in both half-planes of the normalized complex area: a = 0.5*(Clen1 - Clen2); b = 0.5*(Clen1 + Clen2); for mu = 1:n PolyLPolorov(mu) = -a*sin((2*mu - 1)*pi/(2*n))... + 1i*b*cos((2*mu - 1)*pi/(2*n)); PolyPPolorov(mu) = a*sin((2*mu - 1)*pi/(2*n))... + 1i*b*cos((2*mu - 1)*pi/(2*n)); end PolyLPolorov(1:n) PolyPPolorov(1:n)
In case of the transfer function (at odd order), its calculated form is as follows: 0.4913 ˆ NLP (p) ∼ H . (24) = 3 2 p + 0.9883p + 1.2384p + 0.4913 Where p is normalized Laplace operator.
Fig. 2. Normalized low-pass filter: Elliptic placement of three poles of the NLP’s transfer function ˆ NLP (p = jΩ) in the left half-plane of the normalized H complex area (Σ + jΩ).
Fig. 4. Normalized low-pass filter: Normalized radial frequency-group delay characteristics of τg1 (Ω) (obˆ 1 (p = jΩ)) and τg2 (Ω) (obtained by tained by H ˆ 2 (p = jΩ)); the Cartesian coordinates. H In case of the group delay (at odd order and considering the generalized approach, related to the transfer function and its first derivation), its calculated form is as follows (in radians): τg (Ω) ∼ =
0.9883Ω4 − 0.2500Ω2 + 0.6084 . Ω6 − 1.5001Ω4 + 0.5625Ω2 + 0.2414
(28)
Due to MATLAB possibilities and features, it is possible to plot not only normalized radial frequency-magnitude characteristic, but also other characteristics, such as real part and imaginary part of the transfer function, phase characteristic, and group delay in 3D graphs. The MATLAB code can be written like this:
Fig. 3. Normalized low-pass filter: Normalized radial ˆ frequency-magnitude characteristics of H1 (p = jΩ) ˆ (obtained analytically) and H 2 (p = jΩ) (obtained by implemented MATLAB function of cheb(n,Rp)); the Cartesian coordinates and tolerance scheme limits. In case of the normalized radial frequency-magnitude characteristic (at odd order and p = jΩ), its calculated form is as follows: 0.4913 ˆ .(25) HNLP (jΩ) ∼ =√ 6 4 Ω − 1.5001Ω + 0.5625Ω2 + 0.2414
In case of the normalized radial frequency-phase characteristic (at odd order and considering the first quadrant), its calculated form is as follows: Ω3 − 1.2384Ω ϕNLP (Ω) = arctg . (26) 0.4913 − 0.9883Ω2
In case of the characteristic function (based on the characteristic equation), its calculated form is as follows (at odd order n = 3, p = jΩ and j 3 = −j): φNLP (jΩ) ⇒ j 3 T3 (Ω) = 4(jΩ)3 − 3 j 3 Ω . (27)
B0 = 0.4913; A0 = B0; A1 = 1.2384; A2 = 0.9883; A3 = 1; dSigma = 0.01; dOmega = dSigma; [Delta, Sigma] = meshgrid(-1:dSigma:0, -2:dOmega:2); s = Delta + 1i*Sigma; FceZPrenos = B0./(A3*s.^3 + A2*s.^2 + A1*s + A0); FceZSkupinove = - real((-3*s.^2 - 1.9766*s - 1.2384)./... (s.^3 + 0.9883*s.^2 + 1.2384*s + 0.4913)); Imag = imag(FceZPrenos); Real = real(FceZPrenos); mesh(Delta, Sigma, unwrap(atan2(Imag, Real)));
This MATLAB code plots the 3D graph of the normalized radial frequency-phase characteristic at unwrapped mode. For this transfer function form, mentioned in (24), it is also possible to plot it in 3D graph. 2.3 Practical Example for Un-normalized Analogue Filters Due to MATLAB functions, implemented in Signal Processing Toolbox, it is possible to convert the normalized low-pass filter into the un-normalized analogue filters, see summarization in table 1.
Fig. 5. Normalized low-pass filter: 3D graph of the normalized radial frequency-magnitude characteristic, where three poles of the transfer function are shown in the left half-plane of the normalized complex area (Σ + jΩ); the Cartesian coordinates. Where
n Rp b a Wo Bw
approx. order pass-band ripple TF nominator TF denominator cut-off rad. freq. freq. bandwidth
n, Rp = 10−0.05ap , ˆNLP (p), B AˆNLP (p), ω0 , B ≡ Bm .
In case of Wo, this is not only the cut-off radial frequency (in case of the LP and HP filters), but also a central radial frequency (in case of the BP and BS filters), given by this formula (geometric average): p √ ω0 = ωmd ωmh = 2π fmd fmh . (29) For example, let us consider these values of un-normalized radial frequencies (in radians per second):
Fig. 6. Un-normalized low-pass filter: Radial frequencyˆ LP (s = jω); the Cartesian magnitude characteristic H coordinates. ω0 = 2000π rad · s−1 ⇒ f0 = 1000 [Hz] .
Where
HP
Filter Type NLP LP HP
BP
BS
Conversion Formulas Ωp = 1 Ωs ω p = ω 0 Ωp = ω 0 ω s = ω 0 Ωs ωp = ω0 /Ωp = ω0 ωs = ω 0 /Ωs 2 / ω B Ωs = ω02 − ωpd pd m Bm = |ωmh − ωmd | / (2π) Bp = ωph − ωpd / (2π) ωmh < ωph ωmd > ωpd Bm < Bp 2 Ωs = ωpd Bm / ω02 − ωpd Bm = |ωmh − ωmd | / (2π) Bp = ωph − ωpd / (2π) ωmh > ωph ωmd < ωpd Bm > Bp
MATLAB Function che1abp(n,Rp) lp2lp(b,a,Wo) lp2hp(b,a,Wo)
ωp = ω0 Ωp = ω0 = 2000π, ωs = ω0 Ωs = 2.15ω0 = 4300π, ωp = ω0 /Ωp = ω0 = 2000π, ωs = ω0 /Ωs = ω0 /2.15 ∼ = 930.23π.
For designed LP filter, the radial frequency-magnitude characteristic and corresponding MATLAB code are as follows, thus: [NumHDPs, DenHDPs] = lp2lp(NumHs, DenHs, omega0); HDPs = tf(NumHDPs, DenHDPs) [HDP, omegaDP] = freqs(NumHDPs, DenHDPs, omega); plot(omegaDP, abs(HDP))
For the BP filter at ω0 = 2000π rad · s−1 :
For the LP and HP filters: Table 1. Summarization of Conversion Formulas and Corresponding Callings of MATLAB ˇek J. (2009),THE MATHFunctions; Dolec WORKS Inc. (2011)
LP
(30)
ωmd = 1500π, ωmh = 2500π . ωpd = 1200π, ωph = 2800π Where
(31)
Bm = |ωmh − ωmd | / (2π) = 500,
Bp = |ωph − ωpd | / (2π) = 800, 2 Ωp = ω02 − ωmd / (ωmd Bm ) ∼ = 1.17, 2 2 Ωs = ω0 − ωpd / (ωpd Bm ) ∼ = 2.13.
For designed BP filter, the radial frequency-magnitude characteristic and corresponding MATLAB code are as follows, thus: [NumHPPs, DenHPPs] = lp2bp(NumHs, DenHs, omega0, Bomega); HPPs = tf(NumHPPs, DenHPPs)
lp2bp(b,a,Wo,Bw)
lp2bs(b,a,Wo,Bw)
[HPP, omegaPP] = freqs(NumHPPs, DenHPPs, omega); plot(omegaPP, abs(HPP))
Finally, for the BS filter at ω0 = 2000π rad · s−1 : ωpd = 1500π, ωph = 2500π . ωmd = 1200π, ωmh = 2800π
(32)
ˆ z H
−1
−1 2 · π 1 − z ˆ s → . =H −1 f tg fps π 1 + z
(35)
For example, let us consider an un-normalized IIR low−1 pass filter with ω = ω = 2000π rad · s and ωs = 0 p 8000π rad · s−1 . Then radial frequency ratio equals to: pω =
ω0 1 2000π = = 0.25. = ωs 8000π 4
(36)
For this designed filter, the “normalized” radial frequencymagnitude characteristic and corresponding MATLAB code are as follows, thus: For the analogue NLP filter (s-plane): Fig. 7. Un-normalized band-pass filter: Radial frequencyˆ BP (s = jω); the Cartesian magnitude characteristic H coordinates. 3. DESIGNING THE FIRST-ORDER DIGITAL (IIR) CHEBYSHEV FILTER Transfer function of IIR filter (IIR - Infinite Impulse Response) is given by this formula; Dav´ıdek V. et al. (2006):
ˆ z H
−1
=
M P
b0
bm z −m
m=0 N P
1+
n=0
= an z −n
M Q
µ=1 N Q
ν=1
1 − z0µ z −1
. (33)
(1 − z∞ν z −1 )
Designing the IIR filters is based on finding the coefficients ˆ z −1 and denominator Aˆ z −1 of this of nominator B ˆ z −1 transfer function, see (33), or finding zeros of B and poles of Aˆ z −1 to satisfy requirements of some tolerance scheme, whereas this design process is related to the quadrature of magnitude characteristic, and the group delay. Traditional approach is based on acceptance of analogue filters’ approximations, because some methodology of the NLP approximation has been developed yet. Well, this paper tries to demonstrate this idea. In MATLAB, there are two ways, how to approach bilinear transformation, thus: • bilinear(num, den, fs) - prewarped mode is not used to indicate “match” frequency, where num is a nominator of transfer function, den is a denominator of transfer function, and fs is a sampling frequency; THE MATHWORKS Inc. (2011): −1 2 1 − z −1 ˆ z −1 = H ˆ s → 2fs 1 − z = . (34) H 1 + z −1 Ts 1 + z −1 • bilinear(num, den, fs, fp) - fp parameter specifies prewarping, which indicates “match” frequency, for which frequency responses (before and after mapping) match exactly; this parameter is identical with edge of analogue filter’s pass-band; THE MATHWORKS Inc. (2011):
[NumHNDPs, DenHNDPs] = cheby1(n, ap, Omegap, ’low’, ’s’); [HNDPs, omegaNDPs] = freqs(NumHNDPs, DenHNDPs);
For the digital NLP filter (z-plane; IIR; prewarped mode is not considered): [NumHNDPz, DenHNDPz] = cheby1(n, ap,... 2*atan(Omegap/(2*fvz))/pi, ’low’, ’z’); [HNDPz, omegaNDPz] = freqz(NumHNDPz, DenHNDPz);
For the digital LP filter (z-plane; IIR; prewarped mode is not considered): [NumHDP1z, DenHDP1z] = cheby1(n, ap,... 2*atan(omega0*Omegap/(2*fvz))/pi, ’low’, ’z’); [HDP1z, omegaDP1z] = freqz(NumHDP1z, DenHDP1z);
Or (it is an equivalent code; using bilinear transformation, prewarped mode is not considered): [NumHDP2z, DenHDP2z] = bilinear(NumHDPs, DenHDPs, fvz); [HDP2z, omegaDP2z] = freqz(NumHDP2z, DenHDP2z);
For the digital LP filter (z-plane; IIR; prewarped mode is considered): [NumHDP3z, DenHDP3z] = cheby1(n, ap,... 2*atan(omega0*Omegap*tan(pi*f0/fvz)/(2*pi*f0))/pi,... ’low’, ’z’); [HDP3z, omegaDP3z] = freqz(NumHDP3z, DenHDP3z);
Or (it is an equivalent code; using bilinear transformation, prewarped mode is not considered): [NumHDP4z, DenHDP4z] = bilinear(NumHDPs, DenHDPs, fvz, f0); [HDP4z, omegaDP4z] = freqz(NumHDP4z, DenHDP4z);
Fig. 8 shows the “normalized” radial frequency-magnitude characteristics of the designed IIR low-pass filter. For its detailed analysis (e. g. phase delay, group delay, TF nominator, TF denominator, TF poles, and structure), we can use a MATLAB tool called Filter Visualization Tool (fvtool command). In this case, “normalized” means dividing the interval h0, πi by π to get the interval h0, 1i. If we consider the prewarped mode and sampling radial frequency ωs = 2π, then we get (after normalization): ωs = ω0 /pω = 2π ⇒ ω0 /π = pω ωs /π = 0.5.
(37)
ter), s-plane (un-normalized analogue filters), and z-plane (IIR filters; normalized digital low-pass filter and unnormalized digital filters). Designing the NLP filter has five important results (approximation order, characteristic equation, characteristic function, transfer function, and group delay), and other results in form of tolerance scheme, magnitude characteristics, phase characteristics, TF poles, roots, and extrema. All these formulas are used to solve practical examples, which include calculations, 2D graphs, 3D graphs, and some fragments of corresponding MATLAB codes. Due to MATLAB&Simulink possibilities, it is also possible to implement this approach, based on Dav´ıdek V. et al. (2006), not only in Simulink, but also in the field of FIR filters or other kinds of analogue filters. ACKNOWLEDGEMENTS Fig. 8. Un-normalized IIR low-pass filter: “Normalized” radial frequency-magnitude characteristics with the prewarped mode (light curve) and without this mode (dark curve); the Cartesian coordinates. In comparison with the non-prewarped mode, value of ω0 equals to 0.42. So, result of (36) is not considered. 4. CONCLUSIONS This paper deals with using some mathematical background - the Chebyshev approximation, Chebyshev polynomials, selected MATLAB functions from Signal Processing Toolbox, and other tool (e. g. Filter Visualization Tool) to design the first-order (also kind I) analogue Chebyshev filters and IIR Chebyshev filters. Because this paper is mainly centred on the mathematic viewpoint, there are lots of fundamental (maybe essential) formulas (de Moivre, Euler, also Jacobi, Gegenbauer, or Rodrigues). Features of the Chebyshev polynomials are also mentioned and discussed in forms of definition, representations, roots, and extrema. The second part of this paper is focused on designing the first-order Chebyshev filters in p-plane (the NLP fil-
This work is supported by project SP2012/111, named “Data Acquisition and Processing from Large Distributed Systems II” of Student Grant Agency (VSB - Technical University of Ostrava). REFERENCES ˇek M. (2006). AnaDav´ıdek V., Laipert M., and Vlc ˇ logov´e a ˇc´ıslicov´e filtry. CVUT. ˇek J. (2009). Modern´ı uˇcebnice elektroniky 6. d´ıl Dolec Kmitoˇctov´e filtry, gener´ atory sign´ al˚ u a pˇrevodn´ıky dat. BEN - technick´a literatura. ˇas J. et al. (1977). Aplikovan´ Nec a matematika I (A aˇz L). SNTL, Prague. Rektorys K. et al. (2000a). Pˇrehled uˇzit´e matematiky I. Prometheus. Rektorys K. et al. (2000b). Pˇrehled uˇzit´e matematiky II. Prometheus. THE MATHWORKS Inc. (2011). Signal processing toolboxTM : User’s guide. Version 6.16 (Release 2011b), offline.