Implementation

Our implementation provides signal evaluation at an arbitrary time, where time is specified as an unsigned binary fixed-point number in units of the input sampling period (assumed constant). Figure I.11 shows the time register $t$, and Figure I.12 shows an example configuration of the input signal and lowpass filter at a given time. The time register is divided into three fields: The leftmost field gives the number $n$ of samples into the input signal buffer, the middle field is an initial index $l$ into the filter coefficient table $h(l)$, and the rightmost field is interpreted as a number $\epsilon$ between 0 and $1$ for doing linear interpolation between samples $l$ and $l+1$ (initially) of the filter table. The concatenation of $l$ and $\epsilon$ are called $P\in[0,1)$ which is interpreted as the position of the current time between samples $n$ and $n+1$ of the input signal.

Let the three fields have ${n_n}$, $n_l$, and $n_\eta$ bits, respectively. Then the input signal buffer contains $N=2^{n_n}$ samples, and the filter table contains $L=2^n_l$ ``samples per zero-crossing.'' (The term ``zero-crossing'' is precise only for the case of the ideal lowpass; to cover practical cases we generalize ``zero-crossing'' to mean a multiple of time $t_c=1/f_c$, where $f_c$ is the lowpass cutoff frequency.) For example, to use the ideal lowpass filter, the table would contain $h(l)=$sinc$(l/L)$.

Our implementation stores only the ``right wing'' of a symmetric finite-impulse-response (FIR) filter (designed by the window method based on a Kaiser window [340]). It also stores a table of differences $\hbar(l) = h(l+1) - h(l)$ between successive FIR sample values in order to speed up the linear interpolation. The length of each table is then ${\hat N}=L(N_z+1)$.

Consider a sampling-rate conversion by the factor $\rho = f_s^\prime/f_s$. For each output sample, the basic interpolation Eq. (I.4) is performed. The filter table is traversed twice--first to apply the left wing of the FIR filter, and second to apply the right wing. After each output sample is computed, the time register is incremented by $2^{n_l+n_\eta}/\rho$ (i.e., time is incremented by $1/\rho$ in fixed-point format). Suppose the time register $t$ has just been updated, and an interpolated output $y(t)$ is desired. For $\rho\geq1$, output is computed via

$$\begin{eqnarray*} v & \gets & \sum_{i=0}^{\mbox{$h$\ end}} x(n-i) \left[h(l+iL) ... ...$\ end}} x(n+1+i) \left[h(l+iL) + \epsilon \hbar(l+iL)\right], \end{eqnarray*}$$

where $x(n)$ is the current input sample, and $\epsilon \in[0,1)$ is the interpolation factor. When $\rho<1$, the initial $P$ is replaced by $P^\prime=\rho P$, $1-P$ becomes $\rho-P^\prime = \rho(1-P)$, and the step-size through the filter table is reduced to $\rho L$ instead of $L$; this lowers the filter cutoff to avoid aliasing. Note that $\epsilon$ is fixed throughout the computation of an output sample when $\rho\geq1$ but changes when $\rho<1$.

Figure I.11: Time register format.

Figure I.12: Illustration of waveforms and parameters in the interpolator.

When $\rho<1$, more input samples are required to reach the end of the filter table, thus preserving the filtering quality. The number of multiply-adds per second is approximately $(2N_z+1)\max\{f_s,f_s^\prime\}$. Thus the higher sampling rate determines the work rate. Note that for $\rho<1$ there must be $\lceil{N_zf_s/f_s^\prime}\rceil$ extra input samples available before the initial conversion time and after the final conversion time in the input buffer. As $\rho\to 0$, the required extra input data becomes infinite, and some limit must be chosen, thus setting a minimum supported $\rho$. For $\rho\geq1$, only $N_z$ extra input samples are required on the left and right of the data to be resampled, and the upper bound for $\rho$ is determined only by the fixed-point number format, viz., $\rho_{\mbox{max}}= 2^{n_l+n_\eta}$.

As shown below, if $n_c$ denotes the word-length of the stored impulse-response samples, then one may choose $n_l=1+n_c/2$, and $n_\eta=n_c/2$ to obtain $n_c-1$ effective bits of precision in the interpolated impulse response.

Note that rational conversion factors of the form $\rho=L/m$, where $L=2^n_l$ and $m$ is an arbitrary positive integer, do not use the linear interpolation feature (because $\epsilon \equiv0$). In this case our method reduces to the normal type of bandlimited interpolator [91]. With the availability of interpolated lookup, however, the range of conversion factors is boosted to the order of $2^{n_l+n_\eta}/m$. E.g., for $\rho\approx1$, $n_l=9,n_\eta=8$, this is about $5.1$ decimal digits of accuracy in the conversion factor $\rho$. Without interpolation, the number of significant figures in $\rho$ is only about $2.7$.

The number $N_z$ of zero-crossings stored in the table is an independent design parameter. For a given quality specification in terms of aliasing rejection, a trade-off exists between $N_z$ and sacrificed bandwidth. The lost bandwidth is due to the so-called ``transition band'' of the lowpass filter [340]. In general, for a given stop-band specification (such as ``80 dB attenuation''), lowpass filters need approximately twice as many multiply-adds per sample for each halving of the transition band width.

As a practical design example, we use $N_z=13$ in a system designed for high audio quality at $20$% oversampling. Thus, the effective FIR filter is $27$ zero crossings long. The sampling rate in this case would be $50$ kHz.^I.4In the most straightforward filter design, the lowpass filter pass-band would stop and the transition-band would begin at $20$ kHz, and the stop-band would begin (and end) at $25$ kHz. As a further refinement, which reduces the filter design requirements, the transition band is really designed to extend from $20$ kHz to $30$ kHz, so that the half of it between $25$ and $30$ kHz aliases on top of the half between $20$ and $25$ kHz, thereby approximately halving the filter length required. Since the entire transition band lies above the range of human hearing, aliasing within it is not audible.

Using $512$ samples per zero-crossing in the filter table for the above example (which is what we use at CCRMA, and which is somewhat over designed) implies desiging a length $27\times 512 = 13824$ FIR filter having a cut-off frequency near $\pi/512$. It turns out that optimal Chebyshev design procedures such as the Remez multiple exchange algorithm used in the Parks-McLellan software [340] can only handle filter lengths up to a couple hundred or so. It is therefore necessary to use an FIR filter design method which works well at such very high orders, and the window method employed here is one such method.

It is worth noting that a given percentage increase in the original sampling rate (``oversampling'') gives a larger percentage savings in filter computation time, for a given quality specification, because the added bandwidth is a larger percentage of the filter transition bandwidth than it is of the original sampling rate. For example, given a cut-off frequency of $20$ kHz, (ideal for audio work), the transition band available with a sampling rate of $44$ kHz is about $2$ kHz, while a $48$ kHz sampling rate provides a $4$ kHz transition band. Thus, a $10$% increase in sampling rate halves the work per sample in the digital lowpass filter.