charge+B; % set charge for absorptions

function varargout = beam(r,c,beamIntensity)

% Queries the black box or destroys a particle.
%
%  [ROWOUT, COLUMNOUT, SCORE] = BEAM(ROWIN,COLUMNIN) or
%  [ROWOUT, COLUMNOUT, SCORE] = BEAM(ROWIN,COLUMNIN,'low') shines a low
%  intensity light into the black box. If the beam is absorbed by any
%  particle ROWOUT and COLUMNOUT are both 0.
%
%  BEAM(ROWIN,COLUMNIN,'high') shines a high intensity light into the black
%  box, destroying any particle that tries to absorbe it.  
%

% Copyright 2006 The MathWorks, Inc.

persistent A turns charge newd ir ic n beamCtr dc saveinfo

if isempty(beamCtr)
    beamCtr = 0;
end

if nargin==1 % init puzzle and also used to get the number of function calls
    if isempty(r)
        A=[];
        [turns charge] = setCmds(A);
        if saveinfo
            varargout = {beamCtr,dc,0};
        else
            varargout = {beamCtr,0,0};
        end
        return
    else
        if any(r(:)<0)
            saveinfo = true;
            A = -r;
            dc = {};
        else
            saveinfo = false;
            A = r;
        end
        newd  = [0 0 0 0;    % end state     % 1 s, 2 e, 3 n, 4 w, 0 done
            1 2 3 4;    % go straight
            4 3 2 1;    % deflected
            2 1 4 3;    % deflected
            3 4 1 2];   % reflected
        ir = [1 0 -1 0];
        ic = [0 1 0 -1];
        n = length(A);
        [turns charge] = setCmds(A);
        varargout = {0,0,0};
        return
    end
end

if nargin==2
    beamIntensity = 'low';
end

if ~isscalar(r) || ~isscalar(c) || (~r==~c) || rem(r,1) || rem(c,1)
    error('rowIn or columnIn has invalid values.')
end

rrr = r;
ccc = c;

d = (c>0)+2*(r>0)+3*(c<0)+4*(r<0); %1 s, 2 e, 3 n, 4 w
r = abs(r);
c = abs(c);

if  r>n || c>n 
    error('abs(rowIn) or abs(columnIn) has a value larger than n.')
end

r = r+1 + ((n+1)*(d==3));
c = c+1 + ((n+1)*(d==4));
ch = charge(r,c);
rr = r;
cc = c;
if ~ch
    while d
        r = r + ir(d);
        c = c + ic(d);
        rr = [rr r];
        cc = [cc c];
        d = newd(turns(r,c),d);
        ch = ch + charge(r,c);
    end
end
d = (r==n+2)||(c==n+2);
r = rem(r-1,n+1);
c = rem(c-1,n+1);
des = [0;0];
switch beamIntensity
    case 'low'
        beamCtr = beamCtr+1;
        if r&&c
            ch = 0;
            varargout = {0,0,ch};    % beam is absorbed
        elseif d
            varargout = {-r,-c,ch}; % beam leaves at bottom or right
        else
            varargout = {r,c,ch};   % beam leaves at top or left
        end
    case 'high'
        ch = 0;
        if r&&c         % beam is absorbed and particle is destroyed
            beamCtr = beamCtr+A(r,c);   
            A(r,c) = 0;                  % destroy particle
            [turns charge] = setCmds(A); 
            des = [r;c];
        end
        varargout = {0,0,ch};
    otherwise
        error('beamIntensity can only be ''low'' or ''high''.')
end

if saveinfo
    if d
        dc = [dc;{rrr,ccc,rr,cc,ch,isequal(beamIntensity,'high'),des(1),des(2),-r,-c}];
    else
        dc = [dc;{rrr,ccc,rr,cc,ch,isequal(beamIntensity,'high'),des(1),des(2),r,c}];
    end
end

function [turns charge] = setCmds(A)
B = conv2(A,[0 0 0;0 1 0;0 0 0]);
C = conv2(A,[0 1 0;1 0 1;0 1 0])&~B;
turns = conv2(single(A>0),[1 0 2;0 0 0;2 0 1])+2;       % set turns
turns(~conv2(ones(size(A)),[0 0 0;0 1 0;0 0 0])|B) = 1; % set absorption and borders
turns(C) = 2;                                           % remove invalid turns
charge = conv2(A,ones(3)); % set charge for turns
charge(C|B) = 0;           % remove charge at invalid turns
charge = charge+B;         % set charge for absorptions