correct profiles, add Carbon

move to branch IMAS for code camp, so merge this first

matlab/TCV_IMAS/ids2database.m

-function [ids_put_status] = ids2database(shot,run,occurence,ids2put,varargin);
+function [ids_put_status] = ids2database(shot,run,occurence,ids2put,tree_user,tree_tokamak,tree_majorversion,varargin);
 %
 %   [ids_put_status] = ids2database(shot,run,occurence,ids2put,varargin);
 %
 % open/create shot,run and put the ids from ids2put.idsnames (names from fieldnames(ids2put) && within IDS_list)
+% from new ..._env routines, now user/publi database name, tokamak name and major UAL version (3 or 4 at  this stage) need to be specified
+% tree_name: getenv('USER') by default, can be 'public'
+% tree_tokamak: 'tcv' by default
+% tree_majorversion: '4' by default
 %
 % varargin{1}: tbd
 %
 %
 ids_put_status = 0;
 
+tree_user_default = getenv('USER');
+tree_tokamak_default = 'tcv';
+tree_majorversion_default = '4';
+
 % initialize input parser
 p = inputParser;
 % no required inputs here so one can call with empty args and get defaults parameters
@@ -17,6 +25,9 @@ p.addOptional('shot', [], @(x) (isnumeric(x) && isscalar(x) && (x == round(x))))
 p.addOptional('run', [], @(x) (isnumeric(x) && isscalar(x) && (x == round(x)))); % integer
 p.addOptional('occurence', 0, @(x) isempty(x) || (isnumeric(x) && isscalar(x) && (x == round(x)))); % integer
 p.addOptional('ids2put', struct([]), @(x) (isstruct(x)));
+p.addOptional('tree_user', tree_user_default, @(x) (ischar(x)));
+p.addOptional('tree_tokamak', tree_tokamak_default, @(x) (ischar(x)));
+p.addOptional('tree_majorversion', tree_majorversion_default, @(x) (ischar(x)));
 
 p.parse;
 defaults_ids2database = p.Results; % to keep track of defaults
@@ -30,8 +41,20 @@ elseif nargin==2
 elseif nargin==3 && ~ischar(occurence)
   p.parse(shot,run,occurence);
   params = p.Results;
-elseif nargin>=4
-  p.parse(shot,run,occurence,varargin{:});
+elseif nargin==4
+  p.parse(shot,run,occurence,ids2put);
+  params = p.Results;
+elseif nargin==5
+  p.parse(shot,run,occurence,ids2put,tree_user);
+  params = p.Results;
+elseif nargin==6
+  p.parse(shot,run,occurence,ids2put,tree_user,tree_tokamak);
+  params = p.Results;
+elseif nargin==7
+  p.parse(shot,run,occurence,ids2put,tree_user,tree_tokamak,tree_majorversion);
+  params = p.Results;
+elseif nargin>=8
+  p.parse(shot,run,occurence,ids2put,tree_user,tree_tokamak,tree_majorversion,varargin{:});
   params = p.Results;
 else
   p.parse;
@@ -50,7 +73,12 @@ params_ids2database = params;
 
 % check ids_names
 ids_names=fieldnames(ids2put);
-ids_full_list = IDS_list;
+try
+  ids_full_list = IDS_list;
+catch
+  ids_full_list = {'equilibrium', 'magnetics', 'tf', 'pf_active','wall','core_profiles'};
+  warning(['IDS_list not available, quick fix introducing list of ids ready for TCV: ' fprintf('%s ',ids_full_list{:}) char(10)]);
+end
 ids_names_notok = setdiff(ids_names,ids_full_list);
 if ~isempty(ids_names_notok)
   disp(['these subfields are not ids names, so not used: ' sprintf('%s ',ids_names_notok{:})])
@@ -71,9 +99,10 @@ try
     shot_is_new = 0;
   end
   if shot_is_new
-    idx  = imas_create('ids',shot,run,0,0); % 
+    idx  = imas_create_env('ids',shot,run,0,0,params_ids2database.tree_user,params_ids2database.tree_tokamak, ...
+          params_ids2database.tree_majorversion); % 
   else
-    idx  = imas_open('ids',shot,run); % 
+    idx  = imas_open_env('ids',shot,run,params_ids2database.tree_user,params_ids2database.tree_tokamak,params_ids2database.tree_majorversion); % 
   end
 
   %% Put the field

matlab/TCV_IMAS/tcv2ids2database.m

@@ -14,6 +14,9 @@ function [ids_from_tcv,varargout] = tcv2ids2database(shot,run_out,varargin);
 %                'delta' or empty (default): meaning difference in upper is set (standard error_bar
 %                'added': errorbar is added: upper=data+delta and lower=data-delta
 %                'delta_with_lower': as 'delta' but lower also set
+%           'tree_user': getenv('USER') by default, user for database (could be 'public')
+%           'tree_tokamak': tokamak name, default ('tcv')
+%           'tree_majorversion': default '4'
 %
 % example:
 %   ids_from_tcv = tcv2ids2database(62745,9999,'ids_names',{'pf_active'},'error_bar','added'); % to test only one ids
@@ -41,6 +44,9 @@ p.addOptional('run_out', [], @(x) (isnumeric(x) && isscalar(x) && (x == round(x)
 p.addOptional('occurence', 0, @(x) (isnumeric(x) && isscalar(x) && (x == round(x)))); % integer
 p.addOptional('ids_names', {'equilibrium', 'magnetics', 'pf_active','wall', 'tf','core_profiles'}, @(x) isempty(x) | iscell(x)); % char or cell array
 p.addOptional('error_bar', 'delta', @(x) isempty(x) | ischar(x) ); % char or cell array
+p.addOptional('tree_user', getenv('USER'), @(x) isempty(x) | ischar(x) ); % char
+p.addOptional('tree_tokamak', 'tcv', @(x) isempty(x) | ischar(x) ); % char
+p.addOptional('tree_majorversion', '4', @(x) isempty(x) | ischar(x) ); % char
 
 p.parse;
 defaults_tcv2ids2database = p.Results; % to keep track of defaults
@@ -82,7 +88,7 @@ for i=1:length(params.ids_names)
 end
 params.ids_names = ids_names_ok;
 params_tcv2ids2database = params;
-params_not_in_tcv2ids = {'run_out','occurence'};
+params_not_in_tcv2ids = {'run_out','occurence','tree_user','tree_tokamak','tree_majorversion'};
 params_tcv2ids = rmfield(params_tcv2ids2database,params_not_in_tcv2ids);
 [ids_from_tcv,idsok] = tcv2ids(shot,params_tcv2ids);
 ids_from_tcv.params_tcv2ids2database = params_tcv2ids2database;
@@ -97,7 +103,8 @@ if isfield(ids_from_tcv,'tf')
 end
 %%
 
-[ids_put_status] = ids2database(shot,run_out,ids_from_tcv.params_tcv2ids2database.occurence,ids_from_tcv);
+[ids_put_status] = ids2database(shot,run_out,ids_from_tcv.params_tcv2ids2database.occurence,ids_from_tcv, ...
+          params_tcv2ids2database.tree_user,params_tcv2ids2database.tree_tokamak,params_tcv2ids2database.tree_majorversion);
 
 ids_from_tcv.ids_put_status = ids_put_status;
 

matlab/TCV_IMAS/tcv_get_ids_core_profiles.m

@@ -83,6 +83,8 @@ global_quantities.beta_pol = gdat(params_cores_profiles.shot,'betap','machine',m
 global_quantities_desc.beta_pol = 'betap';
 global_quantities.energy_diamagnetic = gdat(params_cores_profiles.shot,'w_mhd','machine',machine);
 global_quantities_desc.energy_diamagnetic = 'w_mhd';
+global_quantities.z_eff_resistive = gdat(params_cores_profiles.shot,'results.conf:z_eff','machine',machine,'fit',1);
+global_quantities_desc.z_eff_resistive = 'results.conf:z_eff';
 
 ids_cores_profiles_description.global_quantities = global_quantities_desc;
 
@@ -129,19 +131,24 @@ for it=1:length(ids_cores_profiles.time)
   ids_cores_profiles.profiles_1d{it}.electrons.temperature_fit.rho_tor_norm = temp_1d.te_rho.grids_1d.rhotornorm(:,it_thom(it));
   ids_cores_profiles.profiles_1d{it}.electrons.temperature_fit.source = {'Thomson, interpos fit'};
   ids_cores_profiles.profiles_1d{it}.electrons.density = temp_1d.fit.ne_rho.data(:,it);
+  ids_cores_profiles.profiles_1d{it}.electrons.density_thermal = ids_cores_profiles.profiles_1d{it}.electrons.density;
   ids_cores_profiles.profiles_1d{it}.electrons.density_fit.measured = temp_1d.ne_rho.data(:,it_thom(it));
   ids_cores_profiles.profiles_1d{it}.electrons.density_fit.time_measurement = temp_1d.ne_rho.t(it_thom(it));
   ids_cores_profiles.profiles_1d{it}.electrons.density_fit.rho_tor_norm = temp_1d.ne_rho.grids_1d.rhotornorm(:,it_thom(it));
   ids_cores_profiles.profiles_1d{it}.electrons.density_fit.source = {'Thomson, interpos fit'};
-  ids_cores_profiles.profiles_1d{it}.electrons.pressure_thermal = 1.6022e-19.*ids_cores_profiles.profiles_1d{it}.electrons.density ...
+  ids_cores_profiles.profiles_1d{it}.electrons.pressure_thermal = 1.6022e-19.*ids_cores_profiles.profiles_1d{it}.electrons.density_thermal ...
       .* ids_cores_profiles.profiles_1d{it}.electrons.temperature;
+  ids_cores_profiles.profiles_1d{it}.zeff = global_quantities.z_eff_resistive.data(it) .* ...
+      ones(size(ids_cores_profiles.profiles_1d{it}.electrons.density));
 end
 
+zeff_error = 0.5;
 switch error_bar
  case 'delta'
   for it=1:length(ids_cores_profiles.time)
     ids_cores_profiles.profiles_1d{it}.electrons.temperature_fit.measured_error_upper = temp_1d.te_rho.error_bar(:,it_thom(it));
     ids_cores_profiles.profiles_1d{it}.electrons.density_fit.measured_error_upper = temp_1d.ne_rho.error_bar(:,it_thom(it));
+    ids_cores_profiles.profiles_1d{it}.zeff_error_upper = zeff_error .* ones(size(ids_cores_profiles.profiles_1d{it}.zeff));
   end
  case 'delta_with_lower'
   for it=1:length(ids_cores_profiles.time)
@@ -151,6 +158,8 @@ switch error_bar
     ids_cores_profiles.profiles_1d{it}.electrons.density_fit.measured_error_upper = temp_1d.ne_rho.error_bar(:,it_thom(it));
     ids_cores_profiles.profiles_1d{it}.electrons.density_fit.measured_error_lower = ...
         ids_cores_profiles.profiles_1d{it}.electrons.density_fit.measured_error_upper;
+    ids_cores_profiles.profiles_1d{it}.zeff_error_upper = zeff_error .* ones(size(ids_cores_profiles.profiles_1d{it}.zeff));
+    ids_cores_profiles.profiles_1d{it}.zeff_error_lower = ids_cores_profiles.profiles_1d{it}.zeff_error_upper;
   end
  case 'added'
   for it=1:length(ids_cores_profiles.time)
@@ -166,6 +175,10 @@ switch error_bar
     ids_cores_profiles.profiles_1d{it}.electrons.density_fit.measured_error_lower = ...
         ids_cores_profiles.profiles_1d{it}.electrons.density_fit.measured ...
         - temp_1d.te_rho.error_bar(:,it_thom(it));
+    ids_cores_profiles.profiles_1d{it}.zeff_error_upper = min(6,ids_cores_profiles.profiles_1d{it}.zeff ...
+        + zeff_error .* ones(size(ids_cores_profiles.profiles_1d{it}.zeff)));
+    ids_cores_profiles.profiles_1d{it}.zeff_error_upper = max(1,ids_cores_profiles.profiles_1d{it}.zeff ...
+        - zeff_error .* ones(size(ids_cores_profiles.profiles_1d{it}.zeff)));
   end
  otherwise
   error(['tcv_ids_bpol_loop: error_bar option not known: ' error_bar])
@@ -175,7 +188,6 @@ end
 temp_1d.cxrs_rho = gdat(params_cores_profiles.shot,'cxrs','machine',machine,'fit',1);
 temp_1d.ti_conf_rho = gdat(params_cores_profiles.shot,'results.conf:ti','machine',machine,'fit',1);
 temp_1d.ni_conf_rho = gdat(params_cores_profiles.shot,'results.conf:ni','machine',machine,'fit',1);
-temp_1d.z_eff_conf = gdat(params_cores_profiles.shot,'results.conf:z_eff','machine',machine,'fit',1);
 data_fullpath_raw = 'Ti(C sometimes B) from cxrs system 1 to 3';
 temp_1d.ti.raw = temp_1d.cxrs_rho.ti.raw;
 temp_1d.ti.raw.shot = temp_1d.cxrs_rho.shot;temp_1d.ti.raw.gdat_params = temp_1d.cxrs_rho.gdat_params;
@@ -191,6 +203,8 @@ temp_1d.ti.fit =get_grids_1d(temp_1d.ti.fit,1,1);
 temp_1d.ni.fit = temp_1d.ni_conf_rho;
 temp_1d.ni.fit =get_grids_1d(temp_1d.ni.fit,1,1);
 it_ti = iround_os(temp_1d.ti.fit.t,ids_cores_profiles.time);
+% assumed 1 impurity with Zp=6
+Zp = 6.;
 for it=1:length(ids_cores_profiles.time)
   ids_cores_profiles.profiles_1d{it}.ion{1}.element{1}.a = 2.;
   ids_cores_profiles.profiles_1d{it}.ion{1}.element{1}.z_n = 1;
@@ -199,16 +213,39 @@ for it=1:length(ids_cores_profiles.time)
   ids_cores_profiles.profiles_1d{it}.ion{1}.element{1}.label = 'D+';
   ids_cores_profiles.profiles_1d{it}.ion{1}.multiple_states_flag = 0;
   ids_cores_profiles.profiles_1d{it}.ion{1}.temperature = temp_1d.ti.fit.data(:,it_ti(it));
-  ids_cores_profiles.profiles_1d{it}.ion{1}.density = temp_1d.ni.fit.data(:,it_ti(it));
-  ids_cores_profiles.profiles_1d{it}.ion{1}.pressure_thermal = 1.6022e-19.*ids_cores_profiles.profiles_1d{it}.ion{1}.density ...
+  ids_cores_profiles.profiles_1d{it}.ion{1}.density = (temp_1d.ni.fit.data(:,it_ti(it)).*Zp-ids_cores_profiles.profiles_1d{it}.electrons.density)./(Zp-1.);
+  ids_cores_profiles.profiles_1d{it}.ion{1}.density_thermal = ids_cores_profiles.profiles_1d{it}.ion{1}.density;
+  ids_cores_profiles.profiles_1d{it}.ion{1}.pressure_thermal = 1.6022e-19.*ids_cores_profiles.profiles_1d{it}.ion{1}.density_thermal ...
       .* ids_cores_profiles.profiles_1d{it}.ion{1}.temperature;
   %
   ids_cores_profiles.profiles_1d{it}.t_i_average = ids_cores_profiles.profiles_1d{it}.ion{1}.temperature;
-  ids_cores_profiles.profiles_1d{it}.n_i_thermal_total = ids_cores_profiles.profiles_1d{it}.ion{1}.density;
+  ids_cores_profiles.profiles_1d{it}.n_i_thermal_total = ids_cores_profiles.profiles_1d{it}.ion{1}.density_thermal;
   ids_cores_profiles.profiles_1d{it}.pressure_ion_total = 1.6022e-19 .* ids_cores_profiles.profiles_1d{it}.n_i_thermal_total ...
       .* ids_cores_profiles.profiles_1d{it}.t_i_average;
   ids_cores_profiles.profiles_1d{it}.pressure_thermal = ids_cores_profiles.profiles_1d{it}.pressure_ion_total ...
       + ids_cores_profiles.profiles_1d{it}.electrons.pressure_thermal;
+  %
+  % C from zeff and above Ti at this stage, should take from cxrs if available but then add something for Zeff matching above
+  %
+  ids_cores_profiles.profiles_1d{it}.ion{2}.element{1}.a = 12.;
+  ids_cores_profiles.profiles_1d{it}.ion{2}.element{1}.z_n = 6.;
+  ids_cores_profiles.profiles_1d{it}.ion{2}.element{1}.atoms_n = 1.;
+  ids_cores_profiles.profiles_1d{it}.ion{2}.element{1}.z_ion = 6.;
+  ids_cores_profiles.profiles_1d{it}.ion{2}.element{1}.label = 'C6+';
+  ids_cores_profiles.profiles_1d{it}.ion{2}.multiple_states_flag = 0;
+  ids_cores_profiles.profiles_1d{it}.ion{2}.temperature = ids_cores_profiles.profiles_1d{it}.ion{1}.temperature;
+  ids_cores_profiles.profiles_1d{it}.ion{2}.density = (ids_cores_profiles.profiles_1d{it}.electrons.density - ids_cores_profiles.profiles_1d{it}.ion{1}.density) ./ ids_cores_profiles.profiles_1d{it}.ion{2}.element{1}.z_ion;
+  ids_cores_profiles.profiles_1d{it}.ion{2}.density_thermal = ids_cores_profiles.profiles_1d{it}.ion{2}.density;
+  ids_cores_profiles.profiles_1d{it}.ion{2}.pressure_thermal = 1.6022e-19.*ids_cores_profiles.profiles_1d{it}.ion{2}.density_thermal ...
+      .* ids_cores_profiles.profiles_1d{it}.ion{2}.temperature;
+  % average/sums
+  ids_cores_profiles.profiles_1d{it}.t_i_average = ids_cores_profiles.profiles_1d{it}.ion{1}.temperature;
+  ids_cores_profiles.profiles_1d{it}.n_i_thermal_total = ids_cores_profiles.profiles_1d{it}.ion{1}.density_thermal + ...
+      ids_cores_profiles.profiles_1d{it}.ion{2}.density_thermal;
+  ids_cores_profiles.profiles_1d{it}.pressure_ion_total = 1.6022e-19 .* ids_cores_profiles.profiles_1d{it}.n_i_thermal_total ...
+      .* ids_cores_profiles.profiles_1d{it}.t_i_average;
+  ids_cores_profiles.profiles_1d{it}.pressure_thermal = ids_cores_profiles.profiles_1d{it}.pressure_ion_total ...
+      + ids_cores_profiles.profiles_1d{it}.electrons.pressure_thermal;  
 end
 if ~isempty(temp_1d.cxrs_rho.ti.fit.data)
   it_raw = iround_os(temp_1d.ti.raw.t,ids_cores_profiles.time);