load_externaldata_TCV_MultiTimes

 This function loads data from TCV database

0001 function external = load_externaldata_TCV_MultiTimes(workdir,opt_gui,style,remote)
0002     %
0003     % This function loads data from TCV database
0004     %
0005     diagnames = {'BT','IP','VLOOP','NE','TE','TI','ZEFF','ZAXIS','ECRH','EQUIL','HXR'};
0006     %
0007     shotdata = iget_datainput_TCV(diagnames,workdir,opt_gui,style);% setup user-defined inputs with GUI or command line
0008     %
0009     % put input in form for luke_get_trdata.m
0010     %
0011     datainput = struct;
0012     datainput.shot = shotdata.shotnum;
0013     if ~isempty(shotdata.times),
0014         datainput.t = shotdata.times;
0015     end
0016     datainput.dtres=[]; % always use input times vector
0017     if ~isempty(shotdata.dtmax),
0018         datainput.dtreslim = shotdata.dtmax;
0019     end
0020     if ~isempty(shotdata.outputdir),
0021         datainput.outputpath = shotdata.outputdir;
0022     end
0023     datainput.flagTPG = shotdata.flagTPG;
0024     datainput.flagFF = shotdata.flagFF;
0025     %
0026     if ~shotdata.flagHXR,
0027         diagnames(strcmp(diagnames,'HXR')) = [];
0028     end
0029     %
0030     if ~shotdata.flagZeff,
0031         diagnames(strcmp(diagnames,'ZEFF')) = [];
0032     end
0033     %
0034     ndiag = length(diagnames);
0035     %
0036     for idiag = 1:ndiag,
0037         %
0038         datainput.diag{idiag}.name = diagnames{idiag};
0039         %
0040         if strcmp(diagnames{idiag},'ECRH'),
0041             if ~isempty(shotdata.isel),
0042                 datainput.diag{idiag}.isel = shotdata.isel;
0043             else
0044                 datainput.diag{idiag}.isel = 1:6;%1:9;%
0045             end
0046         end
0047         %
0048     end
0049     %
0050     if ~isempty(remote),
0051         remote.funcopy = 0;%use file version on remote cluster
0052     end
0053     %
0054     [err,~,dataoutput] = matremote_jd('luke_get_trdata',{datainput},remote);
0055     if ~isempty(err),
0056         disperr_jd(err);
0057         external = '';
0058         return
0059     else
0060         disp(['Data successfully loaded for shot ',num2str(shotdata.shotnum)])
0061     end
0062     %
0063     % fill external structure fields
0064     %
0065     external = struct;
0066     external.basestr = 'TCV';
0067     external.shotnum = num2str(shotdata.shotnum);
0068     %
0069     % get relevant times  - reference is te
0070     %
0071     shotimes = dataoutput.te.time;
0072     %
0073     shotimes = intersect(shotimes,dataoutput.ne.time);
0074     if ~all(isnan(dataoutput.ti.data(:))),
0075         shotimes = intersect(shotimes,dataoutput.ti.time);
0076     else
0077         disp('No ti data');        
0078     end
0079     if shotdata.flagZeff,
0080         shotimes = intersect(shotimes,dataoutput.zeff.time);
0081     end
0082     shotimes = intersect(shotimes,dataoutput.zaxis.time);
0083     shotimes = intersect(shotimes,dataoutput.ip.time);
0084     shotimes = intersect(shotimes,dataoutput.vloop.time);
0085     %
0086     nlaun = length(dataoutput.ecrh.isel);
0087     %
0088     for iisel = 1:nlaun,
0089         shotimes = intersect(shotimes,dataoutput.ecrh.ang{iisel}.time);
0090         shotimes = intersect(shotimes,dataoutput.ecrh.pow{iisel}.time);
0091         if isfield(dataoutput.ecrh,'xfrac'),
0092             shotimes = intersect(shotimes,dataoutput.ecrh.xfrac{iisel}.time);
0093         end
0094         if isfield(dataoutput.ecrh,'xfrac_FF'),
0095             shotimes = intersect(shotimes,dataoutput.ecrh.xfrac_FF{iisel}.time);
0096         end
0097         if isfield(dataoutput.ecrh,'xfrac_TPG'),
0098             shotimes = intersect(shotimes,dataoutput.ecrh.xfrac_TPG{iisel}.time);
0099         end        
0100         if isfield(dataoutput.ecrh,'stray'),
0101             shotimes = intersect(shotimes,dataoutput.ecrh.stray{iisel}.time);
0102         end        
0103     end
0104     %
0105     discnshotimes = setdiff(dataoutput.te.time,shotimes);
0106     if ~isempty(discnshotimes),
0107         disp(['WARNING : ',num2str(length(discnshotimes)),'/',num2str(length(dataoutput.te.time)),' times were discarded as they are missing from some diagnostics']);
0108     end
0109     %
0110     external.shotimes = shotimes;
0111     nt = length(external.shotimes);
0112     %
0113     Xfrac_FF = NaN(nt,nlaun);
0114     Xfrac_TPG = NaN(nt,nlaun);
0115     %
0116     for it = 1:nt,
0117         shotime = external.shotimes(it);
0118         %
0119         external.equil(it).basestr =  external.basestr;
0120         external.equil(it).shotnum =  external.shotnum;
0121         external.equil(it).shotime =  shotime;
0122         %
0123         if isfield(dataoutput.equil,'eqdsk'),
0124             external.equil(it).magnetic.eqdsk = dataoutput.equil.eqdsk{it};% magnetic equilibrium
0125         elseif isfield(dataoutput.equil,'data') && isempty(remote),% eqdsk files created locally, can be read
0126             eqdsk_file = [dataoutput.equil.data,filesep,'eqdsk.',external.shotnum,'t',num2str(shotime,'%6.4f')];
0127             if exist(eqdsk_file,'file'),
0128                 external.equil(it).magnetic.eqdsk = fileread(eqdsk_file);% magnetic equilibrium
0129             else
0130                 error(['File ',eqdsk_file,' could not be found. Aborted.']);
0131             end
0132         else
0133             error('eqdsk data could not be retrieved. Aborted.');
0134         end
0135         %
0136         external.equil(it).prof.te = dataprocess(dataoutput.te,'te',shotime);% te
0137         external.equil(it).prof.ne = dataprocess(dataoutput.ne,'ne',shotime);% ne
0138         %
0139         if all(isnan(dataoutput.ti.data(:))),
0140             external.equil(it).prof.ti = [];
0141         else
0142             external.equil(it).prof.ti = dataprocess(dataoutput.ti,'ti',shotime);% ti
0143         end
0144         %
0145         if shotdata.flagZeff,
0146             external.equil(it).prof.zeff = dataprocess(dataoutput.zeff,'zeff',shotime);% zeff
0147         else
0148             external.equil(it).prof.zeff = NaN;
0149         end
0150         %
0151         external.equil(it).prof.zaxis = dataprocess(dataoutput.zaxis,'zaxis',shotime);% zaxis
0152         external.equil(it).prof.ip = dataprocess(dataoutput.ip,'ip',shotime);% ip
0153         external.equil(it).prof.btsign = dataoutput.btsign.data;% btsign
0154         external.equil(it).prof.ipsign = dataoutput.ipsign.data;% ipsign
0155         %
0156         external.ohm(it).vloop = dataprocess(dataoutput.vloop,'vloop',shotime);% vloop
0157         %
0158         external.wave(it).laun = dataoutput.ecrh.isel;
0159         for iisel = 1:nlaun,
0160             %
0161             angles = selectfields_jd(dataoutput.ecrh.ang{iisel},{'freq','P1_R','P1_theta','P2_R','P2_d','P1_Z','P2_Z','P3_Z'});% mirrors
0162             %
0163             angles.theta = dataprocess(dataoutput.ecrh.ang{iisel},'theta',shotime,'theta');% theta
0164             angles.phi = dataprocess(dataoutput.ecrh.ang{iisel},'phi',shotime,'phi');% phi
0165             angles.theta_eff = dataprocess(dataoutput.ecrh.ang{iisel},'theta_eff',shotime,'theta_eff');% theta_eff
0166             angles.phi_eff = dataprocess(dataoutput.ecrh.ang{iisel},'phi_eff',shotime,'phi_eff');% phi_eff
0167             %
0168             angles.P0 = dataprocess(dataoutput.ecrh.pow{iisel},'pow',shotime);% P0
0169             if isfield(dataoutput.ecrh,'xfrac'),
0170                 angles.Xfrac = dataprocess(dataoutput.ecrh.xfrac{iisel},'xfrac',shotime);% Power fraction in X mode
0171             end
0172             if isfield(dataoutput.ecrh,'xfrac_FF'),
0173                 angles.Xfrac_FF = dataprocess(dataoutput.ecrh.xfrac_FF{iisel},'xfrac_FF',shotime);% Power fraction in X mode
0174                 Xfrac_FF(it,iisel) = angles.Xfrac_FF;
0175             end
0176             if isfield(dataoutput.ecrh,'xfrac_TPG'),
0177                 angles.Xfrac_TPG = dataprocess(dataoutput.ecrh.xfrac_TPG{iisel},'xfrac_TPG',shotime);% Power fraction in X mode
0178                 Xfrac_TPG(it,iisel) = angles.Xfrac_TPG;
0179             end        
0180             if isfield(dataoutput.ecrh,'stray'),
0181                 angles.stray = dataprocess(dataoutput.ecrh.stray{iisel},'stray',shotime);% Power fraction in X mode
0182             end        
0183             %
0184             external.wave(it).angles(iisel) = angles;
0185             %
0186         end
0187         %
0188         if shotdata.flagHXR,
0189             external.hxr = dataoutput.hxr;
0190         end
0191         %
0192     end
0193     %
0194     if ~isfield(angles,'Xfrac'),
0195         Xfrac_list{1} = 'FF code, mean Xfrac';
0196         Xfrac_list{2} = 'TPG code, mean Xfrac';
0197         for iisel = 1:nlaun, 
0198             mask = ~isnan(Xfrac_FF(:,iisel));
0199             if ~any(mask),
0200                 mXfrac_FF = NaN;
0201             else
0202                 mXfrac_FF = mean(Xfrac_FF(mask,iisel));
0203             end
0204             %
0205             mask = ~isnan(Xfrac_TPG(:,iisel));
0206             if ~any(mask),
0207                 mXfrac_TPG = NaN;
0208             else
0209                 mXfrac_TPG = mean(Xfrac_TPG(mask,iisel));
0210             end
0211             %
0212             Xfrac_list{1} = [Xfrac_list{1},' : L',num2str(iisel),' = ',num2str(mXfrac_FF,'%0.2f')];
0213             Xfrac_list{2} = [Xfrac_list{2},' : L',num2str(iisel),' = ',num2str(mXfrac_TPG,'%0.2f')];
0214         end
0215         %
0216         Xfrac_opt = iselect_jd(Xfrac_list,'Which X mode fraction do you want to select ?',opt_gui,style.buttonselectstyle,1);
0217         %
0218         for it = 1:nt,
0219             for iisel = 1:nlaun,
0220                 if Xfrac_opt == 1,
0221                     external.wave(it).angles(iisel).Xfrac = external.wave(it).angles(iisel).Xfrac_FF;
0222                 else
0223                     external.wave(it).angles(iisel).Xfrac = external.wave(it).angles(iisel).Xfrac_TPG;
0224                 end
0225             end
0226         end
0227         %
0228     end
0229     %
0230     external.id = ['TCV_',external.shotnum,'_MultiTimes'];
0231     %
0232 end
0233 %
0234 function dataout = dataprocess(datain,dataname,shotime,fname)
0235     %
0236     if nargin < 4,
0237         fname = 'data';
0238     end
0239     %
0240     itdata = datain.time == shotime;
0241     if any(itdata),
0242         if isfield(datain,'rho'),
0243             dataout = struct;
0244             dataout.rho = datain.rho;
0245             dataout.data = datain.(fname)(:,itdata);
0246         else
0247             dataout = datain.(fname)(itdata);
0248         end
0249     else
0250         error(['Time ',num2str(shotime),' not found in ',dataname,' data. aborted.']);
0251     end
0252     %
0253 end
0254 
0255 
0256 %
0257 %
0258 %
0259 %
0260 %
0261 %
0262 %
0263 %
0264 % if nargin < 7,
0265 %     select = struct;
0266 % elseif ~isstruct(select) && isempty(select),% no user prompt - default options
0267 %     %
0268 %     select.m_opt = 0;
0269 %     select.l_opt = 1;
0270 %     select.laun = 1:6;
0271 %     select.f_opt = 0;
0272 %     select.s_opt = 'n';
0273 %     %
0274 % end
0275 % %
0276 % external = '';
0277 % basestr = 'TCV';
0278 % %
0279 % if ~isfield(select,'m_opt'),
0280 %     select.m_opt = input_dke_yp('Do you want to: (0) load data from database, (1) load data from files',0,[0,1],'',[1,1]);
0281 % end
0282 % %
0283 % eqdsk_path = [workdir,'/'];
0284 % %
0285 % eqdsk_file = ['eqdsk.',shotnum,'t',shotime];
0286 % %
0287 % data_path = [workdir,'/'];
0288 % %
0289 % data_str = ['data_',basestr,'_',shotnum,'_',num2str(t1),'_',num2str(t2)];
0290 % %
0291 % if select.m_opt == 0,
0292 %     %
0293 %     if ~exist('make_eqdsk.m','file'),
0294 %         disp('-----> Access to the TCV database is not granted.');
0295 %         disp('-----> You need the file make_eqdsk.m in the matlab path.');
0296 %         return;
0297 %     end
0298 %     %
0299 %     if ~exist('collect_parameters_for_luke.m','file'),
0300 %         disp('-----> Access to the TCV database is not granted.');
0301 %         disp('-----> You need the file collect_parameters_for_luke.m in the matlab path.');
0302 %         return;
0303 %     end
0304 %     %
0305 %     if ~isfield(select,'l_opt'),
0306 %         select.l_opt = input('Option for liuqe version? (1,2,3) [1] : ');
0307 %     end
0308 %     %
0309 %     if isempty(select.l_opt),
0310 %         select.l_opt = 1;
0311 %     elseif ~isnumeric(select.l_opt) || ~any(select.l_opt == 1:3),
0312 %         return
0313 %     end
0314 %     %
0315 %     make_eqdsk(str2num(shotnum),str2num(shotime),workdir,'TCV',select.l_opt);%call LIUQUE
0316 %     %
0317 %     if ~isfield(select,'laun'),
0318 %         select.laun = input('Vector of EC launchers to account for (ex: [1,3]) : ');
0319 %     end
0320 %     %
0321 %     if isempty(select.laun) || ~isnumeric(select.laun) ,
0322 %         disp('-----> Invalid entry');
0323 %         return
0324 %     end
0325 %     %
0326 %     if ~isfield(select,'f_opt'),
0327 %         select.f_opt = input_dke_yp('Do you want to: (0) load Te from fitted Thomson only, (1) load Te with Xte renormalization',0,[0,1],'',[1,1]);
0328 %     end
0329 %     %
0330 %     launstr = '';
0331 %     for il = 1:length(select.laun),
0332 %         launstr = [launstr,'_',num2str(select.laun(il))];
0333 %     end
0334 %     %
0335 %     data_path = [workdir,'/'];
0336 %     data_file = [data_str,launstr,'.mat'];
0337 %     %
0338 %     collect_parameters_for_luke(str2num(shotnum),t1,t2,select.laun,[data_path,data_file],select.f_opt);%load shot parameters
0339 %     %
0340 % else
0341 %     %
0342 %     [eqdsk_file,eqdsk_path] = igetfile_jd(opt_gui,eqdsk_file,'Please select the EQDSK file',[eqdsk_path,eqdsk_file]);
0343 %     if eqdsk_file == 0,
0344 %         return
0345 %     end
0346 %     %
0347 %     data_file = [data_str,'*.mat'];
0348 %     %
0349 %     [data_file,data_path] = igetfile_jd(opt_gui,data_file,'Please select the data file',[data_path,data_file]);
0350 %     if data_file == 0,
0351 %         return
0352 %     end
0353 %     %
0354 %     i1 = length(data_str);
0355 %     i2 = strfind(data_file,'.mat');
0356 %     %
0357 %     launstr = data_file(i1+1:i2-1);
0358 %     %
0359 % end
0360 % %
0361 % % return structure equil
0362 % %
0363 % external.equil.shotnum = shotnum;
0364 % external.equil.shotime = shotime;
0365 % external.equil.magnetic.eqdsk_path = eqdsk_path;
0366 % external.equil.magnetic.eqdsk_file = eqdsk_file;
0367 % %
0368 % external.equil.prof = load([data_path,data_file],'shot','rho','ne','te','ti','zeff','zaxis','btsign','ipsign','teflag','Ip');
0369 % %
0370 % external.equil.prof.data_path = data_path;
0371 % external.equil.prof.data_file = data_file;
0372 % external.wave = load([data_path,data_file],'shot','laun','angles*');
0373 % external.ohm = load([data_path,data_file],'shot','vloop');
0374 % %
0375 % external.id = ['TCV_',shotnum,'_',shotime,'_',num2str(t1),'_',num2str(t2),'_L',launstr];
0376 % %
0377 % if select.m_opt == 0,
0378 %     %
0379 %     if ~isfield(select,'s_opt'),
0380 %         select.s_opt = input_dke_yp('Do you want to generate an equivalent make_external_TCV script? [y/n]','n',{'y','n'},'',[1,1]);
0381 %     end
0382 %     %
0383 %     if select.s_opt == 'y',
0384 %         %
0385 %         imake_makeexternal_TCV_jd(workdir,external,opt_gui);
0386 %         %
0387 %     end
0388 % end
0389 % %
0390 % if isnan_jd(external,0),
0391 %     disp('Warning : NaNs were found in the structure ''external'' : ')
0392 %     isnan_jd(external,1),
0393 % end
0394 % %
0395 
0396