1 #include <simulation/longlivedneutral/G4LongLivedNeutral.h>
2 #include <simulation/longlivedneutral/G4LongLivedNeutralDecay.h>
4 #include "G4PhysicalConstants.hh"
5 #include "G4SystemOfUnits.hh"
6 #include "G4DynamicParticle.hh"
7 #include "G4PrimaryParticle.hh"
8 #include "G4DecayProducts.hh"
9 #include "G4DecayTable.hh"
10 #include "G4VDecayChannel.hh"
11 #include "G4PhysicsLogVector.hh"
12 #include "G4ParticleChangeForDecay.hh"
19 : G4VRestDiscreteProcess(processName, fDecay),
22 fRemainderLifeTime(-1.0)
26 SetProcessSubType(
static_cast<int>(DECAY));
29 if (GetVerboseLevel() > 1) {
30 G4cout <<
"G4LongLivedNeutralDecay constructor " <<
" Name:" << processName << G4endl;
46 if (!aParticleType.GetPDGStable()) {
55 const G4DynamicParticle* dynamicParticle = aTrack.GetDynamicParticle();
56 const G4PrimaryParticle* primaryParticle = dynamicParticle->GetPrimaryParticle();
57 G4double pTime = primaryParticle->GetProperTime();
65 const G4DynamicParticle* dynamicParticle = aTrack.GetDynamicParticle();
66 const G4PrimaryParticle* primaryParticle = dynamicParticle->GetPrimaryParticle();
67 G4double pTime = primaryParticle->GetProperTime();
68 G4double aMass = primaryParticle->GetMass();
69 G4double aMomentum = primaryParticle->GetTotalMomentum();
70 return c_light * pTime * aMomentum / aMass;
90 const G4DynamicParticle* aParticle = aTrack.GetDynamicParticle();
91 const G4PrimaryParticle* aPrimaryParticle = aParticle->GetPrimaryParticle();
94 const G4DecayProducts* o_products = (aParticle->GetPreAssignedDecayProducts());
95 G4bool isPreAssigned = (o_products !=
nullptr);
96 G4DecayProducts* products =
nullptr;
100 products =
new G4DecayProducts(*o_products);
103 G4ExceptionDescription ed;
104 ed <<
"LongLivedNeutral particle:"
105 <<
" decay probability exist but no pre-assigned decay products defined "
106 <<
"- the particle will be killed;\n";
107 G4Exception(
"G4LongLivedNeutralDecay::DecayIt ",
108 "DECAY101", JustWarning, ed);
115 ClearNumberOfInteractionLengthLeft();
120 G4double ParentEnergy = aPrimaryParticle->GetTotalEnergy();
121 G4double ParentMass = aPrimaryParticle->GetMass();
122 if (ParentEnergy < ParentMass) {
123 G4ExceptionDescription ed;
124 ed <<
"Total Energy is less than its mass - increased the energy"
125 <<
"\n Particle: " <<
"LongLivedNeutral"
126 <<
"\n Energy:" << ParentEnergy / MeV <<
"[MeV]"
127 <<
"\n Mass:" << ParentMass / MeV <<
"[MeV]";
128 G4Exception(
"G4LongLivedNeutralDecay::DecayIt ",
129 "DECAY102", JustWarning, ed);
130 ParentEnergy = ParentMass;
133 G4ThreeVector ParentDirection(aPrimaryParticle->GetMomentumDirection());
136 const G4DynamicParticle dParticle(aPrimaryParticle->GetParticleDefinition(),
137 aPrimaryParticle->GetMomentumDirection(),
138 aPrimaryParticle->GetKineticEnergy(),
139 aPrimaryParticle->GetMass());
141 products->SetParentParticle(dParticle);
143 G4double energyDeposit = 0.0;
144 G4double finalGlobalTime = aTrack.GetGlobalTime();
145 G4double finalLocalTime = aTrack.GetLocalTime();
146 if (aTrack.GetTrackStatus() == fStopButAlive) {
150 energyDeposit += aPrimaryParticle->GetKineticEnergy();
151 if (isPreAssigned) products->Boost(ParentEnergy, ParentDirection);
154 products->Boost(ParentEnergy, ParentDirection);
158 G4int numberOfSecondaries = products->entries();
161 if (GetVerboseLevel() > 1) {
162 G4cout <<
"G4LongLivedNeutralDecay::DoIt : Decay vertex :";
163 G4cout <<
" Time: " << finalGlobalTime / ns <<
"[ns]";
164 G4cout <<
" proper time: " << finalLocalTime / ns <<
"[ns]";
165 G4cout <<
" X:" << (aTrack.GetPosition()).x() / cm <<
"[cm]";
166 G4cout <<
" Y:" << (aTrack.GetPosition()).y() / cm <<
"[cm]";
167 G4cout <<
" Z:" << (aTrack.GetPosition()).z() / cm <<
"[cm]";
169 G4cout <<
"G4LongLivedNeutralDecay::DoIt : decay products in Lab. Frame" << G4endl;
170 products->DumpInfo();
175 G4ThreeVector currentPosition;
176 const G4TouchableHandle thand = aTrack.GetTouchableHandle();
177 for (index = 0; index < numberOfSecondaries; index++) {
179 currentPosition = aTrack.GetPosition();
181 G4Track* secondary =
new G4Track(products->PopProducts(),
185 secondary->SetGoodForTrackingFlag();
186 secondary->SetTouchableHandle(thand);
198 ClearNumberOfInteractionLengthLeft();
206 currentInteractionLength = -1.0;
207 ResetNumberOfInteractionLengthLeft();
216 ClearNumberOfInteractionLengthLeft();
218 currentInteractionLength = -1.0;
223 const G4Track& track,
224 G4double previousStepSize,
225 G4ForceCondition* condition
229 *condition = NotForced;
232 G4double pTime = track.GetDynamicParticle()->GetPrimaryParticle()->GetProperTime();
240 const G4DynamicParticle* dynamicParticle = track.GetDynamicParticle();
241 const G4PrimaryParticle* primaryParticle = dynamicParticle->GetPrimaryParticle();
242 G4double aMass = primaryParticle->GetMass();
243 G4double aMomentum = primaryParticle->GetTotalMomentum();
244 rvalue *= aMomentum / aMass;
249 return previousStepSize * DBL_MAX;
254 const G4Track& track,
255 G4ForceCondition* condition
259 *condition = NotForced;
260 G4double pTime = track.GetDynamicParticle()->GetPrimaryParticle()->GetProperTime();
274 const G4Track& aTrack,
278 if ((aTrack.GetTrackStatus() == fStopButAlive) ||
279 (aTrack.GetTrackStatus() == fStopAndKill)) {