1/*
2Copyright 2018 Google Inc. All Rights Reserved.
3
4Licensed under the Apache License, Version 2.0 (the "License");
5you may not use this file except in compliance with the License.
6You may obtain a copy of the License at
7
8 http://www.apache.org/licenses/LICENSE-2.0
9
10Unless required by applicable law or agreed to in writing, software
11distributed under the License is distributed on an "AS-IS" BASIS,
12WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
13See the License for the specific language governing permissions and
14limitations under the License.
15*/
16
17#include "dsp/utils.h"
18
19#include <algorithm>
20#include <cmath>
21#include <limits>
22#include <random>
23#include <vector>
24
25#include "base/constants_and_types.h"
26#include "base/logging.h"
27#include "base/misc_math.h"
28
29#include "dsp/biquad_filter.h"
30#include "dsp/filter_coefficient_generators.h"
31
32namespace {
33
34// The mean and standard deviation of the normal distribution for bandlimited
35// Gaussian noise.
36const float kMean = 0.0f;
37const float kStandardDeviation = 1.0f;
38
39// Maximum group delay in seconds for each filter. In order to avoid audible
40// distortion, the maximum phase shift of a re-combined stereo sequence should
41// not exceed 5ms at high frequencies. That is why, maximum phase shift of
42// each filter is set to 1/2 of that value.
43const float kMaxGroupDelaySeconds = 0.0025f;
44
45// Phase modulation depth, chosen so that for a given max group delay filters
46// provide the lowest cross-correlation coefficient.
47const float kPhaseModulationDepth = 1.18f;
48
49// Constants used in the generation of uniform random number distributions.
50// https://en.wikipedia.org/wiki/Linear_congruential_generator
51const uint64 kMultiplier = 1664525L;
52const uint64 kIncrement = 1013904223L;
53const float kInt32ToFloat =
54 1.0f / static_cast<float>(std::numeric_limits<uint32>::max());
55
56} // namespace
57
58namespace vraudio {
59
60void GenerateGaussianNoise(float mean, float std_deviation, unsigned seed,
61 AudioBuffer::Channel* noise_channel) {
62 DCHECK(noise_channel);
63 // First generate uniform noise.
64 GenerateUniformNoise(min: 0.0f, max: 1.0f, seed, noise_channel);
65 const size_t length = noise_channel->size();
66
67 // Gaussian distribution with mean and standard deviation in pairs via the
68 // box-muller transform
69 // https://en.wikipedia.org/wiki/Box%E2%80%93Muller_transform.
70 for (size_t i = 0; i < length - 1; i += 2) {
71 const float part_one = std::sqrt(x: -2.0f * std::log(x: (*noise_channel)[i]));
72 const float part_two = kTwoPi * (*noise_channel)[i + 1];
73 const float z0 = part_one * std::cos(x: part_two);
74 const float z1 = part_one * std::sin(x: part_two);
75 (*noise_channel)[i] = std_deviation * z0 + mean;
76 (*noise_channel)[i + 1] = std_deviation * z1 + mean;
77 }
78 // Handle the odd buffer length case cheaply.
79 if (length % 2 > 0) {
80 (*noise_channel)[length - 1] = (*noise_channel)[0];
81 }
82}
83
84void GenerateUniformNoise(float min, float max, unsigned seed,
85 AudioBuffer::Channel* noise_channel) {
86 // Simple random generator to avoid the use of std::uniform_real_distribution
87 // affected by https://gcc.gnu.org/bugzilla/show_bug.cgi?id=56202
88 DCHECK(noise_channel);
89 DCHECK_LT(min, max);
90 const float scaled_conversion_factor = kInt32ToFloat * (max - min);
91 uint32 state = static_cast<uint32>(seed);
92 for (float& sample : *noise_channel) {
93 state = static_cast<uint32>(state * kMultiplier + kIncrement);
94 sample = min + static_cast<float>(state) * scaled_conversion_factor;
95 }
96}
97
98void GenerateBandLimitedGaussianNoise(float center_frequency, int sampling_rate,
99 unsigned seed,
100 AudioBuffer* noise_buffer) {
101
102
103 DCHECK(noise_buffer);
104 DCHECK_GT(sampling_rate, 0);
105 DCHECK_LT(center_frequency, static_cast<float>(sampling_rate) / 2.0f);
106 const size_t num_frames = noise_buffer->num_frames();
107
108 BiquadCoefficients bandpass_coefficients = ComputeBandPassBiquadCoefficients(
109 sample_rate: sampling_rate, centre_frequency: center_frequency, /*bandwidth=*/1);
110 BiquadFilter bandpass_filter(bandpass_coefficients, num_frames);
111
112 for (auto& channel : *noise_buffer) {
113 GenerateGaussianNoise(mean: kMean, std_deviation: kStandardDeviation, seed, noise_channel: &channel);
114 bandpass_filter.Filter(input_channel: channel, output_channel: &channel);
115 bandpass_filter.Clear();
116 }
117}
118
119std::unique_ptr<AudioBuffer> GenerateDecorrelationFilters(int sampling_rate) {
120
121 const int kMaxGroupDelaySamples = static_cast<int>(
122 roundf(x: kMaxGroupDelaySeconds * static_cast<float>(sampling_rate)));
123
124 // Filter coefficients according to:
125 // [1] F. Zotter, M. Frank, "Efficient Phantom Source Widening", Archives of
126 // Acoustics, Vol. 38, No. 1, pp. 27–37 (2013).
127 const float g0 = 1.0f - 0.25f * IntegerPow(base: kPhaseModulationDepth, exp: 2);
128 const float g1 = 0.5f * kPhaseModulationDepth -
129 0.0625f * IntegerPow(base: kPhaseModulationDepth, exp: 3);
130 const float g2 = 0.1250f * IntegerPow(base: kPhaseModulationDepth, exp: 2);
131 std::vector<float> filter1_coefficients{g2, g1, g0, -g1, g2};
132 std::vector<float> filter2_coefficients{g2, -g1, g0, g1, g2};
133
134 const size_t filter_length =
135 filter1_coefficients.size() * kMaxGroupDelaySamples;
136 std::unique_ptr<AudioBuffer> decorrelation_filters(
137 new AudioBuffer(kNumStereoChannels, filter_length));
138 decorrelation_filters->Clear();
139
140 for (size_t coefficient = 0; coefficient < filter1_coefficients.size();
141 ++coefficient) {
142 (*decorrelation_filters)[0][coefficient * kMaxGroupDelaySamples] =
143 filter1_coefficients[coefficient];
144 (*decorrelation_filters)[1][coefficient * kMaxGroupDelaySamples] =
145 filter2_coefficients[coefficient];
146 }
147
148 return decorrelation_filters;
149}
150
151size_t GetNumReverbOctaveBands(int sampling_rate) {
152 DCHECK_GT(sampling_rate, 0);
153
154 const float max_band =
155 log2f(x: 0.5f * static_cast<float>(sampling_rate) / kLowestOctaveBandHz);
156 return std::min(a: kNumReverbOctaveBands, b: static_cast<size_t>(roundf(x: max_band)));
157}
158
159size_t GetNumSamplesFromMilliseconds(float milliseconds, int sampling_rate) {
160 DCHECK_GE(milliseconds, 0.0f);
161 DCHECK_GT(sampling_rate, 0);
162 return static_cast<size_t>(milliseconds * kSecondsFromMilliseconds *
163 static_cast<float>(sampling_rate));
164}
165
166size_t CeilToMultipleOfFramesPerBuffer(size_t size, size_t frames_per_buffer) {
167 DCHECK_NE(frames_per_buffer, 0U);
168 const size_t remainder = size % frames_per_buffer;
169 return remainder == 0 ? std::max(a: size, b: frames_per_buffer)
170 : size + frames_per_buffer - remainder;
171}
172
173void GenerateHannWindow(bool full_window, size_t window_length,
174 AudioBuffer::Channel* buffer) {
175
176 DCHECK(buffer);
177 DCHECK_LE(window_length, buffer->size());
178 const float full_window_scaling_factor =
179 kTwoPi / (static_cast<float>(window_length) - 1.0f);
180 const float half_window_scaling_factor =
181 kTwoPi / (2.0f * static_cast<float>(window_length) - 1.0f);
182 const float scaling_factor =
183 (full_window) ? full_window_scaling_factor : half_window_scaling_factor;
184 for (size_t i = 0; i < window_length; ++i) {
185 (*buffer)[i] =
186 0.5f * (1.0f - std::cos(x: scaling_factor * static_cast<float>(i)));
187 }
188}
189
190} // namespace vraudio
191

source code of qtmultimedia/src/3rdparty/resonance-audio/resonance_audio/dsp/utils.cc