diff --git a/crpptbx/AUG/aug_requests_mapping.m b/crpptbx/AUG/aug_requests_mapping.m index bbb93eb0d1ebcdb936b2b9978f8a9711e3c54423..393c495fddb1f35f2c667632f7109eafc353b233 100644 --- a/crpptbx/AUG/aug_requests_mapping.m +++ b/crpptbx/AUG/aug_requests_mapping.m @@ -175,8 +175,17 @@ switch lower(data_request) case 'ni' mapping.method = 'switchcase'; % especially since might have option fit, etc case 'powers' + mapping.timedim = 1; mapping.label = 'various powers'; mapping.method = 'switchcase'; + case 'psi_axis' + mapping.timedim = 1; + mapping.method = 'switchcase'; % there is psi_axis-psi_edge in FPG but otherwise complicated to get from equil, thus needs swticth case + mapping.label ='psi_\axis' ; + case 'psi_edge' + mapping.timedim = 1; + mapping.method = 'switchcase'; % is set to zero, so not in tree nodes + mapping.label = 'psi\_edge'; case 'q0' mapping.timedim = 1; mapping.label = 'q_0'; @@ -219,9 +228,18 @@ switch lower(data_request) mapping.timedim = 1; mapping.method = 'signal'; mapping.expression = [{'FPG'},{'Raus'},{'AUGD'}]; + case 'rhotor' + mapping.timedim = 2; + mapping.method = 'switchcase'; + mapping.label = 'rhotor\_norm'; + case 'rhotor_edge' + mapping.timedim = 1; + mapping.method = 'switchcase'; + mapping.label = 'rhotor\_edge'; case 'rhovol' + mapping.timedim = 2; mapping.label = 'rhovol\_norm'; - mapping.method = 'switchcase'; % from conf if exist otherwise computes it + mapping.method = 'switchcase'; case 'rmag' mapping.label = 'R\_magaxis'; mapping.timedim = 1; @@ -255,8 +273,12 @@ switch lower(data_request) case 'volume' mapping.label = 'Volume'; mapping.timedim = 1; - mapping.method = 'signal'; + mapping.method = 'switchcase'; mapping.expression = [{'FPG'},{'Vol'},{'AUGD'}]; + case 'volume_rho' + mapping.timedim = 2; + mapping.label = 'volume\_norm'; + mapping.method = 'switchcase'; case 'zeff' mapping.label = 'zeff from cxrs'; mapping.timedim = 1; diff --git a/crpptbx/AUG/gdat_aug.m b/crpptbx/AUG/gdat_aug.m index 33cf62c01b8e34f48c96fa34bf6058d07bc7e7f9..c382962a2211740cbad671b35acbf43907bab467 100644 --- a/crpptbx/AUG/gdat_aug.m +++ b/crpptbx/AUG/gdat_aug.m @@ -699,31 +699,8 @@ elseif strcmp(mapping_for_aug.method,'switchcase') %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% case {'equil'} - exp_name_eff = gdat_data.gdat_params.exp_name; - gdat_data.gdat_params.equil = upper(gdat_data.gdat_params.equil); - DIAG = gdat_data.gdat_params.equil; % 'EQI' by default at this stage, should become EQH? - M_Rmesh_par = sf2par(DIAG,shot,'M','PARMV'); - M_Rmesh = M_Rmesh_par.value + 1; % nb of points - N_Zmesh_par = sf2par(DIAG,shot,'N','PARMV'); - N_Zmesh = N_Zmesh_par.value + 1; % nb of points - NTIME_par = sf2par(DIAG,shot,'NTIME','PARMV'); - NTIME = NTIME_par.value; % nb of points - Lpf_par = rdaAUG_eff(shot,DIAG,'Lpf',exp_name_eff); - % since June, nb of time points in EQ results is not consistent with NTIME and time - % It seems the first NTIME points are correct, so use this explicitely - NTIME_Lpf = length(Lpf_par.value); - if (NTIME < NTIME_Lpf) - if gdat_params.nverbose>=3; disp('WARNING: nb of times points smaller then equil results, use first NTIME points'); end - elseif (NTIME > NTIME_Lpf) - if gdat_params.nverbose >= 1 - disp('ERROR: nb of times points LARGER then equil results') - disp('this is unexpected, so stop there and ask Olivier.Sauter@epfl.ch') - end - return - end - Lpf_tot = Lpf_par.value(1:NTIME); % nb of points: 100000*nb_SOL points + nb_core - Lpf_SOL = fix(Lpf_tot/100000); - Lpf1_t = mod(Lpf_tot,100000)+1; % nb of Lpf points + % get equil params and time array in gdat_data.t + [gdat_data,exp_name_eff,DIAG,NTIME_Lpf,NTIME,Lpf1_t,Lpf_SOL,M_Rmesh,N_Zmesh] = get_EQ_params(gdat_data); % since Lpf depends on time, need to load all first and then loop over time for easier mapping [qpsi,e]=rdaAUG_eff(shot,DIAG,'Qpsi',exp_name_eff); ndimrho = size(qpsi.data,2); @@ -879,9 +856,6 @@ elseif strcmp(mapping_for_aug.method,'switchcase') gdat_data.Rcoils=equil_Rcoil.value; [equil_Zcoil,e]=rdaAUG_eff(shot,DIAG,'Zcl',exp_name_eff); gdat_data.Zcoils=equil_Zcoil.value; - % get time values - [equil_time,e]=rdaAUG_eff(shot,DIAG,'time',exp_name_eff); - gdat_data.t = equil_time.value(1:NTIME); % [IpiPSI,e]=rdaAUG_eff(shot,DIAG,'IpiPSI',exp_name_eff); gdat_data.Ip = IpiPSI.value(1:NTIME); @@ -1093,32 +1067,99 @@ elseif strcmp(mapping_for_aug.method,'switchcase') end %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% - case {'q_rho'} - exp_name_eff = gdat_data.gdat_params.exp_name; - gdat_data.gdat_params.equil = upper(gdat_data.gdat_params.equil); - DIAG = gdat_data.gdat_params.equil; % 'EQI' by default at this stage, should become EQH? - M_Rmesh_par = sf2par(DIAG,shot,'M','PARMV'); - M_Rmesh = M_Rmesh_par.value + 1; % nb of points - N_Zmesh_par = sf2par(DIAG,shot,'N','PARMV'); - N_Zmesh = N_Zmesh_par.value + 1; % nb of points - NTIME_par = sf2par(DIAG,shot,'NTIME','PARMV'); - NTIME = NTIME_par.value; % nb of points - Lpf_par = rdaAUG_eff(shot,DIAG,'Lpf',exp_name_eff); - % since June, nb of time points in EQ results is not consistent with NTIME and time - % It seems the first NTIME points are correct, so use this explicitely - NTIME_Lpf = length(Lpf_par.value); - if (NTIME < NTIME_Lpf) - if gdat_params.nverbose>=3; disp('WARNING: nb of times points smaller then equil results, use first NTIME points'); end - elseif (NTIME > NTIME_Lpf) - if gdat_params.nverbose >= 1 - disp('ERROR: nb of times points LARGER then equil results') - disp('this is unexpected, so stop there and ask Olivier.Sauter@epfl.ch') - end + case {'powers'} + % case total only + params_eff = gdat_data.gdat_params; + % ohmic, use its time-base + params_eff.data_request={'TOT','P_OH'}; + try + gdat_data=gdat_aug(shot,params_eff); + catch + end + gdat_data.ohm.data = gdat_data.data; + gdat_data.ohm.dim = gdat_data.dim; + gdat_data.ohm.t = gdat_data.t; + gdat_data.ohm.label = gdat_data.label; + if isempty(gdat_data.data) || isempty(gdat_data.dim) + if gdat_params.nverbose>=3; disp(['problems with ' params_eff.data_request]); end return end - Lpf_tot = Lpf_par.value(1:NTIME); % nb of points: 100000*nb_SOL points + nb_core - Lpf_SOL = fix(Lpf_tot/100000); - Lpf1_t = mod(Lpf_tot,100000)+1; % nb of Lpf points + mapping_for_aug.timedim = 1; mapping_for_aug.gdat_timedim = 1; + taus = -10; + gdat_data.data = interpos(gdat_data.t,gdat_data.data,5.*taus); + gdat_data.data = reshape(gdat_data.data,length(gdat_data.t),1); + gdat_data.x =[1]; + gdat_data.label={'P_{ohmic}'}; + % nbi + params_eff.data_request={'NIS','PNI'}; + try + nbi=gdat_aug(shot,params_eff); + catch + end + if ~isempty(nbi.data) && ~isempty(nbi.dim) + gdat_data.nbi.data = nbi.data; + gdat_data.nbi.dim = nbi.dim; + gdat_data.nbi.t = nbi.t; + gdat_data.nbi.label = nbi.label; + gdat_data.data(:,end+1) = interpos(gdat_data.nbi.t,gdat_data.nbi.data,gdat_data.t,taus); + gdat_data.x(end+1) =gdat_data.x(end)+1; + gdat_data.label{end+1}='P_{nbi}'; + else + gdat_data.nbi.data = []; + gdat_data.nbi.dim = []; + gdat_data.nbi.t = []; + gdat_data.nbi.label = []; + end + % ic + params_eff.data_request={'ICP','PICRN'}; + try + ic=gdat_aug(shot,params_eff); + catch + end + if ~isempty(ic.data) && ~isempty(ic.dim) + gdat_data.ic.data = ic.data; + gdat_data.ic.dim = ic.dim; + gdat_data.ic.t = ic.t; + gdat_data.ic.label = ic.label; + gdat_data.data(:,end+1) = interpos(gdat_data.ic.t,gdat_data.ic.data,gdat_data.t,taus); + gdat_data.x(end+1) =gdat_data.x(end)+1; + gdat_data.label{end+1}='P_{ic}'; + else + gdat_data.ic.data = []; + gdat_data.ic.dim = []; + gdat_data.ic.t = []; + gdat_data.ic.label = []; + end + % ec + params_eff.data_request={'ECS','PECRH'}; + try + ec=gdat_aug(shot,params_eff); + catch + end + if ~isempty(ec.data) && ~isempty(ec.dim) + gdat_data.ec.data = ec.data; + gdat_data.ec.dim = ec.dim; + gdat_data.ec.t = ec.t; + gdat_data.ec.label = ec.label; + gdat_data.data(:,end+1) = interpos(gdat_data.ec.t,gdat_data.ec.data,gdat_data.t,taus); + gdat_data.x(end+1) =gdat_data.x(end)+1; + gdat_data.label{end+1}='P_{ec}'; + else + gdat_data.ec.data = []; + gdat_data.ec.dim = []; + gdat_data.ec.t = []; + gdat_data.ec.label = []; + end + % add tot power + gdat_data.data(:,end+1) = interpos(gdat_data.t,nansum(gdat_data.data,2),100*taus); + gdat_data.label{end+1}='P_{tot}'; + 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'} + [gdat_data,exp_name_eff,DIAG,NTIME_Lpf,NTIME,Lpf1_t,Lpf_SOL] = get_EQ_params(gdat_data); % since Lpf depends on time, need to load all first and then loop over time for easier mapping [qpsi,e]=rdaAUG_eff(shot,DIAG,'Qpsi',exp_name_eff); ndimrho = size(qpsi.data,2); @@ -1142,6 +1183,15 @@ elseif strcmp(mapping_for_aug.method,'switchcase') max_qValue = 30.; % Note could just put a NaN on LCFS value since ill-defined when diverted for it=1:NTIME Lpf1 = Lpf1_t(it); + % Qpsi and similar data is on (time,radius) with radius being: LCFS..Lpf_points dummy LCFS..SOL part + % change it to (radial,time) and use only Lpf+1 points up to LCFS + ijok=find(qpsi.value(:,1)); % note: eqr fills in only odd points radially + % set NaNs to zeroes + if qpsi.value(ijok(1),1)<0 + gdat_data.qvalue(:,it) = max(qpsi.value(it,Lpf1:-1:1)',-max_qValue); + else + gdat_data.qvalue(:,it) = min(qpsi.value(it,Lpf1:-1:1)',max_qValue); + end % get x values psi_it=psi_tree.value(it,Lpf1:-1:1)'; gdat_data.psi_axis(it)= psi_it(1); @@ -1191,8 +1241,6 @@ elseif strcmp(mapping_for_aug.method,'switchcase') end gdat_data.x = gdat_data.rhopolnorm; % get time values - [equil_time,e]=rdaAUG_eff(shot,DIAG,'time',exp_name_eff); - gdat_data.t = equil_time.value(1:NTIME); gdat_data.data = gdat_data.qvalue; % put q in data gdat_data.units=[]; % not applicable gdat_data.data_fullpath = [DIAG ' from expname: ' gdat_data.gdat_params.exp_name '; q(rhopolnorm,t) in .data(.x,.t)']; @@ -1202,6 +1250,86 @@ elseif strcmp(mapping_for_aug.method,'switchcase') gdat_data.dimunits{1} = 'rhopolnorm'; gdat_data.dimunits{2} = 'time [s]'; + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + case {'psi_axis', 'psi_edge'} + if strcmp(upper(gdat_data.gdat_params.equil),'FPG') + gdat_params_prev = gdat_data.gdat_params; gdat_params_eff = gdat_params_prev; + gdat_params_eff.data_request = [{'FPG'},{'fax-bnd'},{'AUGD'}]; + gdat_data = gdat_aug(gdat_data.shot,gdat_params_eff); + gdat_data.gdat_params = gdat_params_prev; + gdat_data.label = 'psi\_axis-psi\_edge'; + gdat_data.data_fullpath = [gdat_data.data_fullpath ' yields only psi\_axis-psi\_edge from FPG/fax-bnd']; + else + [gdat_data,exp_name_eff,DIAG,NTIME_Lpf,NTIME,Lpf1_t,Lpf_SOL] = get_EQ_params(gdat_data); + % since Lpf depends on time, need to load all first and then loop over time for easier mapping + [psi_tree,e]=rdaAUG_eff(shot,DIAG,'PFL',exp_name_eff); + ndimrho = size(psi_tree.data,2); + if ndimrho==NTIME_Lpf + % data seems to be transposed + ndimrho = size(psi_tree.data,1); + itotransposeback = 1; % seems x,time inverted so transpose and exchange .x and .t + else + itotransposeback = 0; + end + psi_tree=adapt_rda(psi_tree,NTIME,ndimrho,itotransposeback); + ijnan=find(isnan(psi_tree.value)); + psi_tree.value(ijnan)=0; + gdat_data.dim{1} = gdat_data.t; + gdat_data.dimunits{1} = 'time [s]'; + for it=1:NTIME + Lpf1 = Lpf1_t(it); + psi_it=psi_tree.value(it,Lpf1:-1:1)'; + if strcmp(data_request_eff,'psi_axis') + gdat_data.data(it)= psi_it(1); + elseif strcmp(data_request_eff,'psi_edge') + gdat_data.data(it)= psi_it(end); + else + end + end + gdat_data.units = psi_tree.units; + gdat_data.data_fullpath = [DIAG '/PFL extract ' data_request_eff]; + end + + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% + case {'rho_tor', 'rhotor_edge', 'rhovol', 'volume', 'volume_rho'} + if strcmp(upper(gdat_data.gdat_params.equil),'FPG') + gdat_params_prev = gdat_data.gdat_params; gdat_params_eff = gdat_params_prev; + gdat_params_eff.data_request = [{'FPG'},{'fax-bnd'},{'AUGD'}]; + gdat_data = gdat_aug(gdat_data.shot,gdat_params_eff); + gdat_data.gdat_params = gdat_params_prev; + gdat_data.label = 'psi\_axis-psi\_edge'; + gdat_data.data_fullpath = [gdat_data.data_fullpath ' yields only psi\_axis-psi\_edge from FPG/fax-bnd']; + else + [gdat_data,exp_name_eff,DIAG,NTIME_Lpf,NTIME,Lpf1_t,Lpf_SOL] = get_EQ_params(gdat_data); + % since Lpf depends on time, need to load all first and then loop over time for easier mapping + [psi_tree,e]=rdaAUG_eff(shot,DIAG,'PFL',exp_name_eff); + ndimrho = size(psi_tree.data,2); + if ndimrho==NTIME_Lpf + % data seems to be transposed + ndimrho = size(psi_tree.data,1); + itotransposeback = 1; % seems x,time inverted so transpose and exchange .x and .t + else + itotransposeback = 0; + end + psi_tree=adapt_rda(psi_tree,NTIME,ndimrho,itotransposeback); + ijnan=find(isnan(psi_tree.value)); + psi_tree.value(ijnan)=0; + gdat_data.dim{1} = gdat_data.t; + gdat_data.dimunits{1} = 'time [s]'; + for it=1:NTIME + Lpf1 = Lpf1_t(it); + psi_it=psi_tree.value(it,Lpf1:-1:1)'; + if strcmp(data_request_eff,'psi_axis') + gdat_data.data(it)= psi_it(1); + elseif strcmp(data_request_eff,'psi_edge') + gdat_data.data(it)= psi_it(end); + else + end + end + gdat_data.units = psi_tree.units; + gdat_data.data_fullpath = [DIAG '/PFL extract ' data_request_eff]; + end + %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% case {'sxr'} % sxr from B by default or else if 'source','else' is provided @@ -1338,3 +1466,40 @@ error_status=0; return + +%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%% +% functions for portions called several times + +function [gdat_data,exp_name_eff,DIAG,NTIME_Lpf,NTIME,Lpf1_t,Lpf_SOL,M_Rmesh,N_Zmesh] = get_EQ_params(gdat_data); +% get basic params to be able to read results in EQ-like shotfiles +% M_Rmesh,N_Zmesh only needed for equil when 2D quantities are required +% +shot=gdat_data.shot; +exp_name_eff = gdat_data.gdat_params.exp_name; +gdat_data.gdat_params.equil = upper(gdat_data.gdat_params.equil); +DIAG = gdat_data.gdat_params.equil; % 'EQI' by default at this stage, should become EQH? +M_Rmesh_par = sf2par(DIAG,shot,'M','PARMV'); +M_Rmesh = M_Rmesh_par.value + 1; % nb of points +N_Zmesh_par = sf2par(DIAG,shot,'N','PARMV'); +N_Zmesh = N_Zmesh_par.value + 1; % nb of points +NTIME_par = sf2par(DIAG,shot,'NTIME','PARMV'); +NTIME = NTIME_par.value; % nb of points +Lpf_par = rdaAUG_eff(shot,DIAG,'Lpf',exp_name_eff); +% since June, nb of time points in EQ results is not consistent with NTIME and time +% It seems the first NTIME points are correct, so use this explicitely +NTIME_Lpf = length(Lpf_par.value); +if (NTIME < NTIME_Lpf) + if gdat_data.gdat_params.nverbose>=3; disp('WARNING: nb of times points smaller then equil results, use first NTIME points'); end +elseif (NTIME > NTIME_Lpf) + if gdat_data.gdat_params.nverbose >= 1 + disp('ERROR: nb of times points LARGER then equil results') + disp('this is unexpected, so stop there and ask Olivier.Sauter@epfl.ch') + end + return +end +Lpf_tot = Lpf_par.value(1:NTIME); % nb of points: 100000*nb_SOL points + nb_core +Lpf_SOL = fix(Lpf_tot/100000); +Lpf1_t = mod(Lpf_tot,100000)+1; % nb of Lpf points + +[equil_time,e]=rdaAUG_eff(shot,DIAG,'time',exp_name_eff); +gdat_data.t = equil_time.value(1:NTIME); diff --git a/crpptbx/AUG/geteqdskAUG.m b/crpptbx/AUG/geteqdskAUG.m index 7a661a18e9bd4520c294ab915e5051e163e9f4f2..ac5731f4e3467e54dbbae693620b132033df0051 100644 --- a/crpptbx/AUG/geteqdskAUG.m +++ b/crpptbx/AUG/geteqdskAUG.m @@ -91,26 +91,32 @@ hold on psiedge_eff = 0.001*eqdsk.psiaxis + 0.999*eqdsk.psiedge; [hh1 hh2]=contour(eqdsk.rmesh,eqdsk.zmesh,eqdsk.psi',[psiedge_eff psiedge_eff],'k'); axis equal -ij=1; -ij_prev = 0; -nbhh1(ij)=hh1(2,ij_prev+1); -Rbnd{ij}=hh1(1,2:1+nbhh1(ij)); -Zbnd{ij}=hh1(2,2:1+nbhh1(ij)); -ij_prev = ij_prev+1+nbhh1(ij); -while (ij_prev+1<size(hh1,2)) - % next - ij=ij + 1; +if ~isempty(hh1) + ij=1; + ij_prev = 0; nbhh1(ij)=hh1(2,ij_prev+1); - Rbnd{ij}=hh1(1,ij_prev+2:ij_prev+1+nbhh1(ij)); - Zbnd{ij}=hh1(2,ij_prev+2:ij_prev+1+nbhh1(ij)); + Rbnd{ij}=hh1(1,2:1+nbhh1(ij)); + Zbnd{ij}=hh1(2,2:1+nbhh1(ij)); ij_prev = ij_prev+1+nbhh1(ij); + while (ij_prev+1<size(hh1,2)) + % next + ij=ij + 1; + nbhh1(ij)=hh1(2,ij_prev+1); + Rbnd{ij}=hh1(1,ij_prev+2:ij_prev+1+nbhh1(ij)); + Zbnd{ij}=hh1(2,ij_prev+2:ij_prev+1+nbhh1(ij)); + ij_prev = ij_prev+1+nbhh1(ij); + end + % assume LCFS with most points + [zzz irz]=max(nbhh1); + eqdsk.nbbound = nbhh1(irz); + eqdsk.rplas = Rbnd{irz}'; + eqdsk.zplas = Zbnd{irz}'; + plot(eqdsk.rplas,eqdsk.zplas,'k-') +else + eqdsk.nbbound = []; + eqdsk.rplas = []; + eqdsk.zplas = []; end -% assume LCFS with most points -[zzz irz]=max(nbhh1); -eqdsk.nbbound = nbhh1(irz); -eqdsk.rplas = Rbnd{irz}'; -eqdsk.zplas = Zbnd{irz}'; -plot(eqdsk.rplas,eqdsk.zplas,'k-') [aget]=which('geteqdskAUG'); [path1,name2,ext3]=fileparts(aget); diff --git a/crpptbx/gdat_plot.m b/crpptbx/gdat_plot.m index ad09d0063b137b547f3546d6d97527ede0fce985..147f1c46f4852e0c14bc41121c86336ca74b1c14 100644 --- a/crpptbx/gdat_plot.m +++ b/crpptbx/gdat_plot.m @@ -89,7 +89,7 @@ if prod(isfield(gdat_data,{'data','t'})) && ~isempty(gdat_data.data) && ~isempty imagesc(T+tmhdm(1),F/1e3,20*log10(abs(B)));axis xy;colormap jet; ylabel('freq') xlabel(gdat_data.dimunits{1}) - title(gdat_data.label{i}) + title([upper(gdat_data.gdat_params.machine) '#' num2str(gdat_data.shot) ' ' gdat_data.label{i}]) mhd_sum_data = mhd_sum_data + gdat_data.data(:,i); end [B,F,T]=specgram(mhd_sum_data./size(gdat_data.data,2),nfft,1/mean(diff(gdat_data.t)),hanning(nfft),nfft/2); @@ -97,7 +97,7 @@ if prod(isfield(gdat_data,{'data','t'})) && ~isempty(gdat_data.data) && ~isempty imagesc(T+tmhdm(1),F/1e3,20*log10(abs(B)));axis xy;colormap jet; ylabel('freq') xlabel(gdat_data.dimunits{1}) - title(['sum of ' gdat_data.label{:}]) + title([upper(gdat_data.gdat_params.machine) '#' num2str(gdat_data.shot) ' sum of ' gdat_data.label{:}]) end else disp('cannot plot gdat_data, has empty data or t field')