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