File indexing completed on 2024-09-28 07:02:49
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010 #ifndef INCLUDE_EICSMEAR_SMEAR_SMEAR_H_
0011 #define INCLUDE_EICSMEAR_SMEAR_SMEAR_H_
0012
0013 #include <cmath>
0014 #include <fstream>
0015 #include <iostream>
0016 #include <memory>
0017 #include <sstream>
0018
0019 #include <TF1.h>
0020 #include <TF2.h>
0021 #include <TLorentzVector.h>
0022 #include <TMath.h>
0023 #include <TRandom3.h>
0024 #include <TROOT.h>
0025 #include <TString.h>
0026 #include <TSystem.h>
0027 #include <TVector2.h>
0028
0029 #include "eicsmear/erhic/Kinematics.h"
0030 #include "eicsmear/erhic/Particle.h"
0031 #include "eicsmear/erhic/ParticleIdentifier.h"
0032 #include "eicsmear/erhic/VirtualParticle.h"
0033 #include "eicsmear/smear/SmearConstants.h"
0034 #include "eicsmear/smear/ParticleMCS.h"
0035
0036 namespace Smear {
0037
0038
0039
0040
0041
0042
0043
0044
0045
0046
0047
0048
0049 inline int PGenre(const erhic::VirtualParticle& prt) {
0050 int genre(kAll);
0051 const int id = abs(prt.Id());
0052 if (1 == prt.GetStatus()) {
0053 if (id == 11 || id == 22) {
0054 genre = kElectromagnetic;
0055 } else if (id >110) {
0056 genre = kHadronic;
0057 }
0058 }
0059 return genre;
0060 }
0061
0062
0063
0064
0065
0066 inline double FixTheta(double theta) {
0067
0068 while (theta < 0. || theta > TMath::Pi()) {
0069 if (theta < 0.) {
0070 theta = -theta;
0071 }
0072 if (theta > TMath::Pi()) {
0073 theta = TMath::TwoPi() - theta;
0074 }
0075 }
0076 return theta;
0077 }
0078
0079
0080
0081
0082
0083 inline double FixPhi(double phi) {
0084 return TVector2::Phi_0_2pi(phi);
0085 }
0086
0087
0088
0089
0090 inline double GetVariable(const erhic::VirtualParticle& prt, KinType kin) {
0091 double z(0.);
0092 switch (kin) {
0093 case kE:
0094 z = prt.GetE(); break;
0095 case kP:
0096 z = prt.GetP(); break;
0097 case kTheta:
0098 z = prt.GetTheta(); break;
0099 case kPhi:
0100 z = prt.GetPhi(); break;
0101 case kPz:
0102 z = prt.GetPz(); break;
0103 case kPt:
0104 z = prt.GetPt(); break;
0105 default:
0106 break;
0107 }
0108 return z;
0109 }
0110
0111 inline bool IsCoreType(KinType kin) {
0112 if (kin == kE || kin == kP || kin == kTheta || kin == kPhi) return true;
0113 return false;
0114 }
0115
0116 int ParseInputFunction(TString &s, KinType &kin1, KinType &kin2);
0117
0118 }
0119
0120 #endif