Geant4 11.1.1
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
G4FragmentingString Class Reference

#include <G4FragmentingString.hh>

Public Member Functions

 G4FragmentingString (const G4FragmentingString &right)
 
 G4FragmentingString (const G4ExcitedString &excited)
 
 G4FragmentingString (const G4FragmentingString &old, G4ParticleDefinition *newdecay, const G4LorentzVector *momentum)
 
 G4FragmentingString (const G4FragmentingString &old, G4ParticleDefinition *newdecay)
 
 ~G4FragmentingString ()
 
G4FragmentingStringoperator= (const G4FragmentingString &)
 
G4bool operator== (const G4FragmentingString &right) const
 
G4bool operator!= (const G4FragmentingString &right) const
 
G4LorentzVector Get4Momentum () const
 
G4ThreeVector StablePt ()
 
G4ThreeVector DecayPt ()
 
G4double LightConePlus ()
 
G4double LightConeMinus ()
 
G4double LightConeDecay ()
 
G4double Mass () const
 
G4double Mass2 () const
 
G4double MassT2 () const
 
G4ParticleDefinitionGetLeftParton (void) const
 
G4ParticleDefinitionGetRightParton (void) const
 
G4ParticleDefinitionGetStableParton () const
 
G4ParticleDefinitionGetDecayParton () const
 
void SetLeftPartonStable ()
 
void SetRightPartonStable ()
 
G4int GetDecayDirection () const
 
G4bool DecayIsQuark ()
 
G4bool StableIsQuark ()
 
G4bool IsAFourQuarkString (void) const
 
G4LorentzVector GetPstring ()
 
G4LorentzVector GetPleft ()
 
void SetPleft (G4LorentzVector a4momentum)
 
G4LorentzVector GetPright ()
 
void SetPright (G4LorentzVector a4momentum)
 
void LorentzRotate (const G4LorentzRotation &rotation)
 
G4LorentzRotation TransformToCenterOfMass ()
 
G4LorentzRotation TransformToAlignedCms ()
 
void Boost (G4ThreeVector &Velocity)
 

Detailed Description

Definition at line 49 of file G4FragmentingString.hh.

Constructor & Destructor Documentation

◆ G4FragmentingString() [1/4]

G4FragmentingString::G4FragmentingString ( const G4FragmentingString right)

Definition at line 46 of file G4FragmentingString.cc.

47{
48 LeftParton=old.LeftParton;
49 RightParton=old.RightParton;
50 Ptleft=old.Ptleft;
51 Ptright=old.Ptright;
52 Pplus=old.Pplus;
53 Pminus=old.Pminus;
54 theStableParton=old.theStableParton;
55 theDecayParton=old.theDecayParton;
56 decaying=old.decaying;
57 Pstring=old.Pstring;
58 Pleft =old.Pleft;
59 Pright =old.Pright;
60}

◆ G4FragmentingString() [2/4]

G4FragmentingString::G4FragmentingString ( const G4ExcitedString excited)

Definition at line 84 of file G4FragmentingString.cc.

85{
86 LeftParton=excited.GetLeftParton()->GetDefinition();
87 RightParton=excited.GetRightParton()->GetDefinition();
88 Ptleft=excited.GetLeftParton()->Get4Momentum().vect();
89 Ptleft.setZ(0.);
90 Ptright=excited.GetRightParton()->Get4Momentum().vect();
91 Ptright.setZ(0.);
92 theStableParton=0;
93 theDecayParton=0;
94
95 if (excited.GetDirection() > 0) {decaying=Left; }
96 else {decaying=Right;}
97
98 Pleft = excited.GetLeftParton()->Get4Momentum();
99 Pright = excited.GetRightParton()->Get4Momentum();
100 Pstring= Pleft + Pright;
101
102 Pplus = Pstring.plus();
103 Pminus = Pstring.minus();
104}
void setZ(double)
Hep3Vector vect() const
double minus() const
G4int GetDirection(void) const
G4Parton * GetRightParton(void) const
G4Parton * GetLeftParton(void) const
G4ParticleDefinition * GetDefinition()
Definition: G4Parton.hh:161
const G4LorentzVector & Get4Momentum() const
Definition: G4Parton.hh:143

◆ G4FragmentingString() [3/4]

G4FragmentingString::G4FragmentingString ( const G4FragmentingString old,
G4ParticleDefinition newdecay,
const G4LorentzVector momentum 
)

Definition at line 108 of file G4FragmentingString.cc.

111{
112 decaying=None;
113 // Momentum of produced hadron
114 G4LorentzVector Momentum = G4LorentzVector(momentum->vect(),momentum->e());
115
116 if ( old.decaying == Left )
117 {
118 RightParton= old.RightParton;
119 Ptright = old.Ptright;
120 Pright = old.Pright;
121
122 LeftParton = newdecay;
123 Ptleft = old.Ptleft - momentum->vect();
124 Ptleft.setZ(0.);
125 Pleft = old.Pleft - Momentum;
126
127 Pstring = Pleft + Pright;
128 Pplus = Pstring.plus();
129 Pminus = Pstring.minus();
130
131 theDecayParton=GetLeftParton();
132 theStableParton=GetRightParton();
133 decaying = Left;
134 } else if ( old.decaying == Right )
135 {
136 RightParton = newdecay;
137 Ptright = old.Ptright - momentum->vect();
138 Ptright.setZ(0.);
139 Pright = old.Pright - Momentum;
140
141 LeftParton = old.LeftParton;
142 Ptleft = old.Ptleft;
143 Pleft = old.Pleft;
144
145 Pstring = Pleft + Pright;
146 Pplus = Pstring.plus();
147 Pminus = Pstring.minus();
148
149 theDecayParton=GetRightParton();
150 theStableParton=GetLeftParton();
151 decaying = Right;
152 } else
153 {
154 throw G4HadronicException(__FILE__, __LINE__,
155 "G4FragmentingString::G4FragmentingString: no decay Direction defined");
156 }
157}
CLHEP::HepLorentzVector G4LorentzVector
G4ParticleDefinition * GetLeftParton(void) const
G4ParticleDefinition * GetRightParton(void) const

◆ G4FragmentingString() [4/4]

G4FragmentingString::G4FragmentingString ( const G4FragmentingString old,
G4ParticleDefinition newdecay 
)

Definition at line 162 of file G4FragmentingString.cc.

164{
165 decaying=None;
166
167 Ptleft.setX(0.); Ptleft.setY(0.); Ptleft.setZ(0.);
168 Ptright.setX(0.); Ptright.setY(0.); Ptright.setZ(0.);
169 Pplus=0.; Pminus=0.;
170 theStableParton=0; theDecayParton=0;
171
172 Pstring=G4LorentzVector(0.,0.,0.,0.);
173 Pleft =G4LorentzVector(0.,0.,0.,0.);
174 Pright =G4LorentzVector(0.,0.,0.,0.);
175
176 if ( old.decaying == Left )
177 {
178 RightParton= old.RightParton;
179 LeftParton = newdecay;
180 decaying=Left;
181 } else if ( old.decaying == Right )
182 {
183 RightParton = newdecay;
184 LeftParton = old.LeftParton;
185 decaying=Right;
186 } else
187 {
188 throw G4HadronicException(__FILE__, __LINE__,
189 "G4FragmentingString::G4FragmentingString: no decay Direction defined");
190 }
191}
void setY(double)
void setX(double)

◆ ~G4FragmentingString()

G4FragmentingString::~G4FragmentingString ( )

Definition at line 196 of file G4FragmentingString.cc.

197{}

Member Function Documentation

◆ Boost()

void G4FragmentingString::Boost ( G4ThreeVector Velocity)

◆ DecayIsQuark()

G4bool G4FragmentingString::DecayIsQuark ( )

Definition at line 238 of file G4FragmentingString.cc.

239{
240 return theDecayParton->GetParticleSubType()== "quark";
241}
const G4String & GetParticleSubType() const

◆ DecayPt()

G4ThreeVector G4FragmentingString::DecayPt ( )

Definition at line 258 of file G4FragmentingString.cc.

259{
260 if (decaying == Left ) return Ptleft;
261 else if (decaying == Right ) return Ptright;
262 else throw G4HadronicException(__FILE__, __LINE__, "G4FragmentingString::DecayPt: decay side UNdefined!");
263 return G4ThreeVector();
264}
CLHEP::Hep3Vector G4ThreeVector

◆ Get4Momentum()

G4LorentzVector G4FragmentingString::Get4Momentum ( ) const

Definition at line 287 of file G4FragmentingString.cc.

288{
289 return Pstring;
290}

◆ GetDecayDirection()

G4int G4FragmentingString::GetDecayDirection ( ) const

Definition at line 219 of file G4FragmentingString.cc.

220{
221 if (decaying == Left ) return +1;
222 else if (decaying == Right) return -1;
223 else throw G4HadronicException(__FILE__, __LINE__,
224 "G4FragmentingString::GetDecayDirection: decay side UNdefined!");
225 return 0;
226}

Referenced by G4QGSMFragmentation::FragmentString().

◆ GetDecayParton()

G4ParticleDefinition * G4FragmentingString::GetDecayParton ( ) const
inline

Definition at line 137 of file G4FragmentingString.hh.

138{
139 return theDecayParton;
140}

◆ GetLeftParton()

G4ParticleDefinition * G4FragmentingString::GetLeftParton ( void  ) const
inline

◆ GetPleft()

G4LorentzVector G4FragmentingString::GetPleft ( )

Definition at line 312 of file G4FragmentingString.cc.

313{
314 return Pleft;
315}

◆ GetPright()

G4LorentzVector G4FragmentingString::GetPright ( )

Definition at line 317 of file G4FragmentingString.cc.

318{
319 return Pright;
320}

◆ GetPstring()

G4LorentzVector G4FragmentingString::GetPstring ( )

Definition at line 307 of file G4FragmentingString.cc.

308{
309 return Pstring;
310}

◆ GetRightParton()

G4ParticleDefinition * G4FragmentingString::GetRightParton ( void  ) const
inline

◆ GetStableParton()

G4ParticleDefinition * G4FragmentingString::GetStableParton ( ) const
inline

Definition at line 131 of file G4FragmentingString.hh.

132{
133 return theStableParton;
134}

◆ IsAFourQuarkString()

G4bool G4FragmentingString::IsAFourQuarkString ( void  ) const

Definition at line 230 of file G4FragmentingString.cc.

231{
232 return LeftParton->GetParticleSubType()== "di_quark"
233 && RightParton->GetParticleSubType()== "di_quark";
234}

Referenced by G4LundStringFragmentation::FragmentString(), and G4VLongitudinalStringDecay::PossibleHadronMass().

◆ LightConeDecay()

G4double G4FragmentingString::LightConeDecay ( )

Definition at line 278 of file G4FragmentingString.cc.

279{
280 if (decaying == Left ) return Pplus;
281 else if (decaying == Right ) return Pminus;
282 else throw G4HadronicException(__FILE__, __LINE__, "G4FragmentingString::DecayPt: decay side UNdefined!");
283}

◆ LightConeMinus()

G4double G4FragmentingString::LightConeMinus ( )

Definition at line 273 of file G4FragmentingString.cc.

274{
275 return Pminus;
276}

◆ LightConePlus()

G4double G4FragmentingString::LightConePlus ( )

Definition at line 268 of file G4FragmentingString.cc.

269{
270 return Pplus;
271}

◆ LorentzRotate()

void G4FragmentingString::LorentzRotate ( const G4LorentzRotation rotation)
inline

Definition at line 156 of file G4FragmentingString.hh.

157{
158 SetPleft(rotation*Pleft);
159 SetPright(rotation*Pright);
160 Pstring = Pleft+Pright;
161 Ptleft =Pleft.vect(); Ptleft.setZ(0.);
162 Ptright=Pright.vect(); Ptright.setZ(0.);
163 Pplus =Pstring.plus();
164 Pminus=Pstring.minus();
165}
void SetPleft(G4LorentzVector a4momentum)
void SetPright(G4LorentzVector a4momentum)

◆ Mass()

G4double G4FragmentingString::Mass ( ) const

Definition at line 297 of file G4FragmentingString.cc.

298{
299 return Pstring.mag();
300}

Referenced by G4VLongitudinalStringDecay::ProduceOneHadron().

◆ Mass2()

G4double G4FragmentingString::Mass2 ( ) const

Definition at line 292 of file G4FragmentingString.cc.

293{
294 return Pstring.mag2();
295}

◆ MassT2()

G4double G4FragmentingString::MassT2 ( ) const

Definition at line 302 of file G4FragmentingString.cc.

303{
304 return Pplus*Pminus;
305}

◆ operator!=()

G4bool G4FragmentingString::operator!= ( const G4FragmentingString right) const
inline

Definition at line 124 of file G4FragmentingString.hh.

125{
126 return this != &right;
127}

◆ operator=()

G4FragmentingString & G4FragmentingString::operator= ( const G4FragmentingString old)

Definition at line 62 of file G4FragmentingString.cc.

63{
64 if (this != &old)
65 {
66 LeftParton=old.LeftParton;
67 RightParton=old.RightParton;
68 Ptleft=old.Ptleft;
69 Ptright=old.Ptright;
70 Pplus=old.Pplus;
71 Pminus=old.Pminus;
72 theStableParton=old.theStableParton;
73 theDecayParton=old.theDecayParton;
74 decaying=old.decaying;
75 Pstring=old.Pstring;
76 Pleft =old.Pleft;
77 Pright =old.Pright;
78 }
79 return *this;
80}

◆ operator==()

G4bool G4FragmentingString::operator== ( const G4FragmentingString right) const
inline

Definition at line 118 of file G4FragmentingString.hh.

119{
120 return this == &right;
121}

◆ SetLeftPartonStable()

void G4FragmentingString::SetLeftPartonStable ( )

Definition at line 201 of file G4FragmentingString.cc.

202{
203 theStableParton=GetLeftParton();
204 theDecayParton=GetRightParton();
205 decaying=Right;
206}

◆ SetPleft()

void G4FragmentingString::SetPleft ( G4LorentzVector  a4momentum)
inline

Definition at line 184 of file G4FragmentingString.hh.

185{
186 Pleft = a4momentum;
187 Ptleft = Pleft.vect(); Ptleft.setZ(0.);
188 Pstring = Pleft + Pright;
189 Pplus = Pstring.plus();
190 Pminus = Pstring.minus();
191}

Referenced by LorentzRotate().

◆ SetPright()

void G4FragmentingString::SetPright ( G4LorentzVector  a4momentum)
inline

Definition at line 194 of file G4FragmentingString.hh.

195{
196 Pright = a4momentum;
197 Ptright = Pright.vect(); Ptright.setZ(0.);
198 Pstring = Pleft + Pright;
199 Pplus = Pstring.plus();
200 Pminus = Pstring.minus();
201}

Referenced by LorentzRotate().

◆ SetRightPartonStable()

void G4FragmentingString::SetRightPartonStable ( )

Definition at line 210 of file G4FragmentingString.cc.

211{
212 theStableParton=GetRightParton();
213 theDecayParton=GetLeftParton();
214 decaying=Left;
215}

◆ StableIsQuark()

G4bool G4FragmentingString::StableIsQuark ( )

Definition at line 243 of file G4FragmentingString.cc.

244{
245 return theStableParton->GetParticleSubType()== "quark";
246}

◆ StablePt()

G4ThreeVector G4FragmentingString::StablePt ( )

Definition at line 250 of file G4FragmentingString.cc.

251{
252 if (decaying == Left ) return Ptright;
253 else if (decaying == Right ) return Ptleft;
254 else throw G4HadronicException(__FILE__, __LINE__, "G4FragmentingString::DecayPt: decay side UNdefined!");
255 return G4ThreeVector();
256}

◆ TransformToAlignedCms()

G4LorentzRotation G4FragmentingString::TransformToAlignedCms ( )

Definition at line 322 of file G4FragmentingString.cc.

323{
324 G4LorentzVector momentum = Pstring;
325 G4LorentzRotation toAlignedCms(-1*momentum.boostVector());
326 momentum = toAlignedCms * Pleft;
327
328 toAlignedCms.rotateZ(-1*momentum.phi());
329 toAlignedCms.rotateY(-1*momentum.theta());
330
331 Pleft *= toAlignedCms;
332 Pright *= toAlignedCms;
333 Pstring *= toAlignedCms;
334
335 Ptleft = G4ThreeVector(Pleft.vect());
336 Ptleft.setZ(0.);
337 Ptright = G4ThreeVector(Pright.vect());
338 Pplus = Pstring.plus();
339 Pminus = Pstring.minus();
340
341 return toAlignedCms;
342}
double theta() const
Hep3Vector boostVector() const
HepLorentzVector & rotateZ(double)
HepLorentzVector & rotateY(double)

◆ TransformToCenterOfMass()

G4LorentzRotation G4FragmentingString::TransformToCenterOfMass ( )
inline

Definition at line 168 of file G4FragmentingString.hh.

169{
170 G4LorentzVector momentum=Pstring;
171 G4LorentzRotation toCMS(-1*momentum.boostVector());
172
173 Pleft *= toCMS;
174 Pright *= toCMS;
175 Pstring *= toCMS;
176 Ptleft =Pleft.vect(); Ptleft.setZ(0.);
177 Ptright=Pright.vect(); Ptright.setZ(0.);
178 Pplus =Pstring.plus();
179 Pminus=Pstring.minus();
180 return toCMS;
181}

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