function outputsignal= mymultiband1(data, fs, x1,x2,x3,a1,a2,a3,a4,comp1,comp2,comp3,comp4)
% mymultiband 1 will output a monophonic compressed output signal under the
% name 'outputsignal'
% The parameters for this function are listed below:
% data = the monophonic audio signal
% fs = the sample frequency in ?Hz? e.g. 44100 Hz
% x1 = the first crossover frequency in ?Hz? e.g. 150 Hz
% x2 = the second crossover frequency in ?Hz? e.g. 3500 Hz
% x3 = the third crossover frequency in ?Hz? e.g. 12000 Hz
% a1 = tracking filter coefficient for the first band (attack) which is a number <1 where a smaller
% number will react to changes faster than a larger number e.g. 0.2 would be fast whereas 0.9 would be slow.
% a2 = tracking coefficient for frequency band 2
% a3 = tracking coefficient for frequency band 3
% a4 = tracking coefficient for frequency band 4
% comp1 = the amount of compression (gain) added to the signal in band 1. The values must be between -.01 and
% -0.5 where; -0.5 will apply more gain than the latter and may result in audible distortion.
% comp2 = compression (gain) to frequency band 2
% comp3 = compression (gain) to frequency band 3
% comp4 = compression (gain) to frequency band 4
% FL1, FB1, FB2, FH2 are vectors of frequency band edges in pairs, in ascending order
% between 0 and 1. 1 corresponds to the Nyquist frequency or half the
% sampling frequency. At least one frequency band must have a non-zero
% width
% ALP, ABP1, ABP2, AHP are real vectors the same size as F which specifies the
% desired amplitude of the frequency response of the resultant filter B.
x1 = x1; % Crossover frequency 1
xFL1 = x1./(fs./2); % Start band using crossover frequency. Sample frequency is halved for Nyquist.
xFL2 = (x1+40)./(fs./2); % Cutoff band (40Hz further than start band)
FL1 = [0 xFL1 xFL2 1]; % Frequency band vector
ALP = [1 1 0 0]; % Amplitude vector
x2 = x2; % Second crossover frequency
xFB1_1 = x2./(fs./2); % Second start band
xFB1_2 = (x2+40)./(fs./2); % Cutoff band
FB1 = [0 xFL1 xFL2 xFB1_1 xFB1_2 1]; % Frequency band vector
ABP1 = [0 0 1 1 0 0]; % Amplitude vector
x3 = x3; % Third crossover frequency
xFB2_1 = x3./(fs./2); % Third start band
xFB2_2 = (x3+1000)./(fs./2); % Cutoff band
FB2 = [0 xFB1_1 xFB1_2 xFB2_1 xFB2_2 1]; % Frequency band vector
ABP2 = [0 0 1 1 0 0]; % Amplitude vector
FH1 = [0 xFB2_1 xFB2_2 1]; % Frequency band vector
AHP = [0 0 1 1]; % Amplitude vector
[lpb lpa] = yulewalk (30, FL1, ALP); % IIR filter with a filter order of 30
lpsignal = filtfilt (lpb, lpa, data); % Zero phase distortion filter.
% Compressor Function 1
a1= a1;
comp1 = comp1;
h1=filter([(1-a1)^2],[1.0000 -2*a1 a1^2],abs(lpsignal));
h1=h1/max(h1);
h1=h1.^comp1;
Compsignallp =lpsignal.*h1;
[bpb1 bpa1] = yulewalk (30, FB1, ABP1); % IIR filter with a filter order of 30
bpsignal1 = filtfilt (bpb1, bpa1, data); % Zero phase distortion filter.
% Compressor Function 2
a2 = a2;
comp2 = comp2;
h2=filter([(1-a2)^2],[1.0000 -2*a2 a2^2],abs(bpsignal1));
h2=h2/max(h2);
h2=h2.^comp2;
Compsignalbp1=bpsignal1.*h2;
[bpb2 bpa2] = yulewalk (30, FB2, ABP2); % IIR filter with a filter order of 30
bpsignal2 = filtfilt (bpb2, bpa2, data); % Zero phase distortion filter.
% Compressor Function 3
a3 = a3;
comp3 = comp3;
h3=filter([(1-a3)^2],[1.0000 -2*a3 a3^2],abs(bpsignal2));
h3=h3/max(h3);
h3=h3.^comp3;
Compsignalbp2=bpsignal2.*h3;
a4 = a4;
comp4 = comp4;
[hpb hpa] = yulewalk (30, FH1, AHP); % IIR filter with a filter order of 30
hpsignal = filtfilt (hpb, hpa, data); % Zero phase distortion filter.
% Compressor Function 4
h4=filter([(1-a4)^2],[1.0000 -2*a4 a4^2],abs(hpsignal));
h4=h4/max(h4);
h4=h4.^comp4;
Compsignalhp =hpsignal.*h4;
% Re-joining signals to form 1 signal
prenormsig = Compsignallp + Compsignalbp1 + Compsignalbp2 + Compsignalhp;
mymultiband1 = prenormsig.*max(abs(prenormsig))/max(abs(prenormsig)); % Normalise
outputsignal = mymultiband1;