22 ), ( )
( )
( ≤ H e ≤M n=
mn ω n jw n ω , (3-2)
whereHn(ejw)is the frequency response of nth band.
For first octave design, the F37 ~ F39 should meet the ANSI S1.11 37th ~ 39th bands specification. Besides, the passband ripple of the three bands should be small enough because the ripple has effects on H22 ~ H36. For decimation filter and interpolation filters, the alias and image components should be suppressed and the passband ripple should be small enough because it also has effects on H22 ~ H36.
3.2 Filter coefficient design
From the multirate filter bank architecture design, we only need to design F37, F38, F39, D, and I filters. The F37, F38, and F39 are designed as band-pass filters such that the responses meet the ANSI class-2 specification. The D and I filters are designed such that the alias and image distortions are suppressed. Furthermore, the filters’ ripple should take into consideration such that F22 ~ F39 are all meet the ANSI S1.11 class-2 specification. Under all these design constrains, we design filters with minimum complexity.
We design filters by using MATLAB FIR design toolbox. The MATLAB toolbox provides Parks-McClellan optimal equiripple FIR design [18] with firpm function. It needs to specify four parameters. The parameters are band edge frequencies (E) (include lower/upper passband frequencies, lower/upper stopband frequencies), weighting (W), band magnitude (A), and filter order (N).We use some notations to represent them.
z (NFn, EFn, AFn, WFn) : design parameters for band-pass filter Fn, n = 37, 38, 39;
z (ND, ED, AD, WD) : design parameters for low-pass filter D;
z (NI, EI, AI, WI) : design parameters for low-pass filter I.
We specify three band ranges for band-pass filters (F37 ~ F39): one passband = [fpL, fpR] and two stopband = [0, fsL] & [fsR, π] that show in Figure 3-3. Thus the band edge frequencies vector E is in the form of (0, fsL, fpL, fpR, fsR, π). For the low-pass filters, D and I, we specify two band ranges. One is passband = [0, fp], the other is stopband = [fs, π] that show in Figure 3-4. The band edge frequencies vector E is in the form of (0, fp, fs, π).
Figure 3-3 Band ranges of band-pass filter
π freq.
magnitude 1
0
Figure 3-4 Band ranges of low-pass filter
For each filter, the passband is set at 1 (0dB) and the stopbands is set at 0.
However, it is difficult to get the input parameters N and W. We use firpmord function in Matlab that can estimate the parameters of N and W by giving passband ripple and stopband attenuation. Therefore, our design parameters are passband ripple, stopband attenuation, and band edge frequencies of all filters.
The filter coefficients we need to design are finding the optimal values of (N, Eo, Ao, W) for F37, F38, F39, D, and I such that the frequency response of nth filter, Hn(ω), should meet the ANSI S1.11 specification that donates as mn(ω)<Hn(ω)<Mn(ω) for n
= 22 ~ 39. Besides, the computation complexity of overall system is minimized. We define the addition and multiplication numbers as computation complexity in filter bank. The addition numbers are shown in equation (3-3).
where N is the order of the filter. An FIR filter with order N needs N additions and the multirate system reduce the addition numbers since some inputs need not to compute.
Besides, the equation takes the polyphase implementation into consideration already.
Similarly, the multiplication numbers are shown in equation (3-4).
⎡ ⎤ ⎡ ⎤ ⎡ ⎤
(
/2 /2 /2)
1.9375( ⎡
/2⎤ ⎡
/2⎤ )
96875 .
1 × NF37 + NF38 + NF39 + × ND + NI (3-4)
The number of multiplications is almost half to the number of additions because the MATLAB tool designs these filters with linear phase and thus the coefficients are symmetric. From the above two equations, we realize that the complexity is dominated on the filters’ order and every filter has almost the same weighting. So we
should minimize the order of each filter to get the minimum computation complexity.
It is hard to design five filters in the same time since the filter ripple of D and I effect the response of F22 ~ F36. We develop a systematic filter coefficient design flow for our multirate system and try to minimize the complexity of the overall system.
Our design flow can separate into two steps as shown in Figure 3-5. The first step designs F37 ~ F39 independently and then designs D and I filters by assuming F37 ~ F39 are ideal band-pass filters. The second step explores the filters’ ripple to get the optimal design results. Following describes the details of the two steps.
Figure 3-5 Overview of coefficient design flow
z Step 1: Design Fn and of D and I independently.
In this step, we design F37, F38, and F39 filters first by assuming D and I filters are ideal low-pass filters. Then we design D and I filters by assuming F37
~ F39 are ideal band-pass filters.
At first, we design 1/3-octave filters by assuming D and I filters are ideal.
We only need to consider the response of F37 ~ F39 since the D and I filters are ideal. For each band-pass filter Fn, we need to find four design parameters that have mentioned before. These are band edge frequencies (0, fsL, fpL, fpR, fsR, π), weighting on each band, band amplitude, and filter order. The weighting is derived from firpmord function by setting passband ripple and stopband
attenuation. It is time consuming to find the band edge frequencies for a band-pass filter. To reduce the explore time, we translate the problem of finding the band edge frequencies (fsL, fpL, fpR, fsR) into finding the left side transition bandwidth (TBWL), right side transition bandwidth (TBWR), fsL, and fsR, where the fpL = fsL + TBWL and fpR = fsR - TBWR. We develop a algorithm to estimate the TBWL and TBWR and then design the Fn according to the TBWL and TBWR values. Figure 3-6 is the algorithm of designing the F37 ~ F39.
In the algorithm, we first estimate the TBWL and TBWR and then the program goes to a while loop designing filter Fn. In the while loop, we first specify the band edge frequencies by exploring the range of fsL & fsR and then use the firpmord function to estimate weighting and corresponding order with given ripple and attenuation.
After getting the input parameters of firpm function, we design the filter and check the filter’s response. If the filter’s ripple and attenuation do not meet our given ANSI S1.11 specification, we increase the order by one and then re-design the filter. Then we check the designed filter to see if the filter frequency response meets the ANSI S1.11 specification. If the response meets the ANSI S1.11 specification, the program terminates and output the filter coefficients. If the filter response does not meet the specification, we reduce the estimated TBWL and TBWR by 0.95 times and redo the whole program until we find a design that meets ANSI S1.11 specification. From the program, there are two items we do not mention. One is how we estimate the TBWL and TBWR. The other is finding the exploration range of fsL & fsR. For each estimated TBWL and TBWR, we need to explore the points of fsL & fsR such that we can get the all input parameters of the band edge frequencies. Figure 3-7 shows the explore
range of the band edge frequencies. The passband and stopband constrain are shown in Figure 3-7, too.
Estimate TBWL //transition bandwidth at left side Estimate TBWR //transition bandwidth at right side While (1)
{
Find the exploration range of fsL & fsR;
For each (fsL, fsR) pair explored {
Let EF37 =(0, fsL, fsL+TBWL, fsR-TBWR, fsR, pi)
Estimate WF37 and NF37 by ripple 1dB, attenuation -60dB;
Design F37 with (NF37, EF37, WF37)
While (ripple>1dB or attenuation>-60dB) {
NF37 += 1;
re-design F37; }
If(H37 meet the ANSI specification) {
Found = 1;
Break;
} }
TBWL = TBWL*0.95;
TBWR = TBWR*0.95;
}
Figure 3-6 Pseudo code for design F37 (F38, F39)
Since we only need to get fsL & fsR points, the explore process we designed is dividing the search plane for every iteration. In the ith iteration, we divide the fsL’s & fsR’s search plane into i2 regions and check the central points. The equation of the search points is shown in (3-5).
Freq.
Figure 3-7 Explore the band edge frequencies
1
The other item we do not mention is the estimation values of the TBWL and TBWR. Estimating TBWL and TBWR from a band-pass filter is a time consuming work and we simplify the problem by decomposing a band-pass filter into a low-pass filter and a high-pass filter. We use a high-pass filter to find TBWL and a low-pass filter to find TBWR.
Figure 3-8 is the algorithm of finding the TBW of high-pass (low-pass) filter.
We take finding the TBWL for example. The variable fp and fs in this algorithm are the passband frequency and stopband frequency in a high-pass filter. The step is the stepping width. In the beginning, we generate the high-pass (low-pass) specification by cutting the band-pass filter’s specification for checking that shown in Figure 3-9.
Generate the high-pass specification mHP(w) and MHP(w);
TBWL = fm – f1’, step = TBWL/2;
While (step> π/8192) {
For each fs in [f1’, f1] with π/8192 stepping {
Let F = (0, fs, fp = fs + TBWL, pi)
Estimate W & N with ripple=1dB & atten.=-60dB;
Design a high-pass filter H with (N, F, W);
While (H’s ripple>1dB or attenuation>-60dB) {
N +=1; Redesign H with (N, F, A, W);
}
If( mLP < H < MLP(w) ) {
feasible = 1; break;
} }
If (feasible = 1) TBWL += step;
Else TBWL -= step;
step = step / 2;
}
Figure 3-8 Pseudo code of finding transition band width in high-pass filter
MHP( )
Figure 3-9 High-pass and Low-pass filters’ specification
After generating the specification of low-pass filter and high-pass filters, we give the transition band width (TBW) initial value fm – f1. The value fm – f1 is the maximum value that the high-pass filter may meet the specification. Furthermore, we give the stepping width an initial value TBW/2. In the program of finding the TBW, we gradually increase the stopband frequency until we find a solution. The increasing step is the minimum resolution that we used for analysis. We estimate the weighting and the filter order by using firpmord function like design Fn and then design the high-pass filter. We check the response of the designed filter to see if the designed filter’s passband ripple smaller than 1dB and stopband attenuation smaller than -60 dB. The filter order should increase by one and redesign if the passband ripple bigger than 1 dB or the stopband attenuation bigger than -60 dB. If the impulse response of the filter meets the high-pass specification, the next testing TBW will increase by a value of step, else the next testing TBW will decrease by a value of step. The step value is halved every executing iteration. The program runs until the step value samller than FFT resolution (π/8192).
Then, we design decimation filter (D) and interpolation filter (I) with ideal band-pass filter Fn. The ideal Fn means it meet the ANSI S1.11 specification by setting the ripple equal to 0dB and stopband attenuation equal to 60dB. Further, the transitions are infinite sharp. Design D and I for our multirate system should consider the alias and image terms that occur from the down-sample and up-sample operations. D and I filters are low-pass filters so that the band edge frequencies we need to derive are a passband frequency and a stopband frequency. Input parameters for designing D and I are the passband frequencies (fp_D and fp_I), stopband frequencies (fs_D and fs_I), passband ripple, and
stopband attenuation.
From previous discussion, we found the expression of the 36th band’s output is composed of two terms that shown in equation (3-6).
) ( ) ( ) ( ) ( ) ( ) ( ) ( ) ( )
( 39 39 2
36 z X z F z D z I z X z F z D z I z
Y = + − − (3-6)
The term X(-z)F39(z2)D(-z)I(z) is the distortion term because it contains the aliasing component X(-z), so the constraints of D(-z)I(z) is to suppress the distortion term. The term X(z)F39(z2)D(z)I(z) is the desired term, we should make the product of these four components meet the specification.
Figure 3-10 is the specification of F36, and the ideal response of the F39(z2).
From the figure, we can see that the response of D(z)I(z) need to suppress the image term of the F39(z2). Since we want to eliminate the distortion term of the Y36(z), the D(-z)I(z) should suppress the X(-z). The method to suppress the X(-z) is reducing the fs_D and fs_I with small stopband attenuation. Furthermore, D(z)I(z) should not cut the 36th band so that the fp_D and fp_I should bigger than f2 as shown in Figure 3-10. Besides, passband ripple should be small enough for meeting the 1 dB specification. The image of the F39(z2) is between 0.63π to 0.71π that shown in the Figure 3-10, too. To suppress the image of F39(z2), the fs_D and fs_I should small enough (<0.71π).
Figure 3-10 Specification of F36 and the ideal response of F39(z2)
Figure 3-11 Pseudo code of designing D and I
From constrains of designing D and I, the stopband frequencies and ripples of the D and I are the most important items. We should find sufficiently small fs_D and fs_I in our design and then search a small ripple value to meet ANSI S1.11 specification. Figure 3-11 shows the algorithm of designing the D and I.
At first, we set the passband frequencies with the value f2 and try to find the stopband frequencies. We gradually decrease the stopband frequency of the
Set fp_D = fp_I = f2;
For (fs_D = pi, fs_D > f2, fs_D -=step) {
Design D with passband ripple = 1dB, attenuation= -60dB;
For (fs_I = fs_D, fs_I > f2, fs_I -=step) {
Design I with passband ripple = 1dB, attenuation= -60dB;
If (D(-z)I(z) < -60dB)
Break and record the fs_D and fs_I;
} }
rp = 1dB;
While(1)
{ Design D with passband ripple = rp;
Design I with passband ripple = rp;
If (H22~H36 meet the ANSI Spec.) Break and output the coefficients;
rp = rp * 0.95;
}
decimation filter and the interpolation filter from π to f2 by FFT resolution stepping (π/8192). For all stopband frequencies, we design with the ripple value equals to 1 dB and the attenuation value equals to -60 dB and to see if the D(-z)I(z) < -60 dB. If D(-z)I(z) < -60 dB, D and I filters meet our constrain of suppressing the alias term and we record the stopband frequencies. Second, we adjust the ripple by 0.95 times stepping from 1 dB and use the stopband frequencies we recorded to design D and I filters. The second step is executing until there is a design that the responses of 22nd to 39th bands are all meeting the ANSI S1.11 specification while applying ideal F37 ~ F39.
z Step 2: Exploration for optimal ripple distribution
In step 1, we design 1/3-octave filters Fn without considering the ripple of D, I. Besides, we design D and I filters without considering the ripple of Fn we designed. Since F37 ~ F39, D, and I we designed are not ideal, we should take the ripple into consideration. In this step, we explore the ripple of F37 ~ F39 from 0 to 1 dB with 0.1 dB stepping and find the corresponding D and I filters to meet the ANSI S1.11 specification. After designing, we choose the design that has minimum computation complexity as our final result. Figure 3-12 shows the algorithm of the step 2.
Figure 3-12 Pseudo code of step2
For a specify ripple Fn, we get output ripple values of D and I to meet ANSI For (rp_F = 0.1dB, rp_F < 1dB, rp_F += 0.1dB)
{
Design F37, F38, F39 with ripple = rp_F;
Design D and I based on the designed F37, F38, F39;
Calculate and record the complexity associated with rp_F;
}
Find the optimal result that has the minimal complexity
S1.11 specification with minimum complexity. The filter taps and the corresponding complexity results are shown in Table 3-1. For hardware implementation, we adjust the F38, F39, D, and I into even taps or odd taps according to the F37. The adjusting method is just increasing the order by one. For example, if the order of F37 is even, we adjust odd-tap filters of these four filters into even taps. Else we adjust even-tap of these filters into odd taps. The reason we adjust the taps of our design will be mentioned in Chapter 4. We finally choose the 0.9 dB ripple of Fn and 0.1 dB ripple of D and I as our design result since the computation complexity is optimal for our design flow.
Table 3-1 Exploration for optimal ripple distribution
D I F37 F38 F39 MPY ADD
0.0 0.40 29 35 NA NA NA NA NA
0.1 0.20 33 37 57 47 37 176.6 351.3
0.2 0.30 33 37 57 45 37 174.7 347.3
0.3 0.10 35 41 53 57 33 181.5 361.0
0.4 0.30 33 37 53 41 33 162.8 323.7
0.5 0.10 35 41 53 41 33 165.8 329.5
0.6 0.10 35 41 49 39 33 159.8 317.7
0.7 0.10 35 41 47 37 31 153.9 305.9
0.8 0.10 35 41 43 33 27 142.1 282.3
0.9 0.10 35 41 41 33 27 140.2 278.3
1.0 0.02 41 49 41 33 25 145.0 288.0
ripple of F (dB)
ripple of D and I
# tap # of
computations/sample