In order to obtain better practical filter designs, the
interpolation scheme of the previous section is extended to give an
approximation design method [19]. It should be noted at the
outset that the method developed in this section minimizes an
equation-error measure and not the usual frequency-response error
measure.
The number of frequency samples specified, L+1L+1, will be made
larger than the number of filter coefficients, M+N+1M+N+1. This means
that H2H2 is rectangular and, therefore, Equation 11 cannot in
general be satisfied. To formulate an approximation problem, a
length-(L+1)(L+1) error vector εε is introduced in Equation 8
and Equation 9 to give
b
0
=
H
0
a
+
ε
b
0
=
H
0
a
+
ε
(13)Equation Equation 11 becomes
h
1
-
ε
=
-
H
2
a
*
h
1
-
ε
=
-
H
2
a
*
(14)where now H2H2 is rectangular with (L-M)>N(L-M)>N. Using
the same methods as used to derive Equation 11, the error εε is minimized in a least-squared error sense by the
solution of the normal equations [16]
H
2
T
h
1
=
-
H
2
T
H
2
a
*
H
2
T
h
1
=
-
H
2
T
H
2
a
*
(15)If the equations are not singular, the solution is
a
*
=
-
[
H
2
T
H
2
]
-
1
H
2
T
h
1
.
a
*
=
-
[
H
2
T
H
2
]
-
1
H
2
T
h
1
.
(16)If the normal equations are singular, the pseudo-inverse
[16] can be used to obtain a minimum norm or reduced
order solution.
The numerator coefficients are found by the same techniques as
before in Equation 12
which results in the upper M+1M+1 terms in εε
being zero and the total squared equation error being minimum.
As is true for LS-error design of FIR filters, Equation 15 is often
numerically ill-conditioned and Equation 16 should not be used to solve for
a*a*. Special algorithms such as those used by Matlab and
LINPACK [18], [7] should be employed.
The error εε defined in Equation 13 can better be
understood by considering the frequency-domain formulation. Taking
the DFT of Equation 13 gives
B
k
=
H
k
A
k
+
ε
B
k
=
H
k
A
k
+
ε
(18)where εε is the error in trying to satisfy Equation 8 when the
equations are over-specified. This can be reformulated in terms of
EE, the difference between the frequency response samples of
the designed filter and the desired response samples, by dividing
Equation 8 by AkAk to give
E
k
=
B
k
A
k
-
H
k
=
ε
k
A
k
E
k
=
B
k
A
k
-
H
k
=
ε
k
A
k
(19)EE is the error in the solution of the approximation
problem, and εε is the error in the equations defining
the problem. The usual statement of a frequency-domain approximation
problem is in terms of minimizing some measure of EE, but
that results in solving nonlinear equations. The design procedure
developed in this section minimizes the squared error εε,
thus only requiring the solution of linear equations. There is an
important relation between these problems. Equation 19 shows that
minimizing εε is the same as minimizing EE
weighted by AA. However, AA is unknown until after the problem is
solved.
Although this is posed as a frequency-domain design method, the
method of solution for both the interpolation problem and the LS
equation-error problem is the same as the time-domain Prony's
method, discussed in Complex and Minimum Phase Approximation of reference [19].
Numerous modifications and extensions can be made to this method. If the
desired frequency response is close to what can be achieved by an IIR
filter, this method will give a design approximately the same as that of a
true least-squared solution-error method. It can be shown that
ε=0↔E=0ε=0↔E=0. In some cases, improved
results can be obtained by estimating AkAk and using that as a weight on
εε to approximate minimizing EE. There are iterative
methods based on solving Equation 16 and Equation 17 to obtain values
for AkAk. These values are used as weights on εε to solve for
a new set of AkAk used as a new set of weights to solve again for AkAk[19][26]. We found this approach to converge slowly, but a recent paper
using the log-magnitude [14] was more successful. Other
approaches are given in [24], [23], [11]. The solution of
Equation 16 and Equation 17 is sometimes used to obtain starting
values for other iterative optimization algorithms that need good starting
values for convergence.
An interesting iterative design algorithm that can design to approximate
complex or magnitude frequency responses has be recently proposed by
Jackson [11]. A different approach to the same problem was posed
by Soewito [26], [28].
To illustrate this design method a sixth-order lowpass filter was
designed with 41 frequency samples to approximate. The magnitude of
those less than 0.2 Hz is one and of those greater than 0.2 is zero.
The phase was experimentally adjusted to result in a good magnitude
response. The design was performed with Program 9 in the appendix
of [19] and the frequency response is shown in Figure 7-33 of
[19]. Matlab programs have recently been written which are
smaller and easier to understand than those in FORTRAN.
In this section an LS-error approximation method was posed to
design IIR filters. By using an equation-error rather than a
solution-error criterion, a problem resulted that required only the
solution of simultaneous linear equations.
Like the FIR filter version, the IIR frequency sampling design
method and the LS equation-error extension can be used for complex
approximation and, therefore, can design with both magnitude and
phase specifications.
If the desired frequency-response samples are close to what an IIR
filter of the specified order can achieve, this method will produce
a filter very close to what a true least-squared error method would.
However, when the specifications are not consistent with what can be
achieved and the approximating error is large, the results can be
very poor and in some cases, unstable. It is particularly difficult
to set realistic phase response specifications. With this method, it
is even more important to have a design environment that will allow
easy trial-and-error procedure.
Newly published works which will be discussed here are
[14], [13], [22], [10], [20], [17], [12]. Other references can
be found in [19], [14], [6], [26], [28], [5]. The Matlab command invfreqz()
which is an inverse to the freqz()
command gives a
similar or, perhaps, the same result as the method described in this
note but uses a different formulation [15], [25].
Practical problems occur in the design of a filter to
separate signals according to their energy. Because the energy
content of a signal is the integral or sum of the square of the
signal, a mean-squared-error measure is natural. Unfortunately,
for the IIR filter design problem, the optimization procedure is
nonlinear. This was pointed out in the last section where the
equation error was used in order to have a linear problem.
Because of the nonlinear nature of the least-squared-error
minimization, the method of solution becomes dependent on the
desired frequency response, and therefore, there is no single
method for design. The mean-squared error for magnitude
approximation is defined as
q
(
x
)
=
∑
i
=
0
L
|
H
(
ω
i
)
|
-
|
H
d
(
ω
i
)
|
2
q
(
x
)
=
∑
i
=
0
L
|
H
(
ω
i
)
|
-
|
H
d
(
ω
i
)
|
2
(20)where xx is a vector of filter parameters chosen to minimize q,
and the error is sampled at L+1L+1 frequencies ωiωi. Steiglitz
[19] chose the parameter vector x to be the coefficients of a
cascade structure in order to best fit an iterative optimization
scheme. He applied a standard optimization algorithm, the
Fletcher-Powell method, to the minimization of Equation 19. Other
methods which are more directly related to a squared-error measure
can also be used.
Practical difficulties exist in solving this approximation
problem. In some cases, local minima are found rather than the
global minimum. In other cases, convergence of the minimization
algorithm is slow or does not occur at all. Numerical problems can
result from ill-conditioned equations, and there is no guarantee
that the designed filter will be stable.
An important factor is the choice of a desired frequency-
response function Hd(ω)Hd(ω) that does not result in the optimum
approximation having a large error. This often means not having an
abrupt discontinuity between the passband and stopband.
Another factor is the starting of the iterative optimization
algorithm with a set of coefficients in xx that is close to the
optimum. This can be accomplished by using the frequency sampling
method to give a design that can be
used to start a least-squares algorithm. Because the error defined
in Equation 19 is in terms of magnitudes, an unstable design can be
converted to a stable one by moving the unstable pole at a radius of
rr in the zz-plane to a radius of 1/r1/r. This does not change the
magnitude frequency response and does stabilize the effect of that
pole [19].
A generalization of the idea of a squared-error measure is
defined by raising the error to the pp power where pp is a
positive integer. This error is defined by
q
(
x
)
=
∑
i
=
0
L
|
H
(
ω
)
-
H
d
(
ω
)
|
p
q
(
x
)
=
∑
i
=
0
L
|
H
(
ω
)
-
H
d
(
ω
)
|
p
(21)Deczky [8] developed this approach and used the Fletcher-Powell
method to minimize Equation 21. He also applied this method to the
approximation of a desired group-delay function. An important
characteristic of this formulation is that the solution approaches the
Chebyshev or mini-max solution as p becomes large. Initial work shows the
method of iteratively reweighted least squared error (IRLS) as was applied
to the FIR filter design in (Reference) can also be used for LpLp and
constrained least squared error optimal design of IIR filters
[29].
