Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4VLongitudinalStringDecay.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: G4VLongitudinalStringDecay.cc 100828 2016-11-02 15:25:59Z gcosmo $
28 //
29 // -----------------------------------------------------------------------------
30 // GEANT 4 class implementation file
31 //
32 // History: first implementation, Maxim Komogorov, 1-Jul-1998
33 // redesign Gunter Folger, August/September 2001
34 // -----------------------------------------------------------------------------
36 #include "G4PhysicalConstants.hh"
37 #include "G4SystemOfUnits.hh"
38 #include "G4ios.hh"
39 #include "Randomize.hh"
40 #include "G4FragmentingString.hh"
41 
42 #include "G4ParticleDefinition.hh"
43 #include "G4ParticleTypes.hh"
44 #include "G4ParticleChange.hh"
45 #include "G4VShortLivedParticle.hh"
47 #include "G4ParticleTable.hh"
49 #include "G4VDecayChannel.hh"
50 #include "G4DecayTable.hh"
51 
52 #include "G4DiQuarks.hh"
53 #include "G4Quarks.hh"
54 #include "G4Gluons.hh"
55 
56 #include "G4Exp.hh"
57 #include "G4Log.hh"
58 
59 //------------------------debug switches
60 //#define debug_VStringDecay
61 
62 //********************************************************************************
63 // Constructors
64 
66 {
67  MassCut = 0.35*GeV;
68  ClusterMass = 0.15*GeV;
69 
70  SmoothParam = 0.9;
71  StringLoopInterrupt = 1000;
73 
74  // Changable Parameters below.
75  SigmaQT = 0.5 * GeV; // 0.5 0.1
76 
77  StrangeSuppress = 0.44; // 27 % strange quarks produced, ie. u:d:s=1:1:0.27
78  DiquarkSuppress = 0.07;
79  DiquarkBreakProb = 0.1;
80 
81  //... pspin_meson is probability to create pseudo-scalar meson
82  pspin_meson = 0.5;
83 
84  //... pspin_barion is probability to create 1/2 barion
85  pspin_barion = 0.5;
86 
87  //... vectorMesonMix[] is quark mixing parameters for vector mesons (Variable spin = 3)
88  vectorMesonMix.resize(6);
89  vectorMesonMix[0] = 0.0; //AR-20Oct2014 : it was 0.5
90  vectorMesonMix[1] = 0.0;
91  vectorMesonMix[2] = 0.0; //AR-20Oct2014 : it was 0.5
92  vectorMesonMix[3] = 0.0;
93  vectorMesonMix[4] = 1.0;
94  vectorMesonMix[5] = 1.0;
95 
96  //... scalarMesonMix[] is quark mixing parameters for scalar mesons (Variable spin=1)
97  scalarMesonMix.resize(6);
98  scalarMesonMix[0] = 0.5;
99  scalarMesonMix[1] = 0.25;
100  scalarMesonMix[2] = 0.5;
101  scalarMesonMix[3] = 0.25;
102  scalarMesonMix[4] = 1.0;
103  scalarMesonMix[5] = 0.5;
104 
105  // Parameters may be changed until the first fragmentation starts
106  PastInitPhase=false;
108  Kappa = 1.0 * GeV/fermi;
109 }
110 
111 
113 {
114  delete hadronizer;
115 }
116 
117 //=============================================================================
118 
119 // Operators
120 
121 //-----------------------------------------------------------------------------
122 
123 int G4VLongitudinalStringDecay::operator==(const G4VLongitudinalStringDecay &) const
124 {
125  throw G4HadronicException(__FILE__, __LINE__, "G4VLongitudinalStringDecay::operator== forbidden");
126  return false;
127 }
128 
129 //-------------------------------------------------------------------------------------
130 
131 int G4VLongitudinalStringDecay::operator!=(const G4VLongitudinalStringDecay &) const
132 {
133  throw G4HadronicException(__FILE__, __LINE__, "G4VLongitudinalStringDecay::operator!= forbidden");
134  return true;
135 }
136 
137 //***********************************************************************************
138 
139 // For changing Mass Cut used for selection of very small mass strings
141 
142 //-----------------------------------------------------------------------------
143 
144 // For handling a string with very low mass
145 
148 {
149  // Check string decay threshold
150  G4KineticTrackVector * result=0; // return 0 when string exceeds the mass cut
151 
153 
154  G4FragmentingString aString(*string);
155  if ( sqr(FragmentationMass(&aString,0,&hadrons)+MassCut) < aString.Mass2()) {
156  return 0;
157  }
158 
159  // The string mass is very low ---------------------------
160 
161  result=new G4KineticTrackVector;
162 
163  if ( hadrons.second ==0 ) {
164  // Substitute string by light hadron, Note that Energy is not conserved here!
165 
166  #ifdef debug_VStringDecay
167  G4cout << "VlongSF Warning replacing string by single hadron (G4VLongitudinalStringDecay)" <<G4endl;
168  G4cout << hadrons.first->GetParticleName()<<G4endl
169  << "string .. " << string->Get4Momentum() << " "
170  << string->Get4Momentum().m() << G4endl;
171  #endif
172 
173  G4ThreeVector Mom3 = string->Get4Momentum().vect();
174  G4LorentzVector Mom(Mom3, std::sqrt(Mom3.mag2() + sqr(hadrons.first->GetPDGMass())));
175  result->push_back(new G4KineticTrack(hadrons.first, 0, string->GetPosition(), Mom));
176  } else {
177  //... string was qq--qqbar type: Build two stable hadrons,
178 
179  #ifdef debug_VStringDecay
180  G4cout << "VlongSF Warning replacing qq-qqbar string by TWO hadrons (G4VLongitudinalStringDecay)"
181  << hadrons.first->GetParticleName() << " / "
182  << hadrons.second->GetParticleName()
183  << "string .. " << string->Get4Momentum() << " "
184  << string->Get4Momentum().m() << G4endl;
185  #endif
186 
187  G4LorentzVector Mom1, Mom2;
188  Sample4Momentum(&Mom1, hadrons.first->GetPDGMass(), &Mom2,hadrons.second->GetPDGMass(),
189  string->Get4Momentum().mag());
190 
191  result->push_back(new G4KineticTrack(hadrons.first, 0, string->GetPosition(), Mom1));
192  result->push_back(new G4KineticTrack(hadrons.second, 0, string->GetPosition(), Mom2));
193 
194  G4ThreeVector Velocity = string->Get4Momentum().boostVector();
195  result->Boost(Velocity);
196  }
197 
198  return result;
199 }
200 
201 //----------------------------------------------------------------------------------------
202 
204 FragmentationMass( const G4FragmentingString * const string, Pcreate build, pDefPair * pdefs )
205 {
206  G4double mass;
207  static G4ThreadLocal G4bool NeedInit(true);
208  static G4ThreadLocal std::vector<double> *nomix_G4MT_TLS_ = 0 ;
209  if (!nomix_G4MT_TLS_) nomix_G4MT_TLS_ = new std::vector<double>;
210  std::vector<double> &nomix = *nomix_G4MT_TLS_;
211  static G4ThreadLocal G4HadronBuilder * minMassHadronizer;
212  if ( NeedInit )
213  {
214  NeedInit = false;
215  nomix.resize(6);
216  for ( G4int i=0; i<6 ; i++ ) nomix[i]=0;
217 
218  //minMassHadronizer=new G4HadronBuilder(pspin_meson,pspin_barion,nomix,nomix);
219  minMassHadronizer=hadronizer;
220  }
221 
222  if ( build==0 ) build=&G4HadronBuilder::BuildLowSpin;
223 
224  G4ParticleDefinition *Hadron1, *Hadron2=0;
225 
226  if (!string->FourQuarkString() )
227  {
228  // spin 0 meson or spin 1/2 barion will be built
229 
230  Hadron1 = (minMassHadronizer->*build)(string->GetLeftParton(), string->GetRightParton());
231 
232  #ifdef debug_VStringDecay
233  G4cout<<"Quarks at the string ends "<<string->GetLeftParton()->GetParticleName()
234  <<" "<<string->GetRightParton()->GetParticleName()<<G4endl;
235  G4cout<<"(G4VLongitudinalStringDecay) Hadron "<<Hadron1->GetParticleName()<<" "<<Hadron1->GetPDGMass()<<G4endl;
236  #endif
237 
238  mass= (Hadron1)->GetPDGMass();
239  } else {
240  //... string is qq--qqbar: Build two stable hadrons,
241  //... with extra uubar or ddbar quark pair
242  G4int iflc = (G4UniformRand() < 0.5)? 1 : 2;
243  if (string->GetLeftParton()->GetPDGEncoding() < 0) iflc = -iflc;
244 
245  //... theSpin = 4; spin 3/2 baryons will be built
246  Hadron1 = (minMassHadronizer->*build)(string->GetLeftParton(), FindParticle(iflc));
247  Hadron2 = (minMassHadronizer->*build)(string->GetRightParton(), FindParticle(-iflc));
248  mass = (Hadron1)->GetPDGMass() + (Hadron2)->GetPDGMass();
249  }
250 
251  if ( pdefs != 0 )
252  { // need to return hadrons as well....
253  pdefs->first = Hadron1;
254  pdefs->second = Hadron2;
255  }
256 
257  return mass;
258 }
259 
260 //----------------------------------------------------------------------------
261 
263 {
265  if (ptr == NULL)
266  {
267  G4cout << "Particle with encoding "<<Encoding<<" does not exist!!!"<<G4endl;
268  throw G4HadronicException(__FILE__, __LINE__, "Check your particle table");
269  }
270  return ptr;
271 }
272 
273 //*********************************************************************************
274 // For decision on continue or stop string fragmentation
275 // virtual G4bool StopFragmenting(const G4FragmentingString * const string)=0;
276 // virtual G4bool IsFragmentable(const G4FragmentingString * const string)=0;
277 
278 // If a string can not fragment, make last break into 2 hadrons
279 // virtual G4bool SplitLast(G4FragmentingString * string,
280 // G4KineticTrackVector * LeftVector,
281 // G4KineticTrackVector * RightVector)=0;
282 //-----------------------------------------------------------------------------
283 //
284 // If a string fragments, do the following
285 //
286 // For transver of a string to its CMS frame
287 //-----------------------------------------------------------------------------
288 
290 {
291  G4Parton *Left=new G4Parton(*in.GetLeftParton());
292  G4Parton *Right=new G4Parton(*in.GetRightParton());
293  return new G4ExcitedString(Left,Right,in.GetDirection());
294 }
295 
296 //-----------------------------------------------------------------------------
297 
300 {
301 
302  #ifdef debug_VStringDecay
303  G4cout<<G4endl;
304  G4cout<<"Start SplitUP (G4VLongitudinalStringDecay) ========================="<<G4endl;
305  G4cout<<"String partons: " <<string->GetLeftParton()->GetPDGEncoding()<<" "
306  <<string->GetRightParton()->GetPDGEncoding()<<" "
307  <<"Direction " <<string->GetDecayDirection()<<G4endl;
308  #endif
309 
310  //... random choice of string end to use for creating the hadron (decay)
311  G4int SideOfDecay = (G4UniformRand() < 0.5)? 1: -1;
312  if (SideOfDecay < 0)
313  {
314  string->SetLeftPartonStable();
315  } else {
316  string->SetRightPartonStable();
317  }
318 
319  G4ParticleDefinition *newStringEnd;
320  G4ParticleDefinition * HadronDefinition;
321  if (string->DecayIsQuark())
322  {
323  HadronDefinition= QuarkSplitup(string->GetDecayParton(), newStringEnd);
324  } else {
325  HadronDefinition= DiQuarkSplitup(string->GetDecayParton(), newStringEnd);
326  }
327 
328  #ifdef debug_VStringDecay
329  G4cout<<"The parton "<<string->GetDecayParton()->GetPDGEncoding()<<" "
330  <<" produces hadron "<<HadronDefinition->GetParticleName()
331  <<" and is transformed to "<<newStringEnd->GetPDGEncoding()<<G4endl;
332  G4cout<<"The side of the string decay Left/Right (1/-1) "<<SideOfDecay<<G4endl;
333  #endif
334  // create new String from old, ie. keep Left and Right order, but replace decay
335 
336  newString=new G4FragmentingString(*string,newStringEnd); // To store possible quark containt of new string
337 
338  #ifdef debug_VStringDecay
339  G4cout<<"An attempt to determine its energy (SplitEandP)"<<G4endl;
340  #endif
341 
342  G4LorentzVector* HadronMomentum=SplitEandP(HadronDefinition, string, newString);
343 
344  delete newString; newString=0;
345 
346  G4KineticTrack * Hadron =0;
347  if ( HadronMomentum != 0 ) {
348 
349  #ifdef debug_VStringDecay
350  G4cout<<"The attempt was successful"<<G4endl;
351  #endif
352 
354  Hadron = new G4KineticTrack(HadronDefinition, 0,Pos, *HadronMomentum);
355 
356  newString=new G4FragmentingString(*string,newStringEnd,HadronMomentum);
357 
358  delete HadronMomentum;
359  } else {
360 
361  #ifdef debug_VStringDecay
362  G4cout<<"The attempt was not successful !!!"<<G4endl;
363  #endif
364  }
365 
366  #ifdef debug_VStringDecay
367  G4cout<<"End SplitUP (G4VLongitudinalStringDecay) ====================="<<G4endl;
368  #endif
369 
370  return Hadron;
371 }
372 
373 //--------------------------------------------------------------------------------------
374 
377 {
378  // if we have a quark, we need antiquark (or diquark)
379  G4int IsParticle=(decay->GetPDGEncoding()>0) ? -1 : +1;
380 
381  pDefPair QuarkPair = CreatePartonPair(IsParticle);
382  created = QuarkPair.second;
383  return hadronizer->Build(QuarkPair.first, decay);
384 }
385 
386 /* Uzhi June 2014
387 //-----------------------------------------------------------------------------
388 
389 G4ParticleDefinition *G4VLongitudinalStringDecay::
390 DiQuarkSplitup(G4ParticleDefinition* decay, G4ParticleDefinition *&created)
391 {
392  //... can Diquark break or not?
393  if (G4UniformRand() < DiquarkBreakProb ){
394 
395  //... Diquark break
396 
397  G4int stableQuarkEncoding = decay->GetPDGEncoding()/1000;
398  G4int decayQuarkEncoding = (decay->GetPDGEncoding()/100)%10;
399  if (G4UniformRand() < 0.5)
400  {
401  G4int Swap = stableQuarkEncoding;
402  stableQuarkEncoding = decayQuarkEncoding;
403  decayQuarkEncoding = Swap;
404  }
405 
406  G4int IsParticle=(decayQuarkEncoding>0) ? -1 : +1; // if we have a quark, we need antiquark
407 
408  pDefPair QuarkPair = CreatePartonPair(IsParticle,false); // no diquarks wanted
409 
410  //... Build new Diquark
411  G4int QuarkEncoding=QuarkPair.second->GetPDGEncoding();
412  G4int i10 = std::max(std::abs(QuarkEncoding), std::abs(stableQuarkEncoding));
413  G4int i20 = std::min(std::abs(QuarkEncoding), std::abs(stableQuarkEncoding));
414  G4int spin = (i10 != i20 && G4UniformRand() <= 0.5)? 1 : 3;
415  G4int NewDecayEncoding = -1*IsParticle*(i10 * 1000 + i20 * 100 + spin);
416  created = FindParticle(NewDecayEncoding);
417  G4ParticleDefinition * decayQuark=FindParticle(decayQuarkEncoding);
418  G4ParticleDefinition * had=hadronizer->Build(QuarkPair.first, decayQuark);
419  return had;
420  //return hadronizer->Build(QuarkPair.first, decayQuark);
421 
422  } else {
423  //... Diquark does not break
424 
425  G4int IsParticle=(decay->GetPDGEncoding()>0) ? +1 : -1; // if we have a diquark, we need quark
426  pDefPair QuarkPair = CreatePartonPair(IsParticle,false); // no diquarks wanted
427  created = QuarkPair.second;
428 
429  G4ParticleDefinition * had=hadronizer->Build(QuarkPair.first, decay);
430  return had;
431  //return G4ParticleDefinition * had=hadronizer->Build(QuarkPair.first, decay);
432  }
433 }
434 */ // Uzhi June 2014
435 
436 //-----------------------------------------------------------------------------
437 
439 {
440  return (1 + (int)(G4UniformRand()/StrangeSuppress));
441 }
442 
443 //-----------------------------------------------------------------------------
444 
446 CreatePartonPair(G4int NeedParticle,G4bool AllowDiquarks)
447 {
448  // NeedParticle = +1 for Particle, -1 for Antiparticle
449 
450  if ( AllowDiquarks && G4UniformRand() < DiquarkSuppress )
451  {
452  // Create a Diquark - AntiDiquark pair , first in pair is anti to IsParticle
453  G4int q1 = SampleQuarkFlavor();
454  G4int q2 = SampleQuarkFlavor();
455  G4int spin = (q1 != q2 && G4UniformRand() <= 0.5)? 1 : 3;
456  // Convention: quark with higher PDG number is first
457  G4int PDGcode = (std::max(q1,q2) * 1000 + std::min(q1,q2) * 100 + spin) * NeedParticle;
458  return pDefPair (FindParticle(-PDGcode),FindParticle(PDGcode));
459  } else {
460  // Create a Quark - AntiQuark pair, first in pair IsParticle
461  G4int PDGcode=SampleQuarkFlavor()*NeedParticle;
462  return pDefPair (FindParticle(PDGcode),FindParticle(-PDGcode));
463  }
464 
465 }
466 
467 //-----------------------------------------------------------------------------
468 
470 {
471  G4double Pt;
472  if ( ptMax < 0 ) {
473  // sample full gaussian
474  Pt = -G4Log(G4UniformRand());
475  } else {
476  // sample in limited range
477  Pt = -G4Log(G4RandFlat::shoot(G4Exp(-sqr(ptMax)/sqr(SigmaQT)), 1.));
478  }
479  Pt = SigmaQT * std::sqrt(Pt);
480  G4double phi = 2.*pi*G4UniformRand();
481  return G4ThreeVector(Pt * std::cos(phi),Pt * std::sin(phi),0);
482 }
483 
484 //******************************************************************************
485 
488 {
489  // 'yo-yo' formation time
490  //const G4double kappa = 1.0 * GeV/fermi/4.;
492  for(size_t c1 = 0; c1 < Hadrons->size(); c1++)
493  {
494  G4double SumPz = 0;
495  G4double SumE = 0;
496  for(size_t c2 = 0; c2 < c1; c2++)
497  {
498  SumPz += Hadrons->operator[](c2)->Get4Momentum().pz();
499  SumE += Hadrons->operator[](c2)->Get4Momentum().e();
500  }
501  G4double HadronE = Hadrons->operator[](c1)->Get4Momentum().e();
502  G4double HadronPz = Hadrons->operator[](c1)->Get4Momentum().pz();
503  Hadrons->operator[](c1)->SetFormationTime(
504  (theInitialStringMass - 2.*SumPz + HadronE - HadronPz)/(2.*kappa)/c_light );
505 
506  G4ThreeVector aPosition(0, 0, (theInitialStringMass - 2.*SumE - HadronE + HadronPz)/(2.*kappa));
507  Hadrons->operator[](c1)->SetPosition(aPosition);
508 
509  }
510 }
511 
512 //-----------------------------------------------------------------------------
513 
515 {
516  if ( PastInitPhase ) {
517  throw G4HadronicException(__FILE__, __LINE__,
518  "4VLongitudinalStringDecay::SetSigmaTransverseMomentum after FragmentString() not allowed");
519  } else {
520  SigmaQT = aValue;
521  }
522 }
523 
524 //----------------------------------------------------------------------------------------------------------
525 
527 {
528  if ( PastInitPhase ) {
529  throw G4HadronicException(__FILE__, __LINE__,
530  "4VLongitudinalStringDecay::SetStrangenessSuppression after FragmentString() not allowed");
531  } else {
532  StrangeSuppress = aValue;
533  }
534 }
535 
536 //----------------------------------------------------------------------------------------------------------
537 
539 {
540  if ( PastInitPhase ) {
541  throw G4HadronicException(__FILE__, __LINE__,
542  "4VLongitudinalStringDecay::SetDiquarkSuppression after FragmentString() not allowed");
543  } else {
544  DiquarkSuppress = aValue;
545  }
546 }
547 
548 //----------------------------------------------------------------------------------------
549 
551 {
552  if ( PastInitPhase ) {
553  throw G4HadronicException(__FILE__, __LINE__,
554  "4VLongitudinalStringDecay::SetDiquarkBreakProbability after FragmentString() not allowed");
555  } else {
556  DiquarkBreakProb = aValue;
557  }
558 }
559 
560 //----------------------------------------------------------------------------------------------------------
561 
563 {
564  if ( PastInitPhase ) {
565  throw G4HadronicException(__FILE__, __LINE__,
566  "G4VLongitudinalStringDecay::SetVectorMesonProbability after FragmentString() not allowed");
567  } else {
568  pspin_meson = aValue;
569  delete hadronizer;
571  }
572 }
573 
574 //----------------------------------------------------------------------------------------------------------
575 
577 {
578  if ( PastInitPhase ) {
579  throw G4HadronicException(__FILE__, __LINE__,
580  "G4VLongitudinalStringDecay::SetSpinThreeHalfBarionProbability after FragmentString() not allowed");
581  } else {
582  pspin_barion = aValue;
583  delete hadronizer;
585  }
586 }
587 
588 //----------------------------------------------------------------------------------------------------------
589 
590 void G4VLongitudinalStringDecay::SetScalarMesonMixings(std::vector<G4double> aVector)
591 {
592  if ( PastInitPhase ) {
593  throw G4HadronicException(__FILE__, __LINE__,
594  "G4VLongitudinalStringDecay::SetScalarMesonMixings after FragmentString() not allowed");
595  } else {
596  if ( aVector.size() < 6 )
597  throw G4HadronicException(__FILE__, __LINE__,
598  "G4VLongitudinalStringDecay::SetScalarMesonMixings( argument Vector too small");
599  scalarMesonMix[0] = aVector[0];
600  scalarMesonMix[1] = aVector[1];
601  scalarMesonMix[2] = aVector[2];
602  scalarMesonMix[3] = aVector[3];
603  scalarMesonMix[4] = aVector[4];
604  scalarMesonMix[5] = aVector[5];
605  delete hadronizer;
607  }
608 }
609 
610 //----------------------------------------------------------------------------------------------------------
611 
612 void G4VLongitudinalStringDecay::SetVectorMesonMixings(std::vector<G4double> aVector)
613 {
614  if ( PastInitPhase ) {
615  throw G4HadronicException(__FILE__, __LINE__,
616  "G4VLongitudinalStringDecay::SetVectorMesonMixings after FragmentString() not allowed");
617  } else {
618  if ( aVector.size() < 6 )
619  throw G4HadronicException(__FILE__, __LINE__,
620  "G4VLongitudinalStringDecay::SetVectorMesonMixings( argument Vector too small");
621  vectorMesonMix[0] = aVector[0];
622  vectorMesonMix[1] = aVector[1];
623  vectorMesonMix[2] = aVector[2];
624  vectorMesonMix[3] = aVector[3];
625  vectorMesonMix[4] = aVector[4];
626  vectorMesonMix[5] = aVector[5];
627  delete hadronizer;
629  }
630 }
631 
632 //-------------------------------------------------------------------------------------------
634 {
635  Kappa = aValue * GeV/fermi;
636 }
637 //**************************************************************************************
638 
G4double G4ParticleHPJENDLHEData::G4double result
G4ParticleDefinition * GetRightParton(void) const
ThreeVector shoot(const G4int Ap, const G4int Af)
static c2_factory< G4double > c2
G4ParticleDefinition * FindParticle(G4int PDGEncoding)
G4ParticleDefinition * QuarkSplitup(G4ParticleDefinition *decay, G4ParticleDefinition *&created)
CLHEP::Hep3Vector G4ThreeVector
G4ParticleDefinition * BuildLowSpin(G4ParticleDefinition *black, G4ParticleDefinition *white)
G4ExcitedString * CPExcited(const G4ExcitedString &string)
ush Pos
Definition: deflate.h:89
void SetStrangenessSuppression(G4double aValue)
G4ParticleDefinition * Build(G4ParticleDefinition *black, G4ParticleDefinition *white)
G4Parton * GetLeftParton(void) const
G4ThreeVector SampleQuarkPt(G4double ptMax=-1.)
std::vector< G4double > vectorMesonMix
#define G4ThreadLocal
Definition: tls.hh:89
G4double FragmentationMass(const G4FragmentingString *const string, Pcreate build=0, pDefPair *pdefs=0)
virtual void SetMassCut(G4double aValue)
int G4int
Definition: G4Types.hh:78
const G4String & GetParticleName() const
G4ParticleDefinition * GetDecayParton() const
pDefPair CreatePartonPair(G4int NeedParticle, G4bool AllowDiquarks=true)
void SetStringTensionParameter(G4double aValue)
const G4ThreeVector & GetPosition() const
G4ParticleDefinition * GetLeftParton(void) const
G4double Mass2() const
void SetDiquarkBreakProbability(G4double aValue)
#define G4UniformRand()
Definition: Randomize.hh:97
G4GLOB_DLL std::ostream G4cout
ParticleList decay(Cluster *const c)
Carries out a cluster decay.
G4KineticTrackVector * LightFragmentationTest(const G4ExcitedString *const theString)
std::pair< G4ParticleDefinition *, G4ParticleDefinition * > pDefPair
double mag() const
bool G4bool
Definition: G4Types.hh:79
G4ParticleDefinition * FindParticle(G4int Encoding)
virtual void Sample4Momentum(G4LorentzVector *Mom, G4double Mass, G4LorentzVector *AntiMom, G4double AntiMass, G4double InitialMass)=0
void Boost(G4ThreeVector &Velocity)
G4bool FourQuarkString(void) const
std::vector< G4double > scalarMesonMix
G4double G4Log(G4double x)
Definition: G4Log.hh:230
G4double G4Exp(G4double initial_x)
Exponential Function double precision.
Definition: G4Exp.hh:183
void SetVectorMesonProbability(G4double aValue)
G4Parton * GetRightParton(void) const
G4double GetPDGMass() const
static G4ParticleTable * GetParticleTable()
T max(const T t1, const T t2)
brief Return the largest of the two arguments
void SetVectorMesonMixings(std::vector< G4double > aVector)
void SetSpinThreeHalfBarionProbability(G4double aValue)
T min(const T t1, const T t2)
brief Return the smallest of the two arguments
static constexpr double GeV
Definition: G4SIunits.hh:217
double mag2() const
G4int GetDirection(void) const
#define G4endl
Definition: G4ios.hh:61
static constexpr double pi
Definition: G4SIunits.hh:75
T sqr(const T &x)
Definition: templates.hh:145
virtual G4ParticleDefinition * DiQuarkSplitup(G4ParticleDefinition *decay, G4ParticleDefinition *&created)=0
double G4double
Definition: G4Types.hh:76
void SetScalarMesonMixings(std::vector< G4double > aVector)
static constexpr double fermi
Definition: G4SIunits.hh:103
void CalculateHadronTimePosition(G4double theInitialStringMass, G4KineticTrackVector *)
virtual G4LorentzVector * SplitEandP(G4ParticleDefinition *pHadron, G4FragmentingString *string, G4FragmentingString *newString)=0
float c_light
Definition: hepunit.py:257
tuple c1
Definition: plottest35.py:14
G4KineticTrack * Splitup(G4FragmentingString *string, G4FragmentingString *&newString)
void SetDiquarkSuppression(G4double aValue)