function [out_wave] = hammondorganism(input, vibrato, tremolo, freq, depth, output)
% Single Celled Hammond Organism Digital Audio Effect
% Developed by: Hayden P. Tweedley
% Final Revision Date: 08/05/2012
%
% Code created for educational purposes by Hayden P. Tweedley
% extended from code written by S. Disch published in DAFX book 2nd edition,
% pp.75 - 76, copyright Wiley & Sons 2011, available at http://www.dafx.de.
% It may be used for educational purposes and not for commercial
% applications without further permission.
%
% Description:
%
% hammondorganism creates a Leslie-rotary-loudspeaker-like effect by
% simulating some of the fundamental constituent aspects of the Leslie
% sound. These currently include sine modulated vibrato and tremolo.
%
% Input Arguments:
%
% * input: specifies input audio filename of the form 'input.wav'. Inverted
% commas must be included and file must be of wav filetype. File must be
% included in the same directory from which the hammondorganism function is
% called.
%
% * vibrato: specifies depth of the vibrato constituent effect. Argument is
% of the range 0 to 1 where 0 indicates no vibrato effect and 1 indicates
% maximum vibrato effect. For best results, a value of 0.75 or less is
% recommended. 0.4 is close to realistic.
%
% * tremolo: specifies depth of the tremolo constituent effect. Argument is
% of the range 0 to 1 where 0 indicates no tremolo effect and 1 indicates
% maximum tremolo effect. For best results, a value of 0.5 or less is
% recommended. 0.2 is close to realistic.
%
% * freq: specifies the modulation frequency governing both vibrato and
% tremolo in Hz. For best results, a value of between 1 and 10 is
% recommended.
% This spans the modulation frequency between the slow and fast speeds of a
% Leslie rotary loudspeaker. A value of 0 is equivalent to no vibrato
% effect and no tremolo effect. 6.5 Hz is close to realistic for fast
% ("tremolo") Leslie, 3Hz is close to realistic for slow ("chorale") Leslie.
%
% * depth: specifies delay in seconds governing pitch depth of vibrato.
% This delay is responsible for the range of pitch variation heard in the
% vibrato effect. For best results, a value between 0.0005 and 0.001 is
% recommended. 0.0009 is close to realistic.
%
% * output: specifies output audio filename of the form 'output.wav'.
% Inverted commas must be included and file must be of wav filetype
% IE. filename must include .wav extension.
%
% Output Arguments:
%
% * Function returns vector out_wave containing the effected audio data.
% This is duplicated in the output file as specified by the input
% argument: "output". Call function by:
% "x = hammondorganism(input, vibrato, tremolo, freq, flanging, output)"
% (without quotation marks) to write returned out_wave to new vector x,
% else out_wave is written to default vector ans.
%
% Other Outputs:
%
% * An output file is generated, the filename of which is specified by the
% input argument: "output".
%
% * A graph of the input waveform is generated in figure 1
%
% * A graph of the output waveform is generated in figure 2
%
%
% END HELP PAGE
[in_wave, Fs, NBITS ] = wavread(input); % read input waveform,
% and write: data to variable in_wave, samplingfrequency to variable Fs,
% number of bits to variable NBITS.
% variable definitions for ease of debug
vib = vibrato; % percent influence of vibrato - try 0.5 or less
trem = tremolo; % percent influence of tremolo - try 0.5 or less
Modfreq = freq; % modulation frequency specified in Hz - try 10Hz or less
Delay = depth; % basic delay of input sample in sec - try 0.0015 or less
Width = Delay; % Set width equal to delay (cannot be greater than delay)
% for this application, else width just scales delay, IE. same as depth
% input argument except that a minimum basic delay can be hard-coded and the
% width variable can then increase upon this.
% DELAY = Basic delay in # of samples
DELAY=round(Delay*Fs); % = basic delay of input sample in seconds,
% mutliplied by sample rate, rounded to nearest integer
% Modulation width in # of samples (same as DELAY)
WIDTH=round(Width*Fs);
% Error Checking:
if WIDTH>DELAY
error('delay greater than basic delay!');
return;
end
% Modulation frequency in expressed as sample frequency (1/samples)
MODFREQ=Modfreq/Fs; %%IE. frequency is 1/time and FS is sample/time,
% therefore modfreq/fs equals time/sample * 1/time = 1/sample
% Number of samples in WAV-file
LEN=length(in_wave); % Total number of samples in input signal
% Length of the entire delay
L=2+DELAY+WIDTH*2;
% Memory allocation for delay
Delayline=zeros(L,1); % generate a vector of length equal to length of
% entire delay line
% Memory allocation for output vector
out_wave = zeros(length(in_wave),1); % generate a vector of length equal to
% length of entire input signal
% (since output signal will be same length as input signal)
for n=1:(LEN-1) % for n = 1 to the length of in_wave by 1 sample increments
M=MODFREQ; % Modulation frequency expressed as sample frequency
MOD=sin(M*2*pi*n); % Sinusoidal modulation function of frequency M weighted by n
ZEIGER=1+DELAY+WIDTH*(1+vib*MOD); %begin a count 'zeiger' of total delay plus
% total delay modulated by the function MOD weighted by the vibrato intensity coefficient.
i=floor(ZEIGER); % floor rounds to nearest integer towards negative infinity
frac=ZEIGER-i; % frac is difference between zeiger and zeiger rounded down to nearest integer.
Delayline=[in_wave(n);Delayline(1:L-1)]; % construct a column vector delayline.
% Add nth sample of in_wave to the start of the colum vector delayline
% while shifting all existing values down by one row.
%%% Linear Interpolation (not used in this implementation):
% out_wave(n,1)=((1+trem*MOD)*(Delayline(i+1)*frac+Delayline(i)*(1-frac))); % Linear Interpolation.
% write output audio data into vector sample by sample, indexed by nth sample, IE. n by 1 vector.
% output audio data is interpolated linearly between samples, and multiplied
% in amplitude by the modulation function as weighted by the tremolo
% intensity coefficient.
%%%
%%% 3rd order Spline interpolation (superior to linear interpolation):
out_wave(n,1)=(1+trem*MOD)*(Delayline(i+1)*frac^3/6+Delayline(i)*...
((1+frac)^3-4*frac^3)/6+Delayline(i-1)*((2-frac)^3-4*(1-frac)^3)/6+Delayline(i-2)*(1-frac)^3/6);% Linear Interpolation.
% write output audio data into vector sample by sample, indexed by nth sample, IE. n by 1 vector.
% output audio data is interpolated by 3rd order spline between samples, and multiplied
% in amplitude by the modulation function as weighted by the tremolo
% intensity coefficient.
%%%
end
out_wave = out_wave./max(abs(out_wave)); %normalising output audio data
wavwrite(out_wave, Fs, NBITS,output); % write audio output data to output file
% as specified by input parameter.
%plot input and outputs waveforms respectively
figure(1)
plot(in_wave,'r');
title('Input Waveform');
figure(2)
plot(out_wave,'b');
title('Output Waveform');
%%% END