ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/MitPhysics/Mods/src/GeneratorMod.cc
Revision: 1.41
Committed: Tue May 19 11:42:20 2009 UTC (15 years, 11 months ago) by loizides
Content type: text/plain
Branch: MAIN
Changes since 1.40: +3 -3 lines
Log Message: 
Use the new ReqEventObject model.

File Contents

# User Rev Content
1 loizides 1.41 // $Id: GeneratorMod.cc,v 1.40 2009/04/30 08:09:32 loizides Exp $
2 loizides 1.1
3     #include "MitPhysics/Mods/interface/GeneratorMod.h"
4     #include "MitCommon/MathTools/interface/MathUtils.h"
5 loizides 1.8 #include "MitPhysics/Init/interface/ModNames.h"
6 loizides 1.1 #include <TH1D.h>
7     #include <TH2D.h>
8    
9     using namespace mithep;
10    
11     ClassImp(mithep::GeneratorMod)
12    
13     //--------------------------------------------------------------------------------------------------
14 loizides 1.8 GeneratorMod::GeneratorMod(const char *name, const char *title) :
15 loizides 1.1 BaseMod(name,title),
16 ceballos 1.32 fPrintDebug(kFALSE),
17 loizides 1.1 fMCPartName(Names::gkMCPartBrn),
18 ceballos 1.30 fMCMETName(ModNames::gkMCMETName),
19 loizides 1.8 fMCLeptonsName(ModNames::gkMCLeptonsName),
20     fMCAllLeptonsName(ModNames::gkMCAllLeptonsName),
21     fMCTausName(ModNames::gkMCTausName),
22     fMCNeutrinosName(ModNames::gkMCNeutrinosName),
23     fMCQuarksName(ModNames::gkMCQuarksName),
24     fMCqqHsName(ModNames::gkMCqqHsName),
25     fMCBosonsName(ModNames::gkMCBosonsName),
26 ceballos 1.10 fMCPhotonsName(ModNames::gkMCPhotonsName),
27 ceballos 1.32 fMCRadPhotonsName(ModNames::gkMCRadPhotonsName),
28     fMCISRPhotonsName(ModNames::gkMCISRPhotonsName),
29 ceballos 1.14 fPtLeptonMin(0.0),
30     fEtaLeptonMax(5.0),
31     fPtPhotonMin(0.0),
32 loizides 1.15 fEtaPhotonMax(5.0),
33 ceballos 1.32 fPtRadPhotonMin(0.0),
34     fEtaRadPhotonMax(5.0),
35 loizides 1.27 fPdgIdCut(0),
36     fMassMinCut(-FLT_MAX),
37     fMassMaxCut(FLT_MAX),
38 ceballos 1.34 fApplyISRFilter(kFALSE),
39 loizides 1.15 fParticles(0)
40 loizides 1.1 {
41     // Constructor.
42     }
43    
44     //--------------------------------------------------------------------------------------------------
45     void GeneratorMod::Process()
46     {
47 loizides 1.5 // Process entries of the tree.
48 loizides 1.1
49 loizides 1.5 // these arrays will be filled in the loop of particles
50 ceballos 1.30 MetOArr *GenMet = new MetOArr;
51     GenMet->SetName(fMCMETName);
52     GenMet->SetOwner(kTRUE);
53 loizides 1.9 MCParticleOArr *GenLeptons = new MCParticleOArr;
54 loizides 1.13 GenLeptons->SetName(fMCLeptonsName);
55 loizides 1.9 MCParticleOArr *GenAllLeptons = new MCParticleOArr;
56 loizides 1.13 GenAllLeptons->SetName(fMCAllLeptonsName);
57 loizides 1.9 MCParticleOArr *GenTaus = new MCParticleOArr;
58 loizides 1.13 GenTaus->SetName(fMCTausName);
59     GenTaus->SetOwner(kTRUE);
60 loizides 1.9 MCParticleOArr *GenNeutrinos = new MCParticleOArr;
61 loizides 1.13 GenNeutrinos->SetName(fMCNeutrinosName);
62 loizides 1.9 MCParticleOArr *GenQuarks = new MCParticleOArr;
63 loizides 1.13 GenQuarks->SetName(fMCQuarksName);
64 loizides 1.9 MCParticleOArr *GenqqHs = new MCParticleOArr;
65 loizides 1.13 GenqqHs->SetName(fMCqqHsName);
66 loizides 1.9 MCParticleOArr *GenBosons = new MCParticleOArr;
67 loizides 1.13 GenBosons->SetName(fMCBosonsName);
68 ceballos 1.10 MCParticleOArr *GenPhotons = new MCParticleOArr;
69 loizides 1.13 GenPhotons->SetName(fMCPhotonsName);
70 ceballos 1.32 MCParticleOArr *GenRadPhotons = new MCParticleOArr;
71     GenRadPhotons->SetName(fMCRadPhotonsName);
72     MCParticleOArr *GenISRPhotons = new MCParticleOArr;
73     GenISRPhotons->SetName(fMCISRPhotonsName);
74    
75 loizides 1.40 if(fPrintDebug)
76     printf("\n************ Next Event ************\n\n");
77 loizides 1.1
78 loizides 1.5 // load MCParticle branch
79 loizides 1.41 LoadEventObject(fMCPartName, fParticles);
80 loizides 1.5
81 ceballos 1.30 Double_t totalMET[3] = {0.0, 0.0, 0.0};
82 loizides 1.6 Bool_t isqqH = kFALSE;
83 loizides 1.5 for (UInt_t i=0; i<fParticles->GetEntries(); ++i) {
84 loizides 1.7 const MCParticle *p = fParticles->At(i);
85 ceballos 1.32
86 loizides 1.40 if(fPrintDebug)
87     p->Print("l");
88 ceballos 1.32
89 loizides 1.40 // rad photons
90 ceballos 1.34 if(p->Is(MCParticle::kGamma) && p->HasMother() && p->Mother()->Status() == 3 &&
91 ceballos 1.33 (p->Mother()->Is(MCParticle::kEl) || p->Mother()->Is(MCParticle::kMu) ||
92 ceballos 1.34 p->Mother()->Is(MCParticle::kTau)) &&
93 ceballos 1.32 p->Pt() > fPtRadPhotonMin && p->AbsEta() < fEtaRadPhotonMax) {
94     GenRadPhotons->Add(p);
95     }
96    
97     // ISR photons
98     if(p->Is(MCParticle::kGamma) && p->HasMother() && p->Mother()->IsQuark()) {
99     GenISRPhotons->Add(p);
100     }
101    
102 ceballos 1.30 // MET computation at generation level
103     if (p->Status() == 1 && !p->IsNeutrino()) {
104     totalMET[0] = totalMET[0] + p->Px();
105     totalMET[1] = totalMET[1] + p->Py();
106     totalMET[2] = totalMET[2] + p->Pz();
107     }
108 loizides 1.5
109     if (!p->IsGenerated()) continue;
110    
111     // muons/electrons from W/Z decays
112 loizides 1.7 if ((p->Is(MCParticle::kEl) || p->Is(MCParticle::kMu)) && p->Status() == 1) {
113 ceballos 1.14 if (p->Pt() > fPtLeptonMin && p->AbsEta() < fEtaLeptonMax) {
114 loizides 1.5 GenAllLeptons->Add(p);
115     }
116 loizides 1.6 Bool_t isGoodLepton = kFALSE;
117 loizides 1.7 const MCParticle *pm = p;
118 loizides 1.5 while (pm->HasMother() && isGoodLepton == kFALSE) {
119 loizides 1.21 if (pm->PdgId() == 92) // string reached, terminate loop
120     break;
121 loizides 1.29 if (pm->Mother()->Is(MCParticle::kZ) || pm->Mother()->Is(MCParticle::kW) ||
122     pm->Mother()->Is(MCParticle::kZp) || pm->Mother()->Is(MCParticle::kWp) ||
123     pm->Mother()->Is(MCParticle::kH)) {
124 loizides 1.5 GenLeptons->Add(p);
125     isGoodLepton = kTRUE;
126 loizides 1.7 break;
127 loizides 1.8 } else if (pm->Mother()->Is(MCParticle::kPi0) || pm->Mother()->Is(MCParticle::kEta)) {
128 loizides 1.7 // this is fake, but it is a trick to get rid of these cases and abort the loop
129     break;
130     }
131     pm = pm->Mother();
132 loizides 1.1 }
133 loizides 1.5 }
134    
135 loizides 1.7 // hadronic taus
136     else if (p->Is(MCParticle::kTau) && p->Status() == 2) {
137     if (!p->HasDaughter(MCParticle::kEl) && !p->HasDaughter(MCParticle::kMu)) {
138     const MCParticle *tv = p->FindDaughter(MCParticle::kTauNu);
139     if (tv) {
140     MCParticle *pm_f = new MCParticle(*p);
141     pm_f->SetMom(p->Px()-tv->Px(), p->Py()-tv->Py(),
142     p->Pz()-tv->Pz(), p->E()-tv->E());
143     GenTaus->AddOwned(pm_f);
144     } else {
145 loizides 1.18 SendError(kWarning, "Process", "Could not find a tau neutrino!");
146 loizides 1.5 }
147 loizides 1.1 }
148 loizides 1.5 }
149 loizides 1.1
150 loizides 1.5 // neutrinos
151 loizides 1.7 else if (p->Status() == 1 && p->IsNeutrino()) {
152 loizides 1.5 GenNeutrinos->Add(p);
153     }
154 loizides 1.1
155 loizides 1.5 // quarks from W/Z decays or top particles
156 loizides 1.7 else if (p->IsQuark() && p->HasMother()) {
157     if (p->Mother()->Is(MCParticle::kZ) || p->Mother()->Is(MCParticle::kW) ||
158     p->Is(MCParticle::kTop) || p->Mother()->Is(MCParticle::kTop)) {
159 loizides 1.5 GenQuarks->Add(p);
160 loizides 1.1 }
161 loizides 1.5 }
162 loizides 1.1
163 loizides 1.5 // qqH, information about the forward jets
164 loizides 1.8 else if (isqqH == kFALSE && p->Is(MCParticle::kH)) {
165 loizides 1.5 isqqH = kTRUE;
166 loizides 1.17 const MCParticle *pq1 = fParticles->At(i-1);
167     const MCParticle *pq2 = fParticles->At(i-2);
168 loizides 1.7 if (!pq1 || !pq2) {
169     SendError(kWarning, "Process", "Could not find quark pair!");
170 loizides 1.8 } else if (pq1->IsQuark() && pq2->IsQuark() &&
171     pq1->HasMother() && pq2->HasMother() &&
172 ceballos 1.36 pq1->Mother() == pq2->Mother()) {
173 loizides 1.5 GenqqHs->Add(pq1);
174     GenqqHs->Add(pq2);
175 loizides 1.1 }
176 ceballos 1.23
177 loizides 1.7 if (p->Status() == 3)
178 loizides 1.24 GenBosons->Add(p); // take higgs boson in account here rather in next else if
179 loizides 1.5 }
180 loizides 1.1
181 loizides 1.5 // information about bosons: W, Z, h, Z', W', H0, A0, H+
182 loizides 1.7 else if (p->Status() == 3 &&
183     (p->Is(MCParticle::kZ) || p->Is(MCParticle::kW) || p->Is(MCParticle::kH) ||
184     p->Is(MCParticle::kZp) || p->Is(MCParticle::kZpp) ||
185     p->Is(MCParticle::kH0) || p->Is(MCParticle::kA0) || p->Is(MCParticle::kHp))) {
186 loizides 1.5 GenBosons->Add(p);
187 loizides 1.1 }
188 ceballos 1.10
189 ceballos 1.14 // photons
190 ceballos 1.10 else if (p->Status() == 1 && p->Is(MCParticle::kGamma) &&
191 ceballos 1.14 p->Pt() > fPtPhotonMin && p->AbsEta() < fEtaPhotonMax) {
192 ceballos 1.10 GenPhotons->Add(p);
193     }
194    
195 loizides 1.18 // W/Z -> lnu for Madgraph
196 loizides 1.28 if (p->IsParton() && p->NDaughters() >= 2) {
197 ceballos 1.16 CompositeParticle *diBoson = new CompositeParticle();
198 loizides 1.18 if (p->HasDaughter(MCParticle::kMu) && p->HasDaughter(MCParticle::kMuNu)) {
199 phedex 1.39 diBoson->AddDaughter(p->FindDaughter(MCParticle::kMu));
200     diBoson->AddDaughter(p->FindDaughter(MCParticle::kMuNu));
201 loizides 1.20 if (GetFillHist())
202 loizides 1.18 hDVMass[0]->Fill(TMath::Min(diBoson->Mass(),199.999));
203 loizides 1.28 const MCParticle *tmp_mu = p->FindDaughter(MCParticle::kMu);
204     while (tmp_mu->HasDaughter(MCParticle::kMu) &&
205     tmp_mu->FindDaughter(MCParticle::kMu)->IsGenerated())
206     tmp_mu = tmp_mu->FindDaughter(MCParticle::kMu);
207     GenLeptons->Add(tmp_mu);
208 ceballos 1.16 }
209 loizides 1.18 else if (p->HasDaughter(MCParticle::kEl) && p->HasDaughter(MCParticle::kElNu)) {
210 phedex 1.39 diBoson->AddDaughter(p->FindDaughter(MCParticle::kEl));
211     diBoson->AddDaughter(p->FindDaughter(MCParticle::kElNu));
212 loizides 1.20 if (GetFillHist())
213 loizides 1.18 hDVMass[1]->Fill(TMath::Min(diBoson->Mass(),199.999));
214 loizides 1.28 const MCParticle *tmp_e = p->FindDaughter(MCParticle::kEl);
215     while (tmp_e->HasDaughter(MCParticle::kEl) &&
216     tmp_e->FindDaughter(MCParticle::kEl)->IsGenerated())
217     tmp_e = tmp_e->FindDaughter(MCParticle::kEl);
218     GenLeptons->Add(tmp_e);
219 ceballos 1.16 }
220 loizides 1.18 else if (p->HasDaughter(MCParticle::kTau) && p->HasDaughter(MCParticle::kTauNu)) {
221 phedex 1.39 diBoson->AddDaughter(p->FindDaughter(MCParticle::kTau));
222     diBoson->AddDaughter(p->FindDaughter(MCParticle::kTauNu));
223 loizides 1.20 if (GetFillHist())
224 loizides 1.18 hDVMass[2]->Fill(TMath::Min(diBoson->Mass(),199.999));
225     const MCParticle *tau = p->FindDaughter(MCParticle::kTau);
226     if (tau->HasDaughter(MCParticle::kMu))
227     GenLeptons->Add(tau->FindDaughter(MCParticle::kMu));
228     if (tau->HasDaughter(MCParticle::kEl))
229     GenLeptons->Add(tau->FindDaughter(MCParticle::kEl));
230 loizides 1.28 if (tau->HasDaughter(MCParticle::kTau)) {
231     const MCParticle *tau_second = tau->FindDaughter(MCParticle::kTau);
232     if (tau_second->HasDaughter(MCParticle::kMu))
233     GenLeptons->Add(tau_second->FindDaughter(MCParticle::kMu));
234     if (tau_second->HasDaughter(MCParticle::kEl))
235     GenLeptons->Add(tau_second->FindDaughter(MCParticle::kEl));
236     }
237 loizides 1.18 }
238 phedex 1.39 else if (p->HasDaughter(MCParticle::kMu,kTRUE) && p->HasDaughter(-1*MCParticle::kMu,kTRUE)) {
239     diBoson->AddDaughter(p->FindDaughter(MCParticle::kMu,kTRUE));
240     diBoson->AddDaughter(p->FindDaughter(-1*MCParticle::kMu,kTRUE));
241 loizides 1.20 if (GetFillHist())
242 loizides 1.18 hDVMass[3]->Fill(TMath::Min(diBoson->Mass(),199.999));
243 phedex 1.39 const MCParticle *tmp_mu0 = p->FindDaughter(MCParticle::kMu,kTRUE);
244 sixie 1.35 while (tmp_mu0->HasDaughter(MCParticle::kMu) &&
245     tmp_mu0->FindDaughter(MCParticle::kMu)->IsGenerated())
246     tmp_mu0 = tmp_mu0->FindDaughter(MCParticle::kMu);
247 phedex 1.39 const MCParticle *tmp_mu1 = p->FindDaughter(-1*MCParticle::kMu,kTRUE);
248 sixie 1.35 while (tmp_mu1->HasDaughter(MCParticle::kMu) &&
249     tmp_mu1->FindDaughter(MCParticle::kMu)->IsGenerated())
250     tmp_mu1 = tmp_mu1->FindDaughter(MCParticle::kMu);
251     GenLeptons->Add(tmp_mu0);
252     GenLeptons->Add(tmp_mu1);
253 ceballos 1.16 }
254 phedex 1.39 else if (p->HasDaughter(MCParticle::kEl,kTRUE) && p->HasDaughter(-1*MCParticle::kEl,kTRUE)) {
255     diBoson->AddDaughter(p->FindDaughter(MCParticle::kEl,kTRUE));
256     diBoson->AddDaughter(p->FindDaughter(-1*MCParticle::kEl,kTRUE));
257 loizides 1.20 if (GetFillHist())
258 loizides 1.18 hDVMass[4]->Fill(TMath::Min(diBoson->Mass(),199.999));
259 sixie 1.35 const MCParticle *tmp_e0 = p->Daughter(0);
260     while (tmp_e0->HasDaughter(MCParticle::kEl) &&
261     tmp_e0->FindDaughter(MCParticle::kEl)->IsGenerated())
262     tmp_e0 = tmp_e0->FindDaughter(MCParticle::kEl);
263     const MCParticle *tmp_e1 = p->Daughter(1);
264     while (tmp_e1->HasDaughter(MCParticle::kEl) &&
265     tmp_e1->FindDaughter(MCParticle::kEl)->IsGenerated())
266     tmp_e1 = tmp_e1->FindDaughter(MCParticle::kEl);
267     GenLeptons->Add(tmp_e0);
268     GenLeptons->Add(tmp_e1);
269 ceballos 1.16 }
270 phedex 1.39 else if (p->HasDaughter(MCParticle::kTau,kTRUE) && p->HasDaughter(-1*MCParticle::kTau,kTRUE)) {
271     diBoson->AddDaughter(p->FindDaughter(MCParticle::kTau,kTRUE));
272     diBoson->AddDaughter(p->FindDaughter(-1*MCParticle::kTau,kTRUE));
273 loizides 1.20 if (GetFillHist())
274 loizides 1.18 hDVMass[5]->Fill(TMath::Min(diBoson->Mass(),199.999));
275     const MCParticle *tau0 = p->Daughter(0);
276     if (tau0->HasDaughter(MCParticle::kMu))
277     GenLeptons->Add(tau0->FindDaughter(MCParticle::kMu));
278     if (tau0->HasDaughter(MCParticle::kEl))
279     GenLeptons->Add(tau0->FindDaughter(MCParticle::kEl));
280     const MCParticle *tau1 = p->Daughter(1);
281     if (tau1->HasDaughter(MCParticle::kMu))
282     GenLeptons->Add(tau1->FindDaughter(MCParticle::kMu));
283     if (tau1->HasDaughter(MCParticle::kEl))
284     GenLeptons->Add(tau1->FindDaughter(MCParticle::kEl));
285 sixie 1.35 if (tau0->HasDaughter(MCParticle::kTau)) {
286     const MCParticle *tau0_second = tau0->FindDaughter(MCParticle::kTau);
287     if (tau0_second->HasDaughter(MCParticle::kMu))
288     GenLeptons->Add(tau0_second->FindDaughter(MCParticle::kMu));
289     if (tau0_second->HasDaughter(MCParticle::kEl))
290     GenLeptons->Add(tau0_second->FindDaughter(MCParticle::kEl));
291     }
292     if (tau1->HasDaughter(MCParticle::kTau)) {
293     const MCParticle *tau1_second = tau1->FindDaughter(MCParticle::kTau);
294     if (tau1_second->HasDaughter(MCParticle::kMu))
295     GenLeptons->Add(tau1_second->FindDaughter(MCParticle::kMu));
296     if (tau1_second->HasDaughter(MCParticle::kEl))
297     GenLeptons->Add(tau1_second->FindDaughter(MCParticle::kEl));
298     }
299 ceballos 1.16 }
300     delete diBoson;
301     }
302    
303 loizides 1.18 // t -> lnu for Madgraph
304     if (p->Is(MCParticle::kTop)) {
305 ceballos 1.16 CompositeParticle *diBoson = new CompositeParticle();
306 loizides 1.18 if (p->HasDaughter(MCParticle::kMu) && p->HasDaughter(MCParticle::kMuNu)) {
307 ceballos 1.16 diBoson->AddDaughter(p->FindDaughter(MCParticle::kMu));
308     diBoson->AddDaughter(p->FindDaughter(MCParticle::kMuNu));
309 loizides 1.20 if (GetFillHist())
310 loizides 1.18 hDVMass[6]->Fill(TMath::Min(diBoson->Mass(),199.999));
311 ceballos 1.16 GenLeptons->Add(p->FindDaughter(MCParticle::kMu));
312     }
313 loizides 1.18 else if (p->HasDaughter(MCParticle::kEl) && p->HasDaughter(MCParticle::kElNu)) {
314 ceballos 1.16 diBoson->AddDaughter(p->FindDaughter(MCParticle::kEl));
315     diBoson->AddDaughter(p->FindDaughter(MCParticle::kElNu));
316 loizides 1.20 if (GetFillHist())
317 loizides 1.18 hDVMass[7]->Fill(TMath::Min(diBoson->Mass(),199.999));
318 ceballos 1.16 GenLeptons->Add(p->FindDaughter(MCParticle::kEl));
319     }
320 loizides 1.18 else if (p->HasDaughter(MCParticle::kTau) && p->HasDaughter(MCParticle::kTauNu)) {
321 ceballos 1.16 diBoson->AddDaughter(p->FindDaughter(MCParticle::kTau));
322     diBoson->AddDaughter(p->FindDaughter(MCParticle::kTauNu));
323 loizides 1.20 if (GetFillHist())
324 loizides 1.18 hDVMass[8]->Fill(TMath::Min(diBoson->Mass(),199.999));
325     const MCParticle *tau = p->FindDaughter(MCParticle::kTau);
326     if (tau->HasDaughter(MCParticle::kMu))
327     GenLeptons->Add(tau->FindDaughter(MCParticle::kMu));
328     if (tau->HasDaughter(MCParticle::kEl))
329     GenLeptons->Add(tau->FindDaughter(MCParticle::kEl));
330     }
331     else if (!p->HasDaughter(MCParticle::kW)) {
332     for(UInt_t nd=0; nd<p->NDaughters(); ++nd)
333     if (p->Daughter(nd)->IsNot(MCParticle::kBottom) &&
334     p->Daughter(nd)->IsNot(MCParticle::kGamma))
335     diBoson->AddDaughter(p->Daughter(nd));
336 loizides 1.20 if (GetFillHist())
337 loizides 1.18 hDVMass[9]->Fill(TMath::Min(diBoson->Mass(),199.999));
338 ceballos 1.16 }
339     delete diBoson;
340     }
341 ceballos 1.26
342 loizides 1.27 // mass cut for given pid
343     if(fPdgIdCut && p->Is(fPdgIdCut) &&
344 ceballos 1.26 (p->Mass() < fMassMinCut || p->Mass() > fMassMaxCut)) {
345     SkipEvent();
346     return;
347 sixie 1.35 }
348 loizides 1.27 } // end loop of particles
349 loizides 1.1
350 ceballos 1.30 Met *theMET = new Met(totalMET[0], totalMET[1]);
351     theMET->SetElongitudinal(totalMET[2]);
352     GenMet->AddOwned(theMET);
353    
354 loizides 1.19 // sort according to pt
355     GenLeptons->Sort();
356     GenAllLeptons->Sort();
357     GenTaus->Sort();
358     GenNeutrinos->Sort();
359     GenQuarks->Sort();
360     GenqqHs->Sort();
361     GenBosons->Sort();
362     GenPhotons->Sort();
363 ceballos 1.32 GenRadPhotons->Sort();
364     GenISRPhotons->Sort();
365 loizides 1.19
366 loizides 1.8 // add objects to this event for other modules to use
367 ceballos 1.30 AddObjThisEvt(GenMet);
368 loizides 1.13 AddObjThisEvt(GenLeptons);
369     AddObjThisEvt(GenAllLeptons);
370     AddObjThisEvt(GenTaus);
371     AddObjThisEvt(GenNeutrinos);
372     AddObjThisEvt(GenQuarks);
373     AddObjThisEvt(GenqqHs);
374     AddObjThisEvt(GenBosons);
375     AddObjThisEvt(GenPhotons);
376 ceballos 1.32 AddObjThisEvt(GenRadPhotons);
377     AddObjThisEvt(GenISRPhotons);
378    
379 loizides 1.5 // fill histograms if requested
380 loizides 1.20 if (GetFillHist()) {
381 loizides 1.5
382 ceballos 1.30 // MET
383     hDGenMet[0]->Fill(GenMet->At(0)->Pt());
384     hDGenMet[1]->Fill(GenMet->At(0)->Px());
385     hDGenMet[2]->Fill(GenMet->At(0)->Py());
386     hDGenMet[3]->Fill(GenMet->At(0)->Elongitudinal());
387    
388 loizides 1.6 // leptons
389 loizides 1.1 hDGenLeptons[0]->Fill(GenLeptons->GetEntries());
390 loizides 1.5 for(UInt_t i=0; i<GenLeptons->GetEntries(); i++) {
391 loizides 1.1 hDGenLeptons[1]->Fill(GenLeptons->At(i)->Pt());
392 loizides 1.5 hDGenLeptons[2]->Fill(TMath::Abs(GenLeptons->At(i)->Eta()));
393     hDGenLeptons[3]->Fill(GenLeptons->At(i)->PhiDeg());
394     for(UInt_t j=i+1; j<GenLeptons->GetEntries(); j++) {
395 loizides 1.1 CompositeParticle *dilepton = new CompositeParticle();
396     dilepton->AddDaughter(GenLeptons->At(i));
397     dilepton->AddDaughter(GenLeptons->At(j));
398     hDGenLeptons[4]->Fill(dilepton->Mass());
399     delete dilepton;
400     }
401 ceballos 1.22 }
402     // looking at events with two leptons
403     if (GenLeptons->GetEntries() == 2) {
404     hDGenLeptons[5]->Fill(TMath::Min(TMath::Max(TMath::Abs(GenLeptons->At(0)->Eta()),
405     TMath::Abs(GenLeptons->At(1)->Eta())),
406 loizides 1.5 4.999));
407 ceballos 1.22 hDGenLeptons[6]->Fill(TMath::Min(TMath::Min(TMath::Abs(GenLeptons->At(0)->Eta()),
408     TMath::Abs(GenLeptons->At(1)->Eta())),
409 loizides 1.5 4.999));
410 ceballos 1.22 if (TMath::Abs(GenLeptons->At(0)->Eta()) < 2.5 &&
411     TMath::Abs(GenLeptons->At(1)->Eta()) < 2.5) {
412     hDGenLeptons[7]->Fill(TMath::Min(GenLeptons->At(0)->Pt(),199.999));
413     if (GenLeptons->At(0)->Pt() > 20.0) {
414     hDGenLeptons[8]->Fill(TMath::Min(GenLeptons->At(1)->Pt(),199.999));
415     if (GenLeptons->At(1)->Pt() > 10.0) {
416 loizides 1.1 CompositeParticle *dilepton = new CompositeParticle();
417 ceballos 1.22 dilepton->AddDaughter(GenLeptons->At(0));
418     dilepton->AddDaughter(GenLeptons->At(1));
419 loizides 1.1 hDGenLeptons[9]->Fill(TMath::Min(dilepton->Mass(),999.999));
420 ceballos 1.22 if(dilepton->Mass() > 12.0){
421     hDGenLeptons[10]->Fill(MathUtils::DeltaPhi(GenLeptons->At(0)->Phi(),
422     GenLeptons->At(1)->Phi())
423     * 180./ TMath::Pi());
424 loizides 1.40 hDGenLeptons[11]->Fill(MathUtils::DeltaR(*GenLeptons->At(0),
425     *GenLeptons->At(1)));
426 ceballos 1.22 }
427 loizides 1.1 delete dilepton;
428     }
429     }
430     }
431     }
432 ceballos 1.22 // looking at events with three leptons
433     if (GenLeptons->GetEntries() == 3) {
434     if (TMath::Abs(GenLeptons->At(0)->Eta()) < 2.5 &&
435     TMath::Abs(GenLeptons->At(1)->Eta()) < 2.5 &&
436     TMath::Abs(GenLeptons->At(2)->Eta()) < 2.5) {
437     hDGenLeptons[12]->Fill(TMath::Min(GenLeptons->At(0)->Pt(),199.999));
438     if (GenLeptons->At(0)->Pt() > 20.0) {
439     hDGenLeptons[13]->Fill(TMath::Min(GenLeptons->At(1)->Pt(),199.999));
440     hDGenLeptons[14]->Fill(TMath::Min(GenLeptons->At(2)->Pt(),199.999));
441     if (GenLeptons->At(1)->Pt() > 10.0 && GenLeptons->At(2)->Pt() > 10.0) {
442     CompositeParticle *dilepton01 = new CompositeParticle();
443     dilepton01->AddDaughter(GenLeptons->At(0));
444     dilepton01->AddDaughter(GenLeptons->At(1));
445     CompositeParticle *dilepton02 = new CompositeParticle();
446     dilepton02->AddDaughter(GenLeptons->At(0));
447     dilepton02->AddDaughter(GenLeptons->At(2));
448     CompositeParticle *dilepton12 = new CompositeParticle();
449     dilepton12->AddDaughter(GenLeptons->At(1));
450     dilepton12->AddDaughter(GenLeptons->At(2));
451     hDGenLeptons[15]->Fill(TMath::Min(dilepton01->Mass(),999.999));
452     hDGenLeptons[15]->Fill(TMath::Min(dilepton02->Mass(),999.999));
453     hDGenLeptons[15]->Fill(TMath::Min(dilepton12->Mass(),999.999));
454     CompositeParticle *trilepton = new CompositeParticle();
455     trilepton->AddDaughter(GenLeptons->At(0));
456     trilepton->AddDaughter(GenLeptons->At(1));
457     trilepton->AddDaughter(GenLeptons->At(2));
458     hDGenLeptons[16]->Fill(TMath::Min(trilepton->Mass(),999.999));
459 loizides 1.40 Double_t deltaR[3] = {MathUtils::DeltaR(*GenLeptons->At(0),
460     *GenLeptons->At(1)),
461     MathUtils::DeltaR(*GenLeptons->At(0),
462     *GenLeptons->At(2)),
463     MathUtils::DeltaR(*GenLeptons->At(1),
464     *GenLeptons->At(2))};
465 loizides 1.27 Double_t deltaRMin = deltaR[0];
466 loizides 1.40 for(Int_t i=1; i<3; i++)
467     if(deltaRMin > deltaR[i])
468     deltaRMin = deltaR[i];
469 ceballos 1.22 hDGenLeptons[17]->Fill(deltaRMin);
470    
471     delete dilepton01;
472     delete dilepton02;
473     delete dilepton12;
474     delete trilepton;
475     }
476     }
477     }
478     }
479     // looking at events with four leptons
480     if (GenLeptons->GetEntries() == 4) {
481     if (TMath::Abs(GenLeptons->At(0)->Eta()) < 2.5 &&
482     TMath::Abs(GenLeptons->At(1)->Eta()) < 2.5 &&
483     TMath::Abs(GenLeptons->At(2)->Eta()) < 2.5 &&
484     TMath::Abs(GenLeptons->At(3)->Eta()) < 2.5) {
485     hDGenLeptons[18]->Fill(TMath::Min(GenLeptons->At(0)->Pt(),199.999));
486     if (GenLeptons->At(0)->Pt() > 20.0) {
487     hDGenLeptons[19]->Fill(TMath::Min(GenLeptons->At(1)->Pt(),199.999));
488     hDGenLeptons[20]->Fill(TMath::Min(GenLeptons->At(2)->Pt(),199.999));
489     hDGenLeptons[21]->Fill(TMath::Min(GenLeptons->At(3)->Pt(),199.999));
490     if (GenLeptons->At(1)->Pt() > 10.0 && GenLeptons->At(2)->Pt() > 10.0 &&
491     GenLeptons->At(3)->Pt() > 10.0) {
492     CompositeParticle *dilepton01 = new CompositeParticle();
493     dilepton01->AddDaughter(GenLeptons->At(0));
494     dilepton01->AddDaughter(GenLeptons->At(1));
495     CompositeParticle *dilepton02 = new CompositeParticle();
496     dilepton02->AddDaughter(GenLeptons->At(0));
497     dilepton02->AddDaughter(GenLeptons->At(2));
498     CompositeParticle *dilepton03 = new CompositeParticle();
499     dilepton03->AddDaughter(GenLeptons->At(0));
500     dilepton03->AddDaughter(GenLeptons->At(3));
501     CompositeParticle *dilepton12 = new CompositeParticle();
502     dilepton12->AddDaughter(GenLeptons->At(1));
503     dilepton12->AddDaughter(GenLeptons->At(2));
504     CompositeParticle *dilepton13 = new CompositeParticle();
505     dilepton13->AddDaughter(GenLeptons->At(1));
506     dilepton13->AddDaughter(GenLeptons->At(3));
507     CompositeParticle *dilepton23 = new CompositeParticle();
508     dilepton23->AddDaughter(GenLeptons->At(2));
509     dilepton23->AddDaughter(GenLeptons->At(3));
510     hDGenLeptons[22]->Fill(TMath::Min(dilepton01->Mass(),999.999));
511     hDGenLeptons[22]->Fill(TMath::Min(dilepton02->Mass(),999.999));
512     hDGenLeptons[22]->Fill(TMath::Min(dilepton03->Mass(),999.999));
513     hDGenLeptons[22]->Fill(TMath::Min(dilepton12->Mass(),999.999));
514     hDGenLeptons[22]->Fill(TMath::Min(dilepton13->Mass(),999.999));
515     hDGenLeptons[22]->Fill(TMath::Min(dilepton23->Mass(),999.999));
516     CompositeParticle *fourlepton = new CompositeParticle();
517     fourlepton->AddDaughter(GenLeptons->At(0));
518     fourlepton->AddDaughter(GenLeptons->At(1));
519     fourlepton->AddDaughter(GenLeptons->At(2));
520     fourlepton->AddDaughter(GenLeptons->At(3));
521     hDGenLeptons[23]->Fill(TMath::Min(fourlepton->Mass(),999.999));
522 loizides 1.40 Double_t deltaR[6] = {MathUtils::DeltaR(*GenLeptons->At(0),
523     *GenLeptons->At(1)),
524     MathUtils::DeltaR(*GenLeptons->At(0),
525     *GenLeptons->At(2)),
526     MathUtils::DeltaR(*GenLeptons->At(0),
527     *GenLeptons->At(3)),
528     MathUtils::DeltaR(*GenLeptons->At(1),
529     *GenLeptons->At(2)),
530     MathUtils::DeltaR(*GenLeptons->At(1),
531     *GenLeptons->At(3)),
532     MathUtils::DeltaR(*GenLeptons->At(2),
533     *GenLeptons->At(3))};
534 loizides 1.27 Double_t deltaRMin = deltaR[0];
535 loizides 1.40 for(Int_t i=1; i<6; i++)
536     if(deltaRMin > deltaR[i])
537     deltaRMin = deltaR[i];
538 ceballos 1.22 hDGenLeptons[24]->Fill(deltaRMin);
539    
540     delete dilepton01;
541     delete dilepton02;
542     delete dilepton03;
543     delete dilepton12;
544     delete dilepton13;
545     delete dilepton23;
546     delete fourlepton;
547     }
548     }
549     }
550     }
551 loizides 1.1
552 loizides 1.6 // all leptons
553 ceballos 1.3 hDGenAllLeptons[0]->Fill(GenAllLeptons->GetEntries());
554 loizides 1.5 for(UInt_t i=0; i<GenAllLeptons->GetEntries(); i++) {
555 ceballos 1.3 hDGenAllLeptons[1]->Fill(GenAllLeptons->At(i)->Pt());
556     hDGenAllLeptons[2]->Fill(GenAllLeptons->At(i)->Eta());
557 loizides 1.5 hDGenAllLeptons[3]->Fill(GenAllLeptons->At(i)->PhiDeg());
558 ceballos 1.3 }
559    
560 loizides 1.6 // taus
561 loizides 1.1 hDGenTaus[0]->Fill(GenTaus->GetEntries());
562 loizides 1.5 for(UInt_t i=0; i<GenTaus->GetEntries(); i++) {
563 loizides 1.1 hDGenTaus[1]->Fill(GenTaus->At(i)->Pt());
564     hDGenTaus[2]->Fill(GenTaus->At(i)->Eta());
565 loizides 1.5 hDGenTaus[3]->Fill(GenTaus->At(i)->PhiDeg());
566 loizides 1.1 }
567    
568 loizides 1.6 // neutrinos
569 loizides 1.1 hDGenNeutrinos[0]->Fill(GenNeutrinos->GetEntries());
570     CompositeParticle *neutrinoTotal = new CompositeParticle();
571 loizides 1.5 for(UInt_t i=0; i<GenNeutrinos->GetEntries(); i++) {
572     if (GenNeutrinos->At(i)->HasMother())
573 loizides 1.1 neutrinoTotal->AddDaughter(GenNeutrinos->At(i));
574     }
575 loizides 1.5 if (GenNeutrinos->GetEntries() > 0) {
576 loizides 1.1 hDGenNeutrinos[1]->Fill(neutrinoTotal->Pt());
577     hDGenNeutrinos[2]->Fill(neutrinoTotal->Eta());
578 loizides 1.5 hDGenNeutrinos[3]->Fill(neutrinoTotal->PhiDeg());
579 loizides 1.1 }
580     delete neutrinoTotal;
581    
582 loizides 1.6 // quarks
583 loizides 1.1 hDGenQuarks[0]->Fill(GenQuarks->GetEntries());
584 loizides 1.5 for(UInt_t i=0; i<GenQuarks->GetEntries(); i++) {
585     for(UInt_t j=i+1; j<GenQuarks->GetEntries(); j++) {
586 loizides 1.1 CompositeParticle *dijet = new CompositeParticle();
587     dijet->AddDaughter(GenQuarks->At(i));
588     dijet->AddDaughter(GenQuarks->At(j));
589     hDGenQuarks[1]->Fill(dijet->Pt());
590     hDGenQuarks[2]->Fill(dijet->Mass());
591 loizides 1.5 if (TMath::Abs(GenQuarks->At(i)->Eta()) < 2.5 &&
592     TMath::Abs(GenQuarks->At(j)->Eta()) < 2.5) {
593 loizides 1.1 hDGenQuarks[3]->Fill(dijet->Pt());
594     hDGenQuarks[4]->Fill(dijet->Mass());
595     }
596     delete dijet;
597     }
598 ceballos 1.2 // b quark info
599 loizides 1.5 if (GenQuarks->At(i)->AbsPdgId() == 5) {
600 ceballos 1.2 hDGenQuarks[5]->Fill(GenQuarks->At(i)->Pt());
601     hDGenQuarks[6]->Fill(GenQuarks->At(i)->Eta());
602     hDGenQuarks[7]->Fill(GenQuarks->At(i)->Phi());
603 ceballos 1.22 if (GenLeptons->GetEntries() >= 2 &&
604     GenLeptons->At(0)->Pt() > 20 &&
605     GenLeptons->At(1)->Pt() > 15) {
606     if (TMath::Abs(GenLeptons->At(0)->Eta()) < 2.5 &&
607     TMath::Abs(GenLeptons->At(1)->Eta()) < 2.5) {
608 ceballos 1.2 hDGenQuarks[8]->Fill(GenQuarks->At(i)->Pt());
609     hDGenQuarks[9]->Fill(GenQuarks->At(i)->Eta());
610     hDGenQuarks[10]->Fill(GenQuarks->At(i)->Phi());
611     }
612     }
613     }
614     // t quark info
615 loizides 1.5 else if (GenQuarks->At(i)->AbsPdgId() == 6) {
616 ceballos 1.2 hDGenQuarks[11]->Fill(GenQuarks->At(i)->Pt());
617     hDGenQuarks[12]->Fill(GenQuarks->At(i)->Eta());
618     hDGenQuarks[13]->Fill(GenQuarks->At(i)->Phi());
619     }
620     // light quark info
621     else {
622     hDGenQuarks[14]->Fill(GenQuarks->At(i)->Pt());
623     hDGenQuarks[15]->Fill(GenQuarks->At(i)->Eta());
624     hDGenQuarks[16]->Fill(GenQuarks->At(i)->Phi());
625     }
626 loizides 1.1 }
627    
628 loizides 1.6 // wbf
629 loizides 1.5 if (GenqqHs->GetEntries() == 2) {
630 loizides 1.1 hDGenWBF[0]->Fill(MathUtils::DeltaPhi(GenqqHs->At(0)->Phi(),
631     GenqqHs->At(1)->Phi()) * 180./ TMath::Pi());
632 loizides 1.5 hDGenWBF[1]->Fill(TMath::Abs(GenqqHs->At(0)->Eta()-GenqqHs->At(1)->Eta()));
633 loizides 1.1 hDGenWBF[2]->Fill(TMath::Max(GenqqHs->At(0)->Pt(),GenqqHs->At(1)->Pt()));
634     hDGenWBF[3]->Fill(TMath::Min(GenqqHs->At(0)->Pt(),GenqqHs->At(1)->Pt()));
635     CompositeParticle *diqq = new CompositeParticle();
636     diqq->AddDaughter(GenqqHs->At(0));
637     diqq->AddDaughter(GenqqHs->At(1));
638     hDGenWBF[4]->Fill(diqq->Mass());
639     delete diqq;
640     }
641    
642 loizides 1.6 // bosons
643 loizides 1.1 hDGenBosons[0]->Fill(GenBosons->GetEntries());
644 loizides 1.5 for(UInt_t i=0; i<GenBosons->GetEntries(); i++) {
645 loizides 1.1 hDGenBosons[1]->Fill(GenBosons->At(i)->Pt());
646     hDGenBosons[2]->Fill(GenBosons->At(i)->Eta());
647     hDGenBosons[3]->Fill(GenBosons->At(i)->Mass());
648     }
649 ceballos 1.10
650     // photons
651     hDGenPhotons[0]->Fill(GenPhotons->GetEntries());
652     for(UInt_t i=0; i<GenPhotons->GetEntries(); i++) {
653     hDGenPhotons[1]->Fill(GenPhotons->At(i)->Pt());
654     hDGenPhotons[2]->Fill(GenPhotons->At(i)->Eta());
655     }
656 ceballos 1.32
657     // Rad photons
658     hDGenRadPhotons[0]->Fill(GenRadPhotons->GetEntries());
659     for(UInt_t i=0; i<GenRadPhotons->GetEntries(); i++) {
660     hDGenRadPhotons[1]->Fill(TMath::Min(GenRadPhotons->At(i)->Pt(),199.999));
661     hDGenRadPhotons[2]->Fill(TMath::Min(GenRadPhotons->At(i)->AbsEta(),4.999));
662     hDGenRadPhotons[3]->Fill(TMath::Min((double)GenRadPhotons->At(i)->Mother()->Status(),19.499));
663 loizides 1.40 hDGenRadPhotons[4]->Fill(GenRadPhotons->At(i)->IsGenerated() +
664     2*GenRadPhotons->At(i)->IsSimulated());
665 ceballos 1.32 hDGenRadPhotons[5]->Fill(TMath::Min(
666 loizides 1.40 MathUtils::DeltaR(*GenRadPhotons->At(i),
667     *GenRadPhotons->At(i)->Mother()),
668     4.999));
669 ceballos 1.32 Int_t Mother = 0;
670     if(GenRadPhotons->At(i)->Mother()->Is(MCParticle::kMu)) Mother = 1;
671     hDGenRadPhotons[6]->Fill(Mother);
672     }
673    
674     // ISR photons
675     hDGenISRPhotons[0]->Fill(GenISRPhotons->GetEntries());
676     for(UInt_t i=0; i<GenISRPhotons->GetEntries(); i++) {
677     hDGenISRPhotons[1]->Fill(TMath::Min(GenISRPhotons->At(i)->Pt(),199.999));
678     hDGenISRPhotons[2]->Fill(TMath::Min(GenISRPhotons->At(i)->AbsEta(),4.999));
679 loizides 1.40 hDGenISRPhotons[3]->Fill(TMath::Min((Double_t)GenISRPhotons->At(i)->Mother()->Status(),
680     19.499));
681     hDGenISRPhotons[4]->Fill(GenISRPhotons->At(i)->IsGenerated() +
682     2*GenISRPhotons->At(i)->IsSimulated());
683 ceballos 1.32 hDGenISRPhotons[5]->Fill(TMath::Min(
684 loizides 1.40 MathUtils::DeltaR(*GenISRPhotons->At(i),
685     *GenISRPhotons->At(i)->Mother()),4.999));
686 ceballos 1.32 }
687 loizides 1.1 }
688 ceballos 1.34
689     // Apply ISR filter (but filling all histograms)
690 ceballos 1.37 if(fApplyISRFilter == kTRUE && GenISRPhotons->GetEntries() > 0 &&
691 ceballos 1.38 GenISRPhotons->At(0)->Pt() > 15.0){
692 ceballos 1.34 SkipEvent();
693     }
694 loizides 1.1 }
695    
696     //--------------------------------------------------------------------------------------------------
697     void GeneratorMod::SlaveBegin()
698     {
699 loizides 1.5 // Book branch and histograms if wanted.
700    
701 loizides 1.41 ReqEventObject(fMCPartName, fParticles, kTRUE);
702 loizides 1.1
703 loizides 1.5 // fill histograms
704 loizides 1.20 if (GetFillHist()) {
705 loizides 1.5 char sb[1024];
706 ceballos 1.30 // MET
707     sprintf(sb,"hDGenMet_%d", 0); hDGenMet[0] = new TH1D(sb,sb,200,0,200);
708     sprintf(sb,"hDGenMet_%d", 1); hDGenMet[1] = new TH1D(sb,sb,400,-200,200);
709     sprintf(sb,"hDGenMet_%d", 2); hDGenMet[2] = new TH1D(sb,sb,400,-200,200);
710     sprintf(sb,"hDGenMet_%d", 3); hDGenMet[3] = new TH1D(sb,sb,400,-1000,1000);
711     for(Int_t i=0; i<4; i++) AddOutput(hDGenMet[i]);
712    
713 loizides 1.6 // leptons from W
714 loizides 1.1 sprintf(sb,"hDGenLeptons_%d", 0); hDGenLeptons[0] = new TH1D(sb,sb,10,-0.5,9.5);
715     sprintf(sb,"hDGenLeptons_%d", 1); hDGenLeptons[1] = new TH1D(sb,sb,100,0.0,200.0);
716     sprintf(sb,"hDGenLeptons_%d", 2); hDGenLeptons[2] = new TH1D(sb,sb,50,0.0,5.0);
717     sprintf(sb,"hDGenLeptons_%d", 3); hDGenLeptons[3] = new TH1D(sb,sb,90,0.0,180.0);
718     sprintf(sb,"hDGenLeptons_%d", 4); hDGenLeptons[4] = new TH1D(sb,sb,1000,0.0,1000.0);
719     sprintf(sb,"hDGenLeptons_%d", 5); hDGenLeptons[5] = new TH1D(sb,sb,50,0.0,5.0);
720     sprintf(sb,"hDGenLeptons_%d", 6); hDGenLeptons[6] = new TH1D(sb,sb,50,0.0,5.0);
721     sprintf(sb,"hDGenLeptons_%d", 7); hDGenLeptons[7] = new TH1D(sb,sb,100,0.0,200.0);
722     sprintf(sb,"hDGenLeptons_%d", 8); hDGenLeptons[8] = new TH1D(sb,sb,100,0.0,200.0);
723     sprintf(sb,"hDGenLeptons_%d", 9); hDGenLeptons[9] = new TH1D(sb,sb,1000,0.0,1000.0);
724     sprintf(sb,"hDGenLeptons_%d",10); hDGenLeptons[10] = new TH1D(sb,sb,90,0.0,180.0);
725 ceballos 1.22 sprintf(sb,"hDGenLeptons_%d",11); hDGenLeptons[11] = new TH1D(sb,sb,100,0.0,5.0);
726     sprintf(sb,"hDGenLeptons_%d",12); hDGenLeptons[12] = new TH1D(sb,sb,100,0.0,200.0);
727     sprintf(sb,"hDGenLeptons_%d",13); hDGenLeptons[13] = new TH1D(sb,sb,100,0.0,200.0);
728     sprintf(sb,"hDGenLeptons_%d",14); hDGenLeptons[14] = new TH1D(sb,sb,100,0.0,200.0);
729     sprintf(sb,"hDGenLeptons_%d",15); hDGenLeptons[15] = new TH1D(sb,sb,1000,0.0,1000.0);
730     sprintf(sb,"hDGenLeptons_%d",16); hDGenLeptons[16] = new TH1D(sb,sb,1000,0.0,1000.0);
731     sprintf(sb,"hDGenLeptons_%d",17); hDGenLeptons[17] = new TH1D(sb,sb,100,0.0,5.0);
732     sprintf(sb,"hDGenLeptons_%d",18); hDGenLeptons[18] = new TH1D(sb,sb,100,0.0,200.0);
733     sprintf(sb,"hDGenLeptons_%d",19); hDGenLeptons[19] = new TH1D(sb,sb,100,0.0,200.0);
734     sprintf(sb,"hDGenLeptons_%d",20); hDGenLeptons[20] = new TH1D(sb,sb,100,0.0,200.0);
735     sprintf(sb,"hDGenLeptons_%d",21); hDGenLeptons[21] = new TH1D(sb,sb,100,0.0,200.0);
736     sprintf(sb,"hDGenLeptons_%d",22); hDGenLeptons[22] = new TH1D(sb,sb,1000,0.0,1000.0);
737     sprintf(sb,"hDGenLeptons_%d",23); hDGenLeptons[23] = new TH1D(sb,sb,1000,0.0,1000.0);
738     sprintf(sb,"hDGenLeptons_%d",24); hDGenLeptons[24] = new TH1D(sb,sb,100,0.0,5.0);
739     for(Int_t i=0; i<25; i++) AddOutput(hDGenLeptons[i]);
740 loizides 1.1
741 loizides 1.6 // all leptons
742 ceballos 1.3 sprintf(sb,"hDGenAllLeptons_%d", 0); hDGenAllLeptons[0] = new TH1D(sb,sb,10,-0.5,9.5);
743     sprintf(sb,"hDGenAllLeptons_%d", 1); hDGenAllLeptons[1] = new TH1D(sb,sb,100,0.0,200.0);
744 ceballos 1.14 sprintf(sb,"hDGenAllLeptons_%d", 2); hDGenAllLeptons[2] = new TH1D(sb,sb,100,-5.0,5.0);
745 ceballos 1.3 sprintf(sb,"hDGenAllLeptons_%d", 3); hDGenAllLeptons[3] = new TH1D(sb,sb,90,0.0,180.0);
746 loizides 1.5 for(Int_t i=0; i<4; i++) AddOutput(hDGenAllLeptons[i]);
747 ceballos 1.3
748 loizides 1.6 // taus
749 loizides 1.1 sprintf(sb,"hDGenTaus_%d", 0); hDGenTaus[0] = new TH1D(sb,sb,10,-0.5,9.5);
750     sprintf(sb,"hDGenTaus_%d", 1); hDGenTaus[1] = new TH1D(sb,sb,100,0.0,200.0);
751     sprintf(sb,"hDGenTaus_%d", 2); hDGenTaus[2] = new TH1D(sb,sb,100,-5.0,5.0);
752     sprintf(sb,"hDGenTaus_%d", 3); hDGenTaus[3] = new TH1D(sb,sb,90,0.0,180.0);
753 loizides 1.5 for(Int_t i=0; i<4; i++) AddOutput(hDGenTaus[i]);
754 loizides 1.1
755 loizides 1.6 // neutrinos
756 loizides 1.1 sprintf(sb,"hDGenNeutrinos_%d", 0); hDGenNeutrinos[0] = new TH1D(sb,sb,10,-0.5,9.5);
757     sprintf(sb,"hDGenNeutrinos_%d", 1); hDGenNeutrinos[1] = new TH1D(sb,sb,100,0.0,200.0);
758     sprintf(sb,"hDGenNeutrinos_%d", 2); hDGenNeutrinos[2] = new TH1D(sb,sb,100,-5.0,5.0);
759     sprintf(sb,"hDGenNeutrinos_%d", 3); hDGenNeutrinos[3] = new TH1D(sb,sb,90,0.0,180.0);
760 loizides 1.5 for(Int_t i=0; i<4; i++) AddOutput(hDGenNeutrinos[i]);
761 loizides 1.1
762 loizides 1.6 // quarks
763 loizides 1.1 sprintf(sb,"hDGenQuarks_%d", 0); hDGenQuarks[0] = new TH1D(sb,sb,10,-0.5,9.5);
764     sprintf(sb,"hDGenQuarks_%d", 1); hDGenQuarks[1] = new TH1D(sb,sb,200,0.0,400.);
765     sprintf(sb,"hDGenQuarks_%d", 2); hDGenQuarks[2] = new TH1D(sb,sb,2000,0.0,2000.);
766     sprintf(sb,"hDGenQuarks_%d", 3); hDGenQuarks[3] = new TH1D(sb,sb,200,0.0,400.);
767     sprintf(sb,"hDGenQuarks_%d", 4); hDGenQuarks[4] = new TH1D(sb,sb,2000,0.0,2000.);
768 ceballos 1.2 sprintf(sb,"hDGenQuarks_%d", 5); hDGenQuarks[5] = new TH1D(sb,sb,200,0.0,400.);
769     sprintf(sb,"hDGenQuarks_%d", 6); hDGenQuarks[6] = new TH1D(sb,sb,200,-10.0,10.);
770 loizides 1.5 sprintf(sb,"hDGenQuarks_%d", 7); hDGenQuarks[7] = new TH1D(sb,sb,200,-TMath::Pi(),
771     TMath::Pi());
772 ceballos 1.2 sprintf(sb,"hDGenQuarks_%d", 8); hDGenQuarks[8] = new TH1D(sb,sb,200,0.0,400.);
773     sprintf(sb,"hDGenQuarks_%d", 9); hDGenQuarks[9] = new TH1D(sb,sb,200,-10.0,10.);
774 loizides 1.5 sprintf(sb,"hDGenQuarks_%d",10); hDGenQuarks[10] = new TH1D(sb,sb,200,-TMath::Pi(),
775     TMath::Pi());
776 ceballos 1.2 sprintf(sb,"hDGenQuarks_%d",11); hDGenQuarks[11] = new TH1D(sb,sb,200,0.0,400.);
777     sprintf(sb,"hDGenQuarks_%d",12); hDGenQuarks[12] = new TH1D(sb,sb,200,-10.0,10.);
778 loizides 1.5 sprintf(sb,"hDGenQuarks_%d",13); hDGenQuarks[13] = new TH1D(sb,sb,200,-TMath::Pi(),
779     TMath::Pi());
780 ceballos 1.2 sprintf(sb,"hDGenQuarks_%d",14); hDGenQuarks[14] = new TH1D(sb,sb,200,0.0,400.);
781     sprintf(sb,"hDGenQuarks_%d",15); hDGenQuarks[15] = new TH1D(sb,sb,200,-10.0,10.);
782 loizides 1.5 sprintf(sb,"hDGenQuarks_%d",16); hDGenQuarks[16] = new TH1D(sb,sb,200,-TMath::Pi(),
783     TMath::Pi());
784     for(Int_t i=0; i<17; i++) AddOutput(hDGenQuarks[i]);
785 loizides 1.1
786     // qqH
787     sprintf(sb,"hDGenWBF_%d", 0); hDGenWBF[0] = new TH1D(sb,sb,90,0.0,180.);
788     sprintf(sb,"hDGenWBF_%d", 1); hDGenWBF[1] = new TH1D(sb,sb,100,0.0,10.);
789     sprintf(sb,"hDGenWBF_%d", 2); hDGenWBF[2] = new TH1D(sb,sb,200,0.0,400.);
790     sprintf(sb,"hDGenWBF_%d", 3); hDGenWBF[3] = new TH1D(sb,sb,200,0.0,400.);
791     sprintf(sb,"hDGenWBF_%d", 4); hDGenWBF[4] = new TH1D(sb,sb,200,0.0,4000.);
792 loizides 1.5 for(Int_t i=0; i<5; i++) AddOutput(hDGenWBF[i]);
793 loizides 1.1
794 loizides 1.6 // bosons
795 loizides 1.1 sprintf(sb,"hDGenBosons_%d", 0); hDGenBosons[0] = new TH1D(sb,sb,10,-0.5,9.5);
796     sprintf(sb,"hDGenBosons_%d", 1); hDGenBosons[1] = new TH1D(sb,sb,200,0.0,400.0);
797     sprintf(sb,"hDGenBosons_%d", 2); hDGenBosons[2] = new TH1D(sb,sb,100,-5.0,5.0);
798     sprintf(sb,"hDGenBosons_%d", 3); hDGenBosons[3] = new TH1D(sb,sb,2000,0.0,2000.0);
799 loizides 1.5 for(Int_t i=0; i<4; i++) AddOutput(hDGenBosons[i]);
800 ceballos 1.10
801     // photons
802     sprintf(sb,"hDGenPhotons_%d", 0); hDGenPhotons[0] = new TH1D(sb,sb,10,-0.5,9.5);
803     sprintf(sb,"hDGenPhotons_%d", 1); hDGenPhotons[1] = new TH1D(sb,sb,200,0.0,400.0);
804 ceballos 1.14 sprintf(sb,"hDGenPhotons_%d", 2); hDGenPhotons[2] = new TH1D(sb,sb,100,-5.0,5.0);
805 ceballos 1.10 for(Int_t i=0; i<3; i++) AddOutput(hDGenPhotons[i]);
806 ceballos 1.16
807 loizides 1.40 // rad photons
808 ceballos 1.32 sprintf(sb,"hDGenRadPhotons_%d", 0); hDGenRadPhotons[0] = new TH1D(sb,sb,10,-0.5,9.5);
809     sprintf(sb,"hDGenRadPhotons_%d", 1); hDGenRadPhotons[1] = new TH1D(sb,sb,400,0.0,200.0);
810     sprintf(sb,"hDGenRadPhotons_%d", 2); hDGenRadPhotons[2] = new TH1D(sb,sb,100,0.0,5.0);
811     sprintf(sb,"hDGenRadPhotons_%d", 3); hDGenRadPhotons[3] = new TH1D(sb,sb,20,-0.5,19.5);
812     sprintf(sb,"hDGenRadPhotons_%d", 4); hDGenRadPhotons[4] = new TH1D(sb,sb,4,-0.5,3.5);
813     sprintf(sb,"hDGenRadPhotons_%d", 5); hDGenRadPhotons[5] = new TH1D(sb,sb,500,0.0,5.0);
814     sprintf(sb,"hDGenRadPhotons_%d", 6); hDGenRadPhotons[6] = new TH1D(sb,sb,2,-0.5,1.5);
815     for(Int_t i=0; i<7; i++) AddOutput(hDGenRadPhotons[i]);
816    
817     // ISR photons
818     sprintf(sb,"hDGenISRPhotons_%d", 0); hDGenISRPhotons[0] = new TH1D(sb,sb,10,-0.5,9.5);
819     sprintf(sb,"hDGenISRPhotons_%d", 1); hDGenISRPhotons[1] = new TH1D(sb,sb,400,0.0,200.0);
820     sprintf(sb,"hDGenISRPhotons_%d", 2); hDGenISRPhotons[2] = new TH1D(sb,sb,100,0.0,5.0);
821     sprintf(sb,"hDGenISRPhotons_%d", 3); hDGenISRPhotons[3] = new TH1D(sb,sb,20,-0.5,19.5);
822     sprintf(sb,"hDGenISRPhotons_%d", 4); hDGenISRPhotons[4] = new TH1D(sb,sb,4,-0.5,3.5);
823     sprintf(sb,"hDGenISRPhotons_%d", 5); hDGenISRPhotons[5] = new TH1D(sb,sb,500,0.0,5.0);
824     for(Int_t i=0; i<6; i++) AddOutput(hDGenISRPhotons[i]);
825    
826 ceballos 1.16 // auxiliar
827     sprintf(sb,"hDVMass_%d", 0); hDVMass[0] = new TH1D(sb,sb,200,0.,200.);
828     sprintf(sb,"hDVMass_%d", 1); hDVMass[1] = new TH1D(sb,sb,200,0.,200.);
829     sprintf(sb,"hDVMass_%d", 2); hDVMass[2] = new TH1D(sb,sb,200,0.,200.);
830     sprintf(sb,"hDVMass_%d", 3); hDVMass[3] = new TH1D(sb,sb,200,0.,200.);
831     sprintf(sb,"hDVMass_%d", 4); hDVMass[4] = new TH1D(sb,sb,200,0.,200.);
832     sprintf(sb,"hDVMass_%d", 5); hDVMass[5] = new TH1D(sb,sb,200,0.,200.);
833     sprintf(sb,"hDVMass_%d", 6); hDVMass[6] = new TH1D(sb,sb,200,0.,200.);
834     sprintf(sb,"hDVMass_%d", 7); hDVMass[7] = new TH1D(sb,sb,200,0.,200.);
835     sprintf(sb,"hDVMass_%d", 8); hDVMass[8] = new TH1D(sb,sb,200,0.,200.);
836     sprintf(sb,"hDVMass_%d", 9); hDVMass[9] = new TH1D(sb,sb,200,0.,200.);
837     for(Int_t i=0; i<10; i++) AddOutput(hDVMass[i]);
838 loizides 1.1 }
839     }