TkTrackingRegions/plugins/PointSeededTrackingRegionsProducer.h

0001 #ifndef PointSeededTrackingRegionsProducer_h
0002 #define PointSeededTrackingRegionsProducer_h
0003
0004 #include "RecoTracker/TkTrackingRegions/interface/TrackingRegionProducer.h"
0005 #include "RecoTracker/TkTrackingRegions/interface/RectangularEtaPhiTrackingRegion.h"
0006 #include "RecoTracker/MeasurementDet/interface/MeasurementTrackerEvent.h"
0007
0008 #include "FWCore/Framework/interface/Event.h"
0009 #include "FWCore/Framework/interface/EventSetup.h"
0010 #include "FWCore/Utilities/interface/InputTag.h"
0011 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0012 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0013 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0014 #include "FWCore/Framework/interface/ConsumesCollector.h"
0015 #include "DataFormats/Common/interface/Handle.h"
0016 #include "DataFormats/BeamSpot/interface/BeamSpot.h"
0017 #include "DataFormats/VertexReco/interface/Vertex.h"
0018 #include "DataFormats/VertexReco/interface/VertexFwd.h"
0019 #include "DataFormats/Candidate/interface/Candidate.h"
0020 #include "MagneticField/Engine/interface/MagneticField.h"
0021 #include "MagneticField/Records/interface/IdealMagneticFieldRecord.h"
0022 #include "RecoTracker/Record/interface/TrackerMultipleScatteringRecord.h"
0023 #include "RecoTracker/TkMSParametrization/interface/MultipleScatteringParametrisationMaker.h"
0024
0025 #include <TMath.h>
0026 #include <TLorentzVector.h>
0027
0028 /** class PointSeededTrackingRegionsProducer
0029  *
0030  * eta-phi TrackingRegions producer in directions defined by Point-based objects of interest
0031  * from the "input" parameters.
0032  *
0033  * Four operational modes are supported ("mode" parameter):
0034  *
0035  *   BeamSpotFixed:
0036  *     origin is defined by the beam spot
0037  *     z-half-length is defined by a fixed zErrorBeamSpot parameter
0038  *   BeamSpotSigma:
0039  *     origin is defined by the beam spot
0040  *     z-half-length is defined by nSigmaZBeamSpot * beamSpot.sigmaZ
0041  *   VerticesFixed:
0042  *     origins are defined by vertices from VertexCollection (use maximum MaxNVertices of them)
0043  *     z-half-length is defined by a fixed zErrorVetex parameter
0044  *   VerticesSigma:
0045  *     origins are defined by vertices from VertexCollection (use maximum MaxNVertices of them)
0046  *     z-half-length is defined by nSigmaZVertex * vetex.zError
0047  *
0048  *   If, while using one of the "Vertices" modes, there's no vertices in an event, we fall back into
0049  *   either BeamSpotSigma or BeamSpotFixed mode, depending on the positiveness of nSigmaZBeamSpot.
0050  *
0051  */
0052 class PointSeededTrackingRegionsProducer : public TrackingRegionProducer {
0053 public:
0054   typedef enum { BEAM_SPOT_FIXED, BEAM_SPOT_SIGMA, VERTICES_FIXED, VERTICES_SIGMA } Mode;
0055
0056   explicit PointSeededTrackingRegionsProducer(const edm::ParameterSet& conf, edm::ConsumesCollector&& iC)
0057       : token_field(iC.esConsumes()) {
0058     edm::ParameterSet regPSet = conf.getParameter<edm::ParameterSet>("RegionPSet");
0059
0060     // operation mode
0061     std::string modeString = regPSet.getParameter<std::string>("mode");
0062     if (modeString == "BeamSpotFixed")
0063       m_mode = BEAM_SPOT_FIXED;
0064     else if (modeString == "BeamSpotSigma")
0065       m_mode = BEAM_SPOT_SIGMA;
0066     else if (modeString == "VerticesFixed")
0067       m_mode = VERTICES_FIXED;
0068     else if (modeString == "VerticesSigma")
0069       m_mode = VERTICES_SIGMA;
0070     else
0071       edm::LogError("PointSeededTrackingRegionsProducer") << "Unknown mode string: " << modeString;
0072
0073     // basic inputsi
0074     edm::ParameterSet points = regPSet.getParameter<edm::ParameterSet>("points");
0075     etaPoints = points.getParameter<std::vector<double>>("eta");
0076     phiPoints = points.getParameter<std::vector<double>>("phi");
0077     if (!(etaPoints.size() == phiPoints.size()))
0078       throw edm::Exception(edm::errors::Configuration) << "The parameters 'eta' and 'phi' must have the same size";
0079     ;
0080     m_maxNRegions = regPSet.getParameter<int>("maxNRegions");
0081     token_beamSpot = iC.consumes<reco::BeamSpot>(regPSet.getParameter<edm::InputTag>("beamSpot"));
0082     m_maxNVertices = 1;
0083     if (m_mode == VERTICES_FIXED || m_mode == VERTICES_SIGMA) {
0084       token_vertex = iC.consumes<reco::VertexCollection>(regPSet.getParameter<edm::InputTag>("vertexCollection"));
0085       m_maxNVertices = regPSet.getParameter<int>("maxNVertices");
0086     }
0087
0088     // RectangularEtaPhiTrackingRegion parameters:
0089     m_ptMin = regPSet.getParameter<double>("ptMin");
0090     m_originRadius = regPSet.getParameter<double>("originRadius");
0091     m_zErrorBeamSpot = regPSet.getParameter<double>("zErrorBeamSpot");
0092     m_deltaEta = regPSet.getParameter<double>("deltaEta");
0093     m_deltaPhi = regPSet.getParameter<double>("deltaPhi");
0094     m_precise = regPSet.getParameter<bool>("precise");
0095     m_whereToUseMeasurementTracker = RectangularEtaPhiTrackingRegion::stringToUseMeasurementTracker(
0096         regPSet.getParameter<std::string>("whereToUseMeasurementTracker"));
0097     if (m_whereToUseMeasurementTracker != RectangularEtaPhiTrackingRegion::UseMeasurementTracker::kNever) {
0098       token_measurementTracker =
0099           iC.consumes<MeasurementTrackerEvent>(regPSet.getParameter<edm::InputTag>("measurementTrackerName"));
0100     }
0101     m_searchOpt = false;
0102     if (regPSet.exists("searchOpt"))
0103       m_searchOpt = regPSet.getParameter<bool>("searchOpt");
0104
0105     // mode-dependent z-halflength of tracking regions
0106     if (m_mode == VERTICES_SIGMA)
0107       m_nSigmaZVertex = regPSet.getParameter<double>("nSigmaZVertex");
0108     if (m_mode == VERTICES_FIXED)
0109       m_zErrorVetex = regPSet.getParameter<double>("zErrorVetex");
0110     m_nSigmaZBeamSpot = -1.;
0111     if (m_mode == BEAM_SPOT_SIGMA) {
0112       m_nSigmaZBeamSpot = regPSet.getParameter<double>("nSigmaZBeamSpot");
0113       if (m_nSigmaZBeamSpot < 0.)
0114         edm::LogError("PointSeededTrackingRegionsProducer")
0115             << "nSigmaZBeamSpot must be positive for BeamSpotSigma mode!";
0116     }
0117
0118     if (m_precise) {
0119       token_msmaker = iC.esConsumes();
0120     }
0121   }
0122
0123   ~PointSeededTrackingRegionsProducer() override {}
0124
0125   static void fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0126     edm::ParameterSetDescription desc;
0127
0128     edm::ParameterSetDescription descPoints;
0129     descPoints.add<std::vector<double>>("eta", {0.});
0130     descPoints.add<std::vector<double>>("phi", {0.});
0131     desc.add<edm::ParameterSetDescription>("points", descPoints);
0132
0133     desc.add<std::string>("mode", "BeamSpotFixed");
0134     desc.add<int>("maxNRegions", 10);
0135     desc.add<edm::InputTag>("beamSpot", edm::InputTag("hltOnlineBeamSpot"));
0136     desc.add<edm::InputTag>("vertexCollection", edm::InputTag("hltPixelVertices"));
0137     desc.add<int>("maxNVertices", 1);
0138
0139     desc.add<double>("ptMin", 0.9);
0140     desc.add<double>("originRadius", 0.2);
0141     desc.add<double>("zErrorBeamSpot", 24.2);
0142     desc.add<double>("deltaEta", 0.5);
0143     desc.add<double>("deltaPhi", 0.5);
0144     desc.add<bool>("precise", true);
0145
0146     desc.add<double>("nSigmaZVertex", 3.);
0147     desc.add<double>("zErrorVetex", 0.2);
0148     desc.add<double>("nSigmaZBeamSpot", 4.);
0149
0150     desc.add<std::string>("whereToUseMeasurementTracker", "ForSiStrips");
0151     desc.add<edm::InputTag>("measurementTrackerName", edm::InputTag(""));
0152
0153     desc.add<bool>("searchOpt", false);
0154
0155     // Only for backwards-compatibility
0156     edm::ParameterSetDescription descRegion;
0157     descRegion.add<edm::ParameterSetDescription>("RegionPSet", desc);
0158     //edm::ParameterSetDescription descPoint;
0159     //descPoint.add<edm::ParameterSetDescription>("point_input", desc);
0160
0161     descriptions.add("pointSeededTrackingRegion", descRegion);
0162   }
0163
0164   std::vector<std::unique_ptr<TrackingRegion>> regions(const edm::Event& e, const edm::EventSetup& es) const override {
0165     std::vector<std::unique_ptr<TrackingRegion>> result;
0166
0167     // always need the beam spot (as a fall back strategy for vertex modes)
0168     edm::Handle<reco::BeamSpot> bs;
0169     e.getByToken(token_beamSpot, bs);
0170     if (!bs.isValid())
0171       return result;
0172
0173     const auto& field = es.getData(token_field);
0174     const MultipleScatteringParametrisationMaker* msmaker = nullptr;
0175     if (m_precise) {
0176       msmaker = &es.getData(token_msmaker);
0177     }
0178
0179     // this is a default origin for all modes
0180     GlobalPoint default_origin(bs->x0(), bs->y0(), bs->z0());
0181
0182     // vector of origin & halfLength pairs:
0183     std::vector<std::pair<GlobalPoint, float>> origins;
0184
0185     // fill the origins and halfLengths depending on the mode
0186     if (m_mode == BEAM_SPOT_FIXED || m_mode == BEAM_SPOT_SIGMA) {
0187       origins.push_back(std::make_pair(
0188           default_origin, (m_mode == BEAM_SPOT_FIXED) ? m_zErrorBeamSpot : m_nSigmaZBeamSpot * bs->sigmaZ()));
0189     } else if (m_mode == VERTICES_FIXED || m_mode == VERTICES_SIGMA) {
0190       edm::Handle<reco::VertexCollection> vertices;
0191       e.getByToken(token_vertex, vertices);
0192       int n_vert = 0;
0193       for (reco::VertexCollection::const_iterator iv = vertices->begin(), ev = vertices->end();
0194            iv != ev && n_vert < m_maxNVertices;
0195            ++iv) {
0196         if (iv->isFake() || !iv->isValid())
0197           continue;
0198
0199         origins.push_back(std::make_pair(GlobalPoint(iv->x(), iv->y(), iv->z()),
0200                                          (m_mode == VERTICES_FIXED) ? m_zErrorVetex : m_nSigmaZVertex * iv->zError()));
0201         ++n_vert;
0202       }
0203       // no-vertex fall-back case:
0204       if (origins.empty()) {
0205         origins.push_back(std::make_pair(
0206             default_origin, (m_nSigmaZBeamSpot > 0.) ? m_nSigmaZBeamSpot * bs->z0Error() : m_zErrorBeamSpot));
0207       }
0208     }
0209
0210     const MeasurementTrackerEvent* measurementTracker = nullptr;
0211     if (!token_measurementTracker.isUninitialized()) {
0212       edm::Handle<MeasurementTrackerEvent> hmte;
0213       e.getByToken(token_measurementTracker, hmte);
0214       measurementTracker = hmte.product();
0215     }
0216
0217     // create tracking regions (maximum MaxNRegions of them) in directions of the
0218     // points of interest
0219     size_t n_points = etaPoints.size();
0220     int n_regions = 0;
0221     for (size_t i = 0; i < n_points && n_regions < m_maxNRegions; ++i) {
0222       double x = std::cos(phiPoints[i]);
0223       double y = std::sin(phiPoints[i]);
0224       double theta = 2 * std::atan(std::exp(-etaPoints[i]));
0225       double z = 1. / std::tan(theta);
0226
0227       GlobalVector direction(x, y, z);
0228
0229       for (size_t j = 0; j < origins.size() && n_regions < m_maxNRegions; ++j) {
0230         result.push_back(std::make_unique<RectangularEtaPhiTrackingRegion>(direction,         // GlobalVector
0231                                                                            origins[j].first,  // GlobalPoint
0232                                                                            m_ptMin,
0233                                                                            m_originRadius,
0234                                                                            origins[j].second,
0235                                                                            m_deltaEta,
0236                                                                            m_deltaPhi,
0237                                                                            field,
0238                                                                            msmaker,
0239                                                                            m_precise,
0240                                                                            m_whereToUseMeasurementTracker,
0241                                                                            measurementTracker,
0242                                                                            m_searchOpt));
0243         ++n_regions;
0244       }
0245     }
0246     edm::LogInfo("PointSeededTrackingRegionsProducer") << "produced " << n_regions << " regions";
0247
0248     return result;
0249   }
0250
0251 private:
0252   Mode m_mode;
0253
0254   int m_maxNRegions;
0255   edm::EDGetTokenT<reco::VertexCollection> token_vertex;
0256   edm::EDGetTokenT<reco::BeamSpot> token_beamSpot;
0257   int m_maxNVertices;
0258
0259   std::vector<double> etaPoints;
0260   std::vector<double> phiPoints;
0261
0262   float m_ptMin;
0263   float m_originRadius;
0264   float m_zErrorBeamSpot;
0265   float m_deltaEta;
0266   float m_deltaPhi;
0267   bool m_precise;
0268   edm::EDGetTokenT<MeasurementTrackerEvent> token_measurementTracker;
0269   RectangularEtaPhiTrackingRegion::UseMeasurementTracker m_whereToUseMeasurementTracker;
0270   edm::ESGetToken<MagneticField, IdealMagneticFieldRecord> token_field;
0271   edm::ESGetToken<MultipleScatteringParametrisationMaker, TrackerMultipleScatteringRecord> token_msmaker;
0272   bool m_searchOpt;
0273
0274   float m_nSigmaZVertex;
0275   float m_zErrorVetex;
0276   float m_nSigmaZBeamSpot;
0277 };
0278
0279 #endif