ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/MitPhysics/Mods/src/GeneratorMod.cc
Revision: 1.58
Committed: Sat Mar 13 20:51:27 2010 UTC (15 years, 1 month ago) by ceballos
Content type: text/plain
Branch: MAIN
CVS Tags: Mit_014c, Mit_014b, Mit_014a, Mit_014, Mit_014pre3, Mit_014pre2, Mit_014pre1, Mit_013d, Mit_013c, Mit_013b, Mit_013a, Mit_013, Mit_013pre1
Changes since 1.57: +37 -12 lines
Log Message: 
adding some info

File Contents

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