function model = tu1_simulink(num,den,stopTime)
% TU1_SIMULINK Generate a Simulink model for a unit-step experiment.
%   model = tu1_simulink([1 3],[1 6 5],12)
%   Compare the Transfer Fcn and State-Space realizations at x0 = 0.
%   Requires Simulink and Control System Toolbox. Saves a new .slx file.
arguments
    num (1,:) double {mustBeFinite} = [1 3]
    den (1,:) double {mustBeFinite} = [1 6 5]
    stopTime (1,1) double {mustBePositive,mustBeFinite} = 12
end
assert(numel(den) >= 2 && den(1) ~= 0 && numel(num) <= numel(den), ...
    'Use a proper transfer function with a nonconstant denominator.');
[A,B,C,D] = ssdata(ss(tf(num,den)));
stamp = char(datetime('now','Format','yyyyMMdd_HHmmss_SSS'));
model = ['TU1_Comparison_' stamp];
new_system(model);
set_param(model,'Solver','ode45','StopTime',num2str(stopTime),'MaxStep',num2str(stopTime/500));
add_block('simulink/Sources/Step',[model '/Step'], ...
    'Time','0','Before','0','After','1','Position',[40 109 70 141]);
add_block('simulink/Continuous/Transfer Fcn',[model '/Transfer Fcn'], ...
    'Numerator',mat2str(num),'Denominator',mat2str(den),'Position',[170 62 315 102]);
add_block('simulink/Continuous/State-Space',[model '/State-Space'], ...
    'A',mat2str(A),'B',mat2str(B),'C',mat2str(C),'D',mat2str(D), ...
    'X0',mat2str(zeros(size(A,1),1)),'Position',[170 150 315 190]);
add_block('simulink/Signal Routing/Mux',[model '/Mux'], ...
    'Inputs','2','Position',[393 85 399 167]);
add_block('simulink/Sinks/Scope',[model '/Scope'],'Position',[473 108 513 144]);
add_line(model,'Step/1','Transfer Fcn/1','autorouting','on');
add_line(model,'Step/1','State-Space/1','autorouting','on');
add_line(model,'Transfer Fcn/1','Mux/1','autorouting','on');
add_line(model,'State-Space/1','Mux/2','autorouting','on');
add_line(model,'Mux/1','Scope/1','autorouting','on');
save_system(model);
open_system(model);
end