CMS 3D CMS Logo
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
PFRecHitProducerHCAL.cc
Go to the documentation of this file.
1 #include "RecoParticleFlow/PFClusterProducer/plugins/PFRecHitProducerHCAL.h"
2 
3 #include <memory>
4 
5 #include "DataFormats/ParticleFlowReco/interface/PFRecHit.h"
6 
7 #include "DataFormats/HcalRecHit/interface/HFRecHit.h"
8 #include "DataFormats/HcalRecHit/interface/HBHERecHit.h"
9 #include "DataFormats/HcalRecHit/interface/HcalRecHitCollections.h"
10 #include "DataFormats/HcalDetId/interface/HcalSubdetector.h"
11 
12 #include "DataFormats/DetId/interface/DetId.h"
13 #include "DataFormats/EcalDetId/interface/EEDetId.h"
14 #include "DataFormats/EcalDetId/interface/EBDetId.h"
15 
16 #include "FWCore/MessageLogger/interface/MessageLogger.h"
17 #include "FWCore/Framework/interface/ESHandle.h"
18 #include "FWCore/Framework/interface/EventSetup.h"
19 
20 #include "Geometry/CaloGeometry/interface/CaloSubdetectorGeometry.h"
21 #include "Geometry/CaloGeometry/interface/CaloGeometry.h"
22 #include "Geometry/CaloGeometry/interface/CaloCellGeometry.h"
23 #include "Geometry/CaloGeometry/interface/TruncatedPyramid.h"
24 #include "Geometry/Records/interface/CaloGeometryRecord.h"
25 
26 #include "Geometry/CaloTopology/interface/HcalTopology.h"
27 #include "RecoCaloTools/Navigation/interface/CaloNavigator.h"
28 
29 #include "DataFormats/CaloTowers/interface/CaloTowerCollection.h"
30 
31 #include "Geometry/CaloTopology/interface/CaloTowerTopology.h"
32 #include "RecoCaloTools/Navigation/interface/CaloTowerNavigator.h"
33 
34 #include "RecoLocalCalo/HcalRecAlgos/interface/HcalCaloFlagLabels.h"
35 #include "RecoLocalCalo/HcalRecAlgos/interface/HcalSeverityLevelComputer.h"
36 #include "RecoLocalCalo/HcalRecAlgos/interface/HcalSeverityLevelComputerRcd.h"
37 
38 using namespace std;
39 using namespace edm;
40 
41 PFRecHitProducerHCAL::PFRecHitProducerHCAL(const edm::ParameterSet& iConfig)
42  : PFRecHitProducer( iConfig )
43 {
44  // access to the collections of rechits
45 
46  inputTagHcalRecHitsHBHE_ =
47  iConfig.getParameter<InputTag>("hcalRecHitsHBHE");
48 
49  inputTagHcalRecHitsHF_ =
50  iConfig.getParameter<InputTag>("hcalRecHitsHF");
51 
52 
53  inputTagCaloTowers_ =
54  iConfig.getParameter<InputTag>("caloTowers");
55 
56  thresh_HF_ =
57  iConfig.getParameter<double>("thresh_HF");
58  navigation_HF_ =
59  iConfig.getParameter<bool>("navigation_HF");
60  weight_HFem_ =
61  iConfig.getParameter<double>("weight_HFem");
62  weight_HFhad_ =
63  iConfig.getParameter<double>("weight_HFhad");
64 
65  HCAL_Calib_ =
66  iConfig.getParameter<bool>("HCAL_Calib");
67  HF_Calib_ =
68  iConfig.getParameter<bool>("HF_Calib");
69  HCAL_Calib_29 =
70  iConfig.getParameter<double>("HCAL_Calib_29");
71  HF_Calib_29 =
72  iConfig.getParameter<double>("HF_Calib_29");
73 
74  shortFibre_Cut = iConfig.getParameter<double>("ShortFibre_Cut");
75  longFibre_Fraction = iConfig.getParameter<double>("LongFibre_Fraction");
76 
77  longFibre_Cut = iConfig.getParameter<double>("LongFibre_Cut");
78  shortFibre_Fraction = iConfig.getParameter<double>("ShortFibre_Fraction");
79 
80  applyLongShortDPG_ = iConfig.getParameter<bool>("ApplyLongShortDPG");
81 
82  longShortFibre_Cut = iConfig.getParameter<double>("LongShortFibre_Cut");
83  minShortTiming_Cut = iConfig.getParameter<double>("MinShortTiming_Cut");
84  maxShortTiming_Cut = iConfig.getParameter<double>("MaxShortTiming_Cut");
85  minLongTiming_Cut = iConfig.getParameter<double>("MinLongTiming_Cut");
86  maxLongTiming_Cut = iConfig.getParameter<double>("MaxLongTiming_Cut");
87 
88  applyTimeDPG_ = iConfig.getParameter<bool>("ApplyTimeDPG");
89  applyPulseDPG_ = iConfig.getParameter<bool>("ApplyPulseDPG");
90  HcalMaxAllowedHFLongShortSev_ = iConfig.getParameter<int>("HcalMaxAllowedHFLongShortSev");
91  HcalMaxAllowedHFDigiTimeSev_ = iConfig.getParameter<int>("HcalMaxAllowedHFDigiTimeSev");
92  HcalMaxAllowedHFInTimeWindowSev_ = iConfig.getParameter<int>("HcalMaxAllowedHFInTimeWindowSev");
93  HcalMaxAllowedChannelStatusSev_ = iConfig.getParameter<int>("HcalMaxAllowedChannelStatusSev");
94 
95  ECAL_Compensate_ = iConfig.getParameter<bool>("ECAL_Compensate");
96  ECAL_Threshold_ = iConfig.getParameter<double>("ECAL_Threshold");
97  ECAL_Compensation_ = iConfig.getParameter<double>("ECAL_Compensation");
98  ECAL_Dead_Code_ = iConfig.getParameter<unsigned int>("ECAL_Dead_Code");
99 
100  EM_Depth_ = iConfig.getParameter<double>("EM_Depth");
101  HAD_Depth_ = iConfig.getParameter<double>("HAD_Depth");
102 
103  //Get integer values of individual HCAL HF flags
104  hcalHFLongShortFlagValue_=1<<HcalCaloFlagLabels::HFLongShort;
105  hcalHFDigiTimeFlagValue_=1<<HcalCaloFlagLabels::HFDigiTime;
106  hcalHFInTimeWindowFlagValue_=1<<HcalCaloFlagLabels::HFInTimeWindow;
107 
108 
109  //--ab
110  produces<reco::PFRecHitCollection>("HFHAD").setBranchAlias("HFHADRecHits");
111  produces<reco::PFRecHitCollection>("HFEM").setBranchAlias("HFEMRecHits");
112  //--ab
113 }
114 
115 
116 
117 PFRecHitProducerHCAL::~PFRecHitProducerHCAL() {}
118 
119 
120 
121 void PFRecHitProducerHCAL::createRecHits(vector<reco::PFRecHit>& rechits,
122  vector<reco::PFRecHit>& rechitsCleaned,
123  edm::Event& iEvent,
124  const edm::EventSetup& iSetup ) {
125 
126 
127  // this map is necessary to find the rechit neighbours efficiently
128  //C but I should think about using Florian's hashed index to do this.
129  //C in which case the map might not be necessary anymore
130  //C however the hashed index does not seem to be implemented for HCAL
131  //
132  // the key of this map is detId.
133  // the value is the index in the rechits vector
134  map<unsigned, unsigned > idSortedRecHits;
135  map<unsigned, unsigned > idSortedRecHitsHFEM;
136  map<unsigned, unsigned > idSortedRecHitsHFHAD;
137  typedef map<unsigned, unsigned >::iterator IDH;
138 
139 
140  edm::ESHandle<CaloGeometry> geoHandle;
141  iSetup.get<CaloGeometryRecord>().get(geoHandle);
142 
143  // get the hcalBarrel geometry
144  const CaloSubdetectorGeometry *hcalBarrelGeometry =
145  geoHandle->getSubdetectorGeometry(DetId::Hcal, HcalBarrel);
146 
147  // get the endcap geometry
148  const CaloSubdetectorGeometry *hcalEndcapGeometry =
149  geoHandle->getSubdetectorGeometry(DetId::Hcal, HcalEndcap);
150 
151  // Get Hcal Severity Level Computer, so that the severity of each rechit flag/status may be determined
152  edm::ESHandle<HcalSeverityLevelComputer> hcalSevLvlComputerHndl;
153  iSetup.get<HcalSeverityLevelComputerRcd>().get(hcalSevLvlComputerHndl);
154  const HcalSeverityLevelComputer* hcalSevLvlComputer = hcalSevLvlComputerHndl.product();
155 
156  //--ab
157  auto_ptr< vector<reco::PFRecHit> > HFHADRecHits( new vector<reco::PFRecHit> );
158  auto_ptr< vector<reco::PFRecHit> > HFEMRecHits( new vector<reco::PFRecHit> );
159  //--ab
160 
161  // 2 possibilities to make HCAL clustering :
162  // - from the HCAL rechits
163  // - from the CaloTowers.
164  // ultimately, clustering will be done taking CaloTowers as an
165  // input. This possibility is currently under investigation, and
166  // was thus made optional.
167 
168  // in the first step, we will fill the map of PFRecHits hcalrechits
169  // either from CaloTowers or from HCAL rechits.
170 
171  // in the second step, we will perform clustering on this map.
172 
173  if( !(inputTagCaloTowers_ == InputTag()) ) {
174 
175  edm::Handle<CaloTowerCollection> caloTowers;
176  CaloTowerTopology caloTowerTopology;
177  const CaloSubdetectorGeometry *caloTowerGeometry = 0;
178  // = geometry_->getSubdetectorGeometry(id)
179 
180  // get calotowers
181  bool found = iEvent.getByLabel(inputTagCaloTowers_,
182  caloTowers);
183 
184  if(!found) {
185  ostringstream err;
186  err<<"could not find rechits "<<inputTagCaloTowers_;
187  LogError("PFRecHitProducerHCAL")<<err.str()<<endl;
188 
189  throw cms::Exception( "MissingProduct", err.str());
190  }
191  else {
192  assert( caloTowers.isValid() );
193 
194  // get HF rechits
195  edm::Handle<HFRecHitCollection> hfHandle;
196  found = iEvent.getByLabel(inputTagHcalRecHitsHF_,
197  hfHandle);
198 
199  if(!found) {
200  ostringstream err;
201  err<<"could not find HF rechits "<<inputTagHcalRecHitsHF_;
202  LogError("PFRecHitProducerHCAL")<<err.str()<<endl;
203 
204  throw cms::Exception( "MissingProduct", err.str());
205  }
206  else {
207  assert( hfHandle.isValid() );
208  }
209 
210  // get HBHE rechits
211  edm::Handle<HBHERecHitCollection> hbheHandle;
212  found = iEvent.getByLabel(inputTagHcalRecHitsHBHE_,
213  hbheHandle);
214 
215  if(!found) {
216  ostringstream err;
217  err<<"could not find HBHE rechits "<<inputTagHcalRecHitsHBHE_;
218  LogError("PFRecHitProducerHCAL")<<err.str()<<endl;
219 
220  throw cms::Exception( "MissingProduct", err.str());
221  }
222  else {
223  assert( hbheHandle.isValid() );
224  }
225 
226  // create rechits
227  typedef CaloTowerCollection::const_iterator ICT;
228 
229  for(ICT ict=caloTowers->begin(); ict!=caloTowers->end();ict++) {
230 
231  const CaloTower& ct = (*ict);
232 
233  //C
234  if(!caloTowerGeometry)
235  caloTowerGeometry = geoHandle->getSubdetectorGeometry(ct.id());
236 
237 
238  // get the hadronic energy.
239 
240  // Mike: Just ask for the Hadronic part only now!
241  // Patrick : ARGH ! While this is ok for the HCAL, this is
242  // just wrong for the HF (in which em/had are artificially
243  // separated.
244  double energy = ct.hadEnergy();
245  //Auguste: Photons in HF have no hadEnergy in fastsim: -> all RecHit collections are empty with photons.
246  double energyEM = ct.emEnergy(); // For HF !
247  //so test the total energy to deal with the photons in HF:
248  if( (energy+energyEM) < 1e-9 ) continue;
249 
250  assert( ct.constituentsSize() );
251  //Mike: The DetId will be taken by the first Hadronic constituent
252  // of the tower. That is only what we need
253 
254 
255  //get the constituents of the tower
256  const std::vector<DetId>& hits = ct.constituents();
257  const std::vector<DetId>& allConstituents = theTowerConstituentsMap->constituentsOf(ct.id());
258 
259  /*
260  for(unsigned int i=0;i< hits.size();++i) {
261  if(hits[i].det()==DetId::Hcal) {
262  HcalDetId did = hits[i];
263  if ( did.subdet()==HcalEndcap || did.subdet()==HcalForward ) {
264  //double en = hits[i].energy();
265  int ieta = did.ieta();
266  const CaloCellGeometry *thisCell = hcalEndcapGeometry->getGeometry(did);
267  const GlobalPoint& position = thisCell->getPosition();
268  if ( abs(ieta) > 27 && abs(ieta) < 33 && energy > 10. ) {
269  std::cout << "HE/HF hit " << i << " at eta = " << ieta
270  << " with CT energy = " << energy
271  << " at eta, z (hit) = " << position.eta() << " " << position.z()
272  << " at eta, z (cte) = " << ct.emPosition().eta() << " " << ct.emPosition().z()
273  << " at eta, z (cth) = " << ct.hadPosition().eta() << " " << ct.hadPosition().z()
274  << " at eta, z (cto) = " << ct.eta() << " " << ct.vz()
275  << std::endl;
276  }
277  }
278  }
279  }
280  */
281 
282  //Reserve the DetId we are looking for:
283 
284  HcalDetId detid;
285  // EcalDetId edetid;
286  bool foundHCALConstituent = false;
287 
288 
289  //Loop on the calotower constituents and search for HCAL
290  double dead = 0.;
291  double alive = 0.;
292  for(unsigned int i=0;i< hits.size();++i) {
293  if(hits[i].det()==DetId::Hcal) {
294  foundHCALConstituent = true;
295  detid = hits[i];
296  // An HCAL tower was found: Look for dead ECAL channels in the same CaloTower.
297  if ( ECAL_Compensate_ && energy > ECAL_Threshold_ ) {
298  for(unsigned int j=0;j<allConstituents.size();++j) {
299  if ( allConstituents[j].det()==DetId::Ecal ) {
300  alive += 1.;
301  EcalChannelStatus::const_iterator chIt = theEcalChStatus->find(allConstituents[j]);
302  unsigned int dbStatus = chIt != theEcalChStatus->end() ? chIt->getStatusCode() : 0;
303  if ( dbStatus > ECAL_Dead_Code_ ) dead += 1.;
304  }
305  }
306  }
307  // Protection: tower 29 in HF is merged with tower 30.
308  // Just take the position of tower 30 in that case.
309  if ( detid.subdet() == HcalForward && abs(detid.ieta()) == 29 ) continue;
310  break;
311  }
312  }
313 
314  // In case of dead ECAL channel, rescale the HCAL energy...
315  double rescaleFactor = alive > 0. ? 1. + ECAL_Compensation_*dead/alive : 1.;
316 
317  reco::PFRecHit* pfrh = 0;
318  reco::PFRecHit* pfrhCleaned = 0;
319  //---ab: need 2 rechits for the HF:
320  reco::PFRecHit* pfrhHFEM = 0;
321  reco::PFRecHit* pfrhHFHAD = 0;
322  reco::PFRecHit* pfrhHFEMCleaned = 0;
323  reco::PFRecHit* pfrhHFHADCleaned = 0;
324  reco::PFRecHit* pfrhHFEMCleaned29 = 0;
325  reco::PFRecHit* pfrhHFHADCleaned29 = 0;
326 
327  if(foundHCALConstituent)
328  {
329  // std::cout << ", new Energy = " << energy << std::endl;
330  switch( detid.subdet() ) {
331  case HcalBarrel:
332  {
333  if(energy < thresh_Barrel_ ) continue;
334 
335  /*
336  // Check the timing
337  if ( energy > 5. ) {
338  for(unsigned int i=0;i< hits.size();++i) {
339  if( hits[i].det() != DetId::Hcal ) continue;
340  HcalDetId theDetId = hits[i];
341  typedef HBHERecHitCollection::const_iterator iHBHE;
342  iHBHE theHit = hbheHandle->find(theDetId);
343  if ( theHit != hbheHandle->end() )
344  std::cout << "HCAL hit : "
345  << theDetId.ieta() << " " << theDetId.iphi() << " "
346  << theHit->energy() << " " << theHit->time() << std::endl;
347  }
348  }
349  */
350 
351 
352  // if ( HCAL_Calib_ ) energy *= std::min(max_Calib_,myPFCorr->getValues(detid)->getValue());
353  //if ( rescaleFactor > 1. )
354  // std::cout << "Barrel HCAL energy rescaled from = " << energy << " to " << energy*rescaleFactor << std::endl;
355  if ( rescaleFactor > 1. ) {
356  pfrhCleaned = createHcalRecHit( detid,
357  energy,
358  PFLayer::HCAL_BARREL1,
359  hcalBarrelGeometry,
360  ct.id().rawId() );
361  pfrhCleaned->setRescale(rescaleFactor);
362  energy *= rescaleFactor;
363  }
364  pfrh = createHcalRecHit( detid,
365  energy,
366  PFLayer::HCAL_BARREL1,
367  hcalBarrelGeometry,
368  ct.id().rawId() );
369  pfrh->setRescale(rescaleFactor);
370  }
371  break;
372  case HcalEndcap:
373  {
374  if(energy < thresh_Endcap_ ) continue;
375 
376  /*
377  // Check the timing
378  if ( energy > 5. ) {
379  for(unsigned int i=0;i< hits.size();++i) {
380  if( hits[i].det() != DetId::Hcal ) continue;
381  HcalDetId theDetId = hits[i];
382  typedef HBHERecHitCollection::const_iterator iHBHE;
383  iHBHE theHit = hbheHandle->find(theDetId);
384  if ( theHit != hbheHandle->end() )
385  std::cout << "HCAL hit : "
386  << theDetId.ieta() << " " << theDetId.iphi() << " "
387  << theHit->energy() << " " << theHit->time() << std::endl;
388  }
389  }
390  */
391 
392  // Apply tower 29 calibration
393  if ( HCAL_Calib_ && abs(detid.ieta()) == 29 ) energy *= HCAL_Calib_29;
394  //if ( rescaleFactor > 1. )
395  // std::cout << "End-cap HCAL energy rescaled from = " << energy << " to " << energy*rescaleFactor << std::endl;
396  if ( rescaleFactor > 1. ) {
397  pfrhCleaned = createHcalRecHit( detid,
398  energy,
399  PFLayer::HCAL_ENDCAP,
400  hcalEndcapGeometry,
401  ct.id().rawId() );
402  pfrhCleaned->setRescale(rescaleFactor);
403  energy *= rescaleFactor;
404  }
405  pfrh = createHcalRecHit( detid,
406  energy,
407  PFLayer::HCAL_ENDCAP,
408  hcalEndcapGeometry,
409  ct.id().rawId() );
410  pfrh->setRescale(rescaleFactor);
411  }
412  break;
413  case HcalOuter:
414  {
415  }
416  break;
417  case HcalForward:
418  {
419  //---ab: 2 rechits for HF:
420  //double energyemHF = weight_HFem_*ct.emEnergy();
421  //double energyhadHF = weight_HFhad_*ct.hadEnergy();
422  double energyemHF = weight_HFem_ * energyEM;
423  double energyhadHF = weight_HFhad_ * energy;
424  // Some energy in the tower !
425  if((energyemHF+energyhadHF) < thresh_HF_ ) continue;
426 
427  // Some cleaning in the HF
428  double longFibre = energyemHF + energyhadHF/2.;
429  double shortFibre = energyhadHF/2.;
430  int ieta = detid.ieta();
431  int iphi = detid.iphi();
432  HcalDetId theLongDetId (HcalForward, ieta, iphi, 1);
433  HcalDetId theShortDetId (HcalForward, ieta, iphi, 2);
434  typedef HFRecHitCollection::const_iterator iHF;
435  iHF theLongHit = hfHandle->find(theLongDetId);
436  iHF theShortHit = hfHandle->find(theShortDetId);
437  //
438  double theLongHitEnergy = 0.;
439  double theShortHitEnergy = 0.;
440  bool flagShortDPG = false;
441  bool flagLongDPG = false;
442  bool flagShortTimeDPG = false;
443  bool flagLongTimeDPG = false;
444  bool flagShortPulseDPG = false;
445  bool flagLongPulseDPG = false;
446  //
447  if ( theLongHit != hfHandle->end() ) {
448  int theLongFlag = theLongHit->flags();
449  theLongHitEnergy = theLongHit->energy();
450  flagLongDPG = applyLongShortDPG_ && ( hcalSevLvlComputer->getSeverityLevel(theLongDetId, theLongFlag & hcalHFLongShortFlagValue_, 0)> HcalMaxAllowedHFLongShortSev_);
451  flagLongTimeDPG = applyTimeDPG_ && ( hcalSevLvlComputer->getSeverityLevel(theLongDetId, theLongFlag & hcalHFInTimeWindowFlagValue_, 0)> HcalMaxAllowedHFInTimeWindowSev_);
452  flagLongPulseDPG = applyPulseDPG_ && ( hcalSevLvlComputer->getSeverityLevel(theLongDetId, theLongFlag & hcalHFDigiTimeFlagValue_, 0)> HcalMaxAllowedHFDigiTimeSev_);
453 
454  //flagLongDPG = applyLongShortDPG_ && theLongHit->flagField(HcalCaloFlagLabels::HFLongShort)==1;
455  //flagLongTimeDPG = applyTimeDPG_ && theLongHit->flagField(HcalCaloFlagLabels::HFInTimeWindow)==1;
456  //flagLongPulseDPG = applyPulseDPG_ && theLongHit->flagField(HcalCaloFlagLabels::HFDigiTime)==1;
457  }
458  //
459  if ( theShortHit != hfHandle->end() ) {
460  int theShortFlag = theShortHit->flags();
461  theShortHitEnergy = theShortHit->energy();
462  flagShortDPG = applyLongShortDPG_ && ( hcalSevLvlComputer->getSeverityLevel(theShortDetId, theShortFlag & hcalHFLongShortFlagValue_, 0)> HcalMaxAllowedHFLongShortSev_);
463  flagShortTimeDPG = applyTimeDPG_ && ( hcalSevLvlComputer->getSeverityLevel(theShortDetId, theShortFlag & hcalHFInTimeWindowFlagValue_, 0)> HcalMaxAllowedHFInTimeWindowSev_);
464  flagShortPulseDPG = applyPulseDPG_ && ( hcalSevLvlComputer->getSeverityLevel(theShortDetId, theShortFlag & hcalHFDigiTimeFlagValue_, 0)> HcalMaxAllowedHFDigiTimeSev_);
465  //flagShortDPG = applyLongShortDPG_ && theShortHit->flagField(HcalCaloFlagLabels::HFLongShort)==1;
466  //flagShortTimeDPG = applyTimeDPG_ && theShortHit->flagField(HcalCaloFlagLabels::HFInTimeWindow)==1;
467  //flagShortPulseDPG = applyPulseDPG_ && theShortHit->flagField(HcalCaloFlagLabels::HFDigiTime)==1;
468  }
469 
470  // Then check the timing in short and long fibres in all other towers.
471  if ( theShortHitEnergy > longShortFibre_Cut &&
472  ( theShortHit->time() < minShortTiming_Cut ||
473  theShortHit->time() > maxShortTiming_Cut ||
474  flagShortTimeDPG || flagShortPulseDPG ) ) {
475  // rescaleFactor = 0. ;
476  pfrhHFHADCleaned = createHcalRecHit( detid,
477  theShortHitEnergy,
478  PFLayer::HF_HAD,
479  hcalEndcapGeometry,
480  ct.id().rawId() );
481  pfrhHFHADCleaned->setRescale(theShortHit->time());
482  /*
483  std::cout << "ieta/iphi = " << ieta << " " << iphi
484  << ", Energy em/had/long/short = "
485  << energyemHF << " " << energyhadHF << " "
486  << theLongHitEnergy << " " << theShortHitEnergy << " "
487  << ". Time = " << theShortHit->time() << " "
488  << ". The short and long flags : "
489  << flagShortDPG << " " << flagLongDPG << " "
490  << flagShortTimeDPG << " " << flagLongTimeDPG << " "
491  << flagShortPulseDPG << " " << flagLongPulseDPG << " "
492  << ". Short fibres were cleaned." << std::endl;
493  */
494  shortFibre -= theShortHitEnergy;
495  theShortHitEnergy = 0.;
496  }
497 
498 
499  if ( theLongHitEnergy > longShortFibre_Cut &&
500  ( theLongHit->time() < minLongTiming_Cut ||
501  theLongHit->time() > maxLongTiming_Cut ||
502  flagLongTimeDPG || flagLongPulseDPG ) ) {
503  //rescaleFactor = 0. ;
504  pfrhHFEMCleaned = createHcalRecHit( detid,
505  theLongHitEnergy,
506  PFLayer::HF_EM,
507  hcalEndcapGeometry,
508  ct.id().rawId() );
509  pfrhHFEMCleaned->setRescale(theLongHit->time());
510  /*
511  std::cout << "ieta/iphi = " << ieta << " " << iphi
512  << ", Energy em/had/long/short = "
513  << energyemHF << " " << energyhadHF << " "
514  << theLongHitEnergy << " " << theShortHitEnergy << " "
515  << ". Time = " << theLongHit->time() << " "
516  << ". The short and long flags : "
517  << flagShortDPG << " " << flagLongDPG << " "
518  << flagShortTimeDPG << " " << flagLongTimeDPG << " "
519  << flagShortPulseDPG << " " << flagLongPulseDPG << " "
520  << ". Long fibres were cleaned." << std::endl;
521  */
522  longFibre -= theLongHitEnergy;
523  theLongHitEnergy = 0.;
524  }
525 
526  // Some energy must be in the long fibres is there is some energy in the short fibres !
527  if ( theShortHitEnergy > shortFibre_Cut &&
528  ( theLongHitEnergy/theShortHitEnergy < longFibre_Fraction ||
529  flagShortDPG ) ) {
530  // Check if the long-fibre hit was not cleaned already (because hot)
531  // In this case don't apply the cleaning
532  const HcalChannelStatus* theStatus = theHcalChStatus->getValues(theLongDetId);
533  unsigned theStatusValue = theStatus->getValue();
534  int theSeverityLevel = hcalSevLvlComputer->getSeverityLevel(detid, 0, theStatusValue);
535  // The channel is killed
537  if (theSeverityLevel<=HcalMaxAllowedChannelStatusSev_) {
538  // rescaleFactor = 0. ;
539  pfrhHFHADCleaned = createHcalRecHit( detid,
540  theShortHitEnergy,
541  PFLayer::HF_HAD,
542  hcalEndcapGeometry,
543  ct.id().rawId() );
544  pfrhHFHADCleaned->setRescale(theShortHit->time());
545  /*
546  std::cout << "ieta/iphi = " << ieta << " " << iphi
547  << ", Energy em/had/long/short = "
548  << energyemHF << " " << energyhadHF << " "
549  << theLongHitEnergy << " " << theShortHitEnergy << " "
550  << ". Time = " << theShortHit->time() << " "
551  << ". The status value is " << theStatusValue
552  << ". The short and long flags : "
553  << flagShortDPG << " " << flagLongDPG << " "
554  << flagShortTimeDPG << " " << flagLongTimeDPG << " "
555  << flagShortPulseDPG << " " << flagLongPulseDPG << " "
556  << ". Short fibres were cleaned." << std::endl;
557  */
558  shortFibre -= theShortHitEnergy;
559  theShortHitEnergy = 0.;
560  }
561  }
562 
563  if ( theLongHitEnergy > longFibre_Cut &&
564  ( theShortHitEnergy/theLongHitEnergy < shortFibre_Fraction ||
565  flagLongDPG ) ) {
566  // Check if the long-fibre hit was not cleaned already (because hot)
567  // In this case don't apply the cleaning
568  const HcalChannelStatus* theStatus = theHcalChStatus->getValues(theShortDetId);
569  unsigned theStatusValue = theStatus->getValue();
570 
571  int theSeverityLevel = hcalSevLvlComputer->getSeverityLevel(detid, 0, theStatusValue);
572  // The channel is killed
574  if (theSeverityLevel<=HcalMaxAllowedChannelStatusSev_) {
575 
576  //rescaleFactor = 0. ;
577  pfrhHFEMCleaned = createHcalRecHit( detid,
578  theLongHitEnergy,
579  PFLayer::HF_EM,
580  hcalEndcapGeometry,
581  ct.id().rawId() );
582  pfrhHFEMCleaned->setRescale(theLongHit->time());
583  /*
584  std::cout << "ieta/iphi = " << ieta << " " << iphi
585  << ", Energy em/had/long/short = "
586  << energyemHF << " " << energyhadHF << " "
587  << theLongHitEnergy << " " << theShortHitEnergy << " "
588  << ". The status value is " << theStatusValue
589  << ". Time = " << theLongHit->time() << " "
590  << ". The short and long flags : "
591  << flagShortDPG << " " << flagLongDPG << " "
592  << flagShortTimeDPG << " " << flagLongTimeDPG << " "
593  << flagShortPulseDPG << " " << flagLongPulseDPG << " "
594  << ". Long fibres were cleaned." << std::endl;
595  */
596  longFibre -= theLongHitEnergy;
597  theLongHitEnergy = 0.;
598  }
599  }
600 
601  // Special treatment for tower 29
602  // A tower with energy only at ieta = +/- 29 is not physical -> Clean
603  if ( abs(ieta) == 29 ) {
604  // rescaleFactor = 0. ;
605  // Clean long fibres
606  if ( theLongHitEnergy > shortFibre_Cut/2. ) {
607  pfrhHFEMCleaned29 = createHcalRecHit( detid,
608  theLongHitEnergy,
609  PFLayer::HF_EM,
610  hcalEndcapGeometry,
611  ct.id().rawId() );
612  pfrhHFEMCleaned29->setRescale(theLongHit->time());
613  /*
614  std::cout << "ieta/iphi = " << ieta << " " << iphi
615  << ", Energy em/had/long/short = "
616  << energyemHF << " " << energyhadHF << " "
617  << theLongHitEnergy << " " << theShortHitEnergy << " "
618  << ". Long fibres were cleaned." << std::endl;
619  */
620  longFibre -= theLongHitEnergy;
621  theLongHitEnergy = 0.;
622  }
623  // Clean short fibres
624  if ( theShortHitEnergy > shortFibre_Cut/2. ) {
625  pfrhHFHADCleaned29 = createHcalRecHit( detid,
626  theShortHitEnergy,
627  PFLayer::HF_HAD,
628  hcalEndcapGeometry,
629  ct.id().rawId() );
630  pfrhHFHADCleaned29->setRescale(theShortHit->time());
631  /*
632  std::cout << "ieta/iphi = " << ieta << " " << iphi
633  << ", Energy em/had/long/short = "
634  << energyemHF << " " << energyhadHF << " "
635  << theLongHitEnergy << " " << theShortHitEnergy << " "
636  << ". Short fibres were cleaned." << std::endl;
637  */
638  shortFibre -= theShortHitEnergy;
639  theShortHitEnergy = 0.;
640  }
641  }
642  // Check the timing of the long and short fibre rechits
643 
644  // First, check the timing of long and short fibre in eta = 29 if tower 30.
645  else if ( abs(ieta) == 30 ) {
646  int ieta29 = ieta > 0 ? 29 : -29;
647  HcalDetId theLongDetId29 (HcalForward, ieta29, iphi, 1);
648  HcalDetId theShortDetId29 (HcalForward, ieta29, iphi, 2);
649  iHF theLongHit29 = hfHandle->find(theLongDetId29);
650  iHF theShortHit29 = hfHandle->find(theShortDetId29);
651  //
652  double theLongHitEnergy29 = 0.;
653  double theShortHitEnergy29 = 0.;
654  bool flagShortDPG29 = false;
655  bool flagLongDPG29 = false;
656  bool flagShortTimeDPG29 = false;
657  bool flagLongTimeDPG29 = false;
658  bool flagShortPulseDPG29 = false;
659  bool flagLongPulseDPG29 = false;
660  //
661  if ( theLongHit29 != hfHandle->end() ) {
662  int theLongFlag29 = theLongHit29->flags();
663  theLongHitEnergy29 = theLongHit29->energy() ;
664  flagLongDPG29 = applyLongShortDPG_ && ( hcalSevLvlComputer->getSeverityLevel(theLongDetId29, theLongFlag29 & hcalHFLongShortFlagValue_, 0)> HcalMaxAllowedHFLongShortSev_);
665  flagLongTimeDPG29 = applyTimeDPG_ && ( hcalSevLvlComputer->getSeverityLevel(theLongDetId29, theLongFlag29 & hcalHFInTimeWindowFlagValue_, 0)> HcalMaxAllowedHFInTimeWindowSev_);
666  flagLongPulseDPG29 = applyPulseDPG_ && ( hcalSevLvlComputer->getSeverityLevel(theLongDetId29, theLongFlag29 & hcalHFDigiTimeFlagValue_, 0)> HcalMaxAllowedHFDigiTimeSev_);
667 
668  //flagLongDPG29 = applyLongShortDPG_ && theLongHit29->flagField(HcalCaloFlagLabels::HFLongShort)==1;
669  //flagLongTimeDPG29 = applyTimeDPG_ && theLongHit29->flagField(HcalCaloFlagLabels::HFInTimeWindow)==1;
670  //flagLongPulseDPG29 = applyPulseDPG_ && theLongHit29->flagField(HcalCaloFlagLabels::HFDigiTime)==1;
671  }
672  //
673  if ( theShortHit29 != hfHandle->end() ) {
674  int theShortFlag29 = theShortHit29->flags();
675  theShortHitEnergy29 = theShortHit29->energy();
676  flagShortDPG29 = applyLongShortDPG_ && ( hcalSevLvlComputer->getSeverityLevel(theShortDetId29, theShortFlag29 & hcalHFLongShortFlagValue_, 0)> HcalMaxAllowedHFLongShortSev_);
677  flagShortTimeDPG29 = applyTimeDPG_ && ( hcalSevLvlComputer->getSeverityLevel(theShortDetId29, theShortFlag29 & hcalHFInTimeWindowFlagValue_, 0)> HcalMaxAllowedHFInTimeWindowSev_);
678  flagLongPulseDPG29 = applyPulseDPG_ && ( hcalSevLvlComputer->getSeverityLevel(theShortDetId29, theShortFlag29 & hcalHFDigiTimeFlagValue_, 0)> HcalMaxAllowedHFDigiTimeSev_);
679 
680  //flagShortDPG29 = applyLongShortDPG_ && theShortHit29->flagField(HcalCaloFlagLabels::HFLongShort)==1;
681  //flagShortTimeDPG29 = applyTimeDPG_ && theShortHit29->flagField(HcalCaloFlagLabels::HFInTimeWindow)==1;
682  //flagShortPulseDPG29 = applyPulseDPG_ && theShortHit29->flagField(HcalCaloFlagLabels::HFDigiTime)==1;
683  }
684 
685  if ( theLongHitEnergy29 > longShortFibre_Cut &&
686  ( theLongHit29->time() < minLongTiming_Cut ||
687  theLongHit29->time() > maxLongTiming_Cut ||
688  flagLongTimeDPG29 || flagLongPulseDPG29 ) ) {
689  //rescaleFactor = 0. ;
690  pfrhHFEMCleaned29 = createHcalRecHit( detid,
691  theLongHitEnergy29,
692  PFLayer::HF_EM,
693  hcalEndcapGeometry,
694  ct.id().rawId() );
695  pfrhHFEMCleaned29->setRescale(theLongHit29->time());
696  /*
697  std::cout << "ieta/iphi = " << ieta29 << " " << iphi
698  << ", Energy em/had/long/short = "
699  << energyemHF << " " << energyhadHF << " "
700  << theLongHitEnergy29 << " " << theShortHitEnergy29 << " "
701  << ". Time = " << theLongHit29->time() << " "
702  << ". The short and long flags : "
703  << flagShortDPG29 << " " << flagLongDPG29 << " "
704  << flagShortTimeDPG29 << " " << flagLongTimeDPG29 << " "
705  << flagShortPulseDPG29 << " " << flagLongPulseDPG29 << " "
706  << ". Long fibres were cleaned." << std::endl;
707  */
708  longFibre -= theLongHitEnergy29;
709  theLongHitEnergy29 = 0;
710  }
711 
712  if ( theShortHitEnergy29 > longShortFibre_Cut &&
713  ( theShortHit29->time() < minShortTiming_Cut ||
714  theShortHit29->time() > maxShortTiming_Cut ||
715  flagShortTimeDPG29 || flagShortPulseDPG29 ) ) {
716  //rescaleFactor = 0. ;
717  pfrhHFHADCleaned29 = createHcalRecHit( detid,
718  theShortHitEnergy29,
719  PFLayer::HF_HAD,
720  hcalEndcapGeometry,
721  ct.id().rawId() );
722  pfrhHFHADCleaned29->setRescale(theShortHit29->time());
723  /*
724  std::cout << "ieta/iphi = " << ieta29 << " " << iphi
725  << ", Energy em/had/long/short = "
726  << energyemHF << " " << energyhadHF << " "
727  << theLongHitEnergy29 << " " << theShortHitEnergy29 << " "
728  << ". Time = " << theShortHit29->time() << " "
729  << ". The short and long flags : "
730  << flagShortDPG29 << " " << flagLongDPG29 << " "
731  << flagShortTimeDPG29 << " " << flagLongTimeDPG29 << " "
732  << flagShortPulseDPG29 << " " << flagLongPulseDPG29 << " "
733  << ". Short fibres were cleaned." << std::endl;
734  */
735  shortFibre -= theShortHitEnergy29;
736  theShortHitEnergy29 = 0.;
737  }
738 
739  // Some energy must be in the long fibres is there is some energy in the short fibres !
740  if ( theShortHitEnergy29 > shortFibre_Cut &&
741  ( theLongHitEnergy29/theShortHitEnergy29 < 2.*longFibre_Fraction ||
742  flagShortDPG29 ) ) {
743  // Check if the long-fibre hit was not cleaned already (because hot)
744  // In this case don't apply the cleaning
745  const HcalChannelStatus* theStatus = theHcalChStatus->getValues(theLongDetId29);
746  unsigned theStatusValue = theStatus->getValue();
747 
748  int theSeverityLevel = hcalSevLvlComputer->getSeverityLevel(detid, 0, theStatusValue);
749  // The channel is killed
751  if (theSeverityLevel<=HcalMaxAllowedChannelStatusSev_) {
752  //rescaleFactor = 0. ;
753  pfrhHFHADCleaned29 = createHcalRecHit( detid,
754  theShortHitEnergy29,
755  PFLayer::HF_HAD,
756  hcalEndcapGeometry,
757  ct.id().rawId() );
758  pfrhHFHADCleaned29->setRescale(theShortHit29->time());
759  /*
760  std::cout << "ieta/iphi = " << ieta29 << " " << iphi
761  << ", Energy em/had/long/short = "
762  << energyemHF << " " << energyhadHF << " "
763  << theLongHitEnergy29 << " " << theShortHitEnergy29 << " "
764  << ". Time = " << theShortHit29->time() << " "
765  << ". The status value is " << theStatusValue
766  << ". The short and long flags : "
767  << flagShortDPG29 << " " << flagLongDPG29 << " "
768  << flagShortTimeDPG29 << " " << flagLongTimeDPG29 << " "
769  << flagShortPulseDPG29 << " " << flagLongPulseDPG29 << " "
770  << ". Short fibres were cleaned." << std::endl;
771  */
772  shortFibre -= theShortHitEnergy29;
773  theShortHitEnergy29 = 0.;
774  }
775  }
776 
777  // Some energy must be in the short fibres is there is some energy in the long fibres !
778  if ( theLongHitEnergy29 > longFibre_Cut &&
779  ( theShortHitEnergy29/theLongHitEnergy29 < shortFibre_Fraction ||
780  flagLongDPG29 ) ) {
781  // Check if the long-fibre hit was not cleaned already (because hot)
782  // In this case don't apply the cleaning
783  const HcalChannelStatus* theStatus = theHcalChStatus->getValues(theShortDetId29);
784  unsigned theStatusValue = theStatus->getValue();
785  int theSeverityLevel = hcalSevLvlComputer->getSeverityLevel(detid, 0, theStatusValue);
786  // The channel is killed
788  if (theSeverityLevel<=HcalMaxAllowedChannelStatusSev_) {
789 
790  //rescaleFactor = 0. ;
791  pfrhHFEMCleaned29 = createHcalRecHit( detid,
792  theLongHitEnergy29,
793  PFLayer::HF_EM,
794  hcalEndcapGeometry,
795  ct.id().rawId() );
796  pfrhHFEMCleaned29->setRescale(theLongHit29->time());
797  /*
798  std::cout << "ieta/iphi = " << ieta29 << " " << iphi
799  << ", Energy em/had/long/short = "
800  << energyemHF << " " << energyhadHF << " "
801  << theLongHitEnergy29 << " " << theShortHitEnergy29 << " "
802  << ". The status value is " << theStatusValue
803  << ". Time = " << theLongHit29->time() << " "
804  << ". The short and long flags : "
805  << flagShortDPG29 << " " << flagLongDPG29 << " "
806  << flagShortTimeDPG29 << " " << flagLongTimeDPG29 << " "
807  << flagShortPulseDPG29 << " " << flagLongPulseDPG29 << " "
808  << ". Long fibres were cleaned." << std::endl;
809  */
810  longFibre -= theLongHitEnergy29;
811  theLongHitEnergy29 = 0.;
812  }
813  }
814 
815  // Check that the energy in tower 29 is smaller than in tower 30
816  // First in long fibres
817  if ( theLongHitEnergy29 > std::max(theLongHitEnergy,shortFibre_Cut/2) ) {
818  pfrhHFEMCleaned29 = createHcalRecHit( detid,
819  theLongHitEnergy29,
820  PFLayer::HF_EM,
821  hcalEndcapGeometry,
822  ct.id().rawId() );
823  pfrhHFEMCleaned29->setRescale(theLongHit29->time());
824  /*
825  std::cout << "ieta/iphi = " << ieta29 << " " << iphi
826  << ", Energy L29/S29/L30/S30 = "
827  << theLongHitEnergy29 << " " << theShortHitEnergy29 << " "
828  << theLongHitEnergy << " " << theShortHitEnergy << " "
829  << ". Long fibres were cleaned." << std::endl;
830  */
831  longFibre -= theLongHitEnergy29;
832  theLongHitEnergy29 = 0.;
833  }
834  // Second in short fibres
835  if ( theShortHitEnergy29 > std::max(theShortHitEnergy,shortFibre_Cut/2.) ) {
836  pfrhHFHADCleaned29 = createHcalRecHit( detid,
837  theShortHitEnergy29,
838  PFLayer::HF_HAD,
839  hcalEndcapGeometry,
840  ct.id().rawId() );
841  pfrhHFHADCleaned29->setRescale(theShortHit29->time());
842  /*
843  std::cout << "ieta/iphi = " << ieta << " " << iphi
844  << ", Energy L29/S29/L30/S30 = "
845  << theLongHitEnergy29 << " " << theShortHitEnergy29 << " "
846  << theLongHitEnergy << " " << theShortHitEnergy << " "
847  << ". Short fibres were cleaned." << std::endl;
848  */
849  shortFibre -= theShortHitEnergy29;
850  theShortHitEnergy29 = 0.;
851  }
852  }
853 
854 
855  // Determine EM and HAD after cleaning of short and long fibres
856  energyhadHF = 2.*shortFibre;
857  energyemHF = longFibre - shortFibre;
858 
859  // The EM energy might be negative, as it amounts to Long - Short
860  // In that case, put the EM "energy" in the HAD energy
861  // Just to avoid systematic positive bias due to "Short" high fluctuations
862  if ( energyemHF < thresh_HF_ ) {
863  energyhadHF += energyemHF;
864  energyemHF = 0.;
865  }
866 
867  // Apply HCAL calibration factors flor towers close to 29, if requested
868  if ( HF_Calib_ && abs(detid.ieta()) <= 32 ) {
869  energyhadHF *= HF_Calib_29;
870  energyemHF *= HF_Calib_29;
871  }
872 
873  // Create an EM and a HAD rechit if above threshold.
874  if ( energyemHF > thresh_HF_ || energyhadHF > thresh_HF_ ) {
875  pfrhHFEM = createHcalRecHit( detid,
876  energyemHF,
877  PFLayer::HF_EM,
878  hcalEndcapGeometry,
879  ct.id().rawId() );
880  pfrhHFHAD = createHcalRecHit( detid,
881  energyhadHF,
882  PFLayer::HF_HAD,
883  hcalEndcapGeometry,
884  ct.id().rawId() );
885  pfrhHFEM->setEnergyUp(energyhadHF);
886  pfrhHFHAD->setEnergyUp(energyemHF);
887  }
888 
889  }
890  break;
891  default:
892  LogError("PFRecHitProducerHCAL")
893  <<"CaloTower constituent: unknown layer : "
894  <<detid.subdet()<<endl;
895  }
896 
897  if(pfrh) {
898  rechits.push_back( *pfrh );
899  delete pfrh;
900  idSortedRecHits.insert( make_pair(ct.id().rawId(),
901  rechits.size()-1 ) );
902  }
903  if(pfrhCleaned) {
904  rechitsCleaned.push_back( *pfrhCleaned );
905  delete pfrhCleaned;
906  }
907  //---ab: 2 rechits for HF:
908  if(pfrhHFEM) {
909  HFEMRecHits->push_back( *pfrhHFEM );
910  delete pfrhHFEM;
911  idSortedRecHitsHFEM.insert( make_pair(ct.id().rawId(),
912  HFEMRecHits->size()-1 ) );
913  }
914  if(pfrhHFHAD) {
915  HFHADRecHits->push_back( *pfrhHFHAD );
916  delete pfrhHFHAD;
917  idSortedRecHitsHFHAD.insert( make_pair(ct.id().rawId(),
918  HFHADRecHits->size()-1 ) );
919  }
920  //---ab
921  if(pfrhHFEMCleaned) {
922  rechitsCleaned.push_back( *pfrhHFEMCleaned );
923  delete pfrhHFEMCleaned;
924  }
925  if(pfrhHFHADCleaned) {
926  rechitsCleaned.push_back( *pfrhHFHADCleaned );
927  delete pfrhHFHADCleaned;
928  }
929  if(pfrhHFEMCleaned29) {
930  rechitsCleaned.push_back( *pfrhHFEMCleaned29 );
931  delete pfrhHFEMCleaned29;
932  }
933  if(pfrhHFHADCleaned29) {
934  rechitsCleaned.push_back( *pfrhHFHADCleaned29 );
935  delete pfrhHFHADCleaned29;
936  }
937  }
938  }
939  // do navigation
940  for(unsigned i=0; i<rechits.size(); i++ ) {
941  findRecHitNeighboursCT( rechits[i],
942  idSortedRecHits,
943  caloTowerTopology);
944  }
945  for(unsigned i=0; i<HFEMRecHits->size(); i++ ) {
946  findRecHitNeighboursCT( (*HFEMRecHits)[i],
947  idSortedRecHitsHFEM,
948  caloTowerTopology);
949  }
950  for(unsigned i=0; i<HFHADRecHits->size(); i++ ) {
951  findRecHitNeighboursCT( (*HFHADRecHits)[i],
952  idSortedRecHitsHFHAD,
953  caloTowerTopology);
954  }
955  iEvent.put( HFHADRecHits,"HFHAD" );
956  iEvent.put( HFEMRecHits,"HFEM" );
957  }
958  }
959  else if( !(inputTagHcalRecHitsHBHE_ == InputTag()) ) {
960  // clustering is not done on CaloTowers but on HCAL rechits.
961 
962 
963  // get the hcal topology
964  edm::ESHandle<HcalTopology> hcalTopology;
965  iSetup.get<IdealGeometryRecord>().get( hcalTopology );
966 
967  // HCAL rechits
968  // vector<edm::Handle<HBHERecHitCollection> > hcalHandles;
969  edm::Handle<HBHERecHitCollection> hcalHandle;
970 
971 
972  bool found = iEvent.getByLabel(inputTagHcalRecHitsHBHE_,
973  hcalHandle );
974 
975  if(!found) {
976  ostringstream err;
977  err<<"could not find rechits "<<inputTagHcalRecHitsHBHE_;
978  LogError("PFRecHitProducerHCAL")<<err.str()<<endl;
979 
980  throw cms::Exception( "MissingProduct", err.str());
981  }
982  else {
983  assert( hcalHandle.isValid() );
984 
985  const edm::Handle<HBHERecHitCollection>& handle = hcalHandle;
986  for(unsigned irechit=0; irechit<handle->size(); irechit++) {
987  const HBHERecHit& hit = (*handle)[irechit];
988 
989  double energy = hit.energy();
990 
991  reco::PFRecHit* pfrh = 0;
992 
993 
994  const HcalDetId& detid = hit.detid();
995  switch( detid.subdet() ) {
996  case HcalBarrel:
997  {
998  if(energy < thresh_Barrel_ ) continue;
999  pfrh = createHcalRecHit( detid,
1000  energy,
1001  PFLayer::HCAL_BARREL1,
1002  hcalBarrelGeometry );
1003  }
1004  break;
1005  case HcalEndcap:
1006  {
1007  if(energy < thresh_Endcap_ ) continue;
1008  pfrh = createHcalRecHit( detid,
1009  energy,
1010  PFLayer::HCAL_ENDCAP,
1011  hcalEndcapGeometry );
1012  }
1013  break;
1014  case HcalForward:
1015  {
1016  if(energy < thresh_HF_ ) continue;
1017  pfrh = createHcalRecHit( detid,
1018  energy,
1019  PFLayer::HF_HAD,
1020  hcalEndcapGeometry );
1021  }
1022  break;
1023  default:
1024  LogError("PFRecHitProducerHCAL")
1025  <<"HCAL rechit: unknown layer : "<<detid.subdet()<<endl;
1026  continue;
1027  }
1028 
1029  if(pfrh) {
1030  rechits.push_back( *pfrh );
1031  delete pfrh;
1032  idSortedRecHits.insert( make_pair(detid.rawId(),
1033  rechits.size()-1 ) );
1034  }
1035  }
1036 
1037 
1038  // do navigation:
1039  for(unsigned i=0; i<rechits.size(); i++ ) {
1040 
1041  findRecHitNeighbours( rechits[i], idSortedRecHits,
1042  *hcalTopology,
1043  *hcalBarrelGeometry,
1044  *hcalTopology,
1045  *hcalEndcapGeometry);
1046  } // loop for navigation
1047  } // endif hcal rechits were found
1048  } // endif clustering on rechits in hcal
1049 }
1050 
1051 
1052 
1053 
1054 
1055 
1056 reco::PFRecHit*
1057 PFRecHitProducerHCAL::createHcalRecHit( const DetId& detid,
1058  double energy,
1059  PFLayer::Layer layer,
1060  const CaloSubdetectorGeometry* geom,
1061  unsigned newDetId ) {
1062 
1063  const CaloCellGeometry *thisCell = geom->getGeometry(detid);
1064  if(!thisCell) {
1065  edm::LogError("PFRecHitProducerHCAL")
1066  <<"warning detid "<<detid.rawId()<<" not found in layer "
1067  <<layer<<endl;
1068  return 0;
1069  }
1070 
1071  const GlobalPoint& position = thisCell->getPosition();
1072 
1073  double depth_correction = 0.;
1074  switch ( layer ) {
1075  case PFLayer::HF_EM:
1076  depth_correction = position.z() > 0. ? EM_Depth_ : -EM_Depth_;
1077  break;
1078  case PFLayer::HF_HAD:
1079  depth_correction = position.z() > 0. ? HAD_Depth_ : -HAD_Depth_;
1080  break;
1081  default:
1082  break;
1083  }
1084 
1085  unsigned id = detid;
1086  if(newDetId) id = newDetId;
1087  reco::PFRecHit *rh =
1088  new reco::PFRecHit( id, layer, energy,
1089  position.x(), position.y(), position.z()+depth_correction,
1090  0,0,0 );
1091 
1092 
1093 
1094 
1095  // set the corners
1096  const CaloCellGeometry::CornersVec& corners = thisCell->getCorners();
1097 
1098  assert( corners.size() == 8 );
1099 
1100  rh->setNECorner( corners[0].x(), corners[0].y(), corners[0].z()+depth_correction );
1101  rh->setSECorner( corners[1].x(), corners[1].y(), corners[1].z()+depth_correction );
1102  rh->setSWCorner( corners[2].x(), corners[2].y(), corners[2].z()+depth_correction );
1103  rh->setNWCorner( corners[3].x(), corners[3].y(), corners[3].z()+depth_correction );
1104 
1105  return rh;
1106 }
1107 
1108 
1109 
1110 
1111 void
1112 PFRecHitProducerHCAL::findRecHitNeighbours
1113 ( reco::PFRecHit& rh,
1114  const map<unsigned,unsigned >& sortedHits,
1115  const CaloSubdetectorTopology& barrelTopology,
1116  const CaloSubdetectorGeometry& barrelGeometry,
1117  const CaloSubdetectorTopology& endcapTopology,
1118  const CaloSubdetectorGeometry& endcapGeometry ) {
1119 
1120  //cout<<"------PFRecHitProducerHcaL:findRecHitNeighbours navigation value "<<navigation_HF_<<endl;
1121  if(navigation_HF_ == false){
1122  if( rh.layer() == PFLayer::HF_HAD )
1123  return;
1124  if( rh.layer() == PFLayer::HF_EM )
1125  return;
1126  }
1127  DetId detid( rh.detId() );
1128 
1129  const CaloSubdetectorTopology* topology = 0;
1130  const CaloSubdetectorGeometry* geometry = 0;
1131 
1132  switch( rh.layer() ) {
1133  case PFLayer::ECAL_ENDCAP:
1134  topology = &endcapTopology;
1135  geometry = &endcapGeometry;
1136  break;
1137  case PFLayer::ECAL_BARREL:
1138  topology = &barrelTopology;
1139  geometry = &barrelGeometry;
1140  break;
1141  case PFLayer::HCAL_ENDCAP:
1142  topology = &endcapTopology;
1143  geometry = &endcapGeometry;
1144  break;
1145  case PFLayer::HCAL_BARREL1:
1146  topology = &barrelTopology;
1147  geometry = &barrelGeometry;
1148  break;
1149  case PFLayer::PS1:
1150  case PFLayer::PS2:
1151  topology = &barrelTopology;
1152  geometry = &barrelGeometry;
1153  break;
1154  default:
1155  assert(0);
1156  }
1157 
1158  assert( topology && geometry );
1159 
1160  CaloNavigator<DetId> navigator(detid, topology);
1161 
1162  DetId north = navigator.north();
1163 
1164  DetId northeast(0);
1165  if( north != DetId(0) ) {
1166  northeast = navigator.east();
1167  }
1168  navigator.home();
1169 
1170 
1171  DetId south = navigator.south();
1172 
1173 
1174 
1175  DetId southwest(0);
1176  if( south != DetId(0) ) {
1177  southwest = navigator.west();
1178  }
1179  navigator.home();
1180 
1181 
1182  DetId east = navigator.east();
1183  DetId southeast;
1184  if( east != DetId(0) ) {
1185  southeast = navigator.south();
1186  }
1187  navigator.home();
1188  DetId west = navigator.west();
1189  DetId northwest;
1190  if( west != DetId(0) ) {
1191  northwest = navigator.north();
1192  }
1193  navigator.home();
1194 
1195  IDH i = sortedHits.find( north.rawId() );
1196  if(i != sortedHits.end() )
1197  rh.add4Neighbour( i->second );
1198 
1199  i = sortedHits.find( northeast.rawId() );
1200  if(i != sortedHits.end() )
1201  rh.add8Neighbour( i->second );
1202 
1203  i = sortedHits.find( south.rawId() );
1204  if(i != sortedHits.end() )
1205  rh.add4Neighbour( i->second );
1206 
1207  i = sortedHits.find( southwest.rawId() );
1208  if(i != sortedHits.end() )
1209  rh.add8Neighbour( i->second );
1210 
1211  i = sortedHits.find( east.rawId() );
1212  if(i != sortedHits.end() )
1213  rh.add4Neighbour( i->second );
1214 
1215  i = sortedHits.find( southeast.rawId() );
1216  if(i != sortedHits.end() )
1217  rh.add8Neighbour( i->second );
1218 
1219  i = sortedHits.find( west.rawId() );
1220  if(i != sortedHits.end() )
1221  rh.add4Neighbour( i->second );
1222 
1223  i = sortedHits.find( northwest.rawId() );
1224  if(i != sortedHits.end() )
1225  rh.add8Neighbour( i->second );
1226 
1227 
1228 }
1229 
1230 
1231 void
1232 PFRecHitProducerHCAL::findRecHitNeighboursCT
1233 ( reco::PFRecHit& rh,
1234  const map<unsigned, unsigned >& sortedHits,
1235  const CaloSubdetectorTopology& topology ) {
1236  //cout<<"------PFRecHitProducerHcaL:findRecHitNeighboursCT navigation value "<<navigation_HF_<<endl;
1237  // cout<<"----------- rechit print out"<<endl;
1238  // if(( rh.layer() == PFLayer::HF_HAD )||(rh.layer() == PFLayer::HF_EM)) {
1239 
1240  // cout<<rh<<endl;
1241  // }
1242  if(navigation_HF_ == false){
1243  if( rh.layer() == PFLayer::HF_HAD )
1244  return;
1245  if( rh.layer() == PFLayer::HF_EM )
1246  return;
1247  }
1248  CaloTowerDetId ctDetId( rh.detId() );
1249 
1250 
1251  vector<DetId> northids = topology.north(ctDetId);
1252  vector<DetId> westids = topology.west(ctDetId);
1253  vector<DetId> southids = topology.south(ctDetId);
1254  vector<DetId> eastids = topology.east(ctDetId);
1255 
1256 
1257  CaloTowerDetId badId;
1258 
1259  // all the following detids will be CaloTowerDetId
1260  CaloTowerDetId north;
1261  CaloTowerDetId northwest;
1262  CaloTowerDetId northwest2;
1263  CaloTowerDetId west;
1264  CaloTowerDetId west2;
1265  CaloTowerDetId southwest;
1266  CaloTowerDetId southwest2;
1267  CaloTowerDetId south;
1268  CaloTowerDetId southeast;
1269  CaloTowerDetId southeast2;
1270  CaloTowerDetId east;
1271  CaloTowerDetId east2;
1272  CaloTowerDetId northeast;
1273  CaloTowerDetId northeast2;
1274 
1275  // for north and south, there is no ambiguity : 1 or 0 neighbours
1276 
1277  switch( northids.size() ) {
1278  case 0:
1279  break;
1280  case 1:
1281  north = northids[0];
1282  break;
1283  default:
1284  stringstream err("PFRecHitProducerHCAL::findRecHitNeighboursCT : incorrect number of neighbours north: ");
1285  err<<northids.size();
1286  throw( err.str() );
1287  }
1288 
1289  switch( southids.size() ) {
1290  case 0:
1291  break;
1292  case 1:
1293  south = southids[0];
1294  break;
1295  default:
1296  stringstream err("PFRecHitProducerHCAL::findRecHitNeighboursCT : incorrect number of neighbours south: ");
1297  err<<southids.size();
1298  throw( err.str() );
1299  }
1300 
1301  // for east and west, one must take care
1302  // of the pitch change in HCAL endcap.
1303 
1304  switch( eastids.size() ) {
1305  case 0:
1306  break;
1307  case 1:
1308  east = eastids[0];
1309  northeast = getNorth(east, topology);
1310  southeast = getSouth(east, topology);
1311  break;
1312  case 2:
1313  // in this case, 0 is more on the north than 1
1314  east = eastids[0];
1315  east2 = eastids[1];
1316  northeast = getNorth(east, topology );
1317  southeast = getSouth(east2, topology);
1318  northeast2 = getNorth(northeast, topology );
1319  southeast2 = getSouth(southeast, topology);
1320  break;
1321  default:
1322  stringstream err("PFRecHitProducerHCAL::findRecHitNeighboursCT : incorrect number of neighbours eastids: ");
1323  err<<eastids.size();
1324  throw( err.str() );
1325  }
1326 
1327 
1328  switch( westids.size() ) {
1329  case 0:
1330  break;
1331  case 1:
1332  west = westids[0];
1333  northwest = getNorth(west, topology);
1334  southwest = getSouth(west, topology);
1335  break;
1336  case 2:
1337  // in this case, 0 is more on the north than 1
1338  west = westids[0];
1339  west2 = westids[1];
1340  northwest = getNorth(west, topology );
1341  southwest = getSouth(west2, topology );
1342  northwest2 = getNorth(northwest, topology );
1343  southwest2 = getSouth(southwest, topology );
1344  break;
1345  default:
1346  stringstream err("PFRecHitProducerHCAL::findRecHitNeighboursCT : incorrect number of neighbours westids: ");
1347  err<< westids.size();
1348  throw( err.str() );
1349  }
1350 
1351 
1352 
1353 
1354  // find and set neighbours
1355 
1356  IDH i = sortedHits.find( north.rawId() );
1357  if(i != sortedHits.end() )
1358  rh.add4Neighbour( i->second );
1359 
1360  i = sortedHits.find( northeast.rawId() );
1361  if(i != sortedHits.end() )
1362  rh.add8Neighbour( i->second );
1363 
1364  i = sortedHits.find( northeast2.rawId() );
1365  if(i != sortedHits.end() )
1366  rh.add8Neighbour( i->second );
1367 
1368  i = sortedHits.find( south.rawId() );
1369  if(i != sortedHits.end() )
1370  rh.add4Neighbour( i->second );
1371 
1372  i = sortedHits.find( southwest.rawId() );
1373  if(i != sortedHits.end() )
1374  rh.add8Neighbour( i->second );
1375 
1376  i = sortedHits.find( southwest2.rawId() );
1377  if(i != sortedHits.end() )
1378  rh.add8Neighbour( i->second );
1379 
1380  i = sortedHits.find( east.rawId() );
1381  if(i != sortedHits.end() )
1382  rh.add4Neighbour( i->second );
1383 
1384  i = sortedHits.find( east2.rawId() );
1385  if(i != sortedHits.end() )
1386  rh.add4Neighbour( i->second );
1387 
1388  i = sortedHits.find( southeast.rawId() );
1389  if(i != sortedHits.end() )
1390  rh.add8Neighbour( i->second );
1391 
1392  i = sortedHits.find( southeast2.rawId() );
1393  if(i != sortedHits.end() )
1394  rh.add8Neighbour( i->second );
1395 
1396  i = sortedHits.find( west.rawId() );
1397  if(i != sortedHits.end() )
1398  rh.add4Neighbour( i->second );
1399 
1400  i = sortedHits.find( west2.rawId() );
1401  if(i != sortedHits.end() )
1402  rh.add4Neighbour( i->second );
1403 
1404  i = sortedHits.find( northwest.rawId() );
1405  if(i != sortedHits.end() )
1406  rh.add8Neighbour( i->second );
1407 
1408  i = sortedHits.find( northwest2.rawId() );
1409  if(i != sortedHits.end() )
1410  rh.add8Neighbour( i->second );
1411 
1412  // cout<<"----------- rechit print out"<<endl;
1413  // if(( rh.layer() == PFLayer::HF_HAD )||(rh.layer() == PFLayer::HF_EM)) {
1414 
1415  // cout<<rh<<endl;
1416  // }
1417 }
1418 
1419 
1420 
1421 DetId
1422 PFRecHitProducerHCAL::getSouth(const DetId& id,
1423  const CaloSubdetectorTopology& topology) {
1424 
1425  DetId south;
1426  vector<DetId> sids = topology.south(id);
1427  if(sids.size() == 1)
1428  south = sids[0];
1429 
1430  return south;
1431 }
1432 
1433 
1434 
1435 DetId
1436 PFRecHitProducerHCAL::getNorth(const DetId& id,
1437  const CaloSubdetectorTopology& topology) {
1438 
1439  DetId north;
1440  vector<DetId> nids = topology.north(id);
1441  if(nids.size() == 1)
1442  north = nids[0];
1443 
1444  return north;
1445 }
1446 
1447 
reco::PFRecHit::setSECorner
void setSECorner(double posx, double posy, double posz)
Definition: PFRecHit.cc:226
edm::ParameterSet::getParameter
T getParameter(std::string const &) const
reco::PFRecHit::add4Neighbour
void add4Neighbour(unsigned index)
Definition: PFRecHit.cc:176
i
int i
Definition: DBlmapReader.cc:9
reco::PFRecHit::add8Neighbour
void add8Neighbour(unsigned index)
Definition: PFRecHit.cc:181
CaloTower::constituentsSize
size_t constituentsSize() const
Definition: CaloTower.h:74
reco::PFRecHit::setNECorner
void setNECorner(double posx, double posy, double posz)
Definition: PFRecHit.cc:231
HcalDetId::subdet
HcalSubdetector subdet() const
get the subdetector
Definition: HcalDetId.h:32
CaloRecHit::detid
const DetId & detid() const
Definition: CaloRecHit.h:22
reco::PFRecHit::detId
unsigned detId() const
rechit detId
Definition: PFRecHit.h:105
edm::SortedCollection< CaloTower >::const_iterator
std::vector< CaloTower >::const_iterator const_iterator
Definition: SortedCollection.h:80
PV3DBase::y
T y() const
Definition: PV3DBase.h:63
CaloSubdetectorTopology::west
virtual std::vector< DetId > west(const DetId &id) const =0
abs
#define abs(x)
Definition: mlp_lapack.h:159
HcalCondObjectContainer::getValues
const Item * getValues(DetId fId, bool throwOnFail=true) const
Definition: HcalCondObjectContainer.h:123
PFRecHitProducerHCAL::createHcalRecHit
reco::PFRecHit * createHcalRecHit(const DetId &detid, double energy, PFLayer::Layer layer, const CaloSubdetectorGeometry *geom, unsigned newDetId=0)
Definition: PFRecHitProducerHCAL.cc:1057
PFRecHitProducer::theEcalChStatus
const EcalChannelStatus * theEcalChStatus
Definition: PFRecHitProducer.h:70
detailsBasic3DVector::z
float float float z
Definition: extBasic3DVector.h:15
position
static int position[TOTALCHAMBERS][3]
Definition: ReadPGInfo.cc:509
CaloSubdetectorGeometry::getGeometry
virtual 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:43
CaloTowerConstituentsMap::constituentsOf
std::vector< DetId > constituentsOf(const CaloTowerDetId &id) const
Get the constituent detids for this tower id ( not yet implemented )
Definition: CaloTowerConstituentsMap.cc:58
DetId::rawId
uint32_t rawId() const
get the raw id
Definition: DetId.h:45
reco::PFRecHit::layer
PFLayer::Layer layer() const
rechit layer
Definition: PFRecHit.h:108
reco::PFRecHit::setSWCorner
void setSWCorner(double posx, double posy, double posz)
Definition: PFRecHit.cc:221
PFRecHitProducerHCAL::inputTagCaloTowers_
edm::InputTag inputTagCaloTowers_
Definition: PFRecHitProducerHCAL.h:87
iEvent
int iEvent
Definition: GenABIO.cc:243
reco::PFRecHit
Particle flow rechit (rechit + geometry and topology information). See clustering algorithm in PFClus...
Definition: PFRecHit.h:31
PFRecHitProducer
Base producer for particle flow rechits (PFRecHit)
Definition: PFRecHitProducer.h:33
CaloTower::emEnergy
double emEnergy() const
Definition: CaloTower.h:79
max
const T & max(const T &a, const T &b)
Definition: MaterialBudgetTrackerHistos.cc:4
edm::Event::put
OrphanHandle< PROD > put(std::auto_ptr< PROD > product)
Put a new product.
Definition: Event.h:94
CaloRecHit::energy
float energy() const
Definition: CaloRecHit.h:19
reco::PFRecHit::setEnergyUp
void setEnergyUp(double eUp)
Definition: PFRecHit.h:79
PFRecHitProducerHCAL::inputTagHcalRecHitsHF_
edm::InputTag inputTagHcalRecHitsHF_
Definition: PFRecHitProducerHCAL.h:86
PV3DBase::z
T z() const
Definition: PV3DBase.h:64
HcalDetId::ieta
int ieta() const
get the cell ieta
Definition: HcalDetId.h:38
patZpeak.handle
tuple handle
Definition: patZpeak.py:22
j
int j
Definition: DBlmapReader.cc:9
CaloTower::constituents
const std::vector< DetId > & constituents() const
Definition: CaloTower.h:73
PFRecHitProducerHCAL::getSouth
DetId getSouth(const DetId &id, const CaloSubdetectorTopology &topology)
Definition: PFRecHitProducerHCAL.cc:1422
reco::PFRecHit::setRescale
void setRescale(double factor)
Definition: PFRecHit.h:80
edm::HandleBase::isValid
bool isValid() const
Definition: HandleBase.h:76
reco::PFRecHit::setNWCorner
void setNWCorner(double posx, double posy, double posz)
search for pointers to neighbours, using neighbours&#39; DetId.
Definition: PFRecHit.cc:216
edm::Event::getByLabel
bool getByLabel(InputTag const &tag, Handle< PROD > &result) const
Definition: Event.h:361
PFRecHitProducerHCAL::inputTagHcalRecHitsHBHE_
edm::InputTag inputTagHcalRecHitsHBHE_
Definition: PFRecHitProducerHCAL.h:85
PFRecHitProducerHCAL::findRecHitNeighbours
void findRecHitNeighbours(reco::PFRecHit &rh, const std::map< unsigned, unsigned > &sortedHits, const CaloSubdetectorTopology &barrelTopo, const CaloSubdetectorGeometry &barrelGeom, const CaloSubdetectorTopology &endcapTopo, const CaloSubdetectorGeometry &endcapGeom)
Definition: PFRecHitProducerHCAL.cc:1113
PFRecHitProducerHCAL::findRecHitNeighboursCT
void findRecHitNeighboursCT(reco::PFRecHit &rh, const std::map< unsigned, unsigned > &sortedHits, const CaloSubdetectorTopology &topology)
Definition: PFRecHitProducerHCAL.cc:1233
CaloTower::hadEnergy
double hadEnergy() const
Definition: CaloTower.h:80
CaloSubdetectorTopology::north
virtual std::vector< DetId > north(const DetId &id) const =0
PFLayer::Layer
Layer
layer definition
Definition: PFLayer.h:31
HcalDetId::iphi
int iphi() const
get the cell iphi
Definition: HcalDetId.h:40
DetId
Definition: DetId.h:20
CaloTower::id
CaloTowerDetId id() const
Definition: CaloTower.h:72
PFRecHitProducer::theHcalChStatus
const HcalChannelQuality * theHcalChStatus
Definition: PFRecHitProducer.h:69
PFRecHitProducerHCAL::createRecHits
void createRecHits(std::vector< reco::PFRecHit > &rechits, std::vector< reco::PFRecHit > &rechitsCleaned, edm::Event &, const edm::EventSetup &)
Definition: PFRecHitProducerHCAL.cc:121
edm::EventSetup::get
const T & get() const
Definition: EventSetup.h:55
EcalCondObjectContainer::const_iterator
std::vector< Item >::const_iterator const_iterator
Definition: EcalCondObjectContainer.h:15
edm::ESHandle::product
T const * product() const
Definition: ESHandle.h:62
HcalSeverityLevelComputer::getSeverityLevel
int getSeverityLevel(const DetId &myid, const uint32_t &myflag, const uint32_t &mystatus) const
Definition: HcalSeverityLevelComputer.cc:232
PFRecHitProducer::theTowerConstituentsMap
const CaloTowerConstituentsMap * theTowerConstituentsMap
Definition: PFRecHitProducer.h:71
geometry
ESHandle< TrackerGeometry > geometry
Definition: TkLasBeamFitter.cc:176
HI_PhotonSkim_cff.rechits
tuple rechits
Definition: HI_PhotonSkim_cff.py:76
PFRecHitProducer::IDH
std::map< unsigned, unsigned >::const_iterator IDH
Definition: PFRecHitProducer.h:47
CaloSubdetectorTopology::south
virtual std::vector< DetId > south(const DetId &id) const =0
EcalCondObjectContainer::find
const_iterator find(uint32_t rawId) const
Definition: EcalCondObjectContainer.h:58
PFRecHitProducer::thresh_Barrel_
double thresh_Barrel_
rechits with E &lt; threshold will not give rise to a PFRecHit
Definition: PFRecHitProducer.h:63
relval_parameters_module.energy
string energy
Definition: relval_parameters_module.py:29
EcalCondObjectContainer::end
const_iterator end() const
Definition: EcalCondObjectContainer.h:88
PFRecHitProducerHCAL::thresh_HF_
double thresh_HF_
threshold for HF
Definition: PFRecHitProducerHCAL.h:90
PFRecHitProducerHCAL::getNorth
DetId getNorth(const DetId &id, const CaloSubdetectorTopology &topology)
Definition: PFRecHitProducerHCAL.cc:1436
x
Definition: DDAxes.h:10
HiRecoJets_cff.caloTowers
tuple caloTowers
Definition: HiRecoJets_cff.py:9
HcalChannelStatus::getValue
uint32_t getValue() const
Definition: HcalChannelStatus.h:63
CaloCellGeometry::getPosition
const GlobalPoint & getPosition() const
Returns the position of reference for this cell.
Definition: CaloCellGeometry.h:78
PV3DBase::x
T x() const
Definition: PV3DBase.h:62
CaloSubdetectorTopology::east
virtual std::vector< DetId > east(const DetId &id) const =0
PFRecHitProducerHCAL::PFRecHitProducerHCAL
PFRecHitProducerHCAL(const edm::ParameterSet &)
Definition: PFRecHitProducerHCAL.cc:41
CaloCellGeometry::getCorners
virtual const CornersVec & getCorners() const =0
Returns the corner points of this cell&#39;s volume.
Definition: CaloCellGeometry.cc:59
Generated for CMSSW Reference Manual by doxygen 1.8.5