controller_swarm2b(densities, flows, celldata, ts, idx)

function new_flows = controller_swarm2b(densities, flows, celldata, ts, idx)
% CONTROLLER_SWARM2B - implementation of SWARM 2B ramp metering strategy with
%                      modification that takes downstream cell into account.
%
% Call:   new_flows = controller_swarm2b(densities, flows, celldata, ts, idx)
%
% Parameters:
%             densities - vector of densities;
%             flows     - vector of on-ramp flows, same size as 'densities';
%             celldata  - array of freeway cell structures, whose length
%                         must be the same as size of 'densities';
%             ts        - sampling period;
%             idx       - index of the cell with our on-ramp.
%
% Returns:   new_flows - updated vector of on-ramp flows.
%
% Last modified:   11/06/2006.

%
% Alex Kurzhanskiy   <[e-mail]>
%

L = size(densities, 1);
M = size(flows, 1);
N = size(celldata, 2);

if L ~= N
  error('CONTROLLER_SWARM2B: number of densities does not match number of cells.');
end

if M ~= N
  error('CONTROLLER_SWARM2B: number of flows does not match number of cells.');
end

if (idx < 1) | (idx > N)
  error('CONTROLLER_SWARM2B: invalid cell index.');
end

new_flows = flows;  % initialize flow vector with current values

l         = abs(celldata(idx).PMend - celldata(idx).PMstart);  % length of i-th cell
nc        = celldata(idx).FDrhocrit * l;  % critical number of vehicles in i-th cell
n         = densities(idx, 1) * l;  % current number of vehicles in i-th cell
f1        = (nc - n) / ts;  % on-ramp flow that i-th cell allows

%if (idx < N) & (celldata(idx+1).FDrhocrit >= densities(idx+1, 1))
%  if celldata(idx).FRbeta < 1
%    if celldata(idx).FRbeta < 1e-7
%      F = celldata(idx+1).FDfmax;
%    else
%      b = (1 - celldata(idx).FRbeta) / celldata(idx).FRbeta;
%      F = min([celldata(idx+1).FDfmax (b*celldata(idx).FRfmax)]);
%    end
%    b = 1 / (1 - celldata(idx).FRbeta);
%    f2 = b * F  -  celldata(idx).FDfmax;
%  else
%    f2 = celldata(idx).FRfmax - celldata(idx).FDfmax;
%  end
%else
%    f2 = 0;
%end
%f                 = max([0 f1 f2 (f1+f2)]);

f                 = max([0 f1]);
f                 = max([f celldata(idx).ORmlcontroller.Cmin]);
f                 = min([f celldata(idx).ORmlcontroller.Cmax]);

new_flows(idx, 1) = f;  % update flow vector

return;