CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Modules Pages
JetAnalyzer_HeavyIons.cc
Go to the documentation of this file.
1 //
2 // Jet Analyzer class for heavy ion jets. for DQM jet analysis monitoring
3 // For CMSSW_7_4_X, especially reading background subtracted jets
4 // author: Raghav Kunnawalkam Elayavalli, Mohammed Zakaria (co Author)
5 // Jan 12th 2015
6 // Rutgers University, email: raghav.k.e at CERN dot CH
7 // UIC, email: mzakaria @ CERN dot CH
8 
10 
11 using namespace edm;
12 using namespace reco;
13 using namespace std;
14 
15 // declare the constructors:
16 
18  : mInputCollection(iConfig.getParameter<edm::InputTag>("src")),
19  mInputVtxCollection(iConfig.getUntrackedParameter<edm::InputTag>("srcVtx", edm::InputTag("hiSelectedVertex"))),
20  mInputPFCandCollection(iConfig.getParameter<edm::InputTag>("PFcands")),
21  mInputCsCandCollection(iConfig.exists("CScands") ? iConfig.getParameter<edm::InputTag>("CScands")
22  : edm::InputTag()),
23  mOutputFile(iConfig.getUntrackedParameter<std::string>("OutputFile", "")),
24  JetType(iConfig.getUntrackedParameter<std::string>("JetType")),
25  UEAlgo(iConfig.getUntrackedParameter<std::string>("UEAlgo")),
26  mRecoJetPtThreshold(iConfig.getParameter<double>("recoJetPtThreshold")),
27  mReverseEnergyFractionThreshold(iConfig.getParameter<double>("reverseEnergyFractionThreshold")),
28  mRThreshold(iConfig.getParameter<double>("RThreshold")),
29  JetCorrectionService(iConfig.getParameter<std::string>("JetCorrections")) {
30  std::string inputCollectionLabel(mInputCollection.label());
31 
32  isCaloJet = (std::string("calo") == JetType);
33  isJPTJet = (std::string("jpt") == JetType);
34  isPFJet = (std::string("pf") == JetType);
35 
36  //consumes
37  pvToken_ = consumes<std::vector<reco::Vertex>>(edm::InputTag("offlinePrimaryVertices"));
38  caloTowersToken_ = consumes<CaloTowerCollection>(edm::InputTag("towerMaker"));
39  if (isCaloJet)
40  caloJetsToken_ = consumes<reco::CaloJetCollection>(mInputCollection);
41  if (isJPTJet)
42  jptJetsToken_ = consumes<reco::JPTJetCollection>(mInputCollection);
43  if (isPFJet) {
44  if (std::string("Pu") == UEAlgo)
45  basicJetsToken_ = consumes<reco::BasicJetCollection>(mInputCollection);
46  if (std::string("Cs") == UEAlgo)
47  pfJetsToken_ = consumes<reco::PFJetCollection>(mInputCollection);
48  }
49 
50  pfCandToken_ = consumes<reco::PFCandidateCollection>(mInputPFCandCollection);
51  csCandToken_ = mayConsume<reco::PFCandidateCollection>(mInputCsCandCollection);
52  pfCandViewToken_ = consumes<reco::CandidateView>(mInputPFCandCollection);
53  caloCandViewToken_ = consumes<reco::CandidateView>(edm::InputTag("towerMaker"));
54 
55  centralityTag_ = iConfig.getParameter<InputTag>("centralitycollection");
56  centralityToken = consumes<reco::Centrality>(centralityTag_);
57 
58  centralityBinToken = mayConsume<int>(iConfig.exists("centralitybincollection")
59  ? iConfig.getParameter<edm::InputTag>("centralitybincollection")
60  : edm::InputTag());
61 
62  hiVertexToken_ = consumes<std::vector<reco::Vertex>>(mInputVtxCollection);
63 
64  etaToken_ = mayConsume<std::vector<double>>(iConfig.exists("etaMap") ? iConfig.getParameter<edm::InputTag>("etaMap")
65  : edm::InputTag());
66  rhoToken_ = mayConsume<std::vector<double>>(iConfig.exists("rho") ? iConfig.getParameter<edm::InputTag>("rho")
67  : edm::InputTag());
68  rhomToken_ = mayConsume<std::vector<double>>(iConfig.exists("rhom") ? iConfig.getParameter<edm::InputTag>("rhom")
69  : edm::InputTag());
70 
71  // need to initialize the PF cand histograms : which are also event variables
72  if (isPFJet) {
73  mNPFpart = nullptr;
74  mPFPt = nullptr;
75  mPFEta = nullptr;
76  mPFPhi = nullptr;
77 
78  mSumPFPt = nullptr;
79  mSumPFPt_eta = nullptr;
80  mSumSquaredPFPt = nullptr;
81  mSumSquaredPFPt_eta = nullptr;
82  mSumPFPt_HF = nullptr;
83 
84  mPFDeltaR = nullptr;
85  mPFDeltaR_Scaled_R = nullptr;
86 
87  for (int ieta = 0; ieta < etaBins_; ieta++) {
88  mSumPFPtEtaDep[ieta] = nullptr;
89  }
90 
91  //cs-specific histograms
92  mRhoDist_vsEta = nullptr;
93  mRhoMDist_vsEta = nullptr;
94  mRhoDist_vsPt = nullptr;
95  mRhoMDist_vsPt = nullptr;
96 
97  rhoEtaRange = nullptr;
98  for (int ieta = 0; ieta < etaBins_; ieta++) {
99  mCSCandpT_vsPt[ieta] = nullptr;
100  mRhoDist_vsCent[ieta] = nullptr;
101  mRhoMDist_vsCent[ieta] = nullptr;
102  for (int ipt = 0; ipt < ptBins_; ipt++) {
103  mSubtractedEFrac[ipt][ieta] = nullptr;
104  mSubtractedE[ipt][ieta] = nullptr;
105  }
106  }
107 
108  mPFCandpT_vs_eta_Unknown = nullptr; // pf id 0
109  mPFCandpT_vs_eta_ChargedHadron = nullptr; // pf id - 1
110  mPFCandpT_vs_eta_electron = nullptr; // pf id - 2
111  mPFCandpT_vs_eta_muon = nullptr; // pf id - 3
112  mPFCandpT_vs_eta_photon = nullptr; // pf id - 4
113  mPFCandpT_vs_eta_NeutralHadron = nullptr; // pf id - 5
114  mPFCandpT_vs_eta_HadE_inHF = nullptr; // pf id - 6
115  mPFCandpT_vs_eta_EME_inHF = nullptr; // pf id - 7
116 
117  mPFCandpT_Barrel_Unknown = nullptr; // pf id 0
118  mPFCandpT_Barrel_ChargedHadron = nullptr; // pf id - 1
119  mPFCandpT_Barrel_electron = nullptr; // pf id - 2
120  mPFCandpT_Barrel_muon = nullptr; // pf id - 3
121  mPFCandpT_Barrel_photon = nullptr; // pf id - 4
122  mPFCandpT_Barrel_NeutralHadron = nullptr; // pf id - 5
123  mPFCandpT_Barrel_HadE_inHF = nullptr; // pf id - 6
124  mPFCandpT_Barrel_EME_inHF = nullptr; // pf id - 7
125 
126  mPFCandpT_Endcap_Unknown = nullptr; // pf id 0
127  mPFCandpT_Endcap_ChargedHadron = nullptr; // pf id - 1
128  mPFCandpT_Endcap_electron = nullptr; // pf id - 2
129  mPFCandpT_Endcap_muon = nullptr; // pf id - 3
130  mPFCandpT_Endcap_photon = nullptr; // pf id - 4
131  mPFCandpT_Endcap_NeutralHadron = nullptr; // pf id - 5
132  mPFCandpT_Endcap_HadE_inHF = nullptr; // pf id - 6
133  mPFCandpT_Endcap_EME_inHF = nullptr; // pf id - 7
134 
135  mPFCandpT_Forward_Unknown = nullptr; // pf id 0
136  mPFCandpT_Forward_ChargedHadron = nullptr; // pf id - 1
137  mPFCandpT_Forward_electron = nullptr; // pf id - 2
138  mPFCandpT_Forward_muon = nullptr; // pf id - 3
139  mPFCandpT_Forward_photon = nullptr; // pf id - 4
140  mPFCandpT_Forward_NeutralHadron = nullptr; // pf id - 5
141  mPFCandpT_Forward_HadE_inHF = nullptr; // pf id - 6
142  mPFCandpT_Forward_EME_inHF = nullptr; // pf id - 7
143  }
144  if (isCaloJet) {
145  mNCalopart = nullptr;
146  mCaloPt = nullptr;
147  mCaloEta = nullptr;
148  mCaloPhi = nullptr;
149 
150  mSumCaloPt = nullptr;
151  mSumCaloPt_eta = nullptr;
152  mSumSquaredCaloPt = nullptr;
153  mSumSquaredCaloPt_eta = nullptr;
154  mSumCaloPt_HF = nullptr;
155 
156  for (int ieta = 0; ieta < etaBins_; ieta++) {
157  mSumCaloPtEtaDep[ieta] = nullptr;
158  }
159  }
160 
161  mSumpt = nullptr;
162 
163  // Events variables
164  mNvtx = nullptr;
165  mHF = nullptr;
166 
167  // added Jan 12th 2015
168 
169  // Jet parameters
170  mEta = nullptr;
171  mPhi = nullptr;
172  mEnergy = nullptr;
173  mP = nullptr;
174  mPt = nullptr;
175  mMass = nullptr;
176  mConstituents = nullptr;
177  mJetArea = nullptr;
178  mjetpileup = nullptr;
179  mNJets_40 = nullptr;
180  mNJets = nullptr;
181 }
182 
184  edm::Run const &iRun,
185  edm::EventSetup const &iSetup) {
186  ibooker.setCurrentFolder("JetMET/HIJetValidation/" + mInputCollection.label());
187 
188  TH2F *h2D_etabins_vs_pt2 =
189  new TH2F("h2D_etabins_vs_pt2", ";#eta;sum p_{T}^{2}", etaBins_, edge_pseudorapidity, 10000, 0, 10000);
190  TH2F *h2D_etabins_vs_pt =
191  new TH2F("h2D_etabins_vs_pt", ";#eta;sum p_{T}", etaBins_, edge_pseudorapidity, 500, 0, 500);
192  TH2F *h2D_pfcand_etabins_vs_pt =
193  new TH2F("h2D_etabins_vs_pt", ";#eta;sum p_{T}", etaBins_, edge_pseudorapidity, 300, 0, 300);
194 
195  const int nHihfBins = 100;
196  const double hihfBins[nHihfBins + 1] = {
197  0, 11.282305, 11.82962, 12.344717, 13.029054, 13.698554, 14.36821, 15.140326,
198  15.845786, 16.684441, 17.449186, 18.364939, 19.247023, 20.448898, 21.776642, 22.870239,
199  24.405788, 26.366919, 28.340206, 30.661842, 33.657627, 36.656773, 40.028049, 44.274784,
200  48.583706, 52.981358, 56.860199, 61.559853, 66.663689, 72.768196, 78.265915, 84.744431,
201  92.483459, 100.281021, 108.646576, 117.023911, 125.901093, 135.224899, 147.046875, 159.864258,
202  171.06015, 184.76535, 197.687103, 212.873535, 229.276413, 245.175369, 262.498322, 280.54599,
203  299.570801, 317.188446, 336.99881, 357.960144, 374.725922, 400.638367, 426.062103, 453.07251,
204  483.99704, 517.556396, 549.421143, 578.050781, 608.358643, 640.940979, 680.361755, 719.215027,
205  757.798645, 793.882385, 839.83728, 887.268127, 931.233276, 980.856689, 1023.191833, 1080.281494,
206  1138.363892, 1191.303345, 1251.439453, 1305.288818, 1368.290894, 1433.700684, 1501.597412, 1557.918335,
207  1625.636475, 1695.08374, 1761.771484, 1848.941162, 1938.178345, 2027.55603, 2127.364014, 2226.186523,
208  2315.188965, 2399.225342, 2501.608643, 2611.077881, 2726.316162, 2848.74707, 2972.975342, 3096.565674,
209  3219.530762, 3361.178223, 3568.028564, 3765.690186, 50000};
210 
211  TH2F *h2D_etabins_forRho =
212  new TH2F("etabinsForRho", "#rho vs. #eta;#eta;#rho", etaBins_, edge_pseudorapidity, 500, 0, 300);
213  TH2F *h2D_ptBins_forRho = new TH2F("ptBinsForRho", "#rho vs. p_{T};p_{T};#rho", 300, 0, 300, 500, 0, 300);
214  TH2F *h2D_centBins_forRho = new TH2F("centBinsForRho", "dummy;HIHF;#rho", nHihfBins, hihfBins, 500, 0, 300);
215 
216  TH2F *h2D_etabins_forRhoM =
217  new TH2F("etabinsForRho", "#rho_{M} vs. #eta;#eta;#rho_{M}", etaBins_, edge_pseudorapidity, 100, 0, 1.5);
218  TH2F *h2D_ptBins_forRhoM = new TH2F("ptBinsForRho", "#rho_{M} vs. p_{T};p_{T};#rho_{M}", 300, 0, 300, 100, 0, 1.5);
219  TH2F *h2D_centBins_forRhoM = new TH2F("centBinsForRho", "dummy;HIHF;#rho_{M}", nHihfBins, hihfBins, 100, 0, 1.5);
220 
221  if (isPFJet) {
222  mNPFpart = ibooker.book1D("NPFpart", "No of particle flow candidates", 1000, 0, 1000);
223  mPFPt = ibooker.book1D("PFPt", "PF candidate p_{T}", 10000, -500, 500);
224  mPFEta = ibooker.book1D("PFEta", "PF candidate #eta", 120, -6, 6);
225  mPFPhi = ibooker.book1D("PFPhi", "PF candidate #phi", 70, -3.5, 3.5);
226 
227  mPFDeltaR = ibooker.book1D("PFDeltaR", "PF candidate DeltaR", 100, 0, 4); //MZ
229  ibooker.book1D("PFDeltaR_Scaled_R", "PF candidate DeltaR Divided by DeltaR square", 100, 0, 4); //MZ
230 
231  {
232  auto scope = DQMStore::IBooker::UseLumiScope(ibooker);
233  mSumPFPt = ibooker.book1D("SumPFPt", "Sum of initial PF p_{T}", 1000, -10000, 10000);
234  }
235  mSumPFPt_eta = ibooker.book2D("SumPFPt_etaBins", h2D_etabins_vs_pt);
236 
237  mSumSquaredPFPt = ibooker.book1D("SumSquaredPFPt", "Sum of initial PF p_{T} squared", 10000, 0, 10000);
238  mSumSquaredPFPt_eta = ibooker.book2D("SumSquaredPFPt_etaBins", h2D_etabins_vs_pt2);
239 
240  mSumPFPt_HF = ibooker.book2D(
241  "SumPFPt_HF", "HF energy (y axis) vs Sum initial PF p_{T} (x axis)", 1000, -1000, 1000, 1000, 0, 10000);
242 
243  for (int ieta = 0; ieta < etaBins_; ieta++) {
244  int range = 1000;
245  if (ieta < 2 || etaBins_ - ieta <= 2)
246  range = 500;
247  const char *lc = edge_pseudorapidity[ieta] < 0 ? "n" : "p";
248  const char *rc = edge_pseudorapidity[ieta + 1] < 0 ? "n" : "p";
250  Form("mSumCaloPt_%s%.3g_%s%.3g", lc, abs(edge_pseudorapidity[ieta]), rc, abs(edge_pseudorapidity[ieta + 1]));
251  for (int id = 0; id < 2; id++) {
252  if (histoName.find('.') != std::string::npos) {
253  histoName.replace(histoName.find('.'), 1, "p");
254  }
255  }
256  mSumPFPtEtaDep[ieta] = ibooker.book1D(
257  histoName.c_str(),
258  Form("Sum PFPt in the eta range %.3g to %.3g", edge_pseudorapidity[ieta], edge_pseudorapidity[ieta + 1]),
259  500,
260  0,
261  range);
262  }
263 
264  if (std::string("Cs") == UEAlgo) {
265  mRhoDist_vsEta = ibooker.book2D("rhoDist_vsEta", h2D_etabins_forRho);
266  mRhoMDist_vsEta = ibooker.book2D("rhoMDist_vsEta", h2D_etabins_forRhoM);
267  mRhoDist_vsPt = ibooker.book2D("rhoDist_vsPt", h2D_ptBins_forRho);
268  mRhoMDist_vsPt = ibooker.book2D("rhoMDist_vsPt", h2D_ptBins_forRhoM);
269 
270  //this is kind of a janky way to fill the eta since i can't get it from the edm::Event here... - kjung
271  rhoEtaRange = ibooker.book1D("rhoEtaRange", "", 500, -5.5, 5.5);
272  for (int ieta = 0; ieta < etaBins_; ieta++) {
274  ibooker.book1D(Form("csCandPt_etaBin%d", ieta), "CS candidate pt, eta-by-eta", 150, 0, 300);
275 
276  const char *lc = edge_pseudorapidity[ieta] < 0 ? "n" : "p";
277  const char *rc = edge_pseudorapidity[ieta + 1] < 0 ? "n" : "p";
278  std::string histoName = Form(
279  "Dist_vsCent_%s%.3g_%s%.3g", lc, abs(edge_pseudorapidity[ieta]), rc, abs(edge_pseudorapidity[ieta + 1]));
280  for (int id = 0; id < 2; id++) {
281  if (histoName.find('.') != std::string::npos) {
282  histoName.replace(histoName.find('.'), 1, "p");
283  }
284  }
285  std::string rhoName = "rho";
286  rhoName.append(histoName);
287  h2D_centBins_forRho->SetTitle(Form(
288  "#rho vs. HIHF in the range %.3g < #eta < %.3g", edge_pseudorapidity[ieta], edge_pseudorapidity[ieta + 1]));
289  mRhoDist_vsCent[ieta] = ibooker.book2D(rhoName.c_str(), h2D_centBins_forRho);
290  std::string rhoMName = "rhoM";
291  rhoMName.append(histoName);
292  h2D_centBins_forRhoM->SetTitle(Form("#rho_{M} vs. HIHF in the range %.3g < #eta < %.3g",
294  edge_pseudorapidity[ieta + 1]));
295  mRhoMDist_vsCent[ieta] = ibooker.book2D(rhoMName.c_str(), h2D_centBins_forRhoM);
296  for (int ipt = 0; ipt < ptBins_; ipt++) {
297  mSubtractedEFrac[ipt][ieta] =
298  ibooker.book1D(Form("subtractedEFrac_JetPt%d_to_%d_etaBin%d", ptBin[ipt], ptBin[ipt + 1], ieta),
299  "subtracted fraction of CS jet",
300  50,
301  0,
302  1);
303  mSubtractedE[ipt][ieta] =
304  ibooker.book1D(Form("subtractedE_JetPt%d_to_%d_etaBin%d", ptBin[ipt], ptBin[ipt + 1], ieta),
305  "subtracted total of CS jet",
306  300,
307  0,
308  300);
309  }
310  mCSCand_corrPFcand[ieta] = ibooker.book2D(
311  Form("csCandCorrPF%d", ieta), "CS to PF candidate correlation, eta-by-eta", 300, 0, 300, 300, 0, 300);
312  }
313  }
314 
315  mPFCandpT_vs_eta_Unknown = ibooker.book2D("PF_cand_X_unknown", h2D_pfcand_etabins_vs_pt); // pf id 0
316  mPFCandpT_vs_eta_ChargedHadron = ibooker.book2D("PF_cand_chargedHad", h2D_pfcand_etabins_vs_pt); // pf id - 1
317  mPFCandpT_vs_eta_electron = ibooker.book2D("PF_cand_electron", h2D_pfcand_etabins_vs_pt); // pf id - 2
318  mPFCandpT_vs_eta_muon = ibooker.book2D("PF_cand_muon", h2D_pfcand_etabins_vs_pt); // pf id - 3
319  mPFCandpT_vs_eta_photon = ibooker.book2D("PF_cand_photon", h2D_pfcand_etabins_vs_pt); // pf id - 4
320  mPFCandpT_vs_eta_NeutralHadron = ibooker.book2D("PF_cand_neutralHad", h2D_pfcand_etabins_vs_pt); // pf id - 5
321  mPFCandpT_vs_eta_HadE_inHF = ibooker.book2D("PF_cand_HadEner_inHF", h2D_pfcand_etabins_vs_pt); // pf id - 6
322  mPFCandpT_vs_eta_EME_inHF = ibooker.book2D("PF_cand_EMEner_inHF", h2D_pfcand_etabins_vs_pt); // pf id - 7
323 
324  mPFCandpT_Barrel_Unknown = ibooker.book1D("mPFCandpT_Barrel_Unknown",
325  Form(";PF candidate p_{T}, |#eta|<%2.2f; counts", BarrelEta),
326  300,
327  0,
328  300); // pf id - 0
329  mPFCandpT_Barrel_ChargedHadron = ibooker.book1D("mPFCandpT_Barrel_ChargedHadron",
330  Form(";PF candidate p_{T}, |#eta|<%2.2f; counts", BarrelEta),
331  300,
332  0,
333  300); // pf id - 1
334  mPFCandpT_Barrel_electron = ibooker.book1D("mPFCandpT_Barrel_electron",
335  Form(";PF candidate p_{T}, |#eta|<%2.2f; counts", BarrelEta),
336  300,
337  0,
338  300); // pf id - 2
339  mPFCandpT_Barrel_muon = ibooker.book1D("mPFCandpT_Barrel_muon",
340  Form(";PF candidate p_{T}, |#eta|<%2.2f; counts", BarrelEta),
341  300,
342  0,
343  300); // pf id - 3
344  mPFCandpT_Barrel_photon = ibooker.book1D("mPFCandpT_Barrel_photon",
345  Form(";PF candidate p_{T}, |#eta|<%2.2f; counts", BarrelEta),
346  300,
347  0,
348  300); // pf id - 4
349  mPFCandpT_Barrel_NeutralHadron = ibooker.book1D("mPFCandpT_Barrel_NeutralHadron",
350  Form(";PF candidate p_{T}, |#eta|<%2.2f; counts", BarrelEta),
351  300,
352  0,
353  300); // pf id - 5
354  mPFCandpT_Barrel_HadE_inHF = ibooker.book1D("mPFCandpT_Barrel_HadE_inHF",
355  Form(";PF candidate p_{T}, |#eta|<%2.2f; counts", BarrelEta),
356  300,
357  0,
358  300); // pf id - 6
359  mPFCandpT_Barrel_EME_inHF = ibooker.book1D("mPFCandpT_Barrel_EME_inHF",
360  Form(";PF candidate p_{T}, |#eta|<%2.2f; counts", BarrelEta),
361  300,
362  0,
363  300); // pf id - 7
364 
366  ibooker.book1D("mPFCandpT_Endcap_Unknown",
367  Form(";PF candidate p_{T}, %2.2f<|#eta|<%2.2f; counts", BarrelEta, EndcapEta),
368  300,
369  0,
370  300); // pf id - 0
372  ibooker.book1D("mPFCandpT_Endcap_ChargedHadron",
373  Form(";PF candidate p_{T}, %2.2f<|#eta|<%2.2f; counts", BarrelEta, EndcapEta),
374  300,
375  0,
376  300); // pf id - 1
378  ibooker.book1D("mPFCandpT_Endcap_electron",
379  Form(";PF candidate p_{T}, %2.2f<|#eta|<%2.2f; counts", BarrelEta, EndcapEta),
380  300,
381  0,
382  300); // pf id - 2
384  ibooker.book1D("mPFCandpT_Endcap_muon",
385  Form(";PF candidate p_{T}, %2.2f<|#eta|<%2.2f; counts", BarrelEta, EndcapEta),
386  300,
387  0,
388  300); // pf id - 3
390  ibooker.book1D("mPFCandpT_Endcap_photon",
391  Form(";PF candidate p_{T}, %2.2f<|#eta|<%2.2f; counts", BarrelEta, EndcapEta),
392  300,
393  0,
394  300); // pf id - 4
396  ibooker.book1D("mPFCandpT_Endcap_NeutralHadron",
397  Form(";PF candidate p_{T}, %2.2f<|#eta|<%2.2f; counts", BarrelEta, EndcapEta),
398  300,
399  0,
400  300); // pf id - 5
402  ibooker.book1D("mPFCandpT_Endcap_HadE_inHF",
403  Form(";PF candidate p_{T}, %2.2f<|#eta|<%2.2f; counts", BarrelEta, EndcapEta),
404  300,
405  0,
406  300); // pf id - 6
408  ibooker.book1D("mPFCandpT_Endcap_EME_inHF",
409  Form(";PF candidate p_{T}, %2.2f<|#eta|<%2.2f; counts", BarrelEta, EndcapEta),
410  300,
411  0,
412  300); // pf id - 7
413 
415  ibooker.book1D("mPFCandpT_Forward_Unknown",
416  Form(";PF candidate p_{T}, %2.2f<|#eta|<%2.2f; counts", EndcapEta, ForwardEta),
417  300,
418  0,
419  300); // pf id - 0
421  ibooker.book1D("mPFCandpT_Forward_ChargedHadron",
422  Form(";PF candidate p_{T}, %2.2f<|#eta|<%2.2f; counts", EndcapEta, ForwardEta),
423  300,
424  0,
425  300); // pf id - 1
427  ibooker.book1D("mPFCandpT_Forward_electron",
428  Form(";PF candidate p_{T}, %2.2f<|#eta|<%2.2f; counts", EndcapEta, ForwardEta),
429  300,
430  0,
431  300); // pf id - 2
433  ibooker.book1D("mPFCandpT_Forward_muon",
434  Form(";PF candidate p_{T}, %2.2f<|#eta|<%2.2f; counts", EndcapEta, ForwardEta),
435  300,
436  0,
437  300); // pf id - 3
439  ibooker.book1D("mPFCandpT_Forward_photon",
440  Form(";PF candidate p_{T}, %2.2f<|#eta|<%2.2f; counts", EndcapEta, ForwardEta),
441  300,
442  0,
443  300); // pf id - 4
445  ibooker.book1D("mPFCandpT_Forward_NeutralHadron",
446  Form(";PF candidate p_{T}, %2.2f<|#eta|<%2.2f; counts", EndcapEta, ForwardEta),
447  300,
448  0,
449  300); // pf id - 5
451  ibooker.book1D("mPFCandpT_Forward_HadE_inHF",
452  Form(";PF candidate p_{T}, %2.2f<|#eta|<%2.2f; counts", EndcapEta, ForwardEta),
453  300,
454  0,
455  300); // pf id - 6
457  ibooker.book1D("mPFCandpT_Forward_EME_inHF",
458  Form(";PF candidate p_{T}, %2.2f<|#eta|<%2.2f; counts", EndcapEta, ForwardEta),
459  300,
460  0,
461  300); // pf id - 7
462  }
463 
464  if (isCaloJet) {
465  mNCalopart = ibooker.book1D("NCalopart", "No of particle flow candidates", 1000, 0, 10000);
466  mCaloPt = ibooker.book1D("CaloPt", "Calo candidate p_{T}", 1000, -5000, 5000);
467  mCaloEta = ibooker.book1D("CaloEta", "Calo candidate #eta", 120, -6, 6);
468  mCaloPhi = ibooker.book1D("CaloPhi", "Calo candidate #phi", 70, -3.5, 3.5);
469 
470  mSumCaloPt = ibooker.book1D("SumCaloPt", "Sum Calo p_{T}", 1000, -10000, 10000);
471  mSumCaloPt_eta = ibooker.book2D("SumCaloPt_etaBins", h2D_etabins_vs_pt);
472 
473  mSumSquaredCaloPt = ibooker.book1D("SumSquaredCaloPt", "Sum of initial Calo tower p_{T} squared", 10000, 0, 10000);
474  mSumSquaredCaloPt_eta = ibooker.book2D("SumSquaredCaloPt_etaBins", h2D_etabins_vs_pt2);
475 
476  mSumCaloPt_HF =
477  ibooker.book2D("SumCaloPt_HF", "HF Energy (y axis) vs Sum Calo tower p_{T}", 1000, -1000, 1000, 1000, 0, 10000);
478 
479  for (int ieta = 0; ieta < etaBins_; ieta++) {
480  int range = 1000;
481  if (ieta < 2 || etaBins_ - ieta <= 2)
482  range = 5000;
483  const char *lc = edge_pseudorapidity[ieta] < 0 ? "n" : "p";
484  const char *rc = edge_pseudorapidity[ieta + 1] < 0 ? "n" : "p";
486  Form("mSumCaloPt_%s%.3g_%s%.3g", lc, abs(edge_pseudorapidity[ieta]), rc, abs(edge_pseudorapidity[ieta + 1]));
487  for (int id = 0; id < 2; id++) {
488  if (histoName.find('.') != std::string::npos) {
489  histoName.replace(histoName.find('.'), 1, "p");
490  }
491  }
492  mSumCaloPtEtaDep[ieta] = ibooker.book1D(histoName.c_str(),
493  Form("Sum Calo tower Pt in the eta range %.3g to %.3g",
496  1000,
497  -1 * range,
498  range);
499  }
500  }
501 
502  // particle flow variables histograms
503  mSumpt = ibooker.book1D("SumpT", "Sum p_{T} of all the PF candidates per event", 1000, 0, 10000);
504 
505  // Event variables
506  mNvtx = ibooker.book1D("Nvtx", "number of vertices", 60, 0, 60);
507  mHF = ibooker.book1D("HF", "HF energy distribution", 1000, 0, 10000);
508 
509  // Jet parameters
510  mEta = ibooker.book1D("Eta", "Eta", 120, -6, 6);
511  mPhi = ibooker.book1D("Phi", "Phi", 70, -3.5, 3.5);
512  mPt = ibooker.book1D("Pt", "Pt", 1000, 0, 500);
513  mP = ibooker.book1D("P", "P", 100, 0, 1000);
514  mEnergy = ibooker.book1D("Energy", "Energy", 100, 0, 1000);
515  mMass = ibooker.book1D("Mass", "Mass", 100, 0, 200);
516  mConstituents = ibooker.book1D("Constituents", "Constituents", 100, 0, 100);
517  mJetArea = ibooker.book1D("JetArea", "JetArea", 100, 0, 4);
518  mjetpileup = ibooker.book1D("jetPileUp", "jetPileUp", 100, 0, 150);
519  mNJets_40 = ibooker.book1D("NJets_pt_greater_40", "NJets pT > 40 GeV", 50, 0, 50);
520  mNJets = ibooker.book1D("NJets", "NJets", 100, 0, 100);
521 
522  if (mOutputFile.empty())
523  LogInfo("OutputInfo") << " Histograms will NOT be saved";
524  else
525  LogInfo("OutputInfo") << " Histograms will be saved to file:" << mOutputFile;
526 
527  delete h2D_etabins_vs_pt2;
528  delete h2D_etabins_vs_pt;
529  delete h2D_pfcand_etabins_vs_pt;
530  delete h2D_etabins_forRho;
531  delete h2D_ptBins_forRho;
532  delete h2D_centBins_forRho;
533  delete h2D_centBins_forRhoM;
534 }
535 
536 //------------------------------------------------------------------------------
537 // ~JetAnalyzer_HeavyIons
538 //------------------------------------------------------------------------------
540 
541 //------------------------------------------------------------------------------
542 // analyze
543 //------------------------------------------------------------------------------
544 void JetAnalyzer_HeavyIons::analyze(const edm::Event &mEvent, const edm::EventSetup &mSetup) {
545  // switch(mEvent.id().event() == 15296770)
546  // case 1:
547  // break;
548 
549  // Get the primary vertices
551  mEvent.getByToken(pvToken_, pvHandle);
552  reco::Vertex::Point vtx(0, 0, 0);
554 
555  mEvent.getByToken(hiVertexToken_, vtxs);
556  int greatestvtx = 0;
557  int nVertex = vtxs->size();
558 
559  for (unsigned int i = 0; i < vtxs->size(); ++i) {
560  unsigned int daughter = (*vtxs)[i].tracksSize();
561  if (daughter > (*vtxs)[greatestvtx].tracksSize())
562  greatestvtx = i;
563  }
564 
565  if (nVertex <= 0) {
566  vtx = reco::Vertex::Point(0, 0, 0);
567  }
568  vtx = (*vtxs)[greatestvtx].position();
569 
570  int nGoodVertices = 0;
571 
572  if (pvHandle.isValid()) {
573  for (unsigned i = 0; i < pvHandle->size(); i++) {
574  if ((*pvHandle)[i].ndof() > 4 && (fabs((*pvHandle)[i].z()) <= 24) && (fabs((*pvHandle)[i].position().rho()) <= 2))
575  nGoodVertices++;
576  }
577  }
578 
579  mNvtx->Fill(nGoodVertices);
580 
581  // Get the Jet collection
582  //----------------------------------------------------------------------------
583 
584  std::vector<Jet> recoJets;
585  recoJets.clear();
586 
591 
592  // Get the Particle flow candidates and the Voronoi variables
595  edm::Handle<CaloTowerCollection> caloCandidates;
596  edm::Handle<reco::CandidateView> pfcandidates_;
597  edm::Handle<reco::CandidateView> calocandidates_;
598 
599  //Get the new CS stuff
603  if (std::string("Cs") == UEAlgo) {
604  mEvent.getByToken(etaToken_, etaRanges);
605  mEvent.getByToken(rhoToken_, rho);
606  mEvent.getByToken(rhomToken_, rhom);
607  const int rhoSize = (int)etaRanges->size();
608  double rhoRange[rhoSize];
609  for (int irho = 0; irho < rhoSize; irho++) {
610  rhoRange[irho] = etaRanges->at(irho);
611  }
612  double yaxisLimits[501];
613  for (int ibin = 0; ibin < 501; ibin++)
614  yaxisLimits[ibin] = ibin * 2;
615  if (mRhoDist_vsEta->getNbinsX() != rhoSize - 1) {
616  mRhoDist_vsEta->getTH2F()->SetBins(
617  rhoSize - 1, const_cast<double *>(rhoRange), 500, const_cast<double *>(yaxisLimits));
618  mRhoMDist_vsEta->getTH2F()->SetBins(
619  rhoSize - 1, const_cast<double *>(rhoRange), 500, const_cast<double *>(yaxisLimits));
620  }
621  }
622 
623  if (isCaloJet)
625  if (isJPTJet)
626  mEvent.getByToken(jptJetsToken_, jptJets);
627  if (isPFJet) {
628  if (std::string("Pu") == UEAlgo)
629  mEvent.getByToken(basicJetsToken_, basicJets);
630  if (std::string("Cs") == UEAlgo)
631  mEvent.getByToken(pfJetsToken_, pfJets);
632  if (std::string("Vs") == UEAlgo)
633  return; //avoid running Vs jets
634  }
635 
637  mEvent.getByToken(pfCandViewToken_, pfcandidates_);
638  if (std::string("Cs") == UEAlgo)
639  mEvent.getByToken(csCandToken_, csCandidates);
640 
641  mEvent.getByToken(caloTowersToken_, caloCandidates);
642  mEvent.getByToken(caloCandViewToken_, calocandidates_);
643 
644  // get the centrality
646  mEvent.getByToken(centralityToken, cent); //_centralitytag comes from the cfg
647 
648  mHF->Fill(cent->EtHFtowerSum());
649  Float_t HF_energy = cent->EtHFtowerSum();
650 
651  //for later when centrality gets added to RelVal
652  //edm::Handle<int> cbin;
653  //mEvent.getByToken(centralityBinToken, cbin);
654 
655  if (!cent.isValid())
656  return;
657 
658  /*int hibin = -999;
659  if(cbin.isValid()){
660  hibin = *cbin;
661  }*/
662 
663  const reco::PFCandidateCollection *pfCandidateColl = pfCandidates.product();
664 
665  Int_t NPFpart = 0;
666  Int_t NCaloTower = 0;
667  Float_t pfPt = 0;
668  Float_t pfEta = 0;
669  Float_t pfPhi = 0;
670  Int_t pfID = 0;
671  Float_t pfDeltaR = 0;
672  Float_t caloPt = 0;
673  Float_t caloEta = 0;
674  Float_t caloPhi = 0;
675  Float_t SumPt_value = 0;
676 
677  vector<vector<float>> numbers;
678  vector<float> tempVector;
679  numbers.clear();
680  tempVector.clear();
681 
682  if (isCaloJet) {
683  Float_t SumCaloPt[etaBins_];
684  Float_t SumSquaredCaloPt[etaBins_];
685 
686  // Need to set up histograms to get the RMS values for each pT bin
687  TH1F *hSumCaloPt[nedge_pseudorapidity - 1];
688 
689  for (int i = 0; i < etaBins_; ++i) {
690  SumCaloPt[i] = 0;
691  SumSquaredCaloPt[i] = 0;
692  hSumCaloPt[i] = new TH1F(Form("hSumCaloPt_%d", i), "", 10000, -10000, 10000);
693  }
694 
695  for (unsigned icand = 0; icand < caloCandidates->size(); icand++) {
696  const CaloTower &tower = (*caloCandidates)[icand];
697  reco::CandidateViewRef ref(calocandidates_, icand);
698  //10 is tower pT min
699  if (tower.p4(vtx).Et() < 0.1)
700  continue;
701 
702  NCaloTower++;
703 
704  caloPt = tower.p4(vtx).Et();
705  caloEta = tower.p4(vtx).Eta();
706  caloPhi = tower.p4(vtx).Phi();
707 
708  for (size_t k = 0; k < nedge_pseudorapidity - 1; k++) {
709  if (caloEta >= edge_pseudorapidity[k] && caloEta < edge_pseudorapidity[k + 1]) {
710  SumCaloPt[k] = SumCaloPt[k] + caloPt;
711  SumSquaredCaloPt[k] = SumSquaredCaloPt[k] + caloPt * caloPt;
712  break;
713  } // eta selection statement
714 
715  } // eta bin loop
716 
717  SumPt_value = SumPt_value + caloPt;
718 
719  mCaloPt->Fill(caloPt);
720  mCaloEta->Fill(caloEta);
721  mCaloPhi->Fill(caloPhi);
722 
723  } // calo tower candidate loop
724 
725  for (int k = 0; k < nedge_pseudorapidity - 1; k++) {
726  hSumCaloPt[k]->Fill(SumCaloPt[k]);
727 
728  } // eta bin loop
729 
730  Float_t Evt_SumCaloPt = 0;
731  Float_t Evt_SumSquaredCaloPt = 0;
732 
733  for (int ieta = 0; ieta < etaBins_; ieta++) {
734  mSumCaloPtEtaDep[ieta]->Fill(SumCaloPt[ieta]);
735  }
736 
737  for (size_t k = 0; k < nedge_pseudorapidity - 1; k++) {
738  Evt_SumCaloPt = Evt_SumCaloPt + SumCaloPt[k];
739  mSumCaloPt_eta->Fill(edge_pseudorapidity[k], SumCaloPt[k]);
740 
741  Evt_SumSquaredCaloPt = Evt_SumSquaredCaloPt + SumSquaredCaloPt[k];
742  mSumSquaredCaloPt_eta->Fill(edge_pseudorapidity[k], hSumCaloPt[k]->GetRMS(1));
743 
744  delete hSumCaloPt[k];
745 
746  } // eta bin loop
747 
748  mSumCaloPt->Fill(Evt_SumCaloPt);
749  mSumCaloPt_HF->Fill(Evt_SumCaloPt, HF_energy);
750 
751  mSumSquaredCaloPt->Fill(Evt_SumSquaredCaloPt);
752 
753  mNCalopart->Fill(NCaloTower);
754  mSumpt->Fill(SumPt_value);
755 
756  } // is calo jet
757 
758  if (isPFJet) {
759  Float_t SumPFPt[etaBins_];
760 
761  Float_t SumSquaredPFPt[etaBins_];
762 
763  // Need to set up histograms to get the RMS values for each pT bin
764  TH1F *hSumPFPt[nedge_pseudorapidity - 1];
765 
766  for (int i = 0; i < etaBins_; i++) {
767  SumPFPt[i] = 0;
768  SumSquaredPFPt[i] = 0;
769 
770  hSumPFPt[i] = new TH1F(Form("hSumPFPt_%d", i), "", 10000, -10000, 10000);
771  }
772 
773  vector<vector<float>> PF_Space(1, vector<float>(3));
774 
775  if (std::string("Cs") == UEAlgo) {
776  const reco::PFCandidateCollection *csCandidateColl = csCandidates.product();
777 
778  for (unsigned iCScand = 0; iCScand < csCandidateColl->size(); iCScand++) {
779  assert(csCandidateColl->size() <= pfCandidateColl->size());
780  const reco::PFCandidate csCandidate = csCandidateColl->at(iCScand);
781  const reco::PFCandidate pfCandidate = pfCandidateColl->at(iCScand);
782  int ieta = 0;
783  while (csCandidate.eta() > edge_pseudorapidity[ieta] && ieta < etaBins_ - 1)
784  ieta++;
785  mCSCandpT_vsPt[ieta]->Fill(csCandidate.pt());
786  mCSCand_corrPFcand[ieta]->Fill(csCandidate.pt(), pfCandidate.pt());
787  }
788  }
789 
790  for (unsigned icand = 0; icand < pfCandidateColl->size(); icand++) {
791  const reco::PFCandidate pfCandidate = pfCandidateColl->at(icand);
792  reco::CandidateViewRef ref(pfcandidates_, icand);
793  if (pfCandidate.pt() < 5)
794  continue;
795 
796  NPFpart++;
797  pfPt = pfCandidate.pt();
798  pfEta = pfCandidate.eta();
799  pfPhi = pfCandidate.phi();
800  pfID = pfCandidate.particleId();
801 
802  bool isBarrel = false;
803  bool isEndcap = false;
804  bool isForward = false;
805 
806  if (fabs(pfEta) < BarrelEta)
807  isBarrel = true;
808  if (fabs(pfEta) >= BarrelEta && fabs(pfEta) < EndcapEta)
809  isEndcap = true;
810  if (fabs(pfEta) >= EndcapEta && fabs(pfEta) < ForwardEta)
811  isForward = true;
812 
813  switch (pfID) {
814  case 0:
815  mPFCandpT_vs_eta_Unknown->Fill(pfPt, pfEta);
816  if (isBarrel)
818  if (isEndcap)
820  if (isForward)
822  break;
823  case 1:
824  mPFCandpT_vs_eta_ChargedHadron->Fill(pfPt, pfEta);
825  if (isBarrel)
827  if (isEndcap)
829  if (isForward)
831  break;
832  case 2:
833  mPFCandpT_vs_eta_electron->Fill(pfPt, pfEta);
834  if (isBarrel)
836  if (isEndcap)
838  if (isForward)
840  break;
841  case 3:
842  mPFCandpT_vs_eta_muon->Fill(pfPt, pfEta);
843  if (isBarrel)
845  if (isEndcap)
847  if (isForward)
849  break;
850  case 4:
851  mPFCandpT_vs_eta_photon->Fill(pfPt, pfEta);
852  if (isBarrel)
854  if (isEndcap)
856  if (isForward)
858  break;
859  case 5:
860  mPFCandpT_vs_eta_NeutralHadron->Fill(pfPt, pfEta);
861  if (isBarrel)
863  if (isEndcap)
865  if (isForward)
867  break;
868  case 6:
869  mPFCandpT_vs_eta_HadE_inHF->Fill(pfPt, pfEta);
870  if (isBarrel)
872  if (isEndcap)
874  if (isForward)
876  break;
877  case 7:
878  mPFCandpT_vs_eta_EME_inHF->Fill(pfPt, pfEta);
879  if (isBarrel)
881  if (isEndcap)
883  if (isForward)
885  break;
886  }
887 
888  //Fill 2d vector matrix
889  tempVector.push_back(pfPt);
890  tempVector.push_back(pfEta);
891  tempVector.push_back(pfPhi);
892 
893  numbers.push_back(tempVector);
894  tempVector.clear();
895 
896  for (size_t k = 0; k < nedge_pseudorapidity - 1; k++) {
897  if (pfEta >= edge_pseudorapidity[k] && pfEta < edge_pseudorapidity[k + 1]) {
898  SumPFPt[k] = SumPFPt[k] + pfPt;
899  SumSquaredPFPt[k] = SumSquaredPFPt[k] + pfPt * pfPt;
900  break;
901  } // eta selection statement
902 
903  } // eta bin loop
904 
905  SumPt_value = SumPt_value + pfPt;
906 
907  mPFPt->Fill(pfPt);
908  mPFEta->Fill(pfEta);
909  mPFPhi->Fill(pfPhi);
910 
911  } // pf candidate loop
912 
913  for (int k = 0; k < nedge_pseudorapidity - 1; k++) {
914  hSumPFPt[k]->Fill(SumPFPt[k]);
915 
916  } // eta bin loop
917 
918  Float_t Evt_SumPFPt = 0;
919  Float_t Evt_SumSquaredPFPt = 0;
920 
921  for (int ieta = 0; ieta < etaBins_; ieta++) {
922  mSumPFPtEtaDep[ieta]->Fill(SumPFPt[ieta]);
923  }
924 
925  for (size_t k = 0; k < nedge_pseudorapidity - 1; k++) {
926  Evt_SumPFPt = Evt_SumPFPt + SumPFPt[k];
927  mSumPFPt_eta->Fill(edge_pseudorapidity[k], SumPFPt[k]);
928 
929  Evt_SumSquaredPFPt = Evt_SumSquaredPFPt + SumSquaredPFPt[k];
930  mSumSquaredPFPt_eta->Fill(edge_pseudorapidity[k], hSumPFPt[k]->GetRMS(1));
931 
932  delete hSumPFPt[k];
933 
934  } // eta bin loop
935 
936  mSumPFPt->Fill(Evt_SumPFPt);
937  mSumPFPt_HF->Fill(Evt_SumPFPt, HF_energy);
938 
939  mSumSquaredPFPt->Fill(Evt_SumSquaredPFPt);
940 
941  mNPFpart->Fill(NPFpart);
942  mSumpt->Fill(SumPt_value);
943  }
944 
945  if (isCaloJet) {
946  for (unsigned ijet = 0; ijet < caloJets->size(); ijet++) {
947  recoJets.push_back((*caloJets)[ijet]);
948  }
949  }
950 
951  if (isJPTJet) {
952  for (unsigned ijet = 0; ijet < jptJets->size(); ijet++)
953  recoJets.push_back((*jptJets)[ijet]);
954  }
955 
956  if (isPFJet) {
957  if (std::string("Pu") == UEAlgo) {
958  for (unsigned ijet = 0; ijet < basicJets->size(); ijet++) {
959  recoJets.push_back((*basicJets)[ijet]);
960  }
961  }
962  if (std::string("Cs") == UEAlgo) {
963  for (unsigned ijet = 0; ijet < pfJets->size(); ijet++) {
964  recoJets.push_back((*pfJets)[ijet]);
965  }
966  }
967  }
968 
969  if (isCaloJet && !caloJets.isValid()) {
970  return;
971  }
972  if (isJPTJet && !jptJets.isValid()) {
973  return;
974  }
975  if (isPFJet) {
976  if (std::string("Pu") == UEAlgo) {
977  if (!basicJets.isValid())
978  return;
979  }
980  if (std::string("Cs") == UEAlgo) {
981  if (!pfJets.isValid())
982  return;
983  }
984  if (std::string("Vs") == UEAlgo) {
985  return;
986  }
987  }
988 
989  int nJet_40 = 0;
990 
991  mNJets->Fill(recoJets.size());
992 
993  for (unsigned ijet = 0; ijet < recoJets.size(); ijet++) {
994  if (recoJets[ijet].pt() > mRecoJetPtThreshold) {
995  //counting forward and barrel jets
996  // get an idea of no of jets with pT>40 GeV
997  if (recoJets[ijet].pt() > 40)
998  nJet_40++;
999  if (mEta)
1000  mEta->Fill(recoJets[ijet].eta());
1001  if (mjetpileup)
1002  mjetpileup->Fill(recoJets[ijet].pileup());
1003  if (mJetArea)
1004  mJetArea->Fill(recoJets[ijet].jetArea());
1005  if (mPhi)
1006  mPhi->Fill(recoJets[ijet].phi());
1007  if (mEnergy)
1008  mEnergy->Fill(recoJets[ijet].energy());
1009  if (mP)
1010  mP->Fill(recoJets[ijet].p());
1011  if (mPt)
1012  mPt->Fill(recoJets[ijet].pt());
1013  if (mMass)
1014  mMass->Fill(recoJets[ijet].mass());
1015  if (mConstituents)
1017 
1018  if (std::string("Cs") == UEAlgo) {
1019  int ipt = 0, ieta = 0;
1020  while (recoJets[ijet].pt() > ptBin[ipt + 1] && ipt < ptBins_ - 1)
1021  ipt++;
1022  while (recoJets[ijet].eta() > etaRanges->at(ieta + 1) && ieta < (int)(rho->size() - 1))
1023  ieta++;
1024  mSubtractedEFrac[ipt][ieta]->Fill((double)recoJets[ijet].pileup() / (double)recoJets[ijet].energy());
1025  mSubtractedE[ipt][ieta]->Fill(recoJets[ijet].pileup());
1026 
1027  for (unsigned irho = 0; irho < rho->size(); irho++) {
1028  mRhoDist_vsEta->Fill(recoJets[ijet].eta(), rho->at(irho));
1029  mRhoMDist_vsEta->Fill(recoJets[ijet].eta(), rhom->at(irho));
1030  mRhoDist_vsPt->Fill(recoJets[ijet].pt(), rho->at(irho));
1031  mRhoMDist_vsPt->Fill(recoJets[ijet].pt(), rhom->at(irho));
1032  mRhoDist_vsCent[ieta]->Fill(HF_energy, rho->at(irho));
1033  mRhoMDist_vsCent[ieta]->Fill(HF_energy, rhom->at(irho));
1034  }
1035  }
1036 
1037  for (size_t iii = 0; iii < numbers.size(); iii++) {
1038  pfDeltaR = sqrt((numbers[iii][2] - recoJets[ijet].phi()) * (numbers[iii][2] - recoJets[ijet].phi()) +
1039  (numbers[iii][1] - recoJets[ijet].eta()) * (numbers[iii][1] - recoJets[ijet].eta())); //MZ
1040 
1041  mPFDeltaR->Fill(pfDeltaR); //MZ
1042  mPFDeltaR_Scaled_R->Fill(pfDeltaR, 1. / pow(pfDeltaR, 2)); //MZ
1043  }
1044  }
1045  }
1046  if (mNJets_40)
1047  mNJets_40->Fill(nJet_40);
1048 
1049  numbers.clear();
1050 }
MonitorElement * mPFCandpT_Forward_photon
edm::EDGetTokenT< reco::CaloJetCollection > caloJetsToken_
edm::EDGetTokenT< reco::PFCandidateCollection > pfCandToken_
edm::InputTag mInputCsCandCollection
T getParameter(std::string const &) const
Definition: ParameterSet.h:303
MonitorElement * mSumSquaredPFPt
MonitorElement * mPFCandpT_vs_eta_electron
edm::EDGetTokenT< reco::CandidateView > caloCandViewToken_
MonitorElement * mPFCandpT_Forward_EME_inHF
const Double_t EndcapEta
edm::EDGetTokenT< reco::JPTJetCollection > jptJetsToken_
double pt() const final
transverse momentum
double EtHFtowerSum() const
Definition: Centrality.h:24
edm::EDGetTokenT< CaloTowerCollection > caloTowersToken_
virtual void setCurrentFolder(std::string const &fullpath)
Definition: DQMStore.cc:36
MonitorElement * mPFCandpT_Forward_muon
size_type size() const
edm::EDGetTokenT< reco::PFJetCollection > pfJetsToken_
MonitorElement * mCSCandpT_vsPt[etaBins_]
MonitorElement * mPFCandpT_Endcap_NeutralHadron
bool exists(std::string const &parameterName) const
checks if a parameter exists
T const * product() const
Definition: Handle.h:70
MonitorElement * mPFCandpT_vs_eta_EME_inHF
MonitorElement * mSumSquaredCaloPt
bool isForward(DetId const &)
MonitorElement * mSubtractedEFrac[ptBins_][etaBins_]
void analyze(const edm::Event &, const edm::EventSetup &) override
constexpr int pow(int x)
Definition: conifer.h:24
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:536
MonitorElement * mSumPFPtEtaDep[etaBins_]
std::string const & label() const
Definition: InputTag.h:36
MonitorElement * mPFCandpT_Forward_electron
MonitorElement * mRhoMDist_vsEta
edm::EDGetTokenT< std::vector< reco::Vertex > > pvToken_
assert(be >=bs)
MonitorElement * mSumCaloPt_eta
MonitorElement * mPFCandpT_vs_eta_muon
MonitorElement * mCSCand_corrPFcand[etaBins_]
void Fill(long long x)
MonitorElement * mPFCandpT_Endcap_HadE_inHF
edm::EDGetTokenT< int > centralityBinToken
MonitorElement * mPFCandpT_Endcap_ChargedHadron
MonitorElement * mRhoMDist_vsPt
MonitorElement * mPFCandpT_vs_eta_NeutralHadron
MonitorElement * mSumCaloPtEtaDep[etaBins_]
MonitorElement * mSumCaloPt_HF
static constexpr int etaBins_
MonitorElement * mRhoDist_vsEta
MonitorElement * mPFCandpT_vs_eta_photon
edm::EDGetTokenT< reco::Centrality > centralityToken
MonitorElement * mPFCandpT_Forward_ChargedHadron
MonitorElement * mSubtractedE[ptBins_][etaBins_]
T sqrt(T t)
Definition: SSEVec.h:19
MonitorElement * mPFCandpT_Endcap_muon
MonitorElement * mPFCandpT_Endcap_photon
MonitorElement * mPFDeltaR_Scaled_R
edm::EDGetTokenT< std::vector< reco::Vertex > > hiVertexToken_
static constexpr int ptBins_
MonitorElement * mRhoDist_vsPt
MonitorElement * mConstituents
MonitorElement * mPFCandpT_Barrel_photon
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
math::XYZPoint Point
point in the space
Definition: Vertex.h:40
static constexpr int nedge_pseudorapidity
MonitorElement * mSumSquaredCaloPt_eta
MonitorElement * mPFCandpT_Barrel_Unknown
UseScope< MonitorElementData::Scope::LUMI > UseLumiScope
Definition: DQMStore.h:540
edm::EDGetTokenT< std::vector< double > > rhoToken_
edm::EDGetTokenT< std::vector< double > > rhomToken_
MonitorElement * mPFCandpT_Barrel_NeutralHadron
bool isEndcap(GeomDetEnumerators::SubDetector m)
MonitorElement * mPFCandpT_Endcap_Unknown
const Double_t ForwardEta
Log< level::Info, false > LogInfo
MonitorElement * mPFCandpT_vs_eta_ChargedHadron
MonitorElement * mPFCandpT_Forward_HadE_inHF
edm::EDGetTokenT< reco::BasicJetCollection > basicJetsToken_
MonitorElement * mPFCandpT_Barrel_ChargedHadron
MonitorElement * mPFCandpT_Barrel_HadE_inHF
std::vector< reco::PFCandidate > PFCandidateCollection
collection of PFCandidates
edm::EDGetTokenT< reco::PFCandidateCollection > csCandToken_
const double edge_pseudorapidity[etaBins_+1]
edm::InputTag mInputPFCandCollection
MonitorElement * mSumSquaredPFPt_eta
MonitorElement * mPFCandpT_Endcap_EME_inHF
MonitorElement * mPFCandpT_vs_eta_HadE_inHF
MonitorElement * book2D(TString const &name, TString const &title, int nchX, double lowX, double highX, int nchY, double lowY, double highY, FUNC onbooking=NOOP())
Definition: DQMStore.h:212
edm::EDGetTokenT< reco::CandidateView > pfCandViewToken_
void bookHistograms(DQMStore::IBooker &, edm::Run const &, edm::EventSetup const &) override
MonitorElement * mPFCandpT_Forward_Unknown
bool isValid() const
Definition: HandleBase.h:70
Particle reconstructed by the particle flow algorithm.
Definition: PFCandidate.h:41
fixed size matrix
HLT enums.
MonitorElement * mPFCandpT_Barrel_electron
static int position[264][3]
Definition: ReadPGInfo.cc:289
MonitorElement * mPFCandpT_Endcap_electron
edm::EDGetTokenT< std::vector< double > > etaToken_
virtual int getNbinsX() const
get # of bins in X-axis
const Double_t BarrelEta
MonitorElement * mRhoDist_vsCent[etaBins_]
MonitorElement * book1D(TString const &name, TString const &title, int const nchX, double const lowX, double const highX, FUNC onbooking=NOOP())
Definition: DQMStore.h:98
JetAnalyzer_HeavyIons(const edm::ParameterSet &)
const int ptBin[ptBins_+1]
MonitorElement * mPFCandpT_Barrel_EME_inHF
MonitorElement * mRhoMDist_vsCent[etaBins_]
double phi() const final
momentum azimuthal angle
MonitorElement * mPFCandpT_vs_eta_Unknown
MonitorElement * mPFCandpT_Forward_NeutralHadron
Definition: Run.h:45
virtual ParticleType particleId() const
Definition: PFCandidate.h:392
double eta() const final
momentum pseudorapidity
MonitorElement * mPFCandpT_Barrel_muon