function [trace,error,varargout]=loadAUGdata(shot,data_type,varargin)
%
% data_type:
% 'Ip' = current
% 'zmag' = vertical position of the center of the plasma (magnetic axis)
% 'rmag' = radial position of the center of the plasma
% 'sxr' = soft x-ray emission
% 'sxR' = soft x-ray emission with varargout{1} option (requires varargin{5}!)
% 'SXB' = soft x-ray emission from (by default camera J) SXB/J_xx camera (sxb, sxB, etc all work)
% 'SXF' = soft x-ray emission from (by default camera J) SXF/I_xx camera (sxf, sxF, etc all work)
% 'SSX_H' = from SSX, H channel, same for G, I, etc
% 'SSX' : gets the SSX_G by default at this stage
%
% gdat(15133,'MAG/Ipa',1,'AUG') % uses experiment=AUGD shotfiles per default
% gdat(15133,'MAG/Ipa',1) (sufficient at AUG since 'AUG' per defaut, same with gdat(15133,'ip',1)
% gdat(15133,'AUGD/MAG/Ipa',1,'AUG') % to specify experiment explicitely like in:
% gdat(30230,'ECED/RMD/Trad-A',1,'AUG') %
% gdat(30230,'ECED/CEC/Trad-A',1,'AUG') %
% gdat(31053,'MSP/IGP_07Co',1); % N gas opening at -6s (03 as well)
%
% INPUT:
% shot: shot number
% data_type: type of the required data: 'diag_name/sig_name'
%
% examples:
% data_type='SXR/B', 'TOT/beta_N', 'SXB/J_053', 'SXB'
% data_type='POT/ELMa-Han', 'MOD/OddNAmp', 'MOD/EvenNAmp', 'TOT/PNBI_TOT', 'TOT/P_TOT'
%
% Meaning of varargin depends on data_type:
%
% data_type=sxr, sxb, sxf or ece, eced:
% varargin{1}: [i1 i2] : if not empty, assumes need many chords from i1 to i2
% varargin{2}: channel status: 1=unread yet, 0=read
% (for traces with many channel, enables to load additional channels,
% like SXR, ECE, etc.)
% varargin{3}: zmag for varargout{1} computation
% varargin{4}: time range [t1 t2] (to limit data collected)
% varargin{5}: nth to keep only nth time points of traces
% varargin{6}: camera to use: 'B' (default), 'A',
%
% OUTPUT:
% trace.data: data structure
% trace.t: time of reference
% trace.x: space of reference
% .... others related to data
% error: error in loading signal (0=> OK, 1=> error)
%
% Additional Output arguments depending on data_type
%
% data_type=sxR:
% varargout{1}: intersection of the view lines with magnetic axis
%
% functions needed: SF2ML or mdsplus routines
%
% Example:
% [ip,error]=loadAUGdata(shot,'ip');
% [ip,error]=loadAUGdata(shot,'MAG/Ipi');
% [n2,error]=loadAUGdata(shot,'MOD/EvenNAmp');
%
% List of keywords (can be used in adition to 'DIAG/tracenam'), with comments when not obvious:
%
% Ip :
% b0 :
% zmag :
% rmag :
% rgeo :
% zgeo :
% q0 :
% q95 :
% kappa :
% delta :
% deltatop :
% deltabot :
% neint :
% neterho :
% cxrs : CXRS structure from CEZ with vrot, Ti, ...
% cxrs_rho : same as cxrs but project on rho as well (using 'equil' information)
% equil : equilibrium structure, rhopol, rhotor, rhovol, etc from EQI
% equil_fpp : as equil but from FPP
% equil_eqm : as equil but from EQM
% equil_eqr : as equil but from EQR
% equil_eqh : as equil but from EQH
% sxr : from SXR/A or B (from old stuff, not sure still OK)
% sxR : from SXR/A or B adding R of chords (from old stuff, not sure still OK)
% sxb : 'SXB/J' chords
% sxf : 'SXF/I' chords
% transp_EXPNAME : data from TRA using experiment name = "EXPNAME"
% ece :
% ece_rho :
% eced :
% eced_rho :
% eced_rmd :
% Halpha :
% pgyro : for each gyrotrons, power, freq, etc (ask for more)
% powers : power traces for each sources
%
varargout={cell(1,1)};
error=1;
shotfile_exp = 'AUGD';
% To allow multiple ways of writing a specific keyword, use data_type_eff within this routine
data_type_eff=data_type;
if size(data_type,1)==1
% there might be "/" in the name, add backslash to protect it, thus remove in the name afterwards:
inotok=findstr('\/',data_type);
iok=regexp(data_type,'[^\\]/');
iok=iok+1;
if ~isempty(inotok)
% remove \/ and construct index of real separators "/"
for ij=1:length(inotok)
ijk=find(iok>inotok(ij));
if ~isempty(ijk)
iok(ijk) = iok(ijk) - 1;
end
end
end
inobackslash=regexp(data_type,'[^\\]');
data_type=data_type(inobackslash);
index_slash = iok;
% i=findstr('/',data_type);
i = index_slash;
if length(i)==1
% assumes given a la 'MAG/Ipi'
data_type_eff=[{data_type(1:i(1)-1)} ; {data_type(i(1)+1:end)}];
elseif length(i)==2
% assumes given a la 'AUGD/MAG/Ipi' or 'ECED/RMP/Trad-A'
data_type_eff=[{data_type(i(1)+1:i(2)-1)} ; {data_type(i(2)+1:end)}];
shotfile_exp = data_type(1:i(1)-1);
elseif length(i)>2
disp(['more / than expected in tracename: length(i)= ' num2str(length(i))])
data_type
end
elseif isempty(data_type)
data_type_eff = ' ';
end
i_efitm=0;
i_ext=length(data_type_eff)+1;
name_ext='';
if size(data_type_eff,1)==1
data_type_eff_noext=data_type_eff(1:i_ext-1);
if ~isempty(strmatch(data_type_eff_noext,[{'ip'} {'i_p'} {'IP'} {'iP'} {'xip'}],'exact'))
data_type_eff_noext='Ip';
end
if ~isempty(strmatch(data_type_eff_noext,[{'B0'}],'exact'))
data_type_eff_noext='b0';
end
if ~isempty(strmatch(data_type_eff_noext,[{'Te'} {'t_e'} {'TE'} {'T_e'}],'exact'))
data_type_eff_noext='te';
end
if ~isempty(strmatch(data_type_eff_noext,[{'Ne'} {'n_e'} {'NE'} {'N_e'}],'exact'))
data_type_eff_noext='ne';
end
if ~isempty(strmatch(data_type_eff_noext,[{'Terho'}],'exact'))
data_type_eff_noext='terho';
end
if ~isempty(strmatch(lower(data_type_eff_noext),[{'neterho'}],'exact'))
data_type_eff_noext='neterho';
end
if ~isempty(strmatch(lower(data_type_eff_noext),[{'cxrs'} {'vrot'} {'ti'}],'exact'))
data_type_eff_noext='cxrs'; % load full CEZ structure
end
if ~isempty(strmatch(lower(data_type_eff_noext),[{'cxrs_rhos'} {'cxrs_rho'} {'cxrsrho'} {'cxrsrhos'}],'exact'))
data_type_eff_noext='cxrs_rho'; % load full CEZ structure
end
if ~isempty(strmatch(data_type_eff_noext,[{'SXR'}],'exact'))
data_type_eff_noext='sxr';
end
if ~isempty(strmatch(data_type_eff_noext,[{'SXB'} {'sxb'} {'Sxb'} {'sXb'} {'sxB'} {'SXb'}],'exact'))
data_type_eff_noext='sxb';
end
if ~isempty(strmatch(data_type_eff_noext,[{'SXF'} {'sxf'} {'Sxf'} {'sXf'} {'sxf'} {'SXf'}],'exact'))
data_type_eff_noext='sxf';
end
if ~isempty(strmatch(lower(data_type_eff_noext),[{'ssx_g'} {'ssx'}],'exact'))
data_type_eff_noext='ssx_g';
end
if ~isempty(strmatch(lower(data_type_eff_noext),[{'ssx_h'} {'ssx'}],'exact'))
data_type_eff_noext='ssx_h';
end
if ~isempty(strmatch(lower(data_type_eff_noext),[{'ssx_i'} {'ssx'}],'exact'))
data_type_eff_noext='ssx_i';
end
if ~isempty(strmatch(lower(data_type_eff_noext),[{'ssx_j'} {'ssx'}],'exact'))
data_type_eff_noext='ssx_j';
end
if ~isempty(strmatch(lower(data_type_eff_noext),[{'ssx'}],'exact'))
data_type_eff_noext='ssx';
end
if ~isempty(strmatch(lower(data_type_eff_noext(1:min(7,length(data_type_eff_noext)))),[{'transp_'}],'exact'))
shotfile_exp = data_type_eff_noext(8:end)
data_type_eff_noext='transp'
end
if ~isempty(strmatch(data_type_eff_noext,[{'ECE'}],'exact'))
data_type_eff_noext='ece';
end
if ~isempty(strmatch(lower(data_type_eff_noext),[{'ece_rho'}],'exact'))
data_type_eff_noext='ece_rho';
end
if ~isempty(strmatch(upper(data_type_eff_noext),[{'ECED'}],'exact'))
data_type_eff_noext='eced';
end
if ~isempty(strmatch(upper(data_type_eff_noext),[{'ECED_RHO'}],'exact'))
data_type_eff_noext='eced_rho';
end
if ~isempty(strmatch(upper(data_type_eff_noext),[{'ECED_RMD'}],'exact'))
data_type_eff_noext='eced_rmd';
end
if ~isempty(strmatch(data_type_eff_noext,[{'VOL'} {'volume'}],'exact'))
data_type_eff_noext='vol';
end
if ~isempty(strmatch(data_type_eff_noext,[{'q_0'} {'Q0'}],'exact'))
data_type_eff_noext='q0';
end
if ~isempty(strmatch(data_type_eff_noext,[{'q_95'} {'Q95'}],'exact'))
data_type_eff_noext='q95';
end
if ~isempty(strmatch(data_type_eff_noext,[{'elongation'} {'elon'}],'exact'))
data_type_eff_noext='kappa';
end
if ~isempty(strmatch(data_type_eff_noext,[{'triangularity'} {'triang'}],'exact'))
data_type_eff_noext='delta';
end
if ~isempty(strmatch(data_type_eff_noext,[{'deltaup'} {'deltau'} {'triangtop'} {'triangu'} {'triangup'}],'exact'))
data_type_eff_noext='deltatop';
end
if ~isempty(strmatch(data_type_eff_noext,[{'deltalow'} {'deltal'} {'triangbot'} {'triangl'} {'trianglow'}],'exact'))
data_type_eff_noext='deltabot';
end
if ~isempty(strmatch(data_type_eff_noext,[{'qrho_FPP'}],'exact'))
data_type_eff_noext='qrho_fpp';
end
if ~isempty(strmatch(data_type_eff_noext,[{'equil_FPP'}],'exact'))
data_type_eff_noext='equil_fpp';
end
if ~isempty(strmatch(lower(data_type_eff_noext),[{'equil_eqm'}],'exact'))
data_type_eff_noext='equil_eqm';
end
if ~isempty(strmatch(lower(data_type_eff_noext),[{'equil_eqr'}],'exact'))
data_type_eff_noext='equil_eqr';
end
if ~isempty(strmatch(lower(data_type_eff_noext),[{'equil_eqh'}],'exact'))
data_type_eff_noext='equil_eqh';
end
if ~isempty(strmatch(data_type_eff_noext,[{'Rmag'}],'exact'))
data_type_eff_noext='rmag';
end
if ~isempty(strmatch(data_type_eff_noext,[{'Rgeo'}],'exact'))
data_type_eff_noext='rgeo';
end
if ~isempty(strmatch(data_type_eff_noext,[{'Zmag'}],'exact'))
data_type_eff_noext='zmag';
end
if ~isempty(strmatch(data_type_eff_noext,[{'Zgeo'}],'exact'))
data_type_eff_noext='zgeo';
end
% $$$ if ~isempty(strmatch(data_type_eff_noext,[{'Rcont'}],'exact'))
% $$$ data_type_eff_noext='rcont';
% $$$ end
% $$$ if ~isempty(strmatch(data_type_eff_noext,[{'Zcont'}],'exact'))
% $$$ data_type_eff_noext='zcont';
% $$$ end
if ~isempty(strmatch(lower(data_type_eff_noext),[{'betan'} {'betn'} {'beta_n'}],'exact'))
data_type_eff_noext='betan';
end
if ~isempty(strmatch(lower(data_type_eff_noext),[{'betap'} {'betp'} {'betpol'} {'betapol'} {'beta_p'} {'beta_pol'}],'exact'))
data_type_eff_noext='betap';
end
if ~isempty(strmatch(data_type_eff_noext,[{'Ha'} {'ha'} {'Halpha'} {'halpha'}],'exact'))
data_type_eff_noext='Halpha';
end
if ~isempty(strmatch(lower(data_type_eff_noext),[{'pgyro'} {'pec'} {'pech'} {'pecrh'} {'p_ec'} {'p_gyro'}],'exact'))
data_type_eff_noext='pgyro';
end
if ~isempty(strmatch(lower(data_type_eff_noext),[{'powers'} {'ptot'} {'ptots'}],'exact'))
data_type_eff_noext='powers';
end
else
i_ext=length(data_type_eff{2})+1;
name_ext='';
data_type_eff_noext=data_type_eff{2}(1:i_ext-1);
end
% all keywords and corresponding case to run below
AUGkeywrdall=[{'Ip'} {'b0'} {'zmag'} {'rmag'} {'rgeo'} {'zgeo'} {'vol'} {'qrho'} {'qrho_fpp'} {'q0'} {'q95'} {'kappa'} ...
{'delta'} {'deltatop'} {'deltabot'} {'neint'} {'ne'} {'te'} ...
{'nerho'} {'neterho'} {'terho'} {'cxrs'} {'cxrs_rho'} {'equil'} {'equil_fpp'} {'equil_eqm'} ...
{'equil_eqr'} {'equil_eqh'} {'sxr'} {'sxR'} {'sxb'} {'sxf'} {'ssx_g'} {'ssx_h'} {'ssx_i'} {'ssx_j'} {'ssx'} ...
{'transp'} {'ece'} {'ece_rho'} {'eced'} {'eced_rho'} {'eced_rmd'} {'Halpha'} {'pgyro'} {'powers'} ...
{'betan'} {'betap'}]; % rm, not yet implemented: {'rcont'} {'zcont'}
AUGsig.iip=strmatch('Ip',AUGkeywrdall,'exact');
AUGsig.ib0=strmatch('b0',AUGkeywrdall,'exact');
AUGsig.izmag=strmatch('zmag',AUGkeywrdall,'exact');
AUGsig.irmag=strmatch('rmag',AUGkeywrdall,'exact');
AUGsig.irgeo=strmatch('rgeo',AUGkeywrdall,'exact');
AUGsig.izgeo=strmatch('zgeo',AUGkeywrdall,'exact');
% $$$ AUGsig.ircont=strmatch('rcont',AUGkeywrdall,'exact');
% $$$ AUGsig.izcont=strmatch('zcont',AUGkeywrdall,'exact');
AUGsig.ivol=strmatch('vol',AUGkeywrdall,'exact');
AUGsig.iqrho=strmatch('qrho',AUGkeywrdall,'exact');
AUGsig.iqrho_fpp=strmatch('qrho_fpp',AUGkeywrdall,'exact');
AUGsig.iequil=strmatch('equil',AUGkeywrdall,'exact');
AUGsig.iequil_fpp=strmatch('equil_fpp',AUGkeywrdall,'exact');
AUGsig.iequil_eqm=strmatch('equil_eqm',AUGkeywrdall,'exact');
AUGsig.iequil_eqr=strmatch('equil_eqr',AUGkeywrdall,'exact');
AUGsig.iequil_eqh=strmatch('equil_eqh',AUGkeywrdall,'exact');
AUGsig.iq0=strmatch('q0',AUGkeywrdall,'exact');
AUGsig.iq95=strmatch('q95',AUGkeywrdall,'exact');
AUGsig.ikappa=strmatch('kappa',AUGkeywrdall,'exact');
AUGsig.idelta=strmatch('delta',AUGkeywrdall,'exact');
AUGsig.ideltatop=strmatch('deltatop',AUGkeywrdall,'exact');
AUGsig.ideltabot=strmatch('deltabot',AUGkeywrdall,'exact');
AUGsig.ibetan=strmatch('betan',AUGkeywrdall,'exact');
AUGsig.ibetap=strmatch('betap',AUGkeywrdall,'exact');
AUGsig.ineint=strmatch('neint',AUGkeywrdall,'exact');
AUGsig.ine=strmatch('ne',AUGkeywrdall,'exact');
AUGsig.ite=strmatch('te',AUGkeywrdall,'exact');
AUGsig.inerho=strmatch('nerho',AUGkeywrdall,'exact');
AUGsig.iterho=strmatch('terho',AUGkeywrdall,'exact');
AUGsig.ineterho=strmatch('neterho',AUGkeywrdall,'exact');
AUGsig.icxrs=strmatch('cxrs',AUGkeywrdall,'exact');
AUGsig.icxrs_rho=strmatch('cxrs_rho',AUGkeywrdall,'exact');
AUGsig.isxr=strmatch('sxr',AUGkeywrdall,'exact');
AUGsig.isxR=strmatch('sxR',AUGkeywrdall,'exact');
AUGsig.isxb=strmatch('sxb',AUGkeywrdall,'exact');
AUGsig.isxf=strmatch('sxf',AUGkeywrdall,'exact');
AUGsig.issx_g=strmatch('ssx_g',AUGkeywrdall,'exact');
AUGsig.issx_h=strmatch('ssx_h',AUGkeywrdall,'exact');
AUGsig.issx_i=strmatch('ssx_i',AUGkeywrdall,'exact');
AUGsig.issx_j=strmatch('ssx_j',AUGkeywrdall,'exact');
AUGsig.issx=strmatch('ssx',AUGkeywrdall,'exact');
AUGsig.itransp=strmatch('transp',AUGkeywrdall,'exact');
AUGsig.iece=strmatch('ece',AUGkeywrdall,'exact');
AUGsig.ieced=strmatch('eced',AUGkeywrdall,'exact');
AUGsig.iece_rho=strmatch('ece_rho',AUGkeywrdall,'exact');
AUGsig.ieced_rho=strmatch('eced_rho',AUGkeywrdall,'exact');
AUGsig.ieced_rmd=strmatch('eced_rmd',AUGkeywrdall,'exact');
AUGsig.iHalpha=strmatch('Halpha',AUGkeywrdall,'exact');
AUGsig.ipgyro=strmatch('pgyro',AUGkeywrdall,'exact');
AUGsig.ipowers=strmatch('powers',AUGkeywrdall,'exact');
% For each keyword, specify which case to use. As most common is 'simplereaddata', fill in with this and change
% only indices needed. Usually use name of case same as keyword name
AUGkeywrdcase=cell(size(AUGkeywrdall));
AUGkeywrdcase(:)={'simplereaddata'};
AUGkeywrdcase(AUGsig.iqrho)=AUGkeywrdall(AUGsig.iqrho); % special as efit q on psi
AUGkeywrdcase(AUGsig.iqrho_fpp)=AUGkeywrdall(AUGsig.iqrho_fpp); % special as efit q on psi
AUGkeywrdcase(AUGsig.iequil)=AUGkeywrdall(AUGsig.iequil); % special as efit q on psi
AUGkeywrdcase(AUGsig.iequil_fpp)=AUGkeywrdall(AUGsig.iequil_fpp); % special as efit q on psi
AUGkeywrdcase(AUGsig.iequil_eqm)=AUGkeywrdall(AUGsig.iequil_eqm); % special as efit q on psi
AUGkeywrdcase(AUGsig.iequil_eqr)=AUGkeywrdall(AUGsig.iequil_eqr); % special as efit q on psi
AUGkeywrdcase(AUGsig.iequil_eqh)=AUGkeywrdall(AUGsig.iequil_eqh);
%AUGkeywrdcase(AUGsig.idelta)=AUGkeywrdall(AUGsig.idelta); % special as average of triu and tril
AUGkeywrdcase(AUGsig.ine)=AUGkeywrdall(AUGsig.ine); % special as adds error bars
AUGkeywrdcase(AUGsig.ite)=AUGkeywrdall(AUGsig.ite); % idem
AUGkeywrdcase(AUGsig.inerho)=AUGkeywrdall(AUGsig.inerho); % idem
AUGkeywrdcase(AUGsig.ineterho)=AUGkeywrdall(AUGsig.ineterho); % idem
AUGkeywrdcase(AUGsig.iterho)=AUGkeywrdall(AUGsig.iterho); % idem
AUGkeywrdcase(AUGsig.isxr)=AUGkeywrdall(AUGsig.isxr);
AUGkeywrdcase(AUGsig.isxR)=AUGkeywrdall(AUGsig.isxR);
AUGkeywrdcase(AUGsig.isxb)=AUGkeywrdall(AUGsig.isxb);
AUGkeywrdcase(AUGsig.isxf)=AUGkeywrdall(AUGsig.isxf);
AUGkeywrdcase(AUGsig.issx_g)=AUGkeywrdall(AUGsig.issx_g);
AUGkeywrdcase(AUGsig.issx_h)=AUGkeywrdall(AUGsig.issx_h);
AUGkeywrdcase(AUGsig.issx_i)=AUGkeywrdall(AUGsig.issx_i);
AUGkeywrdcase(AUGsig.issx_j)=AUGkeywrdall(AUGsig.issx_j);
AUGkeywrdcase(AUGsig.issx)=AUGkeywrdall(AUGsig.issx);
AUGkeywrdcase(AUGsig.itransp)={data_type_eff_noext};
AUGkeywrdcase(AUGsig.iece)=AUGkeywrdall(AUGsig.iece);
AUGkeywrdcase(AUGsig.ieced)=AUGkeywrdall(AUGsig.ieced);
AUGkeywrdcase(AUGsig.iece_rho)=AUGkeywrdall(AUGsig.iece_rho);
AUGkeywrdcase(AUGsig.ieced_rho)=AUGkeywrdall(AUGsig.ieced_rho);
AUGkeywrdcase(AUGsig.icxrs)=AUGkeywrdall(AUGsig.icxrs);
AUGkeywrdcase(AUGsig.icxrs_rho)=AUGkeywrdall(AUGsig.icxrs_rho);
AUGkeywrdcase(AUGsig.ipgyro)=AUGkeywrdall(AUGsig.ipgyro); % idem
AUGkeywrdcase(AUGsig.ipowers)=AUGkeywrdall(AUGsig.ipowers); % idem
% Information about which dimension has time, always return 2D data as (x,t) array
% as most are 1D arrays with time as first index, fill in with ones and change only those needed
AUGsigtimeindx=ones(size(AUGkeywrdall));
% For the 'simplereaddata' cases, we need the full node in case gdat was called with full location directly
% for the other cases, leave this location empty
AUGsiglocation=cell(2,size(AUGkeywrdall,2));
AUGsiglocation(:)={''};
AUGsiglocation(:,AUGsig.iip)={'MAG'; 'Ipa'};
AUGsiglocation(:,AUGsig.ib0)={'FPC'; 'BTF'}; % at 1.65m (?)
AUGsiglocation(:,AUGsig.izmag)={'FPG'; 'Zmag'};
AUGsiglocation(:,AUGsig.irmag)={'FPG'; 'Rmag'};
AUGsiglocation(:,AUGsig.irgeo)={'FPG'; 'Rgeo'};
AUGsiglocation(:,AUGsig.izgeo)={'FPG'; 'Zgeo'};
% $$$ AUGsiglocation(:,AUGsig.ircont)={'' ; ''}; AUGsigtimeindx(AUGsig.ircont)=2;
% $$$ AUGsiglocation(:,AUGsig.izcont)={'' ; ''}; AUGsigtimeindx(AUGsig.izcont)=2;
AUGsiglocation(:,AUGsig.ivol)={'FPG'; 'Vol'};
AUGsiglocation(:,AUGsig.iq0)={'FPG'; 'q0'};
AUGsiglocation(:,AUGsig.iq95)={'FPG'; 'q95'};
AUGsiglocation(:,AUGsig.ikappa)={'FPG'; 'k'};
AUGsiglocation(:,AUGsig.ideltatop)={'FPG'; 'delRoben'};
AUGsiglocation(:,AUGsig.ideltabot)={'FPG'; 'delRuntn'};
AUGsiglocation(:,AUGsig.ibetan)={'TOT'; 'beta_N'};
AUGsiglocation(:,AUGsig.ibetap)={'FPG'; 'betpol'};
AUGsiglocation(:,AUGsig.ineint)={'DCN'; 'H-1'};
AUGsiglocation(:,AUGsig.iece)={'CEC'; 'Trad-A'};
AUGsiglocation(:,AUGsig.ieced)={'CEC'; 'Trad-A'}; % ECED
AUGsiglocation(:,AUGsig.iece_rho)={'CEC'; 'Trad-A'};
AUGsiglocation(:,AUGsig.ieced_rho)={'CEC'; 'Trad-A'}; % ECED
AUGsiglocation(:,AUGsig.ieced_rmd)={'RMD'; 'Trad-A'}; % ECED_RMD
AUGsiglocation(:,AUGsig.iHalpha)={'POT'; 'ELMa-Han'};
% For the 'simplereaddata' cases, we need the full node in case gdat was called with full location directly
% for the other cases, leave this location empty
AUGexplocation=cell(size(AUGkeywrdall,2),1);
AUGexplocation(:)={'AUGD'};
% cases with a "private" shotfile:
AUGexplocation(AUGsig.ieced)={'ECED'};
AUGexplocation(AUGsig.ieced_rho)={'ECED'};
AUGexplocation(AUGsig.ieced_rmd)={'ECED'};
AUGexplocation(AUGsig.itransp)={shotfile_exp};
if shot==-9
clear trace
for i=1:length(AUGkeywrdall)
fieldname_eff = lower(AUGkeywrdall{i});
keyword_eff = AUGkeywrdall{i};
ij=findstr(fieldname_eff,':');
if ~isempty(ij);
fieldname_eff(ij)='_';
keyword_eff(ij:end+1)=[':i' keyword_eff(ij+1:end)] ;
end
trace.(fieldname_eff) = keyword_eff;
end
% add example for Ip trace full call
trace.ipfullcall = 'MAG/Ipa';
% trace.data=[];
return
end
% initialize order of substructures and allows just a "return" if data empty
trace.data=[];
trace.x=[];
trace.t=[];
trace.dim=[];
trace.dimunits=[];
trace.name=[];
% find index of signal called upon
if size(data_type_eff,1)==2
% in case node name was given in 2 parts directly (as old way)
ii1=strmatch(data_type_eff(1),AUGsiglocation(1,:),'exact');
iiindex=strmatch(data_type_eff_noext,AUGsiglocation(2,ii1),'exact');
if ~isempty(iiindex)
index=ii1(iiindex);
else
index=[];
end
if isempty(index)
% $$$ disp('********************')
% $$$ disp('trace not yet registered.')
% $$$ disp('If standard data, ask andrea.scarabosio@epfl.ch or olivier.sauter@epfl.ch to create a keyqord entry for this data')
% eval(['!mail -s ''' data_type_eff{1} ' ' data_type_eff{2} ' ' num2str(shot) ' ' ...
% getenv('USER') ' AUG'' olivier.sauter@epfl.ch < /dev/null'])
disp('********************')
% temporarily add entry in arrays, so can work below
index=length(AUGkeywrdall)+1;
AUGkeywrdall(end+1)={'new'};
AUGkeywrdcase(end+1)={'simplereaddata'};
AUGsiglocation(1:2,end+1)=[data_type_eff(1) ; {data_type_eff_noext}];
AUGexplocation{end+1}=shotfile_exp;
AUGsigtimeindx(end+1)=0;
AUGkeywrdcase{index}
elseif length(index)>1 || ~strcmp(AUGkeywrdcase{index},'simplereaddata')
% $$$ msgbox(['Problem in loadAUGdata with data_type_eff = ' char(data_type_eff(end)) ...
% $$$ '. Full paths of nodes should only be for case simplereaddata'],'in loadAUGdata','error')
% $$$ error('in loadAUGdata')
% chose a full path which coincides with a keyword
disp('the node path seems to coincide with the following keywords, which might be better to use:')
disp(' ')
disp(AUGkeywrdcase(index))
disp('proceed as is in any case')
index=length(AUGkeywrdall)+1;
AUGkeywrdall(end+1)={'new'};
AUGkeywrdcase(end+1)={'simplereaddata'};
AUGsiglocation(1:2,end+1)=[data_type_eff(1) ; {data_type_eff_noext}];
AUGexplocation{end+1}=shotfile_exp;
AUGsigtimeindx(end+1)=0;
AUGkeywrdcase{index}
end
else
index=strmatch(data_type_eff_noext,AUGkeywrdall,'exact');
if isempty(index)
disp(' ')
disp('********************')
if iscell(data_type_eff)
disp(['no such keyword: ' data_type_eff{1} '/' data_type_eff{2}])
else
disp(['no such keyword: ' data_type_eff])
end
disp(' ')
disp('Available keywords:')
AUGkeywrdall(:)
disp('********************')
return
end
end
disp(' ')
if iscell(data_type_eff)
disp(['loading' ' ' data_type_eff{1} '/' data_type_eff{2} ' from AUG shot #' num2str(shot)]);
else
disp(['loading' ' ' data_type_eff ' from AUG shot #' num2str(shot)]);
end
disp(['case ' AUGkeywrdcase{index}])
disp(' ')
switch AUGkeywrdcase{index}
%&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&
case 'simplereaddata'
ppftype=AUGsiglocation{1,index};
shotfile_exp_eff = AUGexplocation{index};
if i_efitm;
tracename=['eftm' AUGsiglocation{2,index}(5:end) name_ext];
else
tracename=[AUGsiglocation{2,index} name_ext];
end
ij=find(tracename~='''');
tracename=tracename(ij)
[a,e]=rdaAUG_eff(shot,ppftype,tracename,shotfile_exp_eff);
% switch tracename
% special cases if traces do not exist for some shot or other
% end
trace=a;
clear error
error=e;
if isempty(trace.data)
trace.dimunits=[];
elseif length(size(trace.data))==1 | (length(size(trace.data))==2 & size(trace.data,2)==1)
trace.dim=[{trace.t}];
trace.dimunits={'time [s]'};
elseif length(size(trace.data))==2
trace.dim=[{trace.x} ; {trace.t}];
trace.dimunits=[{'R [m] or rho=sqrt(psi_norm)'} ; {'time [s]'}];
else
% disp('how to deal with 3D arrays?')
% find time dimension
ilentime=length(trace.t);
ij=find(size(trace.value)==ilentime);
if ij==1;
% as expected
trace.x = [1:size(trace.value,2)];
trace.dim=[{trace.t} ; {trace.x} ; {[1:size(trace.value,3)]}];
trace.dimunits=[{'time [s]'} ; {''} ; {''}];
else
trace.dim=[{[]} ; {[]} ; {[]}];
trace.dim{ij} = trace.t;
trace.dimunits=[{[]} ; {[]} ; {[]}];
trace.dimunits{ij}='time [s]';
end
end
trace.name=[num2str(shot) '/' ppftype '/' tracename];
%&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&
case {'sxr','sxR'}
% LOAD MULTI CHANNEL DATA
% load AUG soft x-ray data
shotfile_exp_eff = AUGexplocation{index};
if nargin>=3 & ~isempty(varargin{1})
starti=varargin{1}(1);
endi=varargin{1}(2);
else
starti=1;
endi=30;
end
if nargin>=4 & ~isempty(varargin{2})
status=varargin{2};
else
status=ones(endi-starti+1,1);
end
if nargin>=6 & ~isempty(varargin{4})
timerange=varargin{4};
else
timerange=[1 7];
end
if nargin>=7 & ~isempty(varargin{5})
nth=varargin{5};
else
nth=13;
end
if nargin>=8 & ~isempty(varargin{6})
tracename=varargin{6};
else
tracename='B';
end
trace.t=[];
trace.x=[];
ppftype='SXR';
[a,e]=rdaAUG_eff(shot,ppftype,tracename,shotfile_exp_eff,timerange);
trace=a;
trace.dim=[{[starti:endi]'} ; {trace.t}];
trace.x=trace.dim{1};
trace.dimunits=[{'channels'} ; {'time [s]'}];
trace.units='W/m^2';
trace.name=[num2str(shot) '/' ppftype '/' tracename];
% keep only nth points
trace.t=trace.t(1:nth:end);
trace.data=trace.data(:,1:nth:end);
trace.dim{2}=trace.t;
% calculating intersection of the view lines with magnetics axis
if strcmp(data_type_eff_noext,'sxR')
if nargin>=5 & ~isempty(varargin{3})
zmag=varargin{3};
else
zmag=loadAUGdata(shot,'zmag');
end
zmageff=interp1(zmag.t,zmag.data,trace.t);
if strcmp(tracename,'B')
[R_B, Z_B, ang_B,Rsxr]=sxrbgeometry(zmageff);
elseif strcmp(tracename,'A')
[R_A, Z_A, ang_A,Rsxr]=sxrageometry(zmageff);
else
disp(['sxr camera: ' tracename ' not set yet for calculating R projection'])
return
end
radius.data=Rsxr;
radius.t=trace.t;
varargout{1}={radius};
trace.R=radius;
end
%&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&
case {'sxb', 'sxf', 'ssx', 'ssx_g', 'ssx_h', 'ssx_i', 'ssx_j'}
% LOAD MULTI CHANNEL DATA SXB/J_0xx (or other than J camera if specified in varargin{8})
% load AUG soft x-ray data
shotfile_exp_eff = AUGexplocation{index};
if nargin>=3 & ~isempty(varargin{1})
% chords to be loaded
starti=varargin{1}(1);
endi=varargin{1}(2);
elseif strcmp(AUGkeywrdcase{index}(1:3),'ssx')
starti=1;
endi=60;
else
starti=52;
endi=54;
end
if nargin>=4 & ~isempty(varargin{2})
% chords already loaded, 1=not loaded (to load), 0=loaded
status=varargin{2};
else
status=ones(endi-starti+1,1);
end
if nargin>=6 & ~isempty(varargin{4})
timerange=varargin{4};
else
timerange=[0 10];
end
if nargin>=7 & ~isempty(varargin{5})
nth=varargin{5};
else
nth=13;
end
if nargin>=8 & ~isempty(varargin{6})
tracename=varargin{6};
else
if strcmp(AUGkeywrdcase{index},'sxb')
tracename='J';
elseif strcmp(AUGkeywrdcase{index},'sxf')
tracename='I';
elseif strcmp(AUGkeywrdcase{index},'ssx')
tracename='G';
elseif strcmp(AUGkeywrdcase{index}(1:4),'ssx_')
tracename=upper(AUGkeywrdcase{index}(5));
else
disp('should not be here, ask O. Sauter');
end
end
trace.t=[];
trace.x=[];
ppftype=upper(AUGkeywrdcase{index}(1:3));
iok=0;
for ichord=starti:endi
tracename_eff = [tracename '_' num2str(ichord,'%.3d')];
try
[a,e]=rdaAUG_eff(shot,ppftype,tracename_eff,shotfile_exp_eff,timerange);
catch
a = [];
end
if isempty(a) || e~=0
if ~exist('trace_all')
trace_all = struct([]);
else
end
else
if iok==0
trace_all = a;
trace_all = rmfield(trace_all,[{'value'},{'data'}]);
iok = iok+1;
else
iok = iok+1;
end
trace_all.value(ichord,:) = a.value;
trace_all.data(ichord,:) = a.data;
end
end
if ~isempty(trace_all)
trace_all.dim=[{[starti:endi]'} ; {trace.t}];
trace = trace_all;
trace.x=trace.dim{1};
trace.dimunits=[{'channels'} ; {'time [s]'}];
trace.units='W/m^2';
trace.name=[num2str(shot) '/' ppftype '/' tracename];
% keep only nth points
trace.t=trace.t(1:nth:end);
trace.data=trace.data(:,1:nth:end);
trace.dim{2}=trace.t;
trace.value=trace.value(:,1:nth:end);
trace.time_aug.value=trace.time_aug.value(1:nth:end);
else
trace.data = [];
trace.dim = [];
trace.dimunits = [];
trace.x = [];
trace.t = [];
trace.units = [];
trace.name=[num2str(shot) '/' ppftype '/' tracename];
end
%&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&
case {'ece','eced','ece_rho','eced_rho'}
% LOAD MULTI CHANNEL DATA
% load AUG ece data
shotfile_exp_eff = AUGexplocation{index};
if nargin>=3 & ~isempty(varargin{1})
starti=varargin{1}(1);
endi=varargin{1}(2);
else
starti=1;
endi=30;
end
if nargin>=4 & ~isempty(varargin{2})
status=varargin{2};
else
status=ones(endi-starti+1,1);
end
if nargin>=6 & ~isempty(varargin{4})
timerange=varargin{4};
else
timerange=[0 10];
end
if nargin>=7 & ~isempty(varargin{5})
nth=varargin{5};
else
nth=1;
end
trace.t=[];
trace.x=[];
ppftype=AUGsiglocation{1,index};
tracename=AUGsiglocation{2,index};
[a,e]=rdaAUG_eff(shot,ppftype,tracename,shotfile_exp_eff,timerange);
starti=1;
endi=size(a.value,1);
trace=a;
trace.dim=[{[starti:endi]'} ; {trace.t}];
trace.x=trace.dim{1};
trace.dimunits=[{'channels'} ; {'time [s]'}];
trace.units='W/m^2';
trace.name=[num2str(shot) '/' ppftype '/' tracename];
% keep only nth points
trace.t=trace.t(1:nth:end);
trace.data=trace.data(:,1:nth:end);
trace.dim{2}=trace.t;
% get R
try
[aR,e]=rdaAUG_eff(shot,ppftype,'R-A',shotfile_exp_eff,timerange);
catch
end
try
[aZ,e]=rdaAUG_eff(shot,ppftype,'z-A',shotfile_exp_eff,timerange); % (problem now)
%aZ=aR;
%aZ.data=ones(size(aR.data));
catch
disp(['problem with getting z-A in ' ppftype])
end
domatchRtime=0;
if domatchRtime
% interpolate R structure on ece data time array, to ease plot vs R
for i=starti:endi
radius.data(i,:) = interp1(aR.t,aR.data(i,:),trace.t);
zheight.data(i,:) = interp1(aZ.t,aZ.data(i,:),trace.t);
end
radius.t=trace.t;
zheight.t=trace.t;
else
radius.data = aR.data;
radius.t=aR.t;
zheight.data = aZ.data;
zheight.t=aR.t;
end
ij=find(trace.data==0);
trace.data(ij)=NaN;
varargout{1}={radius};
trace.R=radius;
trace.Z=zheight;
if strcmp(AUGkeywrdcase{index},'ece_rho') || strcmp(AUGkeywrdcase{index},'eced_rho')
equil=gdat(shot,'equil',0);
inb_chord_ece=size(trace.R.data,1);
inb_time_ece=size(trace.R.data,2);
psi_out = NaN*ones(inb_chord_ece,inb_time_ece);
rhopsinorm_out = NaN*ones(inb_chord_ece,inb_time_ece);
rhotornorm_out = NaN*ones(inb_chord_ece,inb_time_ece);
rhovolnorm_out = NaN*ones(inb_chord_ece,inb_time_ece);
% constructs intervals within which a given equil is used: [time_equil(i),time_equil(i+1)]
time_equil=[1.5*equil.t(1)-0.5*equil.t(2) 0.5.*(equil.t(1:end-1)+equil.t(2:end)) 1.5*equil.t(end)-0.5*equil.t(end-1)];
iok=find(aR.data(:,1)>0);
for itequil=1:length(time_equil)-1
rr=equil.Rmesh(:,itequil);
zz=equil.Zmesh(:,itequil);
psirz_in = equil.psi2D(:,:,itequil);
it_ece_inequil = find(trace.R.t>=time_equil(itequil) & trace.R.t<=time_equil(itequil+1));
if ~isempty(it_ece_inequil)
rout=trace.R.data(iok,it_ece_inequil);
ijok=find(~isnan(rout));
if ~isempty(ijok)
zout=trace.Z.data(iok,it_ece_inequil);
psi_at_routzout = interpos2Dcartesian(rr,zz,psirz_in,rout,zout);
psi_out(iok,it_ece_inequil) = reshape(psi_at_routzout,length(iok),length(it_ece_inequil));
rhopsinorm_out(iok,it_ece_inequil) = sqrt(abs((psi_out(iok,it_ece_inequil)-equil.psi_axis(itequil))./(equil.psi_lcfs(itequil)-equil.psi_axis(itequil))));
for it_cx=1:length(it_ece_inequil)
rhotornorm_out(iok,it_ece_inequil(it_cx)) = interpos(equil.rhopolnorm(:,itequil),equil.rhotornorm(:,itequil),rhopsinorm_out(iok,it_ece_inequil(it_cx)),-3,[2 2],[0 1]);
rhovolnorm_out(iok,it_ece_inequil(it_cx)) = interpos(equil.rhopolnorm(:,itequil),equil.rhovolnorm(:,itequil),rhopsinorm_out(iok,it_ece_inequil(it_cx)),-3,[2 2],[0 1]);
end
end
end
end
trace.rhos.psi_on_rztime = psi_out;
trace.rhos.rhopsinorm_on_rztime = rhopsinorm_out;
trace.rhos.rhotornorm_on_rztime = rhotornorm_out;
trace.rhos.rhovolnorm_on_rztime = rhovolnorm_out;
trace.rhos.t = trace.R.t;
end
%&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&
case {'cxrs', 'cxrs_rho'}
% LOAD CEZ data
%
shotfile_exp_eff = AUGexplocation{index};
trace.t=[];
trace.x=[];
ppftype='CEZ';
[a,e]=rdaAUG_eff(shot,ppftype,'vrot',shotfile_exp_eff);
if isempty(a) || e~=0
cxrs = struct([]);
else
cxrs.vrot = a;
rmfield(cxrs.vrot,'x'); rmfield(cxrs.vrot,'t'); rmfield(cxrs.vrot,'data');
cxrs.x = a.x;
cxrs.data = a.data;
cxrs.t = a.t;
cxrs.x = a.x;
end
if ~isempty(cxrs)
cxrs.dim=[{cxrs.x} ; {cxrs.t}];
cxrs.dimunits=[{'chord'} ; {'time [s]'}];
cxrs.name=[num2str(shot) '/' ppftype '/vrot;Ti;Ti_c'];
[aerr,e]=rdaAUG_eff(shot,ppftype,'err_vrot',shotfile_exp_eff);
cxrs.vrot.error = aerr.value;
[r_time]=sf2ab(ppftype,shot,'R_time','-exp',shotfile_exp_eff);
cxrs.r_time = r_time.value{1};
[z_time]=sf2ab(ppftype,shot,'z_time','-exp',shotfile_exp_eff);
cxrs.z_time = z_time.value{1};
[a,e]=rdaAUG_eff(shot,ppftype,'Ti',shotfile_exp_eff);
[aerr,e]=rdaAUG_eff(shot,ppftype,'err_Ti',shotfile_exp_eff);
cxrs.ti = a;
cxrs.ti.error = aerr.value;
[a,e]=rdaAUG_eff(shot,ppftype,'Ti_c',shotfile_exp_eff);
[aerr,e]=rdaAUG_eff(shot,ppftype,'err_Ti_c',shotfile_exp_eff);
cxrs.ti_c = a;
cxrs.ti_c.error = aerr.value;
%
if strcmp(AUGkeywrdcase{index},'cxrs_rho')
equil=gdat(shot,'equil',0);
inb_chord_cxrs=size(cxrs.r_time,1);
inb_time_cxrs=size(cxrs.r_time,2);
psi_out = NaN*ones(inb_chord_cxrs,inb_time_cxrs);
rhopsinorm_out = NaN*ones(inb_chord_cxrs,inb_time_cxrs);
rhotornorm_out = NaN*ones(inb_chord_cxrs,inb_time_cxrs);
rhovolnorm_out = NaN*ones(inb_chord_cxrs,inb_time_cxrs);
% constructs intervals within which a given equil is used: [time_equil(i),time_equil(i+1)]
time_equil=[1.5*equil.t(1)-0.5*equil.t(2) 0.5.*(equil.t(1:end-1)+equil.t(2:end)) 1.5*equil.t(end)-0.5*equil.t(end-1)];
iok=find(cxrs.r_time(:,1)>0);
for itequil=1:length(time_equil)-1
rr=equil.Rmesh(:,itequil);
zz=equil.Zmesh(:,itequil);
psirz_in = equil.psi2D(:,:,itequil);
it_cxrs_inequil = find(cxrs.t>=time_equil(itequil) & cxrs.t<=time_equil(itequil+1));
if ~isempty(it_cxrs_inequil)
rout=cxrs.r_time(iok,it_cxrs_inequil);
zout=cxrs.z_time(iok,it_cxrs_inequil);
psi_at_routzout = interpos2Dcartesian(rr,zz,psirz_in,rout,zout);
psi_out(iok,it_cxrs_inequil) = reshape(psi_at_routzout,length(iok),length(it_cxrs_inequil));
rhopsinorm_out(iok,it_cxrs_inequil) = sqrt((psi_out(iok,it_cxrs_inequil)-equil.psi_axis(itequil))./(equil.psi_lcfs(itequil)-equil.psi_axis(itequil)));
for it_cx=1:length(it_cxrs_inequil)
rhotornorm_out(iok,it_cxrs_inequil(it_cx)) = interpos(equil.rhopolnorm(:,itequil),equil.rhotornorm(:,itequil),rhopsinorm_out(iok,it_cxrs_inequil(it_cx)),-3,[2 2],[0 1]);
rhovolnorm_out(iok,it_cxrs_inequil(it_cx)) = interpos(equil.rhopolnorm(:,itequil),equil.rhovolnorm(:,itequil),rhopsinorm_out(iok,it_cxrs_inequil(it_cx)),-3,[2 2],[0 1]);
end
end
end
cxrs.psi_on_rztime = psi_out;
cxrs.rhopsinorm_on_rztime = rhopsinorm_out;
cxrs.rhotornorm_on_rztime = rhotornorm_out;
cxrs.rhovolnorm_on_rztime = rhovolnorm_out;
end
trace = cxrs;
else
trace.data = [];
trace.dim = [];
trace.dimunits = [];
trace.x = [];
trace.t = [];
trace.units = [];
trace.name=[num2str(shot) '/' ppftype '/' 'vrot, Ti_c,...'];
end
%&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&
case {'neterho'} % 'nerho', 'terho',
% LOAD VTA data
% Vertical Thomson core and edge
shotfile_exp_eff = AUGexplocation{index};
trace.t=[];
trace.x=[];
ppftype='VTA';
% if strcmp(AUGkeywrdcase{index},'terho')
[a,e]=rdaAUG_eff(shot,ppftype,'Te_c',shotfile_exp_eff);
if isempty(a) || e~=0
thomson = struct([]);
else
thomson.te_core = a;
rmfield(thomson.te_core,'x'); rmfield(thomson.te_core,'t'); rmfield(thomson.te_core,'data');
thomson.x = a.x;
thomson.data = a.data;
thomson.t = a.t;
thomson.time_core = a.t;
thomson.x = a.x;
end
if ~isempty(thomson)
thomson.dim=[{thomson.x} ; {thomson.t}];
thomson.dimunits=[{'points'} ; {'time [s]'}];
thomson.name=[num2str(shot) '/' ppftype '/te_c;te_e;ne_c;ne_e'];
[alow,e]=rdaAUG_eff(shot,ppftype,'Telow_c',shotfile_exp_eff);
[aup,e]=rdaAUG_eff(shot,ppftype,'Teupp_c',shotfile_exp_eff);
thomson.te_core.error = max(aup.value-thomson.te_core.value,thomson.te_core.value-alow.value);
%
[a,e]=rdaAUG_eff(shot,ppftype,'Ne_c',shotfile_exp_eff);
thomson.ne_core = a;
rmfield(thomson.ne_core,'x'); rmfield(thomson.ne_core,'t'); rmfield(thomson.ne_core,'data');
[alow,e]=rdaAUG_eff(shot,ppftype,'Nelow_c',shotfile_exp_eff);
[aup,e]=rdaAUG_eff(shot,ppftype,'Neupp_c',shotfile_exp_eff);
thomson.ne_core.error = max(aup.value-thomson.ne_core.value,thomson.ne_core.value-alow.value);
%
[r_time]=rdaAUG_eff(shot,ppftype,'R_core',shotfile_exp_eff);
thomson.r_core_time = ones(size(thomson.te_core.value,1),1)*r_time.value;
[z_time]=rdaAUG_eff(shot,ppftype,'Z_core',shotfile_exp_eff);
thomson.z_core_time = z_time.value' * ones(1,size(thomson.te_core.value,2));
%
[a,e]=rdaAUG_eff(shot,ppftype,'Te_e',shotfile_exp_eff);
thomson.te_edge = a;
thomson.time_edge = a.t;
rmfield(thomson.te_edge,'x'); rmfield(thomson.te_edge,'t'); rmfield(thomson.te_edge,'data');
[alow,e]=rdaAUG_eff(shot,ppftype,'Telow_e',shotfile_exp_eff);
[aup,e]=rdaAUG_eff(shot,ppftype,'Teupp_e',shotfile_exp_eff);
thomson.te_edge.error = max(aup.value-thomson.te_edge.value,thomson.te_edge.value-alow.value);
%
[a,e]=rdaAUG_eff(shot,ppftype,'Ne_e',shotfile_exp_eff);
thomson.ne_edge = a;
rmfield(thomson.ne_edge,'x'); rmfield(thomson.ne_edge,'t'); rmfield(thomson.ne_edge,'data');
[alow,e]=rdaAUG_eff(shot,ppftype,'Nelow_e',shotfile_exp_eff);
[aup,e]=rdaAUG_eff(shot,ppftype,'Neupp_e',shotfile_exp_eff);
thomson.ne_edge.error = max(aup.value-thomson.ne_edge.value,thomson.ne_edge.value-alow.value);
%
[r_time]=rdaAUG_eff(shot,ppftype,'R_edge',shotfile_exp_eff);
thomson.r_edge_time = ones(size(thomson.te_edge.value,1),1)*r_time.value;
[z_time]=rdaAUG_eff(shot,ppftype,'Z_edge',shotfile_exp_eff);
thomson.z_edge_time = z_time.value' * ones(1,size(thomson.te_edge.value,2));
%
if strcmp(AUGkeywrdcase{index},'neterho')
equil=gdat(shot,'equil',0);
% core
inb_chord_thomson_core=size(thomson.r_core_time,1);
inb_time_thomson_core=size(thomson.r_core_time,2);
psi_out_core = NaN*ones(inb_chord_thomson_core,inb_time_thomson_core);
rhopsinorm_out_core = NaN*ones(inb_chord_thomson_core,inb_time_thomson_core);
rhotornorm_out_core = NaN*ones(inb_chord_thomson_core,inb_time_thomson_core);
rhovolnorm_out_core = NaN*ones(inb_chord_thomson_core,inb_time_thomson_core);
% edge
inb_chord_thomson_edge=size(thomson.r_edge_time,1);
inb_time_thomson_edge=size(thomson.r_edge_time,2);
psi_out_edge = NaN*ones(inb_chord_thomson_edge,inb_time_thomson_edge);
rhopsinorm_out_edge = NaN*ones(inb_chord_thomson_edge,inb_time_thomson_edge);
rhotornorm_out_edge = NaN*ones(inb_chord_thomson_edge,inb_time_thomson_edge);
rhovolnorm_out_edge = NaN*ones(inb_chord_thomson_edge,inb_time_thomson_edge);
% constructs intervals within which a given equil is used: [time_equil(i),time_equil(i+1)]
time_equil=[1.5*equil.t(1)-0.5*equil.t(2) 0.5.*(equil.t(1:end-1)+equil.t(2:end)) 1.5*equil.t(end)-0.5*equil.t(end-1)];
for itequil=1:length(time_equil)-1
rr=equil.Rmesh(:,itequil);
zz=equil.Zmesh(:,itequil);
psirz_in = equil.psi2D(:,:,itequil);
it_thomson_core_inequil = find(thomson.time_core>=time_equil(itequil) & thomson.time_core<=time_equil(itequil+1));
if ~isempty(it_thomson_core_inequil)
rout_core=thomson.r_core_time(:,it_thomson_core_inequil);
zout_core=thomson.z_core_time(:,it_thomson_core_inequil);
psi_at_routzout = interpos2Dcartesian(rr,zz,psirz_in,rout_core,zout_core);
psi_out_core(:,it_thomson_core_inequil) = reshape(psi_at_routzout,inb_chord_thomson_core,length(it_thomson_core_inequil));
rhopsinorm_out_core(:,it_thomson_core_inequil) = sqrt((psi_out_core(:,it_thomson_core_inequil)-equil.psi_axis(itequil))./(equil.psi_lcfs(itequil)-equil.psi_axis(itequil)));
for it_cx=1:length(it_thomson_core_inequil)
rhotornorm_out_core(:,it_thomson_core_inequil(it_cx)) = ...
interpos(equil.rhopolnorm(:,itequil),equil.rhotornorm(:,itequil),rhopsinorm_out_core(:,it_thomson_core_inequil(it_cx)),-3,[2 2],[0 1]);
rhovolnorm_out_core(:,it_thomson_core_inequil(it_cx)) = ...
interpos(equil.rhopolnorm(:,itequil),equil.rhovolnorm(:,itequil),rhopsinorm_out_core(:,it_thomson_core_inequil(it_cx)),-3,[2 2],[0 1]);
end
end
% edge
it_thomson_edge_inequil = find(thomson.time_edge>=time_equil(itequil) & thomson.time_edge<=time_equil(itequil+1));
if ~isempty(it_thomson_edge_inequil)
rout_edge=thomson.r_edge_time(:,it_thomson_edge_inequil);
zout_edge=thomson.z_edge_time(:,it_thomson_edge_inequil);
psi_at_routzout = interpos2Dcartesian(rr,zz,psirz_in,rout_edge,zout_edge);
psi_out_edge(:,it_thomson_edge_inequil) = reshape(psi_at_routzout,inb_chord_thomson_edge,length(it_thomson_edge_inequil));
rhopsinorm_out_edge(:,it_thomson_edge_inequil) = sqrt((psi_out_edge(:,it_thomson_edge_inequil)-equil.psi_axis(itequil))./(equil.psi_lcfs(itequil)-equil.psi_axis(itequil)));
for it_cx=1:length(it_thomson_edge_inequil)
rhotornorm_out_edge(:,it_thomson_edge_inequil(it_cx)) = ...
interpos(equil.rhopolnorm(:,itequil),equil.rhotornorm(:,itequil),rhopsinorm_out_edge(:,it_thomson_edge_inequil(it_cx)),-3,[2 2],[0 1]);
rhovolnorm_out_edge(:,it_thomson_edge_inequil(it_cx)) = ...
interpos(equil.rhopolnorm(:,itequil),equil.rhovolnorm(:,itequil),rhopsinorm_out_edge(:,it_thomson_edge_inequil(it_cx)),-3,[2 2],[0 1]);
end
end
end
thomson.core_psi_on_rztime = psi_out_core;
thomson.core_rhopsinorm_on_rztime = rhopsinorm_out_core;
thomson.core_rhotornorm_on_rztime = rhotornorm_out_core;
thomson.core_rhovolnorm_on_rztime = rhovolnorm_out_core;
thomson.edge_psi_on_rztime = psi_out_edge;
thomson.edge_rhopsinorm_on_rztime = rhopsinorm_out_edge;
thomson.edge_rhotornorm_on_rztime = rhotornorm_out_edge;
thomson.edge_rhovolnorm_on_rztime = rhovolnorm_out_edge;
end
trace = thomson;
else
trace.data = [];
trace.dim = [];
trace.dimunits = [];
trace.x = [];
trace.t = [];
trace.units = [];
trace.name=[num2str(shot) '/' ppftype '/' tracename];
end
%&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&
case {'te', 'ne'}
shotfile_exp_eff = AUGexplocation{index}
if strcmp(AUGkeywrdcase{index},'te')
[a,e]=rdaAUG_eff(shot,'YPR','Te',shotfile_exp_eff);
else
[a,e]=rdaAUG_eff(shot,'YPR','Ne',shotfile_exp_eff);
end
trace = a;
trace.data = a.value';
trace.dim{1} = trace.x;
trace.dim{2} = trace.t;
trace.units = trace.unit;
trace.dimunits{1} = a.area.name ;
trace.dimunits{2} = a.time_aug.unit;
case {'equil', 'equil_fpp', 'equil_eqm', 'equil_eqr', 'equil_eqh', 'qrho', 'qrho_fpp'}
% In shotfiles time is 1st index in 2D arrays while gdat convention is to have time as last dimension
% so change equil.xxx to (rho,time) on the way...
shotfile_exp_eff = AUGexplocation{index};
if strcmp(AUGkeywrdcase{index},'equil_fpp') || strcmp(AUGkeywrdcase{index},'qrho_fpp')
DIAG = 'FPP';
elseif strcmp(AUGkeywrdcase{index},'equil_eqm')
DIAG = 'EQM';
elseif strcmp(AUGkeywrdcase{index},'equil_eqr')
DIAG = 'EQR';
elseif strcmp(AUGkeywrdcase{index},'equil_eqh')
DIAG = 'EQH';
else
DIAG = 'EQI';
end
M_Rmesh_par = sf2par(DIAG,shot,'M','PARMV');
M_Rmesh = M_Rmesh_par.value + 1; % nb of points
N_Zmesh_par = sf2par(DIAG,shot,'N','PARMV');
N_Zmesh = N_Zmesh_par.value + 1; % nb of points
NTIME_par = sf2par(DIAG,shot,'NTIME','PARMV');
NTIME = NTIME_par.value; % nb of points
Lpf_par = rdaAUG_eff(shot,DIAG,'Lpf',shotfile_exp_eff);
% since June, nb of time points in EQ results is not consistent with NTIME and time
% It seems the first NTIME points are correct, so use this explicitely
NTIME_Lpf = length(Lpf_par.value);
if (NTIME < NTIME_Lpf)
disp('WARNING: nb of times points smaller then equil results, use first NTIME points')
elseif (NTIME > NTIME_Lpf)
disp('ERROR: nb of times points LARGER then equil results')
disp('this is unexpected, so stop there and ask Olivier.Sauter@epfl.ch')
return
end
Lpf_tot = Lpf_par.value(1:NTIME); % nb of points: 100000*nb_SOL points + nb_core
Lpf_SOL = fix(Lpf_tot/100000);
Lpf1_t = mod(Lpf_tot,100000)+1; % nb of Lpf points
% since Lpf depends on time, need to load all first and then loop over time for easier mapping
[qpsi,e]=rdaAUG_eff(shot,DIAG,'Qpsi',shotfile_exp_eff);
ndimrho = size(qpsi.data,2);
if ndimrho==NTIME_Lpf
% data seems to be transposed
ndimrho = size(qpsi.data,1);
itotransposeback = 1; % seems x,time inverted so transpose and exchange .x and .t
else
itotransposeback = 0;
end
qpsi=adapt_rda(qpsi,NTIME,ndimrho,itotransposeback);
ijnan=find(isnan(qpsi.value));
qpsi.value(ijnan)=0;
[psi_tree,e]=rdaAUG_eff(shot,DIAG,'PFL',shotfile_exp);
psi_tree=adapt_rda(psi_tree,NTIME,ndimrho,itotransposeback);
[phi_tree,e]=rdaAUG_eff(shot,DIAG,'TFLx',shotfile_exp);
phi_tree=adapt_rda(phi_tree,NTIME,ndimrho,itotransposeback);
[Vol,e]=rdaAUG_eff(shot,DIAG,'Vol',shotfile_exp);
Vol=adapt_rda(Vol,NTIME,2*ndimrho,itotransposeback);
[Area,e]=rdaAUG_eff(shot,DIAG,'Area',shotfile_exp);
Area=adapt_rda(Area,NTIME,2*ndimrho,itotransposeback);
[Ri,e]=rdaAUG_eff(shot,DIAG,'Ri',shotfile_exp);
Ri=adapt_rda(Ri,NTIME,M_Rmesh,itotransposeback);
[Zj,e]=rdaAUG_eff(shot,DIAG,'Zj',shotfile_exp);
Zj=adapt_rda(Zj,NTIME,N_Zmesh,itotransposeback);
[PFM_tree,e]=rdaAUG_eff(shot,DIAG,'PFM',shotfile_exp);
PFM_tree=adaptPFM_rda(PFM_tree,M_Rmesh,N_Zmesh,NTIME);
[Pres,e]=rdaAUG_eff(shot,DIAG,'Pres',shotfile_exp);
Pres=adapt_rda(Pres,NTIME,2*ndimrho,itotransposeback);
[Jpol,e]=rdaAUG_eff(shot,DIAG,'Jpol',shotfile_exp);
Jpol=adapt_rda(Jpol,NTIME,2*ndimrho,itotransposeback);
[FFP,e]=rdaAUG_eff(shot,DIAG,'FFP',shotfile_exp);
if ~isempty(FFP.value)
FFP=adapt_rda(FFP,NTIME,ndimrho,itotransposeback);
else
FFP.value=NaN*ones(NTIME,max(Lpf1_t));
end
if strcmp(DIAG,'EQI') || strcmp(DIAG,'EQH')
[Rinv,e]=rdaAUG_eff(shot,DIAG,'Rinv',shotfile_exp);
Rinv=adapt_rda(Rinv,NTIME,ndimrho,itotransposeback);
[R2inv,e]=rdaAUG_eff(shot,DIAG,'R2inv',shotfile_exp);
R2inv=adapt_rda(R2inv,NTIME,ndimrho,itotransposeback);
[Bave,e]=rdaAUG_eff(shot,DIAG,'Bave',shotfile_exp);
Bave=adapt_rda(Bave,NTIME,ndimrho,itotransposeback);
[B2ave,e]=rdaAUG_eff(shot,DIAG,'B2ave',shotfile_exp);
B2ave=adapt_rda(B2ave,NTIME,ndimrho,itotransposeback);
[FTRA,e]=rdaAUG_eff(shot,DIAG,'FTRA',shotfile_exp);
FTRA=adapt_rda(FTRA,NTIME,ndimrho,itotransposeback);
else
Rinv.value=[]; R2inv.value=[]; Bave.value=[]; B2ave.value=[]; FTRA.value=[];
end
[LPFx,e]=rdaAUG_eff(shot,DIAG,'LPFx',shotfile_exp);
LPFx.value=LPFx.value(1:NTIME); LPFx.data=LPFx.value; LPFx.t=LPFx.t(1:NTIME);
[PFxx,e]=rdaAUG_eff(shot,DIAG,'PFxx',shotfile_exp);
PFxx=adapt_rda(PFxx,NTIME,max(LPFx.value)+1,itotransposeback);
[RPFx,e]=rdaAUG_eff(shot,DIAG,'RPFx',shotfile_exp);
RPFx=adapt_rda(RPFx,NTIME,max(LPFx.value)+1,itotransposeback);
[zPFx,e]=rdaAUG_eff(shot,DIAG,'zPFx',shotfile_exp);
zPFx=adapt_rda(zPFx,NTIME,max(LPFx.value)+1,itotransposeback);
% seems "LCFS" q-value is far too large, limit to some max (when diverted)
max_qValue = 40.0; % Note could just put a NaN on LCFS value since ill-defined when diverted
for it=1:NTIME
Lpf1 = Lpf1_t(it);
% Qpsi and similar data is on (time,radius) with radius being: LCFS..Lpf_points dummy LCFS..SOL part
% change it to (radial,time) and use only Lpf+1 points up to LCFS
ijok=find(qpsi.value(:,1)); % note: eqr fills in only odd points radially
% set NaNs to zeroes
if qpsi.value(ijok(1),1)<0
equil.qvalue(:,it) = max(qpsi.value(it,Lpf1:-1:1)',-max_qValue);
else
equil.qvalue(:,it) = min(qpsi.value(it,Lpf1:-1:1)',max_qValue);
end
% get x values
equil.psi(:,it)=psi_tree.value(it,Lpf1:-1:1)';
equil.psi_axis(it)= equil.psi(1,it);
equil.psi_lcfs(it)= equil.psi(end,it);
equil.rhopolnorm(:,it) = sqrt(abs((equil.psi(:,it)-equil.psi_axis(it)) ./(equil.psi_lcfs(it)-equil.psi_axis(it))));
if strcmp(DIAG,'EQR');
% q value has only a few values and from center to edge, assume they are from central rhopol values on
% But they are every other point starting from 3rd
ijk=find(equil.qvalue(:,it)~=0);
if length(ijk)>2
% now shots have non-zero axis values in eqr
rhoeff=equil.rhopolnorm(ijk,it);
qeff=equil.qvalue(ijk,it); % radial order was already inverted above
if ijk(1)>1
rhoeff = [0.; rhoeff];
qeff = [qeff(1) ;qeff];
end
ij_nonan=find(~isnan(equil.rhopolnorm(:,it)));
qfit = zeros(size(equil.rhopolnorm(:,it)));
qfit(ij_nonan)=interpos(rhoeff,qeff,equil.rhopolnorm(ij_nonan,it),-0.01,[1 0],[0 0],[300; ones(size(qeff(1:end-1)))]);
else
qfit = zeros(size(equil.rhopolnorm(:,it)));
end
equil.qvalue(:,it) = qfit;
end
% get rhotor values
equil.phi(:,it) = phi_tree.value(it,Lpf1:-1:1)';
equil.rhotornorm(:,it) = sqrt(abs(equil.phi(:,it) ./ equil.phi(end,it)));
% get rhovol values
equil.vol(:,it)=Vol.value(it,2*Lpf1-1:-2:1)';
equil.dvoldpsi(:,it)=Vol.value(it,2*Lpf1:-2:2)'; % 2nd index are dV/dpsi
equil.rhovolnorm(:,it) = sqrt(abs(equil.vol(:,it) ./ equil.vol(end,it)));
equil.area(:,it)=Area.value(it,2*Lpf1-1:-2:1)';
equil.dareadpsi(:,it)=Area.value(it,2*Lpf1:-2:2)'; % 2nd index are dV/dpsi
equil.Rmesh(:,it) = Ri.value(it,1:M_Rmesh);
equil.Zmesh(:,it) = Zj.value(it,1:N_Zmesh);
equil.psi2D(1:M_Rmesh,1:N_Zmesh,it) = PFM_tree.value(1:M_Rmesh,1:N_Zmesh,it);
equil.pressure(:,it)=Pres.value(it,2*Lpf1-1:-2:1)';
equil.dpressuredpsi(:,it)=Pres.value(it,2*Lpf1:-2:2)'; % 2nd index are dV/dpsi
if ~isempty(Jpol.value)
equil.jpol(:,it)=Jpol.value(it,2*Lpf1-1:-2:1)';
equil.djpolpsi(:,it)=Jpol.value(it,2*Lpf1:-2:2)'; % 2nd index are dV/dpsi
else
equil.jpol = [];
equil.djpolpsi = [];
end
equil.ffprime(:,it) = FFP.value(it,Lpf1:-1:1)';
equil.Xpoints.psi(1:LPFx.value(it)+1,it) = PFxx.value(it,1:LPFx.value(it)+1);
equil.Xpoints.Rvalue(1:LPFx.value(it)+1,it) = RPFx.value(it,1:LPFx.value(it)+1);
equil.Xpoints.Zvalue(1:LPFx.value(it)+1,it) = zPFx.value(it,1:LPFx.value(it)+1);
if ~isempty(Rinv.value)
equil.rinv(:,it) = Rinv.value(it,Lpf1:-1:1)';
else
equil.rinv = [];
end
if ~isempty(R2inv.value)
equil.r2inv(:,it) = R2inv.value(it,Lpf1:-1:1)';
else
equil.r2inv = [];
end
if ~isempty(Bave.value)
equil.bave(:,it) = Bave.value(it,Lpf1:-1:1)';
else
equil.bave = [];
end
if ~isempty(B2ave.value)
equil.b2ave(:,it) = B2ave.value(it,Lpf1:-1:1)';
else
equil.b2ave = [];
end
if ~isempty(FTRA.value)
equil.ftra(:,it) = FTRA.value(it,Lpf1:-1:1)';
else
equil.ftra = [];
end
%
end
equil.x = equil.rhopolnorm;
%
[equil_Rcoil,e]=rdaAUG_eff(shot,DIAG,'Rcl',shotfile_exp);
equil.Rcoils=equil_Rcoil.value;
[equil_Zcoil,e]=rdaAUG_eff(shot,DIAG,'Zcl',shotfile_exp);
equil.Zcoils=equil_Zcoil.value;
% get time values
[equil_time,e]=rdaAUG_eff(shot,DIAG,'time',shotfile_exp);
equil.t = equil_time.value(1:NTIME);
%
[IpiPSI,e]=rdaAUG_eff(shot,DIAG,'IpiPSI',shotfile_exp);
equil.Ip = IpiPSI.value(1:NTIME);
%
equil.data = equil.qvalue; % put q in data
equil.unit=[]; % not applicable
equil.name_gdat = [AUGkeywrdcase{index} ' using ' DIAG];
trace = equil;
trace.dim{1} = trace.x;
trace.dim{2} = trace.t;
trace.units = trace.unit;
trace.dimunits{1} = '';
trace.dimunits{2} = 's';
%&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&
case {'pgyro'}
% LOAD MULTI CHANNEL DATA ECS
% powers, frequencies, etc
shotfile_exp_eff = AUGexplocation{index};
trace.t=[];
trace.x=[];
ppftype='ECS';
% pgyro tot in index=9
[a,e]=rdaAUG_eff(shot,ppftype,'PECRH',shotfile_exp_eff);
if isempty(a) || e~=0
trace_all = struct([]);
else
nb_timepoints = length(a.time_aug.value);
trace_all.pgyro = NaN*ones(nb_timepoints,9);
trace_all.freq_ech = NaN*ones(1,8);
trace_all.pgyro(:,9) = reshape(a.data,nb_timepoints,1);
trace_all.t = a.time_aug.value;
end
for i=1:4
% "old" ECRH1 gyrotrons: gyro 1 to 4 in pgyro
tracename_eff = ['PG' num2str(i)];
[a,e]=rdaAUG_eff(shot,ppftype,tracename_eff,shotfile_exp_eff);
if isempty(a) || e~=0
else
trace_all.ecrh1(i) = a;
trace_all.pgyro(:,i) = reshape(a.data,nb_timepoints,1);
end
try
a = sf2par('ECS',shot,'gyr_freq',['P_sy1_g' num2str(i)]);
catch
% gyr_freq not present (old shots for example)
a=[];
end
if isempty(a)
else
trace_all.freq_ecrh1(i) = a;
trace_all.freq_ech(i) = a.value;
end
try
a = sf2par('ECS',shot,'GPolPos',['P_sy1_g' num2str(i)]);
catch
% GPolPos not present
a=[];
end
if isempty(a)
else
trace_all.polpos_ecs(i) = a.value;
end
try
a = sf2par('ECS',shot,'GTorPos',['P_sy1_g' num2str(i)]);
catch
a=[];
end
if isempty(a)
else
trace_all.torpos_ecs(i) = a.value;
end
% "new" ECRH2 gyrotrons: gyro 5 to 8 in pgyro
tracename_eff = ['PG' num2str(i) 'N'];
[a,e]=rdaAUG_eff(shot,ppftype,tracename_eff,shotfile_exp_eff);
if isempty(a) || e~=0
else
trace_all.ecrh2(i) = a;
trace_all.pgyro(:,i+4) = reshape(a.data,nb_timepoints,1);
end
try
a = sf2par('ECS',shot,'gyr_freq',['P_sy2_g' num2str(i)]);
catch
a=[];
end
if isempty(a)
else
trace_all.freq_ecrh2(i) = a;
trace_all.freq_ech(i+4) = a.value;
end
try
a = sf2par('ECS',shot,'GPolPos',['P_sy2_g' num2str(i)]);
catch
a=[];
end
if isempty(a)
else
trace_all.polpos_ecs(i+4) = a.value;
end
try
a = sf2par('ECS',shot,'GTorPos',['P_sy2_g' num2str(i)]);
catch
a=[];
end
if isempty(a)
else
trace_all.torpos_ecs(i+4) = a.value;
end
[a,e]=rdaAUG_eff(shot,'ECN',['G' num2str(i) 'POL'],shotfile_exp_eff);
if isempty(a) || e~=0
else
trace_all.gpol_ecn(i+4) = a;
end
[a,e]=rdaAUG_eff(shot,'ECN',['G' num2str(i) 'TOR'],shotfile_exp_eff);
if isempty(a) || e~=0
else
trace_all.gtor_ecn(i+4) = a;
end
[a,e]=rdaAUG_eff(shot,'ECN',['G' num2str(i) 'PO4'],shotfile_exp_eff);
if isempty(a) || e~=0
else
trace_all.gpo4_ecn(i+4) = a;
end
[a,e]=rdaAUG_eff(shot,'ECN',['G' num2str(i) 'PO8'],shotfile_exp_eff);
if isempty(a) || e~=0
else
trace_all.gpo8_ecn(i+4) = a;
end
end
if ~isempty(trace_all)
trace_all.dim=[{trace_all.t} {[1:9]}];
trace_all.data = trace_all.pgyro;
trace = trace_all;
trace.x=trace.dim{2};
trace.dimunits=[{'time [s]'} {'ECRH1(1:4) ECRH2(1:4) ECtot'}];
trace.units='W';
trace.freq_ech_units = 'GHz';
trace.name=[num2str(shot) '/' ppftype '/' 'PGi and PGiN'];
else
trace.data = [];
trace.dim = [];
trace.dimunits = [];
trace.x = [];
trace.t = [];
trace.units = [];
trace.freq_ech_units =[]';
trace.name=[num2str(shot) '/' ppftype '/' 'PGi and PGiN'];
end
%&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&
case {'powers'}
% load powers from TOT timebase as well
%
shotfile_exp_eff = AUGexplocation{index};
trace.t=[];
trace.x=[];
ppftype='TOT';
% P ohmic
[a,e]=rdaAUG_eff(shot,ppftype,'P_OH',shotfile_exp_eff);
if isempty(a) || e~=0
trace_all = struct([]);
else
nb_timepoints = length(a.time_aug.value);
trace_all.powers = NaN*ones(nb_timepoints,5);
trace_all.powers(:,1) = reshape(a.data,nb_timepoints,1);
trace_all.power_names{1} = [ppftype '/P_OH'];
trace_all.pohmic = a;
trace_all.t = a.time_aug.value;
end
% P NBI
[a,e]=rdaAUG_eff(shot,ppftype,'PNBI_TOT',shotfile_exp_eff);
if isempty(a) || e~=0
else
trace_all.powers(:,2) = reshape(a.data,nb_timepoints,1);
trace_all.power_names{2} = [ppftype '/PNBI_TOT'];
trace_all.pnbi = a;
end
% P ECRH
[a,e]=rdaAUG_eff(shot,ppftype,'PECR_TOT',shotfile_exp_eff);
if isempty(a) || e~=0
else
trace_all.powers(:,3) = reshape(a.data,nb_timepoints,1);
trace_all.power_names{3} = [ppftype '/PECR_TOT'];
trace_all.pecrh = a;
end
% P ICRH
[a,e]=rdaAUG_eff(shot,ppftype,'PICR_TOT',shotfile_exp_eff);
if isempty(a) || e~=0
else
trace_all.powers(:,4) = reshape(a.data,nb_timepoints,1);
trace_all.power_names{4} = [ppftype '/PICR_TOT'];
trace_all.picrh = a;
end
trace_all.powers(:,5) = trace_all.powers(:,1) + trace_all.powers(:,2) + trace_all.powers(:,3) + trace_all.powers(:,4);
% tau_tot
[a,e]=rdaAUG_eff(shot,ppftype,'tau_tot',shotfile_exp_eff);
if isempty(a) || e~=0
else
trace_all.tau_tot = a;
ij=find(~isnan(trace_all.tau_tot.value));
trace_all.tau_tot_spline=interpos(trace_all.t(ij),trace_all.tau_tot.value(ij),trace_all.t,-1e3);
trace_all.tau_tot_spline = max(trace_all.tau_tot_spline,0.);
end
% betaN
[a,e]=rdaAUG_eff(shot,ppftype,'beta_N',shotfile_exp_eff);
if isempty(a) || e~=0
else
trace_all.betan = a;
end
% betaNthermal
[a,e]=rdaAUG_eff(shot,ppftype,'beta_Nth',shotfile_exp_eff);
if isempty(a) || e~=0
else
trace_all.betan_thermal = a;
end
%
if ~isempty(trace_all)
trace_all.dim=[{trace_all.t} {[1:5]}];
trace_all.data = trace_all.powers;
trace = trace_all;
trace.x=trace.dim{2};
trace.dimunits=[{'time [s]'} {'Pohmic Pnbi Pecrh Picrh Ptot'}];
trace.units='W';
trace.name=[num2str(shot) '/' ppftype '/' 'P_OH, PNBI_TOT, PECR_TOT, PICR_TOT'];
else
trace.data = [];
trace.dim = [];
trace.dimunits = [];
trace.x = [];
trace.t = [];
trace.units = [];
trace.freq_ech_units =[]';
trace.name=[num2str(shot) '/' ppftype '/' 'P_OH, PNBI_TOT, PECR_TOT, PICR_TOT'];
end
%&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&
case {'transp'} % for TRANSP
% LOAD TRA data into substructures
shotfile_exp_eff = AUGexplocation{index}
if nargin>=3 & ~isempty(varargin{1})
shotfile_exp_eff = varargin{1};
end
trace.t=[];
trace.x=[];
ppftype='TRA';
TRANSP_signals;
for i=1:size(transp_sig,1)
if strcmp(lower(transp_sig{i,2}),'signal') || strcmp(lower(transp_sig{i,2}),'signal-group')
try
eval(['[trace.' transp_sig{i,1} ',e]=rdaAUG_eff(shot,ppftype,''' transp_sig{i,1} ''',shotfile_exp_eff);']);
catch
eval(['trace.' transp_sig{i,1} '=[];']);
end
elseif strcmp(lower(transp_sig{i,2}),'area-base')
clear adata_area
try
[adata_area]=sf2ab(ppftype,shot,transp_sig{i,1},'-exp',shotfile_exp_eff);
catch
adata_area.value = cell(0);
end
eval(['trace.' transp_sig{i,1} '=adata_area;']);
elseif strcmp(lower(transp_sig{i,2}),'time-base')
clear adata_time
try
[adata_time]=sf2tb(ppftype,shot,transp_sig{i,1},'-exp',shotfile_exp_eff);
catch
adata_time.value = cell(0);
end
eval(['trace.' transp_sig{i,1} '=adata_time;']);
end
end
otherwise
disp('case not yet defined')
end