CommonAlignmentProducer/src/AlignmentGlobalTrackSelector.cc

0001 //Framework
0002 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0003 #include "FWCore/Framework/interface/Event.h"
0004 #include "FWCore/Utilities/interface/EDMException.h"
0005 #include "FWCore/Utilities/interface/InputTag.h"
0006
0007 //DataFormats
0008 #include <DataFormats/Candidate/interface/Particle.h>
0009 #include <DataFormats/Candidate/interface/Candidate.h>
0010 #include <DataFormats/TrackReco/interface/Track.h>
0011 #include <DataFormats/JetReco/interface/CaloJet.h>
0012
0013 #include <DataFormats/MuonReco/interface/Muon.h>
0014
0015 #include <DataFormats/RecoCandidate/interface/RecoCandidate.h>  //for the get<TrackRef>() Call
0016
0017 #include <DataFormats/Math/interface/deltaR.h>
0018
0019 //STL
0020 #include <cmath>
0021
0022 #include "Alignment/CommonAlignmentProducer/interface/AlignmentGlobalTrackSelector.h"
0023
0024 using namespace std;
0025 using namespace edm;
0026
0027 // constructor ----------------------------------------------------------------
0028 AlignmentGlobalTrackSelector::AlignmentGlobalTrackSelector(const edm::ParameterSet& cfg, edm::ConsumesCollector& iC)
0029     : theGMFilterSwitch(cfg.getParameter<bool>("applyGlobalMuonFilter")),
0030       theIsoFilterSwitch(cfg.getParameter<bool>("applyIsolationtest")),
0031       theJetCountFilterSwitch(cfg.getParameter<bool>("applyJetCountFilter")) {
0032   if (theGMFilterSwitch || theIsoFilterSwitch || theJetCountFilterSwitch)
0033     LogDebug("Alignment") << "> applying global Trackfilter ...";
0034
0035   if (theGMFilterSwitch) {
0036     edm::InputTag theMuonSource = cfg.getParameter<InputTag>("muonSource");
0037     theMuonToken = iC.consumes<reco::MuonCollection>(theMuonSource);
0038     theMaxTrackDeltaR = cfg.getParameter<double>("maxTrackDeltaR");
0039     theMinGlobalMuonCount = cfg.getParameter<int>("minGlobalMuonCount");
0040     LogDebug("Alignment") << ">  GlobalMuonFilter : source, maxTrackDeltaR, min. Count       : " << theMuonSource
0041                           << " , " << theMaxTrackDeltaR << " , " << theMinIsolatedCount;
0042   }
0043
0044   if (theIsoFilterSwitch) {
0045     edm::InputTag theJetIsoSource = cfg.getParameter<InputTag>("jetIsoSource");
0046     theJetIsoToken = iC.consumes<reco::CaloJetCollection>(theJetIsoSource);
0047     theMaxJetPt = cfg.getParameter<double>("maxJetPt");
0048     theMinJetDeltaR = cfg.getParameter<double>("minJetDeltaR");
0049     theMinIsolatedCount = cfg.getParameter<int>("minIsolatedCount");
0050     LogDebug("Alignment") << ">  Isolationtest    : source, maxJetPt, minJetDeltaR, min. Count: " << theJetIsoSource
0051                           << " , " << theMaxJetPt << " ," << theMinJetDeltaR << " ," << theMinGlobalMuonCount;
0052   }
0053
0054   if (theJetCountFilterSwitch) {
0055     edm::InputTag theJetCountSource = cfg.getParameter<InputTag>("jetCountSource");
0056     theJetCountToken = iC.consumes<reco::CaloJetCollection>(theJetCountSource);
0057     theMinJetPt = cfg.getParameter<double>("minJetPt");
0058     theMaxJetCount = cfg.getParameter<int>("maxJetCount");
0059     LogDebug("Alignment") << ">  JetCountFilter   : source, minJetPt, maxJetCount             : " << theJetCountSource
0060                           << " , " << theMinJetPt << " ," << theMaxJetCount;
0061   }
0062 }
0063
0064 void AlignmentGlobalTrackSelector::fillPSetDescription(edm::ParameterSetDescription& desc) {
0065   // Global muon finding
0066   desc.add<bool>("applyGlobalMuonFilter", false);
0067   desc.add<edm::InputTag>("muonSource", edm::InputTag("muons"));
0068   desc.add<double>("maxTrackDeltaR", 0.001);
0069   desc.add<int>("minGlobalMuonCount", 1);
0070
0071   // Isolation tests
0072   desc.add<bool>("applyIsolationtest", false);
0073   desc.add<edm::InputTag>("jetIsoSource", edm::InputTag("kt6CaloJets"));
0074   desc.add<double>("maxJetPt", 40.0);  // GeV
0075   desc.add<double>("minJetDeltaR", 0.2);
0076   desc.add<int>("minIsolatedCount", 0);
0077
0078   // Jet count filter
0079   desc.add<bool>("applyJetCountFilter", false);
0080   desc.add<edm::InputTag>("jetCountSource", edm::InputTag("kt6CaloJets"));
0081   desc.add<double>("minJetPt", 40.0);  // GeV
0082   desc.add<int>("maxJetCount", 3);
0083 }
0084
0085 // destructor -----------------------------------------------------------------
0086 AlignmentGlobalTrackSelector::~AlignmentGlobalTrackSelector() {}
0087
0088 ///returns if any of the Filters is used.
0089 bool AlignmentGlobalTrackSelector::useThisFilter() {
0090   return theGMFilterSwitch || theIsoFilterSwitch || theJetCountFilterSwitch;
0091 }
0092
0093 // do selection ---------------------------------------------------------------
0094 AlignmentGlobalTrackSelector::Tracks AlignmentGlobalTrackSelector::select(const Tracks& tracks,
0095                                                                           const edm::Event& iEvent,
0096                                                                           const edm::EventSetup& iSetup) {
0097   Tracks result = tracks;
0098
0099   if (theGMFilterSwitch)
0100     result = findMuons(result, iEvent);
0101   if (theIsoFilterSwitch)
0102     result = checkIsolation(result, iEvent);
0103   if (theJetCountFilterSwitch)
0104     result = checkJetCount(result, iEvent);
0105   LogDebug("Alignment") << ">  Global: tracks all, kept: " << tracks.size() << ", " << result.size();
0106   //  LogDebug("Alignment")<<">  o kept:";
0107   //  printTracks(result);
0108
0109   return result;
0110 }
0111
0112 ///filter for Tracks that match the Track of a global Muon
0113 AlignmentGlobalTrackSelector::Tracks AlignmentGlobalTrackSelector::findMuons(const Tracks& tracks,
0114                                                                              const edm::Event& iEvent) const {
0115   Tracks result;
0116   Tracks globalMuons;
0117
0118   //fill globalMuons with muons
0119   Handle<reco::MuonCollection> muons;
0120   iEvent.getByToken(theMuonToken, muons);
0121
0122   if (muons.isValid()) {
0123     for (reco::MuonCollection::const_iterator itMuon = muons->begin(); itMuon != muons->end(); ++itMuon) {
0124       const reco::Track* muonTrack = (*itMuon).get<reco::TrackRef>().get();
0125       if (!muonTrack) {
0126         LogDebug("Alignment") << "@SUB=AlignmentGlobalTrackSelector::findMuons"
0127                               << "Found muon without track: Standalone Muon!";
0128       } else {
0129         globalMuons.push_back(muonTrack);
0130       }
0131     }
0132   } else {
0133     LogError("Alignment") << "@SUB=AlignmentGlobalTrackSelector::findMuons"
0134                           << ">  could not optain mounCollection!";
0135   }
0136
0137   result = this->matchTracks(tracks, globalMuons);
0138
0139   if (static_cast<int>(result.size()) < theMinGlobalMuonCount)
0140     result.clear();
0141
0142   return result;
0143 }
0144
0145 ///returns only isolated tracks in [cands]
0146 AlignmentGlobalTrackSelector::Tracks AlignmentGlobalTrackSelector::checkIsolation(const Tracks& cands,
0147                                                                                   const edm::Event& iEvent) const {
0148   Tracks result;
0149   result.clear();
0150
0151   Handle<reco::CaloJetCollection> jets;
0152   iEvent.getByToken(theJetIsoToken, jets);
0153
0154   if (jets.isValid()) {
0155     for (Tracks::const_iterator it = cands.begin(); it != cands.end(); ++it) {
0156       bool isolated = true;
0157       for (reco::CaloJetCollection::const_iterator itJet = jets->begin(); itJet != jets->end(); ++itJet)
0158         isolated &= !((*itJet).pt() > theMaxJetPt && deltaR(*(*it), (*itJet)) < theMinJetDeltaR);
0159
0160       if (isolated)
0161         result.push_back(*it);
0162     }
0163     //    LogDebug("Alignment") << "D  Found "<<result.size()<<" isolated of "<< cands.size()<<" Tracks!";
0164
0165   } else
0166     LogError("Alignment") << "@SUB=AlignmentGlobalTrackSelector::checkIsolation"
0167                           << "> could not optain jetCollection!";
0168
0169   if (static_cast<int>(result.size()) < theMinIsolatedCount)
0170     result.clear();
0171
0172   return result;
0173 }
0174
0175 ///returns [tracks] if there are less than theMaxCount Jets with theMinJetPt and an empty set if not
0176 AlignmentGlobalTrackSelector::Tracks AlignmentGlobalTrackSelector::checkJetCount(const Tracks& tracks,
0177                                                                                  const edm::Event& iEvent) const {
0178   Tracks result;
0179   result.clear();
0180
0181   Handle<reco::CaloJetCollection> jets;
0182   iEvent.getByToken(theJetCountToken, jets);
0183
0184   if (jets.isValid()) {
0185     int jetCount = 0;
0186     for (reco::CaloJetCollection::const_iterator itJet = jets->begin(); itJet != jets->end(); ++itJet) {
0187       if ((*itJet).pt() > theMinJetPt)
0188         jetCount++;
0189     }
0190
0191     if (jetCount <= theMaxJetCount)
0192       result = tracks;
0193
0194     LogDebug("Alignment") << "> found " << jetCount << " Jets";
0195   } else
0196     LogError("Alignment") << "@SUB=AlignmentGlobalTrackSelector::checkJetCount"
0197                           << ">  could not optain jetCollection!";
0198
0199   return result;
0200 }
0201
0202 //===================HELPERS===================
0203
0204 ///matches [src] with [comp] returns collection with matching Tracks coming from [src]
0205 AlignmentGlobalTrackSelector::Tracks AlignmentGlobalTrackSelector::matchTracks(const Tracks& src,
0206                                                                                const Tracks& comp) const {
0207   Tracks result;
0208   for (Tracks::const_iterator itComp = comp.begin(); itComp != comp.end(); ++itComp) {
0209     int match = -1;
0210     double min = theMaxTrackDeltaR;
0211     for (unsigned int i = 0; i < src.size(); i++) {
0212       // LogDebug("Alignment") << "> Trackmatch dist: "<<deltaR(src.at(i),*itComp);
0213       if (min > deltaR(*(src.at(i)), *(*itComp))) {
0214         min = deltaR(*(src.at(i)), *(*itComp));
0215         match = static_cast<int>(i);
0216       }
0217     }
0218     if (match > -1)
0219       result.push_back(src.at(match));
0220   }
0221   return result;
0222 }
0223
0224 ///print Information on Track-Collection
0225 void AlignmentGlobalTrackSelector::printTracks(const Tracks& col) const {
0226   int count = 0;
0227   LogDebug("Alignment") << ">......................................";
0228   for (Tracks::const_iterator it = col.begin(); it < col.end(); ++it, ++count) {
0229     LogDebug("Alignment") << ">  Track No. " << count << ": p = (" << (*it)->px() << "," << (*it)->py() << ","
0230                           << (*it)->pz() << ")\n"
0231                           << ">                        pT = " << (*it)->pt() << " eta = " << (*it)->eta()
0232                           << " charge = " << (*it)->charge();
0233   }
0234   LogDebug("Alignment") << ">......................................";
0235 }