[num, den] = pade (0.06, 1);
delay = tf (num, den);
a = ureal('a', 0.5, 'Range', [0 1]);
sys = tf (1, [1 0 -a]);

wk = 0.024 / ellipfilt (4, 1.6, 40, 39);
sigma(wk)
hold on
plot ([50 50 10000], [-40 6 6], 'r--')
hold off
% gain 2 above 50 rad/s
pause

disp ('Perform DK iterations with complex uncertainties')
pause
% gain/phase uncertainty
margin = ucomplex ('margin', 0);

gain = lft (margin, ([1/3 1; 2/3 1]))

g = gain * sys * delay;
gtot = [0 1; wk g];

opts = dksynOptions ('DisplayWhileAutoIter', 'on');
[kc1, clc, gamc] = dksyn (gtot, 1, 1, opts);

pause
disp ('Perform DK iterations with mixed uncertainties')
pause

opts = dksynOptions ('DisplayWhileAutoIter', 'on', ...
                     'MixedMU', 'on');

[kr1, clr, gamr] = dksyn (gtot, 1, 1, opts);

pause

disp ('modified design, corresponds to a (complex) D-scaling on w/z');
pause

gtot = [0 wk; 1 g];

opts = dksynOptions ('DisplayWhileAutoIter','on');
[kc2, clc, gamc] = dksyn (gtot, 1, 1, opts);

pause
disp ('modified design, corresponds to a (real) D-scaling on w/z');
pause

opts = dksynOptions ('DisplayWhileAutoIter','on', 'MixedMU', 'on');
[kr2, clr, gamr] = dksyn (gtot, 1, 1, opts);
pause

step (feedback (sys * delay, -kr2))
pause
disp ('Reduce the order of the controller')
pause

% model reduction of kr
[krb, sig] = balreal (kr2);
sig'
% keep sig > 1
krr = modred (krb, sig < 1);
size(krr)
sigma (kr2, krr)
hold on
plot ([50 50 10000], [40 -6 -6], 'g--')
hold off
