Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2024-04-06 12:26:24

0001 #include <alpaka/alpaka.hpp>
0002 
0003 #include "CondFormats/SiPixelTransient/interface/SiPixelGenError.h"
0004 #include "DataFormats/GeometrySurface/interface/SOARotation.h"
0005 #include "DataFormats/SiPixelClusterSoA/interface/ClusteringConstants.h"
0006 #include "DataFormats/TrackingRecHitSoA/interface/SiPixelHitStatus.h"
0007 #include "Geometry/CommonTopologies/interface/SimplePixelTopology.h"
0008 #include "HeterogeneousCore/AlpakaInterface/interface/CopyToDevice.h"
0009 #include "HeterogeneousCore/AlpakaInterface/interface/config.h"
0010 #include "HeterogeneousCore/AlpakaInterface/interface/memory.h"
0011 #include "RecoLocalTracker/SiPixelRecHits/interface/PixelCPEFastParamsHost.h"
0012 
0013 //-----------------------------------------------------------------------------
0014 //!  The constructor.
0015 //-----------------------------------------------------------------------------
0016 template <typename TrackerTraits>
0017 PixelCPEFastParamsHost<TrackerTraits>::PixelCPEFastParamsHost(edm::ParameterSet const& conf,
0018                                                               const MagneticField* mag,
0019                                                               const TrackerGeometry& geom,
0020                                                               const TrackerTopology& ttopo,
0021                                                               const SiPixelLorentzAngle* lorentzAngle,
0022                                                               const SiPixelGenErrorDBObject* genErrorDBObject,
0023                                                               const SiPixelLorentzAngle* lorentzAngleWidth)
0024     : PixelCPEGenericBase(conf, mag, geom, ttopo, lorentzAngle, genErrorDBObject, lorentzAngleWidth),
0025       buffer_(cms::alpakatools::make_host_buffer<pixelCPEforDevice::ParamsOnDeviceT<TrackerTraits>>()) {
0026   // Use errors from templates or from GenError
0027   if (useErrorsFromTemplates_) {
0028     if (!SiPixelGenError::pushfile(*genErrorDBObject_, this->thePixelGenError_))
0029       throw cms::Exception("InvalidCalibrationLoaded")
0030           << "ERROR: GenErrors not filled correctly. Check the sqlite file. Using SiPixelTemplateDBObject version "
0031           << (*genErrorDBObject_).version();
0032   }
0033 
0034   fillParamsForDevice();
0035 }
0036 
0037 template <typename TrackerTraits>
0038 void PixelCPEFastParamsHost<TrackerTraits>::fillParamsForDevice() {
0039   // this code executes only once per job, computation inefficiency is not an issue
0040   // many code blocks are repeated: better keep the computation local and self consistent as blocks may in future move around, be deleted ...
0041   // It is valid only for Phase1 and the version of GenError in DB used in late 2018 and in 2021
0042 
0043   buffer_->commonParams().theThicknessB = m_DetParams.front().theThickness;
0044   buffer_->commonParams().theThicknessE = m_DetParams.back().theThickness;
0045   buffer_->commonParams().thePitchX = m_DetParams[0].thePitchX;
0046   buffer_->commonParams().thePitchY = m_DetParams[0].thePitchY;
0047 
0048   buffer_->commonParams().numberOfLaddersInBarrel = TrackerTraits::numberOfLaddersInBarrel;
0049 
0050   LogDebug("PixelCPEFastParamsHost") << "pitch & thickness " << buffer_->commonParams().thePitchX << ' '
0051                                      << buffer_->commonParams().thePitchY << "  "
0052                                      << buffer_->commonParams().theThicknessB << ' '
0053                                      << buffer_->commonParams().theThicknessE;
0054 
0055   // zero average geometry
0056   memset(&buffer_->averageGeometry(), 0, sizeof(pixelTopology::AverageGeometryT<TrackerTraits>));
0057 
0058   uint32_t oldLayer = 0;
0059   uint32_t oldLadder = 0;
0060   float rl = 0;
0061   float zl = 0;
0062   float miz = 500, mxz = 0;
0063   float pl = 0;
0064   int nl = 0;
0065 
0066   assert(m_DetParams.size() <= TrackerTraits::numberOfModules);
0067   for (auto i = 0U; i < m_DetParams.size(); ++i) {
0068     auto& p = m_DetParams[i];
0069     auto& g = buffer_->detParams(i);
0070 
0071     g.nRowsRoc = p.theDet->specificTopology().rowsperroc();
0072     g.nColsRoc = p.theDet->specificTopology().colsperroc();
0073     g.nRows = p.theDet->specificTopology().rocsX() * g.nRowsRoc;
0074     g.nCols = p.theDet->specificTopology().rocsY() * g.nColsRoc;
0075 
0076     g.numPixsInModule = g.nRows * g.nCols;
0077 
0078     assert(p.theDet->index() == int(i));
0079     assert(buffer_->commonParams().thePitchY == p.thePitchY);
0080     assert(buffer_->commonParams().thePitchX == p.thePitchX);
0081 
0082     g.isBarrel = GeomDetEnumerators::isBarrel(p.thePart);
0083     g.isPosZ = p.theDet->surface().position().z() > 0;
0084     g.layer = ttopo_.layer(p.theDet->geographicalId());
0085     g.index = i;  // better be!
0086     g.rawId = p.theDet->geographicalId();
0087     auto thickness = g.isBarrel ? buffer_->commonParams().theThicknessB : buffer_->commonParams().theThicknessE;
0088     assert(thickness == p.theThickness);
0089 
0090     auto ladder = ttopo_.pxbLadder(p.theDet->geographicalId());
0091     if (oldLayer != g.layer) {
0092       oldLayer = g.layer;
0093       LogDebug("PixelCPEFastParamsHost") << "new layer at " << i << (g.isBarrel ? " B  " : (g.isPosZ ? " E+ " : " E- "))
0094                                          << g.layer << " starting at " << g.rawId << '\n'
0095                                          << "old layer had " << nl << " ladders";
0096       nl = 0;
0097     }
0098     if (oldLadder != ladder) {
0099       oldLadder = ladder;
0100       LogDebug("PixelCPEFastParamsHost") << "new ladder at " << i
0101                                          << (g.isBarrel ? " B  " : (g.isPosZ ? " E+ " : " E- ")) << ladder
0102                                          << " starting at " << g.rawId << '\n'
0103                                          << "old ladder ave z,r,p mz " << zl / 8.f << " " << rl / 8.f << " " << pl / 8.f
0104                                          << ' ' << miz << ' ' << mxz;
0105       rl = 0;
0106       zl = 0;
0107       pl = 0;
0108       miz = 500;
0109       mxz = 0;
0110       nl++;
0111     }
0112 
0113     g.shiftX = 0.5f * p.lorentzShiftInCmX;
0114     g.shiftY = 0.5f * p.lorentzShiftInCmY;
0115     g.chargeWidthX = p.lorentzShiftInCmX * p.widthLAFractionX;
0116     g.chargeWidthY = p.lorentzShiftInCmY * p.widthLAFractionY;
0117 
0118     g.x0 = p.theOrigin.x();
0119     g.y0 = p.theOrigin.y();
0120     g.z0 = p.theOrigin.z();
0121 
0122     auto vv = p.theDet->surface().position();
0123     auto rr = pixelCPEforDevice::Rotation(p.theDet->surface().rotation());
0124     g.frame = pixelCPEforDevice::Frame(vv.x(), vv.y(), vv.z(), rr);
0125 
0126     zl += vv.z();
0127     miz = std::min(miz, std::abs(vv.z()));
0128     mxz = std::max(mxz, std::abs(vv.z()));
0129     rl += vv.perp();
0130     pl += vv.phi();  // (not obvious)
0131 
0132     // errors .....
0133     ClusterParamGeneric cp;
0134 
0135     cp.with_track_angle = false;
0136 
0137     auto lape = p.theDet->localAlignmentError();
0138     if (lape.invalid())
0139       lape = LocalError();  // zero....
0140 
0141     g.apeXX = lape.xx();
0142     g.apeYY = lape.yy();
0143 
0144     auto toMicron = [&](float x) { return std::min(511, int(x * 1.e4f + 0.5f)); };
0145 
0146     // average angle
0147     auto gvx = p.theOrigin.x() + 40.f * buffer_->commonParams().thePitchX;
0148     auto gvy = p.theOrigin.y();
0149     auto gvz = 1.f / p.theOrigin.z();
0150     //--- Note that the normalization is not required as only the ratio used
0151 
0152     {
0153       // calculate angles (fed into errorFromTemplates)
0154       cp.cotalpha = gvx * gvz;
0155       cp.cotbeta = gvy * gvz;
0156 
0157       errorFromTemplates(p, cp, 20000.);
0158     }
0159 
0160 #ifdef EDM_ML_DEBUG
0161     auto m = 10000.f;
0162     for (float qclus = 15000; qclus < 35000; qclus += 15000) {
0163       errorFromTemplates(p, cp, qclus);
0164       LogDebug("PixelCPEFastParamsHost") << i << ' ' << qclus << ' ' << cp.pixmx << ' ' << m * cp.sigmax << ' '
0165                                          << m * cp.sx1 << ' ' << m * cp.sx2 << ' ' << m * cp.sigmay << ' ' << m * cp.sy1
0166                                          << ' ' << m * cp.sy2;
0167     }
0168     LogDebug("PixelCPEFastParamsHost") << i << ' ' << m * std::sqrt(lape.xx()) << ' ' << m * std::sqrt(lape.yy());
0169 #endif  // EDM_ML_DEBUG
0170 
0171     g.pixmx = std::max(0, cp.pixmx);
0172     g.sx2 = toMicron(cp.sx2);
0173     g.sy1 = std::max(21, toMicron(cp.sy1));  // for some angles sy1 is very small
0174     g.sy2 = std::max(55, toMicron(cp.sy2));  // sometimes sy2 is smaller than others (due to angle?)
0175 
0176     //sample xerr as function of position
0177     // moduleOffsetX is the definition of TrackerTraits::xOffset,
0178     // needs to be calculated because for Phase2 the modules are not uniform
0179     float moduleOffsetX = -(0.5f * float(g.nRows) + TrackerTraits::bigPixXCorrection);
0180     auto const xoff = moduleOffsetX * buffer_->commonParams().thePitchX;
0181 
0182     for (int ix = 0; ix < pixelCPEforDevice::kNumErrorBins; ++ix) {
0183       auto x = xoff * (1.f - (0.5f + float(ix)) / 8.f);
0184       auto gvx = p.theOrigin.x() - x;
0185       auto gvy = p.theOrigin.y();
0186       auto gvz = 1.f / p.theOrigin.z();
0187       cp.cotbeta = gvy * gvz;
0188       cp.cotalpha = gvx * gvz;
0189       errorFromTemplates(p, cp, 20000.f);
0190       g.sigmax[ix] = toMicron(cp.sigmax);
0191       g.sigmax1[ix] = toMicron(cp.sx1);
0192       LogDebug("PixelCPEFastParamsHost") << "sigmax vs x " << i << ' ' << x << ' ' << cp.cotalpha << ' '
0193                                          << int(g.sigmax[ix]) << ' ' << int(g.sigmax1[ix]) << ' ' << 10000.f * cp.sigmay
0194                                          << std::endl;
0195     }
0196 #ifdef EDM_ML_DEBUG
0197     // sample yerr as function of position
0198     // moduleOffsetY is the definition of TrackerTraits::yOffset (removed)
0199     float moduleOffsetY = 0.5f * float(g.nCols) + TrackerTraits::bigPixYCorrection;
0200     auto const yoff = -moduleOffsetY * buffer_->commonParams().thePitchY;
0201 
0202     for (int ix = 0; ix < pixelCPEforDevice::kNumErrorBins; ++ix) {
0203       auto y = yoff * (1.f - (0.5f + float(ix)) / 8.f);
0204       auto gvx = p.theOrigin.x() + 40.f * buffer_->commonParams().thePitchY;
0205       auto gvy = p.theOrigin.y() - y;
0206       auto gvz = 1.f / p.theOrigin.z();
0207       cp.cotbeta = gvy * gvz;
0208       cp.cotalpha = gvx * gvz;
0209       errorFromTemplates(p, cp, 20000.f);
0210       LogDebug("PixelCPEFastParamsHost") << "sigmay vs y " << i << ' ' << y << ' ' << cp.cotbeta << ' '
0211                                          << 10000.f * cp.sigmay << std::endl;
0212     }
0213 #endif  // EDM_ML_DEBUG
0214 
0215     // calculate angles (repeated)
0216     cp.cotalpha = gvx * gvz;
0217     cp.cotbeta = gvy * gvz;
0218     auto aveCB = cp.cotbeta;
0219 
0220     // sample x by charge
0221     int qbin = pixelCPEforDevice::kGenErrorQBins;  // low charge
0222     int k = 0;
0223     for (int qclus = 1000; qclus < 200000; qclus += 1000) {
0224       errorFromTemplates(p, cp, qclus);
0225       if (cp.qBin_ == qbin)
0226         continue;
0227       qbin = cp.qBin_;
0228       g.xfact[k] = cp.sigmax;
0229       g.yfact[k] = cp.sigmay;
0230       g.minCh[k++] = qclus;
0231 #ifdef EDM_ML_DEBUG
0232       LogDebug("PixelCPEFastParamsHost") << i << ' ' << g.rawId << ' ' << cp.cotalpha << ' ' << qclus << ' ' << cp.qBin_
0233                                          << ' ' << cp.pixmx << ' ' << m * cp.sigmax << ' ' << m * cp.sx1 << ' '
0234                                          << m * cp.sx2 << ' ' << m * cp.sigmay << ' ' << m * cp.sy1 << ' ' << m * cp.sy2
0235                                          << std::endl;
0236 #endif  // EDM_ML_DEBUG
0237     }
0238 
0239     assert(k <= pixelCPEforDevice::kGenErrorQBins);
0240 
0241     // fill the rest  (sometimes bin 4 is missing)
0242     for (int kk = k; kk < pixelCPEforDevice::kGenErrorQBins; ++kk) {
0243       g.xfact[kk] = g.xfact[k - 1];
0244       g.yfact[kk] = g.yfact[k - 1];
0245       g.minCh[kk] = g.minCh[k - 1];
0246     }
0247     auto detx = 1.f / g.xfact[0];
0248     auto dety = 1.f / g.yfact[0];
0249     for (int kk = 0; kk < pixelCPEforDevice::kGenErrorQBins; ++kk) {
0250       g.xfact[kk] *= detx;
0251       g.yfact[kk] *= dety;
0252     }
0253     // sample y in "angle"  (estimated from cluster size)
0254     float ys = 8.f - 4.f;  // apperent bias of half pixel (see plot)
0255     // plot: https://indico.cern.ch/event/934821/contributions/3974619/attachments/2091853/3515041/DigilessReco.pdf page 25
0256     // sample yerr as function of "size"
0257     for (int iy = 0; iy < pixelCPEforDevice::kNumErrorBins; ++iy) {
0258       ys += 1.f;  // first bin 0 is for size 9  (and size is in fixed point 2^3)
0259       if (pixelCPEforDevice::kNumErrorBins - 1 == iy)
0260         ys += 8.f;  // last bin for "overflow"
0261       // cp.cotalpha = ys*(buffer_->commonParams().thePitchX/(8.f*thickness));  //  use this to print sampling in "x"  (and comment the line below)
0262       cp.cotbeta = std::copysign(ys * (buffer_->commonParams().thePitchY / (8.f * thickness)), aveCB);
0263       errorFromTemplates(p, cp, 20000.f);
0264       g.sigmay[iy] = toMicron(cp.sigmay);
0265       LogDebug("PixelCPEFastParamsHost") << "sigmax/sigmay " << i << ' ' << (ys + 4.f) / 8.f << ' ' << cp.cotalpha
0266                                          << '/' << cp.cotbeta << ' ' << 10000.f * cp.sigmax << '/' << int(g.sigmay[iy])
0267                                          << std::endl;
0268     }
0269   }  // loop over det
0270 
0271   constexpr int numberOfModulesInLadder = TrackerTraits::numberOfModulesInLadder;
0272   constexpr int numberOfLaddersInBarrel = TrackerTraits::numberOfLaddersInBarrel;
0273   constexpr int numberOfModulesInBarrel = TrackerTraits::numberOfModulesInBarrel;
0274 
0275   constexpr float ladderFactor = 1.f / float(numberOfModulesInLadder);
0276 
0277   constexpr int firstEndcapPos = TrackerTraits::firstEndcapPos;
0278   constexpr int firstEndcapNeg = TrackerTraits::firstEndcapNeg;
0279 
0280   // compute ladder baricenter (only in global z) for the barrel
0281   //
0282   auto& aveGeom = buffer_->averageGeometry();
0283   int il = 0;
0284   for (int im = 0, nm = numberOfModulesInBarrel; im < nm; ++im) {
0285     auto const& g = buffer_->detParams(im);
0286     il = im / numberOfModulesInLadder;
0287     assert(il < int(numberOfLaddersInBarrel));
0288     auto z = g.frame.z();
0289     aveGeom.ladderZ[il] += ladderFactor * z;
0290     aveGeom.ladderMinZ[il] = std::min(aveGeom.ladderMinZ[il], z);
0291     aveGeom.ladderMaxZ[il] = std::max(aveGeom.ladderMaxZ[il], z);
0292     aveGeom.ladderX[il] += ladderFactor * g.frame.x();
0293     aveGeom.ladderY[il] += ladderFactor * g.frame.y();
0294     aveGeom.ladderR[il] += ladderFactor * sqrt(g.frame.x() * g.frame.x() + g.frame.y() * g.frame.y());
0295   }
0296   assert(il + 1 == int(numberOfLaddersInBarrel));
0297   // add half_module and tollerance
0298   constexpr float moduleLength = TrackerTraits::moduleLength;
0299   constexpr float module_tolerance = 0.2f;
0300   for (int il = 0, nl = numberOfLaddersInBarrel; il < nl; ++il) {
0301     aveGeom.ladderMinZ[il] -= (0.5f * moduleLength - module_tolerance);
0302     aveGeom.ladderMaxZ[il] += (0.5f * moduleLength - module_tolerance);
0303   }
0304 
0305   // compute "max z" for first layer in endcap (should we restrict to the outermost ring?)
0306   for (auto im = TrackerTraits::layerStart[firstEndcapPos]; im < TrackerTraits::layerStart[firstEndcapPos + 1]; ++im) {
0307     auto const& g = buffer_->detParams(im);
0308     aveGeom.endCapZ[0] = std::max(aveGeom.endCapZ[0], g.frame.z());
0309   }
0310   for (auto im = TrackerTraits::layerStart[firstEndcapNeg]; im < TrackerTraits::layerStart[firstEndcapNeg + 1]; ++im) {
0311     auto const& g = buffer_->detParams(im);
0312     aveGeom.endCapZ[1] = std::min(aveGeom.endCapZ[1], g.frame.z());
0313   }
0314   // correct for outer ring being closer
0315   aveGeom.endCapZ[0] -= TrackerTraits::endcapCorrection;
0316   aveGeom.endCapZ[1] += TrackerTraits::endcapCorrection;
0317 #ifdef EDM_ML_DEBUG
0318   for (int jl = 0, nl = numberOfLaddersInBarrel; jl < nl; ++jl) {
0319     LogDebug("PixelCPEFastParamsHost") << jl << ':' << aveGeom.ladderR[jl] << '/'
0320                                        << std::sqrt(aveGeom.ladderX[jl] * aveGeom.ladderX[jl] +
0321                                                     aveGeom.ladderY[jl] * aveGeom.ladderY[jl])
0322                                        << ',' << aveGeom.ladderZ[jl] << ',' << aveGeom.ladderMinZ[jl] << ','
0323                                        << aveGeom.ladderMaxZ[jl] << '\n';
0324   }
0325   LogDebug("PixelCPEFastParamsHost") << aveGeom.endCapZ[0] << ' ' << aveGeom.endCapZ[1];
0326 #endif  // EDM_ML_DEBUG
0327 
0328   // fill Layer and ladders geometry
0329   memset(&buffer_->layerGeometry(), 0, sizeof(pixelCPEforDevice::LayerGeometryT<TrackerTraits>));
0330   memcpy(buffer_->layerGeometry().layerStart,
0331          TrackerTraits::layerStart,
0332          sizeof(pixelCPEforDevice::LayerGeometryT<TrackerTraits>::layerStart));
0333   memcpy(buffer_->layerGeometry().layer,
0334          pixelTopology::layer<TrackerTraits>.data(),
0335          pixelTopology::layer<TrackerTraits>.size());
0336   buffer_->layerGeometry().maxModuleStride = pixelTopology::maxModuleStride<TrackerTraits>;
0337 }
0338 
0339 template <typename TrackerTraits>
0340 void PixelCPEFastParamsHost<TrackerTraits>::errorFromTemplates(DetParam const& theDetParam,
0341                                                                ClusterParamGeneric& theClusterParam,
0342                                                                float qclus) const {
0343   float locBz = theDetParam.bz;
0344   float locBx = theDetParam.bx;
0345   LogDebug("PixelCPEFastParamsHost") << "PixelCPEFastParamsHost::localPosition(...) : locBz = " << locBz;
0346 
0347   theClusterParam.pixmx = std::numeric_limits<int>::max();  // max pixel charge for truncation of 2-D cluster
0348 
0349   theClusterParam.sigmay = -999.9;  // CPE Generic y-error for multi-pixel cluster
0350   theClusterParam.sigmax = -999.9;  // CPE Generic x-error for multi-pixel cluster
0351   theClusterParam.sy1 = -999.9;     // CPE Generic y-error for single single-pixel
0352   theClusterParam.sy2 = -999.9;     // CPE Generic y-error for single double-pixel cluster
0353   theClusterParam.sx1 = -999.9;     // CPE Generic x-error for single single-pixel cluster
0354   theClusterParam.sx2 = -999.9;     // CPE Generic x-error for single double-pixel cluster
0355 
0356   float dummy;
0357 
0358   SiPixelGenError gtempl(this->thePixelGenError_);
0359   int gtemplID = theDetParam.detTemplateId;
0360 
0361   theClusterParam.qBin_ = gtempl.qbin(gtemplID,
0362                                       theClusterParam.cotalpha,
0363                                       theClusterParam.cotbeta,
0364                                       locBz,
0365                                       locBx,
0366                                       qclus,
0367                                       false,
0368                                       theClusterParam.pixmx,
0369                                       theClusterParam.sigmay,
0370                                       dummy,
0371                                       theClusterParam.sigmax,
0372                                       dummy,
0373                                       theClusterParam.sy1,
0374                                       dummy,
0375                                       theClusterParam.sy2,
0376                                       dummy,
0377                                       theClusterParam.sx1,
0378                                       dummy,
0379                                       theClusterParam.sx2,
0380                                       dummy);
0381 
0382   theClusterParam.sigmax = theClusterParam.sigmax * pixelCPEforDevice::micronsToCm;
0383   theClusterParam.sx1 = theClusterParam.sx1 * pixelCPEforDevice::micronsToCm;
0384   theClusterParam.sx2 = theClusterParam.sx2 * pixelCPEforDevice::micronsToCm;
0385 
0386   theClusterParam.sigmay = theClusterParam.sigmay * pixelCPEforDevice::micronsToCm;
0387   theClusterParam.sy1 = theClusterParam.sy1 * pixelCPEforDevice::micronsToCm;
0388   theClusterParam.sy2 = theClusterParam.sy2 * pixelCPEforDevice::micronsToCm;
0389 }
0390 
0391 template <>
0392 void PixelCPEFastParamsHost<pixelTopology::Phase2>::errorFromTemplates(DetParam const& theDetParam,
0393                                                                        ClusterParamGeneric& theClusterParam,
0394                                                                        float qclus) const {
0395   theClusterParam.qBin_ = 0.0f;
0396 }
0397 
0398 //-----------------------------------------------------------------------------
0399 //! Hit position in the local frame (in cm).  Unlike other CPE's, this
0400 //! one converts everything from the measurement frame (in channel numbers)
0401 //! into the local frame (in centimeters).
0402 //-----------------------------------------------------------------------------
0403 template <typename TrackerTraits>
0404 LocalPoint PixelCPEFastParamsHost<TrackerTraits>::localPosition(DetParam const& theDetParam,
0405                                                                 ClusterParam& theClusterParamBase) const {
0406   ClusterParamGeneric& theClusterParam = static_cast<ClusterParamGeneric&>(theClusterParamBase);
0407 
0408   if (useErrorsFromTemplates_) {
0409     errorFromTemplates(theDetParam, theClusterParam, theClusterParam.theCluster->charge());
0410   } else {
0411     theClusterParam.qBin_ = 0;
0412   }
0413 
0414   int q_f_X;  //!< Q of the first  pixel  in X
0415   int q_l_X;  //!< Q of the last   pixel  in X
0416   int q_f_Y;  //!< Q of the first  pixel  in Y
0417   int q_l_Y;  //!< Q of the last   pixel  in Y
0418   collect_edge_charges(theClusterParam, q_f_X, q_l_X, q_f_Y, q_l_Y, useErrorsFromTemplates_ && truncatePixelCharge_);
0419 
0420   // do GPU like ...
0421   pixelCPEforDevice::ClusParams cp;
0422 
0423   cp.minRow[0] = theClusterParam.theCluster->minPixelRow();
0424   cp.maxRow[0] = theClusterParam.theCluster->maxPixelRow();
0425   cp.minCol[0] = theClusterParam.theCluster->minPixelCol();
0426   cp.maxCol[0] = theClusterParam.theCluster->maxPixelCol();
0427 
0428   cp.q_f_X[0] = q_f_X;
0429   cp.q_l_X[0] = q_l_X;
0430   cp.q_f_Y[0] = q_f_Y;
0431   cp.q_l_Y[0] = q_l_Y;
0432 
0433   cp.charge[0] = theClusterParam.theCluster->charge();
0434 
0435   auto ind = theDetParam.theDet->index();
0436   pixelCPEforDevice::position<TrackerTraits>(buffer_->commonParams(), buffer_->detParams(ind), cp, 0);
0437   auto xPos = cp.xpos[0];
0438   auto yPos = cp.ypos[0];
0439 
0440   // set the error  (mind ape....)
0441   pixelCPEforDevice::errorFromDB<TrackerTraits>(buffer_->commonParams(), buffer_->detParams(ind), cp, 0);
0442   theClusterParam.sigmax = cp.xerr[0];
0443   theClusterParam.sigmay = cp.yerr[0];
0444 
0445   LogDebug("PixelCPEFastParamsHost") << " in PixelCPEFastParamsHost:localPosition - pos = " << xPos << " " << yPos
0446                                      << " size " << cp.maxRow[0] - cp.minRow[0] << ' ' << cp.maxCol[0] - cp.minCol[0];
0447 
0448   //--- Now put the two together
0449   LocalPoint pos_in_local(xPos, yPos);
0450   return pos_in_local;
0451 }
0452 
0453 //==============  INFLATED ERROR AND ERRORS FROM DB BELOW  ================
0454 
0455 //-------------------------------------------------------------------------
0456 //  Hit error in the local frame
0457 //-------------------------------------------------------------------------
0458 template <typename TrackerTraits>
0459 LocalError PixelCPEFastParamsHost<TrackerTraits>::localError(DetParam const& theDetParam,
0460                                                              ClusterParam& theClusterParamBase) const {
0461   ClusterParamGeneric& theClusterParam = static_cast<ClusterParamGeneric&>(theClusterParamBase);
0462 
0463   auto xerr = theClusterParam.sigmax;
0464   auto yerr = theClusterParam.sigmay;
0465 
0466   LogDebug("PixelCPEFastParamsHost") << " errors  " << xerr << " " << yerr;
0467 
0468   auto xerr_sq = xerr * xerr;
0469   auto yerr_sq = yerr * yerr;
0470 
0471   return LocalError(xerr_sq, 0, yerr_sq);
0472 }
0473 
0474 template <typename TrackerTraits>
0475 void PixelCPEFastParamsHost<TrackerTraits>::fillPSetDescription(edm::ParameterSetDescription& desc) {
0476   // call PixelCPEGenericBase fillPSetDescription to add common rechit errors
0477   PixelCPEGenericBase::fillPSetDescription(desc);
0478 }
0479 
0480 template class PixelCPEFastParamsHost<pixelTopology::Phase1>;
0481 template class PixelCPEFastParamsHost<pixelTopology::Phase2>;
0482 template class PixelCPEFastParamsHost<pixelTopology::HIonPhase1>;