1 /* 2 * Copyright (c) 2012 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_coding/neteq/expand.h" 12 13 #include <assert.h> 14 #include <string.h> // memset 15 16 #include <algorithm> // min, max 17 #include <limits> // numeric_limits<T> 18 19 #include "webrtc/base/safe_conversions.h" 20 #include "webrtc/common_audio/signal_processing/include/signal_processing_library.h" 21 #include "webrtc/modules/audio_coding/neteq/background_noise.h" 22 #include "webrtc/modules/audio_coding/neteq/dsp_helper.h" 23 #include "webrtc/modules/audio_coding/neteq/random_vector.h" 24 #include "webrtc/modules/audio_coding/neteq/statistics_calculator.h" 25 #include "webrtc/modules/audio_coding/neteq/sync_buffer.h" 26 27 namespace webrtc { 28 29 Expand::Expand(BackgroundNoise* background_noise, 30 SyncBuffer* sync_buffer, 31 RandomVector* random_vector, 32 StatisticsCalculator* statistics, 33 int fs, 34 size_t num_channels) 35 : random_vector_(random_vector), 36 sync_buffer_(sync_buffer), 37 first_expand_(true), 38 fs_hz_(fs), 39 num_channels_(num_channels), 40 consecutive_expands_(0), 41 background_noise_(background_noise), 42 statistics_(statistics), 43 overlap_length_(5 * fs / 8000), 44 lag_index_direction_(0), 45 current_lag_index_(0), 46 stop_muting_(false), 47 expand_duration_samples_(0), 48 channel_parameters_(new ChannelParameters[num_channels_]) { 49 assert(fs == 8000 || fs == 16000 || fs == 32000 || fs == 48000); 50 assert(fs <= static_cast<int>(kMaxSampleRate)); // Should not be possible. 51 assert(num_channels_ > 0); 52 memset(expand_lags_, 0, sizeof(expand_lags_)); 53 Reset(); 54 } 55 56 Expand::~Expand() = default; 57 58 void Expand::Reset() { 59 first_expand_ = true; 60 consecutive_expands_ = 0; 61 max_lag_ = 0; 62 for (size_t ix = 0; ix < num_channels_; ++ix) { 63 channel_parameters_[ix].expand_vector0.Clear(); 64 channel_parameters_[ix].expand_vector1.Clear(); 65 } 66 } 67 68 int Expand::Process(AudioMultiVector* output) { 69 int16_t random_vector[kMaxSampleRate / 8000 * 120 + 30]; 70 int16_t scaled_random_vector[kMaxSampleRate / 8000 * 125]; 71 static const int kTempDataSize = 3600; 72 int16_t temp_data[kTempDataSize]; // TODO(hlundin) Remove this. 73 int16_t* voiced_vector_storage = temp_data; 74 int16_t* voiced_vector = &voiced_vector_storage[overlap_length_]; 75 static const size_t kNoiseLpcOrder = BackgroundNoise::kMaxLpcOrder; 76 int16_t unvoiced_array_memory[kNoiseLpcOrder + kMaxSampleRate / 8000 * 125]; 77 int16_t* unvoiced_vector = unvoiced_array_memory + kUnvoicedLpcOrder; 78 int16_t* noise_vector = unvoiced_array_memory + kNoiseLpcOrder; 79 80 int fs_mult = fs_hz_ / 8000; 81 82 if (first_expand_) { 83 // Perform initial setup if this is the first expansion since last reset. 84 AnalyzeSignal(random_vector); 85 first_expand_ = false; 86 expand_duration_samples_ = 0; 87 } else { 88 // This is not the first expansion, parameters are already estimated. 89 // Extract a noise segment. 90 size_t rand_length = max_lag_; 91 // This only applies to SWB where length could be larger than 256. 92 assert(rand_length <= kMaxSampleRate / 8000 * 120 + 30); 93 GenerateRandomVector(2, rand_length, random_vector); 94 } 95 96 97 // Generate signal. 98 UpdateLagIndex(); 99 100 // Voiced part. 101 // Generate a weighted vector with the current lag. 102 size_t expansion_vector_length = max_lag_ + overlap_length_; 103 size_t current_lag = expand_lags_[current_lag_index_]; 104 // Copy lag+overlap data. 105 size_t expansion_vector_position = expansion_vector_length - current_lag - 106 overlap_length_; 107 size_t temp_length = current_lag + overlap_length_; 108 for (size_t channel_ix = 0; channel_ix < num_channels_; ++channel_ix) { 109 ChannelParameters& parameters = channel_parameters_[channel_ix]; 110 if (current_lag_index_ == 0) { 111 // Use only expand_vector0. 112 assert(expansion_vector_position + temp_length <= 113 parameters.expand_vector0.Size()); 114 memcpy(voiced_vector_storage, 115 ¶meters.expand_vector0[expansion_vector_position], 116 sizeof(int16_t) * temp_length); 117 } else if (current_lag_index_ == 1) { 118 // Mix 3/4 of expand_vector0 with 1/4 of expand_vector1. 119 WebRtcSpl_ScaleAndAddVectorsWithRound( 120 ¶meters.expand_vector0[expansion_vector_position], 3, 121 ¶meters.expand_vector1[expansion_vector_position], 1, 2, 122 voiced_vector_storage, temp_length); 123 } else if (current_lag_index_ == 2) { 124 // Mix 1/2 of expand_vector0 with 1/2 of expand_vector1. 125 assert(expansion_vector_position + temp_length <= 126 parameters.expand_vector0.Size()); 127 assert(expansion_vector_position + temp_length <= 128 parameters.expand_vector1.Size()); 129 WebRtcSpl_ScaleAndAddVectorsWithRound( 130 ¶meters.expand_vector0[expansion_vector_position], 1, 131 ¶meters.expand_vector1[expansion_vector_position], 1, 1, 132 voiced_vector_storage, temp_length); 133 } 134 135 // Get tapering window parameters. Values are in Q15. 136 int16_t muting_window, muting_window_increment; 137 int16_t unmuting_window, unmuting_window_increment; 138 if (fs_hz_ == 8000) { 139 muting_window = DspHelper::kMuteFactorStart8kHz; 140 muting_window_increment = DspHelper::kMuteFactorIncrement8kHz; 141 unmuting_window = DspHelper::kUnmuteFactorStart8kHz; 142 unmuting_window_increment = DspHelper::kUnmuteFactorIncrement8kHz; 143 } else if (fs_hz_ == 16000) { 144 muting_window = DspHelper::kMuteFactorStart16kHz; 145 muting_window_increment = DspHelper::kMuteFactorIncrement16kHz; 146 unmuting_window = DspHelper::kUnmuteFactorStart16kHz; 147 unmuting_window_increment = DspHelper::kUnmuteFactorIncrement16kHz; 148 } else if (fs_hz_ == 32000) { 149 muting_window = DspHelper::kMuteFactorStart32kHz; 150 muting_window_increment = DspHelper::kMuteFactorIncrement32kHz; 151 unmuting_window = DspHelper::kUnmuteFactorStart32kHz; 152 unmuting_window_increment = DspHelper::kUnmuteFactorIncrement32kHz; 153 } else { // fs_ == 48000 154 muting_window = DspHelper::kMuteFactorStart48kHz; 155 muting_window_increment = DspHelper::kMuteFactorIncrement48kHz; 156 unmuting_window = DspHelper::kUnmuteFactorStart48kHz; 157 unmuting_window_increment = DspHelper::kUnmuteFactorIncrement48kHz; 158 } 159 160 // Smooth the expanded if it has not been muted to a low amplitude and 161 // |current_voice_mix_factor| is larger than 0.5. 162 if ((parameters.mute_factor > 819) && 163 (parameters.current_voice_mix_factor > 8192)) { 164 size_t start_ix = sync_buffer_->Size() - overlap_length_; 165 for (size_t i = 0; i < overlap_length_; i++) { 166 // Do overlap add between new vector and overlap. 167 (*sync_buffer_)[channel_ix][start_ix + i] = 168 (((*sync_buffer_)[channel_ix][start_ix + i] * muting_window) + 169 (((parameters.mute_factor * voiced_vector_storage[i]) >> 14) * 170 unmuting_window) + 16384) >> 15; 171 muting_window += muting_window_increment; 172 unmuting_window += unmuting_window_increment; 173 } 174 } else if (parameters.mute_factor == 0) { 175 // The expanded signal will consist of only comfort noise if 176 // mute_factor = 0. Set the output length to 15 ms for best noise 177 // production. 178 // TODO(hlundin): This has been disabled since the length of 179 // parameters.expand_vector0 and parameters.expand_vector1 no longer 180 // match with expand_lags_, causing invalid reads and writes. Is it a good 181 // idea to enable this again, and solve the vector size problem? 182 // max_lag_ = fs_mult * 120; 183 // expand_lags_[0] = fs_mult * 120; 184 // expand_lags_[1] = fs_mult * 120; 185 // expand_lags_[2] = fs_mult * 120; 186 } 187 188 // Unvoiced part. 189 // Filter |scaled_random_vector| through |ar_filter_|. 190 memcpy(unvoiced_vector - kUnvoicedLpcOrder, parameters.ar_filter_state, 191 sizeof(int16_t) * kUnvoicedLpcOrder); 192 int32_t add_constant = 0; 193 if (parameters.ar_gain_scale > 0) { 194 add_constant = 1 << (parameters.ar_gain_scale - 1); 195 } 196 WebRtcSpl_AffineTransformVector(scaled_random_vector, random_vector, 197 parameters.ar_gain, add_constant, 198 parameters.ar_gain_scale, 199 current_lag); 200 WebRtcSpl_FilterARFastQ12(scaled_random_vector, unvoiced_vector, 201 parameters.ar_filter, kUnvoicedLpcOrder + 1, 202 current_lag); 203 memcpy(parameters.ar_filter_state, 204 &(unvoiced_vector[current_lag - kUnvoicedLpcOrder]), 205 sizeof(int16_t) * kUnvoicedLpcOrder); 206 207 // Combine voiced and unvoiced contributions. 208 209 // Set a suitable cross-fading slope. 210 // For lag = 211 // <= 31 * fs_mult => go from 1 to 0 in about 8 ms; 212 // (>= 31 .. <= 63) * fs_mult => go from 1 to 0 in about 16 ms; 213 // >= 64 * fs_mult => go from 1 to 0 in about 32 ms. 214 // temp_shift = getbits(max_lag_) - 5. 215 int temp_shift = 216 (31 - WebRtcSpl_NormW32(rtc::checked_cast<int32_t>(max_lag_))) - 5; 217 int16_t mix_factor_increment = 256 >> temp_shift; 218 if (stop_muting_) { 219 mix_factor_increment = 0; 220 } 221 222 // Create combined signal by shifting in more and more of unvoiced part. 223 temp_shift = 8 - temp_shift; // = getbits(mix_factor_increment). 224 size_t temp_length = (parameters.current_voice_mix_factor - 225 parameters.voice_mix_factor) >> temp_shift; 226 temp_length = std::min(temp_length, current_lag); 227 DspHelper::CrossFade(voiced_vector, unvoiced_vector, temp_length, 228 ¶meters.current_voice_mix_factor, 229 mix_factor_increment, temp_data); 230 231 // End of cross-fading period was reached before end of expanded signal 232 // path. Mix the rest with a fixed mixing factor. 233 if (temp_length < current_lag) { 234 if (mix_factor_increment != 0) { 235 parameters.current_voice_mix_factor = parameters.voice_mix_factor; 236 } 237 int16_t temp_scale = 16384 - parameters.current_voice_mix_factor; 238 WebRtcSpl_ScaleAndAddVectorsWithRound( 239 voiced_vector + temp_length, parameters.current_voice_mix_factor, 240 unvoiced_vector + temp_length, temp_scale, 14, 241 temp_data + temp_length, current_lag - temp_length); 242 } 243 244 // Select muting slope depending on how many consecutive expands we have 245 // done. 246 if (consecutive_expands_ == 3) { 247 // Let the mute factor decrease from 1.0 to 0.95 in 6.25 ms. 248 // mute_slope = 0.0010 / fs_mult in Q20. 249 parameters.mute_slope = std::max(parameters.mute_slope, 1049 / fs_mult); 250 } 251 if (consecutive_expands_ == 7) { 252 // Let the mute factor decrease from 1.0 to 0.90 in 6.25 ms. 253 // mute_slope = 0.0020 / fs_mult in Q20. 254 parameters.mute_slope = std::max(parameters.mute_slope, 2097 / fs_mult); 255 } 256 257 // Mute segment according to slope value. 258 if ((consecutive_expands_ != 0) || !parameters.onset) { 259 // Mute to the previous level, then continue with the muting. 260 WebRtcSpl_AffineTransformVector(temp_data, temp_data, 261 parameters.mute_factor, 8192, 262 14, current_lag); 263 264 if (!stop_muting_) { 265 DspHelper::MuteSignal(temp_data, parameters.mute_slope, current_lag); 266 267 // Shift by 6 to go from Q20 to Q14. 268 // TODO(hlundin): Adding 8192 before shifting 6 steps seems wrong. 269 // Legacy. 270 int16_t gain = static_cast<int16_t>(16384 - 271 (((current_lag * parameters.mute_slope) + 8192) >> 6)); 272 gain = ((gain * parameters.mute_factor) + 8192) >> 14; 273 274 // Guard against getting stuck with very small (but sometimes audible) 275 // gain. 276 if ((consecutive_expands_ > 3) && (gain >= parameters.mute_factor)) { 277 parameters.mute_factor = 0; 278 } else { 279 parameters.mute_factor = gain; 280 } 281 } 282 } 283 284 // Background noise part. 285 GenerateBackgroundNoise(random_vector, 286 channel_ix, 287 channel_parameters_[channel_ix].mute_slope, 288 TooManyExpands(), 289 current_lag, 290 unvoiced_array_memory); 291 292 // Add background noise to the combined voiced-unvoiced signal. 293 for (size_t i = 0; i < current_lag; i++) { 294 temp_data[i] = temp_data[i] + noise_vector[i]; 295 } 296 if (channel_ix == 0) { 297 output->AssertSize(current_lag); 298 } else { 299 assert(output->Size() == current_lag); 300 } 301 memcpy(&(*output)[channel_ix][0], temp_data, 302 sizeof(temp_data[0]) * current_lag); 303 } 304 305 // Increase call number and cap it. 306 consecutive_expands_ = consecutive_expands_ >= kMaxConsecutiveExpands ? 307 kMaxConsecutiveExpands : consecutive_expands_ + 1; 308 expand_duration_samples_ += output->Size(); 309 // Clamp the duration counter at 2 seconds. 310 expand_duration_samples_ = 311 std::min(expand_duration_samples_, rtc::checked_cast<size_t>(fs_hz_ * 2)); 312 return 0; 313 } 314 315 void Expand::SetParametersForNormalAfterExpand() { 316 current_lag_index_ = 0; 317 lag_index_direction_ = 0; 318 stop_muting_ = true; // Do not mute signal any more. 319 statistics_->LogDelayedPacketOutageEvent( 320 rtc::checked_cast<int>(expand_duration_samples_) / (fs_hz_ / 1000)); 321 } 322 323 void Expand::SetParametersForMergeAfterExpand() { 324 current_lag_index_ = -1; /* out of the 3 possible ones */ 325 lag_index_direction_ = 1; /* make sure we get the "optimal" lag */ 326 stop_muting_ = true; 327 } 328 329 size_t Expand::overlap_length() const { 330 return overlap_length_; 331 } 332 333 void Expand::InitializeForAnExpandPeriod() { 334 lag_index_direction_ = 1; 335 current_lag_index_ = -1; 336 stop_muting_ = false; 337 random_vector_->set_seed_increment(1); 338 consecutive_expands_ = 0; 339 for (size_t ix = 0; ix < num_channels_; ++ix) { 340 channel_parameters_[ix].current_voice_mix_factor = 16384; // 1.0 in Q14. 341 channel_parameters_[ix].mute_factor = 16384; // 1.0 in Q14. 342 // Start with 0 gain for background noise. 343 background_noise_->SetMuteFactor(ix, 0); 344 } 345 } 346 347 bool Expand::TooManyExpands() { 348 return consecutive_expands_ >= kMaxConsecutiveExpands; 349 } 350 351 void Expand::AnalyzeSignal(int16_t* random_vector) { 352 int32_t auto_correlation[kUnvoicedLpcOrder + 1]; 353 int16_t reflection_coeff[kUnvoicedLpcOrder]; 354 int16_t correlation_vector[kMaxSampleRate / 8000 * 102]; 355 size_t best_correlation_index[kNumCorrelationCandidates]; 356 int16_t best_correlation[kNumCorrelationCandidates]; 357 size_t best_distortion_index[kNumCorrelationCandidates]; 358 int16_t best_distortion[kNumCorrelationCandidates]; 359 int32_t correlation_vector2[(99 * kMaxSampleRate / 8000) + 1]; 360 int32_t best_distortion_w32[kNumCorrelationCandidates]; 361 static const size_t kNoiseLpcOrder = BackgroundNoise::kMaxLpcOrder; 362 int16_t unvoiced_array_memory[kNoiseLpcOrder + kMaxSampleRate / 8000 * 125]; 363 int16_t* unvoiced_vector = unvoiced_array_memory + kUnvoicedLpcOrder; 364 365 int fs_mult = fs_hz_ / 8000; 366 367 // Pre-calculate common multiplications with fs_mult. 368 size_t fs_mult_4 = static_cast<size_t>(fs_mult * 4); 369 size_t fs_mult_20 = static_cast<size_t>(fs_mult * 20); 370 size_t fs_mult_120 = static_cast<size_t>(fs_mult * 120); 371 size_t fs_mult_dist_len = fs_mult * kDistortionLength; 372 size_t fs_mult_lpc_analysis_len = fs_mult * kLpcAnalysisLength; 373 374 const size_t signal_length = static_cast<size_t>(256 * fs_mult); 375 const int16_t* audio_history = 376 &(*sync_buffer_)[0][sync_buffer_->Size() - signal_length]; 377 378 // Initialize. 379 InitializeForAnExpandPeriod(); 380 381 // Calculate correlation in downsampled domain (4 kHz sample rate). 382 int correlation_scale; 383 size_t correlation_length = 51; // TODO(hlundin): Legacy bit-exactness. 384 // If it is decided to break bit-exactness |correlation_length| should be 385 // initialized to the return value of Correlation(). 386 Correlation(audio_history, signal_length, correlation_vector, 387 &correlation_scale); 388 389 // Find peaks in correlation vector. 390 DspHelper::PeakDetection(correlation_vector, correlation_length, 391 kNumCorrelationCandidates, fs_mult, 392 best_correlation_index, best_correlation); 393 394 // Adjust peak locations; cross-correlation lags start at 2.5 ms 395 // (20 * fs_mult samples). 396 best_correlation_index[0] += fs_mult_20; 397 best_correlation_index[1] += fs_mult_20; 398 best_correlation_index[2] += fs_mult_20; 399 400 // Calculate distortion around the |kNumCorrelationCandidates| best lags. 401 int distortion_scale = 0; 402 for (size_t i = 0; i < kNumCorrelationCandidates; i++) { 403 size_t min_index = std::max(fs_mult_20, 404 best_correlation_index[i] - fs_mult_4); 405 size_t max_index = std::min(fs_mult_120 - 1, 406 best_correlation_index[i] + fs_mult_4); 407 best_distortion_index[i] = DspHelper::MinDistortion( 408 &(audio_history[signal_length - fs_mult_dist_len]), min_index, 409 max_index, fs_mult_dist_len, &best_distortion_w32[i]); 410 distortion_scale = std::max(16 - WebRtcSpl_NormW32(best_distortion_w32[i]), 411 distortion_scale); 412 } 413 // Shift the distortion values to fit in 16 bits. 414 WebRtcSpl_VectorBitShiftW32ToW16(best_distortion, kNumCorrelationCandidates, 415 best_distortion_w32, distortion_scale); 416 417 // Find the maximizing index |i| of the cost function 418 // f[i] = best_correlation[i] / best_distortion[i]. 419 int32_t best_ratio = std::numeric_limits<int32_t>::min(); 420 size_t best_index = std::numeric_limits<size_t>::max(); 421 for (size_t i = 0; i < kNumCorrelationCandidates; ++i) { 422 int32_t ratio; 423 if (best_distortion[i] > 0) { 424 ratio = (best_correlation[i] << 16) / best_distortion[i]; 425 } else if (best_correlation[i] == 0) { 426 ratio = 0; // No correlation set result to zero. 427 } else { 428 ratio = std::numeric_limits<int32_t>::max(); // Denominator is zero. 429 } 430 if (ratio > best_ratio) { 431 best_index = i; 432 best_ratio = ratio; 433 } 434 } 435 436 size_t distortion_lag = best_distortion_index[best_index]; 437 size_t correlation_lag = best_correlation_index[best_index]; 438 max_lag_ = std::max(distortion_lag, correlation_lag); 439 440 // Calculate the exact best correlation in the range between 441 // |correlation_lag| and |distortion_lag|. 442 correlation_length = 443 std::max(std::min(distortion_lag + 10, fs_mult_120), 444 static_cast<size_t>(60 * fs_mult)); 445 446 size_t start_index = std::min(distortion_lag, correlation_lag); 447 size_t correlation_lags = static_cast<size_t>( 448 WEBRTC_SPL_ABS_W16((distortion_lag-correlation_lag)) + 1); 449 assert(correlation_lags <= static_cast<size_t>(99 * fs_mult + 1)); 450 451 for (size_t channel_ix = 0; channel_ix < num_channels_; ++channel_ix) { 452 ChannelParameters& parameters = channel_parameters_[channel_ix]; 453 // Calculate suitable scaling. 454 int16_t signal_max = WebRtcSpl_MaxAbsValueW16( 455 &audio_history[signal_length - correlation_length - start_index 456 - correlation_lags], 457 correlation_length + start_index + correlation_lags - 1); 458 correlation_scale = (31 - WebRtcSpl_NormW32(signal_max * signal_max)) + 459 (31 - WebRtcSpl_NormW32(static_cast<int32_t>(correlation_length))) - 31; 460 correlation_scale = std::max(0, correlation_scale); 461 462 // Calculate the correlation, store in |correlation_vector2|. 463 WebRtcSpl_CrossCorrelation( 464 correlation_vector2, 465 &(audio_history[signal_length - correlation_length]), 466 &(audio_history[signal_length - correlation_length - start_index]), 467 correlation_length, correlation_lags, correlation_scale, -1); 468 469 // Find maximizing index. 470 best_index = WebRtcSpl_MaxIndexW32(correlation_vector2, correlation_lags); 471 int32_t max_correlation = correlation_vector2[best_index]; 472 // Compensate index with start offset. 473 best_index = best_index + start_index; 474 475 // Calculate energies. 476 int32_t energy1 = WebRtcSpl_DotProductWithScale( 477 &(audio_history[signal_length - correlation_length]), 478 &(audio_history[signal_length - correlation_length]), 479 correlation_length, correlation_scale); 480 int32_t energy2 = WebRtcSpl_DotProductWithScale( 481 &(audio_history[signal_length - correlation_length - best_index]), 482 &(audio_history[signal_length - correlation_length - best_index]), 483 correlation_length, correlation_scale); 484 485 // Calculate the correlation coefficient between the two portions of the 486 // signal. 487 int32_t corr_coefficient; 488 if ((energy1 > 0) && (energy2 > 0)) { 489 int energy1_scale = std::max(16 - WebRtcSpl_NormW32(energy1), 0); 490 int energy2_scale = std::max(16 - WebRtcSpl_NormW32(energy2), 0); 491 // Make sure total scaling is even (to simplify scale factor after sqrt). 492 if ((energy1_scale + energy2_scale) & 1) { 493 // If sum is odd, add 1 to make it even. 494 energy1_scale += 1; 495 } 496 int32_t scaled_energy1 = energy1 >> energy1_scale; 497 int32_t scaled_energy2 = energy2 >> energy2_scale; 498 int16_t sqrt_energy_product = static_cast<int16_t>( 499 WebRtcSpl_SqrtFloor(scaled_energy1 * scaled_energy2)); 500 // Calculate max_correlation / sqrt(energy1 * energy2) in Q14. 501 int cc_shift = 14 - (energy1_scale + energy2_scale) / 2; 502 max_correlation = WEBRTC_SPL_SHIFT_W32(max_correlation, cc_shift); 503 corr_coefficient = WebRtcSpl_DivW32W16(max_correlation, 504 sqrt_energy_product); 505 // Cap at 1.0 in Q14. 506 corr_coefficient = std::min(16384, corr_coefficient); 507 } else { 508 corr_coefficient = 0; 509 } 510 511 // Extract the two vectors expand_vector0 and expand_vector1 from 512 // |audio_history|. 513 size_t expansion_length = max_lag_ + overlap_length_; 514 const int16_t* vector1 = &(audio_history[signal_length - expansion_length]); 515 const int16_t* vector2 = vector1 - distortion_lag; 516 // Normalize the second vector to the same energy as the first. 517 energy1 = WebRtcSpl_DotProductWithScale(vector1, vector1, expansion_length, 518 correlation_scale); 519 energy2 = WebRtcSpl_DotProductWithScale(vector2, vector2, expansion_length, 520 correlation_scale); 521 // Confirm that amplitude ratio sqrt(energy1 / energy2) is within 0.5 - 2.0, 522 // i.e., energy1 / energy2 is within 0.25 - 4. 523 int16_t amplitude_ratio; 524 if ((energy1 / 4 < energy2) && (energy1 > energy2 / 4)) { 525 // Energy constraint fulfilled. Use both vectors and scale them 526 // accordingly. 527 int32_t scaled_energy2 = std::max(16 - WebRtcSpl_NormW32(energy2), 0); 528 int32_t scaled_energy1 = scaled_energy2 - 13; 529 // Calculate scaled_energy1 / scaled_energy2 in Q13. 530 int32_t energy_ratio = WebRtcSpl_DivW32W16( 531 WEBRTC_SPL_SHIFT_W32(energy1, -scaled_energy1), 532 static_cast<int16_t>(energy2 >> scaled_energy2)); 533 // Calculate sqrt ratio in Q13 (sqrt of en1/en2 in Q26). 534 amplitude_ratio = 535 static_cast<int16_t>(WebRtcSpl_SqrtFloor(energy_ratio << 13)); 536 // Copy the two vectors and give them the same energy. 537 parameters.expand_vector0.Clear(); 538 parameters.expand_vector0.PushBack(vector1, expansion_length); 539 parameters.expand_vector1.Clear(); 540 if (parameters.expand_vector1.Size() < expansion_length) { 541 parameters.expand_vector1.Extend( 542 expansion_length - parameters.expand_vector1.Size()); 543 } 544 WebRtcSpl_AffineTransformVector(¶meters.expand_vector1[0], 545 const_cast<int16_t*>(vector2), 546 amplitude_ratio, 547 4096, 548 13, 549 expansion_length); 550 } else { 551 // Energy change constraint not fulfilled. Only use last vector. 552 parameters.expand_vector0.Clear(); 553 parameters.expand_vector0.PushBack(vector1, expansion_length); 554 // Copy from expand_vector0 to expand_vector1. 555 parameters.expand_vector0.CopyTo(¶meters.expand_vector1); 556 // Set the energy_ratio since it is used by muting slope. 557 if ((energy1 / 4 < energy2) || (energy2 == 0)) { 558 amplitude_ratio = 4096; // 0.5 in Q13. 559 } else { 560 amplitude_ratio = 16384; // 2.0 in Q13. 561 } 562 } 563 564 // Set the 3 lag values. 565 if (distortion_lag == correlation_lag) { 566 expand_lags_[0] = distortion_lag; 567 expand_lags_[1] = distortion_lag; 568 expand_lags_[2] = distortion_lag; 569 } else { 570 // |distortion_lag| and |correlation_lag| are not equal; use different 571 // combinations of the two. 572 // First lag is |distortion_lag| only. 573 expand_lags_[0] = distortion_lag; 574 // Second lag is the average of the two. 575 expand_lags_[1] = (distortion_lag + correlation_lag) / 2; 576 // Third lag is the average again, but rounding towards |correlation_lag|. 577 if (distortion_lag > correlation_lag) { 578 expand_lags_[2] = (distortion_lag + correlation_lag - 1) / 2; 579 } else { 580 expand_lags_[2] = (distortion_lag + correlation_lag + 1) / 2; 581 } 582 } 583 584 // Calculate the LPC and the gain of the filters. 585 // Calculate scale value needed for auto-correlation. 586 correlation_scale = WebRtcSpl_MaxAbsValueW16( 587 &(audio_history[signal_length - fs_mult_lpc_analysis_len]), 588 fs_mult_lpc_analysis_len); 589 590 correlation_scale = std::min(16 - WebRtcSpl_NormW32(correlation_scale), 0); 591 correlation_scale = std::max(correlation_scale * 2 + 7, 0); 592 593 // Calculate kUnvoicedLpcOrder + 1 lags of the auto-correlation function. 594 size_t temp_index = signal_length - fs_mult_lpc_analysis_len - 595 kUnvoicedLpcOrder; 596 // Copy signal to temporary vector to be able to pad with leading zeros. 597 int16_t* temp_signal = new int16_t[fs_mult_lpc_analysis_len 598 + kUnvoicedLpcOrder]; 599 memset(temp_signal, 0, 600 sizeof(int16_t) * (fs_mult_lpc_analysis_len + kUnvoicedLpcOrder)); 601 memcpy(&temp_signal[kUnvoicedLpcOrder], 602 &audio_history[temp_index + kUnvoicedLpcOrder], 603 sizeof(int16_t) * fs_mult_lpc_analysis_len); 604 WebRtcSpl_CrossCorrelation(auto_correlation, 605 &temp_signal[kUnvoicedLpcOrder], 606 &temp_signal[kUnvoicedLpcOrder], 607 fs_mult_lpc_analysis_len, kUnvoicedLpcOrder + 1, 608 correlation_scale, -1); 609 delete [] temp_signal; 610 611 // Verify that variance is positive. 612 if (auto_correlation[0] > 0) { 613 // Estimate AR filter parameters using Levinson-Durbin algorithm; 614 // kUnvoicedLpcOrder + 1 filter coefficients. 615 int16_t stability = WebRtcSpl_LevinsonDurbin(auto_correlation, 616 parameters.ar_filter, 617 reflection_coeff, 618 kUnvoicedLpcOrder); 619 620 // Keep filter parameters only if filter is stable. 621 if (stability != 1) { 622 // Set first coefficient to 4096 (1.0 in Q12). 623 parameters.ar_filter[0] = 4096; 624 // Set remaining |kUnvoicedLpcOrder| coefficients to zero. 625 WebRtcSpl_MemSetW16(parameters.ar_filter + 1, 0, kUnvoicedLpcOrder); 626 } 627 } 628 629 if (channel_ix == 0) { 630 // Extract a noise segment. 631 size_t noise_length; 632 if (distortion_lag < 40) { 633 noise_length = 2 * distortion_lag + 30; 634 } else { 635 noise_length = distortion_lag + 30; 636 } 637 if (noise_length <= RandomVector::kRandomTableSize) { 638 memcpy(random_vector, RandomVector::kRandomTable, 639 sizeof(int16_t) * noise_length); 640 } else { 641 // Only applies to SWB where length could be larger than 642 // |kRandomTableSize|. 643 memcpy(random_vector, RandomVector::kRandomTable, 644 sizeof(int16_t) * RandomVector::kRandomTableSize); 645 assert(noise_length <= kMaxSampleRate / 8000 * 120 + 30); 646 random_vector_->IncreaseSeedIncrement(2); 647 random_vector_->Generate( 648 noise_length - RandomVector::kRandomTableSize, 649 &random_vector[RandomVector::kRandomTableSize]); 650 } 651 } 652 653 // Set up state vector and calculate scale factor for unvoiced filtering. 654 memcpy(parameters.ar_filter_state, 655 &(audio_history[signal_length - kUnvoicedLpcOrder]), 656 sizeof(int16_t) * kUnvoicedLpcOrder); 657 memcpy(unvoiced_vector - kUnvoicedLpcOrder, 658 &(audio_history[signal_length - 128 - kUnvoicedLpcOrder]), 659 sizeof(int16_t) * kUnvoicedLpcOrder); 660 WebRtcSpl_FilterMAFastQ12(&audio_history[signal_length - 128], 661 unvoiced_vector, 662 parameters.ar_filter, 663 kUnvoicedLpcOrder + 1, 664 128); 665 int16_t unvoiced_prescale; 666 if (WebRtcSpl_MaxAbsValueW16(unvoiced_vector, 128) > 4000) { 667 unvoiced_prescale = 4; 668 } else { 669 unvoiced_prescale = 0; 670 } 671 int32_t unvoiced_energy = WebRtcSpl_DotProductWithScale(unvoiced_vector, 672 unvoiced_vector, 673 128, 674 unvoiced_prescale); 675 676 // Normalize |unvoiced_energy| to 28 or 29 bits to preserve sqrt() accuracy. 677 int16_t unvoiced_scale = WebRtcSpl_NormW32(unvoiced_energy) - 3; 678 // Make sure we do an odd number of shifts since we already have 7 shifts 679 // from dividing with 128 earlier. This will make the total scale factor 680 // even, which is suitable for the sqrt. 681 unvoiced_scale += ((unvoiced_scale & 0x1) ^ 0x1); 682 unvoiced_energy = WEBRTC_SPL_SHIFT_W32(unvoiced_energy, unvoiced_scale); 683 int16_t unvoiced_gain = 684 static_cast<int16_t>(WebRtcSpl_SqrtFloor(unvoiced_energy)); 685 parameters.ar_gain_scale = 13 686 + (unvoiced_scale + 7 - unvoiced_prescale) / 2; 687 parameters.ar_gain = unvoiced_gain; 688 689 // Calculate voice_mix_factor from corr_coefficient. 690 // Let x = corr_coefficient. Then, we compute: 691 // if (x > 0.48) 692 // voice_mix_factor = (-5179 + 19931x - 16422x^2 + 5776x^3) / 4096; 693 // else 694 // voice_mix_factor = 0; 695 if (corr_coefficient > 7875) { 696 int16_t x1, x2, x3; 697 // |corr_coefficient| is in Q14. 698 x1 = static_cast<int16_t>(corr_coefficient); 699 x2 = (x1 * x1) >> 14; // Shift 14 to keep result in Q14. 700 x3 = (x1 * x2) >> 14; 701 static const int kCoefficients[4] = { -5179, 19931, -16422, 5776 }; 702 int32_t temp_sum = kCoefficients[0] << 14; 703 temp_sum += kCoefficients[1] * x1; 704 temp_sum += kCoefficients[2] * x2; 705 temp_sum += kCoefficients[3] * x3; 706 parameters.voice_mix_factor = 707 static_cast<int16_t>(std::min(temp_sum / 4096, 16384)); 708 parameters.voice_mix_factor = std::max(parameters.voice_mix_factor, 709 static_cast<int16_t>(0)); 710 } else { 711 parameters.voice_mix_factor = 0; 712 } 713 714 // Calculate muting slope. Reuse value from earlier scaling of 715 // |expand_vector0| and |expand_vector1|. 716 int16_t slope = amplitude_ratio; 717 if (slope > 12288) { 718 // slope > 1.5. 719 // Calculate (1 - (1 / slope)) / distortion_lag = 720 // (slope - 1) / (distortion_lag * slope). 721 // |slope| is in Q13, so 1 corresponds to 8192. Shift up to Q25 before 722 // the division. 723 // Shift the denominator from Q13 to Q5 before the division. The result of 724 // the division will then be in Q20. 725 int temp_ratio = WebRtcSpl_DivW32W16( 726 (slope - 8192) << 12, 727 static_cast<int16_t>((distortion_lag * slope) >> 8)); 728 if (slope > 14746) { 729 // slope > 1.8. 730 // Divide by 2, with proper rounding. 731 parameters.mute_slope = (temp_ratio + 1) / 2; 732 } else { 733 // Divide by 8, with proper rounding. 734 parameters.mute_slope = (temp_ratio + 4) / 8; 735 } 736 parameters.onset = true; 737 } else { 738 // Calculate (1 - slope) / distortion_lag. 739 // Shift |slope| by 7 to Q20 before the division. The result is in Q20. 740 parameters.mute_slope = WebRtcSpl_DivW32W16( 741 (8192 - slope) << 7, static_cast<int16_t>(distortion_lag)); 742 if (parameters.voice_mix_factor <= 13107) { 743 // Make sure the mute factor decreases from 1.0 to 0.9 in no more than 744 // 6.25 ms. 745 // mute_slope >= 0.005 / fs_mult in Q20. 746 parameters.mute_slope = std::max(5243 / fs_mult, parameters.mute_slope); 747 } else if (slope > 8028) { 748 parameters.mute_slope = 0; 749 } 750 parameters.onset = false; 751 } 752 } 753 } 754 755 Expand::ChannelParameters::ChannelParameters() 756 : mute_factor(16384), 757 ar_gain(0), 758 ar_gain_scale(0), 759 voice_mix_factor(0), 760 current_voice_mix_factor(0), 761 onset(false), 762 mute_slope(0) { 763 memset(ar_filter, 0, sizeof(ar_filter)); 764 memset(ar_filter_state, 0, sizeof(ar_filter_state)); 765 } 766 767 void Expand::Correlation(const int16_t* input, 768 size_t input_length, 769 int16_t* output, 770 int* output_scale) const { 771 // Set parameters depending on sample rate. 772 const int16_t* filter_coefficients; 773 size_t num_coefficients; 774 int16_t downsampling_factor; 775 if (fs_hz_ == 8000) { 776 num_coefficients = 3; 777 downsampling_factor = 2; 778 filter_coefficients = DspHelper::kDownsample8kHzTbl; 779 } else if (fs_hz_ == 16000) { 780 num_coefficients = 5; 781 downsampling_factor = 4; 782 filter_coefficients = DspHelper::kDownsample16kHzTbl; 783 } else if (fs_hz_ == 32000) { 784 num_coefficients = 7; 785 downsampling_factor = 8; 786 filter_coefficients = DspHelper::kDownsample32kHzTbl; 787 } else { // fs_hz_ == 48000. 788 num_coefficients = 7; 789 downsampling_factor = 12; 790 filter_coefficients = DspHelper::kDownsample48kHzTbl; 791 } 792 793 // Correlate from lag 10 to lag 60 in downsampled domain. 794 // (Corresponds to 20-120 for narrow-band, 40-240 for wide-band, and so on.) 795 static const size_t kCorrelationStartLag = 10; 796 static const size_t kNumCorrelationLags = 54; 797 static const size_t kCorrelationLength = 60; 798 // Downsample to 4 kHz sample rate. 799 static const size_t kDownsampledLength = kCorrelationStartLag 800 + kNumCorrelationLags + kCorrelationLength; 801 int16_t downsampled_input[kDownsampledLength]; 802 static const size_t kFilterDelay = 0; 803 WebRtcSpl_DownsampleFast( 804 input + input_length - kDownsampledLength * downsampling_factor, 805 kDownsampledLength * downsampling_factor, downsampled_input, 806 kDownsampledLength, filter_coefficients, num_coefficients, 807 downsampling_factor, kFilterDelay); 808 809 // Normalize |downsampled_input| to using all 16 bits. 810 int16_t max_value = WebRtcSpl_MaxAbsValueW16(downsampled_input, 811 kDownsampledLength); 812 int16_t norm_shift = 16 - WebRtcSpl_NormW32(max_value); 813 WebRtcSpl_VectorBitShiftW16(downsampled_input, kDownsampledLength, 814 downsampled_input, norm_shift); 815 816 int32_t correlation[kNumCorrelationLags]; 817 static const int kCorrelationShift = 6; 818 WebRtcSpl_CrossCorrelation( 819 correlation, 820 &downsampled_input[kDownsampledLength - kCorrelationLength], 821 &downsampled_input[kDownsampledLength - kCorrelationLength 822 - kCorrelationStartLag], 823 kCorrelationLength, kNumCorrelationLags, kCorrelationShift, -1); 824 825 // Normalize and move data from 32-bit to 16-bit vector. 826 int32_t max_correlation = WebRtcSpl_MaxAbsValueW32(correlation, 827 kNumCorrelationLags); 828 int16_t norm_shift2 = static_cast<int16_t>( 829 std::max(18 - WebRtcSpl_NormW32(max_correlation), 0)); 830 WebRtcSpl_VectorBitShiftW32ToW16(output, kNumCorrelationLags, correlation, 831 norm_shift2); 832 // Total scale factor (right shifts) of correlation value. 833 *output_scale = 2 * norm_shift + kCorrelationShift + norm_shift2; 834 } 835 836 void Expand::UpdateLagIndex() { 837 current_lag_index_ = current_lag_index_ + lag_index_direction_; 838 // Change direction if needed. 839 if (current_lag_index_ <= 0) { 840 lag_index_direction_ = 1; 841 } 842 if (current_lag_index_ >= kNumLags - 1) { 843 lag_index_direction_ = -1; 844 } 845 } 846 847 Expand* ExpandFactory::Create(BackgroundNoise* background_noise, 848 SyncBuffer* sync_buffer, 849 RandomVector* random_vector, 850 StatisticsCalculator* statistics, 851 int fs, 852 size_t num_channels) const { 853 return new Expand(background_noise, sync_buffer, random_vector, statistics, 854 fs, num_channels); 855 } 856 857 // TODO(turajs): This can be moved to BackgroundNoise class. 858 void Expand::GenerateBackgroundNoise(int16_t* random_vector, 859 size_t channel, 860 int mute_slope, 861 bool too_many_expands, 862 size_t num_noise_samples, 863 int16_t* buffer) { 864 static const size_t kNoiseLpcOrder = BackgroundNoise::kMaxLpcOrder; 865 int16_t scaled_random_vector[kMaxSampleRate / 8000 * 125]; 866 assert(num_noise_samples <= (kMaxSampleRate / 8000 * 125)); 867 int16_t* noise_samples = &buffer[kNoiseLpcOrder]; 868 if (background_noise_->initialized()) { 869 // Use background noise parameters. 870 memcpy(noise_samples - kNoiseLpcOrder, 871 background_noise_->FilterState(channel), 872 sizeof(int16_t) * kNoiseLpcOrder); 873 874 int dc_offset = 0; 875 if (background_noise_->ScaleShift(channel) > 1) { 876 dc_offset = 1 << (background_noise_->ScaleShift(channel) - 1); 877 } 878 879 // Scale random vector to correct energy level. 880 WebRtcSpl_AffineTransformVector( 881 scaled_random_vector, random_vector, 882 background_noise_->Scale(channel), dc_offset, 883 background_noise_->ScaleShift(channel), 884 num_noise_samples); 885 886 WebRtcSpl_FilterARFastQ12(scaled_random_vector, noise_samples, 887 background_noise_->Filter(channel), 888 kNoiseLpcOrder + 1, 889 num_noise_samples); 890 891 background_noise_->SetFilterState( 892 channel, 893 &(noise_samples[num_noise_samples - kNoiseLpcOrder]), 894 kNoiseLpcOrder); 895 896 // Unmute the background noise. 897 int16_t bgn_mute_factor = background_noise_->MuteFactor(channel); 898 NetEq::BackgroundNoiseMode bgn_mode = background_noise_->mode(); 899 if (bgn_mode == NetEq::kBgnFade && too_many_expands && 900 bgn_mute_factor > 0) { 901 // Fade BGN to zero. 902 // Calculate muting slope, approximately -2^18 / fs_hz. 903 int mute_slope; 904 if (fs_hz_ == 8000) { 905 mute_slope = -32; 906 } else if (fs_hz_ == 16000) { 907 mute_slope = -16; 908 } else if (fs_hz_ == 32000) { 909 mute_slope = -8; 910 } else { 911 mute_slope = -5; 912 } 913 // Use UnmuteSignal function with negative slope. 914 // |bgn_mute_factor| is in Q14. |mute_slope| is in Q20. 915 DspHelper::UnmuteSignal(noise_samples, 916 num_noise_samples, 917 &bgn_mute_factor, 918 mute_slope, 919 noise_samples); 920 } else if (bgn_mute_factor < 16384) { 921 // If mode is kBgnOn, or if kBgnFade has started fading, 922 // use regular |mute_slope|. 923 if (!stop_muting_ && bgn_mode != NetEq::kBgnOff && 924 !(bgn_mode == NetEq::kBgnFade && too_many_expands)) { 925 DspHelper::UnmuteSignal(noise_samples, 926 static_cast<int>(num_noise_samples), 927 &bgn_mute_factor, 928 mute_slope, 929 noise_samples); 930 } else { 931 // kBgnOn and stop muting, or 932 // kBgnOff (mute factor is always 0), or 933 // kBgnFade has reached 0. 934 WebRtcSpl_AffineTransformVector(noise_samples, noise_samples, 935 bgn_mute_factor, 8192, 14, 936 num_noise_samples); 937 } 938 } 939 // Update mute_factor in BackgroundNoise class. 940 background_noise_->SetMuteFactor(channel, bgn_mute_factor); 941 } else { 942 // BGN parameters have not been initialized; use zero noise. 943 memset(noise_samples, 0, sizeof(int16_t) * num_noise_samples); 944 } 945 } 946 947 void Expand::GenerateRandomVector(int16_t seed_increment, 948 size_t length, 949 int16_t* random_vector) { 950 // TODO(turajs): According to hlundin The loop should not be needed. Should be 951 // just as good to generate all of the vector in one call. 952 size_t samples_generated = 0; 953 const size_t kMaxRandSamples = RandomVector::kRandomTableSize; 954 while (samples_generated < length) { 955 size_t rand_length = std::min(length - samples_generated, kMaxRandSamples); 956 random_vector_->IncreaseSeedIncrement(seed_increment); 957 random_vector_->Generate(rand_length, &random_vector[samples_generated]); 958 samples_generated += rand_length; 959 } 960 } 961 962 } // namespace webrtc 963
