ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/MitPhysics/Mods/src/GeneratorMod.cc
Revision: 1.64
Committed: Sat Feb 5 05:48:09 2011 UTC (14 years, 3 months ago) by ceballos
Content type: text/plain
Branch: MAIN
CVS Tags: Mit_023, Mit_022a, Mit_022, Mit_020d, TMit_020d, Mit_020c, Mit_021, Mit_021pre2, Mit_021pre1, Mit_020b, Mit_020a, Mit_020, Mit_020pre1
Changes since 1.63: +43 -35 lines
Log Message:
isolation

File Contents

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