function [L,LX,LY] = fgs(tok,shot,t,varargin)
%FGS - FORWARD GRAD-SHAFRANOV SOLVER
% [L,LX,LY] = FGS(TOK,SHOT,T,'PAR',VAL,...)
% Solve the static free-boundary grad shafranov problem. Given conductor
% currents and internal plasma constraints, it returns the full
% equilibrium.
%
% L = FGS('TOK',SHOT,...) returns only the pre-computed quantities as provided by 
% FGSP<TOK>, FGSG<TOK> and FGSC for later use.
% [L,LX] = FGS('TOK',SHOT,...) additionally returns the data in LX
% generated from FGSX<TOK>
%
% For general help,see README.MD and MEQHELP
% For help on LY,  see MEQT and FBTT
% For help on LX,  see FGSX and FGSX<TOK>
% For help on L.P, see MEQP, MEQP<TOK>, FGSP and FGSP<TOK> 
% For help on L.G, see MEQG, MEQG<TOK>, FGSG and FGSG<TOK> 
%
% [+MEQ MatlabEQuilibrium Toolbox+] Swiss Plasma Center EPFL Lausanne 2022. All rights reserved.

assert(ischar(tok),'tok must be a string');
assert(isnumeric(shot),'shot must be numeric');
if nargin>=4, assert(~mod(numel(varargin),2),'varargin must come in ''parameter'',value pairs'); end

% Tokamak specific functions
tok = lower(tok); % tok: 'tcv' 'ids'
fgsptok = str2func(['fgsp' tok]);
fgsgtok = str2func(['fgsg' tok]);
fgsxtok = str2func(['fgsx' tok]);

% Algorithm configuration parameters, geometry and anciliary data
P = fgsptok(shot(1),varargin{:}); % parameters from sources
G = fgsgtok(shot(1),P);  % geometry from sources
P = fgsp(P); 
G = fgsg(G,P); 
L = fgsc(P,G); % Final assembly of all parameters

if nargout >= 2
  assert(isnumeric(t)||isempty(t),'time must be numeric');
  % Get inputs data for simulation
  LX = fgsxtok(t,L);
  LX = fgsx(L,LX); 
end

if nargout >= 3
  % Solve problem
  LY = fgst(L,LX);
end

end
