PixelSeeding/plugins/ThirdHitCorrection.cc

0001 #include "ThirdHitCorrection.h"
0002
0003 #include "TrackingTools/DetLayers/interface/BarrelDetLayer.h"
0004 #include "TrackingTools/DetLayers/interface/ForwardDetLayer.h"
0005 #include "RecoTracker/TkMSParametrization/interface/MultipleScatteringParametrisationMaker.h"
0006
0007 using namespace pixelrecoutilities;
0008
0009 namespace {
0010   template <class T>
0011   inline T sqr(T t) {
0012     return t * t;
0013   }
0014 }  // namespace
0015
0016 using pixelrecoutilities::LongitudinalBendingCorrection;
0017
0018 void ThirdHitCorrection::init(float pt,
0019                               const DetLayer &layer3,
0020                               bool useMultipleScattering,
0021                               const MultipleScatteringParametrisationMaker *msmaker,
0022                               bool useBendingCorrection,
0023                               const MagneticField *bfield) {
0024   theUseMultipleScattering = useMultipleScattering;
0025   theUseBendingCorrection = useBendingCorrection;
0026   if (useBendingCorrection) {
0027     theBendingCorrection.init(pt, *bfield);
0028   }
0029
0030   theMultScattCorrRPhi = 0;
0031   theMScoeff = 0;
0032
0033   theBarrel = layer3.isBarrel();
0034   thePt = pt;
0035
0036   if (theUseMultipleScattering)
0037     sigmaRPhi = msmaker->parametrisation(&layer3);
0038 }
0039
0040 void ThirdHitCorrection::init(float pt,
0041                               const DetLayer &layer1,
0042                               const DetLayer &layer2,
0043                               const DetLayer &layer3,
0044                               bool useMultipleScattering,
0045                               const MultipleScatteringParametrisationMaker *msmaker,
0046                               bool useBendingCorrection,
0047                               const MagneticField *bfield) {
0048   init(pt, layer3, useMultipleScattering, msmaker, useBendingCorrection, bfield);
0049
0050   if (!theUseMultipleScattering)
0051     return;
0052
0053   auto point3 = [&]() -> PixelRecoPointRZ {
0054     if (theBarrel) {
0055       const BarrelDetLayer &bl = static_cast<const BarrelDetLayer &>(layer3);
0056       float rLayer = bl.specificSurface().radius();
0057       auto zmax = 0.5f * layer3.surface().bounds().length();
0058       return PixelRecoPointRZ(rLayer, zmax);
0059     } else {
0060       const ForwardDetLayer &fl = static_cast<const ForwardDetLayer &>(layer3);
0061       auto maxR = fl.specificSurface().outerRadius();
0062       auto layerZ = layer3.position().z();
0063       return PixelRecoPointRZ(maxR, layerZ);
0064     }
0065   };
0066
0067   PixelRecoPointRZ zero(0., 0.);
0068   SimpleLineRZ line(zero, point3());
0069
0070   auto point2 = [&]() -> PixelRecoPointRZ {
0071     if (layer2.isBarrel()) {
0072       const BarrelDetLayer &bl = static_cast<const BarrelDetLayer &>(layer2);
0073       float rLayer = bl.specificSurface().radius();
0074       return PixelRecoPointRZ(rLayer, line.zAtR(rLayer));
0075     } else {
0076       auto layerZ = layer2.position().z();
0077       return PixelRecoPointRZ(line.rAtZ(layerZ), layerZ);
0078     }
0079   };
0080
0081   theMultScattCorrRPhi = 3.f * (*sigmaRPhi)(pt, line.cotLine(), point2(), layer2.seqNum());
0082 }
0083
0084 void ThirdHitCorrection::init(const PixelRecoLineRZ &line, const PixelRecoPointRZ &constraint, int il) {
0085   theLine = line;
0086   if (!theUseMultipleScattering)
0087     return;
0088
0089   // auto newCorr = theMultScattCorrRPhi;
0090   theMultScattCorrRPhi = 3.f * (*sigmaRPhi)(thePt, line.cotLine(), constraint, il);
0091   // std::cout << "ThirdHitCorr " << (theBarrel ? "B " : "F " )<< theMultScattCorrRPhi << ' ' << newCorr << ' ' << newCorr/theMultScattCorrRPhi << std::endl;
0092   float overSinTheta = std::sqrt(1.f + sqr(line.cotLine()));
0093   if (theBarrel) {
0094     theMScoeff = theMultScattCorrRPhi * overSinTheta;
0095   } else {
0096     float overCosTheta = std::abs(line.cotLine()) < 1.e-4f ? 1.e4f : overSinTheta / std::abs(line.cotLine());
0097     theMScoeff = theMultScattCorrRPhi * overCosTheta;
0098   }
0099 }
0100
0101 void ThirdHitCorrection::correctRZRange(Range &range) const {
0102   range.first -= theMScoeff;
0103   range.second += theMScoeff;
0104
0105   if (theUseBendingCorrection) {
0106     if (theBarrel) {
0107       float cotTheta = theLine.cotLine();
0108       if (cotTheta > 0) {
0109         float radius = theLine.rAtZ(range.max());
0110         float corr = theBendingCorrection(radius) * cotTheta;
0111         range.second += corr;
0112       } else {
0113         float radius = theLine.rAtZ(range.min());
0114         float corr = theBendingCorrection(radius) * std::abs(cotTheta);
0115         range.first -= corr;
0116       }
0117     } else {
0118       float radius = range.max();
0119       float corr = theBendingCorrection(radius);
0120       range.first -= corr;
0121     }
0122   }
0123 }