Transcript Lecture 10: LTI design: Analog to Digital transformations

Lecture 10: LTI IIR design: Analog to
Digital and spectral transformations
Instructor:
Dr. Gleb V. Tcheslavski
Contact: [email protected]
Office Hours: Room 2030
Class web site:
http://ee.lamar.edu/gleb/dsp/in
dex.htm
Let’s go! Army will make a human out of you.
ELEN 5346/4304
DSP and Filter Design
Fall 2008
1
2
Preliminary considerations
Techniques for analog filter design are well developed, so it is quite natural to
start designing a digital filter from an analog prototype.
An analog LTI system is stable if all poles of its transfer function H(s) lie in
the left half of the s-plane. Therefore:
1) The j axis of the s-plane should map into the unit circle in the z-plane.
2) The left half of the s-plane should map into the inside of the unit circle in
the z-plane to convert a stable analog filter to a stable digital filter.
A linear phase filter must have a transfer function satisfying:
H ( z)   z  N H ( z 1 )
(10.2.1)
However, in this case, filter would have a reciprocal pole outside the UC and,
therefore, would be unstable. A causal and stable IIR filter cannot have linear
phase.
As a result, when designing IIR filters, only the magnitude response is
specified.
ELEN 5346/4304
DSP and Filter Design
Fall 2008
3
Approximation of derivatives
If the analog filter is specified by the LCC Differential Equation:
d k y(t ) M
d k x(t )
k
  k

k
dt
dt k
k 0
k 0
N
(10.3.1)
the following approximation for the derivative can be used:
y(nT )  y  nT  T  yn  yn1
dy(t )


dt t nT
T
T
(10.3.2)
Where T is the sampling period.
The system function of an analog differentiator:
H (s)  s
(10.3.3)
The system function of a digital differentiator:
1  z 1
H ( z) 
T
ELEN 5346/4304
DSP and Filter Design
Fall 2008
(10.3.4)
4
Approximation of derivatives
Therefore, the frequency-domain equivalent for (10.3.2) is:
1  z 1
s
T
(10.4.1)
Similarly, the kth derivative can be expressed as:
 1 z 
s 

T


1
k
k
(10.4.2)
The system transfer function for the digital IIR filter can be approximated as:
H ( z )  H a ( s) s 1 z1
(10.4.3)
T
Equivalently:
ELEN 5346/4304
DSP and Filter Design
1
z
1  sT
Fall 2008
(10.4.4)
5
Approximation of derivatives
However, it can be shown that this method is
suitable for a quite limited class of filters due
to its mapping property:
There is an attempt to overcome this limitation
by using an alternative mapping:
ynT  kT  ynT kT
dy(t )
1 L
  k
dt t nT T k 1
T
Therefore:
1 L
s    k  z k  z k 
T k 1
By the proper selection of coefficients k, it is possible to map the j axis into
the unit circle. However, this selection is a difficult problem in general.
ELEN 5346/4304
DSP and Filter Design
Fall 2008
(10.5.1)
(10.5.2)
6
Impulse invariance
A digital IIR filter can be obtained by sampling the impulse response of the
analog prototype
gn  h(nT )
(10.6.1)
which may lead to aliasing in the frequency domain  chose small T.
Assuming that the analog filter having N distinct poles is specified by:
N
A
H a ( s)   i
i 1 s   i
coefficients in the partial
fraction expansion
(10.6.2)
poles
N
 it

h
(
t
)

A
e
 i u(t )
l
1
L
(10.6.3)
i 1
Periodical sampling will lead to:
N
 iTn
gn  hl (nT )   Ae
un
i
i 1
ELEN 5346/4304
DSP and Filter Design
Fall 2008
(10.6.4)
7
Impulse invariance
The resulting digital filter will be:

H ( z )   hn z
n 0
n
N

 N
 i nT   n
    Ai e  z   Ai  e iT z 1
n 0  i 1
i 1
n 0




n
(10.7.1)
The inner sum converges since i < 0

e
 iT
n 0
Therefore:
Digital poles at:
z
1

n
1

1  e iT z 1
N
Ai
 iT 1
z
i 1 1  e
H ( z)  
pi  eiT  pi  1BIBO
Not the best method due to aliasing. Also, impulse response can be infinite…
ELEN 5346/4304
DSP and Filter Design
(10.7.2)
Fall 2008
(10.7.3)
(10.7.4)
8
Bilinear transform
The techniques described so far, have severe limitations since they are appropriate
only for LPFs and some BPFs. The bilinear transform does not have such limitations.
The bilinear transform from the s-plane to the z-plane is derived via the trapezoidal
numerical integration of the differential equation describing the analog prototype.
For the given step size T, the BT is given by:
2  1  z 1 
s 

T  1  z 1 
(10.8.1)
This transform is a one-to-one mapping; that is, it maps a single point in the splane to a unique point in the z-plane, and vice versa. The digital transfer function:
G( z)  H a (s) s  2  1 z1 
T  1 z 1 
ELEN 5346/4304
DSP and Filter Design
Fall 2008
(10.8.2)
9
Bilinear transform
We can derive from (10.8.1) that
sT
2
z
sT
1
2
1
In this situation:
T
T
T
  j   1   j
2
2
2  1    jx



T
T
T  1    jx
1



j

1


j


s   j
 0
2
2
2
1
z  1BIBO
Similarly, when  > 0, |z| > 1; therefore, the mapping is correct.
ELEN 5346/4304
DSP and Filter Design
Fall 2008
(10.9.1)
10
Bilinear transform
When
z  e j :
j

2
2  1  z 1  2 1  e j 2 e
s 


1 
 j
T  1 z  T 1 e
T  j 2
e

e j 2  e  j 2 
(10.10.1)
e j 2  e j 2 
2 2 j sin  2 
2
 j tan  2   j
T 2 cos  2 
T
(10.10.2)
Therefore:
Analog freq.
Digital freq.
T

 tan
2
2
We can notice that
(10.10.3)
s s
therefore, a

 p p
“frequency warping” (10.10.3) is needed.
ELEN 5346/4304
DSP and Filter Design
Fall 2008
11
Bilinear transform
The BT destroys the phase response of an
analog filter but preserves linear magnitude.
Therefore, to design a digital filter meeting the
given specs for the magnitude response, we
need to:
1) prewarp the critical bandage frequencies
(p and s) to find their analog equivalents
(p and s) using (10.10.3);
2) design the analog prototype using the
prewarped critical frequencies;
3) apply the BT to obtain the desired digital
transfer function.
Note:
p,a = p,d, s,a = s,d.
analog
digital
(10.11.1)
BT maps the point s =  into the point z = -1.
BT is very good for LPF, HPF, BPF, BSF… – piecewise constant magnitude filters.
ELEN 5346/4304
DSP and Filter Design
Fall 2008
12
Bilinear transform – zero/pole
conversion
Assuming that the analog zeros and poles are known: i.e. the analog transfer
M
function is in the form:
s m 

H  s   ka mN1
(10.12.1)
  s  k 
New zeros
k 1
Applying the BT:
2  T 
H  z   ka
m
m 1
T
N M
N
2   T 
k
k 1
 2   mT 1 
z 
1 

2


T
m 1 
m

N
 2   k T 1 
z 
1 

2


T
k 1 
k

M
M
 1  z

1 N  M
Constant: new gain factor
New poles
N-M zeros at z = -1
(from zeros at )
ELEN 5346/4304
DSP and Filter Design
Fall 2008
(10.12.2)
Bilinear transform – zero/pole
conversion
Therefore:
Digital gain:
M
2  T 
kd  ka
m
m 1
T
N M
N
2   T 
k 1
(10.13.1)
k
Digital poles:
pk 
2  k T
2  k T
(10.13.2)
Digital zeros:
zm 
2   mT
2   mT
(10.13.3)
Need to add N-M zeros at z = -1:
ELEN 5346/4304
13
DSP and Filter Design
Fall 2008
14
“Theoretical” Example
Example 10.1: Design a single-pole digital LPF with a 3-dB bandwidth of 0.2,
by use of the BT applied to the analog filter with the 3-dB bandwidth c :
c
H ( s) 
(10.14.1)
s  c
2
0.65
 c  2

tan
0.1


(10.14.2)
Using the frequency warping: c  tan 



T
T
 2  T
0.65 T
0.65
The analog filter:
(10.14.3)
H ( s) 

s  0.65 T sT  0.65
1
0.245

1

z


0.65
0.325
Applying BT: H ( z ) 


1
1
1
1 z
1  z  0.325  0.325 z
1  0.509 z 1
2
 0.65
1
(10.14.4)
1 z
1  z 1
Which leads to:
ELEN 5346/4304
DSP and Filter Design
H (e j ) 
0.245  1  e j 
1  0.509e
Fall 2008
 j
(10.14.5)
15
“Practical” Example
Example 10.2: Design a digital LPF starting from the elliptic prototype. The
specs for the digital filter are:
p

fp 
 1000
 Hz; f st  s  1500
 Hz; f s  10kHz
2
2
1   p  0.3dB; s  50dB
Note:
T  1 fs  104 s
The critical frequencies normalized for T = 1:
f st
s
 p  2
 0.2 ;s  2
 0.3 ;  1.5
fs
fs
p
fp
1. Frequency pre-warping (for T = 1):


2

tan  p  0.6498rad s ;s  1.01905rad s ; s  1.568
T
2
p
ELEN 5346/4304
DSP and Filter Design
Fall 2008
16
“Practical” Example
2. Table lookup (software):
For the specs (actually, slightly higher than given ones):
1   p  0.28dB; s  50.1dB;N  5;s  p  1.566
we find the poles and zeros for p = 1 rad/s:
1,2  0.09699  j1.0300
5  0.49519
 1,2   j1.6170
 3,4   j 2.4377
Slightly
different
from
Matlab
Imaginary Part
3,4  0.33390  j 0.7177
2
1
0
-1
-2
-3
Stable – all poles are in the left half of the plane.
ELEN 5346/4304
DSP and Filter Design
Fall 2008
-2
-1
0
Real Part
1
2
3
17
“Practical” Example
3. De-normalize the analog design:
Multiply zeros and poles by p:
1,2  1,2  p  0.0630  j 0.6694;3,4  3,4  p  0.2170  j 0.4664;5  5  p  0.3216
1,2  1,2  p   j1.0508; 3,4   3,4  p   j1.5842
4. BT transform of poles and zeros:
pk 
2  k T
2  k T
p1,2  0.7542  j 0.5692
Imaginary Part
p5  0.7228
z1,2  0.5673  j 0.8235
z3,4  0.2289  j 0.9734
0.5
0
-0.5
Since N-M = 1, we must add one zero:
z5  1
DSP and Filter Design
2   mT
2   mT
1
p3,4  0.7278  j 0.3635
ELEN 5346/4304
zm 
-1
-1
Fall 2008
-0.5
0
Real Part
0.5
1
18
“Practical” Example
5. Verification
Need to check whether specs are satisfied.
No: start over (perhaps,
a higher order of the
prototype will help)
Yes: we are done!
In general (especially, for more complicated filters), multiple iterations may be
needed to satisfy specifications (meet or exceed them).
ELEN 5346/4304
DSP and Filter Design
Fall 2008
19
Frequency (spectral) transformation
It is often needed to modify the characteristics of existing digital filter to meet
new specs without starting the design from scratch.
Frequency transformation can convert digital LPF to either another LPF, or HPF,
or BPF, or BSF. The transformation involves replacing the variable z by a rational
function F(ž), while satisfying the following conditions to preserve BIBO:

1
z
1 z 

(10.19.1)
Note: F(z) is an allpass filter:
z 1  l*
z  F ( z )z  F ( z )z  
1
l 1 1   l z
1
1
1
L
Here l  1 to ensure that a stable filter is transferred to another stable filter.
ELEN 5346/4304
DSP and Filter Design
Fall 2008
(10.19.2)
20
Frequency (spectral) transformation
1. LPF to LPF:
1   z z 1  
z  F z 

z   1   z 1
(10.20.1)
e j  

1   e j
(10.20.2)
1
1
e j
Where  is the “old frequency” and  is the “new frequency”.
e j   e j e j  e j    
 j
 j
e  j 2 e  j 2  e  j 2 e j 2  e j 2 e  j 2 
e e
  j 2  j 2 j 2 j 2  j 2  j 2
 j  j
1 e e
e
e
e e  e e 
   
sin 

2 


  
sin 

 2 
We can relate ANY frequency component of the new LPF to the frequency
component of the prototype LPF.
ELEN 5346/4304
DSP and Filter Design
Fall 2008
(10.20.3)
(10.20.4)
21
Frequency (spectral) transformation
   c 
sin  c

2 


   c 
sin  c

 2 
Usually:
(10.21.1)
From (10.20.2):
e j
 j
 j
  e j   e 1   e  2  e j   2e j



2
 j

j

 j 2
1 e
1 e
1 e
Clearly:
ELEN 5346/4304
DSP and Filter Design
(10.21.2)
  2  1 sin  


   tan 1  
 2   2  1 cos  




(10.21.3)
s s

p p
(10.21.4)
Fall 2008
22
Frequency (spectral) transformation
Example of mapping:
c

3
2
"/ 
1
0
=0
 = 0.6
 = -0.9
-1
-2
c
Note: the gain in unchanged
ELEN 5346/4304
DSP and Filter Design
-3
-3
-2
-1
-c
Fall 2008
0
/
1
c
2
3

23
Frequency (spectral) transformation
2. LPF to HPF:
1   z 1
z 1  
z  F  z    1

z 
1   z 1
(10.23.1)
e j  

1   e j
(10.23.2)
1
1
e j
e j   e j  e j  
e  j  e  j

1  e  j e  j
   
cos 

2 


  
cos 

 2 
2



 sin 
1



1

   tan 
 2   2  1 cos  




ELEN 5346/4304
DSP and Filter Design
Fall 2008
(10.23.3)
(10.23.4)
(10.23.5)
24
Frequency (spectral) transformation
Example of mapping:

3
=0
 = 0.6
 = -0.9
2
c
c
"/ 
1
0
-1
-2
-3
-3
-2
-1
-c
ELEN 5346/4304
DSP and Filter Design
Fall 2008
0
/
1
c
2
3

25
Frequency (spectral) transformation
3. LPF to BPF (filter order will be doubled):
 1   1
z 
 1
 1
z 1  
  1 2 2 1
z 
z 1
 1
 1
z 2  2
 c 2  c1 
cos 

2



   c1 
cos  c 2

2


 c 2  c1   c 
 tan  
2

  2 
  cotan 
Here c 2 and c1 are new upper and lower cutoff frequencies.
ELEN 5346/4304
DSP and Filter Design
Fall 2008
(10.25.1)
(10.25.2)
(10.25.3)
26
Frequency (spectral) transformation
A special case: bandwidth preserving
Therefore:
and:
ELEN 5346/4304
c 2  c1  c
(10.26.1)
 1
(10.26.2)
2
1
1

z


z
z
 
1
z 1  


z

1 
 z 1  1
1


z


DSP and Filter Design
Fall 2008
(10.26.3)
27
Frequency (spectral) transformation
4. LPF to BSF (filter order will be doubled):
z 2  2
z 1  

1 
z 1 
1 
1 
1   2 2 1
z 
z 1
1 
1 
 c 2  c1 
cos 

2



   c1 
cos  c 2

2


 c 2  c1   c 
 tan  
2

  2 
  tan 
ELEN 5346/4304
DSP and Filter Design
Fall 2008
(10.27.1)
(10.27.2)
(10.27.3)
28
Spectral transformation: Example
Example 10.3: Design a BW preserving BPF from the LP prototype:
1
20
0.5
15
|G(z)|
Imaginary Part
1  z 1
z 1
G( z) 

1  0.9 z 1 z  0.9
0
-0.5
10
5
-1
-1
-0.5
DC gain:
ELEN 5346/4304
(10.28.1)
DSP and Filter Design
0
0.5
Real Part
0
1
G ( z ) z e j 0 1 
0
0.1
0.2
0.3
Fractional frequency
2
 20  26dB 
1  0.9
Fall 2008
0.4
0.5
(10.28.2)
29
Spectral transformation: Example
Since the filter is BW preserving,  = 1
 j
 z 1   
 
 j
 j  e
z  z 
e  e 
1 
 j 
1


z
1


e




 j0

e
 
 j0
1  e 
For  = 0:
 j0 
1


e


1

1

 1   e  j0 e j0    e  j0 e  j0      e  j0   cos 0
Assuming 0  0.4 - this is where we map DC.   0.309016994
In a general case, both  and  can be found from (10.27.2) and (10.27.2).
 z 1   
z  F (z)  z 
1 
1


z


1
ELEN 5346/4304
1
DSP and Filter Design
1   z 1
z( z   )
F ( z )  1 1

z (z  )
1 z
1
Fall 2008
30
Spectral transformation: Example
Method 1:
1  z 1
GLP ( z ) 
GBP ( z )
1
1  0.9 z
   z 1 
1 z 
1
1
2
1
1 
1


z


z

z
1


z



1


z

 
GBP ( z ) 
1
1
1
1
2




z
1


z
1


z

0.9

z

0.9
z



1  0.9 z 1 
1 
1


z


1
1  z 2

1  1.9 z 1  0.9 z 2
Method 2: SFG replacement…
ELEN 5346/4304
DSP and Filter Design
Fall 2008
31
Spectral transformation: Example
Method 3 (preferred):
z( z   )
 1   z   z 2   z
1 z
 1
zk  F ( zk ) 1  
Zeros:
 z 2  1z1,2
pk  F ( pk )  
Poles:
p( p   )
   p   p 2   p
1 p
 p 2  1.9 p  0.9  0
 p1,2
Therefore:
1.9  1.92  2  4  0.9

 p1,2  0.2936  j 0.9021
2
1  z 2
1  z 2
GBP ( z )  k
 k  1 
1
2
1  1.9 z  0.9 z
1  0.5871z 1  0.9 z 2
Gain factor
The new gain can be determined experimentally…
ELEN 5346/4304
DSP and Filter Design
Fall 2008