diff --git a/macros/CCEScan/deplVphiN_normalized.C b/macros/CCEScan/deplVphiN_normalized.C
index 916b94d..ba3dfdc 100644
--- a/macros/CCEScan/deplVphiN_normalized.C
+++ b/macros/CCEScan/deplVphiN_normalized.C
@@ -317,6 +317,7 @@
                 double flux13_err = v_flux13_err(i_flux);
                 double lumi7     = v_lumi7(i);
                 double lumi8     = v_lumi8(i);
+                double lumi13    = v_lumi13(i);
                 double vdepl_prod = v_vdepl(j);
                 double vdepl_val = v_data(0);
                 double vdepl_err = v_data(1);
@@ -324,8 +325,8 @@
                 double vdepl_tot = vdepl_err;
                 double phiN_val = STTool::FLUKAConvFac;
                 double phiN_err = STTool::FLUKAConvFac;
-                phiN_val *= (flux7_val*lumi7+flux8_val*lumi8);
-                phiN_err *= (flux7_err*lumi7+flux8_err*lumi8);
+                phiN_val *= (flux7_val*lumi7 + flux8_val*lumi8 + flux13_val*lumi13);
+                phiN_err *= (flux7_err*lumi7 + flux8_err*lumi8 + flux13_err*lumi13);
                 // Define maximal and minimal Vdepl values to be drawn
                 double plotMin = 40.0;
                 double plotMax = 290.0;
@@ -344,16 +345,17 @@
                 if (vdepl_prod>vMax || vdepl_prod<=vMin) {
                     continue;
                 }
-                // Add systematic uncertainty to get mean value
-                mean_sys += vdepl_sys;
-                nMeas++;
                 // Remove excluded sectors
                 if (STTool::IsExcluded(sector)) {
                     Warning("deplVphiN","Sector excluded");
                     continue;
                 }
+                // Add systematic uncertainty to get mean value of all uncertainties
+                mean_sys += vdepl_sys;
+                nMeas++;
+                // Print data
                 Info("deplVphiN","Result for %s/%s sector %d, fill %d, radius %7.2f mm:",
-                         det.Data(),lay.Data(),sector,fill,rTT[k]>0.0?rTT[k]:2000.0);
+                      det.Data(),lay.Data(),sector,fill,rTT[k]>0.0?rTT[k]:2000.0);
                 Printf("     phi_N        = (%10.2f+/-%7.2f)*10^12 cm^{-2}",
                              phiN_val/1e12,phiN_err/1e12);
                 Printf("     V_depl       = (%6.2f+/-%6.2f(stat)+/-%6.2f(syst)) V",