手动重新实现双二阶 IIR 滤波器的 vDSP_deq22

Question

我正在将当前使用 Apple 特定（加速）vDSP 函数 vDSP_deq22 的滤波器组移植到 Android（加速不可用）。滤波器组是一组带通滤波器，每个滤波器返回各自频段的 RMS 幅度。目前代码（ObjectiveC++，改编自 NVDSP）如下所示：

- (float) filterContiguousData: (float *)data numFrames:(UInt32)numFrames channel:(UInt32)channel {

    // Init float to store RMS volume
    float rmsVolume = 0.0f;

    // Provide buffer for processing
    float tInputBuffer[numFrames + 2];
    float tOutputBuffer[numFrames + 2];

    // Copy the two frames we stored into the start of the inputBuffer, filling the rest with the current buffer data 
    memcpy(tInputBuffer, gInputKeepBuffer[channel], 2 * sizeof(float));
    memcpy(tOutputBuffer, gOutputKeepBuffer[channel], 2 * sizeof(float));
    memcpy(&(tInputBuffer[2]), data, numFrames * sizeof(float));

    // Do the processing
    vDSP_deq22(tInputBuffer, 1, coefficients, tOutputBuffer, 1, numFrames);
    vDSP_rmsqv(tOutputBuffer, 1, &rmsVolume, numFrames);

    // Copy the last two data points of each array to be put at the start of the next buffer.
    memcpy(gInputKeepBuffer[channel], &(tInputBuffer[numFrames]), 2 * sizeof(float));
    memcpy(gOutputKeepBuffer[channel], &(tOutputBuffer[numFrames]), 2 * sizeof(float));

    return rmsVolume;
}

正如here 所见，deq22 通过递归函数在给定的输入向量上实现双二阶滤波器。这是文档中对函数的描述：

• A =：单精度实数输入向量
  
• IA =：迈向 A.
  
• B =：5 个单精度输入（滤波器系数），步长为 1。
  
• C =：单精度实数输出向量。
  
• IC =：C 的步伐。
  
• N =：要生成的新输出元素的数量。

这是我目前所拥有的（它在 Swift 中，就像我已经在 Android 上运行的其他代码库一样）：

// N is fixed on init to be the same size as buffer.count, below
// 'input' and 'output' are initialised with (N+2) length and filled with 0s

func getFilteredRMSMagnitudeFromBuffer(var buffer: [Float]) -> Float {
    let inputStride = 1 // hardcoded for now
    let outputStride = 1

    input[0] = input[N]
    input[1] = input[N+1]
    output[0] = output[N]
    output[1] = output[N+1]

    // copy the current buffer into input
    input[2 ... N+1] = buffer[0 ..< N]

    // Not sure if this is neccessary, just here to duplicate NVDSP behaviour:
    output[2 ... N+1] = [Float](count: N, repeatedValue: 0)[0 ..< N]

    // Again duplicating NVDSP behaviour, can probably just start at 0:
    var sumOfSquares = (input[0] * input[0]) + (input[1] * input[1])

    for n in (2 ... N+1) {
        let sumG = (0...2).reduce(Float(0)) { total, p in
            return total + input[(n - p) * inputStride] * coefficients[p]
        }

        let sumH = (3...4).reduce(Float(0)) { total, p in
            return total + output[(n - p + 2) * outputStride] * coefficients[p]
        }

        let filteredFrame = sumG - sumH
        output[n] = filteredFrame
        sumOfSquares = filteredFrame * filteredFrame
    }

    let meanSquare = sumOfSquares / Float(N + 2) // we added 2 values by hand, before the loop
    let rootMeanSquare = sqrt(meanSquare)
    return rootMeanSquare
}

该滤波器为 deq22 提供了不同幅度的输出，并且其中似乎有一个循环的颤抖的圆形“噪声”（输入音调恒定，该频率的幅度上下波动）。

我已检查以确保每个实现之间的系数数组相同。每个滤波器实际上似乎都在“工作”，因为它拾取了正确的频率（并且只有那个频率），它只是这种泵浦，并且 RMS 幅度输出比 vDSP 安静得多，通常是数量级：

   Naive    |    vDSP
3.24305e-06   0.000108608
1.57104e-06   5.53645e-05
1.96445e-06   4.33506e-05
2.05422e-06   2.09781e-05
1.44778e-06   1.8729e-05
4.28997e-07   2.72648e-05

有人能看出我的逻辑有问题吗？

编辑：这是一个带有恒定 440Hz 音调的结果的 gif 视频。各种绿色条是单独的过滤带。第 3 个频段（此处显示）是调谐到 440Hz 的频段。

正如预期的那样，NVDSP 版本仅显示与输入音量成比例的恒定（非波动）幅度读数。

Answer 1

好的，sumOfSquares = filteredFrame * filteredFrame 应该是 +=，而不是赋值。所以只计算最后一帧，解释了很多；）

如果您想在 Swift 中进行一些双二阶过滤，请随意使用它。像之前的 NVDSP 一样的 MIT 许可证。