eic-smear  1.0.3
A collection of ROOT classes for Monte Carlo events and a fast-smearing code simulating detector effects for the Electron-Ion Collider task force
Device.cxx
Go to the documentation of this file.
1 
10 #include "eicsmear/smear/Device.h"
11 
12 #include <algorithm>
13 #include <functional>
14 #include <iterator>
15 #include <stdexcept>
16 #include <string>
17 #include <vector>
18 
19 #include <TUUID.h>
20 
23 
24 namespace Smear {
25 
26 bool Device::Init(const TString& kinematicFunction,
27  const TString& resolutionFunction, int genre) {
28  Accept.SetGenre(genre);
29  // Use FormulaString to parse the kinematic function,
30  // though we can't use a FormulaString for the kinematic
31  // function as we need TF1 functionality.
32  FormulaString f(kinematicFunction.Data());
33  // The expression has to have exactly one variable.
34  mSmeared = f.Variables().front();
35  // Use UUID for ROOT function name to avoid instances clashing
36  mKinematicFunction = new TF1(TUUID().AsString(),
37  f.GetString().c_str(), 0., 1.e16);
38  // Set the resolution function.
39  mFormula = new FormulaString(resolutionFunction.Data());
40  return mFormula && mKinematicFunction;
41 }
42 
43 Device::Device(KinType type, const TString& formula, EGenre genre)
44 : mSmeared(type)
45 , mKinematicFunction(NULL)
46 , mFormula(NULL) {
47  Accept.SetGenre(genre);
48  Init(FormulaString::GetKinName(type), formula, genre);
49 }
50 
51 Device::Device(const TString& variable, const TString& resolution,
52  EGenre genre)
53 : mSmeared(kInvalidKinType)
54 , mKinematicFunction(NULL)
55 , mFormula(NULL) {
56  Init(variable, resolution, genre);
57 }
58 
59 Device::Device(const Device& that)
60 : Smearer(that)
61 , mSmeared(that.mSmeared)
62 , mKinematicFunction(NULL)
63 , mFormula(NULL)
64 , mDimensions(that.mDimensions) {
65  if (that.mKinematicFunction) {
66  mKinematicFunction = static_cast<TF1*>(
67  that.mKinematicFunction->Clone(TUUID().AsString()));
68  } // if
69  if (that.mFormula) {
70  mFormula = static_cast<FormulaString*>(that.mFormula->Clone());
71  } // if
72 }
73 
75  if (mFormula) {
76  delete mFormula;
77  mFormula = NULL;
78  } // if
79  if (mKinematicFunction) {
80  delete mKinematicFunction;
81  mKinematicFunction = NULL;
82  } // if
83 }
84 
86  // Test for acceptance and do nothing if it fails.
87  if (!Accept.Is(prt)) {
88  return;
89  } // if
90  // Get each argument for the resolution function from the particle.
91  const std::vector<KinType> vars = mFormula->Variables();
92  std::vector<double> args;
93  for (unsigned i(0); i < vars.size(); ++i) {
94  args.push_back(GetVariable(prt, vars.at(i)));
95  } // for
96  // Evaluate the quantity to smear and the resolution, then throw
97  // a random smeared value.
98  double unsmeared = mKinematicFunction->Eval(GetVariable(prt, mSmeared));
99  double resolution = mFormula->Eval(args);
100  double smeared = mDistribution.Generate(unsmeared, resolution);
101  SetVariable(out, mKinematicFunction->GetX(smeared), mSmeared);
102  // Fix angles to the correct ranges.
103  if (kTheta == mSmeared) {
104  out.theta = FixTheta(out.theta);
105  } else if (kPhi == mSmeared) {
106  out.phi = FixPhi(out.phi);
107  } // else if
108  // Ensure E, p are positive definite
110 }
111 
112 Device* Device::Clone(const char* ) const {
113  return new Device(*this);
114 }
115 
116 void Device::Print(Option_t* /* option */) const {
117  const std::string name = FormulaString::GetKinName(mSmeared);
118  std::cout << "Device smearing " << name << " with sigma(" << name <<
119  ") = " << mFormula->GetInputString() << std::endl;
120 }
121 
122 } // namespace Smear
Smear::FormulaString::Eval
virtual double Eval(const std::vector< double > &) const
Evaluate the formula with the provided arguments.
Definition: FormulaString.cxx:86
Smear::Device::Print
virtual void Print(Option_t *="") const
Print information about this device to standard output.
Definition: Device.cxx:116
Smear::Device::mSmeared
KinType mSmeared
Smeared variable.
Definition: Device.h:118
Smear
Definition: Acceptance.cxx:16
Smear::Acceptance::Is
bool Is(const erhic::VirtualParticle &prt) const
This function determines if the particle provided lies within the acceptance of the detector.
Definition: Acceptance.cxx:46
Smear::FixPhi
double FixPhi(double phi)
Fix an azimuthal angle so that it lies within [0,2*pi).
Definition: Smear.h:101
Smear::HandleBogusValues
void HandleBogusValues(ParticleMCS &prt, KinType kin)
This dictates how the namespace deals with positive definite variables which have been smeared to neg...
Definition: Smear.h:155
ParticleMCS.h
Smear::ParticleMCS::phi
Double32_t phi
Azimuthal angle.
Definition: ParticleMCS.h:180
Smear::Device::~Device
virtual ~Device()
Destructor.
Definition: Device.cxx:74
Smear::Device::mFormula
FormulaString * mFormula
Expression for resolution standard deviation.
Definition: Device.h:120
Smear::Smearer::Accept
Acceptance Accept
Definition: Smearer.h:50
Smear::FormulaString::GetString
virtual std::string GetString() const
Returns the processed formula string with variables substituted.
Definition: FormulaString.cxx:133
Smear::kTheta
@ kTheta
Definition: Smear.h:44
Smear::FormulaString::GetInputString
virtual std::string GetInputString() const
Returns the unprocessed input formula string.
Definition: FormulaString.cxx:141
Smear::SetVariable
void SetVariable(ParticleMCS &prt, double z, KinType kin)
Stores z in the ParticleS.K where K is the kinematic variable associated with kin.
Definition: Smear.h:132
Smear::Device::Device
Device(KinType=kE, const TString &formula="0", EGenre=kAll)
Constructor.
Definition: Device.cxx:43
Smear::Acceptance::SetGenre
void SetGenre(int genre)
Select the class(es) of particles to accept.
Definition: Acceptance.cxx:30
Smear::FixTheta
double FixTheta(double theta)
Fix a polar angle so that it lies within [0,pi].
Definition: Smear.h:85
erhic::VirtualParticle
Abstract base class for a general particle.
Definition: VirtualParticle.h:23
Smear::FormulaString
A formula described by a string with up to four variables.
Definition: FormulaString.h:27
Smear::Device::mDistribution
Distributor mDistribution
Random distribution.
Definition: Device.h:123
Device.h
Smear::kPhi
@ kPhi
Definition: Smear.h:44
FormulaString.h
Smear::Smearer
Abstract base class for objects performing smearing.
Definition: Smearer.h:33
Smear::Device::Clone
virtual Device * Clone(const char *="") const
Returns a dynamically allocated copy of this object.
Definition: Device.cxx:112
Smear::FormulaString::Variables
virtual std::vector< Smear::KinType > Variables() const
Returns a vector of Smear::KinType corresponding to the variables named in the constructor string.
Definition: FormulaString.cxx:98
Smear::ParticleMCS
A smeared Monte Carlo particle.
Definition: ParticleMCS.h:27
Smear::ParticleMCS::theta
Double32_t theta
Polar angle.
Definition: ParticleMCS.h:179
Smear::Device
Performs smearing of a single kinematic variable according to a simple expression defined via a strin...
Definition: Device.h:44
Smear::kInvalidKinType
@ kInvalidKinType
Definition: Smear.h:44
Smear::Device::Init
bool Init(const TString &, const TString &, int)
Definition: Device.cxx:26
Smear::Device::Smear
virtual void Smear(const erhic::VirtualParticle &, ParticleMCS &)
Smear the kinematic value of the input particle and store the result in the ParticleMCS.
Definition: Device.cxx:85
Smear::KinType
KinType
Enumerator listing particle wise kinematic variables.
Definition: Smear.h:43
Smear::GetVariable
double GetVariable(const erhic::VirtualParticle &prt, KinType kin)
Returns the kinematic variable associated with kin from the input particle.
Definition: Smear.h:108
Smear::FormulaString::GetKinName
static std::string GetKinName(KinType)
Returns the name corresponding the a Smear::KinType.
Definition: FormulaString.cxx:153
Smear::EGenre
EGenre
Classes of particles.
Definition: Smear.h:48
Smear::Distributor::Generate
virtual double Generate(double midpoint, double width)
Generate a random value based on a given midpoint and width.
Definition: Distributor.cxx:45
Smear::Device::mKinematicFunction
TF1 * mKinematicFunction
Definition: Device.h:119