L3MuonProducer/plugins/QuarkoniaTrackSelector.cc

0001 #include "RecoMuon/L3MuonProducer/plugins/QuarkoniaTrackSelector.h"
0002
0003 #include "FWCore/Framework/interface/Frameworkfwd.h"
0004
0005 #include "FWCore/Framework/interface/Event.h"
0006 #include "DataFormats/Common/interface/Handle.h"
0007
0008 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0009 #include "FWCore/MessageLogger/interface/MessageDrop.h"
0010
0011 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0012
0013 #include <memory>
0014 #include <iostream>
0015 #include <sstream>
0016
0017 QuarkoniaTrackSelector::QuarkoniaTrackSelector(const edm::ParameterSet& iConfig)
0018     : muonTag_(iConfig.getParameter<edm::InputTag>("muonCandidates")),
0019       trackTag_(iConfig.getParameter<edm::InputTag>("tracks")),
0020       minMasses_(iConfig.getParameter<std::vector<double> >("MinMasses")),
0021       maxMasses_(iConfig.getParameter<std::vector<double> >("MaxMasses")),
0022       checkCharge_(iConfig.getParameter<bool>("checkCharge")),
0023       minTrackPt_(iConfig.getParameter<double>("MinTrackPt")),
0024       minTrackP_(iConfig.getParameter<double>("MinTrackP")),
0025       maxTrackEta_(iConfig.getParameter<double>("MaxTrackEta")) {
0026   muonToken_ = consumes<reco::RecoChargedCandidateCollection>(muonTag_);
0027   trackToken_ = consumes<reco::TrackCollection>(trackTag_);
0028
0029   //register your products
0030   produces<reco::TrackCollection>();
0031   //
0032   // verify mass windows
0033   //
0034   bool massesValid = minMasses_.size() == maxMasses_.size();
0035   if (massesValid) {
0036     for (size_t i = 0; i < minMasses_.size(); ++i) {
0037       if (minMasses_[i] < 0 || maxMasses_[i] < 0 || minMasses_[i] > maxMasses_[i])
0038         massesValid = false;
0039     }
0040   }
0041   if (!massesValid) {
0042     edm::LogError("QuarkoniaTrackSelector") << "Inconsistency in definition of mass windows, "
0043                                             << "no track will be selected";
0044     minMasses_.clear();
0045     maxMasses_.clear();
0046   }
0047
0048   std::ostringstream stream;
0049   stream << "instantiated with parameters\n"
0050          << "  muonTag  = " << muonTag_ << "\n"
0051          << "  trackTag = " << trackTag_ << "\n";
0052   stream << "  mass windows =";
0053   for (size_t i = 0; i < minMasses_.size(); ++i)
0054     stream << " (" << minMasses_[i] << "," << maxMasses_[i] << ")";
0055   stream << "\n";
0056   stream << "  checkCharge  = " << checkCharge_ << "\n";
0057   stream << "  MinTrackPt = " << minTrackPt_ << "\n";
0058   stream << "  MinTrackP = " << minTrackP_ << "\n";
0059   stream << "  MaxTrackEta = " << maxTrackEta_;
0060   LogDebug("QuarkoniaTrackSelector") << stream.str();
0061 }
0062
0063 void QuarkoniaTrackSelector::produce(edm::StreamID, edm::Event& iEvent, const edm::EventSetup& iSetup) const {
0064   //
0065   // the product
0066   //
0067   auto product = std::make_unique<reco::TrackCollection>();
0068   //
0069   // Muons
0070   //
0071   edm::Handle<reco::RecoChargedCandidateCollection> muonHandle;
0072   iEvent.getByToken(muonToken_, muonHandle);
0073   //
0074   // Tracks
0075   //
0076   edm::Handle<reco::TrackCollection> trackHandle;
0077   iEvent.getByToken(trackToken_, trackHandle);
0078   //
0079   // Verification
0080   //
0081   if (!muonHandle.isValid() || !trackHandle.isValid() || minMasses_.empty()) {
0082     iEvent.put(std::move(product));
0083     return;
0084   }
0085   //
0086   // Debug output
0087   //
0088   if (edm::isDebugEnabled()) {
0089     std::ostringstream stream;
0090     stream << "\nInput muons: # / q / pt / p / eta\n";
0091     for (size_t im = 0; im < muonHandle->size(); ++im) {
0092       const reco::RecoChargedCandidate& muon = (*muonHandle)[im];
0093       stream << "   " << im << " " << muon.charge() << " " << muon.pt() << " " << muon.p() << " " << muon.eta() << "\n";
0094     }
0095     stream << "Input tracks: # / q / pt / p / eta\n";
0096     for (size_t it = 0; it < trackHandle->size(); ++it) {
0097       const reco::Track& track = (*trackHandle)[it];
0098       stream << "   " << it << " " << track.charge() << " " << track.pt() << " " << track.p() << " " << track.eta()
0099              << "\n";
0100     }
0101     LogDebug("QuarkoniaTrackSelector") << stream.str();
0102   }
0103   //
0104   // combinations
0105   //
0106   //   std::ostringstream stream;
0107   unsigned int nQ(0);
0108   unsigned int nComb(0);
0109   std::vector<size_t> selectedTrackIndices;
0110   selectedTrackIndices.reserve(muonHandle->size());
0111   reco::Particle::LorentzVector p4Muon;
0112   reco::Particle::LorentzVector p4JPsi;
0113   // muons
0114   for (size_t im = 0; im < muonHandle->size(); ++im) {
0115     const reco::RecoChargedCandidate& muon = (*muonHandle)[im];
0116     int qMuon = muon.charge();
0117     p4Muon = muon.p4();
0118     // tracks
0119     for (size_t it = 0; it < trackHandle->size(); ++it) {
0120       const reco::Track& track = (*trackHandle)[it];
0121       if (track.pt() < minTrackPt_ || track.p() < minTrackP_ || fabs(track.eta()) > maxTrackEta_)
0122         continue;
0123       if (checkCharge_ && track.charge() != -qMuon)
0124         continue;
0125       ++nQ;
0126       reco::Particle::LorentzVector p4Track(
0127           track.px(), track.py(), track.pz(), sqrt(track.p() * track.p() + 0.0111636));
0128       // mass windows
0129       double mass = (p4Muon + p4Track).mass();
0130       //       stream << "Combined mass = " << im << " " << it
0131       //         << " " << mass
0132       //         << " phi " << track.phi() << "\n";
0133       for (size_t j = 0; j < minMasses_.size(); ++j) {
0134         if (mass > minMasses_[j] && mass < maxMasses_[j]) {
0135           ++nComb;
0136           if (find(selectedTrackIndices.begin(), selectedTrackIndices.end(), it) == selectedTrackIndices.end())
0137             selectedTrackIndices.push_back(it);
0138           //      stream << "... adding " << "\n";
0139           break;
0140         }
0141       }
0142     }
0143   }
0144   //   LogDebug("QuarkoniaTrackSelector") << stream.str();
0145   //
0146   // filling of output collection
0147   //
0148   for (size_t i = 0; i < selectedTrackIndices.size(); ++i)
0149     product->push_back((*trackHandle)[selectedTrackIndices[i]]);
0150   //
0151   // debug output
0152   //
0153   if (edm::isDebugEnabled()) {
0154     std::ostringstream stream;
0155     stream << "Total number of combinations = " << muonHandle->size() * trackHandle->size() << " , after charge " << nQ
0156            << " , after mass " << nComb << std::endl;
0157     stream << "Selected " << product->size() << " tracks with # / q / pt / eta\n";
0158     for (size_t i = 0; i < product->size(); ++i) {
0159       const reco::Track& track = (*product)[i];
0160       stream << "  " << i << " " << track.charge() << " " << track.pt() << " " << track.eta() << "\n";
0161     }
0162     LogDebug("QuarkoniaTrackSelector") << stream.str();
0163   }
0164   //
0165   iEvent.put(std::move(product));
0166 }
0167
0168 //define this as a plug-in
0169 #include "FWCore/Framework/interface/MakerMacros.h"
0170 DEFINE_FWK_MODULE(QuarkoniaTrackSelector);