Go to the documentation of this file.00001 #ifndef RecoSelectors_RecoTrackSelector_h
00002 #define RecoSelectors_RecoTrackSelector_h
00003
00004
00005
00006
00007
00008
00009
00010
00011 #include "DataFormats/TrackReco/interface/Track.h"
00012 #include "DataFormats/TrackReco/interface/TrackFwd.h"
00013 #include "DataFormats/BeamSpot/interface/BeamSpot.h"
00014 #include "FWCore/Utilities/interface/InputTag.h"
00015
00016 class RecoTrackSelector {
00017 public:
00018 typedef reco::TrackCollection collection;
00019 typedef std::vector<const reco::Track *> container;
00020 typedef container::const_iterator const_iterator;
00021
00023 RecoTrackSelector() {}
00024 RecoTrackSelector ( const edm::ParameterSet & cfg ) :
00025 ptMin_(cfg.getParameter<double>("ptMin")),
00026 minRapidity_(cfg.getParameter<double>("minRapidity")),
00027 maxRapidity_(cfg.getParameter<double>("maxRapidity")),
00028 tip_(cfg.getParameter<double>("tip")),
00029 lip_(cfg.getParameter<double>("lip")),
00030 minHit_(cfg.getParameter<int>("minHit")),
00031 min3DHit_(cfg.getParameter<int>("min3DHit")),
00032 maxChi2_(cfg.getParameter<double>("maxChi2")),
00033 bsSrc_(cfg.getParameter<edm::InputTag>("beamSpot")),
00034 bs(0)
00035 {
00036 std::vector<std::string> quality = cfg.getParameter<std::vector<std::string> >("quality");
00037 for (unsigned int j=0;j<quality.size();j++) quality_.push_back(reco::TrackBase::qualityByName(quality[j]));
00038 std::vector<std::string> algorithm = cfg.getParameter<std::vector<std::string> >("algorithm");
00039 for (unsigned int j=0;j<algorithm.size();j++) algorithm_.push_back(reco::TrackBase::algoByName(algorithm[j]));
00040 }
00041
00042 RecoTrackSelector ( double ptMin, double minRapidity, double maxRapidity,
00043 double tip, double lip, int minHit, int min3DHit, double maxChi2,
00044 std::vector<std::string> quality , std::vector<std::string> algorithm ) :
00045 ptMin_( ptMin ), minRapidity_( minRapidity ), maxRapidity_( maxRapidity ),
00046 tip_( tip ), lip_( lip ), minHit_( minHit ), min3DHit_( min3DHit), maxChi2_( maxChi2 ), bs(0)
00047 {
00048 for (unsigned int j=0;j<quality.size();j++) quality_.push_back(reco::TrackBase::qualityByName(quality[j]));
00049 for (unsigned int j=0;j<algorithm.size();j++) algorithm_.push_back(reco::TrackBase::algoByName(algorithm[j]));
00050 }
00051
00052 const_iterator begin() const { return selected_.begin(); }
00053 const_iterator end() const { return selected_.end(); }
00054
00055 void select( const edm::Handle<collection>& c, const edm::Event & event, const edm::EventSetup&) {
00056 selected_.clear();
00057 edm::Handle<reco::BeamSpot> beamSpot;
00058 event.getByLabel(bsSrc_,beamSpot);
00059 bs = beamSpot.product();
00060 for( reco::TrackCollection::const_iterator trk = c->begin();
00061 trk != c->end(); ++ trk )
00062 if ( operator()(*trk) ) {
00063 selected_.push_back( & * trk );
00064 }
00065 }
00066
00068 bool operator()( const reco::Track & t, edm::Event& event) {
00069
00070 if ((bs==0)|| (previousEvent != event.id())) {
00071 edm::Handle<reco::BeamSpot> beamSpot;
00072 event.getByLabel(bsSrc_,beamSpot);
00073 bs = beamSpot.product();
00074 previousEvent = event.id();
00075 }
00076 return operator()(t);
00077 }
00078
00079 bool operator()( const reco::Track & t, const reco::BeamSpot* bs_) {
00080 bs = bs_;
00081 return operator()(t);
00082 }
00083
00084 bool operator()( const reco::Track & t) {
00085 bool quality_ok = true;
00086 if (quality_.size()!=0) {
00087 quality_ok = false;
00088 for (unsigned int i = 0; i<quality_.size();++i) {
00089 if (t.quality(quality_[i])){
00090 quality_ok = true;
00091 break;
00092 }
00093 }
00094 }
00095 bool algo_ok = true;
00096 if (algorithm_.size()!=0) {
00097 if (std::find(algorithm_.begin(),algorithm_.end(),t.algo())==algorithm_.end()) algo_ok = false;
00098 }
00099 return
00100 (t.hitPattern().trackerLayersWithMeasurement() >= minHit_ &&
00101 t.hitPattern().pixelLayersWithMeasurement() +
00102 t.hitPattern().numberOfValidStripLayersWithMonoAndStereo() >= min3DHit_ &&
00103 fabs(t.pt()) >= ptMin_ &&
00104 t.eta() >= minRapidity_ && t.eta() <= maxRapidity_ &&
00105 fabs(t.dxy(bs->position())) <= tip_ &&
00106 fabs(t.dsz(bs->position())) <= lip_ &&
00107 t.normalizedChi2()<=maxChi2_ &&
00108 quality_ok &&
00109 algo_ok);
00110 }
00111
00112 size_t size() const { return selected_.size(); }
00113
00114 protected:
00115 double ptMin_;
00116 double minRapidity_;
00117 double maxRapidity_;
00118 double tip_;
00119 double lip_;
00120 int minHit_;
00121 int min3DHit_;
00122 double maxChi2_;
00123 std::vector<reco::TrackBase::TrackQuality> quality_;
00124 std::vector<reco::TrackBase::TrackAlgorithm> algorithm_;
00125 edm::InputTag bsSrc_;
00126 const reco::BeamSpot* bs;
00127 container selected_;
00128 edm::EventID previousEvent;
00129 };
00130
00131 #endif