Geant4  10.01.p02
G4LundStringFragmentation.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 // $Id: G4LundStringFragmentation.cc 88879 2015-03-13 16:13:58Z gcosmo $
28 // GEANT4 tag $Name: $ 1.8
29 //
30 // -----------------------------------------------------------------------------
31 // GEANT 4 class implementation file
32 //
33 // History: first implementation, Maxim Komogorov, 10-Jul-1998
34 // -----------------------------------------------------------------------------
35 #include "G4LundStringFragmentation.hh"
36 #include "G4PhysicalConstants.hh"
37 #include "G4SystemOfUnits.hh"
38 #include "Randomize.hh"
39 #include "G4FragmentingString.hh"
40 #include "G4DiQuarks.hh"
41 #include "G4Quarks.hh"
42 
43 //#define debug_LUNDfragmentation // Uzhi 20.06.2014
44 
45 // Class G4LundStringFragmentation
46 //*************************************************************************************
47 
48 G4LundStringFragmentation::G4LundStringFragmentation()
49 {
50 // ------ For estimation of a minimal string mass ---------------
51  Mass_of_light_quark =140.*MeV;
52  Mass_of_heavy_quark =500.*MeV;
53  Mass_of_string_junction=720.*MeV;
54 // ------ An estimated minimal string mass ----------------------
55  MinimalStringMass = 0.;
56  MinimalStringMass2 = 0.;
57 // ------ Minimal invariant mass used at a string fragmentation -
58  WminLUND = 0.45*GeV; //0.23*GeV; // Uzhi 0.7 -> 0.23 3.8.10 //0.8 1.5
59 // ------ Smooth parameter used at a string fragmentation for ---
60 // ------ smearing sharp mass cut-off ---------------------------
61  SmoothParam = 0.2;
62 
63  SetStringTensionParameter(1.);
64  SetDiquarkBreakProbability(0.05); // 0.05 Uzhi 27.09.2014
65  SetStrangenessSuppression(0.447); // 0.47 Uzhi 27.09.2014
66  SetDiquarkSuppression(0.05); // 0.087 Uzhi 27.09.2014
67 
68 // For treating of small string decays
69  for(G4int i=0; i<3; i++)
70  { for(G4int j=0; j<3; j++)
71  { for(G4int k=0; k<6; k++)
72  {
73  Meson[i][j][k]=0; MesonWeight[i][j][k]=0.;
74  }
75  }
76  }
77 //--------------------------
78  Meson[0][0][0]=111; // dbar-d Pi0
79  MesonWeight[0][0][0]=(1.-pspin_meson)*(1.-scalarMesonMix[0]);
80 
81  Meson[0][0][1]=221; // dbar-d Eta
82  MesonWeight[0][0][1]=(1.-pspin_meson)*(scalarMesonMix[0]-scalarMesonMix[1]);
83 
84  Meson[0][0][2]=331; // dbar-d EtaPrime
85  MesonWeight[0][0][2]=(1.-pspin_meson)*(scalarMesonMix[1]);
86 
87  Meson[0][0][3]=113; // dbar-d Rho0
88  MesonWeight[0][0][3]=pspin_meson*(1.-vectorMesonMix[0]);
89 
90  Meson[0][0][4]=223; // dbar-d Omega
91  MesonWeight[0][0][4]=pspin_meson*(vectorMesonMix[0]);
92 //--------------------------
93 
94  Meson[0][1][0]=211; // dbar-u Pi+
95  MesonWeight[0][1][0]=(1.-pspin_meson);
96 
97  Meson[0][1][1]=213; // dbar-u Rho+
98  MesonWeight[0][1][1]=pspin_meson;
99 //--------------------------
100 
101  Meson[0][2][0]=311; // dbar-s K0bar
102  MesonWeight[0][2][0]=(1.-pspin_meson);
103 
104  Meson[0][2][1]=313; // dbar-s K*0bar
105  MesonWeight[0][2][1]=pspin_meson;
106 //--------------------------
107 //--------------------------
108  Meson[1][0][0]=211; // ubar-d Pi-
109  MesonWeight[1][0][0]=(1.-pspin_meson);
110 
111  Meson[1][0][1]=213; // ubar-d Rho-
112  MesonWeight[1][0][1]=pspin_meson;
113 //--------------------------
114 
115  Meson[1][1][0]=111; // ubar-u Pi0
116  MesonWeight[1][1][0]=(1.-pspin_meson)*(1.-scalarMesonMix[0]);
117 
118  Meson[1][1][1]=221; // ubar-u Eta
119  MesonWeight[1][1][1]=(1.-pspin_meson)*(scalarMesonMix[0]-scalarMesonMix[1]);
120 
121  Meson[1][1][2]=331; // ubar-u EtaPrime
122  MesonWeight[1][1][2]=(1.-pspin_meson)*(scalarMesonMix[1]);
123 
124  Meson[1][1][3]=113; // ubar-u Rho0
125  MesonWeight[1][1][3]=pspin_meson*(1.-vectorMesonMix[0]);
126 
127  Meson[1][1][4]=223; // ubar-u Omega
128  MesonWeight[1][1][4]=pspin_meson*(scalarMesonMix[0]);
129 //--------------------------
130 
131  Meson[1][2][0]=321; // ubar-s K-
132  MesonWeight[1][2][0]=(1.-pspin_meson);
133 
134  Meson[1][2][1]=323; // ubar-s K*-bar -
135  MesonWeight[1][2][1]=pspin_meson;
136 //--------------------------
137 //--------------------------
138 
139  Meson[2][0][0]=311; // sbar-d K0
140  MesonWeight[2][0][0]=(1.-pspin_meson);
141 
142  Meson[2][0][1]=313; // sbar-d K*0
143  MesonWeight[2][0][1]=pspin_meson;
144 //--------------------------
145 
146  Meson[2][1][0]=321; // sbar-u K+
147  MesonWeight[2][1][0]=(1.-pspin_meson);
148 
149  Meson[2][1][1]=323; // sbar-u K*+
150  MesonWeight[2][1][1]=pspin_meson;
151 //--------------------------
152 
153  Meson[2][2][0]=221; // sbar-s Eta
154  MesonWeight[2][2][0]=(1.-pspin_meson)*(1.-scalarMesonMix[5]);
155 
156  Meson[2][2][1]=331; // sbar-s EtaPrime
157  MesonWeight[2][2][1]=(1.-pspin_meson)*(1.-scalarMesonMix[5]);
158 
159  Meson[2][2][3]=333; // sbar-s EtaPrime
160  MesonWeight[2][2][3]=pspin_meson*(vectorMesonMix[5]);
161 //--------------------------
162 
163  for(G4int i=0; i<3; i++)
164  { for(G4int j=0; j<3; j++)
165  { for(G4int k=0; k<3; k++)
166  { for(G4int l=0; l<4; l++)
167  { Baryon[i][j][k][l]=0; BaryonWeight[i][j][k][l]=0.;}
168  }
169  }
170  }
171 
172  G4double pspin_barion_in=pspin_barion;
173  //pspin_barion=0.75;
174 //---------------------------------------
175  Baryon[0][0][0][0]=1114; // Delta-
176  BaryonWeight[0][0][0][0]=1.;
177 
178 //---------------------------------------
179  Baryon[0][0][1][0]=2112; // neutron
180  BaryonWeight[0][0][1][0]=1.-pspin_barion;
181 
182  Baryon[0][0][1][1]=2114; // Delta0
183  BaryonWeight[0][0][1][1]=pspin_barion;
184 
185 //---------------------------------------
186  Baryon[0][0][2][0]=3112; // Sigma-
187  BaryonWeight[0][0][2][0]=1.-pspin_barion;
188 
189  Baryon[0][0][2][1]=3114; // Sigma*-
190  BaryonWeight[0][0][2][1]=pspin_barion;
191 
192 //---------------------------------------
193  Baryon[0][1][0][0]=2112; // neutron
194  BaryonWeight[0][1][0][0]=1.-pspin_barion;
195 
196  Baryon[0][1][0][1]=2114; // Delta0
197  BaryonWeight[0][1][0][1]=pspin_barion;
198 
199 //---------------------------------------
200  Baryon[0][1][1][0]=2212; // proton
201  BaryonWeight[0][1][1][0]=1.-pspin_barion;
202 
203  Baryon[0][1][1][1]=2214; // Delta+
204  BaryonWeight[0][1][1][1]=pspin_barion;
205 
206 //---------------------------------------
207  Baryon[0][1][2][0]=3122; // Lambda
208  BaryonWeight[0][1][2][0]=(1.-pspin_barion)*0.5;
209 
210  Baryon[0][1][2][1]=3212; // Sigma0
211  BaryonWeight[0][1][2][1]=(1.-pspin_barion)*0.5;
212 
213  Baryon[0][1][2][2]=3214; // Sigma*0
214  BaryonWeight[0][1][2][2]=pspin_barion;
215 
216 //---------------------------------------
217  Baryon[0][2][0][0]=3112; // Sigma-
218  BaryonWeight[0][2][0][0]=1.-pspin_barion;
219 
220  Baryon[0][2][0][1]=3114; // Sigma*-
221  BaryonWeight[0][2][0][1]=pspin_barion;
222 
223 //---------------------------------------
224  Baryon[0][2][1][0]=3122; // Lambda
225  BaryonWeight[0][2][1][0]=(1.-pspin_barion)*0.5;
226 
227  Baryon[0][2][1][1]=3212; // Sigma0
228  BaryonWeight[0][2][1][1]=(1.-pspin_barion)*0.5;
229 
230  Baryon[0][2][1][2]=3214; // Sigma*0
231  BaryonWeight[0][2][1][2]=pspin_barion;
232 
233 //---------------------------------------
234  Baryon[0][2][2][0]=3312; // Theta-
235  BaryonWeight[0][2][2][0]=1.-pspin_barion;
236 
237  Baryon[0][2][2][1]=3314; // Theta*-
238  BaryonWeight[0][2][2][1]=pspin_barion;
239 
240 //---------------------------------------
241 //---------------------------------------
242  Baryon[1][0][0][0]=2112; // neutron
243  BaryonWeight[1][0][0][0]=1.-pspin_barion;
244 
245  Baryon[1][0][0][1]=2114; // Delta0
246  BaryonWeight[1][0][0][1]=pspin_barion;
247 
248 //---------------------------------------
249  Baryon[1][0][1][0]=2212; // proton
250  BaryonWeight[1][0][1][0]=1.-pspin_barion;
251 
252  Baryon[1][0][1][1]=2214; // Delta+
253  BaryonWeight[1][0][1][1]=pspin_barion;
254 
255 //---------------------------------------
256  Baryon[1][0][2][0]=3122; // Lambda
257  BaryonWeight[1][0][2][0]=(1.-pspin_barion)*0.5;
258 
259  Baryon[1][0][2][1]=3212; // Sigma0
260  BaryonWeight[1][0][2][1]=(1.-pspin_barion)*0.5;
261 
262  Baryon[1][0][2][2]=3214; // Sigma*0
263  BaryonWeight[1][0][2][2]=pspin_barion;
264 
265 //---------------------------------------
266  Baryon[1][1][0][0]=2212; // proton
267  BaryonWeight[1][1][0][0]=1.-pspin_barion;
268 
269  Baryon[1][1][0][1]=2214; // Delta+
270  BaryonWeight[1][1][0][1]=pspin_barion;
271 
272 //---------------------------------------
273  Baryon[1][1][1][0]=2224; // Delta++
274  BaryonWeight[1][1][1][0]=1.;
275 
276 //---------------------------------------
277  Baryon[1][1][2][0]=3222; // Sigma+
278  BaryonWeight[1][1][2][0]=1.-pspin_barion;
279 
280  Baryon[1][1][2][1]=3224; // Sigma*+
281  BaryonWeight[1][1][2][1]=pspin_barion;
282 
283 //---------------------------------------
284  Baryon[1][2][0][0]=3122; // Lambda
285  BaryonWeight[1][2][0][0]=(1.-pspin_barion)*0.5;
286 
287  Baryon[1][2][0][1]=3212; // Sigma0
288  BaryonWeight[1][2][0][1]=(1.-pspin_barion)*0.5;
289 
290  Baryon[1][2][0][2]=3214; // Sigma*0
291  BaryonWeight[1][2][0][2]=pspin_barion;
292 
293 //---------------------------------------
294  Baryon[1][2][1][0]=3222; // Sigma+
295  BaryonWeight[1][2][1][0]=1.-pspin_barion;
296 
297  Baryon[1][2][1][1]=3224; // Sigma*+
298  BaryonWeight[1][2][1][1]=pspin_barion;
299 
300 //---------------------------------------
301  Baryon[1][2][2][0]=3322; // Theta0
302  BaryonWeight[1][2][2][0]=1.-pspin_barion;
303 
304  Baryon[1][2][2][1]=3324; // Theta*0
305  BaryonWeight[1][2][2][1]=pspin_barion;
306 
307 //---------------------------------------
308 //---------------------------------------
309  Baryon[2][0][0][0]=3112; // Sigma-
310  BaryonWeight[2][0][0][0]=1.-pspin_barion;
311 
312  Baryon[2][0][0][1]=3114; // Sigma*-
313  BaryonWeight[2][0][0][1]=pspin_barion;
314 
315 //---------------------------------------
316  Baryon[2][0][1][0]=3122; // Lambda
317  BaryonWeight[2][0][1][0]=(1.-pspin_barion)*0.5;
318 
319  Baryon[2][0][1][1]=3212; // Sigma0
320  BaryonWeight[2][0][1][1]=(1.-pspin_barion)*0.5;
321 
322  Baryon[2][0][1][2]=3214; // Sigma*0
323  BaryonWeight[2][0][1][2]=pspin_barion;
324 
325 //---------------------------------------
326  Baryon[2][0][2][0]=3312; // Sigma-
327  BaryonWeight[2][0][2][0]=1.-pspin_barion;
328 
329  Baryon[2][0][2][1]=3314; // Sigma*-
330  BaryonWeight[2][0][2][1]=pspin_barion;
331 
332 //---------------------------------------
333  Baryon[2][1][0][0]=3122; // Lambda
334  BaryonWeight[2][1][0][0]=(1.-pspin_barion)*0.5;
335 
336  Baryon[2][1][0][1]=3212; // Sigma0
337  BaryonWeight[2][1][0][1]=(1.-pspin_barion)*0.5;
338 
339  Baryon[2][1][0][2]=3214; // Sigma*0
340  BaryonWeight[2][1][0][2]=pspin_barion;
341 
342 //---------------------------------------
343  Baryon[2][1][1][0]=3222; // Sigma+
344  BaryonWeight[2][1][1][0]=1.-pspin_barion;
345 
346  Baryon[2][1][1][1]=3224; // Sigma*+
347  BaryonWeight[2][1][1][1]=pspin_barion;
348 
349 //---------------------------------------
350  Baryon[2][1][2][0]=3322; // Theta0
351  BaryonWeight[2][1][2][0]=1.-pspin_barion;
352 
353  Baryon[2][1][2][1]=3324; // Theta*0
354  BaryonWeight[2][1][2][1]=pspin_barion;
355 
356 //---------------------------------------
357  Baryon[2][2][0][0]=3312; // Theta-
358  BaryonWeight[2][2][0][0]=1.-pspin_barion;
359 
360  Baryon[2][2][0][1]=3314; // Theta*-
361  BaryonWeight[2][2][0][1]=pspin_barion;
362 
363 //---------------------------------------
364  Baryon[2][2][1][0]=3322; // Theta0
365  BaryonWeight[2][2][1][0]=1.-pspin_barion;
366 
367  Baryon[2][2][1][1]=3324; // Theta*0
368  BaryonWeight[2][2][1][1]=pspin_barion;
369 
370 //---------------------------------------
371  Baryon[2][2][2][0]=3334; // Omega
372  BaryonWeight[2][2][2][0]=1.;
373 
374 //---------------------------------------
375  pspin_barion=pspin_barion_in;
376  /*
377  for(G4int i=0; i<3; i++)
378  { for(G4int j=0; j<3; j++)
379  { for(G4int k=0; k<3; k++)
380  { for(G4int l=0; l<4; l++)
381  { G4cout<<i<<" "<<j<<" "<<k<<" "<<l<<" "<<Baryon[i][j][k][l]<<G4endl;}
382  }
383  }
384  }
385  G4int Uzhi;
386  G4cin>>Uzhi;
387  */
388 
389  Prob_QQbar[0]=StrangeSuppress; // Probability of ddbar production
390  Prob_QQbar[1]=StrangeSuppress; // Probability of uubar production
391  Prob_QQbar[2]=1.0-2.*StrangeSuppress; // Probability of ssbar production Uzhi May 2014
392 
393  //A.R. 25-Jul-2012 : Coverity fix.
394  for ( G4int i=0 ; i<35 ; i++ ) {
395  FS_LeftHadron[i] = 0;
396  FS_RightHadron[i] = 0;
397  FS_Weight[i] = 0.0;
398  }
399  NumberOf_FS = 0;
400 
401 }
402 
403 // --------------------------------------------------------------
404 G4LundStringFragmentation::~G4LundStringFragmentation()
405 {}
406 
407 
408 //--------------------------------------------------------------------------------------
409 void G4LundStringFragmentation::SetMinimalStringMass(const G4FragmentingString * const string)
410 {
411  G4double EstimatedMass=0.;
412  G4int Number_of_quarks=0;
413  G4int Number_of_squarks=0;
414 
415  G4double StringM=string->Get4Momentum().mag();
416 
417  G4int Qleft =std::abs(string->GetLeftParton()->GetPDGEncoding());
418 
419 #ifdef debug_LUNDfragmentation // Uzhi 20.06.2014
420 // G4cout<<"MinStringMass// Input String mass "<<string->Get4Momentum().mag()<<" Qleft "<<Qleft<<G4endl;
421 #endif
422 
423  if( Qleft > 1000)
424  {
425  Number_of_quarks+=2;
426  G4int q1=Qleft/1000;
427  if( q1 < 3) {EstimatedMass +=Mass_of_light_quark;}
428  if( q1 > 2) {EstimatedMass +=Mass_of_heavy_quark; Number_of_squarks++;}
429 
430  G4int q2=(Qleft/100)%10;
431  if( q2 < 3) {EstimatedMass +=Mass_of_light_quark;}
432  if( q2 > 2) {EstimatedMass +=Mass_of_heavy_quark; Number_of_squarks++;}
433 // EstimatedMass +=Mass_of_string_junction;
434 //if((q1 > 2)||(q2 > 2)) EstimatedMass -= 120*MeV; // Uzhi April 2014 ???
435  }
436  else
437  {
438  Number_of_quarks++;
439  if( Qleft < 3) {EstimatedMass +=Mass_of_light_quark;}
440  if( Qleft > 2) {EstimatedMass +=Mass_of_heavy_quark; Number_of_squarks++;}
441  }
442 
443 #ifdef debug_LUNDfragmentation // Uzhi 20.06.2014
444 // G4cout<<"Min mass with Qleft "<<Qleft<<" "<<EstimatedMass<<G4endl;
445 #endif
446  G4int Qright=std::abs(string->GetRightParton()->GetPDGEncoding());
447  if( Qright > 1000)
448  {
449  Number_of_quarks+=2;
450  G4int q1=Qright/1000;
451  if( q1 < 3) {EstimatedMass +=Mass_of_light_quark;}
452  if( q1 > 2) {EstimatedMass +=Mass_of_heavy_quark; Number_of_squarks++;}
453 
454  G4int q2=(Qright/100)%10;
455  if( q2 < 3) {EstimatedMass +=Mass_of_light_quark;}
456  if( q2 > 2) {EstimatedMass +=Mass_of_heavy_quark; Number_of_squarks++;}
457  //EstimatedMass +=Mass_of_string_junction;
458 //if((q1 > 2)||(q2 > 2)) EstimatedMass -= 120*MeV; // Uzhi April 2014 ???
459  }
460  else
461  {
462  Number_of_quarks++;
463  if( Qright < 3) {EstimatedMass +=Mass_of_light_quark;}
464  if( Qright > 2) {EstimatedMass +=Mass_of_heavy_quark; Number_of_squarks++;}
465  }
466 
467 #ifdef debug_LUNDfragmentation // Uzhi 20.06.2014
468 // G4cout<<"Min mass with Qleft and Qright "<<Qright<<" "<<EstimatedMass<<G4endl;
469 #endif
470 
471  if(Number_of_quarks==2){EstimatedMass += 70.*MeV;} //100.*MeV;}
472 // if(Number_of_quarks==3){EstimatedMass += 20.*MeV;}
473  if(Number_of_quarks==3)
474  {
475  if(Number_of_squarks==0) {EstimatedMass += 740.*MeV;} // 700 Uzhi July 2014
476  if(Number_of_squarks==1) {EstimatedMass += 740.*MeV;} // 740 Uzhi Nov 2014
477  if(Number_of_squarks==2) {EstimatedMass += 400.*MeV;}
478  if(Number_of_squarks==3) {EstimatedMass += 382.*MeV;}
479  }
480  if(Number_of_quarks==4)
481  {
482  if((StringM > 1880.) && ( EstimatedMass < 2100)) {EstimatedMass = 2020.;}//1880.;}
483  // if((StringM > 1880.) && ( EstimatedMass < 2100)) {EstimatedMass = 2051.;}
484  else if((StringM > 2232.) && ( EstimatedMass < 2730)){EstimatedMass = 2570.;}
485  else if((StringM > 5130.) && ( EstimatedMass < 3450)){EstimatedMass = 5130.;}
486  else
487  {
488 // VU 30 May 2014 EstimatedMass -=2.*Mass_of_string_junction;
489  if(EstimatedMass <= 1600.*MeV){EstimatedMass-=200.*MeV;}
490  else {EstimatedMass+=100.*MeV;}
491  }
492  }
493 
494 #ifdef debug_LUNDfragmentation // Uzhi 20.06.2014
495 // G4cout<<"EstimatedMass -------------------- "<<EstimatedMass <<G4endl;
496 #endif
497  MinimalStringMass=EstimatedMass;
498  SetMinimalStringMass2(EstimatedMass);
499 }
500 
501 //--------------------------------------------------------------------------------------
502 void G4LundStringFragmentation::SetMinimalStringMass2(const G4double aValue)
503 {
504  MinimalStringMass2=aValue * aValue;
505 }
506 
507 //--------------------------------------------------------------------------------------
508 G4KineticTrackVector* G4LundStringFragmentation::FragmentString(const G4ExcitedString& theString)
509 {
510  // Can no longer modify Parameters for Fragmentation.
511  PastInitPhase=true;
512 
513  SetMassCut(160.*MeV); // For LightFragmentationTest it is required
514  // that no one pi-meson can be produced.
515 
516  G4FragmentingString aString(theString);
517  SetMinimalStringMass(&aString);
518 
519 #ifdef debug_LUNDfragmentation // Uzhi 20.06.2014
520  G4cout<<G4endl<<"LUND StringFragm: String Mass "
521  <<theString.Get4Momentum().mag()<<" Pz "
522  <<theString.Get4Momentum().pz()
523  <<"------------------------------------"<<G4endl;
524  G4cout<<"String ends Direct "<<theString.GetLeftParton()->GetPDGcode()<<" "
525  <<theString.GetRightParton()->GetPDGcode()<<" "
526  <<theString.GetDirection()<< G4endl;
527  G4cout<<"Left mom "<<theString.GetLeftParton()->Get4Momentum()<<G4endl;
528  G4cout<<"Right mom "<<theString.GetRightParton()->Get4Momentum()<<G4endl;
529  G4cout<<"Check for Fragmentation "<<G4endl;
530 #endif
531 
532  G4KineticTrackVector * LeftVector(0);
533 
534  if(!IsFragmentable(&aString)) // produce 1 hadron True ===============
535  {
536 #ifdef debug_LUNDfragmentation // Uzhi 20.06.2014
537  G4cout<<"Non fragmentable - the string is converted to one hadron "<<G4endl;
538 #endif
539  SetMassCut(10000.*MeV);
540  LeftVector=LightFragmentationTest(&theString);
541  SetMassCut(160.*MeV);
542 
543  LeftVector->operator[](0)->SetFormationTime(theString.GetTimeOfCreation());
544  LeftVector->operator[](0)->SetPosition(theString.GetPosition());
545 
546  if(LeftVector->size() > 1)
547  {
548  // 2 hadrons created from qq-qqbar are stored
549  LeftVector->operator[](1)->SetFormationTime(theString.GetTimeOfCreation());
550  LeftVector->operator[](1)->SetPosition(theString.GetPosition());
551  }
552  return LeftVector;
553  } // end of if(!IsFragmentable(&aString))
554 
555 #ifdef debug_LUNDfragmentation // Uzhi 20.06.2014
556  G4cout<<"The string will be fragmented. "<<G4endl;
557 #endif
558 
559  // The string can fragment. At least two particles can be produced.
560  LeftVector =new G4KineticTrackVector;
561  G4KineticTrackVector * RightVector=new G4KineticTrackVector;
562 
563  G4ExcitedString *theStringInCMS=CPExcited(theString);
564 #ifdef debug_LUNDfragmentation // Uzhi 20.06.2014
565  G4cout<<"CMS Left mom "<<theStringInCMS->GetLeftParton()->Get4Momentum()<<G4endl;
566  G4cout<<"CMS Right mom "<<theStringInCMS->GetRightParton()->Get4Momentum()<<G4endl;
567 #endif
568  G4LorentzRotation toCms=theStringInCMS->TransformToAlignedCms();
569 #ifdef debug_LUNDfragmentation // Uzhi 20.06.2014
570  G4cout<<"aligCMS Left mom "<<theStringInCMS->GetLeftParton()->Get4Momentum()<<G4endl;
571  G4cout<<"aligCMS Right mom "<<theStringInCMS->GetRightParton()->Get4Momentum()<<G4endl;
572  G4cout<<G4endl<<"LUND StringFragm: String Mass "
573  <<theStringInCMS->Get4Momentum().mag()<<" Pz Lab "
574  <<theStringInCMS->Get4Momentum().pz()
575  <<"------------------------------------"<<G4endl;
576  G4cout<<"String ends and Direction "<<theStringInCMS->GetLeftParton()->GetPDGcode()<<" "
577  <<theStringInCMS->GetRightParton()->GetPDGcode()<<" "
578  <<theStringInCMS->GetDirection()<< G4endl;
579 #endif
580  G4bool success = Loop_toFragmentString(theStringInCMS, LeftVector, RightVector);
581 
582  delete theStringInCMS;
583 
584  if ( ! success )
585  {
586  std::for_each(LeftVector->begin(), LeftVector->end(), DeleteKineticTrack());
587  LeftVector->clear();
588  std::for_each(RightVector->begin(), RightVector->end(), DeleteKineticTrack());
589  delete RightVector;
590  return LeftVector;
591  }
592 
593  // Join Left- and RightVector into LeftVector in correct order.
594  while(!RightVector->empty())
595  {
596  LeftVector->push_back(RightVector->back());
597  RightVector->erase(RightVector->end()-1);
598  }
599  delete RightVector;
600 
601  CalculateHadronTimePosition(theString.Get4Momentum().mag(), LeftVector);
602 
603  G4LorentzRotation toObserverFrame(toCms.inverse());
604 
605  G4double TimeOftheStringCreation=theString.GetTimeOfCreation();
606  G4ThreeVector PositionOftheStringCreation(theString.GetPosition());
607 
608  for(size_t C1 = 0; C1 < LeftVector->size(); C1++)
609  {
610  G4KineticTrack* Hadron = LeftVector->operator[](C1);
611  G4LorentzVector Momentum = Hadron->Get4Momentum();
612  //G4cout<<"Hadron "<<Hadron->GetDefinition()->GetParticleName()<<" "<<Momentum<<G4endl;
613  Momentum = toObserverFrame*Momentum;
614  Hadron->Set4Momentum(Momentum);
615 
616  G4LorentzVector Coordinate(Hadron->GetPosition(), Hadron->GetFormationTime());
617  Momentum = toObserverFrame*Coordinate;
618  Hadron->SetFormationTime(TimeOftheStringCreation + Momentum.e() - fermi/c_light);
619  G4ThreeVector aPosition(Momentum.vect());
620  Hadron->SetPosition(PositionOftheStringCreation+aPosition);
621  };
622 
623  return LeftVector;
624 }
625 
626 //----------------------------------------------------------------------------------
627 G4bool G4LundStringFragmentation::IsFragmentable(const G4FragmentingString * const string)
628 {
629  SetMinimalStringMass(string);
630 // return sqr(MinimalStringMass + WminLUND) < string->Get4Momentum().mag2();
631 // G4cout<<"MinM StrM "<<MinimalStringMass<<" "<< string->Get4Momentum().mag()<<G4endl;
632  return MinimalStringMass < string->Get4Momentum().mag();
633 }
634 
635 //----------------------------------------------------------------------------------------
636 G4bool G4LundStringFragmentation::StopFragmenting(const G4FragmentingString * const string)
637 {
638  SetMinimalStringMass(string);
639 
640  if (string->FourQuarkString())
641  {
642  return G4UniformRand() < std::exp(-0.0005*(string->Mass() - MinimalStringMass));
643  } else {
644  return G4UniformRand() < std::exp(-0.88e-6*(string->Mass()*string->Mass() -
645  MinimalStringMass*MinimalStringMass));
646  }
647 }
648 
649 //----------------------------------------------------------------------------------------------------------
650 G4bool G4LundStringFragmentation::SplitLast(G4FragmentingString * string,
651  G4KineticTrackVector * LeftVector,
652  G4KineticTrackVector * RightVector)
653 {
654  //... perform last cluster decay
655 #ifdef debug_LUNDfragmentation // Uzhi 20.06.2014
656  G4cout<<"Split last-----------------------------------------"<<G4endl;
657 #endif
658  G4LorentzVector Str4Mom=string->Get4Momentum();
659  G4ThreeVector ClusterVel=string->Get4Momentum().boostVector();
660  G4double StringMass=string->Mass();
661 
662  G4ParticleDefinition * LeftHadron(0), * RightHadron(0);
663 
664  NumberOf_FS=0;
665  for(G4int i=0; i<35; i++) {FS_Weight[i]=0.;}
666 
667 #ifdef debug_LUNDfragmentation // Uzhi 20.06.2014
668  G4cout<<"StrMass "<<StringMass<<" q "
669  <<string->GetLeftParton()->GetParticleName()<<" "
670  <<string->GetRightParton()->GetParticleName()<<G4endl;
671 #endif
672 
673  string->SetLeftPartonStable(); // to query quark contents..
674 
675  if (string->FourQuarkString() )
676  {
677  // The string is qq-qqbar type. Diquarks are on the string ends
678  //G4cout<<"The string is qq-qqbar type. Diquarks are on the string ends"<<G4endl;
679 
680  if(StringMass-MinimalStringMass < 0.)
681  {
682  if (! Diquark_AntiDiquark_belowThreshold_lastSplitting(string, LeftHadron, RightHadron) )
683  {
684  return false;
685  }
686  } else
687  {
688  Diquark_AntiDiquark_aboveThreshold_lastSplitting(string, LeftHadron, RightHadron);
689 
690  if(NumberOf_FS == 0) return false;
691 
692  G4int sampledState = SampleState();
693  if(string->GetLeftParton()->GetPDGEncoding() < 0)
694  {
695  LeftHadron =FS_LeftHadron[sampledState];
696  RightHadron=FS_RightHadron[sampledState];
697  } else
698  {
699  LeftHadron =FS_RightHadron[sampledState];
700  RightHadron=FS_LeftHadron[sampledState];
701  }
702  //G4cout<<"Selected "<<SampledState<<" "<<LeftHadron->GetParticleName()<<" "<<RightHadron->GetParticleName()<<G4endl;
703  }
704  } else
705  {
706  if (string->DecayIsQuark() && string->StableIsQuark() )
707  { //... there are quarks on cluster ends
708 #ifdef debug_LUNDfragmentation // Uzhi 20.06.2014
709  G4cout<<"Q Q string LastSplit"<<G4endl;
710 #endif
711  Quark_AntiQuark_lastSplitting(string, LeftHadron, RightHadron);
712  } else
713  { //... there is a Diquark on one of the cluster ends
714 #ifdef debug_LUNDfragmentation // Uzhi 20.06.2014
715  G4cout<<"DiQ Q string Last Split"<<G4endl;
716 #endif
717 
718  Quark_Diquark_lastSplitting(string, LeftHadron, RightHadron);
719  }
720 
721  if(NumberOf_FS == 0) return false;
722  G4int sampledState = SampleState();
723  LeftHadron =FS_LeftHadron[sampledState];
724  RightHadron=FS_RightHadron[sampledState];
725 
726 #ifdef debug_LUNDfragmentation // Uzhi 20.06.2014
727  G4cout<<"Selected LeftHad RightHad "<<sampledState<<" "
728  <<LeftHadron->GetParticleName()<<" "<<RightHadron->GetParticleName()<<G4endl;
729 #endif
730 
731  } // End of if(!string->FourQuarkString())
732 
733  G4LorentzVector LeftMom, RightMom;
734  G4ThreeVector Pos;
735 
736  Sample4Momentum(&LeftMom, LeftHadron->GetPDGMass(),
737  &RightMom, RightHadron->GetPDGMass(),
738  StringMass);
739 
740  LeftMom.boost(ClusterVel);
741  RightMom.boost(ClusterVel);
742 
743  LeftVector->push_back(new G4KineticTrack(LeftHadron, 0, Pos, LeftMom));
744  RightVector->push_back(new G4KineticTrack(RightHadron, 0, Pos, RightMom));
745 
746  return true;
747 
748 }
749 
750 //----------------------------------------------------------------------------------------------------------
751 void G4LundStringFragmentation::Sample4Momentum(G4LorentzVector* Mom, G4double Mass,
752  G4LorentzVector* AntiMom, G4double AntiMass,
753  G4double InitialMass)
754 {
755  // ------ Sampling of momenta of 2 last produced hadrons --------------------
756  G4ThreeVector Pt;
757  G4double MassMt2, AntiMassMt2;
758  G4double AvailablePz, AvailablePz2;
759 
760 #ifdef debug_LUNDfragmentation // Uzhi 20.06.2014
761  G4cout<<"Sampling of momenta of 2 last produced hadrons ----------------"<<G4endl;
762  G4cout<<"Masses "<<InitialMass<<" "<<Mass<<" "<<AntiMass<<G4endl;
763 #endif
764 
765  if((Mass > 930. || AntiMass > 930.)) //If there is a baryon
766  {
767  // ----------------- Isotropic decay ------------------------------------
768  G4double r_val = sqr(InitialMass*InitialMass - Mass*Mass - AntiMass*AntiMass) -
769  sqr(2.*Mass*AntiMass);
770  G4double Pabs = (r_val > 0.)? std::sqrt(r_val)/(2.*InitialMass) : 0;
771  //G4cout<<"P for isotr decay "<<Pabs<<G4endl;
772 
773  //... sample unit vector
774  G4double pz =1. - 2.*G4UniformRand();
775  G4double st = std::sqrt(1. - pz * pz)*Pabs;
776  G4double phi = 2.*pi*G4UniformRand();
777  G4double px = st*std::cos(phi);
778  G4double py = st*std::sin(phi);
779  pz *= Pabs;
780 
781  Mom->setPx(px); Mom->setPy(py); Mom->setPz(pz);
782  Mom->setE(std::sqrt(Pabs*Pabs + Mass*Mass));
783 
784  AntiMom->setPx(-px); AntiMom->setPy(-py); AntiMom->setPz(-pz);
785  AntiMom->setE (std::sqrt(Pabs*Pabs + AntiMass*AntiMass));
786  //G4int Uzhi; G4cin>>Uzhi;
787  }
788  else
789  {
790  do
791  {
792  // GF 22-May-09, limit sampled pt to allowed range
793 
794  G4double termD = InitialMass*InitialMass -Mass*Mass - AntiMass*AntiMass;
795  G4double termab = 4*sqr(Mass*AntiMass);
796  G4double termN = 2*termD + 4*Mass*Mass + 4*AntiMass*AntiMass;
797  G4double pt2max=(termD*termD - termab )/ termN ;
798  //G4cout<<"Anis "<<pt2max<<" "<<(termD*termD-termab)/(4.*InitialMass*InitialMass)<<G4endl;
799 
800  Pt=SampleQuarkPt(std::sqrt(pt2max)); Pt.setZ(0); G4double Pt2=Pt.mag2();
801  //G4cout<<"Sampl pt2 "<<Pt2<<G4endl;
802  MassMt2 = Mass * Mass + Pt2;
803  AntiMassMt2= AntiMass * AntiMass + Pt2;
804 
805  AvailablePz2= sqr(InitialMass*InitialMass - MassMt2 - AntiMassMt2) -
806  4.*MassMt2*AntiMassMt2;
807  }
808  while(AvailablePz2 < 0.); // GF will occur only for numerical precision problem with limit in sampled pt
809 
810  AvailablePz2 /=(4.*InitialMass*InitialMass);
811  AvailablePz = std::sqrt(AvailablePz2);
812 
813  G4double Px=Pt.getX();
814  G4double Py=Pt.getY();
815 
816  Mom->setPx(Px); Mom->setPy(Py); Mom->setPz(AvailablePz);
817  Mom->setE(std::sqrt(MassMt2+AvailablePz2));
818 
819  AntiMom->setPx(-Px); AntiMom->setPy(-Py); AntiMom->setPz(-AvailablePz);
820  AntiMom->setE (std::sqrt(AntiMassMt2+AvailablePz2));
821  }
822 }
823 
824 //-----------------------------------------------------------------------------
825 G4LorentzVector * G4LundStringFragmentation::SplitEandP(G4ParticleDefinition * pHadron,
826  G4FragmentingString * string,
827  G4FragmentingString * newString)
828 {
829  G4LorentzVector String4Momentum=string->Get4Momentum();
830  G4double StringMT2=string->MassT2();
831  G4double StringMT =std::sqrt(StringMT2);
832 
833  G4double HadronMass = pHadron->GetPDGMass();
834 
835  SetMinimalStringMass(newString);
836 
837 #ifdef debug_LUNDfragmentation // Uzhi 20.06.2014
838  G4cout<<G4endl<<"Start LUND SplitEandP "<<G4endl;
839  G4cout<<"String 4 mom, String M and Mt "<<String4Momentum<<" "<<String4Momentum.mag()<<" "<<std::sqrt(StringMT2)<<G4endl;
840  G4cout<<"Hadron "<<pHadron->GetParticleName()<<G4endl;
841  G4cout<<"HadM MinimalStringMassLeft StringM hM+sM "<<HadronMass<<" "<<MinimalStringMass<<" "
842  <<String4Momentum.mag()<<" "<<HadronMass+MinimalStringMass<<G4endl;
843 #endif
844 
845  if(HadronMass + MinimalStringMass > string->Mass())
846  {
847 #ifdef debug_LUNDfragmentation // Uzhi 20.06.2014
848  G4cout<<"Mass of the string is not sufficient to produce the hadron!"<<G4endl;
849 #endif
850  return 0;
851  }// have to start all over!
852 
853  String4Momentum.setPz(0.);
854  G4ThreeVector StringPt=String4Momentum.vect();
855 
856  // calculate and assign hadron transverse momentum component HadronPx and HadronPy
857  G4ThreeVector HadronPt , RemSysPt;
858  G4double HadronMassT2, ResidualMassT2;
859 
860  //... sample Pt of the hadron
861  G4int attempt=0;
862  do
863  {
864  attempt++; if(attempt > StringLoopInterrupt) return 0;
865 
866  HadronPt =SampleQuarkPt() + string->DecayPt();
867  HadronPt.setZ(0);
868  RemSysPt = StringPt - HadronPt;
869 
870  HadronMassT2 = sqr(HadronMass) + HadronPt.mag2();
871  ResidualMassT2=sqr(MinimalStringMass) + RemSysPt.mag2();
872 
873  } while(std::sqrt(HadronMassT2) + std::sqrt(ResidualMassT2) > StringMT);
874 
875  //... sample z to define hadron longitudinal momentum and energy
876  //... but first check the available phase space
877 
878  G4double Pz2 = (sqr(StringMT2 - HadronMassT2 - ResidualMassT2) -
879  4*HadronMassT2 * ResidualMassT2)/4./StringMT2;
880 
881  if(Pz2 < 0 ) {return 0;} // have to start all over!
882 
883  //... then compute allowed z region z_min <= z <= z_max
884 
885  G4double Pz = std::sqrt(Pz2);
886  G4double zMin = (std::sqrt(HadronMassT2+Pz2) - Pz)/std::sqrt(StringMT2);
887 // G4double zMin = (std::sqrt(HadronMassT2+Pz2) - 0.)/std::sqrt(StringMT2); // Uzhi 2014
888  G4double zMax = (std::sqrt(HadronMassT2+Pz2) + Pz)/std::sqrt(StringMT2);
889 
890  if (zMin >= zMax) return 0; // have to start all over!
891 
892  G4double z = GetLightConeZ(zMin, zMax,
893  string->GetDecayParton()->GetPDGEncoding(), pHadron,
894  HadronPt.x(), HadronPt.y());
895 
896  //... now compute hadron longitudinal momentum and energy
897  // longitudinal hadron momentum component HadronPz
898 
899  HadronPt.setZ(0.5* string->GetDecayDirection() *
900  (z * string->LightConeDecay() -
901  HadronMassT2/(z * string->LightConeDecay())));
902  G4double HadronE = 0.5* (z * string->LightConeDecay() +
903  HadronMassT2/(z * string->LightConeDecay()));
904 
905  G4LorentzVector * a4Momentum= new G4LorentzVector(HadronPt,HadronE);
906 #ifdef debug_LUNDfragmentation // Uzhi 20.06.2014
907  G4cout<<"string->LightConeDecay() "<<string->LightConeDecay()<<G4endl;
908  G4cout<<"HadronPt,HadronE "<<HadronPt<<" "<<HadronE<<G4endl;
909  G4cout<<"String4Momentum "<<String4Momentum<<G4endl;
910 //G4int Uzhi; G4cin>>Uzhi;
911  G4cout<<"Out of LUND SplitEandP "<<G4endl;
912 #endif
913 
914  return a4Momentum;
915 }
916 
917 //-----------------------------------------------------------------------------------------
918 G4double G4LundStringFragmentation::GetLightConeZ(G4double zmin, G4double zmax,
919  G4int PDGEncodingOfDecayParton,
920  G4ParticleDefinition* pHadron,
921  G4double Px, G4double Py)
922 {
923  G4double Mass = pHadron->GetPDGMass();
924 // G4int HadronEncoding=std::abs(pHadron->GetPDGEncoding());
925 
926  G4double Mt2 = Px*Px + Py*Py + Mass*Mass;
927 
928  G4double alund;
929  G4double zOfMaxyf(0.), maxYf(1.), z(0.), yf(1.);
930  if(std::abs(PDGEncodingOfDecayParton) < 1000)
931  { // ---------------- Quark fragmentation ----------------------
932  alund=0.7/GeV/GeV;
933  // If blund get restored, you MUST adapt the calculation of zOfMaxyf.
934  // const G4double blund = 1;
935 
936  zOfMaxyf=alund*Mt2/(alund*Mt2 + 1.);
937  maxYf=(1-zOfMaxyf)/zOfMaxyf * std::exp(-alund*Mt2/zOfMaxyf);
938 
939  do
940  {
941  z = zmin + G4UniformRand()*(zmax-zmin);
942  yf = (1-z)/z * std::exp(-alund*Mt2/z);
943  // yf = std::pow(1.0-z,blund)/z*std::exp(-alund*Mt2/z
944  }
945  while (G4UniformRand()*maxYf > yf);
946  return z;
947  }
948 
949  if(std::abs(PDGEncodingOfDecayParton) > 1000) // Uzhi Sept. 2014
950  {
951 /*
952  if(HadronEncoding < 3000)
953  {
954  maxYf=(zmax-zmin);
955  do
956  {
957  z = zmin + G4UniformRand()*(zmax-zmin);
958  //yf=std::exp(-sqr(z-Zc)/2/sqr(0.28)); // 0.42 0.632 0.28 a'la UrQMD
959  yf =(z-zmin);
960  }
961  while (G4UniformRand()*maxYf > yf);
962  }
963  else
964  { // Strange baryons
965 */
966  z = zmin + G4UniformRand()*(zmax-zmin);
967 // }
968  }
969 
970  return z;
971 
972 }
973 
974 //------------------------------------------------------------------------
975 G4double G4LundStringFragmentation::lambda(G4double S, G4double m1_Sqr, G4double m2_Sqr)
976 {
977  G4double lam = sqr(S - m1_Sqr - m2_Sqr) - 4.*m1_Sqr*m2_Sqr;
978  return lam;
979 }
980 
981 
982 //------------------------------------------------------------------------
983 //------------------------------------------------------------------------
984 // Internal methods introduced to improve the code structure (AR Nov 2011)
985 //------------------------------------------------------------------------
986 //------------------------------------------------------------------------
987 
988 //----------------------------------------------------------------------------------------
989 G4bool G4LundStringFragmentation::Loop_toFragmentString(G4ExcitedString * & theStringInCMS,
990  G4KineticTrackVector * & LeftVector,
991  G4KineticTrackVector * & RightVector)
992 {
993 #ifdef debug_LUNDfragmentation // Uzhi 20.06.2014
994  G4cout<<"Loop_toFrag "<<theStringInCMS->GetLeftParton()->GetPDGcode()<<" "
995  <<theStringInCMS->GetRightParton()->GetPDGcode()<<" "
996  <<theStringInCMS->GetDirection()<< G4endl;
997 #endif
998 
999  G4bool final_success=false;
1000  G4bool inner_success=true;
1001  G4int attempt=0;
1002  while ( ! final_success && attempt++ < StringLoopInterrupt )
1003  { // If the string fragmentation do not be happend, repeat the fragmentation.
1004 
1005  G4FragmentingString *currentString=new G4FragmentingString(*theStringInCMS);
1006  //G4cout<<"Main loop start whilecounter "<<attempt<<G4endl;
1007  // Cleaning up the previously produced hadrons
1008  std::for_each(LeftVector->begin(), LeftVector->end(), DeleteKineticTrack());
1009  LeftVector->clear();
1010  std::for_each(RightVector->begin(), RightVector->end(), DeleteKineticTrack());
1011  RightVector->clear();
1012 
1013  // Main fragmentation loop until the string will not be able to fragment
1014  inner_success=true; // set false on failure.
1015  while (! StopFragmenting(currentString) )
1016  { // Split current string into hadron + new string
1017 #ifdef debug_LUNDfragmentation // Uzhi 20.06.2014
1018  G4cout<<"The string can fragment. "<<G4endl;;
1019 //G4cout<<"1 "<<currentString->GetDecayDirection()<<G4endl;
1020 #endif
1021  G4FragmentingString *newString=0; // used as output from SplitUp.
1022  G4KineticTrack * Hadron=Splitup(currentString,newString);
1023  if ( Hadron != 0 ) // Store the hadron
1024  {
1025 #ifdef debug_LUNDfragmentation // Uzhi 20.06.2014
1026  G4cout<<"Hadron prod at fragm. "<<Hadron->GetDefinition()->GetParticleName()<<G4endl;
1027 //G4cout<<"2 "<<currentString->GetDecayDirection()<<G4endl;
1028 #endif
1029 
1030  if ( currentString->GetDecayDirection() > 0 )
1031  {
1032  LeftVector->push_back(Hadron);
1033  } else
1034  {
1035  RightVector->push_back(Hadron);
1036  }
1037  delete currentString;
1038  currentString=newString;
1039  }
1040  };
1041  // Split remaining string into 2 final hadrons.
1042 #ifdef debug_LUNDfragmentation // Uzhi 20.06.2014
1043  G4cout<<"Split remaining string into 2 final hadrons."<<G4endl;
1044 #endif
1045 // Uzhi Sept 2014
1046  if ( inner_success && SplitLast(currentString, LeftVector, RightVector) )
1047  {
1048  final_success=true;
1049  }
1050 //
1051 //final_success=true; // Uzhi Sept 2014
1052  delete currentString;
1053  } // End of the loop where we try to fragment the string.
1054  return final_success;
1055 }
1056 
1057 //----------------------------------------------------------------------------------------
1058 G4bool G4LundStringFragmentation::
1059 Diquark_AntiDiquark_belowThreshold_lastSplitting(G4FragmentingString * & string,
1060  G4ParticleDefinition * & LeftHadron,
1061  G4ParticleDefinition * & RightHadron)
1062 {
1063  G4double StringMass = string->Mass();
1064  G4int cClusterInterrupt = 0;
1065  do
1066  {
1067  //G4cout<<"cClusterInterrupt "<<cClusterInterrupt<<G4endl;
1068  if (cClusterInterrupt++ >= ClusterLoopInterrupt)
1069  {
1070  return false;
1071  }
1072 
1073  G4int LeftQuark1= string->GetLeftParton()->GetPDGEncoding()/1000;
1074  G4int LeftQuark2=(string->GetLeftParton()->GetPDGEncoding()/100)%10;
1075 
1076  G4int RightQuark1= string->GetRightParton()->GetPDGEncoding()/1000;
1077  G4int RightQuark2=(string->GetRightParton()->GetPDGEncoding()/100)%10;
1078 
1079  if(G4UniformRand()<0.5)
1080  {
1081  LeftHadron =hadronizer->Build(FindParticle( LeftQuark1),
1082  FindParticle(RightQuark1));
1083  RightHadron=hadronizer->Build(FindParticle( LeftQuark2),
1084  FindParticle(RightQuark2));
1085  } else
1086  {
1087  LeftHadron =hadronizer->Build(FindParticle( LeftQuark1),
1088  FindParticle(RightQuark2));
1089  RightHadron=hadronizer->Build(FindParticle( LeftQuark2),
1090  FindParticle(RightQuark1));
1091  }
1092 
1093  //... repeat procedure, if mass of cluster is too low to produce hadrons
1094  //... ClusterMassCut = 0.15*GeV model parameter
1095  }
1096  while ((StringMass <= LeftHadron->GetPDGMass() + RightHadron->GetPDGMass()));
1097 
1098  return true;
1099 }
1100 
1101 //----------------------------------------------------------------------------------------
1102 G4bool G4LundStringFragmentation::
1103 Diquark_AntiDiquark_aboveThreshold_lastSplitting(G4FragmentingString * & string,
1104  G4ParticleDefinition * & LeftHadron,
1105  G4ParticleDefinition * & RightHadron)
1106 {
1107  // StringMass-MinimalStringMass > 0. Creation of 2 baryons is possible ----
1108 
1109  G4double StringMass = string->Mass();
1110  G4double StringMassSqr= sqr(StringMass);
1111  G4ParticleDefinition * Di_Quark;
1112  G4ParticleDefinition * Anti_Di_Quark;
1113 
1114  if(string->GetLeftParton()->GetPDGEncoding() < 0)
1115  {
1116  Anti_Di_Quark =string->GetLeftParton();
1117  Di_Quark=string->GetRightParton();
1118  } else
1119  {
1120  Anti_Di_Quark =string->GetRightParton();
1121  Di_Quark=string->GetLeftParton();
1122  }
1123 
1124  G4int IDAnti_di_quark =Anti_Di_Quark->GetPDGEncoding();
1125  G4int AbsIDAnti_di_quark =std::abs(IDAnti_di_quark);
1126  G4int IDdi_quark =Di_Quark->GetPDGEncoding();
1127  G4int AbsIDdi_quark =std::abs(IDdi_quark);
1128 
1129  G4int ADi_q1=AbsIDAnti_di_quark/1000;
1130  G4int ADi_q2=(AbsIDAnti_di_quark-ADi_q1*1000)/100;
1131 
1132  G4int Di_q1=AbsIDdi_quark/1000;
1133  G4int Di_q2=(AbsIDdi_quark-Di_q1*1000)/100;
1134 
1135  NumberOf_FS=0;
1136  for(G4int ProdQ=1; ProdQ < 4; ProdQ++)
1137  {
1138  G4int StateADiQ=0;
1139  do // while(Meson[AbsIDquark-1][ProdQ-1][StateQ]<>0);
1140  {
1141  LeftHadron=G4ParticleTable::GetParticleTable()->FindParticle(
1142  -Baryon[ADi_q1-1][ADi_q2-1][ProdQ-1][StateADiQ]);
1143  G4double LeftHadronMass=LeftHadron->GetPDGMass();
1144 
1145  //G4cout<<"Anti Bar "<<LeftHadron->GetParticleName()<<G4endl;
1146 
1147  G4int StateDiQ=0;
1148  do // while(Baryon[Di_q1-1][Di_q2-1][ProdQ-1][StateDiQ]<>0);
1149  {
1150  RightHadron=G4ParticleTable::GetParticleTable()->FindParticle(
1151  +Baryon[Di_q1-1][Di_q2-1][ProdQ-1][StateDiQ]);
1152  G4double RightHadronMass=RightHadron->GetPDGMass();
1153 
1154  if(StringMass > LeftHadronMass + RightHadronMass)
1155  {
1156  if ( NumberOf_FS > 34 ) {
1157  G4ExceptionDescription ed;
1158  ed << " NumberOf_FS exceeds its limit: NumberOf_FS=" << NumberOf_FS << G4endl;
1159  G4Exception( "G4LundStringFragmentation::Diquark_AntiDiquark_aboveThreshold_lastSplitting ",
1160  "HAD_LUND_001", JustWarning, ed );
1161  NumberOf_FS = 34;
1162  }
1163 
1164  G4double FS_Psqr=lambda(StringMassSqr,sqr(LeftHadronMass),
1165  sqr(RightHadronMass));
1166  //FS_Psqr=1.;
1167  FS_Weight[NumberOf_FS]=std::sqrt(FS_Psqr)*FS_Psqr*
1168  BaryonWeight[ADi_q1-1][ADi_q2-1][ProdQ-1][StateADiQ]*
1169  BaryonWeight[Di_q1-1][Di_q2-1][ProdQ-1][StateDiQ]*
1170  Prob_QQbar[ProdQ-1];
1171 
1172  FS_LeftHadron[NumberOf_FS] = LeftHadron;
1173  FS_RightHadron[NumberOf_FS]= RightHadron;
1174 
1175  NumberOf_FS++;
1176  } // End of if(StringMass > LeftHadronMass + RightHadronMass)
1177 
1178  StateDiQ++;
1179 
1180  } while(Baryon[Di_q1-1][Di_q2-1][ProdQ-1][StateDiQ]!=0);
1181 
1182  StateADiQ++;
1183  } while(Baryon[ADi_q1-1][ADi_q2-1][ProdQ-1][StateADiQ]!=0);
1184  } // End of for(G4int ProdQ=1; ProdQ < 4; ProdQ++)
1185 
1186  return true;
1187 }
1188 
1189 //----------------------------------------------------------------------------------------
1190 G4bool G4LundStringFragmentation::
1191 Quark_Diquark_lastSplitting(G4FragmentingString * & string,
1192  G4ParticleDefinition * & LeftHadron,
1193  G4ParticleDefinition * & RightHadron)
1194 {
1195  G4double StringMass = string->Mass();
1196  G4double StringMassSqr= sqr(StringMass);
1197 
1198  G4ParticleDefinition * Di_Quark;
1199  G4ParticleDefinition * Quark;
1200 
1201  if(string->GetLeftParton()->GetParticleSubType()== "quark")
1202  {
1203  Quark =string->GetLeftParton();
1204  Di_Quark=string->GetRightParton();
1205  } else
1206  {
1207  Quark =string->GetRightParton();
1208  Di_Quark=string->GetLeftParton();
1209  }
1210 
1211  G4int IDquark =Quark->GetPDGEncoding();
1212  G4int AbsIDquark =std::abs(IDquark);
1213  G4int IDdi_quark =Di_Quark->GetPDGEncoding();
1214  G4int AbsIDdi_quark=std::abs(IDdi_quark);
1215  G4int Di_q1=AbsIDdi_quark/1000;
1216  G4int Di_q2=(AbsIDdi_quark-Di_q1*1000)/100;
1217 
1218  G4int SignDiQ= 1;
1219  if(IDdi_quark < 0) SignDiQ=-1;
1220 
1221  NumberOf_FS=0;
1222  for(G4int ProdQ=1; ProdQ < 4; ProdQ++)
1223  {
1224  G4int SignQ;
1225  if(IDquark > 0)
1226  { SignQ=-1;
1227  if(IDquark == 2) SignQ= 1;
1228  if((IDquark == 1) && (ProdQ == 3)) SignQ= 1; // K0
1229  if((IDquark == 3) && (ProdQ == 1)) SignQ=-1; // K0bar
1230  } else
1231  {
1232  SignQ= 1;
1233  if(IDquark == -2) SignQ=-1;
1234  if((IDquark ==-1) && (ProdQ == 3)) SignQ=-1; // K0bar
1235  if((IDquark ==-3) && (ProdQ == 1)) SignQ= 1; // K0
1236  }
1237 
1238  if(AbsIDquark == ProdQ) SignQ= 1;
1239 
1240  //G4cout<<G4endl;
1241  //G4cout<<"Insert QQbar "<<ProdQ<<" Sign "<<SignQ<<G4endl;
1242 
1243  G4int StateQ=0;
1244  do // while(Meson[AbsIDquark-1][ProdQ-1][StateQ]<>0);
1245  {
1246  LeftHadron=G4ParticleTable::GetParticleTable()->FindParticle(SignQ*
1247  Meson[AbsIDquark-1][ProdQ-1][StateQ]);
1248  G4double LeftHadronMass=LeftHadron->GetPDGMass();
1249 
1250  G4int StateDiQ=0;
1251  do // while(Baryon[Di_q1-1][Di_q2-1][ProdQ-1][StateDiQ]<>0);
1252  {
1253  RightHadron=G4ParticleTable::GetParticleTable()->FindParticle(SignDiQ*
1254  Baryon[Di_q1-1][Di_q2-1][ProdQ-1][StateDiQ]);
1255  G4double RightHadronMass=RightHadron->GetPDGMass();
1256 
1257  if(StringMass > LeftHadronMass + RightHadronMass)
1258  {
1259  if ( NumberOf_FS > 34 ) {
1260  G4ExceptionDescription ed;
1261  ed << " NumberOf_FS exceeds its limit: NumberOf_FS=" << NumberOf_FS << G4endl;
1262  G4Exception( "G4LundStringFragmentation::Quark_Diquark_lastSplitting ",
1263  "HAD_LUND_002", JustWarning, ed );
1264  NumberOf_FS = 34;
1265  }
1266 
1267  G4double FS_Psqr=lambda(StringMassSqr,sqr(LeftHadronMass),
1268  sqr(RightHadronMass));
1269  FS_Weight[NumberOf_FS]=std::sqrt(FS_Psqr)*
1270  MesonWeight[AbsIDquark-1][ProdQ-1][StateQ]*
1271  BaryonWeight[Di_q1-1][Di_q2-1][ProdQ-1][StateDiQ]*
1272  Prob_QQbar[ProdQ-1];
1273 
1274  FS_LeftHadron[NumberOf_FS] = LeftHadron;
1275  FS_RightHadron[NumberOf_FS]= RightHadron;
1276 
1277  NumberOf_FS++;
1278  } // End of if(StringMass > LeftHadronMass + RightHadronMass)
1279 
1280  StateDiQ++;
1281 
1282  } while(Baryon[Di_q1-1][Di_q2-1][ProdQ-1][StateDiQ]!=0);
1283 
1284  StateQ++;
1285  } while(Meson[AbsIDquark-1][ProdQ-1][StateQ]!=0);
1286  } // End of for(G4int ProdQ=1; ProdQ < 4; ProdQ++)
1287 
1288  return true;
1289 }
1290 
1291 //----------------------------------------------------------------------------------------
1292 G4bool G4LundStringFragmentation::
1293 Quark_AntiQuark_lastSplitting(G4FragmentingString * & string,
1294  G4ParticleDefinition * & LeftHadron,
1295  G4ParticleDefinition * & RightHadron)
1296 {
1297  G4double StringMass = string->Mass();
1298  G4double StringMassSqr= sqr(StringMass);
1299 
1300  G4ParticleDefinition * Quark;
1301  G4ParticleDefinition * Anti_Quark;
1302 
1303  if(string->GetLeftParton()->GetPDGEncoding()>0)
1304  {
1305  Quark =string->GetLeftParton();
1306  Anti_Quark=string->GetRightParton();
1307  } else
1308  {
1309  Quark =string->GetRightParton();
1310  Anti_Quark=string->GetLeftParton();
1311  }
1312 
1313  G4int IDquark =Quark->GetPDGEncoding();
1314  G4int AbsIDquark =std::abs(IDquark);
1315  G4int IDanti_quark =Anti_Quark->GetPDGEncoding();
1316  G4int AbsIDanti_quark=std::abs(IDanti_quark);
1317 
1318  NumberOf_FS=0;
1319  for(G4int ProdQ=1; ProdQ < 4; ProdQ++)
1320  {
1321  G4int SignQ=-1;
1322  if(IDquark == 2) SignQ= 1;
1323  if((IDquark == 1) && (ProdQ == 3)) SignQ= 1; // K0
1324  if((IDquark == 3) && (ProdQ == 1)) SignQ=-1; // K0bar
1325  if(IDquark == ProdQ) SignQ= 1;
1326 
1327  G4int SignAQ= 1;
1328  if(IDanti_quark == -2) SignAQ=-1;
1329  if((IDanti_quark ==-1) && (ProdQ == 3)) SignAQ=-1; // K0bar
1330  if((IDanti_quark ==-3) && (ProdQ == 1)) SignAQ= 1; // K0
1331  if(AbsIDanti_quark == ProdQ) SignAQ= 1;
1332 
1333  G4int StateQ=0;
1334  do // while(Meson[AbsIDquark-1][ProdQ-1][StateQ]<>0);
1335  {
1336  LeftHadron=G4ParticleTable::GetParticleTable()->FindParticle(SignQ*
1337  Meson[AbsIDquark-1][ProdQ-1][StateQ]);
1338  G4double LeftHadronMass=LeftHadron->GetPDGMass();
1339 
1340  G4int StateAQ=0;
1341  do // while(Meson[AbsIDanti_quark-1][ProdQ-1][StateAQ]<>0);
1342  {
1343  RightHadron=G4ParticleTable::GetParticleTable()->FindParticle(SignAQ*
1344  Meson[AbsIDanti_quark-1][ProdQ-1][StateAQ]);
1345  G4double RightHadronMass=RightHadron->GetPDGMass();
1346 
1347  if(StringMass > LeftHadronMass + RightHadronMass)
1348  {
1349  if ( NumberOf_FS > 34 ) {
1350  G4ExceptionDescription ed;
1351  ed << " NumberOf_FS exceeds its limit: NumberOf_FS=" << NumberOf_FS << G4endl;
1352  G4Exception( "G4LundStringFragmentation::Quark_AntiQuark_lastSplitting ",
1353  "HAD_LUND_003", JustWarning, ed );
1354  NumberOf_FS = 34;
1355  }
1356 
1357  G4double FS_Psqr=lambda(StringMassSqr,sqr(LeftHadronMass),
1358  sqr(RightHadronMass));
1359  //FS_Psqr=1.;
1360  FS_Weight[NumberOf_FS]=std::sqrt(FS_Psqr)*
1361  MesonWeight[AbsIDquark-1][ProdQ-1][StateQ]*
1362  MesonWeight[AbsIDanti_quark-1][ProdQ-1][StateAQ]*
1363  Prob_QQbar[ProdQ-1];
1364 
1365  if(string->GetLeftParton()->GetPDGEncoding()>0)
1366  {
1367  FS_LeftHadron[NumberOf_FS] = RightHadron;
1368  FS_RightHadron[NumberOf_FS]= LeftHadron;
1369  } else
1370  {
1371  FS_LeftHadron[NumberOf_FS] = LeftHadron;
1372  FS_RightHadron[NumberOf_FS]= RightHadron;
1373  }
1374  NumberOf_FS++;
1375 
1376  } // End of if(StringMass > LeftHadronMass + RightHadronMass)
1377 
1378  StateAQ++;
1379  } while(Meson[AbsIDanti_quark-1][ProdQ-1][StateAQ]!=0);
1380 
1381  StateQ++;
1382  } while(Meson[AbsIDquark-1][ProdQ-1][StateQ]!=0);
1383  } // End of for(G4int ProdQ=1; ProdQ < 4; ProdQ++)
1384 
1385  return true;
1386 }
1387 
1388 //----------------------------------------------------------------------------------------------------------
1389 G4int G4LundStringFragmentation::SampleState(void)
1390 {
1391  if ( NumberOf_FS > 34 ) {
1392  G4ExceptionDescription ed;
1393  ed << " NumberOf_FS exceeds its limit: NumberOf_FS=" << NumberOf_FS << G4endl;
1394  G4Exception( "G4LundStringFragmentation::SampleState ", "HAD_LUND_004", JustWarning, ed );
1395  NumberOf_FS = 34;
1396  }
1397 
1398  G4double SumWeights=0.;
1399 
1400  for(G4int i=0; i<NumberOf_FS; i++) {SumWeights+=FS_Weight[i];}// G4cout<<i<<" "<<FS_Weight[i]<<G4endl;}
1401 
1402  G4double ksi=G4UniformRand();
1403  G4double Sum=0.;
1404  G4int indexPosition = 0;
1405 
1406  for(G4int i=0; i<NumberOf_FS; i++)
1407  {
1408  Sum+=(FS_Weight[i]/SumWeights);
1409  indexPosition=i;
1410  if(Sum >= ksi) break;
1411  }
1412  return indexPosition;
1413 }
1414 // Uzhi June 2014 Insert from G4ExcitedStringDecay.cc
1415 //-----------------------------------------------------------------------------
1416 
1417 G4ParticleDefinition * G4LundStringFragmentation::DiQuarkSplitup(
1418  G4ParticleDefinition* decay,
1419  G4ParticleDefinition *&created)
1420 {
1421  //... can Diquark break or not?
1422  if (G4UniformRand() < DiquarkBreakProb ){
1423  //... Diquark break
1424 
1425  G4int stableQuarkEncoding = decay->GetPDGEncoding()/1000;
1426  G4int decayQuarkEncoding = (decay->GetPDGEncoding()/100)%10;
1427  if (G4UniformRand() < 0.5)
1428  {
1429  G4int Swap = stableQuarkEncoding;
1430  stableQuarkEncoding = decayQuarkEncoding;
1431  decayQuarkEncoding = Swap;
1432  }
1433 
1434  G4int IsParticle=(decayQuarkEncoding>0) ? -1 : +1;
1435  // if we have a quark, we need antiquark)
1436 
1437 //G4cout<<"GetStrangeSuppress() "<<GetStrangeSuppress()<<G4endl;
1438 //G4int Uzhi; G4cin>>Uzhi;
1439 
1440 //G4double StrSup=GetStrangeSuppress();
1441 //StrangeSuppress=0.34;
1442  pDefPair QuarkPair = CreatePartonPair(IsParticle,false); // no diquarks wanted
1443 //StrangeSuppress=StrSup;
1444  //... Build new Diquark
1445  G4int QuarkEncoding=QuarkPair.second->GetPDGEncoding();
1446  G4int i10 = std::max(std::abs(QuarkEncoding), std::abs(stableQuarkEncoding));
1447  G4int i20 = std::min(std::abs(QuarkEncoding), std::abs(stableQuarkEncoding));
1448  G4int spin = (i10 != i20 && G4UniformRand() <= 0.5)? 1 : 3;
1449  G4int NewDecayEncoding = -1*IsParticle*(i10 * 1000 + i20 * 100 + spin);
1450  created = FindParticle(NewDecayEncoding);
1451  G4ParticleDefinition * decayQuark=FindParticle(decayQuarkEncoding);
1452  G4ParticleDefinition * had=hadronizer->Build(QuarkPair.first, decayQuark);
1453  return had;
1454 // return hadronizer->Build(QuarkPair.first, decayQuark);
1455 
1456  } else {
1457  //... Diquark does not break
1458 
1459  G4int IsParticle=(decay->GetPDGEncoding()>0) ? +1 : -1;
1460  // if we have a diquark, we need quark)
1461 
1462 G4double StrSup=GetStrangeSuppress(); // Uzhi Sept. 2014
1463 StrangeSuppress=0.43; //0.42 0.38 // Uzhi Sept. 2014
1464  pDefPair QuarkPair = CreatePartonPair(IsParticle,false); // no diquarks wanted
1465 StrangeSuppress=StrSup; // Uzhi Sept. 2014
1466 
1467  created = QuarkPair.second;
1468 
1469  G4ParticleDefinition * had=hadronizer->Build(QuarkPair.first, decay);
1470  return had;
1471  }
1472 }
1473 // Uzhi June 2014 End of the inserting
G4LundStringFragmentation::FragmentString
virtual G4KineticTrackVector * FragmentString(const G4ExcitedString &theString)
Definition: G4LundStringFragmentation.cc:508
G4FragmentingString::GetRightParton
G4ParticleDefinition * GetRightParton(void) const
Definition: G4FragmentingString.hh:143
G4FragmentingString.hh
MeV
static const double MeV
Definition: G4SIunits.hh:193
G4ParticleTable::FindParticle
G4ParticleDefinition * FindParticle(G4int PDGEncoding)
Definition: G4ParticleTable.cc:544
G4Parton::GetPDGcode
G4int GetPDGcode() const
Definition: G4Parton.hh:124
G4LundStringFragmentation::Quark_AntiQuark_lastSplitting
G4bool Quark_AntiQuark_lastSplitting(G4FragmentingString *&string, G4ParticleDefinition *&LeftHadron, G4ParticleDefinition *&RightHadron)
Definition: G4LundStringFragmentation.cc:1293
G4ExceptionDescription
std::ostringstream G4ExceptionDescription
Definition: globals.hh:76
G4ThreeVector
CLHEP::Hep3Vector G4ThreeVector
Definition: G4ThreeVector.hh:42
G4LundStringFragmentation::SplitEandP
virtual G4LorentzVector * SplitEandP(G4ParticleDefinition *pHadron, G4FragmentingString *string, G4FragmentingString *newString)
Definition: G4LundStringFragmentation.cc:825
G4LorentzRotation
CLHEP::HepLorentzRotation G4LorentzRotation
Definition: G4LorentzRotation.hh:35
G4KineticTrack::GetPosition
const G4ThreeVector & GetPosition() const
Definition: G4KineticTrack.hh:237
G4VLongitudinalStringDecay::CPExcited
G4ExcitedString * CPExcited(const G4ExcitedString &string)
Definition: G4VLongitudinalStringDecay.cc:305
z
G4double z
Definition: TRTMaterials.hh:39
G4LundStringFragmentation::IsFragmentable
virtual G4bool IsFragmentable(const G4FragmentingString *const string)
Definition: G4LundStringFragmentation.cc:627
G4INCL::Math::pi
const G4double pi
Definition: G4INCLGlobals.hh:68
G4KineticTrack::SetFormationTime
void SetFormationTime(G4double aFormationTime)
Definition: G4KineticTrack.hh:230
G4VLongitudinalStringDecay::SetStrangenessSuppression
void SetStrangenessSuppression(G4double aValue)
Definition: G4VLongitudinalStringDecay.cc:553
G4Parton::Get4Momentum
const G4LorentzVector & Get4Momentum() const
Definition: G4Parton.hh:140
G4HadronBuilder::Build
G4ParticleDefinition * Build(G4ParticleDefinition *black, G4ParticleDefinition *white)
Definition: G4HadronBuilder.cc:53
G4LundStringFragmentation::Quark_Diquark_lastSplitting
G4bool Quark_Diquark_lastSplitting(G4FragmentingString *&string, G4ParticleDefinition *&LeftHadron, G4ParticleDefinition *&RightHadron)
Definition: G4LundStringFragmentation.cc:1191
G4ExcitedString::GetLeftParton
G4Parton * GetLeftParton(void) const
Definition: G4ExcitedString.cc:161
G4FragmentingString::Mass
G4double Mass() const
Definition: G4FragmentingString.cc:271
G4VLongitudinalStringDecay::SampleQuarkPt
G4ThreeVector SampleQuarkPt(G4double ptMax=-1.)
Definition: G4VLongitudinalStringDecay.cc:496
G4VLongitudinalStringDecay::vectorMesonMix
std::vector< G4double > vectorMesonMix
Definition: G4VLongitudinalStringDecay.hh:189
G4ParticleDefinition::GetParticleSubType
const G4String & GetParticleSubType() const
Definition: G4ParticleDefinition.hh:181
G4VLongitudinalStringDecay::SetMassCut
virtual void SetMassCut(G4double aValue)
Definition: G4VLongitudinalStringDecay.cc:140
JustWarning
Definition: G4ExceptionSeverity.hh:64
G4int
int G4int
Definition: G4Types.hh:78
G4ParticleDefinition::GetParticleName
const G4String & GetParticleName() const
Definition: G4ParticleDefinition.hh:159
G4FragmentingString::GetDecayParton
G4ParticleDefinition * GetDecayParton() const
Definition: G4FragmentingString.hh:131
G4LundStringFragmentation::GetLightConeZ
virtual G4double GetLightConeZ(G4double zmin, G4double zmax, G4int PartonEncoding, G4ParticleDefinition *pHadron, G4double Px, G4double Py)
Definition: G4LundStringFragmentation.cc:918
G4VLongitudinalStringDecay::CreatePartonPair
pDefPair CreatePartonPair(G4int NeedParticle, G4bool AllowDiquarks=true)
Definition: G4VLongitudinalStringDecay.cc:472
G4VLongitudinalStringDecay::SetStringTensionParameter
void SetStringTensionParameter(G4double aValue)
Definition: G4VLongitudinalStringDecay.cc:656
G4ExcitedString::GetPosition
const G4ThreeVector & GetPosition() const
Definition: G4ExcitedString.hh:149
G4LundStringFragmentation::Sample4Momentum
virtual void Sample4Momentum(G4LorentzVector *Mom, G4double Mass, G4LorentzVector *AntiMom, G4double AntiMass, G4double InitialMass)
Definition: G4LundStringFragmentation.cc:751
G4FragmentingString::GetLeftParton
G4ParticleDefinition * GetLeftParton(void) const
Definition: G4FragmentingString.hh:137
G4VLongitudinalStringDecay::SetDiquarkBreakProbability
void SetDiquarkBreakProbability(G4double aValue)
Definition: G4VLongitudinalStringDecay.cc:575
G4LundStringFragmentation::SetMinimalStringMass2
void SetMinimalStringMass2(const G4double aValue)
Definition: G4LundStringFragmentation.cc:502
G4UniformRand
#define G4UniformRand()
Definition: Randomize.hh:93
G4cout
G4GLOB_DLL std::ostream G4cout
G4INCL::ClusterDecay::decay
ParticleList decay(Cluster *const c)
Carries out a cluster decay.
Definition: G4INCLClusterDecay.cc:553
G4VLongitudinalStringDecay::LightFragmentationTest
G4KineticTrackVector * LightFragmentationTest(const G4ExcitedString *const theString)
Definition: G4VLongitudinalStringDecay.cc:146
G4ExcitedString::Get4Momentum
G4LorentzVector Get4Momentum() const
Definition: G4ExcitedString.hh:161
G4LundStringFragmentation::Loop_toFragmentString
G4bool Loop_toFragmentString(G4ExcitedString *&theStringInCMS, G4KineticTrackVector *&LeftVector, G4KineticTrackVector *&RightVector)
Definition: G4LundStringFragmentation.cc:989
G4KineticTrack::GetFormationTime
G4double GetFormationTime() const
Definition: G4KineticTrack.hh:225
G4VLongitudinalStringDecay::pDefPair
std::pair< G4ParticleDefinition *, G4ParticleDefinition * > pDefPair
Definition: G4VLongitudinalStringDecay.hh:71
G4LundStringFragmentation::SetMinimalStringMass
void SetMinimalStringMass(const G4FragmentingString *const string)
Definition: G4LundStringFragmentation.cc:409
G4bool
bool G4bool
Definition: G4Types.hh:79
G4ExcitedString::TransformToAlignedCms
G4LorentzRotation TransformToAlignedCms()
Definition: G4ExcitedString.hh:237
G4VLongitudinalStringDecay::FindParticle
G4ParticleDefinition * FindParticle(G4int Encoding)
Definition: G4VLongitudinalStringDecay.cc:278
G4KineticTrack::SetPosition
void SetPosition(const G4ThreeVector aPosition)
Definition: G4KineticTrack.hh:242
G4LundStringFragmentation::StopFragmenting
virtual G4bool StopFragmenting(const G4FragmentingString *const string)
Definition: G4LundStringFragmentation.cc:636
GeV
static const double GeV
Definition: G4SIunits.hh:196
Randomize.hh
G4FragmentingString::GetDecayDirection
G4int GetDecayDirection() const
Definition: G4FragmentingString.cc:190
G4KineticTrack::Set4Momentum
void Set4Momentum(const G4LorentzVector &a4Momentum)
Definition: G4KineticTrack.hh:258
G4Quarks.hh
G4LundStringFragmentation.hh
G4LundStringFragmentation::lambda
G4double lambda(G4double s, G4double m1_Sqr, G4double m2_Sqr)
Definition: G4LundStringFragmentation.cc:975
G4FragmentingString::FourQuarkString
G4bool FourQuarkString(void) const
Definition: G4FragmentingString.cc:200
G4VLongitudinalStringDecay::scalarMesonMix
std::vector< G4double > scalarMesonMix
Definition: G4VLongitudinalStringDecay.hh:190
G4LundStringFragmentation::DiQuarkSplitup
virtual G4ParticleDefinition * DiQuarkSplitup(G4ParticleDefinition *decay, G4ParticleDefinition *&created)
Definition: G4LundStringFragmentation.cc:1417
G4Exception
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
G4PhysicalConstants.hh
G4DiQuarks.hh
C1
#define C1
G4ExcitedString::GetRightParton
G4Parton * GetRightParton(void) const
Definition: G4ExcitedString.cc:168
G4ParticleDefinition::GetPDGMass
G4double GetPDGMass() const
Definition: G4ParticleDefinition.hh:161
G4LundStringFragmentation::Diquark_AntiDiquark_aboveThreshold_lastSplitting
G4bool Diquark_AntiDiquark_aboveThreshold_lastSplitting(G4FragmentingString *&string, G4ParticleDefinition *&LeftHadron, G4ParticleDefinition *&RightHadron)
Definition: G4LundStringFragmentation.cc:1103
G4ParticleTable::GetParticleTable
static G4ParticleTable * GetParticleTable()
Definition: G4ParticleTable.cc:96
G4INCL::Math::max
T max(const T t1, const T t2)
brief Return the largest of the two arguments
Definition: G4INCLGlobals.hh:112
G4ExcitedString::GetTimeOfCreation
G4double GetTimeOfCreation() const
Definition: G4ExcitedString.hh:137
G4INCL::Math::min
T min(const T t1, const T t2)
brief Return the smallest of the two arguments
Definition: G4INCLGlobals.hh:117
G4ExcitedString::GetDirection
G4int GetDirection(void) const
Definition: G4ExcitedString.cc:154
G4endl
#define G4endl
Definition: G4ios.hh:61
G4LundStringFragmentation::FS_LeftHadron
G4ParticleDefinition * FS_LeftHadron[35]
Definition: G4LundStringFragmentation.hh:133
G4LundStringFragmentation::Diquark_AntiDiquark_belowThreshold_lastSplitting
G4bool Diquark_AntiDiquark_belowThreshold_lastSplitting(G4FragmentingString *&string, G4ParticleDefinition *&LeftHadron, G4ParticleDefinition *&RightHadron)
Definition: G4LundStringFragmentation.cc:1059
sqr
T sqr(const T &x)
Definition: templates.hh:145
G4LundStringFragmentation::SplitLast
virtual G4bool SplitLast(G4FragmentingString *string, G4KineticTrackVector *LeftVector, G4KineticTrackVector *RightVector)
Definition: G4LundStringFragmentation.cc:650
G4double
double G4double
Definition: G4Types.hh:76
G4SystemOfUnits.hh
G4KineticTrack::Get4Momentum
const G4LorentzVector & Get4Momentum() const
Definition: G4KineticTrack.hh:248
G4VLongitudinalStringDecay::CalculateHadronTimePosition
void CalculateHadronTimePosition(G4double theInitialStringMass, G4KineticTrackVector *)
Definition: G4VLongitudinalStringDecay.cc:513
G4KineticTrack::GetDefinition
const G4ParticleDefinition * GetDefinition() const
Definition: G4KineticTrack.hh:213
fermi
static const double fermi
Definition: G4SIunits.hh:93
G4VLongitudinalStringDecay::Splitup
G4KineticTrack * Splitup(G4FragmentingString *string, G4FragmentingString *&newString)
Definition: G4VLongitudinalStringDecay.cc:314
G4LundStringFragmentation::FS_RightHadron
G4ParticleDefinition * FS_RightHadron[35]
Definition: G4LundStringFragmentation.hh:133
G4VLongitudinalStringDecay::SetDiquarkSuppression
void SetDiquarkSuppression(G4double aValue)
Definition: G4VLongitudinalStringDecay.cc:564
G4LorentzVector
CLHEP::HepLorentzVector G4LorentzVector
Definition: G4LorentzVector.hh:35