1/*
2 * Copyright (C) 2011 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 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 "DynamicsCompressorKernel.h"
34
35#include "AudioUtilities.h"
36#include "DenormalDisabler.h"
37#include <algorithm>
38#include <wtf/MathExtras.h>
39
40namespace WebCore {
41
42using namespace AudioUtilities;
43
44// Metering hits peaks instantly, but releases this fast (in seconds).
45const float meteringReleaseTimeConstant = 0.325f;
46
47const float uninitializedValue = -1;
48
49DynamicsCompressorKernel::DynamicsCompressorKernel(float sampleRate, unsigned numberOfChannels)
50 : m_sampleRate(sampleRate)
51 , m_lastPreDelayFrames(DefaultPreDelayFrames)
52 , m_preDelayReadIndex(0)
53 , m_preDelayWriteIndex(DefaultPreDelayFrames)
54 , m_ratio(uninitializedValue)
55 , m_slope(uninitializedValue)
56 , m_linearThreshold(uninitializedValue)
57 , m_dbThreshold(uninitializedValue)
58 , m_dbKnee(uninitializedValue)
59 , m_kneeThreshold(uninitializedValue)
60 , m_kneeThresholdDb(uninitializedValue)
61 , m_ykneeThresholdDb(uninitializedValue)
62 , m_K(uninitializedValue)
63{
64 setNumberOfChannels(numberOfChannels);
65
66 // Initializes most member variables
67 reset();
68
69 m_meteringReleaseK = static_cast<float>(discreteTimeConstantForSampleRate(meteringReleaseTimeConstant, sampleRate));
70}
71
72void DynamicsCompressorKernel::setNumberOfChannels(unsigned numberOfChannels)
73{
74 if (m_preDelayBuffers.size() == numberOfChannels)
75 return;
76
77 m_preDelayBuffers.clear();
78 for (unsigned i = 0; i < numberOfChannels; ++i)
79 m_preDelayBuffers.append(std::make_unique<AudioFloatArray>(MaxPreDelayFrames));
80}
81
82void DynamicsCompressorKernel::setPreDelayTime(float preDelayTime)
83{
84 // Re-configure look-ahead section pre-delay if delay time has changed.
85 unsigned preDelayFrames = preDelayTime * sampleRate();
86 if (preDelayFrames > MaxPreDelayFrames - 1)
87 preDelayFrames = MaxPreDelayFrames - 1;
88
89 if (m_lastPreDelayFrames != preDelayFrames) {
90 m_lastPreDelayFrames = preDelayFrames;
91 for (unsigned i = 0; i < m_preDelayBuffers.size(); ++i)
92 m_preDelayBuffers[i]->zero();
93
94 m_preDelayReadIndex = 0;
95 m_preDelayWriteIndex = preDelayFrames;
96 }
97}
98
99// Exponential curve for the knee.
100// It is 1st derivative matched at m_linearThreshold and asymptotically approaches the value m_linearThreshold + 1 / k.
101float DynamicsCompressorKernel::kneeCurve(float x, float k)
102{
103 // Linear up to threshold.
104 if (x < m_linearThreshold)
105 return x;
106
107 return m_linearThreshold + (1 - expf(-k * (x - m_linearThreshold))) / k;
108}
109
110// Full compression curve with constant ratio after knee.
111float DynamicsCompressorKernel::saturate(float x, float k)
112{
113 float y;
114
115 if (x < m_kneeThreshold)
116 y = kneeCurve(x, k);
117 else {
118 // Constant ratio after knee.
119 float xDb = linearToDecibels(x);
120 float yDb = m_ykneeThresholdDb + m_slope * (xDb - m_kneeThresholdDb);
121
122 y = decibelsToLinear(yDb);
123 }
124
125 return y;
126}
127
128// Approximate 1st derivative with input and output expressed in dB.
129// This slope is equal to the inverse of the compression "ratio".
130// In other words, a compression ratio of 20 would be a slope of 1/20.
131float DynamicsCompressorKernel::slopeAt(float x, float k)
132{
133 if (x < m_linearThreshold)
134 return 1;
135
136 float x2 = x * 1.001;
137
138 float xDb = linearToDecibels(x);
139 float x2Db = linearToDecibels(x2);
140
141 float yDb = linearToDecibels(kneeCurve(x, k));
142 float y2Db = linearToDecibels(kneeCurve(x2, k));
143
144 float m = (y2Db - yDb) / (x2Db - xDb);
145
146 return m;
147}
148
149float DynamicsCompressorKernel::kAtSlope(float desiredSlope)
150{
151 float xDb = m_dbThreshold + m_dbKnee;
152 float x = decibelsToLinear(xDb);
153
154 // Approximate k given initial values.
155 float minK = 0.1;
156 float maxK = 10000;
157 float k = 5;
158
159 for (int i = 0; i < 15; ++i) {
160 // A high value for k will more quickly asymptotically approach a slope of 0.
161 float slope = slopeAt(x, k);
162
163 if (slope < desiredSlope) {
164 // k is too high.
165 maxK = k;
166 } else {
167 // k is too low.
168 minK = k;
169 }
170
171 // Re-calculate based on geometric mean.
172 k = sqrtf(minK * maxK);
173 }
174
175 return k;
176}
177
178float DynamicsCompressorKernel::updateStaticCurveParameters(float dbThreshold, float dbKnee, float ratio)
179{
180 if (dbThreshold != m_dbThreshold || dbKnee != m_dbKnee || ratio != m_ratio) {
181 // Threshold and knee.
182 m_dbThreshold = dbThreshold;
183 m_linearThreshold = decibelsToLinear(dbThreshold);
184 m_dbKnee = dbKnee;
185
186 // Compute knee parameters.
187 m_ratio = ratio;
188 m_slope = 1 / m_ratio;
189
190 float k = kAtSlope(1 / m_ratio);
191
192 m_kneeThresholdDb = dbThreshold + dbKnee;
193 m_kneeThreshold = decibelsToLinear(m_kneeThresholdDb);
194
195 m_ykneeThresholdDb = linearToDecibels(kneeCurve(m_kneeThreshold, k));
196
197 m_K = k;
198 }
199 return m_K;
200}
201
202void DynamicsCompressorKernel::process(float* sourceChannels[],
203 float* destinationChannels[],
204 unsigned numberOfChannels,
205 unsigned framesToProcess,
206
207 float dbThreshold,
208 float dbKnee,
209 float ratio,
210 float attackTime,
211 float releaseTime,
212 float preDelayTime,
213 float dbPostGain,
214 float effectBlend, /* equal power crossfade */
215
216 float releaseZone1,
217 float releaseZone2,
218 float releaseZone3,
219 float releaseZone4
220 )
221{
222 ASSERT(m_preDelayBuffers.size() == numberOfChannels);
223
224 float sampleRate = this->sampleRate();
225
226 float dryMix = 1 - effectBlend;
227 float wetMix = effectBlend;
228
229 float k = updateStaticCurveParameters(dbThreshold, dbKnee, ratio);
230
231 // Makeup gain.
232 float fullRangeGain = saturate(1, k);
233 float fullRangeMakeupGain = 1 / fullRangeGain;
234
235 // Empirical/perceptual tuning.
236 fullRangeMakeupGain = powf(fullRangeMakeupGain, 0.6f);
237
238 float masterLinearGain = decibelsToLinear(dbPostGain) * fullRangeMakeupGain;
239
240 // Attack parameters.
241 attackTime = std::max(0.001f, attackTime);
242 float attackFrames = attackTime * sampleRate;
243
244 // Release parameters.
245 float releaseFrames = sampleRate * releaseTime;
246
247 // Detector release time.
248 float satReleaseTime = 0.0025f;
249 float satReleaseFrames = satReleaseTime * sampleRate;
250
251 // Create a smooth function which passes through four points.
252
253 // Polynomial of the form
254 // y = a + b*x + c*x^2 + d*x^3 + e*x^4;
255
256 float y1 = releaseFrames * releaseZone1;
257 float y2 = releaseFrames * releaseZone2;
258 float y3 = releaseFrames * releaseZone3;
259 float y4 = releaseFrames * releaseZone4;
260
261 // All of these coefficients were derived for 4th order polynomial curve fitting where the y values
262 // match the evenly spaced x values as follows: (y1 : x == 0, y2 : x == 1, y3 : x == 2, y4 : x == 3)
263 float kA = 0.9999999999999998f*y1 + 1.8432219684323923e-16f*y2 - 1.9373394351676423e-16f*y3 + 8.824516011816245e-18f*y4;
264 float kB = -1.5788320352845888f*y1 + 2.3305837032074286f*y2 - 0.9141194204840429f*y3 + 0.1623677525612032f*y4;
265 float kC = 0.5334142869106424f*y1 - 1.272736789213631f*y2 + 0.9258856042207512f*y3 - 0.18656310191776226f*y4;
266 float kD = 0.08783463138207234f*y1 - 0.1694162967925622f*y2 + 0.08588057951595272f*y3 - 0.00429891410546283f*y4;
267 float kE = -0.042416883008123074f*y1 + 0.1115693827987602f*y2 - 0.09764676325265872f*y3 + 0.028494263462021576f*y4;
268
269 // x ranges from 0 -> 3 0 1 2 3
270 // -15 -10 -5 0db
271
272 // y calculates adaptive release frames depending on the amount of compression.
273
274 setPreDelayTime(preDelayTime);
275
276 const int nDivisionFrames = 32;
277
278 const int nDivisions = framesToProcess / nDivisionFrames;
279
280 unsigned frameIndex = 0;
281 for (int i = 0; i < nDivisions; ++i) {
282 // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
283 // Calculate desired gain
284 // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
285
286 // Fix gremlins.
287 if (std::isnan(m_detectorAverage))
288 m_detectorAverage = 1;
289 if (std::isinf(m_detectorAverage))
290 m_detectorAverage = 1;
291
292 float desiredGain = m_detectorAverage;
293
294 // Pre-warp so we get desiredGain after sin() warp below.
295 float scaledDesiredGain = asinf(desiredGain) / (0.5f * piFloat);
296
297 // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
298 // Deal with envelopes
299 // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
300
301 // envelopeRate is the rate we slew from current compressor level to the desired level.
302 // The exact rate depends on if we're attacking or releasing and by how much.
303 float envelopeRate;
304
305 bool isReleasing = scaledDesiredGain > m_compressorGain;
306
307 // compressionDiffDb is the difference between current compression level and the desired level.
308 float compressionDiffDb = linearToDecibels(m_compressorGain / scaledDesiredGain);
309
310 if (isReleasing) {
311 // Release mode - compressionDiffDb should be negative dB
312 m_maxAttackCompressionDiffDb = -1;
313
314 // Fix gremlins.
315 if (std::isnan(compressionDiffDb))
316 compressionDiffDb = -1;
317 if (std::isinf(compressionDiffDb))
318 compressionDiffDb = -1;
319
320 // Adaptive release - higher compression (lower compressionDiffDb) releases faster.
321
322 // Contain within range: -12 -> 0 then scale to go from 0 -> 3
323 float x = compressionDiffDb;
324 x = std::max(-12.0f, x);
325 x = std::min(0.0f, x);
326 x = 0.25f * (x + 12);
327
328 // Compute adaptive release curve using 4th order polynomial.
329 // Normal values for the polynomial coefficients would create a monotonically increasing function.
330 float x2 = x * x;
331 float x3 = x2 * x;
332 float x4 = x2 * x2;
333 float releaseFrames = kA + kB * x + kC * x2 + kD * x3 + kE * x4;
334
335#define kSpacingDb 5
336 float dbPerFrame = kSpacingDb / releaseFrames;
337
338 envelopeRate = decibelsToLinear(dbPerFrame);
339 } else {
340 // Attack mode - compressionDiffDb should be positive dB
341
342 // Fix gremlins.
343 if (std::isnan(compressionDiffDb))
344 compressionDiffDb = 1;
345 if (std::isinf(compressionDiffDb))
346 compressionDiffDb = 1;
347
348 // As long as we're still in attack mode, use a rate based off
349 // the largest compressionDiffDb we've encountered so far.
350 if (m_maxAttackCompressionDiffDb == -1 || m_maxAttackCompressionDiffDb < compressionDiffDb)
351 m_maxAttackCompressionDiffDb = compressionDiffDb;
352
353 float effAttenDiffDb = std::max(0.5f, m_maxAttackCompressionDiffDb);
354
355 float x = 0.25f / effAttenDiffDb;
356 envelopeRate = 1 - powf(x, 1 / attackFrames);
357 }
358
359 // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
360 // Inner loop - calculate shaped power average - apply compression.
361 // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
362
363 {
364 int preDelayReadIndex = m_preDelayReadIndex;
365 int preDelayWriteIndex = m_preDelayWriteIndex;
366 float detectorAverage = m_detectorAverage;
367 float compressorGain = m_compressorGain;
368
369 int loopFrames = nDivisionFrames;
370 while (loopFrames--) {
371 float compressorInput = 0;
372
373 // Predelay signal, computing compression amount from un-delayed version.
374 for (unsigned i = 0; i < numberOfChannels; ++i) {
375 float* delayBuffer = m_preDelayBuffers[i]->data();
376 float undelayedSource = sourceChannels[i][frameIndex];
377 delayBuffer[preDelayWriteIndex] = undelayedSource;
378
379 float absUndelayedSource = undelayedSource > 0 ? undelayedSource : -undelayedSource;
380 if (compressorInput < absUndelayedSource)
381 compressorInput = absUndelayedSource;
382 }
383
384 // Calculate shaped power on undelayed input.
385
386 float scaledInput = compressorInput;
387 float absInput = scaledInput > 0 ? scaledInput : -scaledInput;
388
389 // Put through shaping curve.
390 // This is linear up to the threshold, then enters a "knee" portion followed by the "ratio" portion.
391 // The transition from the threshold to the knee is smooth (1st derivative matched).
392 // The transition from the knee to the ratio portion is smooth (1st derivative matched).
393 float shapedInput = saturate(absInput, k);
394
395 float attenuation = absInput <= 0.0001f ? 1 : shapedInput / absInput;
396
397 float attenuationDb = -linearToDecibels(attenuation);
398 attenuationDb = std::max(2.0f, attenuationDb);
399
400 float dbPerFrame = attenuationDb / satReleaseFrames;
401
402 float satReleaseRate = decibelsToLinear(dbPerFrame) - 1;
403
404 bool isRelease = (attenuation > detectorAverage);
405 float rate = isRelease ? satReleaseRate : 1;
406
407 detectorAverage += (attenuation - detectorAverage) * rate;
408 detectorAverage = std::min(1.0f, detectorAverage);
409
410 // Fix gremlins.
411 if (std::isnan(detectorAverage))
412 detectorAverage = 1;
413 if (std::isinf(detectorAverage))
414 detectorAverage = 1;
415
416 // Exponential approach to desired gain.
417 if (envelopeRate < 1) {
418 // Attack - reduce gain to desired.
419 compressorGain += (scaledDesiredGain - compressorGain) * envelopeRate;
420 } else {
421 // Release - exponentially increase gain to 1.0
422 compressorGain *= envelopeRate;
423 compressorGain = std::min(1.0f, compressorGain);
424 }
425
426 // Warp pre-compression gain to smooth out sharp exponential transition points.
427 float postWarpCompressorGain = sinf(0.5f * piFloat * compressorGain);
428
429 // Calculate total gain using master gain and effect blend.
430 float totalGain = dryMix + wetMix * masterLinearGain * postWarpCompressorGain;
431
432 // Calculate metering.
433 float dbRealGain = 20 * log10(postWarpCompressorGain);
434 if (dbRealGain < m_meteringGain)
435 m_meteringGain = dbRealGain;
436 else
437 m_meteringGain += (dbRealGain - m_meteringGain) * m_meteringReleaseK;
438
439 // Apply final gain.
440 for (unsigned i = 0; i < numberOfChannels; ++i) {
441 float* delayBuffer = m_preDelayBuffers[i]->data();
442 destinationChannels[i][frameIndex] = delayBuffer[preDelayReadIndex] * totalGain;
443 }
444
445 frameIndex++;
446 preDelayReadIndex = (preDelayReadIndex + 1) & MaxPreDelayFramesMask;
447 preDelayWriteIndex = (preDelayWriteIndex + 1) & MaxPreDelayFramesMask;
448 }
449
450 // Locals back to member variables.
451 m_preDelayReadIndex = preDelayReadIndex;
452 m_preDelayWriteIndex = preDelayWriteIndex;
453 m_detectorAverage = DenormalDisabler::flushDenormalFloatToZero(detectorAverage);
454 m_compressorGain = DenormalDisabler::flushDenormalFloatToZero(compressorGain);
455 }
456 }
457}
458
459void DynamicsCompressorKernel::reset()
460{
461 m_detectorAverage = 0;
462 m_compressorGain = 1;
463 m_meteringGain = 1;
464
465 // Predelay section.
466 for (unsigned i = 0; i < m_preDelayBuffers.size(); ++i)
467 m_preDelayBuffers[i]->zero();
468
469 m_preDelayReadIndex = 0;
470 m_preDelayWriteIndex = DefaultPreDelayFrames;
471
472 m_maxAttackCompressionDiffDb = -1; // uninitialized state
473}
474
475} // namespace WebCore
476
477#endif // ENABLE(WEB_AUDIO)
478