function [rxgp, times_executed, event_log, ssdata, aglog] = coordinate_erg(rxgp,start_time,timestep, run_info, update, eventlog_arraysize, snapshot_interval)

modrxn = rxgp.rxn;
moddepend = rxgp.rxn_dependency_matrix;
numrxns = length(rxgp.rxn);

% because I'm too lazy to change the parameters passed into coordinate_erg
% right now...
if eventlog_arraysize < 0
    evlog_flag = 0;
    eventlog_arraysize = -eventlog_arraysize;
else
    evlog_flag = 1;
end

clear ssdata
clear aglog

% need to fix up the difference between c 1st element (index of 0) and
% matlab first element (index of 1)
for i = 1:numrxns
    modrxn{i}.substrates(:,1) = rxgp.rxn{i}.substrates(:,1) - 1;  
    if isfield(modrxn{i}, 'products')
        modrxn{i}.products(:,1) = rxgp.rxn{i}.products(:,1) - 1;
    end
    
    moddepend{i} = rxgp.rxn_dependency_matrix{i} - 1;
    
    if isfield(modrxn{i}, 'StoP')    
        modrxn{i} = rmfield(modrxn{i},'StoP');
    end

    modrxn{i}.StoP = rateparms_to_vector(rxgp.rxn{i});
end

% also need to deal with indexing issue for compounds to monitor
if ~isempty(rxgp.monitor_comps)
    for i = 1:length(rxgp.monitor_comps)
        mod_monitor(i) = rxgp.monitor_comps(i) - 1;
    end
else
    mod_monitor = [];
end

% create special reaction dependency that is used to update all reactions
% only used if update = 1
moddepend{numrxns+1} = [0:numrxns-1];

% Create structure with run_constants
%         Slava     05/20/07  

run_constants.start_time = start_time;
run_constants.timestep = timestep;
run_constants.volume = run_info.volume;
run_constants.stoch_algor = run_info.stoch_algor;
run_constants.tau_error = run_info.tau_error;
run_constants.update = update;
run_constants.eventlog_arraysize = eventlog_arraysize;
run_constants.randseed = run_info.randseed;
run_constants.snapshot_interval = snapshot_interval;
run_constants.evlog_flag = evlog_flag;

[c_nums,rxnai,rxntime,times_executed,evetime,everxn, ssdata_changes, agevent] = ...
                           gateway_erg(modrxn,moddepend,rxgp.cmpd_nums(:,1),run_constants, mod_monitor);

rxgp.cmpd_nums(:,2) = c_nums;   % only 2nd column changes.  It holds the putative changes from cmpd_nums passed into gateway_erg.

for i = 1:numrxns
    rxgp.rxn{i}.ai = rxnai(i);
    rxgp.rxn{i}.time = rxntime(i);
end

num_events = sum(times_executed);

if num_events > eventlog_arraysize
    event_log = [evetime(1:eventlog_arraysize) everxn(1:eventlog_arraysize)];
else
    event_log = [evetime(1:num_events) everxn(1:num_events)];
end

if snapshot_interval > 0
%    size(agevent)
%    num_snapshots = floor(timestep/snapshot_interval);

    if update ~= 2
        num_snapshots = floor((start_time + timestep)/snapshot_interval) - floor(start_time/snapshot_interval);
        
        if num_snapshots >= 1
            for i = 1:num_snapshots
                ssdata(i,1) = ssdata_changes(1,i);
                aglog(i,1) = agevent(1,i);
                for j = 1:length(rxgp.monitor_comps)
                    ssdata(i,j+1) = rxgp.cmpd_nums(rxgp.monitor_comps(j),1) + ssdata_changes(j+1,i);
                end
               for j = 1:numrxns
                   aglog(i,j+1) = agevent(j+1,i);
              end

            end
        else
            ssdata(1,1) = start_time;
        for j = 1:length(rxgp.monitor_comps)
            ssdata(1,j+1) = rxgp.cmpd_nums(rxgp.monitor_comps(j),1);
        end

            for j = 1:numrxns
              aglog(1,j+1) = 0;
            end
        end
    else    % put initial values into ssdata
        ssdata(1,1) = start_time;
        for j = 1:length(rxgp.monitor_comps)
            ssdata(1,j+1) = rxgp.cmpd_nums(rxgp.monitor_comps(j),1);
        end

        for j = 1:numrxns
              aglog(1,j+1) = 0;
        end

    end
else
    ssdata = [];
    aglog = [];
end

clear evetime
clear everxn
clear ssdata_changes

%aglog

%event_log = [evetime(1:num_events) everxn(1:num_events)];
