Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2025-02-09 23:42:07

0001 //
0002 // Original Author: Felice Pantaleo, CERN
0003 //
0004 
0005 //#define GPU_DEBUG
0006 //#define DUMP_GPU_TK_TUPLES
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   //Common Params
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   //This is needed to have the partial specialization for  isPhase1Topology/isPhase2Topology
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>{// polynomial coefficients for the pT-dependent chi2 cut
0082                                                      {(float)coeff[0], (float)coeff[1], 0.f, 0.f},
0083                                                      // max pT used to determine the chi2 cut
0084                                                      (float)ptMax,
0085                                                      // chi2 scale factor: 8 for broken line fit, ?? for Riemann fit
0086                                                      (float)pset.getParameter<double>("chi2Scale"),
0087                                                      // regional cuts for triplets
0088                                                      {(float)pset.getParameter<double>("tripletMaxTip"),
0089                                                       (float)pset.getParameter<double>("tripletMinPt"),
0090                                                       (float)pset.getParameter<double>("tripletMaxZip")},
0091                                                      // regional cuts for quadruplets
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   //Cell Cuts, as they are the cuts have the same logic for Phase2 and Phase1
0123   //keeping them separate would allow further differentiation in the future
0124   //moving them to topologyCuts and using the same syntax
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<int>("minYsizeB1"),
0134                                     cfg.getParameter<int>("minYsizeB2"),
0135                                     cfg.getParameter<std::vector<int>>("phiCuts")};
0136   }
0137 
0138 }  // namespace
0139 
0140 using namespace std;
0141 
0142 template <typename TrackerTraits>
0143 CAHitNtupletGeneratorOnGPU<TrackerTraits>::CAHitNtupletGeneratorOnGPU(const edm::ParameterSet& cfg,
0144                                                                       edm::ConsumesCollector& iC)
0145     : m_params(makeCommonParams(cfg),
0146                makeCellCuts<TrackerTraits>(cfg),
0147                topologyCuts<TrackerTraits>::makeQualityCuts(cfg.getParameterSet("trackQualityCuts")),
0148                topologyCuts<TrackerTraits>::makeCACuts(cfg)) {
0149 #ifdef DUMP_GPU_TK_TUPLES
0150   printf("TK: %s %s %s %s %s %s %s %s %s %s %s %s %s %s %s %s %s\n",
0151          "tid",
0152          "qual",
0153          "nh",
0154          "nl",
0155          "charge",
0156          "pt",
0157          "eta",
0158          "phi",
0159          "tip",
0160          "zip",
0161          "chi2",
0162          "h1",
0163          "h2",
0164          "h3",
0165          "h4",
0166          "h5",
0167          "hn");
0168 #endif
0169 }
0170 
0171 template <typename TrackerTraits>
0172 void CAHitNtupletGeneratorOnGPU<TrackerTraits>::fillDescriptions(edm::ParameterSetDescription& desc) {
0173   fillDescriptionsCommon(desc);
0174   edm::LogWarning("CAHitNtupletGeneratorOnGPU::fillDescriptions")
0175       << "Note: this fillDescriptions is a dummy one. Most probably you are missing some parameters. \n"
0176          "please implement your TrackerTraits descriptions in CAHitNtupletGeneratorOnGPU. \n";
0177 }
0178 
0179 template <>
0180 void CAHitNtupletGeneratorOnGPU<pixelTopology::Phase1>::fillDescriptions(edm::ParameterSetDescription& desc) {
0181   fillDescriptionsCommon(desc);
0182 
0183   desc.add<bool>("idealConditions", true);
0184   desc.add<bool>("includeJumpingForwardDoublets", false);
0185   desc.add<double>("z0Cut", 12.0);
0186   desc.add<double>("ptCut", 0.5);
0187 
0188   edm::ParameterSetDescription trackQualityCuts;
0189   trackQualityCuts.add<double>("chi2MaxPt", 10.)->setComment("max pT used to determine the pT-dependent chi2 cut");
0190   trackQualityCuts.add<std::vector<double>>("chi2Coeff", {0.9, 1.8})->setComment("chi2 at 1GeV and at ptMax above");
0191   trackQualityCuts.add<double>("chi2Scale", 8.)
0192       ->setComment(
0193           "Factor to multiply the pT-dependent chi2 cut (currently: 8 for the broken line fit, ?? for the Riemann "
0194           "fit)");
0195   trackQualityCuts.add<double>("tripletMinPt", 0.5)->setComment("Min pT for triplets, in GeV");
0196   trackQualityCuts.add<double>("tripletMaxTip", 0.3)->setComment("Max |Tip| for triplets, in cm");
0197   trackQualityCuts.add<double>("tripletMaxZip", 12.)->setComment("Max |Zip| for triplets, in cm");
0198   trackQualityCuts.add<double>("quadrupletMinPt", 0.3)->setComment("Min pT for quadruplets, in GeV");
0199   trackQualityCuts.add<double>("quadrupletMaxTip", 0.5)->setComment("Max |Tip| for quadruplets, in cm");
0200   trackQualityCuts.add<double>("quadrupletMaxZip", 12.)->setComment("Max |Zip| for quadruplets, in cm");
0201 
0202   desc.add<int>("minYsizeB1", 1)->setComment("Min Y cluster size in pixel B1");
0203   desc.add<int>("minYsizeB2", 1)->setComment("Min Y cluster size in pixel B2");
0204 
0205   desc.add<std::vector<int>>(
0206           "phiCuts", std::vector<int>(std::begin(phase1PixelTopology::phicuts), std::end(phase1PixelTopology::phicuts)))
0207       ->setComment("Cuts in phi for cells");
0208 
0209   desc.add<edm::ParameterSetDescription>("trackQualityCuts", trackQualityCuts)
0210       ->setComment(
0211           "Quality cuts based on the results of the track fit:\n  - apply a pT-dependent chi2 cut;\n  - apply \"region "
0212           "cuts\" based on the fit results (pT, Tip, Zip).");
0213 }
0214 
0215 template <>
0216 void CAHitNtupletGeneratorOnGPU<pixelTopology::HIonPhase1>::fillDescriptions(edm::ParameterSetDescription& desc) {
0217   fillDescriptionsCommon(desc);
0218 
0219   desc.add<bool>("idealConditions", false);
0220   desc.add<bool>("includeJumpingForwardDoublets", false);
0221   desc.add<double>("z0Cut", 10.0);
0222   desc.add<double>("ptCut", 0.0);
0223 
0224   edm::ParameterSetDescription trackQualityCuts;
0225   trackQualityCuts.add<double>("chi2MaxPt", 10.)->setComment("max pT used to determine the pT-dependent chi2 cut");
0226   trackQualityCuts.add<std::vector<double>>("chi2Coeff", {0.9, 1.8})->setComment("chi2 at 1GeV and at ptMax above");
0227   trackQualityCuts.add<double>("chi2Scale", 8.)
0228       ->setComment(
0229           "Factor to multiply the pT-dependent chi2 cut (currently: 8 for the broken line fit, ?? for the Riemann "
0230           "fit)");
0231   trackQualityCuts.add<double>("tripletMinPt", 0.0)->setComment("Min pT for triplets, in GeV");
0232   trackQualityCuts.add<double>("tripletMaxTip", 0.1)->setComment("Max |Tip| for triplets, in cm");
0233   trackQualityCuts.add<double>("tripletMaxZip", 6.)->setComment("Max |Zip| for triplets, in cm");
0234   trackQualityCuts.add<double>("quadrupletMinPt", 0.0)->setComment("Min pT for quadruplets, in GeV");
0235   trackQualityCuts.add<double>("quadrupletMaxTip", 0.5)->setComment("Max |Tip| for quadruplets, in cm");
0236   trackQualityCuts.add<double>("quadrupletMaxZip", 6.)->setComment("Max |Zip| for quadruplets, in cm");
0237 
0238   desc.add<int>("minYsizeB1", 36)->setComment("Min Y cluster size in pixel B1");
0239   desc.add<int>("minYsizeB2", 28)->setComment("Min Y cluster size in pixel B2");
0240 
0241   desc.add<std::vector<int>>(
0242           "phiCuts", std::vector<int>(std::begin(phase1PixelTopology::phicuts), std::end(phase1PixelTopology::phicuts)))
0243       ->setComment("Cuts in phi for cells");
0244 
0245   desc.add<edm::ParameterSetDescription>("trackQualityCuts", trackQualityCuts)
0246       ->setComment(
0247           "Quality cuts based on the results of the track fit:\n  - apply a pT-dependent chi2 cut;\n  - apply \"region "
0248           "cuts\" based on the fit results (pT, Tip, Zip).");
0249 }
0250 
0251 template <>
0252 void CAHitNtupletGeneratorOnGPU<pixelTopology::Phase2>::fillDescriptions(edm::ParameterSetDescription& desc) {
0253   fillDescriptionsCommon(desc);
0254 
0255   desc.add<bool>("idealConditions", false);
0256   desc.add<bool>("includeFarForwards", true);
0257   desc.add<bool>("includeJumpingForwardDoublets", true);
0258   desc.add<double>("z0Cut", 7.5);
0259   desc.add<double>("ptCut", 0.85);
0260 
0261   edm::ParameterSetDescription trackQualityCuts;
0262   trackQualityCuts.add<double>("maxChi2", 5.)->setComment("Max normalized chi2");
0263   trackQualityCuts.add<double>("minPt", 0.5)->setComment("Min pT in GeV");
0264   trackQualityCuts.add<double>("maxTip", 0.3)->setComment("Max |Tip| in cm");
0265   trackQualityCuts.add<double>("maxZip", 12.)->setComment("Max |Zip|, in cm");
0266 
0267   desc.add<int>("minYsizeB1", 25)->setComment("Min Y cluster size in pixel B1");
0268   desc.add<int>("minYsizeB2", 15)->setComment("Min Y cluster size in pixel B2");
0269 
0270   desc.add<std::vector<int>>(
0271           "phiCuts", std::vector<int>(std::begin(phase2PixelTopology::phicuts), std::end(phase2PixelTopology::phicuts)))
0272       ->setComment("Cuts in phi for cells");
0273 
0274   desc.add<edm::ParameterSetDescription>("trackQualityCuts", trackQualityCuts)
0275       ->setComment(
0276           "Quality cuts based on the results of the track fit:\n  - apply cuts based on the fit results (pT, Tip, "
0277           "Zip).");
0278 }
0279 
0280 template <typename TrackerTraits>
0281 void CAHitNtupletGeneratorOnGPU<TrackerTraits>::fillDescriptionsCommon(edm::ParameterSetDescription& desc) {
0282   // 87 cm/GeV = 1/(3.8T * 0.3)
0283   // take less than radius given by the hardPtCut and reject everything below
0284   // auto hardCurvCut = 1.f/(0.35 * 87.f);
0285   desc.add<double>("ptmin", 0.9)->setComment("Cut on minimum pt");
0286   desc.add<double>("CAThetaCutBarrel", 0.002)->setComment("Cut on RZ alignement for Barrel");
0287   desc.add<double>("CAThetaCutForward", 0.003)->setComment("Cut on RZ alignment for Forward");
0288   desc.add<double>("hardCurvCut", 1. / (0.35 * 87.))->setComment("Cut on minimum curvature");
0289   desc.add<double>("dcaCutInnerTriplet", 0.15)->setComment("Cut on origin radius when the inner hit is on BPix1");
0290   desc.add<double>("dcaCutOuterTriplet", 0.25)->setComment("Cut on origin radius when the outer hit is on BPix1");
0291   desc.add<bool>("earlyFishbone", true);
0292   desc.add<bool>("lateFishbone", false);
0293   desc.add<bool>("fillStatistics", false);
0294   desc.add<unsigned int>("minHitsPerNtuplet", 4);
0295   desc.add<unsigned int>("maxNumberOfDoublets", TrackerTraits::maxNumberOfDoublets);
0296   desc.add<unsigned int>("minHitsForSharingCut", 10)
0297       ->setComment("Maximum number of hits in a tuple to clean also if the shared hit is on bpx1");
0298 
0299   desc.add<bool>("fitNas4", false)->setComment("fit only 4 hits out of N");
0300   desc.add<bool>("doClusterCut", true);
0301   desc.add<bool>("doZ0Cut", true);
0302   desc.add<bool>("doPtCut", true);
0303   desc.add<bool>("useRiemannFit", false)->setComment("true for Riemann, false for BrokenLine");
0304   desc.add<bool>("doSharedHitCut", true)->setComment("Sharing hit nTuples cleaning");
0305   desc.add<bool>("dupPassThrough", false)->setComment("Do not reject duplicate");
0306   desc.add<bool>("useSimpleTripletCleaner", true)->setComment("use alternate implementation");
0307 }
0308 
0309 template <typename TrackerTraits>
0310 void CAHitNtupletGeneratorOnGPU<TrackerTraits>::beginJob() {
0311   if (m_params.onGPU_) {
0312     // allocate pinned host memory only if CUDA is available
0313     edm::Service<CUDAInterface> cuda;
0314     if (cuda and cuda->enabled()) {
0315       cudaCheck(cudaMalloc(&m_counters, sizeof(Counters)));
0316       cudaCheck(cudaMemset(m_counters, 0, sizeof(Counters)));
0317     }
0318   } else {
0319     m_counters = new Counters();
0320     memset(m_counters, 0, sizeof(Counters));
0321   }
0322 }
0323 
0324 template <typename TrackerTraits>
0325 void CAHitNtupletGeneratorOnGPU<TrackerTraits>::endJob() {
0326   if (m_params.onGPU_) {
0327     // print the gpu statistics and free pinned host memory only if CUDA is available
0328     edm::Service<CUDAInterface> cuda;
0329     if (cuda and cuda->enabled()) {
0330       if (m_params.doStats_) {
0331         // crash on multi-gpu processes
0332         CAHitNtupletGeneratorKernelsGPU<TrackerTraits>::printCounters(m_counters);
0333       }
0334       cudaFree(m_counters);
0335     }
0336   } else {
0337     if (m_params.doStats_) {
0338       CAHitNtupletGeneratorKernelsCPU<TrackerTraits>::printCounters(m_counters);
0339     }
0340     delete m_counters;
0341   }
0342 }
0343 
0344 template <typename TrackerTraits>
0345 TrackSoAHeterogeneousDevice<TrackerTraits> CAHitNtupletGeneratorOnGPU<TrackerTraits>::makeTuplesAsync(
0346     HitsOnDevice const& hits_d, float bfield, cudaStream_t stream) const {
0347   using HelixFitOnGPU = HelixFitOnGPU<TrackerTraits>;
0348   using TrackSoA = TrackSoAHeterogeneousDevice<TrackerTraits>;
0349   using GPUKernels = CAHitNtupletGeneratorKernelsGPU<TrackerTraits>;
0350 
0351   TrackSoA tracks(stream);
0352 
0353   GPUKernels kernels(m_params);
0354   kernels.setCounters(m_counters);
0355   kernels.allocateOnGPU(hits_d.nHits(), stream);
0356 
0357   kernels.buildDoublets(hits_d.view(), hits_d.offsetBPIX2(), stream);
0358 
0359   kernels.launchKernels(hits_d.view(), tracks.view(), stream);
0360 
0361   HelixFitOnGPU fitter(bfield, m_params.fitNas4_);
0362   fitter.allocateOnGPU(kernels.tupleMultiplicity(), tracks.view());
0363   if (m_params.useRiemannFit_) {
0364     fitter.launchRiemannKernels(hits_d.view(), hits_d.nHits(), TrackerTraits::maxNumberOfQuadruplets, stream);
0365   } else {
0366     fitter.launchBrokenLineKernels(hits_d.view(), hits_d.nHits(), TrackerTraits::maxNumberOfQuadruplets, stream);
0367   }
0368   kernels.classifyTuples(hits_d.view(), tracks.view(), stream);
0369 #ifdef GPU_DEBUG
0370   cudaDeviceSynchronize();
0371   cudaCheck(cudaGetLastError());
0372   std::cout << "finished building pixel tracks on GPU" << std::endl;
0373 #endif
0374 
0375   return tracks;
0376 }
0377 
0378 template <typename TrackerTraits>
0379 TrackSoAHeterogeneousHost<TrackerTraits> CAHitNtupletGeneratorOnGPU<TrackerTraits>::makeTuples(HitsOnHost const& hits_h,
0380                                                                                                float bfield) const {
0381   using HelixFitOnGPU = HelixFitOnGPU<TrackerTraits>;
0382   using TrackSoA = TrackSoAHeterogeneousHost<TrackerTraits>;
0383   using CPUKernels = CAHitNtupletGeneratorKernelsCPU<TrackerTraits>;
0384 
0385   TrackSoA tracks;
0386 
0387   CPUKernels kernels(m_params);
0388   kernels.setCounters(m_counters);
0389   kernels.allocateOnGPU(hits_h.nHits(), nullptr);
0390 
0391   kernels.buildDoublets(hits_h.view(), hits_h.offsetBPIX2(), nullptr);
0392   kernels.launchKernels(hits_h.view(), tracks.view(), nullptr);
0393 
0394   if (0 == hits_h.nHits())
0395     return tracks;
0396 
0397   // now fit
0398   HelixFitOnGPU fitter(bfield, m_params.fitNas4_);
0399   fitter.allocateOnGPU(kernels.tupleMultiplicity(), tracks.view());
0400 
0401   if (m_params.useRiemannFit_) {
0402     fitter.launchRiemannKernelsOnCPU(hits_h.view(), hits_h.nHits(), TrackerTraits::maxNumberOfQuadruplets);
0403   } else {
0404     fitter.launchBrokenLineKernelsOnCPU(hits_h.view(), hits_h.nHits(), TrackerTraits::maxNumberOfQuadruplets);
0405   }
0406 
0407   kernels.classifyTuples(hits_h.view(), tracks.view(), nullptr);
0408 
0409 #ifdef GPU_DEBUG
0410   std::cout << "finished building pixel tracks on CPU" << std::endl;
0411 #endif
0412 
0413   // check that the fixed-size SoA does not overflow
0414   auto maxTracks = tracks.view().metadata().size();
0415   auto nTracks = tracks.view().nTracks();
0416   assert(nTracks < maxTracks);
0417   if (nTracks == maxTracks - 1) {
0418     edm::LogWarning("PixelTracks") << "Unsorted reconstructed pixel tracks truncated to " << maxTracks - 1
0419                                    << " candidates";
0420   }
0421 
0422   return tracks;
0423 }
0424 
0425 template class CAHitNtupletGeneratorOnGPU<pixelTopology::Phase1>;
0426 template class CAHitNtupletGeneratorOnGPU<pixelTopology::Phase2>;
0427 template class CAHitNtupletGeneratorOnGPU<pixelTopology::HIonPhase1>;