add shotfile specific output structure

git-svn-id: https://spcsvn.epfl.ch/repos/TCV/gdat/trunk@11977 d63d8f72-b253-0410-a779-e742ad2e26cf

crpptbx/AUG/gdat_aug.m

@@ -1717,8 +1717,10 @@ elseif strcmp(mapping_for_aug.method,'switchcase')
         end
       end
      otherwise
-      disp(['profiles ''fit'' from ' gdat_data.gdat_params.source ' not defined. At this stage only: ' sources_available{:} ...
-            ' are available. Ask O. Suater if need be']);
+      if ~strcmp(gdat_data.gdat_params.source,'VTA')
+        disp(['profiles ''fit'' from ' gdat_data.gdat_params.source ' not defined. At this stage only: ' sources_available{:} ...
+              ' are available. Ask O. Sauter if need be']);
+      end
     end
 
     %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

crpptbx/AUG/torbeam_prepare_inputs_and_run.m

@@ -88,7 +88,11 @@ else
     error('expects inbeam as field of torbeam_in (varargin{1}) from read_ech_py(shot,times_for_torbeam)')
   end
 end
-
+if isempty(times_for_torbeam_eff)
+  disp(char(10))
+  warning(['no times with beams on, return' char(10)])
+  return
+end
 % if torbeam_in provided as input (with its own time input sequence) and times_in provided as well: need to adapt times_in to match available inputs
 % (otherwise would take closest eqdsk)
 if doload==0
@@ -107,6 +111,10 @@ else
     error('expects eqd_all as field of torbeam_in (varargin{1}) from gdat(shot,''eqdsk'',''time'',torbeam_in.inbeam.t_some_beam_on,''write'',0);')
   end
 end
+if length(torbeam_in.eqd_all.eqdsk) < 1
+  error('expects non empty torbeam_in.eqd_all.eqdsk')
+end
+r0exp = torbeam_in.eqd_all.eqdsk{1}.r0;
 
 if doload; disp(['time for reading data: ' num2str(toc)]); end
   
@@ -250,17 +258,40 @@ for it_eff=1:length(times_for_torbeam_eff)
       torbeam_out.jcd(:,ibeam,it_eff) = torbeam_out.t2_new_LIB{ibeam,it_eff}(:,3) * 1e6;
       % structure for shot file, add rho=0 point, use a generic nb_points_shotfile
       x2 = [0.; torbeam_out.pdens_jcd_rho(:,ibeam,it_eff)];
-      pdens2 = [0.; torbeam_out.pdens(:,ibeam,it_eff)];
-      jcd2 = [0.; torbeam_out.jcd(:,ibeam,it_eff)];
+      pdens2 = [torbeam_out.pdens(1,ibeam,it_eff); torbeam_out.pdens(:,ibeam,it_eff)];
+      jcd2 = [torbeam_out.jcd(1,ibeam,it_eff); torbeam_out.jcd(:,ibeam,it_eff)];
+      sig2 = ones(size(x2)); sig2(1)=1e3;
+      sigvol=ones(size(torbeam_out.volumes{ibeam,it_eff}(1:end,1))); sigvol(2) = 1e3; %sigvol(2:3) = 1e3;
       torbeam_out.shotfile.volume(:,it_eff) = ...
-          interpos([0.;torbeam_out.volumes{ibeam,it_eff}(3:end,1)],[0.;torbeam_out.volumes{ibeam,it_eff}(3:end,2)],torbeam_out.shotfile.rhopol);
+          interpos([0.;torbeam_out.volumes{ibeam,it_eff}(2:end,1)],[0.;torbeam_out.volumes{ibeam,it_eff}(2:end,2)].^.5,torbeam_out.shotfile.rhopol,-0.001,[2 0],[0 0],sigvol);
+      torbeam_out.shotfile_vol_ok(it_eff) = 1;
+      if min(diff(torbeam_out.shotfile.volume(:,it_eff))) <=0
+        sigvol(2:3) = 1e3;
+        torbeam_out.shotfile.volume(:,it_eff) = ...
+            interpos([0.;torbeam_out.volumes{ibeam,it_eff}(2:end,1)],[0.;torbeam_out.volumes{ibeam,it_eff}(2:end,2)].^.5,torbeam_out.shotfile.rhopol,-0.001,[2 0],[0 0],sigvol);
+        torbeam_out.shotfile_vol_ok(it_eff) = 2;
+        if min(diff(torbeam_out.shotfile.volume(:,it_eff))) <=0
+          disp(char(10))
+          warning(['non-monotonic volume for it_eff = ' num2str(it_eff) ' (ibeam= ' num2str(ibeam) ')' char(10)])
+          torbeam_out.shotfile_vol_ok(it_eff) = 0;
+        end
+      end
       torbeam_out.shotfile.volume(1,it_eff) = 0.;
-      torbeam_out.shotfile.pdens(:,ibeam,it_eff) = interpos(x2,pdens2,torbeam_out.shotfile.rhopol);
-      torbeam_out.shotfile.jcd(:,ibeam,it_eff) = interpos(x2,jcd2,torbeam_out.shotfile.rhopol);
-      torbeam_out.shotfile.rz_r_ray(1:size(torbeam_out.t1_LIB{ibeam,it_eff},1),1:size(torbeam_out.t1_LIB{ibeam,it_eff},2)/2,ibeam,it_eff) = torbeam_out.t1_LIB{ibeam,it_eff}(:,1:2:end-1);
-      torbeam_out.shotfile.rz_z_ray(1:size(torbeam_out.t1_LIB{ibeam,it_eff},1),1:size(torbeam_out.t1_LIB{ibeam,it_eff},2)/2,ibeam,it_eff) = torbeam_out.t1_LIB{ibeam,it_eff}(:,2:2:end);
-      torbeam_out.shotfile.ry_r_ray(1:size(torbeam_out.t1_LIB{ibeam,it_eff},1),1:size(torbeam_out.t1_LIB{ibeam,it_eff},2)/2,ibeam,it_eff) = torbeam_out.t1tor_LIB{ibeam,it_eff}(:,1:2:end-1);
-      torbeam_out.shotfile.ry_y_ray(1:size(torbeam_out.t1_LIB{ibeam,it_eff},1),1:size(torbeam_out.t1_LIB{ibeam,it_eff},2)/2,ibeam,it_eff) = torbeam_out.t1tor_LIB{ibeam,it_eff}(:,2:2:end);
+      torbeam_out.shotfile.volume(:,it_eff) = torbeam_out.shotfile.volume(:,it_eff).^2;
+      torbeam_out.shotfile.pdens(:,ibeam,it_eff) = interpos(x2,pdens2,torbeam_out.shotfile.rhopol,-0.01,[1 0],[0 0],sig2);
+      torbeam_out.shotfile.jcd(:,ibeam,it_eff) = interpos(x2,jcd2,torbeam_out.shotfile.rhopol,-0.01,[1 0],[0 0],sig2);
+      if torbeam_out.shotfile_vol_ok(it_eff) > 0
+        [~,~,~,torbeam_out.shotfile.pabs_int(:,ibeam,it_eff)] = interpos(torbeam_out.shotfile.volume(:,it_eff),torbeam_out.shotfile.pdens(:,ibeam,it_eff));
+        [~,~,~,torbeam_out.shotfile.icd_int(:,ibeam,it_eff)] = interpos(torbeam_out.shotfile.volume(:,it_eff),torbeam_out.shotfile.jcd(:,ibeam,it_eff));
+        torbeam_out.shotfile.icd_int(:,ibeam,it_eff) = torbeam_out.shotfile.icd_int(:,ibeam,it_eff) ./2 ./ pi ./ r0exp / 1.01; % 1.01 on average while waiting to insert proper metric
+      else
+        torbeam_out.shotfile.pabs_int(:,ibeam,it_eff) = nan(size(torbeam_out.shotfile.volume(:,it_eff)));
+        torbeam_out.shotfile.icd_int(:,ibeam,it_eff) = nan(size(torbeam_out.shotfile.volume(:,it_eff)));
+      end
+      torbeam_out.shotfile.rz_r_ray(1:size(torbeam_out.t1_LIB{ibeam,it_eff},1),1:size(torbeam_out.t1_LIB{ibeam,it_eff},2)/2,ibeam,it_eff) = torbeam_out.t1_LIB{ibeam,it_eff}(:,1:2:end-1)/100.;
+      torbeam_out.shotfile.rz_z_ray(1:size(torbeam_out.t1_LIB{ibeam,it_eff},1),1:size(torbeam_out.t1_LIB{ibeam,it_eff},2)/2,ibeam,it_eff) = torbeam_out.t1_LIB{ibeam,it_eff}(:,2:2:end)/100.;
+      torbeam_out.shotfile.ry_r_ray(1:size(torbeam_out.t1_LIB{ibeam,it_eff},1),1:size(torbeam_out.t1_LIB{ibeam,it_eff},2)/2,ibeam,it_eff) = torbeam_out.t1tor_LIB{ibeam,it_eff}(:,1:2:end-1)/100.;
+      torbeam_out.shotfile.ry_y_ray(1:size(torbeam_out.t1_LIB{ibeam,it_eff},1),1:size(torbeam_out.t1_LIB{ibeam,it_eff},2)/2,ibeam,it_eff) = torbeam_out.t1tor_LIB{ibeam,it_eff}(:,2:2:end)/100.;
       torbeam_out.shotfile.index_abs_ray(1,ibeam,it_eff) = nnn(2) + 1;
       torbeam_out.shotfile.index_abs_ray(2,ibeam,it_eff) = nnn(2) + 1 + nnn(1);
       torbeam_out.shotfile.torbeam_was_run = torbeam_out.run_torbeam;