Line data Source code
1 : /*
2 : * Copyright (c) 2011 The WebRTC project authors. All Rights Reserved.
3 : *
4 : * Use of this source code is governed by a BSD-style license
5 : * that can be found in the LICENSE file in the root of the source
6 : * tree. An additional intellectual property rights grant can be found
7 : * in the file PATENTS. All contributing project authors may
8 : * be found in the AUTHORS file in the root of the source tree.
9 : */
10 :
11 : /* digital_agc.c
12 : *
13 : */
14 :
15 : #include "webrtc/modules/audio_processing/agc/legacy/digital_agc.h"
16 :
17 : #include <string.h>
18 : #ifdef WEBRTC_AGC_DEBUG_DUMP
19 : #include <stdio.h>
20 : #endif
21 :
22 : #include "webrtc/base/checks.h"
23 : #include "webrtc/modules/audio_processing/agc/legacy/gain_control.h"
24 :
25 : // To generate the gaintable, copy&paste the following lines to a Matlab window:
26 : // MaxGain = 6; MinGain = 0; CompRatio = 3; Knee = 1;
27 : // zeros = 0:31; lvl = 2.^(1-zeros);
28 : // A = -10*log10(lvl) * (CompRatio - 1) / CompRatio;
29 : // B = MaxGain - MinGain;
30 : // gains = round(2^16*10.^(0.05 * (MinGain + B * (
31 : // log(exp(-Knee*A)+exp(-Knee*B)) - log(1+exp(-Knee*B)) ) /
32 : // log(1/(1+exp(Knee*B))))));
33 : // fprintf(1, '\t%i, %i, %i, %i,\n', gains);
34 : // % Matlab code for plotting the gain and input/output level characteristic
35 : // (copy/paste the following 3 lines):
36 : // in = 10*log10(lvl); out = 20*log10(gains/65536);
37 : // subplot(121); plot(in, out); axis([-30, 0, -5, 20]); grid on; xlabel('Input
38 : // (dB)'); ylabel('Gain (dB)');
39 : // subplot(122); plot(in, in+out); axis([-30, 0, -30, 5]); grid on;
40 : // xlabel('Input (dB)'); ylabel('Output (dB)');
41 : // zoom on;
42 :
43 : // Generator table for y=log2(1+e^x) in Q8.
44 : enum { kGenFuncTableSize = 128 };
45 : static const uint16_t kGenFuncTable[kGenFuncTableSize] = {
46 : 256, 485, 786, 1126, 1484, 1849, 2217, 2586, 2955, 3324, 3693,
47 : 4063, 4432, 4801, 5171, 5540, 5909, 6279, 6648, 7017, 7387, 7756,
48 : 8125, 8495, 8864, 9233, 9603, 9972, 10341, 10711, 11080, 11449, 11819,
49 : 12188, 12557, 12927, 13296, 13665, 14035, 14404, 14773, 15143, 15512, 15881,
50 : 16251, 16620, 16989, 17359, 17728, 18097, 18466, 18836, 19205, 19574, 19944,
51 : 20313, 20682, 21052, 21421, 21790, 22160, 22529, 22898, 23268, 23637, 24006,
52 : 24376, 24745, 25114, 25484, 25853, 26222, 26592, 26961, 27330, 27700, 28069,
53 : 28438, 28808, 29177, 29546, 29916, 30285, 30654, 31024, 31393, 31762, 32132,
54 : 32501, 32870, 33240, 33609, 33978, 34348, 34717, 35086, 35456, 35825, 36194,
55 : 36564, 36933, 37302, 37672, 38041, 38410, 38780, 39149, 39518, 39888, 40257,
56 : 40626, 40996, 41365, 41734, 42104, 42473, 42842, 43212, 43581, 43950, 44320,
57 : 44689, 45058, 45428, 45797, 46166, 46536, 46905};
58 :
59 : static const int16_t kAvgDecayTime = 250; // frames; < 3000
60 :
61 0 : int32_t WebRtcAgc_CalculateGainTable(int32_t* gainTable, // Q16
62 : int16_t digCompGaindB, // Q0
63 : int16_t targetLevelDbfs, // Q0
64 : uint8_t limiterEnable,
65 : int16_t analogTarget) // Q0
66 : {
67 : // This function generates the compressor gain table used in the fixed digital
68 : // part.
69 : uint32_t tmpU32no1, tmpU32no2, absInLevel, logApprox;
70 : int32_t inLevel, limiterLvl;
71 : int32_t tmp32, tmp32no1, tmp32no2, numFIX, den, y32;
72 0 : const uint16_t kLog10 = 54426; // log2(10) in Q14
73 0 : const uint16_t kLog10_2 = 49321; // 10*log10(2) in Q14
74 0 : const uint16_t kLogE_1 = 23637; // log2(e) in Q14
75 : uint16_t constMaxGain;
76 : uint16_t tmpU16, intPart, fracPart;
77 0 : const int16_t kCompRatio = 3;
78 0 : const int16_t kSoftLimiterLeft = 1;
79 0 : int16_t limiterOffset = 0; // Limiter offset
80 : int16_t limiterIdx, limiterLvlX;
81 : int16_t constLinApprox, zeroGainLvl, maxGain, diffGain;
82 : int16_t i, tmp16, tmp16no1;
83 : int zeros, zerosScale;
84 :
85 : // Constants
86 : // kLogE_1 = 23637; // log2(e) in Q14
87 : // kLog10 = 54426; // log2(10) in Q14
88 : // kLog10_2 = 49321; // 10*log10(2) in Q14
89 :
90 : // Calculate maximum digital gain and zero gain level
91 0 : tmp32no1 = (digCompGaindB - analogTarget) * (kCompRatio - 1);
92 0 : tmp16no1 = analogTarget - targetLevelDbfs;
93 0 : tmp16no1 +=
94 0 : WebRtcSpl_DivW32W16ResW16(tmp32no1 + (kCompRatio >> 1), kCompRatio);
95 0 : maxGain = WEBRTC_SPL_MAX(tmp16no1, (analogTarget - targetLevelDbfs));
96 0 : tmp32no1 = maxGain * kCompRatio;
97 0 : zeroGainLvl = digCompGaindB;
98 0 : zeroGainLvl -= WebRtcSpl_DivW32W16ResW16(tmp32no1 + ((kCompRatio - 1) >> 1),
99 0 : kCompRatio - 1);
100 0 : if ((digCompGaindB <= analogTarget) && (limiterEnable)) {
101 0 : zeroGainLvl += (analogTarget - digCompGaindB + kSoftLimiterLeft);
102 0 : limiterOffset = 0;
103 : }
104 :
105 : // Calculate the difference between maximum gain and gain at 0dB0v:
106 : // diffGain = maxGain + (compRatio-1)*zeroGainLvl/compRatio
107 : // = (compRatio-1)*digCompGaindB/compRatio
108 0 : tmp32no1 = digCompGaindB * (kCompRatio - 1);
109 0 : diffGain =
110 0 : WebRtcSpl_DivW32W16ResW16(tmp32no1 + (kCompRatio >> 1), kCompRatio);
111 0 : if (diffGain < 0 || diffGain >= kGenFuncTableSize) {
112 0 : RTC_DCHECK(0);
113 : return -1;
114 : }
115 :
116 : // Calculate the limiter level and index:
117 : // limiterLvlX = analogTarget - limiterOffset
118 : // limiterLvl = targetLevelDbfs + limiterOffset/compRatio
119 0 : limiterLvlX = analogTarget - limiterOffset;
120 0 : limiterIdx = 2 + WebRtcSpl_DivW32W16ResW16((int32_t)limiterLvlX * (1 << 13),
121 : kLog10_2 / 2);
122 0 : tmp16no1 =
123 0 : WebRtcSpl_DivW32W16ResW16(limiterOffset + (kCompRatio >> 1), kCompRatio);
124 0 : limiterLvl = targetLevelDbfs + tmp16no1;
125 :
126 : // Calculate (through table lookup):
127 : // constMaxGain = log2(1+2^(log2(e)*diffGain)); (in Q8)
128 0 : constMaxGain = kGenFuncTable[diffGain]; // in Q8
129 :
130 : // Calculate a parameter used to approximate the fractional part of 2^x with a
131 : // piecewise linear function in Q14:
132 : // constLinApprox = round(3/2*(4*(3-2*sqrt(2))/(log(2)^2)-0.5)*2^14);
133 0 : constLinApprox = 22817; // in Q14
134 :
135 : // Calculate a denominator used in the exponential part to convert from dB to
136 : // linear scale:
137 : // den = 20*constMaxGain (in Q8)
138 0 : den = WEBRTC_SPL_MUL_16_U16(20, constMaxGain); // in Q8
139 :
140 0 : for (i = 0; i < 32; i++) {
141 : // Calculate scaled input level (compressor):
142 : // inLevel =
143 : // fix((-constLog10_2*(compRatio-1)*(1-i)+fix(compRatio/2))/compRatio)
144 0 : tmp16 = (int16_t)((kCompRatio - 1) * (i - 1)); // Q0
145 0 : tmp32 = WEBRTC_SPL_MUL_16_U16(tmp16, kLog10_2) + 1; // Q14
146 0 : inLevel = WebRtcSpl_DivW32W16(tmp32, kCompRatio); // Q14
147 :
148 : // Calculate diffGain-inLevel, to map using the genFuncTable
149 0 : inLevel = (int32_t)diffGain * (1 << 14) - inLevel; // Q14
150 :
151 : // Make calculations on abs(inLevel) and compensate for the sign afterwards.
152 0 : absInLevel = (uint32_t)WEBRTC_SPL_ABS_W32(inLevel); // Q14
153 :
154 : // LUT with interpolation
155 0 : intPart = (uint16_t)(absInLevel >> 14);
156 0 : fracPart =
157 0 : (uint16_t)(absInLevel & 0x00003FFF); // extract the fractional part
158 0 : tmpU16 = kGenFuncTable[intPart + 1] - kGenFuncTable[intPart]; // Q8
159 0 : tmpU32no1 = tmpU16 * fracPart; // Q22
160 0 : tmpU32no1 += (uint32_t)kGenFuncTable[intPart] << 14; // Q22
161 0 : logApprox = tmpU32no1 >> 8; // Q14
162 : // Compensate for negative exponent using the relation:
163 : // log2(1 + 2^-x) = log2(1 + 2^x) - x
164 0 : if (inLevel < 0) {
165 0 : zeros = WebRtcSpl_NormU32(absInLevel);
166 0 : zerosScale = 0;
167 0 : if (zeros < 15) {
168 : // Not enough space for multiplication
169 0 : tmpU32no2 = absInLevel >> (15 - zeros); // Q(zeros-1)
170 0 : tmpU32no2 = WEBRTC_SPL_UMUL_32_16(tmpU32no2, kLogE_1); // Q(zeros+13)
171 0 : if (zeros < 9) {
172 0 : zerosScale = 9 - zeros;
173 0 : tmpU32no1 >>= zerosScale; // Q(zeros+13)
174 : } else {
175 0 : tmpU32no2 >>= zeros - 9; // Q22
176 : }
177 : } else {
178 0 : tmpU32no2 = WEBRTC_SPL_UMUL_32_16(absInLevel, kLogE_1); // Q28
179 0 : tmpU32no2 >>= 6; // Q22
180 : }
181 0 : logApprox = 0;
182 0 : if (tmpU32no2 < tmpU32no1) {
183 0 : logApprox = (tmpU32no1 - tmpU32no2) >> (8 - zerosScale); // Q14
184 : }
185 : }
186 0 : numFIX = (maxGain * constMaxGain) * (1 << 6); // Q14
187 0 : numFIX -= (int32_t)logApprox * diffGain; // Q14
188 :
189 : // Calculate ratio
190 : // Shift |numFIX| as much as possible.
191 : // Ensure we avoid wrap-around in |den| as well.
192 0 : if (numFIX > (den >> 8) || -numFIX > (den >> 8)) // |den| is Q8.
193 : {
194 0 : zeros = WebRtcSpl_NormW32(numFIX);
195 : } else {
196 0 : zeros = WebRtcSpl_NormW32(den) + 8;
197 : }
198 0 : numFIX *= 1 << zeros; // Q(14+zeros)
199 :
200 : // Shift den so we end up in Qy1
201 0 : tmp32no1 = WEBRTC_SPL_SHIFT_W32(den, zeros - 9); // Q(zeros - 1)
202 0 : y32 = numFIX / tmp32no1; // in Q15
203 : // This is to do rounding in Q14.
204 0 : y32 = y32 >= 0 ? (y32 + 1) >> 1 : -((-y32 + 1) >> 1);
205 :
206 0 : if (limiterEnable && (i < limiterIdx)) {
207 0 : tmp32 = WEBRTC_SPL_MUL_16_U16(i - 1, kLog10_2); // Q14
208 0 : tmp32 -= limiterLvl * (1 << 14); // Q14
209 0 : y32 = WebRtcSpl_DivW32W16(tmp32 + 10, 20);
210 : }
211 0 : if (y32 > 39000) {
212 0 : tmp32 = (y32 >> 1) * kLog10 + 4096; // in Q27
213 0 : tmp32 >>= 13; // In Q14.
214 : } else {
215 0 : tmp32 = y32 * kLog10 + 8192; // in Q28
216 0 : tmp32 >>= 14; // In Q14.
217 : }
218 0 : tmp32 += 16 << 14; // in Q14 (Make sure final output is in Q16)
219 :
220 : // Calculate power
221 0 : if (tmp32 > 0) {
222 0 : intPart = (int16_t)(tmp32 >> 14);
223 0 : fracPart = (uint16_t)(tmp32 & 0x00003FFF); // in Q14
224 0 : if ((fracPart >> 13) != 0) {
225 0 : tmp16 = (2 << 14) - constLinApprox;
226 0 : tmp32no2 = (1 << 14) - fracPart;
227 0 : tmp32no2 *= tmp16;
228 0 : tmp32no2 >>= 13;
229 0 : tmp32no2 = (1 << 14) - tmp32no2;
230 : } else {
231 0 : tmp16 = constLinApprox - (1 << 14);
232 0 : tmp32no2 = (fracPart * tmp16) >> 13;
233 : }
234 0 : fracPart = (uint16_t)tmp32no2;
235 0 : gainTable[i] =
236 0 : (1 << intPart) + WEBRTC_SPL_SHIFT_W32(fracPart, intPart - 14);
237 : } else {
238 0 : gainTable[i] = 0;
239 : }
240 : }
241 :
242 0 : return 0;
243 : }
244 :
245 0 : int32_t WebRtcAgc_InitDigital(DigitalAgc* stt, int16_t agcMode) {
246 0 : if (agcMode == kAgcModeFixedDigital) {
247 : // start at minimum to find correct gain faster
248 0 : stt->capacitorSlow = 0;
249 : } else {
250 : // start out with 0 dB gain
251 0 : stt->capacitorSlow = 134217728; // (int32_t)(0.125f * 32768.0f * 32768.0f);
252 : }
253 0 : stt->capacitorFast = 0;
254 0 : stt->gain = 65536;
255 0 : stt->gatePrevious = 0;
256 0 : stt->agcMode = agcMode;
257 : #ifdef WEBRTC_AGC_DEBUG_DUMP
258 : stt->frameCounter = 0;
259 : #endif
260 :
261 : // initialize VADs
262 0 : WebRtcAgc_InitVad(&stt->vadNearend);
263 0 : WebRtcAgc_InitVad(&stt->vadFarend);
264 :
265 0 : return 0;
266 : }
267 :
268 0 : int32_t WebRtcAgc_AddFarendToDigital(DigitalAgc* stt,
269 : const int16_t* in_far,
270 : size_t nrSamples) {
271 0 : RTC_DCHECK(stt);
272 : // VAD for far end
273 0 : WebRtcAgc_ProcessVad(&stt->vadFarend, in_far, nrSamples);
274 :
275 0 : return 0;
276 : }
277 :
278 0 : int32_t WebRtcAgc_ProcessDigital(DigitalAgc* stt,
279 : const int16_t* const* in_near,
280 : size_t num_bands,
281 : int16_t* const* out,
282 : uint32_t FS,
283 : int16_t lowlevelSignal) {
284 : // array for gains (one value per ms, incl start & end)
285 : int32_t gains[11];
286 :
287 : int32_t out_tmp, tmp32;
288 : int32_t env[10];
289 : int32_t max_nrg;
290 : int32_t cur_level;
291 : int32_t gain32, delta;
292 : int16_t logratio;
293 : int16_t lower_thr, upper_thr;
294 0 : int16_t zeros = 0, zeros_fast, frac = 0;
295 : int16_t decay;
296 : int16_t gate, gain_adj;
297 : int16_t k;
298 : size_t n, i, L;
299 : int16_t L2; // samples/subframe
300 :
301 : // determine number of samples per ms
302 0 : if (FS == 8000) {
303 0 : L = 8;
304 0 : L2 = 3;
305 0 : } else if (FS == 16000 || FS == 32000 || FS == 48000) {
306 0 : L = 16;
307 0 : L2 = 4;
308 : } else {
309 0 : return -1;
310 : }
311 :
312 0 : for (i = 0; i < num_bands; ++i) {
313 0 : if (in_near[i] != out[i]) {
314 : // Only needed if they don't already point to the same place.
315 0 : memcpy(out[i], in_near[i], 10 * L * sizeof(in_near[i][0]));
316 : }
317 : }
318 : // VAD for near end
319 0 : logratio = WebRtcAgc_ProcessVad(&stt->vadNearend, out[0], L * 10);
320 :
321 : // Account for far end VAD
322 0 : if (stt->vadFarend.counter > 10) {
323 0 : tmp32 = 3 * logratio;
324 0 : logratio = (int16_t)((tmp32 - stt->vadFarend.logRatio) >> 2);
325 : }
326 :
327 : // Determine decay factor depending on VAD
328 : // upper_thr = 1.0f;
329 : // lower_thr = 0.25f;
330 0 : upper_thr = 1024; // Q10
331 0 : lower_thr = 0; // Q10
332 0 : if (logratio > upper_thr) {
333 : // decay = -2^17 / DecayTime; -> -65
334 0 : decay = -65;
335 0 : } else if (logratio < lower_thr) {
336 0 : decay = 0;
337 : } else {
338 : // decay = (int16_t)(((lower_thr - logratio)
339 : // * (2^27/(DecayTime*(upper_thr-lower_thr)))) >> 10);
340 : // SUBSTITUTED: 2^27/(DecayTime*(upper_thr-lower_thr)) -> 65
341 0 : tmp32 = (lower_thr - logratio) * 65;
342 0 : decay = (int16_t)(tmp32 >> 10);
343 : }
344 :
345 : // adjust decay factor for long silence (detected as low standard deviation)
346 : // This is only done in the adaptive modes
347 0 : if (stt->agcMode != kAgcModeFixedDigital) {
348 0 : if (stt->vadNearend.stdLongTerm < 4000) {
349 0 : decay = 0;
350 0 : } else if (stt->vadNearend.stdLongTerm < 8096) {
351 : // decay = (int16_t)(((stt->vadNearend.stdLongTerm - 4000) * decay) >>
352 : // 12);
353 0 : tmp32 = (stt->vadNearend.stdLongTerm - 4000) * decay;
354 0 : decay = (int16_t)(tmp32 >> 12);
355 : }
356 :
357 0 : if (lowlevelSignal != 0) {
358 0 : decay = 0;
359 : }
360 : }
361 : #ifdef WEBRTC_AGC_DEBUG_DUMP
362 : stt->frameCounter++;
363 : fprintf(stt->logFile, "%5.2f\t%d\t%d\t%d\t", (float)(stt->frameCounter) / 100,
364 : logratio, decay, stt->vadNearend.stdLongTerm);
365 : #endif
366 : // Find max amplitude per sub frame
367 : // iterate over sub frames
368 0 : for (k = 0; k < 10; k++) {
369 : // iterate over samples
370 0 : max_nrg = 0;
371 0 : for (n = 0; n < L; n++) {
372 0 : int32_t nrg = out[0][k * L + n] * out[0][k * L + n];
373 0 : if (nrg > max_nrg) {
374 0 : max_nrg = nrg;
375 : }
376 : }
377 0 : env[k] = max_nrg;
378 : }
379 :
380 : // Calculate gain per sub frame
381 0 : gains[0] = stt->gain;
382 0 : for (k = 0; k < 10; k++) {
383 : // Fast envelope follower
384 : // decay time = -131000 / -1000 = 131 (ms)
385 0 : stt->capacitorFast =
386 0 : AGC_SCALEDIFF32(-1000, stt->capacitorFast, stt->capacitorFast);
387 0 : if (env[k] > stt->capacitorFast) {
388 0 : stt->capacitorFast = env[k];
389 : }
390 : // Slow envelope follower
391 0 : if (env[k] > stt->capacitorSlow) {
392 : // increase capacitorSlow
393 0 : stt->capacitorSlow = AGC_SCALEDIFF32(500, (env[k] - stt->capacitorSlow),
394 : stt->capacitorSlow);
395 : } else {
396 : // decrease capacitorSlow
397 0 : stt->capacitorSlow =
398 0 : AGC_SCALEDIFF32(decay, stt->capacitorSlow, stt->capacitorSlow);
399 : }
400 :
401 : // use maximum of both capacitors as current level
402 0 : if (stt->capacitorFast > stt->capacitorSlow) {
403 0 : cur_level = stt->capacitorFast;
404 : } else {
405 0 : cur_level = stt->capacitorSlow;
406 : }
407 : // Translate signal level into gain, using a piecewise linear approximation
408 : // find number of leading zeros
409 0 : zeros = WebRtcSpl_NormU32((uint32_t)cur_level);
410 0 : if (cur_level == 0) {
411 0 : zeros = 31;
412 : }
413 0 : tmp32 = (cur_level << zeros) & 0x7FFFFFFF;
414 0 : frac = (int16_t)(tmp32 >> 19); // Q12.
415 0 : tmp32 = (stt->gainTable[zeros - 1] - stt->gainTable[zeros]) * frac;
416 0 : gains[k + 1] = stt->gainTable[zeros] + (tmp32 >> 12);
417 : #ifdef WEBRTC_AGC_DEBUG_DUMP
418 : if (k == 0) {
419 : fprintf(stt->logFile, "%d\t%d\t%d\t%d\t%d\n", env[0], cur_level,
420 : stt->capacitorFast, stt->capacitorSlow, zeros);
421 : }
422 : #endif
423 : }
424 :
425 : // Gate processing (lower gain during absence of speech)
426 0 : zeros = (zeros << 9) - (frac >> 3);
427 : // find number of leading zeros
428 0 : zeros_fast = WebRtcSpl_NormU32((uint32_t)stt->capacitorFast);
429 0 : if (stt->capacitorFast == 0) {
430 0 : zeros_fast = 31;
431 : }
432 0 : tmp32 = (stt->capacitorFast << zeros_fast) & 0x7FFFFFFF;
433 0 : zeros_fast <<= 9;
434 0 : zeros_fast -= (int16_t)(tmp32 >> 22);
435 :
436 0 : gate = 1000 + zeros_fast - zeros - stt->vadNearend.stdShortTerm;
437 :
438 0 : if (gate < 0) {
439 0 : stt->gatePrevious = 0;
440 : } else {
441 0 : tmp32 = stt->gatePrevious * 7;
442 0 : gate = (int16_t)((gate + tmp32) >> 3);
443 0 : stt->gatePrevious = gate;
444 : }
445 : // gate < 0 -> no gate
446 : // gate > 2500 -> max gate
447 0 : if (gate > 0) {
448 0 : if (gate < 2500) {
449 0 : gain_adj = (2500 - gate) >> 5;
450 : } else {
451 0 : gain_adj = 0;
452 : }
453 0 : for (k = 0; k < 10; k++) {
454 0 : if ((gains[k + 1] - stt->gainTable[0]) > 8388608) {
455 : // To prevent wraparound
456 0 : tmp32 = (gains[k + 1] - stt->gainTable[0]) >> 8;
457 0 : tmp32 *= 178 + gain_adj;
458 : } else {
459 0 : tmp32 = (gains[k + 1] - stt->gainTable[0]) * (178 + gain_adj);
460 0 : tmp32 >>= 8;
461 : }
462 0 : gains[k + 1] = stt->gainTable[0] + tmp32;
463 : }
464 : }
465 :
466 : // Limit gain to avoid overload distortion
467 0 : for (k = 0; k < 10; k++) {
468 : // To prevent wrap around
469 0 : zeros = 10;
470 0 : if (gains[k + 1] > 47453132) {
471 0 : zeros = 16 - WebRtcSpl_NormW32(gains[k + 1]);
472 : }
473 0 : gain32 = (gains[k + 1] >> zeros) + 1;
474 0 : gain32 *= gain32;
475 : // check for overflow
476 0 : while (AGC_MUL32((env[k] >> 12) + 1, gain32) >
477 0 : WEBRTC_SPL_SHIFT_W32((int32_t)32767, 2 * (1 - zeros + 10))) {
478 : // multiply by 253/256 ==> -0.1 dB
479 0 : if (gains[k + 1] > 8388607) {
480 : // Prevent wrap around
481 0 : gains[k + 1] = (gains[k + 1] / 256) * 253;
482 : } else {
483 0 : gains[k + 1] = (gains[k + 1] * 253) / 256;
484 : }
485 0 : gain32 = (gains[k + 1] >> zeros) + 1;
486 0 : gain32 *= gain32;
487 : }
488 : }
489 : // gain reductions should be done 1 ms earlier than gain increases
490 0 : for (k = 1; k < 10; k++) {
491 0 : if (gains[k] > gains[k + 1]) {
492 0 : gains[k] = gains[k + 1];
493 : }
494 : }
495 : // save start gain for next frame
496 0 : stt->gain = gains[10];
497 :
498 : // Apply gain
499 : // handle first sub frame separately
500 0 : delta = (gains[1] - gains[0]) * (1 << (4 - L2));
501 0 : gain32 = gains[0] * (1 << 4);
502 : // iterate over samples
503 0 : for (n = 0; n < L; n++) {
504 0 : for (i = 0; i < num_bands; ++i) {
505 0 : tmp32 = out[i][n] * ((gain32 + 127) >> 7);
506 0 : out_tmp = tmp32 >> 16;
507 0 : if (out_tmp > 4095) {
508 0 : out[i][n] = (int16_t)32767;
509 0 : } else if (out_tmp < -4096) {
510 0 : out[i][n] = (int16_t)-32768;
511 : } else {
512 0 : tmp32 = out[i][n] * (gain32 >> 4);
513 0 : out[i][n] = (int16_t)(tmp32 >> 16);
514 : }
515 : }
516 : //
517 :
518 0 : gain32 += delta;
519 : }
520 : // iterate over subframes
521 0 : for (k = 1; k < 10; k++) {
522 0 : delta = (gains[k + 1] - gains[k]) * (1 << (4 - L2));
523 0 : gain32 = gains[k] * (1 << 4);
524 : // iterate over samples
525 0 : for (n = 0; n < L; n++) {
526 0 : for (i = 0; i < num_bands; ++i) {
527 0 : tmp32 = out[i][k * L + n] * (gain32 >> 4);
528 0 : out[i][k * L + n] = (int16_t)(tmp32 >> 16);
529 : }
530 0 : gain32 += delta;
531 : }
532 : }
533 :
534 0 : return 0;
535 : }
536 :
537 0 : void WebRtcAgc_InitVad(AgcVad* state) {
538 : int16_t k;
539 :
540 0 : state->HPstate = 0; // state of high pass filter
541 0 : state->logRatio = 0; // log( P(active) / P(inactive) )
542 : // average input level (Q10)
543 0 : state->meanLongTerm = 15 << 10;
544 :
545 : // variance of input level (Q8)
546 0 : state->varianceLongTerm = 500 << 8;
547 :
548 0 : state->stdLongTerm = 0; // standard deviation of input level in dB
549 : // short-term average input level (Q10)
550 0 : state->meanShortTerm = 15 << 10;
551 :
552 : // short-term variance of input level (Q8)
553 0 : state->varianceShortTerm = 500 << 8;
554 :
555 0 : state->stdShortTerm =
556 : 0; // short-term standard deviation of input level in dB
557 0 : state->counter = 3; // counts updates
558 0 : for (k = 0; k < 8; k++) {
559 : // downsampling filter
560 0 : state->downState[k] = 0;
561 : }
562 0 : }
563 :
564 0 : int16_t WebRtcAgc_ProcessVad(AgcVad* state, // (i) VAD state
565 : const int16_t* in, // (i) Speech signal
566 : size_t nrSamples) // (i) number of samples
567 : {
568 : int32_t out, nrg, tmp32, tmp32b;
569 : uint16_t tmpU16;
570 : int16_t k, subfr, tmp16;
571 : int16_t buf1[8];
572 : int16_t buf2[4];
573 : int16_t HPstate;
574 : int16_t zeros, dB;
575 :
576 : // process in 10 sub frames of 1 ms (to save on memory)
577 0 : nrg = 0;
578 0 : HPstate = state->HPstate;
579 0 : for (subfr = 0; subfr < 10; subfr++) {
580 : // downsample to 4 kHz
581 0 : if (nrSamples == 160) {
582 0 : for (k = 0; k < 8; k++) {
583 0 : tmp32 = (int32_t)in[2 * k] + (int32_t)in[2 * k + 1];
584 0 : tmp32 >>= 1;
585 0 : buf1[k] = (int16_t)tmp32;
586 : }
587 0 : in += 16;
588 :
589 0 : WebRtcSpl_DownsampleBy2(buf1, 8, buf2, state->downState);
590 : } else {
591 0 : WebRtcSpl_DownsampleBy2(in, 8, buf2, state->downState);
592 0 : in += 8;
593 : }
594 :
595 : // high pass filter and compute energy
596 0 : for (k = 0; k < 4; k++) {
597 0 : out = buf2[k] + HPstate;
598 0 : tmp32 = 600 * out;
599 0 : HPstate = (int16_t)((tmp32 >> 10) - buf2[k]);
600 0 : nrg += (out * out) >> 6;
601 : }
602 : }
603 0 : state->HPstate = HPstate;
604 :
605 : // find number of leading zeros
606 0 : if (!(0xFFFF0000 & nrg)) {
607 0 : zeros = 16;
608 : } else {
609 0 : zeros = 0;
610 : }
611 0 : if (!(0xFF000000 & (nrg << zeros))) {
612 0 : zeros += 8;
613 : }
614 0 : if (!(0xF0000000 & (nrg << zeros))) {
615 0 : zeros += 4;
616 : }
617 0 : if (!(0xC0000000 & (nrg << zeros))) {
618 0 : zeros += 2;
619 : }
620 0 : if (!(0x80000000 & (nrg << zeros))) {
621 0 : zeros += 1;
622 : }
623 :
624 : // energy level (range {-32..30}) (Q10)
625 0 : dB = (15 - zeros) << 11;
626 :
627 : // Update statistics
628 :
629 0 : if (state->counter < kAvgDecayTime) {
630 : // decay time = AvgDecTime * 10 ms
631 0 : state->counter++;
632 : }
633 :
634 : // update short-term estimate of mean energy level (Q10)
635 0 : tmp32 = state->meanShortTerm * 15 + dB;
636 0 : state->meanShortTerm = (int16_t)(tmp32 >> 4);
637 :
638 : // update short-term estimate of variance in energy level (Q8)
639 0 : tmp32 = (dB * dB) >> 12;
640 0 : tmp32 += state->varianceShortTerm * 15;
641 0 : state->varianceShortTerm = tmp32 / 16;
642 :
643 : // update short-term estimate of standard deviation in energy level (Q10)
644 0 : tmp32 = state->meanShortTerm * state->meanShortTerm;
645 0 : tmp32 = (state->varianceShortTerm << 12) - tmp32;
646 0 : state->stdShortTerm = (int16_t)WebRtcSpl_Sqrt(tmp32);
647 :
648 : // update long-term estimate of mean energy level (Q10)
649 0 : tmp32 = state->meanLongTerm * state->counter + dB;
650 0 : state->meanLongTerm =
651 0 : WebRtcSpl_DivW32W16ResW16(tmp32, WebRtcSpl_AddSatW16(state->counter, 1));
652 :
653 : // update long-term estimate of variance in energy level (Q8)
654 0 : tmp32 = (dB * dB) >> 12;
655 0 : tmp32 += state->varianceLongTerm * state->counter;
656 0 : state->varianceLongTerm =
657 0 : WebRtcSpl_DivW32W16(tmp32, WebRtcSpl_AddSatW16(state->counter, 1));
658 :
659 : // update long-term estimate of standard deviation in energy level (Q10)
660 0 : tmp32 = state->meanLongTerm * state->meanLongTerm;
661 0 : tmp32 = (state->varianceLongTerm << 12) - tmp32;
662 0 : state->stdLongTerm = (int16_t)WebRtcSpl_Sqrt(tmp32);
663 :
664 : // update voice activity measure (Q10)
665 0 : tmp16 = 3 << 12;
666 : // TODO(bjornv): (dB - state->meanLongTerm) can overflow, e.g., in
667 : // ApmTest.Process unit test. Previously the macro WEBRTC_SPL_MUL_16_16()
668 : // was used, which did an intermediate cast to (int16_t), hence losing
669 : // significant bits. This cause logRatio to max out positive, rather than
670 : // negative. This is a bug, but has very little significance.
671 0 : tmp32 = tmp16 * (int16_t)(dB - state->meanLongTerm);
672 0 : tmp32 = WebRtcSpl_DivW32W16(tmp32, state->stdLongTerm);
673 0 : tmpU16 = (13 << 12);
674 0 : tmp32b = WEBRTC_SPL_MUL_16_U16(state->logRatio, tmpU16);
675 0 : tmp32 += tmp32b >> 10;
676 :
677 0 : state->logRatio = (int16_t)(tmp32 >> 6);
678 :
679 : // limit
680 0 : if (state->logRatio > 2048) {
681 0 : state->logRatio = 2048;
682 : }
683 0 : if (state->logRatio < -2048) {
684 0 : state->logRatio = -2048;
685 : }
686 :
687 0 : return state->logRatio; // Q10
688 : }
|
