ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/MitPhysics/Mods/src/GeneratorMod.cc
Revision: 1.63
Committed: Mon Nov 8 14:20:24 2010 UTC (14 years, 5 months ago) by ceballos
Content type: text/plain
Branch: MAIN
CVS Tags: Mit_018, Mit_017, Mit_017pre3, Mit_017pre2, Mit_017pre1, Mit_016
Changes since 1.62: +8 -6 lines
Log Message:
new cuts

File Contents

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