12 edm::LogInfo(
"OutputInfo") <<
" Hcal RecHit Task histograms will NOT be saved";
48 if (eventype_ ==
"multi")
etype_ = 2;
51 if(sign_ ==
"-")
iz = -1;
52 if(sign_ ==
"*")
iz = 0;
55 if(mc_ ==
"no")
imc = 0;
61 sprintf (histo,
"N_HB" );
63 sprintf (histo,
"N_HE" );
65 sprintf (histo,
"N_HO" );
67 sprintf (histo,
"N_HF" );
73 for (
unsigned int i1 = 0; i1 < 82; i1++) {
74 for (
unsigned int i2 = 0; i2 < 72; i2++) {
75 for (
unsigned int i3 = 0; i3 < 4; i3++) {
76 for (
unsigned int i4 = 0; i4 < 4; i4++) {
85 sprintf (histo,
"ZSmin_map_depth1" );
87 sprintf (histo,
"ZSmin_map_depth2" );
89 sprintf (histo,
"ZSmin_map_depth3" );
91 sprintf (histo,
"ZSmin_map_depth4" );
95 sprintf (histo,
"ZS_Nreco_HB1" );
97 sprintf (histo,
"ZS_Nreco_HB2" );
99 sprintf (histo,
"ZS_Nreco_HE1" );
101 sprintf (histo,
"ZS_Nreco_HE2" );
103 sprintf (histo,
"ZS_Nreco_HE3" );
105 sprintf (histo,
"ZS_Nreco_HO" );
107 sprintf (histo,
"ZS_Nreco_HF1" );
109 sprintf (histo,
"ZS_Nreco_HF2" );
112 sprintf (histo,
"ZSmin_simple1D_HB1" );
114 sprintf (histo,
"ZSmin_simple1D_HB2" );
116 sprintf (histo,
"ZSmin_simple1D_HE1" );
118 sprintf (histo,
"ZSmin_simple1D_HE2" );
120 sprintf (histo,
"ZSmin_simple1D_HE3" );
122 sprintf (histo,
"ZSmin_simple1D_HO" );
124 sprintf (histo,
"ZSmin_simple1D_HF1" );
126 sprintf (histo,
"ZSmin_simple1D_HF2" );
129 sprintf (histo,
"ZSmin_sequential1D_HB1" );
131 sprintf (histo,
"ZSmin_sequential1D_HB2" );
133 sprintf (histo,
"ZSmin_sequential1D_HE1" );
135 sprintf (histo,
"ZSmin_sequential1D_HE2" );
137 sprintf (histo,
"ZSmin_sequential1D_HE3" );
139 sprintf (histo,
"ZSmin_sequential1D_HO" );
141 sprintf (histo,
"ZSmin_sequential1D_HF1" );
143 sprintf (histo,
"ZSmin_sequential1D_HF2" );
151 sprintf (histo,
"emap_depth1" );
153 sprintf (histo,
"emap_depth2" );
155 sprintf (histo,
"emap_depth3" );
157 sprintf (histo,
"emap_depth4" );
162 if (ecalselector_ ==
"yes") {
163 sprintf (histo,
"map_ecal" );
164 map_ecal =
dbe_->
book2D(histo, histo, 70, -3.045, 3.045, 72, -3.1415926536, 3.1415926536);
169 sprintf (histo,
"emean_vs_ieta_HB1" );
171 sprintf (histo,
"emean_vs_ieta_HB2" );
173 sprintf (histo,
"emean_vs_ieta_HE1" );
175 sprintf (histo,
"emean_vs_ieta_HE2" );
177 sprintf (histo,
"emean_vs_ieta_HE3" );
179 sprintf (histo,
"emean_vs_ieta_HO" );
181 sprintf (histo,
"emean_vs_ieta_HF1" );
183 sprintf (histo,
"emean_vs_ieta_HF2" );
187 sprintf (histo,
"RMS_vs_ieta_HB1" );
189 sprintf (histo,
"RMS_vs_ieta_HB2" );
191 sprintf (histo,
"RMS_vs_ieta_HE1" );
193 sprintf (histo,
"RMS_vs_ieta_HE2" );
195 sprintf (histo,
"RMS_vs_ieta_HE3" );
197 sprintf (histo,
"RMS_vs_ieta_HO" );
199 sprintf (histo,
"RMS_vs_ieta_HF1" );
201 sprintf (histo,
"RMS_vs_ieta_HF2" );
205 sprintf (histo,
"emean_seq_HB1" );
207 sprintf (histo,
"emean_seq_HB2" );
209 sprintf (histo,
"emean_seq_HE1" );
211 sprintf (histo,
"emean_seq_HE2" );
213 sprintf (histo,
"emean_seq_HE3" );
215 sprintf (histo,
"emean_seq_HO" );
217 sprintf (histo,
"emean_seq_HF1" );
219 sprintf (histo,
"emean_seq_HF2" );
222 sprintf (histo,
"RMS_seq_HB1" );
224 sprintf (histo,
"RMS_seq_HB2" );
226 sprintf (histo,
"RMS_seq_HE1" );
228 sprintf (histo,
"RMS_seq_HE2" );
230 sprintf (histo,
"RMS_seq_HE3" );
232 sprintf (histo,
"RMS_seq_HO" );
234 sprintf (histo,
"RMS_seq_HF1" );
236 sprintf (histo,
"RMS_seq_HF2" );
242 sprintf (histo,
"occupancy_map_HB1" );
244 sprintf (histo,
"occupancy_map_HB2" );
246 sprintf (histo,
"occupancy_map_HE1" );
248 sprintf (histo,
"occupancy_map_HE2" );
250 sprintf (histo,
"occupancy_map_HE3" );
252 sprintf (histo,
"occupancy_map_HO" );
254 sprintf (histo,
"occupancy_map_HF1" );
256 sprintf (histo,
"occupancy_map_HF2" );
260 sprintf (histo,
"occupancy_vs_ieta_HB1" );
262 sprintf (histo,
"occupancy_vs_ieta_HB2" );
264 sprintf (histo,
"occupancy_vs_ieta_HE1" );
266 sprintf (histo,
"occupancy_vs_ieta_HE2" );
268 sprintf (histo,
"occupancy_vs_ieta_HE3" );
270 sprintf (histo,
"occupancy_vs_ieta_HO" );
272 sprintf (histo,
"occupancy_vs_ieta_HF1" );
274 sprintf (histo,
"occupancy_vs_ieta_HF2" );
279 sprintf (histo,
"occ_sequential1D_HB1" );
281 sprintf (histo,
"occ_sequential1D_HB2" );
283 sprintf (histo,
"occ_sequential1D_HE1" );
285 sprintf (histo,
"occ_sequential1D_HE2" );
287 sprintf (histo,
"occ_sequential1D_HE3" );
289 sprintf (histo,
"occ_sequential1D_HO" );
291 sprintf (histo,
"occ_sequential1D_HF1" );
293 sprintf (histo,
"occ_sequential1D_HF2" );
298 sprintf (histo,
"HcalRecHitTask_RecHit_StatusWord_HB" ) ;
301 sprintf (histo,
"HcalRecHitTask_RecHit_StatusWord_HE" ) ;
304 sprintf (histo,
"HcalRecHitTask_RecHit_StatusWord_HF" ) ;
308 sprintf (histo,
"HcalRecHitTask_RecHit_StatusWord_HF67" ) ;
311 sprintf (histo,
"HcalRecHitTask_RecHit_StatusWord_HO" ) ;
315 sprintf (histo,
"HcalRecHitTask_RecHit_Aux_StatusWord_HB" ) ;
318 sprintf (histo,
"HcalRecHitTask_RecHit_Aux_StatusWord_HE" ) ;
321 sprintf (histo,
"HcalRecHitTask_RecHit_Aux_StatusWord_HF" ) ;
324 sprintf (histo,
"HcalRecHitTask_RecHit_Aux_StatusWord_HO" ) ;
328 if(
imc !=0 && useAllHistos_) {
329 sprintf (histo,
"map_econe_depth1" );
331 dbe_->
book2D(histo, histo, 520, -5.2, 5.2, 72, -3.1415926536, 3.1415926536);
332 sprintf (histo,
"map_econe_depth2" );
334 dbe_->
book2D(histo, histo, 520, -5.2, 5.2, 72, -3.1415926536, 3.1415926536);
335 sprintf (histo,
"map_econe_depth3" );
337 dbe_->
book2D(histo, histo, 520, -5.2, 5.2, 72, -3.1415926536, 3.1415926536);
338 sprintf (histo,
"map_econe_depth4" );
340 dbe_->
book2D(histo, histo, 520, -5.2, 5.2, 72, -3.1415926536, 3.1415926536);
347 if(subdet_ != 0 &&
imc != 0) {
348 sprintf (histo,
"HcalRecHitTask_En_rechits_cone_profile_vs_ieta_all_depths");
351 sprintf (histo,
"HcalRecHitTask_En_rechits_cone_profile_vs_ieta_all_depths_E");
354 sprintf (histo,
"HcalRecHitTask_En_rechits_cone_profile_vs_ieta_all_depths_EH");
359 if(subdet_ != 0 &&
imc != 0) {
363 sprintf (histo,
"HcalRecHitTask_En_rechits_cone_profile_vs_ieta_depth1");
366 sprintf (histo,
"HcalRecHitTask_En_rechits_cone_profile_vs_ieta_depth2");
369 sprintf (histo,
"HcalRecHitTask_En_rechits_cone_profile_vs_ieta_depth3");
372 sprintf (histo,
"HcalRecHitTask_En_rechits_cone_profile_vs_ieta_depth4");
377 if(
etype_ == 1 && subdet_ != 0) {
379 sprintf (histo,
"Delta_phi_cluster-MC");
382 sprintf (histo,
"Delta_eta_cluster-MC");
385 sprintf (histo,
"Delta_phi_simcluster-MC");
388 sprintf (histo,
"Delta_eta_simcluster-MC");
395 sprintf (histo,
"e_hb" ) ;
397 sprintf (histo,
"e_he" ) ;
399 sprintf (histo,
"e_ho" ) ;
401 sprintf (histo,
"e_hfl" ) ;
403 sprintf (histo,
"e_hfs" ) ;
408 if (subdet_ == 1 || subdet_ == 5 ){
412 if(
etype_ == 1 && subdet_ == 1 ) {
414 sprintf (histo,
"HcalRecHitTask_number_of_rechits_in_cone_HB" ) ;
417 sprintf (histo,
"HcalRecHitTask_sum_of_rechits_energy_in_cone_HB" ) ;
421 sprintf (histo,
"HcalRecHitTask_number_of_rechits_above_1GeV_HB");
424 sprintf (histo,
"HcalRecHitTask_sum_of_rechits_energy_HB" ) ;
427 if (ecalselector_ ==
"yes") {
429 sprintf (histo,
"HcalRecHitTask_number_of_ecalrechits_in_cone_HB");
431 sprintf (histo,
"HcalRecHitTask_energy_ecal_plus_hcal_in_cone_HB");
435 sprintf (histo,
"HcalRecHitTask_energy_hcal_vs_ecal_HB");
437 sprintf (histo,
"HcalRecHitTask_energy_ecal_plus_hcal_HB" ) ;
443 sprintf(histo,
"HcalRecHitTask_severityLevel_HB");
446 sprintf (histo,
"HcalRecHitTask_energy_of_rechits_HB" ) ;
449 sprintf (histo,
"HcalRecHitTask_timing_HB" ) ;
453 sprintf (histo,
"HcalRecHitTask_timing_vs_energy_Low_HB" ) ;
456 sprintf (histo,
"HcalRecHitTask_timing_vs_energy_HB" ) ;
459 sprintf (histo,
"HcalRecHitTask_timing_vs_energy_High_HB" ) ;
462 sprintf (histo,
"HcalRecHitTask_timing_vs_energy_profile_Low_HB" ) ;
465 sprintf (histo,
"HcalRecHitTask_timing_vs_energy_profile_HB" ) ;
468 sprintf (histo,
"HcalRecHitTask_timing_vs_energy_profile_High_HB" ) ;
473 sprintf (histo,
"HcalRecHitTask_timing_vs_energy_HB_depth1" ) ;
476 sprintf (histo,
"HcalRecHitTask_timing_vs_energy_HB_depth2" ) ;
480 sprintf (histo,
"HcalRecHitTask_energy_rechits_vs_simhits_HB");
482 sprintf (histo,
"HcalRecHitTask_energy_rechits_vs_simhits_profile_HB");
489 if ( subdet_ == 2 || subdet_ == 5 ){
493 if(
etype_ == 1 && subdet_ == 2 ) {
496 sprintf (histo,
"HcalRecHitTask_number_of_rechits_in_cone_HE" ) ;
499 sprintf (histo,
"HcalRecHitTask_sum_of_rechits_energy_in_cone_HE" ) ;
503 sprintf (histo,
"HcalRecHitTask_number_of_rechits_above_1GeV_HE");
506 sprintf (histo,
"HcalRecHitTask_sum_of_rechits_energy_HE" ) ;
509 if (ecalselector_ ==
"yes") {
510 sprintf (histo,
"HcalRecHitTask_energy_ecal_plus_hcal_HE" ) ;
513 sprintf (histo,
"HcalRecHitTask_energy_hcal_vs_ecal_HE");
516 sprintf (histo,
"HcalRecHitTask_number_of_ecalrechits_in_cone_HE");
518 sprintf (histo,
"HcalRecHitTask_energy_ecal_plus_hcal_in_cone_HE");
526 sprintf(histo,
"HcalRecHitTask_severityLevel_HE");
529 sprintf (histo,
"HcalRecHitTask_energy_of_rechits_HE" ) ;
532 sprintf (histo,
"HcalRecHitTask_timing_HE" ) ;
535 sprintf (histo,
"HcalRecHitTask_timing_vs_energy_Low_HE" ) ;
538 sprintf (histo,
"HcalRecHitTask_timing_vs_energy_HE" ) ;
541 sprintf (histo,
"HcalRecHitTask_timing_vs_energy_profile_Low_HE" ) ;
544 sprintf (histo,
"HcalRecHitTask_timing_vs_energy_profile_HE" ) ;
549 sprintf (histo,
"HcalRecHitTask_timing_vs_energy_HE_depth1" ) ;
552 sprintf (histo,
"HcalRecHitTask_timing_vs_energy_HE_depth2" ) ;
556 sprintf (histo,
"HcalRecHitTask_energy_rechits_vs_simhits_HE");
558 sprintf (histo,
"HcalRecHitTask_energy_rechits_vs_simhits_profile_HE");
566 if ( subdet_ == 3 || subdet_ == 5 ){
570 if(
etype_ == 1 && subdet_ == 3) {
572 sprintf (histo,
"HcalRecHitTask_number_of_rechits_in_cone_HO" ) ;
575 sprintf (histo,
"HcalRecHitTask_sum_of_rechits_energy_in_cone_HO" ) ;
579 sprintf (histo,
"HcalRecHitTask_number_of_rechits_above_1GeV_HO");
582 sprintf (histo,
"HcalRecHitTask_sum_of_rechits_energy_HO" ) ;
587 sprintf(histo,
"HcalRecHitTask_severityLevel_HO");
590 sprintf (histo,
"HcalRecHitTask_energy_of_rechits_HO" ) ;
593 sprintf (histo,
"HcalRecHitTask_timing_HO" ) ;
596 sprintf (histo,
"HcalRecHitTask_timing_vs_energy_HO" ) ;
599 sprintf (histo,
"HcalRecHitTask_timing_vs_energy_High_HO" ) ;
602 sprintf (histo,
"HcalRecHitTask_timing_vs_energy_profile_HO" ) ;
605 sprintf (histo,
"HcalRecHitTask_timing_vs_energy_profile_High_HO" ) ;
611 sprintf (histo,
"HcalRecHitTask_energy_rechits_vs_simhits_HO");
613 sprintf (histo,
"HcalRecHitTask_energy_rechits_vs_simhits_profile_HO");
620 if ( subdet_ == 4 || subdet_ == 5 ){
624 if(
etype_ == 1 && subdet_ == 4) {
627 sprintf (histo,
"HcalRecHitTask_number_of_rechits_in_cone_HF" ) ;
630 sprintf (histo,
"HcalRecHitTask_sum_of_rechits_energy_in_cone_HF" ) ;
633 sprintf (histo,
"HcalRecHitTask_sum_of_rechits_energy_in_cone_HFL" );
636 sprintf (histo,
"HcalRecHitTask_sum_of_rechits_energy_in_cone_HFS");
639 sprintf (histo,
"HcalRecHitTask_sum_of_rechits_energy_HF" ) ;
644 sprintf(histo,
"HcalRecHitTask_severityLevel_HF");
647 sprintf (histo,
"HcalRecHitTask_energy_of_rechits_HF" ) ;
650 sprintf (histo,
"HcalRecHitTask_timing_HF" ) ;
653 sprintf (histo,
"HcalRecHitTask_timing_vs_energy_Low_HF" ) ;
656 sprintf (histo,
"HcalRecHitTask_timing_vs_energy_HF" ) ;
659 sprintf (histo,
"HcalRecHitTask_timing_vs_energy_profile_Low_HF" ) ;
662 sprintf (histo,
"HcalRecHitTask_timing_vs_energy_profile_HF" ) ;
667 sprintf (histo,
"HcalRecHitTask_timing_vs_energy_HFL" ) ;
670 sprintf (histo,
"HcalRecHitTask_timing_vs_energy_HFS" ) ;
674 sprintf (histo,
"HcalRecHitTask_energy_rechits_vs_simhits_HF");
676 sprintf (histo,
"HcalRecHitTask_energy_rechits_vs_simhits_HFL");
678 sprintf (histo,
"HcalRecHitTask_energy_rechits_vs_simhits_HFS");
680 sprintf (histo,
"HcalRecHitTask_energy_rechits_vs_simhits_profile_HF");
682 sprintf (histo,
"HcalRecHitTask_energy_rechits_vs_simhits_profile_HFL");
684 sprintf (histo,
"HcalRecHitTask_energy_rechits_vs_simhits_profile_HFS");
710 double eHcalCone = 0.;
711 double eHcalConeHB = 0.;
712 double eHcalConeHE = 0.;
713 double eHcalConeHO = 0.;
714 double eHcalConeHF = 0.;
715 double eHcalConeHFL = 0.;
716 double eHcalConeHFS = 0.;
719 int nrechitsCone = 0;
720 int nrechitsThresh = 0;
726 double eEcalCone = 0.;
727 int numrechitsEcal = 0;
730 double phi_MC = -999999.;
731 double eta_MC = -999999.;
735 double ehcal_coneMC_1 = 0.;
736 double ehcal_coneMC_2 = 0.;
737 double ehcal_coneMC_3 = 0.;
738 double ehcal_coneMC_4 = 0.;
741 double searchR = 1.0;
746 double etaHot = 99999.;
747 double phiHot = 99999.;
760 std::cout <<
"no HepMCProduct found" << std::endl;
766 double maxPt = -99999.;
768 const HepMC::GenEvent * myGenEvent = evtMC->GetEvent();
769 for ( HepMC::GenEvent::particle_const_iterator
p = myGenEvent->particles_begin();
770 p != myGenEvent->particles_end(); ++
p ) {
771 double phip = (*p)->momentum().phi();
772 double etap = (*p)->momentum().eta();
775 double pt = (*p)->momentum().perp();
776 if(pt > maxPt) { npart++; maxPt = pt; phi_MC = phip; eta_MC = etap; }
835 ev.
getByLabel(
"ecalRecHit",
"EcalRecHitsEB", rhitEB);
840 for (; RecHit != RecHitEnd ; ++RecHit) {
844 geometry->getSubdetectorGeometry (EBid)->getGeometry (EBid) ;
845 double eta = cellGeometry->getPosition ().eta () ;
846 double phi = cellGeometry->getPosition ().phi () ;
847 double en = RecHit->energy();
853 double r =
dR(eta_MC, phi_MC, eta, phi);
863 ev.
getByLabel(
"ecalRecHit",
"EcalRecHitsEE", rhitEE);
865 RecHit = rhitEE.
product()->begin();
866 RecHitEnd = rhitEE.
product()->end();
868 for (; RecHit != RecHitEnd ; ++RecHit) {
872 geometry->getSubdetectorGeometry (EEid)->getGeometry (EEid) ;
873 double eta = cellGeometry->getPosition ().eta () ;
874 double phi = cellGeometry->getPosition ().phi () ;
875 double en = RecHit->energy();
881 double r =
dR(eta_MC, phi_MC, eta, phi);
895 double maxE = -99999.;
906 for (
unsigned int i = 0;
i <
cen.size();
i++) {
919 int index = ieta * 72 + iphi;
928 if( sub == 1 && depth == 1) nhb1++;
929 if( sub == 1 && depth == 2) nhb2++;
930 if( sub == 2 && depth == 1) nhe1++;
931 if( sub == 2 && depth == 2) nhe2++;
932 if( sub == 2 && depth == 3) nhe3++;
933 if( sub == 3 && depth == 4) nho++;
934 if( sub == 4 && depth == 1) nhf1++;
935 if( sub == 4 && depth == 2) nhf2++;
938 if( en <
emap_min[ieta+41][iphi][depth-1][sub-1] )
939 emap_min[ieta+41][iphi][depth-1][sub-1] = en;
943 if(fabs(eta) > 3.) emin = 5.;
945 double r =
dR(eta_MC, phi_MC, eta, phi);
947 if(maxE < en && en > emin) {
968 if( depth == 1 || depth == 2 ) {
971 if (ieta1 < 0) ieta1--;
981 if (depth == 1 && sub == 1 ) {
988 if (depth == 2 && sub == 1) {
995 if (depth == 1 && sub == 2) {
1002 if (depth == 2 && sub == 2) {
1009 if (depth == 3 && sub == 2) {
1023 if (depth == 1 && sub == 4) {
1030 if (depth == 2 && sub == 4) {
1041 if (depth == 1) ehcal_coneMC_1 += en;
1042 if (depth == 2) ehcal_coneMC_2 += en;
1043 if (depth == 3) ehcal_coneMC_3 += en;
1044 if (depth == 4) ehcal_coneMC_4 += en;
1049 unsigned int isw67 = 0;
1050 for (
unsigned int isw = 0; isw < 32; isw++){
1051 statadd = 0x1<<(isw);
1052 if (stwd & statadd){
1058 if (isw == 6) isw67 += 1;
1059 if (isw == 7) isw67 += 2;
1065 for (
unsigned int isw =0; isw < 32; isw++){
1066 statadd = 0x1<<(isw);
1067 if( auxstwd & statadd ){
1091 Nhb->
Fill(
double(nhb1 + nhb2));
1092 Nhe->
Fill(
double(nhe1 + nhe2 + nhe3));
1094 Nhf->
Fill(
double(nhf1 + nhf2));
1111 for (
unsigned int i = 0;
i <
cen.size();
i++) {
1138 double clusEta = 999.;
1139 double clusPhi = 999.;
1142 double HcalCone_d1 = 0.;
1143 double HcalCone_d2 = 0.;
1144 double HcalCone_d3 = 0.;
1145 double HcalCone_d4 = 0.;
1146 double HcalCone = 0.;
1148 int ietaMax1 = 9999;
1149 int ietaMax2 = 9999;
1150 int ietaMax3 = 9999;
1151 int ietaMax4 = 9999;
1153 double enMax1 = -9999.;
1154 double enMax2 = -9999.;
1155 double enMax3 = -9999.;
1156 double enMax4 = -9999.;
1169 for (
unsigned int i = 0;
i <
cen.size();
i++) {
1178 double rhot =
dR(etaHot, phiHot, eta, phi);
1179 if(rhot < partR && en > 1.) {
1180 clusEta = (clusEta * clusEn + eta * en)/(clusEn + en);
1181 clusPhi =
phi12(clusPhi, clusEn, phi, en);
1187 if(en > 1. ) nrechitsThresh++;
1189 double r =
dR(eta_MC, phi_MC, eta, phi);
1191 if(sub == 1) eHcalConeHB += en;
1192 if(sub == 2) eHcalConeHE += en;
1193 if(sub == 3) eHcalConeHO += en;
1196 if (depth == 1) eHcalConeHFL += en;
1197 else eHcalConeHFS += en;
1248 float eta_diff = fabs(eta_MC - eta);
1249 if(eta_diff < etaMax) {
1340 double phidev =
dPhiWsign(clusPhi, phi_MC);
1342 double etadev = clusEta - eta_MC;
1367 if(eHcalConeHF > eps ) {
1434 double maxES = -9999.;
1435 double etaHotS = 1000.;
1436 double phiHotS = 1000.;
1439 ev.
getByLabel(
"g4SimHits",
"HcalHits",hcalHits);
1442 double enSimHits = 0.;
1443 double enSimHitsHB = 0.;
1444 double enSimHitsHE = 0.;
1445 double enSimHitsHO = 0.;
1446 double enSimHitsHF = 0.;
1447 double enSimHitsHFL = 0.;
1448 double enSimHitsHFS = 0.;
1451 for (std::vector<PCaloHit>::const_iterator SimHits = SimHitResult->begin () ; SimHits != SimHitResult->end(); ++SimHits) {
1455 geometry->getSubdetectorGeometry (cell)->getGeometry (cell) ;
1458 double en = SimHits->energy();
1461 if(fabs(etaS) > 3.) emin = 1.;
1463 double r =
dR(eta_MC, phi_MC, etaS, phiS);
1465 if(maxES < en && en > emin ) {
1474 if(sub == 1) enSimHitsHB += en;
1475 if(sub == 2) enSimHitsHE += en;
1476 if(sub == 3) enSimHitsHO += en;
1479 int depth = cell.depth();
1480 if(depth == 1) enSimHitsHFL += en;
1481 else enSimHitsHFS += en;
1489 double clusEta = 999.;
1490 double clusPhi = 999.;
1493 for (std::vector<PCaloHit>::const_iterator SimHits = SimHitResult->begin () ; SimHits != SimHitResult->end(); ++SimHits) {
1497 geometry->getSubdetectorGeometry (cell)->getGeometry (cell) ;
1498 double etaS = cellGeometry->getPosition().eta () ;
1499 double phiS = cellGeometry->getPosition().phi () ;
1500 double en = SimHits->energy();
1503 if(fabs(etaS) > 3.) emin = 1.;
1505 double rhot =
dR(etaHotS, phiHotS, etaS, phiS);
1506 if(rhot < partR && en > emin) {
1507 clusEta = (clusEta * clusEn + etaS * en)/(clusEn + en);
1508 clusPhi =
phi12(clusPhi, clusEn, phiS, en);
1517 double phidev =
dPhiWsign(clusPhi, phi_MC);
1519 double etadev = clusEta - eta_MC;
1525 if(eHcalConeHF > eps) {
1557 using namespace edm;
1577 if( subdet_ == 1 || subdet_ == 2 || subdet_ == 5 || subdet_ == 6 || subdet_ == 0) {
1586 geometry->getSubdetectorGeometry (cell)->getGeometry (cell) ;
1587 double eta = cellGeometry->getPosition().eta () ;
1588 double phi = cellGeometry->getPosition().phi () ;
1589 double zc = cellGeometry->getPosition().z ();
1590 int sub = cell.subdet();
1591 int depth = cell.depth();
1592 int inteta = cell.ieta();
1593 if(inteta > 0) inteta -= 1;
1594 int intphi = cell.iphi()-1;
1595 double en =
j->energy();
1596 double t =
j->time();
1597 int stwd =
j->flags();
1598 int auxstwd =
j->aux();
1607 if((
iz > 0 && eta > 0.) || (
iz < 0 && eta <0.) ||
iz == 0) {
1609 csub.push_back(sub);
1611 ceta.push_back(eta);
1612 cphi.push_back(phi);
1614 cieta.push_back(inteta);
1615 ciphi.push_back(intphi);
1618 cstwd.push_back(stwd);
1625 if( subdet_ == 4 || subdet_ == 5 || subdet_ == 6 || subdet_ == 0) {
1634 geometry->getSubdetectorGeometry (cell)->getGeometry (cell) ;
1635 double eta = cellGeometry->getPosition().eta () ;
1636 double phi = cellGeometry->getPosition().phi () ;
1637 double zc = cellGeometry->getPosition().z ();
1638 int sub = cell.subdet();
1639 int depth = cell.depth();
1640 int inteta = cell.ieta();
1641 if(inteta > 0) inteta -= 1;
1642 int intphi = cell.iphi()-1;
1643 double en =
j->energy();
1644 double t =
j->time();
1645 int stwd =
j->flags();
1646 int auxstwd =
j->aux();
1653 if((
iz > 0 && eta > 0.) || (
iz < 0 && eta <0.) ||
iz == 0) {
1655 csub.push_back(sub);
1657 ceta.push_back(eta);
1658 cphi.push_back(phi);
1660 cieta.push_back(inteta);
1661 ciphi.push_back(intphi);
1664 cstwd.push_back(stwd);
1671 if( subdet_ == 3 || subdet_ == 5 || subdet_ == 6 || subdet_ == 0) {
1679 geometry->getSubdetectorGeometry (cell)->getGeometry (cell) ;
1680 double eta = cellGeometry->getPosition().eta () ;
1681 double phi = cellGeometry->getPosition().phi () ;
1682 double zc = cellGeometry->getPosition().z ();
1683 int sub = cell.subdet();
1684 int depth = cell.depth();
1685 int inteta = cell.ieta();
1686 if(inteta > 0) inteta -= 1;
1687 int intphi = cell.iphi()-1;
1688 double t =
j->time();
1689 double en =
j->energy();
1690 int stwd =
j->flags();
1691 int auxstwd =
j->aux();
1698 if((
iz > 0 && eta > 0.) || (
iz < 0 && eta <0.) ||
iz == 0) {
1699 csub.push_back(sub);
1701 ceta.push_back(eta);
1702 cphi.push_back(phi);
1704 cieta.push_back(inteta);
1705 ciphi.push_back(intphi);
1708 cstwd.push_back(stwd);
1716 double PI = 3.1415926535898;
1717 double deltaphi= phi1 - phi2;
1718 if( phi2 > phi1 ) { deltaphi= phi2 - phi1;}
1719 if(deltaphi > PI) { deltaphi = 2.*PI - deltaphi;}
1720 double deltaeta = eta2 - eta1;
1721 double tmp =
sqrt(deltaeta* deltaeta + deltaphi*deltaphi);
1729 double PI = 3.1415926535898;
1730 double a1 = phi1;
double a2 = phi2;
1732 if( a1 > 0.5*PI && a2 < 0.) a2 += 2*
PI;
1733 if( a2 > 0.5*PI && a1 < 0.) a1 += 2*
PI;
1734 tmp = (a1 * en1 + a2 * en2)/(en1 + en2);
1735 if(tmp > PI) tmp -= 2.*
PI;
1745 double PI = 3.1415926535898;
1746 double a1 = phi1;
double a2 = phi2;
1747 double tmp = a2 - a1;
1749 if(a1 > 0.5 * PI) tmp += 2.*
PI;
1750 if(a2 > 0.5 * PI) tmp -= 2.*
PI;
1760 const uint32_t recHitFlag = hit->
flags();
1765 return severityLevel;
MonitorElement * emean_vs_ieta_HB2
MonitorElement * emap_depth2
MonitorElement * RMS_seq_HB1
MonitorElement * RMS_vs_ieta_HB1
T getUntrackedParameter(std::string const &, T const &) const
MonitorElement * emean_seqHE2
MonitorElement * occupancy_vs_ieta_HO
MonitorElement * meTimeHO
const HcalSeverityLevelComputer * theHcalSevLvlComputer
MonitorElement * meSumRecHitsEnergyConeHE
MonitorElement * occupancy_seqHB1
const HcalChannelQuality * theHcalChStatus
std::vector< PCaloHit > PCaloHitContainer
MonitorElement * occupancy_map_HB1
MonitorElement * occupancy_vs_ieta_HE3
MonitorElement * map_depth2
MonitorElement * meRecHitsEnergyHF
MonitorElement * map_depth1
MonitorElement * meTE_High_HB
MonitorElement * occupancy_vs_ieta_HF1
MonitorElement * RMS_vs_ieta_HF1
HcalSubdetector subdet() const
get the subdetector
MonitorElement * book1D(const char *name, const char *title, int nchX, double lowX, double highX)
Book 1D histogram.
MonitorElement * ZS_seqHE3
double phi12(double phi1, double en1, double phi2, double en2)
MonitorElement * meDeltaPhiS
const DetId & detid() const
std::vector< int > hcalHBSevLvlVec
MonitorElement * RMS_seq_HE2
MonitorElement * meNumRecHitsConeHE
MonitorElement * meSumRecHitsEnergyConeHB
MonitorElement * meTEprofileHE_Low
std::string hcalselector_
#define DEFINE_FWK_MODULE(type)
MonitorElement * emap_depth4
std::vector< uint32_t > cstwd
MonitorElement * RMS_vs_ieta_HE3
MonitorElement * RecHit_Aux_StatusWord_HB
Geom::Phi< T > phi() const
MonitorElement * RMS_vs_ieta_HE1
MonitorElement * emean_vs_ieta_HB1
std::vector< EcalRecHit >::const_iterator const_iterator
MonitorElement * occupancy_seqHO
MonitorElement * meTE_HB2
MonitorElement * ZS_seqHE1
MonitorElement * occupancy_seqHE3
MonitorElement * meEnConeEtaProfile
MonitorElement * RMS_seq_HB2
MonitorElement * occupancy_map_HE3
double dPhiWsign(double phi1, double phi2)
MonitorElement * meRecHitSimHitProfileHF
const Item * getValues(DetId fId, bool throwOnFail=true) const
MonitorElement * meEcalHcalEnergyHB
int hcalSevLvl(const CaloRecHit *hit)
MonitorElement * meRecHitSimHitProfileHFS
MonitorElement * meEnergyHcalVsEcalHE
MonitorElement * occupancy_map_HF1
MonitorElement * meEcalHcalEnergyConeHE
MonitorElement * meRecHitsEnergyHB
MonitorElement * meSumRecHitsEnergyHF
double emap_min[82][72][4][4]
MonitorElement * emap_depth1
MonitorElement * ZS_seqHB2
MonitorElement * meNumRecHitsConeHO
std::vector< double > ceta
MonitorElement * ZS_seqHO
std::string ecalselector_
MonitorElement * occupancy_seqHF2
MonitorElement * meNumRecHitsConeHB
MonitorElement * RMS_vs_ieta_HB2
MonitorElement * meRecHitSimHitHF
MonitorElement * meRecHitsEnergyHE
MonitorElement * occupancy_map_HE1
MonitorElement * meRecHitSimHitProfileHE
MonitorElement * RMS_seq_HE3
MonitorElement * ZS_seqHE2
MonitorElement * meSumRecHitsEnergyHB
MonitorElement * RMS_seq_HO
MonitorElement * meEnConeEtaProfile_E
MonitorElement * meDeltaEtaS
MonitorElement * meRecHitSimHitHB
MonitorElement * meRecHitSimHitHFS
MonitorElement * RMS_vs_ieta_HO
MonitorElement * meSumRecHitsEnergyConeHFS
MonitorElement * RecHit_Aux_StatusWord_HO
MonitorElement * RMS_seq_HE1
std::vector< double > cphi
MonitorElement * meTE_Low_HF
MonitorElement * meNumEcalRecHitsConeHB
MonitorElement * meNumRecHitsThreshHB
MonitorElement * meSumRecHitsEnergyHO
MonitorElement * meTE_Low_HB
std::vector< int > hcalHESevLvlVec
MonitorElement * meEnergyHcalVsEcalHB
MonitorElement * map_depth3
void Fill(HcalDetId &id, double val, std::vector< TH2F > &depth)
MonitorElement * meTE_Low_HE
MonitorElement * RecHit_StatusWord_HF
MonitorElement * RecHit_Aux_StatusWord_HE
MonitorElement * emean_seqHE3
edm::ESHandle< CaloGeometry > geometry
MonitorElement * meSumRecHitsEnergyConeHO
MonitorElement * meTEprofileHB
MonitorElement * meTE_HFL
MonitorElement * ZS_seqHB1
MonitorElement * map_econe_depth1
virtual void fillRecHitsTmp(int subdet_, edm::Event const &ev)
MonitorElement * map_depth4
MonitorElement * map_econe_depth4
double dR(double eta1, double phi1, double eta2, double phi2)
MonitorElement * emean_seqHF1
MonitorElement * emean_vs_ieta_HE3
MonitorElement * meNumRecHitsThreshHO
MonitorElement * meEcalHcalEnergyConeHB
MonitorElement * RMS_seq_HF2
MonitorElement * meNumRecHitsThreshHE
MonitorElement * bookProfile(const char *name, const char *title, int nchX, double lowX, double highX, int nchY, double lowY, double highY, const char *option="s")
edm::InputTag theHFRecHitCollectionLabel
std::vector< int > hcalHFSevLvlVec
MonitorElement * meTimeHF
MonitorElement * occupancy_map_HE2
MonitorElement * occupancy_vs_ieta_HB2
MonitorElement * meSumRecHitsEnergyConeHF
bool getByLabel(InputTag const &tag, Handle< PROD > &result) const
MonitorElement * RecHit_Aux_StatusWord_HF
MonitorElement * meRecHitSimHitHO
MonitorElement * meRecHitSimHitProfileHB
MonitorElement * RecHit_StatusWord_HO
MonitorElement * meEnConeEtaProfile_depth1
MonitorElement * emean_vs_ieta_HF2
virtual void analyze(edm::Event const &ev, edm::EventSetup const &c)
MonitorElement * sevLvl_HE
MonitorElement * meTimeHB
MonitorElement * occupancy_map_HO
edm::InputTag theHORecHitCollectionLabel
MonitorElement * meTE_HE2
MonitorElement * occupancy_vs_ieta_HE2
MonitorElement * meTEprofileHF_Low
MonitorElement * occupancy_map_HF2
MonitorElement * RecHit_StatusWord_HF67
MonitorElement * meTEprofileHO_High
MonitorElement * meTEprofileHE
MonitorElement * emap_depth3
MonitorElement * occupancy_seqHB2
MonitorElement * occupancy_map_HB2
MonitorElement * meTE_HE1
MonitorElement * RMS_vs_ieta_HE2
MonitorElement * meTE_HB1
MonitorElement * meNumRecHitsConeHF
MonitorElement * meTEprofileHF
MonitorElement * meEnConeEtaProfile_depth3
MonitorElement * emean_vs_ieta_HE1
MonitorElement * meRecHitSimHitProfileHO
MonitorElement * emean_seqHO
std::vector< int > hcalHOSevLvlVec
MonitorElement * emean_vs_ieta_HE2
MonitorElement * occupancy_vs_ieta_HF2
MonitorElement * map_ecal
std::vector< double > ctime
MonitorElement * emean_seqHE1
MonitorElement * map_econe_depth2
int getSeverityLevel(const DetId &myid, const uint32_t &myflag, const uint32_t &mystatus) const
MonitorElement * RecHit_StatusWord_HE
MonitorElement * meEnConeEtaProfile_EH
MonitorElement * emean_seqHF2
T const * product() const
MonitorElement * meRecHitSimHitProfileHFL
std::vector< int > cdepth
MonitorElement * meEnConeEtaProfile_depth2
MonitorElement * RMS_vs_ieta_HF2
std::vector< std::vector< double > > tmp
MonitorElement * meRecHitsEnergyHO
std::vector< uint32_t > cauxstwd
edm::InputTag theHBHERecHitCollectionLabel
MonitorElement * sevLvl_HF
MonitorElement * emean_seqHB1
MonitorElement * meTE_HFS
MonitorElement * emean_vs_ieta_HO
MonitorElement * occupancy_seqHF1
MonitorElement * meSumRecHitsEnergyHE
MonitorElement * meRecHitSimHitHFL
MonitorElement * meSumRecHitsEnergyConeHFL
MonitorElement * meNumEcalRecHitsConeHE
MonitorElement * emean_seqHB2
MonitorElement * meTE_High_HO
MonitorElement * meEcalHcalEnergyHE
MonitorElement * emean_vs_ieta_HF1
MonitorElement * meTEprofileHO
MonitorElement * sevLvl_HO
MonitorElement * sevLvl_HB
MonitorElement * meTimeHE
MonitorElement * RMS_seq_HF1
MonitorElement * meDeltaPhi
MonitorElement * meDeltaEta
MonitorElement * meTEprofileHB_Low
MonitorElement * occupancy_vs_ieta_HB1
uint32_t getValue() const
MonitorElement * book2D(const char *name, const char *title, int nchX, double lowX, double highX, int nchY, double lowY, double highY)
Book 2D histogram.
const GlobalPoint & getPosition() const
Returns the position of reference for this cell.
MonitorElement * occupancy_seqHE2
MonitorElement * ZS_seqHF1
MonitorElement * meRecHitSimHitHE
HcalRecHitsValidation(edm::ParameterSet const &conf)
MonitorElement * map_econe_depth3
MonitorElement * occupancy_seqHE1
MonitorElement * meTEprofileHB_High
void setCurrentFolder(const std::string &fullpath)
std::vector< double > cen
MonitorElement * meEnConeEtaProfile_depth4
MonitorElement * occupancy_vs_ieta_HE1
MonitorElement * RecHit_StatusWord_HB
MonitorElement * ZS_seqHF2