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

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
parent 22093e14
No related branches found
No related tags found
No related merge requests found
......@@ -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')
......
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;
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment