From 09d7bcaaff12afad8a01ee96501b690fc7eba873 Mon Sep 17 00:00:00 2001
From: Olivier Sauter <olivier.sauter@epfl.ch>
Date: Wed, 17 Feb 2016 10:12:30 +0000
Subject: [PATCH] add source in cxrs for CMZ, correct description when
expression used
git-svn-id: https://spcsvn.epfl.ch/repos/TCV/gdat/trunk@5478 d63d8f72-b253-0410-a779-e742ad2e26cf
---
crpptbx/AUG/aug_help_parameters.m | 2 +-
crpptbx/AUG/aug_requests_mapping.m | 5 +--
crpptbx/AUG/gdat_aug.m | 38 ++++++++++++++++---
crpptbx/addpsi_fromRZ.m | 60 ++++++++++++++++--------------
4 files changed, 67 insertions(+), 38 deletions(-)
diff --git a/crpptbx/AUG/aug_help_parameters.m b/crpptbx/AUG/aug_help_parameters.m
index f7da4dfa..fef32e41 100644
--- a/crpptbx/AUG/aug_help_parameters.m
+++ b/crpptbx/AUG/aug_help_parameters.m
@@ -37,7 +37,7 @@ help_struct_all.cocos = ['cocos value desired in output, uses eqdsk_cocos_transf
% $$$ help_struct_all.edge = '0 (default), 1 to get edge Thomson values';
% $$$ help_struct_all.fit = '0, no fit profiles, 1 (default) if fit profiles desired as well, relevant for _rho profiles. See also fit_type';
% $$$ help_struct_all.fit_type = 'provenance of fitted profiles ''conf'' (default) from conf nodes, ''avg'' or ''local'' for resp. proffit: nodes';
-help_struct_all.source = 'sxr: ''G'' (default, with ssx), camera name ''J'', ''G'', ...[[F-M], case insensitive';
+help_struct_all.source = 'sxr: ''G'' (default, with ssx), camera name ''J'', ''G'', ...[[F-M], case insensitive;; cxrs: ''CEZ'' (default), ''CMZ''';
help_struct_all.camera = ['[] (default, all), [i1 i2 ...] chord nbs ([1 3 5] if only chords 1, 3 and 5 are desired), ''central'' uses J_049'];
help_struct_all.freq = '''slow'', default (means ssx, 500kHz), lower sampling; ''fast'' full samples 2MHz; integer value nn for downsampling every nn''th points';
%help_struct_all. = '';
diff --git a/crpptbx/AUG/aug_requests_mapping.m b/crpptbx/AUG/aug_requests_mapping.m
index 958f1864..da601b3e 100644
--- a/crpptbx/AUG/aug_requests_mapping.m
+++ b/crpptbx/AUG/aug_requests_mapping.m
@@ -56,8 +56,6 @@ switch lower(data_request)
case 'beta'
mapping.timedim = 1;
mapping.label = '\beta';
- mapping.method = 'signal';
- mapping.expression = [{'TOT'},{'beta'}];
mapping.method = 'expression';
mapping.expression = ['params_eff = gdat_data.gdat_params;params_eff.data_request=''betan'';' ...
'gdat_tmp=gdat_aug(shot,params_eff);' ...
@@ -67,8 +65,7 @@ switch lower(data_request)
'tmp_data_ip=interp1(gdat_tmp2.t,gdat_tmp2.data,gdat_tmp.t,[],NaN);' ...
'tmp_data_b0=interp1(gdat_tmp3.t,gdat_tmp3.data,gdat_tmp.t,[],NaN);' ...
'tmp_data_a=interp1(gdat_tmp4.t,gdat_tmp4.data,gdat_tmp.t,[],NaN);' ...
- 'gdat_tmp.data = 0.01.*abs(gdat_tmp.data.*tmp_data_ip./1e6./tmp_data_a./tmp_data_b0);' ...
- 'gdat_tmp.label=''' mapping.label ''';gdat_tmp.gdat_request=''' data_request ''';'];
+ 'gdat_tmp.data = 0.01.*abs(gdat_tmp.data.*tmp_data_ip./1e6./tmp_data_a./tmp_data_b0);'];
case 'betan'
mapping.timedim = 1;
mapping.label = '\beta_N';
diff --git a/crpptbx/AUG/gdat_aug.m b/crpptbx/AUG/gdat_aug.m
index cb5469d6..816b5244 100644
--- a/crpptbx/AUG/gdat_aug.m
+++ b/crpptbx/AUG/gdat_aug.m
@@ -341,7 +341,12 @@ elseif strcmp(mapping_for_aug.method,'expression')
% assume expression contains an expression to evaluate and which creates a local structure into variable gdat_tmp
% we copy the structure, to make sure default nodes are defined and to avoid if return is an closed object like tdi
% eval_expr = [mapping_for_aug.expression ';'];
+ gdata_data_orig = gdat_data;
+ fields_to_copy = {'gdat_params', 'gdat_request', 'label', 'data_fullpath', 'request_description', 'mapping_for'};
eval([mapping_for_aug.expression ';']);
+ for i=1:length(fields_to_copy)
+ gdat_tmp.(fields_to_copy{i}) = gdata_data_orig.(fields_to_copy{i});
+ end
if isempty(gdat_tmp) || (~isstruct(gdat_tmp) & ~isobject(gdat_tmp))
if (gdat_params.nverbose>=1); warning(['expression does not create a gdat_tmp structure: ' mapping_for_aug.expression]); end
error_status=801;
@@ -388,10 +393,22 @@ 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';
+ if ~isfield(gdat_data.gdat_params,'source') || isempty(gdat_data.gdat_params.source) || strcmp(upper(gdat_data.gdat_params.source),'CEZ')
+ gdat_data.gdat_params.source = 'CEZ';
+ r_node = 'R_time';
+ z_node = 'z_time';
+ elseif strcmp(upper(gdat_data.gdat_params.source),'CMZ')
+ gdat_data.gdat_params.source = 'CMZ';
+ r_node = 'R';
+ z_node = 'z';
+ else
+ warning(['source = ' gdat_data.gdat_params.source ' not expected with data_request= ' data_request_eff]);
+ return
+ end
+ diag_name = gdat_data.gdat_params.source;
% R, Z positions of measurements
try
- [r_time]=sf2ab(diag_name,shot,'R_time','-exp',exp_name_eff);
+ [r_time]=sf2ab(diag_name,shot,r_node,'-exp',exp_name_eff);
catch ME_R_time
% assume no shotfile
getReport(ME_R_time,'basic');
@@ -400,7 +417,7 @@ elseif strcmp(mapping_for_aug.method,'switchcase')
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);
+ [z_time]=sf2ab(diag_name,shot,z_node,'-exp',exp_name_eff);
gdat_data.z = z_time.value{1};
inotok=find(gdat_data.z<=0);
gdat_data.z(inotok) = NaN;
@@ -466,10 +483,19 @@ elseif strcmp(mapping_for_aug.method,'switchcase')
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);
+ if strcmp(upper(gdat_data.gdat_params.source),'CEZ')
+ rout=gdat_data.r(iok,it_cxrs_inequil);
+ zout=gdat_data.z(iok,it_cxrs_inequil);
+ else
+ rout=gdat_data.r(iok);
+ zout=gdat_data.z(iok);
+ end
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));
+ if strcmp(upper(gdat_data.gdat_params.source),'CEZ')
+ psi_out(iok,it_cxrs_inequil) = reshape(psi_at_routzout,length(iok),length(it_cxrs_inequil));
+ else
+ psi_out(iok,it_cxrs_inequil) = repmat(psi_at_routzout,[1,length(it_cxrs_inequil)]);
+ end
rhopolnorm_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),rhopolnorm_out(iok,it_cxrs_inequil(it_cx)),-3,[2 2],[0 1]);
diff --git a/crpptbx/addpsi_fromRZ.m b/crpptbx/addpsi_fromRZ.m
index d3aaa784..cc83bf8d 100644
--- a/crpptbx/addpsi_fromRZ.m
+++ b/crpptbx/addpsi_fromRZ.m
@@ -1,57 +1,63 @@
+function [signal_with_1Dgrids,varargout] = addpsi_fromRZ(main_signal,equil,varargin)
+%
+% [signal_with_1Dgrids,varargout] = addpsi_fromRZ(main_signal,equil,varargin);
+%
+% At this stage tries with cxrs which has cxrs.r(:,time) and cxrs.z(:,time) for the R,Z points
+%
% $$$ shot=30382;
% $$$ time=3.;
% $$$
-% $$$ cxrs=gdat(shot,'cxrs',0);
+% $$$ main_signal=gdat(shot,'main_signal',0);
% $$$ equil=gdat(shot,'equil',0);
if ~exist('equil')
disp('no equil')
return
end
-if ~exist('cxrs')
- disp('no cxrs')
+if ~exist('main_signal')
+ disp('no main_signal')
return
end
-inb_chord_cxrs=size(cxrs.r_time,1);
-inb_time_cxrs=size(cxrs.r_time,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);
+inb_chord_main_signal=size(main_signal.data,1);
+inb_time_main_signal=size(main_signal.data,2);
+psi_out = NaN*ones(inb_chord_main_signal,inb_time_main_signal);
+rhopsinorm_out = NaN*ones(inb_chord_main_signal,inb_time_main_signal);
+rhotornorm_out = NaN*ones(inb_chord_main_signal,inb_time_main_signal);
+rhovolnorm_out = NaN*ones(inb_chord_main_signal,inb_time_main_signal);
% 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(cxrs.r_time(:,1)>0);
+iok=find(main_signal.r(:,1)>0);
tic
for itequil=1:length(time_equil)-1
rr=equil.Rmesh(:,itequil);
zz=equil.Zmesh(:,itequil);
psirz_in = equil.psi2D(:,:,itequil);
- it_cxrs_inequil = find(cxrs.t>=time_equil(itequil) & cxrs.t<=time_equil(itequil+1));
- if ~isempty(it_cxrs_inequil)
- rout=cxrs.r_time(iok,it_cxrs_inequil);
- zout=cxrs.z_time(iok,it_cxrs_inequil);
+ it_main_signal_inequil = find(main_signal.t>=time_equil(itequil) & main_signal.t<=time_equil(itequil+1));
+ if ~isempty(it_main_signal_inequil)
+ rout=main_signal.r(iok,it_main_signal_inequil);
+ zout=main_signal.z(iok,it_main_signal_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]);
+ psi_out(iok,it_main_signal_inequil) = reshape(psi_at_routzout,length(iok),length(it_main_signal_inequil));
+ rhopsinorm_out(iok,it_main_signal_inequil) = sqrt((psi_out(iok,it_main_signal_inequil)-equil.psi_axis(itequil))./(equil.psi_lcfs(itequil)-equil.psi_axis(itequil)));
+ for it_cx=1:length(it_main_signal_inequil)
+ rhotornorm_out(iok,it_main_signal_inequil(it_cx)) = interpos(equil.rhopolnorm(:,itequil),equil.rhotornorm(:,itequil),rhopsinorm_out(iok,it_main_signal_inequil(it_cx)),-3,[2 2],[0 1]);
+ rhovolnorm_out(iok,it_main_signal_inequil(it_cx)) = interpos(equil.rhopolnorm(:,itequil),equil.rhovolnorm(:,itequil),rhopsinorm_out(iok,it_main_signal_inequil(it_cx)),-3,[2 2],[0 1]);
end
end
end
toc
-cxrs_all = cxrs;
-cxrs_all.psi_on_rztime = psi_out;
-cxrs_all.rhopsinorm_on_rztime = rhopsinorm_out;
-cxrs_all.rhotornorm_on_rztime = rhotornorm_out;
-cxrs_all.rhovolnorm_on_rztime = rhovolnorm_out;
+signal_with_1Dgrids = main_signal;
+signal_with_1Dgrids.psi_on_rztime = psi_out;
+signal_with_1Dgrids.rhopsinorm_on_rztime = rhopsinorm_out;
+signal_with_1Dgrids.rhotornorm_on_rztime = rhotornorm_out;
+signal_with_1Dgrids.rhovolnorm_on_rztime = rhovolnorm_out;
-% $$$ iok=find(cxrs.r_time(:,1)>0);
-% $$$ rout=cxrs.r_time(iok,1);
-% $$$ zout=cxrs.z_time(iok,1);
+% $$$ iok=find(main_signal.r(:,1)>0);
+% $$$ rout=main_signal.r(iok,1);
+% $$$ zout=main_signal.z(iok,1);
% $$$
% $$$ % [farray_out,dfarrdx_out,dfarrdy_out,d2farrdx2_out,d2farrdy2_out,d2farrdxdy_out]= ...
% $$$ [farray_out]= interpos2Dcartesian(rr',zz',farray_in,rout,zout);
--
GitLab