ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/MitHzz4l/NonMCBackground/src/SelectionEMU.cc
Revision: 1.4
Committed: Tue Jun 12 22:04:06 2012 UTC (12 years, 11 months ago) by dkralph
Content type: text/plain
Branch: MAIN
Changes since 1.3: +331 -260 lines
Log Message:
*** empty log message ***

File Contents

# User Rev Content
1 dkralph 1.4 //***************************************************************************************************
2     //
3     // selection sync'ed with https://twiki.cern.ch/twiki/bin/viewauth/CMS/HiggsZZ4l2012SummerSync
4     //
5     //***************************************************************************************************
6    
7     // system headers
8     #include <map>
9     #include <utility>
10    
11     // mit headers
12     #include "Vertex.h"
13    
14     // 4L stuff
15 khahn 1.1 #include "SelectionStatus.h"
16     #include "EventData.h"
17     #include "SimpleLepton.h"
18     #include "EfficiencyWeightsInterface.h"
19 dkralph 1.4 #include "ElectronSelection.h"
20     #include "MuonSelection.h"
21 khahn 1.1 #include "IsolationSelection.h"
22     #include "SelectionEMU.h"
23 dkralph 1.4 #include "ReferenceSelection.h"
24     #include "Selection.h"
25     #include "CommonDefs.h"
26 khahn 1.1 #include "SelectionDefs.h"
27 dkralph 1.4 #include "FSR.h"
28     #include "UtilFuncs.h"
29 khahn 1.1
30    
31 dkralph 1.4 extern vector<SimpleLepton> failingLeptons;
32     extern vector<SimpleLepton> passingLeptons;
33 khahn 1.1
34 dkralph 1.4 //--------------------------------------------------------------------------------------------------
35     EventData apply_HZZ4L_EMU_selection(ControlFlags &ctrl, // input control
36     const mithep::EventHeader *info, // input event info
37     const mithep::Array<mithep::Vertex> * vtxArr ,
38     const mithep::Array<mithep::PFCandidate> *pfCandidates,
39     const mithep::Array<mithep::PileupEnergyDensity> *puEnergyDensity,
40     const mithep::Array<mithep::Electron> *electronArr, // input electrons
41     SelectionStatus (*ElectronPreSelector)( ControlFlags &,
42     const mithep::Electron*,
43     const mithep::Vertex *),
44     SelectionStatus (*ElectronIDSelector)( ControlFlags &,
45     const mithep::Electron*,
46     const mithep::Vertex *),
47     SelectionStatus (*ElectronIsoSelector)( ControlFlags &,
48     const mithep::Electron*,
49     const mithep::Vertex *,
50     const mithep::Array<mithep::PFCandidate> *,
51     const mithep::Array<mithep::PileupEnergyDensity> *,
52     mithep::ElectronTools::EElectronEffectiveAreaTarget,
53     vector<const mithep::PFCandidate*>),
54     const mithep::Array<mithep::Muon> *muonArr, // input muons
55     SelectionStatus (*MuonPreSelector)( ControlFlags &,
56     const mithep::Muon*,
57     const mithep::Vertex *,
58     const mithep::Array<mithep::PFCandidate> *),
59     SelectionStatus (*MuonIDSelector)( ControlFlags &,
60     const mithep::Muon*,
61     // const mithep::Vertex &),
62     const mithep::Vertex *,
63     const mithep::Array<mithep::PFCandidate> *),
64     SelectionStatus (*MuonIsoSelector)( ControlFlags &,
65     const mithep::Muon*,
66     const mithep::Vertex *,
67     const mithep::Array<mithep::PFCandidate> *,
68     const mithep::Array<mithep::PileupEnergyDensity> *,
69     mithep::MuonTools::EMuonEffectiveAreaTarget,
70     vector<const mithep::PFCandidate*>)
71     )
72     //--------------------------------------------------------------------------------------------------
73 khahn 1.1 {
74    
75     EventData ret;
76     unsigned evtfail = 0x0;
77     TRandom3 r;
78    
79 dkralph 1.4 failingLeptons.clear();
80     passingLeptons.clear();
81    
82     mithep::MuonTools::EMuonEffectiveAreaTarget eraMu;
83     mithep::ElectronTools::EElectronEffectiveAreaTarget eraEle;
84     getEATargets(ctrl,eraMu,eraEle);
85    
86     const mithep::Vertex * vtx;
87     bool goodVertex = setPV( ctrl, vtxArr, vtx );
88     if(goodVertex) {
89     ret.status.selectionBits.flip(PASS_SKIM2);
90     } else {
91     if(ctrl.debug) cout << "found bad vertex" << endl;
92     ret.status.setStatus(SelectionStatus::FAIL);
93     return ret;
94 khahn 1.1 }
95    
96 dkralph 1.4 //***********************************************************
97     // Lepton Selection
98     //***********************************************************
99 khahn 1.1 vector<SimpleLepton> lepvec;
100 dkralph 1.4 vector<const mithep::PFCandidate*> photonsToVeto;
101    
102 khahn 1.1
103     if( ctrl.debug ) cout << "\tnMuons: " << muonArr->GetEntries() << endl;
104     //----------------------------------------------------
105 dkralph 1.4 for(int i=0; i<muonArr->GetEntries(); i++)
106 khahn 1.1 {
107 dkralph 1.4 const mithep::Muon *mu = (mithep::Muon*)((*muonArr)[i]);
108 khahn 1.1
109     SelectionStatus musel;
110 dkralph 1.4
111     musel |= (*MuonPreSelector)(ctrl,mu,vtx,pfCandidates);
112     if( !(musel.getStatus() & SelectionStatus::PRESELECTION) ) continue;
113    
114     musel |= (*MuonIDSelector)(ctrl,mu,vtx,pfCandidates );
115    
116     if( !(ctrl.doFSR) ) {
117     musel |= (*MuonIsoSelector)(ctrl,mu,vtx,pfCandidates,puEnergyDensity,eraMu,photonsToVeto);
118 khahn 1.1 }
119 dkralph 1.4
120     SimpleLepton tmplep;
121     float pt = mu->Pt();
122     tmplep.vec.SetPtEtaPhiM(pt,
123     mu->Eta(),
124     mu->Phi(),
125     MUON_MASS);
126 khahn 1.1
127 dkralph 1.4 tmplep.type = 13;
128     tmplep.index = i;
129     tmplep.charge = mu->Charge();
130     tmplep.isoTrk = mu->IsoR03SumPt();
131     tmplep.isoEcal = mu->IsoR03EmEt();
132     tmplep.isoHcal = mu->IsoR03HadEt();
133     tmplep.isoPF04 = musel.isoPF04;
134     tmplep.chisoPF04 = musel.chisoPF04;
135     tmplep.gaisoPF04 = musel.gaisoPF04;
136     tmplep.neisoPF04 = musel.neisoPF04;
137     // tmplep.isoPF03 = computePFMuonIso(mu,vtx,pfCandidates,0.3);
138     // tmplep.isoPF04 = computePFMuonIso(mu,vtx,pfCandidates,0.4);
139     tmplep.ip3dSig = mu->Ip3dPVSignificance();
140     tmplep.is4l = false;
141     tmplep.isEB = (fabs(mu->Eta()) < 1.479 ? 1 : 0 );
142     tmplep.isoMVA = musel.isoMVA;
143     tmplep.isTight = musel.tight();
144     tmplep.isLoose = musel.loose();
145     tmplep.status = musel;
146     tmplep.fsrRecoveryAttempted = false;
147     tmplep.tightCutsApplied = false;
148     lepvec.push_back(tmplep);
149     if( ctrl.debug ) cout << endl;
150     }
151    
152     if( ctrl.debug ) { cout << "\tnElectron: " << electronArr->GetEntries() << endl; }
153     // --------------------------------------------------------------------------------
154     for(int i=0; i<electronArr->GetEntries(); i++)
155     {
156     const mithep::Electron *ele = (mithep::Electron*)((*electronArr)[i]);
157 khahn 1.1
158 dkralph 1.4 SelectionStatus elesel;
159    
160     elesel |= (*ElectronPreSelector)(ctrl,ele,vtx);
161     if( !(elesel.getStatus() & SelectionStatus::PRESELECTION) ) continue;
162    
163     elesel |= (*ElectronIDSelector)(ctrl,ele,vtx);
164 khahn 1.1
165 dkralph 1.4 if( !(ctrl.doFSR) ) {
166     elesel |= (*ElectronIsoSelector)(ctrl,ele,vtx,pfCandidates,puEnergyDensity,eraEle,photonsToVeto);
167 khahn 1.1 }
168    
169 dkralph 1.4 SimpleLepton tmplep;
170     float pt = ele->Pt();
171     tmplep.vec.SetPtEtaPhiM( pt,
172     ele->Eta(),
173     ele->Phi(),
174     ELECTRON_MASS );
175 khahn 1.1
176 dkralph 1.4 tmplep.type = 11;
177 khahn 1.1 tmplep.index = i;
178 dkralph 1.4 tmplep.charge = ele->Charge();
179     tmplep.isoTrk = ele->TrackIsolationDr03();
180     tmplep.isoEcal = ele->EcalRecHitIsoDr03();
181     tmplep.isoHcal = ele->HcalTowerSumEtDr03();
182     tmplep.isoPF04 = elesel.isoPF04;
183     tmplep.chisoPF04 = elesel.chisoPF04;
184     tmplep.gaisoPF04 = elesel.gaisoPF04;
185     tmplep.neisoPF04 = elesel.neisoPF04;
186     // tmplep.isoPF03 = computePFEleIso(ele,vtx,pfCandidates,0.3);
187     // tmplep.isoPF04 = computePFEleIso(ele,vtx,pfCandidates,0.4);
188     tmplep.ip3dSig = ele->Ip3dPVSignificance();
189 khahn 1.1 tmplep.is4l = false;
190 dkralph 1.4 tmplep.isEB = ele->IsEB();
191     tmplep.scID = ele->SCluster()->GetUniqueID();
192     tmplep.isTight = elesel.tight();
193     tmplep.isLoose = elesel.loose();
194     tmplep.status = elesel;
195     tmplep.idMVA = elesel.idMVA;
196     tmplep.isoMVA = elesel.isoMVA;
197     tmplep.fsrRecoveryAttempted = false;
198     SelectionStatus tmpstat = electronTagSelection(ele,vtx,pfCandidates);
199     tmplep.tightCutsApplied = tmpstat.tight();
200 khahn 1.1 lepvec.push_back(tmplep);
201 dkralph 1.4 if( ctrl.debug ) cout << endl;
202 khahn 1.1 }
203 khahn 1.3
204 dkralph 1.4 sort( lepvec.begin(), lepvec.end(), SimpleLepton::lep_pt_sort );
205    
206     //********************************************************
207     // Step 2: Lepton Cleaning
208     //********************************************************
209     vector<vector<SimpleLepton>::iterator> electrons_to_erase;
210     for (vector<SimpleLepton>::iterator it1=lepvec.begin(); it1 != lepvec.end(); it1++ ) {
211     if ( abs(it1->type) != 11 ) continue;
212     TVector3 evec = it1->vec.Vect();
213    
214     bool erase_this_electron=false;
215     for (vector<SimpleLepton>::iterator it2=lepvec.begin(); it2 != lepvec.end(); it2++ ) {
216     if ( it2 == it1 ) continue;
217     if ( abs(it2->type) != 13 ) continue;
218     if( !(it2->status.looseIDAndPre()) ) continue;
219     TVector3 mvec = it2->vec.Vect();
220 khahn 1.1
221 dkralph 1.4 if ( evec.DrEtaPhi(mvec) < 0.05 ) {
222     erase_this_electron=true;
223     break;
224 khahn 1.1 }
225 dkralph 1.4 }
226     if( erase_this_electron )
227     electrons_to_erase.push_back(it1);
228     }
229     for( int i=0; i<electrons_to_erase.size(); i++ ) {
230     lepvec.erase(electrons_to_erase[i]);
231     }
232 khahn 1.1
233     //********************************************************
234 dkralph 1.4 // Step 3: Good Leptons
235 khahn 1.1 //********************************************************
236 dkralph 1.4 vector<double> pt_of_leptons_to_erase;
237     for (int i=0; i<lepvec.size(); i++ ) {
238     bool already_pushed=false;
239     if( !(lepvec[i].status.loose()) ) {
240     pt_of_leptons_to_erase.push_back(lepvec[i].vec.Pt());
241     already_pushed = true;
242     failingLeptons.push_back(lepvec[i]); // these should pass preselection
243     } else {
244     passingLeptons.push_back(lepvec[i]);
245     }
246     #ifndef SYNC
247     if( !already_pushed && fabs(lepvec[i].ip3dSig)>4 )
248     pt_of_leptons_to_erase.push_back(lepvec[i].vec.Pt());
249     #endif
250 khahn 1.1 }
251 dkralph 1.4 for( int i=0; i<pt_of_leptons_to_erase.size(); i++ ) {
252     for( vector<SimpleLepton>::iterator it=lepvec.begin(); it != lepvec.end(); it++ ) {
253     SimpleLepton flep = *it;
254     if( flep.vec.Pt() != pt_of_leptons_to_erase[i] ) continue;
255     lepvec.erase(it);
256     break;
257     }
258 khahn 1.1 }
259 dkralph 1.4
260 khahn 1.1 //******************************************************************************
261 dkralph 1.4 // W + (OF SS lepton) Selection
262 khahn 1.1 //******************************************************************************
263 dkralph 1.4 if(ctrl.fakeScheme.Contains("emu-")) {
264     if(has_ssof_lepton(ctrl)) {
265     ret.status.setStatus(SelectionStatus::EVTPASS);
266     ret.Z1leptons.push_back(passingLeptons[0]);
267     ret.Z1leptons.push_back(passingLeptons[0]);
268     ret.Z2leptons.push_back(passingLeptons[0]);
269     ret.Z2leptons.push_back(passingLeptons[0]);
270     } else {
271     ret.status.setStatus(SelectionStatus::FAIL);
272 khahn 1.1 }
273 dkralph 1.4 return ret;
274 khahn 1.1 }
275 dkralph 1.4 }
276 khahn 1.1
277 dkralph 1.4 // //----------------------------------------------------------------------------
278     // //
279     // // Get primary vertices
280     // // Assumes primary vertices are ordered by sum-pT^2 (as should be in CMSSW)
281     // // NOTE: if no PV is found from fitting tracks, the beamspot is used
282     // //
283     // //----------------------------------------------------------------------------
284     // bool setPV(ControlFlags ctrl,
285     // const mithep::Array<mithep::Vertex> * vtxArr,
286     // const mithep::Vertex* &vtx)
287     // //----------------------------------------------------------------------------
288     // {
289    
290     // const mithep::Vertex *bestPV = 0;
291     // float best_sumpt=-1;
292    
293     // // good PV requirements
294     // const UInt_t fMinNTracksFit = 0;
295     // const Double_t fMinNdof = 4;
296     // const Double_t fMaxAbsZ = 24;
297     // const Double_t fMaxRho = 2;
298    
299     // for(int i=0; i<vtxArr->GetEntries(); ++i) {
300     // const mithep::Vertex *pv = (mithep::Vertex*)(vtxArr->At(i));
301     // if( ctrl.debug ) cout << "vertex :: " << i << "\tntrks: " << pv->NTracks() << endl;
302    
303     // // Select best PV for corrected d0; if no PV passing cuts, the first PV in the collection will be used
304     // if(!pv->IsValid()) continue;
305     // if(pv->NTracksFit() < fMinNTracksFit) continue;
306     // if(pv->Ndof() < fMinNdof) continue;
307     // if(fabs(pv->Z()) > fMaxAbsZ) continue;
308     // if(pv->Position().Rho() > fMaxRho) continue;
309    
310     // // take the first ...
311     // bestPV = pv;
312     // break;
313    
314     // // this never reached ...
315     // float tmp_sumpt=0;
316     // for( int t=0; t<pv->NTracks(); t++ )
317     // tmp_sumpt += pv->Trk(t)->Pt();
318    
319     // if( tmp_sumpt > best_sumpt ) {
320     // bestPV = pv;
321     // best_sumpt = tmp_sumpt;
322     // if( ctrl.debug) cout << "new PV set, pt : " << best_sumpt << endl;
323     // }
324     // }
325    
326     // // sync
327     // if(!bestPV)
328     // return false;
329     // else {
330     // vtx = bestPV;
331     // return true;
332     // }
333     // };
334     //----------------------------------------------------------------------------------------
335     bool has_ssof_lepton(ControlFlags &ctrl)
336     {
337     bool has_ssof=false;
338     for(unsigned iw=0; iw<passingLeptons.size(); iw++) {
339     SimpleLepton w_lep = passingLeptons[iw];
340     //????????????????????????????????????????????????????????????????????????????????????????
341     // this is applied in skim (skim also applies ww muon id)
342     // if(abs(w_lep.type) == 11) {
343     // if( !(w_lep.tightCutsApplied) )
344     // continue;
345     // }
346     //????????????????????????????????????????????????????????????????????????????????????????
347     for(unsigned ifake=0; ifake<failingLeptons.size(); ifake++) {
348     SimpleLepton fake_lep = failingLeptons[ifake];
349     if(abs(fake_lep.type) == abs(w_lep.type)) continue;
350     if(fake_lep.charge != w_lep.charge) continue;
351     has_ssof = true;
352     }
353     for(unsigned ipass=0; ipass<passingLeptons.size(); ipass++) {
354     if(ipass == iw) continue;
355     SimpleLepton pass_lep = passingLeptons[ipass];
356     if(abs(pass_lep.type) == abs(w_lep.type)) continue;
357     if(pass_lep.charge != w_lep.charge) continue;
358     has_ssof = true;
359 khahn 1.1 }
360 dkralph 1.4 }
361 khahn 1.1
362 dkralph 1.4 return has_ssof;
363 khahn 1.1 }
364 dkralph 1.4 // //----------------------------------------------------------------------------------------
365     // void getEATargets(ControlFlags &ctrl, mithep::MuonTools::EMuonEffectiveAreaTarget &eraMu, mithep::ElectronTools::EElectronEffectiveAreaTarget &eraEle)
366     // {
367     // if( !ctrl.mc && ctrl.era == 2011 ) {
368     // eraMu = mithep::MuonTools::kMuEAData2011;
369     // eraEle = mithep::ElectronTools::kEleEAData2011;
370     // } else if( !ctrl.mc && ctrl.era == 2012 ) {
371     // eraMu = mithep::MuonTools::kMuEAData2012;
372     // eraEle = mithep::ElectronTools::kEleEAData2012;
373     // } else if( ctrl.mc && ctrl.era == 2011 ) {
374     // eraMu = mithep::MuonTools::kMuEAFall11MC;
375     // eraEle = mithep::ElectronTools::kEleEAFall11MC;
376     // } else if( ctrl.mc && ctrl.era == 2012 ) {
377     // eraMu = mithep::MuonTools::kMuEAData2012;
378     // eraEle = mithep::ElectronTools::kEleEAData2012;
379     // } else {
380     // cerr << "unknown era for effective areas ... quitting." << endl;
381     // exit(1);
382     // }
383     // }