ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/MitPhysics/Mods/src/GeneratorMod.cc
Revision: 1.16
Committed: Thu Dec 4 11:55:19 2008 UTC (16 years, 5 months ago) by ceballos
Content type: text/plain
Branch: MAIN
Changes since 1.15: +87 -2 lines
Log Message: 
adding leptons from W/Z boson with MG

File Contents

# User Rev Content
1 ceballos 1.16 // $Id: GeneratorMod.cc,v 1.15 2008/12/03 10:19:13 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 loizides 1.5 fFillHist(kFALSE),
17 loizides 1.1 fMCPartName(Names::gkMCPartBrn),
18 loizides 1.8 fMCLeptonsName(ModNames::gkMCLeptonsName),
19     fMCAllLeptonsName(ModNames::gkMCAllLeptonsName),
20     fMCTausName(ModNames::gkMCTausName),
21     fMCNeutrinosName(ModNames::gkMCNeutrinosName),
22     fMCQuarksName(ModNames::gkMCQuarksName),
23     fMCqqHsName(ModNames::gkMCqqHsName),
24     fMCBosonsName(ModNames::gkMCBosonsName),
25 ceballos 1.10 fMCPhotonsName(ModNames::gkMCPhotonsName),
26 ceballos 1.14 fPtLeptonMin(0.0),
27     fEtaLeptonMax(5.0),
28     fPtPhotonMin(0.0),
29 loizides 1.15 fEtaPhotonMax(5.0),
30     fParticles(0)
31 loizides 1.1 {
32     // Constructor.
33     }
34    
35     //--------------------------------------------------------------------------------------------------
36     void GeneratorMod::Process()
37     {
38 loizides 1.5 // Process entries of the tree.
39 loizides 1.1
40 loizides 1.5 // these arrays will be filled in the loop of particles
41 loizides 1.9 MCParticleOArr *GenLeptons = new MCParticleOArr;
42 loizides 1.13 GenLeptons->SetName(fMCLeptonsName);
43 loizides 1.9 MCParticleOArr *GenAllLeptons = new MCParticleOArr;
44 loizides 1.13 GenAllLeptons->SetName(fMCAllLeptonsName);
45 loizides 1.9 MCParticleOArr *GenTaus = new MCParticleOArr;
46 loizides 1.13 GenTaus->SetName(fMCTausName);
47     GenTaus->SetOwner(kTRUE);
48 loizides 1.9 MCParticleOArr *GenNeutrinos = new MCParticleOArr;
49 loizides 1.13 GenNeutrinos->SetName(fMCNeutrinosName);
50 loizides 1.9 MCParticleOArr *GenQuarks = new MCParticleOArr;
51 loizides 1.13 GenQuarks->SetName(fMCQuarksName);
52 loizides 1.9 MCParticleOArr *GenqqHs = new MCParticleOArr;
53 loizides 1.13 GenqqHs->SetName(fMCqqHsName);
54 loizides 1.9 MCParticleOArr *GenBosons = new MCParticleOArr;
55 loizides 1.13 GenBosons->SetName(fMCBosonsName);
56 ceballos 1.10 MCParticleOArr *GenPhotons = new MCParticleOArr;
57 loizides 1.13 GenPhotons->SetName(fMCPhotonsName);
58 loizides 1.1
59 loizides 1.5 // load MCParticle branch
60     LoadBranch(fMCPartName);
61    
62 loizides 1.6 Bool_t isqqH = kFALSE;
63 loizides 1.5 for (UInt_t i=0; i<fParticles->GetEntries(); ++i) {
64 loizides 1.7 const MCParticle *p = fParticles->At(i);
65 loizides 1.5
66     if (!p->IsGenerated()) continue;
67    
68     // muons/electrons from W/Z decays
69 loizides 1.7 if ((p->Is(MCParticle::kEl) || p->Is(MCParticle::kMu)) && p->Status() == 1) {
70 ceballos 1.14 if (p->Pt() > fPtLeptonMin && p->AbsEta() < fEtaLeptonMax) {
71 loizides 1.5 GenAllLeptons->Add(p);
72     }
73 loizides 1.6 Bool_t isGoodLepton = kFALSE;
74 loizides 1.7 const MCParticle *pm = p;
75 loizides 1.5 while (pm->HasMother() && isGoodLepton == kFALSE) {
76 loizides 1.7 if (pm->Mother()->Is(MCParticle::kZ) || pm->Mother()->Is(MCParticle::kW)) {
77 loizides 1.5 GenLeptons->Add(p);
78     isGoodLepton = kTRUE;
79 loizides 1.7 break;
80 loizides 1.8 } else if (pm->Mother()->Is(MCParticle::kPi0) || pm->Mother()->Is(MCParticle::kEta)) {
81 loizides 1.7 // this is fake, but it is a trick to get rid of these cases and abort the loop
82 loizides 1.5 isGoodLepton = kTRUE;
83 loizides 1.7 break;
84     }
85     pm = pm->Mother();
86 loizides 1.1 }
87 loizides 1.5 }
88    
89 loizides 1.7 // hadronic taus
90     else if (p->Is(MCParticle::kTau) && p->Status() == 2) {
91     if (!p->HasDaughter(MCParticle::kEl) && !p->HasDaughter(MCParticle::kMu)) {
92     const MCParticle *tv = p->FindDaughter(MCParticle::kTauNu);
93     if (tv) {
94     MCParticle *pm_f = new MCParticle(*p);
95     pm_f->SetMom(p->Px()-tv->Px(), p->Py()-tv->Py(),
96     p->Pz()-tv->Pz(), p->E()-tv->E());
97     GenTaus->AddOwned(pm_f);
98     } else {
99     SendError(kWarning, "Process", "Could not find tau neutrino!");
100 loizides 1.5 }
101 loizides 1.1 }
102 loizides 1.5 }
103 loizides 1.1
104 loizides 1.5 // neutrinos
105 loizides 1.7 else if (p->Status() == 1 && p->IsNeutrino()) {
106 loizides 1.5 GenNeutrinos->Add(p);
107     }
108 loizides 1.1
109 loizides 1.5 // quarks from W/Z decays or top particles
110 loizides 1.7 else if (p->IsQuark() && p->HasMother()) {
111     if (p->Mother()->Is(MCParticle::kZ) || p->Mother()->Is(MCParticle::kW) ||
112     p->Is(MCParticle::kTop) || p->Mother()->Is(MCParticle::kTop)) {
113 loizides 1.5 GenQuarks->Add(p);
114 loizides 1.1 }
115 loizides 1.5 }
116 loizides 1.1
117 loizides 1.5 // qqH, information about the forward jets
118 loizides 1.8 else if (isqqH == kFALSE && p->Is(MCParticle::kH)) {
119 loizides 1.5 isqqH = kTRUE;
120     MCParticle *pq1 = fParticles->At(i-1);
121     MCParticle *pq2 = fParticles->At(i-2);
122 loizides 1.7 if (!pq1 || !pq2) {
123     SendError(kWarning, "Process", "Could not find quark pair!");
124 loizides 1.8 } else if (pq1->IsQuark() && pq2->IsQuark() &&
125     pq1->HasMother() && pq2->HasMother() &&
126     pq1->Mother()->PdgId() == p->Mother()->PdgId() &&
127     pq2->Mother()->PdgId() == p->Mother()->PdgId()) {
128 loizides 1.5 GenqqHs->Add(pq1);
129     GenqqHs->Add(pq2);
130 loizides 1.1 }
131 loizides 1.7 if (p->Status() == 3)
132     GenBosons->Add(p); // take higgs boson in account here rather in next else if
133 loizides 1.5 }
134 loizides 1.1
135 loizides 1.5 // information about bosons: W, Z, h, Z', W', H0, A0, H+
136 loizides 1.7 else if (p->Status() == 3 &&
137     (p->Is(MCParticle::kZ) || p->Is(MCParticle::kW) || p->Is(MCParticle::kH) ||
138     p->Is(MCParticle::kZp) || p->Is(MCParticle::kZpp) ||
139     p->Is(MCParticle::kH0) || p->Is(MCParticle::kA0) || p->Is(MCParticle::kHp))) {
140 loizides 1.5 GenBosons->Add(p);
141 loizides 1.1 }
142 ceballos 1.10
143 ceballos 1.14 // photons
144 ceballos 1.10 else if (p->Status() == 1 && p->Is(MCParticle::kGamma) &&
145 ceballos 1.14 p->Pt() > fPtPhotonMin && p->AbsEta() < fEtaPhotonMax) {
146 ceballos 1.10 GenPhotons->Add(p);
147     }
148    
149 ceballos 1.16 // W/Z -> lnu for Magraph
150     if((p->IsQuark() || p->PdgId() == MCParticle::kGlu) && p->NDaughters() == 2){
151     CompositeParticle *diBoson = new CompositeParticle();
152     diBoson->AddDaughter(p->Daughter(0));
153     diBoson->AddDaughter(p->Daughter(1));
154     if(p->HasDaughter(MCParticle::kMu) && p->HasDaughter(MCParticle::kMuNu)){
155     if (fFillHist) hDVMass[0]->Fill(TMath::Min(diBoson->Mass(),199.999));
156     GenLeptons->Add(p->FindDaughter(MCParticle::kMu));
157     }
158     else if(p->HasDaughter(MCParticle::kEl) && p->HasDaughter(MCParticle::kElNu)){
159     hDVMass[1]->Fill(TMath::Min(diBoson->Mass(),199.999));
160     GenLeptons->Add(p->FindDaughter(MCParticle::kEl));
161     }
162     else if(p->HasDaughter(MCParticle::kTau) && p->HasDaughter(MCParticle::kTauNu)){
163     if (fFillHist) hDVMass[2]->Fill(TMath::Min(diBoson->Mass(),199.999));
164     const MCParticle *Tau = p->FindDaughter(MCParticle::kTau);
165     if (Tau->HasDaughter(MCParticle::kMu)) GenLeptons->Add(Tau->FindDaughter(MCParticle::kMu));
166     if (Tau->HasDaughter(MCParticle::kEl)) GenLeptons->Add(Tau->FindDaughter(MCParticle::kEl));
167     }
168     else if(p->Daughter(0)->AbsPdgId() == MCParticle::kMu && p->Daughter(1)->AbsPdgId() == MCParticle::kMu){
169     if (fFillHist) hDVMass[3]->Fill(TMath::Min(diBoson->Mass(),199.999));
170     GenLeptons->Add(p->Daughter(0));
171     GenLeptons->Add(p->Daughter(1));
172     }
173     else if(p->Daughter(0)->AbsPdgId() == MCParticle::kEl && p->Daughter(1)->AbsPdgId() == MCParticle::kEl){
174     if (fFillHist) hDVMass[4]->Fill(TMath::Min(diBoson->Mass(),199.999));
175     GenLeptons->Add(p->Daughter(0));
176     GenLeptons->Add(p->Daughter(1));
177     }
178     else if(p->Daughter(0)->AbsPdgId() == MCParticle::kTau && p->Daughter(1)->AbsPdgId() == MCParticle::kTau){
179     if (fFillHist) hDVMass[5]->Fill(TMath::Min(diBoson->Mass(),199.999));
180     const MCParticle *Tau0 = p->Daughter(0);
181     if (Tau0->HasDaughter(MCParticle::kMu)) GenLeptons->Add(Tau0->FindDaughter(MCParticle::kMu));
182     if (Tau0->HasDaughter(MCParticle::kEl)) GenLeptons->Add(Tau0->FindDaughter(MCParticle::kEl));
183     const MCParticle *Tau1 = p->Daughter(1);
184     if (Tau1->HasDaughter(MCParticle::kMu)) GenLeptons->Add(Tau1->FindDaughter(MCParticle::kMu));
185     if (Tau1->HasDaughter(MCParticle::kEl)) GenLeptons->Add(Tau1->FindDaughter(MCParticle::kEl));
186     }
187     delete diBoson;
188     }
189    
190     // t -> lnu for Magraph
191     if(p->PdgId() == MCParticle::kTop){
192     CompositeParticle *diBoson = new CompositeParticle();
193     if(p->HasDaughter(MCParticle::kMu) && p->HasDaughter(MCParticle::kMuNu)){
194     diBoson->AddDaughter(p->FindDaughter(MCParticle::kMu));
195     diBoson->AddDaughter(p->FindDaughter(MCParticle::kMuNu));
196     if (fFillHist) hDVMass[6]->Fill(TMath::Min(diBoson->Mass(),199.999));
197     GenLeptons->Add(p->FindDaughter(MCParticle::kMu));
198     }
199     else if(p->HasDaughter(MCParticle::kEl) && p->HasDaughter(MCParticle::kElNu)){
200     diBoson->AddDaughter(p->FindDaughter(MCParticle::kEl));
201     diBoson->AddDaughter(p->FindDaughter(MCParticle::kElNu));
202     if (fFillHist) hDVMass[7]->Fill(TMath::Min(diBoson->Mass(),199.999));
203     GenLeptons->Add(p->FindDaughter(MCParticle::kEl));
204     }
205     else if(p->HasDaughter(MCParticle::kTau) && p->HasDaughter(MCParticle::kTauNu)){
206     diBoson->AddDaughter(p->FindDaughter(MCParticle::kTau));
207     diBoson->AddDaughter(p->FindDaughter(MCParticle::kTauNu));
208     if (fFillHist) hDVMass[8]->Fill(TMath::Min(diBoson->Mass(),199.999));
209     const MCParticle *Tau = p->FindDaughter(MCParticle::kTau);
210     if (Tau->HasDaughter(MCParticle::kMu)) GenLeptons->Add(Tau->FindDaughter(MCParticle::kMu));
211     if (Tau->HasDaughter(MCParticle::kEl)) GenLeptons->Add(Tau->FindDaughter(MCParticle::kEl));
212     }
213     else if(!p->HasDaughter(MCParticle::kW)){
214     for(UInt_t nd=0; nd<p->NDaughters(); nd++)
215     if(p->Daughter(nd)->AbsPdgId() != MCParticle::kBottom &&
216     p->Daughter(nd)->AbsPdgId() != MCParticle::kGamma) diBoson->AddDaughter(p->Daughter(nd));
217     if (fFillHist) hDVMass[9]->Fill(TMath::Min(diBoson->Mass(),199.999));
218     }
219     delete diBoson;
220     }
221 loizides 1.5 }
222 loizides 1.1
223 loizides 1.8 // add objects to this event for other modules to use
224 loizides 1.13 AddObjThisEvt(GenLeptons);
225     AddObjThisEvt(GenAllLeptons);
226     AddObjThisEvt(GenTaus);
227     AddObjThisEvt(GenNeutrinos);
228     AddObjThisEvt(GenQuarks);
229     AddObjThisEvt(GenqqHs);
230     AddObjThisEvt(GenBosons);
231     AddObjThisEvt(GenPhotons);
232 ceballos 1.16
233 loizides 1.5 // fill histograms if requested
234     if (fFillHist) {
235    
236 loizides 1.6 // leptons
237 loizides 1.1 hDGenLeptons[0]->Fill(GenLeptons->GetEntries());
238 loizides 1.5 Int_t idxMaxLep[2] = {-1, -1};
239     Double_t ptMaxLep[2] = {-1, -1};
240     for(UInt_t i=0; i<GenLeptons->GetEntries(); i++) {
241 loizides 1.1 hDGenLeptons[1]->Fill(GenLeptons->At(i)->Pt());
242 loizides 1.5 hDGenLeptons[2]->Fill(TMath::Abs(GenLeptons->At(i)->Eta()));
243     hDGenLeptons[3]->Fill(GenLeptons->At(i)->PhiDeg());
244     for(UInt_t j=i+1; j<GenLeptons->GetEntries(); j++) {
245 loizides 1.1 CompositeParticle *dilepton = new CompositeParticle();
246     dilepton->AddDaughter(GenLeptons->At(i));
247     dilepton->AddDaughter(GenLeptons->At(j));
248     hDGenLeptons[4]->Fill(dilepton->Mass());
249     delete dilepton;
250     }
251 loizides 1.6 // selecting the two highest pt leptons
252 loizides 1.5 if (GenLeptons->At(i)->Pt() > ptMaxLep[0]) {
253 loizides 1.1 ptMaxLep[1] = ptMaxLep[0];
254     idxMaxLep[1] = idxMaxLep[0];
255     ptMaxLep[0] = GenLeptons->At(i)->Pt();
256     idxMaxLep[0] = i;
257     }
258 loizides 1.5 else if (GenLeptons->At(i)->Pt() > ptMaxLep[1]) {
259 loizides 1.1 ptMaxLep[1] = GenLeptons->At(i)->Pt();
260     idxMaxLep[1] = i;
261     }
262     }
263 loizides 1.6 // looking at events with at least two leptons
264 loizides 1.5 if (ptMaxLep[0] > 0 && ptMaxLep[1] > 0) {
265     hDGenLeptons[5]->Fill(TMath::Min(TMath::Max(TMath::Abs(GenLeptons->At(idxMaxLep[0])->Eta()),
266     TMath::Abs(GenLeptons->At(idxMaxLep[1])->Eta())),
267     4.999));
268     hDGenLeptons[6]->Fill(TMath::Min(TMath::Min(TMath::Abs(GenLeptons->At(idxMaxLep[0])->Eta()),
269     TMath::Abs(GenLeptons->At(idxMaxLep[1])->Eta())),
270     4.999));
271     if (TMath::Abs(GenLeptons->At(idxMaxLep[0])->Eta()) < 2.5 &&
272     TMath::Abs(GenLeptons->At(idxMaxLep[1])->Eta()) < 2.5) {
273 loizides 1.1 hDGenLeptons[7]->Fill(TMath::Min(GenLeptons->At(idxMaxLep[0])->Pt(),199.999));
274 loizides 1.5 if (GenLeptons->At(idxMaxLep[0])->Pt() > 20.0) {
275 loizides 1.1 hDGenLeptons[8]->Fill(TMath::Min(GenLeptons->At(idxMaxLep[1])->Pt(),199.999));
276 loizides 1.5 if (GenLeptons->At(idxMaxLep[1])->Pt() > 10.0) {
277 loizides 1.1 CompositeParticle *dilepton = new CompositeParticle();
278     dilepton->AddDaughter(GenLeptons->At(idxMaxLep[0]));
279     dilepton->AddDaughter(GenLeptons->At(idxMaxLep[1]));
280     hDGenLeptons[9]->Fill(TMath::Min(dilepton->Mass(),999.999));
281     hDGenLeptons[10]->Fill(MathUtils::DeltaPhi(GenLeptons->At(idxMaxLep[0])->Phi(),
282     GenLeptons->At(idxMaxLep[1])->Phi())
283     * 180./ TMath::Pi());
284     delete dilepton;
285     }
286     }
287     }
288     }
289    
290 loizides 1.6 // all leptons
291 ceballos 1.3 hDGenAllLeptons[0]->Fill(GenAllLeptons->GetEntries());
292 loizides 1.5 for(UInt_t i=0; i<GenAllLeptons->GetEntries(); i++) {
293 ceballos 1.3 hDGenAllLeptons[1]->Fill(GenAllLeptons->At(i)->Pt());
294     hDGenAllLeptons[2]->Fill(GenAllLeptons->At(i)->Eta());
295 loizides 1.5 hDGenAllLeptons[3]->Fill(GenAllLeptons->At(i)->PhiDeg());
296 ceballos 1.3 }
297    
298 loizides 1.6 // taus
299 loizides 1.1 hDGenTaus[0]->Fill(GenTaus->GetEntries());
300 loizides 1.5 for(UInt_t i=0; i<GenTaus->GetEntries(); i++) {
301 loizides 1.1 hDGenTaus[1]->Fill(GenTaus->At(i)->Pt());
302     hDGenTaus[2]->Fill(GenTaus->At(i)->Eta());
303 loizides 1.5 hDGenTaus[3]->Fill(GenTaus->At(i)->PhiDeg());
304 loizides 1.1 }
305    
306 loizides 1.6 // neutrinos
307 loizides 1.1 hDGenNeutrinos[0]->Fill(GenNeutrinos->GetEntries());
308     CompositeParticle *neutrinoTotal = new CompositeParticle();
309 loizides 1.5 for(UInt_t i=0; i<GenNeutrinos->GetEntries(); i++) {
310     if (GenNeutrinos->At(i)->HasMother())
311 loizides 1.1 neutrinoTotal->AddDaughter(GenNeutrinos->At(i));
312     }
313 loizides 1.5 if (GenNeutrinos->GetEntries() > 0) {
314 loizides 1.1 hDGenNeutrinos[1]->Fill(neutrinoTotal->Pt());
315     hDGenNeutrinos[2]->Fill(neutrinoTotal->Eta());
316 loizides 1.5 hDGenNeutrinos[3]->Fill(neutrinoTotal->PhiDeg());
317 loizides 1.1 }
318     delete neutrinoTotal;
319    
320 loizides 1.6 // quarks
321 loizides 1.1 hDGenQuarks[0]->Fill(GenQuarks->GetEntries());
322 loizides 1.5 for(UInt_t i=0; i<GenQuarks->GetEntries(); i++) {
323     for(UInt_t j=i+1; j<GenQuarks->GetEntries(); j++) {
324 loizides 1.1 CompositeParticle *dijet = new CompositeParticle();
325     dijet->AddDaughter(GenQuarks->At(i));
326     dijet->AddDaughter(GenQuarks->At(j));
327     hDGenQuarks[1]->Fill(dijet->Pt());
328     hDGenQuarks[2]->Fill(dijet->Mass());
329 loizides 1.5 if (TMath::Abs(GenQuarks->At(i)->Eta()) < 2.5 &&
330     TMath::Abs(GenQuarks->At(j)->Eta()) < 2.5) {
331 loizides 1.1 hDGenQuarks[3]->Fill(dijet->Pt());
332     hDGenQuarks[4]->Fill(dijet->Mass());
333     }
334     delete dijet;
335     }
336 ceballos 1.2 // b quark info
337 loizides 1.5 if (GenQuarks->At(i)->AbsPdgId() == 5) {
338 ceballos 1.2 hDGenQuarks[5]->Fill(GenQuarks->At(i)->Pt());
339     hDGenQuarks[6]->Fill(GenQuarks->At(i)->Eta());
340     hDGenQuarks[7]->Fill(GenQuarks->At(i)->Phi());
341 loizides 1.5 if (ptMaxLep[0] > 20 && ptMaxLep[1] > 15) {
342     if (TMath::Abs(GenLeptons->At(idxMaxLep[0])->Eta()) < 2.5 &&
343     TMath::Abs(GenLeptons->At(idxMaxLep[1])->Eta()) < 2.5) {
344 ceballos 1.2 hDGenQuarks[8]->Fill(GenQuarks->At(i)->Pt());
345     hDGenQuarks[9]->Fill(GenQuarks->At(i)->Eta());
346     hDGenQuarks[10]->Fill(GenQuarks->At(i)->Phi());
347     }
348     }
349     }
350     // t quark info
351 loizides 1.5 else if (GenQuarks->At(i)->AbsPdgId() == 6) {
352 ceballos 1.2 hDGenQuarks[11]->Fill(GenQuarks->At(i)->Pt());
353     hDGenQuarks[12]->Fill(GenQuarks->At(i)->Eta());
354     hDGenQuarks[13]->Fill(GenQuarks->At(i)->Phi());
355     }
356     // light quark info
357     else {
358     hDGenQuarks[14]->Fill(GenQuarks->At(i)->Pt());
359     hDGenQuarks[15]->Fill(GenQuarks->At(i)->Eta());
360     hDGenQuarks[16]->Fill(GenQuarks->At(i)->Phi());
361     }
362 loizides 1.1 }
363    
364 loizides 1.6 // wbf
365 loizides 1.5 if (GenqqHs->GetEntries() == 2) {
366 loizides 1.1 hDGenWBF[0]->Fill(MathUtils::DeltaPhi(GenqqHs->At(0)->Phi(),
367     GenqqHs->At(1)->Phi()) * 180./ TMath::Pi());
368 loizides 1.5 hDGenWBF[1]->Fill(TMath::Abs(GenqqHs->At(0)->Eta()-GenqqHs->At(1)->Eta()));
369 loizides 1.1 hDGenWBF[2]->Fill(TMath::Max(GenqqHs->At(0)->Pt(),GenqqHs->At(1)->Pt()));
370     hDGenWBF[3]->Fill(TMath::Min(GenqqHs->At(0)->Pt(),GenqqHs->At(1)->Pt()));
371     CompositeParticle *diqq = new CompositeParticle();
372     diqq->AddDaughter(GenqqHs->At(0));
373     diqq->AddDaughter(GenqqHs->At(1));
374     hDGenWBF[4]->Fill(diqq->Mass());
375     delete diqq;
376     }
377    
378 loizides 1.6 // bosons
379 loizides 1.1 hDGenBosons[0]->Fill(GenBosons->GetEntries());
380 loizides 1.5 for(UInt_t i=0; i<GenBosons->GetEntries(); i++) {
381 loizides 1.1 hDGenBosons[1]->Fill(GenBosons->At(i)->Pt());
382     hDGenBosons[2]->Fill(GenBosons->At(i)->Eta());
383     hDGenBosons[3]->Fill(GenBosons->At(i)->Mass());
384     }
385 ceballos 1.10
386     // photons
387     hDGenPhotons[0]->Fill(GenPhotons->GetEntries());
388     for(UInt_t i=0; i<GenPhotons->GetEntries(); i++) {
389     hDGenPhotons[1]->Fill(GenPhotons->At(i)->Pt());
390     hDGenPhotons[2]->Fill(GenPhotons->At(i)->Eta());
391     }
392 loizides 1.1 }
393     }
394    
395     //--------------------------------------------------------------------------------------------------
396     void GeneratorMod::SlaveBegin()
397     {
398 loizides 1.5 // Book branch and histograms if wanted.
399    
400 loizides 1.1 ReqBranch(fMCPartName, fParticles);
401    
402 loizides 1.5 // fill histograms
403     if (fFillHist == kTRUE) {
404     char sb[1024];
405    
406 loizides 1.6 // leptons from W
407 loizides 1.1 sprintf(sb,"hDGenLeptons_%d", 0); hDGenLeptons[0] = new TH1D(sb,sb,10,-0.5,9.5);
408     sprintf(sb,"hDGenLeptons_%d", 1); hDGenLeptons[1] = new TH1D(sb,sb,100,0.0,200.0);
409     sprintf(sb,"hDGenLeptons_%d", 2); hDGenLeptons[2] = new TH1D(sb,sb,50,0.0,5.0);
410     sprintf(sb,"hDGenLeptons_%d", 3); hDGenLeptons[3] = new TH1D(sb,sb,90,0.0,180.0);
411     sprintf(sb,"hDGenLeptons_%d", 4); hDGenLeptons[4] = new TH1D(sb,sb,1000,0.0,1000.0);
412     sprintf(sb,"hDGenLeptons_%d", 5); hDGenLeptons[5] = new TH1D(sb,sb,50,0.0,5.0);
413     sprintf(sb,"hDGenLeptons_%d", 6); hDGenLeptons[6] = new TH1D(sb,sb,50,0.0,5.0);
414     sprintf(sb,"hDGenLeptons_%d", 7); hDGenLeptons[7] = new TH1D(sb,sb,100,0.0,200.0);
415     sprintf(sb,"hDGenLeptons_%d", 8); hDGenLeptons[8] = new TH1D(sb,sb,100,0.0,200.0);
416     sprintf(sb,"hDGenLeptons_%d", 9); hDGenLeptons[9] = new TH1D(sb,sb,1000,0.0,1000.0);
417     sprintf(sb,"hDGenLeptons_%d",10); hDGenLeptons[10] = new TH1D(sb,sb,90,0.0,180.0);
418 loizides 1.5 for(Int_t i=0; i<11; i++) AddOutput(hDGenLeptons[i]);
419 loizides 1.1
420 loizides 1.6 // all leptons
421 ceballos 1.3 sprintf(sb,"hDGenAllLeptons_%d", 0); hDGenAllLeptons[0] = new TH1D(sb,sb,10,-0.5,9.5);
422     sprintf(sb,"hDGenAllLeptons_%d", 1); hDGenAllLeptons[1] = new TH1D(sb,sb,100,0.0,200.0);
423 ceballos 1.14 sprintf(sb,"hDGenAllLeptons_%d", 2); hDGenAllLeptons[2] = new TH1D(sb,sb,100,-5.0,5.0);
424 ceballos 1.3 sprintf(sb,"hDGenAllLeptons_%d", 3); hDGenAllLeptons[3] = new TH1D(sb,sb,90,0.0,180.0);
425 loizides 1.5 for(Int_t i=0; i<4; i++) AddOutput(hDGenAllLeptons[i]);
426 ceballos 1.3
427 loizides 1.6 // taus
428 loizides 1.1 sprintf(sb,"hDGenTaus_%d", 0); hDGenTaus[0] = new TH1D(sb,sb,10,-0.5,9.5);
429     sprintf(sb,"hDGenTaus_%d", 1); hDGenTaus[1] = new TH1D(sb,sb,100,0.0,200.0);
430     sprintf(sb,"hDGenTaus_%d", 2); hDGenTaus[2] = new TH1D(sb,sb,100,-5.0,5.0);
431     sprintf(sb,"hDGenTaus_%d", 3); hDGenTaus[3] = new TH1D(sb,sb,90,0.0,180.0);
432 loizides 1.5 for(Int_t i=0; i<4; i++) AddOutput(hDGenTaus[i]);
433 loizides 1.1
434 loizides 1.6 // neutrinos
435 loizides 1.1 sprintf(sb,"hDGenNeutrinos_%d", 0); hDGenNeutrinos[0] = new TH1D(sb,sb,10,-0.5,9.5);
436     sprintf(sb,"hDGenNeutrinos_%d", 1); hDGenNeutrinos[1] = new TH1D(sb,sb,100,0.0,200.0);
437     sprintf(sb,"hDGenNeutrinos_%d", 2); hDGenNeutrinos[2] = new TH1D(sb,sb,100,-5.0,5.0);
438     sprintf(sb,"hDGenNeutrinos_%d", 3); hDGenNeutrinos[3] = new TH1D(sb,sb,90,0.0,180.0);
439 loizides 1.5 for(Int_t i=0; i<4; i++) AddOutput(hDGenNeutrinos[i]);
440 loizides 1.1
441 loizides 1.6 // quarks
442 loizides 1.1 sprintf(sb,"hDGenQuarks_%d", 0); hDGenQuarks[0] = new TH1D(sb,sb,10,-0.5,9.5);
443     sprintf(sb,"hDGenQuarks_%d", 1); hDGenQuarks[1] = new TH1D(sb,sb,200,0.0,400.);
444     sprintf(sb,"hDGenQuarks_%d", 2); hDGenQuarks[2] = new TH1D(sb,sb,2000,0.0,2000.);
445     sprintf(sb,"hDGenQuarks_%d", 3); hDGenQuarks[3] = new TH1D(sb,sb,200,0.0,400.);
446     sprintf(sb,"hDGenQuarks_%d", 4); hDGenQuarks[4] = new TH1D(sb,sb,2000,0.0,2000.);
447 ceballos 1.2 sprintf(sb,"hDGenQuarks_%d", 5); hDGenQuarks[5] = new TH1D(sb,sb,200,0.0,400.);
448     sprintf(sb,"hDGenQuarks_%d", 6); hDGenQuarks[6] = new TH1D(sb,sb,200,-10.0,10.);
449 loizides 1.5 sprintf(sb,"hDGenQuarks_%d", 7); hDGenQuarks[7] = new TH1D(sb,sb,200,-TMath::Pi(),
450     TMath::Pi());
451 ceballos 1.2 sprintf(sb,"hDGenQuarks_%d", 8); hDGenQuarks[8] = new TH1D(sb,sb,200,0.0,400.);
452     sprintf(sb,"hDGenQuarks_%d", 9); hDGenQuarks[9] = new TH1D(sb,sb,200,-10.0,10.);
453 loizides 1.5 sprintf(sb,"hDGenQuarks_%d",10); hDGenQuarks[10] = new TH1D(sb,sb,200,-TMath::Pi(),
454     TMath::Pi());
455 ceballos 1.2 sprintf(sb,"hDGenQuarks_%d",11); hDGenQuarks[11] = new TH1D(sb,sb,200,0.0,400.);
456     sprintf(sb,"hDGenQuarks_%d",12); hDGenQuarks[12] = new TH1D(sb,sb,200,-10.0,10.);
457 loizides 1.5 sprintf(sb,"hDGenQuarks_%d",13); hDGenQuarks[13] = new TH1D(sb,sb,200,-TMath::Pi(),
458     TMath::Pi());
459 ceballos 1.2 sprintf(sb,"hDGenQuarks_%d",14); hDGenQuarks[14] = new TH1D(sb,sb,200,0.0,400.);
460     sprintf(sb,"hDGenQuarks_%d",15); hDGenQuarks[15] = new TH1D(sb,sb,200,-10.0,10.);
461 loizides 1.5 sprintf(sb,"hDGenQuarks_%d",16); hDGenQuarks[16] = new TH1D(sb,sb,200,-TMath::Pi(),
462     TMath::Pi());
463     for(Int_t i=0; i<17; i++) AddOutput(hDGenQuarks[i]);
464 loizides 1.1
465     // qqH
466     sprintf(sb,"hDGenWBF_%d", 0); hDGenWBF[0] = new TH1D(sb,sb,90,0.0,180.);
467     sprintf(sb,"hDGenWBF_%d", 1); hDGenWBF[1] = new TH1D(sb,sb,100,0.0,10.);
468     sprintf(sb,"hDGenWBF_%d", 2); hDGenWBF[2] = new TH1D(sb,sb,200,0.0,400.);
469     sprintf(sb,"hDGenWBF_%d", 3); hDGenWBF[3] = new TH1D(sb,sb,200,0.0,400.);
470     sprintf(sb,"hDGenWBF_%d", 4); hDGenWBF[4] = new TH1D(sb,sb,200,0.0,4000.);
471 loizides 1.5 for(Int_t i=0; i<5; i++) AddOutput(hDGenWBF[i]);
472 loizides 1.1
473 loizides 1.6 // bosons
474 loizides 1.1 sprintf(sb,"hDGenBosons_%d", 0); hDGenBosons[0] = new TH1D(sb,sb,10,-0.5,9.5);
475     sprintf(sb,"hDGenBosons_%d", 1); hDGenBosons[1] = new TH1D(sb,sb,200,0.0,400.0);
476     sprintf(sb,"hDGenBosons_%d", 2); hDGenBosons[2] = new TH1D(sb,sb,100,-5.0,5.0);
477     sprintf(sb,"hDGenBosons_%d", 3); hDGenBosons[3] = new TH1D(sb,sb,2000,0.0,2000.0);
478 loizides 1.5 for(Int_t i=0; i<4; i++) AddOutput(hDGenBosons[i]);
479 ceballos 1.10
480     // photons
481     sprintf(sb,"hDGenPhotons_%d", 0); hDGenPhotons[0] = new TH1D(sb,sb,10,-0.5,9.5);
482     sprintf(sb,"hDGenPhotons_%d", 1); hDGenPhotons[1] = new TH1D(sb,sb,200,0.0,400.0);
483 ceballos 1.14 sprintf(sb,"hDGenPhotons_%d", 2); hDGenPhotons[2] = new TH1D(sb,sb,100,-5.0,5.0);
484 ceballos 1.10 for(Int_t i=0; i<3; i++) AddOutput(hDGenPhotons[i]);
485 ceballos 1.16
486     // auxiliar
487     sprintf(sb,"hDVMass_%d", 0); hDVMass[0] = new TH1D(sb,sb,200,0.,200.);
488     sprintf(sb,"hDVMass_%d", 1); hDVMass[1] = new TH1D(sb,sb,200,0.,200.);
489     sprintf(sb,"hDVMass_%d", 2); hDVMass[2] = new TH1D(sb,sb,200,0.,200.);
490     sprintf(sb,"hDVMass_%d", 3); hDVMass[3] = new TH1D(sb,sb,200,0.,200.);
491     sprintf(sb,"hDVMass_%d", 4); hDVMass[4] = new TH1D(sb,sb,200,0.,200.);
492     sprintf(sb,"hDVMass_%d", 5); hDVMass[5] = new TH1D(sb,sb,200,0.,200.);
493     sprintf(sb,"hDVMass_%d", 6); hDVMass[6] = new TH1D(sb,sb,200,0.,200.);
494     sprintf(sb,"hDVMass_%d", 7); hDVMass[7] = new TH1D(sb,sb,200,0.,200.);
495     sprintf(sb,"hDVMass_%d", 8); hDVMass[8] = new TH1D(sb,sb,200,0.,200.);
496     sprintf(sb,"hDVMass_%d", 9); hDVMass[9] = new TH1D(sb,sb,200,0.,200.);
497     for(Int_t i=0; i<10; i++) AddOutput(hDVMass[i]);
498 loizides 1.1 }
499     }