0
$\begingroup$

I am trying to implement the resampling algorithm as described in this paper: https://ccrma.stanford.edu/~jos/resample/resample.pdf

I managed to construct the $h$ (windowed sinc) table, and I have implemented the loop. I am having a hard time to compute the $\eta$ interpolation factor and the $l$ index index for the $h$ table. The "fixed point" implementation of the paper is a bit confusing to me.

Looking at figure 7 in the paper, I think the interpolation factor needs to be computed between the $h$ entries instead of between the original samples. However I don't understand how to compute the starting point ($h(0)$). As in the paper, the table only contains the half wing.

I have just started working with audio code, so I might be misunderstanding something fundamental about the paper.

Here's the main code:

Fs = 44100
Fs_prime = 96000

Nz = 13
L = 128
M = (Nz * 2 + 1) * L

def resample(original, windowed_sinc):
    ratio = Fs_prime / Fs
    t_increment = 1 / Fs_prime
    inteval_orig = 1 / Fs
    Lo = len(original)
    l_half = M // 2

    y = np.zeros(int(ratio * Lo + 0.5))

    t = 0
    l = 0
    for n in range(y.shape[0]):
        n0 = int(t / inteval_orig)
        original_t = n0 * inteval_orig
        eta0 = 1 - (t - original_t ) / inteval_orig
        eta1 = 1 - eta0

        v = 0
        eta = eta0

        for i in range(Nz):
            ni = n0 - i
            if ni < 0 or ni >= Lo:
                continue

            lL = l + (i * L)
            print(lL)

            if lL + 1 < 0 or lL + 1 >= l_half:
                break

            v += original[ni] * ( w[lL] + eta * (w[lL + 1] - w[lL]) )

        eta = eta1
        for i in range(Nz):
            ni = n0 + 1 + i
            if ni < 0 or ni >= Lo:
                continue

            lL = l + (1 + i * L)
            print(lL)

            if lL + 1 >= l_half:
                break

            v += original[ni] * ( w[lL] + eta * (w[lL + 1] - w[lL]) )

        l = (l + n0) % l_half

        y[n] = v

        t += t_increment

    return y
$\endgroup$

1 Answer 1

2
$\begingroup$

You have to first decide how many equally-spaced fractional-sample delays you want. You may have to linearly interpolate between adjacent fractional-sample delays. Call that number $R$. (For extremely high-quality audio resampling, I have found $R=512$ is good. But fewer discrete delays may be acceptable. Like $R=128$.)

Then you have to decide how many non-zero taps you use to compute your fractional-delayed sample. It should be an even number with an equal number of samples preceding the interpolated sample as the number of samples following the interpolated sample. Call this number $N$. I have found that $N=16$ to be good enough for waveform synthesis, but for high-quality audio, I have gone to $N=32$ (so the 32 original samples surrounding the interpolated sample on both sides are used).

Now the number of non-zero coefficients of your ideal brickwall LPF FIR is $N \times R - 1$ (the FIR order is $NR-2$) and the and the target cutoff frequency for the brickwall LPF is

$$ f_\mathrm{cutoff}=\frac{f_\mathrm{Nyquist}}{R} $$

So, first design the impulse response for your FIR filter (sinc or Parks-McClellan or Least-squares). Append a zero preceding the $NR-1$ samples in the FIR. Now you have exactly $NR$ samples.

Now, you have to reshape (that's the MATLAB term) this $1 \times NR$ column matrix into an $R \times N$ matrix. $R$ columns each having $N$ samples. Those are the impulse responses for each of the $R$ fractional-sample delays.

$\endgroup$
2
  • $\begingroup$ Thanks for your reply! And thanks also for dsp.stackexchange.com/a/37715/71569, I used it to implement the sinc window. I now understand the indexing of $h$ in the original paper, if you think of it in terms of a matrix. I guess I still don't get how to determine how to pick the right column. If I have $t_{new}$ that sits between $old_{t0}$ and $old_{t1}$, how do I compute $l$ to index into the matrix? $\endgroup$ Commented Mar 18 at 10:16
  • $\begingroup$ Consider the "equally-spaced fractional-sample delays" as hypothetical samples. If $R=128$ it's like you upsampled by a factor of 128, but you haven't bothered to compute all 128 new samples for each input sample at the original sample rate. Now the fraction that $t_{new}$ sits between $old_{t0}$ and $old_{t1}$ is the same fraction of the 128 samples. You only need to compute two adjacent samples. Then linearly interpolate between those two. , $\endgroup$ Commented Mar 19 at 3:51

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.