function outwave = vibrato(inputwave, fs, Modfreq, Width)

WIDTH = round(Width.*fs./1000); % scaled number of samples in modulation
% width.

MODFREQ = Modfreq./fs; % angular frequency of vibrato oscillations with
% respect to the sampling frequency.
input_length = length(inputwave); % number of samples in input wave.
delay_length = 2 + 2*WIDTH; % length of entire delay.
Delayline = zeros(delay_length,1); % memory allocation for delay.
y = zeros(size(inputwave)); % memory allocation for output vector.

for n = 1:(input_length-1)
    
    Delayline = [inputwave(n); Delayline(1:delay_length-1)]; % Convolves
    % input wave with delay line.
    
    MODULATOR = sin(2*pi*n.*MODFREQ); % Generates sine wave of the specified
    % angular frequency to represent the pitch modulation.
    ZEIGER = 1 + WIDTH.*(1 + MODULATOR); % Alters the amplitude of the sine
    % wave (determines amplitude of pitch modulation).
    i = floor(ZEIGER);
    frac = ZEIGER - i; % Calculates the distance between consecutive samples
    % in preparation for interpolation.

    
    %---Linear Interpolation------------------------------------
    y(n,1) = Delayline(i+1).*frac+Delayline(i).*(1-frac);
end

outwave = y./max(abs(y));