function mask(f,A)
% MASK Performs a simple psychoacoustic masking test by creating
% bandlimited noise around 1000 Hz and a single sinusoid at
% frequency f with amplitude A. It then plays the noise
% alone, and then the noise plus the sinusoid.
%
% f - frequency of sinusoid (0 to 11025)
% A - amplitude of sinusoid (0 to 1)
% Set sampling rate to 22050 Hz
fs = 22050;
% Create a bandpass filter, centered around 1000 Hz. Since the
% sampling rate is 22050, the Nyquist frequency is 11025.
% 1000/11025 is approximately 0.09, hence the freqeuency
% values of 0.08 and 0.1 below. For more info, do 'help butter'.
[b,a] = butter(4,[0.08 0.1]);
% Create a vector of random white noise (equal in all frequencies)
wn = rand(1,22050);
% Filter the white noise with our filter
wf = filter(b,a,wn);
% By filtering, we've reduced the power in the noise, so we normalize:
wf = wf/max(abs(wf));
% Create the sinusoid at frequency f, with amplitude A:
s = A*cos(2*pi*f/fs*[0:fs-1]);
% Play the sounds
sound(wf,22050)
pause(1) % Pause for one second between sounds
sound(wf+s,22050)