CMS 3D CMS Logo
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Groups Pages
MuonTestSummary.cc
Go to the documentation of this file.
1 
2 /*
3  * See header file for a description of this class.
4  *
5  * \author G. Mila - INFN Torino
6  */
7 
8 #include <DQMOffline/Muon/interface/MuonTestSummary.h>
9 
10 // Framework
11 #include <FWCore/Framework/interface/Event.h>
12 #include <FWCore/Framework/interface/EventSetup.h>
13 #include <FWCore/ParameterSet/interface/ParameterSet.h>
14 
15 #include "DQMServices/Core/interface/DQMStore.h"
16 #include "FWCore/MessageLogger/interface/MessageLogger.h"
17 #include "FWCore/Framework/interface/Run.h"
18 
19 #include "DQMOffline/Muon/test/langauFit.C"
20 #include <string>
21 
22 using namespace edm;
23 using namespace std;
24 
25 MuonTestSummary::MuonTestSummary(const edm::ParameterSet &ps) {
26  // parameter initialization for kinematics test
27  etaExpected = ps.getParameter<double>("etaExpected");
28  phiExpected = ps.getParameter<double>("phiExpected");
29  chi2Fraction = ps.getParameter<double>("chi2Fraction");
30  chi2Spread = ps.getParameter<double>("chi2Spread");
31  resEtaSpread_tkGlb = ps.getParameter<double>("resEtaSpread_tkGlb");
32  resEtaSpread_glbSta = ps.getParameter<double>("resEtaSpread_glbSta");
33  resPhiSpread_tkGlb = ps.getParameter<double>("resPhiSpread_tkGlb");
34  resPhiSpread_glbSta = ps.getParameter<double>("resPhiSpread_glbSta");
35  resOneOvPSpread_tkGlb = ps.getParameter<double>("resOneOvPSpread_tkGlb");
36  resOneOvPSpread_glbSta = ps.getParameter<double>("resOneOvPSpread_glbSta");
37  pullEtaSpread = ps.getParameter<double>("pullEtaSpread");
38  pullPhiSpread = ps.getParameter<double>("pullPhiSpread");
39  pullOneOvPSpread = ps.getParameter<double>("pullOneOvPSpread");
40  resChargeLimit_tkGlb = ps.getParameter<double>("resChargeLimit_tkGlb");
41  resChargeLimit_glbSta = ps.getParameter<double>("resChargeLimit_glbSta");
42  resChargeLimit_tkSta = ps.getParameter<double>("resChargeLimit_tkSta");
43  numMatchedExpected_min = ps.getParameter<double>("numMatchedExpected_min");
44  numMatchedExpected_max = ps.getParameter<double>("numMatchedExpected_max");
45  matchesFractionDt_min = ps.getParameter<double>("matchesFractionDt_min");
46  matchesFractionDt_max = ps.getParameter<double>("matchesFractionDt_max");
47  matchesFractionCsc_min = ps.getParameter<double>("matchesFractionCsc_min");
48  matchesFractionCsc_max = ps.getParameter<double>("matchesFractionCsc_max");
49  resSegmTrack_rms_min = ps.getParameter<double>("resSegmTrack_rms_min");
50  resSegmTrack_rms_max = ps.getParameter<double>("resSegmTrack_rms_max");
51  resSegmTrack_mean_min = ps.getParameter<double>("resSegmTrack_mean_min");
52  resSegmTrack_mean_max = ps.getParameter<double>("resSegmTrack_mean_max");
53  expPeakEcalS9_min = ps.getParameter<double>("expPeakEcalS9_min");
54  expPeakEcalS9_max = ps.getParameter<double>("expPeakEcalS9_max");
55  expPeakHadS9_min = ps.getParameter<double>("expPeakHadS9_min");
56  expPeakHadS9_max = ps.getParameter<double>("expPeakHadS9_max");
57  expMultiplicityGlb_max = ps.getParameter<double>("expMultiplicityGlb_max");
58  expMultiplicityTk_max = ps.getParameter<double>("expMultiplicityTk_max");
59  expMultiplicitySta_max = ps.getParameter<double>("expMultiplicitySta_max");
60  expMultiplicityGlb_min = ps.getParameter<double>("expMultiplicityGlb_min");
61  expMultiplicityTk_min = ps.getParameter<double>("expMultiplicityTk_min");
62  expMultiplicitySta_min = ps.getParameter<double>("expMultiplicitySta_min");
63 }
64 
65 MuonTestSummary::~MuonTestSummary() {}
66 
67 void MuonTestSummary::dqmEndJob(DQMStore::IBooker &ibooker, DQMStore::IGetter &igetter) {
68  // BOOKING NEW HISTOGRAMS
69  ibooker.cd();
70  ibooker.setCurrentFolder("Muons/TestSummary");
71 
72  // kinematics test report
73  kinematicsSummaryMap = ibooker.book2D("kinematicsSummaryMap", "Kinematics test summary", 5, 1, 6, 3, 1, 4);
74  kinematicsSummaryMap->setAxisTitle("track monitored", 1);
75  kinematicsSummaryMap->setBinLabel(1, "GLB", 1);
76  kinematicsSummaryMap->setBinLabel(2, "TKfromGLB", 1);
77  kinematicsSummaryMap->setBinLabel(3, "STAfromGLB", 1);
78  kinematicsSummaryMap->setBinLabel(4, "TK", 1);
79  kinematicsSummaryMap->setBinLabel(5, "STA", 1);
80  kinematicsSummaryMap->setAxisTitle("parameter tested", 2);
81  kinematicsSummaryMap->setBinLabel(1, "#chi^{2}", 2);
82  kinematicsSummaryMap->setBinLabel(2, "#eta", 2);
83  kinematicsSummaryMap->setBinLabel(3, "#phi", 2);
84 
85  //chi2 kinematics quality test report
86  chi2TestSummaryMap = ibooker.book2D("chi2TestSummaryMap", "#chi2 quality test summary", 5, 1, 6, 5, 1, 6);
87  chi2TestSummaryMap->setAxisTitle("track monitored", 1);
88  chi2TestSummaryMap->setBinLabel(1, "GLB", 1);
89  chi2TestSummaryMap->setBinLabel(2, "TKfromGLB", 1);
90  chi2TestSummaryMap->setBinLabel(3, "STAfromGLB", 1);
91  chi2TestSummaryMap->setBinLabel(4, "TK", 1);
92  chi2TestSummaryMap->setBinLabel(5, "STA", 1);
93  chi2TestSummaryMap->setAxisTitle("parameter tested", 2);
94  chi2TestSummaryMap->setBinLabel(1, "#chi^{2}", 2);
95  chi2TestSummaryMap->setBinLabel(2, "#eta", 2);
96  chi2TestSummaryMap->setBinLabel(3, "#phi", 2);
97  chi2TestSummaryMap->setBinLabel(4, "#pt", 2);
98  chi2TestSummaryMap->setBinLabel(5, "#q", 2);
99 
100  // residuals test report
101  residualsSummaryMap = ibooker.book2D("residualsSummaryMap", "Residuals test summary", 4, 1, 5, 4, 1, 5);
102  residualsSummaryMap->setAxisTitle("residuals", 1);
103  residualsSummaryMap->setBinLabel(1, "TK-GLB", 1);
104  residualsSummaryMap->setBinLabel(2, "GLB-STA", 1);
105  residualsSummaryMap->setBinLabel(3, "TK-STA", 1);
106  residualsSummaryMap->setBinLabel(4, "TK-STA Pull", 1);
107  residualsSummaryMap->setAxisTitle("parameter tested", 2);
108  residualsSummaryMap->setBinLabel(1, "#eta", 2);
109  residualsSummaryMap->setBinLabel(2, "#phi", 2);
110  residualsSummaryMap->setBinLabel(3, "1/p", 2);
111  residualsSummaryMap->setBinLabel(4, "q", 2);
112 
113  // muonId test report
114  muonIdSummaryMap = ibooker.book2D("muonIdSummaryMap", "muonId test summary", 4, 1, 5, 5, 1, 6);
115  muonIdSummaryMap->setAxisTitle("muons", 1);
116  muonIdSummaryMap->setBinLabel(1, "GLB DT", 1);
117  muonIdSummaryMap->setBinLabel(2, "GLB CSC", 1);
118  muonIdSummaryMap->setBinLabel(3, "TK DT", 1);
119  muonIdSummaryMap->setBinLabel(4, "TK CSC", 1);
120  muonIdSummaryMap->setAxisTitle("tests", 2);
121  muonIdSummaryMap->setBinLabel(1, "#assSeg", 2);
122  muonIdSummaryMap->setBinLabel(2, "x mean", 2);
123  muonIdSummaryMap->setBinLabel(3, "x rms", 2);
124  muonIdSummaryMap->setBinLabel(4, "y mean", 2);
125  muonIdSummaryMap->setBinLabel(5, "y rms", 2);
126 
127  // energy test report
128  energySummaryMap = ibooker.book2D("energySummaryMap", "Energy deposits test summary", 3, 1, 4, 3, 1, 4);
129  energySummaryMap->setAxisTitle("muons", 1);
130  energySummaryMap->setBinLabel(1, "GLB", 1);
131  energySummaryMap->setBinLabel(2, "TK", 1);
132  energySummaryMap->setBinLabel(3, "STA", 1);
133  energySummaryMap->setAxisTitle("calorimeter tested", 2);
134  energySummaryMap->setBinLabel(1, "ECAL", 2);
135  energySummaryMap->setBinLabel(2, "HAD", 2);
136  energySummaryMap->setBinLabel(3, "H0", 2);
137 
138  // multiplicity tests report
139  multiplicitySummaryMap = ibooker.book1D("multiplicitySummaryMap", "muon multiplicity test summary", 3, 1, 4);
140  multiplicitySummaryMap->setAxisTitle("muon");
141  multiplicitySummaryMap->setBinLabel(1, "GLB");
142  multiplicitySummaryMap->setBinLabel(2, "TK");
143  multiplicitySummaryMap->setBinLabel(3, "STA");
144 
145  // summary test report
146  ibooker.setCurrentFolder("Muons/EventInfo");
147  summaryReport = ibooker.bookFloat("reportSummary");
148 
149  summaryReportMap = ibooker.book2D("reportSummaryMap", "Muon Report Summary Map", 3, 1, 4, 7, 1, 8);
150  summaryReportMap->setAxisTitle("muons", 1);
151  summaryReportMap->setBinLabel(1, "GLB", 1);
152  summaryReportMap->setBinLabel(2, "TK", 1);
153  summaryReportMap->setBinLabel(3, "STA", 1);
154  summaryReportMap->setAxisTitle("test", 2);
155  summaryReportMap->setBinLabel(1, "#chi^{2}/Df", 2);
156  summaryReportMap->setBinLabel(2, "#eta", 2);
157  summaryReportMap->setBinLabel(3, "#phi", 2);
158  summaryReportMap->setBinLabel(4, "residuals", 2);
159  summaryReportMap->setBinLabel(5, "muonId", 2);
160  summaryReportMap->setBinLabel(6, "energyDeposits", 2);
161  summaryReportMap->setBinLabel(7, "multiplicity", 2);
162 
163  ibooker.setCurrentFolder("Muons/EventInfo/reportSummaryContents");
164  theSummaryContents.push_back(ibooker.bookFloat("kinematics_GLB"));
165  theSummaryContents.push_back(ibooker.bookFloat("muonId_GLB"));
166  theSummaryContents.push_back(ibooker.bookFloat("residuals_GLB"));
167  theSummaryContents.push_back(ibooker.bookFloat("GLB"));
168  theSummaryContents.push_back(ibooker.bookFloat("kinematics_TK"));
169  theSummaryContents.push_back(ibooker.bookFloat("muonId_TK"));
170  theSummaryContents.push_back(ibooker.bookFloat("residuals_TK"));
171  theSummaryContents.push_back(ibooker.bookFloat("TK"));
172  theSummaryContents.push_back(ibooker.bookFloat("kinematics_STA"));
173  theSummaryContents.push_back(ibooker.bookFloat("residuals_STA"));
174  theSummaryContents.push_back(ibooker.bookFloat("STA"));
175  theSummaryContents.push_back(ibooker.bookFloat("energyDeposits"));
176  theSummaryContents.push_back(ibooker.bookFloat("multiplicity"));
177 
178  // certification report
179  ibooker.setCurrentFolder("Muons/EventInfo");
180  summaryCertification = ibooker.bookFloat("CertificationSummary");
181  summaryCertification->Fill(-1);
182 
183  summaryCertificationMap =
184  ibooker.book2D("CertificationSummaryMap", "Muon Certification Summary Map", 9, 1, 10, 7, 1, 8);
185  summaryCertificationMap->setAxisTitle("muons", 1);
186  summaryCertificationMap->setBinLabel(1, "GLB_Tot", 1);
187  summaryCertificationMap->setBinLabel(2, "TK_Tot", 1);
188  summaryCertificationMap->setBinLabel(3, "STA_tot", 1);
189  summaryCertificationMap->setBinLabel(4, "GLB_B", 1);
190  summaryCertificationMap->setBinLabel(5, "TK_B", 1);
191  summaryCertificationMap->setBinLabel(6, "STA_B", 1);
192  summaryCertificationMap->setBinLabel(7, "GLB_EC", 1);
193  summaryCertificationMap->setBinLabel(8, "TK_EC", 1);
194  summaryCertificationMap->setBinLabel(9, "STA_EC", 1);
195  summaryCertificationMap->setAxisTitle("test", 2);
196  summaryCertificationMap->setBinLabel(1, "#chi^{2}/Df", 2);
197  summaryCertificationMap->setBinLabel(2, "#eta", 2);
198  summaryCertificationMap->setBinLabel(3, "#phi", 2);
199  summaryCertificationMap->setBinLabel(4, "residuals", 2);
200  summaryCertificationMap->setBinLabel(5, "muonId", 2);
201  summaryCertificationMap->setBinLabel(6, "energyDeposits", 2);
202  summaryCertificationMap->setBinLabel(7, "multiplicity", 2);
203 
204  ibooker.setCurrentFolder("Muons/EventInfo/CertificationContents");
205  theCertificationContents.push_back(ibooker.bookFloat("GLB_Tot"));
206  theCertificationContents.push_back(ibooker.bookFloat("STA_Tot"));
207  theCertificationContents.push_back(ibooker.bookFloat("TK_Tot"));
208  theCertificationContents.push_back(ibooker.bookFloat("GLB_B"));
209  theCertificationContents.push_back(ibooker.bookFloat("STA_B"));
210  theCertificationContents.push_back(ibooker.bookFloat("TK_B"));
211  theCertificationContents.push_back(ibooker.bookFloat("GLB_EC"));
212  theCertificationContents.push_back(ibooker.bookFloat("STA_EC"));
213  theCertificationContents.push_back(ibooker.bookFloat("TK_EC"));
214 
215  for (unsigned int icert = 0; icert < theCertificationContents.size(); icert++) {
216  theCertificationContents[icert]->Fill(-1);
217  }
218 
219  // initialisation of histo bins
220  for (int xBin = 1; xBin <= 5; xBin++) {
221  for (int yBin = 1; yBin <= 3; yBin++) {
222  kinematicsSummaryMap->Fill(xBin, yBin, 1);
223  }
224  for (int yBin = 1; yBin <= 5; yBin++) {
225  chi2TestSummaryMap->Fill(xBin, yBin, 1);
226  }
227  }
228  for (int xBin = 1; xBin <= residualsSummaryMap->getNbinsX(); xBin++) {
229  for (int yBin = 1; yBin <= residualsSummaryMap->getNbinsY(); yBin++) {
230  residualsSummaryMap->Fill(xBin, yBin, 1);
231  }
232  }
233  residualsSummaryMap->setBinContent(4, 4, 1); //not used for now
234 
235  for (int xBin = 1; xBin <= muonIdSummaryMap->getNbinsX(); xBin++) {
236  for (int yBin = 1; yBin <= muonIdSummaryMap->getNbinsY(); yBin++) {
237  muonIdSummaryMap->Fill(xBin, yBin, 1);
238  }
239  }
240  for (int xBin = 1; xBin <= 3; xBin++) {
241  for (int yBin = 1; yBin <= 3; yBin++) {
242  energySummaryMap->Fill(xBin, yBin, 1);
243  }
244  }
245  for (int xBin = 1; xBin <= 3; xBin++) {
246  multiplicitySummaryMap->Fill(xBin, 1);
247  }
248  for (int xBin = 1; xBin <= 3; xBin++) {
249  for (int yBin = 1; yBin <= 7; yBin++) {
250  summaryReportMap->Fill(xBin, yBin, 1);
251  }
252  }
253  for (int xBin = 1; xBin <= 9; xBin++) {
254  for (int yBin = 1; yBin <= 7; yBin++) {
255  summaryCertificationMap->Fill(xBin, yBin, 1);
256  }
257  }
260 
261  // fill the kinematics report summary
262  doKinematicsTests(igetter, "GlbMuon_Glb_", 1);
263  doKinematicsTests(igetter, "GlbMuon_Tk_", 2);
264  doKinematicsTests(igetter, "GlbMuon_Sta_", 3);
265  doKinematicsTests(igetter, "TkMuon_", 4);
266  doKinematicsTests(igetter, "StaMuon_", 5);
267 
268  // fill the residuals report summary
269  doResidualsTests(igetter, "TkGlb", "eta", 1);
270  doResidualsTests(igetter, "GlbSta", "eta", 2);
271  doResidualsTests(igetter, "TkSta", "eta", 3);
272  doResidualsTests(igetter, "TkGlb", "phi", 1);
273  doResidualsTests(igetter, "GlbSta", "phi", 2);
274  doResidualsTests(igetter, "TkSta", "phi", 3);
275  doResidualsTests(igetter, "TkGlb", "oneOverp", 1);
276  doResidualsTests(igetter, "GlbSta", "oneOverp", 2);
277  doResidualsTests(igetter, "TkSta", "oneOverp", 3);
278  doResidualsTests(igetter, "GlbMuon", "qComparison", -1);
279 
280  // fill the muonID report summary
281  doMuonIDTests(igetter);
282 
283  // fill the energy report summary
284  doEnergyTests(igetter, "ecalS9PointingMuDepositedEnergy_", "Glb_muons", 1);
285  doEnergyTests(igetter, "hadS9PointingMuDepositedEnergy_", "Glb_muons", 1);
286  doEnergyTests(igetter, "hoS9PointingMuDepositedEnergy_", "Glb_muons", 1);
287  doEnergyTests(igetter, "ecalS9PointingMuDepositedEnergy_", "Tk_muons", 2);
288  doEnergyTests(igetter, "hadS9PointingMuDepositedEnergy_", "Tk_muons", 2);
289  doEnergyTests(igetter, "hoS9PointingMuDepositedEnergy_", "Tk_muons", 2);
290  doEnergyTests(igetter, "ecalS9PointingMuDepositedEnergy_", "Sta_muons", 3);
291  doEnergyTests(igetter, "hadS9PointingMuDepositedEnergy_", "Sta_muons", 3);
292  doEnergyTests(igetter, "hoS9PointingMuDepositedEnergy_", "Sta_muons", 3);
293 
294  // fill the multiplicity test summary
295  doMultiplicityTests(igetter);
296 
297  // fill the final report summary
298 
299  //-- modified GH
300  double residualsSummary = 0;
301  //put the TRK-STA resid & pulls in the first bin ("GLB")
302  //then the GLB-TRK and GLB-STA residuals in the 2nd and 3rd
303  for (int i = 3; i <= residualsSummaryMap->getNbinsX(); i++)
304  for (int j = 1; j <= residualsSummaryMap->getNbinsY(); j++)
305  residualsSummary += residualsSummaryMap->getBinContent(i, j);
306  residualsSummary /= 2 * residualsSummaryMap->getNbinsY();
307  summaryReportMap->setBinContent(1, 4, residualsSummary);
308 
309  residualsSummary = 0;
310  for (int i = 1; i <= 1; i++)
311  for (int j = 1; j <= residualsSummaryMap->getNbinsY(); j++)
312  residualsSummary += residualsSummaryMap->getBinContent(i, j);
313  residualsSummary /= 1 * residualsSummaryMap->getNbinsY();
314  summaryReportMap->setBinContent(2, 4, residualsSummary);
315 
316  residualsSummary = 0;
317  for (int i = 2; i <= 2; i++)
318  for (int j = 1; j <= residualsSummaryMap->getNbinsY(); j++)
319  residualsSummary += residualsSummaryMap->getBinContent(i, j);
320  residualsSummary /= 1 * residualsSummaryMap->getNbinsY();
321  summaryReportMap->setBinContent(3, 4, residualsSummary);
322 
323  //--
324 
325  //-- modified GH
326  float idtest = 0;
327  for (int i = 1; i <= 2; i++)
328  for (int j = 1; j <= 5; j++) {
329  if (j == 3 || j == 5)
330  continue; //ignore pull widths for now
331  idtest += muonIdSummaryMap->getBinContent(i, j);
332  }
333  // idtest/=10.;
334  idtest /= 6.;
335  summaryReportMap->setBinContent(1, 5, idtest);
336  idtest = 0;
337  for (int i = 3; i <= 4; i++)
338  for (int j = 1; j <= 5; j++) {
339  if (j == 3 || j == 5)
340  continue; //ignore pull widths for now
341  idtest += muonIdSummaryMap->getBinContent(i, j);
342  }
343  // idtest/=10.;
344  idtest /= 6.;
345  summaryReportMap->setBinContent(2, 5, idtest);
346  summaryReportMap->setBinContent(3, 5, -1.0 / 6.0);
347  //--
348 
349  summaryReportMap->setBinContent(
350  1, 6, double(energySummaryMap->getBinContent(1, 1) + energySummaryMap->getBinContent(1, 2)) / 2.0);
351  summaryReportMap->setBinContent(
352  2, 6, double(energySummaryMap->getBinContent(2, 1) + energySummaryMap->getBinContent(2, 2)) / 2.0);
353  summaryReportMap->setBinContent(
354  3, 6, double(energySummaryMap->getBinContent(3, 1) + energySummaryMap->getBinContent(3, 2)) / 2.0);
355  summaryReportMap->setBinContent(1, 7, multiplicitySummaryMap->getBinContent(1));
356  summaryReportMap->setBinContent(2, 7, multiplicitySummaryMap->getBinContent(2));
357  summaryReportMap->setBinContent(3, 7, multiplicitySummaryMap->getBinContent(3));
358 
359  double kinematics_GLB = double(summaryReportMap->getBinContent(1, 1) + summaryReportMap->getBinContent(1, 2) +
360  summaryReportMap->getBinContent(1, 3)) /
361  3.0;
362  theSummaryContents[0]->Fill(kinematics_GLB);
363  double muonId_GLB = double(summaryReportMap->getBinContent(1, 5));
364  theSummaryContents[1]->Fill(muonId_GLB);
365  double residuals_GLB = double(summaryReportMap->getBinContent(1, 4));
366  theSummaryContents[2]->Fill(residuals_GLB);
367  double GLB = (kinematics_GLB + muonId_GLB + residuals_GLB) / 3.0;
368  theSummaryContents[3]->Fill(GLB);
369 
370  double kinematics_TK = double(summaryReportMap->getBinContent(2, 1) + summaryReportMap->getBinContent(2, 2) +
371  summaryReportMap->getBinContent(2, 3)) /
372  3.0;
373  theSummaryContents[4]->Fill(kinematics_TK);
374  double muonId_TK = double(summaryReportMap->getBinContent(2, 5));
375  theSummaryContents[5]->Fill(muonId_TK);
376  double residuals_TK = double(summaryReportMap->getBinContent(2, 4));
377  theSummaryContents[6]->Fill(residuals_TK);
378  double TK = double(kinematics_TK + muonId_TK + residuals_TK) / 3.0;
379  theSummaryContents[7]->Fill(TK);
380 
381  double kinematics_STA = double(summaryReportMap->getBinContent(3, 1) + summaryReportMap->getBinContent(3, 2) +
382  summaryReportMap->getBinContent(3, 3)) /
383  3.0;
384  theSummaryContents[8]->Fill(kinematics_STA);
385  double residuals_STA = double(summaryReportMap->getBinContent(3, 4));
386  theSummaryContents[9]->Fill(residuals_STA);
387  double STA = double(kinematics_STA + residuals_STA) / 2.0;
388  theSummaryContents[10]->Fill(STA);
389  double energyDeposits = double(summaryReportMap->getBinContent(1, 6) + summaryReportMap->getBinContent(2, 6) +
390  summaryReportMap->getBinContent(3, 6)) /
391  3.0;
392  theSummaryContents[11]->Fill(energyDeposits);
393  double multiplicity = double(summaryReportMap->getBinContent(1, 7) + summaryReportMap->getBinContent(2, 7) +
394  summaryReportMap->getBinContent(3, 7)) /
395  3.0;
396  theSummaryContents[12]->Fill(multiplicity);
397 
398  summaryReport->Fill((GLB + TK + STA + energyDeposits + multiplicity) / 5.0);
399 
400  //global barrel:
401  float muonIDsummary = 0;
402  // for(int i=2; i<=5; i++)
403  // muonIDsummary += muonIdSummaryMap->getBinContent(2, i);
404  // summaryCertificationMap->setBinContent(4, 5, muonIDsummary/4.);
405  //for now, just report the mean:
406  muonIDsummary += muonIdSummaryMap->getBinContent(1, 2);
407  muonIDsummary += muonIdSummaryMap->getBinContent(1, 4);
408  summaryCertificationMap->setBinContent(4, 5, muonIDsummary / 2.);
409 
410  //global EC:
411  muonIDsummary = 0;
412  // for(int i=2; i<=5; i++)
413  // muonIDsummary += muonIdSummaryMap->getBinContent(2, i);
414  // summaryCertificationMap->setBinContent(7, 5, muonIDsummary/4.);
415  muonIDsummary += muonIdSummaryMap->getBinContent(2, 2);
416  muonIDsummary += muonIdSummaryMap->getBinContent(2, 4);
417  summaryCertificationMap->setBinContent(7, 5, muonIDsummary / 2.);
418 
419  //tracker barrel:
420  muonIDsummary = 0;
421  // for(int i=2; i<=5; i++)
422  // muonIDsummary += muonIdSummaryMap->getBinContent(3, i);
423  // summaryCertificationMap->setBinContent(5, 5, muonIDsummary/4.);
424  muonIDsummary += muonIdSummaryMap->getBinContent(3, 2);
425  muonIDsummary += muonIdSummaryMap->getBinContent(3, 4);
426  summaryCertificationMap->setBinContent(5, 5, muonIDsummary / 2.);
427 
428  //tracker EC:
429  muonIDsummary = 0;
430  // for(int i=2; i<=5; i++)
431  // muonIDsummary += muonIdSummaryMap->getBinContent(4, i);
432  // summaryCertificationMap->setBinContent(8, 5, muonIDsummary/4.);
433  muonIDsummary += muonIdSummaryMap->getBinContent(4, 2);
434  muonIDsummary += muonIdSummaryMap->getBinContent(4, 4);
435  summaryCertificationMap->setBinContent(8, 5, muonIDsummary / 2.);
436 
437  double muonId_GLB_B = double(summaryCertificationMap->getBinContent(4, 5));
438  theCertificationContents[3]->Fill(muonId_GLB_B);
439  double muonId_GLB_EC = double(summaryCertificationMap->getBinContent(7, 5));
440  theCertificationContents[6]->Fill(muonId_GLB_EC);
441 
442  double muonId_TK_B = double(summaryCertificationMap->getBinContent(5, 5));
443  theCertificationContents[5]->Fill(muonId_TK_B);
444  double muonId_TK_EC = double(summaryCertificationMap->getBinContent(8, 5));
445  theCertificationContents[8]->Fill(muonId_TK_EC);
446 }
447 
448 void MuonTestSummary::doKinematicsTests(DQMStore::IGetter &igetter, string muonType, int bin) {
449  // chi2 test
450  string path = "Muons/MuonRecoAnalyzer/" + muonType + "chi2OverDf";
451  MonitorElement *chi2Histo = igetter.get(path);
452 
453  if (chi2Histo) {
454  TH1F *chi2Histo_root = chi2Histo->getTH1F();
455  if (chi2Histo_root->GetEntries() > 20) {
456  //Standard QT based on fraction of events below and above a cut
457  LogTrace(metname) << "chi2 kin test based on fraction for " << muonType << endl;
458  int maxBin = chi2Histo_root->GetMaximumBin();
459  if (chi2Histo_root->Integral(maxBin + 1, chi2Histo_root->GetNbinsX()) != 0) {
460  double fraction = double(chi2Histo_root->Integral(1, maxBin)) /
461  double(chi2Histo_root->Integral(maxBin + 1, chi2Histo_root->GetNbinsX()));
462  LogTrace(metname) << "chi2 fraction for " << muonType << " : " << fraction << " must be within "
463  << chi2Fraction - chi2Spread << "," << chi2Fraction + chi2Spread << endl;
464  if (fraction > (chi2Fraction - chi2Spread) && fraction < (chi2Fraction + chi2Spread))
465  kinematicsSummaryMap->setBinContent(bin, 1, 1);
466  else
467  kinematicsSummaryMap->setBinContent(bin, 1, 0);
468  }
469  } else {
470  LogTrace(metname) << "[MuonTestSummary]: Test of Chi2 Kin not performed for " << muonType
471  << " because # entries < 20 ";
472  }
473  }
474 
475  // pseudorapidity test
476  path = "Muons/MuonRecoAnalyzer/" + muonType + "eta";
477  MonitorElement *etaHisto = igetter.get(path);
478 
479  if (etaHisto) {
480  TH1F *etaHisto_root = etaHisto->getTH1F();
481  if (etaHisto_root->GetEntries() > 20) {
482  //Standard QT based on fraction of events below and above a cut
483  LogTrace(metname) << "eta kin test based on fraction for " << muonType << endl;
484  double binSize =
485  (etaHisto_root->GetXaxis()->GetXmax() - etaHisto_root->GetXaxis()->GetXmin()) / etaHisto_root->GetNbinsX();
486  int binZero = int((0 - etaHisto_root->GetXaxis()->GetXmin()) / binSize);
487  if (etaHisto_root->Integral(1, binZero - 1) != 0 &&
488  etaHisto_root->Integral(binZero, etaHisto_root->GetNbinsX()) != 0) {
489  double symmetryFactor = double(etaHisto_root->Integral(1, binZero - 1)) /
490  double(etaHisto_root->Integral(binZero, etaHisto_root->GetNbinsX()));
491  double errSymmetryFactor =
492  symmetryFactor * sqrt(1.0 / double(etaHisto_root->Integral(1, binZero - 1)) +
493  1.0 / double(etaHisto_root->Integral(binZero, etaHisto_root->GetNbinsX())));
494  LogTrace(metname) << "eta symmetryFactor for " << muonType << " : " << symmetryFactor
495  << " (expected :" << etaExpected << ")" << endl;
496  LogTrace(metname) << "eta errSymmetryFactor for " << muonType << " : " << errSymmetryFactor << endl;
497  double tParameter;
498  if ((symmetryFactor - etaExpected) > 0)
499  tParameter = double(symmetryFactor - etaExpected) / errSymmetryFactor;
500  else
501  tParameter = double(-symmetryFactor + etaExpected) / errSymmetryFactor;
502  LogTrace(metname) << "eta tParameter for " << muonType << " : " << tParameter << " (expected < 1.95)" << endl;
503  if (tParameter < 1.95) //2sigma rejection
504  kinematicsSummaryMap->setBinContent(bin, 2, 1);
505  else
506  kinematicsSummaryMap->setBinContent(bin, 2, 0);
507  }
508 
509  } else {
510  LogTrace(metname) << "[MuonTestSummary]: Test of Eta Kin not performed for " << muonType
511  << " because # entries < 20 ";
512  }
513  }
514 
515  // phi test
516  path = "Muons/MuonRecoAnalyzer/" + muonType + "phi";
517  MonitorElement *phiHisto = igetter.get(path);
518 
519  if (phiHisto) {
520  TH1F *phiHisto_root = phiHisto->getTH1F();
521  if (phiHisto_root->GetEntries() > 20) {
522  //Standard QT based on fraction of events below and above a cut
523  LogTrace(metname) << "phi kin test based on fraction for " << muonType << endl;
524  double binSize =
525  (phiHisto_root->GetXaxis()->GetXmax() - phiHisto_root->GetXaxis()->GetXmin()) / phiHisto_root->GetNbinsX();
526  int binZero = int((0 - phiHisto_root->GetXaxis()->GetXmin()) / binSize);
527  if (phiHisto_root->Integral(binZero + 1, phiHisto_root->GetNbinsX()) != 0 &&
528  phiHisto_root->Integral(1, binZero) != 0) {
529  double symmetryFactor = double(phiHisto_root->Integral(binZero + 1, phiHisto_root->GetNbinsX())) /
530  double(phiHisto_root->Integral(1, binZero));
531  double errSymmetryFactor =
532  symmetryFactor * sqrt(1.0 / double(phiHisto_root->Integral(binZero + 1, phiHisto_root->GetNbinsX())) +
533  1.0 / double(phiHisto_root->Integral(1, binZero)));
534  LogTrace(metname) << "phi symmetryFactor for " << muonType << " : " << symmetryFactor
535  << "(phi expected :" << phiExpected << ")" << endl;
536  LogTrace(metname) << "phi errSymmetryFactor for " << muonType << " : " << errSymmetryFactor << endl;
537  double tParameter;
538  if ((symmetryFactor - phiExpected) > 0)
539  tParameter = double(symmetryFactor - phiExpected) / errSymmetryFactor;
540  else
541  tParameter = double(-symmetryFactor + phiExpected) / errSymmetryFactor;
542  LogTrace(metname) << "phi tParameter for " << muonType << " : " << tParameter << " (expected < 1.95)" << endl;
543  if (tParameter < 1.95) //2sigma rejection
544  kinematicsSummaryMap->setBinContent(bin, 3, 1);
545  else
546  kinematicsSummaryMap->setBinContent(bin, 3, 0);
547  }
548 
549  } else {
550  LogTrace(metname) << "[MuonTestSummary]: Test of Phi Kin not performed for " << muonType
551  << " because # entries < 20 ";
552  }
553  }
554 }
555 //--GH new
556 void MuonTestSummary::GaussFit(
557  string type, string parameter, MonitorElement *Histo, float &mean, float &mean_err, float &sigma, float &sigma_err) {
558  // Gaussian Fit
559  float statMean = Histo->getMean(1);
560  float statSigma = Histo->getRMS(1);
561  TH1F *histo_root = Histo->getTH1F();
562  if (histo_root->GetEntries() > 20) {
563  TF1 *gfit = new TF1("Gaussian", "gaus", (statMean - (2 * statSigma)), (statMean + (2 * statSigma)));
564  try {
565  histo_root->Fit(gfit, "Q0");
566  } catch (cms::Exception &iException) {
567  edm::LogError(metname) << "[MuonTestSummary]: Exception when fitting Res_" << type << "_" << parameter;
568  mean = 1;
569  mean_err = 1;
570  sigma = 1;
571  sigma_err = 1;
572  return;
573  }
574  if (gfit) {
575  mean = gfit->GetParameter(1);
576  mean_err = gfit->GetParErrors()[2];
577  sigma = gfit->GetParameter(2);
578  sigma_err = gfit->GetParErrors()[2];
579  LogTrace(metname) << "Gaussian fit mean: " << mean << " +- " << mean_err << " for " << type << "_" << parameter
580  << endl;
581  LogTrace(metname) << "Gaussina fit sigma: " << sigma << " +- " << sigma_err << " for " << type << "_" << parameter
582  << endl;
583  }
584  } else {
585  LogTrace(metname) << "[MuonTestSummary]: Test of Res_" << type << "_" << parameter
586  << " not performed because # entries < 20 ";
587  //auto-pass if not enough events.
588  mean = 1;
589  mean_err = 1;
590  sigma = 1;
591  sigma_err = 1;
592  }
593 }
594 void MuonTestSummary::doResidualsTests(DQMStore::IGetter &igetter, string type, string parameter, int bin) {
595  // residuals test
596  if (type != "GlbMuon") {
597  string path = "Muons/MuonRecoAnalyzer/Res_" + type + "_" + parameter;
598  MonitorElement *residualsHisto = igetter.get(path);
599 
600  float mean = -1;
601  float mean_err = -1;
602  float sigma = -1;
603  float sigma_err = -1;
604 
605  if (residualsHisto) {
606  LogTrace(metname) << "[MuonTestSummary]: Starting Gaussian fit for Test of Res_" << type << "_" << parameter
607  << endl;
608  GaussFit(type, parameter, residualsHisto, mean, mean_err, sigma, sigma_err);
609 
610  if (sigma != -1 && parameter == "eta" && type == "TkGlb") {
611  if (sigma - sigma_err < resEtaSpread_tkGlb)
612  residualsSummaryMap->setBinContent(bin, 1, 1);
613  else
614  residualsSummaryMap->setBinContent(bin, 1, 0);
615  }
616  if (sigma != -1 && parameter == "eta" && (type == "GlbSta" || type == "TkSta")) {
617  if (sigma - sigma_err < resEtaSpread_glbSta)
618  residualsSummaryMap->setBinContent(bin, 1, 1);
619  else
620  residualsSummaryMap->setBinContent(bin, 1, 0);
621  }
622  if (sigma != -1 && parameter == "phi" && type == "TkGlb") {
623  if (sigma - sigma_err < resPhiSpread_tkGlb)
624  residualsSummaryMap->setBinContent(bin, 2, 1);
625  else
626  residualsSummaryMap->setBinContent(bin, 2, 0);
627  }
628  if (sigma != -1 && parameter == "phi" && (type == "GlbSta" || type == "TkSta")) {
629  if (sigma - sigma_err < resPhiSpread_glbSta)
630  residualsSummaryMap->setBinContent(bin, 2, 1);
631  else
632  residualsSummaryMap->setBinContent(bin, 2, 0);
633  }
634  if (sigma != -1 && parameter == "oneOverp" && type == "TkGlb") {
635  if (sigma - sigma_err < resOneOvPSpread_tkGlb)
636  residualsSummaryMap->setBinContent(bin, 3, 1);
637  else
638  residualsSummaryMap->setBinContent(bin, 3, 0);
639  }
640  if (sigma != -1 && parameter == "oneOverp" && (type == "GlbSta" || type == "TkSta")) {
641  if (sigma - sigma_err < resOneOvPSpread_glbSta)
642  residualsSummaryMap->setBinContent(bin, 3, 1);
643  else
644  residualsSummaryMap->setBinContent(bin, 3, 0);
645  }
646  }
647 
648  //--GH modified
649  if (type == "TkSta") {
650  //look at the pull:
651  string path = "Muons/MuonRecoAnalyzer/Pull_" + type + "_" + parameter;
652  MonitorElement *pullHisto = igetter.get(path);
653 
654  if (pullHisto) {
655  LogTrace(metname) << "[MuonTestSummary]: Starting Gaussian fit for Test of Pull_" << type << "_" << parameter
656  << endl;
657  GaussFit(type, parameter, pullHisto, mean, mean_err, sigma, sigma_err);
658 
659  if (sigma != -1 && parameter == "eta") {
660  if (sigma - sigma_err < pullEtaSpread)
661  residualsSummaryMap->setBinContent(4, 1, 1);
662  else
663  residualsSummaryMap->setBinContent(4, 1, 0);
664  }
665  if (sigma != -1 && parameter == "phi") {
666  if (sigma - sigma_err < pullPhiSpread)
667  residualsSummaryMap->setBinContent(4, 2, 1);
668  else
669  residualsSummaryMap->setBinContent(4, 2, 0);
670  }
671  if (sigma != -1 && parameter == "oneOverp") {
672  if (sigma - sigma_err < pullOneOvPSpread)
673  residualsSummaryMap->setBinContent(4, 3, 1);
674  else
675  residualsSummaryMap->setBinContent(4, 3, 0);
676  }
677 
678  } //have pull histo
679  } //TkSta muons
680  }
681 
682  //this part for Global Muons:
683  else {
684  string path = "Muons/MuonRecoAnalyzer/" + type + "_" + parameter;
685  MonitorElement *residualsHisto = igetter.get(path);
686 
687  if (residualsHisto) {
688  LogTrace(metname) << "[MuonTestSummary]: Test of Charge Comparison " << type << "_" << parameter << endl;
689  if ((residualsHisto->getBinContent(3) + residualsHisto->getBinContent(4)) != 0) {
690  LogTrace(metname) << "Charge comparison TkGlb: "
691  << residualsHisto->getBinContent(4) /
692  double(residualsHisto->getBinContent(3) + residualsHisto->getBinContent(4))
693  << endl;
694  if (residualsHisto->getBinContent(4) /
695  double(residualsHisto->getBinContent(3) + residualsHisto->getBinContent(4)) <
696  resChargeLimit_tkGlb)
697  residualsSummaryMap->setBinContent(1, 4, 1);
698  else
699  residualsSummaryMap->setBinContent(1, 4, 0);
700  }
701  if ((residualsHisto->getBinContent(1) + residualsHisto->getBinContent(2)) != 0) {
702  LogTrace(metname) << "charge comparison GlbSta: "
703  << residualsHisto->getBinContent(2) /
704  double(residualsHisto->getBinContent(1) + residualsHisto->getBinContent(2))
705  << endl;
706  if (residualsHisto->getBinContent(2) /
707  double(residualsHisto->getBinContent(1) + residualsHisto->getBinContent(2)) <
708  resChargeLimit_glbSta)
709  residualsSummaryMap->setBinContent(2, 4, 1);
710  else
711  residualsSummaryMap->setBinContent(2, 4, 0);
712  }
713  if (residualsHisto->getBinContent(5) + residualsHisto->getBinContent(6) != 0) {
714  LogTrace(metname) << "charge comparison TkSta: "
715  << residualsHisto->getBinContent(6) /
716  double(residualsHisto->getBinContent(5) + residualsHisto->getBinContent(6))
717  << endl;
718  if (residualsHisto->getBinContent(6) /
719  double(residualsHisto->getBinContent(5) + residualsHisto->getBinContent(6)) <
720  resChargeLimit_tkSta)
721  residualsSummaryMap->setBinContent(3, 4, 1);
722  else
723  residualsSummaryMap->setBinContent(3, 4, 0);
724  }
725  }
726  }
727 }
728 
729 void MuonTestSummary::doMuonIDTests(DQMStore::IGetter &igetter) {
730  vector<string> muType;
731  muType.push_back("GlobalMuons");
732  muType.push_back("TrackerMuons");
733 
734  for (int i = 0; i <= 1; i++) {
735  // num matches test
736  string path = "Muons/MuonIdDQM/" + muType[i] + "/hNumMatches";
737  MonitorElement *matchesHisto = igetter.get(path);
738 
739  if (matchesHisto) {
740  TH1F *matchesHisto_root = matchesHisto->getTH1F();
741  if (matchesHisto_root->GetMaximumBin() >= numMatchedExpected_min &&
742  matchesHisto_root->GetMaximumBin() <= numMatchedExpected_max)
743  muonIdSummaryMap->setBinContent(i + 1, 1, 1);
744  else
745  muonIdSummaryMap->setBinContent(i + 1, 1, 0);
746  }
747 
748  // num of 0 associated segments
749  double numOneSegm_dt = 0;
750  int numHistos_dt = 0;
751  int numHistos_csc = 0;
752  MonitorElement *DT1Histo = igetter.get("Muons/MuonIdDQM/" + muType[i] + "/hDT1NumSegments");
753  if (DT1Histo) {
754  numHistos_dt++;
755  if (DT1Histo->getEntries() != 0)
756  numOneSegm_dt += double(DT1Histo->getBinContent(2)) / double(DT1Histo->getEntries());
757  }
758  MonitorElement *DT2Histo = igetter.get("Muons/MuonIdDQM/" + muType[i] + "/hDT2NumSegments");
759  if (DT2Histo) {
760  numHistos_dt++;
761  if (DT2Histo->getEntries() != 0)
762  numOneSegm_dt += double(DT2Histo->getBinContent(2)) / double(DT2Histo->getEntries());
763  }
764  MonitorElement *DT3Histo = igetter.get("Muons/MuonIdDQM/" + muType[i] + "/hDT3NumSegments");
765  if (DT3Histo) {
766  numHistos_dt++;
767  if (DT3Histo->getEntries() != 0)
768  numOneSegm_dt += double(DT3Histo->getBinContent(2)) / double(DT3Histo->getEntries());
769  }
770  MonitorElement *DT4Histo = igetter.get("Muons/MuonIdDQM/" + muType[i] + "/hDT4NumSegments");
771  if (DT4Histo) {
772  numHistos_dt++;
773  if (DT4Histo->getEntries() != 0)
774  numOneSegm_dt += double(DT4Histo->getBinContent(2)) / double(DT4Histo->getEntries());
775  }
776  double fraction_dt = 0;
777  if (numOneSegm_dt != 0) {
778  fraction_dt = numOneSegm_dt / double(numHistos_dt);
779  LogTrace(metname) << "fraction_dt: " << fraction_dt << " for " << muType[i] << endl;
780  }
781 
782  double numOneSegm_csc = 0;
783  MonitorElement *CSC1Histo = igetter.get("Muons/MuonIdDQM/" + muType[i] + "/hCSC1NumSegments");
784  if (CSC1Histo) {
785  numHistos_csc++;
786  if (CSC1Histo->getEntries() != 0)
787  numOneSegm_csc += double(CSC1Histo->getBinContent(2)) / double(CSC1Histo->getEntries());
788  }
789  MonitorElement *CSC2Histo = igetter.get("Muons/MuonIdDQM/" + muType[i] + "/hCSC2NumSegments");
790  if (CSC2Histo) {
791  numHistos_csc++;
792  if (CSC2Histo->getEntries() != 0)
793  numOneSegm_csc += double(CSC2Histo->getBinContent(2)) / double(CSC2Histo->getEntries());
794  }
795  MonitorElement *CSC3Histo = igetter.get("Muons/MuonIdDQM/" + muType[i] + "/hCSC3NumSegments");
796  if (CSC3Histo) {
797  numHistos_csc++;
798  if (CSC3Histo->getEntries() != 0)
799  numOneSegm_csc += double(CSC3Histo->getBinContent(2)) / double(CSC3Histo->getEntries());
800  }
801  MonitorElement *CSC4Histo = igetter.get("Muons/MuonIdDQM/" + muType[i] + "/hCSC4NumSegments");
802  if (CSC4Histo) {
803  numHistos_csc++;
804  if (CSC4Histo->getEntries() != 0)
805  numOneSegm_csc += double(CSC4Histo->getBinContent(2)) / double(CSC4Histo->getEntries());
806  }
807  double fraction_csc = 0;
808  if (numOneSegm_csc != 0) {
809  fraction_csc = numOneSegm_csc / double(numHistos_csc);
810  LogTrace(metname) << "fraction_csc: " << fraction_csc << " for " << muType[i] << endl;
811  }
812 
813  //--GH modified
814 
815  if (fraction_dt > matchesFractionDt_min && fraction_dt < matchesFractionDt_max)
816  muonIdSummaryMap->setBinContent(2 * i + 1, 1, 1);
817  else
818  muonIdSummaryMap->setBinContent(2 * i + 1, 1, 0);
819 
820  if (fraction_csc > matchesFractionCsc_min && fraction_csc < matchesFractionCsc_max)
821  muonIdSummaryMap->setBinContent(2 * i + 2, 1, 1);
822  else
823  muonIdSummaryMap->setBinContent(2 * i + 2, 1, 0);
824 
825  //--GH modified
826 
827  // residuals test
828  vector<string> DTXresHistos, DTYresHistos, CSCXresHistos, CSCYresHistos;
829  DTXresHistos.push_back("hDT1Pullx");
830  DTXresHistos.push_back("hDT2Pullx");
831  DTXresHistos.push_back("hDT3Pullx");
832  DTXresHistos.push_back("hDT4Pullx");
833 
834  DTYresHistos.push_back("hDT1Pully");
835  DTYresHistos.push_back("hDT2Pully");
836  DTYresHistos.push_back("hDT3Pully");
837 
838  CSCXresHistos.push_back("hCSC1Pullx");
839  CSCXresHistos.push_back("hCSC2Pullx");
840  CSCXresHistos.push_back("hCSC3Pullx");
841  CSCXresHistos.push_back("hCSC4Pullx");
842 
843  CSCYresHistos.push_back("hCSC1Pully");
844  CSCYresHistos.push_back("hCSC2Pully");
845  CSCYresHistos.push_back("hCSC3Pully");
846  CSCYresHistos.push_back("hCSC4Pully");
847 
848  int numPlot_dtX, numPlot_dtY, numPlot_cscX, numPlot_cscY;
849  double dtSigmaX, dtSigmaY, cscSigmaX, cscSigmaY;
850  double dtSigmaX_err, dtSigmaY_err, cscSigmaX_err, cscSigmaY_err;
851  double dtMeanX, dtMeanY, cscMeanX, cscMeanY;
852  double dtMeanX_err, dtMeanY_err, cscMeanX_err, cscMeanY_err;
853  MuonTestSummary::ResidualCheck(
854  igetter, muType[i], DTXresHistos, numPlot_dtX, dtMeanX, dtMeanX_err, dtSigmaX, dtSigmaX_err);
855  MuonTestSummary::ResidualCheck(
856  igetter, muType[i], DTYresHistos, numPlot_dtY, dtMeanY, dtMeanY_err, dtSigmaY, dtSigmaY_err);
857  MuonTestSummary::ResidualCheck(
858  igetter, muType[i], CSCXresHistos, numPlot_cscX, cscMeanX, cscMeanX_err, cscSigmaX, cscSigmaX_err);
859  MuonTestSummary::ResidualCheck(
860  igetter, muType[i], CSCYresHistos, numPlot_cscY, cscMeanY, cscMeanY_err, cscSigmaY, cscSigmaY_err);
861 
862  LogTrace(metname) << "DT mean must be between: " << resSegmTrack_mean_min << " and " << resSegmTrack_mean_max
863  << endl;
864  LogTrace(metname) << "DT rms must be between: " << resSegmTrack_rms_min << " and " << resSegmTrack_rms_max << endl;
865  LogTrace(metname) << "DT X residual " << muType[i] << " mean: " << dtMeanX << " +- " << dtMeanX_err
866  << ", sigma: " << dtSigmaX << " +- " << dtSigmaX_err << endl;
867  LogTrace(metname) << "DT Y residual " << muType[i] << " mean: " << dtMeanY << " +- " << dtMeanY_err
868  << ", sigma: " << dtSigmaY << " +- " << dtSigmaY_err << endl;
869  LogTrace(metname) << "CSC X residual " << muType[i] << " mean: " << cscMeanX << " +- " << cscMeanX_err
870  << ", sigma: " << cscSigmaX << " +- " << cscSigmaX_err << endl;
871  LogTrace(metname) << "CSC Y residual " << muType[i] << " mean: " << cscMeanY << " +- " << cscMeanY_err
872  << ", sigma: " << cscSigmaY << " +- " << cscSigmaY_err << endl;
873 
874  //require the mean and rms to be within nsig sigma of preferred range;
875  const int nsig = 2;
876  if (numPlot_dtX > 0) {
877  if (dtMeanX + nsig * dtMeanX_err > resSegmTrack_mean_min && dtMeanX - nsig * dtMeanX_err < resSegmTrack_mean_max)
878  muonIdSummaryMap->setBinContent(2 * i + 1, 2, 1);
879  else
880  muonIdSummaryMap->setBinContent(2 * i + 1, 2, 0);
881 
882  if (dtSigmaX + nsig * dtSigmaX_err > resSegmTrack_rms_min &&
883  dtSigmaX - nsig * dtSigmaX_err < resSegmTrack_rms_max)
884  muonIdSummaryMap->setBinContent(2 * i + 1, 3, 1);
885  else
886  muonIdSummaryMap->setBinContent(2 * i + 1, 3, 0);
887  }
888  if (numPlot_dtY > 0) {
889  if (dtMeanY + nsig * dtMeanY_err > resSegmTrack_mean_min && dtMeanY - nsig * dtMeanY_err < resSegmTrack_mean_max)
890  muonIdSummaryMap->setBinContent(2 * i + 1, 4, 1);
891  else
892  muonIdSummaryMap->setBinContent(2 * i + 1, 4, 0);
893 
894  if (dtSigmaY + nsig * dtSigmaY_err > resSegmTrack_rms_min &&
895  dtSigmaY - nsig * dtSigmaX_err < resSegmTrack_rms_max)
896  muonIdSummaryMap->setBinContent(2 * i + 1, 5, 1);
897  else
898  muonIdSummaryMap->setBinContent(2 * i + 1, 5, 0);
899  }
900 
901  if (numPlot_cscX > 0) {
902  if (cscMeanX + nsig * cscMeanX_err > resSegmTrack_mean_min &&
903  cscMeanX - nsig * cscMeanX_err < resSegmTrack_mean_max)
904  muonIdSummaryMap->setBinContent(2 * i + 2, 2, 1);
905  else
906  muonIdSummaryMap->setBinContent(2 * i + 2, 2, 0);
907 
908  if (cscSigmaX + nsig * cscSigmaX_err > resSegmTrack_rms_min &&
909  cscSigmaX - nsig * cscSigmaX_err < resSegmTrack_rms_max)
910  muonIdSummaryMap->setBinContent(2 * i + 2, 3, 1);
911  else
912  muonIdSummaryMap->setBinContent(2 * i + 2, 3, 0);
913  }
914  if (numPlot_cscY > 0) {
915  if (cscMeanY + nsig * cscMeanY_err > resSegmTrack_mean_min &&
916  cscMeanY - nsig * cscMeanY_err < resSegmTrack_mean_max)
917  muonIdSummaryMap->setBinContent(2 * i + 2, 4, 1);
918  else
919  muonIdSummaryMap->setBinContent(2 * i + 2, 4, 0);
920 
921  if (cscSigmaY + nsig * cscSigmaY_err > resSegmTrack_rms_min &&
922  cscSigmaY - nsig * cscSigmaY_err < resSegmTrack_rms_max)
923  muonIdSummaryMap->setBinContent(2 * i + 2, 5, 1);
924  else
925  muonIdSummaryMap->setBinContent(2 * i + 2, 5, 0);
926  }
927 
928  //---- end of modification
929  }
930 }
931 
932 void MuonTestSummary::ResidualCheck(DQMStore::IGetter &igetter,
933  std::string muType,
934  const std::vector<std::string> &resHistos,
935  int &numPlot,
936  double &Mean,
937  double &Mean_err,
938  double &Sigma,
939  double &Sigma_err) {
940  numPlot = 0;
941  Mean = 0;
942  Mean_err = 0;
943  Sigma = 0;
944  Sigma_err = 0;
945  for (uint name = 0; name < resHistos.size(); name++) {
946  MonitorElement *resHisto = igetter.get("Muons/MuonIdDQM/" + muType + "/" + resHistos[name]);
947 
948  if (resHisto) {
949  TH1F *resHisto_root = resHisto->getTH1F();
950  if (resHisto_root->GetEntries() < 20) {
951  LogTrace(metname) << "[MuonTestSummary]: Test of " << muType << " for " << resHistos[name]
952  << " not performed because # entries < 20 ";
953  continue;
954  }
955 
956  //we also want to check if the peak is away from zero.
957  //so, try fitting in 3 sigma around the histogram mean.
958  //alternatively, could use the maximum bin (less sensitive to 1-sided tails).
959  // float mean = resHisto_root->GetMean();
960  float mean = resHisto_root->GetBinLowEdge(resHisto_root->GetMaximumBin());
961  TF1 *gfit = new TF1("Gaussian", "gaus", mean - 3, mean + 3);
962 
963  try {
964  resHisto_root->Fit(gfit, "Q0");
965  } catch (cms::Exception &iException) {
966  edm::LogError(metname) << "[MuonTestSummary]: Exception when fitting " << resHistos[name];
967  continue;
968  }
969  if (gfit) {
970  double mean = gfit->GetParameter(1);
971  double mean_err = gfit->GetParError(1);
972  double sigma = gfit->GetParameter(2);
973  double sigma_err = gfit->GetParError(2);
974  LogTrace(metname) << "meanRes: " << mean << " +- " << mean_err << " for " << resHistos[name] << endl;
975  LogTrace(metname) << "sigmaRes: " << sigma << " +- " << sigma_err << " for " << resHistos[name] << endl;
976 
977  Mean += mean;
978  Mean_err += mean_err * mean_err;
979  Sigma += sigma;
980  Sigma_err += sigma_err * sigma_err;
981  numPlot++;
982  } //if gfit? why would we not have gfit?
983 
984  } //histogram exists...
985  } // loop over residuals histos
986 
987  if (numPlot == 0) { //eg no stats
988  Mean_err = 1;
989  Mean = 1;
990  Sigma_err = 1;
991  Sigma = 1;
992  } else {
993  Mean_err = sqrt(Mean_err);
994  Mean_err /= numPlot;
995  Mean /= numPlot;
996 
997  Sigma_err = sqrt(Sigma_err);
998  Sigma_err /= numPlot;
999  Sigma /= numPlot;
1000  }
1001  return;
1002 }
1003 void MuonTestSummary::doEnergyTests(DQMStore::IGetter &igetter, string histname, string muonType, int binNumber) {
1004  // num matches test
1005  string path = "Muons/MuonEnergyDepositAnalyzer/" + histname + muonType;
1006  MonitorElement *energyHisto = igetter.get(path);
1007  Double_t hPeak = -1, hFWHM = -1;
1008  if (energyHisto) {
1009  TH1F *energyHisto_root = energyHisto->getTH1F();
1010 
1011  // Setting fit range and start values
1012  Double_t fitRange[2];
1013  Double_t startValues[4], parlimitslo[4], parlimitshi[4], fitPar[4], fitParErr[4];
1014 
1015  if (energyHisto->getEntries() > 20) {
1016  if (histname == "ecalS9PointingMuDepositedEnergy_") {
1017  fitRange[0] = 0.04;
1018  fitRange[1] = 3.0;
1019 
1020  startValues[0] = 0.036;
1021  startValues[1] = 0.193;
1022  startValues[2] = 110.0;
1023  startValues[3] = 0.06;
1024  parlimitslo[0] = 0.0;
1025  parlimitslo[1] = 0.;
1026  parlimitslo[2] = 1.0;
1027  parlimitslo[3] = 0.;
1028  parlimitshi[0] = 0.05;
1029  parlimitshi[1] = 0.5;
1030  parlimitshi[2] = 80000.0;
1031  parlimitshi[3] = 0.1;
1032 
1033  Double_t chisqr;
1034  Int_t ndf;
1035  TF1 *fit = langaufit(
1036  energyHisto_root, fitRange, startValues, parlimitslo, parlimitshi, fitPar, fitParErr, &chisqr, &ndf);
1037  if (fit) {
1038  langaupro(fitPar, hPeak, hFWHM);
1039  LogTrace(metname) << "hPeak from langau fit: " << hPeak << " for: " << histname + muonType << endl;
1040  LogTrace(metname) << "hFWHM from langau fit: " << hFWHM << " for: " << histname + muonType << endl;
1041  }
1042  }
1043 
1044  if (histname == "hadS9PointingMuDepositedEnergy_") {
1045  fitRange[0] = 0.0;
1046  fitRange[1] = 7.0;
1047 
1048  startValues[0] = 2.0;
1049  startValues[1] = 2.4;
1050  startValues[2] = 110.0;
1051  startValues[3] = 4.0;
1052  parlimitslo[0] = 0.0;
1053  parlimitslo[1] = 0.;
1054  parlimitslo[2] = 1.0;
1055  parlimitslo[3] = 0.;
1056  parlimitshi[0] = 4.0;
1057  parlimitshi[1] = 4.0;
1058  parlimitshi[2] = 80000.0;
1059  parlimitshi[3] = 8.0;
1060 
1061  Double_t chisqr;
1062  Int_t ndf;
1063  TF1 *fit = langaufit(
1064  energyHisto_root, fitRange, startValues, parlimitslo, parlimitshi, fitPar, fitParErr, &chisqr, &ndf);
1065  if (fit) {
1066  langaupro(fitPar, hPeak, hFWHM);
1067  LogTrace(metname) << "hPeak from langau fit: " << hPeak << " for: " << histname + muonType << endl;
1068  LogTrace(metname) << "hFWHM from langau fit: " << hFWHM << " for: " << histname + muonType << endl;
1069  }
1070  }
1071  } else {
1072  LogTrace(metname) << "[MuonTestSummary]: Test of Energy for " << histname + muonType
1073  << " not performed because # entries < 20 ";
1074  }
1075  }
1076 
1077  if (histname == "ecalS9PointingMuDepositedEnergy_" && hPeak > expPeakEcalS9_min && hPeak < expPeakEcalS9_max)
1078  energySummaryMap->setBinContent(binNumber, 1, 1);
1079  if (histname == "ecalS9PointingMuDepositedEnergy_" && (hPeak != -1) &&
1080  !(hPeak > expPeakEcalS9_min && hPeak < expPeakEcalS9_max))
1081  energySummaryMap->setBinContent(binNumber, 1, 0);
1082 
1083  if (histname == "hadS9PointingMuDepositedEnergy_" && hPeak > expPeakHadS9_min && hPeak < expPeakHadS9_max)
1084  energySummaryMap->setBinContent(binNumber, 2, 1);
1085  if (histname == "hadS9PointingMuDepositedEnergy_" && (hPeak != -1) &&
1086  !(hPeak > expPeakHadS9_min && hPeak < expPeakHadS9_max))
1087  energySummaryMap->setBinContent(binNumber, 2, 0);
1088 
1089  //missing test on ho distributions
1090 }
1091 void MuonTestSummary::doMultiplicityTests(DQMStore::IGetter &igetter) {
1092  MonitorElement *multiplicityHisto = igetter.get("Muons/MuonRecoAnalyzer/muReco");
1093 
1094  if (multiplicityHisto) {
1095  if (multiplicityHisto->getEntries() > 20) {
1096  double multiplicity_GLB = double(multiplicityHisto->getBinContent(1) + multiplicityHisto->getBinContent(2)) /
1097  double(multiplicityHisto->getEntries());
1098  LogTrace(metname) << "multiplicity_GLB: " << multiplicity_GLB << " ExpMultiplicityGlb_min "
1099  << expMultiplicityGlb_min << " ExpMultiplicityGlb_max " << expMultiplicityGlb_max << endl;
1100  double multiplicity_TK = double(multiplicityHisto->getBinContent(3) + multiplicityHisto->getBinContent(4)) /
1101  double(multiplicityHisto->getEntries());
1102  LogTrace(metname) << "multiplicity_TK: " << multiplicity_TK << " ExpMultiplicityTk_min " << expMultiplicityTk_min
1103  << " ExpMultiplicityTk_max " << expMultiplicityTk_max << endl;
1104  double multiplicity_STA = double(multiplicityHisto->getBinContent(3) + multiplicityHisto->getBinContent(5)) /
1105  double(multiplicityHisto->getEntries());
1106  LogTrace(metname) << "multiplicity_STA: " << multiplicity_STA << " ExpMultiplicitySta_min "
1107  << expMultiplicitySta_min << " ExpMultiplicitySta_max " << expMultiplicitySta_max << endl;
1108 
1109  if (multiplicity_GLB > expMultiplicityGlb_min && multiplicity_GLB < expMultiplicityGlb_max)
1110  multiplicitySummaryMap->setBinContent(1, 1);
1111  else
1112  multiplicitySummaryMap->setBinContent(1, 0);
1113 
1114  if (multiplicity_TK > expMultiplicityTk_min && multiplicity_TK < expMultiplicityTk_max)
1115  multiplicitySummaryMap->setBinContent(2, 1);
1116  else
1117  multiplicitySummaryMap->setBinContent(2, 0);
1118 
1119  if (multiplicity_STA > expMultiplicitySta_min && multiplicity_STA < expMultiplicitySta_max)
1120  multiplicitySummaryMap->setBinContent(3, 1);
1121  else
1122  multiplicitySummaryMap->setBinContent(3, 0);
1123  } else {
1124  LogTrace(metname) << "[MuonTestSummary]: Test of Multiplicity not performed because # entries < 20 ";
1125  }
1126  }
1127 }
dqm::implementation::IBooker::bookFloat
MonitorElement * bookFloat(TString const &name, FUNC onbooking=NOOP())
Definition: DQMStore.h:80
mps_fire.i
i
Definition: mps_fire.py:428
dqm::implementation::NavigatorBase::setCurrentFolder
virtual void setCurrentFolder(std::string const &fullpath)
Definition: DQMStore.cc:32
MuonTestSummary::ResidualCheck
void ResidualCheck(DQMStore::IGetter &, std::string muType, const std::vector< std::string > &resHistos, int &numPlot, double &Mean, double &Mean_err, double &Sigma, double &Sigma_err)
Definition: MuonTestSummary.cc:932
cms::cuda::allocator::maxBin
constexpr unsigned int maxBin
Definition: cachingAllocatorCommon.h:16
metname
const std::string metname
Definition: MuonSeedOrcaPatternRecognition.cc:40
cuy.yBin
tuple yBin
Definition: cuy.py:892
dqm::legacy::MonitorElement::getTH1F
virtual TH1F * getTH1F() const
Definition: MonitorElement.h:479
HLT_FULL_cff.fraction
tuple fraction
Definition: HLT_FULL_cff.py:53121
MuonTestSummary::doMultiplicityTests
void doMultiplicityTests(DQMStore::IGetter &)
Definition: MuonTestSummary.cc:1091
edm::LogError
Log< level::Error, false > LogError
Definition: MessageLogger.h:123
type
type
Definition: SiPixelVCal_PayloadInspector.cc:39
AlCaHLTBitMon_QueryRunRegistry.string
string string
Definition: AlCaHLTBitMon_QueryRunRegistry.py:256
LogTrace
#define LogTrace(id)
Definition: MessageLogger.h:234
dqm::impl::MonitorElement::Fill
void Fill(long long x)
Definition: MonitorElement.h:290
MuonTestSummary::doMuonIDTests
void doMuonIDTests(DQMStore::IGetter &)
Definition: MuonTestSummary.cc:729
gpuVertexFinder::fit
auto const & fit
Definition: gpuVertexFinder.cc:25
dqm::impl::MonitorElement::getEntries
virtual double getEntries() const
get # of entries
Definition: MonitorElement.cc:619
mathSSE::sqrt
T sqrt(T t)
Definition: SSEVec.h:19
dqm::implementation::IGetter::get
virtual MonitorElement * get(std::string const &fullpath) const
Definition: DQMStore.cc:673
dqm::impl::MonitorElement::getBinContent
virtual double getBinContent(int binx) const
get content of bin (1-D)
Definition: MonitorElement.cc:583
dqm::impl::MonitorElement::getMean
virtual double getMean(int axis=1) const
get mean value of histogram along x, y or z axis (axis=1, 2, 3 respectively)
Definition: MonitorElement.cc:540
MuonTestSummary::~MuonTestSummary
~MuonTestSummary() override
Destructor.
Definition: MuonTestSummary.cc:65
fed_dqm_sourceclient-live_cfg.path
string path
Definition: fed_dqm_sourceclient-live_cfg.py:47
dqm::implementation::IBooker::book2D
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:177
dqm::impl::MonitorElement::getRMS
virtual double getRMS(int axis=1) const
get RMS of histogram along x, y or z axis (axis=1, 2, 3 respectively)
Definition: MonitorElement.cc:553
edm::ParameterSet::getParameter
T getParameter(std::string const &) const
Definition: ParameterSet.h:303
MuonTestSummary::doKinematicsTests
void doKinematicsTests(DQMStore::IGetter &, std::string, int)
test operations
Definition: MuonTestSummary.cc:448
MuonTestSummary::doEnergyTests
void doEnergyTests(DQMStore::IGetter &, std::string nameHisto, std::string muonType, int bin)
Definition: MuonTestSummary.cc:1003
MuonTestSummary::GaussFit
void GaussFit(std::string type, std::string parameter, MonitorElement *Histo, float &mean, float &mean_err, float &sigma, float &sigma_err)
Definition: MuonTestSummary.cc:556
MuonTestSummary::dqmEndJob
void dqmEndJob(DQMStore::IBooker &, DQMStore::IGetter &) override
Definition: MuonTestSummary.cc:67
dqm::implementation::IBooker::book1D
MonitorElement * book1D(TString const &name, TString const &title, int const nchX, double const lowX, double const highX, FUNC onbooking=NOOP())
Definition: DQMStore.h:98
MuonTestSummary::MuonTestSummary
MuonTestSummary(const edm::ParameterSet &ps)
Constructor.
Definition: MuonTestSummary.cc:25
MuonTestSummary::doResidualsTests
void doResidualsTests(DQMStore::IGetter &, std::string, std::string, int)
Definition: MuonTestSummary.cc:594
dqm::impl::MonitorElement::setAxisTitle
virtual void setAxisTitle(const std::string &title, int axis=1)
set x-, y- or z-axis title (axis=1, 2, 3 respectively)
Definition: MonitorElement.cc:791
Generated for CMSSW Reference Manual by doxygen 1.8.5