frequency_masking.py

"""Frequency masking experiment for COMP28512.

(C) University of Manchester 2015
Author: Andrew Mundy
"""
from __future__ import print_function
import comp28512_utils
from IPython import display
from matplotlib import pyplot as plt
import numpy as np
import os.path
from scipy.io import wavfile
import tempfile


def generate_sound_files(fs=44.1*10**3, full_duration=4.0, signal_duration=2.0,
                         signal_start=2.0, f_mask=1000.0, a_mask=0.5,
                         n_tests=10, path=None):
    """Generate the audio files for the frequency masking experiment.

    Parameters
    ----------
    fs : float
        Sampling frequency.
    full_duration : float
        Full duration of each sample.
    signal_duration : float
        Duration of the non-silent part of each sample.
    signal_start : float
    f_mask : float
        Frequency of the masking signal.
    a_mask : float
        Amplitude of the masking signal.
    n_tests : int
        Number of tests to make.
    path : string or None
        Full path of the directory to write the tests into.  If None then a new
        temporary directory will be created.

    Returns
    -------
    list
        [(frequency, [(AdB, path of wave file)])]
    """
    # Get the directory
    if path is None:
        path = tempfile.mkdtemp(suffix="comp28512")

    # Signal samples
    signal_ns = np.arange(fs * (signal_start + signal_duration),
                          dtype=np.uint32)
    sample_ns = np.arange(fs * signal_start,
                          fs * (signal_start + signal_duration),
                          dtype=np.uint32)

    # Generate the masking signal
    masker = np.zeros(full_duration * fs)
    masker[signal_ns] = a_mask * np.sin(2*np.pi*signal_ns*f_mask/fs)

    # Signal dictionary
    wavs = list()

    # For each sample frequency
    for k in range(1, n_tests + 1):
        # Get the test frequency
        f_test = f_mask + 100*(2*k - 9)  # Magic!

        # Prepare the list of wavs
        f_wavs = list()

        # For various power levels
        for m in range(1, n_tests + 1):
            adb = (m - 10.0)*6
            amplitude = a_mask * 10**(adb / 20)

            # Generate the sound data, first copy the masking signal
            y = np.zeros(masker.size)
            y[:] = masker[:]
            y[sample_ns] += amplitude * np.sin(2*np.pi*f_test*sample_ns/fs)
            y_data = np.int16(y * (2.0**15 - 1))

            # Write to file, add to the list of wav files
            fn = os.path.join(path, "{},{}.wav".format(k, m))
            wavfile.write(fn, fs, y_data)
            f_wavs.append((adb, fn))

        wavs.append((f_test, f_wavs))

    return wavs


def run_experiment(*args, **kwargs):
    """Run the frequency masking experiment in an IPython notebook.

    Parameters
    ----------
    fs : float
        Sampling frequency.
    full_duration : float
        Full duration of each sample.
    signal_duration : float
        Duration of the non-silent part of each sample.
    signal_start : float
    f_mask : float
        Frequency of the masking signal.
    a_mask : float
        Amplitude of the masking signal.
    n_tests : int
        Number of tests to make.
    path : string or None
        Full path of the directory to write the tests into.  If None then a new
        temporary directory will be created.
    """
    display.display_html("<p>Generating audio...</p>", raw=True)

    # Generate the samples given the parameters
    samples = generate_sound_files(*args, **kwargs)

    display.display_html("<p>The following experiment will play a series of "
                        "audio samples, you should indicate whether you hear "
                        "a change occur during a sample. <b>Ensure you are "
                        "wearing headphones.</b></p>", raw=True)

    # Are you sitting comfortably?
    raw_input("Press return to start:")
    # Then I'll begin

    # Run through the experiment
    display.clear_output()
    thresh = np.zeros(len(samples))
    freqs = np.array([f for (f, _) in samples])
    for n, (freq, tests) in enumerate(samples):
        for (adb, wav) in tests:
            display.clear_output()
            print("{:4.0f} Hz at {:4.1f}dB".format(freq, adb))
            play_audio(wav)

            inp = None
            while inp not in ["y", "n"]:
                inp = raw_input("Did you hear any change during the sample? "
                                "(y|n)")

            if inp == "y":
                thresh[n] = adb
                break

    # Plot the results
    fig, ax = plt.subplots(1)
    ax.plot(freqs, thresh)
    ax.set_xlabel("Frequency of masked signal / Hz")
    ax.set_ylabel("Threshold / Amplitude relative to masking signal")
    ax.set_xlim(np.min(freqs), np.max(freqs))
    ax.set_ylim(-55.0, 0)
    ax.grid(True)

    # Return the results
    return freqs, thresh


def play_audio(filename):
    display.display_html(display.HTML(
        "<p><audio autoplay controls>"
        "<source src=\"{}\" type=\"audio/wav\" />"
        "Please listen to {}"
        "</audio></p>".format(comp28512_utils.get_audio_from_file(filename),
                              filename))
    )