ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/MitPhysics/Mods/src/GeneratorMod.cc
Revision: 1.58
Committed: Sat Mar 13 20:51:27 2010 UTC (15 years, 1 month ago) by ceballos
Content type: text/plain
Branch: MAIN
CVS Tags: Mit_014c, Mit_014b, Mit_014a, Mit_014, Mit_014pre3, Mit_014pre2, Mit_014pre1, Mit_013d, Mit_013c, Mit_013b, Mit_013a, Mit_013, Mit_013pre1
Changes since 1.57: +37 -12 lines
Log Message:
adding some info

File Contents

# Content
1 // $Id: GeneratorMod.cc,v 1.57 2009/12/06 14:59:28 ceballos Exp $
2
3 #include "MitPhysics/Mods/interface/GeneratorMod.h"
4 #include "MitCommon/MathTools/interface/MathUtils.h"
5 #include "MitAna/DataTree/interface/MetCol.h"
6 #include "MitAna/DataTree/interface/MCParticleCol.h"
7 #include "MitPhysics/Init/interface/ModNames.h"
8 #include <TH1D.h>
9 #include <TH2D.h>
10
11 using namespace mithep;
12
13 ClassImp(mithep::GeneratorMod)
14
15 //--------------------------------------------------------------------------------------------------
16 GeneratorMod::GeneratorMod(const char *name, const char *title) :
17 BaseMod(name,title),
18 fIsData(kFALSE),
19 fPrintDebug(kFALSE),
20 fCopyArrays(kFALSE),
21 fMCPartName(Names::gkMCPartBrn),
22 fMCMETName(ModNames::gkMCMETName),
23 fMCLeptonsName(ModNames::gkMCLeptonsName),
24 fMCAllLeptonsName(ModNames::gkMCAllLeptonsName),
25 fMCTausName(ModNames::gkMCTausName),
26 fMCNeutrinosName(ModNames::gkMCNeutrinosName),
27 fMCQuarksName(ModNames::gkMCQuarksName),
28 fMCqqHsName(ModNames::gkMCqqHsName),
29 fMCBosonsName(ModNames::gkMCBosonsName),
30 fMCPhotonsName(ModNames::gkMCPhotonsName),
31 fMCRadPhotonsName(ModNames::gkMCRadPhotonsName),
32 fMCISRPhotonsName(ModNames::gkMCISRPhotonsName),
33 fPtLeptonMin(0.0),
34 fEtaLeptonMax(5.0),
35 fPtPhotonMin(0.0),
36 fEtaPhotonMax(5.0),
37 fPtRadPhotonMin(0.0),
38 fEtaRadPhotonMax(5.0),
39 fPdgIdCut(0),
40 fMassMinCut(-FLT_MAX),
41 fMassMaxCut(FLT_MAX),
42 fApplyISRFilter(kFALSE),
43 fParticles(0)
44 {
45 // Constructor
46 fGenLeptons = new MCParticleArr();
47 fGenLeptons->SetName(TString("Pub") + ModNames::gkMCLeptonsName);
48 fGenAllLeptons = new MCParticleArr();
49 fGenAllLeptons->SetName(TString("Pub") + ModNames::gkMCAllLeptonsName);
50 fGenTaus = new MCParticleArr();
51 fGenTaus->SetName(TString("Pub") + ModNames::gkMCTausName);
52 fGenNeutrinos = new MCParticleArr();
53 fGenNeutrinos->SetName(TString("Pub") + ModNames::gkMCNeutrinosName);
54 fGenQuarks = new MCParticleArr();
55 fGenQuarks->SetName(TString("Pub") + ModNames::gkMCQuarksName);
56 fGenqqHs = new MCParticleArr();
57 fGenqqHs->SetName(TString("Pub") + ModNames::gkMCqqHsName);
58 fGenBosons = new MCParticleArr();
59 fGenBosons->SetName(TString("Pub") + ModNames::gkMCBosonsName);
60 fGenPhotons = new MCParticleArr();
61 fGenPhotons->SetName(TString("Pub") + ModNames::gkMCPhotonsName);
62 fGenRadPhotons = new MCParticleArr();
63 fGenRadPhotons->SetName(TString("Pub") + ModNames::gkMCRadPhotonsName);
64 fGenISRPhotons = new MCParticleArr();
65 fGenISRPhotons->SetName(TString("Pub") + ModNames::gkMCISRPhotonsName);
66 }
67
68
69 //--------------------------------------------------------------------------------------------------
70 GeneratorMod::~GeneratorMod()
71 {
72 // Destructor.
73
74 delete fGenLeptons;
75 delete fGenAllLeptons;
76 delete fGenTaus;
77 delete fGenNeutrinos;
78 delete fGenQuarks;
79 delete fGenqqHs;
80 delete fGenBosons;
81 delete fGenPhotons;
82 delete fGenRadPhotons;
83 delete fGenISRPhotons;
84 }
85
86 //--------------------------------------------------------------------------------------------------
87 void GeneratorMod::Process()
88 {
89 // Process entries of the tree.
90
91 // these arrays will be filled in the loop of particles
92 MetOArr *GenMet = new MetOArr;
93 GenMet->SetName(fMCMETName);
94 GenMet->SetOwner(kTRUE);
95 MCParticleOArr *GenLeptons = new MCParticleOArr;
96 GenLeptons->SetName(fMCLeptonsName);
97 MCParticleOArr *GenAllLeptons = new MCParticleOArr;
98 GenAllLeptons->SetName(fMCAllLeptonsName);
99 MCParticleOArr *GenTaus = new MCParticleOArr;
100 GenTaus->SetName(fMCTausName);
101 GenTaus->SetOwner(kTRUE);
102 MCParticleOArr *GenNeutrinos = new MCParticleOArr;
103 GenNeutrinos->SetName(fMCNeutrinosName);
104 MCParticleOArr *GenQuarks = new MCParticleOArr;
105 GenQuarks->SetName(fMCQuarksName);
106 MCParticleOArr *GenqqHs = new MCParticleOArr;
107 GenqqHs->SetName(fMCqqHsName);
108 MCParticleOArr *GenBosons = new MCParticleOArr;
109 GenBosons->SetName(fMCBosonsName);
110 MCParticleOArr *GenPhotons = new MCParticleOArr;
111 GenPhotons->SetName(fMCPhotonsName);
112 MCParticleOArr *GenRadPhotons = new MCParticleOArr;
113 GenRadPhotons->SetName(fMCRadPhotonsName);
114 MCParticleOArr *GenISRPhotons = new MCParticleOArr;
115 GenISRPhotons->SetName(fMCISRPhotonsName);
116
117 MCParticleOArr *GenTempMG0 = new MCParticleOArr;
118
119 Bool_t isOld = kFALSE;
120 Int_t sumV[2] = {0, 0}; // W, Z
121 Int_t sumVVFlavor[9] = {0, 0, 0, 0, 0, 0, 0, 0, 0};
122 Double_t totalMET[3] = {0.0, 0.0, 0.0};
123 Bool_t isqqH = kFALSE;
124 Double_t ptMin = 999999.;
125 if(fIsData == kFALSE){
126 if (fPrintDebug)
127 printf("\n************ Next Event ************\n\n");
128
129 // load MCParticle branch
130 LoadEventObject(fMCPartName, fParticles);
131
132 for (UInt_t i=0; i<fParticles->GetEntries(); ++i) {
133 const MCParticle *p = fParticles->At(i);
134
135 if(fPrintDebug)
136 p->Print("l");
137
138 // rad photons, includes gamma from WWGamma vertex.
139 if( p->Is(MCParticle::kGamma) && p->Status() == 1 &&
140 p->Pt() > fPtRadPhotonMin && p->AbsEta() < fEtaRadPhotonMax &&
141 p->DistinctMother() && p->DistinctMother()->Status() == 3 &&
142 (p->DistinctMother()->Is(MCParticle::kEl) || p->DistinctMother()->Is(MCParticle::kMu) ||
143 p->DistinctMother()->Is(MCParticle::kTau) || p->DistinctMother()->Is(MCParticle::kW))
144 ) {
145 CompositeParticle *object = new CompositeParticle();
146 object->AddDaughter(p);
147 object->AddDaughter(p->DistinctMother());
148 if(object->Mass() > 1.0 || p->DistinctMother()->Is(MCParticle::kW)) GenRadPhotons->Add(p);
149 delete object;
150 }
151
152 // ISR photons
153 if( p->Is(MCParticle::kGamma) && p->Status() == 1 &&
154 p->Pt() > fPtRadPhotonMin && p->AbsEta() < fEtaRadPhotonMax &&
155 p->DistinctMother() && p->DistinctMother()->IsParton()
156 ) {
157 GenISRPhotons->Add(p);
158 }
159
160 // MET computation at generation level
161 if (p->Status() == 1 && !p->IsNeutrino()) {
162 totalMET[0] = totalMET[0] + p->Px();
163 totalMET[1] = totalMET[1] + p->Py();
164 totalMET[2] = totalMET[2] + p->Pz();
165 }
166
167 if (!p->IsGenerated()) continue;
168
169 if ((((p->Is(MCParticle::kEl) || p->Is(MCParticle::kMu)) &&
170 !p->HasMother(MCParticle::kTau, kFALSE)) || p->Is(MCParticle::kTau)) &&
171 p->Status() == 3 &&
172 (p->HasMother(MCParticle::kW, kFALSE) || p->HasMother(MCParticle::kZ, kFALSE))) {
173 if(p->Pt() < ptMin) ptMin = p->Pt();
174 }
175
176 // all muons/electrons
177 if ((p->Is(MCParticle::kEl) || p->Is(MCParticle::kMu)) && p->Status() == 1) {
178 if (p->Pt() > fPtLeptonMin && p->AbsEta() < fEtaLeptonMax) {
179 GenAllLeptons->Add(p);
180 }
181 Bool_t isGoodLepton = kFALSE;
182 const MCParticle *pm = p;
183 while (pm->HasMother() && isGoodLepton == kFALSE) {
184 if (pm->PdgId() == 92) // string reached, terminate loop
185 break;
186 if (pm->Mother()->Is(MCParticle::kZ) || pm->Mother()->Is(MCParticle::kW) ||
187 pm->Mother()->Is(MCParticle::kZp) || pm->Mother()->Is(MCParticle::kWp) ||
188 pm->Mother()->Is(MCParticle::kH)) {
189 GenLeptons->Add(p);
190 isGoodLepton = kTRUE;
191 break;
192 } else if (pm->Mother()->Is(MCParticle::kPi0) || pm->Mother()->Is(MCParticle::kEta)) {
193 // this is fake, but it is a trick to get rid of these cases and abort the loop
194 break;
195 }
196 pm = pm->Mother();
197 }
198 }
199
200 // hadronic taus
201 else if (p->Is(MCParticle::kTau) && p->Status() == 2) {
202 if (!p->HasDaughter(MCParticle::kEl) && !p->HasDaughter(MCParticle::kMu)) {
203 const MCParticle *tv = p->FindDaughter(MCParticle::kTauNu);
204 if (tv) {
205 MCParticle *pm_f = new MCParticle(*p);
206 pm_f->SetMom(p->Px()-tv->Px(), p->Py()-tv->Py(),
207 p->Pz()-tv->Pz(), p->E()-tv->E());
208 GenTaus->AddOwned(pm_f);
209 } else {
210 //SendError(kWarning, "Process", "Could not find a tau neutrino!");
211 }
212 }
213 }
214
215 // neutrinos
216 else if (p->Status() == 1 && p->IsNeutrino()) {
217 GenNeutrinos->Add(p);
218 }
219
220 // quarks from W/Z decays or top particles
221 else if (p->IsQuark() && p->HasMother()) {
222 if (p->Mother()->Is(MCParticle::kZ) || p->Mother()->Is(MCParticle::kW) ||
223 p->Is(MCParticle::kTop) || p->Mother()->Is(MCParticle::kTop)) {
224 GenQuarks->Add(p);
225 }
226 }
227
228 // qqH, information about the forward jets
229 else if (isqqH == kFALSE && p->Is(MCParticle::kH)) {
230 isqqH = kTRUE;
231 const MCParticle *pq1 = fParticles->At(i-1);
232 const MCParticle *pq2 = fParticles->At(i-2);
233 if (!pq1 || !pq2) {
234 SendError(kWarning, "Process", "Could not find quark pair!");
235 } else if (pq1->IsQuark() && pq2->IsQuark() &&
236 pq1->HasMother() && pq2->HasMother() &&
237 pq1->Mother() == pq2->Mother()) {
238 GenqqHs->Add(pq1);
239 GenqqHs->Add(pq2);
240 }
241
242 if (p->Status() == 3)
243 GenBosons->Add(p); // take higgs boson in account here rather in next else if
244 }
245
246 // information about bosons: W, Z, h, Z', W', H0, A0, H+
247 else if ((p->Status() == 3 &&
248 (p->Is(MCParticle::kZ) || p->Is(MCParticle::kW) || p->Is(MCParticle::kH) ||
249 p->Is(MCParticle::kZp) || p->Is(MCParticle::kZpp) ||
250 p->Is(MCParticle::kH0) || p->Is(MCParticle::kA0) || p->Is(MCParticle::kHp))) ||
251 (p->Status() == 2 &&
252 (p->Is(MCParticle::kJPsi) || p->Is(MCParticle::kUpsilon)))) {
253 GenBosons->Add(p);
254 if (p->Is(MCParticle::kW)) sumV[0]++;
255 else if(p->Is(MCParticle::kZ)) sumV[1]++;
256 if (p->Is(MCParticle::kW) && p->HasDaughter(MCParticle::kMu) &&
257 p->HasDaughter(MCParticle::kMuNu))
258 sumVVFlavor[0]++;
259 else if(p->Is(MCParticle::kW) && p->HasDaughter(MCParticle::kEl) &&
260 p->HasDaughter(MCParticle::kElNu))
261 sumVVFlavor[1]++;
262 else if(p->Is(MCParticle::kW) && p->HasDaughter(MCParticle::kTau) &&
263 p->HasDaughter(MCParticle::kTauNu))
264 sumVVFlavor[2]++;
265 else if(p->Is(MCParticle::kZ) && p->HasDaughter(MCParticle::kMu,kTRUE) &&
266 p->HasDaughter(-1*MCParticle::kMu,kTRUE))
267 sumVVFlavor[3]++;
268 else if(p->Is(MCParticle::kZ) && p->HasDaughter(MCParticle::kEl,kTRUE) &&
269 p->HasDaughter(-1*MCParticle::kEl,kTRUE))
270 sumVVFlavor[4]++;
271 else if(p->Is(MCParticle::kZ) && p->HasDaughter(MCParticle::kTau,kTRUE) &&
272 p->HasDaughter(-1*MCParticle::kTau,kTRUE))
273 sumVVFlavor[5]++;
274 else if(p->Is(MCParticle::kZ) && p->HasDaughter(MCParticle::kMuNu,kTRUE) &&
275 p->HasDaughter(-1*MCParticle::kMuNu,kTRUE))
276 sumVVFlavor[6]++;
277 else if(p->Is(MCParticle::kZ) && p->HasDaughter(MCParticle::kElNu,kTRUE) &&
278 p->HasDaughter(-1*MCParticle::kElNu,kTRUE))
279 sumVVFlavor[7]++;
280 else if(p->Is(MCParticle::kZ) && p->HasDaughter(MCParticle::kTauNu,kTRUE) &&
281 p->HasDaughter(-1*MCParticle::kTauNu,kTRUE))
282 sumVVFlavor[8]++;
283 }
284
285 // photons
286 else if (p->Status() == 1 && p->Is(MCParticle::kGamma) &&
287 p->Pt() > fPtPhotonMin && p->AbsEta() < fEtaPhotonMax) {
288 GenPhotons->Add(p);
289 }
290
291 // W/Z -> lnu for Madgraph
292 if (p->IsParton() && p->NDaughters() >= 2) {
293 CompositeParticle *diBoson = new CompositeParticle();
294 if (p->HasDaughter(MCParticle::kMu) && p->HasDaughter(MCParticle::kMuNu)) {
295 isOld = kFALSE;
296 for(UInt_t nl = 0; nl < GenTempMG0->GetEntries(); nl++){
297 if(p->FindDaughter(MCParticle::kMu) == GenTempMG0->At(nl)) {
298 isOld = kTRUE;
299 break;
300 }
301 }
302 if(isOld == kFALSE){
303 GenTempMG0->Add(p->FindDaughter(MCParticle::kMu));
304 diBoson->AddDaughter(p->FindDaughter(MCParticle::kMu));
305 diBoson->AddDaughter(p->FindDaughter(MCParticle::kMuNu));
306 sumV[0]++;
307 sumVVFlavor[0]++;
308 if (GetFillHist())
309 hDVMass[0]->Fill(TMath::Min(diBoson->Mass(),199.999));
310 const MCParticle *tmp_mu = p->FindDaughter(MCParticle::kMu);
311 while (tmp_mu->HasDaughter(MCParticle::kMu) &&
312 tmp_mu->FindDaughter(MCParticle::kMu)->IsGenerated())
313 tmp_mu = tmp_mu->FindDaughter(MCParticle::kMu);
314
315 GenLeptons->Add(tmp_mu);
316 }
317 }
318 if (p->HasDaughter(MCParticle::kEl) && p->HasDaughter(MCParticle::kElNu)) {
319 isOld = kFALSE;
320 for(UInt_t nl = 0; nl < GenTempMG0->GetEntries(); nl++){
321 if(p->FindDaughter(MCParticle::kEl) == GenTempMG0->At(nl)) {
322 isOld = kTRUE;
323 break;
324 }
325 }
326 if(isOld == kFALSE){
327 GenTempMG0->Add(p->FindDaughter(MCParticle::kEl));
328 diBoson->AddDaughter(p->FindDaughter(MCParticle::kEl));
329 diBoson->AddDaughter(p->FindDaughter(MCParticle::kElNu));
330 sumV[0]++;
331 sumVVFlavor[1]++;
332 if (GetFillHist())
333 hDVMass[1]->Fill(TMath::Min(diBoson->Mass(),199.999));
334 const MCParticle *tmp_e = p->FindDaughter(MCParticle::kEl);
335 while (tmp_e->HasDaughter(MCParticle::kEl) &&
336 tmp_e->FindDaughter(MCParticle::kEl)->IsGenerated())
337 tmp_e = tmp_e->FindDaughter(MCParticle::kEl);
338 GenLeptons->Add(tmp_e);
339 }
340 }
341 if (p->HasDaughter(MCParticle::kTau) && p->HasDaughter(MCParticle::kTauNu)) {
342 isOld = kFALSE;
343 for(UInt_t nl = 0; nl < GenTempMG0->GetEntries(); nl++){
344 if(p->FindDaughter(MCParticle::kTau) == GenTempMG0->At(nl)) {
345 isOld = kTRUE;
346 break;
347 }
348 }
349 if(isOld == kFALSE){
350 GenTempMG0->Add(p->FindDaughter(MCParticle::kTau));
351 diBoson->AddDaughter(p->FindDaughter(MCParticle::kTau));
352 diBoson->AddDaughter(p->FindDaughter(MCParticle::kTauNu));
353 sumV[0]++;
354 sumVVFlavor[2]++;
355 if (GetFillHist())
356 hDVMass[2]->Fill(TMath::Min(diBoson->Mass(),199.999));
357 const MCParticle *tau = p->FindDaughter(MCParticle::kTau);
358 if (tau->HasDaughter(MCParticle::kMu))
359 GenLeptons->Add(tau->FindDaughter(MCParticle::kMu));
360 if (tau->HasDaughter(MCParticle::kEl))
361 GenLeptons->Add(tau->FindDaughter(MCParticle::kEl));
362 if (tau->HasDaughter(MCParticle::kTau)) {
363 const MCParticle *tau_second = tau->FindDaughter(MCParticle::kTau);
364 if (tau_second->HasDaughter(MCParticle::kMu))
365 GenLeptons->Add(tau_second->FindDaughter(MCParticle::kMu));
366 if (tau_second->HasDaughter(MCParticle::kEl))
367 GenLeptons->Add(tau_second->FindDaughter(MCParticle::kEl));
368 }
369 }
370 }
371 if (p->HasDaughter(MCParticle::kMu,kTRUE) && p->HasDaughter(-1*MCParticle::kMu,kTRUE)) {
372 isOld = kFALSE;
373 for(UInt_t nl = 0; nl < GenTempMG0->GetEntries(); nl++){
374 if(p->FindDaughter(MCParticle::kMu,kTRUE) == GenTempMG0->At(nl)) {
375 isOld = kTRUE;
376 break;
377 }
378 }
379 if(isOld == kFALSE){
380 GenTempMG0->Add(p->FindDaughter(MCParticle::kMu,kTRUE));
381 diBoson->AddDaughter(p->FindDaughter(MCParticle::kMu,kTRUE));
382 diBoson->AddDaughter(p->FindDaughter(-1*MCParticle::kMu,kTRUE));
383 sumV[1]++;
384 sumVVFlavor[3]++;
385 if (GetFillHist())
386 hDVMass[3]->Fill(TMath::Min(diBoson->Mass(),199.999));
387 const MCParticle *tmp_mu0 = p->FindDaughter(MCParticle::kMu,kTRUE);
388 while (tmp_mu0->HasDaughter(MCParticle::kMu) &&
389 tmp_mu0->FindDaughter(MCParticle::kMu)->IsGenerated())
390 tmp_mu0 = tmp_mu0->FindDaughter(MCParticle::kMu);
391 const MCParticle *tmp_mu1 = p->FindDaughter(-1*MCParticle::kMu,kTRUE);
392 while (tmp_mu1->HasDaughter(MCParticle::kMu) &&
393 tmp_mu1->FindDaughter(MCParticle::kMu)->IsGenerated())
394 tmp_mu1 = tmp_mu1->FindDaughter(MCParticle::kMu);
395 GenLeptons->Add(tmp_mu0);
396 GenLeptons->Add(tmp_mu1);
397 }
398 }
399 if (p->HasDaughter(MCParticle::kEl,kTRUE) && p->HasDaughter(-1*MCParticle::kEl,kTRUE)) {
400 isOld = kFALSE;
401 for(UInt_t nl = 0; nl < GenTempMG0->GetEntries(); nl++){
402 if(p->FindDaughter(MCParticle::kEl,kTRUE) == GenTempMG0->At(nl)) {
403 isOld = kTRUE;
404 break;
405 }
406 }
407 if(isOld == kFALSE){
408 GenTempMG0->Add(p->FindDaughter(MCParticle::kEl,kTRUE));
409 diBoson->AddDaughter(p->FindDaughter(MCParticle::kEl,kTRUE));
410 diBoson->AddDaughter(p->FindDaughter(-1*MCParticle::kEl,kTRUE));
411 sumV[1]++;
412 sumVVFlavor[4]++;
413 if (GetFillHist())
414 hDVMass[4]->Fill(TMath::Min(diBoson->Mass(),199.999));
415 const MCParticle *tmp_e0 = p->Daughter(0);
416 while (tmp_e0->HasDaughter(MCParticle::kEl) &&
417 tmp_e0->FindDaughter(MCParticle::kEl)->IsGenerated())
418 tmp_e0 = tmp_e0->FindDaughter(MCParticle::kEl);
419 const MCParticle *tmp_e1 = p->Daughter(1);
420 while (tmp_e1->HasDaughter(MCParticle::kEl) &&
421 tmp_e1->FindDaughter(MCParticle::kEl)->IsGenerated())
422 tmp_e1 = tmp_e1->FindDaughter(MCParticle::kEl);
423 GenLeptons->Add(tmp_e0);
424 GenLeptons->Add(tmp_e1);
425 }
426 }
427 if (p->HasDaughter(MCParticle::kTau,kTRUE) && p->HasDaughter(-1*MCParticle::kTau,kTRUE)) {
428 isOld = kFALSE;
429 for(UInt_t nl = 0; nl < GenTempMG0->GetEntries(); nl++){
430 if(p->FindDaughter(MCParticle::kTau,kTRUE) == GenTempMG0->At(nl)) {
431 isOld = kTRUE;
432 break;
433 }
434 }
435 if(isOld == kFALSE){
436 GenTempMG0->Add(p->FindDaughter(MCParticle::kTau,kTRUE));
437 diBoson->AddDaughter(p->FindDaughter(MCParticle::kTau,kTRUE));
438 diBoson->AddDaughter(p->FindDaughter(-1*MCParticle::kTau,kTRUE));
439 sumV[1]++;
440 sumVVFlavor[5]++;
441 if (GetFillHist())
442 hDVMass[5]->Fill(TMath::Min(diBoson->Mass(),199.999));
443 const MCParticle *tau0 = p->Daughter(0);
444 if (tau0->HasDaughter(MCParticle::kMu))
445 GenLeptons->Add(tau0->FindDaughter(MCParticle::kMu));
446 if (tau0->HasDaughter(MCParticle::kEl))
447 GenLeptons->Add(tau0->FindDaughter(MCParticle::kEl));
448 const MCParticle *tau1 = p->Daughter(1);
449 if (tau1->HasDaughter(MCParticle::kMu))
450 GenLeptons->Add(tau1->FindDaughter(MCParticle::kMu));
451 if (tau1->HasDaughter(MCParticle::kEl))
452 GenLeptons->Add(tau1->FindDaughter(MCParticle::kEl));
453 if (tau0->HasDaughter(MCParticle::kTau)) {
454 const MCParticle *tau0_second = tau0->FindDaughter(MCParticle::kTau);
455 if (tau0_second->HasDaughter(MCParticle::kMu))
456 GenLeptons->Add(tau0_second->FindDaughter(MCParticle::kMu));
457 if (tau0_second->HasDaughter(MCParticle::kEl))
458 GenLeptons->Add(tau0_second->FindDaughter(MCParticle::kEl));
459 }
460 if (tau1->HasDaughter(MCParticle::kTau)) {
461 const MCParticle *tau1_second = tau1->FindDaughter(MCParticle::kTau);
462 if (tau1_second->HasDaughter(MCParticle::kMu))
463 GenLeptons->Add(tau1_second->FindDaughter(MCParticle::kMu));
464 if (tau1_second->HasDaughter(MCParticle::kEl))
465 GenLeptons->Add(tau1_second->FindDaughter(MCParticle::kEl));
466 }
467 }
468 }
469 if (p->HasDaughter(MCParticle::kMuNu,kTRUE) && p->HasDaughter(-1*MCParticle::kMuNu,kTRUE)) {
470 isOld = kFALSE;
471 for(UInt_t nl = 0; nl < GenTempMG0->GetEntries(); nl++){
472 if(p->FindDaughter(MCParticle::kMuNu,kTRUE) == GenTempMG0->At(nl)) {
473 isOld = kTRUE;
474 break;
475 }
476 }
477 if(isOld == kFALSE){
478 GenTempMG0->Add(p->FindDaughter(MCParticle::kMuNu,kTRUE));
479 diBoson->AddDaughter(p->FindDaughter(MCParticle::kMuNu,kTRUE));
480 diBoson->AddDaughter(p->FindDaughter(-1*MCParticle::kMuNu,kTRUE));
481 sumV[1]++;
482 sumVVFlavor[6]++;
483 if (GetFillHist())
484 hDVMass[6]->Fill(TMath::Min(diBoson->Mass(),199.999));
485 }
486 }
487 if (p->HasDaughter(MCParticle::kElNu,kTRUE) && p->HasDaughter(-1*MCParticle::kElNu,kTRUE)) {
488 isOld = kFALSE;
489 for(UInt_t nl = 0; nl < GenTempMG0->GetEntries(); nl++){
490 if(p->FindDaughter(MCParticle::kElNu,kTRUE) == GenTempMG0->At(nl)) {
491 isOld = kTRUE;
492 break;
493 }
494 }
495 if(isOld == kFALSE){
496 GenTempMG0->Add(p->FindDaughter(MCParticle::kElNu,kTRUE));
497 diBoson->AddDaughter(p->FindDaughter(MCParticle::kElNu,kTRUE));
498 diBoson->AddDaughter(p->FindDaughter(-1*MCParticle::kElNu,kTRUE));
499 sumV[1]++;
500 sumVVFlavor[7]++;
501 if (GetFillHist())
502 hDVMass[7]->Fill(TMath::Min(diBoson->Mass(),199.999));
503 }
504 }
505 if (p->HasDaughter(MCParticle::kTauNu,kTRUE) && p->HasDaughter(-1*MCParticle::kTauNu,kTRUE)) {
506 isOld = kFALSE;
507 for(UInt_t nl = 0; nl < GenTempMG0->GetEntries(); nl++){
508 if(p->FindDaughter(MCParticle::kTauNu,kTRUE) == GenTempMG0->At(nl)) {
509 isOld = kTRUE;
510 break;
511 }
512 }
513 if(isOld == kFALSE){
514 GenTempMG0->Add(p->FindDaughter(MCParticle::kTauNu,kTRUE));
515 diBoson->AddDaughter(p->FindDaughter(MCParticle::kTauNu,kTRUE));
516 diBoson->AddDaughter(p->FindDaughter(-1*MCParticle::kTauNu,kTRUE));
517 sumV[1]++;
518 sumVVFlavor[8]++;
519 if (GetFillHist())
520 hDVMass[8]->Fill(TMath::Min(diBoson->Mass(),199.999));
521 }
522 }
523 delete diBoson;
524 }
525
526 // t -> lnu for Madgraph
527 if (p->Is(MCParticle::kTop)) {
528 CompositeParticle *diBoson = new CompositeParticle();
529 if (p->HasDaughter(MCParticle::kMu) && p->HasDaughter(MCParticle::kMuNu)) {
530 diBoson->AddDaughter(p->FindDaughter(MCParticle::kMu));
531 diBoson->AddDaughter(p->FindDaughter(MCParticle::kMuNu));
532 if (GetFillHist())
533 hDVMass[9]->Fill(TMath::Min(diBoson->Mass(),199.999));
534 GenLeptons->Add(p->FindDaughter(MCParticle::kMu));
535 }
536 else if (p->HasDaughter(MCParticle::kEl) && p->HasDaughter(MCParticle::kElNu)) {
537 diBoson->AddDaughter(p->FindDaughter(MCParticle::kEl));
538 diBoson->AddDaughter(p->FindDaughter(MCParticle::kElNu));
539 if (GetFillHist())
540 hDVMass[10]->Fill(TMath::Min(diBoson->Mass(),199.999));
541 GenLeptons->Add(p->FindDaughter(MCParticle::kEl));
542 }
543 else if (p->HasDaughter(MCParticle::kTau) && p->HasDaughter(MCParticle::kTauNu)) {
544 diBoson->AddDaughter(p->FindDaughter(MCParticle::kTau));
545 diBoson->AddDaughter(p->FindDaughter(MCParticle::kTauNu));
546 if (GetFillHist())
547 hDVMass[11]->Fill(TMath::Min(diBoson->Mass(),199.999));
548 const MCParticle *tau = p->FindDaughter(MCParticle::kTau);
549 if (tau->HasDaughter(MCParticle::kMu))
550 GenLeptons->Add(tau->FindDaughter(MCParticle::kMu));
551 if (tau->HasDaughter(MCParticle::kEl))
552 GenLeptons->Add(tau->FindDaughter(MCParticle::kEl));
553 }
554 else if (!p->HasDaughter(MCParticle::kW)) {
555 for(UInt_t nd=0; nd<p->NDaughters(); ++nd)
556 if (p->Daughter(nd)->IsNot(MCParticle::kBottom) &&
557 p->Daughter(nd)->IsNot(MCParticle::kGamma))
558 diBoson->AddDaughter(p->Daughter(nd));
559 if (GetFillHist())
560 hDVMass[12]->Fill(TMath::Min(diBoson->Mass(),199.999));
561 }
562 delete diBoson;
563 }
564
565 // mass cut for given pid
566 if(fPdgIdCut && p->Is(fPdgIdCut) &&
567 (p->Mass() < fMassMinCut || p->Mass() > fMassMaxCut)) {
568 SkipEvent();
569 return;
570 }
571 } // end loop of particles
572
573 delete GenTempMG0;
574
575 // --------------------------------
576 // Begin special study about VVjets
577 // --------------------------------
578 if(sumV[0] + 4*sumV[1] == 2 || sumV[0] + 4*sumV[1] == 5 || sumV[0] + 4*sumV[1] == 8){
579 MCParticleOArr *GenTempMG1 = new MCParticleOArr;
580 Double_t diBosonMass[2] = {0., 0.};
581 for (UInt_t i=0; i<fParticles->GetEntries(); ++i) {
582 const MCParticle *p = fParticles->At(i);
583
584 if (p->IsParton() && p->NDaughters() >= 2) {
585 CompositeParticle *diBoson = new CompositeParticle();
586 if (p->HasDaughter(MCParticle::kMu) && p->HasDaughter(MCParticle::kMuNu)) {
587 isOld = kFALSE;
588 for(UInt_t nl = 0; nl < GenTempMG1->GetEntries(); nl++){
589 if(p->FindDaughter(MCParticle::kMu) == GenTempMG1->At(nl)) {
590 isOld = kTRUE;
591 break;
592 }
593 }
594 if(isOld == kFALSE){
595 GenTempMG1->Add(p->FindDaughter(MCParticle::kMu));
596 diBoson->AddDaughter(p->FindDaughter(MCParticle::kMu));
597 diBoson->AddDaughter(p->FindDaughter(MCParticle::kMuNu));
598 if (GetFillHist() && sumV[0] + 4*sumV[1] == 2)
599 hDVVMass[0]->Fill(TMath::Min(diBoson->Mass(),199.999));
600 if (GetFillHist() && sumV[0] + 4*sumV[1] == 5)
601 hDVVMass[1]->Fill(TMath::Min(diBoson->Mass(),199.999));
602 }
603 }
604 if (p->HasDaughter(MCParticle::kEl) && p->HasDaughter(MCParticle::kElNu)) {
605 isOld = kFALSE;
606 for(UInt_t nl = 0; nl < GenTempMG1->GetEntries(); nl++){
607 if(p->FindDaughter(MCParticle::kEl) == GenTempMG1->At(nl)) {
608 isOld = kTRUE;
609 break;
610 }
611 }
612 if(isOld == kFALSE){
613 GenTempMG1->Add(p->FindDaughter(MCParticle::kEl));
614 diBoson->AddDaughter(p->FindDaughter(MCParticle::kEl));
615 diBoson->AddDaughter(p->FindDaughter(MCParticle::kElNu));
616 if (GetFillHist() && sumV[0] + 4*sumV[1] == 2)
617 hDVVMass[2]->Fill(TMath::Min(diBoson->Mass(),199.999));
618 if (GetFillHist() && sumV[0] + 4*sumV[1] == 5)
619 hDVVMass[3]->Fill(TMath::Min(diBoson->Mass(),199.999));
620 }
621 }
622 if (p->HasDaughter(MCParticle::kTau) && p->HasDaughter(MCParticle::kTauNu)) {
623 isOld = kFALSE;
624 for(UInt_t nl = 0; nl < GenTempMG1->GetEntries(); nl++){
625 if(p->FindDaughter(MCParticle::kTau) == GenTempMG1->At(nl)) {
626 isOld = kTRUE;
627 break;
628 }
629 }
630 if(isOld == kFALSE){
631 GenTempMG1->Add(p->FindDaughter(MCParticle::kTau));
632 diBoson->AddDaughter(p->FindDaughter(MCParticle::kTau));
633 diBoson->AddDaughter(p->FindDaughter(MCParticle::kTauNu));
634 if (GetFillHist() && sumV[0] + 4*sumV[1] == 2)
635 hDVVMass[4]->Fill(TMath::Min(diBoson->Mass(),199.999));
636 if (GetFillHist() && sumV[0] + 4*sumV[1] == 5)
637 hDVVMass[5]->Fill(TMath::Min(diBoson->Mass(),199.999));
638 }
639 }
640 if (p->HasDaughter(MCParticle::kMu,kTRUE) && p->HasDaughter(-1*MCParticle::kMu,kTRUE)) {
641 isOld = kFALSE;
642 for(UInt_t nl = 0; nl < GenTempMG1->GetEntries(); nl++){
643 if(p->FindDaughter(MCParticle::kMu,kTRUE) == GenTempMG1->At(nl)) {
644 isOld = kTRUE;
645 break;
646 }
647 }
648 if(isOld == kFALSE){
649 GenTempMG1->Add(p->FindDaughter(MCParticle::kMu,kTRUE));
650 diBoson->AddDaughter(p->FindDaughter(MCParticle::kMu,kTRUE));
651 diBoson->AddDaughter(p->FindDaughter(-1*MCParticle::kMu,kTRUE));
652 if (GetFillHist() && sumV[0] + 4*sumV[1] == 5)
653 hDVVMass[6]->Fill(TMath::Min(diBoson->Mass(),199.999));
654 if (GetFillHist() && sumV[0] + 4*sumV[1] == 8)
655 hDVVMass[7]->Fill(TMath::Min(diBoson->Mass(),199.999));
656 }
657 }
658 if (p->HasDaughter(MCParticle::kEl,kTRUE) && p->HasDaughter(-1*MCParticle::kEl,kTRUE)) {
659 isOld = kFALSE;
660 for(UInt_t nl = 0; nl < GenTempMG1->GetEntries(); nl++){
661 if(p->FindDaughter(MCParticle::kEl,kTRUE) == GenTempMG1->At(nl)) {
662 isOld = kTRUE;
663 break;
664 }
665 }
666 if(isOld == kFALSE){
667 GenTempMG1->Add(p->FindDaughter(MCParticle::kEl,kTRUE));
668 diBoson->AddDaughter(p->FindDaughter(MCParticle::kEl,kTRUE));
669 diBoson->AddDaughter(p->FindDaughter(-1*MCParticle::kEl,kTRUE));
670 if (GetFillHist() && sumV[0] + 4*sumV[1] == 5)
671 hDVVMass[8]->Fill(TMath::Min(diBoson->Mass(),199.999));
672 if (GetFillHist() && sumV[0] + 4*sumV[1] == 8)
673 hDVVMass[9]->Fill(TMath::Min(diBoson->Mass(),199.999));
674 }
675 }
676 if (p->HasDaughter(MCParticle::kTau,kTRUE) && p->HasDaughter(-1*MCParticle::kTau,kTRUE)) {
677 isOld = kFALSE;
678 for(UInt_t nl = 0; nl < GenTempMG1->GetEntries(); nl++){
679 if(p->FindDaughter(MCParticle::kTau,kTRUE) == GenTempMG1->At(nl)) {
680 isOld = kTRUE;
681 break;
682 }
683 }
684 if(isOld == kFALSE){
685 GenTempMG1->Add(p->FindDaughter(MCParticle::kTau,kTRUE));
686 diBoson->AddDaughter(p->FindDaughter(MCParticle::kTau,kTRUE));
687 diBoson->AddDaughter(p->FindDaughter(-1*MCParticle::kTau,kTRUE));
688 if (GetFillHist() && sumV[0] + 4*sumV[1] == 5)
689 hDVVMass[10]->Fill(TMath::Min(diBoson->Mass(),199.999));
690 if (GetFillHist() && sumV[0] + 4*sumV[1] == 8)
691 hDVVMass[11]->Fill(TMath::Min(diBoson->Mass(),199.999));
692 }
693 }
694 if (p->HasDaughter(MCParticle::kMuNu,kTRUE) && p->HasDaughter(-1*MCParticle::kMuNu,kTRUE)) {
695 isOld = kFALSE;
696 for(UInt_t nl = 0; nl < GenTempMG1->GetEntries(); nl++){
697 if(p->FindDaughter(MCParticle::kMuNu,kTRUE) == GenTempMG1->At(nl)) {
698 isOld = kTRUE;
699 break;
700 }
701 }
702 if(isOld == kFALSE){
703 GenTempMG1->Add(p->FindDaughter(MCParticle::kMuNu,kTRUE));
704 diBoson->AddDaughter(p->FindDaughter(MCParticle::kMuNu,kTRUE));
705 diBoson->AddDaughter(p->FindDaughter(-1*MCParticle::kMuNu,kTRUE));
706 if (GetFillHist() && sumV[0] + 4*sumV[1] == 5)
707 hDVVMass[12]->Fill(TMath::Min(diBoson->Mass(),199.999));
708 if (GetFillHist() && sumV[0] + 4*sumV[1] == 8)
709 hDVVMass[13]->Fill(TMath::Min(diBoson->Mass(),199.999));
710 }
711 }
712 if (p->HasDaughter(MCParticle::kElNu,kTRUE) && p->HasDaughter(-1*MCParticle::kElNu,kTRUE)) {
713 isOld = kFALSE;
714 for(UInt_t nl = 0; nl < GenTempMG1->GetEntries(); nl++){
715 if(p->FindDaughter(MCParticle::kElNu,kTRUE) == GenTempMG1->At(nl)) {
716 isOld = kTRUE;
717 break;
718 }
719 }
720 if(isOld == kFALSE){
721 GenTempMG1->Add(p->FindDaughter(MCParticle::kElNu,kTRUE));
722 diBoson->AddDaughter(p->FindDaughter(MCParticle::kElNu,kTRUE));
723 diBoson->AddDaughter(p->FindDaughter(-1*MCParticle::kElNu,kTRUE));
724 if (GetFillHist() && sumV[0] + 4*sumV[1] == 5)
725 hDVVMass[14]->Fill(TMath::Min(diBoson->Mass(),199.999));
726 if (GetFillHist() && sumV[0] + 4*sumV[1] == 8)
727 hDVVMass[15]->Fill(TMath::Min(diBoson->Mass(),199.999));
728 }
729 }
730 if (p->HasDaughter(MCParticle::kTauNu,kTRUE) &&
731 p->HasDaughter(-1*MCParticle::kTauNu,kTRUE)) {
732 isOld = kFALSE;
733 for(UInt_t nl = 0; nl < GenTempMG1->GetEntries(); nl++){
734 if(p->FindDaughter(MCParticle::kTauNu,kTRUE) == GenTempMG1->At(nl)) {
735 isOld = kTRUE;
736 break;
737 }
738 }
739 if(isOld == kFALSE){
740 GenTempMG1->Add(p->FindDaughter(MCParticle::kTauNu,kTRUE));
741 diBoson->AddDaughter(p->FindDaughter(MCParticle::kTauNu,kTRUE));
742 diBoson->AddDaughter(p->FindDaughter(-1*MCParticle::kTauNu,kTRUE));
743 if (GetFillHist() && sumV[0] + 4*sumV[1] == 5)
744 hDVVMass[16]->Fill(TMath::Min(diBoson->Mass(),199.999));
745 if (GetFillHist() && sumV[0] + 4*sumV[1] == 8)
746 hDVVMass[17]->Fill(TMath::Min(diBoson->Mass(),199.999));
747 }
748 }
749 if (diBoson && diBosonMass[0] <= 0) diBosonMass[0] = diBoson->Mass();
750 else if(diBoson && diBosonMass[1] <= 0) diBosonMass[1] = diBoson->Mass();
751 delete diBoson;
752 }
753 else if (p->Status() == 3 && (p->Is(MCParticle::kZ) || p->Is(MCParticle::kW))) {
754 if (diBosonMass[0] <= 0) diBosonMass[0] = p->Mass();
755 else if(diBosonMass[1] <= 0) diBosonMass[1] = p->Mass();
756 if (GetFillHist()) {
757 if (sumV[0] + 4*sumV[1] == 2 && p->Is(MCParticle::kW) &&
758 p->HasDaughter(MCParticle::kMu) &&
759 p->HasDaughter(MCParticle::kMuNu))
760 hDVVMass[0]->Fill(TMath::Min(p->Mass(),199.999));
761 else if(sumV[0] + 4*sumV[1] == 2 && p->Is(MCParticle::kW) &&
762 p->HasDaughter(MCParticle::kEl) &&
763 p->HasDaughter(MCParticle::kElNu))
764 hDVVMass[2]->Fill(TMath::Min(p->Mass(),199.999));
765 else if(sumV[0] + 4*sumV[1] == 2 && p->Is(MCParticle::kW) &&
766 p->HasDaughter(MCParticle::kTau) &&
767 p->HasDaughter(MCParticle::kTauNu))
768 hDVVMass[4]->Fill(TMath::Min(p->Mass(),199.999));
769 else if(sumV[0] + 4*sumV[1] == 5 && p->Is(MCParticle::kW) &&
770 p->HasDaughter(MCParticle::kMu) &&
771 p->HasDaughter(MCParticle::kMuNu))
772 hDVVMass[1]->Fill(TMath::Min(p->Mass(),199.999));
773 else if(sumV[0] + 4*sumV[1] == 5 && p->Is(MCParticle::kW) &&
774 p->HasDaughter(MCParticle::kEl) &&
775 p->HasDaughter(MCParticle::kElNu))
776 hDVVMass[3]->Fill(TMath::Min(p->Mass(),199.999));
777 else if(sumV[0] + 4*sumV[1] == 5 && p->Is(MCParticle::kW) &&
778 p->HasDaughter(MCParticle::kTau) &&
779 p->HasDaughter(MCParticle::kTauNu))
780 hDVVMass[5]->Fill(TMath::Min(p->Mass(),199.999));
781 else if(sumV[0] + 4*sumV[1] == 5 && p->Is(MCParticle::kZ) &&
782 p->HasDaughter(MCParticle::kMu,kTRUE) &&
783 p->HasDaughter(-1*MCParticle::kMu,kTRUE))
784 hDVVMass[6]->Fill(TMath::Min(p->Mass(),199.999));
785 else if(sumV[0] + 4*sumV[1] == 5 && p->Is(MCParticle::kZ) &&
786 p->HasDaughter(MCParticle::kEl,kTRUE) &&
787 p->HasDaughter(-1*MCParticle::kEl,kTRUE))
788 hDVVMass[8]->Fill(TMath::Min(p->Mass(),199.999));
789 else if(sumV[0] + 4*sumV[1] == 5 && p->Is(MCParticle::kZ) &&
790 p->HasDaughter(MCParticle::kTau,kTRUE) &&
791 p->HasDaughter(-1*MCParticle::kTau,kTRUE))
792 hDVVMass[10]->Fill(TMath::Min(p->Mass(),199.999));
793 else if(sumV[0] + 4*sumV[1] == 5 && p->Is(MCParticle::kZ) &&
794 p->HasDaughter(MCParticle::kMuNu,kTRUE) &&
795 p->HasDaughter(-1*MCParticle::kMuNu,kTRUE))
796 hDVVMass[12]->Fill(TMath::Min(p->Mass(),199.999));
797 else if(sumV[0] + 4*sumV[1] == 5 && p->Is(MCParticle::kZ) &&
798 p->HasDaughter(MCParticle::kElNu,kTRUE) &&
799 p->HasDaughter(-1*MCParticle::kElNu,kTRUE))
800 hDVVMass[14]->Fill(TMath::Min(p->Mass(),199.999));
801 else if(sumV[0] + 4*sumV[1] == 5 && p->Is(MCParticle::kZ) &&
802 p->HasDaughter(MCParticle::kTauNu,kTRUE) &&
803 p->HasDaughter(-1*MCParticle::kTauNu,kTRUE))
804 hDVVMass[16]->Fill(TMath::Min(p->Mass(),199.999));
805 else if(sumV[0] + 4*sumV[1] == 8 && p->Is(MCParticle::kZ) &&
806 p->HasDaughter(MCParticle::kMu,kTRUE) &&
807 p->HasDaughter(-1*MCParticle::kMu,kTRUE))
808 hDVVMass[7]->Fill(TMath::Min(p->Mass(),199.999));
809 else if(sumV[0] + 4*sumV[1] == 8 && p->Is(MCParticle::kZ) &&
810 p->HasDaughter(MCParticle::kEl,kTRUE) &&
811 p->HasDaughter(-1*MCParticle::kEl,kTRUE))
812 hDVVMass[9]->Fill(TMath::Min(p->Mass(),199.999));
813 else if(sumV[0] + 4*sumV[1] == 8 && p->Is(MCParticle::kZ) &&
814 p->HasDaughter(MCParticle::kTau,kTRUE) &&
815 p->HasDaughter(-1*MCParticle::kTau,kTRUE))
816 hDVVMass[11]->Fill(TMath::Min(p->Mass(),199.999));
817 else if(sumV[0] + 4*sumV[1] == 8 && p->Is(MCParticle::kZ) &&
818 p->HasDaughter(MCParticle::kMuNu,kTRUE) &&
819 p->HasDaughter(-1*MCParticle::kMuNu,kTRUE))
820 hDVVMass[13]->Fill(TMath::Min(p->Mass(),199.999));
821 else if(sumV[0] + 4*sumV[1] == 8 && p->Is(MCParticle::kZ) &&
822 p->HasDaughter(MCParticle::kElNu,kTRUE) &&
823 p->HasDaughter(-1*MCParticle::kElNu,kTRUE))
824 hDVVMass[15]->Fill(TMath::Min(p->Mass(),199.999));
825 else if(sumV[0] + 4*sumV[1] == 8 && p->Is(MCParticle::kZ) &&
826 p->HasDaughter(MCParticle::kTauNu,kTRUE) &&
827 p->HasDaughter(-1*MCParticle::kTauNu,kTRUE))
828 hDVVMass[17]->Fill(TMath::Min(p->Mass(),199.999));
829 }
830 }
831 } // end loop of particles
832 if (GetFillHist()) {
833 if(diBosonMass[0] > 70 && diBosonMass[0] < 110 &&
834 diBosonMass[1] > 70 && diBosonMass[1] < 110){
835 if(sumV[0] + 4*sumV[1] == 2){
836 if (sumVVFlavor[0] == 2) hDVVMass[18]->Fill(0.);
837 else if(sumVVFlavor[1] == 2) hDVVMass[18]->Fill(1.);
838 else if(sumVVFlavor[2] == 2) hDVVMass[18]->Fill(2.);
839 else if(sumVVFlavor[0] == 1 && sumVVFlavor[1] == 1) hDVVMass[18]->Fill(3.);
840 else if(sumVVFlavor[0] == 1 && sumVVFlavor[2] == 1) hDVVMass[18]->Fill(4.);
841 else if(sumVVFlavor[1] == 1 && sumVVFlavor[2] == 1) hDVVMass[18]->Fill(5.);
842 else hDVVMass[18]->Fill(6.);
843 }
844 if(sumV[0] + 4*sumV[1] == 5){
845 if (sumVVFlavor[3] == 1 && sumVVFlavor[0] == 1) hDVVMass[19]->Fill(0.);
846 else if(sumVVFlavor[3] == 1 && sumVVFlavor[1] == 1) hDVVMass[19]->Fill(1.);
847 else if(sumVVFlavor[3] == 1 && sumVVFlavor[2] == 1) hDVVMass[19]->Fill(2.);
848 else if(sumVVFlavor[4] == 1 && sumVVFlavor[0] == 1) hDVVMass[19]->Fill(3.);
849 else if(sumVVFlavor[4] == 1 && sumVVFlavor[1] == 1) hDVVMass[19]->Fill(4.);
850 else if(sumVVFlavor[4] == 1 && sumVVFlavor[2] == 1) hDVVMass[19]->Fill(5.);
851 else if(sumVVFlavor[5] == 1 && sumVVFlavor[0] == 1) hDVVMass[19]->Fill(6.);
852 else if(sumVVFlavor[5] == 1 && sumVVFlavor[1] == 1) hDVVMass[19]->Fill(7.);
853 else if(sumVVFlavor[5] == 1 && sumVVFlavor[2] == 1) hDVVMass[19]->Fill(8.);
854 else { hDVVMass[19]->Fill(9.);
855 /*for(int i=0; i<9; i++) cout << sumVVFlavor[i] << " " ;cout << endl;*/}
856 }
857 if(sumV[0] + 4*sumV[1] == 8 &&
858 sumVVFlavor[3] + sumVVFlavor[4] + sumVVFlavor[5] == 2){
859 if (sumVVFlavor[3] == 2) hDVVMass[20]->Fill(0.);
860 else if(sumVVFlavor[4] == 2) hDVVMass[20]->Fill(1.);
861 else if(sumVVFlavor[5] == 2) hDVVMass[20]->Fill(2.);
862 else if(sumVVFlavor[3] == 1 && sumVVFlavor[4] == 1) hDVVMass[20]->Fill(3.);
863 else if(sumVVFlavor[3] == 1 && sumVVFlavor[5] == 1) hDVVMass[20]->Fill(4.);
864 else if(sumVVFlavor[4] == 1 && sumVVFlavor[5] == 1) hDVVMass[20]->Fill(5.);
865 else hDVVMass[20]->Fill(6.);
866 }
867 else if(sumV[0] + 4*sumV[1] == 8){
868 if (sumVVFlavor[6] == 2) hDVVMass[21]->Fill(0.);
869 else if(sumVVFlavor[7] == 2) hDVVMass[21]->Fill(1.);
870 else if(sumVVFlavor[8] == 2) hDVVMass[21]->Fill(2.);
871 else if(sumVVFlavor[3] == 1 && sumVVFlavor[6] == 1) hDVVMass[21]->Fill(3.);
872 else if(sumVVFlavor[3] == 1 && sumVVFlavor[7] == 1) hDVVMass[21]->Fill(4.);
873 else if(sumVVFlavor[3] == 1 && sumVVFlavor[8] == 1) hDVVMass[21]->Fill(5.);
874 else if(sumVVFlavor[4] == 1 && sumVVFlavor[6] == 1) hDVVMass[21]->Fill(6.);
875 else if(sumVVFlavor[4] == 1 && sumVVFlavor[7] == 1) hDVVMass[21]->Fill(7.);
876 else if(sumVVFlavor[4] == 1 && sumVVFlavor[8] == 1) hDVVMass[21]->Fill(8.);
877 else if(sumVVFlavor[5] == 1 && sumVVFlavor[6] == 1) hDVVMass[21]->Fill(9.);
878 else if(sumVVFlavor[5] == 1 && sumVVFlavor[7] == 1) hDVVMass[21]->Fill(10.);
879 else if(sumVVFlavor[5] == 1 && sumVVFlavor[8] == 1) hDVVMass[21]->Fill(11.);
880 else if(sumVVFlavor[6] == 1 && sumVVFlavor[7] == 1) hDVVMass[21]->Fill(12.);
881 else if(sumVVFlavor[6] == 1 && sumVVFlavor[8] == 1) hDVVMass[21]->Fill(13.);
882 else if(sumVVFlavor[7] == 1 && sumVVFlavor[8] == 1) hDVVMass[21]->Fill(14.);
883 else hDVVMass[21]->Fill(15.);
884 }
885 } // 60<mV1/2<120
886 if(sumV[0] + 4*sumV[1] == 2)
887 hDVVMass[22]->Fill(TMath::Min(TMath::Min(diBosonMass[0],diBosonMass[1]),199.999));
888 if(sumV[0] + 4*sumV[1] == 2)
889 hDVVMass[23]->Fill(TMath::Min(TMath::Max(diBosonMass[0],diBosonMass[1]),199.999));
890 if(sumV[0] + 4*sumV[1] == 5)
891 hDVVMass[24]->Fill(TMath::Min(TMath::Min(diBosonMass[0],diBosonMass[1]),199.999));
892 if(sumV[0] + 4*sumV[1] == 5)
893 hDVVMass[25]->Fill(TMath::Min(TMath::Max(diBosonMass[0],diBosonMass[1]),199.999));
894 if(sumV[0] + 4*sumV[1] == 8)
895 hDVVMass[26]->Fill(TMath::Min(TMath::Min(diBosonMass[0],diBosonMass[1]),199.999));
896 if(sumV[0] + 4*sumV[1] == 8)
897 hDVVMass[27]->Fill(TMath::Min(TMath::Max(diBosonMass[0],diBosonMass[1]),199.999));
898 }
899 delete GenTempMG1;
900 } // WW, WZ or ZZ
901 // --------------------------------
902 // End special study about VVjets
903 // --------------------------------
904
905 Met *theMET = new Met(totalMET[0], totalMET[1]);
906 theMET->SetElongitudinal(totalMET[2]);
907 GenMet->AddOwned(theMET);
908
909 // sort according to pt
910 GenLeptons->Sort();
911 GenAllLeptons->Sort();
912 GenTaus->Sort();
913 GenNeutrinos->Sort();
914 GenQuarks->Sort();
915 GenqqHs->Sort();
916 GenBosons->Sort();
917 GenPhotons->Sort();
918 GenRadPhotons->Sort();
919 GenISRPhotons->Sort();
920 } // Only for Monte Carlo
921
922 // add objects to this event for other modules to use
923 AddObjThisEvt(GenMet);
924 AddObjThisEvt(GenLeptons);
925 AddObjThisEvt(GenAllLeptons);
926 AddObjThisEvt(GenTaus);
927 AddObjThisEvt(GenNeutrinos);
928 AddObjThisEvt(GenQuarks);
929 AddObjThisEvt(GenqqHs);
930 AddObjThisEvt(GenBosons);
931 AddObjThisEvt(GenPhotons);
932 AddObjThisEvt(GenRadPhotons);
933 AddObjThisEvt(GenISRPhotons);
934
935 // --------------------------------
936 // Copy these Collections into the Arrays for Publication for Output Module
937 // --------------------------------
938 fGenLeptons->Delete();
939 fGenAllLeptons->Delete();
940 fGenTaus->Delete();
941 fGenNeutrinos->Delete();
942 fGenQuarks->Delete();
943 fGenqqHs->Delete();
944 fGenBosons->Delete();
945 fGenPhotons->Delete();
946 fGenRadPhotons->Delete();
947 fGenISRPhotons->Delete();
948
949 for (UInt_t i=0; i < GenLeptons->GetEntries(); ++i) {
950 mithep::MCParticle *genParticle = fGenLeptons->Allocate();
951 new (genParticle) mithep::MCParticle(GenLeptons->At(i)->Px(),GenLeptons->At(i)->Py(),
952 GenLeptons->At(i)->Pz(),GenLeptons->At(i)->E(),
953 GenLeptons->At(i)->PdgId(),GenLeptons->At(i)->Status());
954 }
955 for (UInt_t i=0; i < GenAllLeptons->GetEntries(); ++i) {
956 mithep::MCParticle *genParticle = fGenAllLeptons->Allocate();
957 new (genParticle) mithep::MCParticle(GenAllLeptons->At(i)->Px(),GenAllLeptons->At(i)->Py(),
958 GenAllLeptons->At(i)->Pz(),GenAllLeptons->At(i)->E(),
959 GenAllLeptons->At(i)->PdgId(),GenAllLeptons->At(i)->Status());
960 }
961 for (UInt_t i=0; i < GenTaus->GetEntries(); ++i) {
962 mithep::MCParticle *genParticle = fGenTaus->Allocate();
963 new (genParticle) mithep::MCParticle(GenTaus->At(i)->Px(),GenTaus->At(i)->Py(),
964 GenTaus->At(i)->Pz(),GenTaus->At(i)->E(),
965 GenTaus->At(i)->PdgId(),GenTaus->At(i)->Status());
966 }
967 for (UInt_t i=0; i < GenNeutrinos->GetEntries(); ++i) {
968 mithep::MCParticle *genParticle = fGenNeutrinos->Allocate();
969 new (genParticle) mithep::MCParticle(GenNeutrinos->At(i)->Px(),GenNeutrinos->At(i)->Py(),
970 GenNeutrinos->At(i)->Pz(),GenNeutrinos->At(i)->E(),
971 GenNeutrinos->At(i)->PdgId(),GenNeutrinos->At(i)->Status());
972 }
973 for (UInt_t i=0; i < GenQuarks->GetEntries(); ++i) {
974 mithep::MCParticle *genParticle = fGenQuarks->Allocate();
975 new (genParticle) mithep::MCParticle(GenQuarks->At(i)->Px(),GenQuarks->At(i)->Py(),
976 GenQuarks->At(i)->Pz(),GenQuarks->At(i)->E(),
977 GenQuarks->At(i)->PdgId(),GenQuarks->At(i)->Status());
978 }
979 for (UInt_t i=0; i < GenqqHs->GetEntries(); ++i) {
980 mithep::MCParticle *genParticle = fGenqqHs->Allocate();
981 new (genParticle) mithep::MCParticle(GenqqHs->At(i)->Px(),GenqqHs->At(i)->Py(),
982 GenqqHs->At(i)->Pz(),GenqqHs->At(i)->E(),
983 GenqqHs->At(i)->PdgId(),GenqqHs->At(i)->Status());
984 }
985
986 if (fCopyArrays) {
987 // --------------------------------
988 // Copy these Collections into the Arrays for Publication for Output Module
989 // --------------------------------
990 fGenLeptons->Delete();
991 fGenAllLeptons->Delete();
992 fGenTaus->Delete();
993 fGenNeutrinos->Delete();
994 fGenQuarks->Delete();
995 fGenqqHs->Delete();
996 fGenBosons->Delete();
997 fGenPhotons->Delete();
998 fGenRadPhotons->Delete();
999 fGenISRPhotons->Delete();
1000
1001 for (UInt_t i=0; i < GenLeptons->GetEntries(); ++i) {
1002 mithep::MCParticle *genParticle = fGenLeptons->Allocate();
1003 new (genParticle) mithep::MCParticle(GenLeptons->At(i)->Px(),GenLeptons->At(i)->Py(),
1004 GenLeptons->At(i)->Pz(),GenLeptons->At(i)->E(),
1005 GenLeptons->At(i)->PdgId(),GenLeptons->At(i)->Status());
1006 }
1007 for (UInt_t i=0; i < GenAllLeptons->GetEntries(); ++i) {
1008 mithep::MCParticle *genParticle = fGenAllLeptons->Allocate();
1009 new (genParticle) mithep::MCParticle(GenAllLeptons->At(i)->Px(),GenAllLeptons->At(i)->Py(),
1010 GenAllLeptons->At(i)->Pz(),GenAllLeptons->At(i)->E(),
1011 GenAllLeptons->At(i)->PdgId(),
1012 GenAllLeptons->At(i)->Status());
1013 }
1014 for (UInt_t i=0; i < GenTaus->GetEntries(); ++i) {
1015 mithep::MCParticle *genParticle = fGenTaus->Allocate();
1016 new (genParticle) mithep::MCParticle(GenTaus->At(i)->Px(),GenTaus->At(i)->Py(),
1017 GenTaus->At(i)->Pz(),GenTaus->At(i)->E(),
1018 GenTaus->At(i)->PdgId(),
1019 GenTaus->At(i)->Status());
1020 }
1021 for (UInt_t i=0; i < GenNeutrinos->GetEntries(); ++i) {
1022 mithep::MCParticle *genParticle = fGenNeutrinos->Allocate();
1023 new (genParticle) mithep::MCParticle(GenNeutrinos->At(i)->Px(),GenNeutrinos->At(i)->Py(),
1024 GenNeutrinos->At(i)->Pz(),GenNeutrinos->At(i)->E(),
1025 GenNeutrinos->At(i)->PdgId(),
1026 GenNeutrinos->At(i)->Status());
1027 }
1028 for (UInt_t i=0; i < GenQuarks->GetEntries(); ++i) {
1029 mithep::MCParticle *genParticle = fGenQuarks->Allocate();
1030 new (genParticle) mithep::MCParticle(GenQuarks->At(i)->Px(),GenQuarks->At(i)->Py(),
1031 GenQuarks->At(i)->Pz(),GenQuarks->At(i)->E(),
1032 GenQuarks->At(i)->PdgId(),
1033 GenQuarks->At(i)->Status());
1034 }
1035 for (UInt_t i=0; i < GenqqHs->GetEntries(); ++i) {
1036 mithep::MCParticle *genParticle = fGenqqHs->Allocate();
1037 new (genParticle) mithep::MCParticle(GenqqHs->At(i)->Px(),GenqqHs->At(i)->Py(),
1038 GenqqHs->At(i)->Pz(),GenqqHs->At(i)->E(),
1039 GenqqHs->At(i)->PdgId(),
1040 GenqqHs->At(i)->Status());
1041 }
1042 for (UInt_t i=0; i < GenBosons->GetEntries(); ++i) {
1043 mithep::MCParticle *genParticle = fGenBosons->Allocate();
1044 new (genParticle) mithep::MCParticle(GenBosons->At(i)->Px(),GenBosons->At(i)->Py(),
1045 GenBosons->At(i)->Pz(),GenBosons->At(i)->E(),
1046 GenBosons->At(i)->PdgId(),
1047 GenBosons->At(i)->Status());
1048 }
1049 for (UInt_t i=0; i < GenPhotons->GetEntries(); ++i) {
1050 mithep::MCParticle *genParticle = fGenPhotons->Allocate();
1051 new (genParticle) mithep::MCParticle(GenPhotons->At(i)->Px(),GenPhotons->At(i)->Py(),
1052 GenPhotons->At(i)->Pz(),GenPhotons->At(i)->E(),
1053 GenPhotons->At(i)->PdgId(),
1054 GenPhotons->At(i)->Status());
1055 }
1056 for (UInt_t i=0; i < GenRadPhotons->GetEntries(); ++i) {
1057 mithep::MCParticle *genParticle = fGenRadPhotons->Allocate();
1058 new (genParticle) mithep::MCParticle(GenRadPhotons->At(i)->Px(),GenRadPhotons->At(i)->Py(),
1059 GenRadPhotons->At(i)->Pz(),GenRadPhotons->At(i)->E(),
1060 GenRadPhotons->At(i)->PdgId(),
1061 GenRadPhotons->At(i)->Status());
1062 }
1063 for (UInt_t i=0; i < GenISRPhotons->GetEntries(); ++i) {
1064 mithep::MCParticle *genParticle = fGenISRPhotons->Allocate();
1065 new (genParticle) mithep::MCParticle(GenISRPhotons->At(i)->Px(),GenISRPhotons->At(i)->Py(),
1066 GenISRPhotons->At(i)->Pz(),GenISRPhotons->At(i)->E(),
1067 GenISRPhotons->At(i)->PdgId(),
1068 GenISRPhotons->At(i)->Status());
1069 }
1070 }
1071
1072 for (UInt_t i=0; i < GenBosons->GetEntries(); ++i) {
1073 mithep::MCParticle *genParticle = fGenBosons->Allocate();
1074 new (genParticle) mithep::MCParticle(GenBosons->At(i)->Px(),GenBosons->At(i)->Py(),
1075 GenBosons->At(i)->Pz(),GenBosons->At(i)->E(),
1076 GenBosons->At(i)->PdgId(),GenBosons->At(i)->Status());
1077 }
1078 for (UInt_t i=0; i < GenPhotons->GetEntries(); ++i) {
1079 mithep::MCParticle *genParticle = fGenPhotons->Allocate();
1080 new (genParticle) mithep::MCParticle(GenPhotons->At(i)->Px(),GenPhotons->At(i)->Py(),
1081 GenPhotons->At(i)->Pz(),GenPhotons->At(i)->E(),
1082 GenPhotons->At(i)->PdgId(),GenPhotons->At(i)->Status());
1083 }
1084 for (UInt_t i=0; i < GenRadPhotons->GetEntries(); ++i) {
1085 mithep::MCParticle *genParticle = fGenRadPhotons->Allocate();
1086 new (genParticle) mithep::MCParticle(GenRadPhotons->At(i)->Px(),GenRadPhotons->At(i)->Py(),
1087 GenRadPhotons->At(i)->Pz(),GenRadPhotons->At(i)->E(),
1088 GenRadPhotons->At(i)->PdgId(),GenRadPhotons->At(i)->Status());
1089 }
1090 for (UInt_t i=0; i < GenISRPhotons->GetEntries(); ++i) {
1091 mithep::MCParticle *genParticle = fGenISRPhotons->Allocate();
1092 new (genParticle) mithep::MCParticle(GenISRPhotons->At(i)->Px(),GenISRPhotons->At(i)->Py(),
1093 GenISRPhotons->At(i)->Pz(),GenISRPhotons->At(i)->E(),
1094 GenISRPhotons->At(i)->PdgId(),GenISRPhotons->At(i)->Status());
1095 }
1096
1097 //bool theZCut = (GenLeptons->GetEntries() > 0 || GenTaus->GetEntries() > 0) && GenQuarks->GetEntries() == 0;
1098 //if(!theZCut){
1099 // SkipEvent();
1100 // return;
1101 //}
1102
1103 // fill histograms if requested
1104 if (GetFillHist()) {
1105 // MET
1106 hDGenMet[0]->Fill(GenMet->At(0)->Pt());
1107 hDGenMet[1]->Fill(GenMet->At(0)->Px());
1108 hDGenMet[2]->Fill(GenMet->At(0)->Py());
1109 hDGenMet[3]->Fill(GenMet->At(0)->Elongitudinal());
1110
1111 // pt min for leptons from W/Z
1112 hDGenPtMin->Fill(TMath::Min(ptMin, 199.999));
1113 //if(ptMin < 10){
1114 // SkipEvent();
1115 // return;
1116 //}
1117
1118 // leptons
1119 hDGenLeptons[0]->Fill(GenLeptons->GetEntries());
1120 for(UInt_t i=0; i<GenLeptons->GetEntries(); i++) {
1121 hDGenLeptons[1]->Fill(GenLeptons->At(i)->Pt());
1122 hDGenLeptons[2]->Fill(TMath::Abs(GenLeptons->At(i)->Eta()));
1123 hDGenLeptons[3]->Fill(GenLeptons->At(i)->PhiDeg());
1124 for(UInt_t j=i+1; j<GenLeptons->GetEntries(); j++) {
1125 CompositeParticle *dilepton = new CompositeParticle();
1126 dilepton->AddDaughter(GenLeptons->At(i));
1127 dilepton->AddDaughter(GenLeptons->At(j));
1128 hDGenLeptons[4]->Fill(dilepton->Mass());
1129 delete dilepton;
1130 }
1131 }
1132 // looking at events with two leptons
1133 if (GenLeptons->GetEntries() == 2) {
1134 hDGenLeptons[5]->Fill(TMath::Min(TMath::Max(TMath::Abs(GenLeptons->At(0)->Eta()),
1135 TMath::Abs(GenLeptons->At(1)->Eta())),
1136 4.999));
1137 hDGenLeptons[6]->Fill(TMath::Min(TMath::Min(TMath::Abs(GenLeptons->At(0)->Eta()),
1138 TMath::Abs(GenLeptons->At(1)->Eta())),
1139 4.999));
1140 if (TMath::Abs(GenLeptons->At(0)->Eta()) < 2.5 &&
1141 TMath::Abs(GenLeptons->At(1)->Eta()) < 2.5) {
1142 hDGenLeptons[7]->Fill(TMath::Min(GenLeptons->At(0)->Pt(),199.999));
1143 if (GenLeptons->At(0)->Pt() > 20.0) {
1144 hDGenLeptons[8]->Fill(TMath::Min(GenLeptons->At(1)->Pt(),199.999));
1145 if (GenLeptons->At(1)->Pt() > 10.0) {
1146 CompositeParticle *dilepton = new CompositeParticle();
1147 dilepton->AddDaughter(GenLeptons->At(0));
1148 dilepton->AddDaughter(GenLeptons->At(1));
1149 hDGenLeptons[9]->Fill(TMath::Min(dilepton->Mass(),999.999));
1150 if(dilepton->Mass() > 12.0){
1151 hDGenLeptons[10]->Fill(MathUtils::DeltaPhi(GenLeptons->At(0)->Phi(),
1152 GenLeptons->At(1)->Phi())
1153 * 180./ TMath::Pi());
1154 hDGenLeptons[11]->Fill(MathUtils::DeltaR(*GenLeptons->At(0),
1155 *GenLeptons->At(1)));
1156 }
1157 delete dilepton;
1158 }
1159 }
1160 }
1161 }
1162 // looking at events with three leptons
1163 if (GenLeptons->GetEntries() == 3) {
1164 if (TMath::Abs(GenLeptons->At(0)->Eta()) < 2.5 &&
1165 TMath::Abs(GenLeptons->At(1)->Eta()) < 2.5 &&
1166 TMath::Abs(GenLeptons->At(2)->Eta()) < 2.5) {
1167 hDGenLeptons[12]->Fill(TMath::Min(GenLeptons->At(0)->Pt(),199.999));
1168 if (GenLeptons->At(0)->Pt() > 20.0) {
1169 hDGenLeptons[13]->Fill(TMath::Min(GenLeptons->At(1)->Pt(),199.999));
1170 hDGenLeptons[14]->Fill(TMath::Min(GenLeptons->At(2)->Pt(),199.999));
1171 if (GenLeptons->At(1)->Pt() > 10.0 && GenLeptons->At(2)->Pt() > 10.0) {
1172 CompositeParticle *dilepton01 = new CompositeParticle();
1173 dilepton01->AddDaughter(GenLeptons->At(0));
1174 dilepton01->AddDaughter(GenLeptons->At(1));
1175 CompositeParticle *dilepton02 = new CompositeParticle();
1176 dilepton02->AddDaughter(GenLeptons->At(0));
1177 dilepton02->AddDaughter(GenLeptons->At(2));
1178 CompositeParticle *dilepton12 = new CompositeParticle();
1179 dilepton12->AddDaughter(GenLeptons->At(1));
1180 dilepton12->AddDaughter(GenLeptons->At(2));
1181 hDGenLeptons[15]->Fill(TMath::Min(dilepton01->Mass(),999.999));
1182 hDGenLeptons[15]->Fill(TMath::Min(dilepton02->Mass(),999.999));
1183 hDGenLeptons[15]->Fill(TMath::Min(dilepton12->Mass(),999.999));
1184 CompositeParticle *trilepton = new CompositeParticle();
1185 trilepton->AddDaughter(GenLeptons->At(0));
1186 trilepton->AddDaughter(GenLeptons->At(1));
1187 trilepton->AddDaughter(GenLeptons->At(2));
1188 hDGenLeptons[16]->Fill(TMath::Min(trilepton->Mass(),999.999));
1189 Double_t deltaR[3] = {MathUtils::DeltaR(*GenLeptons->At(0),
1190 *GenLeptons->At(1)),
1191 MathUtils::DeltaR(*GenLeptons->At(0),
1192 *GenLeptons->At(2)),
1193 MathUtils::DeltaR(*GenLeptons->At(1),
1194 *GenLeptons->At(2))};
1195 Double_t deltaRMin = deltaR[0];
1196 for(Int_t i=1; i<3; i++)
1197 if(deltaRMin > deltaR[i])
1198 deltaRMin = deltaR[i];
1199 hDGenLeptons[17]->Fill(deltaRMin);
1200
1201 delete dilepton01;
1202 delete dilepton02;
1203 delete dilepton12;
1204 delete trilepton;
1205 }
1206 }
1207 }
1208 }
1209 // looking at events with four leptons
1210 if (GenLeptons->GetEntries() == 4) {
1211 if (TMath::Abs(GenLeptons->At(0)->Eta()) < 2.5 &&
1212 TMath::Abs(GenLeptons->At(1)->Eta()) < 2.5 &&
1213 TMath::Abs(GenLeptons->At(2)->Eta()) < 2.5 &&
1214 TMath::Abs(GenLeptons->At(3)->Eta()) < 2.5) {
1215 hDGenLeptons[18]->Fill(TMath::Min(GenLeptons->At(0)->Pt(),199.999));
1216 if (GenLeptons->At(0)->Pt() > 20.0) {
1217 hDGenLeptons[19]->Fill(TMath::Min(GenLeptons->At(1)->Pt(),199.999));
1218 hDGenLeptons[20]->Fill(TMath::Min(GenLeptons->At(2)->Pt(),199.999));
1219 hDGenLeptons[21]->Fill(TMath::Min(GenLeptons->At(3)->Pt(),199.999));
1220 if (GenLeptons->At(1)->Pt() > 10.0 && GenLeptons->At(2)->Pt() > 10.0 &&
1221 GenLeptons->At(3)->Pt() > 10.0) {
1222 CompositeParticle *dilepton01 = new CompositeParticle();
1223 dilepton01->AddDaughter(GenLeptons->At(0));
1224 dilepton01->AddDaughter(GenLeptons->At(1));
1225 CompositeParticle *dilepton02 = new CompositeParticle();
1226 dilepton02->AddDaughter(GenLeptons->At(0));
1227 dilepton02->AddDaughter(GenLeptons->At(2));
1228 CompositeParticle *dilepton03 = new CompositeParticle();
1229 dilepton03->AddDaughter(GenLeptons->At(0));
1230 dilepton03->AddDaughter(GenLeptons->At(3));
1231 CompositeParticle *dilepton12 = new CompositeParticle();
1232 dilepton12->AddDaughter(GenLeptons->At(1));
1233 dilepton12->AddDaughter(GenLeptons->At(2));
1234 CompositeParticle *dilepton13 = new CompositeParticle();
1235 dilepton13->AddDaughter(GenLeptons->At(1));
1236 dilepton13->AddDaughter(GenLeptons->At(3));
1237 CompositeParticle *dilepton23 = new CompositeParticle();
1238 dilepton23->AddDaughter(GenLeptons->At(2));
1239 dilepton23->AddDaughter(GenLeptons->At(3));
1240 hDGenLeptons[22]->Fill(TMath::Min(dilepton01->Mass(),999.999));
1241 hDGenLeptons[22]->Fill(TMath::Min(dilepton02->Mass(),999.999));
1242 hDGenLeptons[22]->Fill(TMath::Min(dilepton03->Mass(),999.999));
1243 hDGenLeptons[22]->Fill(TMath::Min(dilepton12->Mass(),999.999));
1244 hDGenLeptons[22]->Fill(TMath::Min(dilepton13->Mass(),999.999));
1245 hDGenLeptons[22]->Fill(TMath::Min(dilepton23->Mass(),999.999));
1246 CompositeParticle *fourlepton = new CompositeParticle();
1247 fourlepton->AddDaughter(GenLeptons->At(0));
1248 fourlepton->AddDaughter(GenLeptons->At(1));
1249 fourlepton->AddDaughter(GenLeptons->At(2));
1250 fourlepton->AddDaughter(GenLeptons->At(3));
1251 hDGenLeptons[23]->Fill(TMath::Min(fourlepton->Mass(),999.999));
1252 Double_t deltaR[6] = {MathUtils::DeltaR(*GenLeptons->At(0),
1253 *GenLeptons->At(1)),
1254 MathUtils::DeltaR(*GenLeptons->At(0),
1255 *GenLeptons->At(2)),
1256 MathUtils::DeltaR(*GenLeptons->At(0),
1257 *GenLeptons->At(3)),
1258 MathUtils::DeltaR(*GenLeptons->At(1),
1259 *GenLeptons->At(2)),
1260 MathUtils::DeltaR(*GenLeptons->At(1),
1261 *GenLeptons->At(3)),
1262 MathUtils::DeltaR(*GenLeptons->At(2),
1263 *GenLeptons->At(3))};
1264 Double_t deltaRMin = deltaR[0];
1265 for(Int_t i=1; i<6; i++)
1266 if(deltaRMin > deltaR[i])
1267 deltaRMin = deltaR[i];
1268 hDGenLeptons[24]->Fill(deltaRMin);
1269
1270 delete dilepton01;
1271 delete dilepton02;
1272 delete dilepton03;
1273 delete dilepton12;
1274 delete dilepton13;
1275 delete dilepton23;
1276 delete fourlepton;
1277 }
1278 }
1279 }
1280 }
1281
1282 // all leptons
1283 hDGenAllLeptons[0]->Fill(GenAllLeptons->GetEntries());
1284 for(UInt_t i=0; i<GenAllLeptons->GetEntries(); i++) {
1285 hDGenAllLeptons[1]->Fill(GenAllLeptons->At(i)->Pt());
1286 hDGenAllLeptons[2]->Fill(GenAllLeptons->At(i)->AbsEta());
1287 hDGenAllLeptons[3]->Fill(GenAllLeptons->At(i)->PhiDeg());
1288 }
1289 if(GenAllLeptons->GetEntries() >= 2) hDGenAllLeptons[4]->Fill(GenAllLeptons->At(1)->Pt());
1290
1291 // taus
1292 hDGenTaus[0]->Fill(GenTaus->GetEntries());
1293 for(UInt_t i=0; i<GenTaus->GetEntries(); i++) {
1294 hDGenTaus[1]->Fill(GenTaus->At(i)->Pt());
1295 hDGenTaus[2]->Fill(GenTaus->At(i)->AbsEta());
1296 hDGenTaus[3]->Fill(GenTaus->At(i)->PhiDeg());
1297 }
1298
1299 // neutrinos
1300 hDGenNeutrinos[0]->Fill(GenNeutrinos->GetEntries());
1301 CompositeParticle *neutrinoTotal = new CompositeParticle();
1302 for(UInt_t i=0; i<GenNeutrinos->GetEntries(); i++) {
1303 if (GenNeutrinos->At(i)->HasMother())
1304 neutrinoTotal->AddDaughter(GenNeutrinos->At(i));
1305 }
1306 if (GenNeutrinos->GetEntries() > 0) {
1307 hDGenNeutrinos[1]->Fill(neutrinoTotal->Pt());
1308 hDGenNeutrinos[2]->Fill(neutrinoTotal->AbsEta());
1309 hDGenNeutrinos[3]->Fill(neutrinoTotal->PhiDeg());
1310 }
1311 delete neutrinoTotal;
1312
1313 // quarks
1314 hDGenQuarks[0]->Fill(GenQuarks->GetEntries());
1315 for(UInt_t i=0; i<GenQuarks->GetEntries(); i++) {
1316 for(UInt_t j=i+1; j<GenQuarks->GetEntries(); j++) {
1317 CompositeParticle *dijet = new CompositeParticle();
1318 dijet->AddDaughter(GenQuarks->At(i));
1319 dijet->AddDaughter(GenQuarks->At(j));
1320 hDGenQuarks[1]->Fill(dijet->Pt());
1321 hDGenQuarks[2]->Fill(dijet->Mass());
1322 if (TMath::Abs(GenQuarks->At(i)->Eta()) < 2.5 &&
1323 TMath::Abs(GenQuarks->At(j)->Eta()) < 2.5) {
1324 hDGenQuarks[3]->Fill(dijet->Pt());
1325 hDGenQuarks[4]->Fill(dijet->Mass());
1326 }
1327 delete dijet;
1328 }
1329 // b quark info
1330 if (GenQuarks->At(i)->AbsPdgId() == 5) {
1331 hDGenQuarks[5]->Fill(GenQuarks->At(i)->Pt());
1332 hDGenQuarks[6]->Fill(GenQuarks->At(i)->Eta());
1333 hDGenQuarks[7]->Fill(GenQuarks->At(i)->Phi());
1334 if (GenLeptons->GetEntries() >= 2 &&
1335 GenLeptons->At(0)->Pt() > 20 &&
1336 GenLeptons->At(1)->Pt() > 15) {
1337 if (TMath::Abs(GenLeptons->At(0)->Eta()) < 2.5 &&
1338 TMath::Abs(GenLeptons->At(1)->Eta()) < 2.5) {
1339 hDGenQuarks[8]->Fill(GenQuarks->At(i)->Pt());
1340 hDGenQuarks[9]->Fill(GenQuarks->At(i)->Eta());
1341 hDGenQuarks[10]->Fill(GenQuarks->At(i)->Phi());
1342 }
1343 }
1344 }
1345 // t quark info
1346 else if (GenQuarks->At(i)->AbsPdgId() == 6) {
1347 hDGenQuarks[11]->Fill(GenQuarks->At(i)->Pt());
1348 hDGenQuarks[12]->Fill(GenQuarks->At(i)->Eta());
1349 hDGenQuarks[13]->Fill(GenQuarks->At(i)->Phi());
1350 }
1351 // light quark info
1352 else {
1353 hDGenQuarks[14]->Fill(GenQuarks->At(i)->Pt());
1354 hDGenQuarks[15]->Fill(GenQuarks->At(i)->Eta());
1355 hDGenQuarks[16]->Fill(GenQuarks->At(i)->Phi());
1356 }
1357 }
1358
1359 // wbf
1360 if (GenqqHs->GetEntries() == 2) {
1361 hDGenWBF[0]->Fill(MathUtils::DeltaPhi(GenqqHs->At(0)->Phi(),
1362 GenqqHs->At(1)->Phi()) * 180./ TMath::Pi());
1363 hDGenWBF[1]->Fill(TMath::Abs(GenqqHs->At(0)->Eta()-GenqqHs->At(1)->Eta()));
1364 hDGenWBF[2]->Fill(TMath::Max(GenqqHs->At(0)->Pt(),GenqqHs->At(1)->Pt()));
1365 hDGenWBF[3]->Fill(TMath::Min(GenqqHs->At(0)->Pt(),GenqqHs->At(1)->Pt()));
1366 CompositeParticle *diqq = new CompositeParticle();
1367 diqq->AddDaughter(GenqqHs->At(0));
1368 diqq->AddDaughter(GenqqHs->At(1));
1369 hDGenWBF[4]->Fill(diqq->Mass());
1370 delete diqq;
1371 }
1372
1373 // bosons
1374 hDGenBosons[0]->Fill(GenBosons->GetEntries());
1375 for(UInt_t i=0; i<GenBosons->GetEntries(); i++) {
1376 hDGenBosons[1]->Fill(GenBosons->At(i)->Pt());
1377 hDGenBosons[2]->Fill(GenBosons->At(i)->Eta());
1378 hDGenBosons[3]->Fill(TMath::Min(GenBosons->At(i)->Mass(),1999.999));
1379 hDGenBosons[4]->Fill(TMath::Min(GenBosons->At(i)->Mass(),199.999));
1380 if(GenBosons->At(i)->Is(MCParticle::kW))
1381 hDGenBosons[5]->Fill(TMath::Min(GenBosons->At(i)->Mass(),199.999));
1382 if(GenBosons->At(i)->Is(MCParticle::kZ))
1383 hDGenBosons[6]->Fill(TMath::Min(GenBosons->At(i)->Mass(),199.999));
1384 }
1385 //if(sumV[0] >= 4) printf("More than 3 W bosons (%d)\n",sumV[0]);
1386 //if(sumV[1] >= 4) printf("More than 3 Z bosons (%d)\n",sumV[1]);
1387 hDGenBosons[7]->Fill(TMath::Min((double)(sumV[0] + 4*sumV[1]),12.4999));
1388
1389 // photons
1390 hDGenPhotons[0]->Fill(GenPhotons->GetEntries());
1391 for(UInt_t i=0; i<GenPhotons->GetEntries(); i++) {
1392 hDGenPhotons[1]->Fill(GenPhotons->At(i)->Pt());
1393 hDGenPhotons[2]->Fill(GenPhotons->At(i)->Eta());
1394 }
1395
1396 // Rad photons
1397 hDGenRadPhotons[0]->Fill(GenRadPhotons->GetEntries());
1398 for(UInt_t i=0; i<GenRadPhotons->GetEntries(); i++) {
1399 hDGenRadPhotons[1]->Fill(TMath::Min(GenRadPhotons->At(i)->Pt(),199.999));
1400 hDGenRadPhotons[2]->Fill(TMath::Min(GenRadPhotons->At(i)->AbsEta(),4.999));
1401 hDGenRadPhotons[3]->Fill(TMath::Min((double)GenRadPhotons->At(i)->DistinctMother()->Status(),
1402 19.499));
1403 hDGenRadPhotons[4]->Fill(GenRadPhotons->At(i)->IsGenerated() +
1404 2*GenRadPhotons->At(i)->IsSimulated());
1405 if(GenRadPhotons->At(i)->DistinctMother()){
1406 hDGenRadPhotons[5]->Fill(TMath::Min(
1407 MathUtils::DeltaR(*GenRadPhotons->At(i),
1408 *GenRadPhotons->At(i)->DistinctMother()),
1409 4.999));
1410 CompositeParticle *object = new CompositeParticle();
1411 object->AddDaughter(GenRadPhotons->At(i));
1412 object->AddDaughter(GenRadPhotons->At(i)->DistinctMother());
1413 hDGenRadPhotons[6]->Fill(TMath::Min(object->Mass(),99.999));
1414 delete object;
1415 }
1416 Int_t Mother = 0;
1417 if(GenRadPhotons->At(i)->DistinctMother() &&
1418 GenRadPhotons->At(i)->DistinctMother()->Is(MCParticle::kMu)) Mother = 1;
1419 hDGenRadPhotons[7]->Fill(Mother);
1420 }
1421
1422 // ISR photons
1423 hDGenISRPhotons[0]->Fill(GenISRPhotons->GetEntries());
1424 for(UInt_t i=0; i<GenISRPhotons->GetEntries(); i++) {
1425 hDGenISRPhotons[1]->Fill(TMath::Min(GenISRPhotons->At(i)->Pt(),199.999));
1426 hDGenISRPhotons[2]->Fill(TMath::Min(GenISRPhotons->At(i)->AbsEta(),4.999));
1427 hDGenISRPhotons[3]->Fill(TMath::Min((Double_t)GenISRPhotons->At(i)->DistinctMother()->Status(),
1428 19.499));
1429 hDGenISRPhotons[4]->Fill(GenISRPhotons->At(i)->IsGenerated() +
1430 2*GenISRPhotons->At(i)->IsSimulated());
1431 if(GenISRPhotons->At(i)->DistinctMother()){
1432 hDGenISRPhotons[5]->Fill(TMath::Min(
1433 MathUtils::DeltaR(*GenISRPhotons->At(i),
1434 *GenISRPhotons->At(i)->DistinctMother()),4.999));
1435 CompositeParticle *object = new CompositeParticle();
1436 object->AddDaughter(GenISRPhotons->At(i));
1437 object->AddDaughter(GenISRPhotons->At(i)->DistinctMother());
1438 hDGenISRPhotons[6]->Fill(TMath::Min(object->Mass(),99.999));
1439 delete object;
1440 }
1441 }
1442 } // Fill histograms
1443
1444 // Apply ISR+Rad filter (but filling all histograms)
1445 if (fApplyISRFilter == kTRUE &&
1446 (GenISRPhotons->GetEntries() > 0 || GenRadPhotons->GetEntries() > 0)) {
1447 SkipEvent();
1448 }
1449 }
1450
1451 //--------------------------------------------------------------------------------------------------
1452 void GeneratorMod::SlaveBegin()
1453 {
1454 // Book branch and histograms if wanted.
1455
1456 if(fIsData == kFALSE){
1457 ReqEventObject(fMCPartName, fParticles, kTRUE);
1458 }
1459
1460 // Publish Arrays For the Output Module
1461 PublishObj(fGenLeptons);
1462 PublishObj(fGenAllLeptons);
1463 PublishObj(fGenTaus);
1464 PublishObj(fGenNeutrinos);
1465 PublishObj(fGenQuarks);
1466 PublishObj(fGenqqHs);
1467 PublishObj(fGenBosons);
1468 PublishObj(fGenPhotons);
1469 PublishObj(fGenRadPhotons);
1470 PublishObj(fGenISRPhotons);
1471
1472 // fill histograms
1473 if (GetFillHist()) {
1474 // MET
1475 AddTH1(hDGenMet[0],"hDGenMet_0","Gen MET Pt;p_{t} [GeV];#",200,0,200);
1476 AddTH1(hDGenMet[1],"hDGenMet_1","Gen MET Px;p_{x} [GeV];#",400,-200,200);
1477 AddTH1(hDGenMet[2],"hDGenMet_2","Gen MET Py;p_{y} [GeV];#",400,-200,200);
1478 AddTH1(hDGenMet[3],"hDGenMet_3","Gen MET Pz;p_{z} [GeV];#",400,-1000,1000);
1479
1480 // pt min for leptons from W/Z
1481 AddTH1(hDGenPtMin, "hDGenPtMin","Pt min leptons from W/Z;p_{t} [GeV];#",200,0.0,200.0);
1482
1483 // leptons from W
1484 AddTH1(hDGenLeptons[0], "hDGenLeptons_0",
1485 "Number of leptons from W/Z;N_{leptons};#",10,-0.5,9.5);
1486 AddTH1(hDGenLeptons[1], "hDGenLeptons_1","Pt leptons from W/Z;p_{t} [GeV];#",100,0.0,200.0);
1487 AddTH1(hDGenLeptons[2], "hDGenLeptons_2","Eta leptons from W/Z;#eta;#",50,0.0,5.0);
1488 AddTH1(hDGenLeptons[3], "hDGenLeptons_3","Phi leptons from W/Z;#phi;#",90,0.0,180.0);
1489 AddTH1(hDGenLeptons[4], "hDGenLeptons_4","Dilepton mass from W/Z;m_{ll};#",1000,0.0,1000.0);
1490 AddTH1(hDGenLeptons[5], "hDGenLeptons_5","Eta Max for 2 lepton case;#eta;#",50,0.0,5.0);
1491 AddTH1(hDGenLeptons[6], "hDGenLeptons_6","Eta Min for 2 lepton case;#eta;#",50,0.0,5.0);
1492 AddTH1(hDGenLeptons[7], "hDGenLeptons_7",
1493 "Pt Max for 2 lepton case;p_{t} [GeV];#",100,0.0,200.0);
1494 AddTH1(hDGenLeptons[8], "hDGenLeptons_8",
1495 "Pt Min for 2 lepton case;p_{t} [GeV];#",100,0.0,200.0);
1496 AddTH1(hDGenLeptons[9], "hDGenLeptons_9",
1497 "DiLepton mass for 2 lepton case;p_{t} [GeV];#",1000,0.0,1000.0);
1498 AddTH1(hDGenLeptons[10],"hDGenLeptons_10",
1499 "Delta Phi ll for 2 lepton case;#Delta#phi_{ll};#",90,0.0,180.0);
1500 AddTH1(hDGenLeptons[11],"hDGenLeptons_11","Delta R ll;#Delta R_{ll};#",100,0.0,5.0);
1501 AddTH1(hDGenLeptons[12],"hDGenLeptons_12",
1502 "Pt Max for 3 lepton case;p_{t} [GeV];#",100,0.0,200.0);
1503 AddTH1(hDGenLeptons[13],"hDGenLeptons_13",
1504 "Pt 2nd for 3 lepton case;p_{t} [GeV];#",100,0.0,200.0);
1505 AddTH1(hDGenLeptons[14],"hDGenLeptons_14",
1506 "Pt Min for 3 lepton case;p_{t} [GeV];#",100,0.0,200.0);
1507 AddTH1(hDGenLeptons[15],"hDGenLeptons_15",
1508 "Dilepton mass for 3 lepton case;m_{ll};#",1000,0.0,1000.0);
1509 AddTH1(hDGenLeptons[16],"hDGenLeptons_16",
1510 "Trilepton mass for 3 lepton case;m_{lll};#",1000,0.0,1000.0);
1511 AddTH1(hDGenLeptons[17],"hDGenLeptons_17",
1512 "Delta R Minimum between leptons for 3 lepton case;#Delta R_{ll};#",100,0.0,5.0);
1513 AddTH1(hDGenLeptons[18],"hDGenLeptons_18",
1514 "Pt Max for 4 lepton case;p_{t} [GeV];#",100,0.0,200.0);
1515 AddTH1(hDGenLeptons[19],"hDGenLeptons_19",
1516 "Pt 2nd for 4 lepton case;p_{t} [GeV];#",100,0.0,200.0);
1517 AddTH1(hDGenLeptons[20],"hDGenLeptons_20",
1518 "Pt 3rd for 4 lepton case;p_{t} [GeV];#",100,0.0,200.0);
1519 AddTH1(hDGenLeptons[21],"hDGenLeptons_21",
1520 "Pt 4th for 4 lepton case;#",100,0.0,200.0);
1521 AddTH1(hDGenLeptons[22],"hDGenLeptons_22",
1522 "Dilepton mass for 4 lepton case;m_{ll};#",1000,0.0,1000.0);
1523 AddTH1(hDGenLeptons[23],"hDGenLeptons_23",
1524 "Fourlepton mass for 3 lepton case;m_{llll};#",1000,0.0,1000.0);
1525 AddTH1(hDGenLeptons[24],"hDGenLeptons_24",
1526 "Delta R Minimum between leptons for 4 lepton case;#Delta R_{ll};#",100,0.0,5.0);
1527
1528 // all leptons
1529 AddTH1(hDGenAllLeptons[0], "hDGenAllLeptons_0",
1530 "Number of all leptons;N_{leptons};#",10,-0.5,9.5);
1531 AddTH1(hDGenAllLeptons[1], "hDGenAllLeptons_1","Pt all leptons;p_{t} [GeV];#",400,0.0,200.0);
1532 AddTH1(hDGenAllLeptons[2], "hDGenAllLeptons_2","Eta all leptons;#eta;#",50,0.0,5.0);
1533 AddTH1(hDGenAllLeptons[3], "hDGenAllLeptons_3","Phi all leptons;#phi;#",90,0.0,180.0);
1534 AddTH1(hDGenAllLeptons[4], "hDGenAllLeptons_4","Pt second lepton;p_{t} [GeV];#",400,0.0,200.0);
1535
1536 // taus
1537 AddTH1(hDGenTaus[0], "hDGenTaus_0","Number of taus;N_{tau};#",10,-0.5,9.5);
1538 AddTH1(hDGenTaus[1], "hDGenTaus_1","Pt taus;p_{t} [GeV];#",100,0.0,200.0);
1539 AddTH1(hDGenTaus[2], "hDGenTaus_2","Eta taus;#eta;#",50,0.0,5.0);
1540 AddTH1(hDGenTaus[3], "hDGenTaus_3","Phi taus;#phi;#",90,0.0,180.0);
1541
1542 // neutrinos
1543 AddTH1(hDGenNeutrinos[0], "hDGenNeutrinos_0","Number of neutrinos;N_{#nu};#",10,-0.5,9.5);
1544 AddTH1(hDGenNeutrinos[1], "hDGenNeutrinos_1","Pt neutrinos;p_{t} [GeV];#",100,0.0,200.0);
1545 AddTH1(hDGenNeutrinos[2], "hDGenNeutrinos_2","Eta neutrinos;#eta;#",50,0.0,5.0);
1546 AddTH1(hDGenNeutrinos[3], "hDGenNeutrinos_3","Phi neutrinos;#phi;#",90,0.0,180.0);
1547
1548 // quarks
1549 AddTH1(hDGenQuarks[0], "hDGenQuarks_0", "Number of quarks;N_{quarks};#",10,-0.5,9.5);
1550 AddTH1(hDGenQuarks[1], "hDGenQuarks_1", "dijet pt for quarks;p_{t} [GeV];#",200,0.0,400.);
1551 AddTH1(hDGenQuarks[2], "hDGenQuarks_2", "dijet mass for quarks;m_{jj};#",2000,0.0,2000.);
1552 AddTH1(hDGenQuarks[3], "hDGenQuarks_3",
1553 "dijet pt for quarks with |eta|<2.5;p_{t} [GeV];#",200,0.0,400.);
1554 AddTH1(hDGenQuarks[4], "hDGenQuarks_4",
1555 "dijet mass for quarks with |eta|<2.5;m_{jj};#",2000,0.0,2000.);
1556 AddTH1(hDGenQuarks[5], "hDGenQuarks_5", "Pt for b quarks;p_{t} [GeV];#",200,0.0,400.);
1557 AddTH1(hDGenQuarks[6], "hDGenQuarks_6", "Eta for b quarks;#eta;#",200,-10.0,10.);
1558 AddTH1(hDGenQuarks[7], "hDGenQuarks_7",
1559 "Phi for b quarks;#phi;#",200,-TMath::Pi(),TMath::Pi());
1560 AddTH1(hDGenQuarks[8], "hDGenQuarks_8",
1561 "Pt for b quarks with |eta|<2.5;p_{t} [GeV];#",200,0.0,400.);
1562 AddTH1(hDGenQuarks[9], "hDGenQuarks_9",
1563 "Eta for b quarks with |eta|<2.5;#eta;#",200,-10.0,10.);
1564 AddTH1(hDGenQuarks[10],"hDGenQuarks_10",
1565 "Phi for b quarks with |eta|<2.5;#phi;#",200,-TMath::Pi(),TMath::Pi());
1566 AddTH1(hDGenQuarks[11],"hDGenQuarks_11","Pt for t quarks;p_{t} [GeV];#",200,0.0,400.);
1567 AddTH1(hDGenQuarks[12],"hDGenQuarks_12","Eta for t quarks;#eta;#",200,-10.0,10.);
1568 AddTH1(hDGenQuarks[13],"hDGenQuarks_13",
1569 "Phi for t quarks;#phi;#",200,-TMath::Pi(),TMath::Pi());
1570 AddTH1(hDGenQuarks[14],"hDGenQuarks_14","Pt for light quarks;p_{t} [GeV];#",200,0.0,400.);
1571 AddTH1(hDGenQuarks[15],"hDGenQuarks_15","Eta for light quarks;#eta;#",200,-10.0,10.);
1572 AddTH1(hDGenQuarks[16],"hDGenQuarks_16",
1573 "Phi for light quarks;#phi;#",200,-TMath::Pi(),TMath::Pi());
1574
1575 // qqH
1576 AddTH1(hDGenWBF[0], "hDGenWBF_0",
1577 "Delta Phi jj for WBF quarks;#Delta Phi_{jj};#",90,0.0,180.);
1578 AddTH1(hDGenWBF[1], "hDGenWBF_1",
1579 "Delta Eta jj for WBF quarks;#Delta #eta_{jj};#",100,0.0,10.);
1580 AddTH1(hDGenWBF[2], "hDGenWBF_2",
1581 "Pt max for WBF quarks;p_{t} [GeV];#",200,0.0,400.);
1582 AddTH1(hDGenWBF[3], "hDGenWBF_3",
1583 "Pt min for WBF quarks;p_{t} [GeV];#",200,0.0,400.);
1584 AddTH1(hDGenWBF[4], "hDGenWBF_4",
1585 "dijet mass for WBF quarks;m_{jj};#",200,0.0,4000.);
1586
1587 // bosons
1588 AddTH1(hDGenBosons[0], "hDGenBosons_0", "Number of bosons;N_{bosons};#",10,-0.5,9.5);
1589 AddTH1(hDGenBosons[1], "hDGenBosons_1", "Pt of bosons;p_{t} [GeV];#",200,0.0,400.0);
1590 AddTH1(hDGenBosons[2], "hDGenBosons_2", "Eta of bosons;#eta;#",100,-5.0,5.0);
1591 AddTH1(hDGenBosons[3], "hDGenBosons_3", "Mass of bosons;Mass;#",2000,0.0,2000.0);
1592 AddTH1(hDGenBosons[4], "hDGenBosons_4", "Mass of bosons;m_{V};#",200,0.0,200.0);
1593 AddTH1(hDGenBosons[5], "hDGenBosons_5", "Mass of W bosons;m_{W};#",200,0.0,200.0);
1594 AddTH1(hDGenBosons[6], "hDGenBosons_6", "Mass of Z bosons;m_{Z};#",200,0.0,200.0);
1595 AddTH1(hDGenBosons[7], "hDGenBosons_7",
1596 "Number of W bosons + 4 * Z bosons;Number;#",13,-0.5,12.5);
1597
1598 // photons
1599 AddTH1(hDGenPhotons[0], "hDGenPhotons_0", "Number of photons;N_{photons};#",10,-0.5,9.5);
1600 AddTH1(hDGenPhotons[1], "hDGenPhotons_1", "Pt of photons;p_{t} [GeV];#",200,0.0,400.0);
1601 AddTH1(hDGenPhotons[2], "hDGenPhotons_2", "Eta of photons;#eta;#",100,-5.0,5.0);
1602
1603 // rad photons
1604 AddTH1(hDGenRadPhotons[0], "hDGenRadPhotons_0",
1605 "Number of radiative photons;N_{photons};#",10,-0.5,9.5);
1606 AddTH1(hDGenRadPhotons[1], "hDGenRadPhotons_1",
1607 "Pt of radiative photons;p_{t} [GeV];#",400,0.0,200.0);
1608 AddTH1(hDGenRadPhotons[2], "hDGenRadPhotons_2",
1609 "Eta of radiative photons;#eta;#",100,0.0,5.0);
1610 AddTH1(hDGenRadPhotons[3], "hDGenRadPhotons_3",
1611 "Status of mother of radiative photons;Status;#",20,-0.5,19.5);
1612 AddTH1(hDGenRadPhotons[4], "hDGenRadPhotons_4",
1613 "IsGenerated+2*IsSimulated of radiative photons;IsGenerated+2*IsSimulated;#",
1614 4,-0.5,3.5);
1615 AddTH1(hDGenRadPhotons[5], "hDGenRadPhotons_5",
1616 "Delta R between photon and mother of radiative photons;#Delta R;#",500,0.0,5.0);
1617 AddTH1(hDGenRadPhotons[6], "hDGenRadPhotons_6",
1618 "Mass photon-photon mother;Mass;#",500,0.0,100.0);
1619 AddTH1(hDGenRadPhotons[7], "hDGenRadPhotons_7",
1620 "Number of radiative photon with muon as a mother;Status;#",2,-0.5,1.5);
1621
1622 // ISR photons
1623 AddTH1(hDGenISRPhotons[0], "hDGenISRPhotons_0",
1624 "Number of ISR photons;N_{photons};#",10,-0.5,9.5);
1625 AddTH1(hDGenISRPhotons[1], "hDGenISRPhotons_1",
1626 "Pt of ISR photons;p_{t} [GeV];#",400,0.0,200.0);
1627 AddTH1(hDGenISRPhotons[2], "hDGenISRPhotons_2",
1628 "Eta of ISR photons;#eta;#",100,0.0,5.0);
1629 AddTH1(hDGenISRPhotons[3], "hDGenISRPhotons_3",
1630 "Status of mother of radiative photons;#eta;#",20,-0.5,19.5);
1631 AddTH1(hDGenISRPhotons[4], "hDGenISRPhotons_4",
1632 "IsGenerated+2*IsSimulated of radiative photons;IsGenerated+2*IsSimulated;#",
1633 4,-0.5,3.5);
1634 AddTH1(hDGenISRPhotons[5], "hDGenISRPhotons_5",
1635 "Delta R between photon and mother of ISR photons;#Delta R;#",500,0.0,5.0);
1636 AddTH1(hDGenISRPhotons[6], "hDGenISRPhotons_6",
1637 "Mass photon-photon mother;Mass;#",500,0.0,100.0);
1638
1639 // auxiliar for MG studies
1640 AddTH1(hDVMass[0], "hDVMass_0", "Mass of munu candidates ;Mass [GeV];#",200,0.,200.);
1641 AddTH1(hDVMass[1], "hDVMass_1", "Mass of elnu candidates ;Mass [GeV];#",200,0.,200.);
1642 AddTH1(hDVMass[2], "hDVMass_2", "Mass of taunu candidates ;Mass [GeV];#",200,0.,200.);
1643 AddTH1(hDVMass[3], "hDVMass_3", "Mass of mumu candidates ;Mass [GeV];#",200,0.,200.);
1644 AddTH1(hDVMass[4], "hDVMass_4", "Mass of ee candidates ;Mass [GeV];#",200,0.,200.);
1645 AddTH1(hDVMass[5], "hDVMass_5", "Mass of tautau candidates;Mass [GeV];#",200,0.,200.);
1646 AddTH1(hDVMass[6], "hDVMass_6", "Mass of numunumu candidates;Mass [GeV];#",200,0.,200.);
1647 AddTH1(hDVMass[7], "hDVMass_7", "Mass of nuenue candidates;Mass [GeV];#",200,0.,200.);
1648 AddTH1(hDVMass[8], "hDVMass_8", "Mass of nutaunutau candidates;Mass [GeV];#",200,0.,200.);
1649 AddTH1(hDVMass[9], "hDVMass_9",
1650 "Mass of munu candidates for t events ;Mass [GeV];#",200,0.,200.);
1651 AddTH1(hDVMass[10],"hDVMass_10",
1652 "Mass of elnu candidates for t events ;Mass [GeV];#",200,0.,200.);
1653 AddTH1(hDVMass[11],"hDVMass_11",
1654 "Mass of taunu candidates for t events;Mass [GeV];#",200,0.,200.);
1655 AddTH1(hDVMass[12],"hDVMass_12",
1656 "Mass of qq candidates for t events;Mass [GeV];#",200,0.,200.);
1657
1658 // Special study about VVjets
1659 AddTH1(hDVVMass[0], "hDVVMass_0", "Mass of munu for WW events;Mass [GeV];#",200,0.,200.);
1660 AddTH1(hDVVMass[1], "hDVVMass_1", "Mass of munu WZ events;Mass [GeV];#",200,0.,200.);
1661 AddTH1(hDVVMass[2], "hDVVMass_2", "Mass of elnu WW events;Mass [GeV];#",200,0.,200.);
1662 AddTH1(hDVVMass[3], "hDVVMass_3", "Mass of elnu WZ events;Mass [GeV];#",200,0.,200.);
1663 AddTH1(hDVVMass[4], "hDVVMass_4", "Mass of taunu WW events;Mass [GeV];#",200,0.,200.);
1664 AddTH1(hDVVMass[5], "hDVVMass_5", "Mass of taunu WZ events;Mass [GeV];#",200,0.,200.);
1665 AddTH1(hDVVMass[6], "hDVVMass_6", "Mass of mumu WZ events;Mass [GeV];#",200,0.,200.);
1666 AddTH1(hDVVMass[7], "hDVVMass_7", "Mass of mumu ZZ events;Mass [GeV];#",200,0.,200.);
1667 AddTH1(hDVVMass[8], "hDVVMass_8", "Mass of ee WZ events;Mass [GeV];#",200,0.,200.);
1668 AddTH1(hDVVMass[9], "hDVVMass_9", "Mass of ee ZZ events;Mass [GeV];#",200,0.,200.);
1669 AddTH1(hDVVMass[10],"hDVVMass_10","Mass of tautau WZ events;Mass [GeV];#",200,0.,200.);
1670 AddTH1(hDVVMass[11],"hDVVMass_11","Mass of tautau ZZ events;Mass [GeV];#",200,0.,200.);
1671 AddTH1(hDVVMass[12],"hDVVMass_12","Mass of numunumu WZ events;Mass [GeV];#",200,0.,200.);
1672 AddTH1(hDVVMass[13],"hDVVMass_13","Mass of numunumu ZZ events;Mass [GeV];#",200,0.,200.);
1673 AddTH1(hDVVMass[14],"hDVVMass_14","Mass of nuenue WZ events;Mass [GeV];#",200,0.,200.);
1674 AddTH1(hDVVMass[15],"hDVVMass_15","Mass of nuenue ZZ events;Mass [GeV];#",200,0.,200.);
1675 AddTH1(hDVVMass[16],"hDVVMass_16","Mass of nutaunutau WZ events;Mass [GeV];#",200,0.,200.);
1676 AddTH1(hDVVMass[17],"hDVVMass_17","Mass of nutaunutau ZZ events;Mass [GeV];#",200,0.,200.);
1677 AddTH1(hDVVMass[18],"hDVVMass_18","Ratios for WW events;Type;#",7,-0.5,6.5);
1678 AddTH1(hDVVMass[19],"hDVVMass_19","Ratios for WZ events;Type;#",10,-0.5,9.5);
1679 AddTH1(hDVVMass[20],"hDVVMass_20","Ratios for ZZ2l events;Type;#",7,-0.5,6.5);
1680 AddTH1(hDVVMass[21],"hDVVMass_21","Ratios for ZZ4l events;Type;#",16,-0.5,15.5);
1681 AddTH1(hDVVMass[22],"hDVVMass_22","Maximum mass for WW events;Mass [GeV];#",200,0.,200.);
1682 AddTH1(hDVVMass[23],"hDVVMass_23","Minimum mass for WW events;Mass [GeV];#",200,0.,200.);
1683 AddTH1(hDVVMass[24],"hDVVMass_24","Maximum mass for WZ events;Mass [GeV];#",200,0.,200.);
1684 AddTH1(hDVVMass[25],"hDVVMass_25","Minimum mass for WZ events;Mass [GeV];#",200,0.,200.);
1685 AddTH1(hDVVMass[26],"hDVVMass_26","Maximum mass for ZZ events;Mass [GeV];#",200,0.,200.);
1686 AddTH1(hDVVMass[27],"hDVVMass_27","Minimum mass for ZZ events;Mass [GeV];#",200,0.,200.);
1687 }
1688 }