function y=vibrato(signal,fs,Modfreq,Width)
% This function will result in a flute signal being synthesised into a flute vibrato.
% It will create and plot the window of the input signal first.
% It will then create and plot a Hilbert envelope. This will be the basis
% of the signal to be processed.
% Vibrato will be applied once this enveloped signal is created.
% Resultant y is then plotted.
% Parameters recommended:
% Delay, approx 3 - 10ms. Frequency modulation, 5 - 9 Hz.
% Errors will be given if numbers aren't rounded properly.
[flute, fs] = wavread('filename');
%create and plot the window of the input signal first. this window will be
%used to implement the vibrato function onto.
%this section of code was adapted from the code provided by Luis in our DAS
%tutorials.
integration_time = 0.005; %
integration_time_in_samples = round(integration_time.*fs);
b_value = 1./integration_time_in_samples;
b = repmat(b_value, integration_time_in_samples, 1);
averagesignal = filter(b, 1, abs(flute));
figure (7);
plot(averagesignal);
% create a Hilbert envelope which will be the basis of the signal to be
% processed.
hilbert_env = abs(hilbert(flute)); % create an envelope of the signal (hilbert function)
analytic_flute_response = (flute+1i) .* hilbert_env; % analytic response
flute_envelope = abs(analytic_flute_response); % obtain the envelope of the signal from the magnitude of the Hilbert transform
flute_envelope = flute_envelope./max(flute_envelope); %normalise the resulting envelope
figure (8);
plot(flute_envelope); %plot the created flute envelope.
len = length(flute); %determine the length of the original signal
gausswin_window = (gausswin(1024)); % create a smoothing filter that will be placed over the signal
%a window of approx. 1024 samples is advisable as a good starting point.
smooth_envelope = filter(gausswin_window, len, flute_envelope);
smooth_envelope = smooth_envelope./max(smooth_envelope);
figure (9);
plot(flute_envelope); hold on; plot(smooth_envelope, 'Color', 'red'); hold off
xlabel('number of samples')
%Begin the vibrato processing to the enveloped signal.
Delay = 0.0003; % basic delay of input sample in sec. Approx. 3 - 10 ms recommended.
DELAY = round(Delay*fs); % basic delay in # samples
Width = 0.0003; %(recommended width)
WIDTH = round(Width*fs); % modulation width in # samples
if WIDTH>DELAY
error( 'delay greater than basic delay !!!');
return;
end;
Modfreq = 5; % approx 5 - 9 Hz recommended
MODFREQ=Modfreq/fs; % modulation frequency in # samples
LEN=length(flute_envelope); % # of samples in WAV-file
L=2+DELAY+WIDTH*2; % length of the entire delay
Delayline=zeros(L,1); % memory allocation for delay
y=zeros(size(flute_envelope)); % memory allocation for output vector320
for n=1:(len-1);
M=MODFREQ; %taken from MODFREQ above
MOD=sin(M*2*pi*n); %modulation
ZEIGER=1+DELAY+WIDTH*MOD; %sum of delayed, modulated signal
i=floor(ZEIGER); % floor of the modulated signal
frac=ZEIGER-i; %
Delayline=[flute_envelope(n);Delayline(1:L-1)]; % conjugate flute_envelope with the delayline
%---Linear Interpolation-----------------------------
y(n,1)=Delayline(i+1)*frac+Delayline(i)*(1-frac);
%---Allpass Interpolation------------------------------
%y(n,1)=(Delayline(i+1)+(1-frac)*Delayline(i)-(1-frac)*ya_alt);
%ya_alt=ya(n,1);
end;
figure (10);
plot(y);
end