vDSP_desamp FIR filter sample history continuity

Question

I'm trying to use the vDSP_desamp() routine to perform downsampling of an audio signal that is being constantly captured. According to the Accelerate-framework documentation, this routine performs FIR filtering and decimation.

However I have a conflict, since the API does not provide an input parameter for the circular buffer that to my knowledge, is needed to preserve sample history to perform FIR filtering. This leads me to believe that the sample history is being reset during each function call or perhaps a different algorithm is being applied from the one I know.

Unfortunately the function is a complete black box and no hint is given on how the filtering is performed or how signal continuity should be ensured.

Do you think it would suffice if I pad each frame with the last N samples (N being the number of filter taps) of the previous buffer, to acquire filtering continuity? I would discard those padding-sample outputs after calling vDSP_desamp(), but I guess that the sample history line should be correctly initialized after the padding samples are processed.

Any hints on the innards of vDSP_desamp() are welcome.

Answer 1

Why do you say the function is a complete black box? The documentation shows this pseudo-code for vDSP_desamp(A, I, F, C, N, P):

for (n = 0; n < N; ++n)
{
    sum = 0;
    for (p = 0; p < P; ++p)
        sum += A[n*I+p] * F[p];
    C[n] = sum;
}

You can also look in /System/Library/Frameworks/Accelerate.framework/Frameworks/vecLib.framework/Headers/vDSP.h on any recent OS X system for information on almost all vDSP routines. It shows the same thing in a different form:

for (n = 0; n < N; ++n)
    C[n] = sum(A[n*I+p] * F[p], 0 <= p < P);

From this, it is clear that each output value C[n] is a function of P input values from A[n*I + 0] to A[n*I + P-1]. Supposing that you want the first output value in the next call to continue in the same pattern from the current call, simply consider the last output value in the current call (the last n will be N-1, so the last output value will be C[N-1]) and what would be the next value, C[N]. This C[N] would be a function of input values starting at A[N*I + 0], which is A[N*I].

Then match that to the first output value, C[0], in the next call. C[0] will be a function of input values starting at A[0*I + 0], which is A[0].

Therefore, to make the next call continue from where the current call ends, you would want to copy values from A[N*I] and beyond to A[0] and beyond. (This presumes you are reusing the array A by moving data in it and then appending new data. You could also continue with the data in place by passing A+N*I as the first parameter to vDSP_desamp.)

The number of values to copy would be the number of values you have in A after (and including) A[N*I]. E.g., if you have Number values in A, you can use:

memcpy(A+N*I, A, (Number - N*I) * sizeof *A);

Then put new data at A[Number - N*I] and beyond. (If the source and destination ranges overlapped, you would need to use memmove rather than memcpy, but this would be unusual with vDSP_desamp; usually the data shared by successive calls is a small portion of the overall buffer.)

---

Note: The pseudocode expresses the basic mathematics of the operation. The actual arithmetic in vDSP routines may be arranged differently, for performance reasons. So the actual results may have rounding errors different from the values that would be calculated by using the pseudocode directly.