File indexing completed on 2023-07-28 01:42:34
0001
0002
0003
0004
0005
0006
0007
0008 #include <array>
0009 #include <cassert>
0010 #include <functional>
0011 #include <vector>
0012
0013 #include "CUDADataFormats/Track/interface/TrackSoAHeterogeneousDevice.h"
0014 #include "CUDADataFormats/Track/interface/TrackSoAHeterogeneousHost.h"
0015 #include "CUDADataFormats/TrackingRecHit/interface/TrackingRecHitSoADevice.h"
0016 #include "CUDADataFormats/TrackingRecHit/interface/TrackingRecHitSoAHost.h"
0017 #include "DataFormats/Common/interface/Handle.h"
0018 #include "FWCore/Framework/interface/ConsumesCollector.h"
0019 #include "FWCore/Framework/interface/Event.h"
0020 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0021 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0022 #include "FWCore/ServiceRegistry/interface/Service.h"
0023 #include "FWCore/Utilities/interface/EDMException.h"
0024 #include "FWCore/Utilities/interface/isFinite.h"
0025 #include "HeterogeneousCore/CUDAServices/interface/CUDAInterface.h"
0026 #include "TrackingTools/DetLayers/interface/BarrelDetLayer.h"
0027
0028 #include "CAHitNtupletGeneratorOnGPU.h"
0029
0030 namespace {
0031
0032 using namespace caHitNtupletGenerator;
0033 using namespace gpuPixelDoublets;
0034 using namespace pixelTopology;
0035 using namespace pixelTrack;
0036
0037 template <typename T>
0038 T sqr(T x) {
0039 return x * x;
0040 }
0041
0042
0043 AlgoParams makeCommonParams(edm::ParameterSet const& cfg) {
0044 return AlgoParams({cfg.getParameter<bool>("onGPU"),
0045 cfg.getParameter<unsigned int>("minHitsForSharingCut"),
0046 cfg.getParameter<bool>("useRiemannFit"),
0047 cfg.getParameter<bool>("fitNas4"),
0048 cfg.getParameter<bool>("includeJumpingForwardDoublets"),
0049 cfg.getParameter<bool>("earlyFishbone"),
0050 cfg.getParameter<bool>("lateFishbone"),
0051 cfg.getParameter<bool>("fillStatistics"),
0052 cfg.getParameter<bool>("doSharedHitCut"),
0053 cfg.getParameter<bool>("dupPassThrough"),
0054 cfg.getParameter<bool>("useSimpleTripletCleaner")});
0055 }
0056
0057
0058 template <typename TrackerTraits, typename Enable = void>
0059 struct topologyCuts {};
0060
0061 template <typename TrackerTraits>
0062 struct topologyCuts<TrackerTraits, isPhase1Topology<TrackerTraits>> {
0063 static constexpr CAParamsT<TrackerTraits> makeCACuts(edm::ParameterSet const& cfg) {
0064 return CAParamsT<TrackerTraits>{{
0065 cfg.getParameter<unsigned int>("maxNumberOfDoublets"),
0066 cfg.getParameter<unsigned int>("minHitsPerNtuplet"),
0067 static_cast<float>(cfg.getParameter<double>("ptmin")),
0068 static_cast<float>(cfg.getParameter<double>("CAThetaCutBarrel")),
0069 static_cast<float>(cfg.getParameter<double>("CAThetaCutForward")),
0070 static_cast<float>(cfg.getParameter<double>("hardCurvCut")),
0071 static_cast<float>(cfg.getParameter<double>("dcaCutInnerTriplet")),
0072 static_cast<float>(cfg.getParameter<double>("dcaCutOuterTriplet")),
0073 }};
0074 };
0075
0076 static constexpr pixelTrack::QualityCutsT<TrackerTraits> makeQualityCuts(edm::ParameterSet const& pset) {
0077 auto coeff = pset.getParameter<std::array<double, 2>>("chi2Coeff");
0078 auto ptMax = pset.getParameter<double>("chi2MaxPt");
0079
0080 coeff[1] = (coeff[1] - coeff[0]) / log2(ptMax);
0081 return pixelTrack::QualityCutsT<TrackerTraits>{
0082 {(float)coeff[0], (float)coeff[1], 0.f, 0.f},
0083
0084 (float)ptMax,
0085
0086 (float)pset.getParameter<double>("chi2Scale"),
0087
0088 {(float)pset.getParameter<double>("tripletMaxTip"),
0089 (float)pset.getParameter<double>("tripletMinPt"),
0090 (float)pset.getParameter<double>("tripletMaxZip")},
0091
0092 {(float)pset.getParameter<double>("quadrupletMaxTip"),
0093 (float)pset.getParameter<double>("quadrupletMinPt"),
0094 (float)pset.getParameter<double>("quadrupletMaxZip")}};
0095 }
0096 };
0097
0098 template <typename TrackerTraits>
0099 struct topologyCuts<TrackerTraits, isPhase2Topology<TrackerTraits>> {
0100 static constexpr CAParamsT<TrackerTraits> makeCACuts(edm::ParameterSet const& cfg) {
0101 return CAParamsT<TrackerTraits>{{cfg.getParameter<unsigned int>("maxNumberOfDoublets"),
0102 cfg.getParameter<unsigned int>("minHitsPerNtuplet"),
0103 static_cast<float>(cfg.getParameter<double>("ptmin")),
0104 static_cast<float>(cfg.getParameter<double>("CAThetaCutBarrel")),
0105 static_cast<float>(cfg.getParameter<double>("CAThetaCutForward")),
0106 static_cast<float>(cfg.getParameter<double>("hardCurvCut")),
0107 static_cast<float>(cfg.getParameter<double>("dcaCutInnerTriplet")),
0108 static_cast<float>(cfg.getParameter<double>("dcaCutOuterTriplet"))},
0109 {(bool)cfg.getParameter<bool>("includeFarForwards")}};
0110 }
0111
0112 static constexpr pixelTrack::QualityCutsT<TrackerTraits> makeQualityCuts(edm::ParameterSet const& pset) {
0113 return pixelTrack::QualityCutsT<TrackerTraits>{
0114 static_cast<float>(pset.getParameter<double>("maxChi2")),
0115 static_cast<float>(pset.getParameter<double>("minPt")),
0116 static_cast<float>(pset.getParameter<double>("maxTip")),
0117 static_cast<float>(pset.getParameter<double>("maxZip")),
0118 };
0119 }
0120 };
0121
0122
0123
0124
0125 template <typename TrackerTraits>
0126 CellCutsT<TrackerTraits> makeCellCuts(edm::ParameterSet const& cfg) {
0127 return CellCutsT<TrackerTraits>{cfg.getParameter<bool>("doClusterCut"),
0128 cfg.getParameter<bool>("doZ0Cut"),
0129 cfg.getParameter<bool>("doPtCut"),
0130 cfg.getParameter<bool>("idealConditions"),
0131 (float)cfg.getParameter<double>("z0Cut"),
0132 (float)cfg.getParameter<double>("ptCut"),
0133 cfg.getParameter<std::vector<int>>("phiCuts")};
0134 }
0135
0136 }
0137
0138 using namespace std;
0139
0140 template <typename TrackerTraits>
0141 CAHitNtupletGeneratorOnGPU<TrackerTraits>::CAHitNtupletGeneratorOnGPU(const edm::ParameterSet& cfg,
0142 edm::ConsumesCollector& iC)
0143 : m_params(makeCommonParams(cfg),
0144 makeCellCuts<TrackerTraits>(cfg),
0145 topologyCuts<TrackerTraits>::makeQualityCuts(cfg.getParameterSet("trackQualityCuts")),
0146 topologyCuts<TrackerTraits>::makeCACuts(cfg)) {
0147 #ifdef DUMP_GPU_TK_TUPLES
0148 printf("TK: %s %s %s %s %s %s %s %s %s %s %s %s %s %s %s %s %s\n",
0149 "tid",
0150 "qual",
0151 "nh",
0152 "nl",
0153 "charge",
0154 "pt",
0155 "eta",
0156 "phi",
0157 "tip",
0158 "zip",
0159 "chi2",
0160 "h1",
0161 "h2",
0162 "h3",
0163 "h4",
0164 "h5",
0165 "hn");
0166 #endif
0167 }
0168
0169 template <typename TrackerTraits>
0170 void CAHitNtupletGeneratorOnGPU<TrackerTraits>::fillDescriptions(edm::ParameterSetDescription& desc) {
0171 fillDescriptionsCommon(desc);
0172 edm::LogWarning("CAHitNtupletGeneratorOnGPU::fillDescriptions")
0173 << "Note: this fillDescriptions is a dummy one. Most probably you are missing some parameters. \n"
0174 "please implement your TrackerTraits descriptions in CAHitNtupletGeneratorOnGPU. \n";
0175 }
0176
0177 template <>
0178 void CAHitNtupletGeneratorOnGPU<pixelTopology::Phase1>::fillDescriptions(edm::ParameterSetDescription& desc) {
0179 fillDescriptionsCommon(desc);
0180
0181 desc.add<bool>("idealConditions", true);
0182 desc.add<bool>("includeJumpingForwardDoublets", false);
0183 desc.add<double>("z0Cut", 12.0);
0184 desc.add<double>("ptCut", 0.5);
0185
0186 edm::ParameterSetDescription trackQualityCuts;
0187 trackQualityCuts.add<double>("chi2MaxPt", 10.)->setComment("max pT used to determine the pT-dependent chi2 cut");
0188 trackQualityCuts.add<std::vector<double>>("chi2Coeff", {0.9, 1.8})->setComment("chi2 at 1GeV and at ptMax above");
0189 trackQualityCuts.add<double>("chi2Scale", 8.)
0190 ->setComment(
0191 "Factor to multiply the pT-dependent chi2 cut (currently: 8 for the broken line fit, ?? for the Riemann "
0192 "fit)");
0193 trackQualityCuts.add<double>("tripletMinPt", 0.5)->setComment("Min pT for triplets, in GeV");
0194 trackQualityCuts.add<double>("tripletMaxTip", 0.3)->setComment("Max |Tip| for triplets, in cm");
0195 trackQualityCuts.add<double>("tripletMaxZip", 12.)->setComment("Max |Zip| for triplets, in cm");
0196 trackQualityCuts.add<double>("quadrupletMinPt", 0.3)->setComment("Min pT for quadruplets, in GeV");
0197 trackQualityCuts.add<double>("quadrupletMaxTip", 0.5)->setComment("Max |Tip| for quadruplets, in cm");
0198 trackQualityCuts.add<double>("quadrupletMaxZip", 12.)->setComment("Max |Zip| for quadruplets, in cm");
0199
0200 desc.add<std::vector<int>>(
0201 "phiCuts", std::vector<int>(std::begin(phase1PixelTopology::phicuts), std::end(phase1PixelTopology::phicuts)))
0202 ->setComment("Cuts in phi for cells");
0203
0204 desc.add<edm::ParameterSetDescription>("trackQualityCuts", trackQualityCuts)
0205 ->setComment(
0206 "Quality cuts based on the results of the track fit:\n - apply a pT-dependent chi2 cut;\n - apply \"region "
0207 "cuts\" based on the fit results (pT, Tip, Zip).");
0208 }
0209
0210 template <>
0211 void CAHitNtupletGeneratorOnGPU<pixelTopology::HIonPhase1>::fillDescriptions(edm::ParameterSetDescription& desc) {
0212 fillDescriptionsCommon(desc);
0213
0214 desc.add<bool>("idealConditions", false);
0215 desc.add<bool>("includeJumpingForwardDoublets", false);
0216 desc.add<double>("z0Cut", 10.0);
0217 desc.add<double>("ptCut", 0.0);
0218
0219 edm::ParameterSetDescription trackQualityCuts;
0220 trackQualityCuts.add<double>("chi2MaxPt", 10.)->setComment("max pT used to determine the pT-dependent chi2 cut");
0221 trackQualityCuts.add<std::vector<double>>("chi2Coeff", {0.9, 1.8})->setComment("chi2 at 1GeV and at ptMax above");
0222 trackQualityCuts.add<double>("chi2Scale", 8.)
0223 ->setComment(
0224 "Factor to multiply the pT-dependent chi2 cut (currently: 8 for the broken line fit, ?? for the Riemann "
0225 "fit)");
0226 trackQualityCuts.add<double>("tripletMinPt", 0.0)->setComment("Min pT for triplets, in GeV");
0227 trackQualityCuts.add<double>("tripletMaxTip", 0.1)->setComment("Max |Tip| for triplets, in cm");
0228 trackQualityCuts.add<double>("tripletMaxZip", 6.)->setComment("Max |Zip| for triplets, in cm");
0229 trackQualityCuts.add<double>("quadrupletMinPt", 0.0)->setComment("Min pT for quadruplets, in GeV");
0230 trackQualityCuts.add<double>("quadrupletMaxTip", 0.5)->setComment("Max |Tip| for quadruplets, in cm");
0231 trackQualityCuts.add<double>("quadrupletMaxZip", 6.)->setComment("Max |Zip| for quadruplets, in cm");
0232
0233 desc.add<std::vector<int>>(
0234 "phiCuts", std::vector<int>(std::begin(phase1PixelTopology::phicuts), std::end(phase1PixelTopology::phicuts)))
0235 ->setComment("Cuts in phi for cells");
0236
0237 desc.add<edm::ParameterSetDescription>("trackQualityCuts", trackQualityCuts)
0238 ->setComment(
0239 "Quality cuts based on the results of the track fit:\n - apply a pT-dependent chi2 cut;\n - apply \"region "
0240 "cuts\" based on the fit results (pT, Tip, Zip).");
0241 }
0242
0243 template <>
0244 void CAHitNtupletGeneratorOnGPU<pixelTopology::Phase2>::fillDescriptions(edm::ParameterSetDescription& desc) {
0245 fillDescriptionsCommon(desc);
0246
0247 desc.add<bool>("idealConditions", false);
0248 desc.add<bool>("includeFarForwards", true);
0249 desc.add<bool>("includeJumpingForwardDoublets", true);
0250 desc.add<double>("z0Cut", 7.5);
0251 desc.add<double>("ptCut", 0.85);
0252
0253 edm::ParameterSetDescription trackQualityCuts;
0254 trackQualityCuts.add<double>("maxChi2", 5.)->setComment("Max normalized chi2");
0255 trackQualityCuts.add<double>("minPt", 0.5)->setComment("Min pT in GeV");
0256 trackQualityCuts.add<double>("maxTip", 0.3)->setComment("Max |Tip| in cm");
0257 trackQualityCuts.add<double>("maxZip", 12.)->setComment("Max |Zip|, in cm");
0258
0259 desc.add<std::vector<int>>(
0260 "phiCuts", std::vector<int>(std::begin(phase2PixelTopology::phicuts), std::end(phase2PixelTopology::phicuts)))
0261 ->setComment("Cuts in phi for cells");
0262
0263 desc.add<edm::ParameterSetDescription>("trackQualityCuts", trackQualityCuts)
0264 ->setComment(
0265 "Quality cuts based on the results of the track fit:\n - apply cuts based on the fit results (pT, Tip, "
0266 "Zip).");
0267 }
0268
0269 template <typename TrackerTraits>
0270 void CAHitNtupletGeneratorOnGPU<TrackerTraits>::fillDescriptionsCommon(edm::ParameterSetDescription& desc) {
0271
0272
0273
0274 desc.add<double>("ptmin", 0.9)->setComment("Cut on minimum pt");
0275 desc.add<double>("CAThetaCutBarrel", 0.002)->setComment("Cut on RZ alignement for Barrel");
0276 desc.add<double>("CAThetaCutForward", 0.003)->setComment("Cut on RZ alignment for Forward");
0277 desc.add<double>("hardCurvCut", 1. / (0.35 * 87.))->setComment("Cut on minimum curvature");
0278 desc.add<double>("dcaCutInnerTriplet", 0.15)->setComment("Cut on origin radius when the inner hit is on BPix1");
0279 desc.add<double>("dcaCutOuterTriplet", 0.25)->setComment("Cut on origin radius when the outer hit is on BPix1");
0280 desc.add<bool>("earlyFishbone", true);
0281 desc.add<bool>("lateFishbone", false);
0282 desc.add<bool>("fillStatistics", false);
0283 desc.add<unsigned int>("minHitsPerNtuplet", 4);
0284 desc.add<unsigned int>("maxNumberOfDoublets", TrackerTraits::maxNumberOfDoublets);
0285 desc.add<unsigned int>("minHitsForSharingCut", 10)
0286 ->setComment("Maximum number of hits in a tuple to clean also if the shared hit is on bpx1");
0287
0288 desc.add<bool>("fitNas4", false)->setComment("fit only 4 hits out of N");
0289 desc.add<bool>("doClusterCut", true);
0290 desc.add<bool>("doZ0Cut", true);
0291 desc.add<bool>("doPtCut", true);
0292 desc.add<bool>("useRiemannFit", false)->setComment("true for Riemann, false for BrokenLine");
0293 desc.add<bool>("doSharedHitCut", true)->setComment("Sharing hit nTuples cleaning");
0294 desc.add<bool>("dupPassThrough", false)->setComment("Do not reject duplicate");
0295 desc.add<bool>("useSimpleTripletCleaner", true)->setComment("use alternate implementation");
0296 }
0297
0298 template <typename TrackerTraits>
0299 void CAHitNtupletGeneratorOnGPU<TrackerTraits>::beginJob() {
0300 if (m_params.onGPU_) {
0301
0302 edm::Service<CUDAInterface> cuda;
0303 if (cuda and cuda->enabled()) {
0304 cudaCheck(cudaMalloc(&m_counters, sizeof(Counters)));
0305 cudaCheck(cudaMemset(m_counters, 0, sizeof(Counters)));
0306 }
0307 } else {
0308 m_counters = new Counters();
0309 memset(m_counters, 0, sizeof(Counters));
0310 }
0311 }
0312
0313 template <typename TrackerTraits>
0314 void CAHitNtupletGeneratorOnGPU<TrackerTraits>::endJob() {
0315 if (m_params.onGPU_) {
0316
0317 edm::Service<CUDAInterface> cuda;
0318 if (cuda and cuda->enabled()) {
0319 if (m_params.doStats_) {
0320
0321 CAHitNtupletGeneratorKernelsGPU<TrackerTraits>::printCounters(m_counters);
0322 }
0323 cudaFree(m_counters);
0324 }
0325 } else {
0326 if (m_params.doStats_) {
0327 CAHitNtupletGeneratorKernelsCPU<TrackerTraits>::printCounters(m_counters);
0328 }
0329 delete m_counters;
0330 }
0331 }
0332
0333 template <typename TrackerTraits>
0334 TrackSoAHeterogeneousDevice<TrackerTraits> CAHitNtupletGeneratorOnGPU<TrackerTraits>::makeTuplesAsync(
0335 HitsOnDevice const& hits_d, float bfield, cudaStream_t stream) const {
0336 using HelixFitOnGPU = HelixFitOnGPU<TrackerTraits>;
0337 using TrackSoA = TrackSoAHeterogeneousDevice<TrackerTraits>;
0338 using GPUKernels = CAHitNtupletGeneratorKernelsGPU<TrackerTraits>;
0339
0340 TrackSoA tracks(stream);
0341
0342 GPUKernels kernels(m_params);
0343 kernels.setCounters(m_counters);
0344 kernels.allocateOnGPU(hits_d.nHits(), stream);
0345
0346 kernels.buildDoublets(hits_d.view(), hits_d.offsetBPIX2(), stream);
0347
0348 kernels.launchKernels(hits_d.view(), tracks.view(), stream);
0349
0350 HelixFitOnGPU fitter(bfield, m_params.fitNas4_);
0351 fitter.allocateOnGPU(kernels.tupleMultiplicity(), tracks.view());
0352 if (m_params.useRiemannFit_) {
0353 fitter.launchRiemannKernels(hits_d.view(), hits_d.nHits(), TrackerTraits::maxNumberOfQuadruplets, stream);
0354 } else {
0355 fitter.launchBrokenLineKernels(hits_d.view(), hits_d.nHits(), TrackerTraits::maxNumberOfQuadruplets, stream);
0356 }
0357 kernels.classifyTuples(hits_d.view(), tracks.view(), stream);
0358 #ifdef GPU_DEBUG
0359 cudaDeviceSynchronize();
0360 cudaCheck(cudaGetLastError());
0361 std::cout << "finished building pixel tracks on GPU" << std::endl;
0362 #endif
0363
0364 return tracks;
0365 }
0366
0367 template <typename TrackerTraits>
0368 TrackSoAHeterogeneousHost<TrackerTraits> CAHitNtupletGeneratorOnGPU<TrackerTraits>::makeTuples(HitsOnHost const& hits_h,
0369 float bfield) const {
0370 using HelixFitOnGPU = HelixFitOnGPU<TrackerTraits>;
0371 using TrackSoA = TrackSoAHeterogeneousHost<TrackerTraits>;
0372 using CPUKernels = CAHitNtupletGeneratorKernelsCPU<TrackerTraits>;
0373
0374 TrackSoA tracks;
0375
0376 CPUKernels kernels(m_params);
0377 kernels.setCounters(m_counters);
0378 kernels.allocateOnGPU(hits_h.nHits(), nullptr);
0379
0380 kernels.buildDoublets(hits_h.view(), hits_h.offsetBPIX2(), nullptr);
0381 kernels.launchKernels(hits_h.view(), tracks.view(), nullptr);
0382
0383 if (0 == hits_h.nHits())
0384 return tracks;
0385
0386
0387 HelixFitOnGPU fitter(bfield, m_params.fitNas4_);
0388 fitter.allocateOnGPU(kernels.tupleMultiplicity(), tracks.view());
0389
0390 if (m_params.useRiemannFit_) {
0391 fitter.launchRiemannKernelsOnCPU(hits_h.view(), hits_h.nHits(), TrackerTraits::maxNumberOfQuadruplets);
0392 } else {
0393 fitter.launchBrokenLineKernelsOnCPU(hits_h.view(), hits_h.nHits(), TrackerTraits::maxNumberOfQuadruplets);
0394 }
0395
0396 kernels.classifyTuples(hits_h.view(), tracks.view(), nullptr);
0397
0398 #ifdef GPU_DEBUG
0399 std::cout << "finished building pixel tracks on CPU" << std::endl;
0400 #endif
0401
0402
0403 auto maxTracks = tracks.view().metadata().size();
0404 auto nTracks = tracks.view().nTracks();
0405 assert(nTracks < maxTracks);
0406 if (nTracks == maxTracks - 1) {
0407 edm::LogWarning("PixelTracks") << "Unsorted reconstructed pixel tracks truncated to " << maxTracks - 1
0408 << " candidates";
0409 }
0410
0411 return tracks;
0412 }
0413
0414 template class CAHitNtupletGeneratorOnGPU<pixelTopology::Phase1>;
0415 template class CAHitNtupletGeneratorOnGPU<pixelTopology::Phase2>;
0416 template class CAHitNtupletGeneratorOnGPU<pixelTopology::HIonPhase1>;