%% % saw_criteria module
function varargout = saw_criteria(varargin)
% #codegen

if nargin == 0,
    module_params = struct(...
        's1crit', 0.2,... %critical shear
        'mode', 0, ...%reconnection model: 0 full, 1 incomplete
        'w_norm', 0.); 
    varargout{1} = module_params;
    return %empty call, probably to get default structures
elseif nargin==7;
    x = varargin{1};
    g = varargin{2};
    v = varargin{3};
    u = varargin{4};
    it = varargin{5};
    model  = varargin{6};
    saw_params = varargin{7}; % distribute inputs
else
    error('must call with 0 or 7 inputs');
end


iota = eval_iota(x,g,v,model,true);
shear = eval_shear(x,g,v,model,true);

% bound of of r(q=1)
%%iqle1 = find(iota >= 1., 1, 'last'); % indices of q<=1
%%iqgt1 = find(iota < 1., 1, 'first'); % indices of q>1

% must use for loop: find gives variable size result that codegen does not like
%iq1 = find(iota >= 1., 1, 'last'); % outermost point with iota >=1
for iq1 = numel(iota):-1:1
if iota(iq1)>=1; break; end
end

% s1 = sum(shear([iq1,iq1+1])) / 2. ; % do better (e.g. linear) interpolation later
if iq1 == 1 %loop does not break
    s1 = 0.;
else
    m = (shear(iq1) - shear(iq1-1)) / (iota(iq1) - iota(iq1-1));
    q = shear(iq1-1) - m * iota(iq1-1);
    s1 = m * 1 + q;
end


if s1 >= saw_params.s1crit
    icrash = 4;
else 
    icrash = 0;
end

% s1 debug
% f = fopen('s1.txt', 'a');
% fprintf(f, '%3f \n', s1);
% fclose(f);

varargout{1} = icrash;
return
