Geant4  v4-10.4-release
 모두 클래스 네임스페이스들 파일들 함수 변수 타입정의 열거형 타입 열거형 멤버 Friends 매크로 그룹들 페이지들
G4VLongitudinalStringDecay.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 // $Id: G4VLongitudinalStringDecay.cc 108179 2018-01-18 08:53:39Z 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;
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;// 0.62 0.3 Uzhi Dec. 2017
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.; // Uzhi May 2016 //AR-20Oct2014 : it was 0.5 0.5 Uzhi Dec. 2017
90  vectorMesonMix[1] = 0.375;
91  vectorMesonMix[2] = 0.0; // Uzhi May 2016 //AR-20Oct2014 : it was 0.5 0.5 Uzhi Dec. 2017
92  vectorMesonMix[3] = 0.375;
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;
109  Kappa = 1.0 * GeV/fermi;
110 
111 
112 }
113 
114 
116  {
117  delete hadronizer;
118  }
119 
120 //=============================================================================
121 
122 // Operators
123 
124 //-----------------------------------------------------------------------------
125 
127  {
128  throw G4HadronicException(__FILE__, __LINE__, "G4VLongitudinalStringDecay::operator== forbidden");
129  return false;
130  }
131 
132 //-------------------------------------------------------------------------------------
133 
135  {
136  throw G4HadronicException(__FILE__, __LINE__, "G4VLongitudinalStringDecay::operator!= forbidden");
137  return true;
138  }
139 
140 //***********************************************************************************
141 
142 // For changing Mass Cut used for selection of very small mass strings
145 //-----------------------------------------------------------------------------
146 
147 // For handling a string with very low mass
148 
150  G4ExcitedString * const string)
151 {
152  // Check string decay threshold
153  G4KineticTrackVector * result=0; // return 0 when string exceeds the mass cut
154 
156 
157  G4FragmentingString aString(*string);
158  if ( sqr(FragmentationMass(&aString,0,&hadrons)+MassCut) < aString.Mass2()) {
159  return 0;
160  }
161 
162 // The string mass is very low ---------------------------
163 
164  result=new G4KineticTrackVector;
165 
166  if ( hadrons.second ==0 )
167  {
168 // Substitute string by light hadron, Note that Energy is not conserved here!
169 
170 #ifdef debug_VStringDecay
171  G4cout << "VlongSF Warning replacing string by single hadron (G4VLongitudinalStringDecay)" <<G4endl;
172  G4cout << hadrons.first->GetParticleName()<<G4endl
173  << "string .. " << string->Get4Momentum() << " "
174  << string->Get4Momentum().m() << G4endl;
175 #endif
176 
177  G4ThreeVector Mom3 = string->Get4Momentum().vect();
178  G4LorentzVector Mom(Mom3,
179  std::sqrt(Mom3.mag2() +
180  sqr(hadrons.first->GetPDGMass())));
181  result->push_back(new G4KineticTrack(hadrons.first, 0,
182  string->GetPosition(),
183  Mom));
184  } else
185  {
186 //... string was qq--qqbar type: Build two stable hadrons,
187 
188 #ifdef debug_VStringDecay
189  G4cout << "VlongSF Warning replacing qq-qqbar string by TWO hadrons (G4VLongitudinalStringDecay)"
190  << hadrons.first->GetParticleName() << " / "
191  << hadrons.second->GetParticleName()
192  << "string .. " << string->Get4Momentum() << " "
193  << string->Get4Momentum().m() << G4endl;
194 #endif
195 
196  G4LorentzVector Mom1, Mom2;
197  Sample4Momentum(&Mom1, hadrons.first->GetPDGMass(),
198  &Mom2,hadrons.second->GetPDGMass(),
199  string->Get4Momentum().mag());
200 
201  result->push_back(new G4KineticTrack(hadrons.first, 0,
202  string->GetPosition(),
203  Mom1));
204  result->push_back(new G4KineticTrack(hadrons.second, 0,
205  string->GetPosition(),
206  Mom2));
207 
208  G4ThreeVector Velocity = string->Get4Momentum().boostVector();
209  result->Boost(Velocity);
210  }
211 
212  return result;
213 
214 }
215 
216 //----------------------------------------------------------------------------------------
217 
219  const G4FragmentingString * const string,
220  Pcreate build, pDefPair * pdefs )
221 {
222  G4double mass;
223 
224  if ( build==0 ) build=&G4HadronBuilder::BuildLowSpin;
225 
226  G4ParticleDefinition *Hadron1, *Hadron2=0;
227 
228  if (!string->FourQuarkString() )
229  {
230  // spin 0 meson or spin 1/2 barion will be built
231 
232  Hadron1 = (hadronizer->*build)(string->GetLeftParton(), string->GetRightParton());
233 #ifdef debug_VStringDecay
234  G4cout<<"Quarks at the string ends "<<string->GetLeftParton()->GetParticleName()<<" "<<string->GetRightParton()->GetParticleName()<<G4endl;
235  G4cout<<"(G4VLongitudinalStringDecay) Hadron "<<Hadron1->GetParticleName()<<" "<<Hadron1->GetPDGMass()<<G4endl;
236 #endif
237  mass= (Hadron1)->GetPDGMass();
238  } else
239  {
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 = (hadronizer->*build)(string->GetLeftParton(), FindParticle(iflc));
247  Hadron2 = (hadronizer->*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 can fragment, 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 //-----------------------------------------------------------------------------
299  decay, G4ParticleDefinition *&created)
300 {
301  G4int IsParticle=(decay->GetPDGEncoding()>0) ? -1 : +1; // if we have a quark,
302  // we need antiquark
303  // (or diquark)
304  pDefPair QuarkPair = CreatePartonPair(IsParticle);
305  created = QuarkPair.second;
306  return hadronizer->Build(QuarkPair.first, decay);
307 
308 }
309 
310 //-----------------------------------------------------------------------------
311 
313  {
314  return (1 + (int)(G4UniformRand()/StrangeSuppress));
315  }
316 
317 //-----------------------------------------------------------------------------
318 
320 {
321 // NeedParticle = +1 for Particle, -1 for Antiparticle
322  if ( AllowDiquarks && G4UniformRand() < DiquarkSuppress )
323  {
324  // Create a Diquark - AntiDiquark pair , first in pair is anti to IsParticle
325  G4int q1 = SampleQuarkFlavor();
326  G4int q2 = SampleQuarkFlavor();
327  G4int spin = (q1 != q2 && G4UniformRand() <= 0.5)? 1 : 3;
328  // convention: quark with higher PDG number is first
329  G4int PDGcode = (std::max(q1,q2) * 1000 + std::min(q1,q2) * 100 + spin) * NeedParticle;
330  return pDefPair (FindParticle(-PDGcode),FindParticle(PDGcode));
331 
332 
333  } else {
334  // Create a Quark - AntiQuark pair, first in pair IsParticle
335  G4int PDGcode=SampleQuarkFlavor()*NeedParticle;
336  return pDefPair (FindParticle(PDGcode),FindParticle(-PDGcode));
337  }
338 
339 }
340 
341 //-----------------------------------------------------------------------------
343  {
344  G4double Pt;
345  if ( ptMax < 0 ) {
346  // sample full gaussian
347  Pt = -G4Log(G4UniformRand());
348  } else {
349  // sample in limited range
350  Pt = -G4Log(G4RandFlat::shoot(G4Exp(-sqr(ptMax)/sqr(SigmaQT)), 1.));
351  }
352  Pt = SigmaQT * std::sqrt(Pt);
353  G4double phi = 2.*pi*G4UniformRand();
354  return G4ThreeVector(Pt * std::cos(phi),Pt * std::sin(phi),0);
355  }
356 
357 //******************************************************************************
358 
360  {
361 
362 // `yo-yo` formation time
363 // const G4double kappa = 1.0 * GeV/fermi/4.;
365  for(size_t c1 = 0; c1 < Hadrons->size(); c1++)
366  {
367  G4double SumPz = 0;
368  G4double SumE = 0;
369  for(size_t c2 = 0; c2 < c1; c2++)
370  {
371  SumPz += Hadrons->operator[](c2)->Get4Momentum().pz();
372  SumE += Hadrons->operator[](c2)->Get4Momentum().e();
373  }
374  G4double HadronE = Hadrons->operator[](c1)->Get4Momentum().e();
375  G4double HadronPz = Hadrons->operator[](c1)->Get4Momentum().pz();
376  Hadrons->operator[](c1)->SetFormationTime(
377 (theInitialStringMass - 2.*SumPz + HadronE - HadronPz)/(2.*kappa)/c_light);
378 
379  G4ThreeVector aPosition(0, 0,
380 (theInitialStringMass - 2.*SumE - HadronE + HadronPz)/(2.*kappa));
381  Hadrons->operator[](c1)->SetPosition(aPosition);
382 
383  }
384  }
385 
386 //-----------------------------------------------------------------------------
387 
389 {
390  if ( PastInitPhase ) {
391  throw G4HadronicException(__FILE__, __LINE__, "4VLongitudinalStringDecay::SetSigmaTransverseMomentum after FragmentString() not allowed");
392  } else {
393  SigmaQT = aValue;
394  }
395 }
396 
397 //----------------------------------------------------------------------------------------------------------
398 
400 {
401 // if ( PastInitPhase ) { // Uzhi Oct. 2017
402 // throw G4HadronicException(__FILE__, __LINE__, "4VLongitudinalStringDecay::SetStrangenessSuppression after FragmentString() not allowed");
403 // } else {
404  StrangeSuppress = aValue;
405 // }
406 }
407 
408 //----------------------------------------------------------------------------------------------------------
409 
411 {
412 // if ( PastInitPhase ) { // Uzhi Oct. 2017
413 // throw G4HadronicException(__FILE__, __LINE__, "4VLongitudinalStringDecay::SetDiquarkSuppression after FragmentString() not allowed");
414 // } else {
415  DiquarkSuppress = aValue;
416 // }
417 }
418 
419 //----------------------------------------------------------------------------------------
420 
422 {
423  if ( PastInitPhase ) {
424  throw G4HadronicException(__FILE__, __LINE__, "4VLongitudinalStringDecay::SetDiquarkBreakProbability after FragmentString() not allowed");
425  } else {
426  DiquarkBreakProb = aValue;
427  }
428 }
429 
430 //----------------------------------------------------------------------------------------------------------
431 
433 {
434  if ( PastInitPhase ) {
435  throw G4HadronicException(__FILE__, __LINE__, "G4VLongitudinalStringDecay::SetVectorMesonProbability after FragmentString() not allowed");
436  } else {
437  pspin_meson = aValue;
438  delete hadronizer;
441  }
442 }
443 
444 //----------------------------------------------------------------------------------------------------------
445 
447 {
448  if ( PastInitPhase ) {
449  throw G4HadronicException(__FILE__, __LINE__, "G4VLongitudinalStringDecay::SetSpinThreeHalfBarionProbability after FragmentString() not allowed");
450  } else {
451  pspin_barion = aValue;
452  delete hadronizer;
455  }
456 }
457 
458 //----------------------------------------------------------------------------------------------------------
459 
460 void G4VLongitudinalStringDecay::SetScalarMesonMixings(std::vector<G4double> aVector)
461 {
462  if ( PastInitPhase ) {
463  throw G4HadronicException(__FILE__, __LINE__, "G4VLongitudinalStringDecay::SetScalarMesonMixings after FragmentString() not allowed");
464  } else {
465  if ( aVector.size() < 6 )
466  throw G4HadronicException(__FILE__, __LINE__, "G4VLongitudinalStringDecay::SetScalarMesonMixings( argument Vector too small");
467  scalarMesonMix[0] = aVector[0];
468  scalarMesonMix[1] = aVector[1];
469  scalarMesonMix[2] = aVector[2];
470  scalarMesonMix[3] = aVector[3];
471  scalarMesonMix[4] = aVector[4];
472  scalarMesonMix[5] = aVector[5];
473  delete hadronizer;
476  }
477 }
478 
479 //----------------------------------------------------------------------------------------------------------
480 
481 void G4VLongitudinalStringDecay::SetVectorMesonMixings(std::vector<G4double> aVector)
482 {
483  if ( PastInitPhase ) {
484  throw G4HadronicException(__FILE__, __LINE__, "G4VLongitudinalStringDecay::SetVectorMesonMixings after FragmentString() not allowed");
485  } else {
486  if ( aVector.size() < 6 )
487  throw G4HadronicException(__FILE__, __LINE__, "G4VLongitudinalStringDecay::SetVectorMesonMixings( argument Vector too small");
488  vectorMesonMix[0] = aVector[0];
489  vectorMesonMix[1] = aVector[1];
490  vectorMesonMix[2] = aVector[2];
491  vectorMesonMix[3] = aVector[3];
492  vectorMesonMix[4] = aVector[4];
493  vectorMesonMix[5] = aVector[5];
494  delete hadronizer;
497 
498  }
499 }
500 
501 //-------------------------------------------------------------------------------------------
503 {
504  Kappa = aValue * GeV/fermi;
505 }
506 //**************************************************************************************
507 
G4ParticleDefinition * FindParticle(G4int PDGEncoding)
G4double G4Exp(G4double initial_x)
Exponential Function double precision.
Definition: G4Exp.hh:183
ParticleList decay(Cluster *const c)
Carries out a cluster decay.
T max(const T t1, const T t2)
brief Return the largest of the two arguments
G4Parton * GetRightParton(void) const
G4ParticleDefinition * Build(G4ParticleDefinition *black, G4ParticleDefinition *white)
CLHEP::Hep3Vector G4ThreeVector
const G4ThreeVector & GetPosition() const
static G4ParticleTable * GetParticleTable()
void SetDiquarkSuppression(G4double aValue)
G4ParticleDefinition * BuildLowSpin(G4ParticleDefinition *black, G4ParticleDefinition *white)
G4ParticleDefinition * GetRightParton(void) const
virtual void SetMassCut(G4double aValue)
#define G4endl
Definition: G4ios.hh:61
void Boost(G4ThreeVector &Velocity)
void SetStringTensionParameter(G4double aValue)
const G4String & GetParticleName() const
void SetVectorMesonMixings(std::vector< G4double > aVector)
G4double FragmentationMass(const G4FragmentingString *const string, Pcreate build=0, pDefPair *pdefs=0)
void SetScalarMesonMixings(std::vector< G4double > aVector)
int operator==(const G4VLongitudinalStringDecay &right) const
G4double GetPDGMass() const
void SetSpinThreeHalfBarionProbability(G4double aValue)
void SetDiquarkBreakProbability(G4double aValue)
G4double G4Log(G4double x)
Definition: G4Log.hh:230
std::vector< G4double > vectorMesonMix
G4int GetDirection(void) const
double G4double
Definition: G4Types.hh:76
bool G4bool
Definition: G4Types.hh:79
static constexpr double fermi
Definition: G4SIunits.hh:103
void SetVectorMesonProbability(G4double aValue)
G4ParticleDefinition * FindParticle(G4int Encoding)
G4ThreeVector SampleQuarkPt(G4double ptMax=-1.)
TCanvas * c2
Definition: plot_hist.C:75
void CalculateHadronTimePosition(G4double theInitialStringMass, G4KineticTrackVector *)
G4KineticTrackVector * LightFragmentationTest(const G4ExcitedString *const theString)
#define G4UniformRand()
Definition: Randomize.hh:53
void SetStrangenessSuppression(G4double aValue)
G4Parton * GetLeftParton(void) const
double mag2() const
ThreeVector shoot(const G4int Ap, const G4int Af)
int operator!=(const G4VLongitudinalStringDecay &right) const
G4ParticleDefinition * QuarkSplitup(G4ParticleDefinition *decay, G4ParticleDefinition *&created)
G4double G4ParticleHPJENDLHEData::G4double result
virtual void Sample4Momentum(G4LorentzVector *Mom, G4double Mass, G4LorentzVector *AntiMom, G4double AntiMass, G4double InitialMass)=0
double mag() const
int G4int
Definition: G4Types.hh:78
static constexpr double c_light
G4ParticleDefinition * GetLeftParton(void) const
ifstream in
Definition: comparison.C:7
G4bool FourQuarkString(void) const
pDefPair CreatePartonPair(G4int NeedParticle, G4bool AllowDiquarks=true)
G4GLOB_DLL std::ostream G4cout
G4ExcitedString * CopyExcited(const G4ExcitedString &string)
T sqr(const T &x)
Definition: templates.hh:145
static constexpr double pi
Definition: G4SIunits.hh:75
std::pair< G4ParticleDefinition *, G4ParticleDefinition * > pDefPair
std::vector< G4double > scalarMesonMix
static constexpr double GeV
Definition: G4SIunits.hh:217
G4double Mass2() const
T min(const T t1, const T t2)
brief Return the smallest of the two arguments