DESCRIPTION

This is an implementation in Matlab of auditory illusion known as the Shepard tone. The code is created and tested in Matlab R14SP3 Student Version. The Shepard tone is an auditory illustion in which the resultant tone is percieved to ascend infinitely. This effect is also known as "Shepard-Risset glissando" (Source: http://en.wikipedia.org/wiki/Shepard_Tone). The Shepard tone operates like this: first, a bandwidth of two octaves wide is selected. Let us say, for the sake of explanation, that the bandwidth runs from C4 (low C) to C6 (high C), with C5 (middle C) in the middle. Two notes are created, one positioned at the low C mand another one at the middle C. The middle C tone is loud, while the bottom one is soft. Both tones move in same direction (up) with uniform and constant speed. As both tones move up, the upper tone looses its volume while lower tone gains in volume. By the time the middle C tone reaches the upper C above, it is almost inaudible. The lower note, which started at low C and was almost inaudible, becomes loud by the time it reaches the middle C. As the upper C tone fades away, a new weak tone is generated at the low C frequency, and the process is repeated. The volume of the tone is a Gaussian-shaped envelope. At low C and high C the volume is soft, while at the middle C it is at its highest.

Contents

SETUP

Here we first define our key values.

num_steps=100; %number of steps in one octave equal-tempered interval (the smaller the step, the smoother the transition)
step_ratio=2^(1/num_steps); %actual ratio in frequency between two tones on step apart
lo_c=130; %the low C, which is the base frequency for everything
Fs=44100; %sampling frequency
ampl=0.9; %maximum amplitude
t=0:(1/Fs):0.1; %"official" time-vector for the duration of any given step
t_ext=0:(1/Fs):0.2; %"extended" time-vector we will use to get rid of phase problems
step_dur=size(t,2); %number of samples in the "official" time-vector
step_dur_ext=size(t_ext,2); %number of samples in "extended" time-vector
sound=zeros(1,step_dur_ext*num_steps*2); %to hold sound
snd=inline('ampl*cos(2*pi*f0*t)','f0','t','ampl'); %the function used to generate sound

GENERATING SOUND

Here we sample cosine function to get the sound. We sample sinusoids of different frequency for pre-defined equal periods of time. That means that these sinusoids might not be in phase, which will produce clicking sound. Our job here is to ensure that does not happen. Our approach will be to use the "extended" time vector to create sound thatis longer that we need, then take that sound and cut it at the end of the first complete oscillation after "official" time vector expires.

cut=0; %index to remember where we cut
for k=1:(num_steps*2+1)
    sound(cut+1:cut+step_dur_ext)=snd(lo_c*step_ratio^k,t_ext,ampl); %extended sound
    for m=(cut+step_dur):(cut+step_dur_ext) %the portion remaining after the "official" duration of step has passed
        if sound(m)>=(0.99*ampl) %if we are within 0.95 we are close enough to the end of cosine cycle
            cut=m;
            break;
        end
    end
    sound(cut+1:size(sound,2))=0; %erase everything after cut
end
for k=1:size(sound,2) %now lets remove silence at the end
    if (sound(k)==sound(k+1))&&(sound(k)==sound(k+2)) %no cnahge in waveform - we reached the end of sound
        result=zeros(1,k);%row
        result(1:k)=sound(1:k); %resulting sound
        break;
    end
end
window=0.25*gausswin(size(result,2),3.25);%col
plot(1:size(window,1),window);title('Gaussian volume envelope');
result=result';%row->col
result=result.*window;
result_shifted=circshift(result,[floor(size(result,1)/2) 0]); % the tone that starts at middle C
output=[result+result_shifted;
        result+result_shifted];
wavwrite(output,Fs,'tone.wav');

RESULT

Please listen to the resulting WAV file (3.5MB) at http://www.arsenivanov.com/projects/shepard_tone/tone.wav