From 600655b936366f66b8e4843de9ae342781a3cf76 Mon Sep 17 00:00:00 2001
From: Olivier Sauter <olivier.sauter@epfl.ch>
Date: Thu, 27 Feb 2014 22:05:34 +0000
Subject: [PATCH] update to new sf2sig and add SXB/J series with sxb keyword
git-svn-id: https://spcsvn.epfl.ch/repos/TCV/gdat/trunk@4237 d63d8f72-b253-0410-a779-e742ad2e26cf
---
crpptbx/AUG/loadAUGdata.m | 152 ++++++++++++++++++++++++++++++++++++--
crpptbx/AUG/rdaAUG_eff.m | 87 ++++++++++++++++------
2 files changed, 211 insertions(+), 28 deletions(-)
diff --git a/crpptbx/AUG/loadAUGdata.m b/crpptbx/AUG/loadAUGdata.m
index ca23ae39..8b6d3c4f 100644
--- a/crpptbx/AUG/loadAUGdata.m
+++ b/crpptbx/AUG/loadAUGdata.m
@@ -6,6 +6,7 @@ function [trace,error,varargout]=loadAUGdata(shot,data_type,varargin)
% 'rmag' = radial position of the center of the plasma
% 'sxr' = soft x-ray emission
% 'sxR' = soft x-ray emission with varargout{1} option (requires varargin{5}!)
+% 'SXB' = soft x-ray emission from (by default camera J) SXB/J_xx camera (sxb, sxB, etc all work)
%
% gdat(15133,'MAG/Ipi',1,'AUG')
%
@@ -14,7 +15,7 @@ function [trace,error,varargout]=loadAUGdata(shot,data_type,varargin)
% data_type: type of the required data: 'diag_name/sig_name'
%
% examples:
-% data_type='SXR/B', 'TOT/beta_N'
+% data_type='SXR/B', 'TOT/beta_N', 'SXB/J_053', 'SXB'
% data_type='POT/ELMa-Han', 'MOD/OddNAmp', 'MOD/EvenNAmp', 'TOT/PNBI_TOT'
%
% Meaning of varargin depends on data_type:
@@ -82,6 +83,9 @@ if size(data_type_eff,1)==1
if ~isempty(strmatch(data_type_eff_noext,[{'SXR'}],'exact'))
data_type_eff_noext='sxr';
end
+ if ~isempty(strmatch(data_type_eff_noext,[{'SXB'} {'sxb'} {'Sxb'} {'sXb'} {'sxB'} {'SXb'}],'exact'))
+ data_type_eff_noext='sxb';
+ end
if ~isempty(strmatch(data_type_eff_noext,[{'ECE'}],'exact'))
data_type_eff_noext='ece';
end
@@ -115,6 +119,9 @@ if size(data_type_eff,1)==1
if ~isempty(strmatch(data_type_eff_noext,[{'Zcont'}],'exact'))
data_type_eff_noext='zcont';
end
+ if ~isempty(strmatch(data_type_eff_noext,[{'Ha'} {'ha'} {'Halpha'} {'halpha'}],'exact'))
+ data_type_eff_noext='Halpha';
+ end
else
i_ext=length(data_type_eff{2})+1;
name_ext='';
@@ -122,10 +129,10 @@ else
end
% all keywords and corresponding case to run below
-AUGkeywrdall=[{'Ip'} {'zmag'} {'rmag'} {'rcont'} {'zcont'} {'vol'} {'qrho'} {'q95'} {'kappa'} ...
+AUGkeywrdall=[{'Ip'} {'zmag'} {'rmag'} {'rcont'} {'zcont'} {'vol'} {'qrho'} {'qrho_cliste'} {'q95'} {'kappa'} ...
{'delta'} {'deltatop'} {'deltabot'} {'neint'} ...
{'ne'} {'te'} {'nerho'} {'terho'} ...
- {'sxr'} {'sxR'} {'ece'}];
+ {'sxr'} {'sxR'} {'sxb'} {'ece'} {'Halpha'}];
AUGsig.iip=strmatch('Ip',AUGkeywrdall,'exact');
AUGsig.izmag=strmatch('zmag',AUGkeywrdall,'exact');
AUGsig.irmag=strmatch('rmag',AUGkeywrdall,'exact');
@@ -133,6 +140,7 @@ AUGsig.ircont=strmatch('rcont',AUGkeywrdall,'exact');
AUGsig.izcont=strmatch('zcont',AUGkeywrdall,'exact');
AUGsig.ivol=strmatch('vol',AUGkeywrdall,'exact');
AUGsig.iqrho=strmatch('qrho',AUGkeywrdall,'exact');
+AUGsig.iqrho_cliste=strmatch('qrho_cliste',AUGkeywrdall,'exact');
AUGsig.iq95=strmatch('q95',AUGkeywrdall,'exact');
AUGsig.ikappa=strmatch('kappa',AUGkeywrdall,'exact');
AUGsig.idelta=strmatch('delta',AUGkeywrdall,'exact');
@@ -145,20 +153,24 @@ AUGsig.inerho=strmatch('nerho',AUGkeywrdall,'exact');
AUGsig.iterho=strmatch('terho',AUGkeywrdall,'exact');
AUGsig.isxr=strmatch('sxr',AUGkeywrdall,'exact');
AUGsig.isxR=strmatch('sxR',AUGkeywrdall,'exact');
+AUGsig.isxb=strmatch('sxb',AUGkeywrdall,'exact');
AUGsig.iece=strmatch('ece',AUGkeywrdall,'exact');
+AUGsig.iHalpha=strmatch('Halpha',AUGkeywrdall,'exact');
% For each keyword, specify which case to use. As most common is 'simplereaddata', fill in with this and change
% only indices needed. Usually use name of case same as keyword name
AUGkeywrdcase=cell(size(AUGkeywrdall));
AUGkeywrdcase(:)={'simplereaddata'};
-%AUGkeywrdcase(AUGsig.iqrho)=AUGkeywrdall(AUGsig.iqrho); % special as efit q on psi
+AUGkeywrdcase(AUGsig.iqrho)=AUGkeywrdall(AUGsig.iqrho); % special as efit q on psi
+AUGkeywrdcase(AUGsig.iqrho_cliste)=AUGkeywrdall(AUGsig.iqrho_cliste); % special as efit q on psi
%AUGkeywrdcase(AUGsig.idelta)=AUGkeywrdall(AUGsig.idelta); % special as average of triu and tril
-%AUGkeywrdcase(AUGsig.ine)=AUGkeywrdall(AUGsig.ine); % special as adds error bars
-%AUGkeywrdcase(AUGsig.ite)=AUGkeywrdall(AUGsig.ite); % idem
+AUGkeywrdcase(AUGsig.ine)=AUGkeywrdall(AUGsig.ine); % special as adds error bars
+AUGkeywrdcase(AUGsig.ite)=AUGkeywrdall(AUGsig.ite); % idem
%AUGkeywrdcase(AUGsig.inerho)=AUGkeywrdall(AUGsig.inerho); % idem
%AUGkeywrdcase(AUGsig.iterho)=AUGkeywrdall(AUGsig.iterho); % idem
AUGkeywrdcase(AUGsig.isxr)=AUGkeywrdall(AUGsig.isxr);
AUGkeywrdcase(AUGsig.isxR)=AUGkeywrdall(AUGsig.isxR);
+AUGkeywrdcase(AUGsig.isxb)=AUGkeywrdall(AUGsig.isxb);
%AUGkeywrdcase(AUGsig.iece)=AUGkeywrdall(AUGsig.iece);
% Information about which dimension has time, always return 2D data as (x,t) array
@@ -180,6 +192,8 @@ AUGsiglocation(:,AUGsig.ikappa)={''; ''};
AUGsiglocation(:,AUGsig.ideltatop)={''; ''};
AUGsiglocation(:,AUGsig.ideltabot)={''; ''};
AUGsiglocation(:,AUGsig.ineint)={'DCN'; 'H-1'};
+AUGsiglocation(:,AUGsig.iece)={'CEC'; 'Trad-A'};
+AUGsiglocation(:,AUGsig.iHalpha)={'POT'; 'ELMa-Han'};
% initialize order of substructures and allows just a "return" if data empty
trace.data=[];
@@ -344,6 +358,132 @@ switch AUGkeywrdcase{index}
trace.R=radius;
end
+ %&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&
+ case {'sxb'}
+ % LOAD MULTI CHANNEL DATA SXB/J_0xx (or other than J camera if specified in varargin{8})
+ % load AUG soft x-ray data
+ if nargin>=3 & ~isempty(varargin{1})
+ % chords to be loaded
+ starti=varargin{1}(1);
+ endi=varargin{1}(2);
+ else
+ starti=52;
+ endi=54;
+ end
+ if nargin>=4 & ~isempty(varargin{2})
+ % chords already loaded, 1=not loaded (to load), 0=loaded
+ status=varargin{2};
+ else
+ status=ones(endi-starti+1,1);
+ end
+ if nargin>=6 & ~isempty(varargin{4})
+ timerange=varargin{4};
+ else
+ timerange=[0 10];
+ end
+ if nargin>=7 & ~isempty(varargin{5})
+ nth=varargin{5};
+ else
+ nth=13;
+ end
+ if nargin>=8 & ~isempty(varargin{6})
+ tracename=varargin{6};
+ else
+ tracename='J';
+ end
+ trace.t=[];
+ trace.x=[];
+ ppftype='SXB';
+ iok=0;
+ for ichord=starti:endi
+ tracename_eff = [tracename '_' num2str(ichord,'%.3d')];
+ [a,e]=rdaAUG_eff(shot,ppftype,tracename_eff,timerange);
+ if isempty(a) || e~=0
+ trace_all = struct([]);
+ else
+ if iok==0
+ trace_all = a;
+ trace_all = rmfield(trace_all,[{'value'},{'data'}]);
+ iok = iok+1;
+ else
+ iok = iok+1;
+ end
+ trace_all.value(ichord,:) = a.value;
+ trace_all.data(ichord,:) = a.data;
+ end
+ end
+ trace_all.dim=[{[starti:endi]'} ; {trace.t}];
+ trace = trace_all;
+ trace.x=trace.dim{1};
+ trace.dimunits=[{'channels'} ; {'time [s]'}];
+ trace.units='W/m^2';
+ trace.name=[num2str(shot) '/' ppftype '/' tracename];
+ % keep only nth points
+ trace.t=trace.t(1:nth:end);
+ trace.data=trace.data(:,1:nth:end);
+ trace.dim{2}=trace.t;
+
+ %&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&&
+ case {'te', 'ne'}
+
+ if strcmp(AUGkeywrdcase{index},'te')
+ [a,e]=rdaAUG_eff(shot,'YPR','Te');
+ else
+ [a,e]=rdaAUG_eff(shot,'YPR','Ne');
+ end
+ trace = a;
+ trace.data = a.value';
+ trace.dim{1} = trace.x;
+ trace.dim{2} = trace.t;
+ trace.units = trace.unit;
+ trace.dimunits{1} = a.area.name ;
+ trace.dimunits{2} = a.time_aug.unit;
+
+ case {'qrho', 'qrho_cliste'}
+
+ if strcmp(AUGkeywrdcase{index},'qrho')
+ DIAG = 'FPP';
+ else
+ DIAG = 'EQI';
+ end
+ [a,e]=rdaAUG_eff(shot,DIAG,'Qpsi');
+ % Qpsi has inverted channel/time from CEC
+ a.value = a.value(:,end:-1:1)';
+ a.data = a.value;
+ % get x values
+ [psi,e]=rdaAUG_eff(shot,DIAG,'PFL');
+ psi.value=psi.value(:,end:-1:1)';
+ psi_axis= (ones(size(psi.value,1),1) * psi.value(1,:));
+ a.x = sqrt(abs((psi.value-psi_axis) ./(zeros(size(psi.value))-psi_axis) )); % psi_edge=0 assumed
+ a.psi = psi.value;
+ % get time values
+ [psi_time,e]=rdaAUG_eff(shot,DIAG,'time');
+ if length(psi_time.value) > length(a.x(1,:))
+ a.t = psi_time.value(1:length(a.x(1,:))); % problem with times??
+ elseif length(psi_time.value) < length(a.x(1,:))
+ a.t = psi_time.value; % problem with times??
+ a.t(end+1:length(a.x(1,:))) = psi_time.value(end) + [1:length(a.x(1,:))-length(psi_time.value)].*diff(psi_time.value(end-1:end));
+ else
+ a.t = psi_time.value;
+ end
+ % get rhotor values
+ [phi,e]=rdaAUG_eff(shot,DIAG,'TFLx');
+ phi.value=phi.value(:,end:-1:1)';
+ a.rhotor = sqrt(abs(phi.value ./ (ones(size(phi.value,1),1) * phi.value(end,:))));
+ a.torflux = phi.value;
+ % get rhovol values
+ [Vol,e]=rdaAUG_eff(shot,DIAG,'Vol');
+ Vol.value=Vol.value(:,end-1:-2:1)'; % 2nd index are dV/dpsi
+ a.rhovol = sqrt(abs(Vol.value ./ (ones(size(Vol.value,1),1) * Vol.value(end,:))));
+ a.volume = Vol.value;
+ %
+ trace = a;
+ trace.dim{1} = trace.x;
+ trace.dim{2} = trace.t;
+ trace.units = trace.unit;
+ trace.dimunits{1} = '';
+ trace.dimunits{2} = 's';
+
otherwise
disp('case not yet defined')
diff --git a/crpptbx/AUG/rdaAUG_eff.m b/crpptbx/AUG/rdaAUG_eff.m
index 5086a523..d69299c0 100644
--- a/crpptbx/AUG/rdaAUG_eff.m
+++ b/crpptbx/AUG/rdaAUG_eff.m
@@ -1,6 +1,6 @@
function [adata,error]=rdaAUG_eff(shot,diagname,sigtype,varargin);
%
-% gets data using SF2ML or mdsplus (when mdsplus will be available)
+% gets data using sf2sig or mdsplus (when mdsplus will be available)
% 1D arrays: assumes dimension is time
% 2D arrays: assumes data vs (x,time)
% 3D arrays: assumes data vs (x,time,hsig) (for mdsplus)
@@ -11,10 +11,10 @@ function [adata,error]=rdaAUG_eff(shot,diagname,sigtype,varargin);
% [data,error]=rdaAUG_eff(15133,'MAG','Ipi');
% [data,error]=rdaAUG_eff(15133,'MAG','Ipi',[1 5]);
%
-% set global variable: usemdsplus to decide if SF2ML or mdsplus is used:
+% set global variable: usemdsplus to decide if sf2sig or mdsplus is used:
% >> global usemdsplus
% >> usemdsplus=1 % means use mds to get data
-% >> usemdsplus=0 % means use SF2ML (default if not defined)
+% >> usemdsplus=0 % means use sf2sig (default if not defined)
% if ~exist('usemdsplus'); usemdsplus=0; end
%
@@ -86,30 +86,73 @@ if usemdsplus
end
else
- % use SF2ML
+ % use sf2sig
if isempty(time_int)
- adata=SF2ML(sigtype,diagname,shot);
+ try
+ [adata,adata_time, adata_area]=sf2sig(diagname,shot,sigtype);
+ catch
+ adata = struct([]);
+ adata_time = struct([]);
+ adata_area = [];
+ end
else
- adata=SF2ML(sigtype,diagname,shot,time_int(1),time_int(end));
+ try
+ [adata,adata_time, adata_area]=sf2sig(diagname,shot,sigtype,[time_int(1);time_int(end)]);
+ catch
+ adata = struct([]);
+ adata_time = struct([]);
+ adata_area = [];
+ end
end
- if size(adata.y,1)==1
- adata.x=[];
- adata.y=adata.y';
+
+ if isempty(adata)
+ return
+ end
+
+ % special checks
+ if strcmp(upper(diagname),'SXB')
+ % time missing one point
+ if length(adata.value) == length(adata_time.value)+1
+ adata_time.value=linspace(adata_time.range(1),adata_time.range(2),length(adata.value));
+ adata_time.index(2) = length(adata.value);
+ end
+ end
+
+ adata.time_aug = adata_time;
+
+ adata.area = adata_area;
+
+ if (prod(size(adata.value))==length(adata.value))
+ % only time signal
+ adata.x = [];
+ adata.value=reshape(adata.value,1,length(adata.value));
+ if ~isempty(adata.time_aug)
+ adata.t=adata.time_aug.value;
+ else
+ adata.t=[1:size(adata.value,2)];
+ end
else
- adata.x=[1:size(adata.y,1)];
+ adata.value = adata.value';
+ if ~isempty(adata.time_aug)
+ adata.x=[1:prod(size(adata.value))/length(adata_time.value)];
+ adata.t=adata.time_aug.value;
+ else
+ adata.x=[1:size(adata.value,1)];
+ adata.t=[1:size(adata.value,2)];
+ end
end
- adata.data=adata.y;
- adata.t=adata.t';
-% % transpose data as output in C format, reversed from Fortran and matlab standard
-% ss=size(a);
-% nbofdim=length(ss);
-% if ss(end)==1; nbofdim=nbofdim-1; end
-% nbofdim=max(nbofdim,1);
-% if nbofdim==1
-% data=a;
-% else
-% data=a';
-% end
+ adata.data=adata.value;
+ adata.units = adata.unit;
+ % % transpose data as output in C format, reversed from Fortran and matlab standard
+ % ss=size(a);
+ % nbofdim=length(ss);
+ % if ss(end)==1; nbofdim=nbofdim-1; end
+ % nbofdim=max(nbofdim,1);
+ % if nbofdim==1
+ % data=a;
+ % else
+ % data=a';
+ % end
end
error=0;
--
GitLab