PatExamples/plugins/PatMCMatchingExtended.cc

0001 #include <map>
0002 #include <string>
0003
0004 #include "TH1.h"
0005
0006 #include "FWCore/Framework/interface/Event.h"
0007 #include "FWCore/Framework/interface/one/EDAnalyzer.h"
0008 #include "FWCore/Utilities/interface/InputTag.h"
0009 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0010 #include "FWCore/ServiceRegistry/interface/Service.h"
0011 #include "CommonTools/UtilAlgos/interface/TFileService.h"
0012 #include "Math/VectorUtil.h"
0013
0014 #include "DataFormats/PatCandidates/interface/Muon.h"
0015
0016 class PatMCMatchingExtended : public edm::one::EDAnalyzer<edm::one::SharedResources> {
0017 public:
0018   /// default constructor
0019   explicit PatMCMatchingExtended(const edm::ParameterSet&);
0020   /// default destructor
0021   ~PatMCMatchingExtended() override;
0022
0023 private:
0024   void beginJob() override;
0025   void analyze(const edm::Event&, const edm::EventSetup&) override;
0026   void endJob() override;
0027
0028   // simple map to contain all histograms;
0029   // histograms are booked in the beginJob()
0030   // method
0031   std::map<std::string, TH1F*> histContainer_;
0032
0033   // input tags
0034   edm::EDGetTokenT<edm::View<pat::Muon> > muonSrcToken_;
0035
0036   //counts how often a genParticle with different charge gives a match
0037   unsigned int diffCharge;
0038
0039   //how many muons have no match
0040   unsigned int noMatch;
0041
0042   //how many muons have no status 1 or 3 match, but decay in flight
0043   unsigned int decayInFlight;
0044
0045   //count the number of muons in all events
0046   unsigned int numberMuons;
0047 };
0048
0049 PatMCMatchingExtended::PatMCMatchingExtended(const edm::ParameterSet& iConfig)
0050     : histContainer_(),
0051       muonSrcToken_(consumes<edm::View<pat::Muon> >(iConfig.getUntrackedParameter<edm::InputTag>("muonSrc"))) {
0052   usesResource(TFileService::kSharedResource);
0053 }
0054
0055 PatMCMatchingExtended::~PatMCMatchingExtended() {}
0056
0057 void PatMCMatchingExtended::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0058   // get muon collection
0059   edm::Handle<edm::View<pat::Muon> > muons;
0060   iEvent.getByToken(muonSrcToken_, muons);
0061
0062   for (edm::View<pat::Muon>::const_iterator muon = muons->begin(); muon != muons->end(); ++muon) {
0063     if (muon->genParticleById(0, 1).isNonnull()) {
0064       histContainer_["DR_status1Match"]->Fill(
0065           ROOT::Math::VectorUtil::DeltaR(muon->p4(), (muon->genParticleById(0, 1))->p4()));
0066       histContainer_["DPt_status1Match"]->Fill(muon->pt() - (muon->genParticleById(0, 1))->pt());
0067     }
0068     if (muon->genParticleById(0, 3).isNonnull()) {
0069       histContainer_["DR_status3Match"]->Fill(
0070           ROOT::Math::VectorUtil::DeltaR(muon->p4(), (muon->genParticleById(0, 3))->p4()));
0071       histContainer_["DPt_status3Match"]->Fill(muon->pt() - (muon->genParticleById(0, 3))->pt());
0072     }
0073     if (muon->genParticleById(0, -1).isNonnull()) {
0074       histContainer_["DR_defaultMatch"]->Fill(
0075           ROOT::Math::VectorUtil::DeltaR(muon->p4(), (muon->genParticleById(0, -1))->p4()));
0076       histContainer_["DPt_defaultMatch"]->Fill(muon->pt() - (muon->genParticleById(0, -1))->pt());
0077     }
0078     if (muon->genParticleById(0, 1).isNull() && muon->genParticleById(0, 3).isNull() &&
0079         muon->genParticleById(0, -1).isNull())
0080       noMatch++;
0081     if (muon->genParticleById(0, 1).isNull() && muon->genParticleById(0, 3).isNull() &&
0082         muon->genParticleById(0, -1).isNonnull())
0083       decayInFlight++;
0084
0085     if (muon->genParticleById(-13, 0, 1).isNonnull()) {
0086       diffCharge++;
0087       std::cout << " DIFF CHARGE!!! charge gen: " << muon->genParticleById(-13, 0, true)->charge()
0088                 << " charge reco: " << muon->charge() << std::endl;
0089     }
0090     numberMuons++;
0091   }
0092 }
0093
0094 void PatMCMatchingExtended::beginJob() {
0095   // register to the TFileService
0096   edm::Service<TFileService> fs;
0097
0098   // book histograms:
0099   //DR
0100   histContainer_["DR_defaultMatch"] = fs->make<TH1F>("DR_defaultMatch", "DR_defaultMatch", 100, 0, 0.02);
0101   histContainer_["DR_status1Match"] = fs->make<TH1F>("DR_status1Match", "DR_status1Match", 100, 0, 0.02);
0102   histContainer_["DR_status3Match"] = fs->make<TH1F>("DR_status3Match", "DR_status3Match", 100, 0, 0.02);
0103   //DPT
0104   histContainer_["DPt_defaultMatch"] = fs->make<TH1F>("DPt_defaultMatch", "DPt_defaultMatch", 10, 0, 1.2);
0105   histContainer_["DPt_status1Match"] = fs->make<TH1F>("DPt_status1Match", "DPt_status1Match", 10, 0, 1.2);
0106   histContainer_["DPt_status3Match"] = fs->make<TH1F>("DPt_status3Match", "DPt_status3Match", 10, 0, 1.2);
0107   //some counters
0108   diffCharge = 0;
0109   noMatch = 0;
0110   decayInFlight = 0;
0111   numberMuons = 0;
0112 }
0113
0114 void PatMCMatchingExtended::endJob() {
0115   std::cout << "diffcharge: " << diffCharge << std::endl;
0116   std::cout << "noMatch: " << noMatch << std::endl;
0117   std::cout << "decayInFlight: " << decayInFlight << std::endl;
0118   std::cout << "numberMuons: " << numberMuons << std::endl;
0119 }
0120
0121 #include "FWCore/Framework/interface/MakerMacros.h"
0122 DEFINE_FWK_MODULE(PatMCMatchingExtended);