From 8841acaea3dea573ac7b277406ff900b9dce8917 Mon Sep 17 00:00:00 2001
From: Olivier Sauter <olivier.sauter@epfl.ch>
Date: Tue, 29 Mar 2016 13:28:26 +0000
Subject: [PATCH] add source option to powers and pgyro from within powers in
 .ec substructure

git-svn-id: https://spcsvn.epfl.ch/repos/TCV/gdat/trunk@5594 d63d8f72-b253-0410-a779-e742ad2e26cf
---
 crpptbx/AUG/gdat_aug.m | 191 ++++++++++++++++++++++++++---------------
 crpptbx/TCV/gdat_tcv.m |   5 +-
 2 files changed, 126 insertions(+), 70 deletions(-)

diff --git a/crpptbx/AUG/gdat_aug.m b/crpptbx/AUG/gdat_aug.m
index 82d7afbd..c175a44e 100644
--- a/crpptbx/AUG/gdat_aug.m
+++ b/crpptbx/AUG/gdat_aug.m
@@ -1289,8 +1289,11 @@ elseif strcmp(mapping_for_aug.method,'switchcase')
     % pgyro tot in index=9
     try
       gdat_data=gdat_aug(shot,params_eff);
+      gdat_data.data_request = data_request_eff;
+      gdat_data.gdat_params.data_request = data_request_eff;
     catch
       if gdat_params.nverbose>=3; disp(['problems with ' params_eff.data_request]); end
+      gdat_data.data_request = data_request_eff;
       return
     end
     nb_timepoints = length(gdat_data.t);
@@ -1435,90 +1438,142 @@ elseif strcmp(mapping_for_aug.method,'switchcase')
 
     %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    case {'powers'}
-    % case total only
+    sources_avail = {'ohm','ec','nbi','ic'}; % 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
+      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
+
+    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;
     % ohmic, use its time-base
     params_eff.data_request={'TOT','P_OH'};
     try
-      gdat_data=gdat_aug(shot,params_eff);
+      ohm=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 ~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
       if gdat_params.nverbose>=3; disp(['problems with ' params_eff.data_request]); end
       return
     end
     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);
+    %
+    % 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_{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 = [];
+    gdat_data.label={'P_{ohm}'};
+    gdat_data.units = 'W';
+    %
+    if any(strmatch('ec',gdat_data.gdat_params.source))
+      % ec
+      params_eff.data_request={'ECS','PECRH'};
+      params_eff.data_request='pgyro';
+      try
+	ec=gdat_aug(shot,params_eff);
+      catch
+      end
+      if ~isempty(ec.data) && ~isempty(ec.dim)
+	for i=1:length(fields_to_copy)
+	  % if has pgyro, use not_copy
+	  if isfield(ec,fields_to_copy{i}) && ~any(strmatch(fields_to_not_copy,fields_to_copy{i}))
+	    gdat_data.ec.(fields_to_copy{i}) = ec.(fields_to_copy{i});
+	  end
+	end
+	gdat_data.data(:,end+1) = interpos(-21,gdat_data.ec.t,gdat_data.ec.data(:,end),gdat_data.t);
+	gdat_data.x(end+1) =gdat_data.x(end)+1;
+	gdat_data.label{end+1}='P_{ec}';
+      end
     end
-    % ec
-    params_eff.data_request={'ECS','PECRH'};
-    try
-      ec=gdat_aug(shot,params_eff);
-    catch
+    %
+    if any(strmatch('nb',gdat_data.gdat_params.source))
+      % nbi
+      params_eff.data_request={'NIS','PNI'};
+      try
+	nbi=gdat_aug(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,gdat_data.t);
+	gdat_data.x(end+1) =gdat_data.x(end)+1;
+	gdat_data.label{end+1}='P_{nbi}';
+      end
     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 = [];
+    %
+    if any(strmatch('ic',gdat_data.gdat_params.source))
+      % ic
+      params_eff.data_request={'ICP','PICRN'};
+      try
+	ic=gdat_aug(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
     % add tot power
-    gdat_data.data(:,end+1) = interpos(gdat_data.t,nansum(gdat_data.data,2),100*taus);
+    gdat_data.data(:,end+1) = sum(gdat_data.data,2);
     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} = '';
diff --git a/crpptbx/TCV/gdat_tcv.m b/crpptbx/TCV/gdat_tcv.m
index 96aba49e..dd01515f 100644
--- a/crpptbx/TCV/gdat_tcv.m
+++ b/crpptbx/TCV/gdat_tcv.m
@@ -1108,7 +1108,7 @@ elseif strcmp(mapping_for_tcv.method,'switchcase')
    case {'powers'}
     % note: same time array for all main, ec, ohm, nbi, ...
     % At this stage fill just ech, later add nbi
-    sources_avail = {,'ohm','ec','nbi'}; % note should allow ech, nbh, ohmic in parameter sources
+    sources_avail = {'ohm','ec','nbi'}; % note should allow ech, nbh, ohmic in parameter sources
     % create empty structures for all, so in return one always have same substructres
     for i=1:length(sources_avail)
       gdat_data.(sources_avail{i}).data = [];
@@ -1124,13 +1124,14 @@ elseif strcmp(mapping_for_tcv.method,'switchcase')
       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 = {lower(gdat_data.gdat_params.source)};
+          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
-- 
GitLab