SFrameTools/src/TTbarGen.cxx

#include "include/TTbarGen.h"

using namespace std;

TTbarGen::TTbarGen(  BaseCycleContainer* bcc )
{
  m_pdgId1 = 0;
  m_pdgId2 = 0;
  for(unsigned int i=0; i<bcc->genparticles->size(); ++i)
   {
     GenParticle genp = bcc->genparticles->at(i);

     if ( genp.pdgId() == 6 ) 
       {
         m_Top = genp;
         if (genp.daughter(bcc->genparticles,1) && genp.daughter(bcc->genparticles,2))
           {
             if (genp.daughter(bcc->genparticles,1)->pdgId() == 24) 
               {
                 m_indexW = genp.daughter(bcc->genparticles,1)->index();
                 m_WTop = bcc->genparticles->at(m_indexW);
                 m_indexb = genp.daughter(bcc->genparticles,2)->index();
                 m_bTop = bcc->genparticles->at(m_indexb);
               }
             if (genp.daughter(bcc->genparticles,2)->pdgId() == 24)
               {
                 m_indexW = genp.daughter(bcc->genparticles,2)->index();
                 m_WTop = bcc->genparticles->at(m_indexW);
                 m_indexb = genp.daughter(bcc->genparticles,1)->index();
                 m_bTop = bcc->genparticles->at(m_indexb);
               }
           }
       }
     if ( genp.pdgId() == -6 ) 
       {
         m_Antitop = genp;
         if (genp.daughter(bcc->genparticles,1) && (genp.daughter(bcc->genparticles,2)))
           {
             if (genp.daughter(bcc->genparticles,1)->pdgId() == -24) 
               {
                 m_indexW = genp.daughter(bcc->genparticles,1)->index();
                 m_WAntitop = bcc->genparticles->at(m_indexW);
                 m_indexb = genp.daughter(bcc->genparticles,2)->index();
                 m_bAntitop = bcc->genparticles->at(m_indexb);
               }
             if (genp.daughter(bcc->genparticles,2)->pdgId() == -24)
               {
                 m_indexW = genp.daughter(bcc->genparticles,2)->index();
                 m_WAntitop = bcc->genparticles->at(m_indexW);
                 m_indexb = genp.daughter(bcc->genparticles,1)->index();
                 m_bAntitop = bcc->genparticles->at(m_indexb);
               }
           }
       }
     if ( genp.pdgId() == 24 && genp.mother(bcc->genparticles,1)->pdgId() == 6)
       {
         if (genp.daughter(bcc->genparticles,1) && genp.daughter(bcc->genparticles,2))
           {
             m_index1 = genp.daughter(bcc->genparticles,1)->index();
             m_Wdecay1 = bcc->genparticles->at(m_index1);
             m_pdgId1= m_Wdecay1.pdgId();
             m_index2 = genp.daughter(bcc->genparticles,2)->index();
             m_Wdecay2 = bcc->genparticles->at(m_index2);
             
           }
       }
     if ( genp.pdgId() == -24 && genp.mother(bcc->genparticles,1)->pdgId() == -6)
       {
         if (genp.daughter(bcc->genparticles,1) && (genp.daughter(bcc->genparticles,2)))
           {
             m_index1 = genp.daughter(bcc->genparticles,1)->index();
             m_WMinusdecay1 = bcc->genparticles->at(m_index1);
             m_pdgId2= m_WMinusdecay1.pdgId(); 
             m_index2 = genp.daughter(bcc->genparticles,2)->index();
             m_WMinusdecay2 = bcc->genparticles->at(m_index2);
           }
       }
   }

  // W not linked correctly -> W decay products have top as mother
  if(m_pdgId1==0 || m_pdgId2==0 ){
    //search all particles with top as mother; store that ones which are not W or b (or equivalent light quark)
    bool notfilled1=true;
    bool notfilled2=true;

    for(unsigned int i=0; i<bcc->genparticles->size(); ++i)
      {
        GenParticle genp = bcc->genparticles->at(i);
        

        if(genp.mother(bcc->genparticles,1)){
          if(genp.mother(bcc->genparticles,1)->pdgId()==6 && abs(genp.pdgId())!=24 &&  genp.pdgId()!=1 && genp.pdgId()!=3 && genp.pdgId()!=5 ){
            
            if(notfilled1){
              m_Wdecay1 = genp;
              m_pdgId1 = genp.pdgId();
              notfilled1=false;
            }
            else{
              m_Wdecay2 = genp;
            }
          }
          

          if(genp.mother(bcc->genparticles,1)->pdgId()==-6 && abs(genp.pdgId())!=24 && genp.pdgId()!=-1 && genp.pdgId()!=-3 && genp.pdgId()!=-5 ){
            if(notfilled2){
              m_WMinusdecay1 = genp;
              m_pdgId2 = genp.pdgId();
              notfilled2=false;
            }
            else{
              m_WMinusdecay2 = genp;
            }
          }
          if(genp.mother(bcc->genparticles,1)->pdgId()==6 && abs(genp.pdgId())!=24 &&  (genp.pdgId()==1 || genp.pdgId()==3 || genp.pdgId()==5) ){
            m_bTop = genp;
          }
          if(genp.mother(bcc->genparticles,1)->pdgId()==-6 && abs(genp.pdgId())!=24 &&  (genp.pdgId()==-1 || genp.pdgId()==-3 || genp.pdgId()==-5) ){
            m_bAntitop = genp;
          }
        }

      }
    //fill the missing W bosons
    GenParticle Wplus, Wminus;

    LorentzVector Wp, Wm;
    Wp.SetPxPyPzE( Wdecay1().v4().Px()+ Wdecay2().v4().Px(),  Wdecay1().v4().Py()+ Wdecay2().v4().Py(),  Wdecay1().v4().Pz()+ Wdecay2().v4().Pz(),  Wdecay1().v4().E()+ Wdecay2().v4().E() );
    Wm.SetPxPyPzE( WMinusdecay1().v4().Px()+ WMinusdecay2().v4().Px(),  WMinusdecay1().v4().Py()+ WMinusdecay2().v4().Py(),  WMinusdecay1().v4().Pz()+ WMinusdecay2().v4().Pz(),  WMinusdecay1().v4().E()+ WMinusdecay2().v4().E() );

    Wplus.set_v4( Wp );
    Wplus.set_charge(1.);
    Wplus.set_pdgId(24);
    Wplus.set_status(3);
    Wminus.set_v4( Wm );
    Wminus.set_charge(-1.);
    Wminus.set_pdgId(-24);
    Wminus.set_status(3);
    m_WTop = Wplus;
    m_WAntitop = Wminus;
  }


}   
 

TTbarGen::~TTbarGen()
{
  // default destructor, does nothing
}

GenParticle TTbarGen::Top()
{
  return m_Top;
}

GenParticle TTbarGen::Antitop()
{
  return m_Antitop;
} 

GenParticle TTbarGen::WTop()
{
  return m_WTop;
}

GenParticle TTbarGen::WAntitop()
{
  return m_WAntitop;
}

GenParticle TTbarGen::bTop()
{
  return m_bTop;
}

GenParticle TTbarGen::bAntitop()
{
  return m_bAntitop;
} 

GenParticle TTbarGen::Wdecay1()
{
  return m_Wdecay1;
} 

GenParticle TTbarGen::Wdecay2()
{
  return m_Wdecay2;
} 

GenParticle TTbarGen::WMinusdecay1()
{
  return m_WMinusdecay1;
} 

GenParticle TTbarGen::WMinusdecay2()
{
  return m_WMinusdecay2;
} 

TTbarGen::E_DecayChannel TTbarGen::DecayChannel()
{
  if ((abs(m_pdgId1)==11 || abs(m_pdgId1)==12) && (abs(m_pdgId2)==11 || abs(m_pdgId2)==12)) m_type = e_ee;
  if ((abs(m_pdgId1)==13 || abs(m_pdgId1)==14) && (abs(m_pdgId2)==13 || abs(m_pdgId2)==14)) m_type = e_mumu;
  if ((abs(m_pdgId1)==15 || abs(m_pdgId1)==16) && (abs(m_pdgId2)==15 || abs(m_pdgId2)==16)) m_type = e_tautau;
  if ((abs(m_pdgId1)<=5) && (abs(m_pdgId2)<=5)) m_type = e_had;
  if ((abs(m_pdgId1)==11 || abs(m_pdgId1)==12) && (abs(m_pdgId2)<=5)) m_type = e_ehad;
  if ((abs(m_pdgId2)==11 || abs(m_pdgId2)==12) && (abs(m_pdgId1)<=5)) m_type = e_ehad;
  if ((abs(m_pdgId1)==13 || abs(m_pdgId1)==14) && (abs(m_pdgId2)<=5)) m_type = e_muhad;
  if ((abs(m_pdgId2)==13 || abs(m_pdgId2)==14) && (abs(m_pdgId1)<=5)) m_type = e_muhad;
  if ((abs(m_pdgId1)==15 || abs(m_pdgId1)==16) && (abs(m_pdgId2)<=5)) m_type = e_tauhad;
  if ((abs(m_pdgId2)==15 || abs(m_pdgId2)==16) && (abs(m_pdgId1)<=5)) m_type = e_tauhad;
  if ((abs(m_pdgId1)==11 || abs(m_pdgId1)==12) && (abs(m_pdgId2)==13 || abs(m_pdgId2)==14)) m_type = e_emu;
  if ((abs(m_pdgId2)==11 || abs(m_pdgId2)==12) && (abs(m_pdgId1)==13 || abs(m_pdgId1)==14)) m_type = e_emu;
  if ((abs(m_pdgId1)==11 || abs(m_pdgId1)==12) && (abs(m_pdgId2)==15 || abs(m_pdgId2)==16)) m_type = e_etau;
  if ((abs(m_pdgId2)==11 || abs(m_pdgId2)==12) && (abs(m_pdgId1)==15 || abs(m_pdgId1)==16)) m_type = e_etau;
  if ((abs(m_pdgId1)==13 || abs(m_pdgId1)==14) && (abs(m_pdgId2)==15 || abs(m_pdgId2)==16)) m_type = e_mutau;
  if ((abs(m_pdgId2)==13 || abs(m_pdgId2)==14) && (abs(m_pdgId1)==15 || abs(m_pdgId1)==16)) m_type = e_mutau;
  return m_type;
}

bool TTbarGen::IsTopHadronicDecay()
{
  if (abs(m_pdgId1)<=5) return true;
  else return false;
}

bool TTbarGen::IsAntiTopHadronicDecay()
{
  if (abs(m_pdgId2)<=5) return true;
  else return false;
}

GenParticle TTbarGen::ChargedLepton(){
  GenParticle lepton;  
  if (m_type != e_ehad &&  m_type != e_muhad  && m_type!= e_tauhad){
    std::cerr << "This is no l+jets ttbar event, no charged lepton found" <<std::endl;
    return lepton;
  }

  int nlepton=0;
  if(abs(Wdecay1().pdgId())==11 || abs(Wdecay1().pdgId())==13 || abs(Wdecay1().pdgId())==15){
    lepton = Wdecay1();
    nlepton++;
  }
  if(abs(Wdecay2().pdgId())==11 || abs(Wdecay2().pdgId())==13 || abs(Wdecay2().pdgId())==15){
    lepton = Wdecay2();
    nlepton++;
  }
  if(abs(WMinusdecay1().pdgId())==11 || abs(WMinusdecay1().pdgId())==13 || abs(WMinusdecay1().pdgId())==15){
    lepton = WMinusdecay1();
    nlepton++;
  }
  if(abs(WMinusdecay2().pdgId())==11 || abs(WMinusdecay2().pdgId())==13 || abs(WMinusdecay2().pdgId())==15){
    lepton = WMinusdecay2();
    nlepton++;
  }
  if(nlepton!=1) std::cerr << "Not exactly one lepton found " << nlepton<< std::endl;

  return lepton;

}

GenParticle TTbarGen::Neutrino(){
  GenParticle neutrino;
  
  if (m_type != e_ehad &&  m_type != e_muhad  && m_type!= e_tauhad){
    std::cerr << "This is no l+jets ttbar event, no neutrino found" <<std::endl;
    return neutrino;
  }

  int nneutrino=0;
  if(abs(Wdecay1().pdgId())==12 || abs(Wdecay1().pdgId())==14 || abs(Wdecay1().pdgId())==16){
    neutrino = Wdecay1();
    nneutrino++;
  }
  if(abs(Wdecay2().pdgId())==12 || abs(Wdecay2().pdgId())==14 || abs(Wdecay2().pdgId())==16){
    neutrino = Wdecay2();
    nneutrino++;
  }
  if(abs(WMinusdecay1().pdgId())==12 || abs(WMinusdecay1().pdgId())==14 || abs(WMinusdecay1().pdgId())==16){
    neutrino = WMinusdecay1();
    nneutrino++;
  }
  if(abs(WMinusdecay2().pdgId())==12 || abs(WMinusdecay2().pdgId())==14 || abs(WMinusdecay2().pdgId())==16){
    neutrino = WMinusdecay2();
    nneutrino++;
  }
  if(nneutrino!=1) std::cerr << "Not exactly one neutrino found " << nneutrino<< std::endl;

  return neutrino;

}

GenParticle TTbarGen::TopLep(){
 
  if(ChargedLepton().charge()>0) return Top();
  else return Antitop();

}

GenParticle TTbarGen::TopHad(){

  if(ChargedLepton().charge()<0) return Top();
  else return Antitop();

}


GenParticle TTbarGen::BLep(){
 
  if(ChargedLepton().charge()>0) return bTop();
  else return bAntitop();

}

GenParticle TTbarGen::BHad(){
 
  if(ChargedLepton().charge()<0) return bTop();
  else return bAntitop();

}

GenParticle TTbarGen::WLep(){
 
  if(ChargedLepton().charge()>0) return WTop();
  else return WAntitop();

}

GenParticle TTbarGen::WHad(){
 
  if(ChargedLepton().charge()<0) return WTop();
  else return WAntitop();

}

GenParticle TTbarGen::Q1(){
 
  if(ChargedLepton().charge()>0) return WMinusdecay1();
  else return Wdecay1();

}

GenParticle TTbarGen::Q2(){
 
  if(ChargedLepton().charge()>0) return WMinusdecay2();
  else return Wdecay2();

}

Revision:	1.10
Committed:	Wed Jun 12 12:33:40 2013 UTC (11 years, 10 months ago) by peiffer
Content type:	text/plain
Branch:	MAIN
CVS Tags:	HEAD
Changes since 1.9:	+1 -3 lines
Error occurred while calculating annotation data.
Log Message:	removed ObjectHandler
#	Content
1	#include "include/TTbarGen.h"
2
3	using namespace std;
4
5	TTbarGen::TTbarGen( BaseCycleContainer* bcc )
6	{
7	m_pdgId1 = 0;
8	m_pdgId2 = 0;
9	for(unsigned int i=0; i<bcc->genparticles->size(); ++i)
10	{
11	GenParticle genp = bcc->genparticles->at(i);
12
13	if ( genp.pdgId() == 6 )
14	{
15	m_Top = genp;
16	if (genp.daughter(bcc->genparticles,1) && genp.daughter(bcc->genparticles,2))
17	{
18	if (genp.daughter(bcc->genparticles,1)->pdgId() == 24)
19	{
20	m_indexW = genp.daughter(bcc->genparticles,1)->index();
21	m_WTop = bcc->genparticles->at(m_indexW);
22	m_indexb = genp.daughter(bcc->genparticles,2)->index();
23	m_bTop = bcc->genparticles->at(m_indexb);
24	}
25	if (genp.daughter(bcc->genparticles,2)->pdgId() == 24)
26	{
27	m_indexW = genp.daughter(bcc->genparticles,2)->index();
28	m_WTop = bcc->genparticles->at(m_indexW);
29	m_indexb = genp.daughter(bcc->genparticles,1)->index();
30	m_bTop = bcc->genparticles->at(m_indexb);
31	}
32	}
33	}
34	if ( genp.pdgId() == -6 )
35	{
36	m_Antitop = genp;
37	if (genp.daughter(bcc->genparticles,1) && (genp.daughter(bcc->genparticles,2)))
38	{
39	if (genp.daughter(bcc->genparticles,1)->pdgId() == -24)
40	{
41	m_indexW = genp.daughter(bcc->genparticles,1)->index();
42	m_WAntitop = bcc->genparticles->at(m_indexW);
43	m_indexb = genp.daughter(bcc->genparticles,2)->index();
44	m_bAntitop = bcc->genparticles->at(m_indexb);
45	}
46	if (genp.daughter(bcc->genparticles,2)->pdgId() == -24)
47	{
48	m_indexW = genp.daughter(bcc->genparticles,2)->index();
49	m_WAntitop = bcc->genparticles->at(m_indexW);
50	m_indexb = genp.daughter(bcc->genparticles,1)->index();
51	m_bAntitop = bcc->genparticles->at(m_indexb);
52	}
53	}
54	}
55	if ( genp.pdgId() == 24 && genp.mother(bcc->genparticles,1)->pdgId() == 6)
56	{
57	if (genp.daughter(bcc->genparticles,1) && genp.daughter(bcc->genparticles,2))
58	{
59	m_index1 = genp.daughter(bcc->genparticles,1)->index();
60	m_Wdecay1 = bcc->genparticles->at(m_index1);
61	m_pdgId1= m_Wdecay1.pdgId();
62	m_index2 = genp.daughter(bcc->genparticles,2)->index();
63	m_Wdecay2 = bcc->genparticles->at(m_index2);
64
65	}
66	}
67	if ( genp.pdgId() == -24 && genp.mother(bcc->genparticles,1)->pdgId() == -6)
68	{
69	if (genp.daughter(bcc->genparticles,1) && (genp.daughter(bcc->genparticles,2)))
70	{
71	m_index1 = genp.daughter(bcc->genparticles,1)->index();
72	m_WMinusdecay1 = bcc->genparticles->at(m_index1);
73	m_pdgId2= m_WMinusdecay1.pdgId();
74	m_index2 = genp.daughter(bcc->genparticles,2)->index();
75	m_WMinusdecay2 = bcc->genparticles->at(m_index2);
76	}
77	}
78	}
79
80	// W not linked correctly -> W decay products have top as mother
81	if(m_pdgId1==0 \|\| m_pdgId2==0 ){
82	//search all particles with top as mother; store that ones which are not W or b (or equivalent light quark)
83	bool notfilled1=true;
84	bool notfilled2=true;
85
86	for(unsigned int i=0; i<bcc->genparticles->size(); ++i)
87	{
88	GenParticle genp = bcc->genparticles->at(i);
89
90
91	if(genp.mother(bcc->genparticles,1)){
92	if(genp.mother(bcc->genparticles,1)->pdgId()==6 && abs(genp.pdgId())!=24 && genp.pdgId()!=1 && genp.pdgId()!=3 && genp.pdgId()!=5 ){
93
94	if(notfilled1){
95	m_Wdecay1 = genp;
96	m_pdgId1 = genp.pdgId();
97	notfilled1=false;
98	}
99	else{
100	m_Wdecay2 = genp;
101	}
102	}
103
104
105	if(genp.mother(bcc->genparticles,1)->pdgId()==-6 && abs(genp.pdgId())!=24 && genp.pdgId()!=-1 && genp.pdgId()!=-3 && genp.pdgId()!=-5 ){
106	if(notfilled2){
107	m_WMinusdecay1 = genp;
108	m_pdgId2 = genp.pdgId();
109	notfilled2=false;
110	}
111	else{
112	m_WMinusdecay2 = genp;
113	}
114	}
115	if(genp.mother(bcc->genparticles,1)->pdgId()==6 && abs(genp.pdgId())!=24 && (genp.pdgId()==1 \|\| genp.pdgId()==3 \|\| genp.pdgId()==5) ){
116	m_bTop = genp;
117	}
118	if(genp.mother(bcc->genparticles,1)->pdgId()==-6 && abs(genp.pdgId())!=24 && (genp.pdgId()==-1 \|\| genp.pdgId()==-3 \|\| genp.pdgId()==-5) ){
119	m_bAntitop = genp;
120	}
121	}
122
123	}
124	//fill the missing W bosons
125	GenParticle Wplus, Wminus;
126
127	LorentzVector Wp, Wm;
128	Wp.SetPxPyPzE( Wdecay1().v4().Px()+ Wdecay2().v4().Px(), Wdecay1().v4().Py()+ Wdecay2().v4().Py(), Wdecay1().v4().Pz()+ Wdecay2().v4().Pz(), Wdecay1().v4().E()+ Wdecay2().v4().E() );
129	Wm.SetPxPyPzE( WMinusdecay1().v4().Px()+ WMinusdecay2().v4().Px(), WMinusdecay1().v4().Py()+ WMinusdecay2().v4().Py(), WMinusdecay1().v4().Pz()+ WMinusdecay2().v4().Pz(), WMinusdecay1().v4().E()+ WMinusdecay2().v4().E() );
130
131	Wplus.set_v4( Wp );
132	Wplus.set_charge(1.);
133	Wplus.set_pdgId(24);
134	Wplus.set_status(3);
135	Wminus.set_v4( Wm );
136	Wminus.set_charge(-1.);
137	Wminus.set_pdgId(-24);
138	Wminus.set_status(3);
139	m_WTop = Wplus;
140	m_WAntitop = Wminus;
141	}
142
143
144	}
145
146
147	TTbarGen::~TTbarGen()
148	{
149	// default destructor, does nothing
150	}
151
152	GenParticle TTbarGen::Top()
153	{
154	return m_Top;
155	}
156
157	GenParticle TTbarGen::Antitop()
158	{
159	return m_Antitop;
160	}
161
162	GenParticle TTbarGen::WTop()
163	{
164	return m_WTop;
165	}
166
167	GenParticle TTbarGen::WAntitop()
168	{
169	return m_WAntitop;
170	}
171
172	GenParticle TTbarGen::bTop()
173	{
174	return m_bTop;
175	}
176
177	GenParticle TTbarGen::bAntitop()
178	{
179	return m_bAntitop;
180	}
181
182	GenParticle TTbarGen::Wdecay1()
183	{
184	return m_Wdecay1;
185	}
186
187	GenParticle TTbarGen::Wdecay2()
188	{
189	return m_Wdecay2;
190	}
191
192	GenParticle TTbarGen::WMinusdecay1()
193	{
194	return m_WMinusdecay1;
195	}
196
197	GenParticle TTbarGen::WMinusdecay2()
198	{
199	return m_WMinusdecay2;
200	}
201
202	TTbarGen::E_DecayChannel TTbarGen::DecayChannel()
203	{
204	if ((abs(m_pdgId1)==11 \|\| abs(m_pdgId1)==12) && (abs(m_pdgId2)==11 \|\| abs(m_pdgId2)==12)) m_type = e_ee;
205	if ((abs(m_pdgId1)==13 \|\| abs(m_pdgId1)==14) && (abs(m_pdgId2)==13 \|\| abs(m_pdgId2)==14)) m_type = e_mumu;
206	if ((abs(m_pdgId1)==15 \|\| abs(m_pdgId1)==16) && (abs(m_pdgId2)==15 \|\| abs(m_pdgId2)==16)) m_type = e_tautau;
207	if ((abs(m_pdgId1)<=5) && (abs(m_pdgId2)<=5)) m_type = e_had;
208	if ((abs(m_pdgId1)==11 \|\| abs(m_pdgId1)==12) && (abs(m_pdgId2)<=5)) m_type = e_ehad;
209	if ((abs(m_pdgId2)==11 \|\| abs(m_pdgId2)==12) && (abs(m_pdgId1)<=5)) m_type = e_ehad;
210	if ((abs(m_pdgId1)==13 \|\| abs(m_pdgId1)==14) && (abs(m_pdgId2)<=5)) m_type = e_muhad;
211	if ((abs(m_pdgId2)==13 \|\| abs(m_pdgId2)==14) && (abs(m_pdgId1)<=5)) m_type = e_muhad;
212	if ((abs(m_pdgId1)==15 \|\| abs(m_pdgId1)==16) && (abs(m_pdgId2)<=5)) m_type = e_tauhad;
213	if ((abs(m_pdgId2)==15 \|\| abs(m_pdgId2)==16) && (abs(m_pdgId1)<=5)) m_type = e_tauhad;
214	if ((abs(m_pdgId1)==11 \|\| abs(m_pdgId1)==12) && (abs(m_pdgId2)==13 \|\| abs(m_pdgId2)==14)) m_type = e_emu;
215	if ((abs(m_pdgId2)==11 \|\| abs(m_pdgId2)==12) && (abs(m_pdgId1)==13 \|\| abs(m_pdgId1)==14)) m_type = e_emu;
216	if ((abs(m_pdgId1)==11 \|\| abs(m_pdgId1)==12) && (abs(m_pdgId2)==15 \|\| abs(m_pdgId2)==16)) m_type = e_etau;
217	if ((abs(m_pdgId2)==11 \|\| abs(m_pdgId2)==12) && (abs(m_pdgId1)==15 \|\| abs(m_pdgId1)==16)) m_type = e_etau;
218	if ((abs(m_pdgId1)==13 \|\| abs(m_pdgId1)==14) && (abs(m_pdgId2)==15 \|\| abs(m_pdgId2)==16)) m_type = e_mutau;
219	if ((abs(m_pdgId2)==13 \|\| abs(m_pdgId2)==14) && (abs(m_pdgId1)==15 \|\| abs(m_pdgId1)==16)) m_type = e_mutau;
220	return m_type;
221	}
222
223	bool TTbarGen::IsTopHadronicDecay()
224	{
225	if (abs(m_pdgId1)<=5) return true;
226	else return false;
227	}
228
229	bool TTbarGen::IsAntiTopHadronicDecay()
230	{
231	if (abs(m_pdgId2)<=5) return true;
232	else return false;
233	}
234
235	GenParticle TTbarGen::ChargedLepton(){
236	GenParticle lepton;
237	if (m_type != e_ehad && m_type != e_muhad && m_type!= e_tauhad){
238	std::cerr << "This is no l+jets ttbar event, no charged lepton found" <<std::endl;
239	return lepton;
240	}
241
242	int nlepton=0;
243	if(abs(Wdecay1().pdgId())==11 \|\| abs(Wdecay1().pdgId())==13 \|\| abs(Wdecay1().pdgId())==15){
244	lepton = Wdecay1();
245	nlepton++;
246	}
247	if(abs(Wdecay2().pdgId())==11 \|\| abs(Wdecay2().pdgId())==13 \|\| abs(Wdecay2().pdgId())==15){
248	lepton = Wdecay2();
249	nlepton++;
250	}
251	if(abs(WMinusdecay1().pdgId())==11 \|\| abs(WMinusdecay1().pdgId())==13 \|\| abs(WMinusdecay1().pdgId())==15){
252	lepton = WMinusdecay1();
253	nlepton++;
254	}
255	if(abs(WMinusdecay2().pdgId())==11 \|\| abs(WMinusdecay2().pdgId())==13 \|\| abs(WMinusdecay2().pdgId())==15){
256	lepton = WMinusdecay2();
257	nlepton++;
258	}
259	if(nlepton!=1) std::cerr << "Not exactly one lepton found " << nlepton<< std::endl;
260
261	return lepton;
262
263	}
264
265	GenParticle TTbarGen::Neutrino(){
266	GenParticle neutrino;
267
268	if (m_type != e_ehad && m_type != e_muhad && m_type!= e_tauhad){
269	std::cerr << "This is no l+jets ttbar event, no neutrino found" <<std::endl;
270	return neutrino;
271	}
272
273	int nneutrino=0;
274	if(abs(Wdecay1().pdgId())==12 \|\| abs(Wdecay1().pdgId())==14 \|\| abs(Wdecay1().pdgId())==16){
275	neutrino = Wdecay1();
276	nneutrino++;
277	}
278	if(abs(Wdecay2().pdgId())==12 \|\| abs(Wdecay2().pdgId())==14 \|\| abs(Wdecay2().pdgId())==16){
279	neutrino = Wdecay2();
280	nneutrino++;
281	}
282	if(abs(WMinusdecay1().pdgId())==12 \|\| abs(WMinusdecay1().pdgId())==14 \|\| abs(WMinusdecay1().pdgId())==16){
283	neutrino = WMinusdecay1();
284	nneutrino++;
285	}
286	if(abs(WMinusdecay2().pdgId())==12 \|\| abs(WMinusdecay2().pdgId())==14 \|\| abs(WMinusdecay2().pdgId())==16){
287	neutrino = WMinusdecay2();
288	nneutrino++;
289	}
290	if(nneutrino!=1) std::cerr << "Not exactly one neutrino found " << nneutrino<< std::endl;
291
292	return neutrino;
293
294	}
295
296	GenParticle TTbarGen::TopLep(){
297
298	if(ChargedLepton().charge()>0) return Top();
299	else return Antitop();
300
301	}
302
303	GenParticle TTbarGen::TopHad(){
304
305	if(ChargedLepton().charge()<0) return Top();
306	else return Antitop();
307
308	}
309
310
311	GenParticle TTbarGen::BLep(){
312
313	if(ChargedLepton().charge()>0) return bTop();
314	else return bAntitop();
315
316	}
317
318	GenParticle TTbarGen::BHad(){
319
320	if(ChargedLepton().charge()<0) return bTop();
321	else return bAntitop();
322
323	}
324
325	GenParticle TTbarGen::WLep(){
326
327	if(ChargedLepton().charge()>0) return WTop();
328	else return WAntitop();
329
330	}
331
332	GenParticle TTbarGen::WHad(){
333
334	if(ChargedLepton().charge()<0) return WTop();
335	else return WAntitop();
336
337	}
338
339	GenParticle TTbarGen::Q1(){
340
341	if(ChargedLepton().charge()>0) return WMinusdecay1();
342	else return Wdecay1();
343
344	}
345
346	GenParticle TTbarGen::Q2(){
347
348	if(ChargedLepton().charge()>0) return WMinusdecay2();
349	else return Wdecay2();
350
351	}
352