ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/MitHzz4l/NonMCBackground/src/SelectionEMU.cc
Revision: 1.6
Committed: Wed Jun 13 14:11:22 2012 UTC (12 years, 11 months ago) by dkralph
Content type: text/plain
Branch: MAIN
Changes since 1.5: +2 -1 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 dkralph 1.5 #include "SelectionFuncs.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 dkralph 1.6 SelectionStatus tmpstat = PassWwMuonSel(mu,vtx,pfCandidates);
148     tmplep.tightCutsApplied = tmpstat.tight();
149 dkralph 1.4 lepvec.push_back(tmplep);
150     if( ctrl.debug ) cout << endl;
151     }
152    
153     if( ctrl.debug ) { cout << "\tnElectron: " << electronArr->GetEntries() << endl; }
154     // --------------------------------------------------------------------------------
155     for(int i=0; i<electronArr->GetEntries(); i++)
156     {
157     const mithep::Electron *ele = (mithep::Electron*)((*electronArr)[i]);
158 khahn 1.1
159 dkralph 1.4 SelectionStatus elesel;
160    
161     elesel |= (*ElectronPreSelector)(ctrl,ele,vtx);
162     if( !(elesel.getStatus() & SelectionStatus::PRESELECTION) ) continue;
163    
164     elesel |= (*ElectronIDSelector)(ctrl,ele,vtx);
165 khahn 1.1
166 dkralph 1.4 if( !(ctrl.doFSR) ) {
167     elesel |= (*ElectronIsoSelector)(ctrl,ele,vtx,pfCandidates,puEnergyDensity,eraEle,photonsToVeto);
168 khahn 1.1 }
169    
170 dkralph 1.4 SimpleLepton tmplep;
171     float pt = ele->Pt();
172     tmplep.vec.SetPtEtaPhiM( pt,
173     ele->Eta(),
174     ele->Phi(),
175     ELECTRON_MASS );
176 khahn 1.1
177 dkralph 1.4 tmplep.type = 11;
178 khahn 1.1 tmplep.index = i;
179 dkralph 1.4 tmplep.charge = ele->Charge();
180     tmplep.isoTrk = ele->TrackIsolationDr03();
181     tmplep.isoEcal = ele->EcalRecHitIsoDr03();
182     tmplep.isoHcal = ele->HcalTowerSumEtDr03();
183     tmplep.isoPF04 = elesel.isoPF04;
184     tmplep.chisoPF04 = elesel.chisoPF04;
185     tmplep.gaisoPF04 = elesel.gaisoPF04;
186     tmplep.neisoPF04 = elesel.neisoPF04;
187     // tmplep.isoPF03 = computePFEleIso(ele,vtx,pfCandidates,0.3);
188     // tmplep.isoPF04 = computePFEleIso(ele,vtx,pfCandidates,0.4);
189     tmplep.ip3dSig = ele->Ip3dPVSignificance();
190 khahn 1.1 tmplep.is4l = false;
191 dkralph 1.4 tmplep.isEB = ele->IsEB();
192     tmplep.scID = ele->SCluster()->GetUniqueID();
193     tmplep.isTight = elesel.tight();
194     tmplep.isLoose = elesel.loose();
195     tmplep.status = elesel;
196     tmplep.idMVA = elesel.idMVA;
197     tmplep.isoMVA = elesel.isoMVA;
198     tmplep.fsrRecoveryAttempted = false;
199     SelectionStatus tmpstat = electronTagSelection(ele,vtx,pfCandidates);
200     tmplep.tightCutsApplied = tmpstat.tight();
201 khahn 1.1 lepvec.push_back(tmplep);
202 dkralph 1.4 if( ctrl.debug ) cout << endl;
203 khahn 1.1 }
204 khahn 1.3
205 dkralph 1.4 sort( lepvec.begin(), lepvec.end(), SimpleLepton::lep_pt_sort );
206    
207     //********************************************************
208     // Step 2: Lepton Cleaning
209     //********************************************************
210     vector<vector<SimpleLepton>::iterator> electrons_to_erase;
211     for (vector<SimpleLepton>::iterator it1=lepvec.begin(); it1 != lepvec.end(); it1++ ) {
212     if ( abs(it1->type) != 11 ) continue;
213     TVector3 evec = it1->vec.Vect();
214    
215     bool erase_this_electron=false;
216     for (vector<SimpleLepton>::iterator it2=lepvec.begin(); it2 != lepvec.end(); it2++ ) {
217     if ( it2 == it1 ) continue;
218     if ( abs(it2->type) != 13 ) continue;
219     if( !(it2->status.looseIDAndPre()) ) continue;
220     TVector3 mvec = it2->vec.Vect();
221 khahn 1.1
222 dkralph 1.4 if ( evec.DrEtaPhi(mvec) < 0.05 ) {
223     erase_this_electron=true;
224     break;
225 khahn 1.1 }
226 dkralph 1.4 }
227     if( erase_this_electron )
228     electrons_to_erase.push_back(it1);
229     }
230     for( int i=0; i<electrons_to_erase.size(); i++ ) {
231     lepvec.erase(electrons_to_erase[i]);
232     }
233 khahn 1.1
234     //********************************************************
235 dkralph 1.4 // Step 3: Good Leptons
236 khahn 1.1 //********************************************************
237 dkralph 1.4 vector<double> pt_of_leptons_to_erase;
238     for (int i=0; i<lepvec.size(); i++ ) {
239     bool already_pushed=false;
240     if( !(lepvec[i].status.loose()) ) {
241     pt_of_leptons_to_erase.push_back(lepvec[i].vec.Pt());
242     already_pushed = true;
243     failingLeptons.push_back(lepvec[i]); // these should pass preselection
244     } else {
245     passingLeptons.push_back(lepvec[i]);
246     }
247     #ifndef SYNC
248     if( !already_pushed && fabs(lepvec[i].ip3dSig)>4 )
249     pt_of_leptons_to_erase.push_back(lepvec[i].vec.Pt());
250     #endif
251 khahn 1.1 }
252 dkralph 1.4 for( int i=0; i<pt_of_leptons_to_erase.size(); i++ ) {
253     for( vector<SimpleLepton>::iterator it=lepvec.begin(); it != lepvec.end(); it++ ) {
254     SimpleLepton flep = *it;
255     if( flep.vec.Pt() != pt_of_leptons_to_erase[i] ) continue;
256     lepvec.erase(it);
257     break;
258     }
259 khahn 1.1 }
260 dkralph 1.4
261 khahn 1.1 //******************************************************************************
262 dkralph 1.4 // W + (OF SS lepton) Selection
263 khahn 1.1 //******************************************************************************
264 dkralph 1.4 if(ctrl.fakeScheme.Contains("emu-")) {
265     if(has_ssof_lepton(ctrl)) {
266     ret.status.setStatus(SelectionStatus::EVTPASS);
267     ret.Z1leptons.push_back(passingLeptons[0]);
268     ret.Z1leptons.push_back(passingLeptons[0]);
269     ret.Z2leptons.push_back(passingLeptons[0]);
270     ret.Z2leptons.push_back(passingLeptons[0]);
271     } else {
272     ret.status.setStatus(SelectionStatus::FAIL);
273 khahn 1.1 }
274 dkralph 1.4 return ret;
275 khahn 1.1 }
276 dkralph 1.4 }
277 khahn 1.1
278 dkralph 1.4 // //----------------------------------------------------------------------------
279     // //
280     // // Get primary vertices
281     // // Assumes primary vertices are ordered by sum-pT^2 (as should be in CMSSW)
282     // // NOTE: if no PV is found from fitting tracks, the beamspot is used
283     // //
284     // //----------------------------------------------------------------------------
285     // bool setPV(ControlFlags ctrl,
286     // const mithep::Array<mithep::Vertex> * vtxArr,
287     // const mithep::Vertex* &vtx)
288     // //----------------------------------------------------------------------------
289     // {
290    
291     // const mithep::Vertex *bestPV = 0;
292     // float best_sumpt=-1;
293    
294     // // good PV requirements
295     // const UInt_t fMinNTracksFit = 0;
296     // const Double_t fMinNdof = 4;
297     // const Double_t fMaxAbsZ = 24;
298     // const Double_t fMaxRho = 2;
299    
300     // for(int i=0; i<vtxArr->GetEntries(); ++i) {
301     // const mithep::Vertex *pv = (mithep::Vertex*)(vtxArr->At(i));
302     // if( ctrl.debug ) cout << "vertex :: " << i << "\tntrks: " << pv->NTracks() << endl;
303    
304     // // Select best PV for corrected d0; if no PV passing cuts, the first PV in the collection will be used
305     // if(!pv->IsValid()) continue;
306     // if(pv->NTracksFit() < fMinNTracksFit) continue;
307     // if(pv->Ndof() < fMinNdof) continue;
308     // if(fabs(pv->Z()) > fMaxAbsZ) continue;
309     // if(pv->Position().Rho() > fMaxRho) continue;
310    
311     // // take the first ...
312     // bestPV = pv;
313     // break;
314    
315     // // this never reached ...
316     // float tmp_sumpt=0;
317     // for( int t=0; t<pv->NTracks(); t++ )
318     // tmp_sumpt += pv->Trk(t)->Pt();
319    
320     // if( tmp_sumpt > best_sumpt ) {
321     // bestPV = pv;
322     // best_sumpt = tmp_sumpt;
323     // if( ctrl.debug) cout << "new PV set, pt : " << best_sumpt << endl;
324     // }
325     // }
326    
327     // // sync
328     // if(!bestPV)
329     // return false;
330     // else {
331     // vtx = bestPV;
332     // return true;
333     // }
334     // };
335     //----------------------------------------------------------------------------------------
336     bool has_ssof_lepton(ControlFlags &ctrl)
337     {
338     bool has_ssof=false;
339     for(unsigned iw=0; iw<passingLeptons.size(); iw++) {
340     SimpleLepton w_lep = passingLeptons[iw];
341     //????????????????????????????????????????????????????????????????????????????????????????
342     // this is applied in skim (skim also applies ww muon id)
343     // if(abs(w_lep.type) == 11) {
344     // if( !(w_lep.tightCutsApplied) )
345     // continue;
346     // }
347     //????????????????????????????????????????????????????????????????????????????????????????
348     for(unsigned ifake=0; ifake<failingLeptons.size(); ifake++) {
349     SimpleLepton fake_lep = failingLeptons[ifake];
350     if(abs(fake_lep.type) == abs(w_lep.type)) continue;
351     if(fake_lep.charge != w_lep.charge) continue;
352     has_ssof = true;
353     }
354     for(unsigned ipass=0; ipass<passingLeptons.size(); ipass++) {
355     if(ipass == iw) continue;
356     SimpleLepton pass_lep = passingLeptons[ipass];
357     if(abs(pass_lep.type) == abs(w_lep.type)) continue;
358     if(pass_lep.charge != w_lep.charge) continue;
359     has_ssof = true;
360 khahn 1.1 }
361 dkralph 1.4 }
362 khahn 1.1
363 dkralph 1.4 return has_ssof;
364 khahn 1.1 }
365 dkralph 1.4 // //----------------------------------------------------------------------------------------
366     // void getEATargets(ControlFlags &ctrl, mithep::MuonTools::EMuonEffectiveAreaTarget &eraMu, mithep::ElectronTools::EElectronEffectiveAreaTarget &eraEle)
367     // {
368     // if( !ctrl.mc && ctrl.era == 2011 ) {
369     // eraMu = mithep::MuonTools::kMuEAData2011;
370     // eraEle = mithep::ElectronTools::kEleEAData2011;
371     // } else if( !ctrl.mc && ctrl.era == 2012 ) {
372     // eraMu = mithep::MuonTools::kMuEAData2012;
373     // eraEle = mithep::ElectronTools::kEleEAData2012;
374     // } else if( ctrl.mc && ctrl.era == 2011 ) {
375     // eraMu = mithep::MuonTools::kMuEAFall11MC;
376     // eraEle = mithep::ElectronTools::kEleEAFall11MC;
377     // } else if( ctrl.mc && ctrl.era == 2012 ) {
378     // eraMu = mithep::MuonTools::kMuEAData2012;
379     // eraEle = mithep::ElectronTools::kEleEAData2012;
380     // } else {
381     // cerr << "unknown era for effective areas ... quitting." << endl;
382     // exit(1);
383     // }
384     // }