LeptonSelection/src/MuonSelection.cc

#include <math.h>
#include <iostream>

#include "MuonSelection.h"
#include "ParseArgs.h"
#include "SelectionStatus.h"

#include "Track.h"
#include "MuonTools.h"
#include "MuonQuality.h"
#include "MuonIDMVA.h"

mithep::MuonIDMVA * muIDMVA;
mithep::MuonTools muTools;

extern Float_t computePFMuonIso(const mithep::Muon *mu, 
                        const mithep::Vertex & vtx,
                        const mithep::Array<mithep::PFCandidate> * pfCandidates,
                        const Double_t dRMax);

//--------------------------------------------------------------------------------------------------
SelectionStatus muonDummyVeto(ControlFlags &ctrl, 
                              const mithep::Muon *muon, 
                              const mithep::Vertex &vtx) 
//--------------------------------------------------------------------------------------------------
{ 
  SelectionStatus status;
  status.setStatus(SelectionStatus::PRESELECTION);
  return status;  
}

//--------------------------------------------------------------------------------------------------
SelectionStatus muonCutBasedVeto(ControlFlags &ctrl, 
                                 const mithep::Muon *muon, 
                                 const mithep::Vertex &vtx) 
//--------------------------------------------------------------------------------------------------
{ 
  //
  // Loose cut-based ID for isolation veto
  //
  bool ret = true;
  
  if(!(muon->IsGlobalMuon() || muon->IsTrackerMuon())) ret = false;
  if( muon->NValidHits() < 11 )                        ret = false;
  if( fabs(muon->Ip3dPVSignificance()) >= 4 )          ret = false;

  SelectionStatus status;
  if( ret ) status.setStatus(SelectionStatus::PRESELECTION);
  return status;
}


//--------------------------------------------------------------------------------------------------
SelectionStatus noPreSelection( ControlFlags &ctrl,  const mithep::Muon * mu ) 

//--------------------------------------------------------------------------------------------------
{
        SelectionStatus status;
        status.setStatus(SelectionStatus::PRESELECTION); 
        if(ctrl.debug) cout << "muon presel returning status : " << status.getStatus() << endl;
        return status;
}

//--------------------------------------------------------------------------------------------------
SelectionStatus muonPreSelection( ControlFlags &ctrl,  
                                  const mithep::Muon * mu, 
                                  const mithep::Vertex & vtx, 
                                  const mithep::Array<mithep::PFCandidate> * pfCandidates ) 
//--------------------------------------------------------------------------------------------------
{ 
  bool ret = true; 
  if(ctrl.debug) cout << "inside muonpresel ... " << endl;
  //  bool ret = isMuFO(mu,vtx,pfCandidates);
  if(ctrl.debug) cout << "iS fo? ... " << ret << endl;  
  ret &= ( fabs(mu->Ip3dPVSignificance()) < 4 ); 
  if(ctrl.debug) cout << "and pass IP (" << mu->Ip3dPVSignificance() << ") ? ... " << ret << endl;  
  ret &= ( mu->Pt() >= 5 ); 
  if(ctrl.debug) cout << "and >5 GeV (" << mu->Pt() << ") ? ... " << ret << endl;  
  ret &= ( fabs(mu->Eta()) <= 2.4 ); 
  if(ctrl.debug) cout << "and < 2.4 eta (" << mu->Eta() << ")? ... " << ret << endl;  
  ret &= ( mu->IsTrackerMuon() || mu->IsGlobalMuon() ); 
  if(ctrl.debug) cout << "is Tracker or Global? ... " << ret << endl;  
  if( mu->IsTrackerMuon() && mu->HasTrackerTrk() && !mu->IsGlobalMuon() )  
    ret &=  ( mu->Quality().QualityMask().Mask() & mithep::MuonQuality::AllArbitrated);
  if(ctrl.debug) cout << "if isTrackerMuon, arbitrated ? ... " << ret << endl;  
  ret &= ( fabs(mu->BestTrk()->DzCorrected(vtx)) < 0.1 );   
  if( ctrl.debug ) cout << "elepresel : ret after dZcorr : " << ret <<  endl;


  SelectionStatus status;
  if( ret ) status.setStatus(SelectionStatus::PRESELECTION);
  if(ctrl.debug) cout << "muon presel returning status : " << status.getStatus() << endl;  
  return status;
} 


//--------------------------------------------------------------------------------------------------
SelectionStatus muonReferencePreSelection( ControlFlags &ctrl,  
                                           const mithep::Muon * mu, 
                                           const mithep::Vertex & vtx, 
                                           const mithep::Array<mithep::PFCandidate> * pfCandidates ) 
//--------------------------------------------------------------------------------------------------
{ 
  bool ret = true; 
  if(ctrl.debug) cout << "inside muonpresel ... " << endl;
  ret &= ( fabs(mu->Ip3dPVSignificance()) < 100 ); 
  if(ctrl.debug) cout << "and pass IP (" << mu->Ip3dPVSignificance() << ") ? ... " << ret << endl;  
  ret &= ( mu->Pt() >= 5 ); 
  if(ctrl.debug) cout << "and >5 GeV (" << mu->Pt() << ") ? ... " << ret << endl;  
  ret &= ( fabs(mu->Eta()) <= 2.4 ); 
  if(ctrl.debug) cout << "and < 2.4 eta (" << mu->Eta() << ")? ... " << ret << endl;  
  ret &= ( mu->IsTrackerMuon() || mu->IsGlobalMuon() ); 
  if(ctrl.debug) cout << "is Tracker or Global? ... " << ret << endl;  
  //  if( mu->IsTrackerMuon() && mu->HasTrackerTrk() && !mu->IsGlobalMuon() )  
  //    ret &=  ( mu->Quality().QualityMask().Mask() & mithep::MuonQuality::AllArbitrated);
  //  if(ctrl.debug) cout << "if isTrackerMuon, arbitrated ? ... " << ret << endl;  
  //  ret &= ( mu->BestTrk()->DzCorrected(vtx) < 0.1);   
  //  if( ctrl.debug ) cout << "elepresel : ret after dZcorr : " << ret <<  endl;


  SelectionStatus status;
  if( ret ) status.setStatus(SelectionStatus::PRESELECTION);
  if(ctrl.debug) cout << "muon presel returning status : " << status.getStatus() << endl;  
  return status;
} 


//--------------------------------------------------------------------------------------------------
SelectionStatus muonIDMVASelection(ControlFlags &ctrl, 
                                 const mithep::Muon *mu,
                                 const mithep::Vertex & vtx   )
//--------------------------------------------------------------------------------------------------
{
  double global_rchi2 = (mu->IsGlobalMuon()) ? mu->GlobalTrk()->RChi2() : 0.;

  // 
  // WARNING !!! KH hacked out to sync
  //
  /*
  double mvaval = muIDMVA->MVAValue_ID(mu->Pt(),
                                    mu->Eta(),
                                    mu->IsGlobalMuon(),
                                    mu->IsTrackerMuon(),
                                    mu->TrackerTrk()->RChi2(),
                                    global_rchi2,
                                    (Double_t)(mu->NValidHits()),
                                    (Double_t)(mu->TrackerTrk()->NHits()),
                                    (Double_t)(mu->TrackerTrk()->NPixelHits()),
                                    (Double_t)(mu->NMatches()),
                                    mu->TrkKink(),
                                    muTools.GetSegmentCompatability(mu),
                                    muTools.GetCaloCompatability(mu,true,true),
                                    mu->HadEnergy(),
                                    mu->EmEnergy(),
                                    mu->HadS9Energy(),
                                    mu->EmS9Energy(),
                                    (Bool_t)(ctrl.debug) );
  */
  double mvaval = 0.0;

  SelectionStatus status;
  bool pass = false;
  
  if( mu->IsGlobalMuon() && mu->IsTrackerMuon() 
      && fabs(mu->Eta()) <= 1.5 && mu->Pt() <= 10 && mvaval >= MUON_IDMVA_CUT_BIN0)  pass = true;
  else if( mu->IsGlobalMuon() && mu->IsTrackerMuon()
           && fabs(mu->Eta()) <= 1.5 && mu->Pt() > 10 && mvaval >= MUON_IDMVA_CUT_BIN1)  pass = true;
  else if( mu->IsGlobalMuon() && mu->IsTrackerMuon() 
           && fabs(mu->Eta()) > 1.5 && mu->Pt() <= 10 && mvaval >= MUON_IDMVA_CUT_BIN2)  pass = true;
  else if( mu->IsGlobalMuon() && mu->IsTrackerMuon() 
           && fabs(mu->Eta()) > 1.5 && mu->Pt() > 10 && mvaval >= MUON_IDMVA_CUT_BIN3)  pass = true;
  else if( !(mu->IsGlobalMuon()) && mu->IsTrackerMuon() && mvaval >= MUON_IDMVA_CUT_BIN4 )  pass = true;
  else if( mu->IsGlobalMuon() && !(mu->IsTrackerMuon()) && mvaval >= MUON_IDMVA_CUT_BIN5 )  pass = true;
   

  if( pass ) {
    status.orStatus(SelectionStatus::LOOSEID);
    status.orStatus(SelectionStatus::TIGHTID);
  }

  return status;
}

//--------------------------------------------------------------------------------------------------
void initMuonIDMVA() 
//--------------------------------------------------------------------------------------------------
{ 
  muIDMVA = new mithep::MuonIDMVA();
  vector<string> weightFiles;
  weightFiles.push_back("./data/MuonIDMVAWeights/MuonIDMVA_BDTG_barrel_lowpt_V2.weights.xml");
  weightFiles.push_back("./data/MuonIDMVAWeights/MuonIDMVA_BDTG_barrel_highpt_V2.weights.xml");
  weightFiles.push_back("./data/MuonIDMVAWeights/MuonIDMVA_BDTG_endcap_lowpt_V2.weights.xml");
  weightFiles.push_back("./data/MuonIDMVAWeights/MuonIDMVA_BDTG_endcap_highpt_V2.weights.xml");
  weightFiles.push_back("./data/MuonIDMVAWeights/MuonIDMVA_BDTG_tracker_V2.weights.xml");
  weightFiles.push_back("./data/MuonIDMVAWeights/MuonIDMVA_BDTG_global_V2.weights.xml");
  muIDMVA->Initialize( "MuonIDMVA",
                       mithep::MuonIDMVA::kIDV0,
                       kTRUE, weightFiles);
}

//--------------------------------------------------------------------------------------------------
SelectionStatus muonIDPFSelection(ControlFlags &ctrl, 
                                  const mithep::Muon *mu,
                                  const mithep::Vertex & vtx,
                                  const mithep::Array<mithep::PFCandidate> * pfCandidates )
//--------------------------------------------------------------------------------------------------
{
  bool pass = false;
  SelectionStatus status; // init'ed to FAIL

  // check that it matches to a PF muon
  for( int i=0; i<pfCandidates->GetEntries(); i++ ) {
    const mithep::PFCandidate *pf = (mithep::PFCandidate*)((*pfCandidates)[i]);
    //    if( pf->HasTrackerTrk() ) {
    //      if( (pf->TrackerTrk() == mu->TrackerTrk()) && abs(pf->PFType()) == mithep::PFCandidate::eMuon ) {
    if(  abs(pf->PFType()) == mithep::PFCandidate::eMuon ) {
      pass = true;
      break;
      //    }
    }
  }

  //
  // NB : ipsig here for reference synch ...
  //
  if( pass && fabs(mu->Ip3dPVSignificance()) < 4 ) {
    status.orStatus(SelectionStatus::LOOSEID);
    status.orStatus(SelectionStatus::TIGHTID);
  }

  return status;
}
Revision:	1.17
Committed:	Thu May 10 22:36:12 2012 UTC (13 years ago) by khahn
Content type:	text/plain
Branch:	MAIN
Changes since 1.16:	+6 -5 lines
Log Message:	fscking sync
#	User	Rev	Content
1	khahn	1.1	#include <math.h>
2			#include <iostream>
3
4			#include "MuonSelection.h"
5			#include "ParseArgs.h"
6			#include "SelectionStatus.h"
7
8	khahn	1.4	#include "Track.h"
9			#include "MuonTools.h"
10			#include "MuonQuality.h"
11			#include "MuonIDMVA.h"
12
13			mithep::MuonIDMVA * muIDMVA;
14			mithep::MuonTools muTools;
15
16			extern Float_t computePFMuonIso(const mithep::Muon *mu,
17			const mithep::Vertex & vtx,
18			const mithep::Array<mithep::PFCandidate> * pfCandidates,
19			const Double_t dRMax);
20
21			//--------------------------------------------------------------------------------------------------
22	khahn	1.8	SelectionStatus muonDummyVeto(ControlFlags &ctrl,
23			const mithep::Muon *muon,
24			const mithep::Vertex &vtx)
25			//--------------------------------------------------------------------------------------------------
26			{
27			SelectionStatus status;
28			status.setStatus(SelectionStatus::PRESELECTION);
29			return status;
30			}
31
32			//--------------------------------------------------------------------------------------------------
33			SelectionStatus muonCutBasedVeto(ControlFlags &ctrl,
34			const mithep::Muon *muon,
35			const mithep::Vertex &vtx)
36			//--------------------------------------------------------------------------------------------------
37			{
38			//
39			// Loose cut-based ID for isolation veto
40			//
41			bool ret = true;
42
43			if(!(muon->IsGlobalMuon() \|\| muon->IsTrackerMuon())) ret = false;
44			if( muon->NValidHits() < 11 ) ret = false;
45			if( fabs(muon->Ip3dPVSignificance()) >= 4 ) ret = false;
46
47			SelectionStatus status;
48			if( ret ) status.setStatus(SelectionStatus::PRESELECTION);
49			return status;
50			}
51
52
53			//--------------------------------------------------------------------------------------------------
54	khahn	1.4	SelectionStatus noPreSelection( ControlFlags &ctrl, const mithep::Muon * mu )
55
56			//--------------------------------------------------------------------------------------------------
57			{
58	anlevin	1.3	SelectionStatus status;
59			status.setStatus(SelectionStatus::PRESELECTION);
60			if(ctrl.debug) cout << "muon presel returning status : " << status.getStatus() << endl;
61			return status;
62			}
63	khahn	1.1
64	khahn	1.4	//--------------------------------------------------------------------------------------------------
65			SelectionStatus muonPreSelection( ControlFlags &ctrl,
66			const mithep::Muon * mu,
67			const mithep::Vertex & vtx,
68			const mithep::Array<mithep::PFCandidate> * pfCandidates )
69			//--------------------------------------------------------------------------------------------------
70			{
71			bool ret = true;
72	khahn	1.1	if(ctrl.debug) cout << "inside muonpresel ... " << endl;
73	khahn	1.4	// bool ret = isMuFO(mu,vtx,pfCandidates);
74	khahn	1.1	if(ctrl.debug) cout << "iS fo? ... " << ret << endl;
75	khahn	1.4	ret &= ( fabs(mu->Ip3dPVSignificance()) < 4 );
76			if(ctrl.debug) cout << "and pass IP (" << mu->Ip3dPVSignificance() << ") ? ... " << ret << endl;
77	khahn	1.14	ret &= ( mu->Pt() >= 5 );
78	khahn	1.7	if(ctrl.debug) cout << "and >5 GeV (" << mu->Pt() << ") ? ... " << ret << endl;
79			ret &= ( fabs(mu->Eta()) <= 2.4 );
80			if(ctrl.debug) cout << "and < 2.4 eta (" << mu->Eta() << ")? ... " << ret << endl;
81	khahn	1.11	ret &= ( mu->IsTrackerMuon() \|\| mu->IsGlobalMuon() );
82			if(ctrl.debug) cout << "is Tracker or Global? ... " << ret << endl;
83			if( mu->IsTrackerMuon() && mu->HasTrackerTrk() && !mu->IsGlobalMuon() )
84			ret &= ( mu->Quality().QualityMask().Mask() & mithep::MuonQuality::AllArbitrated);
85			if(ctrl.debug) cout << "if isTrackerMuon, arbitrated ? ... " << ret << endl;
86	khahn	1.15	ret &= ( fabs(mu->BestTrk()->DzCorrected(vtx)) < 0.1 );
87	khahn	1.10	if( ctrl.debug ) cout << "elepresel : ret after dZcorr : " << ret << endl;
88
89	khahn	1.1
90			SelectionStatus status;
91			if( ret ) status.setStatus(SelectionStatus::PRESELECTION);
92			if(ctrl.debug) cout << "muon presel returning status : " << status.getStatus() << endl;
93			return status;
94			}
95
96	khahn	1.4
97	khahn	1.1
98	khahn	1.13	//--------------------------------------------------------------------------------------------------
99			SelectionStatus muonReferencePreSelection( ControlFlags &ctrl,
100			const mithep::Muon * mu,
101			const mithep::Vertex & vtx,
102			const mithep::Array<mithep::PFCandidate> * pfCandidates )
103			//--------------------------------------------------------------------------------------------------
104			{
105			bool ret = true;
106			if(ctrl.debug) cout << "inside muonpresel ... " << endl;
107	khahn	1.16	ret &= ( fabs(mu->Ip3dPVSignificance()) < 100 );
108	khahn	1.13	if(ctrl.debug) cout << "and pass IP (" << mu->Ip3dPVSignificance() << ") ? ... " << ret << endl;
109	khahn	1.14	ret &= ( mu->Pt() >= 5 );
110	khahn	1.13	if(ctrl.debug) cout << "and >5 GeV (" << mu->Pt() << ") ? ... " << ret << endl;
111			ret &= ( fabs(mu->Eta()) <= 2.4 );
112			if(ctrl.debug) cout << "and < 2.4 eta (" << mu->Eta() << ")? ... " << ret << endl;
113			ret &= ( mu->IsTrackerMuon() \|\| mu->IsGlobalMuon() );
114			if(ctrl.debug) cout << "is Tracker or Global? ... " << ret << endl;
115			// if( mu->IsTrackerMuon() && mu->HasTrackerTrk() && !mu->IsGlobalMuon() )
116			// ret &= ( mu->Quality().QualityMask().Mask() & mithep::MuonQuality::AllArbitrated);
117			// if(ctrl.debug) cout << "if isTrackerMuon, arbitrated ? ... " << ret << endl;
118			// ret &= ( mu->BestTrk()->DzCorrected(vtx) < 0.1);
119			// if( ctrl.debug ) cout << "elepresel : ret after dZcorr : " << ret << endl;
120
121
122			SelectionStatus status;
123			if( ret ) status.setStatus(SelectionStatus::PRESELECTION);
124			if(ctrl.debug) cout << "muon presel returning status : " << status.getStatus() << endl;
125			return status;
126			}
127
128
129
130	anlevin	1.3
131	khahn	1.4	//--------------------------------------------------------------------------------------------------
132			SelectionStatus muonIDMVASelection(ControlFlags &ctrl,
133			const mithep::Muon *mu,
134			const mithep::Vertex & vtx )
135			//--------------------------------------------------------------------------------------------------
136			{
137			double global_rchi2 = (mu->IsGlobalMuon()) ? mu->GlobalTrk()->RChi2() : 0.;
138	anlevin	1.3
139	khahn	1.13	//
140			// WARNING !!! KH hacked out to sync
141			//
142			/*
143	khahn	1.4	double mvaval = muIDMVA->MVAValue_ID(mu->Pt(),
144			mu->Eta(),
145			mu->IsGlobalMuon(),
146			mu->IsTrackerMuon(),
147			mu->TrackerTrk()->RChi2(),
148			global_rchi2,
149			(Double_t)(mu->NValidHits()),
150	khahn	1.11	(Double_t)(mu->TrackerTrk()->NHits()),
151	khahn	1.4	(Double_t)(mu->TrackerTrk()->NPixelHits()),
152			(Double_t)(mu->NMatches()),
153			mu->TrkKink(),
154			muTools.GetSegmentCompatability(mu),
155			muTools.GetCaloCompatability(mu,true,true),
156			mu->HadEnergy(),
157			mu->EmEnergy(),
158			mu->HadS9Energy(),
159			mu->EmS9Energy(),
160			(Bool_t)(ctrl.debug) );
161	khahn	1.13	*/
162			double mvaval = 0.0;
163	anlevin	1.3
164	khahn	1.4	SelectionStatus status;
165			bool pass = false;
166
167			if( mu->IsGlobalMuon() && mu->IsTrackerMuon()
168	khahn	1.7	&& fabs(mu->Eta()) <= 1.5 && mu->Pt() <= 10 && mvaval >= MUON_IDMVA_CUT_BIN0) pass = true;
169	khahn	1.4	else if( mu->IsGlobalMuon() && mu->IsTrackerMuon()
170	khahn	1.7	&& fabs(mu->Eta()) <= 1.5 && mu->Pt() > 10 && mvaval >= MUON_IDMVA_CUT_BIN1) pass = true;
171	khahn	1.4	else if( mu->IsGlobalMuon() && mu->IsTrackerMuon()
172	khahn	1.7	&& fabs(mu->Eta()) > 1.5 && mu->Pt() <= 10 && mvaval >= MUON_IDMVA_CUT_BIN2) pass = true;
173	khahn	1.4	else if( mu->IsGlobalMuon() && mu->IsTrackerMuon()
174			&& fabs(mu->Eta()) > 1.5 && mu->Pt() > 10 && mvaval >= MUON_IDMVA_CUT_BIN3) pass = true;
175	khahn	1.11	else if( !(mu->IsGlobalMuon()) && mu->IsTrackerMuon() && mvaval >= MUON_IDMVA_CUT_BIN4 ) pass = true;
176			else if( mu->IsGlobalMuon() && !(mu->IsTrackerMuon()) && mvaval >= MUON_IDMVA_CUT_BIN5 ) pass = true;
177
178	khahn	1.4
179
180			if( pass ) {
181			status.orStatus(SelectionStatus::LOOSEID);
182			status.orStatus(SelectionStatus::TIGHTID);
183	anlevin	1.3	}
184
185	khahn	1.4	return status;
186			}
187	anlevin	1.3
188	khahn	1.4	//--------------------------------------------------------------------------------------------------
189			void initMuonIDMVA()
190			//--------------------------------------------------------------------------------------------------
191			{
192			muIDMVA = new mithep::MuonIDMVA();
193			vector<string> weightFiles;
194	khahn	1.6	weightFiles.push_back("./data/MuonIDMVAWeights/MuonIDMVA_BDTG_barrel_lowpt_V2.weights.xml");
195			weightFiles.push_back("./data/MuonIDMVAWeights/MuonIDMVA_BDTG_barrel_highpt_V2.weights.xml");
196			weightFiles.push_back("./data/MuonIDMVAWeights/MuonIDMVA_BDTG_endcap_lowpt_V2.weights.xml");
197			weightFiles.push_back("./data/MuonIDMVAWeights/MuonIDMVA_BDTG_endcap_highpt_V2.weights.xml");
198			weightFiles.push_back("./data/MuonIDMVAWeights/MuonIDMVA_BDTG_tracker_V2.weights.xml");
199			weightFiles.push_back("./data/MuonIDMVAWeights/MuonIDMVA_BDTG_global_V2.weights.xml");
200	khahn	1.4	muIDMVA->Initialize( "MuonIDMVA",
201			mithep::MuonIDMVA::kIDV0,
202			kTRUE, weightFiles);
203			}
204	khahn	1.9
205			//--------------------------------------------------------------------------------------------------
206			SelectionStatus muonIDPFSelection(ControlFlags &ctrl,
207			const mithep::Muon *mu,
208			const mithep::Vertex & vtx,
209			const mithep::Array<mithep::PFCandidate> * pfCandidates )
210			//--------------------------------------------------------------------------------------------------
211			{
212			bool pass = false;
213			SelectionStatus status; // init'ed to FAIL
214
215			// check that it matches to a PF muon
216			for( int i=0; i<pfCandidates->GetEntries(); i++ ) {
217			const mithep::PFCandidate pf = (mithep::PFCandidate)((*pfCandidates)[i]);
218	khahn	1.17	// if( pf->HasTrackerTrk() ) {
219			// if( (pf->TrackerTrk() == mu->TrackerTrk()) && abs(pf->PFType()) == mithep::PFCandidate::eMuon ) {
220			if( abs(pf->PFType()) == mithep::PFCandidate::eMuon ) {
221			pass = true;
222			break;
223			// }
224	khahn	1.9	}
225			}
226
227	khahn	1.13	//
228			// NB : ipsig here for reference synch ...
229			//
230	khahn	1.16	if( pass && fabs(mu->Ip3dPVSignificance()) < 4 ) {
231	khahn	1.9	status.orStatus(SelectionStatus::LOOSEID);
232			status.orStatus(SelectionStatus::TIGHTID);
233			}
234
235			return status;
236			}