Main Content

acousticFluctuation

Perceived fluctuation strength of acoustic signal

Since R2020b

Description

example

fluctuation = acousticFluctuation(audioIn,fs) returns fluctuation strength in vacil based on Zwicker et al. [1] and ISO 532-1 time-varying loudness [2].

example

fluctuation = acousticFluctuation(audioIn,fs,calibrationFactor) specifies a nondefault microphone calibration factor used to compute loudness.

example

fluctuation = acousticFluctuation(specificLoudnessIn) computes fluctuation using time-varying specific loudness.

example

fluctuation = acousticFluctuation(___,Name,Value) specifies options using one or more Name,Value pair arguments.

Example: fluctuation = acousticFluctuation(audioIn,fs,'SoundField','diffuse') returns fluctuation assuming a diffuse sound field.

[fluctuation,specificFluctuation] = acousticFluctuation(___) also returns specific fluctuation strength.

[fluctuation,specificFluctuation,fMod] = acousticFluctuation(___) also returns the dominant modulation frequency.

example

acousticFluctuation(___) with no output arguments plots fluctuation strength and specific fluctuation strength and displays the modulation frequency textually. If the input is stereo, the 3-D plot shows the sum of both channels.

Examples

collapse all

Measure acoustic fluctuation based on Zwicker et al [2] and ISO 532-1 [1]. Assume the recording level is calibrated such that a 1 kHz tone registers as 100 dB on an SPL meter.

[audioIn,fs] = audioread('WashingMachine-16-44p1-stereo-10secs.wav');

fluctuation = acousticFluctuation(audioIn,fs);

Set up an experiment as indicated by the diagram.

Create an audioDeviceReader object to read from the microphone and an audioDeviceWriter object to write to your speaker.

fs = 48e3;
deviceReader = audioDeviceReader(fs,"SamplesPerFrame",2048);
deviceWriter = audioDeviceWriter(fs);

Create an audioOscillator object to generate a 1 kHz sinusoid.

osc = audioOscillator("sine",1e3,"SampleRate",fs,"SamplesPerFrame",2048);

Create a dsp.AsyncBuffer object to buffer data acquired from the microphone.

dur = 5;
buff = dsp.AsyncBuffer(dur*fs);

For five seconds, play the sinusoid through your speaker and record using your microphone. While the audio streams, note the loudness as reported by your SPL meter. Once complete, read the contents of the buffer object.

numFrames = dur*(fs/osc.SamplesPerFrame);
for ii = 1:numFrames
    audioOut = osc();
    deviceWriter(audioOut);
    
    audioIn = deviceReader();
    write(buff,audioIn);
end

SPLreading = 60.4;

micRecording = read(buff);

To compute the calibration factor for the microphone, use the calibrateMicrophone function.

calibrationFactor = calibrateMicrophone(micRecording(fs+1:end,:),deviceReader.SampleRate,SPLreading);

You can now use the calibration factor you determined to measure the fluctuation of any sound that is acquired through the same microphone recording chain.

Perform the experiment again, this time, add 100% amplitude modulation at 4 Hz. To create the modulation signal, use audioOscillator and specify the amplitude as 0.5 and the DC offset as 0.5 to oscillate between 0 and 1.

mod = audioOscillator("sine",4,"SampleRate",fs, ...
    "Amplitude",0.5,"DCOffset",0.5,"SamplesPerFrame",2048);

dur = 5;
buff = dsp.AsyncBuffer(dur*fs);
numFrames = dur*(fs/osc.SamplesPerFrame);
for ii = 1:numFrames
    audioOut = osc().*mod();
    deviceWriter(audioOut);
    
    audioIn = deviceReader();
    write(buff,audioIn);
end

micRecording = read(buff);

Call acousticFluctuation with the microphone recording, sample rate, and calibration factor. The fluctuation reported from acousticFluctuation uses the true acoustic loudness measurement as specified by 532-1. Display the average fluctuation strength over the 5 seconds.

fluctuation = acousticFluctuation(micRecording,deviceReader.SampleRate,calibrationFactor);
fprintf('Average fluctuation = %d (vacil)',mean(fluctuation(501:end,:)))
Average fluctuation = 1.413824e+00 (vacil)
acousticFluctuation(micRecording,deviceReader.SampleRate,calibrationFactor)

Read in an audio file.

[audioIn,fs] = audioread("Engine-16-44p1-stereo-20sec.wav");

Call acousticLoudness to calculate the specific loudness.

[~,specificLoudness] = acousticLoudness(audioIn,fs,'TimeVarying',true);

Call acousticSharpness without any outputs to plot the acoustic sharpness.

acousticSharpness(specificLoudness,'TimeVarying',true)

Call acousticFluctuation without any outputs to plot the acoustic fluctuation.

acousticFluctuation(specificLoudness)

Generate a pure tone with a 1500 Hz center frequency and approximately 700 Hz frequency deviation at a modulation frequency of 0.25 Hz.

fs = 48e3;

fMod = 0.25;
dur = 20;

numSamples = dur*fs;
t = (0:numSamples-1)/fs;

tone = sin(2*pi*t*fMod)';

fc = 1500;
excursionRatio = 0.47;

excursion = 2*pi*(fc*excursionRatio/fs);
audioIn = modulate(tone,fc,fs,'fm',excursion);

Listen to the first 5 seconds of the audio and plot the spectrogram.

sound(audioIn(1:5*fs),fs)
spectrogram(audioIn(1:5*fs),hann(512,'periodic'),256,1024,fs,'yaxis')

Call acousticFluctuation with no output arguments to plot the acoustic fluctuation strength.

acousticFluctuation(audioIn,fs);

The acousticFluctuation function enables you to specify a known fluctuation frequency. If you do not specify a known fluctuation frequency, the function auto-detects the fluctuation.

Create a dsp.AudioFileReader object to read in an audio signal frame-by-frame. Create an audioOscillator object to create a modulation wave. Apply the modulation wave to the audio file.

fileReader = dsp.AudioFileReader('Engine-16-44p1-stereo-20sec.wav');

fmod = 10.8;
amplitude = 0.15;

osc = audioOscillator('sine',fmod, ...
    "DCOffset",0.5, ...
    "Amplitude",amplitude, ...
    "SampleRate",fileReader.SampleRate, ...
    "SamplesPerFrame",fileReader.SamplesPerFrame);

testSignal = [];
while ~isDone(fileReader)
    x = fileReader();
    testSignal = [testSignal;osc().*fileReader()];
end

Listen to two seconds of the test signal and plot its waveform.

samplesToView = 1:2*fileReader.SampleRate;
sound(testSignal(samplesToView,:),fileReader.SampleRate);

plot(samplesToView/fileReader.SampleRate,testSignal(samplesToView,:))
xlabel('Time (s)')

Plot the acoustic fluctuation. The detected frequency of the modulation is displayed textually.

acousticFluctuation(testSignal,fileReader.SampleRate);

Specify the known modulation frequency and then plot the acoustic fluctuation again.

acousticFluctuation(testSignal,fileReader.SampleRate,'ModulationFrequency',fmod)

Input Arguments

collapse all

Audio input, specified as a column vector (mono) or matrix with two columns (stereo).

Tip

To measure fluctuation strength given any modulation frequency, the recommended minimum signal duration is 10 seconds.

Data Types: single | double

Sample rate in Hz, specified as a positive scalar. The recommended sample rate for new recordings is 48 kHz.

Note

The minimum acceptable sample rate is 8 kHz.

Data Types: single | double

Microphone calibration factor, specified as a positive scalar. The default calibration factor corresponds to a full-scale 1 kHz sine wave with a sound pressure level of 100 dB (SPL). To compute the calibration factor specific to your system, use the calibrateMicrophone function.

Data Types: single | double

Specific loudness in sones/Bark, specified as a T-by-240-by-C array, where:

  • T is one per 2 ms of the time-varying signal.

  • 240 is the number of Bark bins in the domain for specific loudness. The Bark bins are 0.1:0.1:24.

  • C is the number of channels.

You can use the acousticLoudness function to calculate time-varying specific loudness using this syntax:

[~,specificLoudness] = acousticLoudness(audioIn,fs,'TimeVarying',true);

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: acousticFluctuation(audioIn,fs,'ModulationFrequency',50)

Known modulation frequency in Hz, specified either 'auto-detect' or as a scalar or two-element vector with values in the range [0.1,100]. If ModulationFrequency is set to 'auto-detect', then the function limits the search range to between 0.2 and 64 Hz. If the input is mono, then the modulation frequency must be specified as a scalar. If the input is stereo, then the modulation frequency can be specified as either a scalar or two-element vector.

Data Types: single | double | char | string

Sound field of audio recording, specified as 'free' or 'diffuse'.

Data Types: char | string

Reference pressure for dB calculation in pascals, specified as a positive scalar. The default value, 20 micropascals, is the common value of air.

Data Types: single | double

Output Arguments

collapse all

Fluctuation strength in vacil, returned as a K-by-1 column vector or K-by-2 matrix of independent channels. K corresponds to the time dimension.

Data Types: single | double

Specific fluctuation strength in vacil/Bark, returned as a K-by-47 matrix or a K-by-47-by-2 array. The first dimension of specificFluctation, K, corresponds to the time dimension and matches the first dimension of fluctuation. The second dimension of specificFluctation, 47, corresponds to bands on the Bark scale, with centers from 0.5 to 23.5, inclusive, in 0.5 increments. The third dimension of specificFluctation corresponds to the number of channels and matches the second dimension of fluctuation.

Data Types: single | double

Dominant modulation frequency in Hz, returned as a scalar for mono input or a 1-by-2 vector for stereo input.

Data Types: single | double

Algorithms

Acoustic fluctuation strength is a perceptual measurement of slow modulations in amplitude or frequency. The acoustic loudness algorithm is described in [1] and implemented in the acousticLoudness function. The acoustic fluctuation calculation is described in [2]. The algorithm for acoustic fluctuation is outlined as follows.

fluctuation=0.008z=024ΔLdz(fmod4)+(4fmod)

Where fmod is the detected or known modulation frequency and ΔL is the perceived modulation depth. If the modulation frequency is not specified when calling acousticFluctuation, it is auto-detected by peak-picking a frequency-domain representation of the acoustic loudness. The perceived modulation depth, ΔL, is calculated by passing rectified specific loudness bands through ½ octave filters centered around fmod, followed by a lowpass filter to determine the envelope.

References

[1] ISO 532-1:2017(E). "Acoustics – Methods for calculating loudness – Part 1: Zwicker method." International Organization for Standardization.

[2] Zwicker, Eberhard, and H. Fastl. Psychoacoustics: Facts and Models. 2nd updated ed, Springer, 1999.

Extended Capabilities

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

Version History

Introduced in R2020b