Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2025-01-21 01:40:03

0001 #include "RecoLocalTracker/SiStripRecHitConverter/interface/StripCPEfromTrackAngle.h"
0002 #include "Geometry/CommonTopologies/interface/StripTopology.h"
0003 #include "DataFormats/SiStripCluster/interface/SiStripClusterTools.h"
0004 
0005 #include "vdt/vdtMath.h"
0006 
0007 StripCPEfromTrackAngle::StripCPEfromTrackAngle(edm::ParameterSet& conf,
0008                                                const MagneticField& mag,
0009                                                const TrackerGeometry& geom,
0010                                                const SiStripLorentzAngle& lorentz,
0011                                                const SiStripBackPlaneCorrection& backPlaneCorrection,
0012                                                const SiStripConfObject& confObj,
0013                                                const SiStripLatency& latency)
0014     : StripCPE(conf, mag, geom, lorentz, backPlaneCorrection, confObj, latency),
0015       useLegacyError(conf.getParameter<bool>("useLegacyError")),
0016       maxChgOneMIP(conf.getParameter<double>("maxChgOneMIP")),
0017       m_algo(useLegacyError ? Algo::legacy : (maxChgOneMIP < 0 ? Algo::mergeCK : Algo::chargeCK)) {
0018   mLC_P[0] = conf.getParameter<double>("mLC_P0");
0019   mLC_P[1] = conf.getParameter<double>("mLC_P1");
0020   mLC_P[2] = conf.getParameter<double>("mLC_P2");
0021 
0022   mHC_P[SiStripDetId::TIB - 3][0] = conf.getParameter<double>("mTIB_P0");
0023   mHC_P[SiStripDetId::TIB - 3][1] = conf.getParameter<double>("mTIB_P1");
0024   mHC_P[SiStripDetId::TID - 3][0] = conf.getParameter<double>("mTID_P0");
0025   mHC_P[SiStripDetId::TID - 3][1] = conf.getParameter<double>("mTID_P1");
0026   mHC_P[SiStripDetId::TOB - 3][0] = conf.getParameter<double>("mTOB_P0");
0027   mHC_P[SiStripDetId::TOB - 3][1] = conf.getParameter<double>("mTOB_P1");
0028   mHC_P[SiStripDetId::TEC - 3][0] = conf.getParameter<double>("mTEC_P0");
0029   mHC_P[SiStripDetId::TEC - 3][1] = conf.getParameter<double>("mTEC_P1");
0030 }
0031 
0032 void StripCPEfromTrackAngle::fillPSetDescription(edm::ParameterSetDescription& desc) {
0033   desc.add<bool>("useLegacyError", true);
0034   desc.add<double>("maxChgOneMIP", -6000.);
0035   desc.add<double>("mLC_P0", -.326);
0036   desc.add<double>("mLC_P1", .618);
0037   desc.add<double>("mLC_P2", .300);
0038   desc.add<double>("mTIB_P0", -.742);
0039   desc.add<double>("mTIB_P1", .202);
0040   desc.add<double>("mTID_P0", -1.026);
0041   desc.add<double>("mTID_P1", .253);
0042   desc.add<double>("mTOB_P0", -1.427);
0043   desc.add<double>("mTOB_P1", .433);
0044   desc.add<double>("mTEC_P0", -1.885);
0045   desc.add<double>("mTEC_P1", .471);
0046 }
0047 
0048 float StripCPEfromTrackAngle::stripErrorSquared(const unsigned N,
0049                                                 const float uProj,
0050                                                 const SiStripDetId::SubDetector loc) const {
0051   auto fun = [&](float x) -> float { return mLC_P[0] * x * vdt::fast_expf(-x * mLC_P[1]) + mLC_P[2]; };
0052   auto uerr = (N <= 4) ? fun(uProj) : mHC_P[loc - 3][0] + float(N) * mHC_P[loc - 3][1];
0053   return uerr * uerr;
0054 }
0055 
0056 float StripCPEfromTrackAngle::legacyStripErrorSquared(const unsigned N, const float uProj) const {
0057   if UNLIKELY ((float(N) - uProj) > 3.5f)
0058     return float(N * N) / 12.f;
0059   else {
0060     static constexpr float P1 = -0.339;
0061     static constexpr float P2 = 0.90;
0062     static constexpr float P3 = 0.279;
0063     const float uerr = P1 * uProj * vdt::fast_expf(-uProj * P2) + P3;
0064     return uerr * uerr;
0065   }
0066 }
0067 
0068 void StripCPEfromTrackAngle::localParameters(AClusters const& clusters,
0069                                              ALocalValues& retValues,
0070                                              const GeomDetUnit& det,
0071                                              const LocalTrajectoryParameters& ltp) const {
0072   auto const& par = getAlgoParam(det, ltp);
0073   auto const& p = par.p;
0074   auto loc = par.loc;
0075   auto corr = par.corr;
0076   auto afp = par.afullProjection;
0077 
0078   auto fill = [&](unsigned int i, float uerr2) {
0079     const float strip = clusters[i]->barycenter() + corr;
0080     retValues[i].first = p.topology->localPosition(strip, ltp.vector());
0081     retValues[i].second = p.topology->localError(strip, uerr2, ltp.vector());
0082   };
0083 
0084   switch (m_algo) {
0085     case Algo::chargeCK:
0086       for (auto i = 0U; i < clusters.size(); ++i) {
0087         auto dQdx = siStripClusterTools::chargePerCM(*clusters[i], ltp, p.invThickness);
0088         auto N = clusters[i]->amplitudes().size();
0089         auto uerr2 = dQdx > maxChgOneMIP ? legacyStripErrorSquared(N, afp) : stripErrorSquared(N, afp, loc);
0090         fill(i, uerr2);
0091       }
0092       break;
0093     case Algo::legacy:
0094       for (auto i = 0U; i < clusters.size(); ++i) {
0095         auto N = clusters[i]->amplitudes().size();
0096         auto uerr2 = legacyStripErrorSquared(N, afp);
0097         fill(i, uerr2);
0098       }
0099       break;
0100     case Algo::mergeCK:
0101       for (auto i = 0U; i < clusters.size(); ++i) {
0102         auto N = clusters[i]->amplitudes().size();
0103         auto uerr2 = clusters[i]->isMerged() ? legacyStripErrorSquared(N, afp) : stripErrorSquared(N, afp, loc);
0104         fill(i, uerr2);
0105       }
0106       break;
0107   }
0108 }
0109 
0110 StripClusterParameterEstimator::LocalValues StripCPEfromTrackAngle::localParameters(const SiStripCluster& cluster,
0111                                                                                     AlgoParam const& par) const {
0112   auto const& p = par.p;
0113   auto const& ltp = par.ltp;
0114   auto loc = par.loc;
0115   auto corr = par.corr;
0116   auto afp = par.afullProjection;
0117 
0118   float uerr2 = 0;
0119 
0120   auto N = cluster.amplitudes().size();
0121 
0122   switch (m_algo) {
0123     case Algo::chargeCK: {
0124       auto dQdx = siStripClusterTools::chargePerCM(cluster, ltp, p.invThickness);
0125       uerr2 = dQdx > maxChgOneMIP ? legacyStripErrorSquared(N, afp) : stripErrorSquared(N, afp, loc);
0126     } break;
0127     case Algo::legacy:
0128       uerr2 = legacyStripErrorSquared(N, afp);
0129       break;
0130     case Algo::mergeCK:
0131       uerr2 = cluster.isMerged() ? legacyStripErrorSquared(N, afp) : stripErrorSquared(N, afp, loc);
0132       break;
0133   }
0134 
0135   const float strip = cluster.barycenter() + corr;
0136 
0137   return std::make_pair(p.topology->localPosition(strip, ltp.vector()),
0138                         p.topology->localError(strip, uerr2, ltp.vector()));
0139 }
0140 
0141 StripClusterParameterEstimator::LocalValues StripCPEfromTrackAngle::localParameters(
0142     const SiStripCluster& cluster, const GeomDetUnit& det, const LocalTrajectoryParameters& ltp) const {
0143   auto const& par = getAlgoParam(det, ltp);
0144   return localParameters(cluster, par);
0145 }