SFrameTools/src/TTbarGen.cxx

#include "include/TTbarGen.h"

using namespace std;

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