/ElectroWeakAnalysis/ZMuMu/plugins/ZGlobalVsSAIsolationAnalyzer.cc

https://github.com/aivanov-cern/cmssw · C++ · 128 lines · 116 code · 10 blank · 2 comment · 5 complexity · 55afc04c6756e72acbe2a0fe3f32eb6f MD5 · raw file

  1. #include "DataFormats/Common/interface/AssociationVector.h"
  2. #include "FWCore/Framework/interface/EDAnalyzer.h"
  3. #include "FWCore/Utilities/interface/EDMException.h"
  4. #include "DataFormats/Candidate/interface/Candidate.h"
  5. #include "DataFormats/Candidate/interface/Particle.h"
  6. #include "DataFormats/Candidate/interface/CandidateFwd.h"
  7. #include "FWCore/ParameterSet/interface/ParameterSet.h"
  8. #include "FWCore/Framework/interface/Event.h"
  9. #include "FWCore/Framework/interface/EventSetup.h"
  10. #include "FWCore/Utilities/interface/InputTag.h"
  11. #include "FWCore/ServiceRegistry/interface/Service.h"
  12. #include "CommonTools/UtilAlgos/interface/TFileService.h"
  13. #include "DataFormats/Math/interface/deltaR.h"
  14. #include "DataFormats/PatCandidates/interface/Muon.h"
  15. #include "DataFormats/TrackReco/interface/Track.h"
  16. #include "DataFormats/HepMCCandidate/interface/GenParticle.h"
  17. #include "DataFormats/RecoCandidate/interface/IsoDeposit.h"
  18. #include "DataFormats/RecoCandidate/interface/IsoDepositFwd.h"
  19. #include "DataFormats/PatCandidates/interface/Isolation.h"
  20. #include "DataFormats/Common/interface/ValueMap.h"
  21. #include "DataFormats/PatCandidates/interface/GenericParticle.h"
  22. #include "DataFormats/RecoCandidate/interface/IsoDepositVetos.h"
  23. #include "DataFormats/RecoCandidate/interface/IsoDepositDirection.h"
  24. #include <vector>
  25. #include <string>
  26. #include <iostream>
  27. using namespace edm;
  28. using namespace std;
  29. using namespace reco;
  30. using namespace isodeposit;
  31. class ZGlobalVsSAIsolationAnalyzer : public edm::EDAnalyzer {
  32. public:
  33. typedef math::XYZVector Vector;
  34. ZGlobalVsSAIsolationAnalyzer(const edm::ParameterSet& cfg);
  35. private:
  36. virtual void analyze(const edm::Event& event, const edm::EventSetup& setup) override;
  37. virtual void endJob() override;
  38. InputTag src_;
  39. double dRVeto;
  40. double dRTrk, dREcal, dRHcal;
  41. double ptThreshold, etEcalThreshold, etHcalThreshold;
  42. double alpha, beta;
  43. double isoCut_;
  44. unsigned long selGlobal_, selSA_, totGlobal_, totSA_;
  45. bool isolated(const Direction & dir, const pat::IsoDeposit * trkIsoDep,
  46. const pat::IsoDeposit * ecalIsoDep, const pat::IsoDeposit * hcalIsoDep);
  47. void evaluate(const reco::Candidate* dau);
  48. };
  49. ZGlobalVsSAIsolationAnalyzer::ZGlobalVsSAIsolationAnalyzer(const ParameterSet& cfg):
  50. src_(cfg.getParameter<InputTag>("src")),
  51. dRVeto(cfg.getParameter<double>("veto")),
  52. dRTrk(cfg.getParameter<double>("deltaRTrk")),
  53. dREcal(cfg.getParameter<double>("deltaREcal")),
  54. dRHcal(cfg.getParameter<double>("deltaRHcal")),
  55. ptThreshold(cfg.getParameter<double>("ptThreshold")),
  56. etEcalThreshold(cfg.getParameter<double>("etEcalThreshold")),
  57. etHcalThreshold(cfg.getParameter<double>("etHcalThreshold")),
  58. alpha(cfg.getParameter<double>("alpha")),
  59. beta(cfg.getParameter<double>("beta")),
  60. isoCut_(cfg.getParameter<double>("isoCut")),
  61. selGlobal_(0), selSA_(0), totGlobal_(0), totSA_(0) {
  62. }
  63. bool ZGlobalVsSAIsolationAnalyzer::isolated(const Direction & dir, const pat::IsoDeposit * trkIsoDep,
  64. const pat::IsoDeposit * ecalIsoDep, const pat::IsoDeposit * hcalIsoDep) {
  65. IsoDeposit::AbsVetos vetoTrk, vetoEcal, vetoHcal;
  66. vetoTrk.push_back(new ConeVeto(dir, dRVeto));
  67. vetoTrk.push_back(new ThresholdVeto(ptThreshold));
  68. vetoEcal.push_back(new ConeVeto(dir, 0.));
  69. vetoEcal.push_back(new ThresholdVeto(etEcalThreshold));
  70. vetoHcal.push_back(new ConeVeto(dir, 0.));
  71. vetoHcal.push_back(new ThresholdVeto(etHcalThreshold));
  72. double trkIso = trkIsoDep->sumWithin(dir, dRTrk, vetoTrk);
  73. double ecalIso = ecalIsoDep->sumWithin(dir, dREcal, vetoEcal);
  74. double hcalIso = hcalIsoDep->sumWithin(dir, dRHcal, vetoHcal);
  75. double iso = alpha*((0.5*(1+beta)*ecalIso) + (0.5*(1-beta)*hcalIso)) + (1-alpha)*trkIso;
  76. return iso < isoCut_;
  77. }
  78. void ZGlobalVsSAIsolationAnalyzer::evaluate(const reco::Candidate* dau) {
  79. const pat::Muon * mu = dynamic_cast<const pat::Muon *>(&*dau->masterClone());
  80. if(mu == 0) throw Exception(errors::InvalidReference) << "Daughter is not a muon!\n";
  81. const pat::IsoDeposit * trkIsoDep = mu->isoDeposit(pat::TrackIso);
  82. const pat::IsoDeposit * ecalIsoDep = mu->isoDeposit(pat::EcalIso);
  83. const pat::IsoDeposit * hcalIsoDep = mu->isoDeposit(pat::HcalIso);
  84. // global muon
  85. {
  86. Direction dir = Direction(mu->eta(), mu->phi());
  87. if(isolated(dir, trkIsoDep, ecalIsoDep, hcalIsoDep)) selGlobal_++;
  88. totGlobal_++;
  89. }
  90. // stand-alone
  91. {
  92. TrackRef sa = dau->get<TrackRef,reco::StandAloneMuonTag>();
  93. Direction dir = Direction(sa->eta(), sa->phi());
  94. if(isolated(dir, trkIsoDep, ecalIsoDep, hcalIsoDep)) selSA_++;
  95. totSA_++;
  96. }
  97. }
  98. void ZGlobalVsSAIsolationAnalyzer::analyze(const edm::Event& event, const edm::EventSetup& setup) {
  99. Handle<CandidateView> dimuons;
  100. event.getByLabel(src_, dimuons);
  101. for(unsigned int i=0; i< dimuons->size(); ++ i) {
  102. const Candidate & zmm = (* dimuons)[i];
  103. evaluate(zmm.daughter(0));
  104. evaluate(zmm.daughter(1));
  105. }
  106. }
  107. void ZGlobalVsSAIsolationAnalyzer::endJob() {
  108. cout << "Isolation efficiency report:" << endl;
  109. double eff, err;
  110. eff = double(selGlobal_)/double(totGlobal_);
  111. err = sqrt(eff*(1.-eff)/double(totGlobal_));
  112. cout <<"Global: " << selGlobal_ << "/" << totGlobal_ << " = " <<eff <<"+/-" << err<< endl;
  113. eff = double(selSA_)/double(totSA_);
  114. err = sqrt(eff*(1.-eff)/double(totSA_));
  115. cout <<"St.Al.: " << selSA_ << "/" << totSA_ << " = " << eff <<"+/-" << err << endl;
  116. }
  117. #include "FWCore/Framework/interface/MakerMacros.h"
  118. DEFINE_FWK_MODULE(ZGlobalVsSAIsolationAnalyzer);