File indexing completed on 2021-10-01 22:41:05
0001 #ifndef RecoTracker_TkTrackingRegions_GlobalTrackingRegionWithVerticesProducer_H
0002 #define RecoTracker_TkTrackingRegions_GlobalTrackingRegionWithVerticesProducer_H
0003
0004 #include "DataFormats/Common/interface/DetSetVectorNew.h"
0005 #include "DataFormats/SiPixelCluster/interface/SiPixelCluster.h"
0006 #include "RecoTracker/TkTrackingRegions/interface/TrackingRegionProducer.h"
0007 #include "RecoTracker/TkTrackingRegions/interface/GlobalTrackingRegion.h"
0008 #include "FWCore/Framework/interface/Event.h"
0009 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0010 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0011 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0012 #include "DataFormats/VertexReco/interface/Vertex.h"
0013 #include "DataFormats/VertexReco/interface/VertexFwd.h"
0014 #include "DataFormats/BeamSpot/interface/BeamSpot.h"
0015 #include "FWCore/Utilities/interface/InputTag.h"
0016 #include "FWCore/Framework/interface/ConsumesCollector.h"
0017
0018 class GlobalTrackingRegionWithVerticesProducer : public TrackingRegionProducer {
0019 public:
0020 GlobalTrackingRegionWithVerticesProducer(const edm::ParameterSet& cfg, edm::ConsumesCollector&& iC) {
0021 edm::ParameterSet regionPSet = cfg.getParameter<edm::ParameterSet>("RegionPSet");
0022
0023 thePtMin = regionPSet.getParameter<double>("ptMin");
0024 theOriginRadius = regionPSet.getParameter<double>("originRadius");
0025 theNSigmaZ = regionPSet.getParameter<double>("nSigmaZ");
0026 token_beamSpot = iC.consumes<reco::BeamSpot>(regionPSet.getParameter<edm::InputTag>("beamSpot"));
0027 thePrecise = regionPSet.getParameter<bool>("precise");
0028 theUseMS = regionPSet.getParameter<bool>("useMultipleScattering");
0029
0030 theSigmaZVertex = regionPSet.getParameter<double>("sigmaZVertex");
0031 theFixedError = regionPSet.getParameter<double>("fixedError");
0032
0033 theMaxNVertices = regionPSet.getParameter<int>("maxNVertices");
0034
0035 theUseFoundVertices = regionPSet.getParameter<bool>("useFoundVertices");
0036 theUseFakeVertices = regionPSet.getParameter<bool>("useFakeVertices");
0037 theUseFixedError = regionPSet.getParameter<bool>("useFixedError");
0038 token_vertex = iC.consumes<reco::VertexCollection>(regionPSet.getParameter<edm::InputTag>("VertexCollection"));
0039
0040
0041 theOriginRScaling = regionPSet.getParameter<bool>("originRScaling4BigEvts");
0042 thePtMinScaling = regionPSet.getParameter<bool>("ptMinScaling4BigEvts");
0043 theHalfLengthScaling = regionPSet.getParameter<bool>("halfLengthScaling4BigEvts");
0044 theMinOriginR = regionPSet.getParameter<double>("minOriginR");
0045 theMaxPtMin = regionPSet.getParameter<double>("maxPtMin");
0046 theMinHalfLength = regionPSet.getParameter<double>("minHalfLength");
0047 theScalingStart = regionPSet.getParameter<double>("scalingStartNPix");
0048 theScalingEnd = regionPSet.getParameter<double>("scalingEndNPix");
0049 edm::InputTag pixelClustersForScaling = regionPSet.getParameter<edm::InputTag>("pixelClustersForScaling");
0050 if (theOriginRScaling || thePtMinScaling || theHalfLengthScaling)
0051 token_pc = iC.consumes<edmNew::DetSetVector<SiPixelCluster> >(pixelClustersForScaling);
0052
0053 if (theUseMS) {
0054 token_msmaker = iC.esConsumes();
0055 }
0056 }
0057
0058 ~GlobalTrackingRegionWithVerticesProducer() override {}
0059
0060 static void fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0061 edm::ParameterSetDescription desc;
0062
0063 desc.add<bool>("precise", true);
0064 desc.add<bool>("useMultipleScattering", false);
0065 desc.add<edm::InputTag>("beamSpot", edm::InputTag("offlineBeamSpot"));
0066 desc.add<bool>("useFixedError", true);
0067 desc.add<double>("originRadius", 0.2);
0068 desc.add<double>("sigmaZVertex", 3.0);
0069 desc.add<double>("fixedError", 0.2);
0070 desc.add<edm::InputTag>("VertexCollection", edm::InputTag("firstStepPrimaryVertices"));
0071 desc.add<double>("ptMin", 0.9);
0072 desc.add<bool>("useFoundVertices", true);
0073 desc.add<bool>("useFakeVertices", false);
0074 desc.add<int>("maxNVertices", -1)->setComment("-1 for all vertices");
0075 desc.add<double>("nSigmaZ", 4.0);
0076 desc.add<edm::InputTag>("pixelClustersForScaling", edm::InputTag("siPixelClusters"));
0077 desc.add<bool>("originRScaling4BigEvts", false);
0078 desc.add<bool>("ptMinScaling4BigEvts", false);
0079 desc.add<bool>("halfLengthScaling4BigEvts", false);
0080 desc.add<double>("minOriginR", 0);
0081 desc.add<double>("maxPtMin", 1000);
0082 desc.add<double>("minHalfLength", 0);
0083 desc.add<double>("scalingStartNPix", 0.0);
0084 desc.add<double>("scalingEndNPix", 1.0);
0085
0086
0087 edm::ParameterSetDescription descRegion;
0088 descRegion.add<edm::ParameterSetDescription>("RegionPSet", desc);
0089
0090 descriptions.add("globalTrackingRegionWithVertices", descRegion);
0091 }
0092
0093 std::vector<std::unique_ptr<TrackingRegion> > regions(const edm::Event& ev,
0094 const edm::EventSetup& es) const override {
0095 std::vector<std::unique_ptr<TrackingRegion> > result;
0096
0097 GlobalPoint theOrigin;
0098 edm::Handle<reco::BeamSpot> bsHandle;
0099 ev.getByToken(token_beamSpot, bsHandle);
0100 double bsSigmaZ;
0101 if (bsHandle.isValid()) {
0102 const reco::BeamSpot& bs = *bsHandle;
0103 bsSigmaZ = theNSigmaZ * bs.sigmaZ();
0104 theOrigin = GlobalPoint(bs.x0(), bs.y0(), bs.z0());
0105 } else {
0106 throw cms::Exception("Seeding") << "ERROR: input beamSpot is not valid in GlobalTrackingRegionWithVertices";
0107 }
0108
0109 const MultipleScatteringParametrisationMaker* msmaker = nullptr;
0110 if (theUseMS) {
0111 msmaker = &es.getData(token_msmaker);
0112 }
0113
0114 if (theUseFoundVertices) {
0115 edm::Handle<reco::VertexCollection> vertexCollection;
0116 ev.getByToken(token_vertex, vertexCollection);
0117
0118 edm::Handle<edmNew::DetSetVector<SiPixelCluster> > pixelClusterDSV;
0119 bool doScaling = theOriginRScaling || thePtMinScaling || theHalfLengthScaling;
0120 if (doScaling)
0121 ev.getByToken(token_pc, pixelClusterDSV);
0122
0123 for (reco::VertexCollection::const_iterator iV = vertexCollection->begin(); iV != vertexCollection->end(); iV++) {
0124 if (!iV->isValid())
0125 continue;
0126 if (iV->isFake() && !(theUseFakeVertices && theUseFixedError))
0127 continue;
0128 GlobalPoint theOrigin_ = GlobalPoint(iV->x(), iV->y(), iV->z());
0129
0130
0131 if (doScaling) {
0132
0133 double scaledOriginRadius = theOriginRadius;
0134 double scaledHalfLength = theFixedError;
0135 double scaledPtMin = thePtMin;
0136
0137
0138 double nPix = 0;
0139
0140 const edmNew::DetSetVector<SiPixelCluster>& input = *pixelClusterDSV;
0141 nPix = input.dataSize();
0142
0143
0144 if ((nPix > theScalingEnd) || ((theScalingEnd - theScalingStart) <= 0)) {
0145 if (theOriginRScaling)
0146 scaledOriginRadius = theMinOriginR;
0147 if (theHalfLengthScaling)
0148 scaledHalfLength = theMinHalfLength;
0149 if (thePtMinScaling)
0150 scaledPtMin = theMaxPtMin;
0151 }
0152
0153 else if ((nPix <= theScalingEnd) && (nPix > theScalingStart)) {
0154 float slopeFactor = (nPix - theScalingStart) / (theScalingEnd - theScalingStart);
0155 if (theOriginRScaling)
0156 scaledOriginRadius = theOriginRadius - (theOriginRadius - theMinOriginR) * slopeFactor;
0157 if (theHalfLengthScaling)
0158 scaledHalfLength = theFixedError - (theFixedError - theMinHalfLength) * slopeFactor;
0159 if (thePtMinScaling)
0160 scaledPtMin = thePtMin - (thePtMin - theMaxPtMin) * slopeFactor;
0161 }
0162
0163 if (scaledOriginRadius != 0 && scaledHalfLength != 0) {
0164 result.push_back(std::make_unique<GlobalTrackingRegion>(
0165 scaledPtMin, theOrigin_, scaledOriginRadius, scaledHalfLength, thePrecise, theUseMS, msmaker));
0166 }
0167 }
0168
0169 else {
0170 double theOriginHalfLength_ = (theUseFixedError ? theFixedError : (iV->zError()) * theSigmaZVertex);
0171 result.push_back(std::make_unique<GlobalTrackingRegion>(
0172 thePtMin, theOrigin_, theOriginRadius, theOriginHalfLength_, thePrecise, theUseMS, msmaker));
0173 if (theMaxNVertices >= 0 && result.size() >= static_cast<unsigned>(theMaxNVertices))
0174 break;
0175 }
0176 }
0177
0178 if (result.empty() && !(theOriginRScaling || thePtMinScaling || theHalfLengthScaling)) {
0179 result.push_back(std::make_unique<GlobalTrackingRegion>(
0180 thePtMin, theOrigin, theOriginRadius, bsSigmaZ, thePrecise, theUseMS, msmaker));
0181 }
0182 } else {
0183 result.push_back(std::make_unique<GlobalTrackingRegion>(
0184 thePtMin, theOrigin, theOriginRadius, bsSigmaZ, thePrecise, theUseMS, msmaker));
0185 }
0186
0187 return result;
0188 }
0189
0190 private:
0191 double thePtMin;
0192 double theOriginRadius;
0193 double theNSigmaZ;
0194 edm::InputTag theBeamSpotTag;
0195
0196 double theSigmaZVertex;
0197 double theFixedError;
0198 int theMaxNVertices;
0199 bool thePrecise;
0200 bool theUseMS;
0201
0202 bool theUseFoundVertices;
0203 bool theUseFakeVertices;
0204 bool theUseFixedError;
0205 edm::EDGetTokenT<reco::VertexCollection> token_vertex;
0206 edm::EDGetTokenT<reco::BeamSpot> token_beamSpot;
0207 edm::ESGetToken<MultipleScatteringParametrisationMaker, TrackerMultipleScatteringRecord> token_msmaker;
0208
0209
0210 bool theOriginRScaling;
0211 bool thePtMinScaling;
0212 bool theHalfLengthScaling;
0213 double theMinOriginR;
0214 double theMaxPtMin;
0215 double theMinHalfLength;
0216 double theScalingStart;
0217 double theScalingEnd;
0218 edm::EDGetTokenT<edmNew::DetSetVector<SiPixelCluster> > token_pc;
0219 };
0220
0221 #endif