% loop filter with rounded coefficients
a1=1;       k1=1;       b1=1/1024;
a2=1/2;     k2=1;       b2=1/16-1/64;
a3=1/4;     k3=1/2;
a4=1/8;     k4=1/4;
a5=1/8;     k5=1/8;

k = [ k1 k2 k3 k4 k5 ];
a = [ a1 a2 a3 a4 a5 ];
b = [ b1 b2 ];

% transfer function
H = L5_H(k, a, b);
NTF = 1 / (1+H);
M = 64;

figure(1);
f = linspace(eps, 0.5, 300);
[mag,ph,w] = bode(H, 2*pi*f);
mag = 20*log10(mag(1,:));
plot(f, mag);
xlabel('Freqency  [f/fs]');
ylabel('Magnitude  [dB]');
title('Loop Filter H for 5th Order Modulator');
axis([ 0 0.5 -25 160 ]);

figure(2);
f = linspace(1e-4, 1/M, 300);
[mag,ph,w] = bode(H, 2*pi*f);
mag = 20*log10(mag(1,:));
ph  = ph(1,:);
% ph  = ph - 180 * (ph > 180);
f   = w/2/pi;
subplot(2,1,1);
plot(f*M, mag);
ylabel('Magnitude  [dB]');
title('Loop Filter H for 5th Order Modulator');
subplot(2,1,2);
plot(f*M, ph);
xlabel('Freqency  [f/fN]');
ylabel('Phase  [degrees]');

figure(3);
[mag,ph,w] = bode(H);
mag = 20*log10(mag(1,:));
ph  = ph(1,:);
f   = w/2/pi;
subplot(2,1,1);
semilogx(f, mag);
ylabel('Magnitude  [dB]');
title('Loop Filter H for 5th Order Modulator');
subplot(2,1,2);
semilogx(f, ph);
xlabel('Freqency  [f/fs]');
ylabel('Phase  [degrees]');

figure(4);
[mag,ph,w] = bode(NTF);
mag = 20*log10(mag(1,:));
semilogx(w/2/pi, mag);
xlabel('Freqency  [f/fs]');
ylabel('Magnitude  [dB]');
title('NTF for 5th Order Modulator');

figure(5);
rlocus(H);
hold on;
% rlocus(H, 0.45);
