diff --git a/matlab/CHDF/cdf2mat.m b/matlab/CHDF/cdf2mat.m
index c9d7a14e693f33ce81ffa9f00d2d501cb6eb6fcc..f2c1155bd9ac11239d13c77fc09e08f60c2a946f 100644
--- a/matlab/CHDF/cdf2mat.m
+++ b/matlab/CHDF/cdf2mat.m
@@ -1,4 +1,5 @@
-function cdf2mat_out = cdf2mat(pfname)
+function cdf2mat_out = cdf2mat(pfname,varargin)
+% cdf2mat_out = cdf2mat(pfname,varargin)
%
% reads all variables and coordinates from netcdf
% some trials with netcdf, using 50725c01.cdf as test
@@ -7,6 +8,9 @@ function cdf2mat_out = cdf2mat(pfname)
%
% Uses matlab netcdf.open, ncinfo, netcdf.inqVarIDs, netcd.getVar, etc
%
+% varargin{1}: shot number to add (in particular when cdf file does not have it defined
+% shot from file: allinfo=ncinfo(pfname);if ~isempty(allinfo.Attributes),top_attr_names = {allinfo.Attributes(:).Name};allinfo.Attributes(strmatch('shot',top_attr_names,'exact')).Value, else, disp('no Attributes'),end
+%
%
if ~exist(pfname,'file')
@@ -33,6 +37,7 @@ if length(allvarnames) ~= length(allvarids)
end
[varnames_sorted,~]=sort(allvarnames);
+cdf2mat_out.allvarnames_sorted = varnames_sorted;
% to find a variable:
% strmatch('GFUN',allvarnames,'exact')
@@ -56,6 +61,9 @@ for i=1:length(coordnames_sorted)
for ij=1:length(fields_variables_to_copy)
matcdf.coords(i).(fields_variables_to_copy{ij}) = allinfo.Variables(matcdf.coords(i).index_allvarnames).(fields_variables_to_copy{ij});
end
+ % define defaults before fetching the values in Attributes (if provided)
+ matcdf.coords(i).units = '';
+ matcdf.coords(i).long_name = [matcdf.coords(i).name ' empty'];
for jj=1:numel(allinfo.Variables(matcdf.coords(i).index_allvarnames).Attributes)
if strcmp(lower(allinfo.Variables(matcdf.coords(i).index_allvarnames).Attributes(jj).Name),'units')
matcdf.coords(i).units = strtrim(allinfo.Variables(matcdf.coords(i).index_allvarnames).Attributes(jj).Value);
@@ -80,6 +88,7 @@ for i=1:length(varnames_sorted)
matcdf.vars(i).name = varnames_sorted{i};
matcdf.vars(i).index_allvarnames = strmatch(matcdf.vars(i).name,allvarnames,'exact');
matcdf.vars(i).varid = allvarids(matcdf.vars(i).index_allvarnames);
+ % if i==strmatch('ZEFFC',varnames_sorted,'exact'), keyboard; end
if ~isempty(matcdf.vars(i).index_allvarnames)
if strcmp(allinfo.Variables(matcdf.vars(i).index_allvarnames).Datatype,'single')
matcdf.vars(i).data = netcdf.getVar(funnetcdf,matcdf.vars(i).varid,'double');
@@ -89,6 +98,9 @@ for i=1:length(varnames_sorted)
for ij=1:length(fields_variables_to_copy)
matcdf.vars(i).(fields_variables_to_copy{ij}) = allinfo.Variables(matcdf.vars(i).index_allvarnames).(fields_variables_to_copy{ij});
end
+ % define defaults before fetching the values in Attributes (if provided)
+ matcdf.vars(i).units = '';
+ matcdf.vars(i).long_name = [matcdf.vars(i).name ' empty'];
for jj=1:numel(allinfo.Variables(matcdf.vars(i).index_allvarnames).Attributes)
if strcmp(lower(allinfo.Variables(matcdf.vars(i).index_allvarnames).Attributes(jj).Name),'units')
matcdf.vars(i).units = strtrim(allinfo.Variables(matcdf.vars(i).index_allvarnames).Attributes(jj).Value);
@@ -126,9 +138,17 @@ end
if ~isempty(allinfo.Attributes)
top_attr_names = {allinfo.Attributes(:).Name};
ij = strmatch('shot',top_attr_names,'exact');
+else
+ ij = [];
+end
+if ~isempty(ij)
cdf2mat_out.shot = allinfo.Attributes(ij).Value;
else
- cdf2mat_out.shot = NaN;
+ if nargin>1 && isnumeric(varargin{1})
+ cdf2mat_out.shot = varargin{1};
+ else
+ cdf2mat_out.shot = NaN;
+ end
end
cdf2mat_out.fname = pfname;
[a1,a2,a3]=fileparts(pfname);
diff --git a/matlab/D3D/d3d_requests_mapping.m b/matlab/D3D/d3d_requests_mapping.m
index 90695505311c8f4f0cb377308de15cf7d77c9708..a40639a32184a880c4c2f6ec16157a9cfd1ad119 100644
--- a/matlab/D3D/d3d_requests_mapping.m
+++ b/matlab/D3D/d3d_requests_mapping.m
@@ -224,6 +224,11 @@ switch lower(data_request)
mapping.label = 'NEBAR\_R0';
mapping.method = 'signal';
mapping.expression = [{'efit01'},{'\NEBAR_R0'}];
+ case 'nel_v1'
+ mapping.timedim = 1;
+ mapping.label = 'NEBAR\_V1';
+ mapping.method = 'signal';
+ mapping.expression = [{'efit01'},{'\NEBAR_V1'}];
case {'ne', 'ne_rho'}
mapping.timedim = 2;
mapping.label = 'ne';
diff --git a/matlab/D3D/gdat_d3d.m b/matlab/D3D/gdat_d3d.m
index c3b5c90df0453e60fb371d841279af6b0afd3fe1..dd7414586065ab67e4cd6f2ce9d7f9cefeb96c1f 100644
--- a/matlab/D3D/gdat_d3d.m
+++ b/matlab/D3D/gdat_d3d.m
@@ -49,6 +49,7 @@ function [gdat_data,gdat_params,error_status,varargout] = gdat_d3d(shot,data_req
% a5=gdat(32827,'ip'); % standard call
% a6=gdat(32827,'ip','Opt1',123,'Doplot',1,'opt2','Abc'); % standard call with a few options (note all lowercase in output)
%
+% From local: mdsconnect('atlas.gat.com'); mdsvalue('1+3')
% From remote:
% ssh -p 2039 -C -l sautero -L 8002:atlas.gat.com:8000 cybele.gat.com
% And in another session in matlab:
@@ -715,6 +716,12 @@ elseif strcmp(mapping_for_d3d.method,'switchcase')
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
case {'equil'}
+ if isfield(gdat_data.gdat_params,'source') && ~isempty(gdat_data.gdat_params.source) && ischar(gdat_data.gdat_params.source) && strcmp(lower(gdat_data.gdat_params.equil),'efit01')
+ % assume source provided instead of equil for choice of equil, since equil=default and source provided
+ warning('choice of equil should be provided in ''equil'' not in ''source'', moved entry')
+ gdat_data.gdat_params.equil = gdat_data.gdat_params.source;
+ gdat_data.gdat_params = rmfield(gdat_data.gdat_params,'source');
+ end
if ~isfield(gdat_data.gdat_params,'equil') || isempty(gdat_data.gdat_params.equil) || ~ischar(gdat_data.gdat_params.equil)
gdat_data.gdat_params.equil = 'efit03';
else
@@ -844,9 +851,9 @@ elseif strcmp(mapping_for_d3d.method,'switchcase')
if any(strcmp(fieldnames(a),'units'))
gdat_data.(extra_source{i}).units=a.units;
end
- gdat_data.(extra_source{i}).r = mdsdata([main_source 'r']);
- gdat_data.(extra_source{i}).z = mdsdata([main_source 'z']);
- gdat_data.(extra_source{i}).error_bar = mdsdata([nodenameeff '_e'])';
+ gdat_data.(extra_source{i}).r = mdsvalue([main_source 'r']);
+ gdat_data.(extra_source{i}).z = mdsvalue([main_source 'z']);
+ gdat_data.(extra_source{i}).error_bar = mdsvalue([nodenameeff '_e'])';
gdat_data.(extra_source{i}).data_fullpath=[data_request_eff 'from electrons ' nodenameeff];
gdat_data.(extra_source{i}).dim=[{gdat_data.x};{gdat_data.t}];
if strcmp(extra_source{i},'tangential')
@@ -875,7 +882,10 @@ elseif strcmp(mapping_for_d3d.method,'switchcase')
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
case {'neint_ts'}
if ~isfield(gdat_data.gdat_params,'source') || isempty(gdat_data.gdat_params.source) || ~isstruct(gdat_data.gdat_params.source)
- error('requires ''ne_rho'' gdat structure as input in ''source'' option');
+ warning('requires ''ne_rho'' gdat structure as input in ''source'' option. Loading it...');
+ params_eff = gdat_data.gdat_params;
+ params_eff.data_request = 'ne_rho';
+ [gdat_data.gdat_params.source,params_nerho,error_status]=gdat_d3d(shot,params_eff);
end
ne_rho = gdat_data.gdat_params.source;
% Computes R0 and V1 line integrated signal based on ne fit (so inside LCFS), provide also R(psi_ts) and Z(psi_ts) on these lines
@@ -1049,11 +1059,16 @@ elseif strcmp(mapping_for_d3d.method,'switchcase')
return
end
% add rho coordinates
- params_eff.data_request='equil';
- [equil,params_equil,error_status]=gdat_d3d(shot,params_eff);
- if error_status>0
- if gdat_params.nverbose>=3; disp(['problems with ' params_eff.data_request]); end
- return
+ if ~isfield(gdat_data.gdat_params,'equil') || isempty(gdat_data.gdat_params.equil) || ~isstruct(gdat_data.gdat_params.equil)
+ params_eff.data_request='equil';
+ [equil,params_equil,error_status]=gdat_d3d(shot,params_eff);
+ if error_status>0
+ if gdat_params.nverbose>=3; disp(['problems with ' params_eff.data_request]); end
+ return
+ end
+ else
+ params_equil.equil = gdat_data.gdat_params.equil.gdat_params.equil;
+ equil = gdat_data.gdat_params.equil;
end
gdat_data.gdat_params.equil = params_equil.equil;
gdat_data.equil = equil;
@@ -1071,24 +1086,26 @@ elseif strcmp(mapping_for_d3d.method,'switchcase')
rhotornorm_out_source = NaN*ones(inb_chord_source,inb_time_source);
rhovolnorm_out_source = NaN*ones(inb_chord_source,inb_time_source);
% constructs intervals within which a given equil is used: [time_equil(i),time_equil(i+1)]
- time_equil=[min(gdat_data.(extra_source{i}).t(1)-0.1,equil.t(1)-0.1); 0.5.*(equil.t(1:end-1)+equil.t(2:end)) ;max(equil.t(end)+0.1,gdat_data.(extra_source{i}).t(end)+0.1)];
- clear psi_out_core rhopolnorm_out_core rhotornorm_out_core rhovolnorm_out_core
- for itequil=1:numel(time_equil)-1
- rr=equil.Rmesh(:,itequil);
- zz=equil.Zmesh(:,itequil);
- psirz_in = equil.psi2D(:,:,itequil);
- it_core_inequil = find(gdat_data.(extra_source{i}).t>time_equil(itequil) & gdat_data.(extra_source{i}).t<=time_equil(itequil+1));
- if ~isempty(it_core_inequil)
- rout_core=gdat_data.(extra_source{i}).r;
- zout_core=gdat_data.(extra_source{i}).z;
- psi_at_routzout = interpos2Dcartesian(rr,zz,psirz_in,rout_core,zout_core);
- psi_out_core(:,it_core_inequil) = repmat(psi_at_routzout,1,numel(it_core_inequil));
- rhopolnorm_out_core(:,it_core_inequil) = sqrt((psi_out_core(:,it_core_inequil)-equil.psi_axis(itequil))./(equil.psi_lcfs(itequil)-equil.psi_axis(itequil)));
- for it_cx=1:length(it_core_inequil)
- rhotornorm_out_core(:,it_core_inequil(it_cx)) = ...
- interpos(equil.rhopolnorm(:,itequil),equil.rhotornorm(:,itequil),rhopolnorm_out_core(:,it_core_inequil(it_cx)),-3,[2 2],[0 1]);
- rhovolnorm_out_core(:,it_core_inequil(it_cx)) = ...
- interpos(equil.rhopolnorm(:,itequil),equil.rhovolnorm(:,itequil),rhopolnorm_out_core(:,it_core_inequil(it_cx)),-3,[2 2],[0 1]);
+ if ~isempty(gdat_data.(extra_source{i}).t)
+ time_equil=[min(gdat_data.(extra_source{i}).t(1)-0.1,equil.t(1)-0.1); 0.5.*(equil.t(1:end-1)+equil.t(2:end)) ;max(equil.t(end)+0.1,gdat_data.(extra_source{i}).t(end)+0.1)];
+ clear psi_out_core rhopolnorm_out_core rhotornorm_out_core rhovolnorm_out_core
+ for itequil=1:numel(time_equil)-1
+ rr=equil.Rmesh(:,itequil);
+ zz=equil.Zmesh(:,itequil);
+ psirz_in = equil.psi2D(:,:,itequil);
+ it_core_inequil = find(gdat_data.(extra_source{i}).t>time_equil(itequil) & gdat_data.(extra_source{i}).t<=time_equil(itequil+1));
+ if ~isempty(it_core_inequil)
+ rout_core=gdat_data.(extra_source{i}).r;
+ zout_core=gdat_data.(extra_source{i}).z;
+ psi_at_routzout = interpos2Dcartesian(rr,zz,psirz_in,rout_core,zout_core);
+ psi_out_core(:,it_core_inequil) = repmat(psi_at_routzout,1,numel(it_core_inequil));
+ rhopolnorm_out_core(:,it_core_inequil) = sqrt((psi_out_core(:,it_core_inequil)-equil.psi_axis(itequil))./(equil.psi_lcfs(itequil)-equil.psi_axis(itequil)));
+ for it_cx=1:length(it_core_inequil)
+ rhotornorm_out_core(:,it_core_inequil(it_cx)) = ...
+ interpos(equil.rhopolnorm(:,itequil),equil.rhotornorm(:,itequil),rhopolnorm_out_core(:,it_core_inequil(it_cx)),-3,[2 2],[0 1]);
+ rhovolnorm_out_core(:,it_core_inequil(it_cx)) = ...
+ interpos(equil.rhopolnorm(:,itequil),equil.rhovolnorm(:,itequil),rhopolnorm_out_core(:,it_core_inequil(it_cx)),-3,[2 2],[0 1]);
+ end
end
end
end
diff --git a/matlab/IMAS/gdat_imas.m b/matlab/IMAS/gdat_imas.m
index 88b72b2f647f2f8ca0bb8b09601439ac62e11e9d..d2b06fefcbdf8d2fcc5fc3b3289cc933cdc28ddd 100644
--- a/matlab/IMAS/gdat_imas.m
+++ b/matlab/IMAS/gdat_imas.m
@@ -633,7 +633,7 @@ end
function [gdat_data] = get_thomson_raw_data(shot,data_request_eff,gdat_data,doedge,nverbose);
%
try
- time=mdsdata('\results::thomson:times');
+ time=mdsvalue('\results::thomson:times');
catch
gdat_data.error_bar = [];
if strcmp(data_request_eff(1:2),'ne')
@@ -695,7 +695,7 @@ if strcmp(data_request_eff(1:2),'ne')
disp('***********************************************************************')
end
end
-z=mdsdata(['\diagz::thomson_set_up' edge_str_dot ':vertical_pos']);
+z=mdsvalue(['\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;
@@ -716,7 +716,7 @@ function [psiscatvol,psi_max] = psi_scatvol_zshift(gdat_data, zshift, system_str
% 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_z = mdsvalue(['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
diff --git a/matlab/TCV/gdat_tcv.m b/matlab/TCV/gdat_tcv.m
index 17933ceefd2a2c134beca063da92dbee017036f6..8acfc6f38c78d0457b3e330c4ed1c09b7b2be85e 100644
--- a/matlab/TCV/gdat_tcv.m
+++ b/matlab/TCV/gdat_tcv.m
@@ -866,6 +866,94 @@ elseif strcmp(mapping_for_tcv.method,'switchcase')
gdat_data.gdat_params.time_out = [];
end
+ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+ case {'bfields'}
+ params_eff = gdat_data.gdat_params;
+
+ try
+ if params_eff.liuqe == 1
+ % load LIUQE data
+ [L,LY] = mds2meq(gdat_data.shot,'LIUQE.M',[],'ifield',true); % flag to allow field computation
+ LY = meqreprocess(L,LY); % compute fields etc, since not directly loaded from mds2meq
+ datapath = 'Outputs from [L,LY]=mds2meq(gdat_data.shot,''LIUQE.M'',[],''ifield'',true); LY = meqreprocess(L,LY);';
+ Br = permute(LY.Brx,[2,1,3]); % permute for (Z,R)->(R,Z)
+ Bz = permute(LY.Bzx,[2,1,3]);
+ Btor = permute(LY.Btx,[2,1,3]);
+ dim_in = {L.rx,L.zx,LY.t};
+ time_in = LY.t;
+ x_in = {L.rx,L.zx};
+ else % compute B-fields from psi
+ ip_gdat = gdat(gdat_data.shot,'ip');
+ [B,coord,time,psi,fsd] = tcv_mag_field (gdat_data.shot,ip_gdat.t,[],[],[],psitbx_str);
+ datapath ='Outputs from [B,coord,time,psi,fsd] = tcv_mag_field(shot,time,[],[],[],psitbx_str);';
+ Br = B.R;
+ Bz = B.z;
+ Btor = B.phi;
+ dim_in = {coord.R,coord.z,time};
+ time_in = time;
+ x_in = {coord.R,coord.z};
+ end
+
+ gdat_data.data_fullpath = datapath;
+
+ % radial magnetic field
+ gdat_data.Br.data = Br;
+ gdat_data.Br.units = 'T';
+ gdat_data.Br.dim = dim_in;
+ gdat_data.Br.dimunits = {'m','m','s'};
+ gdat_data.Br.t = time_in;
+ gdat_data.Br.x = x_in;
+ gdat_data.Br.data_fullpath = datapath;
+ gdat_data.Br.label = 'Radial magnetic field map in (R,Z)';
+
+ % vertical magnetic field
+ gdat_data.Bz.data = Bz;
+ gdat_data.Bz.units = 'T';
+ gdat_data.Bz.dim = dim_in;
+ gdat_data.Bz.dimunits = {'m','m','s'};
+ gdat_data.Bz.t = time_in;
+ gdat_data.Bz.x = x_in;
+ gdat_data.Bz.data_fullpath = datapath;
+ gdat_data.Bz.label = 'Vertical magnetic field map in (R,Z)';
+
+ % toroidal magnetic field
+ gdat_data.Btor.data = Btor;
+ gdat_data.Btor.units = 'T';
+ gdat_data.Btor.dim = dim_in;
+ gdat_data.Btor.dimunits = {'m','m','s'};
+ gdat_data.Btor.t = time_in;
+ gdat_data.Btor.x = x_in;
+ gdat_data.Btor.data_fullpath = datapath;
+ gdat_data.Btor.label = 'Toroidal magnetic field map in (R,Z)';
+
+ % total magnetic field
+ gdat_data.Btot.data = sqrt(gdat_data.Btor.data.^2 + gdat_data.Bz.data.^2 + gdat_data.Br.data.^2);
+ gdat_data.Btot.units = 'T';
+ gdat_data.Btot.dim = dim_in;
+ gdat_data.Btot.dimunits = {'m','m','s'};
+ gdat_data.Btot.t = time_in;
+ gdat_data.Btot.x = x_in;
+ gdat_data.Btot.data_fullpath = datapath;
+ gdat_data.Btot.label = 'Total magnetic field map in (R,Z)';
+
+ gdat_data.data = gdat_data.Btot.data;
+ gdat_data.units = 'T';
+ gdat_data.dim = dim_in;
+ gdat_data.dimunits = {'m','m','s'};
+ gdat_data.t = time_in;
+ gdat_data.x = x_in;
+ gdat_data.label = 'Total magnetic field map in (R,Z)';
+
+ catch
+
+ warning('Problem obtaining B-fields, check if requested nodes are filled.')
+ gdat_data.Btor = [];
+ gdat_data.Br = [];
+ gdat_data.Bz = [];
+ gdat_data.Btot = [];
+
+ end
+
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
case {'betan', 'beta_tor_norm'}
% 100*beta / |Ip[MA] * B0[T]| * a[m]
@@ -2076,33 +2164,37 @@ elseif strcmp(mapping_for_tcv.method,'switchcase')
n3=tdi('abs(mhdmode("LFS",3,1))');
gdat_data.data_fullpath='abs(mhdmode("LFS",n,1));% for 1, 2, 3';
else
+ % handle source input
if isfield(gdat_data.gdat_params,'source') && ~isempty(gdat_data.gdat_params.source)
if ~any(contains(mhd_source_list,lower(gdat_data.gdat_params.source)))
t1=sprintf('bad source for mhd: ''%s'', should be either ',gdat_data.gdat_params.source);
t2=sprintf('''%s'', ',mhd_source_list{1:end-1});
t3=sprintf('or ''%s''',mhd_source_list{end});
- disp(sprintf('%s %s %s',t1,t2,t3));
+ fprintf('%s %s %s\n',t1,t2,t3);
return
end
else
+ % if no source input decide source based on average z_axis
z_axis=gdat_tcv([],'z_axis');
z_axis_av = 0.;
- if numel(z_axis.data > 10) && isnumeric(z_axis.data)
- z_axis_av = nanmean(z_axis.data([round(numel(z_axis.data)/3):round(numel(z_axis.data)*0.85)]));
+ if numel(z_axis.data) > 10 && isnumeric(z_axis.data)
+ z_axis_av = nanmean(z_axis.data(round(numel(z_axis.data)/3):round(numel(z_axis.data)*0.85)));
end
if z_axis_av > 0.12
- gdat_data.gdat_params.source = '23';
+ gdat_data.gdat_params.source = '23'; % probe array at z = +23cm
elseif z_axis_av < -0.12
- gdat_data.gdat_params.source = '-23';
+ gdat_data.gdat_params.source = '-23'; % probe array at z = -23cm
else
- gdat_data.gdat_params.source = '0';
+ gdat_data.gdat_params.source = '0'; % probe array at z = 0cm
end
t1=sprintf('source set to ''%s'', can be ',gdat_data.gdat_params.source);
t2=sprintf('''%s'', ',mhd_source_list{1:end-1});
t3=sprintf('or ''%s''',mhd_source_list{end});
- disp(sprintf('%s %s %s',t1,t2,t3))
+ fprintf('%s %s %s',t1,t2,t3)
end
- if length(gdat_data.gdat_params.source)>=2 && strcmp(gdat_data.gdat_params.source(1:2),'23')
+
+ % load data dependent on requested source
+ if numel(gdat_data.gdat_params.source)>=2 && strcmp(gdat_data.gdat_params.source(1:2),'23')
aaLFSz23_sect3=tdi('\atlas::DT196_MHD_001:channel_067');
aaLFSz23_sect11=tdi('\atlas::DT196_MHD_001:channel_075');
n1 = aaLFSz23_sect3;
@@ -2135,14 +2227,14 @@ elseif strcmp(mapping_for_tcv.method,'switchcase')
% sect 11 180deg from sec 3
aaHFSz0_sect3=tdi('\atlas::DT196_MHD_002:channel_034');
aaHFSz0_sect11=tdi('\atlas::DT196_MHD_002:channel_038');
- gdat_data.n1_HFS=aaHFSz0_sect11;
+ gdat_data.n1_HFS=aaHFSz0_sect3;
gdat_data.n1_HFS.data = gdat_data.n1_HFS.data - aaHFSz0_sect11.data;
- gdat_data.n2_HFS=aaHFSz0_sect11;
+ gdat_data.n2_HFS=aaHFSz0_sect3;
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
- elseif length(gdat_data.gdat_params.source)>=3 && strcmp(gdat_data.gdat_params.source(1:3),'-23')
+ elseif numel(gdat_data.gdat_params.source)>=3 && strcmp(gdat_data.gdat_params.source(1:3),'-23')
aaLFSzm23_sect3=tdi('\atlas::DT196_MHD_002:channel_003');
aaLFSzm23_sect3.data = - aaLFSzm23_sect3.data; % opposite polarity
aaLFSzm23_sect11=tdi('\atlas::DT196_MHD_002:channel_011');
@@ -2153,7 +2245,7 @@ elseif strcmp(mapping_for_tcv.method,'switchcase')
% n3=n1;
gdat_data.data_fullpath=['\atlas::DT196_MHD_002:channel_003 -+ \atlas::DT196_MHD_002:channel_011 for n=1,2,' ...
'LFS_sect_3/11, z=-23cm'];
- elseif strcmp(lower(gdat_data.gdat_params.source(1:4)),'ltcc')
+ elseif strcmpi(gdat_data.gdat_params.source(1:4),'ltcc')
switch lower(gdat_data.gdat_params.source)
case {'ltcc', 'ltcc-dbpol'}
% since other probes measure dbpol make "ltcc" default to ltcc-dbpol
@@ -2166,8 +2258,8 @@ elseif strcmp(mapping_for_tcv.method,'switchcase')
aaLTCC=tdi('\atlas::dt132_ltcc_001:channel_003');
gdat_data.data_fullpath='ltcc-dbtor: \atlas::dt132_ltcc_001:channel_003';
otherwise
- error(sprintf('case %s not known, choose between ltcc-dbpol(=ltcc), ltcc-dbrad or ltcc-dbtor', ...
- lower(gdat_data.gdat_params.source)));
+ error('case %s not known, choose between ltcc-dbpol(=ltcc), ltcc-dbrad or ltcc-dbtor', ...
+ lower(gdat_data.gdat_params.source));
end
n1 = aaLTCC;
n2.data = [];
diff --git a/matlab/TCV/tcv_requests_mapping.m b/matlab/TCV/tcv_requests_mapping.m
index b68d12da31544a73b3c470c6ccea533a60654b3c..cbae8c287f7114b0b333059ef5fa974a874b3cb8 100644
--- a/matlab/TCV/tcv_requests_mapping.m
+++ b/matlab/TCV/tcv_requests_mapping.m
@@ -69,6 +69,11 @@ switch lower(data_request)
mapping.label = 'B_0';
mapping.method = 'switchcase';
mapping.expression = '';
+ case 'bfields'
+ mapping.timedim = 1;
+ mapping.label = 'bfields';
+ mapping.method = 'switchcase';
+ mapping.expression = '';
case {'beta', 'beta_tor'}
mapping.timedim = 1;
mapping.label = '\beta';
@@ -251,11 +256,21 @@ switch lower(data_request)
mapping.label = 'line integrated el. density';
mapping.method = 'tdi';
mapping.expression = '\results::fir:lin_int_dens';
+ case 'neints'
+ mapping.timedim = 1;
+ mapping.label = 'Array line integrated el. density';
+ mapping.method = 'tdi';
+ mapping.expression = '\tcv_shot::top.results.fir.fir_array:lin_int_dens';
case 'nel'
mapping.timedim = 1;
mapping.label = 'line-averaged el. density';
mapping.method = 'switchcase';
mapping.expression = '\results::fir:n_average';
+ case 'nels'
+ mapping.timedim = 1;
+ mapping.label = 'Array line-averaged el. density';
+ mapping.method = 'tdi';
+ mapping.expression = '\tcv_shot::top.results.fir.fir_array:n_average';
case 'ne_rho'
mapping.timedim = 2;
mapping.label = 'ne';
diff --git a/matlab/TCV_IMAS/tcv_get_ids_bolometer.m b/matlab/TCV_IMAS/tcv_get_ids_bolometer.m
index 3ad56e7d681177661bb1394a35baf7a07b538f48..697e49f28b6a635edf44127e54645d0f09f26d67 100644
--- a/matlab/TCV_IMAS/tcv_get_ids_bolometer.m
+++ b/matlab/TCV_IMAS/tcv_get_ids_bolometer.m
@@ -52,13 +52,11 @@ if status
for ii = 1:nchannel
ids_bolometer.channel{ii}.name = bolo_geom.channel{ii}.name;
-
if strcmp(aa(1),'3')
ids_bolometer.channel{ii}.identifier = bolo_geom.channel{ii}.identifier;
elseif strcmp(aa(1),'4')
ids_bolometer.channel{ii}.description = bolo_geom.channel{ii}.identifier;
end
-
ids_bolometer.channel{ii}.detector.geometry_type= bolo_geom.channel{ii}.detector.geometry_type;
ids_bolometer.channel{ii}.detector.centre.phi= bolo_geom.channel{ii}.detector.centre.phi;
ids_bolometer.channel{ii}.detector.centre.r= bolo_geom.channel{ii}.detector.centre.r;
diff --git a/matlab/TCV_IMAS/tcv_get_ids_ec_launchers.m b/matlab/TCV_IMAS/tcv_get_ids_ec_launchers.m
index e6f3b040e8ce238e03ed89f3baaeb51aebaa838a..56a9a5008965f260f6eb21fbc8ac27d5f3bc67aa 100644
--- a/matlab/TCV_IMAS/tcv_get_ids_ec_launchers.m
+++ b/matlab/TCV_IMAS/tcv_get_ids_ec_launchers.m
@@ -76,9 +76,6 @@ for iant=1:nb_launchers
% non time-dependent quantities, take 1st ok values
ids_ec_launchers.beam{iant}.frequency.time = [pow.ec.t(1) pow.ec.t(end)];
ids_ec_launchers.beam{iant}.frequency.data =[launch_params{iant}{it_ok{iant}(1)}.freq launch_params{iant}{it_ok{iant}(end)}.freq];
- % ids_ec_launchers.beam{iant}.mode.time = [pow.ec.t(1) pow.ec.t(end)];
- % ids_ec_launchers.beam{iant}.mode.data = [-1 -1]; % at this stage assume X mode always, to change when available
- % ids_ec_launchers.beam{iant}.time = [pow.ec.t(1) pow.ec.t(end)];
ids_ec_launchers.beam{iant}.mode = -1; % at this stage assume X mode always, to change when available
for i=1:length(it_ok{iant})
r0 = sqrt(launch_params{iant}{it_ok{iant}(i)}.x0.^2 + launch_params{iant}{it_ok{iant}(i)}.y0.^2) / 100.; % in [m]
@@ -87,15 +84,17 @@ for iant=1:nb_launchers
ids_ec_launchers.beam{iant}.launching_position.z(i) = launch_params{iant}{it_ok{iant}(i)}.z0/100.;
ids_ec_launchers.beam{iant}.launching_position.phi(i) = atan2(launch_params{iant}{it_ok{iant}(i)}.y0/100,r0);
ids_ec_launchers.beam{iant}.time(i) = time_launch;
- kz = cos(launch_params{iant}{it_ok{iant}(i)}.theta_toray*pi/180.);
- kmr = -sin(launch_params{iant}{it_ok{iant}(i)}.theta_toray*pi/180.).*cos(launch_params{iant}{it_ok{iant}(i)}.phi_toray*pi/180.);
- kphi = sin(launch_params{iant}{it_ok{iant}(i)}.theta_toray*pi/180.).*sin(launch_params{iant}{it_ok{iant}(i)}.phi_toray*pi/180.); %*sigma_Rphiz (=+1 for TCV cocos=17)
- if (kz==0 && kmr==0)
- ids_ec_launchers.beam{iant}.steering_angle_pol(i) = 0.;
- else
- ids_ec_launchers.beam{iant}.steering_angle_pol(i) = atan2(-kz,kmr);
- end
- ids_ec_launchers.beam{iant}.steering_angle_tor(i) = asin(kphi);
+
+ % define angles based on TORAY angles, as discussed with F. Poli Oct 2024
+ % Reminder: theta_toray, zero on z-axis, pointing inwards with 90deg
+ % phi_toray, 0deg=outwards radial, 180deg=pointing
+ % towards core, +: clockwise, -: counter-clockwise
+ % IDS angles have same convention as TORBEAM and GRAY
+ pol_angle_ids = launch_params{iant}{it_ok{iant}(i)}.theta_toray-90;
+ tor_angle_ids = launch_params{iant}{it_ok{iant}(i)}.phi_toray-180;
+
+ ids_ec_launchers.beam{iant}.steering_angle_pol(i) = pol_angle_ids*pi/180;
+ ids_ec_launchers.beam{iant}.steering_angle_tor(i) = tor_angle_ids*pi/180;
ids_ec_launchers.beam{iant}.spot.size(1,i) = 0.023;
ids_ec_launchers.beam{iant}.spot.size(2,i) = 0.012;
ids_ec_launchers.beam{iant}.spot.angle(i) = 0.0;
diff --git a/matlab/TCV_IMAS/tcv_get_ids_equilibrium.m b/matlab/TCV_IMAS/tcv_get_ids_equilibrium.m
index b0b22e2d593c954f9091642d07135f68c40c2588..102ff032e284243b7eea96afa8d26d828f8f69c9 100644
--- a/matlab/TCV_IMAS/tcv_get_ids_equilibrium.m
+++ b/matlab/TCV_IMAS/tcv_get_ids_equilibrium.m
@@ -5,7 +5,7 @@ function [ids_equilibrium,ids_equilibrium_description,varargout] = ...
% tcv_get_ids_equilibrium(shot,ids_equil_empty,gdat_params,varargin);
%
%
-% gdat_params: gdat_data.gdat_params to get all params passed from original call,
+% gdat_params: gdat_data.gdat_params to get all params passed from original call,
% in particular error_bar and cocos_out options
%
@@ -527,19 +527,31 @@ profiles_2d.r.data = repmat(repmat(profiles_2d.psi.dim{1},1,numel(profiles_2d.ps
profiles_2d_desc.r = 'from dim{1} of ''psi'' repeated';
profiles_2d.z.data = repmat(repmat(profiles_2d.psi.dim{2}',numel(profiles_2d.psi.dim{1}),1),1,1,numel(profiles_2d.psi.dim{3}));
profiles_2d_desc.z = 'from dim{2} of ''psi'' repeated';
+params_eff.data_request = 'bfields';
+bfields_gdat = gdat(params_equilibrium.shot,params_eff);
+if ~isempty(bfields_gdat.Br)
+ profiles_2d.b_field_r = bfields_gdat.Br;
+ profiles_2d_desc.b_field_r = bfields_gdat.Br.label;
+ profiles_2d.b_field_z= bfields_gdat.Bz;
+ profiles_2d_desc.b_field_z = bfields_gdat.Bz.label;
+ profiles_2d.b_field_tor= bfields_gdat.Btor;
+ profiles_2d_desc.b_field_tor = bfields_gdat.Btor.label;
+else
+ warning('B-fields could not be load from gdat(shot,''bfields''). B-fields only available for LIUQE.M, trial 1.')
+end
% theta = gdat(params_equilibrium.shot,'theta','machine',gdat_params.machine);
profiles_2d_fieldnames = fieldnames(profiles_2d);
special_fields = {'grid', 'grid_type'}; % fields needing non-automatic treatments
for it=1:ntime
- for i=1:numel(profiles_2d_fieldnames)
- if ~any(strcmp(profiles_2d_fieldnames{i},special_fields))
- if ~isstruct(ids_equilibrium.time_slice{it}.profiles_2d{1}.(profiles_2d_fieldnames{i}))
- if ~ischar(profiles_2d.(profiles_2d_fieldnames{i}).data) && ~isempty(profiles_2d.(profiles_2d_fieldnames{i}).data) ...
- && size(profiles_2d.(profiles_2d_fieldnames{i}).data,3)>=it
- ids_equilibrium.time_slice{it}.profiles_2d{1}.(profiles_2d_fieldnames{i}) = ...
- profiles_2d.(profiles_2d_fieldnames{i}).data(:,:,it);
+ for i_name=1:numel(profiles_2d_fieldnames)
+ if ~any(strcmp(profiles_2d_fieldnames{i_name},special_fields))
+ if ~isstruct(ids_equilibrium.time_slice{it}.profiles_2d{1}.(profiles_2d_fieldnames{i_name}))
+ if ~ischar(profiles_2d.(profiles_2d_fieldnames{i_name}).data) && ~isempty(profiles_2d.(profiles_2d_fieldnames{i_name}).data) ...
+ && size(profiles_2d.(profiles_2d_fieldnames{i_name}).data,3)>=it
+ ids_equilibrium.time_slice{it}.profiles_2d{1}.(profiles_2d_fieldnames{i_name}) = ...
+ profiles_2d.(profiles_2d_fieldnames{i_name}).data(:,:,it);
end
else
special_fields{end+1} = profiles_2d_fieldnames{i};
@@ -547,7 +559,6 @@ for it=1:ntime
end
end
end
-
% special cases
for it=1:ntime
ids_equilibrium.time_slice{it}.profiles_2d{1}.grid_type.name = profiles_2d.grid_type.name;
@@ -596,7 +607,6 @@ end
% $$$ ids_equilibrium.time_slice{it}.profiles_2d{2}.psi_error_upper(:,:) = 12.*ones(ldim1,ldim2);
% $$$ ids_equilibrium.time_slice{it}.profiles_2d{2}.psi_error_lower(:,:) = 10.*ones(ldim1,ldim2);
% $$$ end
-
% cocos automatic transform
if ~isempty(which('ids_generic_cocos_nodes_transformation_symbolic'))
[ids_equilibrium,cocoscoeff]=ids_generic_cocos_nodes_transformation_symbolic(ids_equilibrium,'equilibrium',gdat_params.cocos_in, ...
diff --git a/matlab/cell2str.m b/matlab/cell2str.m
new file mode 100644
index 0000000000000000000000000000000000000000..464684efdf629c390bac48e1f3e6df871d274ead
--- /dev/null
+++ b/matlab/cell2str.m
@@ -0,0 +1,50 @@
+function str = cell2str(C,n)
+% function str = cell2str(C,n)
+% Convert cell structure to string such that C = eval(str).
+%
+% [+GenLib General Purpose Library+] Swiss Plasma Center EPFL Lausanne 2022. All rights reserved.
+
+assert(iscell(C),'C must be a cell');
+assert(ndims(C)<=2,'Cell arrays with more than 2 dimensions are not supported');
+
+[nrow,ncol] = size(C);
+
+str = '{';
+for ii=1:nrow
+ for jj=1:ncol
+ celldata = C{ii,jj};
+ if isstruct(celldata)
+ structstrcell = struct2str(celldata,n); % returns cell array of strings
+ datastr = '';
+ for kk = 1:numel(structstrcell)
+ datastr = [datastr,structstrcell{kk},'\n']; %#ok<AGROW>
+ end
+ datastr = sprintf(datastr);
+ % contatenate them to get a single string again
+ elseif iscell(celldata)
+ datastr = cell2str(celldata,n); % recursive call
+ else
+ datastr = field2str(celldata,n);
+ end
+ str = [str,datastr]; %#ok<AGROW>
+ if jj == ncol && ii~=nrow
+ str = [str,';']; %#ok<AGROW>
+ elseif jj~=ncol
+ str = [str,',']; %#ok<AGROW>
+ end
+ end
+end
+str = [str,'}'];
+
+end
+
+function datastr = field2str(data,n)
+
+if isnumeric(data) || islogical(data)
+ datastr = mat2str(data,n,'class');
+elseif ischar(data)
+ datastr = ['''',data,''''];
+else
+ error('%s datatype not supported:',class(data))
+end
+end