#include <MultipleScatteringUpdator.h>

Inheritance diagram for MultipleScatteringUpdator:

Public Member Functions
void	compute (const TrajectoryStateOnSurface &, const PropagationDirection, Effect &effect) const override

	MultipleScatteringUpdator (float mass, float ptMin=-1.)

	~MultipleScatteringUpdator () override
	destructor More...

Public Member Functions inherited from MaterialEffectsUpdator
float	mass () const

	MaterialEffectsUpdator (float mass)

virtual TrajectoryStateOnSurface	updateState (const TrajectoryStateOnSurface &TSoS, const PropagationDirection propDir) const

virtual bool	updateStateInPlace (TrajectoryStateOnSurface &TSoS, const PropagationDirection propDir) const

virtual	~MaterialEffectsUpdator ()

Private Member Functions
MultipleScatteringUpdator *	clone () const override

Private Attributes
float	thePtMin

Additional Inherited Members
Public Types inherited from MaterialEffectsUpdator
typedef materialEffect::Covariance	Covariance

typedef materialEffect::CovIndex	CovIndex

typedef materialEffect::Effect	Effect

Detailed Description

Adds effects from multiple scattering (standard Highland formula) to a trajectory state. Uses radiation length from medium properties. Ported from ORCA.

Author: todorov, cerati

Definition at line 15 of file MultipleScatteringUpdator.h.

Constructor & Destructor Documentation

◆ MultipleScatteringUpdator()

MultipleScatteringUpdator::MultipleScatteringUpdator	(	float	mass,
		float	ptMin = `-1.`
	)

inline

Specify assumed mass of particle for material effects. If ptMin > 0, then the rms muliple scattering angle will be calculated taking into account the uncertainty in the reconstructed track momentum. (By default, it is neglected). However, a lower limit on the possible value of the track Pt will be applied at ptMin, to avoid the rms multiple scattering becoming too big.

Definition at line 23 of file MultipleScatteringUpdator.h.

23 : MaterialEffectsUpdator(mass), thePtMin(ptMin) {}

Referenced by clone().

◆ ~MultipleScatteringUpdator()

MultipleScatteringUpdator::~MultipleScatteringUpdator ( )

inlineoverride

destructor

Definition at line 25 of file MultipleScatteringUpdator.h.

25 {}

Member Function Documentation

◆ clone()

MultipleScatteringUpdator* MultipleScatteringUpdator::clone ( void ) const

inlineoverrideprivatevirtual

Implements MaterialEffectsUpdator.

Definition at line 16 of file MultipleScatteringUpdator.h.

16 { return new MultipleScatteringUpdator(*this); }

References MultipleScatteringUpdator().

◆ compute()

void MultipleScatteringUpdator::compute	(	const TrajectoryStateOnSurface &	TSoS,
		const PropagationDirection	propDir,
		Effect &	effect
	)		const

overridevirtual

Implements MaterialEffectsUpdator.

Definition at line 14 of file MultipleScatteringUpdator.cc.

                                                               {
   //
   // Get surface
   //
   const Surface& surface = TSoS.surface();
   //
   //
   // Now get information on medium
   //
   const MediumProperties& mp = surface.mediumProperties();
   if
     UNLIKELY(mp.radLen() == 0) return;
  
   // Momentum vector
   LocalVector d = TSoS.localMomentum();
   float p2 = d.mag2();
   d *= 1.f / sqrt(p2);
   float xf = 1.f / std::abs(d.z());  // increase of path due to angle of incidence
   // calculate general physics things
   constexpr float amscon = 1.8496e-4;  // (13.6MeV)**2
   const float m2 = mass() * mass();    // use mass hypothesis from constructor
   float e2 = p2 + m2;
   float beta2 = p2 / e2;
   // calculate the multiple scattering angle
   float radLen = mp.radLen() * xf;  // effective rad. length
   float sigt2 = 0.;                 // sigma(alpha)**2
  
   // Calculated rms scattering angle squared.
   float fact = 1.f + 0.038f * unsafe_logf<2>(radLen);
   fact *= fact;
   float a = fact / (beta2 * p2);
   sigt2 = amscon * radLen * a;
  
   if (thePtMin > 0) {
 #ifdef DBG_MSU
     std::cout << "Original rms scattering = " << sqrt(sigt2);
 #endif
     // Inflate estimated rms scattering angle, to take into account
     // that 1/p is not known precisely.
     AlgebraicSymMatrix55 const& covMatrix = TSoS.localError().matrix();
     float error2_QoverP = covMatrix(0, 0);
     // Formula valid for ultra-relativistic particles.
     //      sigt2 *= (1. + p2 * error2_QoverP);
     // Exact formula
     sigt2 *= 1.f + error2_QoverP * (p2 + m2 * beta2 * (5.f + 3.f * p2 * error2_QoverP));
 #ifdef DBG_MSU
     std::cout << " new = " << sqrt(sigt2);
 #endif
     // Convert Pt constraint to P constraint, neglecting uncertainty in
     // track angle.
     float pMin2 = thePtMin * thePtMin * (p2 / TSoS.globalMomentum().perp2());
     // Use P constraint to calculate rms maximum scattering angle.
     float betaMin2 = pMin2 / (pMin2 + m2);
     float a_max = fact / (betaMin2 * pMin2);
     float sigt2_max = amscon * radLen * a_max;
     if (sigt2 > sigt2_max)
       sigt2 = sigt2_max;
 #ifdef DBG_MSU
     std::cout << " after P constraint (" << pMin << ") = " << sqrt(sigt2);
     std::cout << " for track with 1/p=" << 1 / p << "+-" << sqrt(error2_QoverP) << std::endl;
 #endif
   }
  
   float isl2 = 1.f / d.perp2();
   float cl2 = (d.z() * d.z());
   float cf2 = (d.x() * d.x()) * isl2;
   float sf2 = (d.y() * d.y()) * isl2;
   // Create update (transformation of independant variations
   //   on angle in orthogonal planes to local parameters.
   float den = 1.f / (cl2 * cl2);
   using namespace materialEffect;
   effect.deltaCov[msxx] += (den * sigt2) * (sf2 * cl2 + cf2);
   effect.deltaCov[msxy] += (den * sigt2) * (d.x() * d.y());
   effect.deltaCov[msyy] += (den * sigt2) * (cf2 * cl2 + sf2);
  
   /*
     std::cout << "new " <<  theDeltaCov(1,1) << " " <<  theDeltaCov(1,2)  << " " <<  theDeltaCov(2,2) << std::endl;
     oldMUcompute(TSoS,propDir, mass(), thePtMin);
   */
 }