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