src/smear/Tracker.cxx

0001 /**
0002  \file
0003  Implementation of class Smear::Tracker.
0004
0005  \author    Will Foreman
0006  \date      2011-08-19
0007  \copyright 2011 Brookhaven National Lab
0008  */
0009
0010 #include "eicsmear/smear/Tracker.h"
0011
0012 #include <algorithm>
0013 #include <cmath>
0014 #include <limits>
0015 #include <list>
0016
0017 #include <TMath.h>
0018
0019 namespace {
0020
0021 // Functor for testing non-existent intersection points,
0022 struct NoIntersection {
0023   bool operator()(const TVector3& v) const {
0024     return TMath::IsNaN(v.z());
0025   }
0026 };
0027
0028 }  // anonymous namespace
0029
0030 namespace Smear {
0031
0032 Tracker::Tracker(double magneticField, double nRadiationLengths,
0033                  double resolution)
0034 : mFactor(720)  // Assume no vertex constraint.
0035 , mMagField(magneticField)
0036 , mNRadLengths(nRadiationLengths)
0037 , mSigmaRPhi(resolution) {
0038 }
0039
0040 Tracker::~Tracker() {
0041 }
0042
0043 double Tracker::MultipleScatteringContribution(
0044     const erhic::VirtualParticle& p) const {
0045   // Technically should be a factor of particle charge in the numerator
0046   // but this is effectively always one.
0047   double val = 0.016 / 0.3 * p.GetP() * sqrt(mNRadLengths) / LPrime(p) /
0048   p.Get4Vector().Beta() / mMagField;
0049   if (TMath::IsNaN(val)) {
0050     std::cerr << "MS nan!" << std::endl;
0051   }  // if
0052   return val;
0053 }
0054
0055 double Tracker::IntrinsicContribution(
0056     const erhic::VirtualParticle& p) const {
0057   // The factor
0058   //    sqrt(720 * N^3 / ((N-1)(N+1)(N+2)(N+3)))
0059   // is a more exact version of the factor
0060   //    sqrt(720 / (N+5))
0061   // The constant factor in the sqrt depends on whether there is
0062   // a vertex constraint assumed.
0063   double val = sqrt(mFactor * pow(NPoints(p), 3.)) / 0.3 * pow(p.GetP(), 2.)
0064   * mSigmaRPhi / mMagField / pow(LPrime(p), 2.)
0065   / sqrt((NPoints(p)-1)*(NPoints(p)+1)*(NPoints(p)+2)*(NPoints(p)+3));
0066   if (TMath::IsNaN(val)) {
0067     std::cerr << "Intrinsic nan!" << std::endl;
0068   }  // if
0069   return val;
0070 }
0071
0072 double Tracker::Resolution(const erhic::VirtualParticle& p) const {
0073   // Don't compute resolution for very small path length
0074   // as L in the denominator causes it to blow up.
0075   if (!Accepts(p)) {
0076     return 0.;
0077   }  // if
0078   return sqrt(pow(MultipleScatteringContribution(p), 2.) +
0079               pow(IntrinsicContribution(p), 2.));
0080 }
0081
0082 void Tracker::Smear(const erhic::VirtualParticle& pIn,
0083                     ParticleMCS& pOut) {
0084   if (Accepts(pIn) && Accept.Is(pIn)) {
0085     double y = GetVariable(pIn, kP);
0086     // Randomly generate a smeared value from the resolution
0087     // and set it in the smeared particle.
0088     pOut.SetVariable(Distribution.Generate(y, Resolution(pIn)), kP);
0089     // Ensure E, p are positive definite
0090     pOut.HandleBogusValues(kP);
0091     if (pOut.GetP() < 0.) {
0092       std::cerr << "p " << pOut.GetP() << std::endl;
0093     }  // if
0094     if (TMath::IsNaN(pOut.GetP())) {
0095       std::cerr << "p nan" << std::endl;
0096     }  // if
0097   }  // if
0098 }
0099
0100 void Tracker::SetVertexConstraint(bool constrain) {
0101   if (constrain) {
0102     mFactor = 320;
0103   } else {
0104     mFactor = 720;
0105   }  // if
0106 }
0107
0108 }  // namespace Smear