Reduce for loops and add comments

matlab/TCV_IMAS/tcv_get_ids_core_sources.m

@@ -133,97 +133,74 @@ cd_dens       = ec_data.cd_dens.data;
 cd_integrated = ec_data.cd_integrated.data;
 
 p_ec_injected = powers_gdat.ec.data;
-it = iround_os(ec_powers_tgrid,ec_data_tgrid);
-sparse_p_ec_injected = p_ec_injected(it,:); % injected ec power vals corresponding to ec_data_tgrid
+ij = iround_os(ec_powers_tgrid,ec_data_tgrid);
+sparse_p_ec_injected = p_ec_injected(ij,:); % injected ec power vals corresponding to ec_data_tgrid
 n_rho = size(p_dens, 1);
 
 % calculate normalised profiles on ec_data_tgrid grid
-norm_p_dens        = zeros(n_rho,n_active_launchers,nt_ec_data);
-norm_p_integrated  = zeros(n_rho,n_active_launchers,nt_ec_data);
-norm_cd_dens       = zeros(n_rho,n_active_launchers,nt_ec_data);
-norm_cd_integrated = zeros(n_rho,n_active_launchers,nt_ec_data);
-for it = 1:nt_ec_data
-  norm_p_dens_temp        = zeros(n_rho,n_active_launchers);
-  norm_p_integrated_temp  = zeros(n_rho,n_active_launchers);
-  norm_cd_dens_temp       = zeros(n_rho,n_active_launchers);
-  norm_cd_integrated_temp = zeros(n_rho,n_active_launchers);
-  for i_lau = n_active_launchers
-    norm_p_dens_temp(:,i_lau)        = p_dens(:,active_launchers(i_lau),it) ./ sparse_p_ec_injected(it, active_launchers(i_lau));
-    norm_p_integrated_temp(:,i_lau)  = p_integrated(:,active_launchers(i_lau)) ./ sparse_p_ec_injected(it, active_launchers(i_lau));
-    norm_cd_dens_temp(:,i_lau)       = cd_dens(:,active_launchers(i_lau),it) ./ sparse_p_ec_injected(it, i_lau);
-    norm_cd_integrated_temp(:,i_lau) = cd_integrated(:,active_launchers(i_lau),it) ./ sparse_p_ec_injected(it, active_launchers(i_lau));
-  end
-  norm_p_dens(:,:,it)        = norm_p_dens_temp;
-  norm_p_integrated(:,:,it)  = norm_p_integrated_temp;
-  norm_cd_dens(:,:,it)       = norm_cd_dens_temp;
-  norm_cd_integrated(:,:,it) = norm_cd_integrated_temp;
+norm_p_dens        = zeros(n_rho,nb_launchers+1,nt_ec_data);
+norm_p_integrated  = zeros(n_rho,nb_launchers+1,nt_ec_data);
+norm_cd_dens       = zeros(n_rho,nb_launchers+1,nt_ec_data);
+norm_cd_integrated = zeros(n_rho,nb_launchers+1,nt_ec_data);
+for it = 1:nt_ec_data 
+  norm_p_dens(:,:,it)        = p_dens(:,:,it)./sparse_p_ec_injected(it,:);
+  norm_p_integrated(:,:,it)  = p_integrated(:,:,it)./sparse_p_ec_injected(it,:);
+  norm_cd_dens(:,:,it)       = cd_dens(:,:,it)./sparse_p_ec_injected(it,:);
+  norm_cd_integrated(:,:,it) = cd_integrated(:,:,it)./sparse_p_ec_injected(it,:); 
 end
 
 % interpolate normalised p_dens profiles on ec_powers_tgrid 
-interp_norm_p_dens        = zeros(n_rho,n_active_launchers,nt_ec_powers);
-interp_norm_p_integrated  = zeros(n_rho,n_active_launchers,nt_ec_powers);
-interp_norm_cd_dens       = zeros(n_rho,n_active_launchers,nt_ec_powers);
-interp_norm_cd_integrated = zeros(n_rho,n_active_launchers,nt_ec_powers);
-for irho = 1:n_rho
-  for i_lau = 1:n_active_launchers
-    profile_p_dens        = squeeze(norm_p_dens(irho,i_lau,:));
-    profile_p_integrated  = squeeze(norm_p_integrated(irho,i_lau,:));
-    profile_cd_dens       = squeeze(norm_cd_dens(irho,i_lau,:));
-    profile_cd_integrated = squeeze(norm_cd_integrated(irho,i_lau,:));
-    
-    interp_profile_p_dens        = interp1(ec_data_tgrid, profile_p_dens, ec_powers_tgrid);
-    interp_profile_p_integrated  = interp1(ec_data_tgrid, profile_p_integrated, ec_powers_tgrid);
-    interp_profile_cd_dens       = interp1(ec_data_tgrid, profile_cd_dens, ec_powers_tgrid);
-    interp_profile_cd_integrated = interp1(ec_data_tgrid, profile_cd_integrated, ec_powers_tgrid);
-    
-    interp_norm_p_dens(irho,i_lau,:)        = interp_profile_p_dens;
-    interp_norm_p_integrated(irho,i_lau,:)  = interp_profile_p_integrated;
-    interp_norm_cd_dens(irho,i_lau,:)       = interp_profile_cd_dens;
-    interp_norm_cd_integrated(irho,i_lau,:) = interp_profile_cd_integrated;
+interp_norm_p_dens        = zeros(n_rho,nb_launchers+1,nt_ec_powers);
+interp_norm_p_integrated  = zeros(n_rho,nb_launchers+1,nt_ec_powers);
+interp_norm_cd_dens       = zeros(n_rho,nb_launchers+1,nt_ec_powers);
+interp_norm_cd_integrated = zeros(n_rho,nb_launchers+1,nt_ec_powers);
+for i_lau = active_launchers
+  for irho = 1:n_rho  
+    % get power and current density at each rho
+    trace_p_dens        = squeeze(norm_p_dens(irho,i_lau,:));
+    trace_p_integrated  = squeeze(norm_p_integrated(irho,i_lau,:));
+    trace_cd_dens       = squeeze(norm_cd_dens(irho,i_lau,:));
+    trace_cd_integrated = squeeze(norm_cd_integrated(irho,i_lau,:));
+    % interpolate on gdat powers tgrid
+    interp_norm_p_dens(irho,i_lau,:)        = interp1(ec_data_tgrid,trace_p_dens,ec_powers_tgrid);
+    interp_norm_p_integrated(irho,i_lau,:)  = interp1(ec_data_tgrid,trace_p_integrated,ec_powers_tgrid);
+    interp_norm_cd_dens(irho,i_lau,:)       = interp1(ec_data_tgrid,trace_cd_dens,ec_powers_tgrid);
+    interp_norm_cd_integrated(irho,i_lau,:) = interp1(ec_data_tgrid,trace_cd_integrated,ec_powers_tgrid);
   end
 end
 
 % normalised & interpolated profiles * p_ec_injected on ec_powers_tgrid
-interp_p_dens        = zeros(n_rho,n_active_launchers,nt_ec_powers);
-interp_p_integrated  = zeros(n_rho,n_active_launchers,nt_ec_powers);
-interp_cd_dens       = zeros(n_rho,n_active_launchers,nt_ec_powers);
-interp_cd_integrated = zeros(n_rho,n_active_launchers,nt_ec_powers);
-for it = 1:numel(ec_powers_tgrid)
-  unnormalised_p_dens        = zeros(n_rho,n_active_launchers);
-  unnormalised_p_integrated  = zeros(n_rho,n_active_launchers);
-  unnormalised_cd_dens       = zeros(n_rho,n_active_launchers);
-  unnormalised_cd_integrated = zeros(n_rho,n_active_launchers);
-  for i_lau = 1:n_active_launchers
-    unnormalised_p_dens(:,i_lau)        = interp_norm_p_dens(:,i_lau,it).* p_ec_injected(it, active_launchers(i_lau));
-    unnormalised_p_integrated(:,i_lau)  = interp_norm_p_integrated(:,i_lau,it).* p_ec_injected(it, active_launchers(i_lau));
-    unnormalised_cd_dens(:,i_lau)       = interp_norm_cd_dens(:,i_lau,it).* p_ec_injected(it, active_launchers(i_lau));
-    unnormalised_cd_integrated(:,i_lau) = interp_norm_cd_integrated(:,i_lau,it).* p_ec_injected(it, active_launchers(i_lau));
-  end
-  interp_p_dens(:,:,it)        = unnormalised_p_dens;
-  interp_p_integrated(:,:,it)  = unnormalised_p_integrated;
-  interp_cd_dens(:,:,it)       = unnormalised_cd_dens;
-  interp_cd_integrated(:,:,it) = unnormalised_cd_integrated;
+interp_p_dens        = zeros(n_rho,nb_launchers+1,nt_ec_powers);
+interp_p_integrated  = zeros(n_rho,nb_launchers+1,nt_ec_powers);
+interp_cd_dens       = zeros(n_rho,nb_launchers+1,nt_ec_powers);
+interp_cd_integrated = zeros(n_rho,nb_launchers+1,nt_ec_powers);
+for it = 1:nt_ec_powers
+  interp_p_dens(:,:,it)        = interp_norm_p_dens(:,:,it).*p_ec_injected(it,:);
+  interp_p_integrated(:,:,it)  = interp_norm_p_integrated(:,:,it).*p_ec_injected(it,:);
+  interp_cd_dens(:,:,it)       = interp_norm_cd_dens(:,:,it).*p_ec_injected(it,:);
+  interp_cd_integrated(:,:,it) = interp_norm_cd_integrated(:,:,it).*p_ec_injected(it,:);
 end
 
 % fill profiles_1d with interpolated profiles
-for ii = 1:nt_ec_powers
-  
-  for i_lau = active_launchers
-    ids_core_sources.source{last_index+i_lau}.profiles_1d{ii}.time = ec_powers_tgrid(ii); % profiles time 
-    
-    launcher_index = find(active_launchers == i_lau);
+for i_lau = active_launchers
+  for ii = 1:nt_ec_powers  
+    % time
+    ids_core_sources.source{last_index+i_lau}.profiles_1d{ii}.time = ec_powers_tgrid(ii); 
+    % rhotor grid
     ids_core_sources.source{last_index+i_lau}.profiles_1d{ii}.grid.rho_tor_norm = ...
-      ec_data.p_dens.rhotor_norm(1,:);
+      ec_data.p_dens.rhotor_norm(1,:); % constant in time
+    % power density
     ids_core_sources.source{last_index+i_lau}.profiles_1d{ii}.electrons.energy = ...
-      interp_p_dens(:,launcher_index,ii);
+      interp_p_dens(:,i_lau,ii);
+    % integrated power density
     ids_core_sources.source{last_index+i_lau}.profiles_1d{ii}.electrons.power_inside = ...
-      interp_p_integrated(:,launcher_index,ii);
+      interp_p_integrated(:,i_lau,ii);
     % current density (to adapt to <J.B>/B0)
     ids_core_sources.source{last_index+i_lau}.profiles_1d{ii}.j_parallel = ...
-      interp_cd_dens(:,launcher_index,ii);
+      interp_cd_dens(:,i_lau,ii);
     % integrated current density (to adapt to INTEGRAL(<J.B>/B0)*ds_phi)
     ids_core_sources.source{last_index+i_lau}.profiles_1d{ii}.current_parallel_inside = ...
-      interp_cd_integrated(:,launcher_index,ii); 
+      interp_cd_integrated(:,i_lau,ii); 
   end
 end