Consider the filter in (Reference) with transfer function
H
(
z
)
=
B
(
z
)
A
(
z
)
H
(
z
)
=
B
(
z
)
A
(
z
)
(1)We can rewrite it as
B
(
z
)
=
H
(
z
)
A
(
z
)
B
(
z
)
=
H
(
z
)
A
(
z
)
(2)which represents the convolution b(n)=h(n)*a(n)b(n)=h(n)*a(n). This operation can be represented in matrix form [4] as follows,
System (Reference) can be viewed as [1],
If L=M+NL=M+N, equations (Reference) and (Reference) describe the same square systems as equations (Reference) and (Reference), letting ckck in the former systems to become hkhk in the latter ones. Therefore when the number of impulse response samples to be matched is equal to M+NM+N, solving (Reference) for the coefficients akak and bkbk is equivalent to applying Pade's method to matching the first N+MN+M values of the impulse response h(n)h(n). In consequence, this method is known as Pade's method for IIR filter design [4] or simply Pade approximation method [5], [2]. Pade's method is an interpolation algorithm; since the number of samples to be interpolated must be equal to the number of filter parameters, Pade's method is limited to very large filters, which makes it impractical.
The relationship of the above method with Prony's method can be seen by posing Prony's method in the ZZ-domain. Consider the function f(n)f(n) from (Reference),
f
(
n
)
=
∑
k
=
1
N
c
k
e
s
k
n
=
∑
k
=
1
N
c
k
λ
k
n
f
(
n
)
=
∑
k
=
1
N
c
k
e
s
k
n
=
∑
k
=
1
N
c
k
λ
k
n
(6)The ZZ-transform of f(n)f(n) is given by
F
(
z
)
=
Z
{
f
(
n
)
}
=
∑
n
=
-
∞
∞
f
(
n
)
z
-
n
F
(
z
)
=
Z
{
f
(
n
)
}
=
∑
n
=
-
∞
∞
f
(
n
)
z
-
n
(7)where Z{·}Z{·} denotes the ZZ-transform operator. By linearity of the ZZ-transform [5],
F
(
z
)
=
Z
{
f
(
n
)
}
=
∑
k
=
1
N
c
k
Z
{
λ
k
n
}
F
(
z
)
=
Z
{
f
(
n
)
}
=
∑
k
=
1
N
c
k
Z
{
λ
k
n
}
(8)Using
Z
{
λ
k
n
}
=
Z
{
λ
k
n
u
(
n
)
}
=
1
1
-
λ
k
z
-
1
Z
{
λ
k
n
}
=
Z
{
λ
k
n
u
(
n
)
}
=
1
1
-
λ
k
z
-
1
(9)(the left equality is true since f(n)=0f(n)=0 for n<0n<0 by assumption) in Equation 8 we get
F
(
z
)
=
∑
k
=
1
N
c
k
1
-
λ
k
z
-
1
=
c
1
1
-
λ
1
z
-
1
+
c
2
1
-
λ
2
z
-
1
+
⋯
c
N
1
-
λ
N
z
-
1
=
c
1
(
1
-
λ
2
z
-
1
)
⋯
(
1
-
λ
N
z
-
1
)
+
⋯
c
N
(
1
-
λ
1
z
-
1
)
⋯
(
1
-
λ
N
-
1
z
-
1
)
(
1
-
λ
1
z
-
1
)
(
1
-
λ
2
z
-
1
)
⋯
(
1
-
λ
N
z
-
1
)
=
∑
i
=
1
N
c
i
∏
j
=
1
j
≠
i
N
(
1
-
λ
j
z
-
1
)
∏
k
=
1
N
(
1
-
λ
k
z
-
1
)
F
(
z
)
=
∑
k
=
1
N
c
k
1
-
λ
k
z
-
1
=
c
1
1
-
λ
1
z
-
1
+
c
2
1
-
λ
2
z
-
1
+
⋯
c
N
1
-
λ
N
z
-
1
=
c
1
(
1
-
λ
2
z
-
1
)
⋯
(
1
-
λ
N
z
-
1
)
+
⋯
c
N
(
1
-
λ
1
z
-
1
)
⋯
(
1
-
λ
N
-
1
z
-
1
)
(
1
-
λ
1
z
-
1
)
(
1
-
λ
2
z
-
1
)
⋯
(
1
-
λ
N
z
-
1
)
=
∑
i
=
1
N
c
i
∏
j
=
1
j
≠
i
N
(
1
-
λ
j
z
-
1
)
∏
k
=
1
N
(
1
-
λ
k
z
-
1
)
(10)The numerator and denominator in (Reference) are two polynomials in z-1z-1 of degrees N-1N-1 and NN respectively. Expanding both polynomials, equation (Reference) becomes
F
(
z
)
=
b
0
+
b
1
z
-
1
+
⋯
+
b
N
-
1
z
-
(
N
-
1
)
1
+
a
1
z
-
1
+
⋯
+
a
N
z
-
N
F
(
z
)
=
b
0
+
b
1
z
-
1
+
⋯
+
b
N
-
1
z
-
(
N
-
1
)
1
+
a
1
z
-
1
+
⋯
+
a
N
z
-
N
(11)Assuming a0=1a0=1 does not affect the formulation of F(z)F(z). This is equivalent to dividing both the numerator and denominator of (Reference) by a0a0. From Equation 11 it is clear that Prony's method is in fact a particular case of the Pade approximation method described earlier in this section (with M=NM=N). From Equation 10 and Equation 11 we have
c
1
1
-
λ
1
z
-
1
+
c
2
1
-
λ
2
z
-
1
+
⋯
c
N
1
-
λ
N
z
-
1
=
b
0
+
b
1
z
-
1
+
⋯
+
b
N
-
1
z
-
(
N
-
1
)
1
+
a
1
z
-
1
+
⋯
+
a
N
z
-
N
)
c
1
1
-
λ
1
z
-
1
+
c
2
1
-
λ
2
z
-
1
+
⋯
c
N
1
-
λ
N
z
-
1
=
b
0
+
b
1
z
-
1
+
⋯
+
b
N
-
1
z
-
(
N
-
1
)
1
+
a
1
z
-
1
+
⋯
+
a
N
z
-
N
)
(12)Therefore, given 2N2N samples of h(n)=f(n)h(n)=f(n) one can solve for the parameters ckck and λkλk in Equation 6 merely by applying partial fraction expansion [3] on Equation 12.
Consider the systems in (Reference) and (Reference). If L>M+NL>M+N then (Reference) is an overdetermined system and cannot be solved exactly in general. However, it is possible to find a least squares approximation following an approach similar to the one used in the frequency domain design method from (Reference). Given L>M+NL>M+N samples of an impulse response h(n)h(n), we rewrite H^1H^1 as
H
^
1
=
[
h
^
H
2
]
H
^
1
=
[
h
^
H
2
]
(13)where
h
^
=
h
M
+
1
⋮
h
L
and
H
2
=
h
M
h
M
-
1
⋯
h
M
+
1
h
M
⋮
⋱
⋮
h
L
-
1
⋯
h
L
-
N
h
^
=
h
M
+
1
⋮
h
L
and
H
2
=
h
M
h
M
-
1
⋯
h
M
+
1
h
M
⋮
⋱
⋮
h
L
-
1
⋯
h
L
-
N
(14)Following the same formulation of (Reference), the filter coefficients akak and bkbk are found by solving for a^a^ and bb in
a
^
=
-
[
H
2
T
H
2
]
-
1
H
2
T
h
^
a
^
=
-
[
H
2
T
H
2
]
-
1
H
2
T
h
^
(15)where H1H1 is an M×NM×N matrix given by
H
1
=
h
0
0
0
⋯
0
h
1
h
0
0
⋯
0
h
2
h
1
h
0
⋮
⋮
⋱
⋮
h
M
h
M
-
1
⋯
H
1
=
h
0
0
0
⋯
0
h
1
h
0
0
⋯
0
h
2
h
1
h
0
⋮
⋮
⋱
⋮
h
M
h
M
-
1
⋯
(17)The error analysis for these algorithms is equivalent to the one performed in (Reference). In fact, both approaches (frequency sampling versus Pade's method) are quite similar. Among the common properties, both methods use an equation error criterion rather than the more useful solution error. However, it is worth to point out that in the time domain (Pade) approach, samples of the impulse response are used to make the approximation, and the method uses a linear convolution rather than the cyclic one from the frequency domain approach. Also, since the latter method uses a uniform sampling grid within the complete frequency spectrum between ω=0ω=0 and ω=2πω=2π instead of using the first few samples of an infinitely long sequence (h(n)h(n)), the approximation properties of the frequency domain method are superior.
