add source for neint_ts and test if .t exists (TS tangential can be missing

matlab/D3D/gdat_d3d.m

@@ -875,7 +875,10 @@ elseif strcmp(mapping_for_d3d.method,'switchcase')
       %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
     case {'neint_ts'}
       if ~isfield(gdat_data.gdat_params,'source') || isempty(gdat_data.gdat_params.source) || ~isstruct(gdat_data.gdat_params.source)
-        error('requires ''ne_rho'' gdat structure as input in ''source'' option');
+        warning('requires ''ne_rho'' gdat structure as input in ''source'' option. Loading it...');
+        params_eff = gdat_data.gdat_params;
+        params_eff.data_request = 'ne_rho';
+        [gdat_data.gdat_params.source,params_nerho,error_status]=gdat_d3d(shot,params_eff);
       end
       ne_rho = gdat_data.gdat_params.source;
       % Computes R0 and V1 line integrated signal based on ne fit (so inside LCFS), provide also R(psi_ts) and Z(psi_ts) on these lines
@@ -1071,24 +1074,26 @@ elseif strcmp(mapping_for_d3d.method,'switchcase')
           rhotornorm_out_source = NaN*ones(inb_chord_source,inb_time_source);
           rhovolnorm_out_source = NaN*ones(inb_chord_source,inb_time_source);
           % constructs intervals within which a given equil is used: [time_equil(i),time_equil(i+1)]
-          time_equil=[min(gdat_data.(extra_source{i}).t(1)-0.1,equil.t(1)-0.1); 0.5.*(equil.t(1:end-1)+equil.t(2:end)) ;max(equil.t(end)+0.1,gdat_data.(extra_source{i}).t(end)+0.1)];
-          clear psi_out_core rhopolnorm_out_core rhotornorm_out_core rhovolnorm_out_core
-          for itequil=1:numel(time_equil)-1
-            rr=equil.Rmesh(:,itequil);
-            zz=equil.Zmesh(:,itequil);
-            psirz_in = equil.psi2D(:,:,itequil);
-            it_core_inequil = find(gdat_data.(extra_source{i}).t>time_equil(itequil) & gdat_data.(extra_source{i}).t<=time_equil(itequil+1));
-            if ~isempty(it_core_inequil)
-              rout_core=gdat_data.(extra_source{i}).r;
-              zout_core=gdat_data.(extra_source{i}).z;
-              psi_at_routzout = interpos2Dcartesian(rr,zz,psirz_in,rout_core,zout_core);
-              psi_out_core(:,it_core_inequil) = repmat(psi_at_routzout,1,numel(it_core_inequil));
-              rhopolnorm_out_core(:,it_core_inequil) = sqrt((psi_out_core(:,it_core_inequil)-equil.psi_axis(itequil))./(equil.psi_lcfs(itequil)-equil.psi_axis(itequil)));
-              for it_cx=1:length(it_core_inequil)
-                rhotornorm_out_core(:,it_core_inequil(it_cx)) = ...
+          if ~isempty(gdat_data.(extra_source{i}).t)
+            time_equil=[min(gdat_data.(extra_source{i}).t(1)-0.1,equil.t(1)-0.1); 0.5.*(equil.t(1:end-1)+equil.t(2:end)) ;max(equil.t(end)+0.1,gdat_data.(extra_source{i}).t(end)+0.1)];
+            clear psi_out_core rhopolnorm_out_core rhotornorm_out_core rhovolnorm_out_core
+            for itequil=1:numel(time_equil)-1
+              rr=equil.Rmesh(:,itequil);
+              zz=equil.Zmesh(:,itequil);
+              psirz_in = equil.psi2D(:,:,itequil);
+              it_core_inequil = find(gdat_data.(extra_source{i}).t>time_equil(itequil) & gdat_data.(extra_source{i}).t<=time_equil(itequil+1));
+              if ~isempty(it_core_inequil)
+                rout_core=gdat_data.(extra_source{i}).r;
+                zout_core=gdat_data.(extra_source{i}).z;
+                psi_at_routzout = interpos2Dcartesian(rr,zz,psirz_in,rout_core,zout_core);
+                psi_out_core(:,it_core_inequil) = repmat(psi_at_routzout,1,numel(it_core_inequil));
+                rhopolnorm_out_core(:,it_core_inequil) = sqrt((psi_out_core(:,it_core_inequil)-equil.psi_axis(itequil))./(equil.psi_lcfs(itequil)-equil.psi_axis(itequil)));
+                for it_cx=1:length(it_core_inequil)
+                  rhotornorm_out_core(:,it_core_inequil(it_cx)) = ...
                     interpos(equil.rhopolnorm(:,itequil),equil.rhotornorm(:,itequil),rhopolnorm_out_core(:,it_core_inequil(it_cx)),-3,[2 2],[0 1]);
-                rhovolnorm_out_core(:,it_core_inequil(it_cx)) = ...
+                  rhovolnorm_out_core(:,it_core_inequil(it_cx)) = ...
                     interpos(equil.rhopolnorm(:,itequil),equil.rhovolnorm(:,itequil),rhopolnorm_out_core(:,it_core_inequil(it_cx)),-3,[2 2],[0 1]);
+                end
               end
             end
           end