ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/MitHzz4l/LeptonSelection/src/MuonSelection.cc
(Generate patch)

Comparing UserCode/MitHzz4l/LeptonSelection/src/MuonSelection.cc (file contents):
Revision 1.3 by anlevin, Wed Feb 29 10:10:18 2012 UTC vs.
Revision 1.7 by khahn, Mon Apr 30 21:42:16 2012 UTC

# Line 1 | Line 1
1   #include <math.h>
2   #include <iostream>
3  
4 #include "HiggsAnaDefs.hh"
5
4   #include "MuonSelection.h"
5   #include "ParseArgs.h"
6   #include "SelectionStatus.h"
7  
8 < SelectionStatus noPreSelection( ControlFlags &ctrl,  const mithep::TMuon * mu ) {
8 > #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 > SelectionStatus noPreSelection( ControlFlags &ctrl,  const mithep::Muon * mu )
23 >
24 > //--------------------------------------------------------------------------------------------------
25 > {
26          SelectionStatus status;
27          status.setStatus(SelectionStatus::PRESELECTION);
28          if(ctrl.debug) cout << "muon presel returning status : " << status.getStatus() << endl;
29          return status;
30   }
31  
32 < SelectionStatus muonPreSelection( ControlFlags &ctrl,  const mithep::TMuon * mu ) {
32 > //--------------------------------------------------------------------------------------------------
33 > SelectionStatus muonPreSelection( ControlFlags &ctrl,  
34 >                                  const mithep::Muon * mu,
35 >                                  const mithep::Vertex & vtx,
36 >                                  const mithep::Array<mithep::PFCandidate> * pfCandidates )
37 > //--------------------------------------------------------------------------------------------------
38 > {
39 >  bool ret = true;
40    if(ctrl.debug) cout << "inside muonpresel ... " << endl;
41 <  bool ret = isMuFO(mu);
41 >  //  bool ret = isMuFO(mu,vtx,pfCandidates);
42    if(ctrl.debug) cout << "iS fo? ... " << ret << endl;  
43 <  ret &= ( fabs(mu->ip3dSig) < 4 );
44 <  if(ctrl.debug) cout << "and pass IP (" << mu->ip3dSig << ") ? ... " << ret << endl;  
45 <  ret &= ( mu->pt > 5 );
46 <  if(ctrl.debug) cout << "and >5 GeV ? ... " << ret << endl;  
47 <  ret &= ( fabs(mu->eta) < 2.4 );
48 <  if(ctrl.debug) cout << "and < 2.4 eta ? ... " << ret << endl;  
43 >  ret &= ( fabs(mu->Ip3dPVSignificance()) < 4 );
44 >  if(ctrl.debug) cout << "and pass IP (" << mu->Ip3dPVSignificance() << ") ? ... " << ret << endl;  
45 >  ret &= ( mu->Pt() > 5 );
46 >  if(ctrl.debug) cout << "and >5 GeV (" << mu->Pt() << ") ? ... " << ret << endl;  
47 >  ret &= ( fabs(mu->Eta()) <= 2.4 );
48 >  if(ctrl.debug) cout << "and < 2.4 eta (" << mu->Eta() << ")? ... " << ret << endl;  
49 >  ret &=  (mu->IsTrackerMuon() && mu->HasTrackerTrk() &&
50 >           (mu->Quality().QualityMask().Mask() & mithep::MuonQuality::AllArbitrated));
51 >  if(ctrl.debug) cout << "and isTrackerMuon ? ... " << ret << endl;  
52 >  //  ret &=  (mu->IsoR03SumPt()/mu->Pt() < 0.7 );
53 >  ret &=  (mu->IsoR03SumPt()/mu->Pt() < 0.7 );
54 >  if(ctrl.debug) cout << "and loose trk iso ? ... " << ret << endl;  
55  
56    SelectionStatus status;
57    if( ret ) status.setStatus(SelectionStatus::PRESELECTION);
# Line 31 | Line 59 | SelectionStatus muonPreSelection( Contro
59    return status;
60   }
61  
62 < bool isMuFO( const mithep::TMuon * mu ) {
62 >
63 > //--------------------------------------------------------------------------------------------------
64 > bool isMuFO( const mithep::Muon * mu,
65 >             const mithep::Vertex & vtx,
66 >             const mithep::Array<mithep::PFCandidate> * pfCandidates )
67 > //--------------------------------------------------------------------------------------------------
68 > {
69    bool isgood=true;
70 <  float reliso = mu->pfIso03/mu->pt;
70 >  float reliso = computePFMuonIso(mu,vtx,pfCandidates,0.3)/mu->Pt();
71  
72 <  if( mu->pt < 5  )         isgood=false;                          
73 <  if ( fabs(mu->dz) > 0.2 ) isgood=false;   // needed to remove cosmics in HF sample
72 >  if( mu->Pt() < 5  )         isgood=false;                          
73 >  if ( fabs(mu->BestTrk()->DzCorrected(vtx)) > 0.2 ) isgood=false;   // needed to remove cosmics in HF sample
74  
75    //
76  
77    // HF seems to not want tkquality, SS does
78 <  if( mu->pt>20 ) {
79 <    if ( !((mu->typeBits & kGlobal) && mu->nSeg>0  ) ) isgood=false; //&& mu->nValidHits>0
78 >  if( mu->Pt()>20 ) {
79 >    if ( !((mu->IsGlobalMuon()) && mu->NSegments()>0  ) ) isgood=false; //&& mu->nValidHits>0
80    } else {
81 <    if ( !((mu->typeBits & kGlobal) && mu->nSeg>0 ) //&& //&& mu->nValidHits>0
82 <         && !( (mu->typeBits & kTracker) &&  (mu->qualityBits & kTMOneStationLoose)
81 >    if ( !((mu->IsGlobalMuon()) && mu->NSegments()>0 ) //&& //&& mu->nValidHits>0
82 >         && !( (mu->IsTrackerMuon()) &&  (mu->Quality().QualityMask().Mask() & muTools.kTMOneStationLoose)
83                 )
84           )   isgood=false;
85    }
# Line 64 | Line 98 | bool isMuFO( const mithep::TMuon * mu )
98  
99  
100    // comment for more stats for mu fake shape
101 <  if( mu->pt < 20 ) {
101 >  if( mu->Pt() < 20 ) {
102      if( reliso > 3 ) isgood=false;    
103    } else {
104      if( reliso > 5 ) isgood=false;    
# Line 73 | Line 107 | bool isMuFO( const mithep::TMuon * mu )
107    return isgood;
108   };
109  
110 <
111 < SelectionStatus passSoftMuonSelection( ControlFlags &ctrl, const mithep::TMuon * mu ) {
110 > //--------------------------------------------------------------------------------------------------
111 > SelectionStatus passSoftMuonSelection( ControlFlags &ctrl,
112 >                                       const mithep::Muon * mu ,
113 >                                       const mithep::Vertex & vtx)
114 > //--------------------------------------------------------------------------------------------------
115 > {
116  
117    int level=0;
118    unsigned failmask=0x0;
119    
120 <  if(mu->nTkHits < 11 )  
120 >  // Use tracker track when available
121 >  const mithep::Track *muTrk=0;
122 >  if(mu->HasStandaloneTrk()) { muTrk = mu->StandaloneTrk(); }
123 >  else if(mu->HasTrackerTrk())         { muTrk = mu->TrackerTrk(); }
124 >  else if(mu->HasGlobalTrk())     { muTrk = mu->GlobalTrk(); }
125 >  assert(muTrk);                  
126 >
127 >  if(mu->TrackerTrk()->NHits() < 11 )  
128      failmask |= (1<<level);
129    level++;
130  
131 <  if(fabs(mu->d0) > 0.2)
131 >  if(fabs(muTrk->D0Corrected(vtx)) > 0.2)
132      failmask |= (1<<level);
133    level++;
134  
135 <  if(fabs(mu->dz) > 0.1)
135 >  if(fabs(muTrk->DzCorrected(vtx)) > 0.1)
136      failmask |= (1<<level);
137    level++;
138  
139 <  if(!(mu->typeBits & kTracker))
139 >  if(!(mu->IsTrackerMuon()))
140      failmask |= (1<<level);
141    level++;
142  
143 <  if(!(mu->qualityBits & kTMLastStationAngTight))
143 >  //  if(!(mu->Quality().QualityMask().Mask() & muTools.kTMLastStationAngTight))
144 >  if(!(mu->Quality().QualityMask().Mask() & muTools.kTMLastStationTight))
145      failmask |= (1<<level);      
146    level++;
147  
148 <  Double_t iso = (mu->trkIso03 + mu->emIso03 + mu->hadIso03)/mu->pt;
149 <  if(mu->pt>20 && iso<0.1)
148 >  Double_t iso = (mu->IsoR03SumPt() + mu->IsoR03EmEt() + mu->IsoR03HadEt())/mu->Pt();
149 >  if(mu->Pt()>20 && iso<0.1)
150      failmask |= (1<<level);
151  
152    bool passtight  = !(failmask);
# Line 109 | Line 155 | SelectionStatus passSoftMuonSelection( C
155    
156   }
157  
158 <
158 > //--------------------------------------------------------------------------------------------------
159   //unsigned passMuonSelectionZZ( const mithep::TMuon * mu ) {
160 < SelectionStatus passMuonSelectionZZ( ControlFlags &ctrl, const mithep::TMuon * mu ) {
160 > SelectionStatus passMuonSelectionZZ( ControlFlags &ctrl, const mithep::Muon * mu )
161 > //--------------------------------------------------------------------------------------------------
162 > {
163    int level=0;
164    unsigned failmask=0x0;
165  
166 <  if(mu->pt < 5) {
166 >  if(mu->Pt() < 5) {
167      failmask |= (1<<level);
168    }
169  
170    level++;
171 <  if(fabs(mu->eta) > 2.4) {
171 >  if(fabs(mu->Eta()) > 2.4) {
172      failmask |= (1<<level);
173    }
174  
175    level++;
176 <  if(!(mu->typeBits & kGlobal)) {
176 >  if(!(mu->IsGlobalMuon())) {
177      failmask |= (1<<level);
178    }
179  
180    level++;
181 <  if(mu->nTkHits          < 10) {
181 >  if(mu->TrackerTrk()->NHits()  < 10) {
182      failmask |= (1<<level);
183    }
184  
185  
186    level++;  
187 <  if( (mu->trkIso03/mu->pt) > 0.7 ) {
187 >  if( (mu->IsoR03SumPt()/mu->Pt()) > 0.7 ) {
188      failmask |= (1<<level);
189    }
190  
# Line 148 | Line 196 | SelectionStatus passMuonSelectionZZ( Con
196    //  if(fabs(muon->dz)       > 0.1)   return false;  
197    
198    SelectionStatus status;
199 <  if( !failmask ) status.setStatus(SelectionStatus::LOOSEID);
200 <  if( !failmask ) status.setStatus(SelectionStatus::TIGHTID);
199 >  if( !failmask ) status.orStatus(SelectionStatus::LOOSEID);
200 >  if( !failmask ) status.orStatus(SelectionStatus::TIGHTID);
201    return status;
202  
203  
# Line 162 | Line 210 | SelectionStatus passMuonSelectionZZ( Con
210   //
211   // Kevin's WW selection
212   //
213 < SelectionStatus passMuonSelection( ControlFlags &ctrl, const mithep::TMuon * mu ) {
213 > //--------------------------------------------------------------------------------------------------
214 > SelectionStatus passMuonSelection( ControlFlags &ctrl,
215 >                                   const mithep::Muon * mu,
216 >                                   const mithep::Vertex & vtx   )
217 > //--------------------------------------------------------------------------------------------------
218 > {
219 >  SelectionStatus status;
220    int level=0;
221    unsigned failmask=0x0;
222  
223 +  // Use tracker track when available
224 +  const mithep::Track *muTrk=0;
225 +  if(mu->HasTrackerTrk())         { muTrk = mu->TrackerTrk(); }
226 +  else if(mu->HasGlobalTrk())     { muTrk = mu->GlobalTrk(); }
227 +  else if(mu->HasStandaloneTrk()) { muTrk = mu->StandaloneTrk(); }
228 +  assert(muTrk);                  
229 +  
230 +
231  
232    // 0x1
233 <  if(mu->pt < 5) {
233 >  if(muTrk->Pt() < 5) {
234      failmask |= (1<<level);
235    }
236  
237    // 0x2
238    level++;
239 <  if(fabs(mu->eta) > 2.4) {
239 >  if(fabs(muTrk->Eta()) > 2.4) {
240      failmask |= (1<<level);
241    }
242  
243    // 0x4
244    level++;
245 <  if(mu->ptErr/mu->pt > 0.1)   {
245 >  if(muTrk->PtErr()/muTrk->Pt() > 0.1)   {
246      failmask |= (1<<level);
247    }
248  
249    // 0x8
250    level++;
251 <  if( fabs(mu->dz) > 0.1 )   {
251 >  if( fabs(muTrk->DzCorrected(vtx)) > 0.1 )   {
252      failmask |= (1<<level);
253    }
254  
255 +  double muNchi2;
256 +  if(mu->HasStandaloneTrk())      { muNchi2 = mu->StandaloneTrk()->RChi2(); }
257 +  else if(mu->HasGlobalTrk())     { muNchi2 = mu->GlobalTrk()->RChi2();     }
258 +  else if(mu->HasTrackerTrk())    { muNchi2 = mu->TrackerTrk()->RChi2();    }
259  
260 <  Bool_t isGlobal  = (mu->typeBits & kGlobal) && (mu->muNchi2 < 10) && (mu->nMatch > 1) && (mu->nValidHits > 0);
261 <  Bool_t isTracker = (mu->typeBits & kTracker) && (mu->qualityBits & kTMLastStationTight);
260 >  unsigned qualityBits = mu->Quality().QualityMask().Mask();
261 >  
262 >  Bool_t isGlobal  = (mu->IsGlobalMuon()) && (muNchi2 < 10) && (mu->NMatches() > 1) && (mu->NValidHits() > 0);
263 >  Bool_t isTracker = (mu->IsTrackerMuon() ) && (qualityBits & muTools.kTMLastStationTight);
264  
265    // 0x10
266    level++;
267 <  if(!isGlobal && !isTracker) {
267 >  if((!isGlobal && !isTracker) || !(mu->HasTrackerTrk())) {
268      failmask |= (1<<level);
269 +    status.setStatus(SelectionStatus::FAIL);
270 +    return status;    
271    }
272  
273    // 0x20
274    level++;
275 <  if( mu->nTkHits < 10 ) {
275 >  assert(mu->HasTrackerTrk());
276 >  assert(mu->TrackerTrk());
277 >  if( mu->TrackerTrk()->NHits() < 10 ) {
278      failmask |= (1<<level);
279    }
280  
281    level++;
282 <  if(mu->nPixHits         < 1) {
282 >  if(muTrk->NPixelHits()          < 1) {
283      failmask |= (1<<level);
284    }
285  
286  
287    level++;
288 <  if(fabs(mu->d0)>0.02)   {
288 >  if(fabs(muTrk->D0Corrected(vtx))>0.02)   {
289      failmask |= (1<<level);
290    }
291  
# Line 229 | Line 301 | SelectionStatus passMuonSelection( Contr
301    }
302    */
303  
304 <  SelectionStatus status;
305 <  if( !failmask ) status.setStatus(SelectionStatus::LOOSEID);
234 <  if( !failmask ) status.setStatus(SelectionStatus::TIGHTID);
304 >  if( !failmask ) status.orStatus(SelectionStatus::LOOSEID);
305 >  if( !failmask ) status.orStatus(SelectionStatus::TIGHTID);
306    return status;
307  
308  
309   };
310  
240 SelectionStatus passMuonSelectionBackDoor( ControlFlags &ctrl, const mithep::TMuon * mu ) {
241  int level=0;
242  unsigned failmask=0x0;
311  
312 <  if(mu->pt < 5) {
313 <    failmask |= (1<<level);
314 <  }
312 > //--------------------------------------------------------------------------------------------------
313 > SelectionStatus muonIDMVASelection(ControlFlags &ctrl,
314 >                                 const mithep::Muon *mu,
315 >                                 const mithep::Vertex & vtx   )
316 > //--------------------------------------------------------------------------------------------------
317 > {
318 >  assert(mu->IsTrackerMuon()); // must be
319 >  assert(mu->HasTrackerTrk());
320 >  assert(mu->TrackerTrk());
321  
322 <  level++;
323 <  if(fabs(mu->eta) > 2.5) {
324 <    failmask |= (1<<level);
325 <  }
322 >  /*
323 >  cerr << "TrackerTrk: " << mu->TrackerTrk() << endl;
324 >  flush(cerr);
325 >  cerr << "Rchi2: " << mu->TrackerTrk()->RChi2() << endl;
326 >  flush(cerr);
327  
328 <  level++;
329 <  if(!(mu->typeBits & kGlobal) || !(mu->typeBits & kTracker)) {
330 <    failmask |= (1<<level);
331 <  }
328 >  cerr << "calling MVa ... "<< endl;
329 >  flush(cerr);
330 >  */
331 >
332 >  double global_rchi2 = (mu->IsGlobalMuon()) ? mu->GlobalTrk()->RChi2() : 0.;
333 >
334 >  double mvaval = muIDMVA->MVAValue_ID(mu->Pt(),
335 >                                    mu->Eta(),
336 >                                    mu->IsGlobalMuon(),
337 >                                    mu->IsTrackerMuon(),
338 >                                    mu->TrackerTrk()->RChi2(),
339 >                                    global_rchi2,
340 >                                    (Double_t)(mu->NValidHits()),
341 >                                    (Double_t)(mu->TrackerTrk()->NHits()),
342 >                                    (Double_t)(mu->TrackerTrk()->NPixelHits()),
343 >                                    (Double_t)(mu->NMatches()),
344 >                                    mu->TrkKink(),
345 >                                    muTools.GetSegmentCompatability(mu),
346 >                                    muTools.GetCaloCompatability(mu,true,true),
347 >                                    mu->HadEnergy(),
348 >                                    mu->EmEnergy(),
349 >                                    mu->HadS9Energy(),
350 >                                    mu->EmS9Energy(),
351 >                                    (Bool_t)(ctrl.debug) );
352  
258  level++;
259  if(mu->nTkHits          <= 10) {
260    failmask |= (1<<level);
261  }
353  
263        level++;
264        if(mu->ip3dSig > 4) {
265                failmask |= (1<<level);
266        }
267        
268        level++;
269  if(mu->muNchi2          >= 10)  {
270                failmask |= (1<<level);
271        }
272                
273        level++;
274  if(mu->nValidHits       < 1)   {
275                failmask |= (1<<level);
276        }
277        
278        level++;
279  if(mu->nMatch   < 2)   {
280                failmask |= (1<<level);
281        }
282  
354    SelectionStatus status;
355 <        
356 <  if( !failmask ) status.setStatus(SelectionStatus::LOOSEID);
357 <  if( !failmask ) status.setStatus(SelectionStatus::TIGHTID);
358 <        
359 <        return status;
355 >  bool pass = false;
356 >  
357 >  if( mu->IsGlobalMuon() && mu->IsTrackerMuon()
358 >      && fabs(mu->Eta()) <= 1.5 && mu->Pt() <= 10 && mvaval >= MUON_IDMVA_CUT_BIN0)  pass = true;
359 >  else if( mu->IsGlobalMuon() && mu->IsTrackerMuon()
360 >           && fabs(mu->Eta()) <= 1.5 && mu->Pt() > 10 && mvaval >= MUON_IDMVA_CUT_BIN1)  pass = true;
361 >  else if( mu->IsGlobalMuon() && mu->IsTrackerMuon()
362 >           && fabs(mu->Eta()) > 1.5 && mu->Pt() <= 10 && mvaval >= MUON_IDMVA_CUT_BIN2)  pass = true;
363 >  else if( mu->IsGlobalMuon() && mu->IsTrackerMuon()
364 >           && fabs(mu->Eta()) > 1.5 && mu->Pt() > 10 && mvaval >= MUON_IDMVA_CUT_BIN3)  pass = true;
365 >  else if( !(mu->IsGlobalMuon()) && mu->IsTrackerMuon()
366 >           && (mu->Quality().QualityMask().Mask() & mithep::MuonQuality::AllArbitrated) && mvaval >= MUON_IDMVA_CUT_BIN4)
367 >    pass = true;
368 >
369 >
370 >  if( pass ) {
371 >    status.orStatus(SelectionStatus::LOOSEID);
372 >    status.orStatus(SelectionStatus::TIGHTID);
373 >  }
374  
375 < };
375 >  return status;
376 > }
377 >
378 > //--------------------------------------------------------------------------------------------------
379 > void initMuonIDMVA()
380 > //--------------------------------------------------------------------------------------------------
381 > {
382 >  muIDMVA = new mithep::MuonIDMVA();
383 >  vector<string> weightFiles;
384 >  weightFiles.push_back("./data/MuonIDMVAWeights/MuonIDMVA_BDTG_barrel_lowpt_V2.weights.xml");
385 >  weightFiles.push_back("./data/MuonIDMVAWeights/MuonIDMVA_BDTG_barrel_highpt_V2.weights.xml");
386 >  weightFiles.push_back("./data/MuonIDMVAWeights/MuonIDMVA_BDTG_endcap_lowpt_V2.weights.xml");
387 >  weightFiles.push_back("./data/MuonIDMVAWeights/MuonIDMVA_BDTG_endcap_highpt_V2.weights.xml");
388 >  weightFiles.push_back("./data/MuonIDMVAWeights/MuonIDMVA_BDTG_tracker_V2.weights.xml");
389 >  weightFiles.push_back("./data/MuonIDMVAWeights/MuonIDMVA_BDTG_global_V2.weights.xml");
390 >  muIDMVA->Initialize( "MuonIDMVA",
391 >                       mithep::MuonIDMVA::kIDV0,
392 >                       kTRUE, weightFiles);
393 > }

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines