50 arho(0.5), aphi(0.),
an(-0.5), ala(-0.75), aksi(-1.), alft(0.5)
65 #ifdef debug_QGSMfragmentation // Uzhi Oct. 2014
69 <<
"------------------------------------"<<
G4endl;
85 #ifdef debug_QGSMfragmentation // Uzhi Oct. 2014
86 if ( LeftVector != 0 )
G4cout<<
"Non fragmentable - the string is converted to one hadron "<<
G4endl;
89 if ( LeftVector != 0 )
return LeftVector;
91 #ifdef debug_QGSMfragmentation // Uzhi Oct. 2014
102 G4bool success=
false, inner_sucess=
true;
106 #ifdef debug_QGSMfragmentation // Uzhi Oct. 2014
117 RightVector->clear();
120 const G4int maxNumberOfLoops = 1000;
121 G4int loopCounter = -1;
122 while (!
StopFragmenting(currentString) && ++loopCounter < maxNumberOfLoops )
125 #ifdef debug_QGSMfragmentation // Uzhi Oct. 2014
133 #ifdef debug_QGSMfragmentation // Uzhi Oct. 2014
137 LeftVector->push_back(Hadron);
139 RightVector->push_back(Hadron);
141 delete currentString;
142 currentString=newString;
146 #ifdef debug_QGSMfragmentation // Uzhi Oct. 2014
147 G4cout<<
"abandon ... start from the beginning ---------------"<<
G4endl;
151 if (newString)
delete newString;
156 if ( loopCounter >= maxNumberOfLoops ) {
161 #ifdef debug_QGSMfragmentation // Uzhi Oct. 2014
162 G4cout<<
"Split remaining string into 2 final hadrons."<<
G4endl;
166 SplitLast(currentString,LeftVector, RightVector) )
170 delete currentString;
173 delete theStringInCMS;
185 while(!RightVector->empty())
187 LeftVector->push_back(RightVector->back());
188 RightVector->erase(RightVector->end()-1);
196 for(
size_t C1 = 0;
C1 < LeftVector->size();
C1++)
200 Momentum = toObserverFrame*Momentum;
203 Momentum = toObserverFrame*Coordinate;
219 #ifdef debug_QGSMfragmentation // Uzhi Oct. 2014
224 G4double invD1(0.),invD2(0.), r1(0.),r2(0.),r12(0.);
226 G4int absCode = std::abs( PartonEncoding );
230 q1 = absHadronCode/1000; q2 = (absHadronCode % 1000)/100; q3 = (absHadronCode % 100)/10;
232 G4bool StrangeHadron = (q1 == 3) || (q2 == 3) || (q3 == 3);
236 if(absCode == 1 || absCode == 2)
238 if(absHadronCode < 1000)
240 if( !StrangeHadron ) {d1=2.0; d2 = -
arho +
alft;}
244 if( !StrangeHadron ) {d1=0.0; d2 =
arho - 2.0*
an +
alft;}
248 else if(absCode == 3)
250 if(absHadronCode < 1000){d1=1.0 -
aphi; d2 = -
arho +
alft;}
253 }
else throw G4HadronicException(__FILE__, __LINE__,
"Unknown PDGencoding in G4QGSMFragmentation::G4LightConeZ");
255 #ifdef debug_QGSMfragmentation // Uzhi Oct. 2014
261 invD1=1./
d1; invD2=1./
d2;
263 const G4int maxNumberOfLoops = 10000;
264 G4int loopCounter = 0;
271 }
while( ( (r12 > 1.0) || !((zmin <= z)&&(z <= zmax))) && ++loopCounter < maxNumberOfLoops );
272 if ( loopCounter >= maxNumberOfLoops ) {
273 z = 0.5*(zmin + zmax);
280 if(absCode == 1103 || absCode == 2101 ||
281 absCode == 2203 || absCode == 2103)
283 if(absHadronCode < 1000)
285 if( !StrangeHadron ) {d1=1.0; d2=
arho - 2.0*
an +
alft;}
293 #ifdef debug_QGSMfragmentation // Uzhi Oct. 2014
298 invD1=1./
d1; invD2=1./
d2;
300 const G4int maxNumberOfLoops = 10000;
301 G4int loopCounter = 0;
308 }
while( ( (r12 > 1.0) || !((zmin <= z)&&(z <= zmax))) && ++loopCounter < maxNumberOfLoops );
309 if ( loopCounter >= maxNumberOfLoops ) {
310 z = 0.5*(zmin + zmax);
315 else if(absCode == 3101 || absCode == 3103 ||
316 absCode == 3201 || absCode == 3203)
326 const G4int maxNumberOfLoops = 1000;
327 G4int loopCounter = 0;
334 while( (
G4UniformRand() > yf) && ++loopCounter < maxNumberOfLoops );
352 #ifdef debug_QGSMfragmentation // Uzhi Oct. 2014
354 G4cout<<
"String 4 mom, String M "<<
string->Get4Momentum()<<
" "<<
string->Mass()<<
G4endl;
355 G4cout<<
"HadM MinimalStringMassLeft StringM hM+sM "<<HadronMass<<
" "<<MinimalStringMass<<
" "
356 <<
string->Mass()<<
" "<<HadronMass+MinimalStringMass<<
G4endl;
359 if(HadronMass + MinimalStringMass > string->
Mass())
361 #ifdef debug_QGSMfragmentation // Uzhi Oct. 2014
362 G4cout<<
"Mass of the string is not sufficient to produce the hadron!"<<
G4endl;
368 G4double StringMT2 =
string->MassT2();
369 G4double StringMT = std::sqrt(StringMT2);
372 String4Momentum.setPz(0.);
376 G4double HadronMassT2, ResidualMassT2;
386 RemSysPt = StringPt - HadronPt;
388 HadronMassT2 =
sqr(HadronMass) + HadronPt.mag2();
389 ResidualMassT2=
sqr(MinimalStringMass) + RemSysPt.mag2();
391 }
while(std::sqrt(HadronMassT2) + std::sqrt(ResidualMassT2) > StringMT);
396 G4double Pz2 = (
sqr(StringMT2 - HadronMassT2 - ResidualMassT2) -
397 4*HadronMassT2 * ResidualMassT2)/4./StringMT2;
399 if(Pz2 < 0 ) {
return 0;}
404 G4double zMin = (std::sqrt(HadronMassT2+Pz2) - Pz)/std::sqrt(StringMT2);
405 G4double zMax = (std::sqrt(HadronMassT2+Pz2) + Pz)/std::sqrt(StringMT2);
419 if (zMin >= zMax)
return 0;
423 HadronPt.x(), HadronPt.y());
429 (z *
string->LightConeDecay() -
430 HadronMassT2/(z *
string->LightConeDecay())));
431 G4double HadronE = 0.5* (z *
string->LightConeDecay() +
432 HadronMassT2/(z *
string->LightConeDecay()));
436 #ifdef debug_QGSMfragmentation // Uzhi Oct. 2014
437 G4cout<<
"string->GetDecayDirection() string->LightConeDecay() "
438 <<
string->GetDecayDirection()<<
" "<<
string->LightConeDecay()<<
G4endl;
439 G4cout<<
"HadronPt,HadronE "<<HadronPt<<
" "<<HadronE<<
G4endl;
457 G4ThreeVector ClusterVel =
string->Get4Momentum().boostVector();
458 G4double ResidualMass =
string->Mass();
460 #ifdef debug_QGSMfragmentation // Uzhi Oct. 2014
461 G4cout<<
"Split last-----------------------------------------"<<
G4endl;
462 G4cout<<
"StrMass "<<ResidualMass<<
" q's "
463 <<
string->GetLeftParton()->GetParticleName()<<
" "
464 <<
string->GetRightParton()->GetParticleName()<<
G4endl;
468 G4int cClusterInterrupt = 0;
470 const G4int maxNumberOfLoops = 1000;
471 G4int loopCounter = 0;
480 string->SetLeftPartonStable();
490 IsParticle=(
string->GetLeftParton()->GetPDGEncoding()>0) ? -1 : +1;
492 IsParticle=(
string->GetLeftParton()->GetPDGEncoding()>0) ? +1 : -1;
496 quark = QuarkPair.second;
507 while ( (ResidualMass <= LeftHadron->GetPDGMass() + RightHadron->
GetPDGMass() + ClusterMassCut)
508 && ++loopCounter < maxNumberOfLoops );
510 if ( loopCounter >= maxNumberOfLoops ) {
518 &RightMom, RightHadron->
GetPDGMass(), ResidualMass);
519 LeftMom.boost(ClusterVel);
520 RightMom.boost(ClusterVel);
522 #ifdef debug_QGSMfragmentation // Uzhi Oct. 2014
524 G4cout<<
"Left Hadrom P M "<<LeftMom<<
" "<<LeftMom.mag()<<
G4endl;
525 G4cout<<
"Right Hadrom P M "<<RightMom<<
" "<<RightMom.mag()<<
G4endl;
528 LeftVector->push_back(
new G4KineticTrack(LeftHadron, 0, Pos, LeftMom));
529 RightVector->push_back(
new G4KineticTrack(RightHadron, 0, Pos, RightMom));
549 string->Get4Momentum().mag2();
559 G4double r_val =
sqr(InitialMass*InitialMass - Mass*Mass - AntiMass*AntiMass) -
sqr(2.*Mass*AntiMass);
560 G4double Pabs = (r_val > 0.)? std::sqrt(r_val)/(2.*InitialMass) : 0;
564 G4double st = std::sqrt(1. - pz * pz)*Pabs;
570 Mom->setPx(px); Mom->setPy(py); Mom->setPz(pz);
571 Mom->setE(std::sqrt(Pabs*Pabs + Mass*Mass));
573 AntiMom->setPx(-px); AntiMom->setPy(-py); AntiMom->setPz(-pz);
574 AntiMom->setE (std::sqrt(Pabs*Pabs + AntiMass*AntiMass));
595 G4int Swap = stableQuarkEncoding;
596 stableQuarkEncoding = decayQuarkEncoding;
597 decayQuarkEncoding = Swap;
600 G4int IsParticle=(decayQuarkEncoding>0) ? -1 : +1;
609 G4int QuarkEncoding=QuarkPair.second->GetPDGEncoding();
610 G4int i10 =
std::max(std::abs(QuarkEncoding), std::abs(stableQuarkEncoding));
611 G4int i20 =
std::min(std::abs(QuarkEncoding), std::abs(stableQuarkEncoding));
613 G4int NewDecayEncoding = -1*IsParticle*(i10 * 1000 + i20 * 100 + spin);
631 created = QuarkPair.second;
G4ParticleDefinition * BuildHighSpin(G4ParticleDefinition *black, G4ParticleDefinition *white)
static G4Pow * GetInstance()
G4ParticleDefinition * GetRightParton(void) const
G4double powA(G4double A, G4double y) const
G4ParticleDefinition * QuarkSplitup(G4ParticleDefinition *decay, G4ParticleDefinition *&created)
virtual G4bool StopFragmenting(const G4FragmentingString *const string)
CLHEP::Hep3Vector G4ThreeVector
virtual G4bool SplitLast(G4FragmentingString *string, G4KineticTrackVector *LeftVector, G4KineticTrackVector *RightVector)
CLHEP::HepLorentzRotation G4LorentzRotation
const G4ThreeVector & GetPosition() const
G4ExcitedString * CPExcited(const G4ExcitedString &string)
void SetFormationTime(G4double aFormationTime)
void SetStrangenessSuppression(G4double aValue)
const G4LorentzVector & Get4Momentum() const
G4ParticleDefinition * Build(G4ParticleDefinition *black, G4ParticleDefinition *white)
G4int StringLoopInterrupt
virtual G4bool IsFragmentable(const G4FragmentingString *const string)
G4Parton * GetLeftParton(void) const
G4int ClusterLoopInterrupt
G4int GetPDGEncoding() const
G4ThreeVector SampleQuarkPt(G4double ptMax=-1.)
virtual G4KineticTrackVector * FragmentString(const G4ExcitedString &theString)
const G4String & GetParticleSubType() const
G4double FragmentationMass(const G4FragmentingString *const string, Pcreate build=0, pDefPair *pdefs=0)
const G4String & GetParticleName() const
G4ParticleDefinition * GetDecayParton() const
pDefPair CreatePartonPair(G4int NeedParticle, G4bool AllowDiquarks=true)
const G4ThreeVector & GetPosition() const
G4ParticleDefinition * GetLeftParton(void) const
void SetDiquarkBreakProbability(G4double aValue)
virtual G4LorentzVector * SplitEandP(G4ParticleDefinition *pHadron, G4FragmentingString *string, G4FragmentingString *newString)
G4double GetStrangeSuppress()
G4GLOB_DLL std::ostream G4cout
ParticleList decay(Cluster *const c)
Carries out a cluster decay.
G4KineticTrackVector * LightFragmentationTest(const G4ExcitedString *const theString)
G4LorentzVector Get4Momentum() const
G4double GetFormationTime() const
std::pair< G4ParticleDefinition *, G4ParticleDefinition * > pDefPair
G4LorentzRotation TransformToAlignedCms()
G4ParticleDefinition * FindParticle(G4int Encoding)
void SetPosition(const G4ThreeVector aPosition)
G4int GetDecayDirection() const
void Set4Momentum(const G4LorentzVector &a4Momentum)
virtual G4ParticleDefinition * DiQuarkSplitup(G4ParticleDefinition *decay, G4ParticleDefinition *&created)
G4Parton * GetRightParton(void) const
G4HadronBuilder * hadronizer
G4double GetPDGMass() const
virtual void Sample4Momentum(G4LorentzVector *Mom, G4double Mass, G4LorentzVector *AntiMom, G4double AntiMass, G4double InitialMass)
T max(const T t1, const T t2)
brief Return the largest of the two arguments
T min(const T t1, const T t2)
brief Return the smallest of the two arguments
G4double DiquarkBreakProb
G4int GetDirection(void) const
virtual G4double GetLightConeZ(G4double zmin, G4double zmax, G4int PartonEncoding, G4ParticleDefinition *pHadron, G4double Px, G4double Py)
const G4LorentzVector & Get4Momentum() const
void CalculateHadronTimePosition(G4double theInitialStringMass, G4KineticTrackVector *)
const G4ParticleDefinition * GetDefinition() const
G4KineticTrack * Splitup(G4FragmentingString *string, G4FragmentingString *&newString)
void SetDiquarkSuppression(G4double aValue)
CLHEP::HepLorentzVector G4LorentzVector