Main Content

interpolateHRTF

3-D head-related transfer function (HRTF) interpolation

Description

interpolatedHRTF = interpolateHRTF(HRTF,sourcePositions,desiredSourcePositions) returns the interpolated head-related transfer function (HRTF) at the desired position.

example

interpolatedHRTF = interpolateHRTF(sofaObj,desiredSourcePositions) returns the interpolated HRTF using the SOFA file data in sofaObj.

example

interpolatedHRTF = interpolateHRTF(___,Name=Value) specifies options using one or more name-value arguments.

example

interpolateHRTF(___) with no output arguments plots the frequency response of the HRTF at the desired source positions.

example

Examples

collapse all

Modify the 3-D audio image of a sound file by filtering it through a head-related transfer function (HRTF). Set the location of the sound source by specifying the desired azimuth and elevation.

load 'ReferenceHRTF.mat' hrtfData sourcePosition

hrtfData = permute(double(hrtfData),[2,3,1]);

sourcePosition = sourcePosition(:,[1,2]);

Calculate the head-related impulse response (HRIR) using the VBAP algorithm at a desired source position.

desiredAz = 110;
desiredEl = -45;
desiredPosition = [desiredAz desiredEl];

interpolatedIR  = interpolateHRTF(hrtfData,sourcePosition,desiredPosition, ...
                                  Algorithm="vbap");

Create a dsp.AudioFileReader object to read in a file frame by frame. Create an audioDeviceWriter object to play audio to your sound card frame by frame. Create a dsp.FrequencyDomainFIRFilter with the combined left-ear and right-ear filter pair. Since you want to keep the received signals at each ear separate, set SumFilteredOutputs to false.

fileReader = dsp.AudioFileReader("RockDrums-48-stereo-11secs.mp3");
deviceWriter = audioDeviceWriter(SampleRate=fileReader.SampleRate);

spatialFilter = dsp.FrequencyDomainFIRFilter(squeeze(interpolatedIR),SumFilteredOutputs=false);

In an audio stream loop:

  1. Read in a frame of audio data.

  2. Feed the stereo audio data through the filter.

  3. Write the audio to your output device.

while ~isDone(fileReader)
    audioIn = fileReader();

    audioOut = spatialFilter(audioIn);
    deviceWriter(squeeze(audioOut));
end

As a best practice, release your System objects when complete.

release(deviceWriter)
release(fileReader)

Load a SOFA file containing HRTF measurements.

s = sofaread("ReferenceHRTF.sofa");

Specify the desired source positions and then calculate the HRTF at these locations using the interpolateHRTF function.

desiredAz = [-120;-60;0;60;120;0;-120;120];
desiredEl = [-90;90;45;0;-45;0;45;45];
desiredPosition = [desiredAz desiredEl];

interpolatedIR  = interpolateHRTF(s,desiredPosition);

Create an audio file sampled at 48 kHz for compatibility with the SOFA file.

[audio,fs] = audioread("Counting-16-44p1-mono-15secs.wav");
audio = audioresample(audio,InputRate=fs,OutputRate=s.SamplingRate);
audio = audio./max(abs(audio));
audiowrite("Counting-16-48-mono-15secs.wav",audio,s.SamplingRate);

Create a dsp.AudioFileReader object to read in a file frame by frame. Create an audioDeviceWriter object to play audio to your sound card frame by frame. Create a dsp.FrequencyDomainFIRFilter object. Set the Numerator property to the combined left-ear and right-ear filter pair. Since you want to keep the received signals at each ear separate, set SumFilteredOutputs to false.

fileReader = dsp.AudioFileReader("Counting-16-48-mono-15secs.wav");
deviceWriter = audioDeviceWriter(SampleRate=fileReader.SampleRate);

spatialFilter = dsp.FrequencyDomainFIRFilter(squeeze(interpolatedIR(1,:,:)),SumFilteredOutputs=false);

In an audio stream loop:

  1. Read in a frame of audio data.

  2. Feed the audio data through filter.

  3. Write the audio to your output device. If you have a stereo output hardware, such as headphones, you can hear the source shifting position over time.

  4. Modify the desired source position in 2-second intervals by updating the filter coefficients.

durationPerPosition = 2;
samplesPerPosition = durationPerPosition*fileReader.SampleRate;
samplesPerPosition = samplesPerPosition - rem(samplesPerPosition,fileReader.SamplesPerFrame);

sourcePositionIndex = 1;
samplesRead = 0;
while ~isDone(fileReader)
    audioIn = fileReader();
    samplesRead = samplesRead + fileReader.SamplesPerFrame;
    
    audioOut = spatialFilter(audioIn);
    deviceWriter(audioOut);
    
    if mod(samplesRead,samplesPerPosition) == 0
        sourcePositionIndex = sourcePositionIndex + 1;
        spatialFilter.Numerator = squeeze(interpolatedIR(sourcePositionIndex,:,:));
    end
end

As a best practice, release your System objects when complete.

release(deviceWriter)
release(fileReader)

Load a SOFA file containing HRTF measurements.

s = sofaread("ReferenceHRTF.sofa");

Specify the desired source positions and then calculate the HRTF at these locations.

desiredAz = 0:40:180;
desiredEl = -30:15:30;
desiredPosition = [desiredAz', desiredEl'];

interpolatedIR = interpolateHRTF(s, desiredPosition);

Call interpolateHRTF with no output arguments to plot the frequency response of the HRTF.

interpolateHRTF(s, desiredPosition)
legend(Location="southeast")

Figure contains an axes object. The axes object with xlabel Frequency (Hz), ylabel Magnitude (dB) contains 10 objects of type line. These objects represent Measurement 1. Receiver 1, Measurement 1. Receiver 2, Measurement 2. Receiver 1, Measurement 2. Receiver 2, Measurement 3. Receiver 1, Measurement 3. Receiver 2, Measurement 4. Receiver 1, Measurement 4. Receiver 2, Measurement 5. Receiver 1, Measurement 5. Receiver 2.

Input Arguments

collapse all

HRTF values measured at the source positions, specified as a N-by-2-by-M array.

  • N –– Number of known HRTF pairs

  • M –– Number of samples in each known HRTF

If you specify HRTF with real numbers, the function assumes that the input represents an impulse response, and M corresponds to the length of the impulse response. If you specify HRTF with complex numbers, the function assumes that the input represents a transfer function, and M corresponds to the number of bins in the frequency response. The output of the interpolateHRTF function has the same complexity and interpretation as the input.

Data Types: single | double
Complex Number Support: Yes

Source positions corresponding to measured HRTF values, specified as a N-by-2 matrix. N is the number of known HRTF pairs. The two columns correspond to the azimuth and elevation of the source in degrees, respectively.

Azimuth must be in the range [−180,360]. You can use the −180 to 180 convention or the 0 to 360 convention.

Elevation must be in the range [−90,180]. You can use the −90 to 90 convention or the 0 to 180 convention.

Data Types: single | double

Object containing HRTF data from a SOFA file, specified as a SimpleFreeFieldHRIR or SimpleFreeFieldHRTF object. Use the sofaread function to read in data from SOFA files.

Desired source position for HRTF interpolation, specified as a P-by-2 matrix. P is the number of desired source positions. The columns correspond to the desired azimuth and elevation of the source in degrees, respectively.

Azimuth must be in the range [−180,360]. You can use the −180 to 180 convention or the 0 to 360 convention.

Elevation must be in the range [−90,180]. You can use the −90 to 90 convention or the 0 to 180 convention.

Data Types: single | double

Name-Value Arguments

Specify optional pairs of arguments as Name1=Value1,...,NameN=ValueN, where Name is the argument name and Value is the corresponding value. Name-value arguments must appear after other arguments, but the order of the pairs does not matter.

Before R2021a, use commas to separate each name and value, and enclose Name in quotes.

Example: Algorithm="vbap"

Interpolation algorithm, specified as "bilinear", "vbap", or "nearest".

  • "bilinear" –– 3-D bilinear interpolation, as specified by [1].

  • "vbap" –– Vector base amplitude panning interpolation, as specified by [2].

  • "nearest" — Nearest neighbor interpolation, using the Euclidean distance to find the nearest neighbor.

Data Types: char | string

Output Arguments

collapse all

Interpolated HRTF, returned as a P-by-2-by-M array.

  • P –– Number of desired source positions, specified by the number of rows in the desiredSourcePositions input argument.

  • M –– Number of samples in each known HRTF, specified by the number of pages in the HRTF input argument.

interpolatedHRTF has the same complexity and interpretation as the input. If you specify the input, HRTF, with real numbers, the function assumes that the input represents an impulse response. If you specify the input with complex numbers, the function assumes that the input represents a transfer function.

Data Types: single | double
Complex Number Support: Yes

More About

collapse all

SOFA Files

Spatially oriented format for acoustics (SOFA) is a file format for storing spatially oriented acoustic data like HRTFs and binaural or spatial room impulse responses. The AES69 standard [5] defines the SOFA file format.

References

[1] F.P. Freeland, L.W.P. Biscainho and P.S.R. Diniz, "Interpolation of Head-Related Transfer Functions (HRTFS): A multi-source approach." 2004 12th European Signal Processing Conference. Vienna, 2004, pp. 1761–1764.

[2] Pulkki, Ville. "Virtual Sound Source Positioning Using Vector Based Amplitude Panning." Journal of Audio Engineering Society. Vol. 45. Issue 6, pp. 456–466.

[3] Majdak, P., Zotter, F. Brinkmann, F., De Muynke, J., Mihocic, M., and Noisternig, M., “Spatially Oriented Format for Acoustics 2.1: Introduction and Recent Advances.” Journal of the Audio Engineering Society 70, no. 7/8 (July 25, 2022): 565–84. https://doi.org/10.17743/jaes.2022.0026.

[4] Majdak, P., Carpentier, T., Nicol, R., Roginska, A., Suzuki, Y., Watanabe, K., Wierstorf, H., et al., “Spatially Oriented Format for Acoustics: A Data Exchange Format Representing Head-Related Transfer Functions.” In Proceedings of the 134th Convention of the Audio Engineering Society (AES), Paper 8880. Roma, Italy, 2013.

[5] AES69-2022. "AES standard for file exchange - Spatial acoustic data file format." Standard of the Audio Engineering Society. https://www.aes.org/publications/standards/search.cfm?docID=99

Extended Capabilities

C/C++ Code Generation
Generate C and C++ code using MATLAB® Coder™.

Version History

Introduced in R2018b

expand all