Skip to content
Snippets Groups Projects
Commit 7b31b6ec authored by Olivier Sauter's avatar Olivier Sauter
Browse files

Merge branch 'add_IMAS' into 'master' 

add machine='imas' and getting 'ids','source','ids_name' from IMAS database

See merge request spc/tcv/tbx/gdat!36
parents 7ae52b2c a5ea0a15
No related branches found
No related tags found
1 merge request!36add machine='imas' and getting 'ids','source','ids_name' from IMAS database
Checking pipeline status
Hide whitespace changes
Inline Side-by-side
matlab/IMAS/gdat_imas.m 0 → 100644
+ 716
− 0
View file @ 7b31b6ec
function [gdat_data,gdat_params,error_status,varargout] = gdat_imas(shot,data_request,varargin)
%
% function [gdat_data,gdat_params,error_status,varargout] = gdat(shot,data_request,varargin)
%
% Aim: get data from a given machine using full path or keywords.
% data_request are and should be case independent (transformed in lower case in the function and outputs)
%
% If no inputs are provided, return the list of available pre-defined data_request in gdat_data and default parameters gdat_params
%
% IMAS means assume access data from IDS database ala IMAS with tokamakname and dataversion and run_number using imas_open_env and ids_get calls
% Should be the same on ITER hpc and gateway, to be checked
%
% At first only "ids" request implemented, but keywords like ip etc should relate to specific IDS as for any machine
%
% Inputs:
%
% no inputs: return default parameters in a structure form in gdat_params
% shot: shot number (if empty uses 'idx_imas_open_env' in params from previous call or input param)
% data_request: keyword (like 'ip') or trace name or structure containing all parameters but shot number
% varargin{i},varargin{i+1},i=1:nargin-2: additional optional parameters given in pairs: param_name, param_value
% The optional parameters list might depend on the data_request
% examples of optional parameters:
% 'plot',1 (plot is set by default to 0)
% 'machine','IMAS' (the default machine is the local machine)
%
%
% Outputs:
%
% gdat_data: structure containing the data, the time trace if known and other useful information
% gdat_data.t : time trace
% gdat_data.data: requested data values
% gdat_data.dim : values of the various coordinates related to the dimensions of .data(:,:,...)
% note that one of the dim is the time, replicated in .t for clarity
% gdat_data.dimunits : units of the various dimensions, 'dimensionless' if dimensionless
% gdat_data.error_bar : if provided
% gdat_data.gdat_call : list of parameters provided in the gdat call (so can be reproduced)
% gdat_data.shot: shot number
% gdat_data.machine: machine providing the data
% gdat_data.gdat_request: keyword for gdat if relevant
% gdat_data.data_fullpath: full path to the data node if known and relevant, or relevant expression called if relevant
% gdat_data.gdat_params: copy gdat_params for completeness (gdat_params contains a .help structure detailing each parameter)
% gdat_data.xxx: any other relevant information
%
% gdat_params contains the options relevant for the called data_request. It also contains a help structure for each option
% eg.: param1 in gdat_params.param1 and help in gdat_params.help.param1
%
% Note for liuqe parameter: At this time we have liuqe_fortran (which used to be the default) and liuqe_matlab (which is new and now becomes/is the default)
% We still have liuqe1, 2 and 3 nodes for either of them, so you can choose:
% 'liuqe',1 (2 or 3): to get the default liuqe 1, 2 or 3 (which is now the matlab version)
% 'liuqe',11 (12 or 13): to get liuqe_fortran nodes 1, 2 or 3
% 'liuqe',21 (22 or 23): to get liuqe_matlab nodes 1, 2 or 3 (may be needed if we set the default to liuqe_fortran for old shots)
% 'liuqe',-1 : to get FBTE result
%
% Examples:
% (should add working examples for various machines (provides working shot numbers for each machine...))
%
% to get the list of ids: aa=gdat([],'ids','machine','imas'); aa.gdat_params.sources_available
% to get an empty ids: aa=gdat([],'ids','machine','imas','source','core_profiles'); aa.core_profiles
%
% [a1,a2]=gdat;
% a2.data_request = 'ids'; a2.source='equilibrium'; a2.machine = 'imas';
% a3=gdat(48836,a2); % gives input parameters as a structure, allows to call the same for many shots
% a4=gdat('opt1',123,'opt2',[1 2 3],'shot',48832,'data_request','Ip','opt3','aAdB'); % all in pairs
% a5=gdat(48836,'ids','source','core_profiles'); % standard call for a given IDS core_profiles, use {'',''} for a list of IDSs
% a6=gdat(48836,'ip'); % standard call for a specific keyword
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
% Comments for local developer:
% gdat is just a "header routine" calling the gdat for the specific machine gdat_`machine`.m which can be called
% directly, thus which should be able to treat the same type of input arguments
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Prepare some variables, etc
varargout{1}=cell(1,1);
error_status=1;
% construct main default parameters structure
gdat_params.data_request = '';
default_machine = 'imas';
gdat_params.machine=default_machine;
gdat_params.doplot = 0;
gdat_params.nverbose = 1;
gdat_params.run = []; % should be defined by user to make sure
gdat_params.occurence = 0; % default/usual value
% special case need to treat here if params given as structure in data_request:
if isstruct(data_request)
gdat_params = data_request;
end
% construct list of keywords from global set of keywords and specific IMAS set
% get data_request names from centralized table
%%% data_request_names = get_data_request_names; % do not use xlsx anymore but scan respective machine_requests_mapping.m files
data_request_names = get_data_request_names_from_gdat_xxx(default_machine);
% add IMAS specific to all:
if ~isempty(data_request_names.imas)
imas_names = fieldnames(data_request_names.imas);
for i=1:length(imas_names)
data_request_names.all.(imas_names{i}) = data_request_names.imas.(imas_names{i});
end
end
data_request_names_all = sort(fieldnames(data_request_names.all));
% construct default output structure
gdat_data.data = [];
gdat_data.units = [];
gdat_data.dim = [];
gdat_data.dimunits = [];
gdat_data.t = [];
gdat_data.x = [];
gdat_data.shot = [];
gdat_data.gdat_request = [];
gdat_data.gdat_params = gdat_params;
gdat_data.data_fullpath = [];
% Treat inputs:
ivarargin_first_char = 3;
data_request_eff = '';
if nargin>=2
if ischar(data_request);
data_request_eff = data_request;
data_request = data_request; % should not lower until see if keyword
elseif isstruct(data_request)
data_request_eff = data_request.data_request;
end
end
gdat_data.gdat_request = data_request_names_all; % so if return early gives list of possible request names
gdat_data.gdat_params.help = imas_help_parameters(fieldnames(gdat_data.gdat_params));
gdat_params = gdat_data.gdat_params;
% no inputs
if nargin==0
% return defaults and list of keywords
return
end
do_mdsopen_mdsclose = 1;
% treat 1st arg
if nargin>=1
if isempty(shot)
% if shot already open, will use idx from gdat_params, otherwise returns
do_mdsopen_mdsclose = 0;
if ~isfield(gdat_data.gdat_params,'idx_imas_open_env')
if ~strcmp(data_request_eff,'ids'); return; end % empty shot return empty ids so valid command
end
elseif isnumeric(shot)
gdat_data.shot = shot;
elseif ischar(shot)
ivarargin_first_char = 1;
else
if gdat_params.nverbose>=1; warning('type of 1st argument unexpected, should be numeric or char'); end
error_status=2;
return
end
if nargin==1
% Only shot number given. If there is a default data_request set it and continue, otherwise return
return
end
end
% 2nd input argument if not part of pairs
if nargin>=2 && ivarargin_first_char~=1
if isempty(data_request)
return
end
% 2nd arg can be a structure with all options except shot_number, or a string for the pathname or keyword, or the start of pairs string/value for the parameters
if isstruct(data_request)
if ~isfield(data_request,'data_request')
if gdat_params.nverbose>=1; warning('expects field data_request in input parameters structure'); end
error_status=3;
return
end
data_request_eff = data_request.data_request;
% some request are case sensitive because mds names can be like \magnetics::ipol[*,"OH_001"]
data_request.data_request = data_request.data_request;
gdat_params = data_request;
else
% since data_request is char check from nb of args if it is data_request or a start of pairs
if mod(nargin-1,2)==0
ivarargin_first_char = 2;
else
ivarargin_first_char = 3;
data_request_eff = data_request;
end
end
end
if ~isstruct(data_request)
gdat_params.data_request = data_request_eff;
end
% if start pairs from shot or data_request, shift varargin
if ivarargin_first_char==1
varargin_eff{1} = shot;
varargin_eff{2} = data_request;
varargin_eff(3:nargin) = varargin(:);
elseif ivarargin_first_char==2
varargin_eff{1} = data_request;
varargin_eff(2:nargin-1) = varargin(:);
else
varargin_eff(1:nargin-2) = varargin(:);
end
% extract parameters from pairs of varargin:
if (nargin>=ivarargin_first_char)
if mod(nargin-ivarargin_first_char+1,2)==0
for i=1:2:nargin-ivarargin_first_char+1
if ischar(varargin_eff{i})
gdat_params.(lower(varargin_eff{i})) = varargin_eff{i+1}; % cannot enforce lower(params value) if filename or else are case sensitive
else
if gdat_params.nverbose>=1; warning(['input argument nb: ' num2str(i) ' is incorrect, expects a character string']); end
error_status=401;
return
end
end
else
if gdat_params.nverbose>=1; warning('number of input arguments incorrect, cannot make pairs of parameters'); end
error_status=402;
return
end
end
% if it is a request_keyword copy it:
if ischar(data_request_eff) || length(data_request_eff)==1
ij=strmatch(lower(data_request_eff),data_request_names_all,'exact');
else
ij=[];
end
if ~isempty(ij);
gdat_data.gdat_request = data_request_names_all{ij};
gdat_params.data_request = gdat_data.gdat_request;
if isfield(data_request_names.all.(data_request_names_all{ij}),'description') && ~isempty(data_request_names.all.(data_request_names_all{ij}).description)
% copy description of keyword
gdat_data.request_description = data_request_names.all.(data_request_names_all{ij}).description;
end
end
% special treatment if shot and data_request given within pairs
if isfield(gdat_params,'shot')
shot = gdat_params.shot; % should use only gdat_params.shot but change shot to make sure
gdat_data.shot = gdat_params.shot;
gdat_params=rmfield(gdat_params,'shot');
end
if ~isfield(gdat_params,'data_request') || isempty(gdat_params.data_request)
% warning('input for ''data_request'' missing from input arguments') % might be correct, asking for list of requests
error_status=5;
return
end
gdat_data.gdat_params = gdat_params; % This means that before here it is gdat_params which should be updated
clear gdat_params; % to make sure one uses gdat_data.gdat_params now
% re-assign main variables to make sure use the one in gdat_data structure
shot = gdat_data.shot;
data_request_eff = gdat_data.gdat_params.data_request;
if ~isfield(gdat_data.gdat_params,'database_user')
gdat_data.gdat_params.database_user = getenv('USER');
if gdat_data.gdat_params.nverbose >= 3; disp('set database_user to username as default since not provided'); end
end
if ~isfield(gdat_data.gdat_params,'tokamak_name')
gdat_data.gdat_params.tokamak_name = gdat_data.gdat_params.machine;
if gdat_data.gdat_params.nverbose >= 3; disp('set tokamak_name to machine as default since not provided'); end
end
if ~isfield(gdat_data.gdat_params,'data_major_version')
gdat_data.gdat_params.data_major_version = '3';
if gdat_data.gdat_params.nverbose >= 3; disp('set data_major_version to ''3'' as default since not provided'); end
end
error_status = 6; % at least reached this level
gdat_params = gdat_data.gdat_params;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Specifications on how to get the data provided in imas_requests_mapping
mapping_for_imas = imas_requests_mapping(data_request_eff,shot);
gdat_data.label = mapping_for_imas.label;
ishot=NaN;
if do_mdsopen_mdsclose
if ~isfield(gdat_data.gdat_params,'run') || isempty(gdat_data.gdat_params.run)
warning('run number not provided, cannot continue')
return
end
try
gdat_data.gdat_params.idx_imas_open_env = imas_open_env('ids', shot, gdat_data.gdat_params.run,gdat_data.gdat_params.database_user, ...
gdat_data.gdat_params.tokamak_name,gdat_data.gdat_params.data_major_version);
disp('')
catch ME
warning('could not imas_open_env with following gdat_data.gdat_params:');
gdat_data.gdat_params
rethrow(ME)
end
if gdat_data.gdat_params.idx_imas_open_env < 0
if gdat_data.gdat_params.nverbose>=1
warning(['idx < 0, shot/run not opened properly with gdat_data.gdat_params:']);
gdat_data.gdat_params
end
return
end
else
if isfield(gdat_data.gdat_params,'idx_imas_open_env')
if gdat_data.gdat_params.idx_imas_open_env < 0
if gdat_data.gdat_params.nverbose>=1
warning(['idx < 0, shot/run not opened properly with gdat_data.gdat_params:']);
gdat_data.gdat_params
end
return
end
else
if ~isempty(shot)
if gdat_data.gdat_params.nverbose>=1
warning(['no idx in gdat_params, cannot open shot with gdat_data.gdat_params:']);
gdat_data.gdat_params
end
return
end
end
end
gdat_params = gdat_data.gdat_params;
% fill again at end to have full case, but here to have present status in case of early return
gdat_data.gdat_params.help = imas_help_parameters(fieldnames(gdat_data.gdat_params));
gdat_data.mapping_for.imas = mapping_for_imas;
gdat_params = gdat_data.gdat_params;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
if strcmp(mapping_for_imas.method,'expression')
% method="expression"
% assume expression contains an expression to evaluate and which creates a local structure into variable gdat_tmp
% we copy the structure, to make sure default nodes are defined and to avoid if return is an closed object like tdi
% eval_expr = [mapping_for_imas.expression ';'];
eval([mapping_for_imas.expression ';']);
if isempty(gdat_tmp) || (~isstruct(gdat_tmp) & ~isobject(gdat_tmp))
if (gdat_params.nverbose>=1); warning(['expression does not create a gdat_tmp structure: ' mapping_for_imas.expression]); end
error_status=801;
return
end
tmp_fieldnames = setdiff(fieldnames(gdat_tmp),{'gdat_request','label'}); % could/should also remove label in any case
if sum(strcmp(tmp_fieldnames,'data'))==0 % note: cannot do isfield since gdat_tmp might be an object
if (gdat_params.nverbose>=1); warning(['expression does not return a child name ''data'' for ' data_request_eff]); end
end
for i=1:length(tmp_fieldnames)
gdat_data.(tmp_fieldnames{i}) = gdat_tmp.(tmp_fieldnames{i});
end
% add .t and .x in case only dim is provided
% do not allow shifting of timedim since should be treated in the relevant expression
ijdim=find(strcmp(tmp_fieldnames,'dim')==1);
if ~isempty(ijdim)
nbdims = length(gdat_data.dim);
if mapping_for_imas.timedim==-1;
mapping_for_imas.timedim = nbdims;
if (size(gdat_data.data,nbdims)==1 && nbdims>1); mapping_for_imas.timedim = nbdims-1; end
end
dim_nontim = setdiff([1:nbdims],mapping_for_imas.timedim);
ijt=find(strcmp(tmp_fieldnames,'t')==1);
if isempty(ijt)
gdat_data.t = gdat_data.dim{mapping_for_imas.timedim};
end
ijx=find(strcmp(tmp_fieldnames,'x')==1);
if isempty(ijx)
if ~isempty(dim_nontim)
% since most cases have at most 2d, copy as array if data is 2D and as cell if 3D or more
if length(dim_nontim)==1
gdat_data.x = gdat_data.dim{dim_nontim(1)};
else
gdat_data.x = gdat_data.dim(dim_nontim);
end
end
end
gdat_data.data_fullpath=mapping_for_imas.expression;
if isfield(gdat_tmp,'help')
gdat_data.help = gdat_tmp.help;
else
gdat_data.help = [];
end
if isfield(gdat_data.gdat_params,'time_out') && ~isempty(gdat_data.gdat_params.time_out)
[gdat_data] = gdat2time_out(gdat_data,1);
end
end
% end of method "expression"
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
elseif strcmp(mapping_for_imas.method,'switchcase')
switch data_request_eff % not lower(...) since data_request_eff should be lower case already at this stage
% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% First the request names valid for "all" machines:
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
case {'ids'}
ids_empty_path = fullfile(fileparts(mfilename('fullpath')),'..','IMAS_IMAS','ids_empty');
params_eff = gdat_data.gdat_params;
if isfield(params_eff,'source') && ~isempty(params_eff.source)
ids_top_name = params_eff.source;
else
ids_top_name = [];
warning('gdat:EmptyIDSName','Need an ids name in ''source'' parameter\n check substructure gdat_params.sources_available for an ids list');
addpath(ids_empty_path);
assert(~~exist('ids_list','file'),'could not find ids_list.m in %s',ids_empty_path);
gdat_data.gdat_params.sources_available = ids_list;
rmpath(ids_empty_path);
return
end
ids_gen_ok = ~~exist('ids_gen','file');
if ids_gen_ok
ids_empty = ids_gen(ids_top_name); % generate ids from gateway function ids_gen
else
% load empty ids structure from template file
fname = sprintf('ids_empty_%s',ids_top_name);
try
assert(~~exist(ids_empty_path,'dir'),'folder %s not found',ids_empty_path);
addpath(ids_empty_path);
assert(~~exist(fname,'file'),'file %s does not exist in %s',fname,ids_empty_path);
ids_empty = eval(fname);
disp(['use structure in .mat file: ' ids_empty_path '/' fname]);
rmpath(ids_empty_path);
catch ME
fprintf('Could not load empty template for %s\n',ids_top_name);
fprintf('Available templates:\n');
disp(dir(ids_empty_path));
rmpath(ids_empty_path);
rethrow(ME);
end
end
if ~isfield(gdat_data.gdat_params,'error_bar') || isempty(gdat_data.gdat_params.error_bar)
gdat_data.gdat_params.error_bar = 'delta';
end
if ~isfield(gdat_data.gdat_params,'cocos_in') || isempty(gdat_data.gdat_params.cocos_in)
gdat_data.gdat_params.cocos_in = 11;
end
if ~isfield(gdat_data.gdat_params,'cocos_out') || isempty(gdat_data.gdat_params.cocos_out)
gdat_data.gdat_params.cocos_out = 11;
end
if ~isfield(gdat_data.gdat_params,'ipsign_out') || isempty(gdat_data.gdat_params.ipsign_out)
gdat_data.gdat_params.ipsign_out = 0;
end
if ~isfield(gdat_data.gdat_params,'b0sign_out') || isempty(gdat_data.gdat_params.b0sign_out)
gdat_data.gdat_params.b0sign_out = 0;
end
if ~isfield(gdat_data.gdat_params,'ipsign_in') || isempty(gdat_data.gdat_params.ipsign_in)
if gdat_data.gdat_params.ipsign_out~=0
gdat_data.gdat_params.ipsign_in = mdsvalue('\pcs::data:iohfb');
else
gdat_data.gdat_params.ipsign_in = 0;
end
end
if ~isfield(gdat_data.gdat_params,'b0sign_in') || isempty(gdat_data.gdat_params.b0sign_in)
if gdat_data.gdat_params.b0sign_out~=0
gdat_data.gdat_params.b0sign_in = mdsvalue('\pcs::data:if36fb');
else
gdat_data.gdat_params.b0sign_in = 0;
end
end
try
if isfield(gdat_data.gdat_params,'idx_imas_open_env') && gdat_data.gdat_params.idx_imas_open_env >= 0
if gdat_data.gdat_params.occurence == 0
ids_top = ids_get(gdat_data.gdat_params.idx_imas_open_env,ids_top_name);
else
ids_top = ids_get(gdat_data.gdat_params.idx_imas_open_env,[ids_top_name '/' num2str(gdat_data.gdat_params.occurence)]);
end
gdat_data.(ids_top_name) = ids_top;
else
gdat_data.(ids_top_name) = ids_empty;
end
catch ME_imas_ids_get
disp(['there is a problem with: imas_get_ids_' ids_top_name ...
' , may be check if it exists in your path or test it by itself'])
gdat_data.(ids_top_name) = ids_empty;
gdat_data.([ids_top_name '_help']) = getReport(ME_imas_ids_get);
rethrow(ME_imas_ids_get)
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
otherwise
if (gdat_params.nverbose>=1); warning(['switchcase= ' data_request_eff ' not known in gdat_imas']); end
error_status=901;
return
end
else
if (gdat_params.nverbose>=1); warning(['IMAS method=' mapping_for_imas.method ' not known yet, contact Olivier.Sauter@epfl.ch']); end
error_status=602;
return
end
gdat_data.gdat_params.help = imas_help_parameters(fieldnames(gdat_data.gdat_params));
gdat_data.mapping_for.imas = mapping_for_imas;
gdat_params = gdat_data.gdat_params;
error_status=0;
return
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function [firthomratio] = get_fir_thom_rat_data(shot,maintracename,timebase,nverbose);
%
% since depends on shot number for using auto fit and thomson or thomson edge, use tracename and function here
%
% maintracename = 'thomson' or 'thomson_edge
%
% return fir_to_thomson ratio on time=timebase
%
% normally should use "main" thomson one built from auto fit if possible
% take edge one if need be
%
% default: maintracename = "thomson"
%
maintracename_eff = 'thomson';
if exist('maintracename') && ~isempty(maintracename)
maintracename_eff = maintracename;
end
firthomratio = NaN;
if ~exist('timebase') || isempty(timebase)
if (nverbose>=1)
disp('need a timebase in get_fir_thom_rat_data')
end
return
end
firthomratio = NaN*ones(size(timebase));
if strcmp(maintracename_eff,'thomson')
if shot>=23801
tracefirrat=tdi('\results::thomson.profiles.auto:fir_thom_rat'); %time base not same!!
if isempty(tracefirrat.data) || ischar(tracefirrat.data)
if (nverbose>=1); disp('problem with \results::thomson.profiles.auto:fir_thom_rat: empty, use thomson:fir_thom_rat'); end
tracefirrat=tdi('\results::thomson:fir_thom_rat');
end
else
tracefirrat=tdi('\results::thomson:fir_thom_rat');
if isempty(tracefirrat.data) || ischar(tracefirrat.data)
if (nverbose>=1); disp('problem with \results::thomson:fir_thom_rat: empty'); end
end
end
elseif strcmp(maintracename_eff,'thomson_edge')
if shot>=23801
tracefirrat=tdi('\results::thomson.profiles.auto:fir_thom_rat'); %time base not same!!
if isempty(tracefirrat.data) || ischar(tracefirrat.data)
if (nverbose>=1); disp('problem with \results::thomson.profiles.auto:fir_thom_rat: empty, use thomson:fir_thom_rat'); end
tracefirrat=tdi('\results::thomson:fir_thom_rat');
end
else
tracefirrat=tdi('\results::thomson_edge:fir_thom_rat');
if isempty(tracefirrat.data) || ischar(tracefirrat.data)
if (nverbose>=1); disp('problem with \results::thomson_edge:fir_thom_rat: empty'); end
end
end
else
if (nverbose>=1); disp('bad input in get_fir_thom_rat_data'); end
return
end
if ~isempty(tracefirrat.data) && ~ischar(tracefirrat.data) && any(isfinite(tracefirrat.data)) ...
&& ~isempty(tracefirrat.dim) && ~isempty(tracefirrat.dim{1})
firthomratio = interp1(tracefirrat.dim{1},tracefirrat.data,timebase);
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function [gdat_data] = get_thomson_raw_data(shot,data_request_eff,gdat_data,doedge,nverbose);
%
try
time=mdsdata('\results::thomson:times');
catch
gdat_data.error_bar = [];
if strcmp(data_request_eff(1:2),'ne')
tracefirrat_data = [];
gdat_data.firrat=tracefirrat_data;
gdat_data.data_raw = gdat_data.data;
gdat_data.error_bar_raw = gdat_data.error_bar;
end
if (nverbose>=1) && shot<100000
warning('Problems with \results::thomson:times')
warning(['!!!!!!!!!!!!!!!!!!!!!!!!! cannot continue with ' data_request_eff])
end
return
end
if isempty(time) || ischar(time)
thomsontimes=time;
if (nverbose>=1) && shot<100000
warning('!!!!!!!!!!!!!!!!!!!!!!!!!\results::thomson:times is empty? Check')
disp(['!!!!!!!!!!!!!!!!!!!!!!!!! cannot continue with ' data_request_eff])
end
return
end
edge_str_ = '';
edge_str_dot = '';
if doedge
edge_str_ = '_edge';
edge_str_dot = '.edge';
end
nodenameeff = ['\results::thomson' edge_str_dot ':' data_request_eff(1:2)];
tracetdi=tdi(nodenameeff);
if isempty(tracetdi.data) || ischar(tracetdi.data) || isempty(tracetdi.dim)
gdat_data.error_bar = [];
gdat_data.firrat = [];
gdat_data.data_raw = [];
gdat_data.error_bar_raw = [];
return
end
gdat_data.data=tracetdi.data'; % Thomson data as (t,z)
tracestd=tdi(['\results::thomson' edge_str_dot ':' data_request_eff(1:2) ':error_bar']);
gdat_data.error_bar=tracestd.data';
gdat_data.data_fullpath=[nodenameeff '; error_bar'];
% add fir if ne requested
if strcmp(data_request_eff(1:2),'ne')
tracefirrat_data = get_fir_thom_rat_data(shot,['thomson' edge_str_],time,nverbose);
gdat_data.firrat=tracefirrat_data;
gdat_data.data_raw = gdat_data.data;
gdat_data.data = gdat_data.data_raw * diag(tracefirrat_data);
gdat_data.error_bar_raw = gdat_data.error_bar;
gdat_data.error_bar = gdat_data.error_bar_raw * diag(tracefirrat_data);
gdat_data.data_fullpath=[gdat_data.data_fullpath '; fir_thom_rat ; _raw without firrat'];
ij=find(isfinite(tracefirrat_data));
if isempty(ij)
gdat_data.data_fullpath=[gdat_data.data_fullpath ' WARNING use ne_raw since firrat has NaNs thus ne=NaNs; fir_thom_rat ; _raw without firrat'];
disp('***********************************************************************')
disp('WARNING: ne from Thomson has fir_thom_rat with Nans so only ne_raw can be used');
disp('***********************************************************************')
end
end
z=mdsdata(['\diagz::thomson_set_up' edge_str_dot ':vertical_pos']);
gdat_data.dim=[{z};{time}];
gdat_data.dimunits=[{'Z [m]'} ; {'time [s]'}];
gdat_data.x=z;
gdat_data.t=time;
% isfield does not work since tracetdi is not a 'struct' but a tdi object, thus use fieldnames
if any(strcmp(fieldnames(tracetdi),'units'))
gdat_data.units=tracetdi.units;
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
function [psiscatvol,psi_max] = psi_scatvol_zshift(gdat_data, zshift, system_str_dot, psitbx_str, gdat_params )
% Input arguments:
% gdat_data: current structure containing TS data
% zshift: vertical shift of the equilibrium, can be a vector for time-dependent data
% system_str_dot: '' or '.edge' indicating thomson syste,
% psitbx_str: source for the equilibrium map, psitbx-style
% gdat_params: parameters for the gdat call % TODO: Why is it different from gdat_data.gdat_params?
% Use psitbx
t_th = gdat_data.t;
th_z = mdsdata(['dim_of(\thomson' system_str_dot ':te,1)']); % TODO: Isn't this gdat_data.dim{1}?
th_r = 0.9*ones(size(th_z));
ntt = numel(t_th);
nch = numel(th_z); % number of chords
[t_psi,status] = mdsvalue('dim_of(imas_eq("i_p",$))',psitbx_str); % TODO: replace with time_psi when imas_eq supports it
if ~rem(status,2)
warning('problem with getting equilibrium time base in gdat_imas')
display(t_psi);
return
end
zshifteff=zshift;
% set zshifteff time array same as t_psi
switch numel(zshifteff)
case 1
zshifteff=zshifteff * ones(size(t_th));
case numel(t_psi)
zshifteff=interp1(t_psi,zshifteff,t_th);
case numel(t_th)
% ok
otherwise
if (gdat_params.nverbose>=1)
disp(' bad time dimension for zshift')
disp(['it should be 1 or ' num2str(numel(t_th)) ' or ' num2str(numel(t_psi))])
return
end
end
% Get LIUQE times corresponding to TS times
t_eq = t_psi(iround(t_psi,t_th));
% Get PSI map
psi = psitbximas(gdat_data.shot,t_eq,'*0',psitbx_str);
% PSITBXIMAS will have removed duplicate times
i_psi = iround(psi.t,t_eq); % Mapping from psi times to TS times
% Define grid points where to evaluate PSI
if ~any(diff(zshifteff))
% Constant vector
pgrid = psitbxgrid('cylindrical','points',{th_r,th_z - zshifteff(1),0});
else
% We need as many time points in psi as TS times
psi = subs_time(psi,i_psi); % Now has ntt times
i_psi = 1:ntt; % Trivial mapping from psi times to TS times
% Time-varying vector
th_z_eff = repmat(th_z(:),1,ntt) - repmat(reshape(zshifteff,1,ntt),nch,1);
th_r_eff = repmat(th_r(:),1,ntt);
pgrid = psitbxgrid('cylindrical','time-points',{th_r_eff,th_z_eff,0});
end
% Compute psi at TS positions
psi_ts = psitbxf2f(psi,pgrid);
psiscatvol.data = squeeze(psi_ts.x(:,i_psi));
psiscatvol.dim{1} = gdat_data.x;
psiscatvol.dim{2} = t_th;
% NOTE: we should probably not include time points where equilibrium time is far from TS time.
% Compute psi_axis at TS times
% Do not use psi_axis node because:
% - For legacy LIUQE, psi_axis is not computed by psitbx, so we could get a
% different answer and cause complex numbers computing rho_psi
% - For MATLAB LIUQE, psi_axis is only the result of asxymex whereas the
% psitbx adds some Newton iterations so again complex numbers are possible
psi_norm = psitbxp2p(psi,'01');
psi_max.data = psi_norm.psimag(i_psi);
psi_max.dim = {t_th};
matlab/IMAS/imas_help_parameters.m 0 → 100644
+ 60
− 0
View file @ 7b31b6ec
function help_struct = imas_help_parameters(parameter_list)
%
% retrieve from present table the relevant help lines for the parameter_list{:}
% should come from sqlite database at some point...
%
% return the whole help structure if parameter_list empty or not provided
%
% do:
% help_struct = imas_help_parameters(fieldnames(gdat_data.gdat_params));
%
% to get relevant help description
%
% Defaults
help_struct_all = struct(...
'data_request', ['automatically filled in by gdat, name of request used in gdat call.' char(10) ...
'contains current list of keywords if gdat called with no arguments: aa=gdat;' char(10) ...
'Note shot value should not be in params so params can be used to load same data from another shot'] ...
,'machine', 'machine name like ''TCV'', ''AUG'', ''IMAS'' case insensitive' ...
,'doplot', '0 (default), if 1 calls gdat_plot for a new figure, -1 plot over current figure with hold all, see gdat_plot for details' ...
,'nverbose','1 (default) displays warnings, 0: only errors, >=3: displays all extra information' ...
);
% IMAS related
help_struct_all.time = 'eqdsk: time(s) value(s) requested, by default time=1.0s (see time_out for other requests)';
help_struct_all.time_out = ['requested time for output: data points within interval if time_out=[t1 t2], otherwise assumes series of points, uses linear interpolation in that case (default [-Inf Inf])'...
char(10) 'for sxr, mpx: only time interval provided in time_out is relevant'];
help_struct_all.cocos_in_out = ['cocos_in or cocos_out value input or desired in output, uses eqdsk_cocos_transform. Note should use latter if a specific Ip and/or B0 sign' ...
'is wanted. See O. Sauter et al Comput. Phys. Commun. 184 (2013) 293'];
help_struct_all.source = sprintf('%s\n',...
'ids_names for request ''ids'' like magnetics, equilibrium, etc');
help_struct_all.error_bar = sprintf('%s\n','for ids: choice of nodes fill in and how:', ...
'''delta'' (default): only upper fill in with the abs(value) to add or subtract to data to get upper and lower values (symmetric)', ...
'''delta_with_lower'': same as delta but fill in lower node as well (with delta as well, same as upper)', ...
'''added'': add the delta values (old cpo style), so upper=data+error_bar and lower=data+error_bar');
help_struct_all.write = 'eqdsk: write eqdsk while loading data (1, default) or not (0)';
%help_struct_all. = '';
if ~exist('parameter_list') || isempty(parameter_list)
help_struct = help_struct_all;
return
end
if iscell(parameter_list)
parameter_list_eff = parameter_list;
elseif ischar(parameter_list)
% assume only one parameter provided as char string
parameter_list_eff{1} = parameter_list;
else
disp(['unknown type of parameter_list in imas_help_parameters.m'])
parameter_list
return
end
fieldnames_all = fieldnames(help_struct_all);
for i=1:length(parameter_list_eff)
if any(strcmp(fieldnames_all,parameter_list_eff{i}))
help_struct.(parameter_list_eff{i}) = help_struct_all.(parameter_list_eff{i});
end
end
matlab/IMAS/imas_requests_mapping.m 0 → 100644
+ 43
− 0
View file @ 7b31b6ec
function mapping = imas_requests_mapping(data_request,shot)
%
% Information pre-defined for gdat_imas, you can add cases here to match official cases in gdat_imas, allowing backward compatibility
%
% give the shot number in case data origin depends on the shot number, allows to adapt easily
%
% Defaults
mapping = struct(...
'label', '', ...
'method', '', ...
'expression','', ...
'timedim', -1, ... % dim which is the time is the database, to copy in .t, the other dims are in .x (-1 means last dimension)
'gdat_timedim',[], ... % if need to reshape data and dim orders to have timedim as gdat_timedim (shifting time to gdat_timedim)
'min', -inf, ...
'max', inf);
% Note that gdat_timedim is set to timedim at end of this function if empty
% gdat_timedim should have effective index of the time dimension in gdat
if ~exist('data_request') || isempty(data_request)
return
end
% default label: data_request keyword itself
mapping.label = data_request;
%
% label is used for plotting
switch lower(data_request)
case 'ids'
mapping.timedim = 1;
mapping.label = 'ids ala imas';
mapping.method = 'switchcase';
otherwise
mapping.label = data_request;
mapping.method = 'switchcase';
mapping.expression = data_request;
end
if isempty(mapping.gdat_timedim)
mapping.gdat_timedim = mapping.timedim;
end
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment