117{
118
119
120
121
122#ifdef G4VERBOSE
124#endif
125
128
129
131
132
133 const G4int N_DAUGHTER = 3;
135 for (
G4int index = 0; index < N_DAUGHTER; ++index) {
137 }
138
139
142
144 delete parentparticle;
145
146
147 G4double daughtermomentum[N_DAUGHTER];
148
149
150 G4double pmax = (parentmass * parentmass - daughtermass[0] * daughtermass[0]) / 2. / parentmass;
153 const std::size_t MAX_LOOP = 10000;
154 for (std::size_t loop_counter = 0; loop_counter < MAX_LOOP; ++loop_counter) {
155
158 e = std::sqrt(p * p + daughtermass[0] * daughtermass[0]);
159 if (r < spectrum(p, e, parentmass, daughtermass[0])) break;
160 }
161
162
163
164 daughtermomentum[0] = p;
165 G4double costheta, sintheta, phi, sinphi, cosphi;
167 sintheta = std::sqrt((1.0 - costheta) * (1.0 + costheta));
169 sinphi = std::sin(phi);
170 cosphi = std::cos(phi);
171 G4ThreeVector direction0(sintheta * cosphi, sintheta * sinphi, costheta);
172 auto daughterparticle =
174 products->PushProducts(daughterparticle);
175
176
177
179 G4double vmass = std::sqrt((energy2 - daughtermomentum[0]) * (energy2 + daughtermomentum[0]));
180 G4double beta = -1.0 * daughtermomentum[0] / energy2;
182 G4double sinthetan = std::sqrt((1.0 - costhetan) * (1.0 + costhetan));
186
187 G4ThreeVector direction1(sinthetan * cosphin, sinthetan * sinphin, costhetan);
189 auto daughterparticle2 =
191
192
194 p4 = daughterparticle1->Get4Momentum();
195 p4.
boost(direction0.x() * beta, direction0.y() * beta, direction0.z() * beta);
196 daughterparticle1->Set4Momentum(p4);
197 p4 = daughterparticle2->Get4Momentum();
198 p4.
boost(direction0.x() * beta, direction0.y() * beta, direction0.z() * beta);
199 daughterparticle2->Set4Momentum(p4);
200 products->PushProducts(daughterparticle1);
201 products->PushProducts(daughterparticle2);
202 daughtermomentum[1] = daughterparticle1->GetTotalMomentum();
203 daughtermomentum[2] = daughterparticle2->GetTotalMomentum();
204
205
206#ifdef G4VERBOSE
208 G4cout <<
"G4TauLeptonicDecayChannel::DecayIt ";
209 G4cout <<
" create decay products in rest frame " <<
G4endl;
210 products->DumpInfo();
211 }
212#endif
213 return products;
214}
HepLorentzVector & boost(double, double, double)
G4double GetPDGMass() const
G4ParticleDefinition ** G4MT_daughters
void CheckAndFillParent()
G4ParticleDefinition * G4MT_parent
void CheckAndFillDaughters()