112{
113
114
115
116
117#ifdef G4VERBOSE
119#endif
120
123
124
126
127 EMMU = parentmass;
128
129
132 for (
G4int index=0; index<3; index++){
134 sumofdaughtermass += daughtermass[index];
135 }
136
137 EMASS = daughtermass[0];
138
139
142
144 delete parentparticle;
145
146
147
152
154
157
158 G4double W_mue = (EMMU*EMMU+EMASS*EMASS)/(2.*EMMU);
160
162
163
164
165
166
167
168
169
170
171 do{
172
173
174
176
177 x = x0 + rndm*(1.-x0);
178
180
182
183 F_IS = 1./6.*(-2.*x_squared+3.*x-x0_squared);
184 F_AS = 1./6.*std::sqrt(x_squared-x0_squared)*(2.*x-2.+std::sqrt(1.-x0_squared));
185
186 G_IS = 2./9.*(michel_rho-0.75)*(4.*x_squared-3.*x-x0_squared);
187 G_IS = G_IS + michel_eta*(1.-x)*x0;
188
189 G_AS = 3.*(michel_xsi-1.)*(1.-x);
190 G_AS = G_AS+2.*(michel_xsi*michel_delta-0.75)*(4.*x-4.+std::sqrt(1.-x0_squared));
191 G_AS = 1./9.*std::sqrt(x_squared-x0_squared)*G_AS;
192
193 F_IS = F_IS + G_IS;
194 F_AS = F_AS + G_AS;
195
196
197
199
200 G4double F = 6.*F_IS + R_IS/std::sqrt(x_squared-x0_squared);
201
202
203
205
207
208 ctheta = 2.*rndm-1.;
209
210 G4double G = 6.*F_AS - R_AS/std::sqrt(x_squared-x0_squared);
211
212 FG = std::sqrt(x_squared-x0_squared)*F*(1.+(G/F)*ctheta);
213
214 if(FG>FG_max){
215 G4cout<<
"***Problem in Muon Decay *** : FG > FG_max"<<
G4endl;
216 FG_max = FG;
217 }
218
220
221 }while(FG<rndm*FG_max);
222
224
226
228
229 if(energy < EMASS) energy = EMASS;
230
231
233
234 daughtermomentum[0] = std::sqrt(energy*energy - EMASS*EMASS);
235
236 G4double stheta = std::sqrt(1.-ctheta*ctheta);
239
240
241
245
247
248 direction0.rotateUz(parent_polarization);
249
252
254
255
256
257
258 G4double energy2 = parentmass*(1.0 - x/2.0);
259 G4double vmass = std::sqrt((energy2-daughtermomentum[0])*(energy2+daughtermomentum[0]));
260 G4double beta = -1.0*daughtermomentum[0]/energy2;
262 G4double sinthetan = std::sqrt((1.0-costhetan)*(1.0+costhetan));
266
267 G4ThreeVector direction1(sinthetan*cosphin,sinthetan*sinphin,costhetan);
272
273
276 p4.
boost( direction0.x()*beta, direction0.y()*beta, direction0.z()*beta);
279 p4.
boost( direction0.x()*beta, direction0.y()*beta, direction0.z()*beta);
285
286
287#ifdef G4VERBOSE
289 G4cout <<
"G4MuonDecayChannelWithSpin::DecayIt ";
290 G4cout <<
" create decay products in rest frame " <<
G4endl;
292 }
293#endif
294 return products;
295}
G4DLLIMPORT std::ostream G4cout
HepLorentzVector & boost(double, double, double)
G4int PushProducts(G4DynamicParticle *aParticle)
G4LorentzVector Get4Momentum() const
void Set4Momentum(const G4LorentzVector &momentum)
G4double GetTotalMomentum() const
G4double GetPDGMass() const
G4int GetVerboseLevel() const
G4ParticleDefinition * parent
G4ParticleDefinition ** daughters