diff --git a/crpptbx/AUG/aug_requests_mapping.m b/crpptbx/AUG/aug_requests_mapping.m index 17dc0dffcc70c9fe6fb99e004ee763b329a881f0..e92eeaa3edab890fb700b782f3a722d338c3e25f 100644 --- a/crpptbx/AUG/aug_requests_mapping.m +++ b/crpptbx/AUG/aug_requests_mapping.m @@ -101,7 +101,7 @@ switch lower(data_request) mapping.timedim = 1; mapping.method = 'signal'; mapping.expression = [{'FPG'},{'delRuntn'}]; - case 'ece' + case {'ece', 'eced', 'ece_rho', 'eced_rho'} mapping.timedim = 2; mapping.method = 'switchcase'; mapping.expression = ''; diff --git a/crpptbx/AUG/gdat_aug.m b/crpptbx/AUG/gdat_aug.m index 083a28ea4c4fc8e72a890f72f00ec4475bb3433b..eb2bf0cc8071e21ae72fee1d7a68bb6345d16cd9 100644 --- a/crpptbx/AUG/gdat_aug.m +++ b/crpptbx/AUG/gdat_aug.m @@ -382,6 +382,233 @@ elseif strcmp(mapping_for_aug.method,'expression') %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% elseif strcmp(mapping_for_aug.method,'switchcase') switch data_request_eff % not lower(...) since data_request_eff should be lower case already at this stage + case {'cxrs', 'cxrs_rho'} + exp_name_eff = gdat_data.gdat_params.exp_name; + diag_name = 'CEZ'; + % R, Z positions of measurements + [r_time]=sf2ab(diag_name,shot,'R_time','-exp',exp_name_eff); + gdat_data.r = r_time.value{1}; + inotok=find(gdat_data.r<=0); + gdat_data.r(inotok) = NaN; + [z_time]=sf2ab(diag_name,shot,'z_time','-exp',exp_name_eff); + gdat_data.z = z_time.value{1}; + inotok=find(gdat_data.z<=0); + gdat_data.z(inotok) = NaN; + [time]=sf2tb(diag_name,shot,'time','-exp',exp_name_eff); + gdat_data.t = time.value; + gdat_data.dim{1} = {gdat_data.r , gdat_data.z}; + gdat_data.dimunits{1} = 'R, Z [m]'; + gdat_data.dim{2} = gdat_data.t; + gdat_data.dimunits{2} = 't [s]'; + gdat_data.x = gdat_data.dim{1}; + % vrot + [a,error_status]=rdaAUG_eff(shot,diag_name,'vrot',exp_name_eff); + if isempty(a.data) || isempty(a.t) || error_status>0 + if nverbose>=3; + a + disp(['with data_request= ' data_request_eff]) + end + end + gdat_data.vrot.data = a.data; + if any(strcmp(fieldnames(a),'units')); gdat_data.vrot.units=a.units; end + if any(strcmp(fieldnames(a),'name')); gdat_data.vrot.name=a.name; end + gdat_data.vrot.label = 'vrot_tor'; + [aerr,e]=rdaAUG_eff(shot,diag_name,'err_vrot',exp_name_eff); + gdat_data.vrot.error_bar = aerr.data; + % Ti + % [a,e]=rdaAUG_eff(shot,diag_name,'Ti',exp_name_eff); + % [aerr,e]=rdaAUG_eff(shot,diag_name,'err_Ti',exp_name_eff); + [a,e]=rdaAUG_eff(shot,diag_name,'Ti_c',exp_name_eff); + [aerr,e]=rdaAUG_eff(shot,diag_name,'err_Ti_c',exp_name_eff); + gdat_data.ti.data = a.data; + gdat_data.data = a.data; + if any(strcmp(fieldnames(a),'units')); gdat_data.ti.units=a.units; end + if any(strcmp(fieldnames(a),'units')); gdat_data.units=a.units; end + if any(strcmp(fieldnames(a),'name')); gdat_data.ti.name=a.name; end + gdat_data.ti.error_bar = aerr.data; + gdat_data.error_bar = aerr.data; + gdat_data.data_fullpath=[exp_name_eff '/' diag_name '/' a.name ';vrot, Ti in data, {r;z} in .x']; + % + if strcmp(data_request_eff,'cxrs_rho') + params_equil = gdat_data.gdat_params; + params_equil.data_request = 'equil'; + [equil,params_equil,error_status] = gdat_aug(shot,params_equil); + if error_status>0 + if nverbose>=3; disp(['problems with ' params_equil.data_request]); end + return + end + gdat_data.gdat_params.equil = params_equil.equil; + gdat_data.equil = equil; + inb_chord_cxrs=size(gdat_data.r,1); + inb_time_cxrs=size(gdat_data.r,2); + psi_out = NaN*ones(inb_chord_cxrs,inb_time_cxrs); + rhopsinorm_out = NaN*ones(inb_chord_cxrs,inb_time_cxrs); + rhotornorm_out = NaN*ones(inb_chord_cxrs,inb_time_cxrs); + rhovolnorm_out = NaN*ones(inb_chord_cxrs,inb_time_cxrs); + % constructs intervals within which a given equil is used: [time_equil(i),time_equil(i+1)] + time_equil=[1.5*equil.t(1)-0.5*equil.t(2) 0.5.*(equil.t(1:end-1)+equil.t(2:end)) 1.5*equil.t(end)-0.5*equil.t(end-1)]; + iok=find(~isnan(gdat_data.r(:,1))); + for itequil=1:length(time_equil)-1 + rr=equil.Rmesh(:,itequil); + zz=equil.Zmesh(:,itequil); + psirz_in = equil.psi2D(:,:,itequil); + it_cxrs_inequil = find(gdat_data.t>=time_equil(itequil) & gdat_data.t<=time_equil(itequil+1)); + if ~isempty(it_cxrs_inequil) + rout=gdat_data.r(iok,it_cxrs_inequil); + zout=gdat_data.z(iok,it_cxrs_inequil); + psi_at_routzout = interpos2Dcartesian(rr,zz,psirz_in,rout,zout); + psi_out(iok,it_cxrs_inequil) = reshape(psi_at_routzout,length(iok),length(it_cxrs_inequil)); + rhopsinorm_out(iok,it_cxrs_inequil) = sqrt((psi_out(iok,it_cxrs_inequil)-equil.psi_axis(itequil))./(equil.psi_lcfs(itequil)-equil.psi_axis(itequil))); + for it_cx=1:length(it_cxrs_inequil) + rhotornorm_out(iok,it_cxrs_inequil(it_cx)) = interpos(equil.rhopolnorm(:,itequil),equil.rhotornorm(:,itequil),rhopsinorm_out(iok,it_cxrs_inequil(it_cx)),-3,[2 2],[0 1]); + rhovolnorm_out(iok,it_cxrs_inequil(it_cx)) = interpos(equil.rhopolnorm(:,itequil),equil.rhovolnorm(:,itequil),rhopsinorm_out(iok,it_cxrs_inequil(it_cx)),-3,[2 2],[0 1]); + end + end + end + gdat_data.psi = psi_out; + gdat_data.rhopsinorm = rhopsinorm_out; + gdat_data.rhotornorm = rhotornorm_out; + gdat_data.rhovolnorm = rhovolnorm_out; + end + + case {'ece', 'eced', 'ece_rho', 'eced_rho'} + nth_points = 13; + if isfield(gdat_data.gdat_params,'nth_points') && ~isempty(gdat_data.gdat_params.nth_points) + nth_points = gdat_data.gdat_params.nth_points; + else + gdat_data.gdat_params.nth_points = nth_points; + end + channels = -1; + if isfield(gdat_data.gdat_params,'channels') && ~isempty(gdat_data.gdat_params.channels) + channels = gdat_data.gdat_params.channels; + end + if nth_points>=10 + match_rz_to_time = 1; + else + match_rz_to_time = 0; + end + if isfield(gdat_data.gdat_params,'match_rz_to_time') && ~isempty(gdat_data.gdat_params.match_rz_to_time) + match_rz_to_time = gdat_data.gdat_params.match_rz_to_time; + else + gdat_data.gdat_params.match_rz_to_time = match_rz_to_time; + end + time_interval = [-Inf +Inf]; + if isfield(gdat_data.gdat_params,'time_interval') && ~isempty(gdat_data.gdat_params.time_interval) + time_interval = gdat_data.gdat_params.time_interval; + else + gdat_data.gdat_params.time_interval = time_interval; + end + % + if isfield(gdat_data.gdat_params,'diag_name') && ~isempty(gdat_data.gdat_params.diag_name) + diag_name = gdat_data.gdat_params.diag_name; + if strcmp(upper(diag_name),'RMD'); gdat_data.gdat_params.exp_name = 'ECED'; end + else + diag_name = 'CEC'; + gdat_data.gdat_params.diag_name = diag_name; + end + exp_name_eff = gdat_data.gdat_params.exp_name; + if strcmp(data_request_eff,'eced') + exp_name_eff = 'ECED'; + gdat_data.gdat_params.exp_name = exp_name_eff; + end + [a,e]=rdaAUG_eff(shot,diag_name,'Trad-A',exp_name_eff,time_interval); + % [a,e]=rdaAUG_eff(shot,diag_name,'Trad-A',exp_name_eff); + inb_chord = size(a.data,1); + if channels(1)<=0 + channels = [1:inb_chord]; + end + gdat_data.dim{1} = channels; + gdat_data.gdat_params.channels = channels; + if nth_points>1 + gdat_data.data = a.data(channels,1:nth_points:end); + gdat_data.dim{2} = a.t(1:nth_points:end); + else + gdat_data.data = a.data(channels,:); + gdat_data.dim{2} = a.t; + end + gdat_data.x = gdat_data.dim{1}; + gdat_data.t = gdat_data.dim{2}; + gdat_data.dimunits=[{'channels'} ; {'time [s]'}]; + if any(strcmp(fieldnames(a),'units')); gdat_data.units=a.units; end + try + [aR,e]=rdaAUG_eff(shot,diag_name,'R-A',exp_name_eff,time_interval); + catch + end + try + [aZ,e]=rdaAUG_eff(shot,diag_name,'z-A',exp_name_eff,time_interval); + catch + disp(['problem with getting z-A in ' diag_name]) + end + if match_rz_to_time + % interpolate R structure on ece data time array, to ease plot vs R + for i=1:length(channels) + radius.data(i,:) = interp1(aR.t,aR.data(channels(i),:),gdat_data.t); + zheight.data(i,:) = interp1(aZ.t,aZ.data(channels(i),:),gdat_data.t); + end + radiuszheight.t=gdat_data.t; + else + radius.data = aR.data(channels,:); + radiuszheight.t=aR.t; + zheight.data = aZ.data(channels,:); + end + ij=find(gdat_data.data<=0); + gdat_data.data(ij)=NaN; + gdat_data.rz.r=radius.data; + gdat_data.rz.z=zheight.data; + gdat_data.rz.t = radiuszheight.t; + gdat_data.data_fullpath = [exp_name_eff '/' diag_name '/Trad-A with R, Z in .r,.z']; + + if strcmp(data_request_eff,'ece_rho') || strcmp(data_request_eff,'eced_rho') + params_equil = gdat_data.gdat_params; + params_equil.data_request = 'equil'; + [equil,params_equil,error_status] = gdat_aug(shot,params_equil); + if error_status>0 + if nverbose>=3; disp(['problems with ' params_equil.data_request]); end + return + end + gdat_data.gdat_params.equil = params_equil.equil; + gdat_data.equil = equil; + gdat_data.data_fullpath = [exp_name_eff '/' diag_name '/Trad-A with .r,.z projected on equil ' gdat_data.gdat_params.equil ' in .rhos']; + inb_chord_ece=size(gdat_data.rz.r,1); + inb_time_ece=size(gdat_data.rz.r,2); + psi_out = NaN*ones(inb_chord_ece,inb_time_ece); + rhopsinorm_out = NaN*ones(inb_chord_ece,inb_time_ece); + rhotornorm_out = NaN*ones(inb_chord_ece,inb_time_ece); + rhovolnorm_out = NaN*ones(inb_chord_ece,inb_time_ece); + % constructs intervals within which a given equil is used: [time_equil(i),time_equil(i+1)] + time_equil=[1.5*equil.t(1)-0.5*equil.t(2) 0.5.*(equil.t(1:end-1)+equil.t(2:end)) 1.5*equil.t(end)-0.5*equil.t(end-1)]; + for itequil=1:length(time_equil)-1 + rr=equil.Rmesh(:,itequil); + zz=equil.Zmesh(:,itequil); + psirz_in = equil.psi2D(:,:,itequil); + it_ece_inequil = find(gdat_data.rz.t>=time_equil(itequil) & gdat_data.rz.t<=time_equil(itequil+1)); + if ~isempty(it_ece_inequil) + r_it_ece_inequil = gdat_data.rz.r(:,it_ece_inequil); + iok = find(~isnan(r_it_ece_inequil) & r_it_ece_inequil>0); + if ~isempty(iok) + rout=r_it_ece_inequil(iok); + z_it_ece_inequil = gdat_data.rz.z(:,it_ece_inequil); + zout=z_it_ece_inequil(iok); + psi_at_routzout = interpos2Dcartesian(rr,zz,psirz_in,rout,zout); + [ieff,jeff]=ind2sub(size(gdat_data.rz.r(:,it_ece_inequil)),iok); + for ij=1:length(iok) + psi_out(ieff(ij),it_ece_inequil(jeff(ij))) = psi_at_routzout(ij); + end + rhopsinorm_out(:,it_ece_inequil) = sqrt(abs((psi_out(:,it_ece_inequil)-equil.psi_axis(itequil))./(equil.psi_lcfs(itequil)-equil.psi_axis(itequil)))); + for it_cx=1:length(it_ece_inequil) + rhotornorm_out(:,it_ece_inequil(it_cx)) = interpos(equil.rhopolnorm(:,itequil),equil.rhotornorm(:,itequil),rhopsinorm_out(:,it_ece_inequil(it_cx)),-3,[2 2],[0 1]); + rhovolnorm_out(:,it_ece_inequil(it_cx)) = interpos(equil.rhopolnorm(:,itequil),equil.rhovolnorm(:,itequil),rhopsinorm_out(:,it_ece_inequil(it_cx)),-3,[2 2],[0 1]); + end + end + end + end + gdat_data.rhos.psi = psi_out; + gdat_data.rhos.rhopsinorm = rhopsinorm_out; + gdat_data.rhos.rhotornorm = rhotornorm_out; + gdat_data.rhos.rhovolnorm = rhovolnorm_out; + gdat_data.rhos.t = gdat_data.rz.t; + end + case {'eqdsk'} % time_eqdsks=1.; % default time @@ -735,6 +962,7 @@ elseif strcmp(mapping_for_aug.method,'switchcase') return end gdat_data.gdat_params.equil = params_equil.equil; + %gdat_data.equil = equil; % core node_child_nameeff = [upper(data_request_eff(1)) 'e_core']; diff --git a/crpptbx/gdat_plot.m b/crpptbx/gdat_plot.m index 2997459709f6cc0d3f06edc787d175e56cff7003..d75dc2fa6fd31ef3f92cdef3bcd32c8187f8951d 100644 --- a/crpptbx/gdat_plot.m +++ b/crpptbx/gdat_plot.m @@ -7,10 +7,34 @@ function [fighandle]=gdat_plot(gdat_data,varargin); % > 1: create new figure with this number, adding clf % <-1: add to figure number abs(doplot) (with hold all) % -if ~isfield(gdat_data.gdat_params,'doplot') || gdat_data.gdat_params.doplot ==0 + +fighandle = NaN; +gdat_plot_params.dummy=NaN; +if mod(nargin-1,2)==0 + for i=1:2:nargin-1 + if ischar(varargin{i}) + % enforce lower case for any character driven input + gdat_plot_params.(lower(varargin{i})) = varargin{i+1}; + else + warning(['input argument nb: ' num2str(i) ' is incorrect, expects a character string']) + error_status=101; + return + end + end +else + warning('number of input arguments incorrect, cannot make pairs of parameters') + error_status=102; return end +doplot=0; +if isfield(gdat_plot_params,'doplot') && ~isempty(gdat_plot_params.doplot) + doplot = gdat_plot_params.doplot; +elseif isfield(gdat_data.gdat_params,'doplot') || ~isempty(gdat_data.gdat_params.doplot) + doplot = gdat_data.gdat_params.doplot; +end +if doplot==0; return; end + if prod(isfield(gdat_data,{'data','t'})) && ~isempty(gdat_data.data) && ~isempty(gdat_data.t) fighandle = get(0,'CurrentFigure'); if gdat_data.gdat_params.doplot == 1 diff --git a/crpptbx/test_all_requestnames.m b/crpptbx/test_all_requestnames.m new file mode 100644 index 0000000000000000000000000000000000000000..dca619433a0b5e1cc0c8c6eb84391c72f074ca22 --- /dev/null +++ b/crpptbx/test_all_requestnames.m @@ -0,0 +1,13 @@ + +machine='AUG'; +shot=30594; +aa=gdat('machine',machine); +aa=gdat; +all_request_names = aa.gdat_request +%break +istart=1; +for irequest=istart:length(all_request_names) + request=all_request_names{irequest} + ab{irequest} = gdat(shot,request,'machine',machine,'doplot',1); + pause +end