function dataout = mywahwahfunction_2(x, fs, damp, minf, maxf, fw)

% change in centre frequency per sample (Hz)
delta = fw/fs;
%delta=0.1;
% 0.1 => at 44100 samples per second should mean  4.41kHz Fc shift per sec

% create triangle wave of centre frequency values
fc = minf:delta:maxf;

while(length(fc) < length(x) )
    fc = [ fc (maxf:-delta:minf) ];
    fc = [ fc (minf:delta:maxf) ];
end

% trim tri wave to size of input
fc = fc(1:length(x));

% difference equation coefficients
% must be recalculated each time Fc changes
F1 = 2*sin((pi*fc(1))/fs);  
% this dictates size of the pass bands
% the Q value changes 
Q1(n)= 2*damp;               

% create emptly out vectors
yh=zeros(size(x));          
yb=zeros(size(x));
yl=zeros(size(x));

% first sample, to avoid referencing of negative signals
yh(1) = x(1);
yb(1) = F1*yh(1);
yl(1) = F1*yb(1);

% apply difference equation to the sample
for n=2:length(x),    
    yh(n) = x(n) - yl(n-1) - Q1(n)*yb(n-1);
    yb(n) = F1*yh(n) + yb(n-1);
    yl(n) = F1*yb(n) + yl(n-1);
    F1 = 2*sin((pi*fc(n))/fs);
end

% Normalise the signals to create an output data
dataout = yb./max(max(abs(yb)));

