Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-02-14 13:34:46

0001 #ifndef HIProtoTrackSelection_h
0002 #define HIProtoTrackSelection_h
0003 
0004 #include "FWCore/Framework/interface/ConsumesCollector.h"
0005 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0006 #include "FWCore/Framework/interface/Event.h"
0007 
0008 #include "DataFormats/TrackReco/interface/TrackFwd.h"
0009 #include "DataFormats/TrackReco/interface/Track.h"
0010 
0011 #include "DataFormats/VertexReco/interface/Vertex.h"
0012 #include "DataFormats/VertexReco/interface/VertexFwd.h"
0013 
0014 #include "DataFormats/BeamSpot/interface/BeamSpot.h"
0015 #include "FWCore/Utilities/interface/InputTag.h"
0016 
0017 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0018 
0019 #include <algorithm>
0020 #include <iostream>
0021 
0022 /**
0023  Selector to select prototracks that pass certain kinematic cuts based on fast vertex
0024 **/
0025 
0026 class HIProtoTrackSelector {
0027 public:
0028   // input collection type
0029   typedef reco::TrackCollection collection;
0030 
0031   // output collection type
0032   typedef std::vector<const reco::Track*> container;
0033 
0034   // iterator over result collection type.
0035   typedef container::const_iterator const_iterator;
0036 
0037   // constructor from parameter set configurability
0038   HIProtoTrackSelector(const edm::ParameterSet& iConfig, edm::ConsumesCollector&& iC)
0039       : vertexCollection_(iConfig.getParameter<edm::InputTag>("VertexCollection")),
0040         vertexCollectionToken_(iC.consumes<reco::VertexCollection>(vertexCollection_)),
0041         beamSpot_(iConfig.getParameter<edm::InputTag>("beamSpotLabel")),
0042         beamSpotToken_(iC.consumes<reco::BeamSpot>(beamSpot_)),
0043         ptMin_(iConfig.getParameter<double>("ptMin")),
0044         nSigmaZ_(iConfig.getParameter<double>("nSigmaZ")),
0045         minZCut_(iConfig.getParameter<double>("minZCut")),
0046         maxD0Significance_(iConfig.getParameter<double>("maxD0Significance")){};
0047 
0048   // select object from a collection and possibly event content
0049   void select(edm::Handle<reco::TrackCollection>& TCH, const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0050     selected_.clear();
0051 
0052     const collection& c = *(TCH.product());
0053 
0054     // Get fast reco vertex
0055     edm::Handle<reco::VertexCollection> vc;
0056     iEvent.getByToken(vertexCollectionToken_, vc);
0057     const reco::VertexCollection* vertices = vc.product();
0058 
0059     math::XYZPoint vtxPoint(0.0, 0.0, 0.0);
0060     double vzErr = 0.0;
0061 
0062     if (!vertices->empty()) {
0063       vtxPoint = vertices->begin()->position();
0064       vzErr = vertices->begin()->zError();
0065       edm::LogInfo("HeavyIonVertexing") << "Select prototracks compatible with median vertex"
0066                                         << "\n   vz = " << vtxPoint.Z() << "\n   " << nSigmaZ_
0067                                         << " vz sigmas = " << vzErr * nSigmaZ_
0068                                         << "\n   cut at = " << std::max(vzErr * nSigmaZ_, minZCut_);
0069     }
0070     // Supress this warning, since events w/ no vertex are expected
0071     //else {
0072     //edm::LogError("HeavyIonVertexing") << "No vertex found in collection '" << vertexCollection_ << "'";
0073     //}
0074 
0075     // Get beamspot
0076     reco::BeamSpot beamSpot;
0077     edm::Handle<reco::BeamSpot> beamSpotHandle;
0078     iEvent.getByToken(beamSpotToken_, beamSpotHandle);
0079 
0080     math::XYZPoint bsPoint(0.0, 0.0, 0.0);
0081     double bsWidth = 0.0;
0082 
0083     if (beamSpotHandle.isValid()) {
0084       beamSpot = *beamSpotHandle;
0085       bsPoint = beamSpot.position();
0086       bsWidth = sqrt(beamSpot.BeamWidthX() * beamSpot.BeamWidthY());
0087       edm::LogInfo("HeavyIonVertexing") << "Select prototracks compatible with beamspot"
0088                                         << "\n   (x,y,z) = (" << bsPoint.X() << "," << bsPoint.Y() << "," << bsPoint.Z()
0089                                         << ")"
0090                                         << "\n   width = " << bsWidth
0091                                         << "\n   cut at d0/d0sigma = " << maxD0Significance_;
0092     } else {
0093       edm::LogError("HeavyIonVertexing") << "No beam spot available from '" << beamSpot_ << "\n";
0094     }
0095 
0096     // Do selection
0097     int nSelected = 0;
0098     int nRejected = 0;
0099     double d0 = 0.0;
0100     double d0sigma = 0.0;
0101     for (reco::TrackCollection::const_iterator trk = c.begin(); trk != c.end(); ++trk) {
0102       d0 = -1. * trk->dxy(bsPoint);
0103       d0sigma = sqrt(trk->d0Error() * trk->d0Error() + bsWidth * bsWidth);
0104       if (trk->pt() > ptMin_                          // keep only tracks above ptMin
0105           && fabs(d0 / d0sigma) < maxD0Significance_  // keep only tracks with D0 significance less than cut
0106           && fabs(trk->dz(vtxPoint)) < std::max(vzErr * nSigmaZ_, minZCut_)  // within minZCut, nSigmaZ of fast vertex
0107       ) {
0108         nSelected++;
0109         selected_.push_back(&*trk);
0110       } else
0111         nRejected++;
0112     }
0113 
0114     edm::LogInfo("HeavyIonVertexing") << "selected " << nSelected << " prototracks out of " << nRejected + nSelected
0115                                       << "\n";
0116   }
0117 
0118   // iterators over selected objects: collection begin
0119   const_iterator begin() const { return selected_.begin(); }
0120 
0121   // iterators over selected objects: collection end
0122   const_iterator end() const { return selected_.end(); }
0123 
0124   // true if no object has been selected
0125   size_t size() const { return selected_.size(); }
0126 
0127 private:
0128   container selected_;
0129   edm::InputTag vertexCollection_;
0130   edm::EDGetTokenT<reco::VertexCollection> vertexCollectionToken_;
0131   edm::InputTag beamSpot_;
0132   edm::EDGetTokenT<reco::BeamSpot> beamSpotToken_;
0133   double ptMin_;
0134   double nSigmaZ_;
0135   double minZCut_;
0136   double maxD0Significance_;
0137 };
0138 
0139 #endif