diff --git a/crpptbx/JET/get_grids_1d.m b/crpptbx/JET/get_grids_1d.m
new file mode 100644
index 0000000000000000000000000000000000000000..b5cb83b5952f8d7efdd97d18836ef6e82ea72792
--- /dev/null
+++ b/crpptbx/JET/get_grids_1d.m
@@ -0,0 +1,71 @@
+function [gdat_data] = get_grids_1d(gdat_data,nbdim_x,nopt,nverbose);
+%
+% add various rhos in grids_1d
+%
+% assume gdat_data.x=rhopol, gdat_data.t, and gdat_data(x,t)
+%
+% nbdim_x = 1: same x(1:length(x)) (rhopol) for all times
+% = 2: x(:,t), rhopol depends on time like Thomson projection
+%
+% nopt = 0: do not fill in, just make the empty structure
+% = 1: do compute various fields
+%
+% compute psi, rhotor_edge, rhotornorm, volume_edge and rhovol
+%
+% use gdat calls to psi_axis, psi_edge, rhotor, etc with same basic parameters as data.gdat_params
+%
+
+gdat_data.grids_1d.rhopolnorm = gdat_data.x;
+if (nopt == 0) || isempty(gdat_data.x) || isempty(gdat_data.t) || isempty(gdat_data.data) || ischar(gdat_data.data)
+ gdat_data.grids_1d.rhotornorm = [];
+ gdat_data.grids_1d.rhovolnorm = [];
+ gdat_data.grids_1d.psi = [];
+ gdat_data.grids_1d.rhotor_edge = [];
+ gdat_data.grids_1d.volume_edge = [];
+ return
+end
+
+params_eff = gdat_data.gdat_params;params_eff.doplot=0;
+params_eff.data_request='rhotor';
+rhotor = gdat(gdat_data.shot,params_eff);
+params_eff.data_request='rhovol';
+rhovol = gdat(gdat_data.shot,params_eff);
+
+psi0 = interpos(rhotor.t',rhotor.psi_axis,gdat_data.t,-0.01);
+if (nbdim_x == 1)
+ gdat_data.grids_1d.psi = gdat_data.grids_1d.rhopolnorm.^2*reshape(psi0,1,length(psi0));
+elseif (nbdim_x == 2)
+ gdat_data.grids_1d.psi = gdat_data.grids_1d.rhopolnorm.^2.*repmat(reshape(psi0,1,length(psi0)),size(gdat_data.grids_1d.rhopolnorm,1),1);
+else
+ if nverbose>=0; disp(['option: nbdim_x = ' numstr(nbdim_x) ' not implemented, check with O. Sauter']); end
+ return
+end
+gdat_data.grids_1d.rhotornorm = NaN*ones(size(gdat_data.data));
+gdat_data.grids_1d.rhovolnorm = NaN*ones(size(gdat_data.data));
+it_rt = iround_os(rhotor.t,gdat_data.t);
+it_vol = iround_os(rhovol.t,gdat_data.t);
+for it=1:length(gdat_data.t)
+ % do an interpolation on closest point to avoid 2D interp
+ it_rt_eff = it_rt(it);
+ it_vol_eff = it_vol(it);
+ if (nbdim_x == 1)
+ ii=find(~isnan(gdat_data.grids_1d.rhopolnorm));
+ else
+ ii=find(~isnan(gdat_data.grids_1d.rhopolnorm(:,it)));
+ end
+ if (nbdim_x == 1)
+ if length(ii)==length(gdat_data.grids_1d.rhopolnorm)
+ gdat_data.grids_1d.rhotornorm(:,it)=interpos(rhotor.x,rhotor.data(:,it_rt_eff),gdat_data.grids_1d.rhopolnorm);
+ gdat_data.grids_1d.rhovolnorm(:,it)=interpos(rhovol.x,rhovol.data(:,it_vol_eff),gdat_data.grids_1d.rhopolnorm);
+ end
+ else
+ if length(ii)==size(gdat_data.grids_1d.rhopolnorm,1)
+ gdat_data.grids_1d.rhotornorm(:,it)=interpos(rhotor.x,rhotor.data(:,it_rt_eff),gdat_data.grids_1d.rhopolnorm(:,it));
+ gdat_data.grids_1d.rhovolnorm(:,it)=interpos(rhovol.x,rhovol.data(:,it_vol_eff),gdat_data.grids_1d.rhopolnorm(:,it));
+ end
+ end
+end
+gdat_data.grids_1d.rhotor_edge=interpos(rhotor.t',rhotor.rhotor_edge,gdat_data.t',-0.01);
+gdat_data.grids_1d.volume_edge=interpos(rhovol.t',rhovol.volume_edge,gdat_data.t',-0.01);
+
+
diff --git a/crpptbx/README b/crpptbx/README
index 39e8410eef2acdc41524da11006b1428d0cf874e..29112e538cc6cba508082978ba59d79dff5f34ff 100644
--- a/crpptbx/README
+++ b/crpptbx/README
@@ -28,3 +28,6 @@ since then using svntag
2017/07/01: tag 3_10 before changing default to LIUQE.M and having introduced liuqe=1,2,3,11,12,13,21,22,23 or -1 as in GUIprofs and proffit
svn cp https://spcsvn.epfl.ch/repos/TCV/GUIprofs/trunk https://spcsvn.epfl.ch/repos/TCV/GUIprofs/tags/GUIprofs_v3_10 -m"tag version 3_10, before changing to default liuqe to matlab liuqe (and changes in GUIprofs and proffit as well)"
+
+2017/07/01: tag 4_1 change default to LIUQE.M and introduced liuqe=1,2,3,11,12,13,21,22,23 or -1 as in GUIprofs and proffit
+ svn cp https://spcsvn.epfl.ch/repos/TCV/GUIprofs/trunk https://spcsvn.epfl.ch/repos/TCV/GUIprofs/tags/GUIprofs_v4_1 -m"tag version 4_1, changed default liuqe to matlab liuqe (and changes in proffit and gdat as well)"
diff --git a/crpptbx/TCV/gdat_tcv.m b/crpptbx/TCV/gdat_tcv.m
index c636926c3c20515102652e2d32b9a4d6fc7ce828..af483bcfbb5a882b19f932acfaf9b964dd17f13f 100644
--- a/crpptbx/TCV/gdat_tcv.m
+++ b/crpptbx/TCV/gdat_tcv.m
@@ -39,6 +39,13 @@ function [gdat_data,gdat_params,error_status,varargout] = gdat_tcv(shot,data_req
% 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...))
%
@@ -251,16 +258,25 @@ if isfield(gdat_data.gdat_params,'liuqe') && ~isempty(gdat_data.gdat_params.liuq
else
gdat_data.gdat_params.liuqe = liuqe_version;
end
+liuqe_matlab = 1; % now default should be matlab liuqe nodes
+if liuqe_version<0 || (liuqe_version > 10 && liuqe_version < 20)
+ liuqe_matlab = 0;
+end
+liuqe_version_eff = mod(liuqe_version,10);
substr_liuqe = '';
-if liuqe_version==2 || liuqe_version==3
- substr_liuqe = ['_' num2str(liuqe_version)];
+substr_liuqe_tcv_eq = '';
+if liuqe_version_eff==2 || liuqe_version_eff==3 || (liuqe_matlab==1 && liuqe_version_eff==1)
+ substr_liuqe = ['_' num2str(liuqe_version_eff)];
end
-
+if liuqe_version_eff==2 || liuqe_version_eff==3
+ substr_liuqe_tcv_eq = num2str(liuqe_version_eff);
+end
+
% special treatment for model shot=-1 or preparation shot >=100'000
begstr = '';
if shot==-1 || (shot>=100000 && shot < 200000) || liuqe_version==-1
% requires FBTE
- liuqe_version = -1;
+ liuqe_version_eff = -1;
begstr = 'tcv_eq( "';
substr_liuqe = '", "FBTE" )';
end
@@ -268,13 +284,14 @@ end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% Specifications on how to get the data provided in tcv_requests_mapping
-mapping_for_tcv = tcv_requests_mapping(data_request_eff);
+mapping_for_tcv = tcv_requests_mapping(data_request_eff,shot);
gdat_data.label = mapping_for_tcv.label;
ishot=NaN;
if do_mdsopen_mdsclose
% mdsdefaultserver tcv1.epfl.ch; % should be in tcv general path, but set-it in the meantime...
- if liuqe_version==-1
+%%% if liuqe_version_eff==-1
+ if shot==-1
ishot = mdsopen('pcs', shot);
else
if length(data_request_eff)>7 && strcmp(data_request_eff(1:6),'\rtc::')
@@ -301,11 +318,17 @@ gdat_params = gdat_data.gdat_params;
if strcmp(mapping_for_tcv.method(1:3),'tdi')
% need to treat liuqe2, model, etc from options....
substr_tdi = '';
- if strcmp(mapping_for_tcv.method,'tdiliuqe'); substr_tdi = substr_liuqe; end
+ if liuqe_matlab==0 && strcmp(mapping_for_tcv.method,'tdiliuqe'); substr_tdi = substr_liuqe; end
if iscell(mapping_for_tcv.expression)
if length(mapping_for_tcv.expression)>0
% series of arguments for tdi given in each cell
- eval_expr = ['tdi(''' mapping_for_tcv.expression{1} substr_tdi ''''];
+ if liuqe_matlab==1
+ ij = findstr(mapping_for_tcv.expression{1},'equil_');
+ if ~isempty(ij)
+ mapping_for_tcv.expression{1}(ij+5:ij+6) = substr_liuqe;
+ end
+ end
+ eval_expr = ['tdi(''' mapping_for_tcv.expression{1} ''''];
for i=2:length(mapping_for_tcv.expression)
eval_expr = [eval_expr ',''' mapping_for_tcv.expression{i} ''''];
end
@@ -317,14 +340,52 @@ if strcmp(mapping_for_tcv.method(1:3),'tdi')
return
end
else
- if liuqe_version==-1
+ do_liuqem2liuqef = 0;
+ if liuqe_version_eff==-1
mapping_for_tcv_expression_eff = mapping_for_tcv.expression;
- if strcmp(lower(mapping_for_tcv.expression(1:8)),'\results')
+ if length(mapping_for_tcv.expression)>8 && strcmp(lower(mapping_for_tcv.expression(1:8)),'\results')
mapping_for_tcv_expression_eff = mapping_for_tcv.expression(11:end);
+ elseif findstr('tcv_eq',lower(mapping_for_tcv.expression))
+ ij = regexpi(mapping_for_tcv.expression,'''''');
+ if length(ij)<2
+ disp(['expected more double quotes in mapping_for_tcv.expression = ' mapping_for_tcv.expression]);
+ return
+ else
+ mapping_for_tcv_expression_eff = mapping_for_tcv.expression(ij(1)+2:ij(2)-1);
+ end
end
eval_expr = ['tdi(''' begstr mapping_for_tcv_expression_eff substr_liuqe ''');']
else
- eval_expr = ['tdi(''' mapping_for_tcv.expression substr_tdi ''');'];
+ if liuqe_matlab==1
+ ij = findstr(mapping_for_tcv.expression,'equil_');
+ if ~isempty(ij)
+ mapping_for_tcv.expression(ij+5:ij+6) = substr_liuqe;
+ end
+ ij = regexpi(mapping_for_tcv.expression,'LIUQE.*\.M','once');
+ if ~isempty(ij)
+ mapping_for_tcv.expression = [mapping_for_tcv.expression(1:ij+6) substr_liuqe_tcv_eq mapping_for_tcv.expression(ij+7:end)];
+ end
+ else
+ ij = regexpi(mapping_for_tcv.expression,'LIUQE.*\.M','once');
+ if ~isempty(ij)
+ mapping_for_tcv.expression = [mapping_for_tcv.expression(1:ij+4) substr_liuqe_tcv_eq mapping_for_tcv.expression(ij+7:end)];
+ do_liuqem2liuqef = 1;
+ end
+ end
+ eval_expr = ['tdi(''' mapping_for_tcv.expression ''');']
+ end
+ % special cases for matching liuqe fortran and matlab
+ if liuqe_matlab == 0 && do_liuqem2liuqef==1
+ % assume tcv_eq(''keyword'',''LIUQE..'')
+ liuqe_fortran_matlab = liuqefortran2liuqematlab;
+ ij = regexpi(eval_expr,'''''');
+ if numel(ij)>=2
+ liuqe_keyword = eval_expr(ij(1)+2:ij(2)-1);
+ end
+ ij_row=strmatch(liuqe_keyword,liuqe_fortran_matlab(:,2),'exact');
+ if ~isempty(ij_row)
+ eval_expr = [eval_expr(1:ij(1)+1) liuqe_fortran_matlab{ij_row,1} eval_expr(ij(1)+2+length(liuqe_keyword):end)];
+ end
end
aatmp=eval(eval_expr);
end
@@ -336,7 +397,6 @@ if strcmp(mapping_for_tcv.method(1:3),'tdi')
gdat_data.data = aatmp.data;
gdat_data.dim = aatmp.dim;
nbdims = length(gdat_data.dim);
-
if mapping_for_tcv.timedim==-1;
% try to find time dim from units
idim_non1 = []; len_dim = [];
@@ -398,7 +458,7 @@ if strcmp(mapping_for_tcv.method(1:3),'tdi')
else
mapping_for_tcv.gdat_timedim = mapping_for_tcv.timedim;
end
- gdat_data.data_fullpath=[mapping_for_tcv.expression substr_tdi];
+ gdat_data.data_fullpath=[mapping_for_tcv.expression];
% end of method "tdi"
@@ -457,51 +517,113 @@ elseif strcmp(mapping_for_tcv.method,'switchcase')
% First the request names valid for "all" machines:
%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
- case {'a_minor','rgeom'}
+ case {'a_minor','rgeom','r_geom','a_minor_rho','r_geom_rho','rgeom_rho'}
+ if strcmp(data_request_eff,'r_geom'); data_request_eff = 'rgeom'; end
+ if strcmp(data_request_eff,'r_geom_rho'); data_request_eff = 'rgeom_rho'; end
% compute average minor or major radius (on z=zaxis normally)
- nodenameeff=['\results::r_max_psi' substr_liuqe];
- rmaxpsi=tdi(nodenameeff);
- ijnan = find(isnan(rmaxpsi.data));
- if isempty(rmaxpsi.data) || isempty(rmaxpsi.dim) || ischar(rmaxpsi.data) || ...
- ( ~isempty(ijnan) && prod(size(ijnan))==prod(size(rmaxpsi.data)) )
- if (gdat_params.nverbose>=1); warning(['problems loading data for ' nodenameeff ' for data_request= ' data_request_eff]); end
- if (gdat_params.nverbose>=3); disp(['rerun LIUQE?']); end
- return
- end
- nodenameeff2=['\results::r_min_psi' substr_liuqe];
- rminpsi=tdi(nodenameeff2);
- ijnan = find(isnan(rminpsi.data));
- if isempty(rminpsi.data) || isempty(rminpsi.dim) || ...
- ( ~isempty(ijnan) && prod(size(ijnan))==prod(size(rminpsi.data)) )
- if (gdat_params.nverbose>=1); warning(['problems loading data for ' nodenameeff2 ' for data_request= ' data_request_eff]); end
- if (gdat_params.nverbose>=3); disp(['rerun LIUQE?']); end
- return
- end
- ij=find(rmaxpsi.data<0.5 | rmaxpsi.data>1.2);
- if ~isempty(ij); rmaxpsi.data(ij)=NaN; end
- ij=find(rminpsi.data<0.5 | rminpsi.data>1.2);
- if ~isempty(ij); rminpsi.data(ij)=NaN; end
- if strcmp(data_request_eff,'a_minor')
- gdat_data.data=0.5.*(rmaxpsi.data(end,:) - rminpsi.data(end,:));
- gdat_data.data_fullpath=[nodenameeff ' - ' nodenameeff2 ' /2'];
- elseif strcmp(data_request_eff,'rgeom')
- gdat_data.data=0.5.*(rmaxpsi.data(end,:) + rminpsi.data(end,:));
- gdat_data.data_fullpath=[nodenameeff ' + ' nodenameeff2 ' /2'];
+ if liuqe_matlab==0
+ nodenameeff=['tcv_eq(''r_max_psi'',''LIUQE' substr_liuqe_tcv_eq ''')'];
+ rmaxpsi=tdi(nodenameeff);
+ ijnan = find(isnan(rmaxpsi.data));
+ if isempty(rmaxpsi.data) || isempty(rmaxpsi.dim) || ischar(rmaxpsi.data) || ...
+ ( ~isempty(ijnan) && prod(size(ijnan))==prod(size(rmaxpsi.data)) )
+ if (gdat_params.nverbose>=1); warning(['problems loading data for ' nodenameeff ' for data_request= ' data_request_eff]); end
+ if (gdat_params.nverbose>=3); disp(['rerun LIUQE?']); end
+ return
+ end
+ nodenameeff2=['tcv_eq(''r_min_psi'',''LIUQE' substr_liuqe_tcv_eq ''')'];
+ rminpsi=tdi(nodenameeff2);
+ ijnan = find(isnan(rminpsi.data));
+ if isempty(rminpsi.data) || isempty(rminpsi.dim) || ...
+ ( ~isempty(ijnan) && prod(size(ijnan))==prod(size(rminpsi.data)) )
+ if (gdat_params.nverbose>=1); warning(['problems loading data for ' nodenameeff2 ' for data_request= ' data_request_eff]); end
+ if (gdat_params.nverbose>=3); disp(['rerun LIUQE?']); end
+ return
+ end
+ ij=find(rmaxpsi.data<0.5 | rmaxpsi.data>1.2); % points outside TCV vessel
+ if ~isempty(ij); rmaxpsi.data(ij)=NaN; end
+ ij=find(rminpsi.data<0.5 | rminpsi.data>1.2);
+ if ~isempty(ij); rminpsi.data(ij)=NaN; end
+ if strcmp(data_request_eff,'a_minor')
+ gdat_data.data=0.5.*(rmaxpsi.data(end,:) - rminpsi.data(end,:));
+ gdat_data.data_fullpath=[nodenameeff ' - ' nodenameeff2 ' /2'];
+ elseif strcmp(data_request_eff,'a_minor_rho')
+ gdat_data.data=0.5.*(rmaxpsi.data(:,:) - rminpsi.data(:,:));
+ gdat_data.data_fullpath=[nodenameeff ' - ' nodenameeff2 ' /2'];
+ elseif strcmp(data_request_eff,'rgeom')
+ gdat_data.data=0.5.*(rmaxpsi.data(end,:) + rminpsi.data(end,:));
+ gdat_data.data_fullpath=[nodenameeff ' + ' nodenameeff2 ' /2'];
+ elseif strcmp(data_request_eff,'rgeom_rho')
+ gdat_data.data=0.5.*(rmaxpsi.data(:,:) + rminpsi.data(:,:));
+ gdat_data.data_fullpath=[nodenameeff ' + ' nodenameeff2 ' /2'];
+ else
+ if (gdat_params.nverbose>=1); warning(['should not be in this case with data_request_eff = ' data_request_eff]); end
+ return
+ end
+ if strcmp(data_request_eff(end-3:end),'_rho')
+ gdat_data.dim = rmaxpsi.dim;
+ gdat_data.t = gdat_data.dim{mapping_for_tcv.timedim};
+ gdat_data.x = gdat_data.dim{1};
+ gdat_data.dimunits = rmaxpsi.dimunits(2);
+ else
+ gdat_data.dim = rmaxpsi.dim(2);
+ gdat_data.t = gdat_data.dim{mapping_for_tcv.timedim};
+ gdat_data.dimunits = rmaxpsi.dimunits(2);
+ end
+ if any(strcmp(fieldnames(rmaxpsi),'units'))
+ gdat_data.units = rmaxpsi.units;
+ end
else
- if (gdat_params.nverbose>=1); warning(['should not be in this case with data_request_eff = ' data_request_eff]); end
- return
- end
- gdat_data.dim = rmaxpsi.dim(2);
- gdat_data.t = gdat_data.dim{1};
- if any(strcmp(fieldnames(rmaxpsi),'units'))
- gdat_data.units = rmaxpsi.units;
+ if isempty(substr_liuqe); substr_liuqe = '_1'; end
+ if strcmp(data_request_eff,'a_minor') % to update when ready, still buggy, this should be rho,t
+ nodenameeff=['tcv_eq(''a_minor_mid'',''LIUQE.M' substr_liuqe_tcv_eq ''')'];
+ elseif strcmp(data_request_eff,'a_minor_rho')
+ nodenameeff=['tcv_eq(''r_out_mid'',''LIUQE.M' substr_liuqe_tcv_eq ''')'];
+ nodenameeff2=['tcv_eq(''r_in_mid'',''LIUQE.M' substr_liuqe_tcv_eq ''')'];
+ nodenameeff=['tcv_eq(''r_out'',''LIUQE.M' substr_liuqe_tcv_eq ''')'];
+ nodenameeff2=['tcv_eq(''r_in'',''LIUQE.M' substr_liuqe_tcv_eq ''')'];
+ elseif strcmp(data_request_eff,'rgeom')
+ nodenameeff=['tcv_eq(''r_geom_mid'',''LIUQE.M' substr_liuqe_tcv_eq ''')'];
+ elseif strcmp(data_request_eff,'rgeom_rho')
+ nodenameeff=['tcv_eq(''r_out_mid'',''LIUQE.M' substr_liuqe_tcv_eq ''')'];
+ nodenameeff2=['tcv_eq(''r_in_mid'',''LIUQE.M' substr_liuqe_tcv_eq ''')'];
+ nodenameeff=['tcv_eq(''r_out'',''LIUQE.M' substr_liuqe_tcv_eq ''')'];
+ nodenameeff2=['tcv_eq(''r_in'',''LIUQE.M' substr_liuqe_tcv_eq ''')'];
+ else
+ if (gdat_params.nverbose>=1); warning(['should not be in this case with data_request_eff = ' data_request_eff]); end
+ return
+ end
+ aatmp = tdi(nodenameeff);
+ if strcmp(data_request_eff(end-3:end),'_rho')
+ aatmp2 = tdi(nodenameeff2);
+ if strcmp(data_request_eff,'a_minor_rho')
+ gdat_data.data = 0.5.*(aatmp.data-aatmp2.data);
+ elseif strcmp(data_request_eff,'rgeom_rho')
+ gdat_data.data = 0.5.*(aatmp.data+aatmp2.data);
+ end
+ gdat_data.dim = aatmp.dim; % while there is a problem with \tcv_shot::top.results.equil... dim nodes
+ gdat_data.x = gdat_data.dim{1};
+ gdat_data.data_fullpath=[nodenameeff ' and ' nodenameeff2];
+ else
+ gdat_data.data = aatmp.data;
+ gdat_data.dim = aatmp.dim(end); % while there is a problem with \tcv_shot::top.results.equil... dim nodes
+ gdat_data.data_fullpath=[nodenameeff];
+ end
+ gdat_data.t = gdat_data.dim{mapping_for_tcv.timedim};
+ gdat_data.units = aatmp.units(end);
+ gdat_data.dimunits = aatmp.dimunits;
end
- gdat_data.dimunits = rmaxpsi.dimunits(2);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
- case {'zgeom'}
+ case {'zgeom','z_geom'}
% compute average minor or major radius (on z=zaxis normally)
- nodenameeff=['\results::z_contour' substr_liuqe];
+ %
+ if liuqe_matlab==0
+ nodenameeff=['tcv_eq(''z_contour'',''LIUQE.M' substr_liuqe_tcv_eq ''')'];
+ else
+ if isempty(substr_liuqe); substr_liuqe = '_1'; end
+ nodenameeff=['tcv_eq(''z_edge'',''LIUQE.M' substr_liuqe_tcv_eq ''')'];
+ end
zcontour=tdi(nodenameeff);
if isempty(zcontour.data) || isempty(zcontour.dim) % || ischar(zcontour.data) (to add?)
if (gdat_params.nverbose>=1); warning(['problems loading data for ' nodenameeff ' for data_request= ' data_request_eff]); end
@@ -526,7 +648,7 @@ elseif strcmp(mapping_for_tcv.method,'switchcase')
case {'b0'}
% B0 at R0=0.88
R0EXP=0.88;
- if liuqe_version==-1
+ if liuqe_version_eff==-1
nodenameeff = 'tcv_eq("BZERO","FBTE")';
tracetdi=tdi(nodenameeff);
gdat_data.data = tracetdi.data;
@@ -752,19 +874,20 @@ elseif strcmp(mapping_for_tcv.method,'switchcase')
time_eff = time(itime);
% use read_results updated to effectively obtain an eqdsk with sign correct with COCOS=2
[fnames_readresults]=read_results_for_chease(shot,time_eff,liuqe_version,3,[],[],[],zshift,0,1,nrz_out);
+ cocos_read_results_for_chease = 17;
if isempty(fnames_readresults)
if gdat_params.nverbose>=1;
warning(['could not create eqdsk file from read_results_for_chease with data_request= ' data_request_eff]);
end
return
end
- eqdskval=read_eqdsk(fnames_readresults{4},7,0,[],[],1); % LIUQE is 17 but read_results divided psi by 2pi thus 7
+ eqdskval=read_eqdsk(fnames_readresults{4},cocos_read_results_for_chease,0,[],[],1); % read_results now has relevant cocos LIUQE=17
for i=1:length(fnames_readresults)
unix(['rm ' fnames_readresults{i}]);
end
% transform to cocos=2 since read_results originally assumed it was cocos=2
cocos_in = 2;
- [eqdsk_cocos_in, eqdsk_cocosout_IpB0pos,cocos_inout]=eqdsk_cocos_transform(eqdskval,[7 cocos_in]);
+ [eqdsk_cocos_in, eqdsk_cocosout_IpB0pos,cocos_inout]=eqdsk_cocos_transform(eqdskval,[cocos_read_results_for_chease cocos_in]);
fnamefull = fullfile(['/tmp/' getenv('USER')],['EQDSK_' num2str(shot) 't' num2str(time_eff,'%.4f')]);
% We still write COCOS=2 case, since closer to standard (in /tmp)
write_eqdsk(fnamefull,eqdsk_cocos_in,cocos_in,[],[],[],1);
@@ -868,10 +991,7 @@ elseif strcmp(mapping_for_tcv.method,'switchcase')
gdat_data.n2_HFS.data = gdat_data.n2_HFS.data + aaHFSz0_sect11.data;
gdat_data.data_fullpath=['\atlas::DT196_MHD_001:channel_083 -+ \atlas::DT196_MHD_001:channel_091 for n=1,2, LFS_sect3/11, z=0cm; _HFS' ...
' same for HFS: MHD_002:channel_034-+MHD_002:channel_038'];
-
- end
-
-
+ end
else
disp('should not be here in ''mhd'', ask O. Sauter')
return
@@ -921,11 +1041,22 @@ elseif strcmp(mapping_for_tcv.method,'switchcase')
edge_str_ = '_edge';
edge_str_dot = '.edge';
end
- psi_max=gdat([],['\results::thomson' edge_str_dot ':psi_max' substr_liuqe],'nverbose',gdat_params.nverbose);
- psiscatvol=gdat([],['\results::thomson' edge_str_dot ':psiscatvol' substr_liuqe],'nverbose',gdat_params.nverbose);
- if abs(zshift)>1e-5
+ if liuqe_matlab==1 && strcmp(substr_liuqe,'_1'); substr_liuqe = ''; end
+ % psiscatvol obtained from linear interpolation in fact so not quite ok near axis where psi is quadratic
+ recompute_psiscatvol_always = 1;
+ if liuqe_version==-1; recompute_psiscatvol_always = 1; end
+ if abs(zshift)<1e-5 && liuqe_matlab==0 && recompute_psiscatvol_always==0
+ psi_max=gdat([],['\results::thomson' edge_str_dot ':psi_max' substr_liuqe],'nverbose',gdat_params.nverbose);
+ psiscatvol=gdat([],['\results::thomson' edge_str_dot ':psiscatvol' substr_liuqe],'nverbose',gdat_params.nverbose);
+ else
% calculate new psiscatvol
- psitdi=tdi(['\results::psi' substr_liuqe]);
+ if liuqe_matlab==0
+ psitdi = tdi(['tcv_eq(''psi'',''LIUQE' substr_liuqe_tcv_eq ''')']);
+ psiaxis = tdi(['tcv_eq(''psi_axis'',''LIUQE' substr_liuqe_tcv_eq ''')']);
+ else
+ psitdi = tdi(['tcv_eq(''psi'',''LIUQE.M' substr_liuqe_tcv_eq ''')']);
+ psiaxis = tdi(['tcv_eq(''psi_axis'',''LIUQE.M' substr_liuqe_tcv_eq ''')']);
+ end
rmesh=psitdi.dim{1};
zmesh=psitdi.dim{2};
zthom=mdsdata('dim_of(\thomson:te,1)');
@@ -936,7 +1067,7 @@ elseif strcmp(mapping_for_tcv.method,'switchcase')
zeffshift=zeffshift * ones(size(psitdi.dim{3}));
case length(psitdi.dim{3})
% ok
- case length(psiscatvol.dim{1})
+ case length(gdat_data.t)
zeffshift=interp1(psiscatvol.dim{1},zeffshift,psitdi.dim{3});
otherwise
if (gdat_params.nverbose>=1);
@@ -944,16 +1075,25 @@ elseif strcmp(mapping_for_tcv.method,'switchcase')
disp(['it should be 1 or ' num2str(length(psiscatvol.dim{1})) ' or ' num2str(length(psitdi.dim{3}))])
end
end
- for it=1:length(psiscatvol.dim{1})
- itpsitdi=iround(psitdi.dim{3},psiscatvol.dim{1}(it));
+ for it=1:length(gdat_data.t)
+ itpsitdi=iround_os(psitdi.dim{3},gdat_data.t(it));
psirz=psitdi.data(:,:,itpsitdi);
- psiscatvol0=griddata(rmesh,zmesh,psirz',0.9*ones(size(zthom)),zthom-zeffshift(itpsitdi));
+ %psiscatvol0=griddata(rmesh,zmesh,psirz',0.9*ones(size(zthom)),zthom-zeffshift(itpsitdi),'cubic'); % faster with interpos
+ psiscatvol0=interpos2Dcartesian(rmesh,zmesh,psirz,0.9*ones(size(zthom)),zthom-zeffshift(itpsitdi),-0.1,-0.1);
+ psiscatvol0=griddata(rmesh,zmesh,psirz',0.9*ones(size(zthom)),zthom-zeffshift(itpsitdi),'linear');
psiscatvol.data(it,:)=psiscatvol0;
+ % since take closest psi(R,Z) from psitdi, should also take closest for psi_max and not interpolating
+ itpsiaxis = iround_os(psiaxis.dim{1},gdat_data.t(it));
+ psi_max.data(it,1) = psiaxis.data(itpsiaxis);
end
+ psiscatvol.dim{1} = gdat_data.t;
+ psiscatvol.dim{2} = gdat_data.x;
end
if ~isempty(psiscatvol.data) && ~ischar(psiscatvol.data) && ~isempty(psi_max.data) && ~ischar(psi_max.data)
for ir=1:length(psiscatvol.dim{2})
- rho(ir,:)= sqrt(1.-psiscatvol.data(:,ir)./psi_max.data(:))';
+ rho2 = max(1.-psiscatvol.data(:,ir)./psi_max.data,0);
+ rho(ir,:)= sqrt(rho2');
+ % rho(ir,:)= sqrt(1.-psiscatvol.data(:,ir)./psi_max.data(:))';
end
else
if gdat_params.nverbose>=1 && gdat_data.gdat_params.edge==0
@@ -1274,6 +1414,8 @@ elseif strcmp(mapping_for_tcv.method,'switchcase')
gdat_data.nbi.data_fullpath=[nodenameeff ', index ' num2str(index_for_tot_pow) ', CLC NEUTRAL POWER'];
gdat_data.nbi.label='P_{nbi}';
gdat_data.nbi.help = nbh_data_tdi.help;
+ index_for_energy = 33;
+ gdat_data.nbi.energy = nbh_data_tdi.data(:,index_for_energy);
% add to main with linear interpolation and 0 for extrapolated values
gdat_data.data(:,end+1) = interpos(-21,gdat_data.nbi.t,gdat_data.nbi.data,gdat_data.t);
gdat_data.x(end+1) = size(gdat_data.data,2);
@@ -1291,8 +1433,12 @@ elseif strcmp(mapping_for_tcv.method,'switchcase')
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
case {'q_rho'}
% q profile on psi from liuqe
- nodenameeff=['\results::q_psi' substr_liuqe];
- if liuqe_version==-1
+ if liuqe_matlab==0
+ nodenameeff=['tcv_eq(''q_psi'',''LIUQE' substr_liuqe_tcv_eq ''')'];
+ else
+ nodenameeff=['tcv_eq(''q'',''LIUQE.M' substr_liuqe_tcv_eq ''')'];
+ end
+ if liuqe_version_eff==-1
nodenameeff=[begstr 'q_psi' substr_liuqe];
end
tracetdi=tdi(nodenameeff);
@@ -1305,9 +1451,11 @@ elseif strcmp(mapping_for_tcv.method,'switchcase')
gdat_data.dim = tracetdi.dim;
gdat_data.t = gdat_data.dim{2};
gdat_data.data_fullpath=[nodenameeff ' on rhopol'];
- rhopol_eff = ones(size(tracetdi.dim{1}));
- rhopol_eff(:) = sqrt(linspace(0,1,length(tracetdi.dim{1})));
- gdat_data.dim{1} = rhopol_eff;
+ if liuqe_matlab==0
+ rhopol_eff = ones(size(tracetdi.dim{1}));
+ rhopol_eff(:) = sqrt(linspace(0,1,length(tracetdi.dim{1})));
+ gdat_data.dim{1} = rhopol_eff;
+ end
gdat_data.x = gdat_data.dim{1};
gdat_data.dimunits{1} = 'rho_pol~sqrt(\psi_norm)';
gdat_data.dimunits{2} = 's';
@@ -1321,7 +1469,7 @@ elseif strcmp(mapping_for_tcv.method,'switchcase')
% psi at edge, 0 by construction in Liuqe, thus not given
% add surface_psi from surface_flux and d(surface_flux)/dt = vloop
nodenameeff=['\results::psi_axis' substr_liuqe];
- if liuqe_version==-1
+ if liuqe_version_eff==-1
nodenameeff=[begstr 'q_psi' substr_liuqe];
end
tracetdi=tdi(nodenameeff);
@@ -1346,86 +1494,135 @@ elseif strcmp(mapping_for_tcv.method,'switchcase')
gdat_data.vsurf_description = ['time derivative from surface_psi, obtained from \results::surface_flux'];
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
- case {'rhotor_edge','rhotor'}
+ case {'rhotor_edge', 'rhotor', 'rhotor_norm'}
% Phi(LCFS) = int(Bphi dSphi), can use Eq.(11) of "Tokamak coordinate conventions: COCOS" paper:
% O. Sauter, S.Yu. Medvedev, Comput. Phys. Commun. 184 (2013) 293–302
% since cocos=17 for LIUQE we get:
% q = -dPhi/dpsi => Phi = - int(q*dpsi) which should always have the sign of B0
% need to get q_rho but to avoid loop for rhotor in grids_1d, get q_rho explicitely here
- nodenameeff=['\results::q_psi' substr_liuqe];
- if liuqe_version==-1
- nodenameeff=[begstr 'q_psi' substr_liuqe];
- end
- q_rho=tdi(nodenameeff);
- if isempty(q_rho.data) || isempty(q_rho.dim) % || ischar(q_rho.data) (to add?)
- if (gdat_params.nverbose>=1); warning(['problems loading data for q_rho for data_request= ' data_request_eff]); end
- if (gdat_params.nverbose>=3); disp(['rerun LIUQE?']); end
- return
- end
- rhopol_eff = ones(size(q_rho.dim{1}));
- rhopol_eff(:) = sqrt(linspace(0,1,length(q_rho.dim{1})));
- q_rho.dim{1} = rhopol_eff;
- params_eff.data_request='psi_axis'; % psi_edge=0 with LIUQE
- psi_axis=gdat_tcv([],params_eff);
- params_eff.data_request='b0'; % psi_edge=0 with LIUQE
- b0=gdat_tcv([],params_eff);
- b0tpsi = interp1(b0.t,b0.data,psi_axis.t); %q_rho on same time base as psi_axis
- if isempty(psi_axis.data) || isempty(psi_axis.dim) || isempty(q_rho.data) || isempty(q_rho.dim)
- if (gdat_params.nverbose>=1); warning(['problems loading data for ' nodenameeff ' for data_request= ' data_request_eff]); end
- return
- end
- rhoequal = linspace(0,1,length(q_rho.dim{1}));
- if strcmp(data_request_eff,'rhotor_edge')
- gdat_data.data = psi_axis.data; % to have the dimensions correct
- gdat_data.dim = psi_axis.dim;
- gdat_data.t = gdat_data.dim{1};
- gdat_data.data_fullpath='phi from q_rho, psi_axis and integral(-q dpsi)';
- gdat_data.units = 'T m^2';
- gdat_data.dimunits{1} = 's';
- elseif strcmp(data_request_eff,'rhotor')
- gdat_data.data = q_rho.data; % to have the dimensions correct
- gdat_data.dim{1} = ones(size(q_rho.dim{1}));
- gdat_data.dim{1}(:) = rhoequal;
- gdat_data.dim{2} = q_rho.dim{2};
- gdat_data.t = gdat_data.dim{2};
- gdat_data.data_fullpath='sqrt(phitor/phitor_edge), rhotor_edge=sqrt(phitor/B0/pi)';
- gdat_data.units = '';
- gdat_data.dimunits{1} = 'rhopol\_norm';
- gdat_data.dimunits{2} = 's';
- end
- gdat_data.x=gdat_data.dim{1};
- gdat_data.psi_axis = reshape(psi_axis.data,1,length(psi_axis.data));
- if length(gdat_data.psi_axis) ~= length(gdat_data.t)
- disp('problems of time between qrho and psi_axis?')
- return
- end
- for it=1:length(psi_axis.data)
- ij=find(~isnan(q_rho.data(:,it)));
- if ~isempty(ij)
- [qfit,~,~,phi]=interpos(q_rho.dim{1}(ij).^2,q_rho.data(ij,it),rhoequal.^2);
- dataeff = sqrt(phi .* psi_axis.data(it) ./ b0tpsi(it) ./ pi) ; % Delta_psi = -psi_axis
- else
- dataeff = NaN;
+ params_eff = gdat_data.gdat_params;
+ if liuqe_matlab==0
+ nodenameeff=['tcv_eq(''q_psi'',''LIUQE' substr_liuqe_tcv_eq ''')'];
+ if liuqe_version_eff==-1
+ nodenameeff=[begstr 'q_psi' substr_liuqe];
+ end
+ q_rho=tdi(nodenameeff);
+ if isempty(q_rho.data) || isempty(q_rho.dim) % || ischar(q_rho.data) (to add?)
+ if (gdat_params.nverbose>=1); warning(['problems loading data for q_rho for data_request= ' data_request_eff]); end
+ if (gdat_params.nverbose>=3); disp(['rerun LIUQE?']); end
+ return
+ end
+ rhopol_eff = ones(size(q_rho.dim{1}));
+ rhopol_eff(:) = sqrt(linspace(0,1,length(q_rho.dim{1})));
+ q_rho.dim{1} = rhopol_eff;
+ params_eff.data_request='psi_axis'; % psi_edge=0 with LIUQE
+ psi_axis=gdat_tcv(shot,params_eff);
+ params_eff.data_request='b0';
+ b0=gdat_tcv(shot,params_eff);
+ b0tpsi = interp1(b0.t,b0.data,psi_axis.t); %q_rho on same time base as psi_axis
+ if isempty(psi_axis.data) || isempty(psi_axis.dim) || isempty(q_rho.data) || isempty(q_rho.dim)
+ if (gdat_params.nverbose>=1); warning(['problems loading data for ' nodenameeff ' for data_request= ' data_request_eff]); end
+ return
+ end
+ rhoequal = linspace(0,1,length(q_rho.dim{1}));
+ if strcmp(data_request_eff,'rhotor_edge')
+ gdat_data.data = psi_axis.data; % to have the dimensions correct
+ gdat_data.dim = psi_axis.dim;
+ gdat_data.t = gdat_data.dim{1};
+ gdat_data.data_fullpath='phi from q_rho, psi_axis and integral(-q dpsi)';
+ gdat_data.units = 'T m^2';
+ gdat_data.dimunits{1} = 's';
+ elseif strcmp(data_request_eff,'rhotor')
+ gdat_data.data = q_rho.data; % to have the dimensions correct
+ gdat_data.dim{1} = ones(size(q_rho.dim{1}));
+ gdat_data.dim{1}(:) = rhoequal;
+ gdat_data.dim{2} = q_rho.dim{2};
+ gdat_data.t = gdat_data.dim{2};
+ gdat_data.data_fullpath='sqrt(phitor/phitor_edge), rhotor_edge=sqrt(phitor/B0/pi)';
+ gdat_data.units = '';
+ gdat_data.dimunits{1} = 'rhopol\_norm';
+ gdat_data.dimunits{2} = 's';
+ end
+ gdat_data.x=gdat_data.dim{1};
+ gdat_data.psi_axis = reshape(psi_axis.data,1,length(psi_axis.data));
+ if length(gdat_data.psi_axis) ~= length(gdat_data.t)
+ disp('problems of time between qrho and psi_axis?')
+ return
end
+ for it=1:size(q_rho.data,2)
+ ij=find(~isnan(q_rho.data(:,it)));
+ if ~isempty(ij)
+ [qfit,~,~,phi]=interpos(q_rho.dim{1}(ij).^2,q_rho.data(ij,it),rhoequal.^2);
+ dataeff = sqrt(phi .* psi_axis.data(it) ./ b0tpsi(it) ./ pi) ; % Delta_psi = -psi_axis
+ else
+ dataeff = NaN;
+ end
+ if strcmp(data_request_eff,'rhotor_edge')
+ gdat_data.data(it) = dataeff(end);
+ elseif strcmp(data_request_eff,'rhotor')
+ gdat_data.data(:,it) = dataeff./dataeff(end);
+ gdat_data.rhotor_edge(it) = dataeff(end);
+ end
+ gdat_data.b0 = b0tpsi(it);
+ end
+ else
+ params_eff = gdat_data.gdat_params;
+ params_eff.data_request='phi_tor';
+ phi_tor=gdat_tcv([],params_eff);
+ params_eff = gdat_data.gdat_params;
+ params_eff.data_request='b0';
+ b0=gdat_tcv([],params_eff);
+ gdat_data.t = phi_tor.t;
+ ij=find(~isnan(b0.data));
+ gdat_data.b0 = interpos(b0.t(ij),b0.data(ij),gdat_data.t);
if strcmp(data_request_eff,'rhotor_edge')
- gdat_data.data(it) = dataeff(end);
+ gdat_data.data = sqrt(phi_tor.data(end,:)./pi./reshape(gdat_data.b0,1,numel(gdat_data.b0)));
+ gdat_data.dim = phi_tor.dim(2); % while there is a problem with \tcv_shot::top.results.equil... dim nodes
+ gdat_data.dimunits = phi_tor.dimunits{2};
+ gdat_data.units = 'm';
+ gdat_data.data_fullpath=['rhotor_edge=sqrt(Phi(edge)/pi/B0) from phi_tor'];
elseif strcmp(data_request_eff,'rhotor')
- gdat_data.data(:,it) = dataeff./dataeff(end);
- gdat_data.rhotor_edge(it) = dataeff(end);
+ gdat_data.data = sqrt(phi_tor.data./pi./(ones(size(phi_tor.data,1),1)*reshape(gdat_data.b0,1,numel(gdat_data.b0))));
+ gdat_data.rhotor_edge = gdat_data.data(end,:);
+ gdat_data.dim = phi_tor.dim; % while there is a problem with \tcv_shot::top.results.equil... dim nodes
+ gdat_data.x=gdat_data.dim{1};
+ gdat_data.units = 'm';
+ gdat_data.dimunits = phi_tor.dimunits;
+ gdat_data.data_fullpath=['rhotor=sqrt(Phi/pi/B0) from phi_tor'];
+ elseif strcmp(data_request_eff,'rhotor_norm')
+ gdat_data.data = sqrt(phi_tor.data./(ones(size(phi_tor.data,1),1)*reshape(phi_tor.data(end,:),1,size(phi_tor.data,2))));
+ gdat_data.rhotor_edge = sqrt(phi_tor.data(end,:)./pi./reshape(gdat_data.b0,1,numel(gdat_data.b0)));
+ gdat_data.dim = phi_tor.dim; % while there is a problem with \tcv_shot::top.results.equil... dim nodes
+ gdat_data.x=gdat_data.dim{1};
+ gdat_data.units = '';
+ gdat_data.dimunits = phi_tor.dimunits;
+ gdat_data.data_fullpath=['rhotor_norm=sqrt(Phi/Phi(edge)) from phi_tor'];
+ else
+ disp(['data_request_eff = ' data_request_eff ' not defined within rhotor block in gdat_tcv.m']);
+ return
end
- gdat_data.b0 = b0tpsi(it);
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
- case {'rhovol','volume_rho','volume'}
+ case {'rhovol','rho_vol','volume_rho','volume'}
+ if strcmp(data_request_eff,'rho_vol'); data_request_eff = 'rhovol'; end
% volume_rho = vol(rho); volume = vol(LCFS) = vol(rho=1);
% rhovol = sqrt(vol(rho)/vol(rho=1));
- nodenameeff='\results::psitbx:vol';
- if liuqe_version==-1
- nodenameeff=[begstr 'vol' substr_liuqe];
+ if liuqe_matlab==0
+ nodenameeff=['\results::psitbx:vol'];
+ if liuqe_version_eff > 1
+ disp('needs to construct volume')
+ return
+ end
+ else
+ nodenameeff=['tcv_eq(''vol'',''LIUQE.M' substr_liuqe_tcv_eq ''')'];
+ end
+ if liuqe_version_eff==-1
+ data_request_eff = 'volume'; % only LCFS
+ nodenameeff=[begstr 'volume' substr_liuqe];
end
tracetdi=tdi(nodenameeff);
- if isempty(tracetdi.data) || isempty(tracetdi.dim)
+ if (isempty(tracetdi.data) || isempty(tracetdi.dim)) && liuqe_matlab==0
% try to run psitbxput
psitbxput_version = 1.3;
psitbxput(psitbxput_version,shot);
@@ -1437,23 +1634,31 @@ elseif strcmp(mapping_for_tcv.method,'switchcase')
end
gdat_data.units = tracetdi.units;
if strcmp(data_request_eff,'volume')
- gdat_data.data = tracetdi.data(end,:);
- gdat_data.dim{1} = tracetdi.dim{2};
- gdat_data.data_fullpath=['\results::psitbx:vol(end,:)'];
- gdat_data.dimunits{1} = tracetdi.dimunits{2};
+ if length(gdat_data.dim >=2)
+ gdat_data.data = tracetdi.data(end,:);
+ gdat_data.dim{1} = tracetdi.dim{2};
+ gdat_data.dimunits{1} = tracetdi.dimunits{2};
+ else
+ mapping_for_tcv.gdat_timedim = 1;
+ gdat_data.data = tracetdi.data;
+ gdat_data.dim = tracetdi.dim;
+ gdat_data.dimunits = tracetdi.dimunits;
+ gdat_data.volume_edge = gdat_data.data;
+ end
gdat_data.request_description = 'volume(LCFS)=volume(rhopol=1)';
+ gdat_data.data_fullpath=['Volume of LCFS from ' nodenameeff];
else
gdat_data.data = tracetdi.data;
gdat_data.dim = tracetdi.dim;
gdat_data.x = gdat_data.dim{1};
gdat_data.dimunits = tracetdi.dimunits;
if strcmp(data_request_eff,'volume_rho')
- gdat_data.data_fullpath=['\results::psitbx:vol'];
gdat_data.request_description = 'volume(rho)';
+ gdat_data.data_fullpath=['Volume(rho) from ' nodenameeff];
elseif strcmp(data_request_eff,'rhovol')
gdat_data.volume_edge = gdat_data.data(end,:);
gdat_data.data = sqrt(gdat_data.data./repmat(reshape(gdat_data.volume_edge,1,size(gdat_data.data,2)),size(gdat_data.data,1),1));
- gdat_data.data_fullpath='sqrt(\results::psitbx:vol/vol_edge)';
+ gdat_data.data_fullpath=['sqrt(volume/volume(edge) from ' nodenameeff];
gdat_data.request_description = 'sqrt(volume(rho)/volume(edge))';
else
if (gdat_params.nverbose>=1)
@@ -1462,6 +1667,7 @@ elseif strcmp(mapping_for_tcv.method,'switchcase')
return
end
end
+ mapping_for_tcv.gdat_timedim
gdat_data.t = gdat_data.dim{mapping_for_tcv.gdat_timedim};
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
@@ -1516,8 +1722,8 @@ elseif strcmp(mapping_for_tcv.method,'switchcase')
end
% Get data for the source requested
- for ii=1:numel(gdat_data.gdat_params.source)
- switch gdat_data.gdat_params.source{ii}
+ for isource=1:numel(gdat_data.gdat_params.source)
+ switch gdat_data.gdat_params.source{isource}
case 'defined'
DS = take_all_signals(shot); %Read from mds
SDS_DS = define_simulink_signals(DS); %Put them in predifined structure
@@ -1527,33 +1733,34 @@ elseif strcmp(mapping_for_tcv.method,'switchcase')
node = ii;
for jj=1:numel(SDS_DS{ii}) %iter over threads
thread = jj;
- fieldnameslist = fieldnames(SDS_DS{ii}{jj});
-
- if SDS_DS{ii}{1}.conf.hasthreads
- is_with_threads = 1;
- else
- is_with_threads = 0;
- end
- for kk=2:numel(fieldnameslist) %iter over fieldnames (the first one is configuration)
- indices = SDS_DS{ii}{jj}.(fieldnameslist{kk}).ind;
- SDS_DS{ii}{jj}.(fieldnameslist{kk}).data = [];
- for zz=1:numel(indices) %iter over indices
- ind = indices(zz);
- mdsconnect('scd');
- mdsopen('rtc', shot);
- % data expression
- if is_with_threads ==0
- expression = sprintf('\\top.crpprt%.2d.mems.mem_%.2d',node,ind);
- else
- expression = sprintf('\\top.crpprt%.2d.thread%.1d.mems.mem_%.3d',node,thread,ind);
- end
- tmp = tdi(expression);
- if isnumeric(tmp.data) && ~isempty(tmp.data)
- SDS_DS{ii}{jj}.(fieldnameslist{kk}).data =[SDS_DS{ii}{jj}.(fieldnameslist{kk}).data; tmp.data'];
-
- SDS_DS{ii}{jj}.(fieldnameslist{kk}).t = tmp.dim{1};
- else
- fprintf('Warning node: %d thread: %d signal: %s ind %d not available\n', ii,jj,fieldnameslist{kk}, zz );
+ if ~isempty(SDS_DS{ii}{jj})
+ fieldnameslist = fieldnames(SDS_DS{ii}{jj});
+ if SDS_DS{ii}{jj}.conf.hasthreads
+ is_with_threads = 1;
+ else
+ is_with_threads = 0;
+ end
+ for kk=2:numel(fieldnameslist) %iter over fieldnames (the first one is configuration)
+ indices = SDS_DS{ii}{jj}.(fieldnameslist{kk}).ind;
+ SDS_DS{ii}{jj}.(fieldnameslist{kk}).data = [];
+ for zz=1:numel(indices) %iter over indices
+ ind = indices(zz);
+ mdsconnect('scd');
+ mdsopen('rtc', shot);
+ % data expression
+ if is_with_threads ==0
+ expression = sprintf('\\top.crpprt%.2d.mems.mem_%.2d',node,ind);
+ else
+ expression = sprintf('\\top.crpprt%.2d.thread%.1d.mems.mem_%.3d',node,thread,ind);
+ end
+ tmp = tdi(expression);
+ if isnumeric(tmp.data) && ~isempty(tmp.data)
+ SDS_DS{ii}{jj}.(fieldnameslist{kk}).data =[SDS_DS{ii}{jj}.(fieldnameslist{kk}).data; tmp.data'];
+
+ SDS_DS{ii}{jj}.(fieldnameslist{kk}).t = tmp.dim{1};
+ else
+ fprintf('Warning node: %d thread: %d signal: %s ind %d not available\n', ii,jj,fieldnameslist{kk}, zz );
+ end
end
end
end
@@ -1582,29 +1789,33 @@ elseif strcmp(mapping_for_tcv.method,'switchcase')
node = ii;
for jj=1:numel(SDS_AS{ii}) %iter over threads
thread = jj;
- fieldnameslist = fieldnames(SDS_AS{ii}{jj});
- if numel(SDS_AS{ii})>1
- is_with_threads = 1;
- else
- is_with_threads = 0;
- end
- for kk=2:numel(fieldnameslist) %iter over fieldnames (the first one is configuration)
- indices = SDS_AS{ii}{jj}.(fieldnameslist{kk}).ind;
- SDS_AS{ii}{jj}.(fieldnameslist{kk}).data = [];
- for zz=1:numel(indices) %iter over indices
- ind = indices(zz);
- % data expression
- if is_with_threads ==0
- expression = sprintf('\\top.crpprt%.2d.mems.mem_%.2d',node,ind);
- else
- expression = sprintf('\\top.crpprt%.2d.thread%.1d.mems.mem_%.3d',node,thread,ind);
- end
- tmp = tdi(expression);
- if isnumeric(tmp.data)
- SDS_AS{ii}{jj}.(fieldnameslist{kk}).data =[SDS_AS{ii}{jj}.(fieldnameslist{kk}).data; tmp.data'];
- SDS_AS{ii}{jj}.(fieldnameslist{kk}).t = tmp.dim{1};
- else
- fprintf('Warning node: %d thread: %d signal: %s ind %d not available\n', ii,jj,fieldnameslist{kk}, zz );
+ if ~isempty(SDS_AS{ii}{jj})
+ fieldnameslist = fieldnames(SDS_AS{ii}{jj});
+ if numel(SDS_AS{ii})>1
+ is_with_threads = 1;
+ else
+ is_with_threads = 0;
+ end
+ for kk=2:numel(fieldnameslist) %iter over fieldnames (the first one is configuration)
+ indices = SDS_AS{ii}{jj}.(fieldnameslist{kk}).ind;
+ SDS_AS{ii}{jj}.(fieldnameslist{kk}).data = [];
+ for zz=1:numel(indices) %iter over indices
+ ind = indices(zz);
+ % data expression
+ if is_with_threads ==0
+ expression = sprintf('\\top.crpprt%.2d.mems.mem_%.2d',node,ind);
+ else
+ expression = sprintf('\\top.crpprt%.2d.thread%.1d.mems.mem_%.3d',node,thread,ind);
+ end
+ tmp = tdi(expression);
+ if isnumeric(tmp.data)
+ SDS_AS{ii}{jj}.(fieldnameslist{kk}).data =[SDS_AS{ii}{jj}.(fieldnameslist{kk}).data; tmp.data'];
+ SDS_AS{ii}{jj}.(fieldnameslist{kk}).t = tmp.dim{1};
+ else
+ if gdat_params.nverbose>=3
+ fprintf('Warning node: %d thread: %d signal: %s ind %d not available\n', ii,jj,fieldnameslist{kk},zz);
+ end
+ end
end
end
end
@@ -1649,7 +1860,7 @@ elseif strcmp(mapping_for_tcv.method,'switchcase')
else
gdat_data.gdat_params.camera = lower(gdat_data.gdat_params.camera);
end
- if ~any(liuqe_version==[1, 2, 3])
+ if ~any(liuqe_version_eff==[1, 2, 3])
if gdat_data.gdat_params.nverbose>=3
disp(['liuqe_version = ' liuqe_version ' not supported for data_request= ' data_request_eff]);
end
@@ -1666,7 +1877,7 @@ elseif strcmp(mapping_for_tcv.method,'switchcase')
gdat_data.bottom.x = [];
gdat_data.bottom.channel = [];
try
- mpx = mpxdata(shot,'svgr','freq',freq_opt,'liuqe',liuqe_version,'detec',gdat_data.gdat_params.camera, ...
+ mpx = mpxdata(shot,'svgr','freq',freq_opt,'liuqe',liuqe_version_eff,'detec',gdat_data.gdat_params.camera, ...
'time',t_int);
catch
if gdat_data.gdat_params.nverbose>=1
@@ -2015,8 +2226,9 @@ params_eff.data_request='rhotor';
rhotor = gdat(gdat_data.shot,params_eff);
params_eff.data_request='rhovol';
rhovol = gdat(gdat_data.shot,params_eff);
-
-psi0 = interpos(rhotor.t',rhotor.psi_axis,gdat_data.t,-0.01);
+params_eff.data_request='psi_axis';
+psi_axis = gdat(gdat_data.shot,params_eff);
+psi0 = interpos(psi_axis.t,psi_axis.data,gdat_data.t,-0.01);
if (nbdim_x == 1)
gdat_data.grids_1d.psi = gdat_data.grids_1d.rhopolnorm.^2*reshape(psi0,1,length(psi0));
elseif (nbdim_x == 2)
@@ -2041,12 +2253,16 @@ for it=1:length(gdat_data.t)
if (nbdim_x == 1)
if length(ii)==length(gdat_data.grids_1d.rhopolnorm)
gdat_data.grids_1d.rhotornorm(:,it)=interpos(rhotor.x,rhotor.data(:,it_rt_eff),gdat_data.grids_1d.rhopolnorm);
- gdat_data.grids_1d.rhovolnorm(:,it)=interpos(rhovol.x,rhovol.data(:,it_vol_eff),gdat_data.grids_1d.rhopolnorm);
+ if ~isempty(rhovol.x)
+ gdat_data.grids_1d.rhovolnorm(:,it)=interpos(rhovol.x,rhovol.data(:,it_vol_eff),gdat_data.grids_1d.rhopolnorm);
+ end
end
else
if length(ii)==size(gdat_data.grids_1d.rhopolnorm,1)
gdat_data.grids_1d.rhotornorm(:,it)=interpos(rhotor.x,rhotor.data(:,it_rt_eff),gdat_data.grids_1d.rhopolnorm(:,it));
- gdat_data.grids_1d.rhovolnorm(:,it)=interpos(rhovol.x,rhovol.data(:,it_vol_eff),gdat_data.grids_1d.rhopolnorm(:,it));
+ if ~isempty(rhovol.x)
+ gdat_data.grids_1d.rhovolnorm(:,it)=interpos(rhovol.x,rhovol.data(:,it_vol_eff),gdat_data.grids_1d.rhopolnorm(:,it));
+ end
end
end
end
diff --git a/crpptbx/TCV/tcv_help_parameters.m b/crpptbx/TCV/tcv_help_parameters.m
index f2fda0d7cdaaef31fac233d2addf955bc39b00dd..683e1ce1b1b1e9c1ef256d266267110bc6e20420 100644
--- a/crpptbx/TCV/tcv_help_parameters.m
+++ b/crpptbx/TCV/tcv_help_parameters.m
@@ -18,7 +18,7 @@ help_struct_all = struct(...
'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'', 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' ...
- ,'liuqe','liuqe version 1 (default), 2, 3 for LIUQE1, 2, 3 resp. or -1 for model values' ...
+ ,'liuqe','liuqe version 1 (default), 2, 3 for LIUQE1, 2, 3 resp. or -1 for model values (11, 12, 13 for liuqe_fortran)' ...
,'nverbose','1 (default) displays warnings, 0: only errors, >=3: displays all extra information' ...
);
diff --git a/crpptbx/TCV/tcv_requests_mapping.m b/crpptbx/TCV/tcv_requests_mapping.m
index f4f038b7bb172e973c37748bbfdaff83ca6d9171..07d4e7469830176e1c3b60be5eba6f3f55378f1e 100644
--- a/crpptbx/TCV/tcv_requests_mapping.m
+++ b/crpptbx/TCV/tcv_requests_mapping.m
@@ -1,7 +1,9 @@
-function mapping = tcv_requests_mapping(data_request)
+function mapping = tcv_requests_mapping(data_request,shot)
%
% Information pre-defined for gdat_tcv, you can add cases here to match official cases in gdat_tcv, allowing backward compatibility
%
+% give the shot number in case data origin depends on the shot number, allows to adapt easily
+%
% Defaults
mapping = struct(...
@@ -36,6 +38,16 @@ switch lower(data_request)
mapping.label = 'a(LCFS)';
mapping.method = 'switchcase';
mapping.expression = '';
+ case 'a_minor_rho'
+ mapping.timedim = 2;
+ mapping.label = 'a(rho,t)';
+ mapping.method = 'switchcase';
+ mapping.expression = '';
+ case 'area'
+ mapping.timedim = 1;
+ mapping.label = 'area';
+ mapping.method = 'tdiliuqe';
+ mapping.expression = 'tcv_eq(''''area'''',''''LIUQE.M'''')';
case 'b0'
mapping.timedim = 1;
mapping.label = 'B_0';
@@ -46,6 +58,7 @@ switch lower(data_request)
mapping.label = '\beta';
mapping.method = 'tdiliuqe';
mapping.expression = '\results::beta_tor';
+ mapping.expression = 'tcv_eq(''''beta_tor'''',''''LIUQE.M'''')';
case 'betan'
mapping.timedim = 1;
mapping.label = '\beta_N';
@@ -56,6 +69,8 @@ switch lower(data_request)
mapping.label = '\beta_p';
mapping.method = 'tdiliuqe';
mapping.expression = '\results::beta_pol';
+ mapping.expression = '\tcv_shot::top.results.equil_1.results:beta_pol';
+ mapping.expression = 'tcv_eq(''''beta_pol'''',''''LIUQE.M'''')';
case 'cxrs'
mapping.timedim = 2;
mapping.label = 'cxrs';
@@ -70,6 +85,8 @@ switch lower(data_request)
mapping.timedim = 1;
mapping.method = 'tdiliuqe';
mapping.expression = '\results::delta_edge';
+ mapping.expression = '\tcv_shot::top.results.equil_1.results:delta_edge';
+ mapping.expression = 'tcv_eq(''''delta_edge'''',''''LIUQE.M'''')';
% mapping.method = 'expression';
% mapping.expression = ['tdi(''\results::q_psi'');']; % for tests
case 'delta_bottom'
@@ -77,11 +94,20 @@ switch lower(data_request)
mapping.label = 'delta\_bottom';
mapping.method = 'tdiliuqe';
mapping.expression = '\results::delta_ed_bot';
+ mapping.expression = '\tcv_shot::top.results.equil_1.results:delta_bot';
+ mapping.expression = 'tcv_eq(''''delta_bot'''',''''LIUQE.M'''')';
+ case 'delta_rho'
+ mapping.timedim = 1;
+ mapping.method = 'tdiliuqe';
+ mapping.expression = '\tcv_shot::top.results.equil_1.results:delta';
+ mapping.expression = 'tcv_eq(''''delta'''',''''LIUQE.M'''')';
case 'delta_top'
mapping.timedim = 1;
mapping.label = 'delta\_top';
mapping.method = 'tdiliuqe';
mapping.expression = '\results::delta_ed_top';
+ mapping.expression = '\tcv_shot::top.results.equil_1.results:delta_top';
+ mapping.expression = 'tcv_eq(''''delta_top'''',''''LIUQE.M'''')';
case 'ece'
mapping.timedim = 2;
mapping.method = 'switchcase';
@@ -105,19 +131,34 @@ switch lower(data_request)
mapping.label = 'I ohmic transformer';
mapping.method = 'tdi';
mapping.expression = [{'\magnetics::ipol[*,$1]'} {'OH_001'}];
- case 'ip'
+ case 'ip_trapeze'
mapping.timedim = 1;
mapping.label = 'Plasma current';
mapping.method = 'tdi';
mapping.expression = '\magnetics::iplasma:trapeze';
+ case 'ip'
+ mapping.timedim = 1;
+ mapping.label = 'Plasma current';
+ mapping.timedim = 1;
+ mapping.method = 'tdiliuqe';
+ mapping.expression = 'tcv_eq(''''i_pl'''',''''LIUQE.M'''')'; % to be able to get ip consistent with relevant LIUQE value
case 'kappa'
mapping.timedim = 1;
mapping.method = 'tdiliuqe';
mapping.expression = '\results::kappa_edge';
+ mapping.expression = '\tcv_shot::top.results.equil_1.results:kappa_edge';
+ mapping.expression = 'tcv_eq(''''kappa_edge'''',''''LIUQE.M'''')';
+ case 'kappa_rho'
+ mapping.timedim = 1;
+ mapping.method = 'tdiliuqe';
+ mapping.expression = '\tcv_shot::top.results.equil_1.results:kappa';
+ mapping.expression = 'tcv_eq(''''kappa'''',''''LIUQE.M'''')';
case 'li'
mapping.timedim = 1;
mapping.method = 'tdiliuqe';
mapping.expression = '\results::l_i';
+ mapping.expression = '\tcv_shot::top.results.equil_1.results:l_i_3';
+ mapping.expression = 'tcv_eq(''''l_i_3'''',''''LIUQE.M'''')';
case 'mhd'
mapping.timedim = 1;
mapping.label = 'n=1,2, etc';
@@ -175,30 +216,49 @@ switch lower(data_request)
'gdat_tmp2=gdat_aug(shot,params_eff);ij=find(gdat_tmp2.data==0);gdat_tmp2.data(ij)=NaN;' ...
'tmp_data2=interp1(gdat_tmp2.t,gdat_tmp2.data,gdat_tmp.t,[],NaN);' ...
'gdat_tmp.data = gdat_tmp.data./(tmp_data2+1e-5);'];
+ case {'phi_tor', 'phitor', 'toroidal_flux'}
+ mapping.timedim = 2;
+ mapping.method = 'tdiliuqe';
+ mapping.expression = '\results::psi_axis';
+ mapping.expression = '\tcv_shot::top.results.equil_1.results:psi_axis';
+ mapping.expression = 'tcv_eq(''''tor_flux_tot'''',''''LIUQE.M'''')';
case 'powers'
mapping.timedim = 1;
mapping.label = 'various powers';
mapping.method = 'switchcase';
- case 'psi_axis'
+ case {'psi_axis', 'psi_mag'}
mapping.timedim = 1;
mapping.method = 'tdiliuqe';
mapping.expression = '\results::psi_axis';
+ mapping.expression = '\tcv_shot::top.results.equil_1.results:psi_axis';
+ mapping.expression = 'tcv_eq(''''psi_axis'''',''''LIUQE.M'''')';
+ mapping.label = 'psi\_axis with psi_edge=0';
case 'psi_edge'
mapping.timedim = 1;
- mapping.method = 'switchcase'; % is set to zero, so not in tree nodes
+ mapping.method = 'switchcase'; % should be psi_edge(t) so to add to 0 of standard LCFS in LIUQE
mapping.label = 'psi\_edge';
+ mapping.method = 'tdiliuqe';
+ mapping.expression = '\results::surface_flux';
+ mapping.expression = '\tcv_shot::top.results.equil_1.results:psi_surf';
+ mapping.expression = 'tcv_eq(''''psi_surf'''',''''LIUQE.M'''')';
case 'q0'
mapping.timedim = 1;
mapping.method = 'tdiliuqe';
mapping.expression = '\results::q_zero';
+ mapping.expression = '\tcv_shot::top.results.equil_1.results:q_axis';
+ mapping.expression = 'tcv_eq(''''q_axis'''',''''LIUQE.M'''')';
case 'q95'
mapping.timedim = 1;
mapping.method = 'tdiliuqe';
- mapping.expression = '\results::q_95';
+ % mapping.expression = '\results::q_95';
+ mapping.expression = '\tcv_shot::top.results.equil_1.results:q_95';
+ mapping.expression = 'tcv_eq(''''q_95'''',''''LIUQE.M'''')';
case 'qedge'
mapping.timedim = 1;
mapping.method = 'tdiliuqe';
mapping.expression = '\results::q_edge';
+ mapping.expression = '\tcv_shot::top.results.equil_1.results:q_edge';
+ mapping.expression = 'tcv_eq(''''q_edge'''',''''LIUQE.M'''')';
case 'q_rho'
mapping.timedim = 2;
mapping.label = 'q';
@@ -207,27 +267,45 @@ switch lower(data_request)
mapping.timedim = 1;
mapping.method = 'tdiliuqe';
mapping.expression = '\results::r_contour';
- case 'rgeom'
+ mapping.expression = '\tcv_shot::top.results.equil_1.results:r_rho'; % several flux surfaces R coordinates (irho,itheta,t)
+ mapping.expression = 'tcv_eq(''''r_rho'''',''''LIUQE.M'''')';
+ case 'r_contour_edge'
+ mapping.timedim = 1;
+ mapping.method = 'tdiliuqe';
+ mapping.expression = '\results::r_contour';
+ mapping.expression = '\tcv_shot::top.results.equil_1.results:r_surf'; % LCFS R coordinates (r,t)
+ mapping.expression = 'tcv_eq(''''r_edge'''',''''LIUQE.M'''')';
+ case {'rgeom', 'r_geom'}
mapping.timedim = 1;
mapping.label = 'Rgeom';
mapping.method = 'switchcase';
+ case {'rgeom_rho', 'r_geom_rho'}
+ mapping.timedim = 2;
+ mapping.label = 'Rgeom';
+ mapping.method = 'switchcase';
case 'rhotor_edge'
mapping.timedim = 1;
mapping.label = 'rhotor\_edge=sqrt(Phi/pi/B0)';
mapping.method = 'switchcase'; % from conf if exist otherwise computes it
- case 'rhotor'
+ case 'rhotor_norm'
mapping.timedim = 2;
mapping.label = 'rhotor\_norm';
mapping.method = 'switchcase'; % from conf if exist otherwise computes it
- case 'rhovol'
+ case 'rhotor'
+ mapping.timedim = 2;
+ mapping.label = 'rhotor';
+ mapping.method = 'switchcase'; % from conf if exist otherwise computes it
+ case {'rhovol', 'rho_vol'}
mapping.timedim = 2;
mapping.label = 'rhovol\_norm';
mapping.method = 'switchcase'; % from conf if exist otherwise computes it
- case 'rmag'
+ case {'rmag', 'r_axis', 'r_mag'}
mapping.timedim = 1;
mapping.label = 'R\_magaxis';
mapping.method = 'tdiliuqe';
mapping.expression = '\results::r_axis';
+ mapping.expression = '\tcv_shot::top.results.equil_1.results:r_axis';
+ mapping.expression = 'tcv_eq(''''r_axis'''',''''LIUQE.M'''')';
case 'rtc'
mapping.timedim = 1;
mapping.label = 'rtc\_signals';
@@ -272,45 +350,57 @@ switch lower(data_request)
mapping.label = 'Volume(\rho)';
mapping.method = 'switchcase';
% mapping.expression = '\results::psitbx:vol'; (if exists for liuqe2 and 3 as well)
- case 'wmhd'
+ case {'wmhd', 'w_mhd'}
mapping.timedim = 1;
mapping.method = 'tdiliuqe';
mapping.expression = '\results::total_energy';
+ mapping.expression = '\tcv_shot::top.results.equil_1.results:w_mhd';
+ mapping.expression = 'tcv_eq(''''w_mhd'''',''''LIUQE.M'''')';
case 'z_contour'
mapping.timedim = 1;
mapping.method = 'tdiliuqe';
mapping.expression = '\results::z_contour';
+ mapping.expression = '\tcv_shot::top.results.equil_1.results:z_rho'; % several flux surfaces Z coordinates (irho,itheta,t)
+ mapping.expression = 'tcv_eq(''''z_rho'''',''''LIUQE.M'''')';
+ case 'z_contour_edge'
+ mapping.timedim = 1;
+ mapping.method = 'tdiliuqe';
+ mapping.expression = '\results::z_contour';
+ mapping.expression = '\tcv_shot::top.results.equil_1.results:z_surf'; % LCFS Z coordinates (r,t)
+ mapping.expression = 'tcv_eq(''''z_edge'''',''''LIUQE.M'''')';
case 'zeff'
mapping.timedim = 1;
mapping.label = 'zeff from Ip-Ibs';
mapping.method = 'tdi';
mapping.expression = '\results::ibs:z_eff';
- case 'zgeom'
+ case {'zgeom', 'z_geom'}
mapping.timedim = 1;
mapping.label = 'Zgeom';
mapping.method = 'switchcase';
- case 'zmag'
+ case {'zmag', 'z_mag', 'z_axis'}
mapping.timedim = 1;
mapping.label = 'Zmagaxis';
mapping.method = 'tdiliuqe';
mapping.expression = '\results::z_axis';
+ mapping.expression = '\tcv_shot::top.results.equil_1.results:z_axis';
+ mapping.expression = 'tcv_eq(''''z_axis'''',''''LIUQE.M'''')';
%
% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% extra TCV cases (not necessarily in official data_request name list)
% %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
- case '\results::total_energy'
- mapping.timedim = 1;
- mapping.gdat_timedim = 1;
- mapping.method = 'tdiliuqe';
- % mapping.expression = '\results::total_energy:foo';
- mapping.expression = '\results::total_energy';
- case '\results::thomson:psiscatvol'
- mapping.timedim = 1;
- mapping.gdat_timedim = 1;
- mapping.method = 'tdiliuqe';
- % mapping.expression = '\results::thomson:psiscatvol:foo';
- mapping.expression = '\results::thomson:psiscatvol';
+% $$$ case '\results::total_energy'
+% $$$ mapping.timedim = 1;
+% $$$ mapping.gdat_timedim = 1;
+% $$$ mapping.method = 'tdiliuqe';
+% $$$ % mapping.expression = '\results::total_energy:foo';
+% $$$ mapping.expression = '\results::total_energy';
+% $$$ case '\results::thomson:psiscatvol'
+% $$$ mapping.timedim = 1;
+% $$$ mapping.gdat_timedim = 1;
+% $$$ mapping.method = 'tdiliuqe';
+% $$$ % mapping.expression = '\results::thomson:psiscatvol:foo';
+% $$$ mapping.expression = '\results::thomson:psiscatvol';
case 'mpx'
mapping.timedim = 1;
mapping.gdat_timedim = 2;
diff --git a/crpptbx/gdat_plot.m b/crpptbx/gdat_plot.m
index 8f8bc1374728c6fcb9961797b5696a3bc1c3ec23..5cdee4d7cf58be6e25c48980f340caa618fa64c4 100644
--- a/crpptbx/gdat_plot.m
+++ b/crpptbx/gdat_plot.m
@@ -61,7 +61,15 @@ if all(isfield(gdat_data,{'data','t'})) && ~isempty(gdat_data.data) && ~isempty(
axis equal;
title(eval(['gdat_data.eqdsk' endstr '.stitle']));
elseif any(find(size(gdat_data.data)==length(gdat_data.t)))
- plot(gdat_data.t,gdat_data.data);
+ try
+ plot(gdat_data.t,gdat_data.data);
+ catch ME
+ if exist('ME','var')
+ disp('Problem in gdat_plot')
+ getReport(ME)
+ end
+ return
+ end
title([gdat_data.gdat_params.machine ' #' num2str(gdat_data.shot)]);
if isfield(gdat_data,'mapping_for')
xlabel(['time [' gdat_data.dimunits{gdat_data.mapping_for.(gdat_data.gdat_params.machine).gdat_timedim} ']']);
diff --git a/crpptbx/liuqefortran2liuqematlab.m b/crpptbx/liuqefortran2liuqematlab.m
new file mode 100644
index 0000000000000000000000000000000000000000..6d3d783ebd00c9bc8abbafb973f58468e3834455
--- /dev/null
+++ b/crpptbx/liuqefortran2liuqematlab.m
@@ -0,0 +1,86 @@
+function liuqe_fortran_matlab = liuqefortran2liuqematlab(varargin);
+%
+% liuqe_fortran_matlab = liuqefortran2liuqematlab(varargin);
+%
+% returns correspondance of node names, when different between fortran (eq_recon) and matlab (equil_1) mds node names
+%
+% liuqe_fortran_matlab{:,1} give liuqe fortran names
+% liuqe_fortran_matlab{:,2} give liuqe matlab names
+%
+% Note that you can have multiple entries (for liuqe fortran since it has less nodes) but the first one should be the closest match (see r_contour/r_rho/r_surf)
+%
+% varargin: If absent, then return full table
+% varargin{1}: 'node_name_to_convert'
+% varargin{2}: origin of node_name given=varargin{1}: 0: fortran, 1(default): matlab
+% varargin{3}: type of desired output: 0: fortran, 1(default): matlab
+%
+% Thus by default there is no conversion
+%
+% Examples:
+% liuqe_fortran_matlab = liuqefortran2liuqematlab; % => liuqe_fortran_matlab(:,1:2) full cell table
+% liuqe_fortran_matlab = liuqefortran2liuqematlab('q_axis'); => liuqe_fortran_matlab = 'q_axis'
+% liuqe_fortran_matlab = liuqefortran2liuqematlab('q_axis',1,0); => liuqe_fortran_matlab = 'q_zero'
+% liuqe_fortran_matlab = liuqefortran2liuqematlab('q_zero',0,1); => liuqe_fortran_matlab = 'q_axis'
+% liuqe_fortran_matlab = liuqefortran2liuqematlab('r_edge',1,0); => liuqe_fortran_matlab = 'r_contour'
+%
+
+liuqe_fortran_matlab_table = [ ...
+ {'l_i'} , {'l_i_3'} ; ...
+ {'i_p'} , {'i_pl'} ; ...
+ {'surface_flux'} , {'psi_surf'} ; ...
+ {'q_zero'} , {'q_axis'} ; ...
+ {'q_psi'} , {'q'} ; ...
+ {'r_contour'} , {'r_edge'} ; ... % r_rho has all the flux surfaces
+ {'r_min_psi'} , {'r_in'} ; ... % R inboard of rho flux surfaces
+ {'r_max_psi'} , {'r_out'} ; ... % R outboard of rho flux surfaces
+ {'total_energy'} , {'w_mhd'} ; ...
+ {'z_contour'} , {'z_edge'} ; ... % z_rho has all the flux surfaces
+ ];
+
+liuqe_fortran_matlab = [];
+if nargin == 0
+ % return full table
+ liuqe_fortran_matlab = liuqe_fortran_matlab_table;
+ return
+end
+
+
+if isempty(varargin{1})
+ error(['liuqefortran2liuqematlab: unexpected empty 1st argument'])
+ return
+end
+if ~ischar(varargin{1})
+ error(['liuqefortran2liuqematlab: unexpected 1st argument is not a string'])
+ return
+end
+
+liuqe_matlab_in = 1;
+if nargin>=2 && ~isempty(varargin{2})
+ if isnumeric(varargin{2})
+ liuqe_matlab_in = varargin{2};
+ else
+ warning(['liuqefortran2liuqematlab: unexpected 2nd argument type, should be numeric'])
+ varargin{2}
+ return
+ end
+end
+
+liuqe_matlab_out = 1;
+if nargin>=3 && ~isempty(varargin{3})
+ if isnumeric(varargin{3})
+ liuqe_matlab_out = varargin{3};
+ else
+ warning(['liuqefortran2liuqematlab: unexpected 3rd argument type, should be numeric'])
+ varargin{3}
+ return
+ end
+end
+
+% find index of input in corresponding column
+ij = strmatch(varargin{1},liuqe_fortran_matlab_table(:,1+liuqe_matlab_in),'exact');
+if isempty(ij)
+ % assume name is the same
+ liuqe_fortran_matlab = varargin{1};
+else
+ liuqe_fortran_matlab = liuqe_fortran_matlab_table{ij(1),1+liuqe_matlab_out};
+end