ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/MitPhysics/Mods/src/GeneratorMod.cc
Revision: 1.57
Committed: Sun Dec 6 14:59:28 2009 UTC (15 years, 5 months ago) by ceballos
Content type: text/plain
Branch: MAIN
CVS Tags: Mit_012i, Mit_012g, Mit_012f, Mit_012e
Changes since 1.56: +668 -663 lines
Log Message:
small updates

File Contents

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