ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/MitPhysics/Mods/src/MuonIDMod.cc
(Generate patch)

Comparing UserCode/MitPhysics/Mods/src/MuonIDMod.cc (file contents):
Revision 1.55 by sixie, Fri Jul 22 15:45:38 2011 UTC vs.
Revision 1.76 by ceballos, Sat May 5 08:51:06 2012 UTC

# Line 14 | Line 14 | ClassImp(mithep::MuonIDMod)
14   //--------------------------------------------------------------------------------------------------
15    MuonIDMod::MuonIDMod(const char *name, const char *title) :
16    BaseMod(name,title),
17 +  fPrintMVADebugInfo(kFALSE),
18    fMuonBranchName(Names::gkMuonBrn),
19    fCleanMuonsName(ModNames::gkCleanMuonsName),  
20    fNonIsolatedMuonsName("random"),  
# Line 22 | Line 23 | ClassImp(mithep::MuonIDMod)
23    fBeamSpotName(Names::gkBeamSpotBrn),
24    fTrackName(Names::gkTrackBrn),
25    fPFCandidatesName(Names::gkPFCandidatesBrn),
26 <  fMuonIDType("WWMuIdV1"),
26 >  fMuonIDType("WWMuIdV3"),
27    fMuonIsoType("PFIso"),
28    fMuonClassType("Global"),  
29    fTrackIsolationCut(3.0),
# Line 45 | Line 46 | ClassImp(mithep::MuonIDMod)
46    fBeamSpot(0),
47    fTracks(0),
48    fPFCandidates(0),
49 +  fPFNoPileUpCands(0),
50 +  fIntRadius(0.0),
51    fNonIsolatedMuons(0),
52    fNonIsolatedElectrons(0),
53    fPileupEnergyDensityName(Names::gkPileupEnergyDensityBrn),
54 <  fPileupEnergyDensity(0)
54 >  fPileupEnergyDensity(0),
55 >  fMuonTools(0),
56 >  fMuonIDMVA(0),
57 >  fTheRhoType(RhoUtilities::DEFAULT)
58   {
59    // Constructor.
60   }
# Line 67 | Line 73 | void MuonIDMod::Process()
73    LoadEventObject(fBeamSpotName, fBeamSpot);
74    LoadEventObject(fTrackName, fTracks);
75    LoadEventObject(fPFCandidatesName, fPFCandidates);
76 <  if(fMuIsoType == kTrackCaloSliding || fMuIsoType == kCombinedRelativeConeAreaCorrected) {
76 >  if(fMuIsoType == kTrackCaloSliding ||
77 >     fMuIsoType == kCombinedRelativeConeAreaCorrected ||
78 >     fMuIsoType == kPFIsoEffectiveAreaCorrected ||
79 >     fMuIsoType == kMVAIso_BDTG_IDIso ||
80 >     fMuIsoType == kIsoRingsV0_BDTG_Iso ||
81 >     fMuIsoType == kIsoDeltaR
82 >    ) {
83      LoadEventObject(fPileupEnergyDensityName, fPileupEnergyDensity);
84    }
85 +  if(fMuIsoType == kPFRadialIso || fMuIsoType == kIsoDeltaR){
86 +    // Name is hardcoded, can be changed if someone feels to do it
87 +    fPFNoPileUpCands = GetObjThisEvt<PFCandidateCol>("PFNoPileUp");    
88 +  }
89 +
90    MuonOArr *CleanMuons = new MuonOArr;
91    CleanMuons->SetName(fCleanMuonsName);
92  
93    fVertices = GetObjThisEvt<VertexOArr>(fVertexName);
94  
95 <  for (UInt_t i=0; i<fMuons->GetEntries(); ++i) {
95 >  for (UInt_t i=0; i<fMuons->GetEntries() && fVertices->GetEntries() > 0 ; ++i) {
96      const Muon *mu = fMuons->At(i);
97  
98      Bool_t pass = kFALSE;
# Line 156 | Line 173 | void MuonIDMod::Process()
173      if (eta >= fEtaCut)
174        continue;
175  
176 +
177 +    //***********************************************************************************************
178 +    //Debug Info For Lepton MVA
179 +    //***********************************************************************************************
180 +    if( fPrintMVADebugInfo &&
181 +        (fMuIsoType == kMVAIso_BDTG_IDIso || fMuIDType == kMVAID_BDTG_IDIso)
182 +      ) {
183 +      cout << "Event: " << GetEventHeader()->RunNum() << " " << GetEventHeader()->LumiSec() << " "
184 +           << GetEventHeader()->EvtNum() << " : Rho = " << fPileupEnergyDensity->At(0)->Rho()
185 +           << " : Muon " << i << " "
186 +           << endl;
187 +      fMuonIDMVA->MVAValue(mu,fVertices->At(0),fMuonTools,fPFCandidates,fPileupEnergyDensity,kTRUE);
188 +    }
189 +    //***********************************************************************************************
190 +
191 +
192      Double_t RChi2 = 0.0;
193      if     (mu->HasGlobalTrk()) {
194        RChi2 = mu->GlobalTrk()->Chi2()/mu->GlobalTrk()->Ndof();
# Line 211 | Line 244 | void MuonIDMod::Process()
244                   mu->BestTrk()->NPixelHits() > 0 &&
245                   mu->BestTrk()->PtErr()/mu->BestTrk()->Pt() < 0.1;
246          break;
247 <
247 >      case kWWMuIdV3:
248 >        idpass = mu->BestTrk() != 0 &&
249 >                 mu->BestTrk()->NHits() > 10 &&
250 >                 mu->BestTrk()->NPixelHits() > 0 &&
251 >                 mu->BestTrk()->PtErr()/mu->BestTrk()->Pt() < 0.1 &&
252 >                 mu->TrkKink() < 20.0;
253 >        break;
254 >      case kMVAID_BDTG_IDIso:
255 >        {
256 >          Bool_t passDenominatorM2 = (mu->BestTrk() != 0 &&
257 >                                      mu->BestTrk()->NHits() > 10 &&
258 >                                      mu->BestTrk()->NPixelHits() > 0 &&
259 >                                      mu->BestTrk()->PtErr()/mu->BestTrk()->Pt() < 0.1 &&
260 >                                      MuonTools::PassD0Cut(mu, fVertices, 0.20, 0) &&
261 >                                      MuonTools::PassDZCut(mu, fVertices, 0.10, 0) &&
262 >                                      mu->TrkKink() < 20.0
263 >            );  
264 >          idpass =  passDenominatorM2;
265 >          //only evaluate MVA if muon passes M2 denominator to save time
266 >          if (idpass) idpass = PassMuonMVA_BDTG_IdIso(mu, fVertices->At(0), fPileupEnergyDensity);
267 >        }
268 >        break;
269        case kNoId:
270          idpass = kTRUE;
271          break;
# Line 222 | Line 276 | void MuonIDMod::Process()
276      if (!idpass)
277        continue;
278  
279 +    Double_t Rho = 0.;
280 +    if(fPileupEnergyDensity){
281 +      const PileupEnergyDensity *rho =  fPileupEnergyDensity->At(0);
282 +
283 +      switch (fTheRhoType) {
284 +      case RhoUtilities::MIT_RHO_VORONOI_LOW_ETA:
285 +        Rho = rho->RhoLowEta();
286 +        break;
287 +      case RhoUtilities::MIT_RHO_VORONOI_HIGH_ETA:
288 +        Rho = rho->Rho();
289 +        break;
290 +      case RhoUtilities::MIT_RHO_RANDOM_LOW_ETA:
291 +        Rho = rho->RhoRandomLowEta();
292 +        break;
293 +      case RhoUtilities::MIT_RHO_RANDOM_HIGH_ETA:
294 +        Rho = rho->RhoRandom();
295 +        break;
296 +      default:
297 +        Rho = rho->Rho();
298 +      }
299 +
300 +      if ((TMath::IsNaN(fPileupEnergyDensity->At(0)->Rho()) || std::isinf(fPileupEnergyDensity->At(0)->Rho()))) Rho = 0.;
301 +    }
302 +
303      Bool_t isocut = kFALSE;
304      switch (fMuIsoType) {
305        case kTrackCalo:
# Line 235 | Line 313 | void MuonIDMod::Process()
313          break;
314        case kTrackCaloSliding:
315          {
316 <          const PileupEnergyDensity *rho =  fPileupEnergyDensity->At(0);
239 <          Double_t totalIso =  mu->IsoR03SumPt() + mu->IsoR03EmEt() + mu->IsoR03HadEt() - rho->Rho() * TMath::Pi() * 0.3 * 0.3 ;
316 >          Double_t totalIso =  mu->IsoR03SumPt() + mu->IsoR03EmEt() + mu->IsoR03HadEt() - Rho * TMath::Pi() * 0.3 * 0.3 ;
317            // trick to change the signal region cut
318            double theIsoCut = fCombIsolationCut;
319            if(theIsoCut < 0.20){
# Line 262 | Line 339 | void MuonIDMod::Process()
339          break;
340        case kCombinedRelativeConeAreaCorrected:
341          {
342 <          const PileupEnergyDensity *rho =  fPileupEnergyDensity->At(0);
343 <          Double_t totalIso =  mu->IsoR03SumPt() + mu->IsoR03EmEt() + mu->IsoR03HadEt() - rho->Rho() * TMath::Pi() * 0.3 * 0.3 ;
342 >          //const PileupEnergyDensity *rho =  fPileupEnergyDensity->At(0); // Fabian: made Rho customable
343 >          Double_t totalIso =  mu->IsoR03SumPt() + mu->IsoR03EmEt() + mu->IsoR03HadEt() - Rho * TMath::Pi() * 0.3 * 0.3 ;
344            double theIsoCut = fCombRelativeIsolationCut;
345            if (totalIso < (mu->Pt()*theIsoCut)) isocut = kTRUE;
346          }
347          break;          
348 +    case kCombinedRelativeEffectiveAreaCorrected:
349 +      {
350 +        Double_t tmpRho = Rho;   // Fabian: made the Rho type customable.
351 +        //if (!(TMath::IsNaN(fPileupEnergyDensity->At(0)->Rho()) || std::isinf(fPileupEnergyDensity->At(0)->Rho())))
352 +        //tmpRho = fPileupEnergyDensity->At(0)->Rho();
353 +        
354 +          isocut = ( mu->IsoR03SumPt() + mu->IsoR03EmEt() + mu->IsoR03HadEt()
355 +                     -  tmpRho*MuonTools::MuonEffectiveArea(MuonTools::kMuEMIso03, mu->Eta())
356 +                     -  tmpRho*MuonTools::MuonEffectiveArea(MuonTools::kMuHadIso03, mu->Eta())
357 +            ) < (mu->Pt()* 0.40);
358 +        }
359 +        break;          
360        case kPFIso:
361 <      {
361 >        {
362            Double_t pfIsoCutValue = 9999;
363            if(fPFIsolationCut > 0){
364              pfIsoCutValue = fPFIsolationCut;
# Line 288 | Line 377 | void MuonIDMod::Process()
377                }
378              }
379            }
380 <          Double_t totalIso =  IsolationTools::PFMuonIsolation(mu, fPFCandidates, fVertices->At(0), 0.1, 1.0, 0.3, 0.0);
380 >          Double_t totalIso =  IsolationTools::PFMuonIsolation(mu, fPFCandidates, fVertices->At(0), 0.1, 1.0, 0.3, 0.0, fIntRadius);
381 >          if (totalIso < (mu->Pt()*pfIsoCutValue) )
382 >            isocut = kTRUE;
383 >        }
384 >        break;
385 >      case kPFRadialIso:
386 >        {
387 >          Double_t pfIsoCutValue = 9999;
388 >          if(fPFIsolationCut > 0){
389 >            pfIsoCutValue = fPFIsolationCut;
390 >          } else {
391 >            if (mu->Pt() > 20) {
392 >              pfIsoCutValue = 0.10;
393 >            } else {
394 >              pfIsoCutValue = 0.05;
395 >            }
396 >          }
397 >          Double_t totalIso =  IsolationTools::PFRadialMuonIsolation(mu, fPFNoPileUpCands, 1.0, 0.3);
398            if (totalIso < (mu->Pt()*pfIsoCutValue) )
399              isocut = kTRUE;
400          }
401          break;
402 +      case kPFIsoEffectiveAreaCorrected:
403 +        {
404 +          Double_t pfIsoCutValue = 9999;
405 +          if(fPFIsolationCut > 0){
406 +            pfIsoCutValue = fPFIsolationCut;
407 +          } else {
408 +            pfIsoCutValue = fPFIsolationCut; //leave it like this for now
409 +          }
410 +          Double_t EffectiveAreaCorrectedPFIso =  IsolationTools::PFMuonIsolation(mu, fPFCandidates, fVertices->At(0), 0.1, 1.0, 0.3, 0.0, fIntRadius)
411 +            - Rho * MuonTools::MuonEffectiveArea(MuonTools::kMuNeutralIso03, mu->Eta());
412 +          //- fPileupEnergyDensity->At(0)->Rho() * MuonTools::MuonEffectiveArea(MuonTools::kMuNeutralIso03, mu->Eta());  // Fabian: made Rho-type customable
413 +          isocut = EffectiveAreaCorrectedPFIso < (mu->Pt() * pfIsoCutValue);
414 +          break;
415 +        }
416        case kPFIsoNoL:
417          {
418            fNonIsolatedMuons     = GetObjThisEvt<MuonCol>(fNonIsolatedMuonsName);
419            fNonIsolatedElectrons = GetObjThisEvt<ElectronCol>(fNonIsolatedElectronsName);
420  
421 <          Double_t beta = IsolationTools::BetaM(fTracks, mu, fVertices->At(0), 0.0, 0.2, 0.3, 0.02);
422 <          if(beta == 0) beta = 1.0;
423 <          Double_t totalIso =  IsolationTools::PFMuonIsolation(mu, fPFCandidates, fVertices->At(0), fNonIsolatedMuons, fNonIsolatedElectrons, 0.2, 1.0, 0.4, 0.0, 3);
424 <          if (totalIso < (mu->Pt()*fPFIsolationCut) )
421 >          Double_t pfIsoCutValue = 9999;
422 >          if(fPFIsolationCut > 0){
423 >            pfIsoCutValue = fPFIsolationCut;
424 >          } else {
425 >            if (mu->AbsEta() < 1.479) {
426 >              if (mu->Pt() > 20) {
427 >                pfIsoCutValue = 0.13;
428 >              } else {
429 >                pfIsoCutValue = 0.06;
430 >              }
431 >            } else {
432 >              if (mu->Pt() > 20) {
433 >                pfIsoCutValue = 0.09;
434 >              } else {
435 >                pfIsoCutValue = 0.05;
436 >              }
437 >            }
438 >          }
439 >          Double_t totalIso =  IsolationTools::PFMuonIsolation(mu, fPFCandidates, fNonIsolatedMuons, fNonIsolatedElectrons, fVertices->At(0), 0.1, 1.0, 0.3, 0.0, fIntRadius);
440 >          if (totalIso < (mu->Pt()*pfIsoCutValue) )
441              isocut = kTRUE;
442          }
443          break;
444 +      case kMVAIso_BDTG_IDIso:
445 +      {
446 +
447 +        // **************************************************************************
448 +        // Don't use effective area correction denominator. Instead use the old one.
449 +        // **************************************************************************
450 +
451 +        //         Double_t tmpRho = 0;
452 +        //         if (!(TMath::IsNaN(fPileupEnergyDensity->At(0)->Rho()) || isinf(fPileupEnergyDensity->At(0)->Rho())))
453 +        //           tmpRho = fPileupEnergyDensity->At(0)->Rho();
454 +        
455 +        //         isocut = ( mu->IsoR03SumPt() + mu->IsoR03EmEt() + mu->IsoR03HadEt()
456 +        //                    -  tmpRho*MuonTools::MuonEffectiveArea(MuonTools::kMuEMIso03, mu->Eta())
457 +        //                    -  tmpRho*MuonTools::MuonEffectiveArea(MuonTools::kMuHadIso03, mu->Eta())
458 +        //           ) < (mu->Pt()* 0.40);
459 +        
460 +        Double_t totalIso =  IsolationTools::PFMuonIsolation(mu, fPFCandidates, fVertices->At(0), 0.1, 1.0, 0.3, 0.0, fIntRadius);
461 +        isocut = (totalIso < (mu->Pt()*0.4));
462 +
463 +      }
464 +        break;
465 +      case kIsoRingsV0_BDTG_Iso:
466 +      {
467 +        
468 +        isocut = PassMuonIsoRingsV0_BDTG_Iso(mu, fVertices->At(0), fPileupEnergyDensity);
469 +
470 +      }
471 +        break;
472 +      case kIsoDeltaR:
473 +      {
474 +        
475 +        isocut = PassMuonIsoDeltaR(mu, fVertices->At(0), fPileupEnergyDensity);
476 +
477 +      }
478 +        break;
479        case kNoIso:
480          isocut = kTRUE;
481          break;
# Line 361 | Line 532 | void MuonIDMod::SlaveBegin()
532    ReqEventObject(fTrackName, fTracks, kTRUE);
533    ReqEventObject(fPFCandidatesName, fPFCandidates, kTRUE);
534    if (fMuonIsoType.CompareTo("TrackCaloSliding") == 0
535 <      || fMuonIsoType.CompareTo("CombinedRelativeConeAreaCorrected") == 0) {
535 >      || fMuonIsoType.CompareTo("CombinedRelativeConeAreaCorrected") == 0
536 >      || fMuonIsoType.CompareTo("CombinedRelativeEffectiveAreaCorrected") == 0
537 >      || fMuonIsoType.CompareTo("PFIsoEffectiveAreaCorrected") == 0
538 >      || fMuonIsoType.CompareTo("MVA_BDTG_IDIso") == 0
539 >      || fMuonIsoType.CompareTo("IsoRingsV0_BDTG_Iso") == 0
540 >      || fMuonIsoType.CompareTo("IsoDeltaR") == 0
541 >    ) {
542      ReqEventObject(fPileupEnergyDensityName, fPileupEnergyDensity, kTRUE);
543    }
544  
545 +
546    if (fMuonIDType.CompareTo("WMuId") == 0)
547      fMuIDType = kWMuId;
548    else if (fMuonIDType.CompareTo("ZMuId") == 0)
# Line 377 | Line 555 | void MuonIDMod::SlaveBegin()
555      fMuIDType = kWWMuIdV1;
556    else if (fMuonIDType.CompareTo("WWMuIdV2") == 0)
557      fMuIDType = kWWMuIdV2;
558 +  else if (fMuonIDType.CompareTo("WWMuIdV3") == 0)
559 +    fMuIDType = kWWMuIdV3;
560    else if (fMuonIDType.CompareTo("NoId") == 0)
561      fMuIDType = kNoId;
562    else if (fMuonIDType.CompareTo("Custom") == 0) {
563      fMuIDType = kCustomId;
564      SendError(kWarning, "SlaveBegin",
565                "Custom muon identification is not yet implemented.");
566 +  } else if (fMuonIDType.CompareTo("MVA_BDTG_IDIso") == 0) {
567 +    fMuIDType = kMVAID_BDTG_IDIso;
568    } else {
569      SendError(kAbortAnalysis, "SlaveBegin",
570                "The specified muon identification %s is not defined.",
# Line 400 | Line 582 | void MuonIDMod::SlaveBegin()
582      fMuIsoType = kTrackCaloSlidingNoCorrection;
583    else if (fMuonIsoType.CompareTo("CombinedRelativeConeAreaCorrected") == 0)
584      fMuIsoType = kCombinedRelativeConeAreaCorrected;
585 +  else if (fMuonIsoType.CompareTo("CombinedRelativeEffectiveAreaCorrected") == 0)
586 +    fMuIsoType = kCombinedRelativeEffectiveAreaCorrected;
587    else if (fMuonIsoType.CompareTo("PFIso") == 0)
588      fMuIsoType = kPFIso;
589 +  else if (fMuonIsoType.CompareTo("PFRadialIso") == 0)
590 +    fMuIsoType = kPFRadialIso;
591 +  else if (fMuonIsoType.CompareTo("PFIsoEffectiveAreaCorrected") == 0)
592 +    fMuIsoType = kPFIsoEffectiveAreaCorrected;
593    else if (fMuonIsoType.CompareTo("PFIsoNoL") == 0)
594      fMuIsoType = kPFIsoNoL;
595    else if (fMuonIsoType.CompareTo("NoIso") == 0)
# Line 410 | Line 598 | void MuonIDMod::SlaveBegin()
598      fMuIsoType = kCustomIso;
599      SendError(kWarning, "SlaveBegin",
600                "Custom muon isolation is not yet implemented.");
601 +  } else if (fMuonIsoType.CompareTo("MVA_BDTG_IDIso") == 0) {
602 +    fMuIsoType = kMVAIso_BDTG_IDIso;
603 +  } else if (fMuonIsoType.CompareTo("IsoRingsV0_BDTG_Iso") == 0) {
604 +    fMuIsoType = kIsoRingsV0_BDTG_Iso;
605 +  } else if (fMuonIsoType.CompareTo("IsoDeltaR") == 0) {
606 +    fMuIsoType = kIsoDeltaR;
607    } else {
608      SendError(kAbortAnalysis, "SlaveBegin",
609                "The specified muon isolation %s is not defined.",
# Line 437 | Line 631 | void MuonIDMod::SlaveBegin()
631                fMuonClassType.Data());
632      return;
633    }
634 +
635 +
636 +  //If we use MVA ID, need to load MVA weights
637 +  if     (fMuIsoType == kMVAIso_BDTG_IDIso || fMuIDType == kMVAID_BDTG_IDIso) {
638 +    fMuonTools = new MuonTools();
639 +    fMuonIDMVA = new MuonIDMVA();
640 +    fMuonIDMVA->Initialize("BDTG method",
641 +                           string((getenv("CMSSW_BASE")+string("/src/MitPhysics/data/MuonMVAWeights/BarrelPtBin0_IDIsoCombined_BDTG.weights.xml"))),
642 +                           string((getenv("CMSSW_BASE")+string("/src/MitPhysics/data/MuonMVAWeights/EndcapPtBin0_IDIsoCombined_BDTG.weights.xml"))),
643 +                           string((getenv("CMSSW_BASE")+string("/src/MitPhysics/data/MuonMVAWeights/BarrelPtBin1_IDIsoCombined_BDTG.weights.xml"))),
644 +                           string((getenv("CMSSW_BASE")+string("/src/MitPhysics/data/MuonMVAWeights/EndcapPtBin1_IDIsoCombined_BDTG.weights.xml"))),
645 +                           string((getenv("CMSSW_BASE")+string("/src/MitPhysics/data/MuonMVAWeights/BarrelPtBin2_IDIsoCombined_BDTG.weights.xml"))),
646 +                           string((getenv("CMSSW_BASE")+string("/src/MitPhysics/data/MuonMVAWeights/EndcapPtBin2_IDIsoCombined_BDTG.weights.xml"))),
647 +                           MuonIDMVA::kIDIsoCombinedDetIso);
648 +  }
649 +  else if(fMuIsoType == kIsoRingsV0_BDTG_Iso) {
650 +    std::vector<std::string> muonidiso_weightfiles;
651 +    muonidiso_weightfiles.push_back(string((getenv("CMSSW_BASE")+string("/src/MitPhysics/data/MuonMVAWeights/MuonIsoMVA_BDTG_V0_barrel_lowpt.weights.xml"))));
652 +    muonidiso_weightfiles.push_back(string((getenv("CMSSW_BASE")+string("/src/MitPhysics/data/MuonMVAWeights/MuonIsoMVA_BDTG_V0_barrel_highpt.weights.xml"))));
653 +    muonidiso_weightfiles.push_back(string((getenv("CMSSW_BASE")+string("/src/MitPhysics/data/MuonMVAWeights/MuonIsoMVA_BDTG_V0_endcap_lowpt.weights.xml"))));
654 +    muonidiso_weightfiles.push_back(string((getenv("CMSSW_BASE")+string("/src/MitPhysics/data/MuonMVAWeights/MuonIsoMVA_BDTG_V0_endcap_highpt.weights.xml"))));
655 +    muonidiso_weightfiles.push_back(string((getenv("CMSSW_BASE")+string("/src/MitPhysics/data/MuonMVAWeights/MuonIsoMVA_BDTG_V0_tracker.weights.xml"))));
656 +    muonidiso_weightfiles.push_back(string((getenv("CMSSW_BASE")+string("/src/MitPhysics/data/MuonMVAWeights/MuonIsoMVA_BDTG_V0_global.weights.xml"))));
657 +    fMuonTools = new MuonTools();
658 +    fMuonIDMVA = new MuonIDMVA();
659 +    fMuonIDMVA->Initialize("MuonIso_BDTG_IsoRings",
660 +                       MuonIDMVA::kIsoRingsV0,
661 +                       kTRUE,
662 +                       muonidiso_weightfiles);
663 +  }
664 +  else if(fMuIsoType == kIsoDeltaR) {
665 +    std::vector<std::string> muonidiso_weightfiles;
666 +    muonidiso_weightfiles.push_back(string((getenv("CMSSW_BASE")+string("/src/MitPhysics/data/MuonMVAWeights/MuonIsoMVA_santi-V1_LB_BDT.weights.xml"))));
667 +    muonidiso_weightfiles.push_back(string((getenv("CMSSW_BASE")+string("/src/MitPhysics/data/MuonMVAWeights/MuonIsoMVA_santi-V1_LE_BDT.weights.xml"))));
668 +    muonidiso_weightfiles.push_back(string((getenv("CMSSW_BASE")+string("/src/MitPhysics/data/MuonMVAWeights/MuonIsoMVA_santi-V1_HB_BDT.weights.xml"))));
669 +    muonidiso_weightfiles.push_back(string((getenv("CMSSW_BASE")+string("/src/MitPhysics/data/MuonMVAWeights/MuonIsoMVA_santi-V1_HE_BDT.weights.xml"))));
670 +    fMuonTools = new MuonTools();
671 +    fMuonIDMVA = new MuonIDMVA();
672 +    fMuonIDMVA->Initialize("muonHZZ2012IsoDRMVA",
673 +                       MuonIDMVA::kIsoDeltaR,
674 +                       kTRUE,
675 +                       muonidiso_weightfiles);
676 +  }
677 +
678 + }
679 +
680 +
681 + //--------------------------------------------------------------------------------------------------
682 + Bool_t MuonIDMod::PassMuonMVA_BDTG_IdIso(const Muon *mu, const Vertex *vertex,
683 +                                         const PileupEnergyDensityCol *PileupEnergyDensity) const
684 + {
685 +
686 +  const Track *muTrk=0;
687 +  if(mu->HasTrackerTrk())         { muTrk = mu->TrackerTrk();    }
688 +  else if(mu->HasStandaloneTrk()) { muTrk = mu->StandaloneTrk(); }
689 +  
690 +  Double_t MVAValue = fMuonIDMVA->MVAValue(mu,vertex,fMuonTools,fPFCandidates,PileupEnergyDensity);
691 +
692 +  Int_t subdet = 0;
693 +  if (fabs(muTrk->Eta()) < 1.479) subdet = 0;
694 +  else subdet = 1;
695 +  Int_t ptBin = 0;
696 +  if (muTrk->Pt() > 14.5) ptBin = 1;
697 +  if (muTrk->Pt() > 20.0) ptBin = 2;
698 +
699 +  Int_t MVABin = -1;
700 +  if      (subdet == 0 && ptBin == 0) MVABin = 0;
701 +  else if (subdet == 1 && ptBin == 0) MVABin = 1;
702 +  else if (subdet == 0 && ptBin == 1) MVABin = 2;
703 +  else if (subdet == 1 && ptBin == 1) MVABin = 3;
704 +  else if (subdet == 0 && ptBin == 2) MVABin = 4;
705 +  else if (subdet == 1 && ptBin == 2) MVABin = 5;
706 +
707 +  Double_t MVACut = -999;
708 +  if      (MVABin == 0) MVACut = -0.5618;
709 +  else if (MVABin == 1) MVACut = -0.3002;
710 +  else if (MVABin == 2) MVACut = -0.4642;
711 +  else if (MVABin == 3) MVACut = -0.2478;
712 +  else if (MVABin == 4) MVACut =  0.1706;
713 +  else if (MVABin == 5) MVACut =  0.8146;
714 +
715 +  if (MVAValue > MVACut) return kTRUE;
716 +  return kFALSE;
717 + }
718 +
719 + //--------------------------------------------------------------------------------------------------
720 + Bool_t MuonIDMod::PassMuonIsoRingsV0_BDTG_Iso(const Muon *mu, const Vertex *vertex,
721 +                                              const PileupEnergyDensityCol *PileupEnergyDensity) const
722 + {
723 +
724 +  const Track *muTrk=0;
725 +  if(mu->HasTrackerTrk())         { muTrk = mu->TrackerTrk();    }
726 +  else if(mu->HasStandaloneTrk()) { muTrk = mu->StandaloneTrk(); }
727 +  
728 +  ElectronOArr *tempElectrons = new  ElectronOArr;
729 +  MuonOArr     *tempMuons     = new  MuonOArr;
730 +  Double_t MVAValue = fMuonIDMVA->MVAValue(mu,vertex,fMuonTools,fPFCandidates,
731 +                      PileupEnergyDensity,MuonTools::kMuEAFall11MC,tempElectrons,tempMuons,kFALSE);
732 +  delete tempElectrons;
733 +  delete tempMuons;
734 +
735 +  Int_t MVABin = fMuonIDMVA->GetMVABin(muTrk->Eta(), muTrk->Pt(), mu->IsGlobalMuon(), mu->IsTrackerMuon());
736 +
737 +  Double_t MVACut = -999;
738 +  if      (MVABin == 0) MVACut = -0.593;
739 +  else if (MVABin == 1) MVACut =  0.337;
740 +  else if (MVABin == 2) MVACut = -0.767;
741 +  else if (MVABin == 3) MVACut =  0.410;
742 +  else if (MVABin == 4) MVACut = -0.989;
743 +  else if (MVABin == 5) MVACut = -0.995;
744 +
745 +  if (MVAValue > MVACut) return kTRUE;
746 +  return kFALSE;
747 + }
748 +
749 + //--------------------------------------------------------------------------------------------------
750 + Bool_t MuonIDMod::PassMuonIsoDeltaR(const Muon *mu, const Vertex *vertex,
751 +                                    const PileupEnergyDensityCol *PileupEnergyDensity) const
752 + {
753 +
754 +  const Track *muTrk=0;
755 +  if(mu->HasTrackerTrk())         { muTrk = mu->TrackerTrk();    }
756 +  else if(mu->HasStandaloneTrk()) { muTrk = mu->StandaloneTrk(); }
757 +  
758 +  ElectronOArr *tempElectrons = new  ElectronOArr;
759 +  MuonOArr     *tempMuons     = new  MuonOArr;
760 +  Double_t MVAValue = fMuonIDMVA->MVAValue(mu,vertex,fMuonTools,fPFNoPileUpCands,
761 +                      PileupEnergyDensity,MuonTools::kMuEAFall11MC,tempElectrons,tempMuons,kTRUE);
762 +  delete tempElectrons;
763 +  delete tempMuons;
764 +
765 +  Int_t MVABin = fMuonIDMVA->GetMVABin(muTrk->Eta(), muTrk->Pt(), mu->IsGlobalMuon(), mu->IsTrackerMuon());
766 +
767 +  Double_t MVACut = -999;
768 +  if      (MVABin == 0) MVACut =  0.000;
769 +  else if (MVABin == 1) MVACut =  0.000;
770 +  else if (MVABin == 2) MVACut =  0.000;
771 +  else if (MVABin == 3) MVACut =  0.000;
772 +
773 +  if (MVAValue > MVACut) return kTRUE;
774 +  return kFALSE;
775 + }
776 +
777 + //--------------------------------------------------------------------------------------------------
778 + void MuonIDMod::Terminate()
779 + {
780 +  // Run finishing code on the computer (slave) that did the analysis
781 +  delete fMuonIDMVA;
782 +  
783 +  delete fMuonTools;
784   }

Diff Legend

Removed lines
+ Added lines
< Changed lines
> Changed lines