%FL4PMEX  Four point interpolation for the flux on the limiter
%
% This is the MATLAB equivalent implementation of libmeq/fl4p.c.
% A more detailed help is available in FL4PMEX.
%
% [+MEQ MatlabEQuilibrium Toolbox+] Swiss Plasma Center EPFL Lausanne 2022. All rights reserved.
function  Fl = fl4pmexm(Fx,kl,cl,FN,rl,zl,rA,zA)
nz = size(Fx,1);
Fl = Fx(kl+1   ).*cl(1,:)' + ...
     Fx(kl+nz+1).*cl(2,:)' + ...
     Fx(kl+nz+2).*cl(3,:)' + ...
     Fx(kl+2   ).*cl(4,:)';
if nargin >= 4 && FN
  dFl = diff(Fl([end 1:end]));
  Fl(dFl .* dFl([2:end 1]) > 0 | (dFl - dFl([2:end 1])) * FN > 0) = FN;
  
  dRFl = (Fx(kl+nz+1) - Fx(kl   +1)).*(cl(2,:)+cl(1,:))' + ...
         (Fx(kl+nz+2) - Fx(kl   +2)).*(cl(4,:)+cl(3,:))';
  dZFl = (Fx(kl   +2) - Fx(kl   +1)).*(cl(4,:)+cl(1,:))' + ...
         (Fx(kl+nz+2) - Fx(kl+nz+1)).*(cl(3,:)+cl(2,:))';
  dAFl = (dRFl.*(rA-rl) + dZFl.*(zA-zl));
  Fl(dAFl*FN > 0) = FN;
end

end