File indexing completed on 2024-09-05 22:39:24
0001
0002
0003
0004 #include <array>
0005 #include <cassert>
0006 #include <functional>
0007 #include <vector>
0008
0009 #include <alpaka/alpaka.hpp>
0010
0011 #include "DataFormats/TrackSoA/interface/TracksDevice.h"
0012 #include "DataFormats/TrackSoA/interface/TracksHost.h"
0013 #include "DataFormats/TrackSoA/interface/alpaka/TracksSoACollection.h"
0014 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0015 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0016 #include "FWCore/Utilities/interface/Exception.h"
0017 #include "HeterogeneousCore/AlpakaInterface/interface/config.h"
0018
0019 #include "CAHitNtupletGenerator.h"
0020 #include "CAHitNtupletGeneratorKernels.h"
0021 #include "CAPixelDoublets.h"
0022 #include "CAPixelDoubletsAlgos.h"
0023
0024 namespace ALPAKA_ACCELERATOR_NAMESPACE {
0025 namespace {
0026
0027 using namespace caHitNtupletGenerator;
0028 using namespace caPixelDoublets;
0029 using namespace pixelTopology;
0030 using namespace pixelTrack;
0031
0032 template <typename T>
0033 T sqr(T x) {
0034 return x * x;
0035 }
0036
0037
0038 void fillDescriptionsCommon(edm::ParameterSetDescription& desc) {
0039
0040
0041
0042 desc.add<double>("ptmin", 0.9f)->setComment("Cut on minimum pt");
0043 desc.add<double>("CAThetaCutBarrel", 0.002f)->setComment("Cut on RZ alignement for Barrel");
0044 desc.add<double>("CAThetaCutForward", 0.003f)->setComment("Cut on RZ alignment for Forward");
0045 desc.add<double>("hardCurvCut", 1.f / (0.35 * 87.f))
0046 ->setComment("Cut on minimum curvature, used in DCA ntuplet selection");
0047 desc.add<double>("dcaCutInnerTriplet", 0.15f)->setComment("Cut on origin radius when the inner hit is on BPix1");
0048 desc.add<double>("dcaCutOuterTriplet", 0.25f)->setComment("Cut on origin radius when the outer hit is on BPix1");
0049 desc.add<bool>("earlyFishbone", true);
0050 desc.add<bool>("lateFishbone", false);
0051 desc.add<bool>("fillStatistics", false);
0052 desc.add<unsigned int>("minHitsPerNtuplet", 4);
0053 desc.add<unsigned int>("minHitsForSharingCut", 10)
0054 ->setComment("Maximum number of hits in a tuple to clean also if the shared hit is on bpx1");
0055
0056 desc.add<bool>("fitNas4", false)->setComment("fit only 4 hits out of N");
0057 desc.add<bool>("doClusterCut", true);
0058 desc.add<bool>("doZ0Cut", true);
0059 desc.add<bool>("doPtCut", true);
0060 desc.add<bool>("useRiemannFit", false)->setComment("true for Riemann, false for BrokenLine");
0061 desc.add<bool>("doSharedHitCut", true)->setComment("Sharing hit nTuples cleaning");
0062 desc.add<bool>("dupPassThrough", false)->setComment("Do not reject duplicate");
0063 desc.add<bool>("useSimpleTripletCleaner", true)->setComment("use alternate implementation");
0064 }
0065
0066 AlgoParams makeCommonParams(edm::ParameterSet const& cfg) {
0067 return AlgoParams({cfg.getParameter<unsigned int>("minHitsForSharingCut"),
0068 cfg.getParameter<bool>("useRiemannFit"),
0069 cfg.getParameter<bool>("fitNas4"),
0070 cfg.getParameter<bool>("includeJumpingForwardDoublets"),
0071 cfg.getParameter<bool>("earlyFishbone"),
0072 cfg.getParameter<bool>("lateFishbone"),
0073 cfg.getParameter<bool>("fillStatistics"),
0074 cfg.getParameter<bool>("doSharedHitCut"),
0075 cfg.getParameter<bool>("dupPassThrough"),
0076 cfg.getParameter<bool>("useSimpleTripletCleaner")});
0077 }
0078
0079
0080 template <typename TrackerTraits, typename Enable = void>
0081 struct TopologyCuts {};
0082
0083 template <typename TrackerTraits>
0084 struct TopologyCuts<TrackerTraits, isPhase1Topology<TrackerTraits>> {
0085 static constexpr CAParamsT<TrackerTraits> makeCACuts(edm::ParameterSet const& cfg) {
0086 return CAParamsT<TrackerTraits>{{cfg.getParameter<unsigned int>("maxNumberOfDoublets"),
0087 cfg.getParameter<unsigned int>("minHitsPerNtuplet"),
0088 (float)cfg.getParameter<double>("ptmin"),
0089 (float)cfg.getParameter<double>("CAThetaCutBarrel"),
0090 (float)cfg.getParameter<double>("CAThetaCutForward"),
0091 (float)cfg.getParameter<double>("hardCurvCut"),
0092 (float)cfg.getParameter<double>("dcaCutInnerTriplet"),
0093 (float)cfg.getParameter<double>("dcaCutOuterTriplet")}};
0094 };
0095
0096 static constexpr ::pixelTrack::QualityCutsT<TrackerTraits> makeQualityCuts(edm::ParameterSet const& pset) {
0097 auto coeff = pset.getParameter<std::array<double, 2>>("chi2Coeff");
0098 auto ptMax = pset.getParameter<double>("chi2MaxPt");
0099
0100 coeff[1] = (coeff[1] - coeff[0]) / log2(ptMax);
0101 return ::pixelTrack::QualityCutsT<TrackerTraits>{
0102 {(float)coeff[0], (float)coeff[1], 0.f, 0.f},
0103
0104 (float)ptMax,
0105
0106 (float)pset.getParameter<double>("chi2Scale"),
0107
0108 {(float)pset.getParameter<double>("tripletMaxTip"),
0109 (float)pset.getParameter<double>("tripletMinPt"),
0110 (float)pset.getParameter<double>("tripletMaxZip")},
0111
0112 {(float)pset.getParameter<double>("quadrupletMaxTip"),
0113 (float)pset.getParameter<double>("quadrupletMinPt"),
0114 (float)pset.getParameter<double>("quadrupletMaxZip")}};
0115 }
0116 };
0117
0118 template <typename TrackerTraits>
0119 struct TopologyCuts<TrackerTraits, isPhase2Topology<TrackerTraits>> {
0120 static constexpr CAParamsT<TrackerTraits> makeCACuts(edm::ParameterSet const& cfg) {
0121 return CAParamsT<TrackerTraits>{{cfg.getParameter<unsigned int>("maxNumberOfDoublets"),
0122 cfg.getParameter<unsigned int>("minHitsPerNtuplet"),
0123 (float)cfg.getParameter<double>("ptmin"),
0124 (float)cfg.getParameter<double>("CAThetaCutBarrel"),
0125 (float)cfg.getParameter<double>("CAThetaCutForward"),
0126 (float)cfg.getParameter<double>("hardCurvCut"),
0127 (float)cfg.getParameter<double>("dcaCutInnerTriplet"),
0128 (float)cfg.getParameter<double>("dcaCutOuterTriplet")},
0129 {(bool)cfg.getParameter<bool>("includeFarForwards")}};
0130 }
0131
0132 static constexpr ::pixelTrack::QualityCutsT<TrackerTraits> makeQualityCuts(edm::ParameterSet const& pset) {
0133 return ::pixelTrack::QualityCutsT<TrackerTraits>{
0134 static_cast<float>(pset.getParameter<double>("maxChi2")),
0135 static_cast<float>(pset.getParameter<double>("minPt")),
0136 static_cast<float>(pset.getParameter<double>("maxTip")),
0137 static_cast<float>(pset.getParameter<double>("maxZip")),
0138 };
0139 }
0140 };
0141
0142
0143
0144
0145 template <typename TrackerTraits>
0146 CellCutsT<TrackerTraits> makeCellCuts(edm::ParameterSet const& cfg) {
0147 return CellCutsT<TrackerTraits>{cfg.getParameter<bool>("doClusterCut"),
0148 cfg.getParameter<bool>("doZ0Cut"),
0149 cfg.getParameter<bool>("doPtCut"),
0150 cfg.getParameter<bool>("idealConditions"),
0151 (float)cfg.getParameter<double>("cellZ0Cut"),
0152 (float)cfg.getParameter<double>("cellPtCut"),
0153 cfg.getParameter<std::vector<int>>("phiCuts")};
0154 }
0155
0156 }
0157
0158 using namespace std;
0159
0160 template <typename TrackerTraits>
0161 CAHitNtupletGenerator<TrackerTraits>::CAHitNtupletGenerator(const edm::ParameterSet& cfg)
0162 : m_params(makeCommonParams(cfg),
0163 makeCellCuts<TrackerTraits>(cfg),
0164 TopologyCuts<TrackerTraits>::makeQualityCuts(cfg.getParameterSet("trackQualityCuts")),
0165 TopologyCuts<TrackerTraits>::makeCACuts(cfg)) {
0166 #ifdef DUMP_GPU_TK_TUPLES
0167 printf("TK: %s %s %s %s %s %s %s %s %s %s %s %s %s %s %s %s %s\n",
0168 "tid",
0169 "qual",
0170 "nh",
0171 "nl",
0172 "charge",
0173 "pt",
0174 "eta",
0175 "phi",
0176 "tip",
0177 "zip",
0178 "chi2",
0179 "h1",
0180 "h2",
0181 "h3",
0182 "h4",
0183 "h5",
0184 "hn");
0185 #endif
0186 }
0187
0188 template <typename TrackerTraits>
0189 void CAHitNtupletGenerator<TrackerTraits>::fillPSetDescription(edm::ParameterSetDescription& desc) {
0190 static_assert(sizeof(TrackerTraits) == 0,
0191 "Note: this fillPSetDescription is a dummy one. Please specialise it for the correct version of "
0192 "CAHitNtupletGenerator<TrackerTraits>.");
0193 }
0194
0195 template <>
0196 void CAHitNtupletGenerator<pixelTopology::Phase1>::fillPSetDescription(edm::ParameterSetDescription& desc) {
0197 fillDescriptionsCommon(desc);
0198
0199 desc.add<unsigned int>("maxNumberOfDoublets", pixelTopology::Phase1::maxNumberOfDoublets);
0200 desc.add<bool>("idealConditions", true);
0201 desc.add<bool>("includeJumpingForwardDoublets", false);
0202 desc.add<double>("cellZ0Cut", 12.0);
0203 desc.add<double>("cellPtCut", 0.5);
0204
0205 edm::ParameterSetDescription trackQualityCuts;
0206 trackQualityCuts.add<double>("chi2MaxPt", 10.)->setComment("max pT used to determine the pT-dependent chi2 cut");
0207 trackQualityCuts.add<std::vector<double>>("chi2Coeff", {0.9, 1.8})->setComment("chi2 at 1GeV and at ptMax above");
0208 trackQualityCuts.add<double>("chi2Scale", 8.)
0209 ->setComment(
0210 "Factor to multiply the pT-dependent chi2 cut (currently: 8 for the broken line fit, ?? for the Riemann "
0211 "fit)");
0212 trackQualityCuts.add<double>("tripletMinPt", 0.5)->setComment("Min pT for triplets, in GeV");
0213 trackQualityCuts.add<double>("tripletMaxTip", 0.3)->setComment("Max |Tip| for triplets, in cm");
0214 trackQualityCuts.add<double>("tripletMaxZip", 12.)->setComment("Max |Zip| for triplets, in cm");
0215 trackQualityCuts.add<double>("quadrupletMinPt", 0.3)->setComment("Min pT for quadruplets, in GeV");
0216 trackQualityCuts.add<double>("quadrupletMaxTip", 0.5)->setComment("Max |Tip| for quadruplets, in cm");
0217 trackQualityCuts.add<double>("quadrupletMaxZip", 12.)->setComment("Max |Zip| for quadruplets, in cm");
0218 desc.add<edm::ParameterSetDescription>("trackQualityCuts", trackQualityCuts)
0219 ->setComment(
0220 "Quality cuts based on the results of the track fit:\n - apply a pT-dependent chi2 cut;\n - apply "
0221 "\"region "
0222 "cuts\" based on the fit results (pT, Tip, Zip).");
0223
0224 desc.add<std::vector<int>>(
0225 "phiCuts",
0226 std::vector<int>(std::begin(phase1PixelTopology::phicuts), std::end(phase1PixelTopology::phicuts)))
0227 ->setComment("Cuts in phi for cells");
0228 }
0229
0230 template <>
0231 void CAHitNtupletGenerator<pixelTopology::HIonPhase1>::fillPSetDescription(edm::ParameterSetDescription& desc) {
0232 fillDescriptionsCommon(desc);
0233
0234 desc.add<unsigned int>("maxNumberOfDoublets", pixelTopology::HIonPhase1::maxNumberOfDoublets);
0235 desc.add<bool>("idealConditions", false);
0236 desc.add<bool>("includeJumpingForwardDoublets", false);
0237 desc.add<double>("cellZ0Cut", 10.0);
0238 desc.add<double>("cellPtCut", 0.0);
0239
0240 edm::ParameterSetDescription trackQualityCuts;
0241 trackQualityCuts.add<double>("chi2MaxPt", 10.)->setComment("max pT used to determine the pT-dependent chi2 cut");
0242 trackQualityCuts.add<std::vector<double>>("chi2Coeff", {0.9, 1.8})->setComment("chi2 at 1GeV and at ptMax above");
0243 trackQualityCuts.add<double>("chi2Scale", 8.)
0244 ->setComment(
0245 "Factor to multiply the pT-dependent chi2 cut (currently: 8 for the broken line fit, ?? for the Riemann "
0246 "fit)");
0247 trackQualityCuts.add<double>("tripletMinPt", 0.0)->setComment("Min pT for triplets, in GeV");
0248 trackQualityCuts.add<double>("tripletMaxTip", 0.1)->setComment("Max |Tip| for triplets, in cm");
0249 trackQualityCuts.add<double>("tripletMaxZip", 6.)->setComment("Max |Zip| for triplets, in cm");
0250 trackQualityCuts.add<double>("quadrupletMinPt", 0.0)->setComment("Min pT for quadruplets, in GeV");
0251 trackQualityCuts.add<double>("quadrupletMaxTip", 0.5)->setComment("Max |Tip| for quadruplets, in cm");
0252 trackQualityCuts.add<double>("quadrupletMaxZip", 6.)->setComment("Max |Zip| for quadruplets, in cm");
0253
0254 desc.add<edm::ParameterSetDescription>("trackQualityCuts", trackQualityCuts)
0255 ->setComment(
0256 "Quality cuts based on the results of the track fit:\n - apply a pT-dependent chi2 cut;\n - apply "
0257 "\"region "
0258 "cuts\" based on the fit results (pT, Tip, Zip).");
0259
0260 desc.add<std::vector<int>>(
0261 "phiCuts",
0262 std::vector<int>(std::begin(phase1PixelTopology::phicuts), std::end(phase1PixelTopology::phicuts)))
0263 ->setComment("Cuts in phi for cells");
0264 }
0265
0266 template <>
0267 void CAHitNtupletGenerator<pixelTopology::Phase2>::fillPSetDescription(edm::ParameterSetDescription& desc) {
0268 fillDescriptionsCommon(desc);
0269
0270 desc.add<unsigned int>("maxNumberOfDoublets", pixelTopology::Phase2::maxNumberOfDoublets);
0271 desc.add<bool>("idealConditions", false);
0272 desc.add<bool>("includeFarForwards", true);
0273 desc.add<bool>("includeJumpingForwardDoublets", true);
0274 desc.add<double>("cellZ0Cut", 7.5);
0275 desc.add<double>("cellPtCut", 0.85);
0276
0277 edm::ParameterSetDescription trackQualityCuts;
0278 trackQualityCuts.add<double>("maxChi2", 5.)->setComment("Max normalized chi2");
0279 trackQualityCuts.add<double>("minPt", 0.5)->setComment("Min pT in GeV");
0280 trackQualityCuts.add<double>("maxTip", 0.3)->setComment("Max |Tip| in cm");
0281 trackQualityCuts.add<double>("maxZip", 12.)->setComment("Max |Zip|, in cm");
0282 desc.add<edm::ParameterSetDescription>("trackQualityCuts", trackQualityCuts)
0283 ->setComment(
0284 "Quality cuts based on the results of the track fit:\n - apply cuts based on the fit results (pT, Tip, "
0285 "Zip).");
0286
0287 desc.add<std::vector<int>>(
0288 "phiCuts",
0289 std::vector<int>(std::begin(phase2PixelTopology::phicuts), std::end(phase2PixelTopology::phicuts)))
0290 ->setComment("Cuts in phi for cells");
0291 }
0292
0293 template <typename TrackerTraits>
0294 TracksSoACollection<TrackerTraits> CAHitNtupletGenerator<TrackerTraits>::makeTuplesAsync(
0295 HitsOnDevice const& hits_d, ParamsOnDevice const* cpeParams, float bfield, Queue& queue) const {
0296 using HelixFit = HelixFit<TrackerTraits>;
0297 using TrackSoA = TracksSoACollection<TrackerTraits>;
0298 using GPUKernels = CAHitNtupletGeneratorKernels<TrackerTraits>;
0299
0300 TrackSoA tracks(queue);
0301
0302
0303 if (hits_d.view().metadata().size() < 2) {
0304 const auto device = alpaka::getDev(queue);
0305 auto ntracks_d = cms::alpakatools::make_device_view(device, tracks.view().nTracks());
0306 alpaka::memset(queue, ntracks_d, 0);
0307 return tracks;
0308 }
0309 GPUKernels kernels(m_params, hits_d.view().metadata().size(), hits_d.offsetBPIX2(), queue);
0310
0311 kernels.buildDoublets(hits_d.view(), hits_d.offsetBPIX2(), queue);
0312 kernels.launchKernels(hits_d.view(), hits_d.offsetBPIX2(), tracks.view(), queue);
0313
0314 HelixFit fitter(bfield, m_params.fitNas4_);
0315 fitter.allocate(kernels.tupleMultiplicity(), tracks.view());
0316 if (m_params.useRiemannFit_) {
0317 fitter.launchRiemannKernels(
0318 hits_d.view(), cpeParams, hits_d.view().metadata().size(), TrackerTraits::maxNumberOfQuadruplets, queue);
0319 } else {
0320 fitter.launchBrokenLineKernels(
0321 hits_d.view(), cpeParams, hits_d.view().metadata().size(), TrackerTraits::maxNumberOfQuadruplets, queue);
0322 }
0323 kernels.classifyTuples(hits_d.view(), tracks.view(), queue);
0324 #ifdef GPU_DEBUG
0325 alpaka::wait(queue);
0326 std::cout << "finished building pixel tracks on GPU" << std::endl;
0327 #endif
0328
0329 return tracks;
0330 }
0331
0332 template class CAHitNtupletGenerator<pixelTopology::Phase1>;
0333 template class CAHitNtupletGenerator<pixelTopology::Phase2>;
0334 template class CAHitNtupletGenerator<pixelTopology::HIonPhase1>;
0335 }