ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/MitPhysics/Mods/src/GeneratorMod.cc
Revision: 1.68
Committed: Sat Oct 29 14:11:52 2011 UTC (13 years, 6 months ago) by ceballos
Content type: text/plain
Branch: MAIN
CVS Tags: Mit_028a, Mit_028, Mit_027, Mit_027a, Mit_025e, Mit_025d
Changes since 1.67: +23 -1 lines
Log Message:
new filter

File Contents

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