73{
74 auto & tA =
const_cast<G4Track&
>(trackA);
75 auto & tB =
const_cast<G4Track&
>(trackB);
77
79 pChanges->Initialize(trackA, trackB);
80
83
86 if (nbProducts != 0)
87 {
88 const G4double D1 = pMoleculeA->GetDiffusionCoefficient();
89 const G4double D2 = pMoleculeB->GetDiffusionCoefficient();
90 const G4double sqrD1 = D1 == 0. ? 0. : std::sqrt(D1);
91 const G4double sqrD2 = D2 == 0. ? 0. : std::sqrt(D2);
92 const G4double inv_numerator = 1./(sqrD1 + sqrD2);
93 const G4ThreeVector reactionSite = sqrD2 * inv_numerator * tA.GetPosition()
94 + sqrD1 * inv_numerator * tB.GetPosition();
95
97 auto randP = (1-u) * tA.GetPosition() + u * tB.GetPosition();
98
99 for (
G4int j = 0; j < nbProducts; ++j)
100 {
101 auto pProduct =
new G4Molecule(pReactionData->GetProduct(j));
102 auto pProductTrack = pProduct->BuildTrack(trackA.
GetGlobalTime(), (reactionSite + randP)/2);
103 pProductTrack->SetTrackStatus(
fAlive);
105 pChanges->AddSecondary(pProductTrack);
106 }
107 }
108 pChanges->KillParents(true);
109 return pChanges;
110}
G4Molecule * GetMolecule(const G4Track &track)
void UpdatePositionForReaction(G4Track &, G4Track &)
G4int GetNbProducts() const
Data * GetReactionData(Reactant *, Reactant *) const
void Push(G4Track *) override
static G4ITTrackHolder * Instance()
const G4MolecularConfiguration * GetMolecularConfiguration() const
G4double GetGlobalTime() const