1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161 162 163 164 165 166 167 168 169 170 171 172 173 174 175
|
/*
* Copyright (c) 2019 The WebRTC project authors. All Rights Reserved.
*
* Use of this source code is governed by a BSD-style license
* that can be found in the LICENSE file in the root of the source
* tree. An additional intellectual property rights grant can be found
* in the file PATENTS. All contributing project authors may
* be found in the AUTHORS file in the root of the source tree.
*/
#include "modules/audio_processing/ns/signal_model_estimator.h"
#include "modules/audio_processing/ns/fast_math.h"
namespace webrtc {
namespace {
constexpr float kOneByFftSizeBy2Plus1 = 1.f / kFftSizeBy2Plus1;
// Computes the difference measure between input spectrum and a template/learned
// noise spectrum.
float ComputeSpectralDiff(
ArrayView<const float, kFftSizeBy2Plus1> conservative_noise_spectrum,
ArrayView<const float, kFftSizeBy2Plus1> signal_spectrum,
float signal_spectral_sum,
float diff_normalization) {
// spectral_diff = var(signal_spectrum) - cov(signal_spectrum, magnAvgPause)^2
// / var(magnAvgPause)
// Compute average quantities.
float noise_average = 0.f;
for (size_t i = 0; i < kFftSizeBy2Plus1; ++i) {
// Conservative smooth noise spectrum from pause frames.
noise_average += conservative_noise_spectrum[i];
}
noise_average = noise_average * kOneByFftSizeBy2Plus1;
float signal_average = signal_spectral_sum * kOneByFftSizeBy2Plus1;
// Compute variance and covariance quantities.
float covariance = 0.f;
float noise_variance = 0.f;
float signal_variance = 0.f;
for (size_t i = 0; i < kFftSizeBy2Plus1; ++i) {
float signal_diff = signal_spectrum[i] - signal_average;
float noise_diff = conservative_noise_spectrum[i] - noise_average;
covariance += signal_diff * noise_diff;
noise_variance += noise_diff * noise_diff;
signal_variance += signal_diff * signal_diff;
}
covariance *= kOneByFftSizeBy2Plus1;
noise_variance *= kOneByFftSizeBy2Plus1;
signal_variance *= kOneByFftSizeBy2Plus1;
// Update of average magnitude spectrum.
float spectral_diff =
signal_variance - (covariance * covariance) / (noise_variance + 0.0001f);
// Normalize.
return spectral_diff / (diff_normalization + 0.0001f);
}
// Updates the spectral flatness based on the input spectrum.
void UpdateSpectralFlatness(
ArrayView<const float, kFftSizeBy2Plus1> signal_spectrum,
float signal_spectral_sum,
float* spectral_flatness) {
RTC_DCHECK(spectral_flatness);
// Compute log of ratio of the geometric to arithmetic mean (handle the log(0)
// separately).
constexpr float kAveraging = 0.3f;
float avg_spect_flatness_num = 0.f;
for (size_t i = 1; i < kFftSizeBy2Plus1; ++i) {
if (signal_spectrum[i] == 0.f) {
*spectral_flatness -= kAveraging * (*spectral_flatness);
return;
}
}
for (size_t i = 1; i < kFftSizeBy2Plus1; ++i) {
avg_spect_flatness_num += LogApproximation(signal_spectrum[i]);
}
float avg_spect_flatness_denom = signal_spectral_sum - signal_spectrum[0];
avg_spect_flatness_denom = avg_spect_flatness_denom * kOneByFftSizeBy2Plus1;
avg_spect_flatness_num = avg_spect_flatness_num * kOneByFftSizeBy2Plus1;
float spectral_tmp =
ExpApproximation(avg_spect_flatness_num) / avg_spect_flatness_denom;
// Time-avg update of spectral flatness feature.
*spectral_flatness += kAveraging * (spectral_tmp - *spectral_flatness);
}
// Updates the log LRT measures.
void UpdateSpectralLrt(ArrayView<const float, kFftSizeBy2Plus1> prior_snr,
ArrayView<const float, kFftSizeBy2Plus1> post_snr,
ArrayView<float, kFftSizeBy2Plus1> avg_log_lrt,
float* lrt) {
RTC_DCHECK(lrt);
for (size_t i = 0; i < kFftSizeBy2Plus1; ++i) {
float tmp1 = 1.f + 2.f * prior_snr[i];
float tmp2 = 2.f * prior_snr[i] / (tmp1 + 0.0001f);
float bessel_tmp = (post_snr[i] + 1.f) * tmp2;
avg_log_lrt[i] +=
.5f * (bessel_tmp - LogApproximation(tmp1) - avg_log_lrt[i]);
}
float log_lrt_time_avg_k_sum = 0.f;
for (size_t i = 0; i < kFftSizeBy2Plus1; ++i) {
log_lrt_time_avg_k_sum += avg_log_lrt[i];
}
*lrt = log_lrt_time_avg_k_sum * kOneByFftSizeBy2Plus1;
}
} // namespace
SignalModelEstimator::SignalModelEstimator()
: prior_model_estimator_(kLtrFeatureThr) {}
void SignalModelEstimator::AdjustNormalization(int32_t num_analyzed_frames,
float signal_energy) {
diff_normalization_ *= num_analyzed_frames;
diff_normalization_ += signal_energy;
diff_normalization_ /= (num_analyzed_frames + 1);
}
// Update the noise features.
void SignalModelEstimator::Update(
ArrayView<const float, kFftSizeBy2Plus1> prior_snr,
ArrayView<const float, kFftSizeBy2Plus1> post_snr,
ArrayView<const float, kFftSizeBy2Plus1> conservative_noise_spectrum,
ArrayView<const float, kFftSizeBy2Plus1> signal_spectrum,
float signal_spectral_sum,
float signal_energy) {
// Compute spectral flatness on input spectrum.
UpdateSpectralFlatness(signal_spectrum, signal_spectral_sum,
&features_.spectral_flatness);
// Compute difference of input spectrum with learned/estimated noise spectrum.
float spectral_diff =
ComputeSpectralDiff(conservative_noise_spectrum, signal_spectrum,
signal_spectral_sum, diff_normalization_);
// Compute time-avg update of difference feature.
features_.spectral_diff += 0.3f * (spectral_diff - features_.spectral_diff);
signal_energy_sum_ += signal_energy;
// Compute histograms for parameter decisions (thresholds and weights for
// features). Parameters are extracted periodically.
if (--histogram_analysis_counter_ > 0) {
histograms_.Update(features_);
} else {
// Compute model parameters.
prior_model_estimator_.Update(histograms_);
// Clear histograms for next update.
histograms_.Clear();
histogram_analysis_counter_ = kFeatureUpdateWindowSize;
// Update every window:
// Compute normalization for the spectral difference for next estimation.
signal_energy_sum_ = signal_energy_sum_ / kFeatureUpdateWindowSize;
diff_normalization_ = 0.5f * (signal_energy_sum_ + diff_normalization_);
signal_energy_sum_ = 0.f;
}
// Compute the LRT.
UpdateSpectralLrt(prior_snr, post_snr, features_.avg_log_lrt, &features_.lrt);
}
} // namespace webrtc
|