15 edm::LogInfo(
"OutputInfo") <<
" Hcal RecHit Task histograms will NOT be saved";
48 if(sign_ ==
"-")
iz = -1;
49 if(sign_ ==
"*")
iz = 0;
52 if(mc_ ==
"no")
imc = 0;
70 sprintf (histo,
"HcalRecHitTask_En_rechits_cone_profile_vs_ieta_all_depths");
73 sprintf (histo,
"HcalRecHitTask_En_rechits_cone_profile_vs_ieta_all_depths_E");
76 sprintf (histo,
"HcalRecHitTask_En_rechits_cone_profile_vs_ieta_all_depths_EH");
83 sprintf (histo,
"HcalRecHitTask_energy_of_rechits_HB" ) ;
86 sprintf (histo,
"HcalRecHitTask_timing_vs_energy_profile_HB" ) ;
89 sprintf (histo,
"HcalRecHitTask_timing_vs_energy_profile_Low_HB" ) ;
92 sprintf (histo,
"HcalRecHitTask_timing_vs_energy_profile_High_HB" ) ;
96 sprintf (histo,
"HcalRecHitTask_energy_rechits_vs_simhits_HB");
98 sprintf (histo,
"HcalRecHitTask_energy_rechits_vs_simhits_profile_HB");
107 sprintf (histo,
"HcalRecHitTask_energy_of_rechits_HE" ) ;
110 sprintf (histo,
"HcalRecHitTask_timing_vs_energy_profile_Low_HE" ) ;
113 sprintf (histo,
"HcalRecHitTask_timing_vs_energy_profile_HE" ) ;
117 sprintf (histo,
"HcalRecHitTask_energy_rechits_vs_simhits_HE");
119 sprintf (histo,
"HcalRecHitTask_energy_rechits_vs_simhits_profile_HE");
128 sprintf (histo,
"HcalRecHitTask_energy_of_rechits_HO" ) ;
131 sprintf (histo,
"HcalRecHitTask_timing_vs_energy_profile_HO" ) ;
134 sprintf (histo,
"HcalRecHitTask_timing_vs_energy_profile_High_HO" ) ;
138 sprintf (histo,
"HcalRecHitTask_energy_rechits_vs_simhits_HO");
140 sprintf (histo,
"HcalRecHitTask_energy_rechits_vs_simhits_profile_HO");
149 sprintf (histo,
"HcalRecHitTask_energy_of_rechits_HF" ) ;
152 sprintf (histo,
"HcalRecHitTask_timing_vs_energy_profile_Low_HF" ) ;
155 sprintf (histo,
"HcalRecHitTask_timing_vs_energy_profile_HF" ) ;
159 sprintf (histo,
"HcalRecHitTask_energy_rechits_vs_simhits_HF");
161 sprintf (histo,
"HcalRecHitTask_energy_rechits_vs_simhits_HFL");
163 sprintf (histo,
"HcalRecHitTask_energy_rechits_vs_simhits_HFS");
165 sprintf (histo,
"HcalRecHitTask_energy_rechits_vs_simhits_profile_HF");
167 sprintf (histo,
"HcalRecHitTask_energy_rechits_vs_simhits_profile_HFL");
169 sprintf (histo,
"HcalRecHitTask_energy_rechits_vs_simhits_profile_HFS");
190 double eHcalCone = 0.;
191 double eHcalConeHB = 0.;
192 double eHcalConeHE = 0.;
193 double eHcalConeHO = 0.;
194 double eHcalConeHF = 0.;
195 double eHcalConeHFL = 0.;
196 double eHcalConeHFS = 0.;
200 int nrechitsCone = 0;
201 int nrechitsThresh = 0;
207 double eEcalCone = 0.;
208 int numrechitsEcal = 0;
211 double phi_MC = -999999.;
212 double eta_MC = -999999.;
222 edm::LogInfo(
"HcalRecHitsValidation") <<
"no HepMCProduct found";
228 double maxPt = -99999.;
230 const HepMC::GenEvent * myGenEvent = evtMC->
GetEvent();
231 for ( HepMC::GenEvent::particle_const_iterator
p = myGenEvent->particles_begin();
p != myGenEvent->particles_end(); ++
p ) {
232 double phip = (*p)->momentum().phi();
233 double etap = (*p)->momentum().eta();
236 double pt = (*p)->momentum().perp();
237 if(pt > maxPt) { npart++; maxPt =
pt; phi_MC = phip; eta_MC = etap; }
266 for (; RecHit != RecHitEnd ; ++RecHit) {
270 geometry->getSubdetectorGeometry (EBid)->getGeometry (EBid) ;
271 double eta = cellGeometry->getPosition ().eta () ;
272 double phi = cellGeometry->getPosition ().phi () ;
273 double en = RecHit->energy();
277 double r =
dR(eta_MC, phi_MC, eta, phi);
292 for (; RecHit != RecHitEnd ; ++RecHit) {
296 geometry->getSubdetectorGeometry (EEid)->getGeometry (EEid) ;
297 double eta = cellGeometry->getPosition ().eta () ;
298 double phi = cellGeometry->getPosition ().phi () ;
299 double en = RecHit->energy();
303 double r =
dR(eta_MC, phi_MC, eta, phi);
324 double HcalCone = 0.;
332 for (
unsigned int i = 0;
i <
cen.size();
i++) {
343 if(en > 1. ) nrechitsThresh++;
345 double r =
dR(eta_MC, phi_MC, eta, phi);
347 if(sub == 1) eHcalConeHB += en;
348 if(sub == 2) eHcalConeHE += en;
349 if(sub == 3) eHcalConeHO += en;
352 if (depth == 1) eHcalConeHFL += en;
353 else eHcalConeHFS += en;
361 float eta_diff = fabs(eta_MC - eta);
362 if(eta_diff < etaMax) {
419 double enSimHits = 0.;
420 double enSimHitsHB = 0.;
421 double enSimHitsHE = 0.;
422 double enSimHitsHO = 0.;
423 double enSimHitsHF = 0.;
424 double enSimHitsHFL = 0.;
425 double enSimHitsHFS = 0.;
428 for (std::vector<PCaloHit>::const_iterator SimHits = SimHitResult->begin () ; SimHits != SimHitResult->end(); ++SimHits) {
437 depth = cell.depth();
447 double en = SimHits->energy();
449 double r =
dR(eta_MC, phi_MC, etaS, phiS);
454 if(sub == static_cast<int>(
HcalBarrel)) enSimHitsHB += en;
455 if(sub == static_cast<int>(
HcalEndcap)) enSimHitsHE += en;
456 if(sub == static_cast<int>(
HcalOuter)) enSimHitsHO += en;
459 if(depth == 1) enSimHitsHFL += en;
460 else enSimHitsHFS += en;
511 if( subdet_ == 1 || subdet_ == 2 || subdet_ == 5 || subdet_ == 6 || subdet_ == 0) {
525 double eta = cellGeometry->getPosition(cell).eta () ;
526 double phi = cellGeometry->getPosition(cell).phi () ;
527 double zc = cellGeometry->getPosition(cell).z ();
528 int sub = cell.subdet();
529 int depth = cell.depth();
530 int inteta = cell.ieta();
531 if(inteta > 0) inteta -= 1;
532 int intphi = cell.iphi()-1;
533 double en = j->energy();
534 double t = j->time();
536 if((
iz > 0 && eta > 0.) || (
iz < 0 && eta <0.) ||
iz == 0) {
543 cieta.push_back(inteta);
544 ciphi.push_back(intphi);
552 if( subdet_ == 4 || subdet_ == 5 || subdet_ == 6 || subdet_ == 0) {
562 geometry->getSubdetectorGeometry (cell)->getGeometry (cell) ;
564 double eta = cellGeometry->getPosition().eta () ;
565 double phi = cellGeometry->getPosition().phi () ;
566 double zc = cellGeometry->getPosition().z ();
567 int sub = cell.subdet();
568 int depth = cell.depth();
569 int inteta = cell.ieta();
570 if(inteta > 0) inteta -= 1;
571 int intphi = cell.iphi()-1;
572 double en = j->energy();
573 double t = j->time();
575 if((
iz > 0 && eta > 0.) || (
iz < 0 && eta <0.) ||
iz == 0) {
582 cieta.push_back(inteta);
583 ciphi.push_back(intphi);
591 if( subdet_ == 3 || subdet_ == 5 || subdet_ == 6 || subdet_ == 0) {
600 geometry->getSubdetectorGeometry (cell)->getGeometry (cell) ;
602 double eta = cellGeometry->getPosition().eta () ;
603 double phi = cellGeometry->getPosition().phi () ;
604 double zc = cellGeometry->getPosition().z ();
605 int sub = cell.subdet();
606 int depth = cell.depth();
607 int inteta = cell.ieta();
608 if(inteta > 0) inteta -= 1;
609 int intphi = cell.iphi()-1;
610 double t = j->time();
611 double en = j->energy();
613 if((
iz > 0 && eta > 0.) || (
iz < 0 && eta <0.) ||
iz == 0) {
619 cieta.push_back(inteta);
620 ciphi.push_back(intphi);
629 double PI = 3.1415926535898;
630 double deltaphi= phi1 - phi2;
631 if( phi2 > phi1 ) { deltaphi= phi2 - phi1;}
632 if(deltaphi > PI) { deltaphi = 2.*PI - deltaphi;}
633 double deltaeta = eta2 - eta1;
634 double tmp =
sqrt(deltaeta* deltaeta + deltaphi*deltaphi);
642 double PI = 3.1415926535898;
643 double a1 = phi1;
double a2 = phi2;
645 if( a1 > 0.5*PI && a2 < 0.) a2 += 2*
PI;
646 if( a2 > 0.5*PI && a1 < 0.) a1 += 2*
PI;
647 tmp = (a1 * en1 + a2 * en2)/(en1 + en2);
648 if(tmp > PI) tmp -= 2.*
PI;
658 double PI = 3.1415926535898;
659 double a1 = phi1;
double a2 = phi2;
660 double tmp = a2 - a1;
662 if(a1 > 0.5 * PI) tmp += 2.*
PI;
663 if(a2 > 0.5 * PI) tmp -= 2.*
PI;
T getParameter(std::string const &) const
const CaloSubdetectorGeometry * getSubdetectorGeometry(const DetId &id) const
access the subdetector geometry for the given subdetector directly
T getUntrackedParameter(std::string const &, T const &) const
std::vector< PCaloHit > PCaloHitContainer
edm::EDGetTokenT< edm::PCaloHitContainer > tok_hh_
MonitorElement * meRecHitsEnergyHF
HcalSubdetector subdet() const
get the subdetector
double phi12(double phi1, double en1, double phi2, double en2)
MonitorElement * bookProfile(Args &&...args)
bool getByToken(EDGetToken token, Handle< PROD > &result) const
MonitorElement * meTEprofileHE_Low
std::string hcalselector_
#define DEFINE_FWK_MODULE(type)
Geom::Phi< T > phi() const
edm::EDGetTokenT< HFRecHitCollection > tok_hf_
std::vector< EcalRecHit >::const_iterator const_iterator
MonitorElement * meEnConeEtaProfile
double dPhiWsign(double phi1, double phi2)
MonitorElement * meRecHitSimHitProfileHF
MonitorElement * meRecHitSimHitProfileHFS
MonitorElement * meRecHitsEnergyHB
std::vector< double > ceta
std::string ecalselector_
MonitorElement * meRecHitSimHitHF
MonitorElement * meRecHitsEnergyHE
MonitorElement * meRecHitSimHitProfileHE
edm::EDGetTokenT< HORecHitCollection > tok_ho_
MonitorElement * meEnConeEtaProfile_E
MonitorElement * meRecHitSimHitHB
MonitorElement * meRecHitSimHitHFS
std::vector< double > cphi
MonitorElement * meTEprofileHB
MonitorElement * book1D(Args &&...args)
virtual void fillRecHitsTmp(int subdet_, edm::Event const &ev)
double dR(double eta1, double phi1, double eta2, double phi2)
edm::EDGetTokenT< HBHERecHitCollection > tok_hbhe_
MonitorElement * meRecHitSimHitHO
MonitorElement * meRecHitSimHitProfileHB
GlobalPoint getPosition(const DetId &id) const
virtual void analyze(edm::Event const &ev, edm::EventSetup const &c)
virtual void bookHistograms(DQMStore::IBooker &, edm::Run const &, edm::EventSetup const &)
const_iterator end() const
MonitorElement * meTEprofileHF_Low
edm::EDGetTokenT< EBRecHitCollection > tok_EB_
MonitorElement * meTEprofileHO_High
MonitorElement * meTEprofileHE
void setCurrentFolder(const std::string &fullpath)
const HepMC::GenEvent * GetEvent() const
T const * product() const
MonitorElement * book2D(Args &&...args)
MonitorElement * meTEprofileHF
MonitorElement * meRecHitSimHitProfileHO
edm::EDGetTokenT< EERecHitCollection > tok_EE_
std::vector< double > ctime
MonitorElement * meEnConeEtaProfile_EH
MonitorElement * meRecHitSimHitProfileHFL
std::vector< int > cdepth
std::vector< std::vector< double > > tmp
MonitorElement * meRecHitsEnergyHO
MonitorElement * meRecHitSimHitHFL
MonitorElement * meTEprofileHO
edm::EDGetTokenT< edm::HepMCProduct > tok_evt_
DetId relabel(const uint32_t testId) const
MonitorElement * meTEprofileHB_Low
MonitorElement * meRecHitSimHitHE
HcalRecHitsValidation(edm::ParameterSet const &conf)
MonitorElement * meTEprofileHB_High
std::vector< double > cen
const_iterator begin() const