Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2025-01-08 03:35:56

0001 #ifndef CommonTools_RecoAlgos_RecoTrackSelectorBase_h
0002 #define CommonTools_RecoAlgos_RecoTrackSelectorBase_h
0003 
0004 #include "DataFormats/BeamSpot/interface/BeamSpot.h"
0005 #include "DataFormats/Math/interface/deltaPhi.h"
0006 #include "DataFormats/TrackReco/interface/Track.h"
0007 #include "DataFormats/TrackReco/interface/TrackFwd.h"
0008 #include "DataFormats/VertexReco/interface/Vertex.h"
0009 #include "DataFormats/VertexReco/interface/VertexFwd.h"
0010 #include "FWCore/Framework/interface/ConsumesCollector.h"
0011 #include "FWCore/Framework/interface/Event.h"
0012 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0013 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0014 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0015 #include "FWCore/Utilities/interface/InputTag.h"
0016 
0017 class RecoTrackSelectorBase {
0018 public:
0019   RecoTrackSelectorBase() {}
0020   RecoTrackSelectorBase(const edm::ParameterSet& cfg)
0021       : ptMin_(cfg.getParameter<double>("ptMin")),
0022         minRapidity_(cfg.getParameter<double>("minRapidity")),
0023         maxRapidity_(cfg.getParameter<double>("maxRapidity")),
0024         meanPhi_((cfg.getParameter<double>("minPhi") + cfg.getParameter<double>("maxPhi")) / 2.),
0025         rangePhi_((cfg.getParameter<double>("maxPhi") - cfg.getParameter<double>("minPhi")) / 2.),
0026         tip_(cfg.getParameter<double>("tip")),
0027         lip_(cfg.getParameter<double>("lip")),
0028         maxChi2_(cfg.getParameter<double>("maxChi2")),
0029         minHit_(cfg.getParameter<int>("minHit")),
0030         minPixelHit_(cfg.getParameter<int>("minPixelHit")),
0031         minLayer_(cfg.getParameter<int>("minLayer")),
0032         min3DLayer_(cfg.getParameter<int>("min3DLayer")),
0033         usePV_(false),
0034         invertRapidityCut_(cfg.getParameter<bool>("invertRapidityCut")) {
0035     const auto minPhi = cfg.getParameter<double>("minPhi");
0036     const auto maxPhi = cfg.getParameter<double>("maxPhi");
0037     if (minPhi >= maxPhi) {
0038       throw cms::Exception("Configuration")
0039           << "RecoTrackSelectorPhase: minPhi (" << minPhi << ") must be smaller than maxPhi (" << maxPhi
0040           << "). The range is constructed from minPhi to maxPhi around their average.";
0041     }
0042     if (minPhi >= M_PI) {
0043       throw cms::Exception("Configuration")
0044           << "RecoTrackSelectorPhase: minPhi (" << minPhi
0045           << ") must be smaller than PI. The range is constructed from minPhi to maxPhi around their average.";
0046     }
0047     if (maxPhi <= -M_PI) {
0048       throw cms::Exception("Configuration")
0049           << "RecoTrackSelectorPhase: maxPhi (" << maxPhi
0050           << ") must be larger than -PI. The range is constructed from minPhi to maxPhi around their average.";
0051     }
0052 
0053     for (const std::string& quality : cfg.getParameter<std::vector<std::string> >("quality"))
0054       quality_.push_back(reco::TrackBase::qualityByName(quality));
0055     for (const std::string& algorithm : cfg.getParameter<std::vector<std::string> >("algorithm"))
0056       algorithm_.push_back(reco::TrackBase::algoByName(algorithm));
0057     for (const std::string& algorithm : cfg.getParameter<std::vector<std::string> >("originalAlgorithm"))
0058       originalAlgorithm_.push_back(reco::TrackBase::algoByName(algorithm));
0059     for (const std::string& algorithm : cfg.getParameter<std::vector<std::string> >("algorithmMaskContains"))
0060       algorithmMask_.push_back(reco::TrackBase::algoByName(algorithm));
0061   }
0062 
0063   RecoTrackSelectorBase(const edm::ParameterSet& cfg, edm::ConsumesCollector& iC) : RecoTrackSelectorBase(cfg) {
0064     usePV_ = cfg.getParameter<bool>("usePV");
0065     bsSrcToken_ = iC.consumes<reco::BeamSpot>(cfg.getParameter<edm::InputTag>("beamSpot"));
0066     if (usePV_)
0067       vertexToken_ = iC.consumes<reco::VertexCollection>(cfg.getParameter<edm::InputTag>("vertexTag"));
0068   }
0069 
0070   void init(const edm::Event& event, const edm::EventSetup& es) {
0071     edm::Handle<reco::BeamSpot> beamSpot;
0072     event.getByToken(bsSrcToken_, beamSpot);
0073     vertex_ = beamSpot->position();
0074     if (!usePV_)
0075       return;
0076     edm::Handle<reco::VertexCollection> hVtx;
0077     event.getByToken(vertexToken_, hVtx);
0078     if (hVtx->empty())
0079       return;
0080     vertex_ = (*hVtx)[0].position();
0081   }
0082 
0083   static void fillPSetDescription(edm::ParameterSetDescription& desc) {
0084     desc.add<bool>("invertRapidityCut", false);
0085     desc.add<bool>("usePV", false);
0086     desc.add<double>("lip", 300.0);
0087     desc.add<double>("maxChi2", 10000.0);
0088     desc.add<double>("maxPhi", -3.2);
0089     desc.add<double>("maxRapidity", 5.0);
0090     desc.add<double>("minPhi", 3.2);
0091     desc.add<double>("minRapidity", -5.0);
0092     desc.add<double>("ptMin", 0.1);
0093     desc.add<double>("tip", 120.0);
0094     desc.add<edm::InputTag>("beamSpot", edm::InputTag("offlineBeamSpot"));
0095     desc.add<edm::InputTag>("vertexTag", edm::InputTag("offlinePrimaryVertices"));
0096     desc.add<int>("min3DLayer", 0);
0097     desc.add<int>("minHit", 0);
0098     desc.add<int>("minLayer", 3);
0099     desc.add<int>("minPixelHit", 0);
0100     desc.add<std::vector<std::string> >("algorithm", {});
0101     desc.add<std::vector<std::string> >("algorithmMaskContains", {});
0102     desc.add<std::vector<std::string> >("originalAlgorithm", {});
0103     desc.add<std::vector<std::string> >("quality", {});
0104   }
0105 
0106   bool operator()(const reco::TrackRef& tref) const { return (*this)(*tref); }
0107 
0108   bool operator()(const reco::Track& t) const { return (*this)(t, vertex_); }
0109 
0110   bool operator()(const reco::Track& t, const reco::Track::Point& vertex) const {
0111     bool quality_ok = true;
0112     if (!quality_.empty()) {
0113       quality_ok = false;
0114       for (unsigned int i = 0; i < quality_.size(); ++i) {
0115         if (t.quality(quality_[i])) {
0116           quality_ok = true;
0117           break;
0118         }
0119       }
0120     }
0121 
0122     bool algo_ok = true;
0123     if (!algorithm_.empty()) {
0124       if (std::find(algorithm_.begin(), algorithm_.end(), t.algo()) == algorithm_.end())
0125         algo_ok = false;
0126     }
0127     if (!originalAlgorithm_.empty() && algo_ok) {
0128       if (std::find(originalAlgorithm_.begin(), originalAlgorithm_.end(), t.originalAlgo()) == originalAlgorithm_.end())
0129         algo_ok = false;
0130     }
0131     if (!algorithmMask_.empty() && algo_ok) {
0132       if (std::find_if(algorithmMask_.begin(),
0133                        algorithmMask_.end(),
0134                        // for some reason I have to either explicitly give the return type, or use static_cast<bool>()
0135                        [&](reco::TrackBase::TrackAlgorithm algo) -> bool { return t.algoMask()[algo]; }) ==
0136           algorithmMask_.end())
0137         algo_ok = false;
0138     }
0139 
0140     const auto dphi = deltaPhi(t.phi(), meanPhi_);
0141 
0142     auto etaOk = [&](const reco::Track& p) -> bool {
0143       float eta = p.eta();
0144       if (!invertRapidityCut_)
0145         return (eta >= minRapidity_) && (eta <= maxRapidity_);
0146       else
0147         return (eta < minRapidity_ || eta > maxRapidity_);
0148     };
0149 
0150     return ((algo_ok & quality_ok) && t.hitPattern().numberOfValidHits() >= minHit_ &&
0151             t.hitPattern().numberOfValidPixelHits() >= minPixelHit_ &&
0152             t.hitPattern().trackerLayersWithMeasurement() >= minLayer_ &&
0153             t.hitPattern().pixelLayersWithMeasurement() + t.hitPattern().numberOfValidStripLayersWithMonoAndStereo() >=
0154                 min3DLayer_ &&
0155             fabs(t.pt()) >= ptMin_ && etaOk(t) && dphi >= -rangePhi_ && dphi <= rangePhi_ &&
0156             fabs(t.dxy(vertex)) <= tip_ && fabs(t.dsz(vertex)) <= lip_ && t.normalizedChi2() <= maxChi2_);
0157   }
0158 
0159 private:
0160   double ptMin_;
0161   double minRapidity_;
0162   double maxRapidity_;
0163   double meanPhi_;
0164   double rangePhi_;
0165   double tip_;
0166   double lip_;
0167   double maxChi2_;
0168   int minHit_;
0169   int minPixelHit_;
0170   int minLayer_;
0171   int min3DLayer_;
0172   bool usePV_;
0173   bool invertRapidityCut_;
0174 
0175   edm::EDGetTokenT<reco::BeamSpot> bsSrcToken_;
0176   edm::EDGetTokenT<reco::VertexCollection> vertexToken_;
0177 
0178   std::vector<reco::TrackBase::TrackQuality> quality_;
0179   std::vector<reco::TrackBase::TrackAlgorithm> algorithm_;
0180   std::vector<reco::TrackBase::TrackAlgorithm> originalAlgorithm_;
0181   std::vector<reco::TrackBase::TrackAlgorithm> algorithmMask_;
0182 
0183   reco::Track::Point vertex_;
0184 };
0185 
0186 #endif