CMS 3D CMS Logo

List of all members | Public Member Functions | Private Member Functions | Private Attributes
EcalHaloAlgo Class Reference

#include <EcalHaloAlgo.h>

Public Member Functions

reco::EcalHaloData Calculate (const CaloGeometry &TheCaloGeometry, edm::Handle< reco::PhotonCollection > &ThePhotons, edm::Handle< reco::SuperClusterCollection > &TheSuperClusters, edm::Handle< EBRecHitCollection > &TheEBRecHits, edm::Handle< EERecHitCollection > &TheEERecHits, edm::Handle< ESRecHitCollection > &TheESRecHits, edm::Handle< HBHERecHitCollection > &TheHBHERecHits, const edm::EventSetup &TheSetup)
 
bool EBClusterShapeandTimeStudy (reco::HaloClusterCandidateECAL hcand, bool ishlt)
 
 EcalHaloAlgo (edm::ConsumesCollector iC)
 
bool EEClusterShapeandTimeStudy_ITBH (reco::HaloClusterCandidateECAL hcand, bool ishlt)
 
bool EEClusterShapeandTimeStudy_OTBH (reco::HaloClusterCandidateECAL hcand, bool ishlt)
 
float GetAngleCut ()
 
float GetEBRecHitEnergyThreshold ()
 
float GetEERecHitEnergyThreshold ()
 
float GetESRecHitEnergyThreshold ()
 
std::vector< reco::HaloClusterCandidateECALGetHaloClusterCandidateEB (edm::Handle< EcalRecHitCollection > &ecalrechitcoll, edm::Handle< HBHERecHitCollection > &hbherechitcoll, float et_thresh_seedrh)
 
std::vector< reco::HaloClusterCandidateECALGetHaloClusterCandidateEE (edm::Handle< EcalRecHitCollection > &ecalrechitcoll, edm::Handle< HBHERecHitCollection > &hbherechitcoll, float et_thresh_seedrh)
 
float GetPhiWedgeEnergyThreshold ()
 
int GetPhiWedgeNHitsThreshold ()
 
float GetRoundnessCut ()
 
void SetAngleCut (float a=4.)
 
void SetPhiWedgeEnergyThreshold (float SumE)
 
void SetPhiWedgeNHitsThreshold (int nhits)
 
void SetPhiWedgeThresholds (float SumE, int nhits)
 
void SetRecHitEnergyThresholds (float EB, float EE, float ES)
 
void SetRoundnessCut (float r=100.)
 
 ~EcalHaloAlgo ()
 

Private Member Functions

math::XYZPoint getPosition (const DetId &id, reco::Vertex::Point vtx)
 

Private Attributes

float AngleCut
 
float EBRecHitEnergyThreshold
 
float EERecHitEnergyThreshold
 
float ESRecHitEnergyThreshold
 
const CaloGeometrygeo
 
edm::ESGetToken< CaloGeometry, CaloGeometryRecordgeoToken_
 
int NHitsThreshold
 
float RoundnessCut
 
float SumEnergyThreshold
 

Detailed Description

Definition at line 44 of file EcalHaloAlgo.h.

Constructor & Destructor Documentation

◆ EcalHaloAlgo()

EcalHaloAlgo::EcalHaloAlgo ( edm::ConsumesCollector  iC)
explicit

Definition at line 20 of file EcalHaloAlgo.cc.

References AngleCut, EBRecHitEnergyThreshold, EERecHitEnergyThreshold, ESRecHitEnergyThreshold, geo, NHitsThreshold, RoundnessCut, and SumEnergyThreshold.

20  : geoToken_(iC.esConsumes()) {
21  RoundnessCut = 0;
22  AngleCut = 0;
26  SumEnergyThreshold = 0.;
27  NHitsThreshold = 0;
28 
29  geo = nullptr;
30 }
edm::ESGetToken< CaloGeometry, CaloGeometryRecord > geoToken_
Definition: EcalHaloAlgo.h:120
float ESRecHitEnergyThreshold
Definition: EcalHaloAlgo.h:114
float SumEnergyThreshold
Definition: EcalHaloAlgo.h:117
float EERecHitEnergyThreshold
Definition: EcalHaloAlgo.h:113
float EBRecHitEnergyThreshold
Definition: EcalHaloAlgo.h:112
float RoundnessCut
Definition: EcalHaloAlgo.h:107
const CaloGeometry * geo
Definition: EcalHaloAlgo.h:121

◆ ~EcalHaloAlgo()

EcalHaloAlgo::~EcalHaloAlgo ( )
inline

Definition at line 49 of file EcalHaloAlgo.h.

49 {}

Member Function Documentation

◆ Calculate()

EcalHaloData EcalHaloAlgo::Calculate ( const CaloGeometry TheCaloGeometry,
edm::Handle< reco::PhotonCollection > &  ThePhotons,
edm::Handle< reco::SuperClusterCollection > &  TheSuperClusters,
edm::Handle< EBRecHitCollection > &  TheEBRecHits,
edm::Handle< EERecHitCollection > &  TheEERecHits,
edm::Handle< ESRecHitCollection > &  TheESRecHits,
edm::Handle< HBHERecHitCollection > &  TheHBHERecHits,
const edm::EventSetup TheSetup 
)

Definition at line 32 of file EcalHaloAlgo.cc.

References funct::abs(), angle(), edm::SortedCollection< T, SORT >::begin(), CompareTime(), EBRecHitEnergyThreshold, DetId::Ecal, edm::SortedCollection< T, SORT >::end(), geo, geoToken_, GetAngleCut(), edm::EventSetup::getData(), GetHaloClusterCandidateEB(), GetHaloClusterCandidateEE(), reco::EcalHaloData::GetPhiWedges(), GetRoundnessCut(), reco::EcalHaloData::GetShowerShapesAngle(), reco::EcalHaloData::GetShowerShapesRoundness(), CaloGeometry::getSubdetectorGeometry(), reco::EcalHaloData::GetSuperClusters(), mps_fire::i, hit::id, EBDetId::ieta(), edm::helper::Filler< Map >::insert(), EBDetId::iphi(), edm::HandleBase::isValid(), dqmiolumiharvest::j, NHitsThreshold, edm::Handle< T >::product(), edm::RefVector< C, T, F >::push_back(), DetId::rawId(), reco::EcalHaloData::setHaloClusterCandidatesEB(), reco::EcalHaloData::setHaloClusterCandidatesEE(), reco::PhiWedge::SetPlusZOriginConfidence(), jetUpdater_cfi::sort, SumEnergyThreshold, and protons_cff::time.

Referenced by reco::EcalHaloDataProducer::produce().

39  {
40  EcalHaloData TheEcalHaloData;
41 
42  // Store energy sum of rechits as a function of iPhi (iphi goes from 1 to 72)
43  float SumE[361];
44  // Store number of rechits as a function of iPhi
45  int NumHits[361];
46  // Store minimum time of rechit as a function of iPhi
47  float MinTimeHits[361];
48  // Store maximum time of rechit as a function of iPhi
49  float MaxTimeHits[361];
50 
51  // initialize
52  for (int i = 0; i < 361; i++) {
53  SumE[i] = 0.;
54  NumHits[i] = 0;
55  MinTimeHits[i] = 9999.;
56  MaxTimeHits[i] = -9999.;
57  }
58 
59  // Loop over EB RecHits
60  for (EBRecHitCollection::const_iterator hit = TheEBRecHits->begin(); hit != TheEBRecHits->end(); hit++) {
61  // Arbitrary threshold to kill noise (needs to be optimized with data)
62  if (hit->energy() < EBRecHitEnergyThreshold)
63  continue;
64 
65  // Get Det Id of the rechit
66  DetId id = DetId(hit->id());
67 
68  // Get EB geometry
69  const CaloSubdetectorGeometry* TheSubGeometry = TheCaloGeometry.getSubdetectorGeometry(DetId::Ecal, 1);
70  EBDetId EcalID(id.rawId());
71  auto cell = (TheSubGeometry) ? (TheSubGeometry->getGeometry(id)) : nullptr;
72 
73  if (cell) {
74  // GlobalPoint globalpos = cell->getPosition();
75  // float r = TMath::Sqrt ( globalpos.y()*globalpos.y() + globalpos.x()*globalpos.x());
76  int iPhi = EcalID.iphi();
77 
78  if (iPhi < 361) // just to be safe
79  {
80  //iPhi = (iPhi-1)/5 +1; // convert ecal iphi to phiwedge iphi (e.g. there are 5 crystal per phi wedge, as in calotowers )
81  SumE[iPhi] += hit->energy();
82  NumHits[iPhi]++;
83 
84  float time = hit->time();
85  MinTimeHits[iPhi] = time < MinTimeHits[iPhi] ? time : MinTimeHits[iPhi];
86  MaxTimeHits[iPhi] = time > MaxTimeHits[iPhi] ? time : MaxTimeHits[iPhi];
87  }
88  }
89  }
90 
91  //for( int iPhi = 1 ; iPhi < 73; iPhi++ )
92  for (int iPhi = 1; iPhi < 361; iPhi++) {
93  if (SumE[iPhi] >= SumEnergyThreshold && NumHits[iPhi] > NHitsThreshold) {
94  // Build PhiWedge and store to EcalHaloData if energy or #hits pass thresholds
95  PhiWedge wedge(SumE[iPhi], iPhi, NumHits[iPhi], MinTimeHits[iPhi], MaxTimeHits[iPhi]);
96 
97  // Loop over rechits again to calculate direction based on timing info
98 
99  // Loop over EB RecHits
100  std::vector<const EcalRecHit*> Hits;
101  for (EBRecHitCollection::const_iterator hit = TheEBRecHits->begin(); hit != TheEBRecHits->end(); hit++) {
102  if (hit->energy() < EBRecHitEnergyThreshold)
103  continue;
104 
105  // Get Det Id of the rechit
106  DetId id = DetId(hit->id());
107  EBDetId EcalID(id.rawId());
108  int Hit_iPhi = EcalID.iphi();
109  //Hit_iPhi = (Hit_iPhi-1)/5 +1; // convert ecal iphi to phiwedge iphi
110  if (Hit_iPhi != iPhi)
111  continue;
112  Hits.push_back(&(*hit));
113  }
114  std::sort(Hits.begin(), Hits.end(), CompareTime);
115  float MinusToPlus = 0.;
116  float PlusToMinus = 0.;
117  for (unsigned int i = 0; i < Hits.size(); i++) {
118  DetId id_i = DetId(Hits[i]->id());
119  EBDetId EcalID_i(id_i.rawId());
120  int ieta_i = EcalID_i.ieta();
121  for (unsigned int j = (i + 1); j < Hits.size(); j++) {
122  DetId id_j = DetId(Hits[j]->id());
123  EBDetId EcalID_j(id_j.rawId());
124  int ieta_j = EcalID_j.ieta();
125  if (ieta_i > ieta_j)
126  PlusToMinus += TMath::Abs(ieta_j - ieta_i);
127  else
128  MinusToPlus += TMath::Abs(ieta_j - ieta_i);
129  }
130  }
131 
132  float PlusZOriginConfidence = (PlusToMinus + MinusToPlus) ? PlusToMinus / (PlusToMinus + MinusToPlus) : -1.;
133  wedge.SetPlusZOriginConfidence(PlusZOriginConfidence);
134  TheEcalHaloData.GetPhiWedges().push_back(wedge);
135  }
136  }
137 
138  std::vector<float> vShowerShapes_Roundness;
139  std::vector<float> vShowerShapes_Angle;
140  if (TheSuperClusters.isValid()) {
141  for (reco::SuperClusterCollection::const_iterator cluster = TheSuperClusters->begin();
142  cluster != TheSuperClusters->end();
143  cluster++) {
144  if (abs(cluster->eta()) <= 1.48) {
145  vector<float> shapes = EcalClusterTools::roundnessBarrelSuperClusters(*cluster, (*TheEBRecHits.product()));
146  float roundness = shapes[0];
147  float angle = shapes[1];
148 
149  // Check if supercluster belongs to photon and passes the cuts on roundness and angle, if so store the reference to it
150  if ((roundness >= 0 && roundness < GetRoundnessCut()) && angle >= 0 && angle < GetAngleCut()) {
151  edm::Ref<SuperClusterCollection> TheClusterRef(TheSuperClusters, cluster - TheSuperClusters->begin());
152  bool BelongsToPhoton = false;
153  if (ThePhotons.isValid()) {
154  for (reco::PhotonCollection::const_iterator iPhoton = ThePhotons->begin(); iPhoton != ThePhotons->end();
155  iPhoton++) {
156  if (iPhoton->isEB())
157  if (TheClusterRef == iPhoton->superCluster()) {
158  BelongsToPhoton = true;
159  break;
160  }
161  }
162  }
163  //Only store refs to suspicious EB SuperClusters which belong to Photons
164  //Showershape variables are more discriminating for these cases
165  if (BelongsToPhoton) {
166  TheEcalHaloData.GetSuperClusters().push_back(TheClusterRef);
167  }
168  }
169  vShowerShapes_Roundness.push_back(shapes[0]);
170  vShowerShapes_Angle.push_back(shapes[1]);
171  } else {
172  vShowerShapes_Roundness.push_back(-1.);
173  vShowerShapes_Angle.push_back(-1.);
174  }
175  }
176 
177  edm::ValueMap<float>::Filler TheRoundnessFiller(TheEcalHaloData.GetShowerShapesRoundness());
178  TheRoundnessFiller.insert(TheSuperClusters, vShowerShapes_Roundness.begin(), vShowerShapes_Roundness.end());
179  TheRoundnessFiller.fill();
180 
181  edm::ValueMap<float>::Filler TheAngleFiller(TheEcalHaloData.GetShowerShapesAngle());
182  TheAngleFiller.insert(TheSuperClusters, vShowerShapes_Angle.begin(), vShowerShapes_Angle.end());
183  TheAngleFiller.fill();
184  }
185 
186  geo = &TheSetup.getData(geoToken_);
187 
188  //Halo cluster building:
189  //Various clusters are built, depending on the subdetector.
190  //In barrel, one looks for deposits narrow in phi.
191  //In endcaps, one looks for localized deposits (dr condition in EE where r =sqrt(dphi*dphi+deta*deta)
192  //H/E condition is also applied in EB.
193  //The halo cluster building step targets a large efficiency (ideally >99%) for beam halo deposits.
194  //These clusters are used as input for the halo pattern finding methods in EcalHaloAlgo and for the CSC-calo matching methods in GlobalHaloAlgo.
195 
196  //Et threshold hardcoded for now. Might one to get it from config
197  std::vector<HaloClusterCandidateECAL> haloclustercands_EB;
198  haloclustercands_EB = GetHaloClusterCandidateEB(TheEBRecHits, TheHBHERecHits, 5);
199 
200  std::vector<HaloClusterCandidateECAL> haloclustercands_EE;
201  haloclustercands_EE = GetHaloClusterCandidateEE(TheEERecHits, TheHBHERecHits, 10);
202 
203  TheEcalHaloData.setHaloClusterCandidatesEB(haloclustercands_EB);
204  TheEcalHaloData.setHaloClusterCandidatesEE(haloclustercands_EE);
205 
206  return TheEcalHaloData;
207 }
edm::ValueMap< float > & GetShowerShapesRoundness()
Definition: EcalHaloData.h:41
edm::ESGetToken< CaloGeometry, CaloGeometryRecord > geoToken_
Definition: EcalHaloAlgo.h:120
edm::RefVector< reco::SuperClusterCollection > & GetSuperClusters()
Definition: EcalHaloData.h:37
edm::ValueMap< float > & GetShowerShapesAngle()
Definition: EcalHaloData.h:44
T const * product() const
Definition: Handle.h:70
void insert(const H &h, I begin, I end)
Definition: ValueMap.h:53
std::vector< EcalRecHit >::const_iterator const_iterator
int ieta() const
get the crystal ieta
Definition: EBDetId.h:49
float SumEnergyThreshold
Definition: EcalHaloAlgo.h:117
float GetRoundnessCut()
Definition: EcalHaloAlgo.h:92
void setHaloClusterCandidatesEE(const std::vector< HaloClusterCandidateECAL > &x)
Definition: EcalHaloData.h:53
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
void setHaloClusterCandidatesEB(const std::vector< HaloClusterCandidateECAL > &x)
Definition: EcalHaloData.h:52
std::vector< reco::HaloClusterCandidateECAL > GetHaloClusterCandidateEB(edm::Handle< EcalRecHitCollection > &ecalrechitcoll, edm::Handle< HBHERecHitCollection > &hbherechitcoll, float et_thresh_seedrh)
bool getData(T &iHolder) const
Definition: EventSetup.h:122
const_iterator begin() const
float GetAngleCut()
Definition: EcalHaloAlgo.h:94
unsigned int id
const_iterator end() const
std::vector< reco::HaloClusterCandidateECAL > GetHaloClusterCandidateEE(edm::Handle< EcalRecHitCollection > &ecalrechitcoll, edm::Handle< HBHERecHitCollection > &hbherechitcoll, float et_thresh_seedrh)
Definition: DetId.h:17
float EBRecHitEnergyThreshold
Definition: EcalHaloAlgo.h:112
constexpr uint32_t rawId() const
get the raw id
Definition: DetId.h:57
bool isValid() const
Definition: HandleBase.h:70
void push_back(value_type const &ref)
Add a Ref<C, T> to the RefVector.
Definition: RefVector.h:67
bool CompareTime(const EcalRecHit *x, const EcalRecHit *y)
Definition: EcalHaloAlgo.cc:18
const CaloSubdetectorGeometry * getSubdetectorGeometry(const DetId &id) const
access the subdetector geometry for the given subdetector directly
Definition: CaloGeometry.cc:34
const CaloGeometry * geo
Definition: EcalHaloAlgo.h:121
const std::vector< PhiWedge > & GetPhiWedges() const
Definition: EcalHaloData.h:33
T angle(T x1, T y1, T z1, T x2, T y2, T z2)
Definition: angle.h:11

◆ EBClusterShapeandTimeStudy()

bool EcalHaloAlgo::EBClusterShapeandTimeStudy ( reco::HaloClusterCandidateECAL  hcand,
bool  ishlt 
)

Definition at line 484 of file EcalHaloAlgo.cc.

References reco::HaloClusterCandidateECAL::getEtStripIPhiSeedMinus1(), reco::HaloClusterCandidateECAL::getEtStripIPhiSeedPlus1(), reco::HaloClusterCandidateECAL::getNbofCrystalsInEta(), reco::HaloClusterCandidateECAL::getSeedEt(), reco::HaloClusterCandidateECAL::getTimeDiscriminator(), and reco::HaloClusterCandidateECAL::setIsHaloFromPattern().

Referenced by GetHaloClusterCandidateEB().

484  {
485  //Conditions on the central strip size in eta.
486  //For low size, extra conditions on seed et, isolation and cluster timing
487  //The time condition only targets IT beam halo.
488  //EB rechits from OT beam halos are typically too late (around 5 ns or more) and seem therefore already cleaned by the reconstruction.
489 
490  if (hcand.getSeedEt() < 5)
491  return false;
492  if (hcand.getNbofCrystalsInEta() < 4)
493  return false;
494  if (hcand.getNbofCrystalsInEta() == 4 && hcand.getSeedEt() < 10)
495  return false;
496  if (hcand.getNbofCrystalsInEta() == 4 && hcand.getEtStripIPhiSeedPlus1() > 0.1 &&
497  hcand.getEtStripIPhiSeedMinus1() > 0.1)
498  return false;
499  if (hcand.getNbofCrystalsInEta() <= 5 && hcand.getTimeDiscriminator() >= 0.)
500  return false;
501 
502  //For HLT, only use conditions without timing and tighten seed et condition
503  if (ishlt && hcand.getNbofCrystalsInEta() <= 5)
504  return false;
505  if (ishlt && hcand.getSeedEt() < 10)
506  return false;
507 
508  hcand.setIsHaloFromPattern(true);
509 
510  return true;
511 }

◆ EEClusterShapeandTimeStudy_ITBH()

bool EcalHaloAlgo::EEClusterShapeandTimeStudy_ITBH ( reco::HaloClusterCandidateECAL  hcand,
bool  ishlt 
)

Definition at line 532 of file EcalHaloAlgo.cc.

References reco::HaloClusterCandidateECAL::getClusterSize(), reco::HaloClusterCandidateECAL::getSeedEt(), reco::HaloClusterCandidateECAL::getSeedR(), reco::HaloClusterCandidateECAL::getTimeDiscriminator(), and reco::HaloClusterCandidateECAL::setIsHaloFromPattern().

Referenced by GetHaloClusterCandidateEE().

532  {
533  //Separate conditions targeting IT and OT beam halos
534  //For IT beam halos, fakes from collisions are higher => require the cluster size to be small.
535  //Only halos with R>100 cm are considered here.
536  //For lower values, the time difference with particles from collisions is too small
537  //IT outgoing beam halos that interact in EE at low R is probably the most difficult category to deal with:
538  //Their signature is very close to the one of photon from collisions (similar cluster shape and timing)
539  if (hcand.getSeedEt() < 20)
540  return false;
541  if (hcand.getSeedR() < 100)
542  return false;
543  if (hcand.getTimeDiscriminator() < 1)
544  return false;
545  if (hcand.getClusterSize() < 2)
546  return false;
547  if (hcand.getClusterSize() > 4)
548  return false;
549 
550  //The use of time information does not allow this method to work at HLT
551  if (ishlt)
552  return false;
553 
554  hcand.setIsHaloFromPattern(true);
555 
556  return true;
557 }

◆ EEClusterShapeandTimeStudy_OTBH()

bool EcalHaloAlgo::EEClusterShapeandTimeStudy_OTBH ( reco::HaloClusterCandidateECAL  hcand,
bool  ishlt 
)

Definition at line 513 of file EcalHaloAlgo.cc.

References reco::HaloClusterCandidateECAL::getNbEarlyCrystals(), reco::HaloClusterCandidateECAL::getNbLateCrystals(), reco::HaloClusterCandidateECAL::getSeedEt(), reco::HaloClusterCandidateECAL::getSeedTime(), and reco::HaloClusterCandidateECAL::setIsHaloFromPattern().

Referenced by GetHaloClusterCandidateEE().

513  {
514  //Separate conditions targeting IT and OT beam halos
515  //For OT beam halos, just require enough crystals with large T
516  if (hcand.getSeedEt() < 20)
517  return false;
518  if (hcand.getSeedTime() < 0.5)
519  return false;
520  if (hcand.getNbLateCrystals() - hcand.getNbEarlyCrystals() < 2)
521  return false;
522 
523  //The use of time information does not allow this method to work at HLT
524  if (ishlt)
525  return false;
526 
527  hcand.setIsHaloFromPattern(true);
528 
529  return true;
530 }

◆ GetAngleCut()

float EcalHaloAlgo::GetAngleCut ( )
inline

Definition at line 94 of file EcalHaloAlgo.h.

References AngleCut.

Referenced by Calculate().

94 { return AngleCut; }

◆ GetEBRecHitEnergyThreshold()

float EcalHaloAlgo::GetEBRecHitEnergyThreshold ( )
inline

Definition at line 97 of file EcalHaloAlgo.h.

References EBRecHitEnergyThreshold.

97 { return EBRecHitEnergyThreshold; }
float EBRecHitEnergyThreshold
Definition: EcalHaloAlgo.h:112

◆ GetEERecHitEnergyThreshold()

float EcalHaloAlgo::GetEERecHitEnergyThreshold ( )
inline

Definition at line 98 of file EcalHaloAlgo.h.

References EERecHitEnergyThreshold.

98 { return EERecHitEnergyThreshold; }
float EERecHitEnergyThreshold
Definition: EcalHaloAlgo.h:113

◆ GetESRecHitEnergyThreshold()

float EcalHaloAlgo::GetESRecHitEnergyThreshold ( )
inline

Definition at line 99 of file EcalHaloAlgo.h.

References ESRecHitEnergyThreshold.

99 { return ESRecHitEnergyThreshold; }
float ESRecHitEnergyThreshold
Definition: EcalHaloAlgo.h:114

◆ GetHaloClusterCandidateEB()

std::vector< HaloClusterCandidateECAL > EcalHaloAlgo::GetHaloClusterCandidateEB ( edm::Handle< EcalRecHitCollection > &  ecalrechitcoll,
edm::Handle< HBHERecHitCollection > &  hbherechitcoll,
float  et_thresh_seedrh 
)

Definition at line 209 of file EcalHaloAlgo.cc.

References funct::abs(), SiPixelRawToDigiRegional_cfi::deltaPhi, EBClusterShapeandTimeStudy(), CaloRecHit::energy(), EcalRecHit::energy(), PVValHelper::eta, getPosition(), electrons_cff::hoe, HBHERecHit::id(), EcalRecHit::id(), EBDetId::ieta(), phi, edm::RefVector< C, T, F >::push_back(), reco::HaloClusterCandidateECAL::setBeamHaloRecHitsCandidates(), reco::HaloClusterCandidateECAL::setClusterEt(), reco::HaloClusterCandidateECAL::setEtStripIPhiSeedMinus1(), reco::HaloClusterCandidateECAL::setEtStripIPhiSeedPlus1(), reco::HaloClusterCandidateECAL::setHoverE(), reco::HaloClusterCandidateECAL::setIsHaloFromPattern(), reco::HaloClusterCandidateECAL::setIsHaloFromPattern_HLT(), reco::HaloClusterCandidateECAL::setNbofCrystalsInEta(), reco::HaloClusterCandidateECAL::setSeedEt(), reco::HaloClusterCandidateECAL::setSeedEta(), reco::HaloClusterCandidateECAL::setSeedPhi(), reco::HaloClusterCandidateECAL::setSeedR(), reco::HaloClusterCandidateECAL::setSeedTime(), reco::HaloClusterCandidateECAL::setSeedZ(), reco::HaloClusterCandidateECAL::setTimeDiscriminator(), edm::SortedCollection< T, SORT >::size(), mathSSE::sqrt(), EcalRecHit::time(), and extraflags_cff::vtx.

Referenced by Calculate().

212  {
213  std::vector<HaloClusterCandidateECAL> TheHaloClusterCandsEB;
214  reco::Vertex::Point vtx(0, 0, 0);
215 
216  for (size_t ihit = 0; ihit < ecalrechitcoll->size(); ++ihit) {
217  HaloClusterCandidateECAL clustercand;
218 
219  const EcalRecHit& rechit = (*ecalrechitcoll)[ihit];
220  math::XYZPoint rhpos = getPosition(rechit.id(), vtx);
221  //Et condition
222 
223  double rhet = rechit.energy() * sqrt(rhpos.perp2() / rhpos.mag2());
224  if (rhet < et_thresh_seedrh)
225  continue;
226  double eta = rhpos.eta();
227  double phi = rhpos.phi();
228 
229  bool isiso = true;
230  double etcluster(0);
231  int nbcrystalsameeta(0);
232  double timediscriminator(0);
233  double etstrip_iphiseedplus1(0), etstrip_iphiseedminus1(0);
234 
235  //Building the cluster
237  for (size_t jhit = 0; jhit < ecalrechitcoll->size(); ++jhit) {
238  const EcalRecHit& rechitj = (*ecalrechitcoll)[jhit];
239  EcalRecHitRef rhRef(ecalrechitcoll, jhit);
240  math::XYZPoint rhposj = getPosition(rechitj.id(), vtx);
241 
242  double etaj = rhposj.eta();
243  double phij = rhposj.phi();
244 
245  double deta = eta - etaj;
246  double dphi = deltaPhi(phi, phij);
247  if (std::abs(deta) > 0.2)
248  continue; //This means +/-11 crystals in eta
249  if (std::abs(dphi) > 0.08)
250  continue; //This means +/-4 crystals in phi
251 
252  double rhetj = rechitj.energy() * sqrt(rhposj.perp2() / rhposj.mag2());
253  //Rechits with et between 1 and 2 GeV are saved in the rh list but not used in the calculation of the halocluster variables
254  if (rhetj < 1)
255  continue;
256  bhrhcandidates.push_back(rhRef);
257  if (rhetj < 2)
258  continue;
259 
260  if (std::abs(dphi) > 0.03) {
261  isiso = false;
262  break;
263  } //The strip should be isolated
264  if (std::abs(dphi) < 0.01)
265  nbcrystalsameeta++;
266  if (dphi > 0.01)
267  etstrip_iphiseedplus1 += rhetj;
268  if (dphi < -0.01)
269  etstrip_iphiseedminus1 += rhetj;
270  etcluster += rhetj;
271  //Timing discriminator
272  //We assign a weight to the rechit defined as:
273  //Log10(Et)*f(T,R,Z)
274  //where f(T,R,Z) is the separation curve between halo-like and IP-like times.
275  //The time difference between a deposit from a outgoing IT halo and a deposit coming from a particle emitted at the IP is given by:
276  //dt= ( - sqrt(R^2+z^2) + |z| )/c
277  //Here we take R to be 130 cm.
278  //For EB, the function was parametrized as a function of ieta instead of Z.
279  double rhtj = rechitj.time();
280  EBDetId detj = rechitj.id();
281  int rhietaj = detj.ieta();
282  timediscriminator += std::log10(rhetj) *
283  (rhtj + 0.5 * (sqrt(16900 + 9 * rhietaj * rhietaj) - 3 * std::abs(rhietaj)) / c_cm_per_ns);
284  }
285  //Isolation condition
286  if (!isiso)
287  continue;
288 
289  //Calculate H/E
290  double hoe(0);
291  for (size_t jhit = 0; jhit < hbherechitcoll->size(); ++jhit) {
292  const HBHERecHit& rechitj = (*hbherechitcoll)[jhit];
293  math::XYZPoint rhposj = getPosition(rechitj.id(), vtx);
294  double rhetj = rechitj.energy() * sqrt(rhposj.perp2() / rhposj.mag2());
295  if (rhetj < 2)
296  continue;
297  double etaj = rhposj.eta();
298  double phij = rhposj.phi();
299  double deta = eta - etaj;
300  double dphi = deltaPhi(phi, phij);
301  if (std::abs(deta) > 0.2)
302  continue;
303  if (std::abs(dphi) > 0.2)
304  continue;
305  hoe += rhetj / etcluster;
306  }
307  //H/E condition
308  if (hoe > 0.1)
309  continue;
310 
311  clustercand.setClusterEt(etcluster);
312  clustercand.setSeedEt(rhet);
313  clustercand.setSeedEta(eta);
314  clustercand.setSeedPhi(phi);
315  clustercand.setSeedZ(rhpos.Z());
316  clustercand.setSeedR(sqrt(rhpos.perp2()));
317  clustercand.setSeedTime(rechit.time());
318  clustercand.setHoverE(hoe);
319  clustercand.setNbofCrystalsInEta(nbcrystalsameeta);
320  clustercand.setEtStripIPhiSeedPlus1(etstrip_iphiseedplus1);
321  clustercand.setEtStripIPhiSeedMinus1(etstrip_iphiseedminus1);
322  clustercand.setTimeDiscriminator(timediscriminator);
323  clustercand.setBeamHaloRecHitsCandidates(bhrhcandidates);
324 
325  bool isbeamhalofrompattern = EBClusterShapeandTimeStudy(clustercand, false);
326  clustercand.setIsHaloFromPattern(isbeamhalofrompattern);
327 
328  bool isbeamhalofrompattern_hlt = EBClusterShapeandTimeStudy(clustercand, true);
329  clustercand.setIsHaloFromPattern_HLT(isbeamhalofrompattern_hlt);
330 
331  TheHaloClusterCandsEB.push_back(clustercand);
332  }
333 
334  return TheHaloClusterCandsEB;
335 }
float time() const
Definition: EcalRecHit.h:70
size_type size() const
constexpr float energy() const
Definition: CaloRecHit.h:29
int ieta() const
get the crystal ieta
Definition: EBDetId.h:49
T sqrt(T t)
Definition: SSEVec.h:19
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
math::XYZPoint Point
point in the space
Definition: Vertex.h:40
constexpr HcalDetId id() const
get the id
Definition: HBHERecHit.h:41
XYZPointD XYZPoint
point in space with cartesian internal representation
Definition: Point3D.h:12
DetId id() const
get the id
Definition: EcalRecHit.h:77
math::XYZPoint getPosition(const DetId &id, reco::Vertex::Point vtx)
void push_back(value_type const &ref)
Add a Ref<C, T> to the RefVector.
Definition: RefVector.h:67
void setBeamHaloRecHitsCandidates(edm::RefVector< EcalRecHitCollection > x)
bool EBClusterShapeandTimeStudy(reco::HaloClusterCandidateECAL hcand, bool ishlt)
float energy() const
Definition: EcalRecHit.h:68

◆ GetHaloClusterCandidateEE()

std::vector< HaloClusterCandidateECAL > EcalHaloAlgo::GetHaloClusterCandidateEE ( edm::Handle< EcalRecHitCollection > &  ecalrechitcoll,
edm::Handle< HBHERecHitCollection > &  hbherechitcoll,
float  et_thresh_seedrh 
)

Definition at line 337 of file EcalHaloAlgo.cc.

References funct::abs(), SiPixelRawToDigiRegional_cfi::deltaPhi, flavorHistoryFilter_cfi::dr, EEClusterShapeandTimeStudy_ITBH(), EEClusterShapeandTimeStudy_OTBH(), CaloRecHit::energy(), EcalRecHit::energy(), PVValHelper::eta, getPosition(), HBHERecHit::id(), EcalRecHit::id(), phi, funct::pow(), edm::RefVector< C, T, F >::push_back(), reco::HaloClusterCandidateECAL::setBeamHaloRecHitsCandidates(), reco::HaloClusterCandidateECAL::setClusterEt(), reco::HaloClusterCandidateECAL::setClusterSize(), reco::HaloClusterCandidateECAL::setH2overE(), reco::HaloClusterCandidateECAL::setIsHaloFromPattern(), reco::HaloClusterCandidateECAL::setIsHaloFromPattern_HLT(), reco::HaloClusterCandidateECAL::setNbEarlyCrystals(), reco::HaloClusterCandidateECAL::setNbLateCrystals(), reco::HaloClusterCandidateECAL::setSeedEt(), reco::HaloClusterCandidateECAL::setSeedEta(), reco::HaloClusterCandidateECAL::setSeedPhi(), reco::HaloClusterCandidateECAL::setSeedR(), reco::HaloClusterCandidateECAL::setSeedTime(), reco::HaloClusterCandidateECAL::setSeedZ(), reco::HaloClusterCandidateECAL::setTimeDiscriminator(), edm::SortedCollection< T, SORT >::size(), mathSSE::sqrt(), EcalRecHit::time(), and extraflags_cff::vtx.

Referenced by Calculate().

340  {
341  std::vector<HaloClusterCandidateECAL> TheHaloClusterCandsEE;
342 
343  reco::Vertex::Point vtx(0, 0, 0);
344 
345  for (size_t ihit = 0; ihit < ecalrechitcoll->size(); ++ihit) {
346  HaloClusterCandidateECAL clustercand;
347 
348  const EcalRecHit& rechit = (*ecalrechitcoll)[ihit];
349  math::XYZPoint rhpos = getPosition(rechit.id(), vtx);
350  //Et condition
351  double rhet = rechit.energy() * sqrt(rhpos.perp2() / rhpos.mag2());
352  if (rhet < et_thresh_seedrh)
353  continue;
354  double eta = rhpos.eta();
355  double phi = rhpos.phi();
356  double rhr = sqrt(rhpos.perp2());
357 
358  bool isiso = true;
359  double etcluster(0);
360  double timediscriminator(0);
361  int clustersize(0);
362  int nbcrystalssmallt(0);
363  int nbcrystalshight(0);
364  //Building the cluster
366  for (size_t jhit = 0; jhit < ecalrechitcoll->size(); ++jhit) {
367  const EcalRecHit& rechitj = (*ecalrechitcoll)[jhit];
368  EcalRecHitRef rhRef(ecalrechitcoll, jhit);
369  math::XYZPoint rhposj = getPosition(rechitj.id(), vtx);
370 
371  //Ask the hits to be in the same endcap
372  if (rhposj.z() * rhpos.z() < 0)
373  continue;
374 
375  double etaj = rhposj.eta();
376  double phij = rhposj.phi();
377  double dr = sqrt((eta - etaj) * (eta - etaj) + deltaPhi(phi, phij) * deltaPhi(phi, phij));
378 
379  //Outer cone
380  if (dr > 0.3)
381  continue;
382 
383  double rhetj = rechitj.energy() * sqrt(rhposj.perp2() / rhposj.mag2());
384  //Rechits with et between 1 and 2 GeV are saved in the rh list but not used in the calculation of the halocluster variables
385  if (rhetj < 1)
386  continue;
387  bhrhcandidates.push_back(rhRef);
388  if (rhetj < 2)
389  continue;
390 
391  //Isolation between outer and inner cone
392  if (dr > 0.05) {
393  isiso = false;
394  break;
395  } //The deposit should be isolated
396 
397  etcluster += rhetj;
398 
399  //Timing infos:
400  //Here we target both IT and OT beam halo
401  double rhtj = rechitj.time();
402 
403  //Discriminating variables for OT beam halo:
404  if (rhtj > 1)
405  nbcrystalshight++;
406  if (rhtj < 0)
407  nbcrystalssmallt++;
408  //Timing test (likelihood ratio), only for seeds with large R (100 cm) and for crystals with et>5,
409  //This targets IT beam halo (t around - 1ns)
410  if (rhtj > 5) {
411  double corrt_j = rhtj + sqrt(rhposj.x() * rhposj.x() + rhposj.y() * rhposj.y() + 320. * 320.) / c_cm_per_ns -
412  320. / c_cm_per_ns;
413  //BH is modeled by a Gaussian peaking at 0.
414  //Collisions is modeled by a Gaussian peaking at 0.3
415  //The width is similar and taken to be 0.4
416  timediscriminator += 0.5 * (pow((corrt_j - 0.3) / 0.4, 2) - pow((corrt_j - 0.) / 0.4, 2));
417  clustersize++;
418  }
419  }
420  //Isolation condition
421  if (!isiso)
422  continue;
423 
424  //Calculate H2/E
425  //Only second hcal layer is considered as it can happen that a shower initiated in EE reaches HCAL first layer
426  double h2oe(0);
427  for (size_t jhit = 0; jhit < hbherechitcoll->size(); ++jhit) {
428  const HBHERecHit& rechitj = (*hbherechitcoll)[jhit];
429  math::XYZPoint rhposj = getPosition(rechitj.id(), vtx);
430 
431  //Ask the hits to be in the same endcap
432  if (rhposj.z() * rhpos.z() < 0)
433  continue;
434  //Selects only second HCAL layer
435  if (std::abs(rhposj.z()) < 425)
436  continue;
437 
438  double rhetj = rechitj.energy() * sqrt(rhposj.perp2() / rhposj.mag2());
439  if (rhetj < 2)
440  continue;
441 
442  double phij = rhposj.phi();
443  if (std::abs(deltaPhi(phi, phij)) > 0.4)
444  continue;
445 
446  double rhrj = sqrt(rhposj.perp2());
447  if (std::abs(rhr - rhrj) > 50)
448  continue;
449 
450  h2oe += rhetj / etcluster;
451  }
452  //H/E condition
453  if (h2oe > 0.1)
454  continue;
455 
456  clustercand.setClusterEt(etcluster);
457  clustercand.setSeedEt(rhet);
458  clustercand.setSeedEta(eta);
459  clustercand.setSeedPhi(phi);
460  clustercand.setSeedZ(rhpos.Z());
461  clustercand.setSeedR(sqrt(rhpos.perp2()));
462  clustercand.setSeedTime(rechit.time());
463  clustercand.setH2overE(h2oe);
464  clustercand.setNbEarlyCrystals(nbcrystalssmallt);
465  clustercand.setNbLateCrystals(nbcrystalshight);
466  clustercand.setClusterSize(clustersize);
467  clustercand.setTimeDiscriminator(timediscriminator);
468  clustercand.setBeamHaloRecHitsCandidates(bhrhcandidates);
469 
470  bool isbeamhalofrompattern =
471  EEClusterShapeandTimeStudy_ITBH(clustercand, false) || EEClusterShapeandTimeStudy_OTBH(clustercand, false);
472  clustercand.setIsHaloFromPattern(isbeamhalofrompattern);
473 
474  bool isbeamhalofrompattern_hlt =
475  EEClusterShapeandTimeStudy_ITBH(clustercand, true) || EEClusterShapeandTimeStudy_OTBH(clustercand, true);
476  clustercand.setIsHaloFromPattern_HLT(isbeamhalofrompattern_hlt);
477 
478  TheHaloClusterCandsEE.push_back(clustercand);
479  }
480 
481  return TheHaloClusterCandsEE;
482 }
float time() const
Definition: EcalRecHit.h:70
size_type size() const
bool EEClusterShapeandTimeStudy_OTBH(reco::HaloClusterCandidateECAL hcand, bool ishlt)
constexpr float energy() const
Definition: CaloRecHit.h:29
T sqrt(T t)
Definition: SSEVec.h:19
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
math::XYZPoint Point
point in the space
Definition: Vertex.h:40
constexpr HcalDetId id() const
get the id
Definition: HBHERecHit.h:41
XYZPointD XYZPoint
point in space with cartesian internal representation
Definition: Point3D.h:12
bool EEClusterShapeandTimeStudy_ITBH(reco::HaloClusterCandidateECAL hcand, bool ishlt)
DetId id() const
get the id
Definition: EcalRecHit.h:77
math::XYZPoint getPosition(const DetId &id, reco::Vertex::Point vtx)
void push_back(value_type const &ref)
Add a Ref<C, T> to the RefVector.
Definition: RefVector.h:67
void setBeamHaloRecHitsCandidates(edm::RefVector< EcalRecHitCollection > x)
float energy() const
Definition: EcalRecHit.h:68
Power< A, B >::type pow(const A &a, const B &b)
Definition: Power.h:29

◆ GetPhiWedgeEnergyThreshold()

float EcalHaloAlgo::GetPhiWedgeEnergyThreshold ( )
inline

Definition at line 102 of file EcalHaloAlgo.h.

References SumEnergyThreshold.

102 { return SumEnergyThreshold; }
float SumEnergyThreshold
Definition: EcalHaloAlgo.h:117

◆ GetPhiWedgeNHitsThreshold()

int EcalHaloAlgo::GetPhiWedgeNHitsThreshold ( )
inline

Definition at line 103 of file EcalHaloAlgo.h.

References NHitsThreshold.

103 { return NHitsThreshold; }

◆ getPosition()

math::XYZPoint EcalHaloAlgo::getPosition ( const DetId id,
reco::Vertex::Point  vtx 
)
private

Definition at line 559 of file EcalHaloAlgo.cc.

References geo, CaloGeometry::getPosition(), and extraflags_cff::vtx.

Referenced by GetHaloClusterCandidateEB(), and GetHaloClusterCandidateEE().

559  {
560  const GlobalPoint& pos = geo->getPosition(id);
561  math::XYZPoint posV(pos.x() - vtx.x(), pos.y() - vtx.y(), pos.z() - vtx.z());
562  return posV;
563 }
GlobalPoint getPosition(const DetId &id) const
Get the position of a given detector id.
Definition: CaloGeometry.cc:50
XYZPointD XYZPoint
point in space with cartesian internal representation
Definition: Point3D.h:12
const CaloGeometry * geo
Definition: EcalHaloAlgo.h:121

◆ GetRoundnessCut()

float EcalHaloAlgo::GetRoundnessCut ( )
inline

Definition at line 92 of file EcalHaloAlgo.h.

References RoundnessCut.

Referenced by Calculate().

92 { return RoundnessCut; }
float RoundnessCut
Definition: EcalHaloAlgo.h:107

◆ SetAngleCut()

void EcalHaloAlgo::SetAngleCut ( float  a = 4.)
inline

Definition at line 63 of file EcalHaloAlgo.h.

References a, and AngleCut.

Referenced by reco::EcalHaloDataProducer::produce().

63 { AngleCut = a; }
double a
Definition: hdecay.h:119

◆ SetPhiWedgeEnergyThreshold()

void EcalHaloAlgo::SetPhiWedgeEnergyThreshold ( float  SumE)
inline

Definition at line 72 of file EcalHaloAlgo.h.

References SumEnergyThreshold.

72 { SumEnergyThreshold = SumE; }
float SumEnergyThreshold
Definition: EcalHaloAlgo.h:117

◆ SetPhiWedgeNHitsThreshold()

void EcalHaloAlgo::SetPhiWedgeNHitsThreshold ( int  nhits)
inline

Definition at line 73 of file EcalHaloAlgo.h.

References nhits, and NHitsThreshold.

◆ SetPhiWedgeThresholds()

void EcalHaloAlgo::SetPhiWedgeThresholds ( float  SumE,
int  nhits 
)
inline

Definition at line 74 of file EcalHaloAlgo.h.

References nhits, NHitsThreshold, and SumEnergyThreshold.

Referenced by reco::EcalHaloDataProducer::produce().

74  {
75  SumEnergyThreshold = SumE;
77  }
float SumEnergyThreshold
Definition: EcalHaloAlgo.h:117

◆ SetRecHitEnergyThresholds()

void EcalHaloAlgo::SetRecHitEnergyThresholds ( float  EB,
float  EE,
float  ES 
)
inline

◆ SetRoundnessCut()

void EcalHaloAlgo::SetRoundnessCut ( float  r = 100.)
inline

Definition at line 61 of file EcalHaloAlgo.h.

References alignCSCRings::r, and RoundnessCut.

Referenced by reco::EcalHaloDataProducer::produce().

Member Data Documentation

◆ AngleCut

float EcalHaloAlgo::AngleCut
private

Definition at line 109 of file EcalHaloAlgo.h.

Referenced by EcalHaloAlgo(), GetAngleCut(), and SetAngleCut().

◆ EBRecHitEnergyThreshold

float EcalHaloAlgo::EBRecHitEnergyThreshold
private

◆ EERecHitEnergyThreshold

float EcalHaloAlgo::EERecHitEnergyThreshold
private

◆ ESRecHitEnergyThreshold

float EcalHaloAlgo::ESRecHitEnergyThreshold
private

◆ geo

const CaloGeometry* EcalHaloAlgo::geo
private

Definition at line 121 of file EcalHaloAlgo.h.

Referenced by Calculate(), EcalHaloAlgo(), and getPosition().

◆ geoToken_

edm::ESGetToken<CaloGeometry, CaloGeometryRecord> EcalHaloAlgo::geoToken_
private

Definition at line 120 of file EcalHaloAlgo.h.

Referenced by Calculate().

◆ NHitsThreshold

int EcalHaloAlgo::NHitsThreshold
private

◆ RoundnessCut

float EcalHaloAlgo::RoundnessCut
private

Definition at line 107 of file EcalHaloAlgo.h.

Referenced by EcalHaloAlgo(), GetRoundnessCut(), and SetRoundnessCut().

◆ SumEnergyThreshold

float EcalHaloAlgo::SumEnergyThreshold
private