CMS 3D CMS Logo
List of all members | Public Types | Public Member Functions | Private Attributes
PreshowerClusterAlgo Class Reference

#include <PreshowerClusterAlgo.h>

Public Types

typedef std::set< DetId > HitsID
 
typedef math::XYZPoint Point
 
typedef std::map< DetId, EcalRecHit > RecHitsMap
 

Public Member Functions

void findRoad (ESDetId strip, EcalPreshowerNavigator theESNav, int plane)
 
bool goodStrip (RecHitsMap::iterator candidate_it)
 
reco::PreshowerCluster makeOneCluster (ESDetId strip, HitsID *used_strips, RecHitsMap *rechits_map, const CaloSubdetectorGeometry *geometry_p, const CaloSubdetectorTopology *topology_p)
 
 PreshowerClusterAlgo ()
 
 PreshowerClusterAlgo (double stripEnergyCut, double clusterEnergyCut, int nStripCut)
 
 ~PreshowerClusterAlgo ()
 

Private Attributes

double preshClusterEnergyCut_
 
int preshSeededNstr_
 
double preshStripEnergyCut_
 
RecHitsMap * rechits_map
 
std::vector< ESDetId > road_2d
 
HitsID * used_s
 

Detailed Description

Definition at line 22 of file PreshowerClusterAlgo.h.

Member Typedef Documentation

typedef std::set<DetId> PreshowerClusterAlgo::HitsID

Definition at line 27 of file PreshowerClusterAlgo.h.

typedef math::XYZPoint PreshowerClusterAlgo::Point

Definition at line 24 of file PreshowerClusterAlgo.h.

typedef std::map<DetId, EcalRecHit> PreshowerClusterAlgo::RecHitsMap

Definition at line 26 of file PreshowerClusterAlgo.h.

Constructor & Destructor Documentation

PreshowerClusterAlgo::PreshowerClusterAlgo ( )
inline

Definition at line 29 of file PreshowerClusterAlgo.h.

29 : preshStripEnergyCut_(0.), preshClusterEnergyCut_(0.), preshSeededNstr_(15) {}
PreshowerClusterAlgo::PreshowerClusterAlgo ( double  stripEnergyCut,
double  clusterEnergyCut,
int  nStripCut 
)
inline

Definition at line 31 of file PreshowerClusterAlgo.h.

32  : preshStripEnergyCut_(stripEnergyCut), preshClusterEnergyCut_(clusterEnergyCut), preshSeededNstr_(nStripCut) {}
PreshowerClusterAlgo::~PreshowerClusterAlgo ( )
inline

Definition at line 34 of file PreshowerClusterAlgo.h.

References findRoad(), goodStrip(), makeOneCluster(), rechits_map, and digitizers_cfi::strip.

34 {};

Member Function Documentation

void PreshowerClusterAlgo::findRoad ( ESDetId  strip,
EcalPreshowerNavigator  theESNav,
int  plane 
)

Definition at line 269 of file PreshowerClusterAlgo.cc.

References CaloNavigator< T, TOPO >::east(), CaloNavigator< T, TOPO >::home(), LogTrace, GetRecoTauVFromDQM_MC_cff::next, CaloNavigator< T, TOPO >::north(), preshSeededNstr_, road_2d, CaloNavigator< T, TOPO >::setHome(), CaloNavigator< T, TOPO >::south(), digitizers_cfi::strip, and CaloNavigator< T, TOPO >::west().

Referenced by makeOneCluster(), and ~PreshowerClusterAlgo().

269  {
270  if (strip == ESDetId(0))
271  return;
272 
273  ESDetId next;
274  theESNav.setHome(strip);
275  // First, add a central strip to the road
276  road_2d.push_back(strip);
277  LogTrace("PreShowerClusterAlgo") << "findRoad starts from strip" << strip;
278 
279  if (plane == 1) {
280  // east road
281  int n_east = 0;
282  LogTrace("PreShowerClusterAlgo") << " Go to the East ";
283 
284  while (((next = theESNav.east()) != ESDetId(0) && next != strip)) {
285  LogTrace("PreShowerClusterAlgo") << "East:" << n_east << "current strip is" << next;
286 
287  road_2d.push_back(next);
288  ++n_east;
289  if (n_east == preshSeededNstr_)
290  break;
291  }
292  // west road
293  int n_west = 0;
294  LogTrace("PreShowerClusterAlgo") << "Go to the West";
295 
296  theESNav.home();
297  while (((next = theESNav.west()) != ESDetId(0) && next != strip)) {
298  LogTrace("PreShowerClusterAlgo") << "West: " << n_west << "current strip is" << next;
299 
300  road_2d.push_back(next);
301  ++n_west;
302  if (n_west == preshSeededNstr_)
303  break;
304  }
305  LogTrace("PreShowerClusterAlgo") << "Total number of strips found in the road at 1-st plane is" << n_east + n_west;
306 
307  } else if (plane == 2) {
308  // north road
309  int n_north = 0;
310  LogTrace("PreShowerClusterAlgo") << "Go to the North";
311 
312  while (((next = theESNav.north()) != ESDetId(0) && next != strip)) {
313  LogTrace("PreShowerClusterAlgo") << "North:" << n_north << "current strip is" << next;
314 
315  road_2d.push_back(next);
316  ++n_north;
317  if (n_north == preshSeededNstr_)
318  break;
319  }
320  // south road
321  int n_south = 0;
322  LogTrace("PreShowerClusterAlgo") << "Go to the South";
323 
324  theESNav.home();
325  while (((next = theESNav.south()) != ESDetId(0) && next != strip)) {
326  LogTrace("PreShowerClusterAlgo") << "South:" << n_south << "current strip is" << next;
327 
328  road_2d.push_back(next);
329  ++n_south;
330  if (n_south == preshSeededNstr_)
331  break;
332  }
333 
334  LogTrace("PreShowerClusterAlgo") << "Total number of strips found in the road at 2-nd plane is"
335  << n_south + n_north;
336 
337  } else {
338  LogTrace("PreShowerClusterAlgo") << " Wrong plane number, null cluster will be returned!";
339 
340  } // end of if
341 
342  theESNav.home();
343 }
GetRecoTauVFromDQM_MC_cff.next
next
Definition: GetRecoTauVFromDQM_MC_cff.py:31
CaloNavigator::setHome
void setHome(const T &startingPoint)
set the starting position
Definition: CaloNavigator.h:91
CaloNavigator::west
T west() const
move the navigator west
Definition: CaloNavigator.h:49
CaloNavigator::south
T south() const
move the navigator south
Definition: CaloNavigator.h:37
LogTrace
#define LogTrace(id)
Definition: MessageLogger.h:671
PreshowerClusterAlgo::road_2d
std::vector< ESDetId > road_2d
Definition: PreshowerClusterAlgo.h:51
CaloNavigator::east
T east() const
move the navigator east
Definition: CaloNavigator.h:43
CaloNavigator::home
void home() const
move the navigator back to the starting point
Definition: CaloNavigator.h:97
digitizers_cfi.strip
strip
Definition: digitizers_cfi.py:19
CaloNavigator::north
T north() const
move the navigator north
Definition: CaloNavigator.h:31
bool PreshowerClusterAlgo::goodStrip ( RecHitsMap::iterator  candidate_it)

Definition at line 250 of file PreshowerClusterAlgo.cc.

References LogTrace, preshStripEnergyCut_, rechits_map, and used_s.

Referenced by makeOneCluster(), and ~PreshowerClusterAlgo().

250  {
251  if (used_s->find(candidate_it->first) != used_s->end())
252  LogTrace("PreShowerClusterAlgo") << " This strip is in use";
253  if (candidate_it == rechits_map->end())
254  LogTrace("PreShowerClusterAlgo") << " No such a strip in rechits_map";
255  if (candidate_it->second.energy() <= preshStripEnergyCut_)
256  LogTrace("PreShowerClusterAlgo") << "Strip energy" << candidate_it->second.energy() << "is below threshold";
257 
258  // crystal should not be included...
259  if ((used_s->find(candidate_it->first) != used_s->end()) || // ...if it already belongs to a cluster
260  (candidate_it == rechits_map->end()) || // ...if it corresponds to a hit
261  (candidate_it->second.energy() <= preshStripEnergyCut_)) // ...if it has a negative or zero energy
262  {
263  return false;
264  }
265  return true;
266 }
LogTrace
#define LogTrace(id)
Definition: MessageLogger.h:671
reco::PreshowerCluster PreshowerClusterAlgo::makeOneCluster ( ESDetId  strip,
HitsID *  used_strips,
RecHitsMap *  rechits_map,
const CaloSubdetectorGeometry *  geometry_p,
const CaloSubdetectorTopology *  topology_p 
)

Definition at line 13 of file PreshowerClusterAlgo.cc.

References edm::SortedCollection< T, SORT >::begin(), SimDataFormats::CaloAnalysis::cp, CaloNavigator< T, TOPO >::east(), edm::SortedCollection< T, SORT >::end(), CaloParticle::energy(), reco::CaloCluster::energy(), reco::CaloCluster::eta(), findRoad(), newFWLiteAna::found, CaloSubdetectorGeometry::getGeometry(), goodStrip(), CaloNavigator< T, TOPO >::home(), LogTrace, HLT_2018_cff::navigator, GetRecoTauVFromDQM_MC_cff::next, reco::PreshowerCluster::nhits(), CaloNavigator< T, TOPO >::north(), reco::CaloCluster::phi(), ESDetId::plane(), reco::CaloCluster::position(), position, preshClusterEnergyCut_, edm::SortedCollection< T, SORT >::push_back(), rechits_map, road_2d, CaloNavigator< T, TOPO >::setHome(), edm::SortedCollection< T, SORT >::size(), CaloNavigator< T, TOPO >::south(), ESDetId::strip(), used_s, CaloNavigator< T, TOPO >::west(), PV3DBase< T, PVType, FrameType >::x(), reco::CaloCluster::x(), PV3DBase< T, PVType, FrameType >::y(), reco::CaloCluster::y(), PV3DBase< T, PVType, FrameType >::z(), and reco::CaloCluster::z().

Referenced by ~PreshowerClusterAlgo().

17  {
18  road_2d.clear();
19 
20  rechits_map = the_rechitsMap_p;
21 
22  used_s = used_strips;
23 
24  int plane = strip.plane();
25 
26  LogTrace("PreShowerClusterAlgo") << "Preshower Seeded Algorithm - looking for clusters";
27  LogTrace("PreShowerClusterAlgo") << "Preshower is intersected at strip" << strip.strip() << ",at plane" << plane;
28 
29  // create null-cluster
30  std::vector<std::pair<DetId, float> > dummy;
31  Point posi(0, 0, 0);
32  LogTrace("PreShowerClusterAlgo") << " Creating a null-cluster";
33 
34  reco::PreshowerCluster nullcluster = reco::PreshowerCluster(0., posi, dummy, plane);
35 
36  if (strip == ESDetId(0))
37  return nullcluster; //works in case of no intersected strip found (e.g. in the Barrel)
38 
39  // Collection of cluster strips
40  EcalRecHitCollection clusterRecHits;
41  // Map of strips for position calculation
42  RecHitsMap recHits_pos;
43 
44  //Make a navigator, and set it to the strip cell.
45  EcalPreshowerNavigator navigator(strip, topology_p);
46  navigator.setHome(strip);
47  //search for neighbours in the central road
48  findRoad(strip, navigator, plane);
49  LogTrace("PreShowerClusterAlgo") << "Total number of strips in the central road:" << road_2d.size();
50 
51  if (plane == 1) {
52  ESDetId strip_north = navigator.north();
53  findRoad(strip_north, navigator, plane);
54  navigator.home();
55  ESDetId strip_south = navigator.south();
56  findRoad(strip_south, navigator, plane);
57  navigator.home();
58  }
59  if (plane == 2) {
60  ESDetId strip_east = navigator.east();
61  findRoad(strip_east, navigator, plane);
62  navigator.home();
63  ESDetId strip_west = navigator.west();
64  findRoad(strip_west, navigator, plane);
65  navigator.home();
66  }
67  LogTrace("PreShowerClusterAlgo") << "Total number of strips in all three roads:" << road_2d.size();
68 
69  // Start clustering from strip with max Energy in the road
70  float E_max = 0.;
71  bool found = false;
72  RecHitsMap::iterator max_it;
73  // Loop over strips:
74  std::vector<ESDetId>::iterator itID;
75  for (itID = road_2d.begin(); itID != road_2d.end(); itID++) {
76  LogTrace("PreShowerClusterAlgo") << "ID =" << *itID;
77 
78  RecHitsMap::iterator strip_it = rechits_map->find(*itID);
79  //if ( strip_it->second.energy() < 0 ) std::cout << " ##### E = " << strip_it->second.energy() << std::endl;
80  if (strip_it == rechits_map->end() || !goodStrip(strip_it))
81  continue;
82  LogTrace("PreShowerClusterAlgo") << " strip is " << ESDetId(strip_it->first) << "E =" << strip_it->second.energy();
83 
84  float E = strip_it->second.energy();
85  if (E > E_max) {
86  E_max = E;
87  found = true;
88  max_it = strip_it;
89  }
90  }
91  // no hottest strip found ==> null cluster will be returned
92  if (!found)
93  return nullcluster;
94 
95  // First, save the hottest strip
96  clusterRecHits.push_back(max_it->second);
97  recHits_pos.insert(std::make_pair(max_it->first, max_it->second));
98  used_s->insert(max_it->first);
99  LogTrace("PreShowerClusterAlgo") << "Central hottest strip" << ESDetId(max_it->first)
100  << "is saved with energy E =" << E_max;
101 
102  // Find positions of adjacent strips:
103  ESDetId next, strip_1, strip_2;
104  navigator.setHome(max_it->first);
105  ESDetId startES = max_it->first;
106 
107  if (plane == 1) {
108  // Save two neighbouring strips to the east
109  int nadjacents_east = 0;
110  while ((next = navigator.east()) != ESDetId(0) && next != startES && nadjacents_east < 2) {
111  ++nadjacents_east;
112  LogTrace("PreShowerClusterAlgo") << "Adjacent east #" << nadjacents_east << ":" << next;
113 
114  RecHitsMap::iterator strip_it = rechits_map->find(next);
115 
116  if (strip_it == rechits_map->end() || !goodStrip(strip_it))
117  continue;
118  // Save strip for clustering if it exists, not already in use, and satisfies an energy threshold
119  clusterRecHits.push_back(strip_it->second);
120  // save strip for position calculation
121  if (nadjacents_east == 1)
122  strip_1 = next;
123  used_s->insert(strip_it->first);
124  LogTrace("PreShowerClusterAlgo") << "East adjacent strip #" << nadjacents_east
125  << "is saved with energy E =" << strip_it->second.energy();
126  }
127  // Save two neighbouring strips to the west
128  navigator.home();
129  int nadjacents_west = 0;
130  while ((next = navigator.west()) != ESDetId(0) && next != startES && nadjacents_west < 2) {
131  ++nadjacents_west;
132  LogTrace("PreShowerClusterAlgo") << "Adjacent west #" << nadjacents_west << ":" << next;
133 
134  RecHitsMap::iterator strip_it = rechits_map->find(next);
135  if (strip_it == rechits_map->end() || !goodStrip(strip_it))
136  continue;
137  clusterRecHits.push_back(strip_it->second);
138  if (nadjacents_west == 1)
139  strip_2 = next;
140  used_s->insert(strip_it->first);
141  LogTrace("PreShowerClusterAlgo") << "West adjacent strip #" << nadjacents_west
142  << "is saved with energy E =" << strip_it->second.energy();
143  }
144  } else if (plane == 2) {
145  // Save two neighbouring strips to the north
146  int nadjacents_north = 0;
147  while ((next = navigator.north()) != ESDetId(0) && next != startES && nadjacents_north < 2) {
148  ++nadjacents_north;
149  LogTrace("PreShowerClusterAlgo") << "Adjacent north #" << nadjacents_north << ":" << next;
150 
151  RecHitsMap::iterator strip_it = rechits_map->find(next);
152  if (strip_it == rechits_map->end() || !goodStrip(strip_it))
153  continue;
154  clusterRecHits.push_back(strip_it->second);
155  if (nadjacents_north == 1)
156  strip_1 = next;
157  used_s->insert(strip_it->first);
158  LogTrace("PreShowerClusterAlgo") << "North adjacent strip #" << nadjacents_north
159  << "is saved with energy E =" << strip_it->second.energy();
160  }
161  // Save two neighbouring strips to the south
162  navigator.home();
163  int nadjacents_south = 0;
164  while ((next = navigator.south()) != ESDetId(0) && next != startES && nadjacents_south < 2) {
165  ++nadjacents_south;
166  LogTrace("PreShowerClusterAlgo") << "Adjacent south #" << nadjacents_south << ":" << next;
167 
168  RecHitsMap::iterator strip_it = rechits_map->find(next);
169  if (strip_it == rechits_map->end() || !goodStrip(strip_it))
170  continue;
171  clusterRecHits.push_back(strip_it->second);
172  if (nadjacents_south == 1)
173  strip_2 = next;
174  used_s->insert(strip_it->first);
175  LogTrace("PreShowerClusterAlgo") << "South adjacent strip #" << nadjacents_south
176  << "is saved with energy E =" << strip_it->second.energy();
177  }
178  } else {
179  LogTrace("PreShowerClusterAlgo") << " Wrong plane number" << plane << ", null cluster will be returned! "
180  << std::endl;
181  return nullcluster;
182  } // end of if
183 
184  LogTrace("PreShowerClusterAlgo") << "Total size of clusterRecHits is" << clusterRecHits.size();
185  LogTrace("PreShowerClusterAlgo") << "Two adjacent strips for position calculation are:" << strip_1 << "and"
186  << strip_2;
187 
188  // strips for position calculation
189  RecHitsMap::iterator strip_it1, strip_it2;
190  if (strip_1 != ESDetId(0)) {
191  strip_it1 = rechits_map->find(strip_1);
192  recHits_pos.insert(std::make_pair(strip_it1->first, strip_it1->second));
193  }
194  if (strip_2 != ESDetId(0)) {
195  strip_it2 = rechits_map->find(strip_2);
196  recHits_pos.insert(std::make_pair(strip_it2->first, strip_it2->second));
197  }
198 
199  RecHitsMap::iterator cp;
200  double energy_pos = 0;
201  double x_pos = 0;
202  double y_pos = 0;
203  double z_pos = 0;
204  for (cp = recHits_pos.begin(); cp != recHits_pos.end(); cp++) {
205  double E = cp->second.energy();
206  energy_pos += E;
207  auto thisCell = geometry_p->getGeometry(cp->first);
208  const GlobalPoint& position = thisCell->getPosition();
209  x_pos += E * position.x();
210  y_pos += E * position.y();
211  z_pos += E * position.z();
212  }
213  if (energy_pos > 0.) {
214  x_pos /= energy_pos;
215  y_pos /= energy_pos;
216  z_pos /= energy_pos;
217  }
218  Point pos(x_pos, y_pos, z_pos);
219  LogTrace("PreShowerClusterAlgo") << "ES Cluster position ="
220  << "(" << x_pos << "," << y_pos << "," << z_pos << ")";
221 
222  EcalRecHitCollection::iterator it;
223  double Eclust = 0;
224 
225  std::vector<std::pair<DetId, float> > usedHits;
226  for (it = clusterRecHits.begin(); it != clusterRecHits.end(); it++) {
227  Eclust += it->energy();
228  usedHits.push_back(std::pair<DetId, float>(it->id(), 1.));
229  }
230 
231  // ES cluster is created from vector clusterRecHits
232  reco::PreshowerCluster cluster = reco::PreshowerCluster(Eclust, pos, usedHits, plane);
233  LogTrace("PreShowerClusterAlgo") << " ES Cluster is created with:";
234  LogTrace("PreShowerClusterAlgo") << " energy =" << cluster.energy();
235  LogTrace("PreShowerClusterAlgo") << " (eta,phi) ="
236  << "(" << cluster.eta() << "," << cluster.phi() << ")";
237  LogTrace("PreShowerClusterAlgo") << " nhits =" << cluster.nhits();
238  LogTrace("PreShowerClusterAlgo") << " radius =" << cluster.position().r();
239  LogTrace("PreShowerClusterAlgo") << " (x,y,z) ="
240  << "(" << cluster.x() << ", " << cluster.y() << "," << cluster.z() << ")";
241 
242  // return the cluster if its energy is greater a threshold
243  if (cluster.energy() > preshClusterEnergyCut_)
244  return cluster;
245  else
246  return nullcluster;
247 }
reco::CaloCluster::position
const math::XYZPoint & position() const
cluster centroid position
Definition: CaloCluster.h:153
ESDetId::strip
int strip() const
Definition: ESDetId.h:47
PreshowerClusterAlgo::findRoad
void findRoad(ESDetId strip, EcalPreshowerNavigator theESNav, int plane)
Definition: PreshowerClusterAlgo.cc:269
PV3DBase::y
T y() const
Definition: PV3DBase.h:60
edm::SortedCollection::push_back
void push_back(T const &t)
Definition: SortedCollection.h:190
GetRecoTauVFromDQM_MC_cff.next
next
Definition: GetRecoTauVFromDQM_MC_cff.py:31
HLT_2018_cff.navigator
navigator
Definition: HLT_2018_cff.py:11719
reco::CaloCluster::z
double z() const
z coordinate of cluster centroid
Definition: CaloCluster.h:177
reco::CaloCluster::eta
double eta() const
pseudorapidity of cluster centroid
Definition: CaloCluster.h:180
PV3DBase::z
T z() const
Definition: PV3DBase.h:61
CaloParticle::energy
float energy() const
Energy. Note this is taken from the first SimTrack only.
Definition: CaloParticle.h:98
PreshowerClusterAlgo::goodStrip
bool goodStrip(RecHitsMap::iterator candidate_it)
Definition: PreshowerClusterAlgo.cc:250
reco::CaloCluster::x
double x() const
x coordinate of cluster centroid
Definition: CaloCluster.h:171
reco::CaloCluster::energy
double energy() const
cluster energy
Definition: CaloCluster.h:148
reco::PreshowerCluster::nhits
int nhits() const
Number of RecHits the cluster.
Definition: PreshowerCluster.h:35
LogTrace
#define LogTrace(id)
Definition: MessageLogger.h:671
PreshowerClusterAlgo::road_2d
std::vector< ESDetId > road_2d
Definition: PreshowerClusterAlgo.h:51
edm::SortedCollection< EcalRecHit >::iterator
std::vector< EcalRecHit >::iterator iterator
Definition: SortedCollection.h:83
edm::SortedCollection::end
const_iterator end() const
Definition: SortedCollection.h:269
CaloSubdetectorGeometry::getGeometry
virtual std::shared_ptr< const CaloCellGeometry > getGeometry(const DetId &id) const
Get the cell geometry of a given detector id. Should return false if not found.
Definition: CaloSubdetectorGeometry.cc:36
RecHitsMap
std::map< DetId, EcalRecHit > RecHitsMap
Definition: CosmicClusterAlgo.h:25
edm::SortedCollection::size
size_type size() const
Definition: SortedCollection.h:217
Point
Structure Point Contains parameters of Gaussian fits to DMRs.
Definition: DMRtrends.cc:57
position
static int position[264][3]
Definition: ReadPGInfo.cc:289
ESDetId::plane
int plane() const
Definition: ESDetId.h:41
reco::CaloCluster::y
double y() const
y coordinate of cluster centroid
Definition: CaloCluster.h:174
reco::CaloCluster::phi
double phi() const
azimuthal angle of cluster centroid
Definition: CaloCluster.h:183
PV3DBase::x
T x() const
Definition: PV3DBase.h:59
edm::SortedCollection::begin
const_iterator begin() const
Definition: SortedCollection.h:264

Member Data Documentation

double PreshowerClusterAlgo::preshClusterEnergyCut_
private

Definition at line 48 of file PreshowerClusterAlgo.h.

Referenced by makeOneCluster().

int PreshowerClusterAlgo::preshSeededNstr_
private

Definition at line 49 of file PreshowerClusterAlgo.h.

Referenced by findRoad().

double PreshowerClusterAlgo::preshStripEnergyCut_
private

Definition at line 47 of file PreshowerClusterAlgo.h.

Referenced by goodStrip().

RecHitsMap* PreshowerClusterAlgo::rechits_map
private

Definition at line 54 of file PreshowerClusterAlgo.h.

Referenced by goodStrip(), makeOneCluster(), and ~PreshowerClusterAlgo().

std::vector<ESDetId> PreshowerClusterAlgo::road_2d
private

Definition at line 51 of file PreshowerClusterAlgo.h.

Referenced by findRoad(), and makeOneCluster().

HitsID* PreshowerClusterAlgo::used_s
private

Definition at line 57 of file PreshowerClusterAlgo.h.

Referenced by goodStrip(), and makeOneCluster().

Generated for CMSSW Reference Manual by doxygen 1.8.11