-
Olivier Sauter authoredgit-svn-id: https://spcsvn.epfl.ch/repos/TCV/gdat/trunk@4385 d63d8f72-b253-0410-a779-e742ad2e26cf
Olivier Sauter authoredgit-svn-id: https://spcsvn.epfl.ch/repos/TCV/gdat/trunk@4385 d63d8f72-b253-0410-a779-e742ad2e26cf
loadAUGdata.m 42.51 KiB
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)
%
% 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') %
%
% 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 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
% 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
% ece :
% eced :
% 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_eff,1)==1
i=findstr('/',data_type);
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
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,[{'ECE'}],'exact'))
data_type_eff_noext='ece';
end
if ~isempty(strmatch(upper(data_type_eff_noext),[{'ECED'}],'exact'))
data_type_eff_noext='eced';
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(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(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'} {'rcont'} {'zcont'} {'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'} {'sxr'} {'sxR'} {'sxb'} {'ece'} {'eced'} {'Halpha'} {'pgyro'} {'powers'}];
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.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.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.iece=strmatch('ece',AUGkeywrdall,'exact');
AUGsig.ieced=strmatch('eced',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.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.iece)=AUGkeywrdall(AUGsig.iece);
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)={''; ''};
AUGsiglocation(:,AUGsig.iq0)={'FPG'; 'q0'};
AUGsiglocation(:,AUGsig.iq95)={'FPG'; 'q95'};
AUGsiglocation(:,AUGsig.ikappa)={''; ''};
AUGsiglocation(:,AUGsig.ideltatop)={''; ''};
AUGsiglocation(:,AUGsig.ideltabot)={''; ''};
AUGsiglocation(:,AUGsig.ineint)={'DCN'; 'H-1'};
AUGsiglocation(:,AUGsig.iece)={'CEC'; 'Trad-A'};
AUGsiglocation(:,AUGsig.ieced)={'RMD'; 'Trad-A'}; % ECED
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'};
% 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;
elseif ~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')
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?')
trace.dim=[{trace.x} ; {trace.t} ; {d}];
trace.dimunits=[{'R [m] or rho=sqrt(psi_norm)'} ; {'time [s]'} ; {'d'}];
trace.d=d;
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'}
% 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);
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
tracename='J';
end
trace.t=[];
trace.x=[];
ppftype='SXB';
iok=0;
for ichord=starti:endi
tracename_eff = [tracename '_' num2str(ichord,'%.3d')];
[a,e]=rdaAUG_eff(shot,ppftype,tracename_eff,shotfile_exp_eff,timerange);
if isempty(a) || e~=0
trace_all = struct([]);
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 {'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 '/' tracename];
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', 'qrho', 'qrho_fpp'}
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';
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);
Lpf_tot = Lpf_par.value; % nb of points: 100000*nb_SOL points + nb_core
Lpf_SOL = fix(Lpf_tot(1:NTIME)/100000);
Lpf1_t = mod(Lpf_tot(1:NTIME),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);
[psi_tree,e]=rdaAUG_eff(shot,DIAG,'PFL',shotfile_exp);
[phi_tree,e]=rdaAUG_eff(shot,DIAG,'TFLx',shotfile_exp);
[Vol,e]=rdaAUG_eff(shot,DIAG,'Vol',shotfile_exp);
[Area,e]=rdaAUG_eff(shot,DIAG,'Area',shotfile_exp);
[Ri,e]=rdaAUG_eff(shot,DIAG,'Ri',shotfile_exp);
[Zj,e]=rdaAUG_eff(shot,DIAG,'Zj',shotfile_exp);
[PFM_tree,e]=rdaAUG_eff(shot,DIAG,'PFM',shotfile_exp);
[Pres,e]=rdaAUG_eff(shot,DIAG,'Pres',shotfile_exp);
[FFP,e]=rdaAUG_eff(shot,DIAG,'FFP',shotfile_exp); if isempty(FFP.value); FFP.value=NaN*ones(NTIME,max(Lpf1_t)); end
[LPFx,e]=rdaAUG_eff(shot,DIAG,'LPFx',shotfile_exp);
[PFxx,e]=rdaAUG_eff(shot,DIAG,'PFxx',shotfile_exp);
[RPFx,e]=rdaAUG_eff(shot,DIAG,'RPFx',shotfile_exp);
[zPFx,e]=rdaAUG_eff(shot,DIAG,'zPFx',shotfile_exp);
% seems "LCFS" q-value is far too large, limit to some max (when diverted)
max_qValue = 40; % 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
if qpsi.value(1,2)<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))));
% 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
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);
%
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';
%
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
a = sf2par('ECS',shot,'gyr_freq',['P_sy1_g' num2str(i)]);
if isempty(a)
else
trace_all.freq_ecrh1(i) = a;
trace_all.freq_ech(i) = a.value;
end
a = sf2par('ECS',shot,'GPolPos',['P_sy1_g' num2str(i)]);
if isempty(a)
else
trace_all.polpos_ecs(i) = a.value;
end
a = sf2par('ECS',shot,'GTorPos',['P_sy1_g' num2str(i)]);
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
a = sf2par('ECS',shot,'gyr_freq',['P_sy2_g' num2str(i)]);
if isempty(a)
else
trace_all.freq_ecrh2(i) = a;
trace_all.freq_ech(i+4) = a.value;
end
a = sf2par('ECS',shot,'GPolPos',['P_sy2_g' num2str(i)]);
if isempty(a)
else
trace_all.polpos_ecs(i+4) = a.value;
end
a = sf2par('ECS',shot,'GTorPos',['P_sy2_g' num2str(i)]);
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
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
otherwise
disp('case not yet defined')
end