Bug 1265408 - Import IIRFilter from blink; r=padenot
Imported from git revision 57f70919a0a3da5ba002b896778b580986343e08.
MozReview-Commit-ID: 8QF0wWEHI8
new file mode 100644
--- /dev/null
+++ b/dom/media/webaudio/blink/IIRFilter.cpp
@@ -0,0 +1,133 @@
+// Copyright 2016 The Chromium Authors. All rights reserved.
+// Use of this source code is governed by a BSD-style license that can be
+// found in the LICENSE file.
+
+#include "platform/audio/IIRFilter.h"
+
+#include "wtf/MathExtras.h"
+#include <complex>
+
+namespace blink {
+
+// The length of the memory buffers for the IIR filter. This MUST be a power of two and must be
+// greater than the possible length of the filter coefficients.
+const int kBufferLength = 32;
+static_assert(kBufferLength >= IIRFilter::kMaxOrder + 1,
+ "Internal IIR buffer length must be greater than maximum IIR Filter order.");
+
+IIRFilter::IIRFilter(const AudioDoubleArray* feedforward, const AudioDoubleArray* feedback)
+ : m_bufferIndex(0)
+ , m_feedback(feedback)
+ , m_feedforward(feedforward)
+{
+ // These are guaranteed to be zero-initialized.
+ m_xBuffer.allocate(kBufferLength);
+ m_yBuffer.allocate(kBufferLength);
+}
+
+IIRFilter::~IIRFilter()
+{
+}
+
+void IIRFilter::reset()
+{
+ m_xBuffer.zero();
+ m_yBuffer.zero();
+}
+
+static std::complex<double> evaluatePolynomial(const double* coef, std::complex<double> z, int order)
+{
+ // Use Horner's method to evaluate the polynomial P(z) = sum(coef[k]*z^k, k, 0, order);
+ std::complex<double> result = 0;
+
+ for (int k = order; k >= 0; --k)
+ result = result * z + std::complex<double>(coef[k]);
+
+ return result;
+}
+
+void IIRFilter::process(const float* sourceP, float* destP, size_t framesToProcess)
+{
+ // Compute
+ //
+ // y[n] = sum(b[k] * x[n - k], k = 0, M) - sum(a[k] * y[n - k], k = 1, N)
+ //
+ // where b[k] are the feedforward coefficients and a[k] are the feedback coefficients of the
+ // filter.
+
+ // This is a Direct Form I implementation of an IIR Filter. Should we consider doing a
+ // different implementation such as Transposed Direct Form II?
+ const double* feedback = m_feedback->data();
+ const double* feedforward = m_feedforward->data();
+
+ ASSERT(feedback);
+ ASSERT(feedforward);
+
+ // Sanity check to see if the feedback coefficients have been scaled appropriately. It must
+ // be EXACTLY 1!
+ ASSERT(feedback[0] == 1);
+
+ int feedbackLength = m_feedback->size();
+ int feedforwardLength = m_feedforward->size();
+ int minLength = std::min(feedbackLength, feedforwardLength);
+
+ double* xBuffer = m_xBuffer.data();
+ double* yBuffer = m_yBuffer.data();
+
+ for (size_t n = 0; n < framesToProcess; ++n) {
+ // To help minimize roundoff, we compute using double's, even though the filter coefficients
+ // only have single precision values.
+ double yn = feedforward[0] * sourceP[n];
+
+ // Run both the feedforward and feedback terms together, when possible.
+ for (int k = 1; k < minLength; ++k) {
+ int n = (m_bufferIndex - k) & (kBufferLength - 1);
+ yn += feedforward[k] * xBuffer[n];
+ yn -= feedback[k] * yBuffer[n];
+ }
+
+ // Handle any remaining feedforward or feedback terms.
+ for (int k = minLength; k < feedforwardLength; ++k)
+ yn += feedforward[k] * xBuffer[(m_bufferIndex - k) & (kBufferLength - 1)];
+
+ for (int k = minLength; k < feedbackLength; ++k)
+ yn -= feedback[k] * yBuffer[(m_bufferIndex - k) & (kBufferLength - 1)];
+
+ // Save the current input and output values in the memory buffers for the next output.
+ m_xBuffer[m_bufferIndex] = sourceP[n];
+ m_yBuffer[m_bufferIndex] = yn;
+
+ m_bufferIndex = (m_bufferIndex + 1) & (kBufferLength - 1);
+
+ destP[n] = yn;
+ }
+}
+
+void IIRFilter::getFrequencyResponse(int nFrequencies, const float* frequency, float* magResponse, float* phaseResponse)
+{
+ // Evaluate the z-transform of the filter at the given normalized frequencies from 0 to 1. (One
+ // corresponds to the Nyquist frequency.)
+ //
+ // The z-tranform of the filter is
+ //
+ // H(z) = sum(b[k]*z^(-k), k, 0, M) / sum(a[k]*z^(-k), k, 0, N);
+ //
+ // The desired frequency response is H(exp(j*omega)), where omega is in [0, 1).
+ //
+ // Let P(x) = sum(c[k]*x^k, k, 0, P) be a polynomial of order P. Then each of the sums in H(z)
+ // is equivalent to evaluating a polynomial at the point 1/z.
+
+ for (int k = 0; k < nFrequencies; ++k) {
+ // zRecip = 1/z = exp(-j*frequency)
+ double omega = -piDouble * frequency[k];
+ std::complex<double> zRecip = std::complex<double>(cos(omega), sin(omega));
+
+ std::complex<double> numerator = evaluatePolynomial(m_feedforward->data(), zRecip, m_feedforward->size() - 1);
+ std::complex<double> denominator = evaluatePolynomial(m_feedback->data(), zRecip, m_feedback->size() - 1);
+ std::complex<double> response = numerator / denominator;
+ magResponse[k] = static_cast<float>(abs(response));
+ phaseResponse[k] = static_cast<float>(atan2(imag(response), real(response)));
+ }
+}
+
+} // namespace blink
new file mode 100644
--- /dev/null
+++ b/dom/media/webaudio/blink/IIRFilter.h
@@ -0,0 +1,59 @@
+// Copyright 2016 The Chromium Authors. All rights reserved.
+// Use of this source code is governed by a BSD-style license that can be
+// found in the LICENSE file.
+
+#ifndef IIRFilter_h
+#define IIRFilter_h
+
+#include "platform/PlatformExport.h"
+#include "platform/audio/AudioArray.h"
+#include "wtf/Vector.h"
+
+namespace blink {
+
+class PLATFORM_EXPORT IIRFilter final {
+public:
+ // The maximum IIR filter order. This also limits the number of feedforward coefficients. The
+ // maximum number of coefficients is 20 according to the spec.
+ const static size_t kMaxOrder = 19;
+ IIRFilter(const AudioDoubleArray* feedforwardCoef, const AudioDoubleArray* feedbackCoef);
+ ~IIRFilter();
+
+ void process(const float* sourceP, float* destP, size_t framesToProcess);
+
+ void reset();
+
+ void getFrequencyResponse(int nFrequencies,
+ const float* frequency,
+ float* magResponse,
+ float* phaseResponse);
+
+private:
+ // Filter memory
+ //
+ // For simplicity, we assume |m_xBuffer| and |m_yBuffer| have the same length, and the length is
+ // a power of two. Since the number of coefficients has a fixed upper length, the size of
+ // xBuffer and yBuffer is fixed. |m_xBuffer| holds the old input values and |m_yBuffer| holds
+ // the old output values needed to compute the new output value.
+ //
+ // m_yBuffer[m_bufferIndex] holds the most recent output value, say, y[n]. Then
+ // m_yBuffer[m_bufferIndex - k] is y[n - k]. Similarly for m_xBuffer.
+ //
+ // To minimize roundoff, these arrays are double's instead of floats.
+ AudioDoubleArray m_xBuffer;
+ AudioDoubleArray m_yBuffer;
+
+ // Index into the xBuffer and yBuffer arrays where the most current x and y values should be
+ // stored. xBuffer[bufferIndex] corresponds to x[n], the current x input value and
+ // yBuffer[bufferIndex] is where y[n], the current output value.
+ int m_bufferIndex;
+
+ // Coefficients of the IIR filter. To minimize storage, these point to the arrays given in the
+ // constructor.
+ const AudioDoubleArray* m_feedback;
+ const AudioDoubleArray* m_feedforward;
+};
+
+} // namespace blink
+
+#endif // IIRFilter_h