MuMinusCaptureFix Class Reference

List of all members.


Detailed Description

class MuMinusCaptureFix will fixup the G4MuonMinusCaptureAtRest process.

The constructor replaces muMinusCaptureAtRest with this class, so all that is needed is to create an instance of this class with the desired arguments. That must occur AFTER the physics processes are completely set up but before any tracking.

This class adds additional neutrons to mu- capture. The neutrons are added with a Poisson distribution with mean=neutronMeanNumber, and have an exponential distribution in kinetic energy: (1/neutronMeanEnergy)*exp(-KE/neutronMeanEnergy) As the muonic atom cascades to its ground state it forgets the incident mu- direction, so neutrons are generated isotropically in the lab.

The extra neutrons are added only to those captures that are hadronic (i.e. not decay in orbit). The value of neutronMeanNumber should reflect this. The default values correspond to Al.

Author: Tom Roberts Date: September 8, 2010 This class is completely independent of G4beamline and can be copied into other code (keep the above include-s and definitions).

Public Member Functions

 MuMinusCaptureFix (double _neutronMeanNumber=DEFAULT_NEUTRONMEANNUMBER, double _neutronMeanEnergy=DEFAULT_NEUTRONMEANENERGY)
 Constructor.
virtual G4VParticleChange * AtRestDoIt (const G4Track &, const G4Step &)
 AtRestDoit() from G4MuonMinusCaptureAtRest.
G4Track * createExtraNeutron (const G4ThreeVector &position, double time)
 createExtraNeutron() creates an extra neutron with the given KE distribution, isotropically at position.

Private Attributes

double neutronMeanNumber
double neutronMeanEnergy
G4ParticleChange change


Constructor & Destructor Documentation

MuMinusCaptureFix::MuMinusCaptureFix ( double  _neutronMeanNumber = DEFAULT_NEUTRONMEANNUMBER,
double  _neutronMeanEnergy = DEFAULT_NEUTRONMEANENERGY 
)

Constructor.

References neutronMeanEnergy, and neutronMeanNumber.

00079                         : G4MuonMinusCaptureAtRest("muMinusCaptureAtRestFIX"),
00080                           change()
00081 {
00082         neutronMeanNumber = _neutronMeanNumber;
00083         neutronMeanEnergy = _neutronMeanEnergy;
00084 
00085         // replace G4MuonMinusCaptureAtRest with this MuMinusCaptureFix.
00086         G4ParticleDefinition *mu = G4MuonMinus::MuonMinus();
00087         G4ProcessManager *pm = mu->GetProcessManager();
00088 again:  G4ProcessVector *pv = pm->GetProcessList();
00089         int n = pv->size();
00090         for(int i=0; i<n; ++i) {
00091                 G4VProcess *process = (*pv)[i];
00092                 if(process->GetProcessName() == "muMinusCaptureAtRest") {
00093                         pm->RemoveProcess(process);
00094                         goto again;
00095                 }
00096         }
00097         pm->AddProcess(this);
00098         pm->SetProcessOrdering(this,idxAtRest);
00099 }


Member Function Documentation

G4VParticleChange * MuMinusCaptureFix::AtRestDoIt ( const G4Track &  track,
const G4Step &  step 
) [virtual]

AtRestDoit() from G4MuonMinusCaptureAtRest.

References change, createExtraNeutron(), and neutronMeanNumber.

00104 {
00105         int nExtraNeutrons = CLHEP::RandPoissonQ::shoot(neutronMeanNumber);
00106 
00107         // localtime will become the earliest time of any secondary hadron
00108         // from the original process; extra neutrons are not added to
00109         // decay in orbit captures (events with no secondary hadrons).
00110         // Note the capture gives prompt gammas which must not be considered.
00111         double localtime=DBL_MAX;
00112         bool haveHadron=false;
00113 
00114         // call original process's AtRestDoit()
00115         G4VParticleChange *p=G4MuonMinusCaptureAtRest::AtRestDoIt(track,step);
00116         int n0=p->GetNumberOfSecondaries();
00117         G4ThreeVector position=track.GetPosition();
00118 
00119         // setup a ParticleChange with room for the extra neutrons (this is a
00120         // copy of all modifications to aParticleChange in the original code)
00121         change.Initialize(track);
00122         change.SetNumberOfSecondaries(n0+nExtraNeutrons);
00123         for(int i=0; i<n0; ++i) {
00124                 G4Track *t=p->GetSecondary(i);
00125                 change.AddSecondary(new G4Track(*t)); // avoid double deletes
00126                 G4String type=t->GetDefinition()->GetParticleType();
00127                 if(type == "baryon" || type == "meson" || type == "nucleus") {
00128                         haveHadron = true;
00129                         double time=t->GetGlobalTime();
00130                         if(time < localtime) localtime = time;
00131                 }
00132         }
00133         change.ProposeLocalEnergyDeposit(0.0);
00134         change.ProposeTrackStatus(fStopAndKill);
00135         p->Clear();
00136 
00137         // store extra neutrons
00138         if(haveHadron) {
00139             for(int i=0; i<nExtraNeutrons; ++i) {
00140                 G4Track *neutron = createExtraNeutron(position,localtime);
00141                 change.AddSecondary(neutron);
00142             }
00143         }
00144 
00145         return &change;
00146 }

G4Track * MuMinusCaptureFix::createExtraNeutron ( const G4ThreeVector &  position,
double  time 
)

createExtraNeutron() creates an extra neutron with the given KE distribution, isotropically at position.

References KE, M_PI, and neutronMeanEnergy.

Referenced by AtRestDoIt().

00150 {
00151         G4ParticleDefinition *neutron=G4Neutron::Neutron();
00152         double KE=CLHEP::RandExponential::shoot(neutronMeanEnergy);
00153         double costheta=2.0*(G4UniformRand()-0.5);
00154         double sintheta=sqrt(1.0-costheta*costheta);
00155         double phi=2.0*M_PI*G4UniformRand();
00156         G4ThreeVector direction(sintheta*cos(phi),sintheta*sin(phi),costheta);
00157         G4DynamicParticle *dyn = new G4DynamicParticle(neutron,direction,KE);
00158         G4Track *track = new G4Track(dyn,time,position);
00159         return track;
00160 }


Member Data Documentation

Referenced by AtRestDoIt(), and MuMinusCaptureFix().

G4ParticleChange MuMinusCaptureFix::change [private]

Referenced by AtRestDoIt().


The documentation for this class was generated from the following file:
g4beamline