%% Script to perform radial scan using FGS
%
% [+MEQ MatlabEQuilibrium Toolbox+] Swiss Plasma Center EPFL Lausanne 2022. All rights reserved.

clear;

% [L,LX0] = fgs('tcv',65668,0.3,'fancyplot',true,'insrc','liu','debug',1,'iterq',10);
[L,LX0] = fgs('tcv',65279,0.6,'fancyplot',true,'insrc','liu','debug',1,'iterq',10);

%%
factgrid = linspace(0.95,1.08,8);

nr=2; nc = ceil(numel(factgrid)/nr);
figure(1); clf;
LY = LX0;
for ii=1:numel(factgrid)
  LX = LY; % previous as init
  LX.Ip = factgrid(ii)*LX0.Ip; % scale Ip
  LX.qA = LX0.qA; LX.bp = LX0.bp;
  LY = fgst(L,LX);
  
  subplot(nr,nc,ii);
  meqplott(L,LY);
  drawnow;
  
  R = LY.rgeom(end); a = LY.aminor(end); %(R,a of LCFS)
  titstr = sprintf('bp=%2.2f\nli=%2.2f\nR/a=%2.2f',LY.bp,LY.li,R/a);
  
  title(titstr,'fontsize',8)
  bpgrid(ii) = LY.bp;
  ragrid(ii) = R/a;
  rgrid(ii) = R;
end

%%
figure(2)
subplot(211)
plot(rgrid,ragrid)
xlabel('R'),ylabel('R/a');
subplot(212)
plot(rgrid,log(ragrid))
xlabel('R'),ylabel('ln(R/a)');
