From 0c3f4fcb41880df4fae3b32c1dba27d242b66ea3 Mon Sep 17 00:00:00 2001 From: Olivier Sauter <olivier.sauter@epfl.ch> Date: Tue, 4 Jul 2017 13:07:06 +0000 Subject: [PATCH] rm geteqdskJET.m and add to CHEASEgui subfolder for consistency git-svn-id: https://spcsvn.epfl.ch/repos/TCV/gdat/trunk@7752 d63d8f72-b253-0410-a779-e742ad2e26cf --- crpptbx/JET/geteqdskJET.m | 412 -------------------------------------- 1 file changed, 412 deletions(-) delete mode 100644 crpptbx/JET/geteqdskJET.m diff --git a/crpptbx/JET/geteqdskJET.m b/crpptbx/JET/geteqdskJET.m deleted file mode 100644 index f4ed9d20..00000000 --- a/crpptbx/JET/geteqdskJET.m +++ /dev/null @@ -1,412 +0,0 @@ -function [efitdata,eqd]=geteqdskJET(shot,times_in,nrg,nzg,savedir,deltaz,efitlab,uid,seqd,varargin); -% -% function [efitdata]=geteqdskJET(shot,times_in,nrg,nzg,savedir,deltaz,efitlab,uid,seqd,varargin); -% -% save eqdsk file to savedir/JET_EQDSK_shot_t(times_in(:)) -% returns data used in structure efitdata (efitdata.shot, efitdata.Rbnd, efitdataZbnd....) -% -% examples: -% efitdata=geteqdskJET(shot,times_in,33,65,'/tmp') -% efitdata=geteqdskJET(shot,times_in,33,65,'/tmp',[],[],[],[],1,efitdata); % gives data so no need to load again -% -% INPUTS : -% shot : shot number -% times_in : single time slice of analysis -% -% OPTIONAL -% nrg, nzg: nb of R and Z points for R,Z grid: default: 129, -129 -% if nzg<0 force symmetric box aorund z=deltaz -% savedir: directory to save eqdsk file (default: '/tmp/username') -% deltaz: shift equilibrium vertically -% = 0: (default) no shift -% = -99: shift so that zmag=0 (thus deltaz=-zmag) -% otherwise: shift Z position by this amount -% efitlab : 'efit' or 'eftm' (default: 'efit') -% uid{i}: user id if different from main database (default: 'jetppf') -% seq(i): sequence number if not last one required (default: 0) -% i: 1 (efit), 2 (chain2: lid2), 3 (equi) -% -% varargin{1}: plot option: 0: do not plot contours, >0 plot contour with varargin{1} nb of contours (60 is good) -% varargin{2}: efitdata in input (no need to load again if shot=efitdata.shot) -% -% OUTPUTS -% efitdata: structure containing all data vs time so that one can call geteqdsk for another time quicker for same shot -% -%~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~ - -if ~exist('nrg') || isempty(nrg) - nrg = 129; -end -if ~exist('nzg') || isempty(nzg) - nzg = -129; -end -if nargin<=4 | isempty(savedir) - savedir=['/tmp/' getenv('USER')]; - if ~exist(savedir,'dir') - unix(['mkdir ' savedir]); - end -end -if unix(['test -d ' savedir]) - disp(['Problems in geteqdskJET, savedir=' savedir ' is not a directory']) - return -end -if nargin<=5 | isempty(deltaz) - deltaz=0; -end -if nargin<=6 | isempty(efitlab) - efitlab='efit'; -end -iexefit=1; -iexchain2=2; -iexequi=3; -uid0{iexefit}='jetppf'; -uid0{iexchain2}='jetppf'; -uid0{iexequi}='jetppf'; -seqd0(1:iexequi)=0; -if nargin>=8 & ~isempty(uid) - for i=length(uid)+1:length(uid0) - uid{i}=uid0{i} - end -else - uid=uid0; -end -if nargin>=9 & ~isempty(seqd) - for i=length(seqd)+1:length(seqd0) - seqd(i)=seqd0(i); - end -else - seqd=seqd0; -end -for i=1:length(uid) - s_extra{i}=['?uid=' uid{i} '+seq=' num2str(seqd(i))]; -end - -ncont=0; -if nargin>=10 & ~isempty(varargin{1}) - ncont=varargin{1}; -end - -iread=1; -if nargin>=11 & ~isempty(varargin{2}) - efitdata=varargin{2}; - if shot==efitdata.shot - iread=0; - else - disp(['shot=' num2str(shot) ' different from efitdata.shot=' num2str(efitdata.shot)]) - warning('reload data with new shot') - end -end - -do_plot=0; -exp_machine='jet'; -if iread==1 - efitdata.shot=shot; - % get data needed - efitdata.Rbnd=gdat(shot,[{'ppf'}, {efitlab}, {['RBND' s_extra{iexefit}]}],'machine',exp_machine,'doplot',do_plot); - efitdata.tefit=efitdata.Rbnd.t; - efitdata.Zbnd=gdat(shot,[{'ppf'}, {efitlab}, {['ZBND' s_extra{iexefit}]}],'machine',exp_machine,'doplot',do_plot); - efitdata.R0=gdat(shot,[{'ppf'}, {efitlab}, {['RGEO' s_extra{iexefit}]}],'machine',exp_machine,'doplot',do_plot); - efitdata.a=gdat(shot,[{'ppf'}, {efitlab}, {['CR0' s_extra{iexefit}]}],'machine',exp_machine,'doplot',do_plot); - efitdata.Z0=gdat(shot,[{'ppf'}, {efitlab}, {['ZGEO' s_extra{iexefit}]}],'machine',exp_machine,'doplot',do_plot); - efitdata.rmag=gdat(shot,[{'ppf'}, {efitlab}, {['rmag' s_extra{iexefit}]}],'machine',exp_machine,'doplot',do_plot); - efitdata.zmag=gdat(shot,[{'ppf'}, {efitlab}, {['zmag' s_extra{iexefit}]}],'machine',exp_machine,'doplot',do_plot); - efitdata.faxs=gdat(shot,[{'ppf'}, {efitlab}, {['faxs' s_extra{iexefit}]}],'machine',exp_machine,'doplot',do_plot); - efitdata.fbnd=gdat(shot,[{'ppf'}, {efitlab}, {['fbnd' s_extra{iexefit}]}],'machine',exp_machine,'doplot',do_plot); - efitdata.bvac=gdat(shot,[{'ppf'}, {efitlab}, {['bvac' s_extra{iexefit}]}],'machine',exp_machine,'doplot',do_plot); - efitdata.ip=gdat(shot,[{'ppf'}, {efitlab}, {['xip' s_extra{iexefit}]}],'machine',exp_machine,'doplot',do_plot); - efitdata.F=gdat(shot,[{'ppf'}, {efitlab}, {['F' s_extra{iexefit}]}],'machine',exp_machine,'doplot',do_plot); - efitdata.psin=efitdata.F.x; - efitdata.P=gdat(shot,[{'ppf'}, {efitlab}, {['P' s_extra{iexefit}]}],'machine',exp_machine,'doplot',do_plot); - efitdata.Q=gdat(shot,[{'ppf'}, {efitlab}, {['Q' s_extra{iexefit}]}],'machine',exp_machine,'doplot',do_plot); - efitdata.kappa=gdat(shot,[{'ppf'}, {efitlab}, {['ELON' s_extra{iexefit}]}],'machine',exp_machine,'doplot',do_plot); - efitdata.q95=gdat(shot,[{'ppf'}, {efitlab}, {['q95' s_extra{iexefit}]}],'machine',exp_machine,'doplot',do_plot); - efitdata.btpd=gdat(shot,[{'ppf'}, {efitlab}, {['btpd' s_extra{iexefit}]}],'machine',exp_machine,'doplot',do_plot); - efitdata.bttd=gdat(shot,[{'ppf'}, {efitlab}, {['bttd' s_extra{iexefit}]}],'machine',exp_machine,'doplot',do_plot); - efitdata.btnd=gdat(shot,[{'ppf'}, {efitlab}, {['btnd' s_extra{iexefit}]}],'machine',exp_machine,'doplot',do_plot); - efitdata.btpm=gdat(shot,[{'ppf'}, {efitlab}, {['btpm' s_extra{iexefit}]}],'machine',exp_machine,'doplot',do_plot); - efitdata.bttm=gdat(shot,[{'ppf'}, {efitlab}, {['bttm' s_extra{iexefit}]}],'machine',exp_machine,'doplot',do_plot); - efitdata.btnm=gdat(shot,[{'ppf'}, {efitlab}, {['btnm' s_extra{iexefit}]}],'machine',exp_machine,'doplot',do_plot); - efitdata.xli=gdat(shot,[{'ppf'}, {efitlab}, {['li3m' s_extra{iexefit}]}],'machine',exp_machine,'doplot',do_plot); - efitdata.sspr=gdat(shot,[{'ppf'}, {efitlab}, {['sspr' s_extra{iexefit}]}],'machine',exp_machine,'doplot',do_plot); - efitdata.sspi=gdat(shot,[{'ppf'}, {efitlab}, {['sspi' s_extra{iexefit}]}],'machine',exp_machine,'doplot',do_plot); -% $$$ % add for profiles -% $$$ efitdata.ti=gdat(shot,[{'ppf'},{'TION'},{['TI' s_extra{iexchain2}]}],'machine',exp_machine,'doplot',do_plot); -% $$$ efitdata.p_tion=gdat(shot,[{'ppf'},{'TION'},{['p' s_extra{iexchain2}]}],'machine',exp_machine,'doplot',do_plot); -% $$$ efitdata.pi=gdat(shot,[{'ppf'},{'NION'},{['DD' s_extra{iexchain2}]}],'machine',exp_machine,'doplot',do_plot); -% $$$ efitdata.zef=gdat(shot,[{'ppf'},{'NION'},{['ZEF' s_extra{iexchain2}]}],'machine',exp_machine,'doplot',do_plot); -% $$$ % add for calculating NTM parameters -% $$$ efitdata.bpol=gdat(shot,[{'ppf'},{'equi'},{['bpol' s_extra{iexequi}]}],'machine',exp_machine,'doplot',do_plot); -% $$$ efitdata.bpo2=gdat(shot,[{'ppf'},{'equi'},{['bpo2' s_extra{iexequi}]}],'machine',exp_machine,'doplot',do_plot); -% $$$ efitdata.qmag=gdat(shot,[{'ppf'},{'' efitlab ''},{['qmag' s_extra{iexefit}]}],'machine',exp_machine,'doplot',do_plot); -% $$$ efitdata.lidrpe=gdat(shot,[{'ppf'},{'lidr'},{['pe']}],'machine',exp_machine,'doplot',do_plot); -% $$$ efitdata.nexav=gdat(shot,[{'ppf'},{'nex'},{['av']}],'machine',exp_machine,'doplot',do_plot); -% $$$ efitdata.nbi=gdat(shot,[{'ppf'},{'nbi'},{['ptot']}],'machine',exp_machine,'doplot',do_plot); -% $$$ efitdata.icrh=gdat(shot,[{'ppf'},{'icrh'},{['ptot']}],'machine',exp_machine,'doplot',do_plot); -% $$$ % old, not working anymore as of June 2017? -% $$$ % $$$ efitdata.ptot=gdat(shot,[{'ppf'},{'mg3'},{['yto']}],'machine',exp_machine,'doplot',do_plot); -% $$$ % $$$ efitdata.halpha=gdat(shot,[{'jpf'},{'dd'},{['s3-ad35']}],'machine',exp_machine,'doplot',do_plot); -% $$$ % $$$ efitdata.n1=gdat(shot,[{'jpf'},{'da'},{['c1-g101']}],'machine',exp_machine,'doplot',do_plot); -% $$$ % $$$ efitdata.n2=gdat(shot,[{'jpf'},{'da'},{['c1-g102']}],'machine',exp_machine,'doplot',do_plot); -end -tefit=efitdata.tefit; -index_efit = iround_os(tefit,times_in); -for it=1:numel(index_efit) - index_efit_eff=index_efit(it); - time_efit=tefit(index_efit_eff); - - Rbnd=efitdata.Rbnd.data(:,index_efit_eff); - Zbnd=efitdata.Zbnd.data(:,index_efit_eff); - R0=efitdata.R0.data(index_efit_eff); - Z0=efitdata.Z0.data(index_efit_eff); - rmag=efitdata.rmag.data(index_efit_eff); - zmag=efitdata.zmag.data(index_efit_eff); - if deltaz==-99 - deltaz=-zmag; - end - faxs=efitdata.faxs.data(index_efit_eff); - fbnd=efitdata.fbnd.data(index_efit_eff); - bvac=efitdata.bvac.data(index_efit_eff); - ip=efitdata.ip.data(index_efit_eff); - psi_efit=efitdata.psin; - F=efitdata.F.data(:,index_efit_eff); - P=efitdata.P.data(:,index_efit_eff); - Q=efitdata.Q.data(:,index_efit_eff); - kappa=efitdata.kappa.data(index_efit_eff); - q95=efitdata.q95.data(index_efit_eff); - if index_efit_eff<=length(efitdata.btpd.data); btpd=efitdata.btpd.data(index_efit_eff); else; btpd=0; end - if index_efit_eff<=length(efitdata.bttd.data); bttd=efitdata.bttd.data(index_efit_eff); else; bttd=0; end - if index_efit_eff<=length(efitdata.btnd.data); btnd=efitdata.btnd.data(index_efit_eff); else; btnd=0; end - if index_efit_eff<=length(efitdata.btpm.data); btpm=efitdata.btpm.data(index_efit_eff); else; btpm=0; end - if index_efit_eff<=length(efitdata.bttm.data); bttm=efitdata.bttm.data(index_efit_eff); else; bttm=0; end - if index_efit_eff<=length(efitdata.btnm.data); btnm=efitdata.btnm.data(index_efit_eff); else; btnm=0; end - xli=efitdata.xli.data(index_efit_eff); - [R,Z,psinrz]=psinrzjet(shot,time_efit,nrg,nzg,efitlab,uid{iexefit},seqd(iexefit),ncont,efitdata.sspr,efitdata.sspi,deltaz); - nzg=abs(nzg); - - % define file name - s = sprintf('%.6f',time_efit); - fname=['EQDSK.' num2str(shot) 't' s]; - fnamefull=fullfile(savedir,fname); - - % fill in eqdsk structures - eqd{it}.fnamefull = fnamefull; - eqd{it}.fname = fname; - eqd{it}.cocos = 7; - % 1st eqdsk line: 48 characters for file description and then, 3, nr, nz - tdate=date; - ss=['JET #' num2str(shot) ' t= ' num2str(time_efit) ' ' efitlab '?uid=' uid{iexefit} '+seq=' num2str(seqd(iexefit)) ', ' tdate]; - if length(ss)<48 - ss(end:48)=' '; - else - ss=ss(1:48); - end - eqd{it}.stitle = ss; - eqd{it}.nr = nrg; - eqd{it}.nz = nzg; - rmin=min(R); - rmax=max(R); - zmin=min(Z); - zmax=max(Z); - eqd{it}.rboxlen = rmax-rmin; - eqd{it}.zboxlen = zmax-zmin; - eqd{it}.r0 = R0; - eqd{it}.rboxleft = rmin; - eqd{it}.zmid = 0.5*(zmin+zmax)+deltaz; - eqd{it}.rmesh=linspace(eqd{it}.rboxleft,eqd{it}.rboxleft+eqd{it}.rboxlen,eqd{it}.nr)'; - eqd{it}.zmesh=linspace(eqd{it}.zmid-eqd{it}.zboxlen/2,eqd{it}.zmid+eqd{it}.zboxlen/2,eqd{it}.nz)'; - eqd{it}.raxis = rmag; - eqd{it}.zaxis = zmag+deltaz; - eqd{it}.psiaxis = faxs; - eqd{it}.psiedge = fbnd; - eqd{it}.b0 = bvac; - eqd{it}.ip = ip; - rdel=+0.8*(R(2)-R(1)); - zdel=+0.8*(Z(2)-Z(1)); - eqd{it}.rlim=[eqd{it}.rboxleft-rdel ; eqd{it}.rboxleft+eqd{it}.rboxlen+rdel ; ... - eqd{it}.rboxleft+eqd{it}.rboxlen+rdel ; eqd{it}.rboxleft-rdel ; ... - eqd{it}.rboxleft-rdel]; - eqd{it}.zlim=[eqd{it}.zmid-eqd{it}.zboxlen/2-zdel ; eqd{it}.zmid-eqd{it}.zboxlen/2-zdel ; ... - eqd{it}.zmid+eqd{it}.zboxlen/2+zdel ; eqd{it}.zmid+eqd{it}.zboxlen/2+zdel ; ... - eqd{it}.zmid-eqd{it}.zboxlen/2-zdel]; - eqd{it}.nblim = 5; - psieq=[0:1/(nrg-1):1]; - eqd{it}.psimesh = psieq; - eqd{it}.rhopsi = sqrt(eqd{it}.psimesh); - psi_efit_eff=faxs+psi_efit.*(fbnd-faxs); - psieq_eff=faxs+psieq.*(fbnd-faxs); - [G Gprime]=interpos(13,psi_efit_eff,F,psieq_eff); - eqd{it}.F = G; - eqd{it}.FFprime = G.*Gprime; - [press pressprime]=interpos(13,psi_efit_eff,P,psieq_eff); - eqd{it}.p = press; - eqd{it}.pprime = pressprime; - psirz=faxs+psinrz.*(fbnd-faxs); - eqd{it}.psi = psirz; - eqd{it}.psirz = reshape(eqd{it}.psi,prod(size(eqd{it}.psi)),1); - y=interpos(13,psi_efit,Q,psieq,-0.01); - eqd{it}.q = y; - npts=length(Rbnd); - eqd{it}.nbbound = npts; - eqd{it}.rplas = Rbnd; - eqd{it}.zplas = Zbnd+deltaz; - % Some useful data to compare with recomputed equilibria - eqd{it}.extralines = []; - eqd{it}.extralines{end+1} = [' ' num2str(fbnd-faxs) ' psiedge-psiax ']; - eqd{it}.extralines{end+1} = [' ' num2str(rmag) ' r-magaxe ']; - eqd{it}.extralines{end+1} = [' ' num2str(zmag+deltaz) ' z-magaxe ']; - eqd{it}.extralines{end+1} = [' ' num2str(0.5.*(min(Zbnd)+max(Zbnd))+deltaz) ' z0 (zaver) ']; - eqd{it}.extralines{end+1} = [' ' num2str(R0) ' r-major ']; - eqd{it}.extralines{end+1} = [' ' num2str(bvac) ' B0 ']; - eqd{it}.extralines{end+1} = [' ' num2str(ip/4e-7/pi) ' CURRT -> I-p [A]: ' num2str(ip/4e-7/pi,ip)]; - eqd{it}.extralines{end+1} = [' ' num2str(kappa) ' kappa ']; - eqd{it}.extralines{end+1} = [' ' num2str(Q(1)) ' q_0 ']; - y=interpos(13,psi_efit,Q,[0:0.01:1],1e-6); - eqd{it}.extralines{end+1} = [' ' num2str(Q(end)) ' q_edge, from int: ' num2str(y(end))]; - eqd{it}.extralines{end+1} = [' ' num2str(q95) ' q_95, from int: ' num2str(y(96))]; - [qmin ind]=min(y); - eqd{it}.extralines{end+1} = [' ' num2str(qmin) ' q_min, at psinorm= ' num2str(y(ind))]; - eqd{it}.extralines{end+1} = [' ' num2str(btpd) ' beta_pol(wdia) ']; - eqd{it}.extralines{end+1} = [' ' num2str(btpm) ' beta_pol(efit) ']; - eqd{it}.extralines{end+1} = [' ' num2str(bttd) ' beta_tor(wdia) ']; - eqd{it}.extralines{end+1} = [' ' num2str(bttm) ' beta_tor(efit) ']; - eqd{it}.extralines{end+1} = [' ' num2str(btnd) ' beta_N(wdia) ']; - eqd{it}.extralines{end+1} = [' ' num2str(btnm) ' beta_N(efit) ']; - eqd{it}.extralines{end+1} = [' ' num2str(xli) ' li ']; - - eqd{it}.extralines{end+1} = [' ' num2str(time_efit) 'time ']; - eqd{it}.extralines{end+1} = [' ' num2str(shot) 'shot number ']; - eqd{it}.extralines{end+1} = ['efit version : %s ' num2str(efitlab)]; - eqd{it}.extralines{end+1} = ['uid : ' num2str(uid{iexefit})]; - eqd{it}.extralines{end+1} = ['seq : ' num2str(seqd(iexefit))]; - % write to file (it write the positive Ip, B0 as well automatically - write_eqdsk(eqd{it}.fnamefull,eqd{it}); % cocos contained in eqdsk - -% $$$ fid=fopen(fnamefull,'w'); -% $$$ -% $$$ % 1st eqdsk line: 48 characters for file description and then, 3, nr, nz -% $$$ tdate=date; -% $$$ ss=['JET #' num2str(shot) ' t= ' num2str(time_efit) ' ' efitlab '?uid=' uid{iexefit} '+seq=' num2str(seqd(iexefit)) ', ' tdate]; -% $$$ if length(ss)<48 -% $$$ ss(end:48)=' '; -% $$$ else -% $$$ ss=ss(1:48); -% $$$ end -% $$$ -% $$$ fprintf(fid,'%s%4d%4d%4d\n',ss,3,nrg,nzg); -% $$$ -% $$$ % 2nd line: rboxlen, zboxlen, r0, rboxlft, zboxmid -% $$$ fprintf(fid,'%16.9e%16.9e%16.9e%16.9e%16.9e\n',rmax-rmin,zmax-zmin,R0,rmin,0.5*(zmin+zmax)+deltaz); -% $$$ -% $$$ % 3rd line: rmag, zmag, psimag, psiedge, B0 -% $$$ fprintf(fid,'%16.9e%16.9e%16.9e%16.9e%16.9e\n',rmag,zmag+deltaz,faxs,fbnd,bvac); -% $$$ -% $$$ % 4th line: Ip, psiax1, psiax2, raxis1, raxis2 -% $$$ fprintf(fid,'%16.9e%16.9e%16.9e%16.9e%16.9e\n',ip,faxs,0.,rmag,0.); -% $$$ -% $$$ % 5th line: zaxis1, zaxis2, psi_sep, R_xpoint, Z_xpoint -% $$$ fprintf(fid,'%16.9e%16.9e%16.9e%16.9e%16.9e\n',zmag+deltaz,0.0,fbnd,-1,-1); -% $$$ -% $$$ % 6th entry: F(psi) on nr equidistant psi mesh -% $$$ psieq=[0:1/(nrg-1):1]; -% $$$ psi_efit_eff=faxs+psi_efit.*(fbnd-faxs); -% $$$ psieq_eff=faxs+psieq.*(fbnd-faxs); -% $$$ [G Gprime]=interpos(13,psi_efit_eff,F,psieq_eff); -% $$$ fprintf(fid,'%16.9e%16.9e%16.9e%16.9e%16.9e\n',G(1:end-1)); -% $$$ fprintf(fid,'%16.9e\n',G(end)); -% $$$ -% $$$ % 7th entry: p(psi) on nr equidistant psi mesh -% $$$ [press pressprime]=interpos(13,psi_efit_eff,P,psieq_eff); -% $$$ fprintf(fid,'%16.9e%16.9e%16.9e%16.9e%16.9e\n',press(1:end-1)); -% $$$ fprintf(fid,'%16.9e\n',press(end)); -% $$$ -% $$$ % 8th entry: FF'(psi) on nr equidistant psi mesh -% $$$ y=G.*Gprime; -% $$$ fprintf(fid,'%16.9e%16.9e%16.9e%16.9e%16.9e\n',y(1:end-1)); -% $$$ fprintf(fid,'%16.9e\n',y(end)); -% $$$ -% $$$ % 9th entry: p'(psi) on nr equidistant psi mesh (in MKSA) -% $$$ fprintf(fid,'%16.9e%16.9e%16.9e%16.9e%16.9e\n',pressprime(1:end-1)); -% $$$ fprintf(fid,'%16.9e\n',pressprime(end)); -% $$$ -% $$$ % 10th entry: psi(i,j) -% $$$ psirz=faxs+psinrz.*(fbnd-faxs); -% $$$ fprintf(fid,'%16.9e%16.9e%16.9e%16.9e%16.9e\n',psirz); -% $$$ if mod(nrg*nzg,5)~=0 -% $$$ fprintf(fid,'\n'); -% $$$ end -% $$$ -% $$$ % 11th entry: q profile on nr equidistant psi mesh -% $$$ y=interpos(13,psi_efit,Q,psieq,1e-6); -% $$$ fprintf(fid,'%16.9e%16.9e%16.9e%16.9e%16.9e\n',y(1:end-1)); -% $$$ fprintf(fid,'%16.9e\n',y(end)); -% $$$ -% $$$ % 12th entry: (R,Z) plasma boundary and wall position -% $$$ npts=length(Rbnd); -% $$$ fprintf(fid,'%5d%5d\n',npts,5); -% $$$ fprintf(fid,'%16.9e%16.9e%16.9e%16.9e%16.9e\n',[Rbnd Zbnd+deltaz]'); -% $$$ if mod(2*npts,5) ~= 0 -% $$$ fprintf(fid,'\n'); -% $$$ end -% $$$ rdel=+0.8*(R(2)-R(1)); -% $$$ zdel=+0.8*(Z(2)-Z(1)); -% $$$ fprintf(fid,'%16.9e%16.9e%16.9e%16.9e%16.9e\n',min(Rbnd)-rdel,min(Zbnd)+deltaz-zdel, ... -% $$$ max(Rbnd)+rdel,min(Zbnd)+deltaz-zdel,max(Rbnd)+rdel,max(Zbnd)+deltaz+zdel, ... -% $$$ min(Rbnd)-rdel,max(Zbnd)+deltaz+zdel,min(Rbnd)-rdel,min(Zbnd)+deltaz-zdel); -% $$$ -% $$$ % Some useful data to compare with recomputed equilibria -% $$$ fprintf(fid,'%18.8e psiedge-psiax\n',fbnd-faxs); -% $$$ fprintf(fid,'%18.8e r-magaxe\n',rmag); -% $$$ fprintf(fid,'%18.8e z-magaxe\n',zmag+deltaz); -% $$$ fprintf(fid,'%18.8e z0 (zaver)\n',0.5.*(min(Zbnd)+max(Zbnd))+deltaz); -% $$$ fprintf(fid,'%18.8e r-major\n',R0); -% $$$ fprintf(fid,'%18.8e B0\n',bvac); -% $$$ fprintf(fid,'%18.8e CURRT -> I-p [A]: %18.8e \n',ip/4e-7/pi,ip); -% $$$ fprintf(fid,'%18.8e kappa\n',kappa); -% $$$ fprintf(fid,'%18.8e q_0\n',Q(1)); -% $$$ y=interpos(13,psi_efit,Q,[0:0.01:1],1e-6); -% $$$ fprintf(fid,'%18.8e q_edge, from int: %18.8e\n',Q(end),y(end)); -% $$$ fprintf(fid,'%18.8e q_95, from int: %18.8e\n',q95,y(96)); -% $$$ [qmin ind]=min(y); -% $$$ fprintf(fid,'%18.8e q_min, at psinorm= %18.8e\n',qmin,y(ind)); -% $$$ fprintf(fid,'%18.8e beta_pol(wdia)\n',btpd); -% $$$ fprintf(fid,'%18.8e beta_pol(efit)\n',btpm); -% $$$ fprintf(fid,'%18.8e beta_tor(wdia)\n',bttd); -% $$$ fprintf(fid,'%18.8e beta_tor(efit)\n',bttm); -% $$$ fprintf(fid,'%18.8e beta_N(wdia)\n',btnd); -% $$$ fprintf(fid,'%18.8e beta_N(efit)\n',btnm); -% $$$ fprintf(fid,'%18.8e li\n',xli); -% $$$ -% $$$ fprintf(fid,'\n%18.8e time\n',time_efit); -% $$$ fprintf(fid,' %d shot number\n',shot); -% $$$ fprintf(fid,' efit version : %s\n',efitlab); -% $$$ fprintf(fid,' uid : %s, %s, %s\n',uid{1},uid{2},uid{3}); -% $$$ fprintf(fid,' seq : %d, %d, %d\n',seqd(1:3)); -% $$$ -% $$$ fclose(fid); -% $$$ disp(['wrote ',fnamefull]); - - if ncont>0 - figure; - pos=get(gcf,'position'); - set(gcf,'position',[pos(1)+0.5*pos(3) 0.8*abs(pos(2)-pos(4)) pos(3) 2*pos(4)]) - subplot(3,1,1) - plot(psi_efit,P/P(1),'-') - hold on - plot(psi_efit,F/F(end),'r-') - ss=sprintf('%.4f',time_efit); - title(['JET #' num2str(shot) ' t= ' ss]) - legend('P/P(0)','F/F(edge)',3) - subplot(3,1,2) - plot(psieq,pressprime/abs(pressprime(1)),'-') - hold on - plot(psieq,G.*Gprime/abs(G(1).*Gprime(1)),'r-') - legend('Pprime/|Pprime(0)|','F*Fprime/|F*Fprime(0)|',2) - subplot(3,1,3) - plot(psi_efit,Q,'-') - hold on - aa=axis; - plot([0.95 0.95],[aa(3) aa(4)],'k--') - grid on - xlabel('\psi/\psi_{edge}') - end -end -- GitLab