SND@LHC Software
Loading...
Searching...
No Matches
ShipMCTrack Class Reference

#include <ShipMCTrack.h>

Inheritance diagram for ShipMCTrack:
Collaboration diagram for ShipMCTrack:

Public Member Functions

 ShipMCTrack ()
 
 ShipMCTrack (Int_t pdgCode, Int_t motherID, Double_t px, Double_t py, Double_t pz, Double_t E, Double_t x, Double_t y, Double_t z, Double_t t, Int_t nPoints, Double_t w)
 
 ShipMCTrack (const ShipMCTrack &track)
 
 ShipMCTrack (TParticle *particle)
 
virtual ~ShipMCTrack ()
 
void Print (Int_t iTrack=0) const
 
Int_t GetPdgCode () const
 
Int_t GetMotherId () const
 
Double_t GetPx () const
 
Double_t GetPy () const
 
Double_t GetPz () const
 
Double_t GetStartX () const
 
Double_t GetStartY () const
 
Double_t GetStartZ () const
 
Double_t GetStartT () const
 
void SetProcID (Int_t i)
 
Int_t GetProcID () const
 
TString GetProcName () const
 
Double_t GetMass () const
 
Double_t GetEnergy () const
 
Double_t GetPt () const
 
Double_t GetP () const
 
Double_t GetRapidity () const
 
void MultiplyWeight (Double_t w)
 
void SetWeight (Double_t w)
 
Double_t GetWeight () const
 
void GetMomentum (TVector3 &momentum)
 
void Get4Momentum (TLorentzVector &momentum)
 
void GetStartVertex (TVector3 &vertex)
 
Int_t GetNPoints (DetectorId detId) const
 
void SetMotherId (Int_t id)
 
void SetNPoints (Int_t iDet, Int_t np)
 

Private Member Functions

 ClassDef (ShipMCTrack, 8)
 

Private Attributes

Int_t fPdgCode
 
Int_t fMotherId
 
Double32_t fPx
 
Double32_t fPy
 
Double32_t fPz
 
Double32_t fM
 
Double32_t fStartX
 
Double32_t fStartY
 
Double32_t fStartZ
 
Double32_t fStartT
 
Double32_t fW
 
Int_t fProcID
 
Int_t fNPoints
 

Detailed Description

Definition at line 26 of file ShipMCTrack.h.

Constructor & Destructor Documentation

◆ ShipMCTrack() [1/4]

ShipMCTrack::ShipMCTrack ( )

Default constructor

Definition at line 13 of file ShipMCTrack.cxx.

14 : TObject(),
15 fPdgCode(0),
16 fMotherId(-1),
17 fPx(0.),
18 fPy(0.),
19 fPz(0.),
20 fM(-1.),
21 fStartX(0.),
22 fStartY(0.),
23 fStartZ(0.),
24 fStartT(0.),
25 fNPoints(0),
26 fW(1.),
27 fProcID(44)
28{
29}
Double32_t fPy
Int_t fPdgCode
Definition ShipMCTrack.h:96
Double32_t fStartX
Double32_t fW
Double32_t fStartZ
Double32_t fM
Double32_t fStartY
Double32_t fStartT
Int_t fMotherId
Definition ShipMCTrack.h:99
Double32_t fPz
Double32_t fPx

◆ ShipMCTrack() [2/4]

ShipMCTrack::ShipMCTrack ( Int_t  pdgCode,
Int_t  motherID,
Double_t  px,
Double_t  py,
Double_t  pz,
Double_t  E,
Double_t  x,
Double_t  y,
Double_t  z,
Double_t  t,
Int_t  nPoints,
Double_t  w 
)

Standard constructor

Definition at line 35 of file ShipMCTrack.cxx.

38 : TObject(),
39 fPdgCode(pdgCode),
40 fMotherId(motherId),
41 fPx(px),
42 fPy(py),
43 fPz(pz),
44 fM(M),
45 fStartX(x),
46 fStartY(y),
47 fStartZ(z),
48 fStartT(t),
49 fNPoints(nPoints),
50 fW(w)
51{
52}

◆ ShipMCTrack() [3/4]

ShipMCTrack::ShipMCTrack ( const ShipMCTrack track)

Copy constructor

Definition at line 58 of file ShipMCTrack.cxx.

59 : TObject(track),
60 fPdgCode(track.fPdgCode),
61 fMotherId(track.fMotherId),
62 fPx(track.fPx),
63 fPy(track.fPy),
64 fPz(track.fPz),
65 fM(track.fM),
66 fStartX(track.fStartX),
67 fStartY(track.fStartY),
68 fStartZ(track.fStartZ),
69 fStartT(track.fStartT),
70 fNPoints(track.fNPoints),
71 fProcID(track.GetProcID()),
72 fW(track.GetWeight())
73{
74}

◆ ShipMCTrack() [4/4]

ShipMCTrack::ShipMCTrack ( TParticle *  particle)

Constructor from TParticle

Definition at line 80 of file ShipMCTrack.cxx.

81 : TObject(),
82 fPdgCode(part->GetPdgCode()),
83 fMotherId(part->GetMother(0)),
84 fPx(part->Px()),
85 fPy(part->Py()),
86 fPz(part->Pz()),
87 fM(TMath::Sqrt( part->Energy()*part->Energy()-part->P()*part->P() )),
88 fStartX(part->Vx()),
89 fStartY(part->Vy()),
90 fStartZ(part->Vz()),
91 fStartT(part->T()*1e09),
92 fNPoints(0),
93 fProcID(part->GetUniqueID()),
94 fW(part->GetWeight())
95{
96}

◆ ~ShipMCTrack()

ShipMCTrack::~ShipMCTrack ( )
virtual

Destructor

Definition at line 102 of file ShipMCTrack.cxx.

102{ }

Member Function Documentation

◆ ClassDef()

ShipMCTrack::ClassDef ( ShipMCTrack  ,
 
)
private

◆ Get4Momentum()

void ShipMCTrack::Get4Momentum ( TLorentzVector &  momentum)
inline

Definition at line 146 of file ShipMCTrack.h.

147{
148 momentum.SetXYZT(fPx,fPy,fPz,GetEnergy());
149}
Double_t GetEnergy() const

◆ GetEnergy()

Double_t ShipMCTrack::GetEnergy ( ) const

Definition at line 121 of file ShipMCTrack.cxx.

122{
123 if (fM<0){
124// older data, mass not made persistent
125 Double_t mass = GetMass();
126 return TMath::Sqrt(mass*mass + fPx*fPx + fPy*fPy + fPz*fPz );
127 }else{
128 if (fPdgCode==22){return TMath::Sqrt(fPx*fPx + fPy*fPy + fPz*fPz );}
129 return TMath::Sqrt(fM*fM + fPx*fPx + fPy*fPy + fPz*fPz );
130 }
131}
Double_t GetMass() const
mass(particle)
Definition hnl.py:47

◆ GetMass()

Double_t ShipMCTrack::GetMass ( ) const

Definition at line 133 of file ShipMCTrack.cxx.

134{
135 if (fM<0){
136// older data, mass not made persistent
137 if ( TDatabasePDG::Instance() ) {
138 TParticlePDG* particle = TDatabasePDG::Instance()->GetParticle(fPdgCode);
139 if ( particle ) { return particle->Mass(); }
140 else { return 0.; }
141 }
142 }
143 return fM;
144}

◆ GetMomentum()

void ShipMCTrack::GetMomentum ( TVector3 &  momentum)
inline

Definition at line 140 of file ShipMCTrack.h.

141{
142 momentum.SetXYZ(fPx,fPy,fPz);
143}

◆ GetMotherId()

Int_t ShipMCTrack::GetMotherId ( ) const
inline

Definition at line 59 of file ShipMCTrack.h.

59{ return fMotherId; }

◆ GetNPoints()

Int_t ShipMCTrack::GetNPoints ( DetectorId  detId) const

Accessors to the number of MCPoints in the detectors

Definition at line 169 of file ShipMCTrack.cxx.

170{
171/* // TODO: Where does this come from
172 if ( detId == kREF ) { return ( fNPoints & 1); }
173 else if ( detId == kTutDet ) { return ( (fNPoints & ( 7 << 1) ) >> 1); }
174 else if ( detId == kFairRutherford ) { return ( (fNPoints & (31 << 4) ) >> 4); }
175 else {
176 LOG(ERROR) << "Unknown detector ID " << detId ;
177 return 0;
178 }
179*/
180 return 0;
181}

◆ GetP()

Double_t ShipMCTrack::GetP ( ) const
inline

Definition at line 73 of file ShipMCTrack.h.

73{ return TMath::Sqrt(fPx*fPx+fPy*fPy+fPz*fPz); }

◆ GetPdgCode()

Int_t ShipMCTrack::GetPdgCode ( ) const
inline

Accessors

Definition at line 58 of file ShipMCTrack.h.

58{ return fPdgCode; }

◆ GetProcID()

Int_t ShipMCTrack::GetProcID ( ) const
inline

Definition at line 68 of file ShipMCTrack.h.

68{ return fProcID; }

◆ GetProcName()

TString ShipMCTrack::GetProcName ( ) const
inline

Definition at line 69 of file ShipMCTrack.h.

69{ return TMCProcessName[fProcID]; }

◆ GetPt()

Double_t ShipMCTrack::GetPt ( ) const
inline

Definition at line 72 of file ShipMCTrack.h.

72{ return TMath::Sqrt(fPx*fPx+fPy*fPy); }

◆ GetPx()

Double_t ShipMCTrack::GetPx ( ) const
inline

Definition at line 60 of file ShipMCTrack.h.

60{ return fPx; }

◆ GetPy()

Double_t ShipMCTrack::GetPy ( ) const
inline

Definition at line 61 of file ShipMCTrack.h.

61{ return fPy; }

◆ GetPz()

Double_t ShipMCTrack::GetPz ( ) const
inline

Definition at line 62 of file ShipMCTrack.h.

62{ return fPz; }

◆ GetRapidity()

Double_t ShipMCTrack::GetRapidity ( ) const

Definition at line 157 of file ShipMCTrack.cxx.

158{
159 Double_t e = GetEnergy();
160 Double_t y = 0.5 * TMath::Log( (e+fPz) / (e-fPz) );
161 return y;
162}

◆ GetStartT()

Double_t ShipMCTrack::GetStartT ( ) const
inline

Definition at line 66 of file ShipMCTrack.h.

66{ return fStartT; }

◆ GetStartVertex()

void ShipMCTrack::GetStartVertex ( TVector3 &  vertex)
inline

Definition at line 152 of file ShipMCTrack.h.

153{
154 vertex.SetXYZ(fStartX,fStartY,fStartZ);
155}

◆ GetStartX()

Double_t ShipMCTrack::GetStartX ( ) const
inline

Definition at line 63 of file ShipMCTrack.h.

63{ return fStartX; }

◆ GetStartY()

Double_t ShipMCTrack::GetStartY ( ) const
inline

Definition at line 64 of file ShipMCTrack.h.

64{ return fStartY; }

◆ GetStartZ()

Double_t ShipMCTrack::GetStartZ ( ) const
inline

Definition at line 65 of file ShipMCTrack.h.

65{ return fStartZ; }

◆ GetWeight()

Double_t ShipMCTrack::GetWeight ( ) const

Definition at line 149 of file ShipMCTrack.cxx.

150{
151 return fW;
152}

◆ MultiplyWeight()

void ShipMCTrack::MultiplyWeight ( Double_t  w)
inline

Definition at line 75 of file ShipMCTrack.h.

75{fW = fW*w;}

◆ Print()

void ShipMCTrack::Print ( Int_t  iTrack = 0) const

Output to screen

Definition at line 108 of file ShipMCTrack.cxx.

109{
110 LOG(DEBUG) << "Track " << trackId << ", mother : " << fMotherId << ", Type "
111 << fPdgCode << ", momentum (" << fPx << ", " << fPy << ", "
112 << fPz << ") GeV" ;
113 /* LOG(DEBUG2) << " Ref " << GetNPoints(kREF)
114 << ", TutDet " << GetNPoints(kTutDet)
115 << ", Rutherford " << GetNPoints(kFairRutherford) ;
116*/
117}

◆ SetMotherId()

void ShipMCTrack::SetMotherId ( Int_t  id)
inline

Modifiers

Definition at line 88 of file ShipMCTrack.h.

◆ SetNPoints()

void ShipMCTrack::SetNPoints ( Int_t  iDet,
Int_t  np 
)

Definition at line 187 of file ShipMCTrack.cxx.

188{
189/*
190 if ( iDet == kREF ) {
191 if ( nPoints < 0 ) { nPoints = 0; }
192 else if ( nPoints > 1 ) { nPoints = 1; }
193 fNPoints = ( fNPoints & ( ~ 1 ) ) | nPoints;
194 }
195
196 else if ( iDet == kTutDet ) {
197 if ( nPoints < 0 ) { nPoints = 0; }
198 else if ( nPoints > 7 ) { nPoints = 7; }
199 fNPoints = ( fNPoints & ( ~ ( 7 << 1 ) ) ) | ( nPoints << 1 );
200 }
201
202 else if ( iDet == kFairRutherford ) {
203 if ( nPoints < 0 ) { nPoints = 0; }
204 else if ( nPoints > 31 ) { nPoints = 31; }
205 fNPoints = ( fNPoints & ( ~ ( 31 << 4 ) ) ) | ( nPoints << 4 );
206 }
207
208 else LOG(ERROR) << "Unknown detector ID " << iDet ;
209*/
210}

◆ SetProcID()

void ShipMCTrack::SetProcID ( Int_t  i)
inline

Definition at line 67 of file ShipMCTrack.h.

67{ fProcID = i; }
int i
Definition ShipAna.py:86

◆ SetWeight()

void ShipMCTrack::SetWeight ( Double_t  w)
inline

Definition at line 76 of file ShipMCTrack.h.

76{fW = w;}

Member Data Documentation

◆ fM

Double32_t ShipMCTrack::fM
private

Definition at line 102 of file ShipMCTrack.h.

◆ fMotherId

Int_t ShipMCTrack::fMotherId
private

Index of mother track. -1 for primary particles.

Definition at line 99 of file ShipMCTrack.h.

◆ fNPoints

Int_t ShipMCTrack::fNPoints
private

Bitvector representing the number of MCPoints for this track in each subdetector. The detectors are represented by REF: Bit 0 (1 bit, max. value 1) MVD: Bit 1 - 3 (3 bits, max. value 7) STS: Bit 4 - 8 (5 bits, max. value 31) RICH: Bit 9 (1 bit, max. value 1) MUCH: Bit 10 - 14 (5 bits, max. value 31) TRD: Bit 15 - 19 (5 bits, max. value 31) TOF: Bit 20 - 23 (4 bits, max. value 15) ECAL: Bit 24 (1 bit, max. value 1) ZDC: Bit 25 (1 bit, max. value 1) The respective point numbers can be accessed and modified with the inline functions. Bits 26-31 are spare for potential additional detectors.

Definition at line 128 of file ShipMCTrack.h.

◆ fPdgCode

Int_t ShipMCTrack::fPdgCode
private

PDG particle code

Definition at line 96 of file ShipMCTrack.h.

◆ fProcID

Int_t ShipMCTrack::fProcID
private

Geant4 process ID which created the particle

Definition at line 111 of file ShipMCTrack.h.

◆ fPx

Double32_t ShipMCTrack::fPx
private

Momentum components at start vertex [GeV]

Definition at line 102 of file ShipMCTrack.h.

◆ fPy

Double32_t ShipMCTrack::fPy
private

Definition at line 102 of file ShipMCTrack.h.

◆ fPz

Double32_t ShipMCTrack::fPz
private

Definition at line 102 of file ShipMCTrack.h.

◆ fStartT

Double32_t ShipMCTrack::fStartT
private

Definition at line 105 of file ShipMCTrack.h.

◆ fStartX

Double32_t ShipMCTrack::fStartX
private

Coordinates of start vertex [cm, ns]

Definition at line 105 of file ShipMCTrack.h.

◆ fStartY

Double32_t ShipMCTrack::fStartY
private

Definition at line 105 of file ShipMCTrack.h.

◆ fStartZ

Double32_t ShipMCTrack::fStartZ
private

Definition at line 105 of file ShipMCTrack.h.

◆ fW

Double32_t ShipMCTrack::fW
private

weight

Definition at line 108 of file ShipMCTrack.h.


The documentation for this class was generated from the following files: