CMS 3D CMS Logo

MuonRecoAnalyzer.cc
Go to the documentation of this file.
2 
10 
12 
13 #include <string>
14 #include "TMath.h"
15 using namespace std;
16 using namespace edm;
17 
19  parameters = pSet;
20 
21  // Input booleans
22  IsminiAOD = parameters.getParameter<bool>("IsminiAOD");
23  doMVA = parameters.getParameter<bool>("doMVA");
24  // the services:
25  theService = new MuonServiceProxy(parameters.getParameter<ParameterSet>("ServiceParameters"));
26  theMuonCollectionLabel_ = consumes<edm::View<reco::Muon> >(parameters.getParameter<edm::InputTag>("MuonCollection"));
27  theVertexLabel_ = consumes<reco::VertexCollection>(pSet.getParameter<InputTag>("inputTagVertex"));
28  theBeamSpotLabel_ = consumes<reco::BeamSpot>(pSet.getParameter<InputTag>("inputTagBeamSpot"));
29  dcsStatusCollection_ =
30  consumes<DcsStatusCollection>(pSet.getUntrackedParameter<std::string>("dcsStatusCollection", "scalersRawToDigi"));
31 
32  ptBin = parameters.getParameter<int>("ptBin");
33  ptMin = parameters.getParameter<double>("ptMin");
34  ptMax = parameters.getParameter<double>("ptMax");
35  pResBin = parameters.getParameter<int>("pResBin");
36  pResMin = parameters.getParameter<double>("pResMin");
37  pResMax = parameters.getParameter<double>("pResMax");
38  rhBin = parameters.getParameter<int>("rhBin");
39  rhMin = parameters.getParameter<double>("rhMin");
40  rhMax = parameters.getParameter<double>("rhMax");
41  pBin = parameters.getParameter<int>("pBin");
42  pMin = parameters.getParameter<double>("pMin");
43  pMax = parameters.getParameter<double>("pMax");
44  chi2Bin = parameters.getParameter<int>("chi2Bin");
45  chi2Min = parameters.getParameter<double>("chi2Min");
46  chi2Max = parameters.getParameter<double>("chi2Max");
47  phiBin = parameters.getParameter<int>("phiBin");
48  phiMin = parameters.getParameter<double>("phiMin");
49  phiMax = parameters.getParameter<double>("phiMax");
50  tunePBin = parameters.getParameter<int>("tunePBin");
51  tunePMax = parameters.getParameter<double>("tunePMax");
52  tunePMin = parameters.getParameter<double>("tunePMin");
53  thetaBin = parameters.getParameter<int>("thetaBin");
54  thetaMin = parameters.getParameter<double>("thetaMin");
55  thetaMax = parameters.getParameter<double>("thetaMax");
56  etaBin = parameters.getParameter<int>("etaBin");
57  etaMin = parameters.getParameter<double>("etaMin");
58  etaMax = parameters.getParameter<double>("etaMax");
59 
60  theFolder = parameters.getParameter<string>("folder");
61 }
62 
63 MuonRecoAnalyzer::~MuonRecoAnalyzer() { delete theService; }
65  edm::Run const& /*iRun*/,
66  edm::EventSetup const& /* iSetup */) {
67  ibooker.cd();
68  ibooker.setCurrentFolder(theFolder);
69 
70  muReco = ibooker.book1D("muReco", "muon reconstructed tracks", 6, 1, 7);
71  muReco->setBinLabel(1, "glb+tk+sta");
72  muReco->setBinLabel(2, "glb+sta");
73  muReco->setBinLabel(3, "tk+sta");
74  muReco->setBinLabel(4, "tk");
75  muReco->setBinLabel(5, "sta");
76  muReco->setBinLabel(6, "calo");
77 
78  int binFactor = 4;
79 
81  // monitoring of eta parameter
83  std::string histname = "GlbMuon_";
84  etaGlbTrack.push_back(ibooker.book1D(histname + "Glb_eta", "#eta_{GLB}", etaBin, etaMin, etaMax));
85  etaGlbTrack.push_back(ibooker.book1D(histname + "Tk_eta", "#eta_{TKfromGLB}", etaBin, etaMin, etaMax));
86  etaGlbTrack.push_back(ibooker.book1D(histname + "Sta_eta", "#eta_{STAfromGLB}", etaBin, etaMin, etaMax));
87  etaResolution.push_back(ibooker.book1D(
88  "Res_TkGlb_eta", "#eta_{TKfromGLB} - #eta_{GLB}", etaBin * binFactor, etaMin / 3000, etaMax / 3000));
89  etaResolution.push_back(ibooker.book1D(
90  "Res_GlbSta_eta", "#eta_{GLB} - #eta_{STAfromGLB}", etaBin * binFactor, etaMin / 100, etaMax / 100));
91  etaResolution.push_back(ibooker.book1D(
92  "Res_TkSta_eta", "#eta_{TKfromGLB} - #eta_{STAfromGLB}", etaBin * binFactor, etaMin / 100, etaMax / 100));
93  etaResolution.push_back(ibooker.book2D("ResVsEta_TkGlb_eta",
94  "(#eta_{TKfromGLB} - #eta_{GLB}) vs #eta_{GLB}",
95  etaBin,
96  etaMin,
97  etaMax,
98  etaBin * binFactor,
99  etaMin / 3000,
100  etaMax / 3000));
101  etaResolution.push_back(ibooker.book2D("ResVsEta_GlbSta_eta",
102  "(#eta_{GLB} - #eta_{STAfromGLB}) vs #eta_{GLB}",
103  etaBin,
104  etaMin,
105  etaMax,
106  etaBin * binFactor,
107  etaMin / 100,
108  etaMax / 100));
109  etaResolution.push_back(ibooker.book2D("ResVsEta_TkSta_eta",
110  "(#eta_{TKfromGLB} - #eta_{STAfromGLB}) vs #eta_{TKfromGLB}",
111  etaBin,
112  etaMin,
113  etaMax,
114  etaBin * binFactor,
115  etaMin / 100,
116  etaMax / 100));
117  etaPull = ibooker.book1D("Pull_TkSta_eta", "#eta_{TKfromGLB} - #eta_{GLB} / error", 100, -10, 10);
118  etaTrack = ibooker.book1D("TkMuon_eta", "#eta_{TK}", etaBin, etaMin, etaMax);
119  etaStaTrack = ibooker.book1D("StaMuon_eta", "#eta_{STA}", etaBin, etaMin, etaMax);
120  etaEfficiency.push_back(ibooker.book1D("StaEta", "#eta_{STAfromGLB}", etaBin, etaMin, etaMax));
121  etaEfficiency.push_back(
122  ibooker.book1D("StaEta_ifCombinedAlso", "#eta_{STAfromGLB} if isGlb=true", etaBin, etaMin, etaMax));
123 
125  // monitoring of theta parameter
127  thetaGlbTrack.push_back(ibooker.book1D(histname + "Glb_theta", "#theta_{GLB}", thetaBin, thetaMin, thetaMax));
128  thetaGlbTrack[0]->setAxisTitle("rad");
129  thetaGlbTrack.push_back(ibooker.book1D(histname + "Tk_theta", "#theta_{TKfromGLB}", thetaBin, thetaMin, thetaMax));
130  thetaGlbTrack[1]->setAxisTitle("rad");
131  thetaGlbTrack.push_back(ibooker.book1D(histname + "Sta_theta", "#theta_{STAfromGLB}", thetaBin, thetaMin, thetaMax));
132  thetaGlbTrack[2]->setAxisTitle("rad");
133  thetaResolution.push_back(ibooker.book1D("Res_TkGlb_theta",
134  "#theta_{TKfromGLB} - #theta_{GLB}",
135  thetaBin * binFactor,
136  -(thetaMax / 3000),
137  thetaMax / 3000));
138  thetaResolution[0]->setAxisTitle("rad");
139  thetaResolution.push_back(ibooker.book1D("Res_GlbSta_theta",
140  "#theta_{GLB} - #theta_{STAfromGLB}",
141  thetaBin * binFactor,
142  -(thetaMax / 100),
143  thetaMax / 100));
144  thetaResolution[1]->setAxisTitle("rad");
145  thetaResolution.push_back(ibooker.book1D("Res_TkSta_theta",
146  "#theta_{TKfromGLB} - #theta_{STAfromGLB}",
147  thetaBin * binFactor,
148  -(thetaMax / 100),
149  thetaMax / 100));
150  thetaResolution[2]->setAxisTitle("rad");
151  thetaResolution.push_back(ibooker.book2D("ResVsTheta_TkGlb_theta",
152  "(#theta_{TKfromGLB} - #theta_{GLB}) vs #theta_{GLB}",
153  thetaBin,
154  thetaMin,
155  thetaMax,
156  thetaBin * binFactor,
157  -(thetaMax / 3000),
158  thetaMax / 3000));
159  thetaResolution[3]->setAxisTitle("rad", 1);
160  thetaResolution[3]->setAxisTitle("rad", 2);
161  thetaResolution.push_back(ibooker.book2D("ResVsTheta_GlbSta_theta",
162  "(#theta_{GLB} - #theta_{STAfromGLB}) vs #theta_{GLB}",
163  thetaBin,
164  thetaMin,
165  thetaMax,
166  thetaBin * binFactor,
167  -(thetaMax / 100),
168  thetaMax / 100));
169  thetaResolution[4]->setAxisTitle("rad", 1);
170  thetaResolution[4]->setAxisTitle("rad", 2);
171  thetaResolution.push_back(ibooker.book2D("ResVsTheta_TkSta_theta",
172  "(#theta_{TKfromGLB} - #theta_{STAfromGLB}) vs #theta_{TKfromGLB}",
173  thetaBin,
174  thetaMin,
175  thetaMax,
176  thetaBin * binFactor,
177  -(thetaMax / 100),
178  thetaMax / 100));
179  thetaResolution[5]->setAxisTitle("rad", 1);
180  thetaResolution[5]->setAxisTitle("rad", 2);
181  thetaPull = ibooker.book1D("Pull_TkSta_theta", "#theta_{TKfromGLB} - #theta_{STAfromGLB} / error", 100, -10, 10);
182  thetaTrack = ibooker.book1D("TkMuon_theta", "#theta_{TK}", thetaBin, thetaMin, thetaMax);
183  thetaTrack->setAxisTitle("rad");
184  thetaStaTrack = ibooker.book1D("StaMuon_theta", "#theta_{STA}", thetaBin, thetaMin, thetaMax);
185  thetaStaTrack->setAxisTitle("rad");
186 
187  // monitoring tunePMuonBestTrack Pt
188  tunePResolution = ibooker.book1D(
189  "Res_TuneP_pt", "Pt_{MuonBestTrack}-Pt_{tunePMuonBestTrack}/Pt_{MuonBestTrack}", tunePBin, tunePMin, tunePMax);
190 
191  // monitoring of phi paramater
192  phiGlbTrack.push_back(ibooker.book1D(histname + "Glb_phi", "#phi_{GLB}", phiBin, phiMin, phiMax));
193  phiGlbTrack[0]->setAxisTitle("rad");
194  phiGlbTrack.push_back(ibooker.book1D(histname + "Tk_phi", "#phi_{TKfromGLB}", phiBin, phiMin, phiMax));
195  phiGlbTrack[1]->setAxisTitle("rad");
196  phiGlbTrack.push_back(ibooker.book1D(histname + "Sta_phi", "#phi_{STAfromGLB}", phiBin, phiMin, phiMax));
197  phiGlbTrack[2]->setAxisTitle("rad");
198  phiResolution.push_back(ibooker.book1D(
199  "Res_TkGlb_phi", "#phi_{TKfromGLB} - #phi_{GLB}", phiBin * binFactor, phiMin / 3000, phiMax / 3000));
200  phiResolution[0]->setAxisTitle("rad");
201  phiResolution.push_back(ibooker.book1D(
202  "Res_GlbSta_phi", "#phi_{GLB} - #phi_{STAfromGLB}", phiBin * binFactor, phiMin / 100, phiMax / 100));
203  phiResolution[1]->setAxisTitle("rad");
204  phiResolution.push_back(ibooker.book1D(
205  "Res_TkSta_phi", "#phi_{TKfromGLB} - #phi_{STAfromGLB}", phiBin * binFactor, phiMin / 100, phiMax / 100));
206  phiResolution[2]->setAxisTitle("rad");
207  phiResolution.push_back(ibooker.book2D("ResVsPhi_TkGlb_phi",
208  "(#phi_{TKfromGLB} - #phi_{GLB}) vs #phi_GLB",
209  phiBin,
210  phiMin,
211  phiMax,
212  phiBin * binFactor,
213  phiMin / 3000,
214  phiMax / 3000));
215  phiResolution[3]->setAxisTitle("rad", 1);
216  phiResolution[3]->setAxisTitle("rad", 2);
217  phiResolution.push_back(ibooker.book2D("ResVsPhi_GlbSta_phi",
218  "(#phi_{GLB} - #phi_{STAfromGLB}) vs #phi_{GLB}",
219  phiBin,
220  phiMin,
221  phiMax,
222  phiBin * binFactor,
223  phiMin / 100,
224  phiMax / 100));
225  phiResolution[4]->setAxisTitle("rad", 1);
226  phiResolution[4]->setAxisTitle("rad", 2);
227  phiResolution.push_back(ibooker.book2D("ResVsPhi_TkSta_phi",
228  "(#phi_{TKfromGLB} - #phi_{STAfromGLB}) vs #phi_{TKfromGLB}",
229  phiBin,
230  phiMin,
231  phiMax,
232  phiBin * binFactor,
233  phiMin / 100,
234  phiMax / 100));
235  phiResolution[5]->setAxisTitle("rad", 1);
236  phiResolution[5]->setAxisTitle("rad", 2);
237  phiPull = ibooker.book1D("Pull_TkSta_phi", "#phi_{TKfromGLB} - #phi_{STAfromGLB} / error", 100, -10, 10);
238  phiTrack = ibooker.book1D("TkMuon_phi", "#phi_{TK}", phiBin, phiMin, phiMax);
239  phiTrack->setAxisTitle("rad");
240  phiStaTrack = ibooker.book1D("StaMuon_phi", "#phi_{STA}", phiBin, phiMin, phiMax);
241  phiStaTrack->setAxisTitle("rad");
242  phiEfficiency.push_back(ibooker.book1D("StaPhi", "#phi_{STAfromGLB}", phiBin, phiMin, phiMax));
243  phiEfficiency[0]->setAxisTitle("rad");
244  phiEfficiency.push_back(
245  ibooker.book1D("StaPhi_ifCombinedAlso", "#phi_{STAfromGLB} if the isGlb=true", phiBin, phiMin, phiMax));
246  phiEfficiency[1]->setAxisTitle("rad");
247 
248  // monitoring of the chi2 parameter
249  chi2OvDFGlbTrack.push_back(
250  ibooker.book1D(histname + "Glb_chi2OverDf", "#chi_{2}OverDF_{GLB}", chi2Bin, chi2Min, chi2Max));
251  chi2OvDFGlbTrack.push_back(
252  ibooker.book1D(histname + "Tk_chi2OverDf", "#chi_{2}OverDF_{TKfromGLB}", phiBin, chi2Min, chi2Max));
253  chi2OvDFGlbTrack.push_back(
254  ibooker.book1D(histname + "Sta_chi2OverDf", "#chi_{2}OverDF_{STAfromGLB}", chi2Bin, chi2Min, chi2Max));
255  chi2OvDFTrack = ibooker.book1D("TkMuon_chi2OverDf", "#chi_{2}OverDF_{TK}", chi2Bin, chi2Min, chi2Max);
256  chi2OvDFStaTrack = ibooker.book1D("StaMuon_chi2OverDf", "#chi_{2}OverDF_{STA}", chi2Bin, chi2Min, chi2Max);
257  //--------------------------
258  probchi2GlbTrack.push_back(ibooker.book1D(histname + "Glb_probchi", "Prob #chi_{GLB}", 120, chi2Min, 1.20));
259  probchi2GlbTrack.push_back(ibooker.book1D(histname + "Tk_probchi", "Prob #chi_{TKfromGLB}", 120, chi2Min, 1.20));
260  probchi2GlbTrack.push_back(ibooker.book1D(histname + "Sta_probchi", "Prob #chi_{STAfromGLB}", 120, chi2Min, 1.20));
261  probchi2Track = ibooker.book1D("TkMuon_probchi", "Prob #chi_{TK}", 120, chi2Min, 1.20);
262  probchi2StaTrack = ibooker.book1D("StaMuon_probchi", "Prob #chi_{STA}", 120, chi2Min, 1.20);
263 
264  // monitoring of the momentum
265  pGlbTrack.push_back(ibooker.book1D(histname + "Glb_p", "p_{GLB}", pBin, pMin, pMax));
266  pGlbTrack[0]->setAxisTitle("GeV");
267  pGlbTrack.push_back(ibooker.book1D(histname + "Tk_p", "p_{TKfromGLB}", pBin, pMin, pMax));
268  pGlbTrack[1]->setAxisTitle("GeV");
269  pGlbTrack.push_back(ibooker.book1D(histname + "Sta_p", "p_{STAfromGLB}", pBin, pMin, pMax));
270  pGlbTrack[2]->setAxisTitle("GeV");
271  pTrack = ibooker.book1D("TkMuon_p", "p_{TK}", pBin, pMin, pMax);
272  pTrack->setAxisTitle("GeV");
273  pStaTrack = ibooker.book1D("StaMuon_p", "p_{STA}", pBin, pMin, pMax);
274  pStaTrack->setAxisTitle("GeV");
275 
276  // monitoring of the transverse momentum
277  ptGlbTrack.push_back(ibooker.book1D(histname + "Glb_pt", "pt_{GLB}", ptBin, ptMin, ptMax));
278  ptGlbTrack[0]->setAxisTitle("GeV");
279  ptGlbTrack.push_back(ibooker.book1D(histname + "Tk_pt", "pt_{TKfromGLB}", ptBin, ptMin, ptMax));
280  ptGlbTrack[1]->setAxisTitle("GeV");
281  ptGlbTrack.push_back(ibooker.book1D(histname + "Sta_pt", "pt_{STAfromGLB}", ptBin, ptMin, ptMax));
282  ptGlbTrack[2]->setAxisTitle("GeV");
283  ptTrack = ibooker.book1D("TkMuon_pt", "pt_{TK}", ptBin, ptMin, ptMax);
284  ptTrack->setAxisTitle("GeV");
285  ptStaTrack = ibooker.book1D("StaMuon_pt", "pt_{STA}", ptBin, ptMin, pMax);
286  ptStaTrack->setAxisTitle("GeV");
287 
288  //monitoring of variables needed by the MVA soft muon
289 
290  ptSoftMuonMVA = ibooker.book1D("ptSoftMuonMVA", "pt_{SoftMuon}", 50, 0, 50);
291  deltaRSoftMuonMVA = ibooker.book1D("deltaRSoftMuonMVA", "#Delta R", 50, 0, 5);
292  gNchi2SoftMuonMVA = ibooker.book1D("gNchi2SoftMuonMVA", "gNchi2", 50, 0, 3);
293  vMuHitsSoftMuonMVA = ibooker.book1D("vMuHitsSoftMuonMVA", "vMuHits", 50, 0, 50);
294  mNuStationsSoftMuonMVA = ibooker.book1D("mNuStationsSoftMuonMVA", "mNuStations", 6, 0, 6);
295  dxyRefSoftMuonMVA = ibooker.book1D("dxyRefSoftMuonMVA", "dxyRef", 50, -0.1, 0.1);
296  dzRefSoftMuonMVA = ibooker.book1D("dzRefSoftMuonMVA", "dzRef", 50, -0.1, 0.1);
297  LWHSoftMuonMVA = ibooker.book1D("LWHSoftMuonMVA", "LWH", 20, 0, 20);
298  valPixHitsSoftMuonMVA = ibooker.book1D("valPixHitsSoftMuonMVA", "valPixHits", 8, 0, 8);
299  innerChi2SoftMuonMVA = ibooker.book1D("innerChi2SoftMuonMVA", "innerChi2", 50, 0, 3);
300  outerChi2SoftMuonMVA = ibooker.book1D("outerChi2SoftMuonMVA", "outerChi2", 50, 0, 4);
301  iValFracSoftMuonMVA = ibooker.book1D("iValFracSoftMuonMVA", "iValFrac", 50, 0.5, 1.0);
302  segCompSoftMuonMVA = ibooker.book1D("segCompSoftMuonMVA", "segComp", 50, 0, 1.2);
303  chi2LocMomSoftMuonMVA = ibooker.book1D("chi2LocMomSoftMuonMVA", "chi2LocMom", 50, 0, 40);
304  chi2LocPosSoftMuonMVA = ibooker.book1D("chi2LocPosSoftMuonMVA", "chi2LocPos", 0, 0, 8);
305  glbTrackTailProbSoftMuonMVA = ibooker.book1D("glbTrackTailProbSoftMuonMVA", "glbTrackTailProb", 50, 0, 8);
306  NTrkVHitsSoftMuonMVA = ibooker.book1D("NTrkVHitsSoftMuonMVA", "NTrkVHits", 50, 0, 35);
307  kinkFinderSoftMuonMVA = ibooker.book1D("kinkFinderSoftMuonMVA", "kinkFinder", 50, 0, 30);
308  vRPChitsSoftMuonMVA = ibooker.book1D("vRPChitsSoftMuonMVA", "vRPChits", 50, 0, 50);
309  glbKinkFinderSoftMuonMVA = ibooker.book1D("glbKinkFinderSoftMuonMVA", "glbKinkFinder", 50, 0, 50);
310  glbKinkFinderLogSoftMuonMVA = ibooker.book1D("glbKinkFinderLogSoftMuonMVA", "glbKinkFinderLog", 50, 0, 50);
311  staRelChi2SoftMuonMVA = ibooker.book1D("staRelChi2SoftMuonMVA", "staRelChi2", 50, 0, 2);
312  glbDeltaEtaPhiSoftMuonMVA = ibooker.book1D("glbDeltaEtaPhiSoftMuonMVA", "glbDeltaEtaPhi", 50, 0, 0.15);
313  trkRelChi2SoftMuonMVA = ibooker.book1D("trkRelChi2SoftMuonMVA", "trkRelChi2", 50, 0, 1.2);
314  vDThitsSoftMuonMVA = ibooker.book1D("vDThitsSoftMuonMVA", "vDThits", 50, 0, 50);
315  vCSChitsSoftMuonMVA = ibooker.book1D("vCSChitsSoftMuonMVA", "vCSChits", 50, 0, 50);
316  timeAtIpInOutSoftMuonMVA = ibooker.book1D("timeAtIpInOutSoftMuonMVA", "timeAtIpInOut", 50, -10.0, 10.0);
317  timeAtIpInOutErrSoftMuonMVA = ibooker.book1D("timeAtIpInOutErrSoftMuonMVA", "timeAtIpInOutErr", 50, 0, 3.5);
318  getMuonHitsPerStationSoftMuonMVA =
319  ibooker.book1D("getMuonHitsPerStationSoftMuonMVA", "getMuonHitsPerStation", 6, 0, 6);
320  QprodSoftMuonMVA = ibooker.book1D("QprodSoftMuonMVA", "Qprod", 4, -2, 2);
321 
322  // monitoring of the muon charge
323  qGlbTrack.push_back(ibooker.book1D(histname + "Glb_q", "q_{GLB}", 5, -2.5, 2.5));
324  qGlbTrack.push_back(ibooker.book1D(histname + "Tk_q", "q_{TKfromGLB}", 5, -2.5, 2.5));
325  qGlbTrack.push_back(ibooker.book1D(histname + "Sta_q", "q_{STAformGLB}", 5, -2.5, 2.5));
326  qGlbTrack.push_back(ibooker.book1D(
327  histname + "qComparison", "comparison between q_{GLB} and q_{TKfromGLB}, q_{STAfromGLB}", 8, 0.5, 8.5));
328  qGlbTrack[3]->setBinLabel(1, "qGlb=qSta");
329  qGlbTrack[3]->setBinLabel(2, "qGlb!=qSta");
330  qGlbTrack[3]->setBinLabel(3, "qGlb=qTk");
331  qGlbTrack[3]->setBinLabel(4, "qGlb!=qTk");
332  qGlbTrack[3]->setBinLabel(5, "qSta=qTk");
333  qGlbTrack[3]->setBinLabel(6, "qSta!=qTk");
334  qGlbTrack[3]->setBinLabel(7, "qGlb!=qSta,qGlb!=Tk");
335  qGlbTrack[3]->setBinLabel(8, "qGlb=qSta,qGlb=Tk");
336  qTrack = ibooker.book1D("TkMuon_q", "q_{TK}", 5, -2.5, 2.5);
337  qStaTrack = ibooker.book1D("StaMuon_q", "q_{STA}", 5, -2.5, 2.5);
338 
340  // monitoring of the momentum resolution
341  qOverpResolution.push_back(ibooker.book1D(
342  "Res_TkGlb_qOverp", "(q/p)_{TKfromGLB} - (q/p)_{GLB}", pResBin * binFactor * 2, pResMin / 10, pResMax / 10));
343  qOverpResolution[0]->setAxisTitle("GeV^{-1}");
344  qOverpResolution.push_back(
345  ibooker.book1D("Res_GlbSta_qOverp", "(q/p)_{GLB} - (q/p)_{STAfromGLB}", pResBin * binFactor, pResMin, pResMax));
346  qOverpResolution[1]->setAxisTitle("GeV^{-1}");
347  qOverpResolution.push_back(ibooker.book1D(
348  "Res_TkSta_qOverp", "(q/p)_{TKfromGLB} - (q/p)_{STAfromGLB}", pResBin * binFactor, pResMin, pResMax));
349  qOverpResolution[2]->setAxisTitle("GeV^{-1}");
350  qOverpPull = ibooker.book1D("Pull_TkSta_qOverp", "(q/p)_{TKfromGLB} - (q/p)_{STAfromGLB} / error", 100, -10, 10);
351 
352  oneOverpResolution.push_back(ibooker.book1D(
353  "Res_TkGlb_oneOverp", "(1/p)_{TKfromGLB} - (1/p)_{GLB}", pResBin * binFactor * 2, pResMin / 10, pResMax / 10));
354  oneOverpResolution[0]->setAxisTitle("GeV^{-1}");
355  oneOverpResolution.push_back(
356  ibooker.book1D("Res_GlbSta_oneOverp", "(1/p)_{GLB} - (1/p)_{STAfromGLB}", pResBin * binFactor, pResMin, pResMax));
357  oneOverpResolution[1]->setAxisTitle("GeV^{-1}");
358  oneOverpResolution.push_back(ibooker.book1D(
359  "Res_TkSta_oneOverp", "(q/p)_{TKfromGLB} - (q/p)_{STAfromGLB}", pResBin * binFactor, pResMin, pResMax));
360  oneOverpResolution[2]->setAxisTitle("GeV^{-1}");
361  oneOverpPull = ibooker.book1D("Pull_TkSta_oneOverp", "(1/p)_{TKfromGLB} - (1/p)_{STAfromGLB} / error", 100, -10, 10);
362 
363  qOverptResolution.push_back(ibooker.book1D("Res_TkGlb_qOverpt",
364  "(q/p_{t})_{TKfromGLB} - (q/p_{t})_{GLB}",
365  pResBin * binFactor * 2,
366  pResMin / 10,
367  pResMax / 10));
368  qOverptResolution[0]->setAxisTitle("GeV^{-1}");
369  qOverptResolution.push_back(ibooker.book1D(
370  "Res_GlbSta_qOverpt", "(q/p_{t})_{GLB} - (q/p_{t})_{STAfromGLB}", pResBin * binFactor, pResMin, pResMax));
371  qOverptResolution[1]->setAxisTitle("GeV^{-1}");
372  qOverptResolution.push_back(ibooker.book1D(
373  "Res_TkSta_qOverpt", "(q/p_{t})_{TKfromGLB} - (q/p_{t})_{STAfromGLB}", pResBin * binFactor, pResMin, pResMax));
374  qOverptResolution[2]->setAxisTitle("GeV^{-1}");
375  qOverptPull = ibooker.book1D("Pull_TkSta_qOverpt", "(q/pt)_{TKfromGLB} - (q/pt)_{STAfromGLB} / error", 100, -10, 10);
376 
377  oneOverptResolution.push_back(ibooker.book1D("Res_TkGlb_oneOverpt",
378  "(1/p_{t})_{TKfromGLB} - (1/p_{t})_{GLB}",
379  pResBin * binFactor * 2,
380  pResMin / 10,
381  pResMax / 10));
382  oneOverptResolution[0]->setAxisTitle("GeV^{-1}");
383  oneOverptResolution.push_back(ibooker.book1D(
384  "Res_GlbSta_oneOverpt", "(1/p_{t})_{GLB} - (1/p_{t})_{STAfromGLB}", pResBin * binFactor, pResMin, pResMax));
385  oneOverptResolution[1]->setAxisTitle("GeV^{-1}");
386  oneOverptResolution.push_back(ibooker.book1D(
387  "Res_TkSta_oneOverpt", "(1/p_{t})_{TKfromGLB} - (1/p_{t})_{STAfromGLB}", pResBin * binFactor, pResMin, pResMax));
388  oneOverptResolution[2]->setAxisTitle("GeV^{-1}");
389  oneOverptResolution.push_back(ibooker.book2D("ResVsEta_TkGlb_oneOverpt",
390  "(#eta_{TKfromGLB} - #eta_{GLB}) vs (1/p_{t})_{GLB}",
391  etaBin,
392  etaMin,
393  etaMax,
394  pResBin * binFactor * 2,
395  pResMin / 10,
396  pResMax / 10));
397  oneOverptResolution[3]->setAxisTitle("GeV^{-1}", 2);
398  oneOverptResolution.push_back(ibooker.book2D("ResVsEta_GlbSta_oneOverpt",
399  "(#eta_{GLB} - #eta_{STAfromGLB} vs (1/p_{t})_{GLB}",
400  etaBin,
401  etaMin,
402  etaMax,
403  pResBin * binFactor,
404  pResMin,
405  pResMax));
406  oneOverptResolution[4]->setAxisTitle("GeV^{-1}", 2);
407  oneOverptResolution.push_back(ibooker.book2D("ResVsEta_TkSta_oneOverpt",
408  "(#eta_{TKfromGLB} - #eta_{STAfromGLB}) vs (1/p_{t})_{TKfromGLB}",
409  etaBin,
410  etaMin,
411  etaMax,
412  pResBin * binFactor,
413  pResMin,
414  pResMax));
415  oneOverptResolution[5]->setAxisTitle("GeV^{-1}", 2);
416  oneOverptResolution.push_back(ibooker.book2D("ResVsPhi_TkGlb_oneOverpt",
417  "(#phi_{TKfromGLB} - #phi_{GLB}) vs (1/p_{t})_{GLB}",
418  phiBin,
419  phiMin,
420  phiMax,
421  pResBin * binFactor * 2,
422  pResMin / 10,
423  pResMax / 10));
424  oneOverptResolution[6]->setAxisTitle("rad", 1);
425  oneOverptResolution[6]->setAxisTitle("GeV^{-1}", 2);
426  oneOverptResolution.push_back(ibooker.book2D("ResVsPhi_GlbSta_oneOverpt",
427  "(#phi_{GLB} - #phi_{STAfromGLB} vs (1/p_{t})_{GLB}",
428  phiBin,
429  phiMin,
430  phiMax,
431  pResBin * binFactor,
432  pResMin,
433  pResMax));
434  oneOverptResolution[7]->setAxisTitle("rad", 1);
435  oneOverptResolution[7]->setAxisTitle("GeV^{-1}", 2);
436  oneOverptResolution.push_back(ibooker.book2D("ResVsPhi_TkSta_oneOverpt",
437  "(#phi_{TKfromGLB} - #phi_{STAfromGLB}) vs (1/p_{t})_{TKfromGLB}",
438  phiBin,
439  phiMin,
440  phiMax,
441  pResBin * binFactor,
442  pResMin,
443  pResMax));
444  oneOverptResolution[8]->setAxisTitle("rad", 1);
445  oneOverptResolution[8]->setAxisTitle("GeV^{-1}", 2);
446  oneOverptResolution.push_back(ibooker.book2D("ResVsPt_TkGlb_oneOverpt",
447  "((1/p_{t})_{TKfromGLB} - (1/p_{t})_{GLB}) vs (1/p_{t})_{GLB}",
448  ptBin / 5,
449  ptMin,
450  ptMax / 100,
451  pResBin * binFactor * 2,
452  pResMin / 10,
453  pResMax / 10));
454  oneOverptResolution[9]->setAxisTitle("GeV^{-1}", 1);
455  oneOverptResolution[9]->setAxisTitle("GeV^{-1}", 2);
456  oneOverptResolution.push_back(ibooker.book2D("ResVsPt_GlbSta_oneOverpt",
457  "((1/p_{t})_{GLB} - (1/p_{t})_{STAfromGLB} vs (1/p_{t})_{GLB}",
458  ptBin / 5,
459  ptMin,
460  ptMax / 100,
461  pResBin * binFactor,
462  pResMin,
463  pResMax));
464  oneOverptResolution[10]->setAxisTitle("GeV^{-1}", 1);
465  oneOverptResolution[10]->setAxisTitle("GeV^{-1}", 2);
466  oneOverptResolution.push_back(
467  ibooker.book2D("ResVsPt_TkSta_oneOverpt",
468  "((1/p_{t})_{TKfromGLB} - (1/p_{t})_{STAfromGLB}) vs (1/p_{t})_{TKfromGLB}",
469  ptBin / 5,
470  ptMin,
471  ptMax / 100,
472  pResBin * binFactor,
473  pResMin,
474  pResMax));
475  oneOverptResolution[11]->setAxisTitle("GeV^{-1}", 1);
476  oneOverptResolution[11]->setAxisTitle("GeV^{-1}", 2);
477  oneOverptPull =
478  ibooker.book1D("Pull_TkSta_oneOverpt", "(1/pt)_{TKfromGLB} - (1/pt)_{STAfromGLB} / error", 100, -10, 10);
479 
481  // monitoring of the phi-eta
482  phiVsetaGlbTrack.push_back(ibooker.book2D(
483  histname + "Glb_phiVSeta", "#phi vs #eta (GLB)", etaBin / 2, etaMin, etaMax, phiBin / 2, phiMin, phiMax));
484  phiVsetaGlbTrack.push_back(ibooker.book2D(
485  histname + "Tk_phiVSeta", "#phi vs #eta (TKfromGLB)", etaBin / 2, etaMin, etaMax, phiBin / 2, phiMin, phiMax));
486  phiVsetaGlbTrack.push_back(ibooker.book2D(
487  histname + "Sta_phiVseta", "#phi vs #eta (STAfromGLB)", etaBin / 2, etaMin, etaMax, phiBin / 2, phiMin, phiMax));
488 
489  phiVsetaGlbTrack_badlumi.push_back(ibooker.book2D(
490  histname + "Glb_phiVSeta_badlumi", "#phi vs #eta (GLB)", etaBin / 2, etaMin, etaMax, phiBin / 2, phiMin, phiMax));
491  phiVsetaGlbTrack_badlumi.push_back(ibooker.book2D(histname + "Tk_phiVSeta_badlumi",
492  "#phi vs #eta (TKfromGLB)",
493  etaBin / 2,
494  etaMin,
495  etaMax,
496  phiBin / 2,
497  phiMin,
498  phiMax));
499  phiVsetaGlbTrack_badlumi.push_back(ibooker.book2D(histname + "Sta_phiVseta_badlumi",
500  "#phi vs #eta (STAfromGLB)",
501  etaBin / 2,
502  etaMin,
503  etaMax,
504  phiBin / 2,
505  phiMin,
506  phiMax));
507 
509  // monitoring of the recHits provenance
510  rhAnalysis.push_back(ibooker.book1D("StaRh_Frac_inGlb", "recHits_{STAinGLB} / recHits_{GLB}", rhBin, rhMin, rhMax));
511  rhAnalysis.push_back(ibooker.book1D("TkRh_Frac_inGlb", "recHits_{TKinGLB} / recHits_{GLB}", rhBin, rhMin, rhMax));
512  rhAnalysis.push_back(
513  ibooker.book1D("StaRh_inGlb_Div_RhAssoSta", "recHits_{STAinGLB} / recHits_{STAfromGLB}", rhBin, rhMin, rhMax));
514  rhAnalysis.push_back(
515  ibooker.book1D("TkRh_inGlb_Div_RhAssoTk", "recHits_{TKinGLB} / recHits_{TKfromGLB}", rhBin, rhMin, rhMax));
516  rhAnalysis.push_back(ibooker.book1D(
517  "GlbRh_Div_RhAssoStaTk", "recHits_{GLB} / (recHits_{TKfromGLB}+recHits_{STAfromGLB})", rhBin, rhMin, rhMax));
518  rhAnalysis.push_back(ibooker.book1D("invalidRh_Frac_inTk", "Invalid recHits / rechits_{GLB}", rhBin, rhMin, rhMax));
519 
521  // monitoring of the muon system rotation w.r.t. tracker
522  muVStkSytemRotation.push_back(ibooker.book2D(
523  "muVStkSytemRotation_posMu", "pT_{TK} / pT_{GLB} vs pT_{GLB} for #mu^{+}", 50, 0, 200, 100, 0.8, 1.2));
524  muVStkSytemRotation.push_back(ibooker.book2D(
525  "muVStkSytemRotation_negMu", "pT_{TK} / pT_{GLB} vs pT_{GLB} for #mu^{-}", 50, 0, 200, 100, 0.8, 1.2));
526 }
527 
528 void MuonRecoAnalyzer::GetRes(reco::TrackRef t1, reco::TrackRef t2, string par, float& res, float& pull) {
529  float p1 = 0, p2 = 0, p1e = 1, p2e = 1;
530 
531  if (par == "eta") {
532  p1 = t1->eta();
533  p1e = t1->etaError();
534  p2 = t2->eta();
535  p2e = t2->etaError();
536  } else if (par == "theta") {
537  p1 = t1->theta();
538  p1e = t1->thetaError();
539  p2 = t2->theta();
540  p2e = t2->thetaError();
541  } else if (par == "phi") {
542  p1 = t1->phi();
543  p1e = t1->phiError();
544  p2 = t2->phi();
545  p2e = t2->phiError();
546  } else if (par == "qOverp") {
547  p1 = t1->charge() / t1->p();
548  p1e = t1->qoverpError();
549  p2 = t2->charge() / t2->p();
550  p2e = t2->qoverpError();
551  } else if (par == "oneOverp") {
552  p1 = 1. / t1->p();
553  p1e = t1->qoverpError();
554  p2 = 1. / t2->p();
555  p2e = t2->qoverpError();
556  } else if (par == "qOverpt") {
557  p1 = t1->charge() / t1->pt();
558  p1e = t1->ptError() * p1 * p1;
559  p2 = t2->charge() / t2->pt();
560  p2e = t2->ptError() * p2 * p2;
561  } else if (par == "oneOverpt") {
562  p1 = 1. / t1->pt();
563  p1e = t1->ptError() * p1 * p1;
564  p2 = 1. / t2->pt();
565  p2e = t2->ptError() * p2 * p2;
566  }
567 
568  res = p1 - p2;
569  if (p1e != 0 || p2e != 0)
570  pull = res / sqrt(p1e * p1e + p2e * p2e);
571  else
572  pull = -99;
573  return;
574 }
575 
577  LogTrace(metname) << "[MuonRecoAnalyzer] Analyze the mu";
578  theService->update(iSetup);
579 
580  // Take the muon container
582  iEvent.getByToken(theMuonCollectionLabel_, muons);
583 
586  if (doMVA) {
587  iEvent.getByToken(theBeamSpotLabel_, beamSpot);
588  if (!beamSpot.isValid()) {
589  edm::LogInfo("MuonRecoAnalyzer") << "Error: Can't get the beamspot" << endl;
590  doMVA = false;
591  }
592  iEvent.getByToken(theVertexLabel_, vertex);
593  if (!vertex.isValid()) {
594  edm::LogInfo("MuonRecoAnalyzer") << "Error: Can't get the vertex collection" << endl;
595  doMVA = false;
596  }
597  }
598 
599  //In this part we determine if we want to fill the plots for events where the DCS flag was set to bad
601  bool fillBadLumi = false;
602  if (iEvent.getByToken(dcsStatusCollection_, dcsStatus) && dcsStatus.isValid()) {
603  for (auto const& dcsStatusItr : *dcsStatus) {
604  if (!dcsStatusItr.ready(DcsStatus::CSCp))
605  fillBadLumi = true;
606  if (!dcsStatusItr.ready(DcsStatus::CSCm))
607  fillBadLumi = true;
608  if (!dcsStatusItr.ready(DcsStatus::DT0))
609  fillBadLumi = true;
610  if (!dcsStatusItr.ready(DcsStatus::DTp))
611  fillBadLumi = true;
612  if (!dcsStatusItr.ready(DcsStatus::DTm))
613  fillBadLumi = true;
614  if (!dcsStatusItr.ready(DcsStatus::EBp))
615  fillBadLumi = true;
616  if (!dcsStatusItr.ready(DcsStatus::EBm))
617  fillBadLumi = true;
618  if (!dcsStatusItr.ready(DcsStatus::EEp))
619  fillBadLumi = true;
620  if (!dcsStatusItr.ready(DcsStatus::EEm))
621  fillBadLumi = true;
622  if (!dcsStatusItr.ready(DcsStatus::ESp))
623  fillBadLumi = true;
624  if (!dcsStatusItr.ready(DcsStatus::ESm))
625  fillBadLumi = true;
626  if (!dcsStatusItr.ready(DcsStatus::HBHEa))
627  fillBadLumi = true;
628  if (!dcsStatusItr.ready(DcsStatus::HBHEb))
629  fillBadLumi = true;
630  if (!dcsStatusItr.ready(DcsStatus::HBHEc))
631  fillBadLumi = true;
632  if (!dcsStatusItr.ready(DcsStatus::HF))
633  fillBadLumi = true;
634  if (!dcsStatusItr.ready(DcsStatus::HO))
635  fillBadLumi = true;
636  if (!dcsStatusItr.ready(DcsStatus::BPIX))
637  fillBadLumi = true;
638  if (!dcsStatusItr.ready(DcsStatus::FPIX))
639  fillBadLumi = true;
640  if (!dcsStatusItr.ready(DcsStatus::RPC))
641  fillBadLumi = true;
642  if (!dcsStatusItr.ready(DcsStatus::TIBTID))
643  fillBadLumi = true;
644  if (!dcsStatusItr.ready(DcsStatus::TOB))
645  fillBadLumi = true;
646  if (!dcsStatusItr.ready(DcsStatus::TECp))
647  fillBadLumi = true;
648  if (!dcsStatusItr.ready(DcsStatus::TECm))
649  fillBadLumi = true;
650  //if (!dcsStatusItr.ready(DcsStatus::CASTOR)) fillBadLumi = true;
651  }
652  }
653 
654  float res = 0, pull = 0;
655  if (!muons.isValid())
656  return;
657 
658  for (edm::View<reco::Muon>::const_iterator muon = muons->begin(); muon != muons->end(); ++muon) {
659  //Needed for MVA soft muon
660 
661  reco::TrackRef gTrack = muon->globalTrack();
662  reco::TrackRef iTrack = muon->innerTrack();
663  reco::TrackRef oTrack = muon->outerTrack();
664  if (iTrack.isNonnull() && oTrack.isNonnull() && gTrack.isNonnull()) {
665  const reco::HitPattern gHits = gTrack->hitPattern();
666  const reco::HitPattern iHits = iTrack->hitPattern();
667  const reco::MuonQuality muonQuality = muon->combinedQuality();
668  int pvIndex = 0;
669  math::XYZPoint refPoint;
670  if (doMVA) {
671  pvIndex = getPv(iTrack.index(), &(*vertex)); //HFDumpUtitilies
672  if (pvIndex > -1) {
673  refPoint = vertex->at(pvIndex).position();
674  } else {
675  if (beamSpot.isValid()) {
676  refPoint = beamSpot->position();
677  } else {
678  edm::LogInfo("MuonRecoAnalyzer") << "ERROR: No beam sport found!" << endl;
679  }
680  }
681  }
682  ptSoftMuonMVA->Fill(iTrack->eta());
683  deltaRSoftMuonMVA->Fill(getDeltaR(*iTrack, *oTrack));
684  gNchi2SoftMuonMVA->Fill(gTrack->normalizedChi2());
685  vMuHitsSoftMuonMVA->Fill(gHits.numberOfValidMuonHits());
686  mNuStationsSoftMuonMVA->Fill(muon->numberOfMatchedStations());
687  if (doMVA) {
688  dxyRefSoftMuonMVA->Fill(iTrack->dxy(refPoint));
689  dzRefSoftMuonMVA->Fill(iTrack->dz(refPoint));
690  }
691  LWHSoftMuonMVA->Fill(iHits.trackerLayersWithMeasurement());
692  valPixHitsSoftMuonMVA->Fill(iHits.numberOfValidPixelHits());
693  innerChi2SoftMuonMVA->Fill(iTrack->normalizedChi2());
694  outerChi2SoftMuonMVA->Fill(oTrack->normalizedChi2());
695  iValFracSoftMuonMVA->Fill(iTrack->validFraction());
696  //segCompSoftMuonMVA->Fill(reco::Muon::segmentCompatibility(*muon));
697  chi2LocMomSoftMuonMVA->Fill(muonQuality.chi2LocalMomentum);
698  chi2LocPosSoftMuonMVA->Fill(muonQuality.chi2LocalPosition);
699  glbTrackTailProbSoftMuonMVA->Fill(muonQuality.glbTrackProbability);
700  NTrkVHitsSoftMuonMVA->Fill(iHits.numberOfValidTrackerHits());
701  kinkFinderSoftMuonMVA->Fill(muonQuality.trkKink);
702  vRPChitsSoftMuonMVA->Fill(gHits.numberOfValidMuonRPCHits());
703  glbKinkFinderSoftMuonMVA->Fill(muonQuality.glbKink);
704  glbKinkFinderLogSoftMuonMVA->Fill(TMath::Log(2 + muonQuality.glbKink));
705  staRelChi2SoftMuonMVA->Fill(muonQuality.staRelChi2);
706  glbDeltaEtaPhiSoftMuonMVA->Fill(muonQuality.globalDeltaEtaPhi);
707  trkRelChi2SoftMuonMVA->Fill(muonQuality.trkRelChi2);
708  vDThitsSoftMuonMVA->Fill(gHits.numberOfValidMuonDTHits());
709  vCSChitsSoftMuonMVA->Fill(gHits.numberOfValidMuonCSCHits());
710  timeAtIpInOutSoftMuonMVA->Fill(muon->time().timeAtIpInOut);
711  timeAtIpInOutErrSoftMuonMVA->Fill(muon->time().timeAtIpInOutErr);
712  //getMuonHitsPerStationSoftMuonMVA->Fill(gTrack);
713  QprodSoftMuonMVA->Fill((iTrack->charge() * oTrack->charge()));
714  }
715 
716  if (muon->isGlobalMuon()) {
717  LogTrace(metname) << "[MuonRecoAnalyzer] The mu is global - filling the histos";
718  if (muon->isTrackerMuon() && muon->isStandAloneMuon())
719  muReco->Fill(1);
720  if (!(muon->isTrackerMuon()) && muon->isStandAloneMuon())
721  muReco->Fill(2);
722  if (!muon->isStandAloneMuon())
723  LogTrace(metname) << "[MuonRecoAnalyzer] ERROR: the mu is global but not standalone!";
724  // get the track combinig the information from both the Tracker and the Spectrometer
725  reco::TrackRef recoCombinedGlbTrack = muon->combinedMuon();
726 
727  // get the track using only the tracker data
728  reco::TrackRef recoTkGlbTrack = muon->track();
729  // get the track using only the mu spectrometer data
730  reco::TrackRef recoStaGlbTrack = muon->standAloneMuon();
731  etaGlbTrack[0]->Fill(recoCombinedGlbTrack->eta());
732  etaGlbTrack[1]->Fill(recoTkGlbTrack->eta());
733  etaGlbTrack[2]->Fill(recoStaGlbTrack->eta());
734 
735  phiVsetaGlbTrack[0]->Fill(recoCombinedGlbTrack->eta(), recoCombinedGlbTrack->phi());
736  phiVsetaGlbTrack[1]->Fill(recoTkGlbTrack->eta(), recoTkGlbTrack->phi());
737  phiVsetaGlbTrack[2]->Fill(recoStaGlbTrack->eta(), recoStaGlbTrack->phi());
738 
739  if (fillBadLumi) {
740  phiVsetaGlbTrack_badlumi[0]->Fill(recoCombinedGlbTrack->eta(), recoCombinedGlbTrack->phi());
741  phiVsetaGlbTrack_badlumi[1]->Fill(recoTkGlbTrack->eta(), recoTkGlbTrack->phi());
742  phiVsetaGlbTrack_badlumi[2]->Fill(recoStaGlbTrack->eta(), recoStaGlbTrack->phi());
743  }
744 
745  GetRes(recoTkGlbTrack, recoCombinedGlbTrack, "eta", res, pull);
746  etaResolution[0]->Fill(res);
747  GetRes(recoCombinedGlbTrack, recoStaGlbTrack, "eta", res, pull);
748  etaResolution[1]->Fill(res);
749  GetRes(recoTkGlbTrack, recoStaGlbTrack, "eta", res, pull);
750  etaResolution[2]->Fill(res);
751  etaPull->Fill(pull);
752  etaResolution[3]->Fill(recoCombinedGlbTrack->eta(), recoTkGlbTrack->eta() - recoCombinedGlbTrack->eta());
753  etaResolution[4]->Fill(recoCombinedGlbTrack->eta(), -recoStaGlbTrack->eta() + recoCombinedGlbTrack->eta());
754  etaResolution[5]->Fill(recoCombinedGlbTrack->eta(), recoTkGlbTrack->eta() - recoStaGlbTrack->eta());
755 
756  thetaGlbTrack[0]->Fill(recoCombinedGlbTrack->theta());
757  thetaGlbTrack[1]->Fill(recoTkGlbTrack->theta());
758  thetaGlbTrack[2]->Fill(recoStaGlbTrack->theta());
759  GetRes(recoTkGlbTrack, recoCombinedGlbTrack, "theta", res, pull);
760  thetaResolution[0]->Fill(res);
761  GetRes(recoCombinedGlbTrack, recoStaGlbTrack, "theta", res, pull);
762  thetaResolution[1]->Fill(res);
763 
764  GetRes(recoTkGlbTrack, recoStaGlbTrack, "theta", res, pull);
765  thetaResolution[2]->Fill(res);
766  thetaPull->Fill(pull);
767  thetaResolution[3]->Fill(recoCombinedGlbTrack->theta(), recoTkGlbTrack->theta() - recoCombinedGlbTrack->theta());
768  thetaResolution[4]->Fill(recoCombinedGlbTrack->theta(),
769  -recoStaGlbTrack->theta() + recoCombinedGlbTrack->theta());
770  thetaResolution[5]->Fill(recoCombinedGlbTrack->theta(), recoTkGlbTrack->theta() - recoStaGlbTrack->theta());
771 
772  phiGlbTrack[0]->Fill(recoCombinedGlbTrack->phi());
773  phiGlbTrack[1]->Fill(recoTkGlbTrack->phi());
774  phiGlbTrack[2]->Fill(recoStaGlbTrack->phi());
775  GetRes(recoTkGlbTrack, recoCombinedGlbTrack, "phi", res, pull);
776  phiResolution[0]->Fill(res);
777  GetRes(recoCombinedGlbTrack, recoStaGlbTrack, "phi", res, pull);
778  phiResolution[1]->Fill(res);
779  GetRes(recoTkGlbTrack, recoStaGlbTrack, "phi", res, pull);
780  phiResolution[2]->Fill(res);
781  phiPull->Fill(pull);
782  phiResolution[3]->Fill(recoCombinedGlbTrack->phi(), recoTkGlbTrack->phi() - recoCombinedGlbTrack->phi());
783  phiResolution[4]->Fill(recoCombinedGlbTrack->phi(), -recoStaGlbTrack->phi() + recoCombinedGlbTrack->phi());
784  phiResolution[5]->Fill(recoCombinedGlbTrack->phi(), recoTkGlbTrack->phi() - recoStaGlbTrack->phi());
785 
786  chi2OvDFGlbTrack[0]->Fill(recoCombinedGlbTrack->normalizedChi2());
787  chi2OvDFGlbTrack[1]->Fill(recoTkGlbTrack->normalizedChi2());
788  chi2OvDFGlbTrack[2]->Fill(recoStaGlbTrack->normalizedChi2());
789  //-------------------------
790  // double probchi = TMath::Prob(recoCombinedGlbTrack->normalizedChi2(),recoCombinedGlbTrack->ndof());
791  // cout << "rellenando histos."<<endl;
792  probchi2GlbTrack[0]->Fill(TMath::Prob(recoCombinedGlbTrack->chi2(), recoCombinedGlbTrack->ndof()));
793  probchi2GlbTrack[1]->Fill(TMath::Prob(recoTkGlbTrack->chi2(), recoTkGlbTrack->ndof()));
794  probchi2GlbTrack[2]->Fill(TMath::Prob(recoStaGlbTrack->chi2(), recoStaGlbTrack->ndof()));
795  // cout << "rellenados histos."<<endl;
796  //-------------------------
797 
798  pGlbTrack[0]->Fill(recoCombinedGlbTrack->p());
799  pGlbTrack[1]->Fill(recoTkGlbTrack->p());
800  pGlbTrack[2]->Fill(recoStaGlbTrack->p());
801 
802  ptGlbTrack[0]->Fill(recoCombinedGlbTrack->pt());
803  ptGlbTrack[1]->Fill(recoTkGlbTrack->pt());
804  ptGlbTrack[2]->Fill(recoStaGlbTrack->pt());
805 
806  qGlbTrack[0]->Fill(recoCombinedGlbTrack->charge());
807  qGlbTrack[1]->Fill(recoTkGlbTrack->charge());
808  qGlbTrack[2]->Fill(recoStaGlbTrack->charge());
809  if (recoCombinedGlbTrack->charge() == recoStaGlbTrack->charge())
810  qGlbTrack[3]->Fill(1);
811  else
812  qGlbTrack[3]->Fill(2);
813  if (recoCombinedGlbTrack->charge() == recoTkGlbTrack->charge())
814  qGlbTrack[3]->Fill(3);
815  else
816  qGlbTrack[3]->Fill(4);
817  if (recoStaGlbTrack->charge() == recoTkGlbTrack->charge())
818  qGlbTrack[3]->Fill(5);
819  else
820  qGlbTrack[3]->Fill(6);
821  if (recoCombinedGlbTrack->charge() != recoStaGlbTrack->charge() &&
822  recoCombinedGlbTrack->charge() != recoTkGlbTrack->charge())
823  qGlbTrack[3]->Fill(7);
824  if (recoCombinedGlbTrack->charge() == recoStaGlbTrack->charge() &&
825  recoCombinedGlbTrack->charge() == recoTkGlbTrack->charge())
826  qGlbTrack[3]->Fill(8);
827 
828  GetRes(recoTkGlbTrack, recoCombinedGlbTrack, "qOverp", res, pull);
829  qOverpResolution[0]->Fill(res);
830  GetRes(recoCombinedGlbTrack, recoStaGlbTrack, "qOverp", res, pull);
831  qOverpResolution[1]->Fill(res);
832  GetRes(recoTkGlbTrack, recoStaGlbTrack, "qOverp", res, pull);
833  qOverpResolution[2]->Fill(res);
834  qOverpPull->Fill(pull);
835 
836  GetRes(recoTkGlbTrack, recoCombinedGlbTrack, "oneOverp", res, pull);
837  oneOverpResolution[0]->Fill(res);
838  GetRes(recoCombinedGlbTrack, recoStaGlbTrack, "oneOverp", res, pull);
839  oneOverpResolution[1]->Fill(res);
840  GetRes(recoTkGlbTrack, recoStaGlbTrack, "oneOverp", res, pull);
841  oneOverpResolution[2]->Fill(res);
842  oneOverpPull->Fill(pull);
843 
844  GetRes(recoTkGlbTrack, recoCombinedGlbTrack, "qOverpt", res, pull);
845  qOverptResolution[0]->Fill(res);
846  GetRes(recoCombinedGlbTrack, recoStaGlbTrack, "qOverpt", res, pull);
847  qOverptResolution[1]->Fill(res);
848  GetRes(recoTkGlbTrack, recoStaGlbTrack, "qOverpt", res, pull);
849  qOverptResolution[2]->Fill(res);
850  qOverptPull->Fill(pull);
851 
852  GetRes(recoTkGlbTrack, recoCombinedGlbTrack, "oneOverpt", res, pull);
853  oneOverptResolution[0]->Fill(res);
854  GetRes(recoCombinedGlbTrack, recoStaGlbTrack, "oneOverpt", res, pull);
855  oneOverptResolution[1]->Fill(res);
856  GetRes(recoTkGlbTrack, recoStaGlbTrack, "oneOverpt", res, pull);
857  oneOverptResolution[2]->Fill(res);
858  oneOverptPull->Fill(pull);
859 
860  // //--- Test new tunePMuonBestTrack() method from Muon.h
861 
862  reco::TrackRef recoBestTrack = muon->muonBestTrack();
863 
864  reco::TrackRef recoTunePBestTrack = muon->tunePMuonBestTrack();
865 
866  double bestTrackPt = recoBestTrack->pt();
867 
868  double tunePBestTrackPt = recoTunePBestTrack->pt();
869 
870  double tunePBestTrackRes = (bestTrackPt - tunePBestTrackPt) / bestTrackPt;
871 
872  tunePResolution->Fill(tunePBestTrackRes);
873 
874  oneOverptResolution[3]->Fill(recoCombinedGlbTrack->eta(),
875  (1 / recoTkGlbTrack->pt()) - (1 / recoCombinedGlbTrack->pt()));
876  oneOverptResolution[4]->Fill(recoCombinedGlbTrack->eta(),
877  -(1 / recoStaGlbTrack->pt()) + (1 / recoCombinedGlbTrack->pt()));
878  oneOverptResolution[5]->Fill(recoCombinedGlbTrack->eta(),
879  (1 / recoTkGlbTrack->pt()) - (1 / recoStaGlbTrack->pt()));
880  oneOverptResolution[6]->Fill(recoCombinedGlbTrack->phi(),
881  (1 / recoTkGlbTrack->pt()) - (1 / recoCombinedGlbTrack->pt()));
882  oneOverptResolution[7]->Fill(recoCombinedGlbTrack->phi(),
883  -(1 / recoStaGlbTrack->pt()) + (1 / recoCombinedGlbTrack->pt()));
884  oneOverptResolution[8]->Fill(recoCombinedGlbTrack->phi(),
885  (1 / recoTkGlbTrack->pt()) - (1 / recoStaGlbTrack->pt()));
886  oneOverptResolution[9]->Fill(recoCombinedGlbTrack->pt(),
887  (1 / recoTkGlbTrack->pt()) - (1 / recoCombinedGlbTrack->pt()));
888  oneOverptResolution[10]->Fill(recoCombinedGlbTrack->pt(),
889  -(1 / recoStaGlbTrack->pt()) + (1 / recoCombinedGlbTrack->pt()));
890  oneOverptResolution[11]->Fill(recoCombinedGlbTrack->pt(),
891  (1 / recoTkGlbTrack->pt()) - (1 / recoStaGlbTrack->pt()));
892 
893  if (!IsminiAOD) {
894  // valid hits Glb track
895  double rhGlb = recoCombinedGlbTrack->found();
896  // valid hits Glb track from Tracker
897  double rhGlb_StaProvenance = 0;
898  // valid hits Glb track from Sta system
899  double rhGlb_TkProvenance = 0;
900 
901  for (trackingRecHit_iterator recHit = recoCombinedGlbTrack->recHitsBegin();
902  recHit != recoCombinedGlbTrack->recHitsEnd();
903  ++recHit) {
904  if ((*recHit)->isValid()) {
905  DetId id = (*recHit)->geographicalId();
906  if (id.det() == DetId::Muon)
907  rhGlb_StaProvenance++;
908  if (id.det() == DetId::Tracker)
909  rhGlb_TkProvenance++;
910  }
911  }
912  // valid hits Sta track associated to Glb track
913  double rhStaGlb = recoStaGlbTrack->recHitsSize();
914  // valid hits Traker track associated to Glb track
915  double rhTkGlb = recoTkGlbTrack->found();
916  // invalid hits Traker track associated to Glb track
917  double rhTkGlb_notValid = recoTkGlbTrack->lost();
918 
919  // fill the histos
920  rhAnalysis[0]->Fill(rhGlb_StaProvenance / rhGlb);
921  rhAnalysis[1]->Fill(rhGlb_TkProvenance / rhGlb);
922  rhAnalysis[2]->Fill(rhGlb_StaProvenance / rhStaGlb);
923  rhAnalysis[3]->Fill(rhGlb_TkProvenance / rhTkGlb);
924  rhAnalysis[4]->Fill(rhGlb / (rhStaGlb + rhTkGlb));
925  rhAnalysis[5]->Fill(rhTkGlb_notValid / rhGlb);
926  }
927  // aligment plots (mu system w.r.t. tracker rotation)
928  if (recoCombinedGlbTrack->charge() > 0)
929  muVStkSytemRotation[0]->Fill(recoCombinedGlbTrack->pt(), recoTkGlbTrack->pt() / recoCombinedGlbTrack->pt());
930  else
931  muVStkSytemRotation[1]->Fill(recoCombinedGlbTrack->pt(), recoTkGlbTrack->pt() / recoCombinedGlbTrack->pt());
932  }
933 
934  if (muon->isTrackerMuon() && !(muon->isGlobalMuon())) {
935  LogTrace(metname) << "[MuonRecoAnalyzer] The mu is tracker only - filling the histos";
936  if (muon->isStandAloneMuon())
937  muReco->Fill(3);
938  if (!(muon->isStandAloneMuon()))
939  muReco->Fill(4);
940 
941  // get the track using only the tracker data
942  reco::TrackRef recoTrack = muon->track();
943 
944  etaTrack->Fill(recoTrack->eta());
945  thetaTrack->Fill(recoTrack->theta());
946  phiTrack->Fill(recoTrack->phi());
947  chi2OvDFTrack->Fill(recoTrack->normalizedChi2());
948  probchi2Track->Fill(TMath::Prob(recoTrack->chi2(), recoTrack->ndof()));
949  pTrack->Fill(recoTrack->p());
950  ptTrack->Fill(recoTrack->pt());
951  qTrack->Fill(recoTrack->charge());
952  }
953 
954  if (muon->isStandAloneMuon() && !(muon->isGlobalMuon())) {
955  LogTrace(metname) << "[MuonRecoAnalyzer] The mu is STA only - filling the histos";
956  if (!(muon->isTrackerMuon()))
957  muReco->Fill(5);
958 
959  // get the track using only the mu spectrometer data
960  reco::TrackRef recoStaTrack = muon->standAloneMuon();
961 
962  etaStaTrack->Fill(recoStaTrack->eta());
963  thetaStaTrack->Fill(recoStaTrack->theta());
964  phiStaTrack->Fill(recoStaTrack->phi());
965  chi2OvDFStaTrack->Fill(recoStaTrack->normalizedChi2());
966  probchi2StaTrack->Fill(TMath::Prob(recoStaTrack->chi2(), recoStaTrack->ndof()));
967  pStaTrack->Fill(recoStaTrack->p());
968  ptStaTrack->Fill(recoStaTrack->pt());
969  qStaTrack->Fill(recoStaTrack->charge());
970  }
971 
972  if (muon->isCaloMuon() && !(muon->isGlobalMuon()) && !(muon->isTrackerMuon()) && !(muon->isStandAloneMuon()))
973  muReco->Fill(6);
974 
975  //efficiency plots
976 
977  // get the track using only the mu spectrometer data
978  reco::TrackRef recoStaGlbTrack = muon->standAloneMuon();
979 
980  if (muon->isStandAloneMuon()) {
981  etaEfficiency[0]->Fill(recoStaGlbTrack->eta());
982  phiEfficiency[0]->Fill(recoStaGlbTrack->phi());
983  }
984  if (muon->isStandAloneMuon() && muon->isGlobalMuon()) {
985  etaEfficiency[1]->Fill(recoStaGlbTrack->eta());
986  phiEfficiency[1]->Fill(recoStaGlbTrack->phi());
987  }
988  }
989 }
990 
991 //Needed by MVA Soft Muon
993  double dphi = acos(cos(track1.phi() - track2.phi()));
994  double deta = track1.eta() - track2.eta();
995  return sqrt(dphi * dphi + deta * deta);
996 }
997 
998 // ----------------------------------------------------------------------
1000  if (vc) {
1001  for (unsigned int i = 0; i < vc->size(); ++i) {
1002  reco::Vertex::trackRef_iterator v1TrackIter;
1003  reco::Vertex::trackRef_iterator v1TrackBegin = vc->at(i).tracks_begin();
1004  reco::Vertex::trackRef_iterator v1TrackEnd = vc->at(i).tracks_end();
1005  for (v1TrackIter = v1TrackBegin; v1TrackIter != v1TrackEnd; v1TrackIter++) {
1006  if (static_cast<unsigned int>(tidx) == v1TrackIter->key())
1007  return i;
1008  }
1009  }
1010  }
1011  return -1;
1012 }
float chi2LocalPosition
chi2 value for the STA-TK matching of local position
Definition: MuonQuality.h:19
MonitorElement * book1D(TString const &name, TString const &title, int const nchX, double const lowX, double const highX)
Definition: DQMStore.cc:239
T getParameter(std::string const &) const
T getUntrackedParameter(std::string const &, T const &) const
bool isNonnull() const
Checks for non-null.
Definition: Ref.h:238
float chi2LocalMomentum
chi2 value for the STA-TK matching of local momentum
Definition: MuonQuality.h:21
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:525
const std::string metname
float glbKink
value of the kink algorithm applied to the global track
Definition: MuonQuality.h:13
void setCurrentFolder(std::string const &fullpath)
Definition: DQMStore.cc:418
float glbTrackProbability
the tail probability (-ln(P)) of the global fit
Definition: MuonQuality.h:29
key_type index() const
Definition: Ref.h:253
double phi() const
azimuthal angle of momentum vector
Definition: TrackBase.h:614
std::vector< Vertex > VertexCollection
collection of Vertex objects
Definition: VertexFwd.h:9
float trkKink
value of the kink algorithm applied to the inner track stub
Definition: MuonQuality.h:11
float trkRelChi2
chi2 value for the inner track stub with respect to the global track
Definition: MuonQuality.h:15
int trackerLayersWithMeasurement() const
Definition: HitPattern.cc:513
Definition: Electron.h:6
void analyze(const edm::Event &, const edm::EventSetup &) override
Inizialize parameters for histo binning.
int numberOfValidMuonCSCHits() const
Definition: HitPattern.h:837
int numberOfValidMuonRPCHits() const
Definition: HitPattern.h:841
int iEvent
Definition: GenABIO.cc:224
double getDeltaR(reco::Track track1, reco::Track track2)
double eta() const
pseudorapidity of momentum vector
Definition: TrackBase.h:617
void Fill(HcalDetId &id, double val, std::vector< TH2F > &depth)
T sqrt(T t)
Definition: SSEVec.h:19
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
float globalDeltaEtaPhi
global delta-Eta-Phi of STA-TK matching
Definition: MuonQuality.h:25
bool isValid() const
Definition: HandleBase.h:70
virtual void setBinLabel(int bin, const std::string &label, int axis=1)
set bin label for x, y or z axis (axis=1, 2, 3 respectively)
double p2[4]
Definition: TauolaWrapper.h:90
#define LogTrace(id)
void GetRes(reco::TrackRef t1, reco::TrackRef t2, std::string par, float &res, float &pull)
void bookHistograms(DQMStore::IBooker &, edm::Run const &, edm::EventSetup const &) override
float staRelChi2
chi2 value for the outer track stub with respect to the global track
Definition: MuonQuality.h:17
Definition: DetId.h:17
int getPv(int tidx, const reco::VertexCollection *vc)
XYZPointD XYZPoint
point in space with cartesian internal representation
Definition: Point3D.h:12
int numberOfValidTrackerHits() const
Definition: HitPattern.h:789
HLT enums.
double p1[4]
Definition: TauolaWrapper.h:89
boost::indirect_iterator< typename seq_t::const_iterator > const_iterator
Definition: View.h:86
MonitorElement * book2D(TString const &name, TString const &title, int nchX, double lowX, double highX, int nchY, double lowY, double highY)
Definition: DQMStore.cc:266
int numberOfValidMuonDTHits() const
Definition: HitPattern.h:833
int numberOfValidPixelHits() const
Definition: HitPattern.h:801
const Point & position() const
position
Definition: BeamSpot.h:59
int numberOfValidMuonHits() const
Definition: HitPattern.h:793
~MuonRecoAnalyzer() override
Destructor.
std::vector< TrackBaseRef >::const_iterator trackRef_iterator
The iteratator for the vector<TrackRef>
Definition: Vertex.h:37
Definition: Run.h:45
MuonRecoAnalyzer(const edm::ParameterSet &)
Constructor.
virtual void setAxisTitle(const std::string &title, int axis=1)
set x-, y- or z-axis title (axis=1, 2, 3 respectively)