ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/FirstData/Macros/ProjectQQ.C
Revision: 1.1
Committed: Thu Dec 17 18:09:59 2009 UTC (15 years, 4 months ago) by hwoehri
Content type: text/plain
Branch: MAIN
CVS Tags: HEAD
Log Message:
*** empty log message ***

File Contents

# Content
1 #define ProjectQQ_cxx
2 #include "ProjectQQ.h"
3 #include <TH2.h>
4 #include <TStyle.h>
5 #include <TCanvas.h>
6 #include <TLorentzVector.h>
7 #include "TClonesArray.h"
8
9 #include "commonVar.inc"
10
11 //dimuon related histograms
12 TH1F *hReco_mass[kNbCharge][kNbCath][kNbCuts];
13 TH1F *hReco_pT[kNbDimuonSet][kNbCharge][kNbCath][kNbCuts];
14 TH1F *hReco_rap[kNbDimuonSet][kNbCharge][kNbCath][kNbCuts];
15 TH1F *hRecoBest_mass[kNbCharge][kNbCath];
16 TH1F *hRecoBest_pT[kNbDimuonSet][kNbCharge][kNbCath];
17 TH1F *hRecoBest_rap[kNbDimuonSet][kNbCharge][kNbCath];
18 TH1F *hMult_QQ[kNbCharge][kNbCath][kNbCuts];
19 TH2F *hReco_eta1_eta2[kNbDimuonSet][kNbCharge][kNbCath+1][kNbCuts];
20 TH1F *hRecoNoBest;
21
22 //single muon histos for two dimuon mass windows:
23 TH1F *hMuon_eta[kNbDimuonSet][kNbCharge][kNbCath];
24 TH1F *hMuon_pT[kNbDimuonSet][kNbCharge][kNbCath];
25 TH2F *hMuon_pT_eta[kNbDimuonSet][kNbCharge][kNbMuCath];
26 // quality histogams
27 TH1F *hChi2GlobalFit[kNbDimuonSet][kNbCharge];
28 TH1F *hChi2TrackerFit[kNbDimuonSet][kNbCharge][kNbMuCath];
29 TH1F *hTM2DCompatibilityTight[kNbDimuonSet][kNbCharge];
30 TH1F *hTMLastStationOptimizedLowPtLoose[kNbDimuonSet][kNbCharge];
31 TH1F *hTrackerMuonArbitrated[kNbDimuonSet][kNbCharge];//tracker only
32 TH1F *hCaloCompatibility[kNbDimuonSet][kNbCharge][kNbMuCath];
33 TH1F *hNHitsSilicon[kNbDimuonSet][kNbCharge][kNbMuCath];
34 TH1F *hDimuVtxProb[kNbDimuonSet][kNbCharge][kNbCath];
35
36 //=======================================
37 void ProjectQQ::Loop(Bool_t removeQQ, Bool_t matchMC)
38 {
39 if (fChain == 0) return;
40
41 Long64_t nentries = fChain->GetEntries();
42 // nentries = 10000;
43
44 Long64_t nb = 0;
45
46 Int_t countAllRecoQQEvents = 0, countNotFoundBestQQEvents = 0;
47 for (Long64_t jentry=0; jentry<nentries;jentry++) {
48
49 if(jentry % 1000 == 0) printf("event %d/%d\n", (Int_t) jentry, (Int_t) nentries);
50
51 Long64_t ientry = LoadTree(jentry);
52 if (ientry < 0) break;
53 nb = fChain->GetEntry(jentry);
54
55 if(removeQQ && Mc_QQ_size > 0) continue; //Mc_QQ_size contains the nb. of MC J/psi's
56
57 Int_t chargeID;
58 Int_t nComb[kNbCharge][kNbCath][kNbCuts];
59 for(int iCharge = 0; iCharge < kNbCharge; iCharge++)
60 for(int iCat = 0; iCat < kNbCath; iCat++)
61 for(int iCut = 0; iCut < kNbCuts; iCut++)
62 nComb[iCharge][iCat][iCut] = 0;
63
64 //loop over all quarkonia
65 for(int iQQ = 0; iQQ < Reco_QQ_size; iQQ++){
66
67 chargeID = 0;
68 if(Reco_QQ_sign[iQQ] < 0) chargeID = 1;
69 else if(Reco_QQ_sign[iQQ] > 0) chargeID = 2;
70
71 TLorentzVector *Reco_QQ = (TLorentzVector *) Reco_QQ_4mom->At(iQQ);
72 Double_t mass_Reco_QQ = Reco_QQ->M();
73 Double_t pT_Reco_QQ = Reco_QQ->Pt();
74 Double_t rap_Reco_QQ = Reco_QQ->Rapidity();
75
76 //get the individual muons:
77 TLorentzVector *Reco_QQ_posMu, *Reco_QQ_negMu;
78 Bool_t etaCorr = kFALSE;
79 if(Reco_QQ_type[iQQ] == 0){ //global-global
80 if(Reco_QQ_mupl[iQQ] < Reco_mu_glb_size && Reco_QQ_mumi[iQQ] < Reco_mu_glb_size){
81 Reco_QQ_posMu = (TLorentzVector *) Reco_mu_glb_4mom->At(Reco_QQ_mupl[iQQ]);
82 Reco_QQ_negMu = (TLorentzVector *) Reco_mu_glb_4mom->At(Reco_QQ_mumi[iQQ]);
83 etaCorr = kTRUE;
84 }
85 }
86 else if(Reco_QQ_type[iQQ] == 2){ //tracker-tracker
87 if(Reco_QQ_mupl[iQQ] < Reco_mu_trk_size && Reco_QQ_mumi[iQQ] < Reco_mu_trk_size){
88 Reco_QQ_posMu = (TLorentzVector *) Reco_mu_trk_4mom->At(Reco_QQ_mupl[iQQ]);
89 Reco_QQ_negMu = (TLorentzVector *) Reco_mu_trk_4mom->At(Reco_QQ_mumi[iQQ]);
90 etaCorr = kTRUE;
91 }
92 }
93
94 Bool_t isMatched[2];//only relevant for MC
95 for (int iCut=0; iCut<kNbCuts; iCut++)
96 {
97 //start with the MC matching:
98 // if (iCut == 0 && matchMC){
99 if (matchMC){
100
101 if(Mc_QQmupl_indx[0] > 1 || Mc_QQmumi_indx[0] > 1) continue; //sanity check
102
103 isMatched[0] = kFALSE;
104 isMatched[1] = kFALSE;
105
106 switch (Reco_QQ_type[iQQ]) {
107 case 0: // global-global
108 //first MC muon
109 if(Reco_QQ_muhpt[iQQ] < Reco_mu_glb_size){//sanity check
110 if(deltaR((TLorentzVector*)Mc_mu_4mom->At(Mc_QQmupl_indx[0]),
111 (TLorentzVector*)Reco_mu_glb_4mom->At(Reco_QQ_muhpt[iQQ])) < MAX_deltaR ||
112 deltaR((TLorentzVector*)Mc_mu_4mom->At(Mc_QQmumi_indx[0]),
113 (TLorentzVector*)Reco_mu_glb_4mom->At(Reco_QQ_muhpt[iQQ])) < MAX_deltaR)
114 isMatched[0] = kTRUE;
115 }
116 //second MC muon
117 if(Reco_QQ_mulpt[iQQ] < Reco_mu_glb_size){//sanity check
118 if(deltaR((TLorentzVector*)Mc_mu_4mom->At(Mc_QQmupl_indx[0]),
119 (TLorentzVector*)Reco_mu_glb_4mom->At(Reco_QQ_mulpt[iQQ])) < MAX_deltaR ||
120 deltaR((TLorentzVector*)Mc_mu_4mom->At(Mc_QQmumi_indx[0]),
121 (TLorentzVector*)Reco_mu_glb_4mom->At(Reco_QQ_mulpt[iQQ])) < MAX_deltaR)
122 isMatched[1] = kTRUE;
123 }
124 break;
125 case 1: // global-tracker
126 //first MC muon
127 if(Reco_QQ_muhpt[iQQ] < Reco_mu_glb_size){//sanity check
128 if(deltaR((TLorentzVector*)Mc_mu_4mom->At(Mc_QQmupl_indx[0]),
129 (TLorentzVector*)Reco_mu_glb_4mom->At(Reco_QQ_muhpt[iQQ])) < MAX_deltaR ||
130 deltaR((TLorentzVector*)Mc_mu_4mom->At(Mc_QQmumi_indx[0]),
131 (TLorentzVector*)Reco_mu_glb_4mom->At(Reco_QQ_muhpt[iQQ])) < MAX_deltaR)
132 isMatched[0] = kTRUE;
133 }
134 //second MC muon
135 if(Reco_QQ_mulpt[iQQ] < Reco_mu_trk_size){//sanity check
136 if(deltaR((TLorentzVector*)Mc_mu_4mom->At(Mc_QQmupl_indx[0]),
137 (TLorentzVector*)Reco_mu_trk_4mom->At(Reco_QQ_mulpt[iQQ])) < MAX_deltaR ||
138 deltaR((TLorentzVector*)Mc_mu_4mom->At(Mc_QQmumi_indx[0]),
139 (TLorentzVector*)Reco_mu_trk_4mom->At(Reco_QQ_mulpt[iQQ])) < MAX_deltaR)
140 isMatched[1] = kTRUE;
141 }
142 break;
143 case 2: // tracker-tracker
144 //first MC muon
145 if(Reco_QQ_muhpt[iQQ] < Reco_mu_trk_size){//sanity check
146 if(deltaR((TLorentzVector*)Mc_mu_4mom->At(Mc_QQmupl_indx[0]),
147 (TLorentzVector*)Reco_mu_trk_4mom->At(Reco_QQ_muhpt[iQQ])) < MAX_deltaR ||
148 deltaR((TLorentzVector*)Mc_mu_4mom->At(Mc_QQmumi_indx[0]),
149 (TLorentzVector*)Reco_mu_trk_4mom->At(Reco_QQ_muhpt[iQQ])) < MAX_deltaR)
150 isMatched[0] = kTRUE;
151 }
152 //second MC muon
153 if(Reco_QQ_mulpt[iQQ] < Reco_mu_trk_size){//sanity check
154 if(deltaR((TLorentzVector*)Mc_mu_4mom->At(Mc_QQmupl_indx[0]),
155 (TLorentzVector*)Reco_mu_trk_4mom->At(Reco_QQ_mulpt[iQQ])) < MAX_deltaR ||
156 deltaR((TLorentzVector*)Mc_mu_4mom->At(Mc_QQmumi_indx[0]),
157 (TLorentzVector*)Reco_mu_trk_4mom->At(Reco_QQ_mulpt[iQQ])) < MAX_deltaR)
158 isMatched[1] = kTRUE;
159 }
160 break;
161 case 3: // global-calo
162 //first MC muon
163 if(Reco_QQ_muhpt[iQQ] < Reco_mu_glb_size){//sanity check
164 if(deltaR((TLorentzVector*)Mc_mu_4mom->At(Mc_QQmupl_indx[0]),
165 (TLorentzVector*)Reco_mu_glb_4mom->At(Reco_QQ_muhpt[iQQ])) < MAX_deltaR ||
166 deltaR((TLorentzVector*)Mc_mu_4mom->At(Mc_QQmumi_indx[0]),
167 (TLorentzVector*)Reco_mu_glb_4mom->At(Reco_QQ_muhpt[iQQ])) < MAX_deltaR)
168 isMatched[0] = kTRUE;
169 }
170 //second MC muon
171 if(Reco_QQ_mulpt[iQQ] < Reco_mu_cal_size){//sanity check
172 if(deltaR((TLorentzVector*)Mc_mu_4mom->At(Mc_QQmupl_indx[0]),
173 (TLorentzVector*)Reco_mu_cal_4mom->At(Reco_QQ_mulpt[iQQ])) < MAX_deltaR ||
174 deltaR((TLorentzVector*)Mc_mu_4mom->At(Mc_QQmumi_indx[0]),
175 (TLorentzVector*)Reco_mu_cal_4mom->At(Reco_QQ_mulpt[iQQ])) < MAX_deltaR)
176 isMatched[1] = kTRUE;
177 }
178 break;
179 case 4: // tracker-calo
180 //first MC muon
181 if(Reco_QQ_muhpt[iQQ] < Reco_mu_trk_size){//sanity check
182 if(deltaR((TLorentzVector*)Mc_mu_4mom->At(Mc_QQmupl_indx[0]),
183 (TLorentzVector*)Reco_mu_trk_4mom->At(Reco_QQ_muhpt[iQQ])) < MAX_deltaR ||
184 deltaR((TLorentzVector*)Mc_mu_4mom->At(Mc_QQmumi_indx[0]),
185 (TLorentzVector*)Reco_mu_trk_4mom->At(Reco_QQ_muhpt[iQQ])) < MAX_deltaR)
186 isMatched[0] = kTRUE;
187 }
188 //second MC muon
189 if(Reco_QQ_mulpt[iQQ] < Reco_mu_cal_size){//sanity check
190 if(deltaR((TLorentzVector*)Mc_mu_4mom->At(Mc_QQmupl_indx[0]),
191 (TLorentzVector*)Reco_mu_cal_4mom->At(Reco_QQ_mulpt[iQQ])) < MAX_deltaR ||
192 deltaR((TLorentzVector*)Mc_mu_4mom->At(Mc_QQmumi_indx[0]),
193 (TLorentzVector*)Reco_mu_cal_4mom->At(Reco_QQ_mulpt[iQQ])) < MAX_deltaR)
194 isMatched[1] = kTRUE;
195 }
196 break;
197 case 5: // calo-calo
198 //first MC muon
199 if(Reco_QQ_muhpt[iQQ] < Reco_mu_cal_size){//sanity check
200 if(deltaR((TLorentzVector*)Mc_mu_4mom->At(Mc_QQmupl_indx[0]),
201 (TLorentzVector*)Reco_mu_cal_4mom->At(Reco_QQ_muhpt[iQQ])) < MAX_deltaR ||
202 deltaR((TLorentzVector*)Mc_mu_4mom->At(Mc_QQmumi_indx[0]),
203 (TLorentzVector*)Reco_mu_cal_4mom->At(Reco_QQ_muhpt[iQQ])) < MAX_deltaR)
204 isMatched[0] = kTRUE;
205 }
206 //second MC muon
207 if(Reco_QQ_mulpt[iQQ] < Reco_mu_cal_size){//sanity check
208 if(deltaR((TLorentzVector*)Mc_mu_4mom->At(Mc_QQmupl_indx[0]),
209 (TLorentzVector*)Reco_mu_cal_4mom->At(Reco_QQ_mulpt[iQQ])) < MAX_deltaR ||
210 deltaR((TLorentzVector*)Mc_mu_4mom->At(Mc_QQmumi_indx[0]),
211 (TLorentzVector*)Reco_mu_cal_4mom->At(Reco_QQ_mulpt[iQQ])) < MAX_deltaR)
212 isMatched[1] = kTRUE;
213 }
214 break;
215 default:
216 break;
217 }
218 //check the MC matching before continuing...
219 if(isMatched[0] == kFALSE || isMatched[1] == kFALSE)
220 continue;
221 }
222
223 if (iCut>0) {
224 // distance cut
225 switch (Reco_QQ_type[iQQ]) {
226 case 0: // global-global
227 if ( !(fabs(Reco_mu_glb_d0[Reco_QQ_muhpt[iQQ]]) < MAX_d0_trk &&
228 fabs(Reco_mu_glb_dz[Reco_QQ_muhpt[iQQ]]) < MAX_dz_trk &&
229 fabs(Reco_mu_glb_d0[Reco_QQ_mulpt[iQQ]]) < MAX_d0_trk &&
230 fabs(Reco_mu_glb_dz[Reco_QQ_mulpt[iQQ]]) < MAX_dz_trk) )
231 continue;
232 break;
233 case 1: // global-tracker
234 if ( !(fabs(Reco_mu_glb_d0[Reco_QQ_muhpt[iQQ]]) < MAX_d0_trk &&
235 fabs(Reco_mu_glb_dz[Reco_QQ_muhpt[iQQ]]) < MAX_dz_trk &&
236 fabs(Reco_mu_trk_d0[Reco_QQ_mulpt[iQQ]]) < MAX_d0_trk &&
237 fabs(Reco_mu_trk_dz[Reco_QQ_mulpt[iQQ]]) < MAX_dz_trk) )
238 continue;
239 break;
240 case 2: // tracker-tracker
241 if ( !(fabs(Reco_mu_trk_d0[Reco_QQ_muhpt[iQQ]]) < MAX_d0_trk &&
242 fabs(Reco_mu_trk_dz[Reco_QQ_muhpt[iQQ]]) < MAX_dz_trk &&
243 fabs(Reco_mu_trk_d0[Reco_QQ_mulpt[iQQ]]) < MAX_d0_trk &&
244 fabs(Reco_mu_trk_dz[Reco_QQ_mulpt[iQQ]]) < MAX_dz_trk) )
245 continue;
246 break;
247 case 3: // global-calo
248 if ( !(fabs(Reco_mu_glb_d0[Reco_QQ_muhpt[iQQ]]) < MAX_d0_trk &&
249 fabs(Reco_mu_glb_dz[Reco_QQ_muhpt[iQQ]]) < MAX_dz_trk &&
250 fabs(Reco_mu_cal_d0[Reco_QQ_mulpt[iQQ]]) < MAX_d0_trk &&
251 fabs(Reco_mu_cal_dz[Reco_QQ_mulpt[iQQ]]) < MAX_dz_trk) )
252 continue;
253 break;
254 case 4: // tracker-calo
255 if ( !(fabs(Reco_mu_trk_d0[Reco_QQ_muhpt[iQQ]]) < MAX_d0_trk &&
256 fabs(Reco_mu_trk_dz[Reco_QQ_muhpt[iQQ]]) < MAX_dz_trk &&
257 fabs(Reco_mu_cal_d0[Reco_QQ_mulpt[iQQ]]) < MAX_d0_trk &&
258 fabs(Reco_mu_cal_dz[Reco_QQ_mulpt[iQQ]]) < MAX_dz_trk) )
259 continue;
260 break;
261 case 5: // calo-calo
262 if ( !(fabs(Reco_mu_cal_d0[Reco_QQ_muhpt[iQQ]]) < MAX_d0_trk &&
263 fabs(Reco_mu_cal_dz[Reco_QQ_muhpt[iQQ]]) < MAX_dz_trk &&
264 fabs(Reco_mu_cal_d0[Reco_QQ_mulpt[iQQ]]) < MAX_d0_trk &&
265 fabs(Reco_mu_cal_dz[Reco_QQ_mulpt[iQQ]]) < MAX_dz_trk) )
266 continue;
267 break;
268 default:
269 break;
270 }
271 }
272
273 if (iCut>1) {
274 // for tracker muons: TM2DCompatibilityTight || TMLastStationOptimizedLowPtTight
275 // for calo muons: calo compatibility
276 switch (Reco_QQ_type[iQQ]) {
277 case 0: // global-global
278 break;
279 case 1: // global-tracker
280 if ((Reco_mu_trk_PIDmask[Reco_QQ_mulpt[iQQ]] & (int) pow(2,5))/(int) pow(2,5) == 0 &&
281 (Reco_mu_trk_PIDmask[Reco_QQ_mulpt[iQQ]] & (int) pow(2,8))/(int) pow(2,8) == 0)
282 continue;
283 break;
284 case 2: // tracker-tracker
285 if ( ((Reco_mu_trk_PIDmask[Reco_QQ_muhpt[iQQ]] & (int) pow(2,5))/(int) pow(2,5) == 0 &&
286 (Reco_mu_trk_PIDmask[Reco_QQ_muhpt[iQQ]] & (int) pow(2,8))/(int) pow(2,8) == 0) ||
287 ((Reco_mu_trk_PIDmask[Reco_QQ_mulpt[iQQ]] & (int) pow(2,5))/(int) pow(2,5) == 0 &&
288 (Reco_mu_trk_PIDmask[Reco_QQ_mulpt[iQQ]] & (int) pow(2,8))/(int) pow(2,8) == 0) )
289 continue;
290 break;
291 case 3: // global-calo
292 if ( !(Reco_mu_cal_caloComp[Reco_QQ_mulpt[iQQ]] > MIN_caloComp) )
293 continue;
294 break;
295 case 4: // tracker-calo
296 if ( ((Reco_mu_trk_PIDmask[Reco_QQ_muhpt[iQQ]] & (int) pow(2,5))/(int) pow(2,5) == 0 &&
297 (Reco_mu_trk_PIDmask[Reco_QQ_muhpt[iQQ]] & (int) pow(2,8))/(int) pow(2,8) == 0) ||
298 Reco_mu_cal_caloComp[Reco_QQ_mulpt[iQQ]] < MIN_caloComp)
299 continue;
300 break;
301 case 5: // calo-calo
302 if ( !(Reco_mu_cal_caloComp[Reco_QQ_muhpt[iQQ]] > MIN_caloComp &&
303 Reco_mu_cal_caloComp[Reco_QQ_mulpt[iQQ]] > MIN_caloComp) )
304 continue;
305 break;
306 default:
307 break;
308 }
309 }
310
311 if (iCut>2) {
312 // for global: chi2/ndf (global fit)
313 // for tracker and calo: chi2/ndf (track fit) && #hits
314 switch (Reco_QQ_type[iQQ]) {
315 case 0: // global-global
316 if ( !(Reco_mu_glb_normChi2[Reco_QQ_muhpt[iQQ]] < MAX_normchi2_glb &&
317 Reco_mu_glb_normChi2[Reco_QQ_mulpt[iQQ]] < MAX_normchi2_glb) )
318 continue;
319 break;
320 case 1: // global-tracker
321 if ( !(Reco_mu_glb_normChi2[Reco_QQ_muhpt[iQQ]] < MAX_normchi2_glb &&
322 Reco_mu_trk_normChi2[Reco_QQ_mulpt[iQQ]] < MAX_normchi2_trk &&
323 Reco_mu_trk_nhitstrack[Reco_QQ_mulpt[iQQ]] > MIN_nhits_trk) )
324 continue;
325 break;
326 case 2: // tracker-tracker
327 if ( !(Reco_mu_trk_normChi2[Reco_QQ_muhpt[iQQ]] < MAX_normchi2_trk &&
328 Reco_mu_trk_nhitstrack[Reco_QQ_muhpt[iQQ]] > MIN_nhits_trk &&
329 Reco_mu_trk_normChi2[Reco_QQ_mulpt[iQQ]] < MAX_normchi2_trk &&
330 Reco_mu_trk_nhitstrack[Reco_QQ_mulpt[iQQ]] > MIN_nhits_trk) )
331 continue;
332 break;
333 case 3: // global-calo
334 if ( !(Reco_mu_glb_normChi2[Reco_QQ_muhpt[iQQ]] < MAX_normchi2_glb &&
335 Reco_mu_cal_normChi2[Reco_QQ_mulpt[iQQ]] < MAX_normchi2_cal &&
336 Reco_mu_cal_nhitstrack[Reco_QQ_mulpt[iQQ]] > MIN_nhits_trk) )
337 continue;
338 break;
339 case 4: // tracker-calo
340 if ( !(Reco_mu_trk_normChi2[Reco_QQ_muhpt[iQQ]] < MAX_normchi2_trk &&
341 Reco_mu_trk_nhitstrack[Reco_QQ_muhpt[iQQ]] > MIN_nhits_trk &&
342 Reco_mu_cal_normChi2[Reco_QQ_mulpt[iQQ]] < MAX_normchi2_cal &&
343 Reco_mu_cal_nhitstrack[Reco_QQ_mulpt[iQQ]] > MIN_nhits_trk) )
344 continue;
345 break;
346 case 5: // calo-calo
347 if ( !(Reco_mu_cal_normChi2[Reco_QQ_muhpt[iQQ]] < MAX_normchi2_cal &&
348 Reco_mu_cal_nhitstrack[Reco_QQ_muhpt[iQQ]] > MIN_nhits_trk &&
349 Reco_mu_cal_normChi2[Reco_QQ_mulpt[iQQ]] < MAX_normchi2_cal &&
350 Reco_mu_cal_nhitstrack[Reco_QQ_mulpt[iQQ]] > MIN_nhits_trk) )
351 continue;
352 break;
353 default:
354 break;
355 }
356 }
357
358 if (iCut>3) {
359 // dimuon vertex probability
360 if (Reco_QQ_probChi2[iQQ]<= MIN_vtxprob_jpsi)
361 continue;
362 }
363
364 nComb[chargeID][Reco_QQ_type[iQQ]][iCut]++;
365
366 //mass histograms
367 hReco_mass[chargeID][Reco_QQ_type[iQQ]][iCut]->Fill(mass_Reco_QQ);
368 if(chargeID == 1 || chargeID == 2) //LS histo
369 hReco_mass[3][Reco_QQ_type[iQQ]][iCut]->Fill(mass_Reco_QQ);
370
371 //fill pT and rap histos only for the signal:
372 for(int iMassW = 0; iMassW < kNbDimuonSet; iMassW++){
373 if(mass_Reco_QQ > massMIN[iMassW] && mass_Reco_QQ < massMAX[iMassW]){
374 //pT histograms
375 hReco_pT[iMassW][chargeID][Reco_QQ_type[iQQ]][iCut]->Fill(pT_Reco_QQ);
376 if(chargeID == 1 || chargeID == 2) //LS histo
377 hReco_pT[iMassW][3][Reco_QQ_type[iQQ]][iCut]->Fill(pT_Reco_QQ);
378 //rap histograms
379 hReco_rap[iMassW][chargeID][Reco_QQ_type[iQQ]][iCut]->Fill(rap_Reco_QQ);
380 if(chargeID == 1 || chargeID == 2) //LS histo
381 hReco_rap[iMassW][3][Reco_QQ_type[iQQ]][iCut]->Fill(rap_Reco_QQ);
382
383 //correlating the single muons:
384 if(etaCorr)
385 hReco_eta1_eta2[iMassW][chargeID][Reco_QQ_type[iQQ]][iCut]->Fill(Reco_QQ_posMu->Eta(),Reco_QQ_negMu->Eta());
386
387 // prepare single muon indices to avoid double counting
388 int glMuID[10000], trMuID[10000],caMuID[10000];
389 for (int i=0; i<10000; i++){
390 glMuID[i] = -1;
391 trMuID[i] = -1;
392 caMuID[i] = -1;
393 }
394
395 // fill single muon quantities
396 // if (iCut==0) { // d0 and dz
397 // switch (Reco_QQ_type[iQQ]) {
398 // case 0: // global+global
399 // hChi2GlobalFit[0][chargeID]->Fill(Reco_mu_glb_normChi2[Reco_QQ_muhpt[iQQ]]);
400 // hChi2GlobalFit[0][chargeID]->Fill(Reco_mu_glb_normChi2[Reco_QQ_mulpt[iQQ]]);
401 // break;
402 // case 1: // global+tracker
403 // hChi2GlobalFit[0][chargeID]->Fill(Reco_mu_glb_normChi2[Reco_QQ_muhpt[iQQ]]);
404 // hChi2TrackerFit[0][chargeID]->Fill(Reco_mu_trk_normChi2[Reco_QQ_mulpt[iQQ]]);
405 // break;
406 // case 2: // tracker+tracker
407 // hChi2TrackerFit[0][chargeID]->Fill(Reco_mu_trk_normChi2[Reco_QQ_muhpt[iQQ]]);
408 // hChi2TrackerFit[0][chargeID]->Fill(Reco_mu_trk_normChi2[Reco_QQ_mulpt[iQQ]]);
409 // break;
410 // case 3: // global+calo
411 // hChi2GlobalFit[0][chargeID]->Fill(Reco_mu_glb_normChi2[Reco_QQ_muhpt[iQQ]]);
412 // hChi2TrackerFit[0][chargeID]->Fill(Reco_mu_cal_normChi2[Reco_QQ_mulpt[iQQ]]);
413 // break;
414 // case 4: // tracker+calo
415 // hChi2TrackerFit[0][chargeID]->Fill(Reco_mu_trk_normChi2[Reco_QQ_muhpt[iQQ]]);
416 // hChi2TrackerFit[0][chargeID]->Fill(Reco_mu_cal_normChi2[Reco_QQ_mulpt[iQQ]]);
417 // break;
418 // case 5: // calo+calo
419 // hChi2TrackerFit[0][chargeID]->Fill(Reco_mu_cal_normChi2[Reco_QQ_muhpt[iQQ]]);
420 // hChi2TrackerFit[0][chargeID]->Fill(Reco_mu_cal_normChi2[Reco_QQ_mulpt[iQQ]]);
421 // break;
422 // default:
423 // break;
424 // }
425 // }
426
427 if (iCut==1) { // comp, after d0 and dz
428 switch (Reco_QQ_type[iQQ]) {
429 case 0: // global+global
430 break;
431 case 1: // global+tracker
432 if (trMuID[Reco_QQ_mulpt[iQQ]]!=1) {
433 hTM2DCompatibilityTight[iMassW][chargeID]->Fill((Reco_mu_trk_PIDmask[Reco_QQ_mulpt[iQQ]] & (int) pow(2,5))/(int) pow(2,5));
434 hTMLastStationOptimizedLowPtLoose[iMassW][chargeID]->Fill((Reco_mu_trk_PIDmask[Reco_QQ_mulpt[iQQ]] & (int) pow(2,8))/(int) pow(2,8));
435 hTrackerMuonArbitrated[iMassW][chargeID]->Fill((Reco_mu_trk_PIDmask[Reco_QQ_mulpt[iQQ]] & (int)pow(2,1))/(int)pow(2,1));
436 hCaloCompatibility[iMassW][chargeID][1]->Fill(Reco_mu_trk_caloComp[Reco_QQ_mulpt[iQQ]]);
437 trMuID[Reco_QQ_mulpt[iQQ]]=1;
438 }
439 break;
440 case 2: // tracker+tracker
441 if (trMuID[Reco_QQ_muhpt[iQQ]]!=1) {
442 hTM2DCompatibilityTight[iMassW][chargeID]->Fill((Reco_mu_trk_PIDmask[Reco_QQ_muhpt[iQQ]] & (int) pow(2,5))/(int) pow(2,5));
443 hTMLastStationOptimizedLowPtLoose[iMassW][chargeID]->Fill((Reco_mu_trk_PIDmask[Reco_QQ_muhpt[iQQ]] & (int) pow(2,8))/(int) pow(2,8));
444 hTrackerMuonArbitrated[iMassW][chargeID]->Fill((Reco_mu_trk_PIDmask[Reco_QQ_muhpt[iQQ]] & (int)pow(2,1))/(int)pow(2,1));
445
446 hCaloCompatibility[iMassW][chargeID][1]->Fill(Reco_mu_trk_caloComp[Reco_QQ_muhpt[iQQ]]);
447 trMuID[Reco_QQ_muhpt[iQQ]]=1;
448 }
449 if (trMuID[Reco_QQ_mulpt[iQQ]]!=1) {
450 hTM2DCompatibilityTight[iMassW][chargeID]->Fill((Reco_mu_trk_PIDmask[Reco_QQ_mulpt[iQQ]] & (int) pow(2,5))/(int) pow(2,5));
451 hTMLastStationOptimizedLowPtLoose[iMassW][chargeID]->Fill((Reco_mu_trk_PIDmask[Reco_QQ_mulpt[iQQ]] & (int) pow(2,8))/(int) pow(2,8));
452 hTrackerMuonArbitrated[iMassW][chargeID]->Fill((Reco_mu_trk_PIDmask[Reco_QQ_mulpt[iQQ]] & (int)pow(2,1))/(int)pow(2,1));
453
454 hCaloCompatibility[iMassW][chargeID][1]->Fill(Reco_mu_trk_caloComp[Reco_QQ_mulpt[iQQ]]);
455 trMuID[Reco_QQ_mulpt[iQQ]]=1;
456 }
457 break;
458 case 3: // global+calo
459 if (caMuID[Reco_QQ_mulpt[iQQ]]!=1) {
460 if (Reco_mu_cal_normChi2[Reco_QQ_mulpt[iQQ]] < MAX_normchi2_cal &&
461 Reco_mu_cal_nhitstrack[Reco_QQ_mulpt[iQQ]] > MIN_nhits_trk)
462 hCaloCompatibility[iMassW][chargeID][2]->Fill(Reco_mu_cal_caloComp[Reco_QQ_mulpt[iQQ]]);
463 caMuID[Reco_QQ_mulpt[iQQ]]=1;
464 }
465 break;
466 case 4: // tracker+calo
467 if (trMuID[Reco_QQ_muhpt[iQQ]]!=1) {
468 hTM2DCompatibilityTight[iMassW][chargeID]->Fill((Reco_mu_trk_PIDmask[Reco_QQ_muhpt[iQQ]] & (int) pow(2,5))/(int) pow(2,5));
469 hTMLastStationOptimizedLowPtLoose[iMassW][chargeID]->Fill((Reco_mu_trk_PIDmask[Reco_QQ_muhpt[iQQ]] & (int) pow(2,8))/(int) pow(2,8));
470 hTrackerMuonArbitrated[iMassW][chargeID]->Fill((Reco_mu_trk_PIDmask[Reco_QQ_muhpt[iQQ]] & (int)pow(2,1))/(int)pow(2,1));
471 hCaloCompatibility[iMassW][chargeID][1]->Fill(Reco_mu_trk_caloComp[Reco_QQ_muhpt[iQQ]]);
472 trMuID[Reco_QQ_muhpt[iQQ]]=1;
473 }
474 if (caMuID[Reco_QQ_mulpt[iQQ]]!=1) {
475 if (Reco_mu_cal_normChi2[Reco_QQ_mulpt[iQQ]] < MAX_normchi2_cal &&
476 Reco_mu_cal_nhitstrack[Reco_QQ_mulpt[iQQ]] > MIN_nhits_trk)
477 hCaloCompatibility[iMassW][chargeID][2]->Fill(Reco_mu_cal_caloComp[Reco_QQ_mulpt[iQQ]]);
478 caMuID[Reco_QQ_mulpt[iQQ]]=1;
479 }
480 break;
481 case 5: // calo+calo
482 if (caMuID[Reco_QQ_muhpt[iQQ]]!=1) {
483 if (Reco_mu_cal_normChi2[Reco_QQ_muhpt[iQQ]] < MAX_normchi2_cal &&
484 Reco_mu_cal_nhitstrack[Reco_QQ_muhpt[iQQ]] > MIN_nhits_trk) {
485 hCaloCompatibility[iMassW][chargeID][2]->Fill(Reco_mu_cal_caloComp[Reco_QQ_muhpt[iQQ]]);
486 caMuID[Reco_QQ_muhpt[iQQ]]=1;
487 }
488 }
489 if (caMuID[Reco_QQ_mulpt[iQQ]]!=1) {
490 if (Reco_mu_cal_normChi2[Reco_QQ_mulpt[iQQ]] < MAX_normchi2_cal &&
491 Reco_mu_cal_nhitstrack[Reco_QQ_mulpt[iQQ]] > MIN_nhits_trk) {
492 hCaloCompatibility[iMassW][chargeID][2]->Fill(Reco_mu_cal_caloComp[Reco_QQ_mulpt[iQQ]]);
493 caMuID[Reco_QQ_mulpt[iQQ]]=1;
494 }
495 }
496 break;
497 default:
498 break;
499 }
500 }
501
502 if (iCut==2) { // chi2, after comp
503 switch (Reco_QQ_type[iQQ]) {
504 case 0: // global+global
505 if (glMuID[Reco_QQ_muhpt[iQQ]]!=1) {
506 hChi2GlobalFit[iMassW][chargeID]->Fill(Reco_mu_glb_normChi2[Reco_QQ_muhpt[iQQ]]);
507 hNHitsSilicon[iMassW][chargeID][0]->Fill(Reco_mu_glb_nhitstrack[Reco_QQ_muhpt[iQQ]]);
508 glMuID[Reco_QQ_muhpt[iQQ]]=1;
509 }
510 if (glMuID[Reco_QQ_mulpt[iQQ]]!=1) {
511 hChi2GlobalFit[iMassW][chargeID]->Fill(Reco_mu_glb_normChi2[Reco_QQ_mulpt[iQQ]]);
512 hNHitsSilicon[iMassW][chargeID][0]->Fill(Reco_mu_glb_nhitstrack[Reco_QQ_mulpt[iQQ]]);
513 glMuID[Reco_QQ_mulpt[iQQ]]=1;
514 }
515 break;
516 case 1: // global+tracker
517 if (glMuID[Reco_QQ_muhpt[iQQ]]!=1) {
518 hChi2GlobalFit[iMassW][chargeID]->Fill(Reco_mu_glb_normChi2[Reco_QQ_muhpt[iQQ]]);
519 hNHitsSilicon[iMassW][chargeID][0]->Fill(Reco_mu_glb_nhitstrack[Reco_QQ_muhpt[iQQ]]);
520 glMuID[Reco_QQ_mulpt[iQQ]]=1;
521 }
522 if (trMuID[Reco_QQ_mulpt[iQQ]]!=1) {
523 hChi2TrackerFit[iMassW][chargeID][1]->Fill(Reco_mu_trk_normChi2[Reco_QQ_mulpt[iQQ]]);
524 hNHitsSilicon[iMassW][chargeID][1]->Fill(Reco_mu_trk_nhitstrack[Reco_QQ_mulpt[iQQ]]);
525 trMuID[Reco_QQ_mulpt[iQQ]]=1;
526 }
527 break;
528 case 2: // tracker+tracker
529 if (trMuID[Reco_QQ_muhpt[iQQ]]!=1) {
530 hChi2TrackerFit[iMassW][chargeID][1]->Fill(Reco_mu_trk_normChi2[Reco_QQ_muhpt[iQQ]]);
531 hNHitsSilicon[iMassW][chargeID][1]->Fill(Reco_mu_trk_nhitstrack[Reco_QQ_muhpt[iQQ]]);
532 trMuID[Reco_QQ_muhpt[iQQ]]=1;
533 }
534 if (trMuID[Reco_QQ_mulpt[iQQ]]!=1) {
535 hChi2TrackerFit[iMassW][chargeID][1]->Fill(Reco_mu_trk_normChi2[Reco_QQ_mulpt[iQQ]]);
536 hNHitsSilicon[iMassW][chargeID][1]->Fill(Reco_mu_trk_nhitstrack[Reco_QQ_mulpt[iQQ]]);
537 trMuID[Reco_QQ_mulpt[iQQ]]=1;
538 }
539 break;
540 case 3: // global+calo
541 if (glMuID[Reco_QQ_muhpt[iQQ]]!=1) {
542 hChi2GlobalFit[iMassW][chargeID]->Fill(Reco_mu_glb_normChi2[Reco_QQ_muhpt[iQQ]]);
543 hNHitsSilicon[iMassW][chargeID][0]->Fill(Reco_mu_glb_nhitstrack[Reco_QQ_muhpt[iQQ]]);
544 glMuID[Reco_QQ_muhpt[iQQ]]=1;
545 }
546 if (caMuID[Reco_QQ_mulpt[iQQ]]!=1) {
547 hChi2TrackerFit[iMassW][chargeID][2]->Fill(Reco_mu_cal_normChi2[Reco_QQ_mulpt[iQQ]]);
548 hNHitsSilicon[iMassW][chargeID][2]->Fill(Reco_mu_cal_nhitstrack[Reco_QQ_mulpt[iQQ]]);
549 caMuID[Reco_QQ_mulpt[iQQ]]=1;
550 }
551 break;
552 case 4: // tracker+calo
553 if (trMuID[Reco_QQ_muhpt[iQQ]]!=1) {
554 hChi2TrackerFit[iMassW][chargeID][1]->Fill(Reco_mu_trk_normChi2[Reco_QQ_muhpt[iQQ]]);
555 hNHitsSilicon[iMassW][chargeID][1]->Fill(Reco_mu_trk_nhitstrack[Reco_QQ_muhpt[iQQ]]);
556 trMuID[Reco_QQ_muhpt[iQQ]]=1;
557 }
558 if (caMuID[Reco_QQ_mulpt[iQQ]]!=1) {
559 hChi2TrackerFit[iMassW][chargeID][2]->Fill(Reco_mu_cal_normChi2[Reco_QQ_mulpt[iQQ]]);
560 hNHitsSilicon[iMassW][chargeID][2]->Fill(Reco_mu_cal_nhitstrack[Reco_QQ_mulpt[iQQ]]);
561 caMuID[Reco_QQ_mulpt[iQQ]]=1;
562 }
563 break;
564 case 5: // calo+calo
565 if (caMuID[Reco_QQ_muhpt[iQQ]]!=1) {
566 hChi2TrackerFit[iMassW][chargeID][2]->Fill(Reco_mu_cal_normChi2[Reco_QQ_muhpt[iQQ]]);
567 hNHitsSilicon[iMassW][chargeID][2]->Fill(Reco_mu_cal_nhitstrack[Reco_QQ_muhpt[iQQ]]);
568 caMuID[Reco_QQ_muhpt[iQQ]]=1;
569 }
570 if (caMuID[Reco_QQ_mulpt[iQQ]]!=1) {
571 hChi2TrackerFit[iMassW][chargeID][2]->Fill(Reco_mu_cal_normChi2[Reco_QQ_mulpt[iQQ]]);
572 hNHitsSilicon[iMassW][chargeID][2]->Fill(Reco_mu_cal_nhitstrack[Reco_QQ_mulpt[iQQ]]);
573 caMuID[Reco_QQ_mulpt[iQQ]]=1;
574 }
575 break;
576 default:
577 break;
578 }
579 }
580
581 if (iCut==3) { // dimuon vertex prob, after comp
582 hDimuVtxProb[iMassW][chargeID][Reco_QQ_type[iQQ]]->Fill(Reco_QQ_probChi2[iQQ]);
583 if(chargeID == 1 || chargeID == 2) //LS histo
584 hDimuVtxProb[iMassW][3][Reco_QQ_type[iQQ]]->Fill(Reco_QQ_probChi2[iQQ]);
585 }
586 } // fill for signal
587 }//different mass windows
588 }//cuts
589 }//#Reco_QQ combinations
590
591 for(int iCharge = 0; iCharge < kNbCharge; iCharge++)
592 for(int iCat = 0; iCat < kNbCath; iCat++)
593 for(int iCut = 0; iCut < kNbCuts; iCut++)
594 if(nComb[iCharge][iCat][iCut] > 0)
595 hMult_QQ[chargeID][iCat][iCut]->Fill(nComb[iCharge][iCat][iCut]);
596 for(int iCat = 0; iCat < kNbCath; iCat++)
597 for(int iCut = 0; iCut < kNbCuts; iCut++)
598 hMult_QQ[3][iCat][iCut]->Fill(nComb[1][iCat][iCut]+nComb[2][iCat][iCut]);
599
600 //now fill the histos for the "best" QQ pair only
601 Int_t iDBestQQ = -1;
602 if(Reco_QQ_size > 0){
603 countAllRecoQQEvents++;
604 iDBestQQ = theBestQQ();
605 // if(iDBestQQ >= 0){
606 // printf("best QQ is %d (out of %d) and is of type %s\n",
607 // iDBestQQ, Reco_QQ_size, oniaCatName[Reco_QQ_type[iDBestQQ]]);
608 // }
609 // else{
610 // printf("best QQ not identified... %d total QQbar\n", Reco_QQ_size);
611 // countNotFoundBestQQEvents++;
612 // }
613 if(iDBestQQ < 0){
614 countNotFoundBestQQEvents++;
615 for(int iRec = 0; iRec < Reco_QQ_size; iRec++){
616
617 if(Reco_QQ_sign[iRec] != 0) continue;
618 //we count ALL reco QQpairs in the event:
619 hRecoNoBest->Fill(Reco_QQ_type[iRec]);
620 }
621 }
622 else{
623
624 TLorentzVector *Reco_QQ = (TLorentzVector *) Reco_QQ_4mom->At(iDBestQQ);
625 Double_t mass_Reco_QQ = Reco_QQ->M();
626 Double_t pT_Reco_QQ = Reco_QQ->Pt();
627 Double_t rap_Reco_QQ = Reco_QQ->Rapidity();
628
629 chargeID = 0;
630 if(Reco_QQ_sign[iDBestQQ] < 0) chargeID = 1;
631 else if(Reco_QQ_sign[iDBestQQ] > 0) chargeID = 2;
632
633 Bool_t isMatched[2] = {kTRUE, kTRUE};//only relevant for MC; for data must be set to kTRUE
634 if (matchMC){
635
636 isMatched[0] = kFALSE;
637 isMatched[1] = kFALSE;
638
639 if(Mc_QQmupl_indx[0] < 2 && Mc_QQmumi_indx[0] < 2){ //sanity check
640
641 switch (Reco_QQ_type[iDBestQQ]) {
642 case 0: // global-global
643 //first MC muon
644 if(Reco_QQ_muhpt[iDBestQQ] < Reco_mu_glb_size){//sanity check
645 if(deltaR((TLorentzVector*)Mc_mu_4mom->At(Mc_QQmupl_indx[0]),
646 (TLorentzVector*)Reco_mu_glb_4mom->At(Reco_QQ_muhpt[iDBestQQ])) < MAX_deltaR ||
647 deltaR((TLorentzVector*)Mc_mu_4mom->At(Mc_QQmumi_indx[0]),
648 (TLorentzVector*)Reco_mu_glb_4mom->At(Reco_QQ_muhpt[iDBestQQ])) < MAX_deltaR)
649 isMatched[0] = kTRUE;
650 }
651 //second MC muon
652 if(Reco_QQ_mulpt[iDBestQQ] < Reco_mu_glb_size){//sanity check
653 if(deltaR((TLorentzVector*)Mc_mu_4mom->At(Mc_QQmupl_indx[0]),
654 (TLorentzVector*)Reco_mu_glb_4mom->At(Reco_QQ_mulpt[iDBestQQ])) < MAX_deltaR ||
655 deltaR((TLorentzVector*)Mc_mu_4mom->At(Mc_QQmumi_indx[0]),
656 (TLorentzVector*)Reco_mu_glb_4mom->At(Reco_QQ_mulpt[iDBestQQ])) < MAX_deltaR)
657 isMatched[1] = kTRUE;
658 }
659 break;
660 case 1: // global-tracker
661 //first MC muon
662 if(Reco_QQ_muhpt[iDBestQQ] < Reco_mu_glb_size){//sanity check
663 if(deltaR((TLorentzVector*)Mc_mu_4mom->At(Mc_QQmupl_indx[0]),
664 (TLorentzVector*)Reco_mu_glb_4mom->At(Reco_QQ_muhpt[iDBestQQ])) < MAX_deltaR ||
665 deltaR((TLorentzVector*)Mc_mu_4mom->At(Mc_QQmumi_indx[0]),
666 (TLorentzVector*)Reco_mu_glb_4mom->At(Reco_QQ_muhpt[iDBestQQ])) < MAX_deltaR)
667 isMatched[0] = kTRUE;
668 }
669 //second MC muon
670 if(Reco_QQ_mulpt[iDBestQQ] < Reco_mu_trk_size){//sanity check
671 if(deltaR((TLorentzVector*)Mc_mu_4mom->At(Mc_QQmupl_indx[0]),
672 (TLorentzVector*)Reco_mu_trk_4mom->At(Reco_QQ_mulpt[iDBestQQ])) < MAX_deltaR ||
673 deltaR((TLorentzVector*)Mc_mu_4mom->At(Mc_QQmumi_indx[0]),
674 (TLorentzVector*)Reco_mu_trk_4mom->At(Reco_QQ_mulpt[iDBestQQ])) < MAX_deltaR)
675 isMatched[1] = kTRUE;
676 }
677 break;
678 case 2: // tracker-tracker
679 //first MC muon
680 if(Reco_QQ_muhpt[iDBestQQ] < Reco_mu_trk_size){//sanity check
681 if(deltaR((TLorentzVector*)Mc_mu_4mom->At(Mc_QQmupl_indx[0]),
682 (TLorentzVector*)Reco_mu_trk_4mom->At(Reco_QQ_muhpt[iDBestQQ])) < MAX_deltaR ||
683 deltaR((TLorentzVector*)Mc_mu_4mom->At(Mc_QQmumi_indx[0]),
684 (TLorentzVector*)Reco_mu_trk_4mom->At(Reco_QQ_muhpt[iDBestQQ])) < MAX_deltaR)
685 isMatched[0] = kTRUE;
686 }
687 //second MC muon
688 if(Reco_QQ_mulpt[iDBestQQ] < Reco_mu_trk_size){//sanity check
689 if(deltaR((TLorentzVector*)Mc_mu_4mom->At(Mc_QQmupl_indx[0]),
690 (TLorentzVector*)Reco_mu_trk_4mom->At(Reco_QQ_mulpt[iDBestQQ])) < MAX_deltaR ||
691 deltaR((TLorentzVector*)Mc_mu_4mom->At(Mc_QQmumi_indx[0]),
692 (TLorentzVector*)Reco_mu_trk_4mom->At(Reco_QQ_mulpt[iDBestQQ])) < MAX_deltaR)
693 isMatched[1] = kTRUE;
694 }
695 break;
696 case 3: // global-calo
697 //first MC muon
698 if(Reco_QQ_muhpt[iDBestQQ] < Reco_mu_glb_size){//sanity check
699 if(deltaR((TLorentzVector*)Mc_mu_4mom->At(Mc_QQmupl_indx[0]),
700 (TLorentzVector*)Reco_mu_glb_4mom->At(Reco_QQ_muhpt[iDBestQQ])) < MAX_deltaR ||
701 deltaR((TLorentzVector*)Mc_mu_4mom->At(Mc_QQmumi_indx[0]),
702 (TLorentzVector*)Reco_mu_glb_4mom->At(Reco_QQ_muhpt[iDBestQQ])) < MAX_deltaR)
703 isMatched[0] = kTRUE;
704 }
705 //second MC muon
706 if(Reco_QQ_mulpt[iDBestQQ] < Reco_mu_cal_size){//sanity check
707 if(deltaR((TLorentzVector*)Mc_mu_4mom->At(Mc_QQmupl_indx[0]),
708 (TLorentzVector*)Reco_mu_cal_4mom->At(Reco_QQ_mulpt[iDBestQQ])) < MAX_deltaR ||
709 deltaR((TLorentzVector*)Mc_mu_4mom->At(Mc_QQmumi_indx[0]),
710 (TLorentzVector*)Reco_mu_cal_4mom->At(Reco_QQ_mulpt[iDBestQQ])) < MAX_deltaR)
711 isMatched[1] = kTRUE;
712 }
713 break;
714 case 4: // tracker-calo
715 //first MC muon
716 if(Reco_QQ_muhpt[iDBestQQ] < Reco_mu_trk_size){//sanity check
717 if(deltaR((TLorentzVector*)Mc_mu_4mom->At(Mc_QQmupl_indx[0]),
718 (TLorentzVector*)Reco_mu_trk_4mom->At(Reco_QQ_muhpt[iDBestQQ])) < MAX_deltaR ||
719 deltaR((TLorentzVector*)Mc_mu_4mom->At(Mc_QQmumi_indx[0]),
720 (TLorentzVector*)Reco_mu_trk_4mom->At(Reco_QQ_muhpt[iDBestQQ])) < MAX_deltaR)
721 isMatched[0] = kTRUE;
722 }
723 //second MC muon
724 if(Reco_QQ_mulpt[iDBestQQ] < Reco_mu_cal_size){//sanity check
725 if(deltaR((TLorentzVector*)Mc_mu_4mom->At(Mc_QQmupl_indx[0]),
726 (TLorentzVector*)Reco_mu_cal_4mom->At(Reco_QQ_mulpt[iDBestQQ])) < MAX_deltaR ||
727 deltaR((TLorentzVector*)Mc_mu_4mom->At(Mc_QQmumi_indx[0]),
728 (TLorentzVector*)Reco_mu_cal_4mom->At(Reco_QQ_mulpt[iDBestQQ])) < MAX_deltaR)
729 isMatched[1] = kTRUE;
730 }
731 break;
732 case 5: // calo-calo
733 //first MC muon
734 if(Reco_QQ_muhpt[iDBestQQ] < Reco_mu_cal_size){//sanity check
735 if(deltaR((TLorentzVector*)Mc_mu_4mom->At(Mc_QQmupl_indx[0]),
736 (TLorentzVector*)Reco_mu_cal_4mom->At(Reco_QQ_muhpt[iDBestQQ])) < MAX_deltaR ||
737 deltaR((TLorentzVector*)Mc_mu_4mom->At(Mc_QQmumi_indx[0]),
738 (TLorentzVector*)Reco_mu_cal_4mom->At(Reco_QQ_muhpt[iDBestQQ])) < MAX_deltaR)
739 isMatched[0] = kTRUE;
740 }
741 //second MC muon
742 if(Reco_QQ_mulpt[iDBestQQ] < Reco_mu_cal_size){//sanity check
743 if(deltaR((TLorentzVector*)Mc_mu_4mom->At(Mc_QQmupl_indx[0]),
744 (TLorentzVector*)Reco_mu_cal_4mom->At(Reco_QQ_mulpt[iDBestQQ])) < MAX_deltaR ||
745 deltaR((TLorentzVector*)Mc_mu_4mom->At(Mc_QQmumi_indx[0]),
746 (TLorentzVector*)Reco_mu_cal_4mom->At(Reco_QQ_mulpt[iDBestQQ])) < MAX_deltaR)
747 isMatched[1] = kTRUE;
748 }
749 break;
750 default:
751 break;
752 }
753 }
754 }
755 //check the MC matching before continuing...
756 if(isMatched[0] == kTRUE && isMatched[1] == kTRUE){
757
758 //mass histograms
759 hRecoBest_mass[chargeID][Reco_QQ_type[iDBestQQ]]->Fill(mass_Reco_QQ);
760 if(chargeID == 1 || chargeID == 2) //LS histo
761 hRecoBest_mass[3][Reco_QQ_type[iDBestQQ]]->Fill(mass_Reco_QQ);
762
763 //fill pT and rap histos only for the signal:
764 for(int iMassW = 0; iMassW < kNbDimuonSet; iMassW++){
765 if(mass_Reco_QQ > massMIN[iMassW] && mass_Reco_QQ < massMAX[iMassW]){
766 //pT histograms
767 hRecoBest_pT[iMassW][chargeID][Reco_QQ_type[iDBestQQ]]->Fill(pT_Reco_QQ);
768 if(chargeID == 1 || chargeID == 2) //LS histo
769 hRecoBest_pT[iMassW][3][Reco_QQ_type[iDBestQQ]]->Fill(pT_Reco_QQ);
770 //rap histograms
771 hRecoBest_rap[iMassW][chargeID][Reco_QQ_type[iDBestQQ]]->Fill(rap_Reco_QQ);
772 if(chargeID == 1 || chargeID == 2) //LS histo
773 hRecoBest_rap[iMassW][3][Reco_QQ_type[iDBestQQ]]->Fill(rap_Reco_QQ);
774 }
775 //===============================================================
776 //single muon histograms for best dimuon within 3.0 < M < 3.2 GeV
777 //and for dimuons with M > 2 GeV
778 //===============================================================
779 TLorentzVector *Reco_mu_1, *Reco_mu_2;
780 Int_t iMuCat[2] = {-1,-1};
781 if(Reco_QQ_type[iDBestQQ] == 0){//gl-gl
782 Reco_mu_1 = (TLorentzVector *) Reco_mu_glb_4mom->At(Reco_QQ_muhpt[iDBestQQ]);
783 Reco_mu_2 = (TLorentzVector *) Reco_mu_glb_4mom->At(Reco_QQ_mulpt[iDBestQQ]);
784 iMuCat[0] = 0; iMuCat[1] = 0;
785 }
786 else if(Reco_QQ_type[iDBestQQ] == 1){//gl-tr
787 Reco_mu_1 = (TLorentzVector *) Reco_mu_glb_4mom->At(Reco_QQ_muhpt[iDBestQQ]);
788 Reco_mu_2 = (TLorentzVector *) Reco_mu_trk_4mom->At(Reco_QQ_mulpt[iDBestQQ]);
789 iMuCat[0] = 0; iMuCat[1] = 1;
790 }
791 else if(Reco_QQ_type[iDBestQQ] == 2){//tr-tr
792 Reco_mu_1 = (TLorentzVector *) Reco_mu_trk_4mom->At(Reco_QQ_muhpt[iDBestQQ]);
793 Reco_mu_2 = (TLorentzVector *) Reco_mu_trk_4mom->At(Reco_QQ_mulpt[iDBestQQ]);
794 iMuCat[0] = 1; iMuCat[1] = 1;
795 }
796 else if(Reco_QQ_type[iDBestQQ] == 3){//gl-calo
797 Reco_mu_1 = (TLorentzVector *) Reco_mu_glb_4mom->At(Reco_QQ_muhpt[iDBestQQ]);
798 Reco_mu_2 = (TLorentzVector *) Reco_mu_cal_4mom->At(Reco_QQ_mulpt[iDBestQQ]);
799 iMuCat[0] = 0; iMuCat[1] = 2;
800 }
801 else if(Reco_QQ_type[iDBestQQ] == 4){//tr-calo
802 Reco_mu_1 = (TLorentzVector *) Reco_mu_trk_4mom->At(Reco_QQ_muhpt[iDBestQQ]);
803 Reco_mu_2 = (TLorentzVector *) Reco_mu_cal_4mom->At(Reco_QQ_mulpt[iDBestQQ]);
804 iMuCat[0] = 1; iMuCat[1] = 2;
805 }
806 else if(Reco_QQ_type[iDBestQQ] == 5){//calo-calo
807 Reco_mu_1 = (TLorentzVector *) Reco_mu_cal_4mom->At(Reco_QQ_muhpt[iDBestQQ]);
808 Reco_mu_2 = (TLorentzVector *) Reco_mu_cal_4mom->At(Reco_QQ_mulpt[iDBestQQ]);
809 iMuCat[0] = 2; iMuCat[1] = 2;
810 }
811
812 Double_t etaMu[2] = {Reco_mu_1->Eta(), Reco_mu_2->Eta()};
813 Double_t pTMu[2] = {Reco_mu_1->Pt(), Reco_mu_2->Pt()};
814 for(int iMu = 0; iMu < 2; iMu++){
815 hMuon_eta[iMassW][chargeID][iMuCat[iMu]]->Fill(etaMu[iMu]);
816 hMuon_pT[iMassW][chargeID][iMuCat[iMu]]->Fill(pTMu[iMu]);
817 hMuon_pT_eta[iMassW][chargeID][iMuCat[iMu]]->Fill(etaMu[iMu], pTMu[iMu]);
818 }
819 }//mass windows
820 }//matched: relevant for MC
821 }//best QQ
822 }//Reco_QQ_size > 0
823 }//events
824
825 if(countAllRecoQQEvents > 0){
826 hRecoNoBest->Scale(1./countAllRecoQQEvents);
827 printf("\n\n The number of events w/o identification of the best QQbar is %d out of %d, i.e. %1.1f %%\n", countNotFoundBestQQEvents, countAllRecoQQEvents, 100.*countNotFoundBestQQEvents/countAllRecoQQEvents);
828 }
829 }
830
831 //=============================
832 int ProjectQQ::theBestQQ() {
833
834 int theBest = -1;
835 float thehighestPt = -1.;
836
837 for (int iqq=0; iqq<Reco_QQ_size; iqq++) {
838
839 if (Reco_QQ_sign[iqq] == 0 && Reco_QQ_type[iqq] == 0 ) {
840
841 int thehptMu = Reco_QQ_muhpt[iqq]; if (thehptMu >= Reco_mu_glb_size) continue;
842 int thelptMu = Reco_QQ_mulpt[iqq]; if (thelptMu >= Reco_mu_glb_size) continue;
843 if (Reco_QQ_probChi2[iqq] > MIN_vtxprob_jpsi &&
844 //Reco_mu_glb_nhitstrack[thehptMu] > MIN_nhits_trk &&
845 Reco_mu_glb_normChi2[thehptMu] < MAX_normchi2_glb &&
846 //(((Reco_mu_glb_nhitsPixB[thehptMu] + Reco_mu_glb_nhitsPixE[thehptMu]) > MIN_nhits_pixel) || ((Reco_mu_glb_nhitsPixB[thehptMu] + Reco_mu_glb_nhitsPixE[thehptMu]) > MIN_nhits_pixel-1 && Reco_mu_glb_nhitsPix1Hit[thehptMu] == 1)) &&
847 fabs(Reco_mu_glb_d0[thehptMu]) < MAX_d0_trk &&
848 fabs(Reco_mu_glb_dz[thehptMu]) < MAX_dz_trk &&
849 //Reco_mu_glb_nhitstrack[thelptMu] > MIN_nhits_trk &&
850 Reco_mu_glb_normChi2[thelptMu] < MAX_normchi2_glb &&
851 //(((Reco_mu_glb_nhitsPixB[thelptMu] + Reco_mu_glb_nhitsPixE[thelptMu]) > MIN_nhits_pixel) || ((Reco_mu_glb_nhitsPixB[thelptMu] + Reco_mu_glb_nhitsPixE[thelptMu]) > MIN_nhits_pixel-1 && Reco_mu_glb_nhitsPix1Hit[thelptMu] == 1)) &&
852 fabs(Reco_mu_glb_d0[thelptMu]) < MAX_d0_trk &&
853 fabs(Reco_mu_glb_dz[thelptMu]) < MAX_dz_trk
854 ) {
855 return iqq;
856 }
857 }
858 }
859
860 for (int iqq=0; iqq<Reco_QQ_size; iqq++) {
861
862 if (Reco_QQ_sign[iqq] == 0 && Reco_QQ_type[iqq] == 1 ) {
863
864 int thehptMu = Reco_QQ_muhpt[iqq]; if (thehptMu >= Reco_mu_glb_size) continue;
865 int thelptMu = Reco_QQ_mulpt[iqq]; if (thelptMu >= Reco_mu_trk_size) continue;
866
867 if ( Reco_QQ_probChi2[iqq] > MIN_vtxprob_jpsi &&
868 //Reco_mu_glb_nhitstrack[thehptMu] > MIN_nhits_trk &&
869 Reco_mu_glb_normChi2[thehptMu] < MAX_normchi2_glb &&
870 //(((Reco_mu_glb_nhitsPixB[thehptMu] + Reco_mu_glb_nhitsPixE[thehptMu]) > MIN_nhits_pixel) || ((Reco_mu_glb_nhitsPixB[thehptMu] + Reco_mu_glb_nhitsPixE[thehptMu]) > MIN_nhits_pixel-1 && Reco_mu_glb_nhitsPix1Hit[thehptMu] == 1)) &&
871 fabs(Reco_mu_glb_d0[thehptMu]) < MAX_d0_trk &&
872 fabs(Reco_mu_glb_dz[thehptMu]) < MAX_dz_trk &&
873 Reco_mu_trk_nhitstrack[thelptMu] > MIN_nhits_trk &&
874 ((Reco_mu_trk_PIDmask[thelptMu] & (int)pow(2,5))/(int)pow(2,5) > 0 || (Reco_mu_trk_PIDmask[thelptMu] & (int)pow(2,8))/(int)pow(2,8) > 0) &&
875 //(Reco_mu_trk_nhitsPixB[thelptMu] + Reco_mu_trk_nhitsPixE[thelptMu]) > MIN_nhits_pixel &&
876 Reco_mu_trk_normChi2[thelptMu] < MAX_normchi2_trk &&
877 fabs(Reco_mu_trk_d0[thelptMu]) < MAX_d0_trk &&
878 fabs(Reco_mu_trk_dz[thelptMu]) < MAX_dz_trk) {
879
880 TLorentzVector *theTrMumom = (TLorentzVector*)Reco_mu_trk_4mom->At(thelptMu);
881 if (theTrMumom->Perp() > thehighestPt) {
882 thehighestPt = theTrMumom->Perp();
883 theBest = iqq;
884 }
885 }
886 }
887 }
888
889 if (theBest >= 0) return theBest;
890
891 for (int iqq=0; iqq<Reco_QQ_size; iqq++) {
892
893 if (Reco_QQ_sign[iqq] == 0 && Reco_QQ_type[iqq] == 2 ) {
894
895 int thehptMu = Reco_QQ_muhpt[iqq]; if (thehptMu >= Reco_mu_trk_size) continue;
896 int thelptMu = Reco_QQ_mulpt[iqq]; if (thelptMu >= Reco_mu_trk_size) continue;
897
898 if ( Reco_QQ_probChi2[iqq] > MIN_vtxprob_jpsi &&
899 Reco_mu_trk_nhitstrack[thehptMu] > MIN_nhits_trk &&
900 ((Reco_mu_trk_PIDmask[thehptMu] & (int)pow(2,5))/(int)pow(2,5) > 0 || (Reco_mu_trk_PIDmask[thehptMu] & (int)pow(2,8))/(int)pow(2,8) > 0) &&
901 //(Reco_mu_trk_nhitsPixB[thehptMu] + Reco_mu_trk_nhitsPixE[thehptMu]) > MIN_nhits_pixel &&
902 Reco_mu_trk_normChi2[thehptMu] < MAX_normchi2_trk &&
903 fabs(Reco_mu_trk_d0[thehptMu]) < MAX_d0_trk &&
904 fabs(Reco_mu_trk_dz[thehptMu]) < MAX_dz_trk &&
905 Reco_mu_trk_nhitstrack[thelptMu] > MIN_nhits_trk &&
906 ((Reco_mu_trk_PIDmask[thelptMu] & (int)pow(2,5))/(int)pow(2,5) > 0 || (Reco_mu_trk_PIDmask[thelptMu] & (int)pow(2,8))/(int)pow(2,8) > 0) &&
907 //(Reco_mu_trk_nhitsPixB[thelptMu] + Reco_mu_trk_nhitsPixE[thelptMu]) > MIN_nhits_pixel &&
908 Reco_mu_trk_normChi2[thelptMu] < MAX_normchi2_trk &&
909 fabs(Reco_mu_trk_d0[thelptMu]) < MAX_d0_trk &&
910 fabs(Reco_mu_trk_dz[thelptMu]) < MAX_dz_trk) {
911
912 TLorentzVector *theTrMumom = (TLorentzVector*)Reco_mu_trk_4mom->At(thehptMu);
913 if (theTrMumom->Perp() > thehighestPt) {
914 thehighestPt = theTrMumom->Perp();
915 theBest = iqq;
916 }
917 }
918 }
919 }
920
921 if (theBest >= 0) return theBest;
922
923 for (int iqq=0; iqq<Reco_QQ_size; iqq++) {
924
925 if (Reco_QQ_sign[iqq] == 0 && Reco_QQ_type[iqq] == 3 ) {
926
927 int thehptMu = Reco_QQ_muhpt[iqq];
928 int thelptMu = Reco_QQ_mulpt[iqq];
929 if (thelptMu >= Reco_mu_cal_size) {
930 // cout << "Non deve succedere! cmIndex = " << thelptMu+1 << " cmSize = " << Reco_mu_cal_size << endl;
931 continue;
932 }
933
934 if ( //Reco_mu_glb_nhitstrack[thehptMu] > MIN_nhits_trk &&
935 Reco_mu_glb_normChi2[thehptMu] < MAX_normchi2_glb &&
936 Reco_mu_cal_nhitstrack[thelptMu] > MIN_nhits_trk &&
937 Reco_mu_cal_normChi2[thelptMu] < MAX_normchi2_cal &&
938 Reco_mu_cal_caloComp[thelptMu] > MIN_caloComp) {
939
940 TLorentzVector *theCaMumom = (TLorentzVector*)Reco_mu_cal_4mom->At(thelptMu);
941 if (theCaMumom->Perp() > thehighestPt) {
942 thehighestPt = theCaMumom->Perp();
943 theBest = iqq;
944 }
945 }
946 }
947 }
948
949 return theBest;
950 }
951 //=========================================
952 Double_t ProjectQQ::deltaR(TLorentzVector* t, TLorentzVector* u){
953
954 return sqrt(pow(t->Eta()-u->Eta(),2) +pow(PhiInRange(t->Phi()-u->Phi()),2));
955 }
956 //=========================================
957 double ProjectQQ::PhiInRange(double phi){
958
959 double phiout = phi;
960
961 if( phiout > 2*M_PI || phiout < -2*M_PI) {
962 phiout = fmod( phiout, 2*M_PI);
963 }
964 if (phiout <= -M_PI) phiout += 2*M_PI;
965 else if (phiout > M_PI) phiout -= 2*M_PI;
966
967 return phiout;
968 }