|
1 /* |
|
2 * Copyright (C) 2010 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 "HRTFElevation.h" |
|
30 |
|
31 #include "speex/speex_resampler.h" |
|
32 #include "mozilla/PodOperations.h" |
|
33 #include "AudioSampleFormat.h" |
|
34 |
|
35 #include "IRC_Composite_C_R0195-incl.cpp" |
|
36 |
|
37 using namespace std; |
|
38 using namespace mozilla; |
|
39 |
|
40 namespace WebCore { |
|
41 |
|
42 const int elevationSpacing = irc_composite_c_r0195_elevation_interval; |
|
43 const int firstElevation = irc_composite_c_r0195_first_elevation; |
|
44 const int numberOfElevations = MOZ_ARRAY_LENGTH(irc_composite_c_r0195); |
|
45 |
|
46 const unsigned HRTFElevation::NumberOfTotalAzimuths = 360 / 15 * 8; |
|
47 |
|
48 const int rawSampleRate = irc_composite_c_r0195_sample_rate; |
|
49 |
|
50 // Number of frames in an individual impulse response. |
|
51 const size_t ResponseFrameSize = 256; |
|
52 |
|
53 size_t HRTFElevation::sizeOfIncludingThis(mozilla::MallocSizeOf aMallocSizeOf) const |
|
54 { |
|
55 size_t amount = aMallocSizeOf(this); |
|
56 |
|
57 amount += m_kernelListL.SizeOfExcludingThis(aMallocSizeOf); |
|
58 for (size_t i = 0; i < m_kernelListL.Length(); i++) { |
|
59 amount += m_kernelListL[i]->sizeOfIncludingThis(aMallocSizeOf); |
|
60 } |
|
61 |
|
62 return amount; |
|
63 } |
|
64 |
|
65 size_t HRTFElevation::fftSizeForSampleRate(float sampleRate) |
|
66 { |
|
67 // The IRCAM HRTF impulse responses were 512 sample-frames @44.1KHz, |
|
68 // but these have been truncated to 256 samples. |
|
69 // An FFT-size of twice impulse response size is used (for convolution). |
|
70 // So for sample rates of 44.1KHz an FFT size of 512 is good. |
|
71 // We double the FFT-size only for sample rates at least double this. |
|
72 // If the FFT size is too large then the impulse response will be padded |
|
73 // with zeros without the fade-out provided by HRTFKernel. |
|
74 MOZ_ASSERT(sampleRate > 1.0 && sampleRate < 1048576.0); |
|
75 |
|
76 // This is the size if we were to use all raw response samples. |
|
77 unsigned resampledLength = |
|
78 floorf(ResponseFrameSize * sampleRate / rawSampleRate); |
|
79 // Keep things semi-sane, with max FFT size of 1024 and minimum of 4. |
|
80 // "size |= 3" ensures a minimum of 4 (with the size++ below) and sets the |
|
81 // 2 least significant bits for rounding up to the next power of 2 below. |
|
82 unsigned size = min(resampledLength, 1023U); |
|
83 size |= 3; |
|
84 // Round up to the next power of 2, making the FFT size no more than twice |
|
85 // the impulse response length. This doubles size for values that are |
|
86 // already powers of 2. This works by filling in 7 bits to right of the |
|
87 // most significant bit. The most significant bit is no greater than |
|
88 // 1 << 9, and the least significant 2 bits were already set above. |
|
89 size |= (size >> 1); |
|
90 size |= (size >> 2); |
|
91 size |= (size >> 4); |
|
92 size++; |
|
93 MOZ_ASSERT((size & (size - 1)) == 0); |
|
94 |
|
95 return size; |
|
96 } |
|
97 |
|
98 nsReturnRef<HRTFKernel> HRTFElevation::calculateKernelForAzimuthElevation(int azimuth, int elevation, SpeexResamplerState* resampler, float sampleRate) |
|
99 { |
|
100 int elevationIndex = (elevation - firstElevation) / elevationSpacing; |
|
101 MOZ_ASSERT(elevationIndex >= 0 && elevationIndex <= numberOfElevations); |
|
102 |
|
103 int numberOfAzimuths = irc_composite_c_r0195[elevationIndex].count; |
|
104 int azimuthSpacing = 360 / numberOfAzimuths; |
|
105 MOZ_ASSERT(numberOfAzimuths * azimuthSpacing == 360); |
|
106 |
|
107 int azimuthIndex = azimuth / azimuthSpacing; |
|
108 MOZ_ASSERT(azimuthIndex * azimuthSpacing == azimuth); |
|
109 |
|
110 const int16_t (&impulse_response_data)[ResponseFrameSize] = |
|
111 irc_composite_c_r0195[elevationIndex].azimuths[azimuthIndex]; |
|
112 |
|
113 // When libspeex_resampler is compiled with FIXED_POINT, samples in |
|
114 // speex_resampler_process_float are rounded directly to int16_t, which |
|
115 // only works well if the floats are in the range +/-32767. On such |
|
116 // platforms it's better to resample before converting to float anyway. |
|
117 #ifdef MOZ_SAMPLE_TYPE_S16 |
|
118 # define RESAMPLER_PROCESS speex_resampler_process_int |
|
119 const int16_t* response = impulse_response_data; |
|
120 const int16_t* resampledResponse; |
|
121 #else |
|
122 # define RESAMPLER_PROCESS speex_resampler_process_float |
|
123 float response[ResponseFrameSize]; |
|
124 ConvertAudioSamples(impulse_response_data, response, ResponseFrameSize); |
|
125 float* resampledResponse; |
|
126 #endif |
|
127 |
|
128 // Note that depending on the fftSize returned by the panner, we may be truncating the impulse response. |
|
129 const size_t resampledResponseLength = fftSizeForSampleRate(sampleRate) / 2; |
|
130 |
|
131 nsAutoTArray<AudioDataValue, 2 * ResponseFrameSize> resampled; |
|
132 if (sampleRate == rawSampleRate) { |
|
133 resampledResponse = response; |
|
134 MOZ_ASSERT(resampledResponseLength == ResponseFrameSize); |
|
135 } else { |
|
136 resampled.SetLength(resampledResponseLength); |
|
137 resampledResponse = resampled.Elements(); |
|
138 speex_resampler_skip_zeros(resampler); |
|
139 |
|
140 // Feed the input buffer into the resampler. |
|
141 spx_uint32_t in_len = ResponseFrameSize; |
|
142 spx_uint32_t out_len = resampled.Length(); |
|
143 RESAMPLER_PROCESS(resampler, 0, response, &in_len, |
|
144 resampled.Elements(), &out_len); |
|
145 |
|
146 if (out_len < resampled.Length()) { |
|
147 // The input should have all been processed. |
|
148 MOZ_ASSERT(in_len == ResponseFrameSize); |
|
149 // Feed in zeros get the data remaining in the resampler. |
|
150 spx_uint32_t out_index = out_len; |
|
151 in_len = speex_resampler_get_input_latency(resampler); |
|
152 out_len = resampled.Length() - out_index; |
|
153 RESAMPLER_PROCESS(resampler, 0, nullptr, &in_len, |
|
154 resampled.Elements() + out_index, &out_len); |
|
155 out_index += out_len; |
|
156 // There may be some uninitialized samples remaining for very low |
|
157 // sample rates. |
|
158 PodZero(resampled.Elements() + out_index, |
|
159 resampled.Length() - out_index); |
|
160 } |
|
161 |
|
162 speex_resampler_reset_mem(resampler); |
|
163 } |
|
164 |
|
165 #ifdef MOZ_SAMPLE_TYPE_S16 |
|
166 nsAutoTArray<float, 2 * ResponseFrameSize> floatArray; |
|
167 floatArray.SetLength(resampledResponseLength); |
|
168 float *floatResponse = floatArray.Elements(); |
|
169 ConvertAudioSamples(resampledResponse, |
|
170 floatResponse, resampledResponseLength); |
|
171 #else |
|
172 float *floatResponse = resampledResponse; |
|
173 #endif |
|
174 #undef RESAMPLER_PROCESS |
|
175 |
|
176 return HRTFKernel::create(floatResponse, resampledResponseLength, sampleRate); |
|
177 } |
|
178 |
|
179 // The range of elevations for the IRCAM impulse responses varies depending on azimuth, but the minimum elevation appears to always be -45. |
|
180 // |
|
181 // Here's how it goes: |
|
182 static int maxElevations[] = { |
|
183 // Azimuth |
|
184 // |
|
185 90, // 0 |
|
186 45, // 15 |
|
187 60, // 30 |
|
188 45, // 45 |
|
189 75, // 60 |
|
190 45, // 75 |
|
191 60, // 90 |
|
192 45, // 105 |
|
193 75, // 120 |
|
194 45, // 135 |
|
195 60, // 150 |
|
196 45, // 165 |
|
197 75, // 180 |
|
198 45, // 195 |
|
199 60, // 210 |
|
200 45, // 225 |
|
201 75, // 240 |
|
202 45, // 255 |
|
203 60, // 270 |
|
204 45, // 285 |
|
205 75, // 300 |
|
206 45, // 315 |
|
207 60, // 330 |
|
208 45 // 345 |
|
209 }; |
|
210 |
|
211 nsReturnRef<HRTFElevation> HRTFElevation::createBuiltin(int elevation, float sampleRate) |
|
212 { |
|
213 if (elevation < firstElevation || |
|
214 elevation > firstElevation + numberOfElevations * elevationSpacing || |
|
215 (elevation / elevationSpacing) * elevationSpacing != elevation) |
|
216 return nsReturnRef<HRTFElevation>(); |
|
217 |
|
218 // Spacing, in degrees, between every azimuth loaded from resource. |
|
219 // Some elevations do not have data for all these intervals. |
|
220 // See maxElevations. |
|
221 static const unsigned AzimuthSpacing = 15; |
|
222 static const unsigned NumberOfRawAzimuths = 360 / AzimuthSpacing; |
|
223 static_assert(AzimuthSpacing * NumberOfRawAzimuths == 360, |
|
224 "Not a multiple"); |
|
225 static const unsigned InterpolationFactor = |
|
226 NumberOfTotalAzimuths / NumberOfRawAzimuths; |
|
227 static_assert(NumberOfTotalAzimuths == |
|
228 NumberOfRawAzimuths * InterpolationFactor, "Not a multiple"); |
|
229 |
|
230 HRTFKernelList kernelListL; |
|
231 kernelListL.SetLength(NumberOfTotalAzimuths); |
|
232 |
|
233 SpeexResamplerState* resampler = sampleRate == rawSampleRate ? nullptr : |
|
234 speex_resampler_init(1, rawSampleRate, sampleRate, |
|
235 SPEEX_RESAMPLER_QUALITY_DEFAULT, nullptr); |
|
236 |
|
237 // Load convolution kernels from HRTF files. |
|
238 int interpolatedIndex = 0; |
|
239 for (unsigned rawIndex = 0; rawIndex < NumberOfRawAzimuths; ++rawIndex) { |
|
240 // Don't let elevation exceed maximum for this azimuth. |
|
241 int maxElevation = maxElevations[rawIndex]; |
|
242 int actualElevation = min(elevation, maxElevation); |
|
243 |
|
244 kernelListL[interpolatedIndex] = calculateKernelForAzimuthElevation(rawIndex * AzimuthSpacing, actualElevation, resampler, sampleRate); |
|
245 |
|
246 interpolatedIndex += InterpolationFactor; |
|
247 } |
|
248 |
|
249 if (resampler) |
|
250 speex_resampler_destroy(resampler); |
|
251 |
|
252 // Now go back and interpolate intermediate azimuth values. |
|
253 for (unsigned i = 0; i < NumberOfTotalAzimuths; i += InterpolationFactor) { |
|
254 int j = (i + InterpolationFactor) % NumberOfTotalAzimuths; |
|
255 |
|
256 // Create the interpolated convolution kernels and delays. |
|
257 for (unsigned jj = 1; jj < InterpolationFactor; ++jj) { |
|
258 float x = float(jj) / float(InterpolationFactor); // interpolate from 0 -> 1 |
|
259 |
|
260 kernelListL[i + jj] = HRTFKernel::createInterpolatedKernel(kernelListL[i], kernelListL[j], x); |
|
261 } |
|
262 } |
|
263 |
|
264 return nsReturnRef<HRTFElevation>(new HRTFElevation(&kernelListL, elevation, sampleRate)); |
|
265 } |
|
266 |
|
267 nsReturnRef<HRTFElevation> HRTFElevation::createByInterpolatingSlices(HRTFElevation* hrtfElevation1, HRTFElevation* hrtfElevation2, float x, float sampleRate) |
|
268 { |
|
269 MOZ_ASSERT(hrtfElevation1 && hrtfElevation2); |
|
270 if (!hrtfElevation1 || !hrtfElevation2) |
|
271 return nsReturnRef<HRTFElevation>(); |
|
272 |
|
273 MOZ_ASSERT(x >= 0.0 && x < 1.0); |
|
274 |
|
275 HRTFKernelList kernelListL; |
|
276 kernelListL.SetLength(NumberOfTotalAzimuths); |
|
277 |
|
278 const HRTFKernelList& kernelListL1 = hrtfElevation1->kernelListL(); |
|
279 const HRTFKernelList& kernelListL2 = hrtfElevation2->kernelListL(); |
|
280 |
|
281 // Interpolate kernels of corresponding azimuths of the two elevations. |
|
282 for (unsigned i = 0; i < NumberOfTotalAzimuths; ++i) { |
|
283 kernelListL[i] = HRTFKernel::createInterpolatedKernel(kernelListL1[i], kernelListL2[i], x); |
|
284 } |
|
285 |
|
286 // Interpolate elevation angle. |
|
287 double angle = (1.0 - x) * hrtfElevation1->elevationAngle() + x * hrtfElevation2->elevationAngle(); |
|
288 |
|
289 return nsReturnRef<HRTFElevation>(new HRTFElevation(&kernelListL, static_cast<int>(angle), sampleRate)); |
|
290 } |
|
291 |
|
292 void HRTFElevation::getKernelsFromAzimuth(double azimuthBlend, unsigned azimuthIndex, HRTFKernel* &kernelL, HRTFKernel* &kernelR, double& frameDelayL, double& frameDelayR) |
|
293 { |
|
294 bool checkAzimuthBlend = azimuthBlend >= 0.0 && azimuthBlend < 1.0; |
|
295 MOZ_ASSERT(checkAzimuthBlend); |
|
296 if (!checkAzimuthBlend) |
|
297 azimuthBlend = 0.0; |
|
298 |
|
299 unsigned numKernels = m_kernelListL.Length(); |
|
300 |
|
301 bool isIndexGood = azimuthIndex < numKernels; |
|
302 MOZ_ASSERT(isIndexGood); |
|
303 if (!isIndexGood) { |
|
304 kernelL = 0; |
|
305 kernelR = 0; |
|
306 return; |
|
307 } |
|
308 |
|
309 // Return the left and right kernels, |
|
310 // using symmetry to produce the right kernel. |
|
311 kernelL = m_kernelListL[azimuthIndex]; |
|
312 int azimuthIndexR = (numKernels - azimuthIndex) % numKernels; |
|
313 kernelR = m_kernelListL[azimuthIndexR]; |
|
314 |
|
315 frameDelayL = kernelL->frameDelay(); |
|
316 frameDelayR = kernelR->frameDelay(); |
|
317 |
|
318 int azimuthIndex2L = (azimuthIndex + 1) % numKernels; |
|
319 double frameDelay2L = m_kernelListL[azimuthIndex2L]->frameDelay(); |
|
320 int azimuthIndex2R = (numKernels - azimuthIndex2L) % numKernels; |
|
321 double frameDelay2R = m_kernelListL[azimuthIndex2R]->frameDelay(); |
|
322 |
|
323 // Linearly interpolate delays. |
|
324 frameDelayL = (1.0 - azimuthBlend) * frameDelayL + azimuthBlend * frameDelay2L; |
|
325 frameDelayR = (1.0 - azimuthBlend) * frameDelayR + azimuthBlend * frameDelay2R; |
|
326 } |
|
327 |
|
328 } // namespace WebCore |