| OLD | NEW |
|---|
| (Empty) |
| 1 /* | |
| 2 * Copyright (c) 2013 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 #include "webrtc/modules/audio_processing/aecm/aecm_core.h" | |
| 12 | |
| 13 #include <assert.h> | |
| 14 #include <stddef.h> | |
| 15 #include <stdlib.h> | |
| 16 | |
| 17 #include "webrtc/common_audio/ring_buffer.h" | |
| 18 #include "webrtc/common_audio/signal_processing/include/real_fft.h" | |
| 19 #include "webrtc/modules/audio_processing/aecm/echo_control_mobile.h" | |
| 20 #include "webrtc/modules/audio_processing/utility/delay_estimator_wrapper.h" | |
| 21 #include "webrtc/system_wrappers/include/compile_assert_c.h" | |
| 22 #include "webrtc/system_wrappers/include/cpu_features_wrapper.h" | |
| 23 #include "webrtc/typedefs.h" | |
| 24 | |
| 25 // Square root of Hanning window in Q14. | |
| 26 #if defined(WEBRTC_DETECT_NEON) || defined(WEBRTC_HAS_NEON) | |
| 27 // Table is defined in an ARM assembly file. | |
| 28 extern const ALIGN8_BEG int16_t WebRtcAecm_kSqrtHanning[] ALIGN8_END; | |
| 29 #else | |
| 30 static const ALIGN8_BEG int16_t WebRtcAecm_kSqrtHanning[] ALIGN8_END = { | |
| 31 0, 399, 798, 1196, 1594, 1990, 2386, 2780, 3172, | |
| 32 3562, 3951, 4337, 4720, 5101, 5478, 5853, 6224, | |
| 33 6591, 6954, 7313, 7668, 8019, 8364, 8705, 9040, | |
| 34 9370, 9695, 10013, 10326, 10633, 10933, 11227, 11514, | |
| 35 11795, 12068, 12335, 12594, 12845, 13089, 13325, 13553, | |
| 36 13773, 13985, 14189, 14384, 14571, 14749, 14918, 15079, | |
| 37 15231, 15373, 15506, 15631, 15746, 15851, 15947, 16034, | |
| 38 16111, 16179, 16237, 16286, 16325, 16354, 16373, 16384 | |
| 39 }; | |
| 40 #endif | |
| 41 | |
| 42 #ifdef AECM_WITH_ABS_APPROX | |
| 43 //Q15 alpha = 0.99439986968132 const Factor for magnitude approximation | |
| 44 static const uint16_t kAlpha1 = 32584; | |
| 45 //Q15 beta = 0.12967166976970 const Factor for magnitude approximation | |
| 46 static const uint16_t kBeta1 = 4249; | |
| 47 //Q15 alpha = 0.94234827210087 const Factor for magnitude approximation | |
| 48 static const uint16_t kAlpha2 = 30879; | |
| 49 //Q15 beta = 0.33787806009150 const Factor for magnitude approximation | |
| 50 static const uint16_t kBeta2 = 11072; | |
| 51 //Q15 alpha = 0.82247698684306 const Factor for magnitude approximation | |
| 52 static const uint16_t kAlpha3 = 26951; | |
| 53 //Q15 beta = 0.57762063060713 const Factor for magnitude approximation | |
| 54 static const uint16_t kBeta3 = 18927; | |
| 55 #endif | |
| 56 | |
| 57 static const int16_t kNoiseEstQDomain = 15; | |
| 58 static const int16_t kNoiseEstIncCount = 5; | |
| 59 | |
| 60 static void ComfortNoise(AecmCore* aecm, | |
| 61 const uint16_t* dfa, | |
| 62 ComplexInt16* out, | |
| 63 const int16_t* lambda); | |
| 64 | |
| 65 static void WindowAndFFT(AecmCore* aecm, | |
| 66 int16_t* fft, | |
| 67 const int16_t* time_signal, | |
| 68 ComplexInt16* freq_signal, | |
| 69 int time_signal_scaling) { | |
| 70 int i = 0; | |
| 71 | |
| 72 // FFT of signal | |
| 73 for (i = 0; i < PART_LEN; i++) { | |
| 74 // Window time domain signal and insert into real part of | |
| 75 // transformation array |fft| | |
| 76 int16_t scaled_time_signal = time_signal[i] << time_signal_scaling; | |
| 77 fft[i] = (int16_t)((scaled_time_signal * WebRtcAecm_kSqrtHanning[i]) >> 14); | |
| 78 scaled_time_signal = time_signal[i + PART_LEN] << time_signal_scaling; | |
| 79 fft[PART_LEN + i] = (int16_t)(( | |
| 80 scaled_time_signal * WebRtcAecm_kSqrtHanning[PART_LEN - i]) >> 14); | |
| 81 } | |
| 82 | |
| 83 // Do forward FFT, then take only the first PART_LEN complex samples, | |
| 84 // and change signs of the imaginary parts. | |
| 85 WebRtcSpl_RealForwardFFT(aecm->real_fft, fft, (int16_t*)freq_signal); | |
| 86 for (i = 0; i < PART_LEN; i++) { | |
| 87 freq_signal[i].imag = -freq_signal[i].imag; | |
| 88 } | |
| 89 } | |
| 90 | |
| 91 static void InverseFFTAndWindow(AecmCore* aecm, | |
| 92 int16_t* fft, | |
| 93 ComplexInt16* efw, | |
| 94 int16_t* output, | |
| 95 const int16_t* nearendClean) { | |
| 96 int i, j, outCFFT; | |
| 97 int32_t tmp32no1; | |
| 98 // Reuse |efw| for the inverse FFT output after transferring | |
| 99 // the contents to |fft|. | |
| 100 int16_t* ifft_out = (int16_t*)efw; | |
| 101 | |
| 102 // Synthesis | |
| 103 for (i = 1, j = 2; i < PART_LEN; i += 1, j += 2) { | |
| 104 fft[j] = efw[i].real; | |
| 105 fft[j + 1] = -efw[i].imag; | |
| 106 } | |
| 107 fft[0] = efw[0].real; | |
| 108 fft[1] = -efw[0].imag; | |
| 109 | |
| 110 fft[PART_LEN2] = efw[PART_LEN].real; | |
| 111 fft[PART_LEN2 + 1] = -efw[PART_LEN].imag; | |
| 112 | |
| 113 // Inverse FFT. Keep outCFFT to scale the samples in the next block. | |
| 114 outCFFT = WebRtcSpl_RealInverseFFT(aecm->real_fft, fft, ifft_out); | |
| 115 for (i = 0; i < PART_LEN; i++) { | |
| 116 ifft_out[i] = (int16_t)WEBRTC_SPL_MUL_16_16_RSFT_WITH_ROUND( | |
| 117 ifft_out[i], WebRtcAecm_kSqrtHanning[i], 14); | |
| 118 tmp32no1 = WEBRTC_SPL_SHIFT_W32((int32_t)ifft_out[i], | |
| 119 outCFFT - aecm->dfaCleanQDomain); | |
| 120 output[i] = (int16_t)WEBRTC_SPL_SAT(WEBRTC_SPL_WORD16_MAX, | |
| 121 tmp32no1 + aecm->outBuf[i], | |
| 122 WEBRTC_SPL_WORD16_MIN); | |
| 123 | |
| 124 tmp32no1 = (ifft_out[PART_LEN + i] * | |
| 125 WebRtcAecm_kSqrtHanning[PART_LEN - i]) >> 14; | |
| 126 tmp32no1 = WEBRTC_SPL_SHIFT_W32(tmp32no1, | |
| 127 outCFFT - aecm->dfaCleanQDomain); | |
| 128 aecm->outBuf[i] = (int16_t)WEBRTC_SPL_SAT(WEBRTC_SPL_WORD16_MAX, | |
| 129 tmp32no1, | |
| 130 WEBRTC_SPL_WORD16_MIN); | |
| 131 } | |
| 132 | |
| 133 // Copy the current block to the old position | |
| 134 // (aecm->outBuf is shifted elsewhere) | |
| 135 memcpy(aecm->xBuf, aecm->xBuf + PART_LEN, sizeof(int16_t) * PART_LEN); | |
| 136 memcpy(aecm->dBufNoisy, | |
| 137 aecm->dBufNoisy + PART_LEN, | |
| 138 sizeof(int16_t) * PART_LEN); | |
| 139 if (nearendClean != NULL) | |
| 140 { | |
| 141 memcpy(aecm->dBufClean, | |
| 142 aecm->dBufClean + PART_LEN, | |
| 143 sizeof(int16_t) * PART_LEN); | |
| 144 } | |
| 145 } | |
| 146 | |
| 147 // Transforms a time domain signal into the frequency domain, outputting the | |
| 148 // complex valued signal, absolute value and sum of absolute values. | |
| 149 // | |
| 150 // time_signal [in] Pointer to time domain signal | |
| 151 // freq_signal_real [out] Pointer to real part of frequency domain array | |
| 152 // freq_signal_imag [out] Pointer to imaginary part of frequency domain | |
| 153 // array | |
| 154 // freq_signal_abs [out] Pointer to absolute value of frequency domain | |
| 155 // array | |
| 156 // freq_signal_sum_abs [out] Pointer to the sum of all absolute values in | |
| 157 // the frequency domain array | |
| 158 // return value The Q-domain of current frequency values | |
| 159 // | |
| 160 static int TimeToFrequencyDomain(AecmCore* aecm, | |
| 161 const int16_t* time_signal, | |
| 162 ComplexInt16* freq_signal, | |
| 163 uint16_t* freq_signal_abs, | |
| 164 uint32_t* freq_signal_sum_abs) { | |
| 165 int i = 0; | |
| 166 int time_signal_scaling = 0; | |
| 167 | |
| 168 int32_t tmp32no1 = 0; | |
| 169 int32_t tmp32no2 = 0; | |
| 170 | |
| 171 // In fft_buf, +16 for 32-byte alignment. | |
| 172 int16_t fft_buf[PART_LEN4 + 16]; | |
| 173 int16_t *fft = (int16_t *) (((uintptr_t) fft_buf + 31) & ~31); | |
| 174 | |
| 175 int16_t tmp16no1; | |
| 176 #ifndef WEBRTC_ARCH_ARM_V7 | |
| 177 int16_t tmp16no2; | |
| 178 #endif | |
| 179 #ifdef AECM_WITH_ABS_APPROX | |
| 180 int16_t max_value = 0; | |
| 181 int16_t min_value = 0; | |
| 182 uint16_t alpha = 0; | |
| 183 uint16_t beta = 0; | |
| 184 #endif | |
| 185 | |
| 186 #ifdef AECM_DYNAMIC_Q | |
| 187 tmp16no1 = WebRtcSpl_MaxAbsValueW16(time_signal, PART_LEN2); | |
| 188 time_signal_scaling = WebRtcSpl_NormW16(tmp16no1); | |
| 189 #endif | |
| 190 | |
| 191 WindowAndFFT(aecm, fft, time_signal, freq_signal, time_signal_scaling); | |
| 192 | |
| 193 // Extract imaginary and real part, calculate the magnitude for | |
| 194 // all frequency bins | |
| 195 freq_signal[0].imag = 0; | |
| 196 freq_signal[PART_LEN].imag = 0; | |
| 197 freq_signal_abs[0] = (uint16_t)WEBRTC_SPL_ABS_W16(freq_signal[0].real); | |
| 198 freq_signal_abs[PART_LEN] = (uint16_t)WEBRTC_SPL_ABS_W16( | |
| 199 freq_signal[PART_LEN].real); | |
| 200 (*freq_signal_sum_abs) = (uint32_t)(freq_signal_abs[0]) + | |
| 201 (uint32_t)(freq_signal_abs[PART_LEN]); | |
| 202 | |
| 203 for (i = 1; i < PART_LEN; i++) | |
| 204 { | |
| 205 if (freq_signal[i].real == 0) | |
| 206 { | |
| 207 freq_signal_abs[i] = (uint16_t)WEBRTC_SPL_ABS_W16(freq_signal[i].imag); | |
| 208 } | |
| 209 else if (freq_signal[i].imag == 0) | |
| 210 { | |
| 211 freq_signal_abs[i] = (uint16_t)WEBRTC_SPL_ABS_W16(freq_signal[i].real); | |
| 212 } | |
| 213 else | |
| 214 { | |
| 215 // Approximation for magnitude of complex fft output | |
| 216 // magn = sqrt(real^2 + imag^2) | |
| 217 // magn ~= alpha * max(|imag|,|real|) + beta * min(|imag|,|real|) | |
| 218 // | |
| 219 // The parameters alpha and beta are stored in Q15 | |
| 220 | |
| 221 #ifdef AECM_WITH_ABS_APPROX | |
| 222 tmp16no1 = WEBRTC_SPL_ABS_W16(freq_signal[i].real); | |
| 223 tmp16no2 = WEBRTC_SPL_ABS_W16(freq_signal[i].imag); | |
| 224 | |
| 225 if(tmp16no1 > tmp16no2) | |
| 226 { | |
| 227 max_value = tmp16no1; | |
| 228 min_value = tmp16no2; | |
| 229 } else | |
| 230 { | |
| 231 max_value = tmp16no2; | |
| 232 min_value = tmp16no1; | |
| 233 } | |
| 234 | |
| 235 // Magnitude in Q(-6) | |
| 236 if ((max_value >> 2) > min_value) | |
| 237 { | |
| 238 alpha = kAlpha1; | |
| 239 beta = kBeta1; | |
| 240 } else if ((max_value >> 1) > min_value) | |
| 241 { | |
| 242 alpha = kAlpha2; | |
| 243 beta = kBeta2; | |
| 244 } else | |
| 245 { | |
| 246 alpha = kAlpha3; | |
| 247 beta = kBeta3; | |
| 248 } | |
| 249 tmp16no1 = (int16_t)((max_value * alpha) >> 15); | |
| 250 tmp16no2 = (int16_t)((min_value * beta) >> 15); | |
| 251 freq_signal_abs[i] = (uint16_t)tmp16no1 + (uint16_t)tmp16no2; | |
| 252 #else | |
| 253 #ifdef WEBRTC_ARCH_ARM_V7 | |
| 254 __asm __volatile( | |
| 255 "smulbb %[tmp32no1], %[real], %[real]\n\t" | |
| 256 "smlabb %[tmp32no2], %[imag], %[imag], %[tmp32no1]\n\t" | |
| 257 :[tmp32no1]"+&r"(tmp32no1), | |
| 258 [tmp32no2]"=r"(tmp32no2) | |
| 259 :[real]"r"(freq_signal[i].real), | |
| 260 [imag]"r"(freq_signal[i].imag) | |
| 261 ); | |
| 262 #else | |
| 263 tmp16no1 = WEBRTC_SPL_ABS_W16(freq_signal[i].real); | |
| 264 tmp16no2 = WEBRTC_SPL_ABS_W16(freq_signal[i].imag); | |
| 265 tmp32no1 = tmp16no1 * tmp16no1; | |
| 266 tmp32no2 = tmp16no2 * tmp16no2; | |
| 267 tmp32no2 = WebRtcSpl_AddSatW32(tmp32no1, tmp32no2); | |
| 268 #endif // WEBRTC_ARCH_ARM_V7 | |
| 269 tmp32no1 = WebRtcSpl_SqrtFloor(tmp32no2); | |
| 270 | |
| 271 freq_signal_abs[i] = (uint16_t)tmp32no1; | |
| 272 #endif // AECM_WITH_ABS_APPROX | |
| 273 } | |
| 274 (*freq_signal_sum_abs) += (uint32_t)freq_signal_abs[i]; | |
| 275 } | |
| 276 | |
| 277 return time_signal_scaling; | |
| 278 } | |
| 279 | |
| 280 int WebRtcAecm_ProcessBlock(AecmCore* aecm, | |
| 281 const int16_t* farend, | |
| 282 const int16_t* nearendNoisy, | |
| 283 const int16_t* nearendClean, | |
| 284 int16_t* output) { | |
| 285 int i; | |
| 286 | |
| 287 uint32_t xfaSum; | |
| 288 uint32_t dfaNoisySum; | |
| 289 uint32_t dfaCleanSum; | |
| 290 uint32_t echoEst32Gained; | |
| 291 uint32_t tmpU32; | |
| 292 | |
| 293 int32_t tmp32no1; | |
| 294 | |
| 295 uint16_t xfa[PART_LEN1]; | |
| 296 uint16_t dfaNoisy[PART_LEN1]; | |
| 297 uint16_t dfaClean[PART_LEN1]; | |
| 298 uint16_t* ptrDfaClean = dfaClean; | |
| 299 const uint16_t* far_spectrum_ptr = NULL; | |
| 300 | |
| 301 // 32 byte aligned buffers (with +8 or +16). | |
| 302 // TODO(kma): define fft with ComplexInt16. | |
| 303 int16_t fft_buf[PART_LEN4 + 2 + 16]; // +2 to make a loop safe. | |
| 304 int32_t echoEst32_buf[PART_LEN1 + 8]; | |
| 305 int32_t dfw_buf[PART_LEN2 + 8]; | |
| 306 int32_t efw_buf[PART_LEN2 + 8]; | |
| 307 | |
| 308 int16_t* fft = (int16_t*) (((uintptr_t) fft_buf + 31) & ~ 31); | |
| 309 int32_t* echoEst32 = (int32_t*) (((uintptr_t) echoEst32_buf + 31) & ~ 31); | |
| 310 ComplexInt16* dfw = (ComplexInt16*)(((uintptr_t)dfw_buf + 31) & ~31); | |
| 311 ComplexInt16* efw = (ComplexInt16*)(((uintptr_t)efw_buf + 31) & ~31); | |
| 312 | |
| 313 int16_t hnl[PART_LEN1]; | |
| 314 int16_t numPosCoef = 0; | |
| 315 int16_t nlpGain = ONE_Q14; | |
| 316 int delay; | |
| 317 int16_t tmp16no1; | |
| 318 int16_t tmp16no2; | |
| 319 int16_t mu; | |
| 320 int16_t supGain; | |
| 321 int16_t zeros32, zeros16; | |
| 322 int16_t zerosDBufNoisy, zerosDBufClean, zerosXBuf; | |
| 323 int far_q; | |
| 324 int16_t resolutionDiff, qDomainDiff, dfa_clean_q_domain_diff; | |
| 325 | |
| 326 const int kMinPrefBand = 4; | |
| 327 const int kMaxPrefBand = 24; | |
| 328 int32_t avgHnl32 = 0; | |
| 329 | |
| 330 // Determine startup state. There are three states: | |
| 331 // (0) the first CONV_LEN blocks | |
| 332 // (1) another CONV_LEN blocks | |
| 333 // (2) the rest | |
| 334 | |
| 335 if (aecm->startupState < 2) | |
| 336 { | |
| 337 aecm->startupState = (aecm->totCount >= CONV_LEN) + | |
| 338 (aecm->totCount >= CONV_LEN2); | |
| 339 } | |
| 340 // END: Determine startup state | |
| 341 | |
| 342 // Buffer near and far end signals | |
| 343 memcpy(aecm->xBuf + PART_LEN, farend, sizeof(int16_t) * PART_LEN); | |
| 344 memcpy(aecm->dBufNoisy + PART_LEN, nearendNoisy, sizeof(int16_t) * PART_LEN); | |
| 345 if (nearendClean != NULL) | |
| 346 { | |
| 347 memcpy(aecm->dBufClean + PART_LEN, | |
| 348 nearendClean, | |
| 349 sizeof(int16_t) * PART_LEN); | |
| 350 } | |
| 351 | |
| 352 // Transform far end signal from time domain to frequency domain. | |
| 353 far_q = TimeToFrequencyDomain(aecm, | |
| 354 aecm->xBuf, | |
| 355 dfw, | |
| 356 xfa, | |
| 357 &xfaSum); | |
| 358 | |
| 359 // Transform noisy near end signal from time domain to frequency domain. | |
| 360 zerosDBufNoisy = TimeToFrequencyDomain(aecm, | |
| 361 aecm->dBufNoisy, | |
| 362 dfw, | |
| 363 dfaNoisy, | |
| 364 &dfaNoisySum); | |
| 365 aecm->dfaNoisyQDomainOld = aecm->dfaNoisyQDomain; | |
| 366 aecm->dfaNoisyQDomain = (int16_t)zerosDBufNoisy; | |
| 367 | |
| 368 | |
| 369 if (nearendClean == NULL) | |
| 370 { | |
| 371 ptrDfaClean = dfaNoisy; | |
| 372 aecm->dfaCleanQDomainOld = aecm->dfaNoisyQDomainOld; | |
| 373 aecm->dfaCleanQDomain = aecm->dfaNoisyQDomain; | |
| 374 dfaCleanSum = dfaNoisySum; | |
| 375 } else | |
| 376 { | |
| 377 // Transform clean near end signal from time domain to frequency domain. | |
| 378 zerosDBufClean = TimeToFrequencyDomain(aecm, | |
| 379 aecm->dBufClean, | |
| 380 dfw, | |
| 381 dfaClean, | |
| 382 &dfaCleanSum); | |
| 383 aecm->dfaCleanQDomainOld = aecm->dfaCleanQDomain; | |
| 384 aecm->dfaCleanQDomain = (int16_t)zerosDBufClean; | |
| 385 } | |
| 386 | |
| 387 // Get the delay | |
| 388 // Save far-end history and estimate delay | |
| 389 WebRtcAecm_UpdateFarHistory(aecm, xfa, far_q); | |
| 390 if (WebRtc_AddFarSpectrumFix(aecm->delay_estimator_farend, | |
| 391 xfa, | |
| 392 PART_LEN1, | |
| 393 far_q) == -1) { | |
| 394 return -1; | |
| 395 } | |
| 396 delay = WebRtc_DelayEstimatorProcessFix(aecm->delay_estimator, | |
| 397 dfaNoisy, | |
| 398 PART_LEN1, | |
| 399 zerosDBufNoisy); | |
| 400 if (delay == -1) | |
| 401 { | |
| 402 return -1; | |
| 403 } | |
| 404 else if (delay == -2) | |
| 405 { | |
| 406 // If the delay is unknown, we assume zero. | |
| 407 // NOTE: this will have to be adjusted if we ever add lookahead. | |
| 408 delay = 0; | |
| 409 } | |
| 410 | |
| 411 if (aecm->fixedDelay >= 0) | |
| 412 { | |
| 413 // Use fixed delay | |
| 414 delay = aecm->fixedDelay; | |
| 415 } | |
| 416 | |
| 417 // Get aligned far end spectrum | |
| 418 far_spectrum_ptr = WebRtcAecm_AlignedFarend(aecm, &far_q, delay); | |
| 419 zerosXBuf = (int16_t) far_q; | |
| 420 if (far_spectrum_ptr == NULL) | |
| 421 { | |
| 422 return -1; | |
| 423 } | |
| 424 | |
| 425 // Calculate log(energy) and update energy threshold levels | |
| 426 WebRtcAecm_CalcEnergies(aecm, | |
| 427 far_spectrum_ptr, | |
| 428 zerosXBuf, | |
| 429 dfaNoisySum, | |
| 430 echoEst32); | |
| 431 | |
| 432 // Calculate stepsize | |
| 433 mu = WebRtcAecm_CalcStepSize(aecm); | |
| 434 | |
| 435 // Update counters | |
| 436 aecm->totCount++; | |
| 437 | |
| 438 // This is the channel estimation algorithm. | |
| 439 // It is base on NLMS but has a variable step length, | |
| 440 // which was calculated above. | |
| 441 WebRtcAecm_UpdateChannel(aecm, | |
| 442 far_spectrum_ptr, | |
| 443 zerosXBuf, | |
| 444 dfaNoisy, | |
| 445 mu, | |
| 446 echoEst32); | |
| 447 supGain = WebRtcAecm_CalcSuppressionGain(aecm); | |
| 448 | |
| 449 | |
| 450 // Calculate Wiener filter hnl[] | |
| 451 for (i = 0; i < PART_LEN1; i++) | |
| 452 { | |
| 453 // Far end signal through channel estimate in Q8 | |
| 454 // How much can we shift right to preserve resolution | |
| 455 tmp32no1 = echoEst32[i] - aecm->echoFilt[i]; | |
| 456 aecm->echoFilt[i] += (tmp32no1 * 50) >> 8; | |
| 457 | |
| 458 zeros32 = WebRtcSpl_NormW32(aecm->echoFilt[i]) + 1; | |
| 459 zeros16 = WebRtcSpl_NormW16(supGain) + 1; | |
| 460 if (zeros32 + zeros16 > 16) | |
| 461 { | |
| 462 // Multiplication is safe | |
| 463 // Result in | |
| 464 // Q(RESOLUTION_CHANNEL+RESOLUTION_SUPGAIN+ | |
| 465 // aecm->xfaQDomainBuf[diff]) | |
| 466 echoEst32Gained = WEBRTC_SPL_UMUL_32_16((uint32_t)aecm->echoFilt[i], | |
| 467 (uint16_t)supGain); | |
| 468 resolutionDiff = 14 - RESOLUTION_CHANNEL16 - RESOLUTION_SUPGAIN; | |
| 469 resolutionDiff += (aecm->dfaCleanQDomain - zerosXBuf); | |
| 470 } else | |
| 471 { | |
| 472 tmp16no1 = 17 - zeros32 - zeros16; | |
| 473 resolutionDiff = 14 + tmp16no1 - RESOLUTION_CHANNEL16 - | |
| 474 RESOLUTION_SUPGAIN; | |
| 475 resolutionDiff += (aecm->dfaCleanQDomain - zerosXBuf); | |
| 476 if (zeros32 > tmp16no1) | |
| 477 { | |
| 478 echoEst32Gained = WEBRTC_SPL_UMUL_32_16((uint32_t)aecm->echoFilt[i], | |
| 479 supGain >> tmp16no1); | |
| 480 } else | |
| 481 { | |
| 482 // Result in Q-(RESOLUTION_CHANNEL+RESOLUTION_SUPGAIN-16) | |
| 483 echoEst32Gained = (aecm->echoFilt[i] >> tmp16no1) * supGain; | |
| 484 } | |
| 485 } | |
| 486 | |
| 487 zeros16 = WebRtcSpl_NormW16(aecm->nearFilt[i]); | |
| 488 assert(zeros16 >= 0); // |zeros16| is a norm, hence non-negative. | |
| 489 dfa_clean_q_domain_diff = aecm->dfaCleanQDomain - aecm->dfaCleanQDomainOld; | |
| 490 if (zeros16 < dfa_clean_q_domain_diff && aecm->nearFilt[i]) { | |
| 491 tmp16no1 = aecm->nearFilt[i] << zeros16; | |
| 492 qDomainDiff = zeros16 - dfa_clean_q_domain_diff; | |
| 493 tmp16no2 = ptrDfaClean[i] >> -qDomainDiff; | |
| 494 } else { | |
| 495 tmp16no1 = dfa_clean_q_domain_diff < 0 | |
| 496 ? aecm->nearFilt[i] >> -dfa_clean_q_domain_diff | |
| 497 : aecm->nearFilt[i] << dfa_clean_q_domain_diff; | |
| 498 qDomainDiff = 0; | |
| 499 tmp16no2 = ptrDfaClean[i]; | |
| 500 } | |
| 501 tmp32no1 = (int32_t)(tmp16no2 - tmp16no1); | |
| 502 tmp16no2 = (int16_t)(tmp32no1 >> 4); | |
| 503 tmp16no2 += tmp16no1; | |
| 504 zeros16 = WebRtcSpl_NormW16(tmp16no2); | |
| 505 if ((tmp16no2) & (-qDomainDiff > zeros16)) { | |
| 506 aecm->nearFilt[i] = WEBRTC_SPL_WORD16_MAX; | |
| 507 } else { | |
| 508 aecm->nearFilt[i] = qDomainDiff < 0 ? tmp16no2 << -qDomainDiff | |
| 509 : tmp16no2 >> qDomainDiff; | |
| 510 } | |
| 511 | |
| 512 // Wiener filter coefficients, resulting hnl in Q14 | |
| 513 if (echoEst32Gained == 0) | |
| 514 { | |
| 515 hnl[i] = ONE_Q14; | |
| 516 } else if (aecm->nearFilt[i] == 0) | |
| 517 { | |
| 518 hnl[i] = 0; | |
| 519 } else | |
| 520 { | |
| 521 // Multiply the suppression gain | |
| 522 // Rounding | |
| 523 echoEst32Gained += (uint32_t)(aecm->nearFilt[i] >> 1); | |
| 524 tmpU32 = WebRtcSpl_DivU32U16(echoEst32Gained, | |
| 525 (uint16_t)aecm->nearFilt[i]); | |
| 526 | |
| 527 // Current resolution is | |
| 528 // Q-(RESOLUTION_CHANNEL+RESOLUTION_SUPGAIN- max(0,17-zeros16- zeros32)) | |
| 529 // Make sure we are in Q14 | |
| 530 tmp32no1 = (int32_t)WEBRTC_SPL_SHIFT_W32(tmpU32, resolutionDiff); | |
| 531 if (tmp32no1 > ONE_Q14) | |
| 532 { | |
| 533 hnl[i] = 0; | |
| 534 } else if (tmp32no1 < 0) | |
| 535 { | |
| 536 hnl[i] = ONE_Q14; | |
| 537 } else | |
| 538 { | |
| 539 // 1-echoEst/dfa | |
| 540 hnl[i] = ONE_Q14 - (int16_t)tmp32no1; | |
| 541 if (hnl[i] < 0) | |
| 542 { | |
| 543 hnl[i] = 0; | |
| 544 } | |
| 545 } | |
| 546 } | |
| 547 if (hnl[i]) | |
| 548 { | |
| 549 numPosCoef++; | |
| 550 } | |
| 551 } | |
| 552 // Only in wideband. Prevent the gain in upper band from being larger than | |
| 553 // in lower band. | |
| 554 if (aecm->mult == 2) | |
| 555 { | |
| 556 // TODO(bjornv): Investigate if the scaling of hnl[i] below can cause | |
| 557 // speech distortion in double-talk. | |
| 558 for (i = 0; i < PART_LEN1; i++) | |
| 559 { | |
| 560 hnl[i] = (int16_t)((hnl[i] * hnl[i]) >> 14); | |
| 561 } | |
| 562 | |
| 563 for (i = kMinPrefBand; i <= kMaxPrefBand; i++) | |
| 564 { | |
| 565 avgHnl32 += (int32_t)hnl[i]; | |
| 566 } | |
| 567 assert(kMaxPrefBand - kMinPrefBand + 1 > 0); | |
| 568 avgHnl32 /= (kMaxPrefBand - kMinPrefBand + 1); | |
| 569 | |
| 570 for (i = kMaxPrefBand; i < PART_LEN1; i++) | |
| 571 { | |
| 572 if (hnl[i] > (int16_t)avgHnl32) | |
| 573 { | |
| 574 hnl[i] = (int16_t)avgHnl32; | |
| 575 } | |
| 576 } | |
| 577 } | |
| 578 | |
| 579 // Calculate NLP gain, result is in Q14 | |
| 580 if (aecm->nlpFlag) | |
| 581 { | |
| 582 for (i = 0; i < PART_LEN1; i++) | |
| 583 { | |
| 584 // Truncate values close to zero and one. | |
| 585 if (hnl[i] > NLP_COMP_HIGH) | |
| 586 { | |
| 587 hnl[i] = ONE_Q14; | |
| 588 } else if (hnl[i] < NLP_COMP_LOW) | |
| 589 { | |
| 590 hnl[i] = 0; | |
| 591 } | |
| 592 | |
| 593 // Remove outliers | |
| 594 if (numPosCoef < 3) | |
| 595 { | |
| 596 nlpGain = 0; | |
| 597 } else | |
| 598 { | |
| 599 nlpGain = ONE_Q14; | |
| 600 } | |
| 601 | |
| 602 // NLP | |
| 603 if ((hnl[i] == ONE_Q14) && (nlpGain == ONE_Q14)) | |
| 604 { | |
| 605 hnl[i] = ONE_Q14; | |
| 606 } else | |
| 607 { | |
| 608 hnl[i] = (int16_t)((hnl[i] * nlpGain) >> 14); | |
| 609 } | |
| 610 | |
| 611 // multiply with Wiener coefficients | |
| 612 efw[i].real = (int16_t)(WEBRTC_SPL_MUL_16_16_RSFT_WITH_ROUND(dfw[i].real, | |
| 613 hnl[i], 14)); | |
| 614 efw[i].imag = (int16_t)(WEBRTC_SPL_MUL_16_16_RSFT_WITH_ROUND(dfw[i].imag, | |
| 615 hnl[i], 14)); | |
| 616 } | |
| 617 } | |
| 618 else | |
| 619 { | |
| 620 // multiply with Wiener coefficients | |
| 621 for (i = 0; i < PART_LEN1; i++) | |
| 622 { | |
| 623 efw[i].real = (int16_t)(WEBRTC_SPL_MUL_16_16_RSFT_WITH_ROUND(dfw[i].real, | |
| 624 hnl[i], 14)); | |
| 625 efw[i].imag = (int16_t)(WEBRTC_SPL_MUL_16_16_RSFT_WITH_ROUND(dfw[i].imag, | |
| 626 hnl[i], 14)); | |
| 627 } | |
| 628 } | |
| 629 | |
| 630 if (aecm->cngMode == AecmTrue) | |
| 631 { | |
| 632 ComfortNoise(aecm, ptrDfaClean, efw, hnl); | |
| 633 } | |
| 634 | |
| 635 InverseFFTAndWindow(aecm, fft, efw, output, nearendClean); | |
| 636 | |
| 637 return 0; | |
| 638 } | |
| 639 | |
| 640 static void ComfortNoise(AecmCore* aecm, | |
| 641 const uint16_t* dfa, | |
| 642 ComplexInt16* out, | |
| 643 const int16_t* lambda) { | |
| 644 int16_t i; | |
| 645 int16_t tmp16; | |
| 646 int32_t tmp32; | |
| 647 | |
| 648 int16_t randW16[PART_LEN]; | |
| 649 int16_t uReal[PART_LEN1]; | |
| 650 int16_t uImag[PART_LEN1]; | |
| 651 int32_t outLShift32; | |
| 652 int16_t noiseRShift16[PART_LEN1]; | |
| 653 | |
| 654 int16_t shiftFromNearToNoise = kNoiseEstQDomain - aecm->dfaCleanQDomain; | |
| 655 int16_t minTrackShift; | |
| 656 | |
| 657 assert(shiftFromNearToNoise >= 0); | |
| 658 assert(shiftFromNearToNoise < 16); | |
| 659 | |
| 660 if (aecm->noiseEstCtr < 100) | |
| 661 { | |
| 662 // Track the minimum more quickly initially. | |
| 663 aecm->noiseEstCtr++; | |
| 664 minTrackShift = 6; | |
| 665 } else | |
| 666 { | |
| 667 minTrackShift = 9; | |
| 668 } | |
| 669 | |
| 670 // Estimate noise power. | |
| 671 for (i = 0; i < PART_LEN1; i++) | |
| 672 { | |
| 673 // Shift to the noise domain. | |
| 674 tmp32 = (int32_t)dfa[i]; | |
| 675 outLShift32 = tmp32 << shiftFromNearToNoise; | |
| 676 | |
| 677 if (outLShift32 < aecm->noiseEst[i]) | |
| 678 { | |
| 679 // Reset "too low" counter | |
| 680 aecm->noiseEstTooLowCtr[i] = 0; | |
| 681 // Track the minimum. | |
| 682 if (aecm->noiseEst[i] < (1 << minTrackShift)) | |
| 683 { | |
| 684 // For small values, decrease noiseEst[i] every | |
| 685 // |kNoiseEstIncCount| block. The regular approach below can not | |
| 686 // go further down due to truncation. | |
| 687 aecm->noiseEstTooHighCtr[i]++; | |
| 688 if (aecm->noiseEstTooHighCtr[i] >= kNoiseEstIncCount) | |
| 689 { | |
| 690 aecm->noiseEst[i]--; | |
| 691 aecm->noiseEstTooHighCtr[i] = 0; // Reset the counter | |
| 692 } | |
| 693 } | |
| 694 else | |
| 695 { | |
| 696 aecm->noiseEst[i] -= ((aecm->noiseEst[i] - outLShift32) | |
| 697 >> minTrackShift); | |
| 698 } | |
| 699 } else | |
| 700 { | |
| 701 // Reset "too high" counter | |
| 702 aecm->noiseEstTooHighCtr[i] = 0; | |
| 703 // Ramp slowly upwards until we hit the minimum again. | |
| 704 if ((aecm->noiseEst[i] >> 19) > 0) | |
| 705 { | |
| 706 // Avoid overflow. | |
| 707 // Multiplication with 2049 will cause wrap around. Scale | |
| 708 // down first and then multiply | |
| 709 aecm->noiseEst[i] >>= 11; | |
| 710 aecm->noiseEst[i] *= 2049; | |
| 711 } | |
| 712 else if ((aecm->noiseEst[i] >> 11) > 0) | |
| 713 { | |
| 714 // Large enough for relative increase | |
| 715 aecm->noiseEst[i] *= 2049; | |
| 716 aecm->noiseEst[i] >>= 11; | |
| 717 } | |
| 718 else | |
| 719 { | |
| 720 // Make incremental increases based on size every | |
| 721 // |kNoiseEstIncCount| block | |
| 722 aecm->noiseEstTooLowCtr[i]++; | |
| 723 if (aecm->noiseEstTooLowCtr[i] >= kNoiseEstIncCount) | |
| 724 { | |
| 725 aecm->noiseEst[i] += (aecm->noiseEst[i] >> 9) + 1; | |
| 726 aecm->noiseEstTooLowCtr[i] = 0; // Reset counter | |
| 727 } | |
| 728 } | |
| 729 } | |
| 730 } | |
| 731 | |
| 732 for (i = 0; i < PART_LEN1; i++) | |
| 733 { | |
| 734 tmp32 = aecm->noiseEst[i] >> shiftFromNearToNoise; | |
| 735 if (tmp32 > 32767) | |
| 736 { | |
| 737 tmp32 = 32767; | |
| 738 aecm->noiseEst[i] = tmp32 << shiftFromNearToNoise; | |
| 739 } | |
| 740 noiseRShift16[i] = (int16_t)tmp32; | |
| 741 | |
| 742 tmp16 = ONE_Q14 - lambda[i]; | |
| 743 noiseRShift16[i] = (int16_t)((tmp16 * noiseRShift16[i]) >> 14); | |
| 744 } | |
| 745 | |
| 746 // Generate a uniform random array on [0 2^15-1]. | |
| 747 WebRtcSpl_RandUArray(randW16, PART_LEN, &aecm->seed); | |
| 748 | |
| 749 // Generate noise according to estimated energy. | |
| 750 uReal[0] = 0; // Reject LF noise. | |
| 751 uImag[0] = 0; | |
| 752 for (i = 1; i < PART_LEN1; i++) | |
| 753 { | |
| 754 // Get a random index for the cos and sin tables over [0 359]. | |
| 755 tmp16 = (int16_t)((359 * randW16[i - 1]) >> 15); | |
| 756 | |
| 757 // Tables are in Q13. | |
| 758 uReal[i] = (int16_t)((noiseRShift16[i] * WebRtcAecm_kCosTable[tmp16]) >> | |
| 759 13); | |
| 760 uImag[i] = (int16_t)((-noiseRShift16[i] * WebRtcAecm_kSinTable[tmp16]) >> | |
| 761 13); | |
| 762 } | |
| 763 uImag[PART_LEN] = 0; | |
| 764 | |
| 765 for (i = 0; i < PART_LEN1; i++) | |
| 766 { | |
| 767 out[i].real = WebRtcSpl_AddSatW16(out[i].real, uReal[i]); | |
| 768 out[i].imag = WebRtcSpl_AddSatW16(out[i].imag, uImag[i]); | |
| 769 } | |
| 770 } | |
| 771 | |
| OLD | NEW |
|---|
Chromium Code Reviews
Issue 1857153002:
Changed AECM to be built using C++ (Closed)
Base URL: https://chromium.googlesource.com/external/webrtc.git@master