Truncated frequency sampling design using the inverse FFT or IDCT

The design problem is posed by defining an error measure qq as a sum of the squared differences between the actual and the desired frequency response over a set of LL frequency samples. This error function is defined as

where Hd(ωk)Hd(ωk) are the LL samples of the desired response. This problem is easier to formulate and solve if the frequency samples are equally spaced as in ((Reference)) which gives

and the problem is restricted to linear-phase filters where the real-valued amplitude A(ω)A(ω) can be approximated rather than the complex frequency response H(ω)H(ω). For approximations to a complex response, see "Complex L 2 and Minimum Phase Approximation".

Linear phase and equally spaced samples cause (Equation 4) to become

or with a simpler notation

A very powerful property of the Fourier transform allows a straightforward design of least-squared-error FIR filters. Parseval's Theorem, which is based on the orthogonality of the DFT, states that the error defined by (Equation 7) in the frequency domain can also be calculated in the time domain by

where hd(n)hd(n) is the length-L symmetric FIR filter that has the LL frequency response amplitude samples AdkAdk. This may be calculated by the frequency sampling method of (Reference) using the special formulas such as ((Reference)) for length LL or the inverse DFT. The filter to be designed has a length-N symmetric impulse response h(n)h(n) with LL frequency response samples AkAk.

Because the filter h(n)h(n) is of length-N and symmetric, the error equation (Equation 8) can be split into two sums

where h^(n)h^(n) and h^d(n)h^d(n) are the inverse DTFTs of AkAk and AdkAdk respectively, which means they are the h(n)h(n) and hd(n)hd(n) shifted to be symmetic about n=0n=0. This requires the number of frequency samples LL must be odd.

Equation Equation 9 clearly shows that to minimize qq, the NN values of h(n)h(n) are chosen to be equal to the equivalent NN values of hd(n)hd(n) making the first sum equal zero. In other words, h(n)h(n) is obtained by symmetrically truncating hd(n)hd(n). The residual error is then given by the second summation above. An examination of the residual error as a function of NN may aid in the choice of the filter length NN.

For the Type 1 linear-phase FIR filter (described in (Reference)) which has an odd length NN and an even-symmetric impulse response, the LL equally spaced samples of the frequency response from ((Reference)) gives

for k=0,1,2,....,L-1k=0,1,2,....,L-1, where M=(N-1)/2M=(N-1)/2. This formula was derived as a special case of the DFT applied to the Type 1 real, even-symmetric FIR filter coefficients to calculate the sampled amplitude of the frequency response (perhaps better posed using a(n)a(n)). It was noted in (Reference) that it is also a cosine transform and it can be shown that this transformation is orthogonal over the independent values of AkAk, just as the DFT is.

The desired ideal amplitude gives the ideal impulse response hd(n)hd(n) from ((Reference)0 by

for n=0,1,⋯,L-1n=0,1,⋯,L-1. This is used in (Equation 9), and is the ideal impulse response that is truncated and shifted to give a causal, symmetric h(n)h(n).

Use of the alternative equally-spaced sampling in ((Reference)), which has no sample at zero frequency, requires hd(n)hd(n) be calculated from ((Reference)) and ((Reference)). The Type 2 filters with even NN are developed in a similar way and use the design formulas ((Reference)) and ((Reference)). These methods are summarized by:

The filter design procedure for an odd-length Type 1 filter is to first design an odd-length-L FIR filter by the frequency sampling method from ((Reference)) or ((Reference)) or the IDFT, then to symmetrically truncate it to the desired odd-length N and shift it to make h(n)h(n) causal. To design an even-length Type 2 filter , start with an even-length-L frequency-sampling design from ((Reference)) or ((Reference)) or the IDFT and symmetrically truncate and shift. The resulting length-N FIR filters are optimal LS-error approximations to the desired frequency response over the LL frequency samples.

This approach can also be applied to the general arbitrary phase FIR filter design problem.

Weighted, Unevenly Sampled Discrete Least Squared Error Filter Design by Solving Simultaneous Equations

It is sometimes desirable to formulate the least squared error design problem using unequally-spaced frequency samples and/or a weighting function on the error. This is not possible using the IDFT or derived formulas above and requires a different approach to the solution.

Samples of the amplitude response derived for NN odd in ((Reference)) are given by

for k=0,1,⋯,L-1k=0,1,⋯,L-1. This relates the LL frequency samples A(ωk)A(ωk) to the M+1 independent values of the symmetric length-N impulse response h(n). In the design problem where the AkAk are given and the values for h(n) are to be found, this represents LL equations with M+1 unknowns. Because of the symmetries of A(ω)A(ω) shown in (Reference), only half of the LL values of AkAk are independent; however, in some cases, to have proper weights on all LL samples, all must be calculated.

Equation 12 sampled at LL arbitrary frequencies can be written as a matrix equation

where aa is an M+1M+1 length vector with elements which are the first half of h(n)h(n). CC is an LL by (M+1)(M+1) matrix of the cosine terms from (Equation 12), and AA is a length-L vector of the frequency samples A(ωk)A(ωk).

If the formula for the calculation of LL values of the frequency response of a length-N FIR filter in (Equation 10) is used to define an error vector of differences as defined in (Equation 3) and the result is written in the matrix formulation of ((Reference)), the error becomes

or

where ee is a vector of differences between the actual and desired samples of the frequency response. The error measure defined in (Equation 4) becomes the quadratic form

For L>NL>N, equation (Equation 13) is over determined and cannot, in general, be solved for aa. The filter design error measure is the norm of ee, as given in (Equation 16). This error measure is minimized by making ee orthogonal to the columns of CC in (Equation 15). Multiplying both sides of (Equation 14 by the transpose of CC gives

In order for qq to be minimum, ee must be orthogonal to the columns of CC and, therefore, CTeCTe must be zero. Hence, the optimal aa must satisfy the “normal equations" Entry 10, Entry 20, Entry 9 which are

and which can be rewritten in terms of the pseudo-inverse Entry 10, Entry 9 as

If L=NL=N, this becomes the regular frequency-sampling problem and can be solved with zero error. For the case of interest in this section, where L>NL>N, there are still only M+1M+1 equations to be solved. For L>>NL>>N, equation (Equation 10) may be ill-conditioned, and (Equation 19) should not be used to solve them. Special methods will be necessary to avoid serious numerical problems Entry 9.

If a weighted error function is desired, (Equation 4) is modified to give

The normal equations of (Equation 18) become

where WW is a positive-definite matrix of the weights. If zero weights are desired, the effect is be achieved by removing those frequencies from the set of LL frequencies, not by using a zero value weight which would violate the vector-space conditions of a well-posed minimization problem.

Although developed here for the linear-phase filter, (Equation 21) is a very general design approach for the FIR filter that allows arbitrary phase, as well as uneven frequency sampling and a weighting function in the error definition. For the arbitrary phase case, a complex FF is obtained from sampling ((Reference)) and the full h(n)h(n) is used. For the special case of the equally-spaced frequency samples and linear- phase filter with unity weighting, the solution of (Equation 18) or (Equation 21) is the same as given by the frequency sampling design formulas.

One of the important uses of the unequally spaced frequency samples is to create a transition band between the pass and stopbands where there are no samples. This “don't care" band does not contribute to the error measure q and allows better approximation to occur over the pass and stopbands.

Of the many ways to solve (Equation 18) or (Equation 21), one of the easiest and most reliable is the use of Matlab , which has a special command to solve this least-mean-squared error problem. Equation (Equation 19) should not be solved directly. For large LL, it is ill-conditioned and a direct solution will probably have large errors. Matlab uses special algorithms to minimize these numerical errors.

This approach was applied to the same problems that were solved by frequency sampling in the previous section. For N=LN=L, the same results are obtained, thus verifying the theoretical prediction. As LL becomes larger compared to NN, more control is exerted over the behavior between the original sample points. As LL becomes large compared to NN, the solution approaches the same results as obtained where the error is defined as a continuous function of frequency and the integral of the squared error is minimized. Although the solution of the normal equations is a powerful and flexible technique, it can be slow, have numerical problems, and require large amounts of computer memory.

Examples of Discrete Least Squared Error Filter Design

Here we will give examples of several least squared error designs of FIR filters.

Figure 1: Frequency Response of Length-15 FIR Filter Designed by Least Squared Error

As for the frequency sampling design, we see a good lowpass filter frequency response with the actual amplitude interpolating the desired values at different points from the frequency sampling example in Figure 1 even though the length and band edge are the same. Notice there is less over shoot but more ripple near f=0f=0. The Gibbs phenomenon is the same as for the Fourier series.

If a transition band is introduces in the ideal amplitude response between f=0.4f=0.4 and f=0.6f=0.6 with a straight line, the overshoot is reduced significantly but with a slightly slower transition from the pass to stop band. This is illustrated in Figure 2.

Figure 2: Frequency Response of Length-15 FIR Filter with a Transition Band Designed by Least Squared Error

The Unweighted Least Integral Squared Error Approximation

If the error measure is the unweighted integral squared error defined in (Equation 22), Parseval's theorem gives the equivalent time-domain formulation for the error to be

In general, this ideal response is infinite in duration and, therefore, cannot be realized exactly by an actual FIR filter.

As was done in the case of the discrete error measure, we break the infin