9 #define BOOK1D(name,title,nbinsx,lowx,highx) \ 10 h##name = dbe_ ? dbe_->book1D(#name,title,nbinsx,lowx,highx)->getTH1F() \ 11 : new TH1F(#name,title,nbinsx,lowx,highx) 14 #define BOOK2D(name,title,nbinsx,lowx,highx,nbinsy,lowy,highy) \ 15 h##name = dbe_ ? dbe_->book2D(#name,title,nbinsx,lowx,highx,nbinsy,lowy,highy)->getTH2F() \ 16 : new TH2F(#name,title,nbinsx,lowx,highx,nbinsy,lowy,highy) 19 #define DBOOK1D(name,title,nbinsx,lowx,highx) \ 20 BOOK1D(B##name,"Barrel "#title,nbinsx,lowx,highx); BOOK1D(E##name,"Endcap "#title,nbinsx,lowx,highx); BOOK1D(F##name,"Forward "#title,nbinsx,lowx,highx); 21 #define DBOOK2D(name,title,nbinsx,lowx,highx,nbinsy,lowy,highy) \ 22 BOOK2D(B##name,"Barrel "#title,nbinsx,lowx,highx,nbinsy,lowy,highy); BOOK2D(E##name,"Endcap "#title,nbinsx,lowx,highx,nbinsy,lowy,highy); BOOK2D(F##name,"Forward "#title,nbinsx,lowx,highx,nbinsy,lowy,highy); 26 #define SETAXES(name,xtitle,ytitle) \ 27 h##name->GetXaxis()->SetTitle(xtitle); h##name->GetYaxis()->SetTitle(ytitle) 30 #define DSETAXES(name,xtitle,ytitle) \ 31 SETAXES(B##name,xtitle,ytitle);SETAXES(E##name,xtitle,ytitle);SETAXES(F##name,xtitle,ytitle); 36 #define PT (plotAgainstReco_)?"reconstructed P_{T}" :"generated P_{T}" 37 #define P (plotAgainstReco_)?"generated P" :"generated P" 56 cout <<
"Benchmark output written to file " <<
outputFile_.c_str() << endl;
61 cout <<
"No output file specified ("<<
outputFile_<<
"). Results will not be saved!" << endl;
85 cout<<
"PFJetBenchmark Setup parameters =============================================="<<endl;
86 cout <<
"Filename to write histograms " << Filename<<endl;
87 cout <<
"PFJetBenchmark debug " <<
debug_<< endl;
90 cout <<
"deltaRMax " << deltaRMax << endl;
91 cout <<
"benchmarkLabel " << benchmarkLabel_ << endl;
98 string path =
"PFTask/Benchmarks/"+ benchmarkLabel_ +
"/";
99 if (plotAgainstReco) path +=
"Reco";
else path +=
"Gen";
101 dbe_->setCurrentFolder(path);
107 cout <<
"Info: DQM is not available to provide data storage service. Using TFile to save histograms. "<<endl;
111 sprintf(cutString,
"Jet multiplicity P_{T}>%4.1f GeV",
recPt_cut);
114 BOOK1D(jetsPt,
"Jets P_{T} Distribution",100, 0, 500);
116 sprintf(cutString,
"Jets #eta Distribution |#eta|<%4.1f",
maxEta_cut);
117 BOOK1D(jetsEta,cutString,100, -5, 5);
119 BOOK2D(RPtvsEta,
"#DeltaP_{T}/P_{T} vs #eta",200, -5., 5., 200,-1,1);
120 BOOK2D(RNeutvsEta,
"R_{Neutral} vs #eta",200, -5., 5., 200,-1,1);
121 BOOK2D(RNEUTvsEta,
"R_{HCAL+ECAL} vs #eta",200, -5., 5., 200,-1,1);
122 BOOK2D(RNONLvsEta,
"R_{HCAL+ECAL - Hcal Only} vs #eta",200, -5., 5., 200,-1,1);
123 BOOK2D(RHCALvsEta,
"R_{HCAL} vs #eta",200, -5., 5., 200,-1,1);
124 BOOK2D(RHONLvsEta,
"R_{HCAL only} vs #eta",200, -5., 5., 200,-1,1);
125 BOOK2D(RCHEvsEta,
"R_{Charged} vs #eta",200, -5., 5., 200,-1,1);
126 BOOK2D(NCHvsEta,
"N_{Charged} vs #eta",200, -5., 5., 200,0.,200.);
127 BOOK2D(NCH0vsEta,
"N_{Charged} vs #eta, iter 0",200, -5., 5., 200,0.,200.);
128 BOOK2D(NCH1vsEta,
"N_{Charged} vs #eta, iter 1",200, -5., 5., 200,0.,200.);
129 BOOK2D(NCH2vsEta,
"N_{Charged} vs #eta, iter 2",200, -5., 5., 200,0.,200.);
130 BOOK2D(NCH3vsEta,
"N_{Charged} vs #eta, iter 3",200, -5., 5., 200,0.,200.);
131 BOOK2D(NCH4vsEta,
"N_{Charged} vs #eta, iter 4",200, -5., 5., 200,0.,200.);
132 BOOK2D(NCH5vsEta,
"N_{Charged} vs #eta, iter 5",200, -5., 5., 200,0.,200.);
133 BOOK2D(NCH6vsEta,
"N_{Charged} vs #eta, iter 6",200, -5., 5., 200,0.,200.);
136 DBOOK1D(RCHE,#DeltaE/E (charged had),80,-2,2);
137 DBOOK1D(RNHE,#DeltaE/E (neutral had),80,-2,2);
138 DBOOK1D(RNEE,#DeltaE/E (neutral em),80,-2,2);
139 DBOOK1D(Rneut,#DeltaE/E (neutral),80,-2,2);
140 DBOOK1D(NCH, #N_{charged},200,0.,200.);
141 DBOOK2D(RPtvsPt,#DeltaP_{
T}/P_{
T} vs P_{
T},250, 0, 500, 200,-1,1);
142 DBOOK2D(RCHEvsPt,#DeltaE/E (charged had) vs P_{
T},250, 0, 500, 120,-1,2);
143 DBOOK2D(RNHEvsPt,#DeltaE/E (neutral had) vs P_{
T},250, 0, 500, 120,-1,2);
144 DBOOK2D(RNEEvsPt,#DeltaE/E (neutral em) vs P_{
T},250, 0, 500, 120,-1,2);
145 DBOOK2D(RneutvsPt,#DeltaE/E (neutral) vs P_{
T},250, 0, 500, 120,-1,2);
146 DBOOK2D(NCHvsPt,N_{charged} vs P_{
T},250,0,500,200,0.,200.);
147 DBOOK2D(NCH0vsPt, N_{charged} vs P_{
T} iter 0,250,0,500,200,0.,200.);
148 DBOOK2D(NCH1vsPt, N_{charged} vs P_{
T} iter 1,250,0,500,200,0.,200.);
149 DBOOK2D(NCH2vsPt, N_{charged} vs P_{
T} iter 2,250,0,500,200,0.,200.);
150 DBOOK2D(NCH3vsPt, N_{charged} vs P_{
T} iter 3,250,0,500,200,0.,200.);
151 DBOOK2D(NCH4vsPt, N_{charged} vs P_{
T} iter 4,250,0,500,200,0.,200.);
152 DBOOK2D(NCH5vsPt, N_{charged} vs P_{
T} iter 5,250,0,500,200,0.,200.);
153 DBOOK2D(NCH6vsPt, N_{charged} vs P_{
T} iter 6,250,0,500,200,0.,200.);
158 DBOOK2D(RNONLvsP,#DeltaE/E (
ECAL+
HCAL-only) vs P,250, 0, 1000, 150,-1.5,1.5);
159 DBOOK2D(RHCALvsP,#DeltaE/E (
HCAL) vs P,250, 0, 1000, 150,-1.5,1.5);
160 DBOOK2D(RHONLvsP,#DeltaE/E (
HCAL only) vs P,250, 0, 1000, 150,-1.5,1.5);
161 DBOOK1D(RPt20_40,#DeltaP_{
T}/P_{
T},80,-1,1);
162 DBOOK1D(RPt40_60,#DeltaP_{
T}/P_{
T},80,-1,1);
163 DBOOK1D(RPt60_80,#DeltaP_{
T}/P_{
T},80,-1,1);
164 DBOOK1D(RPt80_100,#DeltaP_{
T}/P_{
T},80,-1,1);
165 DBOOK1D(RPt100_150,#DeltaP_{
T}/P_{
T},80,-1,1);
166 DBOOK1D(RPt150_200,#DeltaP_{
T}/P_{
T},80,-1,1);
167 DBOOK1D(RPt200_250,#DeltaP_{
T}/P_{
T},80,-1,1);
168 DBOOK1D(RPt250_300,#DeltaP_{
T}/P_{
T},80,-1,1);
169 DBOOK1D(RPt300_400,#DeltaP_{
T}/P_{
T},160,-1,1);
170 DBOOK1D(RPt400_500,#DeltaP_{
T}/P_{
T},160,-1,1);
171 DBOOK1D(RPt500_750,#DeltaP_{
T}/P_{
T},160,-1,1);
172 DBOOK1D(RPt750_1250,#DeltaP_{
T}/P_{
T},160,-1,1);
173 DBOOK1D(RPt1250_2000,#DeltaP_{
T}/P_{
T},160,-1,1);
174 DBOOK1D(RPt2000_5000,#DeltaP_{
T}/P_{
T},160,-1,1);
176 DBOOK2D(DEtavsPt,#Delta#
eta vs P_{
T},1000,0,2000,500,-0.5,0.5);
177 DBOOK2D(DPhivsPt,#Delta#
phi vs P_{
T},1000,0,2000,500,-0.5,0.5);
178 BOOK2D(DEtavsEta,
"#Delta#eta vs P_{T}",1000,-5.,+5.,500,-0.5,0.5);
179 BOOK2D(DPhivsEta,
"#Delta#phi vs P_{T}",1000,-5.,+5.,500,-0.5,0.5);
186 SETAXES(jetsEta,
"#eta",
"Number of Events");
187 SETAXES(RNeutvsEta,
"#eta",
"#DeltaE/E (Neutral)");
188 SETAXES(RNEUTvsEta,
"#eta",
"#DeltaE/E (ECAL+HCAL)");
189 SETAXES(RNONLvsEta,
"#eta",
"#DeltaE/E (ECAL+HCAL-only)");
190 SETAXES(RHCALvsEta,
"#eta",
"#DeltaE/E (HCAL)");
191 SETAXES(RHONLvsEta,
"#eta",
"#DeltaE/E (HCAL Only)");
192 SETAXES(RCHEvsEta,
"#eta",
"#DeltaE/E (Charged)");
193 SETAXES(RPtvsEta,
"#eta",
"#DeltaP_{T}/P_{T}");
194 SETAXES(DEtavsEta,
"#eta",
"#Delta#eta");
195 SETAXES(DPhivsEta,
"#eta",
"#Delta#phi");
197 DSETAXES(RPt,
"#DeltaP_{T}/P_{T}",
"Events");
198 DSETAXES(RPt20_40,
"#DeltaP_{T}/P_{T}",
"Events");
199 DSETAXES(RPt40_60,
"#DeltaP_{T}/P_{T}",
"Events");
200 DSETAXES(RPt60_80,
"#DeltaP_{T}/P_{T}",
"Events");
201 DSETAXES(RPt80_100,
"#DeltaP_{T}/P_{T}",
"Events");
202 DSETAXES(RPt100_150,
"#DeltaP_{T}/P_{T}",
"Events");
203 DSETAXES(RPt150_200,
"#DeltaP_{T}/P_{T}",
"Events");
204 DSETAXES(RPt200_250,
"#DeltaP_{T}/P_{T}",
"Events");
205 DSETAXES(RPt250_300,
"#DeltaP_{T}/P_{T}",
"Events");
206 DSETAXES(RPt300_400,
"#DeltaP_{T}/P_{T}",
"Events");
207 DSETAXES(RPt400_500,
"#DeltaP_{T}/P_{T}",
"Events");
208 DSETAXES(RPt500_750,
"#DeltaP_{T}/P_{T}",
"Events");
209 DSETAXES(RPt750_1250,
"#DeltaP_{T}/P_{T}",
"Events");
210 DSETAXES(RPt1250_2000,
"#DeltaP_{T}/P_{T}",
"Events");
211 DSETAXES(RPt2000_5000,
"#DeltaP_{T}/P_{T}",
"Events");
212 DSETAXES(RCHE,
"#DeltaE/E(charged had)",
"Events");
213 DSETAXES(RNHE,
"#DeltaE/E(neutral had)",
"Events");
214 DSETAXES(RNEE,
"#DeltaE/E(neutral em)",
"Events");
215 DSETAXES(Rneut,
"#DeltaE/E(neutral)",
"Events");
217 DSETAXES(RCHEvsPt,
PT,
"#DeltaE/E(charged had)");
218 DSETAXES(RNHEvsPt,
PT,
"#DeltaE/E(neutral had)");
219 DSETAXES(RNEEvsPt,
PT,
"#DeltaE/E(neutral em)");
220 DSETAXES(RneutvsPt,
PT,
"#DeltaE/E(neutral)");
221 DSETAXES(RHCALvsP, P,
"#DeltaE/E(HCAL)");
222 DSETAXES(RHONLvsP, P,
"#DeltaE/E(HCAL-only)");
223 DSETAXES(RNEUTvsP, P,
"#DeltaE/E(ECAL+HCAL)");
224 DSETAXES(RNONLvsP, P,
"#DeltaE/E(ECAL+HCAL-only)");
239 for(
unsigned i=0;
i<pfJets.size();
i++) {
242 unsigned highJets = 0;
243 for(
unsigned j=0; j<pfJets.size(); j++) {
244 if ( j !=
i && pfJets[j].
pt() > pfJets[
i].pt() ) highJets++;
250 double rec_pt = pfj.
pt();
251 double rec_eta = pfj.
eta();
252 double rec_phi = pfj.
phi();
269 bool Forward =
false;
270 if (
std::abs(rec_eta) < 1.4 ) Barrel =
true;
288 double true_E = truth->
p();
289 double true_pt = truth->
pt();
290 double true_eta = truth->
eta();
291 double true_phi = truth->
phi();
294 else {pt_denom = true_pt;}
296 double true_ChargedHadEnergy;
297 double true_NeutralHadEnergy;
298 double true_NeutralEmEnergy;
299 gettrue (truth, true_ChargedHadEnergy, true_NeutralHadEnergy, true_NeutralEmEnergy);
300 double true_NeutralEnergy = true_NeutralHadEnergy + true_NeutralEmEnergy;
304 double rec_NeutralEnergy = rec_NeutralHadEnergy + rec_NeutralEmEnergy;
307 std::vector <unsigned int> chMult(9, static_cast<unsigned int>(0));
308 for (
unsigned ic = 0; ic < constituents.size (); ++ic) {
309 if ( constituents[ic]->
particleId() > 3 )
continue;
311 if ( trackRef.
isNull() ) {
317 unsigned int iter = trackRef->algo()>6 ? 6: trackRef->algo();
328 double cut1 = 0.0001;
329 double cut2 = 0.0001;
330 double cut3 = 0.0001;
331 double cut4 = 0.0001;
332 double cut5 = 0.0001;
333 double cut6 = 0.0001;
334 double cut7 = 0.0001;
336 double resChargedHadEnergy= 0.;
337 double resNeutralHadEnergy= 0.;
338 double resNeutralEmEnergy= 0.;
339 double resNeutralEnergy= 0.;
341 double resHCALEnergy = 0.;
342 double resNEUTEnergy = 0.;
343 if ( rec_NeutralHadEnergy > cut6 && rec_ChargedHadEnergy < cut1 ) {
344 double true_NEUTEnergy = true_NeutralHadEnergy + true_NeutralEmEnergy;
345 double true_HCALEnergy = true_NEUTEnergy - rec_NeutralEmEnergy;
346 double rec_NEUTEnergy = rec_NeutralHadEnergy+rec_NeutralEmEnergy;
347 double rec_HCALEnergy = rec_NeutralHadEnergy;
348 resHCALEnergy = (rec_HCALEnergy-true_HCALEnergy)/rec_HCALEnergy;
349 resNEUTEnergy = (rec_NEUTEnergy-true_NEUTEnergy)/rec_NEUTEnergy;
350 if ( rec_NeutralEmEnergy > cut7 ) {
358 if (true_pt > 0.0001){
359 resPt = (rec_pt -true_pt)/true_pt ;
361 if (true_ChargedHadEnergy > cut1){
362 resChargedHadEnergy = (rec_ChargedHadEnergy- true_ChargedHadEnergy)/true_ChargedHadEnergy;
364 if (true_NeutralHadEnergy > cut2){
365 resNeutralHadEnergy = (rec_NeutralHadEnergy- true_NeutralHadEnergy)/true_NeutralHadEnergy;
368 if (rec_NeutralHadEnergy > cut3){
369 resNeutralHadEnergy = (rec_NeutralHadEnergy- true_NeutralHadEnergy)/rec_NeutralHadEnergy;
371 if (true_NeutralEmEnergy > cut4){
372 resNeutralEmEnergy = (rec_NeutralEmEnergy- true_NeutralEmEnergy)/true_NeutralEmEnergy;
374 if (true_NeutralEnergy > cut5){
375 resNeutralEnergy = (rec_NeutralEnergy- true_NeutralEnergy)/true_NeutralEnergy;
382 if ( ( resPt > 0.2 && true_pt > 100. ) ||
383 ( resPt < -0.5 && true_pt > 100. ) ) {
388 <<
" resPt = " << resPt
389 <<
" resCharged " << resChargedHadEnergy
390 <<
" resNeutralHad " << resNeutralHadEnergy
391 <<
" resNeutralEm " << resNeutralEmEnergy
392 <<
" pT (T/R) " << true_pt <<
"/" << rec_pt
393 <<
" Eta (T/R) " << truth->
eta() <<
"/" << rec_eta
394 <<
" Phi (T/R) " << truth->
phi() <<
"/" << rec_phi
400 for(
unsigned j=0; j<pfJets.size(); j++) {
402 if ( j !=
i &&
algo_->
deltaR(&pfj,&pfo) < dRo && pfo.
pt() > 0.25*pfj.
pt()) {
412 for(
unsigned j=0; j<genJets.size(); j++) {
414 if ( gjo != truth &&
algo_->
deltaR(truth,gjo) < dRgo && gjo->
pt() > 0.25*truth->
pt() ) {
421 std::cout <<
"Excess probably due to overlapping jets (DR = " <<
algo_->
deltaR(genoj,pfoj) <<
")," 422 <<
" at DeltaR(T/R) = " << dRgo <<
"/" << dRo
423 <<
" with pT(T/R) " << genoj->
pt() <<
"/" << pfoj->
pt()
424 <<
" and Eta (T/R) " << genoj->
eta() <<
"/" << pfoj->
eta()
425 <<
" and Phi (T/R) " << genoj->
phi() <<
"/" << pfoj->
phi()
434 cout <<
i <<
" =========PFJet Pt "<< rec_pt
435 <<
" eta " << rec_eta
436 <<
" phi " << rec_phi
437 <<
" Charged Had Energy " << rec_ChargedHadEnergy
438 <<
" Neutral Had Energy " << rec_NeutralHadEnergy
439 <<
" Neutral elm Energy " << rec_NeutralEmEnergy << endl;
440 cout <<
" matching Gen Jet Pt " << true_pt
441 <<
" eta " << truth->
eta()
442 <<
" phi " << truth->
phi()
443 <<
" Charged Had Energy " << true_ChargedHadEnergy
444 <<
" Neutral Had Energy " << true_NeutralHadEnergy
445 <<
" Neutral elm Energy " << true_NeutralEmEnergy << endl;
451 cout <<
"==============deltaR " << deltaR <<
" resPt " << resPt
452 <<
" resChargedHadEnergy " << resChargedHadEnergy
453 <<
" resNeutralHadEnergy " << resNeutralHadEnergy
454 <<
" resNeutralEmEnergy " << resNeutralEmEnergy
467 hNCHvsEta->Fill(true_eta, rec_ChargedMultiplicity);
477 if(plot2)
hRCHEvsEta->Fill(true_eta, resChargedHadEnergy);
478 if(plot5)
hRNeutvsEta->Fill(true_eta, resNeutralEnergy);
493 if ( pt_denom > 20. && pt_denom < 40. )
hBRPt20_40->Fill (resPt);
494 if ( pt_denom > 40. && pt_denom < 60. )
hBRPt40_60->Fill (resPt);
495 if ( pt_denom > 60. && pt_denom < 80. )
hBRPt60_80->Fill (resPt);
496 if ( pt_denom > 80. && pt_denom < 100. )
hBRPt80_100->Fill (resPt);
497 if ( pt_denom > 100. && pt_denom < 150. )
hBRPt100_150->Fill (resPt);
498 if ( pt_denom > 150. && pt_denom < 200. )
hBRPt150_200->Fill (resPt);
499 if ( pt_denom > 200. && pt_denom < 250. )
hBRPt200_250->Fill (resPt);
500 if ( pt_denom > 250. && pt_denom < 300. )
hBRPt250_300->Fill (resPt);
501 if ( pt_denom > 300. && pt_denom < 400. )
hBRPt300_400->Fill (resPt);
502 if ( pt_denom > 400. && pt_denom < 500. )
hBRPt400_500->Fill (resPt);
503 if ( pt_denom > 500. && pt_denom < 750. )
hBRPt500_750->Fill (resPt);
504 if ( pt_denom > 750. && pt_denom < 1250. )
hBRPt750_1250->Fill (resPt);
505 if ( pt_denom > 1250. && pt_denom < 2000. )
hBRPt1250_2000->Fill (resPt);
506 if ( pt_denom > 2000. && pt_denom < 5000. )
hBRPt2000_5000->Fill (resPt);
507 hBNCH->Fill(rec_ChargedMultiplicity);
516 hBNCHvsPt->Fill(pt_denom,rec_ChargedMultiplicity);
523 if(plot2)
hBRCHE->Fill(resChargedHadEnergy);
524 if(plot3)
hBRNHE->Fill(resNeutralHadEnergy);
525 if(plot4)
hBRNEE->Fill(resNeutralEmEnergy);
526 if(plot5)
hBRneut->Fill(resNeutralEnergy);
527 if(plot1)
hBRPtvsPt->Fill(pt_denom, resPt);
528 if(plot2)
hBRCHEvsPt->Fill(pt_denom, resChargedHadEnergy);
529 if(plot3)
hBRNHEvsPt->Fill(pt_denom, resNeutralHadEnergy);
530 if(plot4)
hBRNEEvsPt->Fill(pt_denom, resNeutralEmEnergy);
531 if(plot5)
hBRneutvsPt->Fill(pt_denom, resNeutralEnergy);
546 if ( pt_denom > 20. && pt_denom < 40. )
hERPt20_40->Fill (resPt);
547 if ( pt_denom > 40. && pt_denom < 60. )
hERPt40_60->Fill (resPt);
548 if ( pt_denom > 60. && pt_denom < 80. )
hERPt60_80->Fill (resPt);
549 if ( pt_denom > 80. && pt_denom < 100. )
hERPt80_100->Fill (resPt);
550 if ( pt_denom > 100. && pt_denom < 150. )
hERPt100_150->Fill (resPt);
551 if ( pt_denom > 150. && pt_denom < 200. )
hERPt150_200->Fill (resPt);
552 if ( pt_denom > 200. && pt_denom < 250. )
hERPt200_250->Fill (resPt);
553 if ( pt_denom > 250. && pt_denom < 300. )
hERPt250_300->Fill (resPt);
554 if ( pt_denom > 300. && pt_denom < 400. )
hERPt300_400->Fill (resPt);
555 if ( pt_denom > 400. && pt_denom < 500. )
hERPt400_500->Fill (resPt);
556 if ( pt_denom > 500. && pt_denom < 750. )
hERPt500_750->Fill (resPt);
557 if ( pt_denom > 750. && pt_denom < 1250. )
hERPt750_1250->Fill (resPt);
558 if ( pt_denom > 1250. && pt_denom < 2000. )
hERPt1250_2000->Fill (resPt);
559 if ( pt_denom > 2000. && pt_denom < 5000. )
hERPt2000_5000->Fill (resPt);
560 hENCH->Fill(rec_ChargedMultiplicity);
561 hENCHvsPt->Fill(pt_denom,rec_ChargedMultiplicity);
576 if(plot2)
hERCHE->Fill(resChargedHadEnergy);
577 if(plot3)
hERNHE->Fill(resNeutralHadEnergy);
578 if(plot4)
hERNEE->Fill(resNeutralEmEnergy);
579 if(plot5)
hERneut->Fill(resNeutralEnergy);
580 if(plot1)
hERPtvsPt->Fill(pt_denom, resPt);
581 if(plot2)
hERCHEvsPt->Fill(pt_denom, resChargedHadEnergy);
582 if(plot3)
hERNHEvsPt->Fill(pt_denom, resNeutralHadEnergy);
583 if(plot4)
hERNEEvsPt->Fill(pt_denom, resNeutralEmEnergy);
584 if(plot5)
hERneutvsPt->Fill(pt_denom, resNeutralEnergy);
598 if ( pt_denom > 20. && pt_denom < 40. )
hFRPt20_40->Fill (resPt);
599 if ( pt_denom > 40. && pt_denom < 60. )
hFRPt40_60->Fill (resPt);
600 if ( pt_denom > 60. && pt_denom < 80. )
hFRPt60_80->Fill (resPt);
601 if ( pt_denom > 80. && pt_denom < 100. )
hFRPt80_100->Fill (resPt);
602 if ( pt_denom > 100. && pt_denom < 150. )
hFRPt100_150->Fill (resPt);
603 if ( pt_denom > 150. && pt_denom < 200. )
hFRPt150_200->Fill (resPt);
604 if ( pt_denom > 200. && pt_denom < 250. )
hFRPt200_250->Fill (resPt);
605 if ( pt_denom > 250. && pt_denom < 300. )
hFRPt250_300->Fill (resPt);
606 if ( pt_denom > 300. && pt_denom < 400. )
hFRPt300_400->Fill (resPt);
607 if ( pt_denom > 400. && pt_denom < 500. )
hFRPt400_500->Fill (resPt);
608 if ( pt_denom > 500. && pt_denom < 750. )
hFRPt500_750->Fill (resPt);
609 if ( pt_denom > 750. && pt_denom < 1250. )
hFRPt750_1250->Fill (resPt);
610 if ( pt_denom > 1250. && pt_denom < 2000. )
hFRPt1250_2000->Fill (resPt);
611 if ( pt_denom > 2000. && pt_denom < 5000. )
hFRPt2000_5000->Fill (resPt);
618 if(plot2)
hFRCHE->Fill(resChargedHadEnergy);
619 if(plot3)
hFRNHE->Fill(resNeutralHadEnergy);
620 if(plot4)
hFRNEE->Fill(resNeutralEmEnergy);
621 if(plot5)
hFRneut->Fill(resNeutralEnergy);
622 if(plot1)
hFRPtvsPt->Fill(pt_denom, resPt);
623 if(plot2)
hFRCHEvsPt->Fill(pt_denom, resChargedHadEnergy);
624 if(plot3)
hFRNHEvsPt->Fill(pt_denom, resNeutralHadEnergy);
625 if(plot4)
hFRNEEvsPt->Fill(pt_denom, resNeutralEmEnergy);
626 if(plot5)
hFRneutvsPt->Fill(pt_denom, resNeutralEnergy);
645 double& true_NeutralHadEnergy,
double& true_NeutralEmEnergy){
647 true_NeutralEmEnergy = 0.;
648 true_ChargedHadEnergy = 0.;
649 true_NeutralHadEnergy = 0.;
651 for (
unsigned i = 0;
i < mcparts.size ();
i++) {
657 true_NeutralEmEnergy +=
e;
663 true_ChargedHadEnergy +=
e;
669 true_NeutralHadEnergy +=
e;
677 cout<<setiosflags(ios::right);
679 cout<<setprecision(3);
681 cout <<
"PFJet p/px/py/pz/pt: " << pfj->
p() <<
"/" << pfj->
px ()
682 <<
"/" << pfj->
py() <<
"/" << pfj->
pz() <<
"/" << pfj->
pt() << endl
683 <<
" eta/phi: " << pfj->
eta () <<
"/" << pfj->
phi () << endl
684 <<
" PFJet specific:" << std::endl
699 cout <<
"GenJet p/px/py/pz/pt: " << truth->
p() <<
'/' << truth->
px ()
700 <<
'/' << truth->
py() <<
'/' << truth->
pz() <<
'/' << truth->
pt() << endl
701 <<
" eta/phi: " << truth->
eta () <<
'/' << truth->
phi () << endl
702 <<
" # of constituents: " << mcparts.size() << endl;
703 cout <<
" constituents:" << endl;
704 for (
unsigned i = 0;
i < mcparts.size ();
i++) {
706 cout <<
" #" <<
i <<
" PDG code:" << mcpart->
pdgId()
707 <<
", p/pt/eta/phi: " << mcpart->
p() <<
'/' << mcpart->
pt()
708 <<
'/' << mcpart->
eta() <<
'/' << mcpart->
phi() << endl;
int pdgId() const final
PDG identifier.
#define BOOK1D(name, title, nbinsx, lowx, highx)
double eta() const final
momentum pseudorapidity
float chargedEmEnergy() const
chargedEmEnergy
double px() const final
x coordinate of momentum vector
std::vector< GenJet > GenJetCollection
collection of GenJet objects
#define SETAXES(name, xtitle, ytitle)
double pt() const final
transverse momentum
int chargedMultiplicity() const
chargedMultiplicity
Jets made from PFObjects.
float neutralEmEnergy() const
neutralEmEnergy
virtual ~PFJetBenchmark()
void gettrue(const reco::GenJet *truth, double &true_ChargedHadEnergy, double &true_NeutralHadEnergy, double &true_NeutralEmEnergy)
XYZTLorentzVectorD XYZTLorentzVector
Lorentz vector with cylindrical internal representation using pseudorapidity.
#define BOOK2D(name, title, nbinsx, lowx, highx, nbinsy, lowy, highy)
double pz() const final
z coordinate of momentum vector
static double deltaR(const T *, const U *)
double energy() const final
energy
#define DBOOK2D(name, title, nbinsx, lowx, highx, nbinsy, lowy, highy)
Abs< T >::type abs(const T &t)
Jets made from MC generator particles.
void setup(std::string Filename, bool debug, bool plotAgainstReco=false, bool onlyTwoJets=true, double deltaRMax=0.1, std::string benchmarkLabel_="ParticleFlow", double recPt=-1, double maxEta=-1, DQMStore *dbe_store=0)
virtual std::vector< const GenParticle * > getGenConstituents() const
get all constituents
bool isNull() const
Checks for null.
int neutralMultiplicity() const
neutralMultiplicity
void process(const reco::PFJetCollection &, const reco::GenJetCollection &)
std::string print() const override
Print object in details.
double p() const final
magnitude of momentum vector
double resChargedHadEnergyMax_
static const Collection::value_type * matchByDeltaR(const T *, const Collection *)
double py() const final
y coordinate of momentum vector
std::pair< OmniClusterRef, TrackingParticleRef > P
std::vector< PFJet > PFJetCollection
collection of PFJet objects
double resNeutralHadEnergyMax_
void printPFJet(const reco::PFJet *)
virtual std::vector< reco::PFCandidatePtr > getPFConstituents() const
get all constituents
float neutralHadronEnergy() const
neutralHadronEnergy
float chargedMuEnergy() const
chargedMuEnergy
double resNeutralEmEnergyMax_
#define DSETAXES(name, xtitle, ytitle)
double phi() const final
momentum azimuthal angle
void printGenJet(const reco::GenJet *)
#define DBOOK1D(name, title, nbinsx, lowx, highx)
float chargedHadronEnergy() const
chargedHadronEnergy