ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/MitPhysics/Mods/src/GeneratorMod.cc
Revision: 1.62
Committed: Mon Sep 27 15:49:24 2010 UTC (14 years, 7 months ago) by ceballos
Content type: text/plain
Branch: MAIN
CVS Tags: Mit_015b, Mit_015a, Mit_015
Changes since 1.61: +14 -6 lines
Log Message:
small fix

File Contents

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