Fix nhit calculation and debug output in TFCSLateralShapeParametrizationHitChain

Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/src/TFCSLateralShapeParametrizationHitChain.cxx

 /*
-  Copyright (C) 2002-2019 CERN for the benefit of the ATLAS collaboration
+  Copyright (C) 2002-2020 CERN for the benefit of the ATLAS collaboration
 */
 
 #include "ISF_FastCaloSimEvent/TFCSLateralShapeParametrizationHitChain.h"
 #include "ISF_FastCaloSimEvent/FastCaloSim_CaloCell_ID.h"
 
 #include "ISF_FastCaloSimEvent/TFCSSimulationState.h"
+#include "TMath.h"
 
 //=============================================
 //======= TFCSLateralShapeParametrizationHitChain =========
@@ -66,7 +67,7 @@ float TFCSLateralShapeParametrizationHitChain::get_E_hit(TFCSSimulationState& si
 {  
   const int nhits = get_number_of_hits(simulstate,truth,extrapol);
   const int sample = calosample();
-  if(nhits<=0 || sample<0) return -1.;
+  if(nhits<=0 || sample<0) return simulstate.E(sample);
   const float maxWeight = getMaxWeight();// maxWeight = -1 if  shapeWeight class is not in m_chain  
   
   if(maxWeight>0) return simulstate.E(sample)/(maxWeight*nhits); // maxWeight is used only if shapeWeight class is in m_chain
@@ -89,17 +90,21 @@ float TFCSLateralShapeParametrizationHitChain::get_sigma2_fluctuation(TFCSSimula
 
 FCSReturnCode TFCSLateralShapeParametrizationHitChain::simulate(TFCSSimulationState& simulstate,const TFCSTruthState* truth, const TFCSExtrapolationState* extrapol) const
 {
-  // Call get_number_of_hits() only once, as it could contain a random number
-  int nhit = get_number_of_hits(simulstate, truth, extrapol);
-  if (nhit <= 0) {
-    ATH_MSG_ERROR("TFCSLateralShapeParametrizationHitChain::simulate(): number of hits could not be calculated");
+  const float Elayer=simulstate.E(calosample());
+  const float Ehit=get_E_hit(simulstate,truth,extrapol);
+  if (Ehit <= 0) {
+    ATH_MSG_ERROR("TFCSLateralShapeParametrizationHitChain::simulate(): Ehit negative Ehit="<<Ehit);
     return FCSFatal;
   }
 
-  const float Elayer=simulstate.E(calosample());
-  const float Ehit=get_E_hit(simulstate,truth,extrapol);
+  // Call get_number_of_hits() only once inside get_E_hit(...), 
+  // as it could contain a random number
+  int nhit = TMath::Nint(Elayer/Ehit);
+  if(nhit<1) nhit=1;
+
   float sumEhit=0;
 
+  MSG::Level old_level=level();
   const bool debug = msgLvl(MSG::DEBUG);
   if (debug) {
     ATH_MSG_DEBUG("E("<<calosample()<<")="<<simulstate.E(calosample())<<" #hits~"<<nhit);
@@ -113,7 +118,7 @@ FCSReturnCode TFCSLateralShapeParametrizationHitChain::simulate(TFCSSimulationSt
     hit.E()=Ehit;
     for(TFCSLateralShapeParametrizationHitBase* hitsim : m_chain) {
       if (debug) {
-        if (ihit < 2) hitsim->setLevel(MSG::DEBUG);
+        if (ihit < 2) hitsim->setLevel(old_level);
         else hitsim->setLevel(MSG::INFO);
       }
 
@@ -127,7 +132,14 @@ FCSReturnCode TFCSLateralShapeParametrizationHitChain::simulate(TFCSSimulationSt
           //if(sumEhit+hit.E()>Elayer) hit.E()=Elayer-sumEhit;//sum of all hit energies needs to be Elayer: correct last hit accordingly
           break;
         } else {
-          if (status == FCSFatal) return FCSFatal;
+          if (status == FCSFatal) {
+            if (debug) {
+              for(TFCSLateralShapeParametrizationHitBase* reset : m_chain) {
+                reset->setLevel(old_level);
+              }
+            }  
+            return FCSFatal;
+          }  
         }    
 
         if (i == FCS_RETRY_COUNT) {
@@ -147,6 +159,11 @@ FCSReturnCode TFCSLateralShapeParametrizationHitChain::simulate(TFCSSimulationSt
     }  
   } while (sumEhit<Elayer);
 
+  if (debug) {
+    for(TFCSLateralShapeParametrizationHitBase* reset : m_chain) {
+      reset->setLevel(old_level);
+    }
+  }  
   return FCSSuccess;
 }