Geant4  10.02.p02
G4FragmentingString.cc
Go to the documentation of this file.
1 //
2 // ********************************************************************
3 // * License and Disclaimer *
4 // * *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
10 // * *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
17 // * *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 //
26 //
27 
28 
29 // ------------------------------------------------------------
30 // GEANT 4 class implementation file
31 //
32 // ---------------- G4FragmentingString ----------------
33 // by Gunter Folger, September 2001.
34 // class for an excited string used in Fragmention
35 // ------------------------------------------------------------
36 
37 
38 // G4FragmentingString
39 #include "G4FragmentingString.hh"
40 #include "G4ExcitedString.hh"
41 
42 //---------------------------------------------------------------------------------
43 
44 //---------------------------------------------------------------------------------
45 
47 {
50  Ptleft=old.Ptleft;
51  Ptright=old.Ptright;
52  Pplus=old.Pplus;
53  Pminus=old.Pminus;
56  decaying=old.decaying;
57 }
58 
60 {
61  if (this != &old)
62  {
65  Ptleft=old.Ptleft;
66  Ptright=old.Ptright;
67  Pplus=old.Pplus;
68  Pminus=old.Pminus;
71  decaying=old.decaying;
72  }
73  return *this;
74 }
75 
76 //---------------------------------------------------------------------------------
77 
79 {
82  Ptleft=excited.GetLeftParton()->Get4Momentum().vect();
83  Ptleft.setZ(0.);
84  Ptright=excited.GetRightParton()->Get4Momentum().vect();
85  Ptright.setZ(0.);
86  G4LorentzVector P=excited.Get4Momentum();
87  Pplus =P.e() + P.pz();
88  Pminus=P.e() - P.pz();
91 // decaying=None; // Uzhi 19.06.2014
92  if(excited.GetDirection() > 0) {decaying=Left; } // Uzhi 20.06.2014
93  else {decaying=Right;} // Uzhi 20.06.2014
94 }
95 
96 //---------------------------------------------------------------------------------
97 
99  G4ParticleDefinition * newdecay,
100  const G4LorentzVector *momentum)
101 {
102  decaying=None;
103  if ( old.decaying == Left )
104  {
106  Ptright = old.Ptright;
107  LeftParton = newdecay;
108  Ptleft = old.Ptleft - momentum->vect();
109  Ptleft.setZ(0.);
112  decaying=Left; // Uzhi 19.06.2014
113  } else if ( old.decaying == Right )
114  {
115  RightParton = newdecay;
116  Ptright = old.Ptright - momentum->vect();
117  Ptright.setZ(0.);
118  LeftParton = old.LeftParton;
119  Ptleft = old.Ptleft;
122  decaying=Right; // Uzhi 19.06.2014
123  } else
124  {
125  throw G4HadronicException(__FILE__, __LINE__, "G4FragmentingString::G4FragmentingString: no decay Direction defined");
126  }
127  Pplus = old.Pplus - (momentum->e() + momentum->pz());
128  Pminus = old.Pminus - (momentum->e() - momentum->pz());
129 
130  //G4double Eold=0.5 * (old.Pplus + old.Pminus);
131  //G4double Enew=0.5 * (Pplus + Pminus);
132 }
133 
134 
135 //---------------------------------------------------------------------------------
136 
138  G4ParticleDefinition * newdecay)
139 {
140  decaying=None;
141 
142  Ptleft.setX(0.); Ptleft.setY(0.); Ptleft.setZ(0.);
143  Ptright.setX(0.); Ptright.setY(0.); Ptright.setZ(0.);
144  Pplus=0.; Pminus=0.;
146 
147  if ( old.decaying == Left )
148  {
149  RightParton= old.RightParton;
150  LeftParton = newdecay;
151  decaying=Left; // Uzhi 19.06.2014
152  } else if ( old.decaying == Right )
153  {
154  RightParton = newdecay;
155  LeftParton = old.LeftParton;
156  decaying=Right; // Uzhi 19.06.2014
157  } else
158  {
159  throw G4HadronicException(__FILE__, __LINE__, "G4FragmentingString::G4FragmentingString: no decay Direction defined");
160  }
161 }
162 
163 
164 //---------------------------------------------------------------------------------
165 
167 {}
168 
169 
170 //---------------------------------------------------------------------------------
171 
173 {
176  decaying=Right;
177 }
178 
179 //---------------------------------------------------------------------------------
180 
182 {
185  decaying=Left;
186 }
187 
188 //---------------------------------------------------------------------------------
189 
191 {
192  if (decaying == Left ) return +1;
193  else if (decaying == Right) return -1;
194  else throw G4HadronicException(__FILE__, __LINE__, "G4FragmentingString::GetDecayDirection: decay side UNdefined!");
195  return 0;
196 }
197 
198 //---------------------------------------------------------------------------------
199 
201 {
202  return LeftParton->GetParticleSubType()== "di_quark"
203  && RightParton->GetParticleSubType()== "di_quark";
204 }
205 
206 //---------------------------------------------------------------------------------
207 
209 {
210  return theDecayParton->GetParticleSubType()== "quark";
211 }
212 
214 {
215  return theStableParton->GetParticleSubType()== "quark";
216 }
217 
218 //---------------------------------------------------------------------------------
219 
221 {
222  if (decaying == Left ) return Ptright;
223  else if (decaying == Right ) return Ptleft;
224  else throw G4HadronicException(__FILE__, __LINE__, "G4FragmentingString::DecayPt: decay side UNdefined!");
225  return G4ThreeVector();
226 }
227 
229 {
230  if (decaying == Left ) return Ptleft;
231  else if (decaying == Right ) return Ptright;
232  else throw G4HadronicException(__FILE__, __LINE__, "G4FragmentingString::DecayPt: decay side UNdefined!");
233  return G4ThreeVector();
234 }
235 
236 //---------------------------------------------------------------------------------
237 
239 {
240  return Pplus;
241 }
242 
244 {
245  return Pminus;
246 }
247 
249 {
250  if (decaying == Left ) return Pplus;
251  else if (decaying == Right ) return Pminus;
252  else throw G4HadronicException(__FILE__, __LINE__, "G4FragmentingString::DecayPt: decay side UNdefined!");
253 }
254 
255 //---------------------------------------------------------------------------------
256 
258 {
259  G4LorentzVector momentum(Ptleft+Ptright,0);
260  momentum.setPz(0.5*(Pplus-Pminus));
261  momentum.setE(0.5*(Pplus+Pminus));
262  return momentum;
263 }
264 
266 {
267  return Pplus*Pminus - (Ptleft+Ptright).mag2();
268 }
269 
271 {
272  return std::sqrt(this->Mass2());
273 }
274 
276 {
277  return Pplus*Pminus;
278 }
G4ParticleDefinition * RightParton
G4FragmentingString(const G4FragmentingString &right)
G4ParticleDefinition * GetRightParton(void) const
CLHEP::Hep3Vector G4ThreeVector
const G4LorentzVector & Get4Momentum() const
Definition: G4Parton.hh:140
G4Parton * GetLeftParton(void) const
G4double Mass() const
G4ParticleDefinition * GetDefinition()
Definition: G4Parton.hh:158
const G4String & GetParticleSubType() const
int G4int
Definition: G4Types.hh:78
static double P[]
G4ParticleDefinition * GetLeftParton(void) const
G4double Mass2() const
G4ParticleDefinition * theStableParton
G4LorentzVector Get4Momentum() const
bool G4bool
Definition: G4Types.hh:79
G4int GetDecayDirection() const
G4bool FourQuarkString(void) const
G4ParticleDefinition * LeftParton
G4LorentzVector Get4Momentum() const
G4Parton * GetRightParton(void) const
G4ParticleDefinition * theDecayParton
G4FragmentingString & operator=(const G4FragmentingString &)
G4double MassT2() const
G4int GetDirection(void) const
double G4double
Definition: G4Types.hh:76
CLHEP::HepLorentzVector G4LorentzVector