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 // destructor -----------------------------------------------------------------
0065 AlignmentGlobalTrackSelector::~AlignmentGlobalTrackSelector() {}
0066
0067 ///returns if any of the Filters is used.
0068 bool AlignmentGlobalTrackSelector::useThisFilter() {
0069   return theGMFilterSwitch || theIsoFilterSwitch || theJetCountFilterSwitch;
0070 }
0071
0072 // do selection ---------------------------------------------------------------
0073 AlignmentGlobalTrackSelector::Tracks AlignmentGlobalTrackSelector::select(const Tracks& tracks,
0074                                                                           const edm::Event& iEvent,
0075                                                                           const edm::EventSetup& iSetup) {
0076   Tracks result = tracks;
0077
0078   if (theGMFilterSwitch)
0079     result = findMuons(result, iEvent);
0080   if (theIsoFilterSwitch)
0081     result = checkIsolation(result, iEvent);
0082   if (theJetCountFilterSwitch)
0083     result = checkJetCount(result, iEvent);
0084   LogDebug("Alignment") << ">  Global: tracks all, kept: " << tracks.size() << ", " << result.size();
0085   //  LogDebug("Alignment")<<">  o kept:";
0086   //  printTracks(result);
0087
0088   return result;
0089 }
0090
0091 ///filter for Tracks that match the Track of a global Muon
0092 AlignmentGlobalTrackSelector::Tracks AlignmentGlobalTrackSelector::findMuons(const Tracks& tracks,
0093                                                                              const edm::Event& iEvent) const {
0094   Tracks result;
0095   Tracks globalMuons;
0096
0097   //fill globalMuons with muons
0098   Handle<reco::MuonCollection> muons;
0099   iEvent.getByToken(theMuonToken, muons);
0100
0101   if (muons.isValid()) {
0102     for (reco::MuonCollection::const_iterator itMuon = muons->begin(); itMuon != muons->end(); ++itMuon) {
0103       const reco::Track* muonTrack = (*itMuon).get<reco::TrackRef>().get();
0104       if (!muonTrack) {
0105         LogDebug("Alignment") << "@SUB=AlignmentGlobalTrackSelector::findMuons"
0106                               << "Found muon without track: Standalone Muon!";
0107       } else {
0108         globalMuons.push_back(muonTrack);
0109       }
0110     }
0111   } else {
0112     LogError("Alignment") << "@SUB=AlignmentGlobalTrackSelector::findMuons"
0113                           << ">  could not optain mounCollection!";
0114   }
0115
0116   result = this->matchTracks(tracks, globalMuons);
0117
0118   if (static_cast<int>(result.size()) < theMinGlobalMuonCount)
0119     result.clear();
0120
0121   return result;
0122 }
0123
0124 ///returns only isolated tracks in [cands]
0125 AlignmentGlobalTrackSelector::Tracks AlignmentGlobalTrackSelector::checkIsolation(const Tracks& cands,
0126                                                                                   const edm::Event& iEvent) const {
0127   Tracks result;
0128   result.clear();
0129
0130   Handle<reco::CaloJetCollection> jets;
0131   iEvent.getByToken(theJetIsoToken, jets);
0132
0133   if (jets.isValid()) {
0134     for (Tracks::const_iterator it = cands.begin(); it != cands.end(); ++it) {
0135       bool isolated = true;
0136       for (reco::CaloJetCollection::const_iterator itJet = jets->begin(); itJet != jets->end(); ++itJet)
0137         isolated &= !((*itJet).pt() > theMaxJetPt && deltaR(*(*it), (*itJet)) < theMinJetDeltaR);
0138
0139       if (isolated)
0140         result.push_back(*it);
0141     }
0142     //    LogDebug("Alignment") << "D  Found "<<result.size()<<" isolated of "<< cands.size()<<" Tracks!";
0143
0144   } else
0145     LogError("Alignment") << "@SUB=AlignmentGlobalTrackSelector::checkIsolation"
0146                           << "> could not optain jetCollection!";
0147
0148   if (static_cast<int>(result.size()) < theMinIsolatedCount)
0149     result.clear();
0150
0151   return result;
0152 }
0153
0154 ///returns [tracks] if there are less than theMaxCount Jets with theMinJetPt and an empty set if not
0155 AlignmentGlobalTrackSelector::Tracks AlignmentGlobalTrackSelector::checkJetCount(const Tracks& tracks,
0156                                                                                  const edm::Event& iEvent) const {
0157   Tracks result;
0158   result.clear();
0159
0160   Handle<reco::CaloJetCollection> jets;
0161   iEvent.getByToken(theJetCountToken, jets);
0162
0163   if (jets.isValid()) {
0164     int jetCount = 0;
0165     for (reco::CaloJetCollection::const_iterator itJet = jets->begin(); itJet != jets->end(); ++itJet) {
0166       if ((*itJet).pt() > theMinJetPt)
0167         jetCount++;
0168     }
0169
0170     if (jetCount <= theMaxJetCount)
0171       result = tracks;
0172
0173     LogDebug("Alignment") << "> found " << jetCount << " Jets";
0174   } else
0175     LogError("Alignment") << "@SUB=AlignmentGlobalTrackSelector::checkJetCount"
0176                           << ">  could not optain jetCollection!";
0177
0178   return result;
0179 }
0180
0181 //===================HELPERS===================
0182
0183 ///matches [src] with [comp] returns collection with matching Tracks coming from [src]
0184 AlignmentGlobalTrackSelector::Tracks AlignmentGlobalTrackSelector::matchTracks(const Tracks& src,
0185                                                                                const Tracks& comp) const {
0186   Tracks result;
0187   for (Tracks::const_iterator itComp = comp.begin(); itComp != comp.end(); ++itComp) {
0188     int match = -1;
0189     double min = theMaxTrackDeltaR;
0190     for (unsigned int i = 0; i < src.size(); i++) {
0191       // LogDebug("Alignment") << "> Trackmatch dist: "<<deltaR(src.at(i),*itComp);
0192       if (min > deltaR(*(src.at(i)), *(*itComp))) {
0193         min = deltaR(*(src.at(i)), *(*itComp));
0194         match = static_cast<int>(i);
0195       }
0196     }
0197     if (match > -1)
0198       result.push_back(src.at(match));
0199   }
0200   return result;
0201 }
0202
0203 ///print Information on Track-Collection
0204 void AlignmentGlobalTrackSelector::printTracks(const Tracks& col) const {
0205   int count = 0;
0206   LogDebug("Alignment") << ">......................................";
0207   for (Tracks::const_iterator it = col.begin(); it < col.end(); ++it, ++count) {
0208     LogDebug("Alignment") << ">  Track No. " << count << ": p = (" << (*it)->px() << "," << (*it)->py() << ","
0209                           << (*it)->pz() << ")\n"
0210                           << ">                        pT = " << (*it)->pt() << " eta = " << (*it)->eta()
0211                           << " charge = " << (*it)->charge();
0212   }
0213   LogDebug("Alignment") << ">......................................";
0214 }