diff --git a/crpptbx/TCV/gdat_tcv.m b/crpptbx/TCV/gdat_tcv.m
index 555bff0a1cbb434e10ee2d139d5f886ee307b6c3..106b268f0ac0e130852d6b10902d8f268765adbb 100644
--- a/crpptbx/TCV/gdat_tcv.m
+++ b/crpptbx/TCV/gdat_tcv.m
@@ -936,6 +936,31 @@ elseif strcmp(mapping_for_tcv.method,'switchcase')
gdat_data.dim{1}=rho;
gdat_data.dimunits=[{'sqrt(psi_norm)'} ; {'time [s]'}];
gdat_data.x=rho;
+ % add various rhos in grids_1d
+ 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);
+ gdat_data.grids_1d.rhopolnorm = gdat_data.x;
+ psi0 = interpos(rhotor.t',rhotor.psi_axis,gdat_data.t,-0.01);
+ 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);
+ 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);
+ ii=find(~isnan(gdat_data.grids_1d.rhopolnorm(:,it)));
+ 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(:,it));
+ gdat_data.grids_1d.rhovolnorm(:,it)=interpos(rhovol.x,rhovol.data(:,it_vol_eff),gdat_data.grids_1d.rhopolnorm(:,it));
+ end
+ end
+ gdat_data.grids_1d.rhotor_edge=interpos(rhotor.t',rhotor.rhotor_edge,gdat_data.t',-0.1);
+
%%%%%%%%%%%
% if nete_rho, copy data as .ne, get .te and put pe=e ne Te in data:
if strcmp(data_request_eff(1:4),'nete')
@@ -1282,6 +1307,27 @@ elseif strcmp(mapping_for_tcv.method,'switchcase')
gdat_data.dimunits{2} = 's';
gdat_data.units = '';
gdat_data.request_description = 'q(rhopol\_norm)';
+ % add various rhos in grids_1d
+ 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);
+ gdat_data.grids_1d.rhopolnorm = gdat_data.x;
+ psi0 = interpos(rhotor.t',rhotor.psi_axis,gdat_data.t,-0.01);
+ gdat_data.grids_1d.psi = gdat_data.grids_1d.rhopolnorm.^2*reshape(psi0,1,length(psi0));
+ 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);
+ 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
+ gdat_data.grids_1d.rhotor_edge=interpos(rhotor.t',rhotor.rhotor_edge,gdat_data.t',-0.1);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
case {'psi_edge'}
@@ -1310,14 +1356,20 @@ elseif strcmp(mapping_for_tcv.method,'switchcase')
% 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
- params_eff = gdat_data.gdat_params;
- params_eff.data_request='q_rho';
- q_rho=gdat_tcv([],params_eff);
+ % 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
@@ -1328,14 +1380,14 @@ elseif strcmp(mapping_for_tcv.method,'switchcase')
return
end
rhoequal = linspace(0,1,length(q_rho.dim{1}));
- if strcmp(data_request,'rhotor_edge')
+ 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,'rhotor')
+ 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;
@@ -1347,17 +1399,22 @@ elseif strcmp(mapping_for_tcv.method,'switchcase')
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.x(ij).^2,q_rho.data(ij,it),rhoequal.^2);
+ [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,'rhotor_edge')
+ if strcmp(data_request_eff,'rhotor_edge')
gdat_data.data(it) = dataeff(end);
- elseif strcmp(data_request,'rhotor')
+ elseif strcmp(data_request_eff,'rhotor')
gdat_data.data(:,it) = dataeff./dataeff(end);
gdat_data.rhotor_edge(it) = dataeff(end);
end
@@ -1384,7 +1441,7 @@ elseif strcmp(mapping_for_tcv.method,'switchcase')
end
end
gdat_data.units = tracetdi.units;
- if strcmp(data_request,'volume')
+ 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,:)'];
@@ -1393,11 +1450,12 @@ elseif strcmp(mapping_for_tcv.method,'switchcase')
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,'volume_rho')
+ if strcmp(data_request_eff,'volume_rho')
gdat_data.data_fullpath=['\results::psitbx:vol'];
gdat_data.request_description = 'volume(rho)';
- elseif strcmp(data_request,'rhovol')
+ 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)';
@@ -1409,7 +1467,7 @@ elseif strcmp(mapping_for_tcv.method,'switchcase')
return
end
end
- gdat_data.t = gdat_data.dim{mapping_for_tcv.gdat_timedim};
+ gdat_data.t = gdat_data.dim{mapping_for_tcv.gdat_timedim};
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
case {'sxr', 'mpx'}