Design Examples Using Least Square IIR Filter Optimization Program

File LangIIR.zip consist of all files needed for this demo.
I recreated all scripts from Mr. Lang's dissertation as well. You can find it under the link diss_lang.zip .

Contents

Example 1: Low Order Lowpass Filter with Approximately Linear Passband Phase (M=4, N=4)

M=4;
N=4;
tau=5;
om=pi*[linspace(0,0.2,20),linspace(0.4,1,60)];
D=[exp(-1i*om(1:20)*tau),zeros(1,60)];
W=ones(1,80);
[b,a,e]=mpiir_l2(M,N,om,D,W,0.98);
   ERROR    MAX.RADIUS    STEP       SLOPE
 9.280e+000 7.036e-001 1.353e+000 -1.532e+000 
 3.499e+000 7.302e-001 1.162e+000 -1.675e+000 
 1.066e+000 7.495e-001 8.979e-001 -7.679e-001 
 7.201e-001 7.236e-001 9.646e-001 -7.300e-001 
 4.508e-001 8.181e-001 8.228e-001 -1.374e+000 
 1.851e-001 7.819e-001 3.298e-001 -7.704e-001 
 1.270e-001 7.878e-001 4.324e-001 -3.952e-001 
 1.054e-001 7.897e-001 5.743e-001 -2.697e-001 
 9.155e-002 7.875e-001 6.806e-001 -3.394e-001 
 7.479e-002 7.871e-001 7.235e-001 -6.559e-001 
 4.806e-002 8.382e-001 5.452e-001 -1.179e+000 
 2.384e-002 8.574e-001 8.139e-002 -7.582e-001 
 1.665e-002 8.678e-001 3.464e-002 -3.804e-001 
 1.413e-002 8.742e-001 2.118e-002 -1.572e-001 
 1.325e-002 8.782e-001 1.250e-002 -5.872e-002 
 1.294e-002 8.805e-001 7.406e-003 -2.110e-002 
 1.283e-002 8.819e-001 4.402e-003 -7.492e-003 
 1.279e-002 8.827e-001 2.623e-003 -2.653e-003 
 1.278e-002 8.832e-001 1.565e-003 -9.401e-004 
 1.278e-002 8.835e-001 9.341e-004 -3.335e-004 
 1.277e-002 8.837e-001 5.579e-004 -1.185e-004 
 1.277e-002 8.838e-001 3.333e-004 -4.214e-005 
 1.277e-002 8.839e-001 1.991e-004 -1.500e-005 
 1.277e-002 8.839e-001 1.189e-004 -5.341e-006 
 1.277e-002 8.839e-001 7.106e-005 -1.903e-006 

Filter calculated successfully

h=fvtool(b,a);
axis([0 1 -60 5])
title('Magnitude of Frequency Response (dB)')

Figure 5.2(a)

h=fvtool(b,a);
title('Passband Magnitude Response (dB)')
axis([0 0.2 -0.4 0.4])

Figure 5.2(b)

N=8192;
w=pi*linspace(0,1,N);
H=freqz(b,a,w);
plot(w/pi,phase(H)-phase(exp(-1i*pi/N*(0:N-1)*tau)));
grid on;
title('Passband Phase Error (radians)');
xlabel('w/pi');
axis([0 0.2 -0.04 0.04]);

Figure 5.2(c)

h=fvtool(b,a);
set(h,'Analysis','grpdelay')
title('Passband Group Delay (samples)')
axis([0 0.2 4.5 6.5])

Figure 5.2(d)

h=fvtool(b,a);
set(h,'Analysis','polezero')
title('Poles (x) and Zeros (o)')

Figure 5.3(a)

Example 2: Lowpass Filter with Approximately Linear Passband Phase and Magnitude Constraints (M=15, N=15)

M=15;
N=15;
tau=15;
om=pi*[linspace(0,0.4,40),linspace(0.56,1,44)];
D=[exp(-1i*om(1:40)*tau),zeros(1,44)];
W=[ones(1,40),100*ones(1,44)];
[b,a,e]=mpiir_l2(M,N,om,D,W,0.8263);
   ERROR    MAX.RADIUS    STEP       SLOPE
 1.974e+001 6.747e-001 1.164e+000 -1.297e+000 
 8.797e+000 7.132e-001 8.404e-001 -1.482e+000 
 3.585e+000 7.510e-001 6.443e-001 -1.400e+000 
 1.518e+000 7.803e-001 5.471e-001 -1.342e+000 
 6.750e-001 8.023e-001 4.900e-001 -1.237e+000 
 3.245e-001 8.138e-001 4.496e-001 -1.138e+000 
 1.675e-001 8.194e-001 4.177e-001 -1.051e+000 
 9.196e-002 8.230e-001 4.043e-001 -9.621e-001 
 5.367e-002 8.251e-001 5.276e-001 -8.558e-001 
 3.329e-002 8.258e-001 5.472e-001 -7.480e-001 
 2.167e-002 8.261e-001 3.317e-001 -6.324e-001 
 1.527e-002 8.260e-001 1.708e-001 -4.855e-001 
 1.183e-002 8.258e-001 9.125e-002 -3.434e-001 
 9.933e-003 8.256e-001 5.687e-002 -2.345e-001 
 8.838e-003 8.253e-001 4.383e-002 -1.611e-001 
 8.161e-003 8.252e-001 4.270e-002 -1.155e-001 
 7.707e-003 8.257e-001 5.023e-002 -8.971e-002 
 7.370e-003 8.260e-001 6.600e-002 -7.821e-002 
 7.085e-003 8.261e-001 9.127e-002 -7.822e-002 
 6.810e-003 8.262e-001 1.285e-001 -8.849e-002 
 6.513e-003 8.262e-001 1.808e-001 -1.079e-001 
 6.174e-003 8.263e-001 2.548e-001 -1.315e-001 
 5.802e-003 8.263e-001 3.812e-001 -1.358e-001 
 5.458e-003 8.263e-001 4.426e-001 -1.165e-001 
 5.135e-003 8.263e-001 4.631e-001 -1.418e-001 
 4.747e-003 8.263e-001 4.161e-001 -1.799e-001 
 4.302e-003 8.263e-001 3.435e-001 -2.126e-001 
 3.846e-003 8.263e-001 2.523e-001 -2.331e-001 
 3.409e-003 8.263e-001 1.431e-001 -2.559e-001 
 2.984e-003 8.262e-001 5.680e-002 -2.879e-001 
 2.573e-003 8.262e-001 1.292e-002 -3.518e-001 
 2.156e-003 8.262e-001 6.537e-002 -4.199e-001 
 1.768e-003 8.262e-001 1.171e-001 -4.212e-001 
 1.497e-003 8.262e-001 4.402e-002 -1.293e-001 
 1.409e-003 8.263e-001 8.785e-002 -8.741e-002 
 1.348e-003 8.263e-001 1.263e-001 -7.747e-002 
 1.296e-003 8.263e-001 1.762e-001 -8.091e-002 
 1.244e-003 8.263e-001 2.414e-001 -1.002e-001 
 1.180e-003 8.263e-001 3.204e-001 -1.427e-001 
 1.094e-003 8.263e-001 4.063e-001 -2.264e-001 
 9.658e-004 8.263e-001 4.887e-001 -3.967e-001 
 7.682e-004 8.263e-001 3.504e-001 -4.515e-001 
 5.986e-004 8.263e-001 1.664e-001 -3.132e-001 
 5.085e-004 8.263e-001 7.748e-002 -1.763e-001 
 4.652e-004 8.263e-001 4.328e-002 -9.154e-002 
 4.443e-004 8.263e-001 2.496e-002 -4.650e-002 
 4.341e-004 8.263e-001 1.932e-002 -2.356e-002 
 4.290e-004 8.263e-001 2.129e-002 -1.209e-002 
 4.264e-004 8.263e-001 2.863e-002 -6.470e-003 
 4.250e-004 8.263e-001 4.118e-002 -3.829e-003 
 4.242e-004 8.263e-001 6.032e-002 -2.752e-003 
 4.236e-004 8.263e-001 8.835e-002 -2.539e-003 
 4.231e-004 8.263e-001 1.281e-001 -2.826e-003 
 4.225e-004 8.263e-001 1.364e-001 -2.542e-003 
 4.220e-004 8.263e-001 1.109e-001 -1.759e-003 
 4.217e-004 8.263e-001 8.961e-002 -1.193e-003 
 4.215e-004 8.263e-001 7.204e-002 -7.915e-004 
 4.213e-004 8.263e-001 5.760e-002 -5.170e-004 
 4.212e-004 8.263e-001 4.593e-002 -3.358e-004 
 4.212e-004 8.263e-001 3.652e-002 -2.158e-004 
 4.211e-004 8.263e-001 2.896e-002 -1.392e-004 
 4.211e-004 8.263e-001 2.296e-002 -9.015e-005 
 4.211e-004 8.263e-001 1.819e-002 -5.917e-005 

Filter calculated successfully

h=fvtool(b,a);
axis([0 1 -90 5])
title('Magnitude of Frequency Response (dB)')

Figure 5.4(a)

h=fvtool(b,a);
title('Passband Magnitude Response (dB)')
axis([0 0.4 -0.1 0.1])

Figure 5.4(b)

N=8192;
w=pi*linspace(0,1,N);
H=freqz(b,a,w);
plot(w/pi,phase(H)-phase(exp(-1i*pi/N*(0:N-1)*tau)));
grid on;
title('Passband Phase Error (radians)');
xlabel('w/pi');
axis([0 0.4 -0.015 0.015]);

Figure 5.4(c)

h=fvtool(b,a);
set(h,'Analysis','grpdelay')
title('Passband Group Delay (samples)')
axis([0 0.4 14.8 15.5])

Figure 5.4(d)

h=fvtool(b,a);
set(h,'Analysis','polezero')
title('Poles (x) and Zeros (o)')

Figure 5.5(a)

Example 3: Highpass Filter with Approximately Linear Passband Phase (M=14, N=6)

M=14;
N=6;
tau=12;
om=pi*[linspace(0,0.475,50),linspace(0.525,1,50)];
D=[zeros(1,50),exp(-1i*om(51:100)*tau)];
W=ones(1,100);
[b,a,e]=mpiir_l2(M,N,om,D,W,0.9276);
   ERROR    MAX.RADIUS    STEP       SLOPE
 1.718e+000 7.097e-001 4.997e-001 -1.565e+000 
 7.082e-001 7.469e-001 7.790e-001 -9.926e-001 
 4.437e-001 7.494e-001 4.708e-001 -5.530e-001 
 3.456e-001 7.965e-001 5.654e-001 -3.702e-001 
 2.794e-001 8.745e-001 8.098e-001 -5.245e-001 
 2.049e-001 8.655e-001 6.442e-001 -6.088e-001 
 1.546e-001 8.552e-001 2.368e-001 -3.722e-001 
 1.302e-001 8.685e-001 2.773e-001 -2.147e-001 
 1.174e-001 8.769e-001 3.462e-001 -1.482e-001 
 1.089e-001 8.839e-001 4.610e-001 -1.380e-001 
 1.012e-001 8.911e-001 5.714e-001 -1.783e-001 
 9.151e-002 9.010e-001 6.456e-001 -2.902e-001 
 7.697e-002 9.112e-001 5.712e-001 -3.926e-001 
 6.183e-002 9.188e-001 2.386e-001 -2.786e-001 
 5.359e-002 9.231e-001 8.965e-002 -1.500e-001 
 4.970e-002 9.253e-001 3.429e-002 -7.564e-002 
 4.786e-002 9.264e-001 1.320e-002 -3.813e-002 
 4.696e-002 9.270e-001 4.937e-003 -1.926e-002 
 4.651e-002 9.273e-001 1.797e-003 -9.742e-003 
 4.628e-002 9.274e-001 6.476e-004 -4.925e-003 
 4.617e-002 9.275e-001 2.488e-004 -2.489e-003 
 4.611e-002 9.276e-001 1.260e-004 -1.259e-003 
 4.608e-002 9.276e-001 8.095e-005 -6.366e-004 
 4.607e-002 9.276e-001 4.814e-005 -3.228e-004 
 4.606e-002 9.276e-001 2.892e-005 -1.633e-004 
 4.606e-002 9.276e-001 1.709e-005 -8.263e-005 

Filter calculated successfully

h=fvtool(b,a);
axis([0 1 -50 5])
title('Magnitude of Frequency Response (dB)')

Figure 5.6(a)

h=fvtool(b,a);
title('Passband Magnitude Response (dB)')
axis([0.53 1 -0.2 0.2])

Figure 5.6(b)

N=8192;
w=pi*linspace(0,1,N);
H=freqz(b,a,w);
plot(w/pi,mod(phase(H)-phase(exp(-1i*pi/N*(0:N-1)*tau))+pi,2*pi)-pi);
grid on;
title('Passband Phase Error (radians)');
xlabel('w/pi');
axis([0.5 1 -0.06 0.06]);

Figure 5.6(c)

h=fvtool(b,a);
set(h,'Analysis','grpdelay')
title('Passband Group Delay (samples)')
axis([0.53 1 10 14])

Figure 5.6(d)

h=fvtool(b,a);
set(h,'Analysis','polezero')
title('Poles (x) and Zeros (o)')
axis([-1.5 1.5 -1.5 1.5])

Figure 5.7(a)

Example 4: Narrow Bandpass Filter with Approximately Linear Passband Phase (M=20, N=8)

M=20;
N=8;
tau=20;
om=pi*[linspace(0,0.36,36),linspace(0.4,0.5,10),linspace(0.54,1,46)];
D=[zeros(1,36),exp(-1i*om(37:46)*tau),zeros(1,46)];
W=[100*ones(1,36),ones(1,10),100*ones(1,46)];
[b,a,e]=mpiir_l2(M,N,om,D,W,0.98);
   ERROR    MAX.RADIUS    STEP       SLOPE
 9.070e+000 8.877e-001 1.979e+000 -2.872e-001 
 7.212e+000 9.403e-001 8.745e-001 -5.000e-001 
 5.251e+000 9.577e-001 8.267e-001 -4.071e-001 
 4.342e+000 9.505e-001 1.184e+000 -1.605e-001 
 3.969e+000 9.400e-001 1.164e+000 -2.535e-001 
 3.280e+000 9.216e-001 7.838e-001 -5.230e-001 
 2.153e+000 9.528e-001 3.441e-001 -7.251e-001 
 1.508e+000 9.423e-001 2.206e-001 -2.737e-001 
 1.350e+000 9.436e-001 2.837e-001 -1.233e-001 
 1.282e+000 9.455e-001 3.225e-001 -6.556e-002 
 1.245e+000 9.454e-001 4.222e-001 -4.978e-002 
 1.215e+000 9.457e-001 5.883e-001 -6.161e-002 
 1.174e+000 9.442e-001 7.407e-001 -1.057e-001 
 1.098e+000 9.547e-001 7.887e-001 -2.360e-001 
 9.064e-001 9.693e-001 4.477e-001 -4.362e-001 
 6.435e-001 9.750e-001 2.132e-001 -6.106e-001 
 4.383e-001 9.710e-001 5.675e-001 -4.661e-001 
 3.507e-001 9.690e-001 5.604e-001 -3.447e-001 
 3.067e-001 9.699e-001 4.295e-001 -2.198e-001 
 2.814e-001 9.710e-001 2.610e-001 -1.130e-001 
 2.694e-001 9.727e-001 1.663e-001 -5.685e-002 
 2.634e-001 9.733e-001 1.107e-001 -2.854e-002 
 2.603e-001 9.740e-001 1.134e-001 -1.851e-002 
 2.581e-001 9.743e-001 1.644e-001 -1.727e-002 
 2.559e-001 9.746e-001 2.393e-001 -2.226e-002 
 2.530e-001 9.748e-001 3.458e-001 -3.366e-002 
 2.484e-001 9.750e-001 4.671e-001 -5.595e-002 
 2.408e-001 9.753e-001 5.751e-001 -1.014e-001 
 2.265e-001 9.759e-001 6.326e-001 -2.086e-001 
 1.962e-001 9.763e-001 4.803e-001 -3.611e-001 
 1.535e-001 9.773e-001 2.031e-001 -3.467e-001 
 1.256e-001 9.782e-001 8.498e-002 -2.277e-001 
 1.113e-001 9.789e-001 3.959e-002 -1.332e-001 
 1.039e-001 9.794e-001 1.972e-002 -7.519e-002 
 1.000e-001 9.796e-001 1.030e-002 -4.193e-002 
 9.795e-002 9.798e-001 5.462e-003 -2.325e-002 
 9.681e-002 9.799e-001 2.995e-003 -1.285e-002 
 9.619e-002 9.799e-001 1.642e-003 -7.083e-003 
 9.585e-002 9.800e-001 9.394e-004 -3.919e-003 
 9.566e-002 9.800e-001 5.351e-004 -2.152e-003 
 9.556e-002 9.800e-001 3.037e-004 -1.191e-003 
 9.550e-002 9.800e-001 1.822e-004 -6.597e-004 
 9.547e-002 9.800e-001 1.039e-004 -3.575e-004 
 9.546e-002 9.800e-001 6.150e-005 -2.025e-004 
 9.545e-002 9.800e-001 3.701e-005 -1.147e-004 
 9.544e-002 9.800e-001 2.236e-005 -6.450e-005 

Filter calculated successfully

h=fvtool(b,a);
axis([0 1 -80 5])
title('Magnitude of Frequency Response (dB)')

Figure 5.9(a)

h=fvtool(b,a);
title('Passband Magnitude Response (dB)')
axis([0.53 1 -0.2 0.2])

Figure 5.9(b)

N=8192;
w=pi*linspace(0,1,N);
H=freqz(b,a,w);
plot(w/pi,mod(phase(H)-phase(exp(-1i*pi/N*(0:N-1)*tau))+pi,2*pi)-pi);
grid on;
title('Passband Phase Error (radians)');
xlabel('w/pi');
axis([0.4 0.5 -0.15 0.15]);

Figure 5.9(c)

Example 5: IIR Equalization and Anti-Aliasing Filter (M=20, N=4)

In this example we want an IIR filter that equalizes circuit representing analog filter of order 3.

tau=35;M=20;N=4;r=0.97;
om=pi*[linspace(0,1/16,20),linspace(3/16,1,260)];
Dt=[exp(-1i*om(1:20)*tau),zeros(1,260)];
Wt=[ones(1,20),15*ones(1,260)];
z=[];p=[-0.6493,1/(-0.3246-1i*1.0325),1/(-0.3246+1i*1.0325)];k=0.7606;
[num,den]=zp2tf(z,p,k);
Ha=freqs(num,den,om*16/pi);
D=Dt./Ha;W=Wt.*abs(Ha).^2;
[b,a]=mpiir_l2(M,N,om,D,W,r);

% plot graphs
N=8192;
w=pi*linspace(0,1,N);
Ha=freqs(num,den,w*16/pi);
Hd=freqz(b,a,w);
H=Ha.*Hd;
Hlog=20*log10(abs(H));
plot(w/pi,Hlog);
grid on;
title('Magnitude of Frequency Response (dB)');
xlabel('w/pi');
axis([0 1 -100 5]);
   ERROR    MAX.RADIUS    STEP       SLOPE
 5.272e-001 7.944e-001 2.035e+000 -1.027e+000 
 2.548e-001 8.940e-001 2.714e+000 -1.409e+000 
 9.807e-002 9.315e-001 8.653e-001 -6.334e-001 
 7.564e-002 9.505e-001 2.588e-001 -1.315e-001 
 7.208e-002 9.600e-001 1.282e-001 -2.714e-002 
 7.137e-002 9.546e-001 1.203e-001 -5.948e-003 
 7.121e-002 9.462e-001 1.757e-001 -1.397e-003 
 7.117e-002 9.432e-001 2.678e-001 -4.214e-004 
 7.116e-002 9.565e-001 4.004e-001 -2.592e-004 
 7.115e-002 9.635e-001 4.992e-001 -3.540e-004 
 7.113e-002 9.669e-001 5.267e-001 -7.970e-004 
 7.109e-002 9.656e-001 4.864e-001 -2.824e-003 
 7.086e-002 9.677e-001 4.705e-001 -1.797e-002 
 6.901e-002 9.683e-001 8.731e-001 -1.447e-001 
 5.953e-002 9.671e-001 4.124e-001 -4.421e-001 
 4.383e-002 9.688e-001 5.379e-001 -3.342e-001 
 3.695e-002 9.695e-001 7.251e-001 -2.316e-001 
 3.148e-002 9.698e-001 5.768e-001 -3.323e-001 
 2.370e-002 9.699e-001 4.071e-001 -5.273e-001 
 1.599e-002 9.700e-001 2.501e-001 -6.054e-001 
 1.101e-002 9.699e-001 1.769e-001 -4.578e-001 
 8.701e-003 9.699e-001 1.765e-001 -2.475e-001 
 7.752e-003 9.699e-001 9.819e-002 -1.414e-001 
 7.269e-003 9.698e-001 6.326e-002 -7.730e-002 
 7.023e-003 9.697e-001 4.177e-002 -4.192e-002 
 6.894e-003 9.696e-001 2.825e-002 -2.275e-002 
 6.826e-003 9.693e-001 1.949e-002 -1.240e-002 
 6.789e-003 9.690e-001 1.367e-002 -6.806e-003 
 6.768e-003 9.685e-001 9.726e-003 -3.799e-003 
 6.757e-003 9.678e-001 7.115e-003 -2.204e-003 
 6.751e-003 9.667e-001 5.701e-003 -1.398e-003 
 6.746e-003 9.650e-001 5.728e-003 -1.052e-003 
 6.743e-003 9.625e-001 7.408e-003 -1.006e-003 
 6.740e-003 9.590e-001 1.076e-002 -1.192e-003 
 6.736e-003 9.590e-001 1.608e-002 -1.610e-003 
 6.731e-003 9.590e-001 2.415e-002 -2.311e-003 
 6.723e-003 9.590e-001 3.628e-002 -3.403e-003 
 6.711e-003 9.590e-001 5.459e-002 -5.075e-003 
 6.694e-003 9.590e-001 8.233e-002 -7.648e-003 
 6.668e-003 9.589e-001 1.245e-001 -1.167e-002 
 6.629e-003 9.588e-001 1.883e-001 -1.809e-002 
 6.567e-003 9.586e-001 2.823e-001 -2.876e-002 
 6.468e-003 9.582e-001 4.092e-001 -4.743e-002 
 6.303e-003 9.575e-001 5.497e-001 -8.313e-002 
 6.006e-003 9.560e-001 6.614e-001 -1.624e-001 
 5.398e-003 9.506e-001 7.314e-001 -3.908e-001 
 3.795e-003 9.367e-001 3.537e-001 -7.723e-001 
 1.960e-003 9.239e-001 1.574e-001 -1.160e+000 
 7.937e-004 9.474e-001 1.554e-001 -1.394e+000 
 2.947e-004 9.588e-001 1.346e-001 -1.397e+000 
 1.251e-004 9.643e-001 8.828e-002 -1.006e+000 
 7.678e-005 9.665e-001 4.355e-002 -3.690e-001 
 6.600e-005 9.675e-001 2.204e-002 -1.079e-001 
 6.329e-005 9.680e-001 1.148e-002 -2.971e-002 
 6.258e-005 9.683e-001 5.929e-003 -8.151e-003 
 6.238e-005 9.684e-001 3.148e-003 -2.254e-003 
 6.233e-005 9.684e-001 1.710e-003 -6.286e-004 
 6.231e-005 9.685e-001 9.052e-004 -1.763e-004 
 6.231e-005 9.685e-001 4.846e-004 -4.969e-005 
 6.231e-005 9.685e-001 2.611e-004 -1.404e-005 
 6.231e-005 9.685e-001 1.410e-004 -3.977e-006 
 6.231e-005 9.685e-001 7.479e-005 -1.128e-006 

Filter calculated successfully

Figure 5.11(a)

plot(w/pi,Hlog);
grid on;
title('Passband Magnitude (dB)');
xlabel('w/pi');
axis([0 1/16 -0.03 0.03]);

Figure 5.11(b)

plot(w/pi,phase(H)-phase(exp(-1i*pi/N*(0:N-1)*tau)));
grid on;
title('Passband Phase Error (radians)');
xlabel('w/pi');
axis([0 1/16 -0.004 0.004]);

Figure 5.11(c)

plot(w(1:end-1)/pi,-diff(phase(H))*N/pi);
grid on;
title('Passband Group Delay (samples)');
xlabel('w/pi');
axis([0 1/16 34.55 35.55]);

Figure 5.11(d)

Example 6: Comparison to a cascade of Cauer Lowpass Filter and and Phase Equalizer (M=10, N=6)

In this example we design Lowpass filter with approximately linear phase in passband (instead of combination of Cauer Lowpass Filter and and Phase Equalizer)

M=10;N=6;
tau=8.5;
om=pi*[linspace(0,0.5,50),linspace(0.6,1,40)];
D=0.5*(1+10^(-0.5/20))*[exp(-1i*om(1:50)*tau) zeros(1,40)];
W=[ones(1,50) 10*ones(1,40)];
[b,a]=mpiir_l2(M,N,om,D,W,0.92);
   ERROR    MAX.RADIUS    STEP       SLOPE
 2.304e+000 6.960e-001 6.943e-001 -1.601e+000 
 9.306e-001 7.448e-001 7.267e-001 -9.340e-001 
 5.458e-001 8.497e-001 8.356e-001 -1.021e+000 
 3.130e-001 8.205e-001 2.937e-001 -6.235e-001 
 2.285e-001 8.224e-001 2.930e-001 -3.595e-001 
 1.922e-001 8.306e-001 3.820e-001 -2.313e-001 
 1.714e-001 8.363e-001 4.963e-001 -2.155e-001 
 1.528e-001 8.422e-001 5.988e-001 -2.955e-001 
 1.292e-001 8.547e-001 6.603e-001 -5.051e-001 
 9.606e-002 8.796e-001 5.853e-001 -7.766e-001 
 6.357e-002 8.985e-001 2.099e-001 -6.035e-001 
 4.764e-002 9.072e-001 3.368e-002 -3.132e-001 
 4.125e-002 9.129e-001 2.672e-002 -1.513e-001 
 3.852e-002 9.161e-001 1.385e-002 -6.679e-002 
 3.737e-002 9.179e-001 7.138e-003 -3.017e-002 
 3.684e-002 9.188e-001 3.689e-003 -1.429e-002 
 3.659e-002 9.194e-001 1.947e-003 -7.084e-003 
 3.647e-002 9.197e-001 1.017e-003 -3.634e-003 
 3.640e-002 9.198e-001 5.543e-004 -1.910e-003 
 3.637e-002 9.199e-001 2.925e-004 -1.015e-003 
 3.635e-002 9.199e-001 1.601e-004 -5.452e-004 
 3.634e-002 9.200e-001 8.702e-005 -2.934e-004 
 3.634e-002 9.200e-001 4.606e-005 -1.581e-004 
 3.633e-002 9.200e-001 2.408e-005 -8.557e-005 

Filter calculated successfully

h=fvtool(b,a);
axis([0 1 -60 5])
title('Magnitude of Frequency Response (dB)')

Figure 5.12(a)

h=fvtool(b,a);
title('Passband Magnitude Response (dB)')
axis([0 0.55 -0.55 0.05])

Figure 5.12(b)

N=8192;
w=pi*linspace(0,1,N);
H=freqz(b,a,w);
plot(w/pi,phase(H)-phase(exp(-1i*pi/N*(0:N-1)*tau)));
grid on;
title('Passband Phase Error (radians)');
xlabel('w/pi');
axis([0 0.5 -0.06 0.06]);

Figure 5.12(c)

h=fvtool(b,a);
set(h,'Analysis','grpdelay')
title('Passband Group Delay (dB)')
axis([0 0.5 7 17])

Figure 5.12(d)

h=fvtool(b,a);
set(h,'Analysis','polezero')
title('Poles (x) and Zeros (o)')

Figure 5.13(a)

Example 7: Comparison to a Coupled Allpass Filters (M=40, N=6)

tau=37;M=40;N=6;r=0.98;
om=pi*[linspace(0,0.2,20),linspace(0.23,1,77)];
D=[exp(-1i*om(1:20)*tau),zeros(1,77)];
W=[ones(1,20),1000*ones(1,77)];
[b,a]=mpiir_l2(M,N,om,D,W,r);
   ERROR    MAX.RADIUS    STEP       SLOPE
 1.363e+001 7.775e-001 1.693e+000 -9.622e-001 
 8.554e+000 8.297e-001 1.145e+000 -8.074e-001 
 5.415e+000 8.506e-001 1.075e+000 -3.995e-001 
 4.387e+000 8.336e-001 1.159e+000 -3.344e-001 
 3.486e+000 9.026e-001 9.075e-001 -5.945e-001 
 2.389e+000 8.835e-001 4.135e-001 -4.233e-001 
 1.907e+000 8.881e-001 5.040e-001 -2.415e-001 
 1.683e+000 8.934e-001 6.900e-001 -1.832e-001 
 1.514e+000 8.956e-001 8.093e-001 -2.291e-001 
 1.292e+000 9.175e-001 7.254e-001 -4.261e-001 
 9.432e-001 9.303e-001 3.827e-001 -5.274e-001 
 7.177e-001 9.277e-001 2.539e-001 -3.545e-001 
 6.166e-001 9.353e-001 3.707e-001 -2.163e-001 
 5.515e-001 9.386e-001 5.067e-001 -1.516e-001 
 5.072e-001 9.417e-001 6.387e-001 -1.593e-001 
 4.592e-001 9.456e-001 7.111e-001 -2.679e-001 
 3.751e-001 9.525e-001 4.784e-001 -4.547e-001 
 2.706e-001 9.623e-001 2.306e-001 -5.589e-001 
 1.959e-001 9.644e-001 1.129e-001 -4.442e-001 
 1.588e-001 9.690e-001 1.719e-001 -3.072e-001 
 1.357e-001 9.716e-001 2.406e-001 -1.905e-001 
 1.236e-001 9.735e-001 3.456e-001 -1.251e-001 
 1.162e-001 9.749e-001 4.707e-001 -1.032e-001 
 1.100e-001 9.760e-001 5.841e-001 -1.252e-001 
 1.023e-001 9.771e-001 6.549e-001 -2.139e-001 
 8.891e-002 9.773e-001 6.519e-001 -4.174e-001 
 6.599e-002 9.781e-001 2.406e-001 -4.085e-001 
 5.274e-002 9.789e-001 6.837e-002 -2.142e-001 
 4.744e-002 9.793e-001 2.391e-002 -1.079e-001 
 4.499e-002 9.796e-001 9.436e-003 -5.776e-002 
 4.373e-002 9.798e-001 3.837e-003 -3.187e-002 
 4.304e-002 9.799e-001 1.535e-003 -1.780e-002 
 4.266e-002 9.799e-001 5.682e-004 -1.003e-002 
 4.245e-002 9.800e-001 2.549e-004 -5.789e-003 
 4.233e-002 9.800e-001 1.102e-004 -3.309e-003 
 4.226e-002 9.800e-001 3.729e-005 -1.873e-003 
 4.222e-002 9.800e-001 2.365e-005 -1.081e-003 
 4.219e-002 9.800e-001 8.230e-006 -5.939e-004 
 4.218e-002 9.800e-001 4.498e-006 -3.248e-004 
 4.218e-002 9.800e-001 3.297e-006 -1.811e-004 
 4.217e-002 9.800e-001 2.054e-006 -8.711e-005 

Filter calculated successfully

h=fvtool(b,a);
axis([0 1 -80 5])
title('Magnitude of Frequency Response (dB)')

Figure 5.14(a)

h=fvtool(b,a);
title('Passband Magnitude Response (dB)')
axis([0 0.22 -1.5 0.5])

Figure 5.14(b)

N=8192;
w=pi*linspace(0,1,N);
H=freqz(b,a,w);
plot(w/pi,phase(H)-phase(exp(-1i*pi/N*(0:N-1)*tau)));
grid on;
title('Passband Phase Error (radians)');
xlabel('w/pi');
axis([0 0.22 -0.2 0.6]);

Figure 5.12(c)

h=fvtool(b,a);
set(h,'Analysis','grpdelay')
title('Passband Group Delay (dB)')
axis([0 0.2 20 50])

Figure 5.14(d)

h=fvtool(b,a);
set(h,'Analysis','polezero')
title('Poles (x) and Zeros (o)')

Figure 5.15(a)