ECCE @ EIC Software
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4DiffractiveSplitableHadron.cc
Go to the documentation of this file. Or view the newest version in sPHENIX GitHub for file G4DiffractiveSplitableHadron.cc
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 // ---------------- G4DiffractiveSplitableHadron----------------
33 // by Gunter Folger, August 1998.
34 // class splitting an interacting particle. Used by FTF String Model.
35 // ------------------------------------------------------------
36 
38 
39 #include "G4ParticleDefinition.hh"
40 #include "Randomize.hh"
41 
42 
43 //============================================================================
44 
46 {
47  PartonIndex = -1;
48  G4LorentzVector tmp=G4LorentzVector(0.,0.,0.,0.);
49  Parton[0] = new G4Parton( 1 );
50  Parton[1] = new G4Parton(-1 );
51 
52  Parton[0]->Set4Momentum(tmp); Parton[1]->Set4Momentum(tmp);
53 }
54 
55 
56 //============================================================================
57 
59  G4VSplitableHadron( aPrimary )
60 {
61  PartonIndex = -2;
62  Parton[0] = NULL;
63 }
64 
65 
66 //============================================================================
67 
69  G4VSplitableHadron( aNucleon )
70 {
71  PartonIndex = -2;
72  Parton[0] = NULL;
73 }
74 
75 
76 //============================================================================
77 
79  G4VSplitableHadron( aNucleon )
80 {
81  PartonIndex = -2;
82  Parton[0] = NULL;
83 }
84 
85 
86 //============================================================================
87 
89 
90 
91 //============================================================================
92 
94 
95  if ( IsSplit() ) return;
96  Splitting();
97  // Split once only...
98  if ( Parton[0] != NULL ) return;
99 
100  // flavours of quark ends
101  G4int PDGcode = GetDefinition()->GetPDGEncoding();
102  G4int stringStart, stringEnd;
103  ChooseStringEnds( PDGcode, &stringStart, &stringEnd );
104 
105  Parton[0] = new G4Parton( stringStart );
106  Parton[1] = new G4Parton( stringEnd );
107 
108  G4LorentzVector tmp=G4LorentzVector(0.,0.,0.,0.);
109  Parton[0]->Set4Momentum(tmp); Parton[1]->Set4Momentum(tmp);
110 
111  /* // Inversion of a string
112  if ( G4UniformRand() < 1.75 ) { //0.75
113  Parton[0] = new G4Parton( stringStart );
114  Parton[1] = new G4Parton( stringEnd );
115  } else {
116  Parton[0] = new G4Parton( stringEnd );
117  Parton[1] = new G4Parton( stringStart );
118  }
119  */
120 
121  PartonIndex = -1;
122 }
123 
124 
125 //============================================================================
126 
128  ++PartonIndex;
129  if ( PartonIndex > 1 || PartonIndex < 0 ) return NULL;
130  G4int PartonInd( PartonIndex );
131  if ( PartonIndex == 1 ) PartonIndex = -1;
132  return Parton[ PartonInd ];
133 }
134 
135 
136 //============================================================================
137 
139  ++PartonIndex;
140  if ( PartonIndex > 1 || PartonIndex < 0 ) return NULL;
141  G4int PartonInd( PartonIndex );
142  if ( PartonIndex == 1 ) PartonIndex = -1;
143  return Parton[ PartonInd ];
144 }
145 
146 
147 //============================================================================
148 
150  delete Parton[0];
151  Parton[0] = new G4Parton( PDGcode );
152  G4LorentzVector tmp=G4LorentzVector(0.,0.,0.,0.);
153  Parton[0]->Set4Momentum(tmp);
154 }
155 
156 
157 //============================================================================
158 
160  delete Parton[1];
161  Parton[1] = new G4Parton( PDGcode );
162  G4LorentzVector tmp=G4LorentzVector(0.,0.,0.,0.);
163  Parton[1]->Set4Momentum(tmp);
164 }
165 
166 
167 //============================================================================
168 
170  G4int* bEnd ) const {
171  G4int absPDGcode = std::abs( PDGcode );
172 
173  if ( absPDGcode < 1000 ) { //-------------------- Meson -------------
174  G4int heavy(0), light(0);
175  if (!((absPDGcode == 111)||(absPDGcode == 221)||(absPDGcode == 331)))
176  { // Ordinary mesons =======================
177  heavy = absPDGcode/100;
178  light = (absPDGcode % 100)/10;
179  //G4int anti = std::pow( -1 , std::max( heavy, light ) );
180  G4int anti = 1 - 2*( std::max( heavy, light ) % 2 );
181  if (PDGcode < 0 ) anti *= -1;
182  heavy *= anti;
183  light *= -1 * anti;
184  }
185  else
186  { // Pi0, Eta, Eta' =======================
187  if ( G4UniformRand() < 0.5 ) {heavy = 1; light = -1;}
188  else {heavy = 2; light = -2;}
189  }
190  if ( G4UniformRand() < 0.5 ) {
191  *aEnd = heavy;
192  *bEnd = light;
193  } else {
194  *aEnd = light;
195  *bEnd = heavy;
196  }
197  } else { //-------------------- Baryon --------------
198  G4int j1000 = PDGcode/1000;
199  G4int j100 = (PDGcode % 1000)/100;
200  G4int j10 = (PDGcode % 100)/10;
201 
202  G4double SuppresUUDDSS=1.0/2.0;
203  if ((j1000 == j100) && (j1000 == j10)) SuppresUUDDSS=1.;
204 
205  const G4int maxNumberOfLoops = 1000;
206  G4int loopCounter = 0;
207  do
208  {
209  G4double random = G4UniformRand();
210 
211  if (random < 0.33333)
212  {
213  if (( j100 == j10 ) && ( G4UniformRand() > SuppresUUDDSS )) continue;
214  *aEnd = j1000;
215  if ( j100 == j10 ) {*bEnd = Diquark( j100, j10, 1 );}
216  else
217  if ( G4UniformRand() > 0.25) {*bEnd = Diquark( j100, j10, 0 );}
218  else {*bEnd = Diquark( j100, j10, 1 );}
219  break;
220  }
221  else if (random < 0.66667)
222  {
223  if (( j1000 == j10 ) && ( G4UniformRand() > SuppresUUDDSS )) continue;
224  *aEnd = j100;
225  if ( j1000 == j10 ) {*bEnd = Diquark( j1000, j10, 1 );}
226  else
227  if ( G4UniformRand() > 0.25) {*bEnd = Diquark( j1000, j10, 0 );}
228  else {*bEnd = Diquark( j1000, j10, 1 );}
229  break;
230  }
231  else
232  {
233  if (( j1000 == j100 ) && ( G4UniformRand() > SuppresUUDDSS )) continue;
234  *aEnd = j10;
235  if ( j1000 == j100 ) {*bEnd = Diquark( j1000, j100, 1 );}
236  else
237  if ( G4UniformRand() > 0.25) {*bEnd = Diquark( j1000, j100, 0 );}
238  else {*bEnd = Diquark( j1000, j100, 1 );}
239  break;
240  }
241  } while ( (true) &&
242  ++loopCounter < maxNumberOfLoops ); /* Loop checking, 10.08.2015, A.Ribon */
243  if ( loopCounter >= maxNumberOfLoops ) {
244  *aEnd = j10; *bEnd = Diquark( j1000, j100, 1 ); // Just something acceptable, without any physics consideration.
245  }
246 
247  }
248 }
249 
250 
251 //============================================================================
252 
254  G4int diquarkPDG = std::max( std::abs( aquark ), std::abs( bquark ) ) * 1000 +
255  std::min( std::abs( aquark ), std::abs( bquark ) ) * 100 +
256  2*Spin + 1;
257  return ( aquark > 0 && bquark > 0 ) ? diquarkPDG : -1*diquarkPDG;
258 }
259