1 /* 2 * Copyright (C) 2010 Google Inc. All rights reserved. 3 * 4 * Redistribution and use in source and binary forms, with or without 5 * modification, are permitted provided that the following conditions 6 * are met: 7 * 8 * 1. Redistributions of source code must retain the above copyright 9 * notice, this list of conditions and the following disclaimer. 10 * 2. Redistributions in binary form must reproduce the above copyright 11 * notice, this list of conditions and the following disclaimer in the 12 * documentation and/or other materials provided with the distribution. 13 * 3. Neither the name of Apple Computer, Inc. ("Apple") nor the names of 14 * its contributors may be used to endorse or promote products derived 15 * from this software without specific prior written permission. 16 * 17 * THIS SOFTWARE IS PROVIDED BY APPLE AND ITS CONTRIBUTORS "AS IS" AND ANY 18 * EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED 19 * WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE 20 * DISCLAIMED. IN NO EVENT SHALL APPLE OR ITS CONTRIBUTORS BE LIABLE FOR ANY 21 * DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES 22 * (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; 23 * LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND 24 * ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT 25 * (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF 26 * THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE. 27 */ 28 29 #include "config.h" 30 31 #if ENABLE(WEB_AUDIO) 32 33 #include "Biquad.h" 34 35 #include <algorithm> 36 #include <stdio.h> 37 #include <wtf/MathExtras.h> 38 39 #if OS(DARWIN) 40 #include <Accelerate/Accelerate.h> 41 #endif 42 43 namespace WebCore { 44 45 const int kBufferSize = 1024; 46 47 Biquad::Biquad() 48 { 49 #if OS(DARWIN) 50 // Allocate two samples more for filter history 51 m_inputBuffer.resize(kBufferSize + 2); 52 m_outputBuffer.resize(kBufferSize + 2); 53 #endif 54 55 // Initialize as pass-thru (straight-wire, no filter effect) 56 m_a0 = 1.0; 57 m_a1 = 0.0; 58 m_a2 = 0.0; 59 m_b1 = 0.0; 60 m_b2 = 0.0; 61 62 m_g = 1.0; 63 64 reset(); // clear filter memory 65 } 66 67 void Biquad::process(const float* sourceP, float* destP, size_t framesToProcess) 68 { 69 #if OS(DARWIN) 70 // Use vecLib if available 71 processFast(sourceP, destP, framesToProcess); 72 #else 73 int n = framesToProcess; 74 75 // Create local copies of member variables 76 double x1 = m_x1; 77 double x2 = m_x2; 78 double y1 = m_y1; 79 double y2 = m_y2; 80 81 double a0 = m_a0; 82 double a1 = m_a1; 83 double a2 = m_a2; 84 double b1 = m_b1; 85 double b2 = m_b2; 86 87 while (n--) { 88 // FIXME: this can be optimized by pipelining the multiply adds... 89 float x = *sourceP++; 90 float y = a0*x + a1*x1 + a2*x2 - b1*y1 - b2*y2; 91 92 y *= m_g; 93 94 *destP++ = y; 95 96 // Update state variables 97 x2 = x1; 98 x1 = x; 99 y2 = y1; 100 y1 = y; 101 } 102 103 // Local variables back to member 104 m_x1 = x1; 105 m_x2 = x2; 106 m_y1 = y1; 107 m_y2 = y2; 108 109 m_a0 = a0; 110 m_a1 = a1; 111 m_a2 = a2; 112 m_b1 = b1; 113 m_b2 = b2; 114 #endif 115 } 116 117 #if OS(DARWIN) 118 119 // Here we have optimized version using Accelerate.framework 120 121 void Biquad::processFast(const float* sourceP, float* destP, size_t framesToProcess) 122 { 123 // Filter coefficients 124 double B[5]; 125 B[0] = m_a0; 126 B[1] = m_a1; 127 B[2] = m_a2; 128 B[3] = m_b1; 129 B[4] = m_b2; 130 131 double* inputP = m_inputBuffer.data(); 132 double* outputP = m_outputBuffer.data(); 133 134 double* input2P = inputP + 2; 135 double* output2P = outputP + 2; 136 137 // Break up processing into smaller slices (kBufferSize) if necessary. 138 139 int n = framesToProcess; 140 141 while (n > 0) { 142 int framesThisTime = n < kBufferSize ? n : kBufferSize; 143 144 // Copy input to input buffer 145 for (int i = 0; i < framesThisTime; ++i) 146 input2P[i] = *sourceP++; 147 148 processSliceFast(inputP, outputP, B, framesThisTime); 149 150 // Copy output buffer to output (converts float -> double). 151 for (int i = 0; i < framesThisTime; ++i) 152 *destP++ = static_cast<float>(output2P[i]); 153 154 n -= framesThisTime; 155 } 156 } 157 158 void Biquad::processSliceFast(double* sourceP, double* destP, double* coefficientsP, size_t framesToProcess) 159 { 160 // Use double-precision for filter stability 161 vDSP_deq22D(sourceP, 1, coefficientsP, destP, 1, framesToProcess); 162 163 // Save history. Note that sourceP and destP reference m_inputBuffer and m_outputBuffer respectively. 164 // These buffers are allocated (in the constructor) with space for two extra samples so it's OK to access 165 // array values two beyond framesToProcess. 166 sourceP[0] = sourceP[framesToProcess - 2 + 2]; 167 sourceP[1] = sourceP[framesToProcess - 1 + 2]; 168 destP[0] = destP[framesToProcess - 2 + 2]; 169 destP[1] = destP[framesToProcess - 1 + 2]; 170 } 171 172 #endif // OS(DARWIN) 173 174 175 void Biquad::reset() 176 { 177 m_x1 = m_x2 = m_y1 = m_y2 = 0.0; 178 179 #if OS(DARWIN) 180 // Two extra samples for filter history 181 double* inputP = m_inputBuffer.data(); 182 inputP[0] = 0.0; 183 inputP[1] = 0.0; 184 185 double* outputP = m_outputBuffer.data(); 186 outputP[0] = 0.0; 187 outputP[1] = 0.0; 188 #endif 189 } 190 191 void Biquad::setLowpassParams(double cutoff, double resonance) 192 { 193 resonance = std::max(0.0, resonance); // can't go negative 194 195 double g = pow(10.0, 0.05 * resonance); 196 double d = sqrt((4.0 - sqrt(16.0 - 16.0 / (g * g))) / 2.0); 197 198 // Compute biquad coefficients for lopass filter 199 double theta = piDouble * cutoff; 200 double sn = 0.5 * d * sin(theta); 201 double beta = 0.5 * (1.0 - sn) / (1.0 + sn); 202 double gamma = (0.5 + beta) * cos(theta); 203 double alpha = 0.25 * (0.5 + beta - gamma); 204 205 m_a0 = 2.0 * alpha; 206 m_a1 = 2.0 * 2.0*alpha; 207 m_a2 = 2.0 * alpha; 208 m_b1 = 2.0 * -gamma; 209 m_b2 = 2.0 * beta; 210 } 211 212 void Biquad::setHighpassParams(double cutoff, double resonance) 213 { 214 resonance = std::max(0.0, resonance); // can't go negative 215 216 double g = pow(10.0, 0.05 * resonance); 217 double d = sqrt((4.0 - sqrt(16.0 - 16.0 / (g * g))) / 2.0); 218 219 // Compute biquad coefficients for highpass filter 220 double theta = piDouble * cutoff; 221 double sn = 0.5 * d * sin(theta); 222 double beta = 0.5 * (1.0 - sn) / (1.0 + sn); 223 double gamma = (0.5 + beta) * cos(theta); 224 double alpha = 0.25 * (0.5 + beta + gamma); 225 226 m_a0 = 2.0 * alpha; 227 m_a1 = 2.0 * -2.0*alpha; 228 m_a2 = 2.0 * alpha; 229 m_b1 = 2.0 * -gamma; 230 m_b2 = 2.0 * beta; 231 } 232 233 void Biquad::setLowShelfParams(double cutoff, double dbGain) 234 { 235 double theta = piDouble * cutoff; 236 237 double A = pow(10.0, dbGain / 40.0); 238 double S = 1.0; // filter slope (1.0 is max value) 239 double alpha = 0.5 * sin(theta) * sqrt((A + 1.0 / A) * (1.0 / S - 1.0) + 2.0); 240 241 double k = cos(theta); 242 double k2 = 2.0 * sqrt(A) * alpha; 243 244 double b0 = A * ((A + 1.0) - (A - 1.0) * k + k2); 245 double b1 = 2.0 * A * ((A - 1.0) - (A + 1.0) * k); 246 double b2 = A * ((A + 1.0) - (A - 1.0) * k - k2); 247 double a0 = (A + 1.0) + (A - 1.0) * k + k2; 248 double a1 = -2.0 * ((A - 1.0) + (A + 1.0) * k); 249 double a2 = (A + 1.0) + (A - 1.0) * k - k2; 250 251 double a0Inverse = 1.0 / a0; 252 253 m_a0 = b0 * a0Inverse; 254 m_a1 = b1 * a0Inverse; 255 m_a2 = b2 * a0Inverse; 256 m_b1 = a1 * a0Inverse; 257 m_b2 = a2 * a0Inverse; 258 } 259 260 void Biquad::setZeroPolePairs(const Complex &zero, const Complex &pole) 261 { 262 m_a0 = 1.0; 263 m_a1 = -2.0 * zero.real(); 264 265 double zeroMag = abs(zero); 266 m_a2 = zeroMag * zeroMag; 267 268 m_b1 = -2.0 * pole.real(); 269 270 double poleMag = abs(pole); 271 m_b2 = poleMag * poleMag; 272 } 273 274 void Biquad::setAllpassPole(const Complex &pole) 275 { 276 Complex zero = Complex(1.0, 0.0) / pole; 277 setZeroPolePairs(zero, pole); 278 } 279 280 } // namespace WebCore 281 282 #endif // ENABLE(WEB_AUDIO) 283
