81 theProjectilePotential(0)
142 theDefinition(aDefinition),
143 theFormationTime(aFormationTime),
144 thePosition(aPosition),
145 the4Momentum(a4Momentum),
146 theFermi3Momentum(0),
147 theTotal4Momentum(a4Momentum),
150 theProjectilePotential(0)
170 if (theDecayTable != 0)
194 G4bool * theDaughterIsShortLived = 0;
200 for (index =
nChannels - 1; index >= 0; index--)
205 if (nDaughters == 2 || nDaughters == 3)
213 theDaughterIsShortLived =
new G4bool[nDaughters];
215 for (n = 0; n < nDaughters; n++)
232 if ( !theDaughterIsShortLived[0] && !theDaughterIsShortLived[1] )
252 else if ( !theDaughterIsShortLived[0] && theDaughterIsShortLived[1] )
274 else if ( theDaughterIsShortLived[0] && !theDaughterIsShortLived[1] )
304 else if ( theDaughterIsShortLived[0] && theDaughterIsShortLived[1] )
331 G4int nShortLived = 0;
332 if ( theDaughterIsShortLived[0] )
336 if ( theDaughterIsShortLived[1] )
342 if ( theDaughterIsShortLived[2] )
348 if ( nShortLived == 0 )
357 else if ( nShortLived >= 1 )
375 G4double theMomRatio = theActualMom / thePoleMom;
387 delete [] theDaughterIsShortLived;
388 theDaughterIsShortLived = 0;
419 : theDefinition(nucleon->GetDefinition()),
421 thePosition(aPosition),
422 the4Momentum(a4Momentum),
423 theFermi3Momentum(nucleon->GetMomentum()),
426 theActualMass(nucleon->GetDefinition()->GetPDGMass()),
431 theProjectilePotential(0)
474 return (
this == & right);
481 return (
this != & right);
500 G4cerr <<
"Error condition encountered in G4KineticTrack::Decay()"<<
G4endl;
501 G4cerr <<
" track has no particle definition associated."<<
G4endl;
507 G4cerr <<
"Error condition encountered in G4KineticTrack::Decay()"<<
G4endl;
508 G4cerr <<
" particle definiton has no decay table associated."<<
G4endl;
517 if (theTotalActualWidth !=0)
522 G4bool isChannelBelowThreshold =
true;
523 const G4int maxNumberOfLoops = 10000;
524 G4int loopCounter = 0;
530 G4int theNumberOfDaughters;
542 for (index =
nChannels - 1; index >= 0; index--)
545 theCumActualWidth[index] = theSumActualWidth;
553 for (index =
nChannels - 1; index >= 0; index--)
555 if (r < theCumActualWidth[index])
564 delete [] theCumActualWidth;
568 G4cerr <<
"Error condition encountered in G4KineticTrack::Decay()"<<
G4endl;
569 G4cerr <<
" decay channel has 0x0 channel associated."<<
G4endl;
579 theDaughtersName1 =
"";
580 theDaughtersName2 =
"";
581 theDaughtersName3 =
"";
582 theDaughtersName4 =
"";
584 for (
G4int i=0; i < 4; i++) masses[i]=0.;
585 G4int shortlivedDaughters[4];
586 G4int numberOfShortliveds(0);
588 for (
G4int aD=0; aD < theNumberOfDaughters ; aD++)
594 shortlivedDaughters[numberOfShortliveds]=aD;
595 numberOfShortliveds++;
601 switch (theNumberOfDaughters)
607 theDaughtersName2 =
"";
608 theDaughtersName3 =
"";
609 theDaughtersName4 =
"";
614 theDaughtersName3 =
"";
615 theDaughtersName4 =
"";
616 if ( numberOfShortliveds == 1)
618 G4double massmax=theParentMass - SumLongLivedMass;
620 masses[shortlivedDaughters[0]]= aSampler.
SampleMass(aDaughter,massmax);
621 }
else if ( numberOfShortliveds == 2) {
628 masses[shortlivedDaughters[zero]]=aSampler.
SampleMass(aDaughter,massmax);
629 massmax=theParentMass - masses[shortlivedDaughters[zero]];
630 aDaughter=theDecayChannel->
GetDaughter(shortlivedDaughters[one]);
631 masses[shortlivedDaughters[one]]=aSampler.
SampleMass(aDaughter,massmax);
638 theDaughtersName4 =
"";
639 if ( numberOfShortliveds == 1)
641 G4double massmax=theParentMass - SumLongLivedMass;
643 masses[shortlivedDaughters[0]]= aSampler.
SampleMass(aDaughter,massmax);
651 if ( numberOfShortliveds == 1)
653 G4double massmax=theParentMass - SumLongLivedMass;
655 masses[shortlivedDaughters[0]]= aSampler.
SampleMass(aDaughter,massmax);
657 if ( theNumberOfDaughters > 4 ) {
659 ed <<
"More than 4 decay daughters: kept only the first 4" <<
G4endl;
669 for (
G4int i=0; i < 4; i++) sumDaughterMasses += masses[i];
670 if ( theParentMass - sumDaughterMasses > 0.0 ) isChannelBelowThreshold =
false;
672 }
while ( isChannelBelowThreshold && ++loopCounter < maxNumberOfLoops );
681 theNumberOfDaughters,
688 if(!theDecayProducts)
691 ed <<
"Error condition encountered: phase-space decay failed." <<
G4endl
693 <<
"\t the channel index is: "<< chosench <<
" of "<<
nChannels <<
"channels" <<
G4endl
694 <<
"\t " << theNumberOfDaughters <<
" daughter particles: "
695 << theDaughtersName1 <<
" " << theDaughtersName2 <<
" " << theDaughtersName3 <<
" "
696 << theDaughtersName4 <<
G4endl;
713 for (
G4int i=dEntries; i > 0; i--)
715 theDynamicParticle = theDecayProducts->
PopProducts();
719 momentumBalanceCMS += theDynamicParticle->
Get4Momentum();
721 energyMomentumBalance -= momentum;
726 delete theDynamicParticle;
728 delete theDecayProducts;
729 if(getenv(
"DecayEnergyBalanceCheck"))
730 std::cout <<
"DEBUGGING energy balance in cms and lab, charge baryon balance : "
731 << momentumBalanceCMS <<
" "
732 <<energyMomentumBalance <<
" "
736 return theDecayProductList;
752 std::sqrt(
std::max((((mass * mass) - (mass1 + xmass) * (mass1 + xmass)) *
753 ((mass * mass) - (mass1 - xmass) * (mass1 - xmass))),0.0)) *
754 BrWig(gamma2, mass2, xmass);
765 std::sqrt(
std::max((((mass * mass) - (mass1 + xmass) * (mass1 + xmass)) *
766 ((mass * mass) - (mass1 - xmass) * (mass1 - xmass))),0.0)) *
767 BrWig(gamma2, mass2, xmass);
781 BrWig(gamma2, mass2, xmass);
795 const G4double theUpperLimit = mass - xmass;
796 const G4int nIterations = 100;
808 const G4int nIterations = 100;
810 if (theLowerLimit>=theUpperLimit)
return 0.0;
814 theLowerLimit, theUpperLimit, nIterations);
815 return theIntegralOverMass2;
821 const G4int nIterations = 100;
823 if (theLowerLimit>=theUpperLimit)
return 0.0;
827 theLowerLimit, theUpperLimit, nIterations);
828 return theIntegralOverMass2;
836 const G4int nIterations = 100;
838 if (theLowerLimit>=theUpperLimit)
return 0.0;
842 theLowerLimit, theUpperLimit, nIterations);
843 return theIntegralOverMass2;
T max(const T t1, const T t2)
brief Return the largest of the two arguments
G4double GetPDGWidth() const
G4int operator==(const G4KineticTrack &right) const
void G4SwapObj(T *a, T *b)
std::ostringstream G4ExceptionDescription
G4DecayTable * GetDecayTable() const
const G4LorentzVector & Get4Momentum() const
static G4KaonZero * KaonZero()
G4LorentzVector the4Momentum
G4int GetNumberOfDaughters() const
const G4LorentzVector & GetTrackingMomentum() const
G4double IntegrateCMMomentum(const G4double lowerLimit) const
G4double * theDaughterMass
const G4String & GetParticleName() const
G4bool nucleon(G4int ityp)
const G4ParticleDefinition * GetDefinition() const
G4double GetPDGCharge() const
G4double theProjectilePotential
G4double SampleMass(const G4double poleMass, const G4double gamma, const G4double minMass, const G4double maxMass) const
G4double * theActualWidth
G4DynamicParticle * PopProducts()
G4double GetPDGMass() const
const G4ParticleDefinition * theDefinition
G4double IntegrandFunction1(G4double xmass) const
CascadeState theStateToNucleus
static G4AntiKaonZero * AntiKaonZero()
const G4String & GetParentName() const
G4double EvaluateCMMomentum(const G4double mass, const G4double *m_ij) const
static G4KaonZeroLong * KaonZeroLong()
G4ThreeVector thePosition
G4double IntegrandFunction3(G4double xmass) const
G4ParticleDefinition * GetDefinition() const
G4double IntegrandFunction4(G4double xmass) const
G4int GetBaryonNumber() const
G4int operator!=(const G4KineticTrack &right) const
static G4ThreadLocal G4double G4KineticTrack_Gmass
G4double theFormationTime
G4VDecayChannel * GetDecayChannel(G4int index) const
G4double GetActualMass() const
G4LorentzVector Get4Momentum() const
G4double IntegrateCMMomentum2() const
G4LorentzVector theTotal4Momentum
G4ParticleDefinition * GetDaughter(G4int anIndex)
const G4String & GetDaughterName(G4int anIndex) const
G4GLOB_DLL std::ostream G4cerr
G4double Simpson(T &typeT, F f, G4double a, G4double b, G4int n)
G4double G4ParticleHPJENDLHEData::G4double result
static G4KaonZeroShort * KaonZeroShort()
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
G4double GetFormationTime() const
G4KineticTrack & operator=(const G4KineticTrack &right)
G4int GetnChannels() const
G4double BrWig(const G4double Gamma, const G4double rmass, const G4double mass) const
G4LorentzVector theFermi3Momentum
G4bool IsShortLived() const
G4double IntegrandFunction2(G4double xmass) const
const G4ThreeVector & GetPosition() const
virtual G4DecayProducts * DecayIt(G4double mass=0.0)
G4double EvaluateTotalActualWidth()
void Set4Momentum(const G4LorentzVector &a4Momentum)
G4double * theDaughterWidth
G4KineticTrackVector * Decay()
static G4ThreadLocal G4double G4KineticTrack_xmass1
G4double GetMinimumMass(const G4ParticleDefinition *p) const