% Set up the parameters
imacmic_numpar;
% Range of damping for the macro controller
D_1 = [0:0.1:2];
% Compute the damping factor d and save in arrray DD
% -- uses zeros of system
% $$$ DD = [];
% $$$ for d_1 = D_1
% $$$ [A,B,C,D] = macmic_dm(m_1,m_2,b_1,b_2,k_2,k_t,b_t,b_3,m_3,p_1,d_1);
% $$$ zz = ss2zp(A,B,C,D,1);
% $$$ z = (zz(:,2));
% $$$ re = sort(real(z));
% $$$ im = sort(imag(z));
% $$$ d = -re./sqrt(re.^2 + im.^2);
% $$$ d = sort(d);
% $$$ DD = [DD d];
% $$$ end;
% $$$
% $$$ plot(D_1,min(DD)); grid
% $$$ xlabel('Derivative gain');
% $$$ ylabel('Damping factor');
% $$$
% $$$
% Compute the damping factor d and save in arrray DD
% -- uses poles of inverse system
DD1 = [];
for d_1 = D_1
[A] = imacmic_sm;
z = eig(A);
re = (real(z));
im = (imag(z));
d = -re./sqrt(re.^2 + im.^2);
d = sort(d);
DD1 = [DD1 d];
end;
minDD1 = min(DD1);
plot(D_1,minDD1); grid
xlabel('Derivative gain');
ylabel('Damping factor');
Maximum_Damping = max(minDD1)
index = minDD1 == Maximum_Damping*ones(size(minDD1));
Opt_gain = index*D_1'
% Plot on disk
set term postscript
set output "damping.ps"
plot(D_1,minDD1); grid