function recwhistle2
%RECWHISTLE2 records sound continuously and calls the program improvise
% when a sound is detected.

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% TITLE: Musical Interaction
% NAME: recwhistle2.m
% AUTHOR: Patrick Ramer
% Software: Matlab 5.3
% Modifications:
% Conclusions:
% Errors:
%
% Version: 0.2
% Used files: improvise7.m
% Organization: LAMI EPFL
% Date: 20.8.1999
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% CONSTANTS
samplerate = 8000; % sampling rate
buffer_ms = 200;  % buffer in ms
rec_unit = 1000;  % unit to record a whole melody in ms
nrg_threshold = 0.05; % if signal energy is higher than threshold recording begins

% ------------------------------------------------------------------------------
% MAIN PROGRAM
% ------------------------------------------------------------------------------

% number of samples
n_buffer = samplerate*buffer_ms/1000; % number of samples in buffer
n_rec = samplerate*rec_unit/1000;  %number of samples in record unit

% endless loop
loop = 0;
while loop == 0
   disp('Whistle now...');
 % wait for whistle
 buffer = wavrecord(n_buffer,samplerate);
 energy = max(abs(buffer));
 while energy < nrg_threshold
    buffer = wavrecord(n_buffer,samplerate);
    energy = max(abs(buffer));
 end

 % record whistle until no more whistle
 i = 1;
 y(:,i) = wavrecord(n_rec,samplerate);
 energy = max(abs(y(:,i)));
 while energy > nrg_threshold
    i = i + 1;
    y(:,i) = wavrecord(n_rec,samplerate);
    energy = max(abs(y(:,i)));
 end
 n = i;

 % create one sound vector
 s = buffer;
 for i = 1:n
    first = n_buffer + (i-1)*n_rec + 1;
    last = n_buffer + i*n_rec;
    s(first:last,1) = y(:,i);
 end

 % improvise
 disp('Please wait...');
 improvise7(s,samplerate);
end

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
 

function improvise7(signal,sr)
%IMPROVISE7(signal) calculates and filters the spectrogram of a melody.
%  The notes must always be seperated by a short rest.
%   It saves the melody in a format which contains only note-on's and pitches.
%   No quantisation is made. Module to recognize the key of the melody.
%  Adds a second voice to the melody. Implementation of a sound playing module.
%   The melody is passed in the vector signal.
%  sr: sampling rate of the signal

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% TITLE: Musical Interaction
% NAME: improvise7.m
% AUTHOR: Patrick Ramer
% Software: Matlab 5.3
% Modifications: is now called by rec_whistle
%                error handling
%
% Conclusions:
% Errors: If second voice is not contained in defined tones.
%
% Version: 0.7
% Used files:
% Organization: LAMI EPFL
% Date: 20.8.1999
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

% length(time) = fix(N/n_interval - overlap/n_interval)
% length(freq) = n_window/2 + 1 (fréquence 0 Hz : +1)

% GLOBAL CONSTANTS
global SEMITONE NYQUIST
SEMITONE = 2^(1/12);
NYQUIST = 0.5;

% CONSTANTS
window = 30; %window is the width of the time windows in ms
interval = 10; %interval is the distance between the centerpoints of windows in ms
n_tones = 37; % number of tones including rest
disp_pause = 1; % pause (in sec) needed to display message before leaving program
 

% ------------------------------------------------------------------------------
% MAIN PROGRAM
% ------------------------------------------------------------------------------

% Load tones
tones = tonefreq(n_tones);

% melody
% returns the melody in terms of frequencies, one frequency for each time window
[melody,time,freq] = extract_melody(signal,sr,window,interval,tones);
n_time = length(time);

% format music in terms of frequencies
music_freq = music_format(melody,time);

% number of tones including rest
n_music2 = size(music_freq);
n_music = n_music2(2); % n_music = length(music)

% if there's a melody
if n_music > 1
 % classification in pitches (index of tones)
 music(1,:) = music_freq(1,:);
 music(2,:) = classificate(music_freq(2,:),freq,tones);

 % get key
 key = get_key(music(2,1:n_music-1),n_music-1);

 % create improvisation
 improvised = create_improvisation(music,n_music,key);

 % play music
   play_2voices(improvised,tones);
else
   disp('No melody');
   pause(disp_pause);   % wait to display message
end
 

% -------------------------------------------------------------------------
% FUNCTIONS
% -------------------------------------------------------------------------
 

% -------------------------------------------------------------------------
function impro = create_improvisation(music,n,key);
% CREATE_IMPROVISATION adds an improvisation to the melody in music.
% Last column of music is a rest.
% key is the key of the melody.
% music: double row vector containing note on time information in the
% first row and pitch (number of tone) in the second row.
% n: number of tones including rest (last element)

impro = add_voice2(music,n,key);
 

% -------------------------------------------------------------------------
function voices = add_voice2(music,n,key);
% ADD_VOICE2 adds a second voice to the melody in music.
% Last column of music is a rest.
% key is the key of the melody.
% Error occurs in function play_2voices if 2nd voice is not in tones.
% music: double row vector containing note on time information in the
% first row and pitch (number of tone) in the second row.
% n: number of tones including rest (last element)

% CONSTANTS
% Attention: the definition of these constants depends on the definitions in
% tonefreq
A4 = 2;
C5 = 5;
B5 = 16;

octave = 12;
stones_oct = 7; % scale tones per octave
chromatic_degree = [1 3 5 6 8 10 12]; % number of the degrees in the chromatic scale
 

% Normalization of the melody
% All tones must be between 1 (for C) and 12 (for B)
melody = music(2,:);
for i = 1:n-1
   % too high tones
   while melody(i) < C5
      melody(i) = melody(i) + octave;
   end

   % too low tones
   while melody(i) > B5
      melody(i) = melody(i) - octave;
   end
   melody(i) = melody(i) - C5 + 1;
end

% Determination of the degree
for i = 1:n-1
   m = melody(i);
   chromatic_number = m - key + 1;
   if chromatic_number < 1
      chromatic_number = chromatic_number + octave;
   end
   degree(i) = close_index(chromatic_degree,chromatic_number);
end

% add second voice
for i = 1:n-1
   switch degree(i)
   case 1
      music(3,i) = music(2,i) - 5;
   case {2,4,5}
      music(3,i) = music(2,i) - 3;
   case {3,6,7}
      music(3,i) = music(2,i) - 4;
   end
end

% add rest in second voice
music(3,n) = 1;

% test for tones that are not contained in table
% test in both voices
for i = 1:n-1
   if music(2,i) < A4
      music(2,i) = music(3,i) + octave;
   end
   if music(3,i) < A4
      music(3,i) = music(3,i) + octave;
   end
end

% return
voices = music;
 

% -------------------------------------------------------------------------
function key = get_key(melody,n)
% GET_KEY returns the key of the melody which consists of the tones in melody.
% Recognizes only major keys. Minor keys are interpreted as corresponding major keys.
% It applies criterions as long as there are more than 1 possible keys.
% n: number of tones in melody

% CONSTANTS
% Attention: the definition of these constants depends on the definitions in
% tonefreq
C5 = 5;
B5 = 16;
octave = 12;
stones_oct = 7; % scale tones per octave
third = 8;   % depends on scales
dominant = 9;  % depends on scales

% initialization
n_diff_tones = 0;

% Normalization of the melody
% All tones must be between 1 (for C) and 12 (for B)
% Each tone must appear only once
for i = 1:n
   % too high tones
   while melody(i) < C5
      melody(i) = melody(i) + octave;
   end

   % too low tones
   while melody(i) > B5
      melody(i) = melody(i) - octave;
   end
   melody(i) = melody(i) - C5 + 1;

   % add only the first time the tone appears
   flag = 1;
   for j = 1:n_diff_tones
      if diff_mel(j) == melody(i)
         flag = 0;
      end
   end
   if flag == 1
      n_diff_tones = n_diff_tones + 1;
      diff_mel(n_diff_tones) = melody(i);
   end
end

% create scales
scales = create_scales;

% test for melody tones in all scales.
% all scales which contain most tones of the melody
% will be in the variable key.
% 12*7*nb_of_tones operations: could be optimised
%test every scale
for i = 1:octave
   % test every tone of the melody
   n_tones_scale(i) = 0; % number of tones that are found in scale
   for j = 1:n_diff_tones
      % test every tone of the scale
      flag_tone = 0;
      for k = 1:stones_oct
         if scales(i,k) == diff_mel(j)
            % tone of the melody is found in scale
            flag_tone = 1;
         end
      end
      if flag_tone == 1
         % if tone is found in the scale
     n_tones_scale(i) = n_tones_scale(i) + 1;
      end
   end
end
% find keys with maximum match
[max_scale key(1)] = max(n_tones_scale);
n_keys = 1; % number of possible keys
for i = key(1) + 1:octave
   if n_tones_scale(i) == max_scale
      n_keys = n_keys + 1;
      key(n_keys) = i;
   end
end

% beginning criterion of the melody (1st)
if n_keys > 1
   n_keys_loop = 0; % number of possible keys, used in loops
   for i = 1:n_keys
      k = key(i);
  if (melody(1) == k)
         n_keys_loop = n_keys_loop + 1;
         key_loop(n_keys_loop) = k;
      end
   end
   % accept only if there's at least one possible key left
   if n_keys_loop ~= 0
    key = key_loop;
    n_keys = n_keys_loop;
 end
end

% beginning criterion of the melody (5th)
if n_keys > 1
   n_keys_loop = 0; % number of possible keys, used in loops
   for i = 1:n_keys
      k = key(i);
  if (melody(1) == scales(k,dominant))
         n_keys_loop = n_keys_loop + 1;
         key_loop(n_keys_loop) = k;
      end
   end
   % accept only if there's at least one possible key left
   if n_keys_loop ~= 0
    key = key_loop;
    n_keys = n_keys_loop;
 end
end

% beginning criterion of the melody (3rd)
if n_keys > 1
   n_keys_loop = 0; % number of possible keys, used in loops
   for i = 1:n_keys
      k = key(i);
  if melody(1) == scales(k,third)
         n_keys_loop = n_keys_loop + 1;
         key_loop(n_keys_loop) = k;
      end
   end
   % accept only if there's at least one possible key left
   if n_keys_loop ~= 0
    key = key_loop;
    n_keys = n_keys_loop;
 end
end

% end criterion of the melody (1st, 3rd, 5th)
if n_keys > 1
   n_keys_loop = 0; % number of possible keys, used in loops
   key_loop = [];
   for i = 1:n_keys
      k = key(i);
  if (melody(n) == k) | (melody(n) == scales(k,third)) | (melody(n) == scales(k,dominant))
         n_keys_loop = n_keys_loop + 1;
         key_loop(n_keys_loop) = k;
      end
   end
   % accept only if there's at least one possible key left
   if n_keys_loop ~= 0
    key = key_loop;
    n_keys = n_keys_loop;
 end
end

% do random key if there's more than one possible key left
if n_keys > 1
 key = key(fix(rand * n_keys) + 1);
end
 

% -------------------------------------------------------------------------
function scales = create_scales;
% CREATE_SCALES returns the number of the tones in the chromatic scale
% for each major scale.
% Returns a 12 by 9 matrix.
% Each row corresponds to a scale consisting of 7 different tones
% 8. column are thirds, 9. are dominants
scales = [1 3 5 6 8 10 12 5 8; % C
   1 2 4 6 7 9 11 6 9;    % Csharp
   2 3 5 7 8 10 12 7 10;   % D
   1 3 4 6 8 9 11 8 11;    % Dsharp
   2 4 5 7 9 10 12 9 12;   % E
   1 3 5 6 8 10 11 10 1;   % F
   2 4 6 7 9 11 12 11 2;   % Fsharp
   1 3 5 7 8 10 12 12 3;   % G
   1 2 4 6 8 9 11 1 4;    % Gsharp
   2 3 5 7 9 10 12 2 5;    % A
   1 3 4 6 8 10 11 3 6;    % Asharp
   2 4 5 7 9 11 12 4 7];   % B
 

% -------------------------------------------------------------------------
function play_2voices(music,tones)
% PLAY_1VOICE plays a melody consisting of 1 voice.
% music: three row vector containing note-on time information in the
% first row and pitch (number of tone) in the following rows.
% One row for each voice

% CONSTANTS
transpose = 1/4; % transpose two octaves lower
sr_play = 8000; % sampling rate for playback
rest = 10; % short rest between tones in ms
pi = 3.1415926535897;

n2 = size(music);
n = n2(2); %length(music)

% create time vector
time_play = 0:1/sr_play:music(1,n);

% initialization
last = 0;

% create vector to play
for i = 1:n-1
   first = last + 1;
   last = close_index(time_play,music(1,i+1)) - 1;
   tone_end = last - fix(rest * sr_play/1000);
   voice1(first:tone_end) = sin(2 * pi * tones(music(2,i)) * transpose * time_play(first:tone_end));
   voice2(first:tone_end) = sin(2 * pi * tones(music(3,i)) * transpose * time_play(first:tone_end));
   voice1(tone_end+1:last) = 0; %rest between tones
   voice2(tone_end+1:last) = 0; %rest between tones
end
play = voice1 + voice2;

% play the scaled sound
soundsc(play,sr_play);

% wait for melody to end
duration = ceil(music(1,n)) + 1;
pause(duration);

% -------------------------------------------------------------------------
function m = classificate(music_freq,freq,tones);
% CLASSIFICATE does a relative classification (intervals) of the melody
% in music_freq, which is presented in terms of frequencies.
% Last element in music_freq is a rest.
% freq: frequency vector
% tones: tone frequencies including rest

% CONSTANTS
global SEMITONE

% Load tones
n_tones = length(tones);

n2 = size(music_freq);
n = n2(2); %length(music_freq)

if n > 1
   % take every tone as reference
   for i = 1:n-1
      ref = i;
      % determine first tone
  % if there is no tone, then put a rest
    mel(ref,ref) = close_index(tones,music_freq(ref));
    f1 = music_freq(ref);

  % calculate n-1 intervals
  for j = 1:n-1
     f2 = music_freq(j);
     interval = round(log(f2/f1)/log(SEMITONE));
     mel(ref,j) = mel(ref,ref) + interval;
      end
   end
   m = round(sum(mel)/(n-1));
else
   m(1) = 1;
end

% add a rest at the end
if n > 1
   m(n) = 1;
end
 

% -------------------------------------------------------------------------
function m = music_format(mel,t)
% MUSIC_FORMAT saves the melody in a musical format.
% It saves the melody in a format which contains only note-on's and frequencies.
% Nevertheless the last element is a note off.
% m is a double row vector containing 'note on' time information in the
% first row and frequency in the second row.
% A 'note on' is created when there is a change in frequency. So every note can
% be preceded by a rest or another note.
% No quantisation is done. The 'note on' information is absolute time.
% If there is no tone in melody, the program terminates with an error
% mel: vector of fundamental frequency for every instant
% t: time vector
n_time = length(t);
j = 0; % number of tone in the melody

% first element
if (mel(1) ~= 0)
 j = 1;
   music(:,j) = [t(1) ; mel(1)];
end

% elements 2 to n_time
% recognizes note on and note off
for i = 2:n_time
   if mel(i) ~= mel(i-1)
      j = j + 1;
      music(:,j) = [t(i) ; mel(i)];
   end
end

% last element must be a note off
if (mel(n_time) ~= 0)
   j = j + 1;
   music(:,j) = [t(n_time) ; 0];
end

% if there's no tone
if j == 0
   music = [0 ; 0];
end

% delete note off's except last one
n2 = size(music);
n = n2(2);    % n = length(music)
for i = 1:2:n-1
   m(:,(i+1)/2) = music(:,i);
end
m(:,fix(n/2 + 1)) = music(:,n);

% cut off the first rest before the first tone
if m(1,1) ~= 0
 m(1,:) = m(1,:) - m(1,1);
end
 

% -------------------------------------------------------------------------
function [m,time,freq] = extract_melody(signal,samplerate,window,interval,tones)
% EXTRACT_MELODY returns the melody in terms of frequencies,
% one frequency for each time window
% executes Short Term Fourier Transform and determines the fundamental frequency
% at every instant in time. The melody is filtered.
% The melody is returned in row vector containing the fundamental frequency
% for every instant.
% signal: sound signal vector
% window: width of the window in ms.
% interval: distance between the centerpoints of windows in ms.
% tones: tone frequencies including rest

% CONSTANTS
global NYQUIST
noise_threshold = 20; % signal to average rate (frequency)
low_cutoff_dist = 0.7; % percentage of lowest frequency being cut off

% Parameters
N = length(signal);
n_window = fix(window * samplerate/1000);
n_overlap = fix(n_window - interval * samplerate/1000);

% Calculates and plots the spectrogram
[sp,freq,time] = specgram(signal,n_window,samplerate,n_window,n_overlap);

n_freq = length(freq);
n_time = length(time);

%Spectrum
spectrum = abs(sp(1:n_freq,:));     % number of columns = length(time)

% Low frequency noise filtering
% All amplitudes for frequencies below low_cutoff are set to zero
low_cutoff = close_index(freq,low_cutoff_dist * tones(2));
spectrum(1:low_cutoff,:) = 0;

% Calculation of fundamental frequencies
[maxvalue indexfreq] = max(spectrum); % 2 row vector

% Filtering of noise
% If the signal is not higher than a constant times the average,
% then it is considered being a rest.
average = sum(spectrum)/n_freq;
for i = 1:n_time
 if (maxvalue(i) < noise_threshold * average(i))
    indexfreq(i) = 1;
 end
end

% the melody in terms of frequencies
freq_melody = freq(indexfreq);

% filter melody to have static tones
filt_freq = filter_melody(freq_melody,interval);

% melody
m = filt_freq;
 

% -------------------------------------------------------------------------
function y = filter_melody(freq_melody,interval)
% FILTER_MELODY filters the tones to be static.
% It looks for the transitions from a rest to a tone and from a tone to a rest
% and sets all the tones in-between to the average tone.
% It prohibits ton durations below min_duration
% freq_melody: row vector containing the fundamental frequency of each time window
% interval: interval between time windows in ms

% CONSTANTS
min_duration = 50; % a tone must last for 50 ms at least

l = fix(min_duration/interval); % A vote is made over l elements

l2 = fix(l/2);
n = length(freq_melody);

% initializations
y = freq_melody; % output
flag = 0; % flag = 1 when frequency not equal to 0
s = 0; % Sum of fundamental frequencies of all time windows of one tone

% first element
if (freq_melody(1) ~= 0)
   trans12 = 1; % beginning of a tone, frequency passes from 0 to another value
   flag = 1;  % flag = 1 when frequency not equal to 0
   s = freq_melody(1);
end

% elements 2 to n (n-1)
for i = 2:n
   % rising
   if (freq_melody(i) ~= 0) & (freq_melody(i-1) == 0)
      trans12 = i;
      flag = 1;
   end
   % descending
   if (freq_melody(i) == 0) & (freq_melody(i-1) ~= 0)
      if i - trans12 < l % minimal duration not respected
       y(trans12:i-1) = 0;    % frequency 0 Hz, representing a rest
      else % minimal duration respected -> vote
         vote = round(s/(i-trans12));
         y(trans12:i-1) = vote;
      end
      s = 0;
      flag = 0;
   end
   % static
   if flag == 1
      s = s + freq_melody(i);
   end
end

% last element
if (freq_melody(n) ~= 0)
   for j = trans12:n
    y(j) = round(s/(n-trans12+1));
   end
end

% -----------------------------------------------------------------------
function f = tonefreq(n)
% TONEFREQ(n) defines the frequencies corresponding to the tones.
% First frequency is zero, representing a rest. Second frequency is A4.
% All the following are a consecutive semitone higher until number n.
% The frequencies are stored in a one row matrice.
global SEMITONE
A4 = 440;
rest = 0;
f(1) = rest;
f(2) = A4;
for i = 3:n
   f(i) = f(i-1) * SEMITONE;
end

% -----------------------------------------------------------------------
function i = close_index(v,s)
% CLOSE_INDEX(v,s) returns the index of the value of vector v
% which is the closest to the value s
d = v - s;
[mindist i] = min(abs(d));

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%