Line data Source code
1 : /*
2 : * Copyright (C) 2012 Google Inc. All rights reserved.
3 : *
4 : * Redistribution and use in source and binary forms, with or without
5 : * modification, are permitted provided that the following conditions
6 : * are met:
7 : *
8 : * 1. Redistributions of source code must retain the above copyright
9 : * notice, this list of conditions and the following disclaimer.
10 : * 2. Redistributions in binary form must reproduce the above copyright
11 : * notice, this list of conditions and the following disclaimer in the
12 : * documentation and/or other materials provided with the distribution.
13 : * 3. Neither the name of Apple Computer, Inc. ("Apple") nor the names of
14 : * its contributors may be used to endorse or promote products derived
15 : * from this software without specific prior written permission.
16 : *
17 : * THIS SOFTWARE IS PROVIDED BY APPLE AND ITS CONTRIBUTORS "AS IS" AND ANY
18 : * EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE IMPLIED
19 : * WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR PURPOSE ARE
20 : * DISCLAIMED. IN NO EVENT SHALL APPLE OR ITS CONTRIBUTORS BE LIABLE FOR ANY
21 : * DIRECT, INDIRECT, INCIDENTAL, SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES
22 : * (INCLUDING, BUT NOT LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES;
23 : * LOSS OF USE, DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND
24 : * ON ANY THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
25 : * (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE OF
26 : * THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
27 : */
28 :
29 : #include "PeriodicWave.h"
30 : #include <algorithm>
31 : #include <cmath>
32 : #include <limits>
33 : #include "mozilla/FFTBlock.h"
34 :
35 : const unsigned MinPeriodicWaveSize = 4096; // This must be a power of two.
36 : const unsigned MaxPeriodicWaveSize = 8192; // This must be a power of two.
37 : const float CentsPerRange = 1200 / 3; // 1/3 Octave.
38 :
39 : using namespace mozilla;
40 : using mozilla::dom::OscillatorType;
41 :
42 : namespace WebCore {
43 :
44 : already_AddRefed<PeriodicWave>
45 0 : PeriodicWave::create(float sampleRate,
46 : const float* real,
47 : const float* imag,
48 : size_t numberOfComponents,
49 : bool disableNormalization)
50 : {
51 0 : bool isGood = real && imag && numberOfComponents > 0;
52 0 : MOZ_ASSERT(isGood);
53 0 : if (isGood) {
54 : RefPtr<PeriodicWave> periodicWave =
55 : new PeriodicWave(sampleRate, numberOfComponents,
56 0 : disableNormalization);
57 :
58 : // Limit the number of components used to those for frequencies below the
59 : // Nyquist of the fixed length inverse FFT.
60 0 : size_t halfSize = periodicWave->m_periodicWaveSize / 2;
61 0 : numberOfComponents = std::min(numberOfComponents, halfSize);
62 0 : periodicWave->m_numberOfComponents = numberOfComponents;
63 0 : periodicWave->m_realComponents = new AudioFloatArray(numberOfComponents);
64 0 : periodicWave->m_imagComponents = new AudioFloatArray(numberOfComponents);
65 0 : memcpy(periodicWave->m_realComponents->Elements(), real,
66 0 : numberOfComponents * sizeof(float));
67 0 : memcpy(periodicWave->m_imagComponents->Elements(), imag,
68 0 : numberOfComponents * sizeof(float));
69 :
70 0 : return periodicWave.forget();
71 : }
72 0 : return nullptr;
73 : }
74 :
75 : already_AddRefed<PeriodicWave>
76 0 : PeriodicWave::createSine(float sampleRate)
77 : {
78 : RefPtr<PeriodicWave> periodicWave =
79 0 : new PeriodicWave(sampleRate, MinPeriodicWaveSize, false);
80 0 : periodicWave->generateBasicWaveform(OscillatorType::Sine);
81 0 : return periodicWave.forget();
82 : }
83 :
84 : already_AddRefed<PeriodicWave>
85 0 : PeriodicWave::createSquare(float sampleRate)
86 : {
87 : RefPtr<PeriodicWave> periodicWave =
88 0 : new PeriodicWave(sampleRate, MinPeriodicWaveSize, false);
89 0 : periodicWave->generateBasicWaveform(OscillatorType::Square);
90 0 : return periodicWave.forget();
91 : }
92 :
93 : already_AddRefed<PeriodicWave>
94 0 : PeriodicWave::createSawtooth(float sampleRate)
95 : {
96 : RefPtr<PeriodicWave> periodicWave =
97 0 : new PeriodicWave(sampleRate, MinPeriodicWaveSize, false);
98 0 : periodicWave->generateBasicWaveform(OscillatorType::Sawtooth);
99 0 : return periodicWave.forget();
100 : }
101 :
102 : already_AddRefed<PeriodicWave>
103 0 : PeriodicWave::createTriangle(float sampleRate)
104 : {
105 : RefPtr<PeriodicWave> periodicWave =
106 0 : new PeriodicWave(sampleRate, MinPeriodicWaveSize, false);
107 0 : periodicWave->generateBasicWaveform(OscillatorType::Triangle);
108 0 : return periodicWave.forget();
109 : }
110 :
111 0 : PeriodicWave::PeriodicWave(float sampleRate, size_t numberOfComponents, bool disableNormalization)
112 : : m_sampleRate(sampleRate)
113 : , m_centsPerRange(CentsPerRange)
114 : , m_maxPartialsInBandLimitedTable(0)
115 : , m_normalizationScale(1.0f)
116 0 : , m_disableNormalization(disableNormalization)
117 : {
118 0 : float nyquist = 0.5 * m_sampleRate;
119 :
120 0 : if (numberOfComponents <= MinPeriodicWaveSize) {
121 0 : m_periodicWaveSize = MinPeriodicWaveSize;
122 : } else {
123 0 : unsigned npow2 = powf(2.0f, floorf(logf(numberOfComponents - 1.0)/logf(2.0f) + 1.0f));
124 0 : m_periodicWaveSize = std::min(MaxPeriodicWaveSize, npow2);
125 : }
126 :
127 0 : m_numberOfRanges = (unsigned)(3.0f*logf(m_periodicWaveSize)/logf(2.0f));
128 0 : m_bandLimitedTables.SetLength(m_numberOfRanges);
129 0 : m_lowestFundamentalFrequency = nyquist / maxNumberOfPartials();
130 0 : m_rateScale = m_periodicWaveSize / m_sampleRate;
131 0 : }
132 :
133 0 : size_t PeriodicWave::sizeOfIncludingThis(mozilla::MallocSizeOf aMallocSizeOf) const
134 : {
135 0 : size_t amount = aMallocSizeOf(this);
136 :
137 0 : amount += m_bandLimitedTables.ShallowSizeOfExcludingThis(aMallocSizeOf);
138 0 : for (size_t i = 0; i < m_bandLimitedTables.Length(); i++) {
139 0 : if (m_bandLimitedTables[i]) {
140 0 : amount += m_bandLimitedTables[i]->ShallowSizeOfIncludingThis(aMallocSizeOf);
141 : }
142 : }
143 :
144 0 : return amount;
145 : }
146 :
147 0 : void PeriodicWave::waveDataForFundamentalFrequency(float fundamentalFrequency, float* &lowerWaveData, float* &higherWaveData, float& tableInterpolationFactor)
148 : {
149 :
150 : // Negative frequencies are allowed, in which case we alias
151 : // to the positive frequency.
152 0 : fundamentalFrequency = fabsf(fundamentalFrequency);
153 :
154 : // We only need to rebuild to the tables if the new fundamental
155 : // frequency is low enough to allow for more partials below the
156 : // Nyquist frequency.
157 0 : unsigned numberOfPartials = numberOfPartialsForRange(0);
158 0 : float nyquist = 0.5 * m_sampleRate;
159 0 : if (fundamentalFrequency != 0.0) {
160 0 : numberOfPartials = std::min(numberOfPartials, (unsigned)(nyquist / fundamentalFrequency));
161 : }
162 0 : if (numberOfPartials > m_maxPartialsInBandLimitedTable) {
163 0 : for (unsigned rangeIndex = 0; rangeIndex < m_numberOfRanges; ++rangeIndex) {
164 0 : m_bandLimitedTables[rangeIndex] = 0;
165 : }
166 :
167 : // We need to create the first table to determine the normalization
168 : // constant.
169 0 : createBandLimitedTables(fundamentalFrequency, 0);
170 0 : m_maxPartialsInBandLimitedTable = numberOfPartials;
171 : }
172 :
173 : // Calculate the pitch range.
174 0 : float ratio = fundamentalFrequency > 0 ? fundamentalFrequency / m_lowestFundamentalFrequency : 0.5;
175 0 : float centsAboveLowestFrequency = logf(ratio)/logf(2.0f) * 1200;
176 :
177 : // Add one to round-up to the next range just in time to truncate
178 : // partials before aliasing occurs.
179 0 : float pitchRange = 1 + centsAboveLowestFrequency / m_centsPerRange;
180 :
181 0 : pitchRange = std::max(pitchRange, 0.0f);
182 0 : pitchRange = std::min(pitchRange, static_cast<float>(m_numberOfRanges - 1));
183 :
184 : // The words "lower" and "higher" refer to the table data having
185 : // the lower and higher numbers of partials. It's a little confusing
186 : // since the range index gets larger the more partials we cull out.
187 : // So the lower table data will have a larger range index.
188 0 : unsigned rangeIndex1 = static_cast<unsigned>(pitchRange);
189 0 : unsigned rangeIndex2 = rangeIndex1 < m_numberOfRanges - 1 ? rangeIndex1 + 1 : rangeIndex1;
190 :
191 0 : if (!m_bandLimitedTables[rangeIndex1].get())
192 0 : createBandLimitedTables(fundamentalFrequency, rangeIndex1);
193 :
194 0 : if (!m_bandLimitedTables[rangeIndex2].get())
195 0 : createBandLimitedTables(fundamentalFrequency, rangeIndex2);
196 :
197 0 : lowerWaveData = m_bandLimitedTables[rangeIndex2]->Elements();
198 0 : higherWaveData = m_bandLimitedTables[rangeIndex1]->Elements();
199 :
200 : // Ranges from 0 -> 1 to interpolate between lower -> higher.
201 0 : tableInterpolationFactor = rangeIndex2 - pitchRange;
202 0 : }
203 :
204 0 : unsigned PeriodicWave::maxNumberOfPartials() const
205 : {
206 0 : return m_periodicWaveSize / 2;
207 : }
208 :
209 0 : unsigned PeriodicWave::numberOfPartialsForRange(unsigned rangeIndex) const
210 : {
211 : // Number of cents below nyquist where we cull partials.
212 0 : float centsToCull = rangeIndex * m_centsPerRange;
213 :
214 : // A value from 0 -> 1 representing what fraction of the partials to keep.
215 0 : float cullingScale = pow(2, -centsToCull / 1200);
216 :
217 : // The very top range will have all the partials culled.
218 0 : unsigned numberOfPartials = cullingScale * maxNumberOfPartials();
219 :
220 0 : return numberOfPartials;
221 : }
222 :
223 : // Convert into time-domain wave buffers.
224 : // One table is created for each range for non-aliasing playback
225 : // at different playback rates. Thus, higher ranges have more
226 : // high-frequency partials culled out.
227 0 : void PeriodicWave::createBandLimitedTables(float fundamentalFrequency,
228 : unsigned rangeIndex)
229 : {
230 0 : unsigned fftSize = m_periodicWaveSize;
231 : unsigned i;
232 :
233 0 : const float *realData = m_realComponents->Elements();
234 0 : const float *imagData = m_imagComponents->Elements();
235 :
236 : // This FFTBlock is used to cull partials (represented by frequency bins).
237 0 : FFTBlock frame(fftSize);
238 :
239 : // Find the starting bin where we should start culling the aliasing
240 : // partials for this pitch range. We need to clear out the highest
241 : // frequencies to band-limit the waveform.
242 0 : unsigned numberOfPartials = numberOfPartialsForRange(rangeIndex);
243 : // Also limit to the number of components that are provided.
244 0 : numberOfPartials = std::min(numberOfPartials, m_numberOfComponents - 1);
245 :
246 : // Limit number of partials to those below Nyquist frequency
247 0 : float nyquist = 0.5 * m_sampleRate;
248 0 : if (fundamentalFrequency != 0.0) {
249 0 : numberOfPartials = std::min(numberOfPartials,
250 0 : (unsigned)(nyquist / fundamentalFrequency));
251 : }
252 :
253 : // Copy from loaded frequency data and generate complex conjugate
254 : // because of the way the inverse FFT is defined.
255 : // The coefficients of higher partials remain zero, as initialized in
256 : // the FFTBlock constructor.
257 0 : for (i = 0; i < numberOfPartials + 1; ++i) {
258 0 : frame.RealData(i) = realData[i];
259 0 : frame.ImagData(i) = -imagData[i];
260 : }
261 :
262 : // Clear any DC-offset.
263 0 : frame.RealData(0) = 0;
264 : // Clear value which has no effect.
265 0 : frame.ImagData(0) = 0;
266 :
267 : // Create the band-limited table.
268 0 : AlignedAudioFloatArray* table = new AlignedAudioFloatArray(m_periodicWaveSize);
269 0 : m_bandLimitedTables[rangeIndex] = table;
270 :
271 : // Apply an inverse FFT to generate the time-domain table data.
272 0 : float* data = m_bandLimitedTables[rangeIndex]->Elements();
273 0 : frame.GetInverseWithoutScaling(data);
274 :
275 : // For the first range (which has the highest power), calculate
276 : // its peak value then compute normalization scale.
277 0 : if (!m_disableNormalization && !rangeIndex) {
278 : float maxValue;
279 0 : maxValue = AudioBufferPeakValue(data, m_periodicWaveSize);
280 :
281 0 : if (maxValue)
282 0 : m_normalizationScale = 1.0f / maxValue;
283 : }
284 :
285 : // Apply normalization scale.
286 0 : if (!m_disableNormalization) {
287 0 : AudioBufferInPlaceScale(data, m_normalizationScale, m_periodicWaveSize);
288 : }
289 0 : }
290 :
291 0 : void PeriodicWave::generateBasicWaveform(OscillatorType shape)
292 : {
293 0 : const float piFloat = float(M_PI);
294 0 : unsigned fftSize = periodicWaveSize();
295 0 : unsigned halfSize = fftSize / 2;
296 :
297 0 : m_numberOfComponents = halfSize;
298 0 : m_realComponents = new AudioFloatArray(halfSize);
299 0 : m_imagComponents = new AudioFloatArray(halfSize);
300 0 : float* realP = m_realComponents->Elements();
301 0 : float* imagP = m_imagComponents->Elements();
302 :
303 : // Clear DC and imag value which is ignored.
304 0 : realP[0] = 0;
305 0 : imagP[0] = 0;
306 :
307 0 : for (unsigned n = 1; n < halfSize; ++n) {
308 0 : float omega = 2 * piFloat * n;
309 0 : float invOmega = 1 / omega;
310 :
311 : // Fourier coefficients according to standard definition.
312 : float a; // Coefficient for cos().
313 : float b; // Coefficient for sin().
314 :
315 : // Calculate Fourier coefficients depending on the shape.
316 : // Note that the overall scaling (magnitude) of the waveforms
317 : // is normalized in createBandLimitedTables().
318 0 : switch (shape) {
319 : case OscillatorType::Sine:
320 : // Standard sine wave function.
321 0 : a = 0;
322 0 : b = (n == 1) ? 1 : 0;
323 0 : break;
324 : case OscillatorType::Square:
325 : // Square-shaped waveform with the first half its maximum value
326 : // and the second half its minimum value.
327 0 : a = 0;
328 0 : b = invOmega * ((n & 1) ? 2 : 0);
329 0 : break;
330 : case OscillatorType::Sawtooth:
331 : // Sawtooth-shaped waveform with the first half ramping from
332 : // zero to maximum and the second half from minimum to zero.
333 0 : a = 0;
334 0 : b = -invOmega * cos(0.5 * omega);
335 0 : break;
336 : case OscillatorType::Triangle:
337 : // Triangle-shaped waveform going from its maximum value to
338 : // its minimum value then back to the maximum value.
339 0 : a = 0;
340 0 : if (n & 1) {
341 0 : b = 2 * (2 / (n * piFloat) * 2 / (n * piFloat)) * ((((n - 1) >> 1) & 1) ? -1 : 1);
342 : } else {
343 0 : b = 0;
344 : }
345 0 : break;
346 : default:
347 0 : NS_NOTREACHED("invalid oscillator type");
348 0 : a = 0;
349 0 : b = 0;
350 0 : break;
351 : }
352 :
353 0 : realP[n] = a;
354 0 : imagP[n] = b;
355 : }
356 0 : }
357 :
358 : } // namespace WebCore
|
