NonMCBackground/src/SelectionEMU.cc

#include "SelectionStatus.h" 
#include "EventData.h" 
#include "SimpleLepton.h" 
#include "EfficiencyWeightsInterface.h"

#include "HZZBDTElectronSelection.h"
#include "IsolationSelection.h"
#include "PassHLT.h"
#include "SelectionEMU.h"

#include "ExternData.h"
#include "SelectionDefs.h"


// +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
// +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
EventData apply_EMU_selection(ControlFlags &ctrl,           // input control 
                              mithep::TEventInfo *info,     // input event info
                              TClonesArray *electronArr,    // input electrons
                              SelectionStatus (*ElectronPreSelector)( ControlFlags &, const mithep::TElectron*),
                              SelectionStatus (*ElectronIDSelector)( ControlFlags &, const mithep::TElectron*),
                              SelectionStatus (*ElectronIsoSelector)( ControlFlags &, const mithep::TElectron*),
                              TClonesArray *muonArr,         // input muons
                              SelectionStatus (*MuonPreSelector)( ControlFlags &, const mithep::TMuon*),
                              SelectionStatus (*MuonIDSelector)( ControlFlags &, const mithep::TMuon*),
                              SelectionStatus (*MuonIsoSelector)( ControlFlags &, const mithep::TMuon*) )
// +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
// +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
{       

  EventData ret;
  unsigned evtfail = 0x0;
  TRandom3 r;
  
  if( ctrl.debug ) {
    cout << "Run: " << info->runNum 
         << "\tEvt: " << info->evtNum
         << "\tLumi: " << info->lumiSec
         << endl;
  }

  if( !ctrl.mc ) {
    // not accounting for overlap atm
    RunLumiRangeMap::RunLumiPairType rl(info->runNum, info->lumiSec);      
    if( !(rlrm.HasRunLumi(rl)) )  {
      if( ctrl.debug ) cout << "\tfails JSON" << endl;
      ret.status.setStatus(0);
      return ret;
    }
  }
  
  
  //********************************************************
  // Trigger
  //********************************************************
//   if( !ctrl.mc ) { 
//     if( !(passHLTSingleMuon(info->triggerBits) )  ) { 
//       if( ctrl.debug ) cout << "\tfails trigger" << endl;
//       evtfail |= (1<<EVTFAIL_TRIGGER);
//       ret.status.setStatus(0);
//       return ret;
//     }        
//   } 
  if( ctrl.debug ) { 
    cout << "presel nlep: " << muonArr->GetEntries() + electronArr->GetEntries()
         << "\tnmuon: "    << muonArr->GetEntries()
         << "\tnelectron: " << electronArr->GetEntries()
         << endl;
  }

  //********************************************************
  // Lepton Selection
  //********************************************************
  vector<SimpleLepton> lepvec;
  if( muonArr->GetEntries() != 1 && electronArr->GetEntries() != 1 ) { 
    ret.status.setStatus(0);
    return ret;
  }

  //    
  if( ctrl.debug ) cout << "\tnMuons: " << muonArr->GetEntries() << endl;
  //----------------------------------------------------
  for(Int_t i=0; i<muonArr->GetEntries(); i++) 
    {
      const mithep::TMuon *mu = (mithep::TMuon*)((*muonArr)[i]);      
      
      SelectionStatus musel;
      if(ctrl.debug) cout << "musel.status  before anything: " << musel.getStatus() << endl;
      musel |= (*MuonPreSelector)(ctrl,mu);
      if(ctrl.debug) cout << "musel.status  after presel: " << musel.getStatus() << endl;
      musel |= (*MuonIDSelector)(ctrl,mu);
      if(ctrl.debug) cout << "musel.status  after ID: " << musel.getStatus() << endl;
      musel |= (*MuonIsoSelector)(ctrl,mu);
      if(ctrl.debug) cout << "musel.status  after iso: " << musel.getStatus() << endl;

      if( ctrl.debug ) { 
        cout << "muon:: pt: " << mu->pt 
             << "\teta: " << mu->eta 
             << "\tisorel: " <<  mu->pfIso03/mu->pt
             << "\tstatus: " << hex << musel.getStatus() << dec
             << endl;
      }
      
      if ( musel.pass() ) {

        SimpleLepton tmplep;
        float pt = mu->pt;
        tmplep.vecorig->SetPtEtaPhiM(pt, 
                                    mu->eta, 
                                    mu->phi,
                                    MUON_MASS);
        
        if( ctrl.do_escale_up ) { 
          pt=scale_smear_muon_Up(pt, 1,  r);
        }
        if( ctrl.do_escale_down ) { 
          pt=scale_smear_muon_Down(pt, 1,  r);
        }
        
        tmplep.vec->SetPtEtaPhiM(pt, 
                                mu->eta, 
                                mu->phi,
                                MUON_MASS);
        
        tmplep.type    = 13;
        tmplep.index   = i;
        tmplep.charge  = mu->q;
        tmplep.isoTrk  = mu->trkIso03;
        tmplep.isoEcal = mu->emIso03;
        tmplep.isoHcal = mu->hadIso03;
        tmplep.isoPF03 = mu->pfIso03;
        tmplep.isoPF04 = mu->pfIso04;
        tmplep.ip3dSig = mu->ip3dSig;
        tmplep.is4l    = false;
        tmplep.isEB    = (fabs(mu->eta) < 1.479 ? 1 : 0 ); 
        tmplep.isTight = musel.tight();
        tmplep.isLoose = musel.loose();
        lepvec.push_back(tmplep);
        if( ctrl.debug ) { cout << "muon passes ... " << endl;}
      }
      //  }
    }    
  
  
    //
    if( ctrl.debug ) { cout << "\tnElectron: " << electronArr->GetEntries() << endl; }
    // --------------------------------------------------------------------------------
    for(Int_t i=0; i<electronArr->GetEntries(); i++) 
      {
        const mithep::TElectron *ele = (mithep::TElectron*)((*electronArr)[i]);
        
        Bool_t isMuonOverlap = kFALSE;
        for (int k=0; k<lepvec.size(); ++k) {
          TVector3 tmplep;
          tmplep.SetPtEtaPhi(ele->pt, ele->eta, ele->phi);
          if ( abs(lepvec[k].type == 13) && lepvec[k].vec->Vect().DrEtaPhi(tmplep) < 0.15 ) {
            if( ctrl.debug ) cout << "-----> isMuonOverlap! " << endl;
            isMuonOverlap = kTRUE; 
            break;
          }
        }
        
        SelectionStatus elesel;
        elesel |= (*ElectronPreSelector)(ctrl,ele);
        elesel |= (*ElectronIDSelector)(ctrl,ele);
        elesel |= (*ElectronIsoSelector)(ctrl,ele);
        if( elesel.getStatus() & SelectionStatus::PRESELECTION ) 
          elesel.setStatus(SelectionStatus::LOOSESELECTION);

        if( ctrl.debug ){ 
          cout << "\tscEt: " << ele->scEt 
               << "\tscEta: " << ele->scEta
               << "\tncluster: " << ele->ncluster
               << "\tisorel: " <<  ele->pfIso04/ele->pt
               << "\tstatus: " << hex << elesel.getStatus() << dec
               << endl;
        }

        if ( elesel.pass() && !isMuonOverlap ) 
          {
            SimpleLepton tmplep;
            
            float pt = ele->pt;
            tmplep.vecorig->SetPtEtaPhiM( pt,
                                          ele->eta,
                                          ele->phi,
                                          ELECTRON_MASS );
            
            if( ctrl.do_escale ) { 
              pt=scale_smear_electron(pt, ele->isEB, r);
            }
            if( ctrl.do_escale_up ) { 
              pt=scale_smear_electron_Up(pt, ele->isEB,  r);
            }
            if( ctrl.do_escale_down ) { 
              pt=scale_smear_electron_Down(pt, ele->isEB,  r);
            }

            
            tmplep.vec->SetPtEtaPhiM( pt,
                                      ele->eta,
                                      ele->phi,
                                      ELECTRON_MASS );
            
            tmplep.type    = 11;
            tmplep.index   = i;
            tmplep.charge  = ele->q;
            tmplep.isoTrk  = ele->trkIso03;
            tmplep.isoEcal = ele->emIso03;
            tmplep.isoHcal = ele->hadIso03;
            tmplep.isoPF03 = ele->pfIso03;
            tmplep.isoPF04 = ele->pfIso04;
            tmplep.ip3dSig = ele->ip3dSig;
            tmplep.is4l    = false;
            tmplep.isEB    = ele->isEB;
            tmplep.scID    = ele->scID;
            tmplep.isTight = elesel.tight();
            tmplep.isLoose = elesel.pass();
            lepvec.push_back(tmplep);
            if( ctrl.debug ) { cout << "\telectron passes ... " << endl; }
          }
      }
    
    
    //********************************************************
    // Dump Stuff
    //********************************************************
    sort( lepvec.begin(), lepvec.end(), SimpleLepton::lep_pt_sort );
    int nmu=0, nele=0;
    for( int i=0; i<lepvec.size(); i++ ) {
      if(ctrl.debug) cout << "lepvec :: index: " << i 
                          << "\tpt: " << lepvec[i].vec->Pt()
                          << "\ttype: " << lepvec[i].type 
                          << endl;
      if( abs(lepvec[i].type) == 11 ) nele++;
      else nmu++;
    }
    if( ctrl.debug ) { 
      cout << "postsel nlep: " << lepvec.size()
           << "\tnmuon: " << nmu
           << "\tnelectron: " << nele
           << endl;
    }
    
    
    //******************************************************************************
    // Selection
    //******************************************************************************
    float bestMass=-1; int best_mu_index=-1, best_ele_index=-1;
    for(int i = 0; i<lepvec.size(); i++) {      // get a tight muon
      if( abs(lepvec[i].type) != 13 ) continue; 
      if( !(lepvec[i].isTight) ) continue; 
      //      if( lepvec[i].vec->Pt() < 35 ) continue; 
      if( ctrl.debug ) cout << "got a muon, index: " << i << endl;

      for(int j = 0; j<lepvec.size(); j++) {     // get a loose electron
        if( j == i ) continue; 
        if( abs(lepvec[j].type) != 11 ) continue; 
        if( !(lepvec[j].isLoose) ) continue; 
        if (lepvec[i].charge == lepvec[j].charge) continue;          
        //      if (lepvec[i].charge != lepvec[j].charge) continue;          
        if( ctrl.debug ) cout << "got a electron, index: " << j << endl;

        float tmpMass = ( *(lepvec[i].vec) + *(lepvec[j].vec) ).M();
        if( ctrl.debug ) cout << "tmp vs best :: " << tmpMass << "," << bestMass << endl;
        if( tmpMass > bestMass ) { 
          bestMass = tmpMass;
          best_mu_index=i;
          best_ele_index=j;
        }

      }
    }

    if( bestMass > 0 ) { 
      if( ctrl.debug ) cout << "EMU candidate, mass : " << bestMass << endl;
      if( lepvec[best_ele_index].isTight ) 
        ret.status.setStatus(SelectionStatus::TIGHTSELECTION);
      else
        ret.status.setStatus(SelectionStatus::LOOSESELECTION);
      ret.Z1leptons.push_back(lepvec[best_mu_index]);
      ret.Z1leptons.push_back(lepvec[best_ele_index]);
      ret.Z2leptons.push_back(lepvec[best_mu_index]);
      ret.Z2leptons.push_back(lepvec[best_ele_index]);
      return ret;
    }

    ret.status.setStatus(0);
    return ret;

}


Revision:	1.2
Committed:	Tue Feb 14 21:29:18 2012 UTC (13 years, 3 months ago) by khahn
Content type:	text/plain
Branch:	MAIN
Changes since 1.1:	+10 -11 lines
Log Message:	various fixes
#	User	Rev	Content
1	khahn	1.1	#include "SelectionStatus.h"
2			#include "EventData.h"
3			#include "SimpleLepton.h"
4			#include "EfficiencyWeightsInterface.h"
5
6			#include "HZZBDTElectronSelection.h"
7			#include "IsolationSelection.h"
8			#include "PassHLT.h"
9			#include "SelectionEMU.h"
10
11			#include "ExternData.h"
12			#include "SelectionDefs.h"
13
14
15
16			// +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
17			// +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
18			EventData apply_EMU_selection(ControlFlags &ctrl, // input control
19			mithep::TEventInfo *info, // input event info
20			TClonesArray *electronArr, // input electrons
21			SelectionStatus (ElectronPreSelector)( ControlFlags &, const mithep::TElectron),
22			SelectionStatus (ElectronIDSelector)( ControlFlags &, const mithep::TElectron),
23			SelectionStatus (ElectronIsoSelector)( ControlFlags &, const mithep::TElectron),
24			TClonesArray *muonArr, // input muons
25			SelectionStatus (MuonPreSelector)( ControlFlags &, const mithep::TMuon),
26			SelectionStatus (MuonIDSelector)( ControlFlags &, const mithep::TMuon),
27			SelectionStatus (MuonIsoSelector)( ControlFlags &, const mithep::TMuon) )
28			// +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
29			// +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
30			{
31
32			EventData ret;
33			unsigned evtfail = 0x0;
34			TRandom3 r;
35
36			if( ctrl.debug ) {
37			cout << "Run: " << info->runNum
38			<< "\tEvt: " << info->evtNum
39			<< "\tLumi: " << info->lumiSec
40			<< endl;
41			}
42
43			if( !ctrl.mc ) {
44			// not accounting for overlap atm
45			RunLumiRangeMap::RunLumiPairType rl(info->runNum, info->lumiSec);
46			if( !(rlrm.HasRunLumi(rl)) ) {
47			if( ctrl.debug ) cout << "\tfails JSON" << endl;
48			ret.status.setStatus(0);
49			return ret;
50			}
51			}
52
53
54			//********************************************************
55			// Trigger
56			//********************************************************
57	khahn	1.2	// if( !ctrl.mc ) {
58			// if( !(passHLTSingleMuon(info->triggerBits) ) ) {
59	khahn	1.1	// if( ctrl.debug ) cout << "\tfails trigger" << endl;
60			// evtfail \|= (1<<EVTFAIL_TRIGGER);
61			// ret.status.setStatus(0);
62			// return ret;
63			// }
64	khahn	1.2	// }
65	khahn	1.1	if( ctrl.debug ) {
66			cout << "presel nlep: " << muonArr->GetEntries() + electronArr->GetEntries()
67			<< "\tnmuon: " << muonArr->GetEntries()
68			<< "\tnelectron: " << electronArr->GetEntries()
69			<< endl;
70			}
71
72			//********************************************************
73			// Lepton Selection
74			//********************************************************
75			vector<SimpleLepton> lepvec;
76			if( muonArr->GetEntries() != 1 && electronArr->GetEntries() != 1 ) {
77			ret.status.setStatus(0);
78			return ret;
79			}
80
81			//
82			if( ctrl.debug ) cout << "\tnMuons: " << muonArr->GetEntries() << endl;
83			//----------------------------------------------------
84			for(Int_t i=0; i<muonArr->GetEntries(); i++)
85			{
86			const mithep::TMuon mu = (mithep::TMuon)((*muonArr)[i]);
87
88			SelectionStatus musel;
89			if(ctrl.debug) cout << "musel.status before anything: " << musel.getStatus() << endl;
90			musel \|= (*MuonPreSelector)(ctrl,mu);
91			if(ctrl.debug) cout << "musel.status after presel: " << musel.getStatus() << endl;
92			musel \|= (*MuonIDSelector)(ctrl,mu);
93			if(ctrl.debug) cout << "musel.status after ID: " << musel.getStatus() << endl;
94			musel \|= (*MuonIsoSelector)(ctrl,mu);
95			if(ctrl.debug) cout << "musel.status after iso: " << musel.getStatus() << endl;
96
97			if( ctrl.debug ) {
98			cout << "muon:: pt: " << mu->pt
99			<< "\teta: " << mu->eta
100			<< "\tisorel: " << mu->pfIso03/mu->pt
101			<< "\tstatus: " << hex << musel.getStatus() << dec
102			<< endl;
103			}
104
105			if ( musel.pass() ) {
106
107			SimpleLepton tmplep;
108			float pt = mu->pt;
109			tmplep.vecorig->SetPtEtaPhiM(pt,
110			mu->eta,
111			mu->phi,
112			MUON_MASS);
113
114			if( ctrl.do_escale_up ) {
115			pt=scale_smear_muon_Up(pt, 1, r);
116			}
117			if( ctrl.do_escale_down ) {
118			pt=scale_smear_muon_Down(pt, 1, r);
119			}
120
121			tmplep.vec->SetPtEtaPhiM(pt,
122			mu->eta,
123			mu->phi,
124			MUON_MASS);
125
126			tmplep.type = 13;
127			tmplep.index = i;
128			tmplep.charge = mu->q;
129			tmplep.isoTrk = mu->trkIso03;
130			tmplep.isoEcal = mu->emIso03;
131			tmplep.isoHcal = mu->hadIso03;
132			tmplep.isoPF03 = mu->pfIso03;
133			tmplep.isoPF04 = mu->pfIso04;
134			tmplep.ip3dSig = mu->ip3dSig;
135			tmplep.is4l = false;
136			tmplep.isEB = (fabs(mu->eta) < 1.479 ? 1 : 0 );
137			tmplep.isTight = musel.tight();
138			tmplep.isLoose = musel.loose();
139			lepvec.push_back(tmplep);
140			if( ctrl.debug ) { cout << "muon passes ... " << endl;}
141			}
142			// }
143			}
144
145
146
147			//
148			if( ctrl.debug ) { cout << "\tnElectron: " << electronArr->GetEntries() << endl; }
149			// --------------------------------------------------------------------------------
150			for(Int_t i=0; i<electronArr->GetEntries(); i++)
151			{
152			const mithep::TElectron ele = (mithep::TElectron)((*electronArr)[i]);
153
154			Bool_t isMuonOverlap = kFALSE;
155			for (int k=0; k<lepvec.size(); ++k) {
156			TVector3 tmplep;
157			tmplep.SetPtEtaPhi(ele->pt, ele->eta, ele->phi);
158	khahn	1.2	if ( abs(lepvec[k].type == 13) && lepvec[k].vec->Vect().DrEtaPhi(tmplep) < 0.15 ) {
159	khahn	1.1	if( ctrl.debug ) cout << "-----> isMuonOverlap! " << endl;
160			isMuonOverlap = kTRUE;
161			break;
162			}
163			}
164
165			SelectionStatus elesel;
166			elesel \|= (*ElectronPreSelector)(ctrl,ele);
167			elesel \|= (*ElectronIDSelector)(ctrl,ele);
168			elesel \|= (*ElectronIsoSelector)(ctrl,ele);
169	khahn	1.2	if( elesel.getStatus() & SelectionStatus::PRESELECTION )
170	khahn	1.1	elesel.setStatus(SelectionStatus::LOOSESELECTION);
171
172			if( ctrl.debug ){
173			cout << "\tscEt: " << ele->scEt
174			<< "\tscEta: " << ele->scEta
175			<< "\tncluster: " << ele->ncluster
176			<< "\tisorel: " << ele->pfIso04/ele->pt
177			<< "\tstatus: " << hex << elesel.getStatus() << dec
178			<< endl;
179			}
180
181			if ( elesel.pass() && !isMuonOverlap )
182			{
183			SimpleLepton tmplep;
184
185			float pt = ele->pt;
186			tmplep.vecorig->SetPtEtaPhiM( pt,
187			ele->eta,
188			ele->phi,
189			ELECTRON_MASS );
190
191			if( ctrl.do_escale ) {
192			pt=scale_smear_electron(pt, ele->isEB, r);
193			}
194			if( ctrl.do_escale_up ) {
195			pt=scale_smear_electron_Up(pt, ele->isEB, r);
196			}
197			if( ctrl.do_escale_down ) {
198			pt=scale_smear_electron_Down(pt, ele->isEB, r);
199			}
200
201
202			tmplep.vec->SetPtEtaPhiM( pt,
203			ele->eta,
204			ele->phi,
205			ELECTRON_MASS );
206
207			tmplep.type = 11;
208			tmplep.index = i;
209			tmplep.charge = ele->q;
210			tmplep.isoTrk = ele->trkIso03;
211			tmplep.isoEcal = ele->emIso03;
212			tmplep.isoHcal = ele->hadIso03;
213			tmplep.isoPF03 = ele->pfIso03;
214			tmplep.isoPF04 = ele->pfIso04;
215			tmplep.ip3dSig = ele->ip3dSig;
216			tmplep.is4l = false;
217			tmplep.isEB = ele->isEB;
218			tmplep.scID = ele->scID;
219			tmplep.isTight = elesel.tight();
220			tmplep.isLoose = elesel.pass();
221			lepvec.push_back(tmplep);
222			if( ctrl.debug ) { cout << "\telectron passes ... " << endl; }
223			}
224			}
225
226
227			//********************************************************
228			// Dump Stuff
229			//********************************************************
230			sort( lepvec.begin(), lepvec.end(), SimpleLepton::lep_pt_sort );
231			int nmu=0, nele=0;
232			for( int i=0; i<lepvec.size(); i++ ) {
233			if(ctrl.debug) cout << "lepvec :: index: " << i
234			<< "\tpt: " << lepvec[i].vec->Pt()
235			<< "\ttype: " << lepvec[i].type
236			<< endl;
237			if( abs(lepvec[i].type) == 11 ) nele++;
238			else nmu++;
239			}
240			if( ctrl.debug ) {
241			cout << "postsel nlep: " << lepvec.size()
242			<< "\tnmuon: " << nmu
243			<< "\tnelectron: " << nele
244			<< endl;
245			}
246
247
248			//******************************************************************************
249			// Selection
250			//******************************************************************************
251			float bestMass=-1; int best_mu_index=-1, best_ele_index=-1;
252	khahn	1.2	for(int i = 0; i<lepvec.size(); i++) { // get a tight muon
253	khahn	1.1	if( abs(lepvec[i].type) != 13 ) continue;
254			if( !(lepvec[i].isTight) ) continue;
255	khahn	1.2	// if( lepvec[i].vec->Pt() < 35 ) continue;
256	khahn	1.1	if( ctrl.debug ) cout << "got a muon, index: " << i << endl;
257
258	khahn	1.2	for(int j = 0; j<lepvec.size(); j++) { // get a loose electron
259	khahn	1.1	if( j == i ) continue;
260			if( abs(lepvec[j].type) != 11 ) continue;
261			if( !(lepvec[j].isLoose) ) continue;
262	khahn	1.2	if (lepvec[i].charge == lepvec[j].charge) continue;
263			// if (lepvec[i].charge != lepvec[j].charge) continue;
264	khahn	1.1	if( ctrl.debug ) cout << "got a electron, index: " << j << endl;
265
266			float tmpMass = ( (lepvec[i].vec) + (lepvec[j].vec) ).M();
267			if( ctrl.debug ) cout << "tmp vs best :: " << tmpMass << "," << bestMass << endl;
268			if( tmpMass > bestMass ) {
269			bestMass = tmpMass;
270			best_mu_index=i;
271			best_ele_index=j;
272			}
273
274			}
275			}
276
277			if( bestMass > 0 ) {
278			if( ctrl.debug ) cout << "EMU candidate, mass : " << bestMass << endl;
279			if( lepvec[best_ele_index].isTight )
280			ret.status.setStatus(SelectionStatus::TIGHTSELECTION);
281			else
282			ret.status.setStatus(SelectionStatus::LOOSESELECTION);
283			ret.Z1leptons.push_back(lepvec[best_mu_index]);
284			ret.Z1leptons.push_back(lepvec[best_ele_index]);
285			ret.Z2leptons.push_back(lepvec[best_mu_index]);
286			ret.Z2leptons.push_back(lepvec[best_ele_index]);
287			return ret;
288			}
289
290			ret.status.setStatus(0);
291			return ret;
292
293			}
294
295