Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-06-12 03:09:37

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.existsAs<bool>("useLegacyError") ? conf.getParameter<bool>("useLegacyError") : true),
0016       maxChgOneMIP(conf.existsAs<float>("maxChgOneMIP") ? conf.getParameter<double>("maxChgOneMIP") : -6000.),
0017       m_algo(useLegacyError ? Algo::legacy : (maxChgOneMIP < 0 ? Algo::mergeCK : Algo::chargeCK)) {
0018   mLC_P[0] = conf.existsAs<double>("mLC_P0") ? conf.getParameter<double>("mLC_P0") : -.326;
0019   mLC_P[1] = conf.existsAs<double>("mLC_P1") ? conf.getParameter<double>("mLC_P1") : .618;
0020   mLC_P[2] = conf.existsAs<double>("mLC_P2") ? conf.getParameter<double>("mLC_P2") : .300;
0021 
0022   mHC_P[SiStripDetId::TIB - 3][0] = conf.existsAs<double>("mTIB_P0") ? conf.getParameter<double>("mTIB_P0") : -.742;
0023   mHC_P[SiStripDetId::TIB - 3][1] = conf.existsAs<double>("mTIB_P1") ? conf.getParameter<double>("mTIB_P1") : .202;
0024   mHC_P[SiStripDetId::TID - 3][0] = conf.existsAs<double>("mTID_P0") ? conf.getParameter<double>("mTID_P0") : -1.026;
0025   mHC_P[SiStripDetId::TID - 3][1] = conf.existsAs<double>("mTID_P1") ? conf.getParameter<double>("mTID_P1") : .253;
0026   mHC_P[SiStripDetId::TOB - 3][0] = conf.existsAs<double>("mTOB_P0") ? conf.getParameter<double>("mTOB_P0") : -1.427;
0027   mHC_P[SiStripDetId::TOB - 3][1] = conf.existsAs<double>("mTOB_P1") ? conf.getParameter<double>("mTOB_P1") : .433;
0028   mHC_P[SiStripDetId::TEC - 3][0] = conf.existsAs<double>("mTEC_P0") ? conf.getParameter<double>("mTEC_P0") : -1.885;
0029   mHC_P[SiStripDetId::TEC - 3][1] = conf.existsAs<double>("mTEC_P1") ? conf.getParameter<double>("mTEC_P1") : .471;
0030 }
0031 
0032 float StripCPEfromTrackAngle::stripErrorSquared(const unsigned N,
0033                                                 const float uProj,
0034                                                 const SiStripDetId::SubDetector loc) const {
0035   auto fun = [&](float x) -> float { return mLC_P[0] * x * vdt::fast_expf(-x * mLC_P[1]) + mLC_P[2]; };
0036   auto uerr = (N <= 4) ? fun(uProj) : mHC_P[loc - 3][0] + float(N) * mHC_P[loc - 3][1];
0037   return uerr * uerr;
0038 }
0039 
0040 float StripCPEfromTrackAngle::legacyStripErrorSquared(const unsigned N, const float uProj) const {
0041   if UNLIKELY ((float(N) - uProj) > 3.5f)
0042     return float(N * N) / 12.f;
0043   else {
0044     static constexpr float P1 = -0.339;
0045     static constexpr float P2 = 0.90;
0046     static constexpr float P3 = 0.279;
0047     const float uerr = P1 * uProj * vdt::fast_expf(-uProj * P2) + P3;
0048     return uerr * uerr;
0049   }
0050 }
0051 
0052 void StripCPEfromTrackAngle::localParameters(AClusters const& clusters,
0053                                              ALocalValues& retValues,
0054                                              const GeomDetUnit& det,
0055                                              const LocalTrajectoryParameters& ltp) const {
0056   auto const& par = getAlgoParam(det, ltp);
0057   auto const& p = par.p;
0058   auto loc = par.loc;
0059   auto corr = par.corr;
0060   auto afp = par.afullProjection;
0061 
0062   auto fill = [&](unsigned int i, float uerr2) {
0063     const float strip = clusters[i]->barycenter() + corr;
0064     retValues[i].first = p.topology->localPosition(strip, ltp.vector());
0065     retValues[i].second = p.topology->localError(strip, uerr2, ltp.vector());
0066   };
0067 
0068   switch (m_algo) {
0069     case Algo::chargeCK:
0070       for (auto i = 0U; i < clusters.size(); ++i) {
0071         auto dQdx = siStripClusterTools::chargePerCM(*clusters[i], ltp, p.invThickness);
0072         auto N = clusters[i]->amplitudes().size();
0073         auto uerr2 = dQdx > maxChgOneMIP ? legacyStripErrorSquared(N, afp) : stripErrorSquared(N, afp, loc);
0074         fill(i, uerr2);
0075       }
0076       break;
0077     case Algo::legacy:
0078       for (auto i = 0U; i < clusters.size(); ++i) {
0079         auto N = clusters[i]->amplitudes().size();
0080         auto uerr2 = legacyStripErrorSquared(N, afp);
0081         fill(i, uerr2);
0082       }
0083       break;
0084     case Algo::mergeCK:
0085       for (auto i = 0U; i < clusters.size(); ++i) {
0086         auto N = clusters[i]->amplitudes().size();
0087         auto uerr2 = clusters[i]->isMerged() ? legacyStripErrorSquared(N, afp) : stripErrorSquared(N, afp, loc);
0088         fill(i, uerr2);
0089       }
0090       break;
0091   }
0092 }
0093 
0094 StripClusterParameterEstimator::LocalValues StripCPEfromTrackAngle::localParameters(const SiStripCluster& cluster,
0095                                                                                     AlgoParam const& par) const {
0096   auto const& p = par.p;
0097   auto const& ltp = par.ltp;
0098   auto loc = par.loc;
0099   auto corr = par.corr;
0100   auto afp = par.afullProjection;
0101 
0102   float uerr2 = 0;
0103 
0104   auto N = cluster.amplitudes().size();
0105 
0106   switch (m_algo) {
0107     case Algo::chargeCK: {
0108       auto dQdx = siStripClusterTools::chargePerCM(cluster, ltp, p.invThickness);
0109       uerr2 = dQdx > maxChgOneMIP ? legacyStripErrorSquared(N, afp) : stripErrorSquared(N, afp, loc);
0110     } break;
0111     case Algo::legacy:
0112       uerr2 = legacyStripErrorSquared(N, afp);
0113       break;
0114     case Algo::mergeCK:
0115       uerr2 = cluster.isMerged() ? legacyStripErrorSquared(N, afp) : stripErrorSquared(N, afp, loc);
0116       break;
0117   }
0118 
0119   const float strip = cluster.barycenter() + corr;
0120 
0121   return std::make_pair(p.topology->localPosition(strip, ltp.vector()),
0122                         p.topology->localError(strip, uerr2, ltp.vector()));
0123 }
0124 
0125 StripClusterParameterEstimator::LocalValues StripCPEfromTrackAngle::localParameters(
0126     const SiStripCluster& cluster, const GeomDetUnit& det, const LocalTrajectoryParameters& ltp) const {
0127   auto const& par = getAlgoParam(det, ltp);
0128   return localParameters(cluster, par);
0129 }