blob: 905274fa97bb0da56f45d4c949dd18f0bd992d56 [file] [log] [blame]
andrew@webrtc.orge03cafa2013-11-11 20:10:011/*
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
Mirko Bonadei92ea95e2017-09-15 04:47:3111#include "modules/audio_processing/aecm/aecm_core.h"
andrew@webrtc.orge03cafa2013-11-11 20:10:0112
andrew@webrtc.orge03cafa2013-11-11 20:10:0113#include <stddef.h>
14#include <stdlib.h>
15
peah27045122016-04-11 05:38:1416extern "C" {
Mirko Bonadei92ea95e2017-09-15 04:47:3117#include "common_audio/ring_buffer.h"
18#include "common_audio/signal_processing/include/real_fft.h"
peah27045122016-04-11 05:38:1419}
Mirko Bonadei92ea95e2017-09-15 04:47:3120#include "modules/audio_processing/aecm/echo_control_mobile.h"
21#include "modules/audio_processing/utility/delay_estimator_wrapper.h"
peahbdb7af62016-04-12 21:47:4022extern "C" {
Mirko Bonadei92ea95e2017-09-15 04:47:3123#include "system_wrappers/include/cpu_features_wrapper.h"
peah27045122016-04-11 05:38:1424}
kwiberg9e2be5f2016-09-14 12:23:2225
Mirko Bonadei92ea95e2017-09-15 04:47:3126#include "rtc_base/checks.h"
Karl Wiberge40468b2017-11-22 09:42:2627#include "rtc_base/numerics/safe_conversions.h"
Mirko Bonadei92ea95e2017-09-15 04:47:3128#include "rtc_base/sanitizer.h"
andrew@webrtc.orge03cafa2013-11-11 20:10:0129
30// Square root of Hanning window in Q14.
andrew@webrtc.orge03cafa2013-11-11 20:10:0131static const ALIGN8_BEG int16_t WebRtcAecm_kSqrtHanning[] ALIGN8_END = {
Yves Gerey665174f2018-06-19 13:03:0532 0, 399, 798, 1196, 1594, 1990, 2386, 2780, 3172, 3562, 3951,
33 4337, 4720, 5101, 5478, 5853, 6224, 6591, 6954, 7313, 7668, 8019,
34 8364, 8705, 9040, 9370, 9695, 10013, 10326, 10633, 10933, 11227, 11514,
35 11795, 12068, 12335, 12594, 12845, 13089, 13325, 13553, 13773, 13985, 14189,
36 14384, 14571, 14749, 14918, 15079, 15231, 15373, 15506, 15631, 15746, 15851,
37 15947, 16034, 16111, 16179, 16237, 16286, 16325, 16354, 16373, 16384};
andrew@webrtc.orge03cafa2013-11-11 20:10:0138
39#ifdef AECM_WITH_ABS_APPROX
Yves Gerey665174f2018-06-19 13:03:0540// Q15 alpha = 0.99439986968132 const Factor for magnitude approximation
andrew@webrtc.orge03cafa2013-11-11 20:10:0141static const uint16_t kAlpha1 = 32584;
Yves Gerey665174f2018-06-19 13:03:0542// Q15 beta = 0.12967166976970 const Factor for magnitude approximation
andrew@webrtc.orge03cafa2013-11-11 20:10:0143static const uint16_t kBeta1 = 4249;
Yves Gerey665174f2018-06-19 13:03:0544// Q15 alpha = 0.94234827210087 const Factor for magnitude approximation
andrew@webrtc.orge03cafa2013-11-11 20:10:0145static const uint16_t kAlpha2 = 30879;
Yves Gerey665174f2018-06-19 13:03:0546// Q15 beta = 0.33787806009150 const Factor for magnitude approximation
andrew@webrtc.orge03cafa2013-11-11 20:10:0147static const uint16_t kBeta2 = 11072;
Yves Gerey665174f2018-06-19 13:03:0548// Q15 alpha = 0.82247698684306 const Factor for magnitude approximation
andrew@webrtc.orge03cafa2013-11-11 20:10:0149static const uint16_t kAlpha3 = 26951;
Yves Gerey665174f2018-06-19 13:03:0550// Q15 beta = 0.57762063060713 const Factor for magnitude approximation
andrew@webrtc.orge03cafa2013-11-11 20:10:0151static const uint16_t kBeta3 = 18927;
52#endif
53
54static const int16_t kNoiseEstQDomain = 15;
55static const int16_t kNoiseEstIncCount = 5;
56
pbos@webrtc.orge468bc92014-12-18 09:11:3357static void ComfortNoise(AecmCore* aecm,
andrew@webrtc.orge03cafa2013-11-11 20:10:0158 const uint16_t* dfa,
pbos@webrtc.orge468bc92014-12-18 09:11:3359 ComplexInt16* out,
andrew@webrtc.orge03cafa2013-11-11 20:10:0160 const int16_t* lambda);
61
pbos@webrtc.orge468bc92014-12-18 09:11:3362static void WindowAndFFT(AecmCore* aecm,
63 int16_t* fft,
64 const int16_t* time_signal,
65 ComplexInt16* freq_signal,
66 int time_signal_scaling) {
andrew@webrtc.orge03cafa2013-11-11 20:10:0167 int i = 0;
68
69 // FFT of signal
70 for (i = 0; i < PART_LEN; i++) {
71 // Window time domain signal and insert into real part of
72 // transformation array |fft|
Alex Loikoee67ca32018-01-12 12:48:0673 int16_t scaled_time_signal = time_signal[i] * (1 << time_signal_scaling);
bjornv@webrtc.orgb38b0092015-03-10 06:40:0274 fft[i] = (int16_t)((scaled_time_signal * WebRtcAecm_kSqrtHanning[i]) >> 14);
Alex Loikoee67ca32018-01-12 12:48:0675 scaled_time_signal = time_signal[i + PART_LEN] * (1 << time_signal_scaling);
Yves Gerey665174f2018-06-19 13:03:0576 fft[PART_LEN + i] = (int16_t)(
77 (scaled_time_signal * WebRtcAecm_kSqrtHanning[PART_LEN - i]) >> 14);
andrew@webrtc.orge03cafa2013-11-11 20:10:0178 }
79
80 // Do forward FFT, then take only the first PART_LEN complex samples,
81 // and change signs of the imaginary parts.
82 WebRtcSpl_RealForwardFFT(aecm->real_fft, fft, (int16_t*)freq_signal);
83 for (i = 0; i < PART_LEN; i++) {
84 freq_signal[i].imag = -freq_signal[i].imag;
85 }
86}
87
pbos@webrtc.orge468bc92014-12-18 09:11:3388static void InverseFFTAndWindow(AecmCore* aecm,
andrew@webrtc.orge03cafa2013-11-11 20:10:0189 int16_t* fft,
pbos@webrtc.orge468bc92014-12-18 09:11:3390 ComplexInt16* efw,
andrew@webrtc.orge03cafa2013-11-11 20:10:0191 int16_t* output,
pbos@webrtc.orge468bc92014-12-18 09:11:3392 const int16_t* nearendClean) {
andrew@webrtc.orge03cafa2013-11-11 20:10:0193 int i, j, outCFFT;
94 int32_t tmp32no1;
95 // Reuse |efw| for the inverse FFT output after transferring
96 // the contents to |fft|.
97 int16_t* ifft_out = (int16_t*)efw;
98
99 // Synthesis
100 for (i = 1, j = 2; i < PART_LEN; i += 1, j += 2) {
101 fft[j] = efw[i].real;
102 fft[j + 1] = -efw[i].imag;
103 }
104 fft[0] = efw[0].real;
105 fft[1] = -efw[0].imag;
106
107 fft[PART_LEN2] = efw[PART_LEN].real;
108 fft[PART_LEN2 + 1] = -efw[PART_LEN].imag;
109
110 // Inverse FFT. Keep outCFFT to scale the samples in the next block.
111 outCFFT = WebRtcSpl_RealInverseFFT(aecm->real_fft, fft, ifft_out);
112 for (i = 0; i < PART_LEN; i++) {
113 ifft_out[i] = (int16_t)WEBRTC_SPL_MUL_16_16_RSFT_WITH_ROUND(
Yves Gerey665174f2018-06-19 13:03:05114 ifft_out[i], WebRtcAecm_kSqrtHanning[i], 14);
andrew@webrtc.orge03cafa2013-11-11 20:10:01115 tmp32no1 = WEBRTC_SPL_SHIFT_W32((int32_t)ifft_out[i],
Yves Gerey665174f2018-06-19 13:03:05116 outCFFT - aecm->dfaCleanQDomain);
andrew@webrtc.orge03cafa2013-11-11 20:10:01117 output[i] = (int16_t)WEBRTC_SPL_SAT(WEBRTC_SPL_WORD16_MAX,
118 tmp32no1 + aecm->outBuf[i],
119 WEBRTC_SPL_WORD16_MIN);
120
Yves Gerey665174f2018-06-19 13:03:05121 tmp32no1 =
122 (ifft_out[PART_LEN + i] * WebRtcAecm_kSqrtHanning[PART_LEN - i]) >> 14;
123 tmp32no1 = WEBRTC_SPL_SHIFT_W32(tmp32no1, outCFFT - aecm->dfaCleanQDomain);
124 aecm->outBuf[i] = (int16_t)WEBRTC_SPL_SAT(WEBRTC_SPL_WORD16_MAX, tmp32no1,
125 WEBRTC_SPL_WORD16_MIN);
andrew@webrtc.orge03cafa2013-11-11 20:10:01126 }
127
128 // Copy the current block to the old position
129 // (aecm->outBuf is shifted elsewhere)
130 memcpy(aecm->xBuf, aecm->xBuf + PART_LEN, sizeof(int16_t) * PART_LEN);
Yves Gerey665174f2018-06-19 13:03:05131 memcpy(aecm->dBufNoisy, aecm->dBufNoisy + PART_LEN,
andrew@webrtc.orge03cafa2013-11-11 20:10:01132 sizeof(int16_t) * PART_LEN);
Yves Gerey665174f2018-06-19 13:03:05133 if (nearendClean != NULL) {
134 memcpy(aecm->dBufClean, aecm->dBufClean + PART_LEN,
andrew@webrtc.orge03cafa2013-11-11 20:10:01135 sizeof(int16_t) * PART_LEN);
136 }
137}
138
139// Transforms a time domain signal into the frequency domain, outputting the
140// complex valued signal, absolute value and sum of absolute values.
141//
142// time_signal [in] Pointer to time domain signal
143// freq_signal_real [out] Pointer to real part of frequency domain array
144// freq_signal_imag [out] Pointer to imaginary part of frequency domain
145// array
146// freq_signal_abs [out] Pointer to absolute value of frequency domain
147// array
148// freq_signal_sum_abs [out] Pointer to the sum of all absolute values in
149// the frequency domain array
150// return value The Q-domain of current frequency values
151//
pbos@webrtc.orge468bc92014-12-18 09:11:33152static int TimeToFrequencyDomain(AecmCore* aecm,
andrew@webrtc.orge03cafa2013-11-11 20:10:01153 const int16_t* time_signal,
pbos@webrtc.orge468bc92014-12-18 09:11:33154 ComplexInt16* freq_signal,
andrew@webrtc.orge03cafa2013-11-11 20:10:01155 uint16_t* freq_signal_abs,
pbos@webrtc.orge468bc92014-12-18 09:11:33156 uint32_t* freq_signal_sum_abs) {
andrew@webrtc.orge03cafa2013-11-11 20:10:01157 int i = 0;
158 int time_signal_scaling = 0;
159
160 int32_t tmp32no1 = 0;
161 int32_t tmp32no2 = 0;
162
163 // In fft_buf, +16 for 32-byte alignment.
164 int16_t fft_buf[PART_LEN4 + 16];
Yves Gerey665174f2018-06-19 13:03:05165 int16_t* fft = (int16_t*)(((uintptr_t)fft_buf + 31) & ~31);
andrew@webrtc.orge03cafa2013-11-11 20:10:01166
167 int16_t tmp16no1;
168#ifndef WEBRTC_ARCH_ARM_V7
169 int16_t tmp16no2;
170#endif
171#ifdef AECM_WITH_ABS_APPROX
172 int16_t max_value = 0;
173 int16_t min_value = 0;
174 uint16_t alpha = 0;
175 uint16_t beta = 0;
176#endif
177
178#ifdef AECM_DYNAMIC_Q
179 tmp16no1 = WebRtcSpl_MaxAbsValueW16(time_signal, PART_LEN2);
180 time_signal_scaling = WebRtcSpl_NormW16(tmp16no1);
181#endif
182
183 WindowAndFFT(aecm, fft, time_signal, freq_signal, time_signal_scaling);
184
185 // Extract imaginary and real part, calculate the magnitude for
186 // all frequency bins
187 freq_signal[0].imag = 0;
188 freq_signal[PART_LEN].imag = 0;
189 freq_signal_abs[0] = (uint16_t)WEBRTC_SPL_ABS_W16(freq_signal[0].real);
Yves Gerey665174f2018-06-19 13:03:05190 freq_signal_abs[PART_LEN] =
191 (uint16_t)WEBRTC_SPL_ABS_W16(freq_signal[PART_LEN].real);
192 (*freq_signal_sum_abs) =
193 (uint32_t)(freq_signal_abs[0]) + (uint32_t)(freq_signal_abs[PART_LEN]);
andrew@webrtc.orge03cafa2013-11-11 20:10:01194
Yves Gerey665174f2018-06-19 13:03:05195 for (i = 1; i < PART_LEN; i++) {
196 if (freq_signal[i].real == 0) {
andrew@webrtc.orge03cafa2013-11-11 20:10:01197 freq_signal_abs[i] = (uint16_t)WEBRTC_SPL_ABS_W16(freq_signal[i].imag);
Yves Gerey665174f2018-06-19 13:03:05198 } else if (freq_signal[i].imag == 0) {
andrew@webrtc.orge03cafa2013-11-11 20:10:01199 freq_signal_abs[i] = (uint16_t)WEBRTC_SPL_ABS_W16(freq_signal[i].real);
Yves Gerey665174f2018-06-19 13:03:05200 } else {
201// Approximation for magnitude of complex fft output
202// magn = sqrt(real^2 + imag^2)
203// magn ~= alpha * max(|imag|,|real|) + beta * min(|imag|,|real|)
204//
205// The parameters alpha and beta are stored in Q15
andrew@webrtc.orge03cafa2013-11-11 20:10:01206
207#ifdef AECM_WITH_ABS_APPROX
208 tmp16no1 = WEBRTC_SPL_ABS_W16(freq_signal[i].real);
209 tmp16no2 = WEBRTC_SPL_ABS_W16(freq_signal[i].imag);
210
Yves Gerey665174f2018-06-19 13:03:05211 if (tmp16no1 > tmp16no2) {
andrew@webrtc.orge03cafa2013-11-11 20:10:01212 max_value = tmp16no1;
213 min_value = tmp16no2;
Yves Gerey665174f2018-06-19 13:03:05214 } else {
andrew@webrtc.orge03cafa2013-11-11 20:10:01215 max_value = tmp16no2;
216 min_value = tmp16no1;
217 }
218
219 // Magnitude in Q(-6)
Yves Gerey665174f2018-06-19 13:03:05220 if ((max_value >> 2) > min_value) {
andrew@webrtc.orge03cafa2013-11-11 20:10:01221 alpha = kAlpha1;
222 beta = kBeta1;
Yves Gerey665174f2018-06-19 13:03:05223 } else if ((max_value >> 1) > min_value) {
andrew@webrtc.orge03cafa2013-11-11 20:10:01224 alpha = kAlpha2;
225 beta = kBeta2;
Yves Gerey665174f2018-06-19 13:03:05226 } else {
andrew@webrtc.orge03cafa2013-11-11 20:10:01227 alpha = kAlpha3;
228 beta = kBeta3;
229 }
bjornv@webrtc.orgb38b0092015-03-10 06:40:02230 tmp16no1 = (int16_t)((max_value * alpha) >> 15);
231 tmp16no2 = (int16_t)((min_value * beta) >> 15);
andrew@webrtc.orge03cafa2013-11-11 20:10:01232 freq_signal_abs[i] = (uint16_t)tmp16no1 + (uint16_t)tmp16no2;
233#else
234#ifdef WEBRTC_ARCH_ARM_V7
235 __asm __volatile(
Yves Gerey665174f2018-06-19 13:03:05236 "smulbb %[tmp32no1], %[real], %[real]\n\t"
237 "smlabb %[tmp32no2], %[imag], %[imag], %[tmp32no1]\n\t"
238 : [tmp32no1] "+&r"(tmp32no1), [tmp32no2] "=r"(tmp32no2)
239 : [real] "r"(freq_signal[i].real), [imag] "r"(freq_signal[i].imag));
andrew@webrtc.orge03cafa2013-11-11 20:10:01240#else
241 tmp16no1 = WEBRTC_SPL_ABS_W16(freq_signal[i].real);
242 tmp16no2 = WEBRTC_SPL_ABS_W16(freq_signal[i].imag);
bjornv@webrtc.org758d6d42015-01-08 17:52:56243 tmp32no1 = tmp16no1 * tmp16no1;
244 tmp32no2 = tmp16no2 * tmp16no2;
bjornv@webrtc.org6e71d172014-08-25 07:44:52245 tmp32no2 = WebRtcSpl_AddSatW32(tmp32no1, tmp32no2);
Yves Gerey665174f2018-06-19 13:03:05246#endif // WEBRTC_ARCH_ARM_V7
andrew@webrtc.orge03cafa2013-11-11 20:10:01247 tmp32no1 = WebRtcSpl_SqrtFloor(tmp32no2);
248
249 freq_signal_abs[i] = (uint16_t)tmp32no1;
Yves Gerey665174f2018-06-19 13:03:05250#endif // AECM_WITH_ABS_APPROX
andrew@webrtc.orge03cafa2013-11-11 20:10:01251 }
252 (*freq_signal_sum_abs) += (uint32_t)freq_signal_abs[i];
253 }
254
255 return time_signal_scaling;
256}
257
oprypin8ad0e582017-09-05 10:00:37258int RTC_NO_SANITIZE("signed-integer-overflow") // bugs.webrtc.org/8200
Yves Gerey665174f2018-06-19 13:03:05259 WebRtcAecm_ProcessBlock(AecmCore* aecm,
260 const int16_t* farend,
261 const int16_t* nearendNoisy,
262 const int16_t* nearendClean,
263 int16_t* output) {
andrew@webrtc.orge03cafa2013-11-11 20:10:01264 int i;
265
266 uint32_t xfaSum;
267 uint32_t dfaNoisySum;
268 uint32_t dfaCleanSum;
269 uint32_t echoEst32Gained;
270 uint32_t tmpU32;
271
272 int32_t tmp32no1;
273
274 uint16_t xfa[PART_LEN1];
275 uint16_t dfaNoisy[PART_LEN1];
276 uint16_t dfaClean[PART_LEN1];
277 uint16_t* ptrDfaClean = dfaClean;
278 const uint16_t* far_spectrum_ptr = NULL;
279
280 // 32 byte aligned buffers (with +8 or +16).
pbos@webrtc.orge468bc92014-12-18 09:11:33281 // TODO(kma): define fft with ComplexInt16.
Yves Gerey665174f2018-06-19 13:03:05282 int16_t fft_buf[PART_LEN4 + 2 + 16]; // +2 to make a loop safe.
andrew@webrtc.orge03cafa2013-11-11 20:10:01283 int32_t echoEst32_buf[PART_LEN1 + 8];
284 int32_t dfw_buf[PART_LEN2 + 8];
285 int32_t efw_buf[PART_LEN2 + 8];
286
Yves Gerey665174f2018-06-19 13:03:05287 int16_t* fft = (int16_t*)(((uintptr_t)fft_buf + 31) & ~31);
288 int32_t* echoEst32 = (int32_t*)(((uintptr_t)echoEst32_buf + 31) & ~31);
pbos@webrtc.orge468bc92014-12-18 09:11:33289 ComplexInt16* dfw = (ComplexInt16*)(((uintptr_t)dfw_buf + 31) & ~31);
290 ComplexInt16* efw = (ComplexInt16*)(((uintptr_t)efw_buf + 31) & ~31);
andrew@webrtc.orge03cafa2013-11-11 20:10:01291
292 int16_t hnl[PART_LEN1];
293 int16_t numPosCoef = 0;
294 int16_t nlpGain = ONE_Q14;
295 int delay;
296 int16_t tmp16no1;
297 int16_t tmp16no2;
298 int16_t mu;
299 int16_t supGain;
300 int16_t zeros32, zeros16;
301 int16_t zerosDBufNoisy, zerosDBufClean, zerosXBuf;
302 int far_q;
bjornv@webrtc.orgc0ba4392014-07-03 13:38:53303 int16_t resolutionDiff, qDomainDiff, dfa_clean_q_domain_diff;
andrew@webrtc.orge03cafa2013-11-11 20:10:01304
305 const int kMinPrefBand = 4;
306 const int kMaxPrefBand = 24;
307 int32_t avgHnl32 = 0;
308
309 // Determine startup state. There are three states:
310 // (0) the first CONV_LEN blocks
311 // (1) another CONV_LEN blocks
312 // (2) the rest
313
Yves Gerey665174f2018-06-19 13:03:05314 if (aecm->startupState < 2) {
315 aecm->startupState =
316 (aecm->totCount >= CONV_LEN) + (aecm->totCount >= CONV_LEN2);
andrew@webrtc.orge03cafa2013-11-11 20:10:01317 }
318 // END: Determine startup state
319
320 // Buffer near and far end signals
321 memcpy(aecm->xBuf + PART_LEN, farend, sizeof(int16_t) * PART_LEN);
322 memcpy(aecm->dBufNoisy + PART_LEN, nearendNoisy, sizeof(int16_t) * PART_LEN);
Yves Gerey665174f2018-06-19 13:03:05323 if (nearendClean != NULL) {
324 memcpy(aecm->dBufClean + PART_LEN, nearendClean,
andrew@webrtc.orge03cafa2013-11-11 20:10:01325 sizeof(int16_t) * PART_LEN);
326 }
327
328 // Transform far end signal from time domain to frequency domain.
Yves Gerey665174f2018-06-19 13:03:05329 far_q = TimeToFrequencyDomain(aecm, aecm->xBuf, dfw, xfa, &xfaSum);
andrew@webrtc.orge03cafa2013-11-11 20:10:01330
331 // Transform noisy near end signal from time domain to frequency domain.
Yves Gerey665174f2018-06-19 13:03:05332 zerosDBufNoisy =
333 TimeToFrequencyDomain(aecm, aecm->dBufNoisy, dfw, dfaNoisy, &dfaNoisySum);
andrew@webrtc.orge03cafa2013-11-11 20:10:01334 aecm->dfaNoisyQDomainOld = aecm->dfaNoisyQDomain;
335 aecm->dfaNoisyQDomain = (int16_t)zerosDBufNoisy;
336
Yves Gerey665174f2018-06-19 13:03:05337 if (nearendClean == NULL) {
andrew@webrtc.orge03cafa2013-11-11 20:10:01338 ptrDfaClean = dfaNoisy;
339 aecm->dfaCleanQDomainOld = aecm->dfaNoisyQDomainOld;
340 aecm->dfaCleanQDomain = aecm->dfaNoisyQDomain;
341 dfaCleanSum = dfaNoisySum;
Yves Gerey665174f2018-06-19 13:03:05342 } else {
andrew@webrtc.orge03cafa2013-11-11 20:10:01343 // Transform clean near end signal from time domain to frequency domain.
Yves Gerey665174f2018-06-19 13:03:05344 zerosDBufClean = TimeToFrequencyDomain(aecm, aecm->dBufClean, dfw, dfaClean,
andrew@webrtc.orge03cafa2013-11-11 20:10:01345 &dfaCleanSum);
346 aecm->dfaCleanQDomainOld = aecm->dfaCleanQDomain;
347 aecm->dfaCleanQDomain = (int16_t)zerosDBufClean;
348 }
349
350 // Get the delay
351 // Save far-end history and estimate delay
352 WebRtcAecm_UpdateFarHistory(aecm, xfa, far_q);
Yves Gerey665174f2018-06-19 13:03:05353 if (WebRtc_AddFarSpectrumFix(aecm->delay_estimator_farend, xfa, PART_LEN1,
andrew@webrtc.orge03cafa2013-11-11 20:10:01354 far_q) == -1) {
355 return -1;
356 }
Yves Gerey665174f2018-06-19 13:03:05357 delay = WebRtc_DelayEstimatorProcessFix(aecm->delay_estimator, dfaNoisy,
358 PART_LEN1, zerosDBufNoisy);
359 if (delay == -1) {
andrew@webrtc.orge03cafa2013-11-11 20:10:01360 return -1;
Yves Gerey665174f2018-06-19 13:03:05361 } else if (delay == -2) {
andrew@webrtc.orge03cafa2013-11-11 20:10:01362 // If the delay is unknown, we assume zero.
363 // NOTE: this will have to be adjusted if we ever add lookahead.
364 delay = 0;
365 }
366
Yves Gerey665174f2018-06-19 13:03:05367 if (aecm->fixedDelay >= 0) {
andrew@webrtc.orge03cafa2013-11-11 20:10:01368 // Use fixed delay
369 delay = aecm->fixedDelay;
370 }
371
372 // Get aligned far end spectrum
373 far_spectrum_ptr = WebRtcAecm_AlignedFarend(aecm, &far_q, delay);
Yves Gerey665174f2018-06-19 13:03:05374 zerosXBuf = (int16_t)far_q;
375 if (far_spectrum_ptr == NULL) {
andrew@webrtc.orge03cafa2013-11-11 20:10:01376 return -1;
377 }
378
379 // Calculate log(energy) and update energy threshold levels
Yves Gerey665174f2018-06-19 13:03:05380 WebRtcAecm_CalcEnergies(aecm, far_spectrum_ptr, zerosXBuf, dfaNoisySum,
andrew@webrtc.orge03cafa2013-11-11 20:10:01381 echoEst32);
382
383 // Calculate stepsize
384 mu = WebRtcAecm_CalcStepSize(aecm);
385
386 // Update counters
387 aecm->totCount++;
388
389 // This is the channel estimation algorithm.
390 // It is base on NLMS but has a variable step length,
391 // which was calculated above.
Yves Gerey665174f2018-06-19 13:03:05392 WebRtcAecm_UpdateChannel(aecm, far_spectrum_ptr, zerosXBuf, dfaNoisy, mu,
andrew@webrtc.orge03cafa2013-11-11 20:10:01393 echoEst32);
394 supGain = WebRtcAecm_CalcSuppressionGain(aecm);
395
andrew@webrtc.orge03cafa2013-11-11 20:10:01396 // Calculate Wiener filter hnl[]
Yves Gerey665174f2018-06-19 13:03:05397 for (i = 0; i < PART_LEN1; i++) {
andrew@webrtc.orge03cafa2013-11-11 20:10:01398 // Far end signal through channel estimate in Q8
399 // How much can we shift right to preserve resolution
400 tmp32no1 = echoEst32[i] - aecm->echoFilt[i];
Alessio Bazzicaba68aab2017-10-16 11:45:58401 aecm->echoFilt[i] +=
402 rtc::dchecked_cast<int32_t>((int64_t{tmp32no1} * 50) >> 8);
andrew@webrtc.orge03cafa2013-11-11 20:10:01403
404 zeros32 = WebRtcSpl_NormW32(aecm->echoFilt[i]) + 1;
405 zeros16 = WebRtcSpl_NormW16(supGain) + 1;
Yves Gerey665174f2018-06-19 13:03:05406 if (zeros32 + zeros16 > 16) {
andrew@webrtc.orge03cafa2013-11-11 20:10:01407 // Multiplication is safe
408 // Result in
409 // Q(RESOLUTION_CHANNEL+RESOLUTION_SUPGAIN+
410 // aecm->xfaQDomainBuf[diff])
Yves Gerey665174f2018-06-19 13:03:05411 echoEst32Gained =
412 WEBRTC_SPL_UMUL_32_16((uint32_t)aecm->echoFilt[i], (uint16_t)supGain);
andrew@webrtc.orge03cafa2013-11-11 20:10:01413 resolutionDiff = 14 - RESOLUTION_CHANNEL16 - RESOLUTION_SUPGAIN;
414 resolutionDiff += (aecm->dfaCleanQDomain - zerosXBuf);
Yves Gerey665174f2018-06-19 13:03:05415 } else {
andrew@webrtc.orge03cafa2013-11-11 20:10:01416 tmp16no1 = 17 - zeros32 - zeros16;
Yves Gerey665174f2018-06-19 13:03:05417 resolutionDiff =
418 14 + tmp16no1 - RESOLUTION_CHANNEL16 - RESOLUTION_SUPGAIN;
andrew@webrtc.orge03cafa2013-11-11 20:10:01419 resolutionDiff += (aecm->dfaCleanQDomain - zerosXBuf);
Yves Gerey665174f2018-06-19 13:03:05420 if (zeros32 > tmp16no1) {
andrew@webrtc.orge03cafa2013-11-11 20:10:01421 echoEst32Gained = WEBRTC_SPL_UMUL_32_16((uint32_t)aecm->echoFilt[i],
bjornv@webrtc.orgd4fe8242014-10-13 13:01:13422 supGain >> tmp16no1);
Yves Gerey665174f2018-06-19 13:03:05423 } else {
andrew@webrtc.orge03cafa2013-11-11 20:10:01424 // Result in Q-(RESOLUTION_CHANNEL+RESOLUTION_SUPGAIN-16)
bjornv@webrtc.orgf87c0af2014-10-15 12:51:23425 echoEst32Gained = (aecm->echoFilt[i] >> tmp16no1) * supGain;
andrew@webrtc.orge03cafa2013-11-11 20:10:01426 }
427 }
428
429 zeros16 = WebRtcSpl_NormW16(aecm->nearFilt[i]);
kwiberg9e2be5f2016-09-14 12:23:22430 RTC_DCHECK_GE(zeros16, 0); // |zeros16| is a norm, hence non-negative.
bjornv@webrtc.orgc0ba4392014-07-03 13:38:53431 dfa_clean_q_domain_diff = aecm->dfaCleanQDomain - aecm->dfaCleanQDomainOld;
432 if (zeros16 < dfa_clean_q_domain_diff && aecm->nearFilt[i]) {
Alex Loiko600bdb42018-01-25 10:42:43433 tmp16no1 = aecm->nearFilt[i] * (1 << zeros16);
bjornv@webrtc.orgc0ba4392014-07-03 13:38:53434 qDomainDiff = zeros16 - dfa_clean_q_domain_diff;
435 tmp16no2 = ptrDfaClean[i] >> -qDomainDiff;
436 } else {
437 tmp16no1 = dfa_clean_q_domain_diff < 0
Alex Loiko736d2f72018-01-22 10:52:31438 ? aecm->nearFilt[i] >> -dfa_clean_q_domain_diff
439 : aecm->nearFilt[i] * (1 << dfa_clean_q_domain_diff);
andrew@webrtc.orge03cafa2013-11-11 20:10:01440 qDomainDiff = 0;
bjornv@webrtc.orgc0ba4392014-07-03 13:38:53441 tmp16no2 = ptrDfaClean[i];
andrew@webrtc.orge03cafa2013-11-11 20:10:01442 }
andrew@webrtc.orge03cafa2013-11-11 20:10:01443 tmp32no1 = (int32_t)(tmp16no2 - tmp16no1);
bjornv@webrtc.orgf87c0af2014-10-15 12:51:23444 tmp16no2 = (int16_t)(tmp32no1 >> 4);
andrew@webrtc.orge03cafa2013-11-11 20:10:01445 tmp16no2 += tmp16no1;
446 zeros16 = WebRtcSpl_NormW16(tmp16no2);
bjornv@webrtc.orgc0ba4392014-07-03 13:38:53447 if ((tmp16no2) & (-qDomainDiff > zeros16)) {
andrew@webrtc.orge03cafa2013-11-11 20:10:01448 aecm->nearFilt[i] = WEBRTC_SPL_WORD16_MAX;
bjornv@webrtc.orgc0ba4392014-07-03 13:38:53449 } else {
Alex Loiko736d2f72018-01-22 10:52:31450 aecm->nearFilt[i] = qDomainDiff < 0 ? tmp16no2 * (1 << -qDomainDiff)
bjornv@webrtc.orgc0ba4392014-07-03 13:38:53451 : tmp16no2 >> qDomainDiff;
andrew@webrtc.orge03cafa2013-11-11 20:10:01452 }
453
454 // Wiener filter coefficients, resulting hnl in Q14
Yves Gerey665174f2018-06-19 13:03:05455 if (echoEst32Gained == 0) {
andrew@webrtc.orge03cafa2013-11-11 20:10:01456 hnl[i] = ONE_Q14;
Yves Gerey665174f2018-06-19 13:03:05457 } else if (aecm->nearFilt[i] == 0) {
andrew@webrtc.orge03cafa2013-11-11 20:10:01458 hnl[i] = 0;
Yves Gerey665174f2018-06-19 13:03:05459 } else {
andrew@webrtc.orge03cafa2013-11-11 20:10:01460 // Multiply the suppression gain
461 // Rounding
462 echoEst32Gained += (uint32_t)(aecm->nearFilt[i] >> 1);
Yves Gerey665174f2018-06-19 13:03:05463 tmpU32 =
464 WebRtcSpl_DivU32U16(echoEst32Gained, (uint16_t)aecm->nearFilt[i]);
andrew@webrtc.orge03cafa2013-11-11 20:10:01465
466 // Current resolution is
467 // Q-(RESOLUTION_CHANNEL+RESOLUTION_SUPGAIN- max(0,17-zeros16- zeros32))
468 // Make sure we are in Q14
469 tmp32no1 = (int32_t)WEBRTC_SPL_SHIFT_W32(tmpU32, resolutionDiff);
Yves Gerey665174f2018-06-19 13:03:05470 if (tmp32no1 > ONE_Q14) {
andrew@webrtc.orge03cafa2013-11-11 20:10:01471 hnl[i] = 0;
Yves Gerey665174f2018-06-19 13:03:05472 } else if (tmp32no1 < 0) {
andrew@webrtc.orge03cafa2013-11-11 20:10:01473 hnl[i] = ONE_Q14;
Yves Gerey665174f2018-06-19 13:03:05474 } else {
andrew@webrtc.orge03cafa2013-11-11 20:10:01475 // 1-echoEst/dfa
476 hnl[i] = ONE_Q14 - (int16_t)tmp32no1;
Yves Gerey665174f2018-06-19 13:03:05477 if (hnl[i] < 0) {
andrew@webrtc.orge03cafa2013-11-11 20:10:01478 hnl[i] = 0;
479 }
480 }
481 }
Yves Gerey665174f2018-06-19 13:03:05482 if (hnl[i]) {
andrew@webrtc.orge03cafa2013-11-11 20:10:01483 numPosCoef++;
484 }
485 }
486 // Only in wideband. Prevent the gain in upper band from being larger than
487 // in lower band.
Yves Gerey665174f2018-06-19 13:03:05488 if (aecm->mult == 2) {
andrew@webrtc.orge03cafa2013-11-11 20:10:01489 // TODO(bjornv): Investigate if the scaling of hnl[i] below can cause
490 // speech distortion in double-talk.
Yves Gerey665174f2018-06-19 13:03:05491 for (i = 0; i < PART_LEN1; i++) {
bjornv@webrtc.orgb38b0092015-03-10 06:40:02492 hnl[i] = (int16_t)((hnl[i] * hnl[i]) >> 14);
andrew@webrtc.orge03cafa2013-11-11 20:10:01493 }
494
Yves Gerey665174f2018-06-19 13:03:05495 for (i = kMinPrefBand; i <= kMaxPrefBand; i++) {
andrew@webrtc.orge03cafa2013-11-11 20:10:01496 avgHnl32 += (int32_t)hnl[i];
497 }
kwiberg9e2be5f2016-09-14 12:23:22498 RTC_DCHECK_GT(kMaxPrefBand - kMinPrefBand + 1, 0);
andrew@webrtc.orge03cafa2013-11-11 20:10:01499 avgHnl32 /= (kMaxPrefBand - kMinPrefBand + 1);
500
Yves Gerey665174f2018-06-19 13:03:05501 for (i = kMaxPrefBand; i < PART_LEN1; i++) {
502 if (hnl[i] > (int16_t)avgHnl32) {
andrew@webrtc.orge03cafa2013-11-11 20:10:01503 hnl[i] = (int16_t)avgHnl32;
504 }
505 }
506 }
507
508 // Calculate NLP gain, result is in Q14
Yves Gerey665174f2018-06-19 13:03:05509 if (aecm->nlpFlag) {
510 for (i = 0; i < PART_LEN1; i++) {
andrew@webrtc.orge03cafa2013-11-11 20:10:01511 // Truncate values close to zero and one.
Yves Gerey665174f2018-06-19 13:03:05512 if (hnl[i] > NLP_COMP_HIGH) {
andrew@webrtc.orge03cafa2013-11-11 20:10:01513 hnl[i] = ONE_Q14;
Yves Gerey665174f2018-06-19 13:03:05514 } else if (hnl[i] < NLP_COMP_LOW) {
andrew@webrtc.orge03cafa2013-11-11 20:10:01515 hnl[i] = 0;
516 }
517
518 // Remove outliers
Yves Gerey665174f2018-06-19 13:03:05519 if (numPosCoef < 3) {
andrew@webrtc.orge03cafa2013-11-11 20:10:01520 nlpGain = 0;
Yves Gerey665174f2018-06-19 13:03:05521 } else {
andrew@webrtc.orge03cafa2013-11-11 20:10:01522 nlpGain = ONE_Q14;
523 }
524
525 // NLP
Yves Gerey665174f2018-06-19 13:03:05526 if ((hnl[i] == ONE_Q14) && (nlpGain == ONE_Q14)) {
andrew@webrtc.orge03cafa2013-11-11 20:10:01527 hnl[i] = ONE_Q14;
Yves Gerey665174f2018-06-19 13:03:05528 } else {
bjornv@webrtc.orgb38b0092015-03-10 06:40:02529 hnl[i] = (int16_t)((hnl[i] * nlpGain) >> 14);
andrew@webrtc.orge03cafa2013-11-11 20:10:01530 }
531
532 // multiply with Wiener coefficients
Yves Gerey665174f2018-06-19 13:03:05533 efw[i].real = (int16_t)(
534 WEBRTC_SPL_MUL_16_16_RSFT_WITH_ROUND(dfw[i].real, hnl[i], 14));
535 efw[i].imag = (int16_t)(
536 WEBRTC_SPL_MUL_16_16_RSFT_WITH_ROUND(dfw[i].imag, hnl[i], 14));
andrew@webrtc.orge03cafa2013-11-11 20:10:01537 }
Yves Gerey665174f2018-06-19 13:03:05538 } else {
andrew@webrtc.orge03cafa2013-11-11 20:10:01539 // multiply with Wiener coefficients
Yves Gerey665174f2018-06-19 13:03:05540 for (i = 0; i < PART_LEN1; i++) {
541 efw[i].real = (int16_t)(
542 WEBRTC_SPL_MUL_16_16_RSFT_WITH_ROUND(dfw[i].real, hnl[i], 14));
543 efw[i].imag = (int16_t)(
544 WEBRTC_SPL_MUL_16_16_RSFT_WITH_ROUND(dfw[i].imag, hnl[i], 14));
andrew@webrtc.orge03cafa2013-11-11 20:10:01545 }
546 }
547
Yves Gerey665174f2018-06-19 13:03:05548 if (aecm->cngMode == AecmTrue) {
andrew@webrtc.orge03cafa2013-11-11 20:10:01549 ComfortNoise(aecm, ptrDfaClean, efw, hnl);
550 }
551
552 InverseFFTAndWindow(aecm, fft, efw, output, nearendClean);
553
554 return 0;
555}
556
pbos@webrtc.orge468bc92014-12-18 09:11:33557static void ComfortNoise(AecmCore* aecm,
andrew@webrtc.orge03cafa2013-11-11 20:10:01558 const uint16_t* dfa,
pbos@webrtc.orge468bc92014-12-18 09:11:33559 ComplexInt16* out,
560 const int16_t* lambda) {
andrew@webrtc.orge03cafa2013-11-11 20:10:01561 int16_t i;
562 int16_t tmp16;
563 int32_t tmp32;
564
565 int16_t randW16[PART_LEN];
566 int16_t uReal[PART_LEN1];
567 int16_t uImag[PART_LEN1];
568 int32_t outLShift32;
569 int16_t noiseRShift16[PART_LEN1];
570
571 int16_t shiftFromNearToNoise = kNoiseEstQDomain - aecm->dfaCleanQDomain;
572 int16_t minTrackShift;
573
kwiberg9e2be5f2016-09-14 12:23:22574 RTC_DCHECK_GE(shiftFromNearToNoise, 0);
575 RTC_DCHECK_LT(shiftFromNearToNoise, 16);
andrew@webrtc.orge03cafa2013-11-11 20:10:01576
Yves Gerey665174f2018-06-19 13:03:05577 if (aecm->noiseEstCtr < 100) {
andrew@webrtc.orge03cafa2013-11-11 20:10:01578 // Track the minimum more quickly initially.
579 aecm->noiseEstCtr++;
580 minTrackShift = 6;
Yves Gerey665174f2018-06-19 13:03:05581 } else {
andrew@webrtc.orge03cafa2013-11-11 20:10:01582 minTrackShift = 9;
583 }
584
585 // Estimate noise power.
Yves Gerey665174f2018-06-19 13:03:05586 for (i = 0; i < PART_LEN1; i++) {
andrew@webrtc.orge03cafa2013-11-11 20:10:01587 // Shift to the noise domain.
588 tmp32 = (int32_t)dfa[i];
bjornv@webrtc.org750423c2014-09-30 09:26:36589 outLShift32 = tmp32 << shiftFromNearToNoise;
andrew@webrtc.orge03cafa2013-11-11 20:10:01590
Yves Gerey665174f2018-06-19 13:03:05591 if (outLShift32 < aecm->noiseEst[i]) {
andrew@webrtc.orge03cafa2013-11-11 20:10:01592 // Reset "too low" counter
593 aecm->noiseEstTooLowCtr[i] = 0;
594 // Track the minimum.
Yves Gerey665174f2018-06-19 13:03:05595 if (aecm->noiseEst[i] < (1 << minTrackShift)) {
andrew@webrtc.orge03cafa2013-11-11 20:10:01596 // For small values, decrease noiseEst[i] every
597 // |kNoiseEstIncCount| block. The regular approach below can not
598 // go further down due to truncation.
599 aecm->noiseEstTooHighCtr[i]++;
Yves Gerey665174f2018-06-19 13:03:05600 if (aecm->noiseEstTooHighCtr[i] >= kNoiseEstIncCount) {
andrew@webrtc.orge03cafa2013-11-11 20:10:01601 aecm->noiseEst[i]--;
Yves Gerey665174f2018-06-19 13:03:05602 aecm->noiseEstTooHighCtr[i] = 0; // Reset the counter
andrew@webrtc.orge03cafa2013-11-11 20:10:01603 }
Yves Gerey665174f2018-06-19 13:03:05604 } else {
605 aecm->noiseEst[i] -=
606 ((aecm->noiseEst[i] - outLShift32) >> minTrackShift);
andrew@webrtc.orge03cafa2013-11-11 20:10:01607 }
Yves Gerey665174f2018-06-19 13:03:05608 } else {
andrew@webrtc.orge03cafa2013-11-11 20:10:01609 // Reset "too high" counter
610 aecm->noiseEstTooHighCtr[i] = 0;
611 // Ramp slowly upwards until we hit the minimum again.
Yves Gerey665174f2018-06-19 13:03:05612 if ((aecm->noiseEst[i] >> 19) > 0) {
andrew@webrtc.orge03cafa2013-11-11 20:10:01613 // Avoid overflow.
614 // Multiplication with 2049 will cause wrap around. Scale
615 // down first and then multiply
616 aecm->noiseEst[i] >>= 11;
617 aecm->noiseEst[i] *= 2049;
Yves Gerey665174f2018-06-19 13:03:05618 } else if ((aecm->noiseEst[i] >> 11) > 0) {
andrew@webrtc.orge03cafa2013-11-11 20:10:01619 // Large enough for relative increase
620 aecm->noiseEst[i] *= 2049;
621 aecm->noiseEst[i] >>= 11;
Yves Gerey665174f2018-06-19 13:03:05622 } else {
andrew@webrtc.orge03cafa2013-11-11 20:10:01623 // Make incremental increases based on size every
624 // |kNoiseEstIncCount| block
625 aecm->noiseEstTooLowCtr[i]++;
Yves Gerey665174f2018-06-19 13:03:05626 if (aecm->noiseEstTooLowCtr[i] >= kNoiseEstIncCount) {
andrew@webrtc.orge03cafa2013-11-11 20:10:01627 aecm->noiseEst[i] += (aecm->noiseEst[i] >> 9) + 1;
Yves Gerey665174f2018-06-19 13:03:05628 aecm->noiseEstTooLowCtr[i] = 0; // Reset counter
andrew@webrtc.orge03cafa2013-11-11 20:10:01629 }
630 }
631 }
632 }
633
Yves Gerey665174f2018-06-19 13:03:05634 for (i = 0; i < PART_LEN1; i++) {
bjornv@webrtc.orgf87c0af2014-10-15 12:51:23635 tmp32 = aecm->noiseEst[i] >> shiftFromNearToNoise;
Yves Gerey665174f2018-06-19 13:03:05636 if (tmp32 > 32767) {
andrew@webrtc.orge03cafa2013-11-11 20:10:01637 tmp32 = 32767;
bjornv@webrtc.org750423c2014-09-30 09:26:36638 aecm->noiseEst[i] = tmp32 << shiftFromNearToNoise;
andrew@webrtc.orge03cafa2013-11-11 20:10:01639 }
640 noiseRShift16[i] = (int16_t)tmp32;
641
642 tmp16 = ONE_Q14 - lambda[i];
bjornv@webrtc.orgb38b0092015-03-10 06:40:02643 noiseRShift16[i] = (int16_t)((tmp16 * noiseRShift16[i]) >> 14);
andrew@webrtc.orge03cafa2013-11-11 20:10:01644 }
645
646 // Generate a uniform random array on [0 2^15-1].
647 WebRtcSpl_RandUArray(randW16, PART_LEN, &aecm->seed);
648
649 // Generate noise according to estimated energy.
Yves Gerey665174f2018-06-19 13:03:05650 uReal[0] = 0; // Reject LF noise.
andrew@webrtc.orge03cafa2013-11-11 20:10:01651 uImag[0] = 0;
Yves Gerey665174f2018-06-19 13:03:05652 for (i = 1; i < PART_LEN1; i++) {
andrew@webrtc.orge03cafa2013-11-11 20:10:01653 // Get a random index for the cos and sin tables over [0 359].
bjornv@webrtc.orgb38b0092015-03-10 06:40:02654 tmp16 = (int16_t)((359 * randW16[i - 1]) >> 15);
andrew@webrtc.orge03cafa2013-11-11 20:10:01655
656 // Tables are in Q13.
Yves Gerey665174f2018-06-19 13:03:05657 uReal[i] =
658 (int16_t)((noiseRShift16[i] * WebRtcAecm_kCosTable[tmp16]) >> 13);
659 uImag[i] =
660 (int16_t)((-noiseRShift16[i] * WebRtcAecm_kSinTable[tmp16]) >> 13);
andrew@webrtc.orge03cafa2013-11-11 20:10:01661 }
662 uImag[PART_LEN] = 0;
663
Yves Gerey665174f2018-06-19 13:03:05664 for (i = 0; i < PART_LEN1; i++) {
bjornv@webrtc.org6e71d172014-08-25 07:44:52665 out[i].real = WebRtcSpl_AddSatW16(out[i].real, uReal[i]);
666 out[i].imag = WebRtcSpl_AddSatW16(out[i].imag, uImag[i]);
andrew@webrtc.orge03cafa2013-11-11 20:10:01667 }
668}