File indexing completed on 2023-04-15 01:47:33
0001 #ifndef RecoTracker_PixelSeeding_CAHitQuadrupletGenerator_h
0002 #define RecoTracker_PixelSeeding_CAHitQuadrupletGenerator_h
0003
0004 #include "RecoTracker/TkSeedingLayers/interface/SeedComparitorFactory.h"
0005 #include "RecoTracker/TkSeedingLayers/interface/SeedComparitor.h"
0006 #include "RecoTracker/PixelTrackFitting/interface/RZLine.h"
0007 #include "RecoTracker/TkSeedGenerator/interface/FastCircleFit.h"
0008 #include "RecoTracker/TkMSParametrization/interface/PixelRecoUtilities.h"
0009 #include "RecoTracker/TkMSParametrization/interface/LongitudinalBendingCorrection.h"
0010 #include "DataFormats/SiPixelDetId/interface/PixelSubdetector.h"
0011 #include "MagneticField/Records/interface/IdealMagneticFieldRecord.h"
0012
0013 #include "RecoTracker/TkHitPairs/interface/HitPairGeneratorFromLayerPair.h"
0014 #include "RecoTracker/TkHitPairs/interface/LayerHitMapCache.h"
0015 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0016 #include "FWCore/Utilities/interface/EDGetToken.h"
0017
0018 #include "RecoTracker/TkHitPairs/interface/IntermediateHitDoublets.h"
0019 #include "RecoTracker/PixelSeeding/interface/OrderedHitSeeds.h"
0020 #include "RecoTracker/PixelSeeding/interface/CACut.h"
0021
0022 class TrackingRegion;
0023 class SeedingLayerSetsHits;
0024
0025 namespace edm {
0026 class Event;
0027 class EventSetup;
0028 class ParameterSetDescription;
0029 }
0030
0031 class CAHitQuadrupletGenerator {
0032 public:
0033 typedef LayerHitMapCache LayerCacheType;
0034
0035 static constexpr unsigned int minLayers = 4;
0036 typedef OrderedHitSeeds ResultType;
0037
0038 public:
0039 CAHitQuadrupletGenerator(const edm::ParameterSet& cfg, edm::ConsumesCollector&& iC)
0040 : CAHitQuadrupletGenerator(cfg, iC) {}
0041 CAHitQuadrupletGenerator(const edm::ParameterSet& cfg, edm::ConsumesCollector& iC);
0042
0043 ~CAHitQuadrupletGenerator() = default;
0044
0045 static void fillDescriptions(edm::ParameterSetDescription& desc);
0046 static const char* fillDescriptionsLabel() { return "caHitQuadrupletDefault"; }
0047
0048 void initEvent(const edm::Event& ev, const edm::EventSetup& es);
0049
0050 void hitNtuplets(const IntermediateHitDoublets& regionDoublets,
0051 std::vector<OrderedHitSeeds>& result,
0052 const SeedingLayerSetsHits& layers);
0053
0054 private:
0055 LayerCacheType theLayerCache;
0056
0057 std::unique_ptr<SeedComparitor> theComparitor;
0058
0059 const edm::ESGetToken<MagneticField, IdealMagneticFieldRecord> theFieldToken;
0060 const MagneticField* theField = nullptr;
0061
0062 class QuantityDependsPtEval {
0063 public:
0064 QuantityDependsPtEval(float v1, float v2, float c1, float c2)
0065 : value1_(v1), value2_(v2), curvature1_(c1), curvature2_(c2) {}
0066
0067 float value(float curvature) const {
0068 if (value1_ == value2_)
0069 return value1_;
0070
0071 if (curvature1_ < curvature)
0072 return value1_;
0073 if (curvature2_ < curvature && curvature <= curvature1_)
0074 return value2_ + (curvature - curvature2_) / (curvature1_ - curvature2_) * (value1_ - value2_);
0075 return value2_;
0076 }
0077
0078 private:
0079 const float value1_;
0080 const float value2_;
0081 const float curvature1_;
0082 const float curvature2_;
0083 };
0084
0085
0086
0087
0088
0089
0090 class QuantityDependsPt {
0091 public:
0092 explicit QuantityDependsPt(const edm::ParameterSet& pset)
0093 : value1_(pset.getParameter<double>("value1")),
0094 value2_(pset.getParameter<double>("value2")),
0095 pt1_(pset.getParameter<double>("pt1")),
0096 pt2_(pset.getParameter<double>("pt2")),
0097 enabled_(pset.getParameter<bool>("enabled")) {
0098 if (enabled_ && pt1_ >= pt2_)
0099 throw cms::Exception("Configuration") << "PixelQuadrupletGenerator::QuantityDependsPt: pt1 (" << pt1_
0100 << ") needs to be smaller than pt2 (" << pt2_ << ")";
0101 if (pt1_ <= 0)
0102 throw cms::Exception("Configuration")
0103 << "PixelQuadrupletGenerator::QuantityDependsPt: pt1 needs to be > 0; is " << pt1_;
0104 if (pt2_ <= 0)
0105 throw cms::Exception("Configuration")
0106 << "PixelQuadrupletGenerator::QuantityDependsPt: pt2 needs to be > 0; is " << pt2_;
0107 }
0108
0109 QuantityDependsPtEval evaluator(const MagneticField& field) const {
0110 if (enabled_) {
0111 return QuantityDependsPtEval(value1_,
0112 value2_,
0113 PixelRecoUtilities::curvature(1.f / pt1_, field),
0114 PixelRecoUtilities::curvature(1.f / pt2_, field));
0115 }
0116 return QuantityDependsPtEval(value2_, value2_, 0.f, 0.f);
0117 }
0118
0119 private:
0120 const float value1_;
0121 const float value2_;
0122 const float pt1_;
0123 const float pt2_;
0124 const bool enabled_;
0125 };
0126
0127 const float extraHitRPhitolerance;
0128
0129 const QuantityDependsPt maxChi2;
0130 const bool fitFastCircle;
0131 const bool fitFastCircleChi2Cut;
0132 const bool useBendingCorrection;
0133
0134 CACut caThetaCut;
0135 CACut caPhiCut;
0136 const float caHardPtCut = 0.f;
0137 };
0138 #endif