From c68ffb9e7daf22b0f302a2c5cd71e2a810c08b39 Mon Sep 17 00:00:00 2001 From: Michael Duehrssen <michael.duehrssen@cern.ch> Date: Mon, 4 May 2020 17:29:56 +0200 Subject: [PATCH] Allow possibility to initialize the hit of the HitChain only once --- .../TFCSLateralShapeParametrizationHitChain.h | 11 +++- ...SLateralShapeParametrizationFluctChain.cxx | 53 ++++++++++++---- ...FCSLateralShapeParametrizationHitChain.cxx | 63 +++++++++++++------ 3 files changed, 93 insertions(+), 34 deletions(-) diff --git a/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/ISF_FastCaloSimEvent/TFCSLateralShapeParametrizationHitChain.h b/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/ISF_FastCaloSimEvent/TFCSLateralShapeParametrizationHitChain.h index 4cd1818b2011..c36c08acca82 100644 --- a/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/ISF_FastCaloSimEvent/TFCSLateralShapeParametrizationHitChain.h +++ b/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/ISF_FastCaloSimEvent/TFCSLateralShapeParametrizationHitChain.h @@ -24,8 +24,9 @@ public: const Chain_t& chain() const {return m_chain;}; Chain_t& chain() {return m_chain;}; void push_back( const Chain_t::value_type& value ) {m_chain.push_back(value);}; - //TODO: add generic functionality to determine the number of hits or center position only once - // and not for every iteration of the hit chain + + unsigned int get_nr_of_init() const {return m_ninit;}; + void set_nr_of_init(unsigned int ninit) {m_ninit=ninit;}; /// set which instance should determine the number of hits virtual void set_number_of_hits_simul(TFCSLateralShapeParametrizationHitBase* sim) {m_number_of_hits_simul=sim;}; @@ -48,11 +49,15 @@ public: void Print(Option_t *option = "") const override; protected: + void PropagateMSGLevel(MSG::Level level) const; + Chain_t m_chain; + int m_ninit=0; private: TFCSLateralShapeParametrizationHitBase* m_number_of_hits_simul; - ClassDefOverride(TFCSLateralShapeParametrizationHitChain,1) //TFCSLateralShapeParametrizationHitChain + + ClassDefOverride(TFCSLateralShapeParametrizationHitChain,2) //TFCSLateralShapeParametrizationHitChain }; #endif diff --git a/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/src/TFCSLateralShapeParametrizationFluctChain.cxx b/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/src/TFCSLateralShapeParametrizationFluctChain.cxx index 0cf78123233c..07ca7e48cae0 100644 --- a/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/src/TFCSLateralShapeParametrizationFluctChain.cxx +++ b/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/src/TFCSLateralShapeParametrizationFluctChain.cxx @@ -35,16 +35,43 @@ float TFCSLateralShapeParametrizationFluctChain::get_E_hit(TFCSSimulationState& FCSReturnCode TFCSLateralShapeParametrizationFluctChain::simulate(TFCSSimulationState& simulstate,const TFCSTruthState* truth, const TFCSExtrapolationState* extrapol) const { + MSG::Level old_level=level(); + const bool debug = msgLvl(MSG::DEBUG); + + //Execute the first get_nr_of_init() simulate calls only once. Used for example to initialize the center position + TFCSLateralShapeParametrizationHitBase::Hit hit; + hit.reset_center(); + auto hitloopstart=m_chain.begin(); + if(get_nr_of_init()>0) { + if (debug) { + PropagateMSGLevel(old_level); + ATH_MSG_DEBUG("E("<<calosample()<<")="<<simulstate.E(calosample())<<" before init"); + } + + hitloopstart+=get_nr_of_init(); + for(auto hititr=m_chain.begin(); hititr!=hitloopstart; ++hititr) { + TFCSLateralShapeParametrizationHitBase* hitsim=*hititr; + + FCSReturnCode status = hitsim->simulate_hit(hit, simulstate, truth, extrapol); + + if (status != FCSSuccess) { + ATH_MSG_ERROR("TFCSLateralShapeParametrizationFluctChain::simulate(): simulate_hit init call failed"); + return FCSFatal; + } + } + } + + //Initialize hit energy only now, as init loop above might change the layer energy const float Elayer=simulstate.E(calosample()); if (Elayer == 0) { ATH_MSG_VERBOSE("Elayer=0, nothing to do"); return FCSSuccess; } - // Call get_number_of_hits() only once, as it could contain a random number + // Call get_sigma2_fluctuation only once, as it could contain a random number float sigma2 = get_sigma2_fluctuation(simulstate, truth, extrapol); if (sigma2 >= s_max_sigma2_fluctuation) { - ATH_MSG_ERROR("TFCSLateralShapeParametrizationHitChain::simulate(): fluctuation of hits could not be calculated"); + ATH_MSG_ERROR("TFCSLateralShapeParametrizationFluctChain::simulate(): fluctuation of hits could not be calculated"); return FCSFatal; } @@ -58,16 +85,14 @@ FCSReturnCode TFCSLateralShapeParametrizationFluctChain::simulate(TFCSSimulation float error2_sumEhit=0; float error2=2*s_max_sigma2_fluctuation; - const bool debug = msgLvl(MSG::DEBUG); if (debug) { + PropagateMSGLevel(old_level); ATH_MSG_DEBUG("E("<<calosample()<<")="<<Elayer<<" sigma2="<<sigma2); } int ihit=0; int ifail=0; int itotalfail=0; - TFCSLateralShapeParametrizationHitBase::Hit hit; - hit.reset_center(); do { hit.reset(); //hit.E()=Eavghit; @@ -75,16 +100,20 @@ FCSReturnCode TFCSLateralShapeParametrizationFluctChain::simulate(TFCSSimulation hit.E()=CLHEP::RandGauss::shoot(simulstate.randomEngine(), Eavghit, m_RMS*Eavghit); } while (std::abs(hit.E())<absEavghit_tenth); bool failed=false; - for(TFCSLateralShapeParametrizationHitBase* hitsim : m_chain) { - if (debug) { - if (ihit < 2) hitsim->setLevel(MSG::DEBUG); - else hitsim->setLevel(MSG::INFO); - } + if(debug) if(ihit==2) { + //Switch debug output back to INFO to avoid huge logs + PropagateMSGLevel(MSG::INFO); + } + for(auto hititr=hitloopstart; hititr!=m_chain.end(); ++hititr) { + TFCSLateralShapeParametrizationHitBase* hitsim=*hititr; FCSReturnCode status = hitsim->simulate_hit(hit, simulstate, truth, extrapol); if (status == FCSSuccess) continue; - if (status == FCSFatal) return FCSFatal; + if (status == FCSFatal) { + if (debug) PropagateMSGLevel(old_level); + return FCSFatal; + } failed=true; ++ifail; ++itotalfail; @@ -103,12 +132,14 @@ FCSReturnCode TFCSLateralShapeParametrizationFluctChain::simulate(TFCSSimulation ATH_MSG_ERROR("TFCSLateralShapeParametrizationFluctChain::simulate(): simulate_hit call failed after " << FCS_RETRY_COUNT << "retries"); } if(ifail >= FCS_RETRY_COUNT*FCS_RETRY_COUNT) { + if (debug) PropagateMSGLevel(old_level); return FCSFatal; } } } while (error2>sigma2); if (debug) { + PropagateMSGLevel(old_level); ATH_MSG_DEBUG("E("<<calosample()<<")="<<Elayer<<" sumE="<<sumEhit<<"+-"<<TMath::Sqrt(error2_sumEhit)<<" ~ "<<TMath::Sqrt(error2_sumEhit)/sumEhit*100<<"% rel error^2="<<error2<<" sigma^2="<<sigma2<<" ~ "<<TMath::Sqrt(sigma2)*100<<"% hits="<<ihit<<" fail="<<itotalfail); } diff --git a/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/src/TFCSLateralShapeParametrizationHitChain.cxx b/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/src/TFCSLateralShapeParametrizationHitChain.cxx index 0543c4b891e3..d6bd27a1bbfe 100644 --- a/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/src/TFCSLateralShapeParametrizationHitChain.cxx +++ b/Simulation/ISF/ISF_FastCaloSim/ISF_FastCaloSimEvent/src/TFCSLateralShapeParametrizationHitChain.cxx @@ -126,8 +126,42 @@ float TFCSLateralShapeParametrizationHitChain::get_sigma2_fluctuation(TFCSSimula return s_max_sigma2_fluctuation; } +void TFCSLateralShapeParametrizationHitChain::PropagateMSGLevel(MSG::Level level) const +{ + for(TFCSLateralShapeParametrizationHitBase* reset : m_chain) { + reset->setLevel(level); + } +} + FCSReturnCode TFCSLateralShapeParametrizationHitChain::simulate(TFCSSimulationState& simulstate,const TFCSTruthState* truth, const TFCSExtrapolationState* extrapol) const { + MSG::Level old_level=level(); + const bool debug = msgLvl(MSG::DEBUG); + + //Execute the first get_nr_of_init() simulate calls only once. Used for example to initialize the center position + TFCSLateralShapeParametrizationHitBase::Hit hit; + hit.reset_center(); + auto hitloopstart=m_chain.begin(); + if(get_nr_of_init()>0) { + if (debug) { + PropagateMSGLevel(old_level); + ATH_MSG_DEBUG("E("<<calosample()<<")="<<simulstate.E(calosample())<<" before init"); + } + + hitloopstart+=get_nr_of_init(); + for(auto hititr=m_chain.begin(); hititr!=hitloopstart; ++hititr) { + TFCSLateralShapeParametrizationHitBase* hitsim=*hititr; + + FCSReturnCode status = hitsim->simulate_hit(hit, simulstate, truth, extrapol); + + if (status != FCSSuccess) { + ATH_MSG_ERROR("TFCSLateralShapeParametrizationHitChain::simulate(): simulate_hit init call failed"); + return FCSFatal; + } + } + } + + //Initialize hit energy only now, as init loop above might change the layer energy const float Elayer=simulstate.E(calosample()); if (Elayer == 0) { ATH_MSG_VERBOSE("Elayer=0, nothing to do"); @@ -146,23 +180,21 @@ FCSReturnCode TFCSLateralShapeParametrizationHitChain::simulate(TFCSSimulationSt float sumEhit=0; - MSG::Level old_level=level(); - const bool debug = msgLvl(MSG::DEBUG); if (debug) { + PropagateMSGLevel(old_level); ATH_MSG_DEBUG("E("<<calosample()<<")="<<simulstate.E(calosample())<<" #hits~"<<nhit); } int ihit=0; - TFCSLateralShapeParametrizationHitBase::Hit hit; - hit.reset_center(); do { hit.reset(); hit.E()=Ehit; - for(TFCSLateralShapeParametrizationHitBase* hitsim : m_chain) { - if (debug) { - if (ihit < 2) hitsim->setLevel(old_level); - else hitsim->setLevel(MSG::INFO); - } + if(debug) if(ihit==2) { + //Switch debug output back to INFO to avoid huge logs + PropagateMSGLevel(MSG::INFO); + } + for(auto hititr=hitloopstart; hititr!=m_chain.end(); ++hititr) { + TFCSLateralShapeParametrizationHitBase* hitsim=*hititr; for (int i = 0; i <= FCS_RETRY_COUNT; i++) { //TODO: potentially change logic in case of a retry to redo the whole hit chain from an empty hit instead of just redoing one step in the hit chain @@ -171,15 +203,10 @@ FCSReturnCode TFCSLateralShapeParametrizationHitChain::simulate(TFCSSimulationSt FCSReturnCode status = hitsim->simulate_hit(hit, simulstate, truth, extrapol); if (status == FCSSuccess) { - //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) { - if (debug) { - for(TFCSLateralShapeParametrizationHitBase* reset : m_chain) { - reset->setLevel(old_level); - } - } + if (debug) PropagateMSGLevel(old_level); return FCSFatal; } } @@ -201,11 +228,7 @@ FCSReturnCode TFCSLateralShapeParametrizationHitChain::simulate(TFCSSimulationSt } } while (std::abs(sumEhit)<std::abs(Elayer)); - if (debug) { - for(TFCSLateralShapeParametrizationHitBase* reset : m_chain) { - reset->setLevel(old_level); - } - } + if (debug) PropagateMSGLevel(old_level); return FCSSuccess; } -- GitLab