Skip to content
Snippets Groups Projects
Commit f8cdc498 authored by Olivier Sauter's avatar Olivier Sauter
Browse files

added powers to jet and radiated from bolo 

git-svn-id: https://spcsvn.epfl.ch/repos/TCV/gdat/trunk@7760 d63d8f72-b253-0410-a779-e742ad2e26cf
parent 14f0bd30
No related branches found
No related tags found
No related merge requests found
Hide whitespace changes
Inline Side-by-side
crpptbx/JET/gdat_jet.m
+ 288
− 76
View file @ f8cdc498
......@@ -782,12 +782,123 @@ elseif strcmp(mapping_for_jet.method,'switchcase')
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
case {'ne','te'}
nodenameeff_rho = 'rho';
params_eff = gdat_data.gdat_params;
if isfield(params_eff,'source') && ~isempty(params_eff.source)
if strcmp(lower(params_eff.source),'chain2')
params_eff.source = 'hrt2';
end
if strcmp(lower(params_eff.source),'hrt2')
nodenameeff_rho = []; % profiles already on rho
end
else
params_eff.source = 'hrts';
end
gdat_data.gdat_params = params_eff;
params_eff.doplot = 0;
% ne or Te from Thomson data on raw z mesh vs (z,t)
if ~(isfield(gdat_data.gdat_params,'edge') && ~isempty(gdat_data.gdat_params.edge) && ...
gdat_data.gdat_params.edge>0)
gdat_data.gdat_params.edge = 0;
params_eff.data_request = {'ppf',params_eff.source,data_request_eff};
aa = gdat_jet(shot,params_eff);
if isempty(aa.data) || isempty(aa.t)
if gdat_params.nverbose>=3;
aa
disp(['with data_request= ' data_request_eff])
end
return
end
gdat_data.data = aa.data;
gdat_data.dim = aa.dim;
gdat_data.t = aa.dim{2};
gdat_data.x = aa.dim{1};
gdat_data.dimunits=[{'R [m]'} ; {'time [s]'}];
gdat_data.data_fullpath=[data_request_eff ' from ' params_eff.source];
gdat_data.(data_request_eff).data = aa.data;
gdat_data.(data_request_eff).t = aa.t;
gdat_data.(data_request_eff).r = aa.x;
if any(strcmp(fieldnames(aa),'units')) && ~isempty(aa.units)
gdat_data.units=aa.units;
else
if strcmp(data_request_eff,'ne')
gdat_data.units = 'm^{-3}';
elseif strcmp(data_request_eff,'te')
gdat_data.units = 'eV';
end
end
params_eff.data_request = {'ppf',params_eff.source,['d' data_request_eff]};
aaerr = gdat_jet(shot,params_eff);
gdat_data.error_bar = aaerr.data;
gdat_data.(data_request_eff).error_bar = aaerr.data;
if ~isempty(nodenameeff_rho)
params_eff.data_request = {'ppf',params_eff.source,'rho'};
aarho = gdat_jet(shot,params_eff);
gdat_data.rhopol = abs(aarho.data);
gdat_data.(data_request_eff).rhopol = abs(aarho.data); % keep rho >0
else
gdat_data.rhopol = gdat_data.x;
gdat_data.(data_request_eff).rhopol = gdat_data.x;
gdat_data.dimunits{1} = 'rhopol';
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
case {'ni','ti'}
nodenameeff_rho = 'rho';
data_type_eff = data_request_eff;
if strcmp(data_request_eff,'ni'); data_type_eff = 'dd'; end
params_eff = gdat_data.gdat_params;
if isfield(params_eff,'source') && ~isempty(params_eff.source)
if strcmp(lower(params_eff.source),'chain2')
params_eff.source = [data_request_eff(1) 'ion'];
end
if strcmp(lower(params_eff.source(2:end)),'ion')
nodenameeff_rho = []; % profiles already on rho
end
else
params_eff.source = [datar_equest_eff(1) 'ion']; % only chain2 at this stage
end
gdat_data.gdat_params = params_eff;
params_eff.doplot = 0;
% ne or Te from Thomson data on raw z mesh vs (z,t)
params_eff.data_request = {'ppf',params_eff.source,data_type_eff};
aa = gdat_jet(shot,params_eff);
if isempty(aa.data) || isempty(aa.t)
if gdat_params.nverbose>=3;
aa
disp(['with data_request= ' data_request_eff])
end
return
end
gdat_data.data = aa.data;
gdat_data.dim = aa.dim;
gdat_data.t = aa.dim{2};
gdat_data.x = aa.dim{1};
gdat_data.dimunits=[{'R [m]'} ; {'time [s]'}];
gdat_data.data_fullpath=[data_request_eff ' from ' params_eff.source];
gdat_data.(data_request_eff).data = aa.data;
gdat_data.(data_request_eff).t = aa.t;
gdat_data.(data_request_eff).r = aa.x;
if any(strcmp(fieldnames(aa),'units')) && ~isempty(aa.units)
gdat_data.units=aa.units;
else
if strcmp(data_request_eff,'ni')
gdat_data.units = 'm^{-3}';
elseif strcmp(data_request_eff,'ti')
gdat_data.units = 'eV';
end
end
params_eff.data_request = {'ppf',params_eff.source,['e' data_type_eff]}; % not given for dd but ok will be empty
aaerr = gdat_jet(shot,params_eff);
gdat_data.error_bar = aaerr.data;
gdat_data.(data_request_eff).error_bar = aaerr.data;
if ~isempty(nodenameeff_rho)
params_eff.data_request = {'ppf',params_eff.source,'rho'};
aarho = gdat_jet(shot,params_eff);
gdat_data.rhopol = abs(aarho.data);
gdat_data.(data_request_eff).rhopol = abs(aarho.data); % keep rho >0
else
gdat_data.rhopol = gdat_data.x;
gdat_data.(data_request_eff).rhopol = gdat_data.x;
gdat_data.dimunits{1} = 'rhopol';
end
[gdat_data] = get_thomson_raw_data(shot,data_request_eff,gdat_data,gdat_data.gdat_params.edge,gdat_params.nverbose);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
case {'ne_rho', 'te_rho', 'nete_rho'}
......@@ -1023,83 +1134,151 @@ elseif strcmp(mapping_for_jet.method,'switchcase')
case {'powers'}
% note: same time array for all main, ec, ohm, nbi, ...
% At this stage fill just ech, later add nbi
% EC
nodenameeff='\results::toray.input:p_gyro';
tracetdi=tdi(nodenameeff);
if isempty(tracetdi.data) || isempty(tracetdi.dim) % || ischar(tracetdi.data) (to add?)
if (gdat_params.nverbose>=1); warning(['problems loading data for ' nodenameeff ' for data_request= ' data_request_eff]); end
gdat_data.ec.data = [];
gdat_data.ec.units = [];
gdat_data.ec.dim=[];
gdat_data.ec.dimunits=[];
gdat_data.ec.t=[];
gdat_data.ec.x=[];
gdat_data.ec.data_fullpath=[];
gdat_data.ec.label='';
sources_avail = {'ohm','nbi','ic','rad'}; % note should allow ech, nbh, ohmic in parameter sources
for i=1:length(sources_avail)
gdat_data.(sources_avail{i}).data = [];
gdat_data.(sources_avail{i}).units = [];
gdat_data.(sources_avail{i}).dim=[];
gdat_data.(sources_avail{i}).dimunits=[];
gdat_data.(sources_avail{i}).t=[];
gdat_data.(sources_avail{i}).x=[];
gdat_data.(sources_avail{i}).data_fullpath=[];
gdat_data.(sources_avail{i}).label=[];
end
if ~isfield(gdat_data.gdat_params,'source') || isempty(gdat_data.gdat_params.source)
gdat_data.gdat_params.source = sources_avail;
elseif ~iscell(gdat_data.gdat_params.source)
if ischar(gdat_data.gdat_params.source)
gdat_data.gdat_params.source = lower(gdat_data.gdat_params.source);
if ~any(strmatch(gdat_data.gdat_params.source,lower(sources_avail)))
if (gdat_params.nverbose>=1)
warning(['source= ' gdat_data.gdat_params.source ' is not part of the available sources: ' sprintf('''%s'' ',sources_avail{:})]);
end
return
else
gdat_data.gdat_params.source = {gdat_data.gdat_params.source};
end
else
if (gdat_params.nverbose>=1); warning([' source parameter not compatible with: ' sprintf('''%s'' ',sources_avail{:})]); end
return
end
else
gdat_data.ec.data = tracetdi.data*1e3; % at this stage p_gyro is in kW'
gdat_data.ec.units = 'W';
gdat_data.ec.dim=tracetdi.dim;
gdat_data.ec.dimunits=tracetdi.dimunits;
gdat_data.ec.t=tracetdi.dim{1};
gdat_data.ec.x=tracetdi.dim{2};
gdat_data.ec.data_fullpath=[nodenameeff];
gdat_data.ec.label='P_{EC}';
% set ec time as reference
gdat_data.t = gdat_data.ec.t;
gdat_data.dim{1} = gdat_data.t;
gdat_data.dimunits{1} = 's';
gdat_data.units = 'W';
for i=1:length(gdat_data.gdat_params.source)
gdat_data.gdat_params.source{i} = lower(gdat_data.gdat_params.source{i});
if ~any(strmatch(gdat_data.gdat_params.source{i},lower(sources_avail)))
if gdat_data.gdat_params.nverbose>=1
warning(['source = ' gdat_data.gdat_params.source{i} ' not expected with data_request= ' data_request_eff])
end
end
end
end
% always start from ohmic so can use this time as base time since should yield full shot
% Ohmic
gdat_data.ohm.data = [];
gdat_data.ohm.units = '';
gdat_data.ohm.dim=gdat_data.dim;
gdat_data.ohm.dimunits=gdat_data.dimunits;
gdat_data.ohm.t=[];
gdat_data.ohm.x=[];
gdat_data.ohm.data_fullpath=[];
gdat_data.ohm.label='';
% get ohmic power simply from vloop*Ip (minus sign for JET)
ip=gdat([],'ip');
vloop=gdat([],'vloop');
tension = -1e5;
if isempty(gdat_data.t)
gdat_data.t = vloop.t;
gdat_data.dim{1} = gdat_data.t;
gdat_data.ohm.data = vloop.data;
gdat_data.dimunits{1} = 's';
gdat_data.units = 'W';
else
vloop_smooth=interpos(vloop.t,vloop.data,gdat_data.t,tension);
ip_t = interp1(ip.t,ip.data,gdat_data.t);
gdat_data.ohm.data = -vloop_smooth.*ip_t;
end
gdat_data.ohm.units = gdat_data.units;
gdat_data.ohm.dim=gdat_data.dim;
gdat_data.ohm.dimunits=gdat_data.dimunits;
gdat_data.ohm.t=gdat_data.t;
gdat_data.ohm.x=[];
gdat_data.ohm.data_fullpath=['-vloop(tens=' num2str(tension,'%.0e') ')*ip, from gdat'];
gdat_data.ohm.label='P_{OHM}';
% total power from each and total
gdat_data.data(:,1) = gdat_data.ohm.data;
if ~isempty(gdat_data.ec.data)
gdat_data.data(:,2) = gdat_data.ec.data(:,10);
gdat_data.data(:,3) = gdat_data.ec.data(:,10) + gdat_data.ohm.data;
gdat_data.dim{2} = [1:3];
gdat_data.dimunits{2} = 'Pohm;Pec;Ptot';
gdat_data.data_fullpath=['tot power from EC and ohm'];
gdat_data.label = 'P_{ohm};P_{EC};P_{tot}';
fields_to_copy = {'data','units','dim','dimunits','t','x','data_fullpath','label','help','gdat_params'};
fields_to_not_copy = {'shot','gdat_request'};
% total of each source in .data, but full data in subfield like pgyro in .ec, to check for nbi
params_eff = gdat_data.gdat_params;
params_eff.source=[]; % use default for individual calls
% ohmic, use its time-base
params_eff.data_request='p_ohmic';
try
ohm=gdat_jet(shot,params_eff);
catch
ohm.data = [];
ohm.dim = [];
end
if ~isempty(ohm.data) && ~isempty(ohm.dim)
for i=1:length(fields_to_copy)
if isfield(ohm,fields_to_copy{i})
gdat_data.ohm.(fields_to_copy{i}) = ohm.(fields_to_copy{i});
end
end
gdat_data.ohm.raw_data = gdat_data.ohm.data;
else
gdat_data.dim{2} = [1];
gdat_data.dimunits{2} = 'Pohm=Ptot';
gdat_data.data_fullpath=['tot power from ohm'];
gdat_data.label = 'P_{tot}=P_{ohm}';
if gdat_params.nverbose>=3; disp(['problems with ' params_eff.data_request]); end
return
end
mapping_for_jet.timedim = 1; mapping_for_jet.gdat_timedim = 1;
taus = -10;
%
% add each source in main.data, on ohm time array
gdat_data.t = gdat_data.ohm.t;
gdat_data.dim{1} = gdat_data.t;
gdat_data.dimunits{1} = 's';
gdat_data.ohm.data = interpos(gdat_data.t,gdat_data.ohm.raw_data,5.*taus);
gdat_data.data = reshape(gdat_data.ohm.data,length(gdat_data.t),1);
gdat_data.ohm.tension = 5.*taus;
gdat_data.x =[1];
gdat_data.label={'P_{ohm}'};
gdat_data.units = 'W';
%
if any(strmatch('nb',gdat_data.gdat_params.source))
% nbi
params_eff.data_request={'ppf','nbi','ptot'};
try
nbi=gdat_jet(shot,params_eff);
catch
end
if ~isempty(nbi.data) && ~isempty(nbi.dim)
for i=1:length(fields_to_copy)
if isfield(nbi,fields_to_copy{i})
gdat_data.nbi.(fields_to_copy{i}) = nbi.(fields_to_copy{i});
end
end
% 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(:,end),gdat_data.t);
gdat_data.x(end+1) =gdat_data.x(end)+1;
gdat_data.label{end+1}='P_{nbi}';
end
end
%
if any(strmatch('ic',gdat_data.gdat_params.source))
% ic
params_eff.data_request={'ppf','icrh','ptot'};
try
ic=gdat_jet(shot,params_eff);
catch
end
if ~isempty(ic.data) && ~isempty(ic.dim)
for i=1:length(fields_to_copy)
if isfield(ic,fields_to_copy{i})
gdat_data.ic.(fields_to_copy{i}) = ic.(fields_to_copy{i});
end
end
gdat_data.data(:,end+1) = interpos(-21,gdat_data.ic.t,gdat_data.ic.data,gdat_data.t);
gdat_data.x(end+1) =gdat_data.x(end)+1;
gdat_data.label{end+1}='P_{ic}';
end
end
if any(strmatch('rad',gdat_data.gdat_params.source))
% radiated power
params_eff.data_request='p_rad';
try
rad=gdat_jet(shot,params_eff);
catch
end
if ~isempty(rad.data) && ~isempty(rad.dim)
for i=1:length(fields_to_copy)
if isfield(rad,fields_to_copy{i})
gdat_data.rad.(fields_to_copy{i}) = rad.(fields_to_copy{i});
end
end
% add to main with linear interpolation and 0 for extrapolated values
gdat_data.data(:,end+1) = interpos(-21,gdat_data.rad.t,gdat_data.rad.data,gdat_data.t);
gdat_data.x(end+1) =gdat_data.x(end)+1;
gdat_data.label{end+1}='P_{rad}';
end
end
% add tot power
sources_coeff = [ 1, 1, 1, 0];
aa=sum(gdat_data.data.*repmat(reshape(sources_coeff,1,numel(sources_coeff)),size(gdat_data.data,1),1),2);
gdat_data.data(:,end+1) = aa;
% gdat_data.data(:,end+1) = sum(gdat_data.data,2);
gdat_data.label{end+1}='P_{tot} without Prad';
gdat_data.x(end+1) =gdat_data.x(end)+1;
gdat_data.dim{2} = gdat_data.x; gdat_data.dimunits{2} = '';
gdat_data.data_fullpath = 'tot powers from each sources, and total power in .data(:,end)';
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
case {'q_rho'}
......@@ -1128,6 +1307,39 @@ elseif strcmp(mapping_for_jet.method,'switchcase')
% add grids_1d to have rhotor, etc
gdat_data = get_grids_1d(gdat_data,1,1,gdat_params.nverbose);
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
case {'p_rad', 'prad'}
params_eff = gdat_data.gdat_params;
if ~isfield(params_eff,'source') || isempty(params_eff.source)
params_eff.source = 'bolo';
end
gdat_data.gdat_params = params_eff;
switch params_eff.source
case 'prad'
params_eff.data_request = {'ppf','prad','prad'};
rad=gdat_jet(shot,params_eff);
gdat_data.label={'P_{rad} from prad'};
otherwise
% 'bolo' by default
gdat_data.gdat_params.source = 'bolo';
params_eff.data_request = {'ppf','bolo','topi'};
rad=gdat_jet(shot,params_eff);
gdat_data.label={'P_{rad} from bolo/topi'};
params_eff.data_request = {'ppf','bolo','tobh'};
rad_bulk1=gdat_jet(shot,params_eff);
gdat_data.rad_bulk_h = rad_bulk1.data;
params_eff.data_request = {'ppf','bolo','tobu'};
rad_bulk2=gdat_jet(shot,params_eff);
gdat_data.rad_bulk_u = rad_bulk2.data;
gdat_data.rad_bulk_avg = 0.5.*(rad_bulk1.data+rad_bulk2.data);
end
gdat_data.t = rad.t;
gdat_data.dim{1} = gdat_data.t;
gdat_data.dimunits{1} = 's';
gdat_data.data = rad.data;
gdat_data.units = 'W';
gdat_data.data_fullpath = params_eff.data_request;
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
case {'psi_edge'}
% psi at edge, 0 by construction in Liuqe, thus not given
......
crpptbx/JET/jet_requests_mapping.m
+ 5
− 4
View file @ f8cdc498
......@@ -129,7 +129,7 @@ switch lower(data_request)
case 'mhd'
mapping.timedim = 1;
mapping.label = 'n=1 and n=2 amplitude';
mapping.method = 'expression';
mapping.method = 'expression'; % {'JPF','DA','C1M-H306'} {'JPF','DA','C1M-T009'} (big) {'JPF','DA','C1M-MIT27'}
mapping.expression = ['params_eff = gdat_data.gdat_params;params_eff.data_request={''JPF'',''DA'',''C1H-I301''}; ' ...
'gdat_tmp=gdat_jet(shot,params_eff);gdat_tmp.data=reshape(gdat_tmp.data,length(gdat_tmp.data),1 );' ...
'gdat_tmp.dim{1}=gdat_tmp.t;gdat_tmp.dim{2}=[1 2];gdat_tmp.x=gdat_tmp.dim{2};' ...
......@@ -203,7 +203,7 @@ switch lower(data_request)
mapping.timedim = 1;
mapping.method = 'signal';
mapping.expression = [{'PPF'},{'EQUI'},{'PHIT'}]; % note this is only chain2, so should check with efit...
case 'p_ohmic'
case {'p_ohmic', 'p_ohm', 'pohm'}
mapping.label = 'p\_ohmic';
mapping.timedim = 1;
mapping.method = 'signal';
......@@ -211,8 +211,9 @@ switch lower(data_request)
case {'p_rad', 'prad'}
mapping.label = 'p\_rad';
mapping.timedim = 1;
mapping.method = 'signal';
mapping.expression = [{'PPF'},{'PRAD'},{'PRAD'}];
% $$$ mapping.method = 'signal';
% $$$ mapping.expression = [{'PPF'},{'bolo'},{'topi'}];
mapping.method = 'switchcase';
case 'powers'
mapping.timedim = 1;
mapping.label = 'various powers';
......
crpptbx/gdat_plot.m
+ 1
− 1
View file @ f8cdc498
......@@ -91,7 +91,7 @@ if all(isfield(gdat_data,{'data','t'})) && ~isempty(gdat_data.data) && ~isempty(
if strcmp(gdat_data.gdat_request,'mhd')
% special case, add legend instead of ylabel
delete(hylabel);
legend(gdat_data.label,2);
legend(gdat_data.label);
% add spectrogram per signal
mhd_sum_data = 0.;
nfft=1024;
......
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment