diff --git a/crpptbx/IMAS/get_ids_equilibrium_fixed_boundary.m b/crpptbx/IMAS/get_ids_equilibrium_fixed_boundary.m
index 60ca92c45b08b24afdc29eb551b8c791f412c46c..9223e759b26af21768108c62a74defe9535904a9 100644
--- a/crpptbx/IMAS/get_ids_equilibrium_fixed_boundary.m
+++ b/crpptbx/IMAS/get_ids_equilibrium_fixed_boundary.m
@@ -46,26 +46,226 @@ if nargout>=3; varargout{1} = rmfield(p.Results,'shot'); end
params_equilibrium = params;
ids_equilibrium = params_equilibrium.empty_struct;
+clear ids_equilibrium_description % so can add new subfields
%
% ids_properties
%
-ids_properties.comment = params_equilibrium.comment;
-ids_properties.homogeneous_time = 1;
-ids_properties.source = ['TCV mds for shot = ' num2str(params_equilibrium.shot)];
-ids_properties.provider = 'get_ids_equilibrium_fixed_boundary';
-ids_properties.creation_date = date;
-
-ids_equilibrium.ids_properties = ids_properties;
+ids_equilibrium.ids_properties.comment = params_equilibrium.comment;
+ids_equilibrium.ids_properties.homogeneous_time = 1;
+ids_equilibrium.ids_properties.source = ['TCV mds for shot = ' num2str(params_equilibrium.shot)];
+ids_equilibrium.ids_properties.provider = 'get_ids_equilibrium_fixed_boundary';
+ids_equilibrium.ids_properties.creation_date = date;
+% As a general rule, for a new substructure under the main ids, construct a local structure like:
+% "global_quantities" with subfields being the relevant data to get and a local structure:
+% "global_quantities_desc" which contains the same subfields themselves containing the gdat string aftre shot used
%
% vacuum_toroidal_field and time, using homogeneous
%
-b0=gdat(params_equilibrium.shot,'b0','source','liuqe'); % to get on liuqe time array
-vacuum_toroidal_field.r0 = b0.r0;
-vacuum_toroidal_field.b0 = b0.data;
+vacuum_toroidal_field.b0=gdat(params_equilibrium.shot,'b0','source','liuqe'); % to get on liuqe time array
+vacuum_toroidal_field_desc.b0 = '''b0'',''source'',''liuqe''';
+vacuum_toroidal_field_desc.r0 = '.r0 subfield from: [''b0'',''source'',''liuqe'']';
+ids_equilibrium.vacuum_toroidal_field.r0 = vacuum_toroidal_field.b0.r0;
+ids_equilibrium.vacuum_toroidal_field.b0 = vacuum_toroidal_field.b0.data;
+ids_equilibrium_description.vacuum_toroidal_field = vacuum_toroidal_field_desc;
+
+ids_equilibrium.time = vacuum_toroidal_field.b0.t;
+ids_equilibrium_description.time = '.t subfield from: [''b0'',''source'',''liuqe'']';
+
+ids_equilibrium.time = [0.6666 0.68333];
+
+ids_equilibrium.time_slice(1:length(ids_equilibrium.time)) = ids_equilibrium.time_slice(1);
+
+% load time array data to copy to time_slices
+
+% global_quantities data into local global_quantities.* structure with correct end names and global_quantities_desc.* with description. Use temp.* and temp_desc.* structures for temporary data
+
+% brute force solution load all eqdsks
+% $$$ for it=1:length(ids_equilibrium.time)
+% $$$ ids_equilibrium.time(it)
+% $$$ temp.eqdsks{it}=gdat(params_equilibrium.shot,'eqdsk','time',ids_equilibrium.time(it),'write',0);
+% $$$ end
+% $$$ temp_desc.eqdsks{1} = '''eqdsk'',''time'',ids_equilibrium.time(it)';
+
+global_quantities.area = gdat(params_equilibrium.shot,'area_edge');
+global_quantities_desc.area = 'area_edge';
+global_quantities.beta_normal = gdat(params_equilibrium.shot,'betan');
+global_quantities_desc.beta_normal = 'betan';
+global_quantities.beta_pol = gdat(params_equilibrium.shot,'betap');
+global_quantities_desc.beta_pol = 'betap';
+global_quantities.beta_tor = gdat(params_equilibrium.shot,'beta');
+global_quantities_desc.beta_tor = 'beta';
+global_quantities.energy_mhd = gdat(params_equilibrium.shot,'w_mhd');
+global_quantities_desc.energy_mhd = 'w_mhd';
+global_quantities.ip = gdat(params_equilibrium.shot,'ip');
+global_quantities_desc.ip = 'ip';
+% length_pol = gdat(params_equilibrium.shot,'length_pol'); % to be added
+global_quantities.li_3 = gdat(params_equilibrium.shot,'li');
+global_quantities_desc.li_3 = 'li';
+temp.r_magnetic_axis = gdat(params_equilibrium.shot,'r_axis');
+temp_desc.r_magnetic_axis = 'r_axis';
+temp.z_magnetic_axis = gdat(params_equilibrium.shot,'z_axis');
+temp_desc.z_magnetic_axis = 'z_axis';
+temp.psi_axis = gdat(params_equilibrium.shot,'psi_axis'); % needs to add psi_edge sincepsi_axis liuqe assuming 0 dege value
+temp_desc.psi_axis = 'psi_axis';
+global_quantities.psi_boundary = gdat(params_equilibrium.shot,'psi_edge');
+global_quantities_desc.psi_boundary = 'psi_edge';
+global_quantities.q_95 = gdat(params_equilibrium.shot,'q95');
+global_quantities_desc.q_95 = 'q95';
+global_quantities.q_axis = gdat(params_equilibrium.shot,'q0'); % will be checked with q_rho?
+global_quantities_desc.q_axis = 'q0';
+temp.q_rho = gdat(params_equilibrium.shot,'q_rho');
+temp_desc.q_rho = 'q_rho';
+% surface = gdat(params_equilibrium.shot,'surface'); % to be added
+global_quantities.volume = gdat(params_equilibrium.shot,'volume');
+global_quantities_desc.volume = 'volume';
+global_quantities.w_mhd = gdat(params_equilibrium.shot,'w_mhd');
+global_quantities_desc.w_mhd = 'w_mhd';
+
+global_quantities_fieldnames = fieldnames(global_quantities);
+special_fields = {'magnetic_axis', 'psi_axis', 'q_min'}; % fields needing non-automatic treatments
+for it=1:length(ids_equilibrium.time)
+ for i=1:length(global_quantities_fieldnames)
+ if ~any(strcmp(global_quantities_fieldnames{i},special_fields))
+ if ~isstruct(ids_equilibrium.time_slice{it}.global_quantities.(global_quantities_fieldnames{i}))
+ ids_equilibrium.time_slice{it}.global_quantities.(global_quantities_fieldnames{i}) = ...
+ global_quantities.(global_quantities_fieldnames{i}).data(it);
+ else
+ special_fields{end+1} = global_quantities_fieldnames{i};
+ end
+ end
+ end
+end
-ids_equilibrium.vacuum_toroidal_field = vacuum_toroidal_field;
-ids_equilibrium.time = b0.t;
+% special case
+for it=1:length(ids_equilibrium.time)
+ ids_equilibrium.time_slice{it}.global_quantities.magnetic_axis.r = temp.r_magnetic_axis.data(it);
+ ids_equilibrium.time_slice{it}.global_quantities.magnetic_axis.z = temp.z_magnetic_axis.data(it);
+ ids_equilibrium.time_slice{it}.global_quantities.psi_axis = temp.psi_axis.data(it) + ...
+ ids_equilibrium.time_slice{it}.global_quantities.psi_boundary;
+ [zz,izz]=min(temp.q_rho.data(:,it));
+ ids_equilibrium.time_slice{it}.global_quantities.q_min.value = zz;
+ ids_equilibrium.time_slice{it}.global_quantities.q_min.rho_tor_norm = temp.q_rho.grids_1d.rhotornorm(izz);
+end
+
+% for boundary in addition to lcfs
+% active_limiter_point = gdat(params_equilibrium.shot,'active_limiter_point');
+boundary.elongation = gdat(params_equilibrium.shot,'kappa');
+boundary_desc.elongation = 'kappa';
+% elongation_lower = gdat(params_equilibrium.shot,'elongation_lower');
+% elongation_upper = gdat(params_equilibrium.shot,'elongation_upper');
+boundary.minor_radius = gdat(params_equilibrium.shot,'a_minor');
+boundary_desc.minor_radius = 'a_minor';
+% squareness_lower_inner = gdat(params_equilibrium.shot,'squareness_lower_inner');
+% squareness_lower_outer = gdat(params_equilibrium.shot,'squareness_lower_outer');
+% squareness_upper_inner = gdat(params_equilibrium.shot,'squareness_upper_inner');
+% squareness_upper_outer = gdat(params_equilibrium.shot,'squareness_upper_outer');
+% strike_point = gdat(params_equilibrium.shot,'strike_point');
+boundary.triangularity = gdat(params_equilibrium.shot,'delta');
+boundary_desc.triangularity = 'delta';
+boundary.triangularity_lower = gdat(params_equilibrium.shot,'delta_bottom');
+boundary_desc.triangularity_lower = 'delta_bottom';
+boundary.triangularity_upper = gdat(params_equilibrium.shot,'delta_top');
+boundary_desc.triangularity_upper = 'delta_top';
+temp.n_x_point = gdat(params_equilibrium.shot,'tcv_eq(''''n_xpts'''',''''liuqe.m'''')');
+temp_desc.n_x_point = '''tcv_eq(''''n_xpts'''',''''liuqe.m'''')''';
+temp.r_x_point = gdat(params_equilibrium.shot,'tcv_eq(''''r_xpts'''',''''liuqe.m'''')');
+temp_desc.r_x_point = '''tcv_eq(''''r_xpts'''',''''liuqe.m'''')''';
+temp.z_x_point = gdat(params_equilibrium.shot,'tcv_eq(''''z_xpts'''',''''liuqe.m'''')');
+temp_desc.z_x_point = '''tcv_eq(''''z_xpts'''',''''liuqe.m'''')''';
+temp.rgeom = gdat(params_equilibrium.shot,'rgeom');
+temp_desc.rgeom = 'rgeom';
+temp.zgeom = gdat(params_equilibrium.shot,'zgeom');
+temp_desc.zgeom = 'zgeom';
+temp.r_lcfs = gdat(params_equilibrium.shot,'r_contour_edge');
+temp_desc.r_lcfs = 'r_contour_edge';
+temp.z_lcfs = gdat(params_equilibrium.shot,'z_contour_edge');
+temp_desc.z_lcfs = 'z_contour_edge';
+
+boundary_fieldnames = fieldnames(boundary);
+special_fields = {'lcfs', 'geometric_axis', 'x_point'}; % fields needing non-automatic treatments
+for it=1:length(ids_equilibrium.time)
+ for i=1:length(boundary_fieldnames)
+ if ~any(strcmp(boundary_fieldnames{i},special_fields))
+ if ~isstruct(ids_equilibrium.time_slice{it}.boundary.(boundary_fieldnames{i}))
+ ids_equilibrium.time_slice{it}.boundary.(boundary_fieldnames{i}) = ...
+ boundary.(boundary_fieldnames{i}).data(it);
+ else
+ special_fields{end+1} = boundary_fieldnames{i};
+ end
+ end
+ end
+end
+
+% special cases
+for it=1:length(ids_equilibrium.time)
+ ids_equilibrium.time_slice{it}.boundary.outline.r = temp.r_lcfs.data(it);
+ ids_equilibrium.time_slice{it}.boundary.outline.z = temp.z_lcfs.data(it);
+ ids_equilibrium.time_slice{it}.boundary.lcfs.r = ids_equilibrium.time_slice{it}.boundary.outline.r;
+ ids_equilibrium.time_slice{it}.boundary.lcfs.z = ids_equilibrium.time_slice{it}.boundary.outline.z;
+ ids_equilibrium.time_slice{it}.boundary.geometric_axis.r = temp.rgeom.data(it);
+ ids_equilibrium.time_slice{it}.boundary.geometric_axis.z = temp.zgeom.data(it);
+ if temp.n_x_point.data(it) > 0
+ ids_equilibrium.time_slice{it}.boundary.x_point(1:temp.n_x_point.data(it)) = ids_equilibrium.time_slice{it}.boundary.x_point(1);
+ for i=1:length(temp.n_x_point.data(it))
+ ids_equilibrium.time_slice{it}.boundary.x_point{i}.r = temp.r_x_point.data(i,it);
+ ids_equilibrium.time_slice{it}.boundary.x_point{i}.z = temp.z_x_point.data(i,it);
+ end
+ end
+end
+
+%
+%% profiles_1d (cannot use eqdsk since not same radial mesh)
+%
+% area = gdat(params_equilibrium.shot,'area');
+% b_average = gdat(params_equilibrium.shot,'b_average');
+% beta_pol = gdat(params_equilibrium.shot,'beta_pol');
+% b_field_average = gdat(params_equilibrium.shot,'b_field_average');
+% b_field_max = gdat(params_equilibrium.shot,'b_field_max');
+% b_field_min = gdat(params_equilibrium.shot,'b_field_min');
+% b_max = gdat(params_equilibrium.shot,'b_max');
+% b_min = gdat(params_equilibrium.shot,'b_min');
+% darea_dpsi = gdat(params_equilibrium.shot,'darea_dpsi');
+% darea_drho_tor = gdat(params_equilibrium.shot,'darea_drho_tor');
+profiles_1d.dpressure_dpsi = gdat(params_equilibrium.shot,'pprime');
+% dpsi_drho_tor = gdat(params_equilibrium.shot,'dpsi_drho_tor');
+% dvolume_dpsi = gdat(params_equilibrium.shot,'dvolume_dpsi');
+% dvolume_drho_tor = gdat(params_equilibrium.shot,'dvolume_drho_tor');
+% elongation = gdat(params_equilibrium.shot,'elongation');
+profiles_1d.f_df_dpsi = gdat(params_equilibrium.shot,'ffprime');
+profiles_1d.f = gdat(params_equilibrium.shot,'rbphi_rho');
+% geometric_axis = gdat(params_equilibrium.shot,'geometric_axis');
+% gm1 = gdat(params_equilibrium.shot,'gm1');
+% gm2 = gdat(params_equilibrium.shot,'gm2');
+% gm3 = gdat(params_equilibrium.shot,'gm3');
+% gm4 = gdat(params_equilibrium.shot,'gm4');
+% gm5 = gdat(params_equilibrium.shot,'gm5');
+% gm6 = gdat(params_equilibrium.shot,'gm6');
+% gm7 = gdat(params_equilibrium.shot,'gm7');
+% gm8 = gdat(params_equilibrium.shot,'gm8');
+% gm9 = gdat(params_equilibrium.shot,'gm9');
+% j_parallel = gdat(params_equilibrium.shot,'j_parallel');
+% j_tor = gdat(params_equilibrium.shot,'j_tor');
+% magnetic_shear = gdat(params_equilibrium.shot,'magnetic_shear');
+% mass_density = gdat(params_equilibrium.shot,'mass_density');
+profiles_1d.phi = gdat(params_equilibrium.shot,'phi_tor');
+profiles_1d.pressure = gdat(params_equilibrium.shot,'pressure');
+% psi = gdat(params_equilibrium.shot,'psi_rho'); % (could take from .x of any like rhotor and psi_axis, psi_edge from global_quantities)
+profiles_1d.q = gdat(params_equilibrium.shot,'q_rho');
+profiles_1d.rho_tor = gdat(params_equilibrium.shot,'rhotor');
+%rho_tor_norm = gdat(params_equilibrium.shot,'rhotor_norm'); % from rho_tor
+profiles_1d.rho_volume_norm = gdat(params_equilibrium.shot,'rhovol');
+% r_inboard = gdat(params_equilibrium.shot,'r_inboard');
+% r_outboard = gdat(params_equilibrium.shot,'r_outboard');
+% squareness_lower_inner = gdat(params_equilibrium.shot,'squareness_lower_inner');
+% squareness_lower_outer = gdat(params_equilibrium.shot,'squareness_lower_outer');
+% squareness_upper_inner = gdat(params_equilibrium.shot,'squareness_upper_inner');
+% squareness_upper_outer = gdat(params_equilibrium.shot,'squareness_upper_outer');
+% surface = gdat(params_equilibrium.shot,'surface');
+% trapped_fraction = gdat(params_equilibrium.shot,'trapped_fraction');
+% triangularity_lower = gdat(params_equilibrium.shot,'triangularity_lower');
+% triangularity_upper = gdat(params_equilibrium.shot,'triangularity_upper');
+profiles_1d.volume = gdat(params_equilibrium.shot,'volume_rho');
-time_slice(1:length(ids_equilibrium.time)) = ids_equilibrium.time_slice(1);
+keyboard
diff --git a/crpptbx/TCV/gdat_tcv.m b/crpptbx/TCV/gdat_tcv.m
index adfb5aa4911b60565684051d9cc13884a0a41858..d3f2b77db13a8ac974f9d916e98a880dcc7cb19c 100644
--- a/crpptbx/TCV/gdat_tcv.m
+++ b/crpptbx/TCV/gdat_tcv.m
@@ -932,6 +932,34 @@ elseif strcmp(mapping_for_tcv.method,'switchcase')
gdat_data.gdat_params.time = time;
if (gdat_params.nverbose>=3); disp(['"time" is expected as an option, choose default time = ' num2str(time)]); end
end
+ default_write_eqdsk = 1;
+ if isfield(gdat_data.gdat_params,'write') && ~isempty(gdat_data.gdat_params.write)
+ if ischar(gdat_data.gdat_params.write)
+ if strcmp(lower(gdat_data.gdat_params.write),'no')
+ gdat_data.gdat_params.write = 0;
+ elseif strcmp(lower(gdat_data.gdat_params.write),'yes')
+ gdat_data.gdat_params.write = 1;
+ else
+ if (gdat_params.nverbose>=3); disp(['expects 0 or 1, event yes or no for write option, not: ' ...
+ gdat_data.gdat_params.write ', use default 1']);
+ end
+ gdat_data.gdat_params.write = [];
+ end
+ end
+ if gdat_data.gdat_params.write~=1 && gdat_data.gdat_params.write~=0
+ gdat_data.gdat_params.write = default_write_eqdsk;
+ end
+ else
+ gdat_data.gdat_params.write = default_write_eqdsk;
+ end
+ default_map_eqdsk_psirz = 0;
+ if isfield(gdat_data.gdat_params,'map_eqdsk_psirz') && ~isempty(gdat_data.gdat_params.map_eqdsk_psirz)
+ if ~isnumeric(gdat_data.gdat_params.map_eqdsk_psirz) || (gdat_data.gdat_params.map_eqdsk_psirz~=1 && gdat_data.gdat_params.map_eqdsk_psirz~=0)
+ gdat_data.gdat_params.map_eqdsk_psirz = default_map_eqdsk_psirz_eqdsk;
+ end
+ else
+ gdat_data.gdat_params.map_eqdsk_psirz = default_map_eqdsk_psirz;
+ end
gdat_data.gdat_params.time = time;
gdat_data.t = time;
zshift = 0.;
@@ -959,7 +987,13 @@ elseif strcmp(mapping_for_tcv.method,'switchcase')
end
return
end
- eqdskval=read_eqdsk(fnames_readresults{4},cocos_read_results_for_chease,0,[],[],1); % read_results now has relevant cocos LIUQE=17
+ ii = regexpi(fnames_readresults{4},'eqdsksigns');
+ if ~isempty(ii)
+ eqdskval=read_eqdsk(fnames_readresults{4},cocos_read_results_for_chease,0,[],[],1); % read_results now has relevant cocos LIUQE= 17
+ else
+ disp('wrong name check order of eqdsk outputs')
+ return
+ end
for i=1:length(fnames_readresults)
unix(['rm ' fnames_readresults{i}]);
end
@@ -968,7 +1002,9 @@ elseif strcmp(mapping_for_tcv.method,'switchcase')
[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);
+ if gdat_data.gdat_params.write==1
+ write_eqdsk(fnamefull,eqdsk_cocos_in,cocos_in,[],[],[],1);
+ end
% Now gdat_tcv should return the convention from LIUQE which is COCOS=17, except if specified in option
% create standard filename name from shot, time_eff (cocos will be added by write_eqdsk)
cocos_out = 17;
@@ -982,20 +1018,32 @@ elseif strcmp(mapping_for_tcv.method,'switchcase')
% cannot use it for psi(R,Z) in .data and .x since R, Z might be different at different times,
% so project psi(R,Z) on Rmesh, Zmesh of 1st time
if length(time) > 1
- gdat_data.eqdsk{itime} = write_eqdsk(fnamefull,eqdsk_cocosout,cocos_out);
- if itime==1
- gdat_data.data(:,:,itime) = gdat_data.eqdsk{itime}.psi;
- gdat_data.dim{1} = gdat_data.eqdsk{itime}.rmesh;
- gdat_data.dim{2} = gdat_data.eqdsk{itime}.zmesh;
+ if gdat_data.gdat_params.write==1
+ gdat_data.eqdsk{itime} = write_eqdsk(fnamefull,eqdsk_cocosout,cocos_out);
else
- xx=repmat(reshape(gdat_data.dim{1},length(gdat_data.dim{1}),1),1,size(gdat_data.eqdsk{itime}.psi,2));
- yy=repmat(reshape(gdat_data.dim{2},1,length(gdat_data.dim{2})),size(gdat_data.eqdsk{itime}.psi,1),1);
- aa = interpos2Dcartesian(gdat_data.eqdsk{itime}.rmesh,gdat_data.eqdsk{itime}.zmesh ...
+ gdat_data.eqdsk{itime} = eqdsk_cocosout;
+ end
+ if gdat_data.gdat_params.map_eqdsk_psirz==1
+ if itime==1
+ gdat_data.data(:,:,itime) = gdat_data.eqdsk{itime}.psi;
+ gdat_data.dim{1} = gdat_data.eqdsk{itime}.rmesh;
+ gdat_data.dim{2} = gdat_data.eqdsk{itime}.zmesh;
+ else
+ xx=repmat(reshape(gdat_data.dim{1},length(gdat_data.dim{1}),1),1,size(gdat_data.eqdsk{itime}.psi,2));
+ yy=repmat(reshape(gdat_data.dim{2},1,length(gdat_data.dim{2})),size(gdat_data.eqdsk{itime}.psi,1),1);
+ aa = interpos2Dcartesian(gdat_data.eqdsk{itime}.rmesh,gdat_data.eqdsk{itime}.zmesh ...
,gdat_data.eqdsk{itime}.psi,xx,yy,-1,-1);
- gdat_data.data(:,:,itime) = aa;
+ gdat_data.data(:,:,itime) = aa;
+ end
+ else
+ % do not map all psi(r,z) onto same mesh and leave .data, .dim empty (much faster)
end
else
- gdat_data.eqdsk = write_eqdsk(fnamefull,eqdsk_cocosout,cocos_out);
+ if gdat_data.gdat_params.write==1
+ gdat_data.eqdsk = write_eqdsk(fnamefull,eqdsk_cocosout,cocos_out);
+ else
+ gdat_data.eqdsk = eqdsk_cocosout;
+ end
gdat_data.data = gdat_data.eqdsk.psi;
gdat_data.dim{1} = gdat_data.eqdsk.rmesh;
gdat_data.dim{2} = gdat_data.eqdsk.zmesh;
@@ -1003,7 +1051,11 @@ elseif strcmp(mapping_for_tcv.method,'switchcase')
end
gdat_data.dim{3} = gdat_data.t;
gdat_data.x = gdat_data.dim(1:2);
- gdat_data.data_fullpath=['psi(R,Z) and eqdsk from read_eqdsk from LIUQE' num2str(liuqe_version) ';zshift=' num2str(zshift)];
+ if gdat_data.gdat_params.map_eqdsk_psirz==1
+ gdat_data.data_fullpath=['psi(R,Z,t) on same R,Zmesh in .data and eqdsk{itime} from read_eqdsk from LIUQE' num2str(liuqe_version) ';zshift=' num2str(zshift)];
+ else
+ gdat_data.data_fullpath=['eqdsk{itime} from read_eqdsk from LIUQE' num2str(liuqe_version) ';zshift=' num2str(zshift)];
+ end
gdat_data.units = 'T m^2';
gdat_data.dimunits = {'m','m','s'};
gdat_data.request_description = ['data=psi, x=(R,Z), eqdsk contains eqdsk structure with which ' ...
@@ -1842,6 +1894,26 @@ elseif strcmp(mapping_for_tcv.method,'switchcase')
if (length(gdat_data.dim)>=mapping_for_tcv.gdat_timedim); gdat_data.t = gdat_data.dim{mapping_for_tcv.gdat_timedim}; end
if (length(gdat_data.dim)>0); gdat_data.x = gdat_data.dim{1}; end
+ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+ case {'pprime', 'pprime_rho'}
+ if liuqe_matlab==0
+ nodenameeff = ['tcv_eq(''' liuqefortran2liuqematlab('pprime_rho',1,0) ''',''LIUQE' substr_liuqe_tcv_eq ''')'];
+ tracetdi=tdi(nodenameeff);
+ tracetdi.dim{1} = sqrt(tracetdi.dim{1}); % correct x-axis psi_norm to rhopol
+ tracetdi.data = tracetdi.data ./2 ./pi; % correct node assumption
+ else
+ nodenameeff = ['tcv_eq(''pprime_rho'',''LIUQE.M' substr_liuqe_tcv_eq ''')'];
+ tracetdi=tdi(nodenameeff);
+ end
+ gdat_data.data = tracetdi.data;
+ gdat_data.dim = tracetdi.dim;
+ gdat_data.t = gdat_data.dim{mapping_for_tcv.gdat_timedim};
+ gdat_data.x = gdat_data.dim{setdiff([1 2],mapping_for_tcv.gdat_timedim)};
+ gdat_data.data_fullpath=nodenameeff;
+ gdat_data.dimunits = tracetdi.dimunits;
+ gdat_data.units = tracetdi.units;
+ gdat_data.request_description = 'pprime=dp/dpsi';
+
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
case {'psi_edge'}
% psi at edge, 0 by construction in Liuqe, thus not given
@@ -1872,6 +1944,25 @@ elseif strcmp(mapping_for_tcv.method,'switchcase')
gdat_data.vsurf = vsurf;
gdat_data.vsurf_description = ['time derivative from surface_psi, obtained from \results::surface_flux'];
+ %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
+ case {'r_contour_edge', 'z_contour_edge'}
+ if liuqe_matlab==0
+ nodenameeff=['\results::' data_request_eff(1) '_contour' substr_liuqe];
+ % npts_contour = tdi(['\results::npts_contour' substr_liuqe]);
+ tracetdi=tdi(nodenameeff);
+ gdat_data.request_description = 'NaNs when smaller nb of boundary points at given time, can use \results::npts_contour';
+ else
+ nodenameeff = ['tcv_eq(''' data_request_eff(1) '_edge'',''LIUQE.M' substr_liuqe_tcv_eq ''')'];
+ tracetdi=tdi(nodenameeff);
+ gdat_data.request_description = 'lcfs on same nb of points for all times';
+ end
+ gdat_data.data = tracetdi.data;
+ gdat_data.dim = tracetdi.dim;
+ gdat_data.t = gdat_data.dim{mapping_for_tcv.gdat_timedim};
+ gdat_data.data_fullpath=nodenameeff;
+ gdat_data.dimunits = tracetdi.dimunits;
+ gdat_data.units = tracetdi.units;
+
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
case {'rhotor_edge', 'rhotor', 'rhotor_norm'}
% Phi(LCFS) = int(Bphi dSphi), can use Eq.(11) of "Tokamak coordinate conventions: COCOS" paper:
diff --git a/crpptbx/TCV/tcv_help_parameters.m b/crpptbx/TCV/tcv_help_parameters.m
index 5a7d2f5f561e05e6da6bf0d8f9c25c533d51c5bf..db22c4ebb7cd51caf11b25d1dc1b7f226370c921 100644
--- a/crpptbx/TCV/tcv_help_parameters.m
+++ b/crpptbx/TCV/tcv_help_parameters.m
@@ -52,7 +52,8 @@ help_struct_all.camera = ['sxr: for MPX: ''central'', ''top'' (default), ''botto
help_struct_all.freq = '''slow'', default, lower sampling; ''fast'' full samples for both mpx and xtomo';
help_struct_all.max_adcs = 'rtc: source=''adcs'' maximum nb of adc channels loaded for each board in each active node';
help_struct_all.nfft = '512 (default) changes time resolution in spectrogram in gdat_plot for ''mhd'' request';
-
+help_struct_all.map_eqdsk_psirz = 'eqdsk: if time array provided, maps all psi(R,Z,t) on same R,Zmesh in .data (1) or not (0, default)';
+help_struct_all.write = 'eqdsk: write eqdsk while loading data (1, default) or not (0)';
%help_struct_all. = '';
diff --git a/crpptbx/TCV/tcv_requests_mapping.m b/crpptbx/TCV/tcv_requests_mapping.m
index da7cecd65acfd5995aa74eb2993d7a7bc987dc3a..93b32b12741719fd3d59a45bd0f03c7884d6e370 100644
--- a/crpptbx/TCV/tcv_requests_mapping.m
+++ b/crpptbx/TCV/tcv_requests_mapping.m
@@ -236,6 +236,10 @@ switch lower(data_request)
mapping.timedim = 1;
mapping.label = 'various powers';
mapping.method = 'switchcase';
+ case {'pprime', 'pprime_rho'}
+ mapping.timedim = 2;
+ mapping.label = 'pprime';
+ mapping.method = 'switchcase';
case {'psi_axis', 'psi_mag'}
mapping.timedim = 1;
mapping.method = 'tdiliuqe';
@@ -289,6 +293,10 @@ switch lower(data_request)
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'''')';
+
+ mapping.label = 'R\_lcfs';
+ mapping.method = 'switchcase';
+
% $$$ mapping.method = 'expression';
% $$$ mapping.expression = ['params_eff = gdat_data.gdat_params;params_eff.data_request=''tcv_eq(''''''''r_edge'''''''',''''''''LIUQE.M'''''''')''; ' ...
% $$$ 'gdat_tmp=gdat_tcv(shot,params_eff);gdat_tmp.data=gdat_tmp.data(:,:,end);' ...
@@ -354,6 +362,10 @@ switch lower(data_request)
case 'transp'
mapping.label = 'transp output';
mapping.method = 'switchcase';
+ case {'ttprime', 'ttprime_rho'}
+ mapping.timedim = 2;
+ mapping.label = 'pprime';
+ mapping.method = 'switchcase';
case 'vloop'
mapping.timedim = 1;
mapping.label = 'Vloop';
@@ -398,6 +410,10 @@ switch lower(data_request)
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'''')';
+
+ mapping.label = 'Z\_lcfs';
+ mapping.method = 'switchcase';
+
% $$$ mapping.method = 'expression';
% $$$ mapping.expression = ['params_eff = gdat_data.gdat_params;params_eff.data_request=''tcv_eq(''''''''z_edge'''''''',''''''''LIUQE.M'''''''')''; ' ...
% $$$ 'gdat_tmp=gdat_tcv(shot,params_eff);gdat_tmp.data=gdat_tmp.data(:,:,end);' ...
diff --git a/crpptbx/liuqefortran2liuqematlab.m b/crpptbx/liuqefortran2liuqematlab.m
index 6d3d783ebd00c9bc8abbafb973f58468e3834455..b5effd13341d8dfd9bcfeef771b9e475dea27806 100644
--- a/crpptbx/liuqefortran2liuqematlab.m
+++ b/crpptbx/liuqefortran2liuqematlab.m
@@ -27,6 +27,7 @@ function liuqe_fortran_matlab = liuqefortran2liuqematlab(varargin);
liuqe_fortran_matlab_table = [ ...
{'l_i'} , {'l_i_3'} ; ...
{'i_p'} , {'i_pl'} ; ...
+ {'pprime_psi'} , {'pprime_rho'} ; ... % warning on different x-mesh
{'surface_flux'} , {'psi_surf'} ; ...
{'q_zero'} , {'q_axis'} ; ...
{'q_psi'} , {'q'} ; ...