VHbbDataFormats/interface/VHbbCandidateTools.h

#ifndef VHBBCANDIDATETOOLS_H
#define VHBBCANDIDATETOOLS_H

#include "VHbbAnalysis/VHbbDataFormats/interface/VHbbCandidate.h"

#include <iostream>

struct CompareJetPtMuons {
  bool operator()( const VHbbEvent::MuonInfo& j1, const  VHbbEvent::MuonInfo& j2 ) const {
    return j1.p4.Pt() > j2.p4.Pt();
  }
};
struct CompareJetPtElectrons {
  bool operator()( const VHbbEvent::ElectronInfo& j1, const  VHbbEvent::ElectronInfo& j2 ) const {
    return j1.p4.Pt() > j2.p4.Pt();
  }
};


class VHbbCandidateTools {
 public:

 VHbbCandidateTools(bool verbose = false): verbose_(verbose){}

  float deltaPhi(float in2, float in1){
    float dphi = in2-in1;
    if ( dphi > M_PI ) {
      dphi -= 2.0*M_PI;
    } else if ( dphi <= -M_PI ) {
      dphi += 2.0*M_PI;
    }
    return dphi;
  }
  
  
  VHbbCandidate getHZtaumuCandidate(const VHbbCandidate & in, bool & ok, std::vector<unsigned int>& muPos, std::vector<unsigned int>& tauPos){
    if (verbose_) std::cout <<" getHZtaumuCandidate input mu "<<in.V.muons.size()<<" e "<<in.V.electrons.size()<< " tau " << in.V.taus.size() << std::endl;
    ok = false;
    VHbbCandidate temp=in;
    // Require exactly one tau and muon, and no electrons
          std::vector<VHbbEvent::MuonInfo> Nmuons_ = temp.V.muons;
          std::vector<VHbbEvent::ElectronInfo> Nelectrons_ = temp.V.electrons;
int Nmuons = 0;
int Nelectrons = 0;
int postition = 0;
          for(unsigned int i = 0; i<temp.V.muons.size();i++){
                  if (Nmuons_[i].p4.Pt()>20) {
                  Nmuons++;
        postition = i;
                  }       }
          for(unsigned int j = 0; j<temp.V.electrons.size();j++){
                  if (Nelectrons_[j].p4.Pt()>20) Nelectrons++;
          }
          
    if (temp.V.taus.size()!=1) return in;
    if (Nmuons!=1) return in ;
    if (Nelectrons!=0) return in ;
    temp.V.p4 = temp.V.taus[0].p4 + temp.V.muons[postition].p4;
    temp.V.firstLepton = muPos[postition];
    temp.V.secondLepton = tauPos[0];
    
    ok = true;
    return temp;
    
  }
  VHbbCandidate getHWtaunCandidate(const VHbbCandidate & in, bool & ok , std::vector<unsigned int>& pos){
    
    if (verbose_){
      std::cout <<" getHWtaunCandidate input mu "<<in.V.muons.size()<<" e "<<in.V.electrons.size()<< " tau " << in.V.taus.size() << std::endl;
      std::cout << " pos.size()=" << pos.size() << std::endl;
    }
    
    ok = false;
    VHbbCandidate temp=in;
    // require a tau and no electrons or muons
          std::vector<VHbbEvent::MuonInfo> Nmuons_ = temp.V.muons;
          std::vector<VHbbEvent::ElectronInfo> Nelectrons_ = temp.V.electrons;
          int Nmuons = 0;
          int Nelectrons = 0;
          int postition = 0;
          for(unsigned int i = 0; i<temp.V.muons.size();i++){
                  if (Nmuons_[i].p4.Pt()>20) {
                          Nmuons++;
                          postition = i;
                  }       }
          for(unsigned int j = 0; j<temp.V.electrons.size();j++){
                  if (Nelectrons_[j].p4.Pt()>20) Nelectrons++;
          }
          
    if (temp.V.taus.size()!=1) return in;
    if (Nmuons!=0) return in ;
    if (Nelectrons!=0) return in ;
    if (temp.V.mets.size()<1) return in ;
    temp.V.p4 = temp.V.taus[0].p4+temp.V.mets[0].p4;
    temp.V.firstLepton=pos[0];
    
    if (verbose_){
      std::cout << "Done with getHWtaunCandidate" << std::endl;
    }
    
    ok=true;
    return temp;
  }
 
        VHbbCandidate getHZemuCandidate(const VHbbCandidate & in, bool & ok, std::vector<unsigned int>& pos){
                if (verbose_){
                        std::cout <<" getHZemuCandidate input mu "<<in.V.muons.size()<<" e "<<in.V.electrons.size()<<std::endl;
                }
                ok = false;
                VHbbCandidate temp=in;
                
                //
                // change: allow for additional leptons; by definition 
                //
                if (temp.V.muons.size()<1) return in ;
                if (temp.V.electrons.size()<1) return in ;
                if (temp.V.mets.size()<1) return in;
                
                //    if (temp.V.electrons.size()!=0) return in ;
                std::vector<VHbbEvent::MuonInfo> muons_ = temp.V.muons;
                std::vector<VHbbEvent::ElectronInfo> electrons_ = temp.V.electrons;
                
                //orthogonality with Zmumu and Zee
                if (temp.V.muons.size()>1){
                        if (muons_[0].p4.Pt()>20 && muons_[1].p4.Pt()>20 && (muons_[0].charge*muons_[1].charge)< 0)return in;
                }
                if (temp.V.electrons.size()>1){
                        if (electrons_[0].p4.Pt()>20 && electrons_[1].p4.Pt()>20&&(electrons_[0].charge*electrons_[1].charge)< 0)return in;
                }
                
                
                // beware: assumes already sorted!!!!
                
                //    CompareJetPtMuons ptComparator;
                //    std::sort(muons_.begin(), muons_.end(), ptComparator);
                if (muons_[0].p4.Pt()<10 || electrons_[0].p4.Pt()<10 ) return in;
                
                //
                // now I need to ask also for the charge
                //
                int tauEle=-99, tauMu =-99;
                float Tau_pT = -99.99;
                int tmpEle =0, tmpMu =0;
                int N_pairs = 0;
                //
                // i need to find a proper pair
                //
                
                for (unsigned int ele_it=0; ele_it< electrons_.size(); ++ele_it){
                        for (unsigned int mu_it=0; mu_it< muons_.size(); ++mu_it){
                                if (muons_[mu_it].p4.Pt()>10 && electrons_[ele_it].p4.Pt()>10 ){
                                        if (  electrons_[ele_it].charge * muons_[mu_it].charge < 0) {
                                                tmpEle = ele_it;
                                                tmpMu = mu_it;
                                                if (N_pairs ==0){
                                                        tauEle = ele_it;
                                                        tauMu = mu_it;
                                                }
                                                N_pairs++;
                                                temp.V.p4 = muons_[mu_it].p4+electrons_[ele_it].p4;
                                                if (temp.V.p4.Pt() > Tau_pT) {
                                                        tauEle = ele_it;
                                                        tauMu = mu_it;
                                                        Tau_pT = temp.V.p4.Pt();
                                                }
                                        } // if opposite charge
                                }//leption pt requirement
                        }  // muon for loop
                } //electron for loop
                if (tauEle==-99) return in;
                if (N_pairs>1) std::cout <<" Number emu pairs: "<<N_pairs<<std::endl;
                
                
                temp.V.p4 = muons_[tauMu].p4+electrons_[tauEle].p4;
                std::vector<VHbbEvent::MuonInfo> muons2_;
                for (std::vector<VHbbEvent::MuonInfo>::const_iterator it = muons_.begin(); it!= muons_.end(); ++it){
                        if (it->p4.Pt()>10) muons2_.push_back(*it);
                }
                temp.V.muons = muons2_;
                
                // beware; assumes already sorted
                
                //    CompareJetPtElectrons ptComparator2;
                //    std::sort(electrons_.begin(), electrons_.end(), ptComparator2);
                std::vector<VHbbEvent::ElectronInfo> electrons2_;
                for (std::vector<VHbbEvent::ElectronInfo>::const_iterator it = electrons_.begin(); it!= electrons_.end(); ++it){
                        if (it->p4.Pt()>10) electrons2_.push_back(*it);
                }
                temp.V.electrons = electrons2_;
                
                //
                // consider all 
                //
                
                
                //    if (temp.V.Pt()<150 ) return in;
                //    if (temp.H.Pt()<150) return in;
                //    if (temp.H.firstJet().csv< 0.9) return in;
                //    if (temp.H.secondJet().csv<0.5) return in;
                //    if (deltaPhi(temp.V.Phi(),temp.H.Phi())<2.7) return in;
                //    if (temp.V.FourMomentum.Mass()<75 || temp.V.FourMomentum.Mass()>105) return in; 
                //    if (temp.additionalJets.size()>0) return in;
                //    if (std::Abs(deltaTheta) ????
                
                temp.V.firstLepton = pos[tauEle];
                temp.V.secondLepton = pos[tauMu];
                ok = true;
                
                return temp;
        }  
        

  VHbbCandidate getHZmumuCandidate(const VHbbCandidate & in, bool & ok, std::vector<unsigned int>& pos){
    if (verbose_){
      std::cout <<" getHZmumuCandidate input mu "<<in.V.muons.size()<<" e "<<in.V.electrons.size()<<std::endl;
    }
    ok = false;
    VHbbCandidate temp=in;
    
    //
    // change: allow for additional leptons; by definition 
    //
          std::vector<VHbbEvent::MuonInfo> Nmuons_ = temp.V.muons;
          std::vector<VHbbEvent::ElectronInfo> Nelectrons_ = temp.V.electrons;
          int Nmuons = 0;
          int Nelectrons = 0;
          int postition = 0;
          for(unsigned int i = 0; i<temp.V.muons.size();i++){
                  if (Nmuons_[i].p4.Pt()>20) {
                          Nmuons++;
                          postition = i;
                  }       }
          for(unsigned int j = 0; j<temp.V.electrons.size();j++){
                  if (Nelectrons_[j].p4.Pt()>20) Nelectrons++;
          }
          
    if (Nmuons<2) return in ;
    //    if (temp.V.electrons.size()!=0) return in ;
    std::vector<VHbbEvent::MuonInfo> muons_ = temp.V.muons;

    // beware: assumes already sorted!!!!

    //    CompareJetPtMuons ptComparator;
    //    std::sort(muons_.begin(), muons_.end(), ptComparator);
    if (muons_[0].p4.Pt()<20 || muons_[1].p4.Pt()<20 ) return in;
    
    //
    // now I need to ask also for the charge
    //
    int selectMu2=1;
    
    if (muons_[0].charge* muons_[selectMu2].charge > 0){
      if (muons_.size() ==2) return in;
      //
      // i need to find a proper pair
      //
      
      for (unsigned int it=2; it< muons_.size(); ++it){
        if (  muons_[it].charge * muons_[0].charge < 0) {
          selectMu2 = it;
          break;
        }
        if (selectMu2 == 1) return in;
      }  
    }
    temp.V.p4 = muons_[0].p4+muons_[selectMu2].p4;
    std::vector<VHbbEvent::MuonInfo> muons2_;
    for (std::vector<VHbbEvent::MuonInfo>::const_iterator it = muons_.begin(); it!= muons_.end(); ++it){
      if (it->p4.Pt()>20) muons2_.push_back(*it);
    }
    temp.V.muons = muons2_;
    
    // the same for electrons
    std::vector<VHbbEvent::ElectronInfo> electrons_ = temp.V.electrons;

    // beware; assumes already sorted

    //    CompareJetPtElectrons ptComparator2;
    //    std::sort(electrons_.begin(), electrons_.end(), ptComparator2);
    std::vector<VHbbEvent::ElectronInfo> electrons2_;
    for (std::vector<VHbbEvent::ElectronInfo>::const_iterator it = electrons_.begin(); it!= electrons_.end(); ++it){
      if (it->p4.Pt()>20) electrons2_.push_back(*it);
    }
    temp.V.electrons = electrons2_;
    
    //
    // consider all 
    //
    
    
    //    if (temp.V.Pt()<150 ) return in;
    //    if (temp.H.Pt()<150) return in;
    //    if (temp.H.firstJet().csv< 0.9) return in;
    //    if (temp.H.secondJet().csv<0.5) return in;
    //    if (deltaPhi(temp.V.Phi(),temp.H.Phi())<2.7) return in;
    //    if (temp.V.FourMomentum.Mass()<75 || temp.V.FourMomentum.Mass()>105) return in; 
    //    if (temp.additionalJets.size()>0) return in;
    //    if (std::Abs(deltaTheta) ????

    temp.V.firstLepton = pos[0];
    temp.V.secondLepton = pos[selectMu2];
   ok = true;
    return temp;
  }  
  VHbbCandidate getHZeeCandidate(const VHbbCandidate & in, bool & ok, std::vector<unsigned int>& pos){
    if (verbose_){
      std::cout <<" getHZeeCandidate input mu "<<in.V.electrons.size()<<" e "<<in.V.muons.size()<<std::endl;
    }
    ok = false;
    VHbbCandidate temp=in;
    
    //
    // change: allow for additional leptons; by definition 
    //
          std::vector<VHbbEvent::MuonInfo> Nmuons_ = temp.V.muons;
          std::vector<VHbbEvent::ElectronInfo> Nelectrons_ = temp.V.electrons;
          int Nmuons = 0;
          int Nelectrons = 0;
          int postition = 0;
          for(unsigned int i = 0; i<temp.V.muons.size();i++){
                  if (Nmuons_[i].p4.Pt()>20) {
                          Nmuons++;
                          postition = i;
                  }       }
          for(unsigned int j = 0; j<temp.V.electrons.size();j++){
                  if (Nelectrons_[j].p4.Pt()>20) Nelectrons++;
          }
          
    if (Nelectrons<2) return in ;
    //    if (temp.V.electrons.size()!=0) return in ;
    std::vector<VHbbEvent::ElectronInfo> electrons_ = temp.V.electrons;

    // beware assumes already sorted

    //    CompareJetPtElectrons ptComparator;
    //    std::sort(electrons_.begin(), electrons_.end(), ptComparator);
    if (electrons_[0].p4.Pt()<20 || electrons_[1].p4.Pt()<20 ) return in;
    //
    // now I need to ask also for the charge
    //
    int selectE2=1;
    
    if (electrons_[0].charge* electrons_[selectE2].charge > 0){
      if (electrons_.size() ==2) return in;
      //
      // i need to find a proper pair
      //
      
      for (unsigned int it=2; it< electrons_.size(); ++it){
        if (  electrons_[it].charge * electrons_[0].charge < 0) {
          selectE2 = it;
          break;
        }
        if (selectE2 == 1) return in;
      }  
    }
   temp.V.p4 = electrons_[0].p4+electrons_[selectE2].p4;
    std::vector<VHbbEvent::ElectronInfo> electrons2_;
    for (std::vector<VHbbEvent::ElectronInfo>::const_iterator it = electrons_.begin(); it!= electrons_.end(); ++it){
      if (it->p4.Pt()>20) electrons2_.push_back(*it);
    }
    temp.V.electrons = electrons2_;
  
    // the same for muonss
    std::vector<VHbbEvent::MuonInfo> muons_ = temp.V.muons;

    // beware assumes already sorted

    //    CompareJetPtMuons ptComparator2;
    //    std::sort(muons_.begin(), muons_.end(), ptComparator2);
    std::vector<VHbbEvent::MuonInfo> muons2_;
    for (std::vector<VHbbEvent::MuonInfo>::const_iterator it = muons_.begin(); it!= muons_.end(); ++it){
      if (it->p4.Pt()>20) muons2_.push_back(*it);
    }
    temp.V.muons = muons2_;

    temp.V.firstLepton = pos[0];
    temp.V.secondLepton = pos[selectE2];

    
   ok = true;
    return temp;
  }  
  VHbbCandidate getHZnnCandidate(const VHbbCandidate & in, bool & ok){
    if (verbose_){
      std::cout <<" getHZnnCandidate input mu "<<in.V.muons.size()<<" e "<<in.V.electrons.size()<<" met "<<in.V.mets.size()<<std::endl;
    }
    
    ok = false;
    VHbbCandidate temp=in;
    if (temp.V.mets.size()!=1) return in;
//always build a NuNu candidate, exclusion come from if/else series on the caller side
// this allow to still build candidates for mu+e  and same sign dilept
//    if (temp.V.muons.size()!=0) return in ;
//    if (temp.V.electrons.size()!=0) return in ;
    
    temp.V.p4 = temp.V.mets[0].p4;
    if (verbose_) {
      std::cout <<" debug met "<< temp.V.mets[0].metSig << " " <<  temp.V.mets[0].sumEt<< std::endl;
    }   
//    if (temp.V.mets[0].metSig<5) return in;
    if (temp.V.mets[0].p4.Pt()<80) return in;
    //    if (temp.H.p4.Pt()<150)return in;
    //    if (temp.H.firstJet().csv< 0.9) return in;
    //    if (temp.H.secondJet().csv<0.5) return in;
    //    if (deltaPhi(temp.V.p4.Phi(),temp.H.p4.Phi())<2.95) return in;
    //    if (temp.V.electrons.size()>0 || temp.V.muons.size()>0 ) return in;
    //    if (temp.additionalJets.size()>0) return in;
    //    if (std::Abs(deltaTheta) ????
    
    ok = true;
    return temp;
  }

  VHbbCandidate getHWmunCandidate(const VHbbCandidate & in, bool & ok , std::vector<unsigned int>& pos){
    ok = false;
    VHbbCandidate temp=in;
    // require a muon and no electrons
          std::vector<VHbbEvent::MuonInfo> Nmuons_ = temp.V.muons;
          std::vector<VHbbEvent::ElectronInfo> Nelectrons_ = temp.V.electrons;
          int Nmuons = 0;
          int Nelectrons = 0;
          int postition = 0;
          for(unsigned int i = 0; i<temp.V.muons.size();i++){
                  if (Nmuons_[i].p4.Pt()>20) {
                          Nmuons++;
                          postition = i;
                  }       }
          for(unsigned int j = 0; j<temp.V.electrons.size();j++){
                  if (Nelectrons_[j].p4.Pt()>20) Nelectrons++;
          }
          
    if (Nmuons!=1) return in ;
    if (Nelectrons!=0) return in ;
    if (temp.V.mets.size()<1) return in ;
        std::vector<VHbbEvent::MuonInfo> muons_ = temp.V.muons;
        if (muons_[postition].p4.Pt()<20) return in;
    //
    /*pT(W) > 150 GeV (pt(W) computed using lepton px,y and PF MET x and y components)
      pT(H) > 150 GeV
      best btag, CSV > 0.90
      second-best btag, CSV > 0.50
      Dphi(W,H) > 2.95
      no additional isolated leptons (pT > 15 GeV)
      
      same lepton definition as in the physics objects section 
      
      No additional ak5PFjets (pT > 30 GeV; |eta| < 2.4)
      MET>35. for the electron BDT analysis
      |cos(theta * )| (TBO)
      color flow pull angle (TBO)
      We don't cut on the transverse mass (for boosted objects cutting on the transverse mass introduces an inefficiency due to the angle between the MET and the lepton being close to 0.) 
    */
    
    temp.V.p4 = temp.V.muons[postition].p4+temp.V.mets[0].p4;
    temp.V.firstLepton=pos[0];

    ok=true;
    return temp;
  }

  VHbbCandidate getHWenCandidate(const VHbbCandidate & in, bool & ok, std::vector<unsigned int>& pos){
    ok = false;
    VHbbCandidate temp=in;
          std::vector<VHbbEvent::MuonInfo> Nmuons_ = temp.V.muons;
          std::vector<VHbbEvent::ElectronInfo> Nelectrons_ = temp.V.electrons;
          int Nmuons = 0;
          int Nelectrons = 0;
          int postition = 0;
          for(unsigned int i = 0; i<temp.V.muons.size();i++){
                  if (Nmuons_[i].p4.Pt()>20) {
                          Nmuons++;
                         
                  }       }
          for(unsigned int j = 0; j<temp.V.electrons.size();j++){
                  if (Nelectrons_[j].p4.Pt()>20) {
                  Nelectrons++;
                   postition = j;
                  } }
          
    if (Nelectrons!=1) return in ;
    if (Nmuons!=0) return in ;
    if (temp.V.mets.size()<1) return in ;
        std::vector<VHbbEvent::ElectronInfo> electrons_ = temp.V.electrons;
        if (electrons_[postition].p4.Pt()<20) return in;
    //
    /*pT(W) > 150 GeV (pt(W) computed using lepton px,y and PF MET x and y components)
      pT(H) > 150 GeV
      best btag, CSV > 0.90
      second-best btag, CSV > 0.50
      Dphi(W,H) > 2.95
      no additional isolated leptons (pT > 15 GeV)
      
      same lepton definition as in the physics objects section 
      
      No additional ak5PFjets (pT > 30 GeV; |eta| < 2.4)
      MET>35. for the electron BDT analysis
      |cos(theta * )| (TBO)
      color flow pull angle (TBO)
      We don't cut on the transverse mass (for boosted objects cutting on the transverse mass introduces an inefficiency due to the angle between the MET and the lepton being close to 0.) 
    */
    
    ok=true;
    temp.V.firstLepton=pos[0];
    temp.V.p4 = temp.V.electrons[0].p4+temp.V.mets[0].p4;
    return temp;
  }

 public: 
  bool verbose_;
  
 public:
  static float getDeltaTheta( const VHbbEvent::SimpleJet & j1, const VHbbEvent::SimpleJet & j2 ) {
    double deltaTheta = 1e10;
    TLorentzVector pi(0,0,0,0);
    TLorentzVector v_j1 = j1.chargedTracksFourMomentum;
    TLorentzVector v_j2 = j2.chargedTracksFourMomentum;
    
    if( v_j2.Mag() == 0 
        || v_j1.Mag() == 0 )
      return deltaTheta = 1e10;
    
    //use j1 to calculate the pull vector
    TVector2 t = j1.tVector;
    
    if( t.Mod() == 0 )
      return deltaTheta = 1e10;
    
    Double_t dphi =  v_j2.Phi()- v_j1.Phi();
    if ( dphi > M_PI ) {
      dphi -= 2.0*M_PI;
    } else if ( dphi <= -M_PI ) {
      dphi += 2.0*M_PI;
    }
    Double_t deltaeta = v_j2.Rapidity() - v_j1.Rapidity();
    TVector2 BBdir( deltaeta, dphi );
    
    deltaTheta = t.DeltaPhi(BBdir);
    
    return deltaTheta;
    
  }
  
  float getHelicity(const VHbbEvent::SimpleJet& j, TVector3 boost) const {
    double hel = 1e10;
    TLorentzVector jet = j.p4;
    jet.Boost( -boost );
    hel = TMath::Cos( jet.Vect().Angle( boost ) );
    return hel;
  }


};


#endif


Revision:	1.13.2.1
Committed:	Fri Jul 27 15:40:08 2012 UTC (12 years, 9 months ago) by wilken
Content type:	text/plain
Branch:	V21emuCandidate
Changes since 1.13:	+210 -13 lines
Log Message:	adding emu candidate
#	User	Rev	Content
1	tboccali	1.1	#ifndef VHBBCANDIDATETOOLS_H
2			#define VHBBCANDIDATETOOLS_H
3
4			#include "VHbbAnalysis/VHbbDataFormats/interface/VHbbCandidate.h"
5
6			#include <iostream>
7
8	tboccali	1.7	struct CompareJetPtMuons {
9			bool operator()( const VHbbEvent::MuonInfo& j1, const VHbbEvent::MuonInfo& j2 ) const {
10			return j1.p4.Pt() > j2.p4.Pt();
11			}
12			};
13			struct CompareJetPtElectrons {
14			bool operator()( const VHbbEvent::ElectronInfo& j1, const VHbbEvent::ElectronInfo& j2 ) const {
15			return j1.p4.Pt() > j2.p4.Pt();
16			}
17			};
18
19
20
21	tboccali	1.1	class VHbbCandidateTools {
22			public:
23
24			VHbbCandidateTools(bool verbose = false): verbose_(verbose){}
25
26			float deltaPhi(float in2, float in1){
27			float dphi = in2-in1;
28			if ( dphi > M_PI ) {
29			dphi -= 2.0*M_PI;
30			} else if ( dphi <= -M_PI ) {
31			dphi += 2.0*M_PI;
32			}
33			return dphi;
34			}
35	tboccali	1.9
36
37	sethzenz	1.13	VHbbCandidate getHZtaumuCandidate(const VHbbCandidate & in, bool & ok, std::vector<unsigned int>& muPos, std::vector<unsigned int>& tauPos){
38			if (verbose_) std::cout <<" getHZtaumuCandidate input mu "<<in.V.muons.size()<<" e "<<in.V.electrons.size()<< " tau " << in.V.taus.size() << std::endl;
39			ok = false;
40			VHbbCandidate temp=in;
41			// Require exactly one tau and muon, and no electrons
42	wilken	1.13.2.1	std::vector<VHbbEvent::MuonInfo> Nmuons_ = temp.V.muons;
43			std::vector<VHbbEvent::ElectronInfo> Nelectrons_ = temp.V.electrons;
44			int Nmuons = 0;
45			int Nelectrons = 0;
46			int postition = 0;
47			for(unsigned int i = 0; i<temp.V.muons.size();i++){
48			if (Nmuons_[i].p4.Pt()>20) {
49			Nmuons++;
50			postition = i;
51			} }
52			for(unsigned int j = 0; j<temp.V.electrons.size();j++){
53			if (Nelectrons_[j].p4.Pt()>20) Nelectrons++;
54			}
55
56	sethzenz	1.13	if (temp.V.taus.size()!=1) return in;
57	wilken	1.13.2.1	if (Nmuons!=1) return in ;
58			if (Nelectrons!=0) return in ;
59			temp.V.p4 = temp.V.taus[0].p4 + temp.V.muons[postition].p4;
60			temp.V.firstLepton = muPos[postition];
61	sethzenz	1.13	temp.V.secondLepton = tauPos[0];
62
63			ok = true;
64			return temp;
65
66			}
67			VHbbCandidate getHWtaunCandidate(const VHbbCandidate & in, bool & ok , std::vector<unsigned int>& pos){
68
69			if (verbose_){
70			std::cout <<" getHWtaunCandidate input mu "<<in.V.muons.size()<<" e "<<in.V.electrons.size()<< " tau " << in.V.taus.size() << std::endl;
71			std::cout << " pos.size()=" << pos.size() << std::endl;
72			}
73
74			ok = false;
75			VHbbCandidate temp=in;
76			// require a tau and no electrons or muons
77	wilken	1.13.2.1	std::vector<VHbbEvent::MuonInfo> Nmuons_ = temp.V.muons;
78			std::vector<VHbbEvent::ElectronInfo> Nelectrons_ = temp.V.electrons;
79			int Nmuons = 0;
80			int Nelectrons = 0;
81			int postition = 0;
82			for(unsigned int i = 0; i<temp.V.muons.size();i++){
83			if (Nmuons_[i].p4.Pt()>20) {
84			Nmuons++;
85			postition = i;
86			} }
87			for(unsigned int j = 0; j<temp.V.electrons.size();j++){
88			if (Nelectrons_[j].p4.Pt()>20) Nelectrons++;
89			}
90
91	sethzenz	1.13	if (temp.V.taus.size()!=1) return in;
92	wilken	1.13.2.1	if (Nmuons!=0) return in ;
93			if (Nelectrons!=0) return in ;
94	sethzenz	1.13	if (temp.V.mets.size()<1) return in ;
95			temp.V.p4 = temp.V.taus[0].p4+temp.V.mets[0].p4;
96			temp.V.firstLepton=pos[0];
97
98			if (verbose_){
99			std::cout << "Done with getHWtaunCandidate" << std::endl;
100			}
101
102			ok=true;
103			return temp;
104			}
105
106	wilken	1.13.2.1	VHbbCandidate getHZemuCandidate(const VHbbCandidate & in, bool & ok, std::vector<unsigned int>& pos){
107			if (verbose_){
108			std::cout <<" getHZemuCandidate input mu "<<in.V.muons.size()<<" e "<<in.V.electrons.size()<<std::endl;
109			}
110			ok = false;
111			VHbbCandidate temp=in;
112
113			//
114			// change: allow for additional leptons; by definition
115			//
116			if (temp.V.muons.size()<1) return in ;
117			if (temp.V.electrons.size()<1) return in ;
118			if (temp.V.mets.size()<1) return in;
119
120			// if (temp.V.electrons.size()!=0) return in ;
121			std::vector<VHbbEvent::MuonInfo> muons_ = temp.V.muons;
122			std::vector<VHbbEvent::ElectronInfo> electrons_ = temp.V.electrons;
123
124			//orthogonality with Zmumu and Zee
125			if (temp.V.muons.size()>1){
126			if (muons_[0].p4.Pt()>20 && muons_[1].p4.Pt()>20 && (muons_[0].charge*muons_[1].charge)< 0)return in;
127			}
128			if (temp.V.electrons.size()>1){
129			if (electrons_[0].p4.Pt()>20 && electrons_[1].p4.Pt()>20&&(electrons_[0].charge*electrons_[1].charge)< 0)return in;
130			}
131
132
133			// beware: assumes already sorted!!!!
134
135			// CompareJetPtMuons ptComparator;
136			// std::sort(muons_.begin(), muons_.end(), ptComparator);
137			if (muons_[0].p4.Pt()<10 \|\| electrons_[0].p4.Pt()<10 ) return in;
138
139			//
140			// now I need to ask also for the charge
141			//
142			int tauEle=-99, tauMu =-99;
143			float Tau_pT = -99.99;
144			int tmpEle =0, tmpMu =0;
145			int N_pairs = 0;
146			//
147			// i need to find a proper pair
148			//
149
150			for (unsigned int ele_it=0; ele_it< electrons_.size(); ++ele_it){
151			for (unsigned int mu_it=0; mu_it< muons_.size(); ++mu_it){
152			if (muons_[mu_it].p4.Pt()>10 && electrons_[ele_it].p4.Pt()>10 ){
153			if ( electrons_[ele_it].charge * muons_[mu_it].charge < 0) {
154			tmpEle = ele_it;
155			tmpMu = mu_it;
156			if (N_pairs ==0){
157			tauEle = ele_it;
158			tauMu = mu_it;
159			}
160			N_pairs++;
161			temp.V.p4 = muons_[mu_it].p4+electrons_[ele_it].p4;
162			if (temp.V.p4.Pt() > Tau_pT) {
163			tauEle = ele_it;
164			tauMu = mu_it;
165			Tau_pT = temp.V.p4.Pt();
166			}
167			} // if opposite charge
168			}//leption pt requirement
169			} // muon for loop
170			} //electron for loop
171			if (tauEle==-99) return in;
172			if (N_pairs>1) std::cout <<" Number emu pairs: "<<N_pairs<<std::endl;
173
174
175			temp.V.p4 = muons_[tauMu].p4+electrons_[tauEle].p4;
176			std::vector<VHbbEvent::MuonInfo> muons2_;
177			for (std::vector<VHbbEvent::MuonInfo>::const_iterator it = muons_.begin(); it!= muons_.end(); ++it){
178			if (it->p4.Pt()>10) muons2_.push_back(*it);
179			}
180			temp.V.muons = muons2_;
181
182			// beware; assumes already sorted
183
184			// CompareJetPtElectrons ptComparator2;
185			// std::sort(electrons_.begin(), electrons_.end(), ptComparator2);
186			std::vector<VHbbEvent::ElectronInfo> electrons2_;
187			for (std::vector<VHbbEvent::ElectronInfo>::const_iterator it = electrons_.begin(); it!= electrons_.end(); ++it){
188			if (it->p4.Pt()>10) electrons2_.push_back(*it);
189			}
190			temp.V.electrons = electrons2_;
191
192			//
193			// consider all
194			//
195
196
197			// if (temp.V.Pt()<150 ) return in;
198			// if (temp.H.Pt()<150) return in;
199			// if (temp.H.firstJet().csv< 0.9) return in;
200			// if (temp.H.secondJet().csv<0.5) return in;
201			// if (deltaPhi(temp.V.Phi(),temp.H.Phi())<2.7) return in;
202			// if (temp.V.FourMomentum.Mass()<75 \|\| temp.V.FourMomentum.Mass()>105) return in;
203			// if (temp.additionalJets.size()>0) return in;
204			// if (std::Abs(deltaTheta) ????
205
206			temp.V.firstLepton = pos[tauEle];
207			temp.V.secondLepton = pos[tauMu];
208			ok = true;
209
210			return temp;
211			}
212
213	sethzenz	1.13
214	tboccali	1.9
215			VHbbCandidate getHZmumuCandidate(const VHbbCandidate & in, bool & ok, std::vector<unsigned int>& pos){
216	tboccali	1.1	if (verbose_){
217			std::cout <<" getHZmumuCandidate input mu "<<in.V.muons.size()<<" e "<<in.V.electrons.size()<<std::endl;
218			}
219			ok = false;
220			VHbbCandidate temp=in;
221	tboccali	1.9
222	tboccali	1.7	//
223			// change: allow for additional leptons; by definition
224			//
225	wilken	1.13.2.1	std::vector<VHbbEvent::MuonInfo> Nmuons_ = temp.V.muons;
226			std::vector<VHbbEvent::ElectronInfo> Nelectrons_ = temp.V.electrons;
227			int Nmuons = 0;
228			int Nelectrons = 0;
229			int postition = 0;
230			for(unsigned int i = 0; i<temp.V.muons.size();i++){
231			if (Nmuons_[i].p4.Pt()>20) {
232			Nmuons++;
233			postition = i;
234			} }
235			for(unsigned int j = 0; j<temp.V.electrons.size();j++){
236			if (Nelectrons_[j].p4.Pt()>20) Nelectrons++;
237			}
238
239			if (Nmuons<2) return in ;
240	tboccali	1.7	// if (temp.V.electrons.size()!=0) return in ;
241			std::vector<VHbbEvent::MuonInfo> muons_ = temp.V.muons;
242	tboccali	1.9
243			// beware: assumes already sorted!!!!
244
245			// CompareJetPtMuons ptComparator;
246			// std::sort(muons_.begin(), muons_.end(), ptComparator);
247	tboccali	1.7	if (muons_[0].p4.Pt()<20 \|\| muons_[1].p4.Pt()<20 ) return in;
248	tboccali	1.9
249			//
250			// now I need to ask also for the charge
251			//
252			int selectMu2=1;
253
254			if (muons_[0].charge* muons_[selectMu2].charge > 0){
255			if (muons_.size() ==2) return in;
256			//
257			// i need to find a proper pair
258			//
259
260			for (unsigned int it=2; it< muons_.size(); ++it){
261			if ( muons_[it].charge * muons_[0].charge < 0) {
262			selectMu2 = it;
263			break;
264			}
265			if (selectMu2 == 1) return in;
266			}
267			}
268			temp.V.p4 = muons_[0].p4+muons_[selectMu2].p4;
269	tboccali	1.7	std::vector<VHbbEvent::MuonInfo> muons2_;
270			for (std::vector<VHbbEvent::MuonInfo>::const_iterator it = muons_.begin(); it!= muons_.end(); ++it){
271			if (it->p4.Pt()>20) muons2_.push_back(*it);
272			}
273			temp.V.muons = muons2_;
274	tboccali	1.9
275	tboccali	1.7	// the same for electrons
276			std::vector<VHbbEvent::ElectronInfo> electrons_ = temp.V.electrons;
277	tboccali	1.9
278			// beware; assumes already sorted
279
280			// CompareJetPtElectrons ptComparator2;
281			// std::sort(electrons_.begin(), electrons_.end(), ptComparator2);
282	tboccali	1.7	std::vector<VHbbEvent::ElectronInfo> electrons2_;
283			for (std::vector<VHbbEvent::ElectronInfo>::const_iterator it = electrons_.begin(); it!= electrons_.end(); ++it){
284			if (it->p4.Pt()>20) electrons2_.push_back(*it);
285			}
286			temp.V.electrons = electrons2_;
287
288	tboccali	1.9	//
289	tboccali	1.7	// consider all
290			//
291
292	tboccali	1.1
293			// if (temp.V.Pt()<150 ) return in;
294			// if (temp.H.Pt()<150) return in;
295			// if (temp.H.firstJet().csv< 0.9) return in;
296			// if (temp.H.secondJet().csv<0.5) return in;
297			// if (deltaPhi(temp.V.Phi(),temp.H.Phi())<2.7) return in;
298			// if (temp.V.FourMomentum.Mass()<75 \|\| temp.V.FourMomentum.Mass()>105) return in;
299			// if (temp.additionalJets.size()>0) return in;
300			// if (std::Abs(deltaTheta) ????
301	tboccali	1.9
302			temp.V.firstLepton = pos[0];
303			temp.V.secondLepton = pos[selectMu2];
304			ok = true;
305	tboccali	1.1	return temp;
306			}
307	tboccali	1.9	VHbbCandidate getHZeeCandidate(const VHbbCandidate & in, bool & ok, std::vector<unsigned int>& pos){
308	tboccali	1.1	if (verbose_){
309	tboccali	1.7	std::cout <<" getHZeeCandidate input mu "<<in.V.electrons.size()<<" e "<<in.V.muons.size()<<std::endl;
310	tboccali	1.1	}
311			ok = false;
312			VHbbCandidate temp=in;
313	tboccali	1.9
314	tboccali	1.1	//
315	tboccali	1.7	// change: allow for additional leptons; by definition
316	tboccali	1.1	//
317	wilken	1.13.2.1	std::vector<VHbbEvent::MuonInfo> Nmuons_ = temp.V.muons;
318			std::vector<VHbbEvent::ElectronInfo> Nelectrons_ = temp.V.electrons;
319			int Nmuons = 0;
320			int Nelectrons = 0;
321			int postition = 0;
322			for(unsigned int i = 0; i<temp.V.muons.size();i++){
323			if (Nmuons_[i].p4.Pt()>20) {
324			Nmuons++;
325			postition = i;
326			} }
327			for(unsigned int j = 0; j<temp.V.electrons.size();j++){
328			if (Nelectrons_[j].p4.Pt()>20) Nelectrons++;
329			}
330
331			if (Nelectrons<2) return in ;
332	tboccali	1.7	// if (temp.V.electrons.size()!=0) return in ;
333			std::vector<VHbbEvent::ElectronInfo> electrons_ = temp.V.electrons;
334	tboccali	1.9
335			// beware assumes already sorted
336
337			// CompareJetPtElectrons ptComparator;
338			// std::sort(electrons_.begin(), electrons_.end(), ptComparator);
339	tboccali	1.7	if (electrons_[0].p4.Pt()<20 \|\| electrons_[1].p4.Pt()<20 ) return in;
340	tboccali	1.9	//
341			// now I need to ask also for the charge
342			//
343			int selectE2=1;
344
345			if (electrons_[0].charge* electrons_[selectE2].charge > 0){
346			if (electrons_.size() ==2) return in;
347			//
348			// i need to find a proper pair
349			//
350
351			for (unsigned int it=2; it< electrons_.size(); ++it){
352			if ( electrons_[it].charge * electrons_[0].charge < 0) {
353			selectE2 = it;
354			break;
355			}
356			if (selectE2 == 1) return in;
357			}
358			}
359			temp.V.p4 = electrons_[0].p4+electrons_[selectE2].p4;
360	tboccali	1.7	std::vector<VHbbEvent::ElectronInfo> electrons2_;
361			for (std::vector<VHbbEvent::ElectronInfo>::const_iterator it = electrons_.begin(); it!= electrons_.end(); ++it){
362			if (it->p4.Pt()>20) electrons2_.push_back(*it);
363			}
364			temp.V.electrons = electrons2_;
365
366			// the same for muonss
367			std::vector<VHbbEvent::MuonInfo> muons_ = temp.V.muons;
368	tboccali	1.9
369			// beware assumes already sorted
370
371			// CompareJetPtMuons ptComparator2;
372			// std::sort(muons_.begin(), muons_.end(), ptComparator2);
373	tboccali	1.7	std::vector<VHbbEvent::MuonInfo> muons2_;
374			for (std::vector<VHbbEvent::MuonInfo>::const_iterator it = muons_.begin(); it!= muons_.end(); ++it){
375			if (it->p4.Pt()>20) muons2_.push_back(*it);
376			}
377			temp.V.muons = muons2_;
378	tboccali	1.9
379			temp.V.firstLepton = pos[0];
380			temp.V.secondLepton = pos[selectE2];
381
382	tboccali	1.7
383			ok = true;
384	tboccali	1.1	return temp;
385			}
386			VHbbCandidate getHZnnCandidate(const VHbbCandidate & in, bool & ok){
387			if (verbose_){
388			std::cout <<" getHZnnCandidate input mu "<<in.V.muons.size()<<" e "<<in.V.electrons.size()<<" met "<<in.V.mets.size()<<std::endl;
389			}
390
391			ok = false;
392			VHbbCandidate temp=in;
393			if (temp.V.mets.size()!=1) return in;
394	arizzi	1.10	//always build a NuNu candidate, exclusion come from if/else series on the caller side
395			// this allow to still build candidates for mu+e and same sign dilept
396			// if (temp.V.muons.size()!=0) return in ;
397			// if (temp.V.electrons.size()!=0) return in ;
398	tboccali	1.1
399	tboccali	1.6	temp.V.p4 = temp.V.mets[0].p4;
400	tboccali	1.1	if (verbose_) {
401			std::cout <<" debug met "<< temp.V.mets[0].metSig << " " << temp.V.mets[0].sumEt<< std::endl;
402			}
403	arizzi	1.12	// if (temp.V.mets[0].metSig<5) return in;
404	tboccali	1.11	if (temp.V.mets[0].p4.Pt()<80) return in;
405	tboccali	1.6	// if (temp.H.p4.Pt()<150)return in;
406	tboccali	1.1	// if (temp.H.firstJet().csv< 0.9) return in;
407			// if (temp.H.secondJet().csv<0.5) return in;
408	tboccali	1.6	// if (deltaPhi(temp.V.p4.Phi(),temp.H.p4.Phi())<2.95) return in;
409	tboccali	1.1	// if (temp.V.electrons.size()>0 \|\| temp.V.muons.size()>0 ) return in;
410			// if (temp.additionalJets.size()>0) return in;
411			// if (std::Abs(deltaTheta) ????
412
413			ok = true;
414			return temp;
415			}
416
417	tboccali	1.9	VHbbCandidate getHWmunCandidate(const VHbbCandidate & in, bool & ok , std::vector<unsigned int>& pos){
418	tboccali	1.1	ok = false;
419			VHbbCandidate temp=in;
420	tboccali	1.2	// require a muon and no electrons
421	wilken	1.13.2.1	std::vector<VHbbEvent::MuonInfo> Nmuons_ = temp.V.muons;
422			std::vector<VHbbEvent::ElectronInfo> Nelectrons_ = temp.V.electrons;
423			int Nmuons = 0;
424			int Nelectrons = 0;
425			int postition = 0;
426			for(unsigned int i = 0; i<temp.V.muons.size();i++){
427			if (Nmuons_[i].p4.Pt()>20) {
428			Nmuons++;
429			postition = i;
430			} }
431			for(unsigned int j = 0; j<temp.V.electrons.size();j++){
432			if (Nelectrons_[j].p4.Pt()>20) Nelectrons++;
433			}
434
435			if (Nmuons!=1) return in ;
436			if (Nelectrons!=0) return in ;
437	tboccali	1.3	if (temp.V.mets.size()<1) return in ;
438	wilken	1.13.2.1	std::vector<VHbbEvent::MuonInfo> muons_ = temp.V.muons;
439			if (muons_[postition].p4.Pt()<20) return in;
440	tboccali	1.1	//
441	tboccali	1.2	/*pT(W) > 150 GeV (pt(W) computed using lepton px,y and PF MET x and y components)
442			pT(H) > 150 GeV
443			best btag, CSV > 0.90
444			second-best btag, CSV > 0.50
445			Dphi(W,H) > 2.95
446			no additional isolated leptons (pT > 15 GeV)
447
448			same lepton definition as in the physics objects section
449
450			No additional ak5PFjets (pT > 30 GeV; \|eta\| < 2.4)
451			MET>35. for the electron BDT analysis
452			\|cos(theta * )\| (TBO)
453			color flow pull angle (TBO)
454			We don't cut on the transverse mass (for boosted objects cutting on the transverse mass introduces an inefficiency due to the angle between the MET and the lepton being close to 0.)
455			*/
456
457	wilken	1.13.2.1	temp.V.p4 = temp.V.muons[postition].p4+temp.V.mets[0].p4;
458	tboccali	1.9	temp.V.firstLepton=pos[0];
459	tboccali	1.2
460			ok=true;
461	tboccali	1.3	return temp;
462	tboccali	1.1	}
463
464	tboccali	1.9	VHbbCandidate getHWenCandidate(const VHbbCandidate & in, bool & ok, std::vector<unsigned int>& pos){
465	tboccali	1.1	ok = false;
466			VHbbCandidate temp=in;
467	wilken	1.13.2.1	std::vector<VHbbEvent::MuonInfo> Nmuons_ = temp.V.muons;
468			std::vector<VHbbEvent::ElectronInfo> Nelectrons_ = temp.V.electrons;
469			int Nmuons = 0;
470			int Nelectrons = 0;
471			int postition = 0;
472			for(unsigned int i = 0; i<temp.V.muons.size();i++){
473			if (Nmuons_[i].p4.Pt()>20) {
474			Nmuons++;
475
476			} }
477			for(unsigned int j = 0; j<temp.V.electrons.size();j++){
478			if (Nelectrons_[j].p4.Pt()>20) {
479			Nelectrons++;
480			postition = j;
481			} }
482
483			if (Nelectrons!=1) return in ;
484			if (Nmuons!=0) return in ;
485	tboccali	1.3	if (temp.V.mets.size()<1) return in ;
486	wilken	1.13.2.1	std::vector<VHbbEvent::ElectronInfo> electrons_ = temp.V.electrons;
487			if (electrons_[postition].p4.Pt()<20) return in;
488	tboccali	1.1	//
489	tboccali	1.2	/*pT(W) > 150 GeV (pt(W) computed using lepton px,y and PF MET x and y components)
490			pT(H) > 150 GeV
491			best btag, CSV > 0.90
492			second-best btag, CSV > 0.50
493			Dphi(W,H) > 2.95
494			no additional isolated leptons (pT > 15 GeV)
495
496			same lepton definition as in the physics objects section
497
498			No additional ak5PFjets (pT > 30 GeV; \|eta\| < 2.4)
499			MET>35. for the electron BDT analysis
500			\|cos(theta * )\| (TBO)
501			color flow pull angle (TBO)
502			We don't cut on the transverse mass (for boosted objects cutting on the transverse mass introduces an inefficiency due to the angle between the MET and the lepton being close to 0.)
503			*/
504
505			ok=true;
506	tboccali	1.9	temp.V.firstLepton=pos[0];
507	tboccali	1.6	temp.V.p4 = temp.V.electrons[0].p4+temp.V.mets[0].p4;
508	tboccali	1.3	return temp;
509	tboccali	1.1	}
510
511			public:
512			bool verbose_;
513	tboccali	1.2
514			public:
515	tboccali	1.9	static float getDeltaTheta( const VHbbEvent::SimpleJet & j1, const VHbbEvent::SimpleJet & j2 ) {
516	tboccali	1.5	double deltaTheta = 1e10;
517			TLorentzVector pi(0,0,0,0);
518			TLorentzVector v_j1 = j1.chargedTracksFourMomentum;
519			TLorentzVector v_j2 = j2.chargedTracksFourMomentum;
520
521			if( v_j2.Mag() == 0
522			\|\| v_j1.Mag() == 0 )
523			return deltaTheta = 1e10;
524
525			//use j1 to calculate the pull vector
526			TVector2 t = j1.tVector;
527
528			if( t.Mod() == 0 )
529			return deltaTheta = 1e10;
530
531			Double_t dphi = v_j2.Phi()- v_j1.Phi();
532			if ( dphi > M_PI ) {
533			dphi -= 2.0*M_PI;
534			} else if ( dphi <= -M_PI ) {
535			dphi += 2.0*M_PI;
536			}
537			Double_t deltaeta = v_j2.Rapidity() - v_j1.Rapidity();
538			TVector2 BBdir( deltaeta, dphi );
539
540			deltaTheta = t.DeltaPhi(BBdir);
541
542			return deltaTheta;
543
544			}
545
546	tboccali	1.4	float getHelicity(const VHbbEvent::SimpleJet& j, TVector3 boost) const {
547			double hel = 1e10;
548	tboccali	1.6	TLorentzVector jet = j.p4;
549	tboccali	1.4	jet.Boost( -boost );
550			hel = TMath::Cos( jet.Vect().Angle( boost ) );
551			return hel;
552			}
553
554	tboccali	1.1
555			};
556
557
558
559			#endif
560
561
562
563
564
565