38 EBuncalibrechitCollection_Token_ = consumes<EBUncalibratedRecHitCollection>(ps.
getParameter<
edm::InputTag>(
"EBuncalibrechitCollection"));
39 EEuncalibrechitCollection_Token_ = consumes<EEUncalibratedRecHitCollection>(ps.
getParameter<
edm::InputTag>(
"EEuncalibrechitCollection"));
48 if ( outputFile_.size() != 0 ) {
49 LogInfo(
"OutputInfo") <<
" Ecal RecHits Task histograms will be saved to '" << outputFile_.c_str() <<
"'";
51 LogInfo(
"OutputInfo") <<
" Ecal RecHits Task histograms will NOT be saved";
62 meEBRecHitSimHitRatio_ = 0;
63 meEERecHitSimHitRatio_ = 0;
64 meESRecHitSimHitRatio_ = 0;
65 meEBRecHitSimHitRatio1011_ = 0;
66 meEERecHitSimHitRatio1011_ = 0;
67 meEBRecHitSimHitRatio12_ = 0;
68 meEERecHitSimHitRatio12_ = 0;
69 meEBRecHitSimHitRatio13_ = 0;
70 meEERecHitSimHitRatio13_ = 0;
71 meEBRecHitSimHitRatioGt35_ = 0;
72 meEERecHitSimHitRatioGt35_ = 0;
73 meEBUnRecHitSimHitRatio_ = 0;
74 meEEUnRecHitSimHitRatio_ = 0;
75 meEBUnRecHitSimHitRatioGt35_ = 0;
76 meEEUnRecHitSimHitRatioGt35_ = 0;
78 meEBe5x5OverSimHits_ = 0;
81 meEEe5x5OverSimHits_ = 0;
84 meEBRecHitLog10Energy_ = 0;
85 meEERecHitLog10Energy_ = 0;
86 meESRecHitLog10Energy_ = 0;
87 meEBRecHitLog10EnergyContr_ = 0;
88 meEERecHitLog10EnergyContr_ = 0;
89 meESRecHitLog10EnergyContr_ = 0;
90 meEBRecHitLog10Energy5x5Contr_ = 0;
91 meEERecHitLog10Energy5x5Contr_ = 0;
93 meEBRecHitsOccupancyFlag5_6_ = 0;
94 meEBRecHitsOccupancyFlag8_9_ = 0;
95 meEERecHitsOccupancyPlusFlag5_6_ = 0;
96 meEERecHitsOccupancyMinusFlag5_6_ = 0;
97 meEERecHitsOccupancyPlusFlag8_9_ = 0;
98 meEERecHitsOccupancyMinusFlag8_9_ = 0;
100 meEBRecHitFlags_ = 0;
101 meEBRecHitSimHitvsSimHitFlag5_6_ = 0;
102 meEBRecHitSimHitFlag6_ = 0;
103 meEBRecHitSimHitFlag7_ = 0;
104 meEB5x5RecHitSimHitvsSimHitFlag8_ = 0;
106 meEERecHitFlags_ = 0;
107 meEERecHitSimHitvsSimHitFlag5_6_ = 0;
108 meEERecHitSimHitFlag6_ = 0;
109 meEERecHitSimHitFlag7_ = 0;
123 histo =
"EcalRecHitsTask Gun Momentum";
124 meGunEnergy_ = ibooker.
book1D(histo.c_str(), histo.c_str(), 100, 0., 1000.);
126 histo =
"EcalRecHitsTask Gun Eta";
127 meGunEta_ = ibooker.
book1D(histo.c_str(), histo.c_str(), 700, -3.5, 3.5);
129 histo =
"EcalRecHitsTask Gun Phi";
130 meGunPhi_ = ibooker.
book1D(histo.c_str(), histo.c_str(), 360, 0., 360.);
132 histo =
"EcalRecHitsTask Barrel RecSimHit Ratio";
133 meEBRecHitSimHitRatio_ = ibooker.
book1D(histo.c_str(), histo.c_str(), 80, 0., 2.);
135 histo =
"EcalRecHitsTask Endcap RecSimHit Ratio";
136 meEERecHitSimHitRatio_ = ibooker.
book1D(histo.c_str(), histo.c_str(), 80, 0., 2.);
138 histo =
"EcalRecHitsTask Preshower RecSimHit Ratio";
139 meESRecHitSimHitRatio_ = ibooker.
book1D(histo.c_str(), histo.c_str(), 80, 0., 2.);
141 histo =
"EcalRecHitsTask Barrel RecSimHit Ratio gt 3p5 GeV";
142 meEBRecHitSimHitRatioGt35_ = ibooker.
book1D(histo.c_str(), histo.c_str(), 80, 0.9, 1.1);
144 histo =
"EcalRecHitsTask Endcap RecSimHit Ratio gt 3p5 GeV";
145 meEERecHitSimHitRatioGt35_ = ibooker.
book1D(histo.c_str(), histo.c_str(), 80, 0.9, 1.1);
147 histo =
"EcalRecHitsTask Barrel Unc RecSimHit Ratio";
148 meEBUnRecHitSimHitRatio_ = ibooker.
book1D(histo.c_str(), histo.c_str(), 80, 0., 2.);
150 histo =
"EcalRecHitsTask Endcap Unc RecSimHit Ratio";
151 meEEUnRecHitSimHitRatio_ = ibooker.
book1D(histo.c_str(), histo.c_str(), 80, 0., 2.);
153 histo =
"EcalRecHitsTask Barrel RecSimHit Ratio Channel Status=10 11";
154 meEBRecHitSimHitRatio1011_ = ibooker.
book1D(histo.c_str(), histo.c_str(), 80, 0., 2.);
156 histo =
"EcalRecHitsTask Endcap RecSimHit Ratio Channel Status=10 11";
157 meEERecHitSimHitRatio1011_ = ibooker.
book1D(histo.c_str(), histo.c_str(), 80, 0., 2.);
159 histo =
"EcalRecHitsTask Barrel RecSimHit Ratio Channel Status=12";
160 meEBRecHitSimHitRatio12_ = ibooker.
book1D(histo.c_str(), histo.c_str(), 80, 0., 2.);
162 histo =
"EcalRecHitsTask Endcap RecSimHit Ratio Channel Status=12";
163 meEERecHitSimHitRatio12_ = ibooker.
book1D(histo.c_str(), histo.c_str(), 80, 0., 2.);
165 histo =
"EcalRecHitsTask Barrel RecSimHit Ratio Channel Status=13";
166 meEBRecHitSimHitRatio13_ = ibooker.
book1D(histo.c_str(), histo.c_str(), 80, 0., 2.);
168 histo =
"EcalRecHitsTask Endcap RecSimHit Ratio Channel Status=13";
169 meEERecHitSimHitRatio13_ = ibooker.
book1D(histo.c_str(), histo.c_str(), 80, 0., 2.);
171 histo =
"EcalRecHitsTask Barrel Unc RecSimHit Ratio gt 3p5 GeV";
172 meEBUnRecHitSimHitRatioGt35_ = ibooker.
book1D(histo.c_str(), histo.c_str(), 80, 0.9, 1.1);
174 histo =
"EcalRecHitsTask Endcap Unc RecSimHit Ratio gt 3p5 GeV";
175 meEEUnRecHitSimHitRatioGt35_ = ibooker.
book1D(histo.c_str(), histo.c_str(), 80, 0.9, 1.1);
177 histo =
"EcalRecHitsTask Barrel Rec E5x5";
178 meEBe5x5_ = ibooker.
book1D(histo.c_str(), histo.c_str(), 4000, 0., 400.);
180 histo =
"EcalRecHitsTask Barrel Rec E5x5 over Sim E5x5";
181 meEBe5x5OverSimHits_ = ibooker.
book1D(histo.c_str(), histo.c_str(), 80, 0.9, 1.1);
183 histo =
"EcalRecHitsTask Barrel Rec E5x5 over gun energy";
184 meEBe5x5OverGun_ = ibooker.
book1D(histo.c_str(), histo.c_str(), 80, 0.9, 1.1);
186 histo =
"EcalRecHitsTask Endcap Rec E5x5";
187 meEEe5x5_ = ibooker.
book1D(histo.c_str(), histo.c_str(), 4000, 0., 400.);
189 histo =
"EcalRecHitsTask Endcap Rec E5x5 over Sim E5x5";
190 meEEe5x5OverSimHits_ = ibooker.
book1D(histo.c_str(), histo.c_str(), 80, 0.9, 1.1);
192 histo =
"EcalRecHitsTask Endcap Rec E5x5 over gun energy";
193 meEEe5x5OverGun_ = ibooker.
book1D(histo.c_str(), histo.c_str(), 80, 0.9, 1.1);
195 meEBRecHitLog10Energy_ = ibooker.
book1D(
"EcalRecHitsTask Barrel Log10 Energy",
"EcalRecHitsTask Barrel Log10 Energy", 90, -5., 4. );
196 meEERecHitLog10Energy_ = ibooker.
book1D(
"EcalRecHitsTask Endcap Log10 Energy",
"EcalRecHitsTask Endcap Log10 Energy", 90, -5., 4. );
197 meESRecHitLog10Energy_ = ibooker.
book1D(
"EcalRecHitsTask Preshower Log10 Energy",
"EcalRecHitsTask Preshower Log10 Energy", 90, -5., 4. );
198 meEBRecHitLog10EnergyContr_ = ibooker.
bookProfile(
"EcalRecHits Barrel Log10En vs Hit Contribution",
"EcalRecHits Barrel Log10En vs Hit Contribution", 90, -5., 4., 100, 0., 1. );
199 meEERecHitLog10EnergyContr_ = ibooker.
bookProfile(
"EcalRecHits Endcap Log10En vs Hit Contribution",
"EcalRecHits Endcap Log10En vs Hit Contribution", 90, -5., 4., 100, 0., 1. );
200 meESRecHitLog10EnergyContr_ = ibooker.
bookProfile(
"EcalRecHits Preshower Log10En vs Hit Contribution",
"EcalRecHits Preshower Log10En vs Hit Contribution", 90, -5., 4., 100, 0., 1. );
201 meEBRecHitLog10Energy5x5Contr_ = ibooker.
bookProfile(
"EcalRecHits Barrel Log10En5x5 vs Hit Contribution",
"EcalRecHits Barrel Log10En5x5 vs Hit Contribution", 90, -5., 4., 100, 0., 1. );
202 meEERecHitLog10Energy5x5Contr_ = ibooker.
bookProfile(
"EcalRecHits Endcap Log10En5x5 vs Hit Contribution",
"EcalRecHits Endcap Log10En5x5 vs Hit Contribution", 90, -5., 4., 100, 0., 1. );
205 histo =
"EB Occupancy Flag=5 6";
206 meEBRecHitsOccupancyFlag5_6_ = ibooker.
book2D(histo, histo, 170, -85., 85., 360, 0., 360.);
207 histo =
"EB Occupancy Flag=8 9";
208 meEBRecHitsOccupancyFlag8_9_ = ibooker.
book2D(histo, histo, 170, -85., 85., 360, 0., 360.);
210 histo =
"EE+ Occupancy Flag=5 6";
211 meEERecHitsOccupancyPlusFlag5_6_ = ibooker.
book2D(histo, histo, 100, 0., 100., 100, 0., 100.);
212 histo =
"EE- Occupancy Flag=5 6";
213 meEERecHitsOccupancyMinusFlag5_6_ = ibooker.
book2D(histo, histo, 100, 0., 100., 100, 0., 100.);
214 histo =
"EE+ Occupancy Flag=8 9";
215 meEERecHitsOccupancyPlusFlag8_9_ = ibooker.
book2D(histo, histo, 100, 0., 100., 100, 0., 100.);
216 histo =
"EE- Occupancy Flag=8 9";
217 meEERecHitsOccupancyMinusFlag8_9_ = ibooker.
book2D(histo, histo, 100, 0., 100., 100, 0., 100.);
220 histo =
"EcalRecHitsTask Barrel Reco Flags";
221 meEBRecHitFlags_ = ibooker.
book1D(histo.c_str(), histo.c_str(), 10, 0., 10.);
222 histo =
"EcalRecHitsTask Endcap Reco Flags";
223 meEERecHitFlags_ = ibooker.
book1D(histo.c_str(), histo.c_str(), 10, 0., 10.);
224 histo =
"EcalRecHitsTask Barrel RecSimHit Ratio vs SimHit Flag=5 6";
225 meEBRecHitSimHitvsSimHitFlag5_6_ = ibooker.
book2D(histo.c_str(), histo.c_str(), 80, 0., 2., 4000, 0., 400. );
226 histo =
"EcalRecHitsTask Endcap RecSimHit Ratio vs SimHit Flag=5 6";
227 meEERecHitSimHitvsSimHitFlag5_6_ = ibooker.
book2D(histo.c_str(), histo.c_str(), 80, 0., 2., 4000, 0., 400. );
228 histo =
"EcalRecHitsTask Barrel RecSimHit Ratio Flag=6";
229 meEBRecHitSimHitFlag6_ = ibooker.
book1D(histo.c_str(), histo.c_str(), 80, 0., 2.);
230 histo =
"EcalRecHitsTask Endcap RecSimHit Ratio Flag=6";
231 meEERecHitSimHitFlag6_ = ibooker.
book1D(histo.c_str(), histo.c_str(), 80, 0., 2.);
232 histo =
"EcalRecHitsTask Barrel RecSimHit Ratio Flag=7";
233 meEBRecHitSimHitFlag7_ = ibooker.
book1D(histo.c_str(), histo.c_str(), 80, 0., 2.);
234 histo =
"EcalRecHitsTask Endcap RecSimHit Ratio Flag=7";
235 meEERecHitSimHitFlag7_ = ibooker.
book1D(histo.c_str(), histo.c_str(), 80, 0., 2.);
236 histo =
"EcalRecHitsTask Barrel 5x5 RecSimHit Ratio vs SimHit Flag=8";
237 meEB5x5RecHitSimHitvsSimHitFlag8_ = ibooker.
book2D(histo.c_str(), histo.c_str(), 80, 0., 2., 4000, 0., 400. );
248 LogInfo(
"EcalRecHitsTask, EventInfo: ") <<
" Run = " << e.
id().
run() <<
" Event = " << e.
id().
event();
254 const double barrelADCtoGeV_ = agc->
getEBValue();
255 const double endcapADCtoGeV_ = agc->
getEEValue();
260 if (!MCEvt.
isValid()) { skipMC =
true; }
264 bool skipBarrel =
false;
267 e.
getByToken( EBuncalibrechitCollection_Token_, EcalUncalibRecHitEB);
268 if (EcalUncalibRecHitEB.
isValid()) {
269 EBUncalibRecHit = EcalUncalibRecHitEB.
product() ;
274 bool skipEndcap =
false;
277 e.
getByToken( EEuncalibrechitCollection_Token_, EcalUncalibRecHitEE);
278 if (EcalUncalibRecHitEE.
isValid()){
279 EEUncalibRecHit = EcalUncalibRecHitEE.
product () ;
286 e.
getByToken( EBrechitCollection_Token_, EcalRecHitEB);
288 EBRecHit = EcalRecHitEB.
product();
295 e.
getByToken( EErechitCollection_Token_, EcalRecHitEE);
297 EERecHit = EcalRecHitEE.
product ();
302 bool skipPreshower =
false;
305 e.
getByToken( ESrechitCollection_Token_, EcalRecHitES);
307 ESRecHit = EcalRecHitES.
product ();
309 skipPreshower =
true;
317 for ( HepMC::GenEvent::particle_const_iterator
p = MCEvt->GetEvent()->particles_begin();
p != MCEvt->GetEvent()->particles_end(); ++
p ) {
318 double htheta = (*p)->momentum().theta();
319 double heta = -99999.;
320 if(
tan(htheta * 0.5) > 0 ) {
321 heta = -
log(
tan(htheta * 0.5));
323 double hphi = (*p)->momentum().phi();
324 hphi = (hphi>=0) ? hphi : hphi+2*
M_PI;
325 hphi = hphi /
M_PI * 180.;
327 LogDebug(
"EventInfo") <<
"EcalRecHitsTask: Particle gun type form MC = " <<
abs((*p)->pdg_id())
328 <<
"\n" <<
"Energy = "<< (*p)->momentum().e()
329 <<
"\n" <<
"Eta = " << heta
330 <<
"\n" <<
"Phi = " << hphi;
332 if ( (*p)->momentum().e() > eGun ) eGun = (*p)->momentum().e();
334 if (meGunEnergy_) meGunEnergy_->Fill((*p)->momentum().e());
335 if (meGunEta_) meGunEta_ ->Fill(heta);
336 if (meGunPhi_) meGunPhi_ ->Fill(hphi);
347 std::auto_ptr<MixCollection<PCaloHit> >
352 const int ebcSize = 90;
353 double ebcontr[ebcSize];
354 double ebcontr25[ebcSize];
355 for(
int i=0;
i<ebcSize;
i++ ) { ebcontr[
i] = 0.0; ebcontr25[
i] = 0.0; }
362 <<
"CaloHit " << hitItr->getName() <<
" DetID = " << hitItr->id() <<
"\n"
363 <<
"Energy = " << hitItr->energy() <<
" Time = " << hitItr->time() <<
"\n"
364 <<
"EBDetId = " << ebid.
ieta() <<
" " << ebid.
iphi();
366 uint32_t crystid = ebid.
rawId();
367 ebSimMap[crystid] += hitItr->energy();
378 if( myRecHit == EBRecHit->
end() )
continue;
379 ebRecMap[EBid.
rawId()] += myRecHit->energy();
382 ebtotal += myRecHit->energy();
383 if( myRecHit->energy() > 0 ) {
384 if( meEBRecHitLog10Energy_ ) meEBRecHitLog10Energy_->Fill( log10( myRecHit->energy() ) );
385 int log10i = int( ( log10( myRecHit->energy() ) + 5. ) * 10. );
386 if( log10i >=0 and log10i < ebcSize ) ebcontr[ log10i ] += myRecHit->energy();
390 if ( ebSimMap[EBid.
rawId()] != 0. ) {
391 double uncEnergy = uncalibRecHit->amplitude()*barrelADCtoGeV_;
392 if (meEBUnRecHitSimHitRatio_) {meEBUnRecHitSimHitRatio_ ->Fill(uncEnergy/ebSimMap[EBid.
rawId()]);}
393 if (meEBUnRecHitSimHitRatioGt35_ && (myRecHit->energy()>3.5)){meEBUnRecHitSimHitRatioGt35_->Fill(uncEnergy/ebSimMap[EBid.
rawId()]);}
396 if (myRecHit != EBRecHit->
end()) {
397 if ( ebSimMap[EBid.
rawId()] != 0. ) {
398 if (meEBRecHitSimHitRatio_) {meEBRecHitSimHitRatio_ ->Fill(myRecHit->energy()/ebSimMap[EBid.
rawId()]);}
399 if (meEBRecHitSimHitRatioGt35_ && (myRecHit->energy()>3.5)){meEBRecHitSimHitRatioGt35_->Fill(myRecHit->energy()/ebSimMap[EBid.
rawId()]);}
408 if( csmi != ecs->
end() ) csc = *csmi;
412 if( meEBRecHitSimHitRatio1011_ != 0 &&
413 ( sc == 10 || sc == 11 ) ) { meEBRecHitSimHitRatio1011_->Fill(myRecHit->energy()/ebSimMap[EBid.
rawId()]); }
414 if( meEBRecHitSimHitRatio12_ != 0 && sc == 12 ) { meEBRecHitSimHitRatio12_->Fill(myRecHit->energy()/ebSimMap[EBid.
rawId()]); }
420 double ttSimEnergy = 0;
424 for( std::vector<DetId>::const_iterator dit = vid.begin(); dit != vid.end(); dit++ ) {
426 ttSimEnergy += ebSimMap[ttEBid.
rawId()];
428 if( vid.size() != 0 ) ttSimEnergy = ttSimEnergy / vid.size();
430 if( meEBRecHitSimHitRatio13_ != 0 && sc == 13 && ttSimEnergy != 0 )
431 meEBRecHitSimHitRatio13_->Fill(myRecHit->energy()/ttSimEnergy);
433 int flag = myRecHit->recoFlag();
434 if( meEBRecHitFlags_ != 0 ) meEBRecHitFlags_->Fill( flag );
436 meEBRecHitSimHitvsSimHitFlag5_6_->Fill( myRecHit->energy()/ebSimMap[EBid.
rawId()], ebSimMap[EBid.
rawId()] );
438 meEBRecHitSimHitFlag6_->Fill( myRecHit->energy()/ebSimMap[EBid.
rawId()] );
440 meEBRecHitSimHitFlag6_->Fill( myRecHit->energy()/ebSimMap[EBid.
rawId()] );
442 meEB5x5RecHitSimHitvsSimHitFlag8_->Fill( myRecHit->energy()/ttSimEnergy, ttSimEnergy );
445 meEBRecHitsOccupancyFlag5_6_ ->
Fill(EBid.
ieta(), EBid.
iphi());
447 meEBRecHitsOccupancyFlag8_9_ ->
Fill(EBid.
ieta(), EBid.
iphi());
456 uint32_t ebcenterid = getUnitWithMaxEnergy(ebRecMap);
459 int by = myEBid.
iphi();
460 int bz = myEBid.
zside();
461 findBarrelMatrix(5,5,bx,by,bz,ebRecMap);
464 for (
unsigned int i = 0;
i < crystalMatrix.size();
i++ ) {
465 e5x5rec += ebRecMap[crystalMatrix[
i]];
466 e5x5sim += ebSimMap[crystalMatrix[
i]];
467 if( ebRecMap[crystalMatrix[
i]] > 0 ) {
468 int log10i25 = int( ( log10( ebRecMap[crystalMatrix[
i]] ) + 5. ) * 10. );
469 if( log10i25 >=0 && log10i25 < ebcSize ) ebcontr25[ log10i25 ] += ebRecMap[crystalMatrix[
i]];
473 if( meEBe5x5_ ) meEBe5x5_->Fill(e5x5rec);
474 if ( e5x5sim > 0. && meEBe5x5OverSimHits_ ) meEBe5x5OverSimHits_->Fill(e5x5rec/e5x5sim);
475 if ( eGun > 0. && meEBe5x5OverGun_ ) meEBe5x5OverGun_->Fill(e5x5rec/eGun);
477 if( meEBRecHitLog10EnergyContr_ && ebtotal != 0 ) {
478 for(
int i=0;
i<ebcSize;
i++ ) {
479 meEBRecHitLog10EnergyContr_->Fill( -5.+(
float(
i)+0.5)/10., ebcontr[
i]/ebtotal );
483 if( meEBRecHitLog10Energy5x5Contr_ && e5x5rec != 0 ) {
484 for(
int i=0;
i<ebcSize;
i++ ) {
485 meEBRecHitLog10Energy5x5Contr_->Fill( -5.+(
float(
i)+0.5)/10., ebcontr25[
i]/e5x5rec );
493 if ( ! skipEndcap ) {
497 std::auto_ptr<MixCollection<PCaloHit> >
502 const int eecSize = 90;
503 double eecontr[eecSize];
504 double eecontr25[eecSize];
505 for(
int i=0;
i<eecSize;
i++ ) { eecontr[
i] = 0.0; eecontr25[
i] = 0.0; }
512 <<
" CaloHit " << hitItr->getName() <<
" DetID = " << hitItr->id() <<
"\n"
513 <<
"Energy = " << hitItr->energy() <<
" Time = " << hitItr->time() <<
"\n"
514 <<
"EEDetId side " << eeid.
zside() <<
" = " << eeid.
ix() <<
" " << eeid.
iy();
516 uint32_t crystid = eeid.
rawId();
517 eeSimMap[crystid] += hitItr->energy();
528 if( myRecHit == EERecHit->
end() )
continue;
529 eeRecMap[EEid.
rawId()] += myRecHit->energy();
532 eetotal += myRecHit->energy();
533 if( myRecHit->energy() > 0 ) {
534 if( meEERecHitLog10Energy_ ) meEERecHitLog10Energy_->Fill( log10( myRecHit->energy() ) );
535 int log10i = int( ( log10( myRecHit->energy() ) + 5. ) * 10. );
536 if( log10i >=0 and log10i < eecSize ) eecontr[ log10i ] += myRecHit->energy();
540 if ( eeSimMap[EEid.
rawId()] != 0. ) {
541 double uncEnergy = uncalibRecHit->amplitude()*endcapADCtoGeV_;
542 if (meEEUnRecHitSimHitRatio_) {meEEUnRecHitSimHitRatio_ ->Fill(uncEnergy/eeSimMap[EEid.
rawId()]);}
543 if (meEEUnRecHitSimHitRatioGt35_ && (myRecHit->energy()>3.5)){meEEUnRecHitSimHitRatioGt35_->Fill(uncEnergy/eeSimMap[EEid.
rawId()]);}
546 if (myRecHit != EERecHit->
end()) {
547 if ( eeSimMap[EEid.
rawId()] != 0. ) {
548 if (meEERecHitSimHitRatio_) {meEERecHitSimHitRatio_ ->Fill(myRecHit->energy()/eeSimMap[EEid.
rawId()]); }
549 if (meEERecHitSimHitRatioGt35_ && (myRecHit->energy()>3.5)){meEERecHitSimHitRatioGt35_->Fill(myRecHit->energy()/eeSimMap[EEid.
rawId()]); }
558 if( csmi != ecs->
end() ) csc = *csmi;
560 if( meEERecHitSimHitRatio1011_ != 0 &&
561 ( sc == 10 || sc == 11 ) ) { meEERecHitSimHitRatio1011_->Fill(myRecHit->energy()/eeSimMap[EEid.
rawId()]); }
562 if( meEERecHitSimHitRatio12_ != 0 && sc == 12 ) { meEERecHitSimHitRatio12_->Fill(myRecHit->energy()/eeSimMap[EEid.
rawId()]); }
563 if( meEERecHitSimHitRatio13_ != 0 && sc == 13 ) { meEERecHitSimHitRatio13_->Fill(myRecHit->energy()/eeSimMap[EEid.
rawId()]); }
566 int flag = myRecHit->recoFlag();
567 if( meEERecHitFlags_ != 0 ) meEERecHitFlags_->Fill( flag );
569 meEERecHitSimHitvsSimHitFlag5_6_->Fill( myRecHit->energy()/eeSimMap[EEid.
rawId()], eeSimMap[EEid.
rawId()] );
571 meEERecHitSimHitFlag6_->Fill( myRecHit->energy()/eeSimMap[EEid.
rawId()] );
573 meEERecHitSimHitFlag6_->Fill( myRecHit->energy()/eeSimMap[EEid.
rawId()] );
575 if (EEid.
zside() > 0) {
577 meEERecHitsOccupancyPlusFlag5_6_ ->
Fill(EEid.
ix(), EEid.
iy());
579 meEERecHitsOccupancyPlusFlag8_9_ ->
Fill(EEid.
ix(), EEid.
iy());
581 if (EEid.
zside() < 0) {
583 meEERecHitsOccupancyMinusFlag5_6_ ->
Fill(EEid.
ix(), EEid.
iy());
585 meEERecHitsOccupancyMinusFlag8_9_ ->
Fill(EEid.
ix(), EEid.
iy());
594 uint32_t eecenterid = getUnitWithMaxEnergy(eeRecMap);
596 int bx = myEEid.
ix();
597 int by = myEEid.
iy();
598 int bz = myEEid.
zside();
599 findEndcapMatrix(5,5,bx,by,bz,eeRecMap);
602 for (
unsigned int i = 0;
i < crystalMatrix.size();
i++ ) {
603 e5x5rec += eeRecMap[crystalMatrix[
i]];
604 e5x5sim += eeSimMap[crystalMatrix[
i]];
605 if( eeRecMap[crystalMatrix[
i]] > 0 ) {
606 int log10i25 = int( ( log10( eeRecMap[crystalMatrix[
i]] ) + 5. ) * 10. );
607 if( log10i25 >=0 && log10i25 < eecSize ) eecontr25[ log10i25 ] += eeRecMap[crystalMatrix[
i]];
611 if( meEEe5x5_ ) meEEe5x5_->Fill(e5x5rec);
612 if ( e5x5sim > 0. && meEEe5x5OverSimHits_ ) meEEe5x5OverSimHits_->Fill(e5x5rec/e5x5sim);
613 if ( eGun > 0. && meEEe5x5OverGun_ ) meEEe5x5OverGun_->Fill(e5x5rec/eGun);
616 if( meEERecHitLog10EnergyContr_ && eetotal != 0 ) {
617 for(
int i=0;
i<eecSize;
i++ ) {
618 meEERecHitLog10EnergyContr_->Fill( -5.+(
float(
i)+0.5)/10., eecontr[
i]/eetotal );
622 if( meEERecHitLog10Energy5x5Contr_ && e5x5rec != 0 ) {
623 for(
int i=0;
i<eecSize;
i++ ) {
624 meEERecHitLog10Energy5x5Contr_->Fill( -5.+(
float(
i)+0.5)/10., eecontr25[
i]/e5x5rec );
632 if ( ! skipPreshower ) {
636 std::auto_ptr<MixCollection<PCaloHit> >
640 const int escSize = 90;
641 double escontr[escSize];
642 for(
int i=0;
i<escSize;
i++ ) { escontr[
i] = 0.0; }
650 <<
" CaloHit " << hitItr->getName() <<
" DetID = " << hitItr->id() <<
"\n"
651 <<
"Energy = " << hitItr->energy() <<
" Time = " << hitItr->time() <<
"\n"
652 <<
"ESDetId strip " << esid.
strip() <<
" = " << esid.
six() <<
" " << esid.
siy();
654 uint32_t crystid = esid.
rawId();
655 esSimMap[crystid] += hitItr->energy();
662 if ( esSimMap[ESid.
rawId()] != 0. ) {
665 estotal += recHit->energy();
666 if( recHit->energy() > 0 ) {
667 if( meESRecHitLog10Energy_ ) meESRecHitLog10Energy_->Fill( log10( recHit->energy() ) );
668 int log10i = int( ( log10( recHit->energy() ) + 5. ) * 10. );
669 if( log10i >=0 and log10i < escSize ) escontr[ log10i ] += recHit->energy();
672 if (meESRecHitSimHitRatio_) {
673 meESRecHitSimHitRatio_ ->Fill(recHit->energy()/esSimMap[ESid.
rawId()]);
680 if( meESRecHitLog10EnergyContr_ && estotal != 0 ) {
681 for(
int i=0;
i<escSize;
i++ ) {
682 meESRecHitLog10EnergyContr_->Fill( -5.+(
float(
i)+0.5)/10., escontr[
i]/estotal );
695 uint32_t unitWithMaxEnergy = 0;
696 float maxEnergy = 0.;
698 MapType::iterator iter;
699 for (iter = themap.begin(); iter != themap.end(); iter++) {
701 if (maxEnergy < (*iter).second) {
702 maxEnergy = (*iter).second;
703 unitWithMaxEnergy = (*iter).first;
707 return unitWithMaxEnergy;
711 int CentralEta,
int CentralPhi,
int CentralZ,
714 int goBackInEta = nCellInEta/2;
715 int goBackInPhi = nCellInPhi/2;
716 int matrixSize = nCellInEta*nCellInPhi;
717 crystalMatrix.clear();
718 crystalMatrix.resize(matrixSize);
720 int startEta = CentralZ*CentralEta - goBackInEta;
721 int startPhi = CentralPhi - goBackInPhi;
724 for (
int ieta = startEta; ieta < startEta+nCellInEta; ieta ++ ) {
725 for(
int iphi = startPhi; iphi < startPhi + nCellInPhi; iphi++ ) {
727 if (
abs(ieta) > 85 ||
abs(ieta)<1 ) {
continue; }
728 if (iphi< 1) { index =
EBDetId(ieta,iphi+360).
rawId(); }
729 else if(iphi>360) { index =
EBDetId(ieta,iphi-360).
rawId(); }
731 crystalMatrix[i++] =
index;
738 int CentralX,
int CentralY,
int CentralZ,
740 int goBackInX = nCellInX/2;
741 int goBackInY = nCellInY/2;
742 crystalMatrix.clear();
744 int startX = CentralX - goBackInX;
745 int startY = CentralY - goBackInY;
747 for (
int ix = startX; ix < startX+nCellInX; ix ++ ) {
749 for(
int iy = startY; iy < startY + nCellInY; iy++ ) {
757 crystalMatrix.push_back(index);
T getParameter(std::string const &) const
EventNumber_t event() const
T getUntrackedParameter(std::string const &, T const &) const
~EcalRecHitsValidation()
Destructor.
MonitorElement * bookProfile(Args &&...args)
EcalRecHitsValidation(const edm::ParameterSet &ps)
Constructor.
bool getByToken(EDGetToken token, Handle< PROD > &result) const
void analyze(const edm::Event &e, const edm::EventSetup &c)
Analyze.
Code getStatusCode() const
return decoded status
std::vector< EcalUncalibratedRecHit >::const_iterator const_iterator
uint32_t getUnitWithMaxEnergy(MapType &themap)
int iphi() const
get the crystal iphi
uint32_t rawId() const
get the raw id
void findEndcapMatrix(int nCellInX, int nCellInY, int CentralX, int CentralY, int CentralZ, MapType &themap)
void Fill(HcalDetId &id, double val, std::vector< TH2F > &depth)
MonitorElement * book1D(Args &&...args)
Tan< T >::type tan(const T &t)
Abs< T >::type abs(const T &t)
EcalTrigTowerDetId tower() const
get the HCAL/trigger iphi of this crystal
int ieta() const
get the crystal ieta
std::vector< DetId > constituentsOf(const EcalTrigTowerDetId &id) const
Get the constituent detids for this tower id.
void bookHistograms(DQMStore::IBooker &i, edm::Run const &, edm::EventSetup const &) override
const_iterator end() const
void setCurrentFolder(const std::string &fullpath)
T const * product() const
MonitorElement * book2D(Args &&...args)
static bool validDetId(int crystal_ix, int crystal_iy, int iz)
std::vector< Item >::const_iterator const_iterator
T const * product() const
iterator find(key_type k)
void findBarrelMatrix(int nCellInEta, int nCellInPhi, int CentralEta, int CentralPhi, int CentralZ, MapType &themap)
const_iterator find(uint32_t rawId) const
const_iterator end() const
std::map< uint32_t, float, std::less< uint32_t > > MapType
int ietaAbs() const
get the absolute value of the crystal ieta
const_iterator begin() const
int zside() const
get the z-side of the crystal (1/-1)