Skip to content
Commits on Source (9)
......@@ -2291,8 +2291,8 @@ StatusCode RpcDigitizationTool::DetectionEfficiency(const Identifier* IdEtaRpcSt
OnlyEtaEff = OnlyEtaEff * eff_sf;
OnlyPhiEff = OnlyPhiEff * eff_sf;
ATH_MSG_INFO ( "Scale factor applied to fractional charged particle : " << eff_sf );
}
//std::cout<<"lqy test: qcharge = "<<qcharge<<" , qbetagamma = "<<qbetagamma<<" , scale factor = "<<eff_sf<<std::endl;
}
......
......@@ -127,7 +127,7 @@ bool DerivationFramework::VHLowTrackJetFilterTool::eventPassesFilter() const
for (auto electron : *electrons){
if(electron->pt()<m_electronPtCut) continue;
if( (fabs(electron->caloCluster()->etaBE(2))>1.37 && fabs(electron->caloCluster()->etaBE(2)<1.52))
if( (fabs(electron->caloCluster()->etaBE(2))>1.37 && fabs(electron->caloCluster()->etaBE(2))<1.52)
|| fabs(electron->caloCluster()->etaBE(2))>2.47 ) continue;
if(electron->isolation(xAOD::Iso::topoetcone20)/electron->pt()>0.2) continue;
......
21.0.105.1
\ No newline at end of file
21.0.106.1
\ No newline at end of file
21.0.106
\ No newline at end of file
21.0.107
\ No newline at end of file
......@@ -21,12 +21,19 @@ public:
bool Initialize(TH1* hist);
virtual void Print(Option_t *option = "") const override;
virtual void setMinWeight(float minWeight) {m_minWeight=minWeight;}
virtual void setMaxWeight(float maxWeight) {m_maxWeight=maxWeight;}
virtual float getMinWeight() const override;
virtual float getMaxWeight() const override;
protected:
/// Histogram to be used for the shape simulation
/// The histogram x-axis should be in dR^2=deta^2+dphi^2
TH1* m_hist{nullptr};
float m_minWeight{-1.};
float m_maxWeight{-1.};
ClassDefOverride(TFCSHistoLateralShapeWeight,1) //TFCSHistoLateralShapeWeight
ClassDefOverride(TFCSHistoLateralShapeWeight,2) //TFCSHistoLateralShapeWeight
};
#if defined(__ROOTCLING__) && defined(__FastCaloSimStandAlone__)
......
......@@ -15,6 +15,9 @@ public:
// set and get the amount of hit energy fluctation around E/n for n hits
void set_hit_RMS(float RMS);
float get_hit_RMS() const {return m_RMS;};
///Get hit energy from layer energy and number of hits
virtual float get_E_hit(TFCSSimulationState& simulstate,const TFCSTruthState* truth, const TFCSExtrapolationState* extrapol) const override;
virtual FCSReturnCode simulate(TFCSSimulationState& simulstate,const TFCSTruthState* truth, const TFCSExtrapolationState* extrapol) const override;
......
......@@ -18,6 +18,13 @@ public:
///Call get_number_of_hits() only once per shower simulation, as it could be calculated with random numbers and give different results each time. Return a value of -1 if this instance can't determine
virtual int get_number_of_hits(TFCSSimulationState& simulstate,const TFCSTruthState* truth, const TFCSExtrapolationState* extrapol) const;
///Get hit energy from layer energy and number of hits
virtual float get_E_hit(TFCSSimulationState& simulstate,const TFCSTruthState* truth, const TFCSExtrapolationState* extrapol) const;
///Get minimum and maximum value of weight for hit energy reweighting
virtual float getMinWeight() const;
virtual float getMaxWeight() const;
class Hit
{
......
......@@ -33,9 +33,16 @@ public:
/// Call get_number_of_hits() only once, as it could contain a random number
virtual int get_number_of_hits(TFCSSimulationState& simulstate,const TFCSTruthState* truth, const TFCSExtrapolationState* extrapol) const;
///Get hit energy from layer energy and number of hits
virtual float get_E_hit(TFCSSimulationState& simulstate,const TFCSTruthState* truth, const TFCSExtrapolationState* extrapol) const;
///Give the effective size sigma^2 of the fluctuations that should be generated.
virtual float get_sigma2_fluctuation(TFCSSimulationState& simulstate,const TFCSTruthState* truth, const TFCSExtrapolationState* extrapol) const;
///Get minimum and maximum value of weight for hit energy reweighting
virtual float getMinWeight() const;
virtual float getMaxWeight() const;
static constexpr float s_max_sigma2_fluctuation=1000;//! Do not persistify!
......
......@@ -18,14 +18,23 @@
TFCSHistoLateralShapeWeight::TFCSHistoLateralShapeWeight(const char* name, const char* title) :
TFCSLateralShapeParametrizationHitBase(name,title)
{
}
{}
TFCSHistoLateralShapeWeight::~TFCSHistoLateralShapeWeight()
{
if(m_hist) delete m_hist;
}
float TFCSHistoLateralShapeWeight::getMinWeight() const
{
return m_minWeight;
}
float TFCSHistoLateralShapeWeight::getMaxWeight() const
{
return m_maxWeight;
}
FCSReturnCode TFCSHistoLateralShapeWeight::simulate_hit(Hit& hit,TFCSSimulationState& simulstate,const TFCSTruthState* /*truth*/, const TFCSExtrapolationState* /*extrapol*/)
{
if (!simulstate.randomEngine()) {
......@@ -82,7 +91,7 @@ void TFCSHistoLateralShapeWeight::Print(Option_t *option) const
TFCSLateralShapeParametrizationHitBase::Print(option);
if(longprint) {
if(m_hist) ATH_MSG_INFO(optprint <<" Histogram: "<<m_hist->GetNbinsX()<<" bins ["<<m_hist->GetXaxis()->GetXmin()<<","<<m_hist->GetXaxis()->GetXmax()<<"]");
if(m_hist) ATH_MSG_INFO(optprint <<" Histogram: "<<m_hist->GetNbinsX()<<" bins ["<<m_hist->GetXaxis()->GetXmin()<<","<<m_hist->GetXaxis()->GetXmax()<<"]" << " min weight: " << m_minWeight << " max weight: " << m_maxWeight);
else ATH_MSG_INFO(optprint <<" no Histogram");
}
}
......
......@@ -55,20 +55,37 @@ FCSReturnCode TFCSHistoLateralShapeWeightHitAndMiss::simulate_hit(Hit& hit,TFCSS
if(RMS>0) {
weight=CLHEP::RandGauss::shoot(simulstate.randomEngine(), meanweight, RMS);
}
if(meanweight<=1) {
//if meanweight<=1, give lower energy to hit.
//TFCSLateralShapeParametrizationHitChain needs to be able to generate more hits in this case
hit.E()*=weight;
} else {
//if meanweight>1, accept only 1/meanweight events, but give them a higher energy increased by weight.
//This leads to larger fluctuations, while keeping the shape unchanged.
float prob=1.0/meanweight;
/* -------------------------------------------------------------------
* Weight is used to scale hit energy.
*
* if (meanweight > m_minWeight and meanweight < m_maxWeight)
* Hit is accecpted with probability of m_minWeight/meanweight.
* If not accepted, weight is set to zero (this is
* equivalent to not accept the hit).
*
* if (meanweight<=m_minWeight)
* Give lower energy to hit with probability 1.
* TFCSLateralShapeParametrizationHitChain needs to be able
* to generate more hits in this case
*
* if (meanweight >= m_maxWeight)
* meanweight is above upper threshold. It is set to m_maxWeight.
* -------------------------------------------------------------------
*/
if (meanweight >= m_maxWeight) {
weight = m_maxWeight;
}
else if (meanweight > m_minWeight){
float prob=m_minWeight/meanweight;
float rnd=CLHEP::RandFlat::shoot(simulstate.randomEngine());
if(rnd<prob) hit.E()*=weight;
else hit.E()=0;
if(rnd>=prob) weight = 0.;
}
hit.E()*=weight;
ATH_MSG_DEBUG("HIT: E="<<hit.E()<<" dR_mm="<<delta_r_mm<<" meanweight="<<meanweight<<" weight="<<weight);
return FCSSuccess;
}
......
......@@ -21,6 +21,18 @@ TFCSLateralShapeParametrizationFluctChain::TFCSLateralShapeParametrizationFluctC
{
}
float TFCSLateralShapeParametrizationFluctChain::get_E_hit(TFCSSimulationState& simulstate,const TFCSTruthState* truth, const TFCSExtrapolationState* extrapol) const
{
const float sigma2 = get_sigma2_fluctuation(simulstate,truth,extrapol);
const int sample = calosample();
if(sigma2<=0 || sample<0) return -1.;
const float maxWeight = getMaxWeight();
if(maxWeight>0) return simulstate.E(sample)*sigma2/maxWeight;
else return simulstate.E(sample)*sigma2;
}
FCSReturnCode TFCSLateralShapeParametrizationFluctChain::simulate(TFCSSimulationState& simulstate,const TFCSTruthState* truth, const TFCSExtrapolationState* extrapol) const
{
// Call get_number_of_hits() only once, as it could contain a random number
......@@ -34,10 +46,9 @@ FCSReturnCode TFCSLateralShapeParametrizationFluctChain::simulate(TFCSSimulation
if(sigma2<1e-8) sigma2=1e-8;
const float Elayer=simulstate.E(calosample());
const int nhit = 1.0 / sigma2;
//Make a good guess of the needed hit energy, assuming all hits would have the same energy
const float Eavghit=Elayer/nhit;
const float Eavghit=get_E_hit(simulstate,truth,extrapol);
const float Eavghit_tenth=Eavghit/10;
float sumEhit=0;
float error2_sumEhit=0;
......
......@@ -26,6 +26,24 @@ int TFCSLateralShapeParametrizationHitBase::get_number_of_hits(TFCSSimulationSta
return -1;
}
float TFCSLateralShapeParametrizationHitBase::get_E_hit(TFCSSimulationState& simulstate,const TFCSTruthState* truth, const TFCSExtrapolationState* extrapol) const
{
const int nhits = get_number_of_hits(simulstate,truth,extrapol);
const int sample = calosample();
if(nhits<=0 || sample<0) return -1.;
else return simulstate.E(sample)/nhits;
}
float TFCSLateralShapeParametrizationHitBase::getMinWeight() const
{
return -1.;
}
float TFCSLateralShapeParametrizationHitBase::getMaxWeight() const
{
return -1.;
}
FCSReturnCode TFCSLateralShapeParametrizationHitBase::simulate_hit(Hit& hit,TFCSSimulationState& /*simulstate*/,const TFCSTruthState* /*truth*/, const TFCSExtrapolationState* extrapol)
{
int cs=calosample();
......
......@@ -43,6 +43,36 @@ int TFCSLateralShapeParametrizationHitChain::get_number_of_hits(TFCSSimulationSt
return 1;
}
float TFCSLateralShapeParametrizationHitChain::getMinWeight() const
{
for(TFCSLateralShapeParametrizationHitBase* hitsim : m_chain) {
float weight=hitsim->getMinWeight();
if(weight>0.) return weight;
}
return -1.;
}
float TFCSLateralShapeParametrizationHitChain::getMaxWeight() const
{
for(TFCSLateralShapeParametrizationHitBase* hitsim : m_chain) {
float weight=hitsim->getMaxWeight();
if(weight>0.) return weight;
}
return -1.;
}
float TFCSLateralShapeParametrizationHitChain::get_E_hit(TFCSSimulationState& simulstate,const TFCSTruthState* truth, const TFCSExtrapolationState* extrapol) const
{
const int nhits = get_number_of_hits(simulstate,truth,extrapol);
const int sample = calosample();
if(nhits<=0 || sample<0) return -1.;
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
else return simulstate.E(sample)/nhits; // Otherwise, old definition of E_hit is used
}
float TFCSLateralShapeParametrizationHitChain::get_sigma2_fluctuation(TFCSSimulationState& simulstate,const TFCSTruthState* truth, const TFCSExtrapolationState* extrapol) const
{
if(m_number_of_hits_simul) {
......@@ -67,7 +97,7 @@ FCSReturnCode TFCSLateralShapeParametrizationHitChain::simulate(TFCSSimulationSt
}
const float Elayer=simulstate.E(calosample());
const float Ehit=Elayer/nhit;
const float Ehit=get_E_hit(simulstate,truth,extrapol);
float sumEhit=0;
const bool debug = msgLvl(MSG::DEBUG);
......@@ -107,8 +137,12 @@ FCSReturnCode TFCSLateralShapeParametrizationHitChain::simulate(TFCSSimulationSt
}
sumEhit+=hit.E();
++ihit;
if(ihit>100*nhit && ihit>1000) {
ATH_MSG_WARNING("TFCSLateralShapeParametrizationHitChain::simulate(): aborting hit chain, iterated " << 100*nhit << " times, expected " << nhit<<" times. Deposited E("<<calosample()<<")="<<sumEhit<<" expected E="<<Elayer);
if( (ihit==10*nhit) || (ihit==100*nhit) ) {
ATH_MSG_WARNING("TFCSLateralShapeParametrizationHitChain::simulate(): Iterated " << ihit << " times, expected " << nhit <<" times. Deposited E("<<calosample()<<")="<<sumEhit);
}
if(ihit>1000*nhit && ihit>1000) {
ATH_MSG_WARNING("TFCSLateralShapeParametrizationHitChain::simulate(): Aborting hit chain, iterated " << 1000*nhit << " times, expected " << nhit <<" times. Deposited E("<<calosample()<<")="<<sumEhit<<" expected E="<<Elayer);
break;
}
} while (sumEhit<Elayer);
......