% Butterworth filters
freq = 10.^(-3:0.001:3);

% Begin with a clean graphical slate...
for fig=1:6 %bodeplot creates 6 figures.
   figure(fig)
   hold off
endfor


% 2-pole Butterworth
Q = [0.707];
fn = [1];
Butter2 = multistage(freq, Q, fn);
bodeplot(freq, Butter2, 'r')

% 4-pole Butterworth
Q = [0.541, 1.307];
fn = [1,1];
Butter4 = multistage(freq, Q, fn);
bodeplot(freq, Butter4, 'g')

% 6-pole Butterworth
Q = [0.518, 0.707, 1.934];
fn = [1,1,1];
Butter6 = multistage(freq, Q, fn);
bodeplot(freq, Butter6, 'b')

% 8-pole Butterworth
Q = [0.510, 0.601, 0.900, 2.564];
fn = [1,1,1,1];
Butter8 = multistage(freq, Q, fn);
bodeplot(freq, Butter8, 'm')

legendtext=['2-pole';'4-pole';'6-pole';...
   '8-pole'];

figure(1)
title('Gain of Butterworth filters','fontsize',18)
xlabel('Frequency (Hz)','fontsize',16)
ylabel('Gain (dB)','fontsize',16)
legend(legendtext)
axis([0.01,1000,-160,2])
print('ButterworthGain.pdf')
axis([0.1,10,-20,5])
print('ButterworthGainZoom.pdf')

figure(2)
title('Phase of Butterworth filters','fontsize',18)
xlabel('Frequency (Hz)','fontsize',16)
ylabel('Phase (degrees)','fontsize',16)
legend(legendtext)
print('ButterworthPhase.pdf')

figure(3)
title('Time delay of Butterworth filters','fontsize',18)
xlabel('Frequency (Hz)','fontsize',16)
ylabel('Time delay (fc = 1)','fontsize',16)
legend(legendtext)
print('ButterworthTimeDelay.pdf')

figure(4)
title('Time delay vs. attenuation of Butterworth filters','fontsize',18)
xlabel('Attenuation (dB)','fontsize',16)
ylabel('Time delay (s)','fontsize',16)
legend(legendtext)
axis([-3,10])
print('ButterworthTimeDelay-Attenuation.pdf')

figure(5)
title('Step response of Butterworth filters','fontsize',18)
xlabel('Time (s)','fontsize',16)
ylabel('Response','fontsize',16)
legend(legendtext)
axis([-6,16])
print('ButterworthStepResponse.pdf')

figure(6)
title('Impulse response of Butterworth filters','fontsize',18)
xlabel('Time (s)','fontsize',16)
ylabel('Response','fontsize',16)
legend(legendtext)
axis([-6,16])
print('ButterworthImpulseResponse.pdf')
