63 if (theParentName ==
"mu+") {
70 }
else if (theParentName ==
"mu-") {
80 G4cout <<
"G4MuonDecayChannel:: constructor :";
81 G4cout <<
" parent particle is not muon but ";
138 const int N_DAUGHTER=3;
143 for (
G4int index=0; index<N_DAUGHTER; index++){
145 sumofdaughtermass += daughtermass[index];
153 delete parentparticle;
156 G4double daughtermomentum[N_DAUGHTER];
158 G4double xmax = (1.0+daughtermass[0]*daughtermass[0]/parentmass/parentmass);
164 G4double EMax=parentmass/2-daughtermass[0];
166 const size_t MAX_LOOP=1000;
168 for (
size_t loop1=0; loop1 <MAX_LOOP; ++loop1){
170 for (
size_t loop2 =0; loop2<MAX_LOOP; ++loop2){
173 if (gam <= x*(1.-x))
break;
177 if ( Ene >= (1.-Ee))
break;
185 G4double costheta,sintheta,rphi,rtheta,rpsi;
186 costheta= 1.-2./Ee-2./Ene+2./Ene/Ee;
187 sintheta=std::sqrt(1.-costheta*costheta);
195 rot.
set(rphi,rtheta,rpsi);
198 daughtermomentum[0]=std::sqrt(Ee*Ee*EMax*EMax+2.0*Ee*EMax * daughtermass[0]);
209 daughtermomentum[1]=std::sqrt(Ene*Ene*EMax*EMax+2.0*Ene*EMax * daughtermass[1]);
219 daughtermomentum[2]=std::sqrt(Enm*Enm*EMax*EMax +2.0*Enm*EMax*daughtermass[2]);
220 G4ThreeVector direction2(-Ene/Enm*sintheta,0,-Ee/Enm-Ene/Enm*costheta);
225 direction2 * daughtermomentum[2]);
231 G4cout <<
"G4MuonDecayChannel::DecayIt ";
232 G4cout <<
" create decay products in rest frame " <<
G4endl;