function out = autowah(in, fc, fb, fw, mix, fs)

% autowah           - auto wah-wah effect
%
% out = autowah(in, fc, fb, fw, mix, fs)
%
% input parameters:
% - in:   input
% - fc:   center frequency [Hz]
% - fb:   bandwidth [Hz]
% - fw:   wah-wah frequency of LFO [Hz]
% - mix:  ratio of dry to wet
% - fs:   sampling frequency

len=length(in)+2;
y=zeros(len,1);
yh=zeros(len,1);
x=zeros(len,1);
x(3:len)=in;

b=(tan(pi*fb/fs)-1)/(tan(pi*fb/fs)+1);
dTH=2*pi*fw/fs;
TH=0:dTH:2*pi*fw/fs*(len-1);
c=-cos(2*pi*(fc+fb/2*cos(TH))/fs);  % modulation

  % non-recursive: XX = b*x(n)+c*(1-b)*x(n-1)+x(n-2)
XX(3:len)=-b.*x(3:len)+c(3:len)'.*(1-b).*x(2:len-1)+x(1:len-2); % non-recursiv part
  % recursiv part
cc=c*(1-b);
for n=3:len
    yh(n) = XX(n)-cc(n)*yh(n-1)+b*yh(n-2);
end
y=x-yh;     % bandpass

out=(1-mix)*x(3:end)+mix*y(3:end);  % mix
out=out/max(abs(out));   % normalise