function runsimulationcallback

% CC Violeta Ivanova, violeta@mit.edu.
% RUNSIMULATIONCALLBACK is a callback function for the Run Simulation push button
% in the SYSTEMS123 program. The callback starts the Monte Carlo simulations 
% of three supermarket checkout systems.
% Developed to support course material in 1.010 Uncertainty in Engineering.
% Faculty: Daniele Veneziano. TA: Andreas Langousis.

fig1 = openfig('simulationsGUI.fig', 'reuse');
figure(fig1)

disp('Running Monte Carlo simulations ...');

N_counters_popup      = findobj( fig1, 'Tag', 'N_counters_popup');
N_states_popup        = findobj( fig1, 'Tag', 'N_states_popup');
T_arrival_text        = findobj( fig1, 'Tag', 'T_arrival_text');
T_service_text        = findobj( fig1, 'Tag', 'T_service_text');
accuracy_popup        = findobj( fig1, 'Tag', 'accuracy_popup');

% Get arrival and service times from GUI
arrivalstr            = get( T_arrival_text,  'String');
servicestr            = get( T_service_text,  'String');

% Set mean arrival time and mean service time for the simulations
t_arrival         = str2num ( arrivalstr );
t_service         = str2num ( servicestr );

% Get number of open counters for simulations from GUI
N_counters           = get( N_counters_popup,   'Value' );   

% Get number of simulation states from GUI
N_states              = get( N_states_popup,     'Value' ); 

% Set number of simulation states N_sim for computations

switch N_states
    case 1, N_sim = 1500;
    case 2, N_sim = 5000;
    case 3, N_sim = 10000;
    case 4, N_sim = 15000;
    case 5, N_sim = 20000;
    case 6, N_sim = 25000;
    case 7, N_sim = 30000;
    case 8, N_sim = 40000;
    case 9, N_sim = 50000;
end

% Get accuracy from GUI
accuracy             = get( accuracy_popup,     'Value' ); 

% Set accuracy accu for computations

switch accuracy
    case 1, accu = 0.01;
    case 2, accu = 0.005;
    case 3, accu = 0.001;
    case 4, accu = 0.0005;
    case 5, accu = 0.0001;
end

[aveNwait1, lastNwait1, avwaitT1] = system1( N_counters, t_arrival, t_service, N_sim, accu);
[aveNwait2, lastNwait2, avwaitT2] = system2( N_counters, t_arrival, t_service, N_sim, accu);
[aveNwait3, lastNwait3, avwaitT3] = system3( N_counters, t_arrival, t_service, N_sim, accu);

% Plotting simulation results
disp('Plotting simulation results, please wait ...');
figure(2);
set(2, 'Position', [15 35 1244 800]); 
titlestr = ['Number of Open Counters: ' num2str(N_counters) ...
            ', Mean Customer Interarrival Time: ' num2str(t_arrival) ...
            ' min, Mean Service Time per Customer: ' num2str(t_service) ...
            ' min, Max. Number of Simulation Steps: ' num2str(N_sim), ...
            ', c = ', num2str(N_counters*t_arrival/t_service)];
set(2, 'Name', titlestr, 'Color', [0.95 0.87 0.73], 'NumberTitle', 'off');

% System 1 
s1 = subplot(2, 3, 1);
     p1 = plot(avwaitT1(:,2), avwaitT1(:,3), ...
                'Color', [0. 0.5 0.], ...
                'LineWidth', [2.0]);
     title('SYSTEM 1: Random Allocation', ...
                'FontWeight', 'bold', ...
                'FontSize', 16);
     xlabel('Simulation Time [min]', 'FontSize', 12); 
     ylabel('Mean Customer Waiting Time [min]', 'FontSize', 12);
     axes1 = gca;
     ylim1 = get(axes1, 'YLim');
     xlim1 = get(axes1, 'XLim');

% System 2
s2 = subplot(2, 3, 2);
     p2 = plot(avwaitT2(:,2), avwaitT2(:,3), ...
                'Color', [0. 0.5 0.], ...
                'LineWidth', [2.0]);
     title('SYSTEM 2: Shortest Queue', ...
                'FontWeight', 'bold', ...
                'FontSize', 16);
     xlabel('Simulation Time [min]', 'FontSize', 12); 
     ylabel('Mean Customer Waiting Time [min]', 'FontSize', 12);
     axes2 = gca;
     ylim2 = get(axes2, 'YLim');
     xlim2 = get(axes2, 'XLim');
     
% System 3
s3 = subplot(2, 3, 3);
     p3 = plot(avwaitT3(:,2), avwaitT3(:,3), ...
                'Color', [0. 0.5 0.], ...
                'LineWidth', [2.0]);
     title('SYSTEM 3: Single Queue', ...
                'FontWeight', 'bold', ...
                'FontSize', 16);
     xlabel('Simulation Time [min]', 'FontSize', 12); 
     ylabel('Mean Customer Waiting Time [min]', 'FontSize', 12);
     axes3 = gca;
     ylim3 = get(axes3, 'YLim');
     xlim3 = get(axes3, 'XLim');
     
% Make axes equal for subplots 1, 2, 3
     xmax = max([xlim1 xlim2 xlim3]);
     ymax = max([ylim1 ylim2 ylim3]);
     set (axes1, 'XLim', [0 xmax], 'YLim', [0 ymax]);
     set (axes2, 'XLim', [0 xmax], 'YLim', [0 ymax]);
     set (axes3, 'XLim', [0 xmax], 'YLim', [0 ymax]);
     
% System 1 
s4 = subplot(2, 3, 4);
     customers1 = [aveNwait1 lastNwait1];
     p4 = barh(customers1);
     title('System 1: Mean and St. Dev. of Waiting Customers', ...
                'FontWeight', 'bold', ...
                'FontSize', 12);
     ylabel('Counter #', 'FontSize', 12); 
     xlabel('Number of Customers (n)', 'FontSize', 12);
     % uncomment next line to get legend on subplot 4 (same as subplot 6)
     % legend(p4, 'N average', 'N last state');
     axes4 = gca;
     xlim4 = get(axes4, 'XLim');
     
% System 2 
s5 = subplot(2, 3, 5);
     customers2 = [aveNwait2 lastNwait2];
     p5 = barh(customers2);
     title('System 2: Mean and St. Dev. of Waiting Customers', ...
                'FontWeight', 'bold', ...
                'FontSize', 12);
     ylabel('Counter #', 'FontSize', 12); 
     xlabel('Number of Customers (n)', 'FontSize', 12);
     % uncomment next line to get legend on subplot 5 (same as subplot 6)
     % legend(p5, 'N average', 'N last state');
     axes5 = gca;
     xlim5 = get(axes5, 'XLim');
     
% System 3 
s6 = subplot(2, 3, 6);
     customers3 = [aveNwait3 lastNwait3];
     p6 = barh(customers3);
      title('System 3: Mean and St. Dev. of Waiting Customers', ...
                'FontWeight', 'bold', ...
                'FontSize', 12);
     ylabel('Single Queue', 'FontSize', 12); 
     xlabel('Number of Customers (n)', 'FontSize', 12);
     % same legend for subplots 4, 5, 6
     legend(p6, 'Mean', 'St.Dev.');
     axes6 = gca;
     xlim6 = get(axes6, 'XLim');
     set(axes6, 'YTick', []);
     
% Make axes equal for subplots 4, 5, 6
     
     xlim = max([xlim4 xlim5 xlim6]);
     set (axes4, 'XLim', [0 xlim]);
     set (axes5, 'XLim', [0 xlim]);
     set (axes6, 'XLim', [0 xlim]);
     
     xticks = [0 : xlim];
     set (axes4, 'XTick', xticks);
     set (axes5, 'XTick', xticks);
     set (axes6, 'XTick', xticks);
     
disp('Simulation completed.');

return