Line data    Source code

       1             : #ifndef CRPROPA_NUCLEARDECAY_H
       2             : #define CRPROPA_NUCLEARDECAY_H
       3             : 
       4             : #include "crpropa/Module.h"
       5             : 
       6             : #include <vector>
       7             : 
       8             : namespace crpropa {
       9             : /**
      10             :  * \addtogroup EnergyLosses
      11             :  * @{
      12             :  */
      13             : 
      14             : /**
      15             :  @class NuclearDecay
      16             :  @brief Nuclear decay of unstable nuclei.
      17             : 
      18             :  This module simulates the nuclear decay of unstable nuclei using data from NuDat2.
      19             :  All decay modes are considered: alpha, beta+- and gamma decay, as well as proton- and neutron dripping.
      20             :  The resulting non-hadronic secondary particles (e+, e-, neutrinos, gamma) can optionally be created.
      21             : 
      22             :  For details on the preprocessing of the NuDat2 data refer to "CRPropa3-data/calc_decay.py".
      23             :  */
      24             : class NuclearDecay: public Module {
      25             : private:
      26             :         double limit;
      27             :         bool haveElectrons;
      28             :         bool havePhotons;
      29             :         bool haveNeutrinos;
      30        1800 :         struct DecayMode {
      31             :                 int channel; // (#beta- #beta+ #alpha #proton #neutron)
      32             :                 double rate; // decay rate in [1/m]
      33             :                 std::vector<double> energy; // photon energies of ensuing gamma decays
      34             :                 std::vector<double> intensity; // probabilities of ensuing gamma decays
      35             :         };
      36             :         std::vector<std::vector<DecayMode> > decayTable; // decayTable[Z * 31 + N] = vector<DecayMode>
      37             :         std::string interactionTag = "ND";
      38             : 
      39             : public:
      40             :         /** Constructor.
      41             :          @param electrons               if true, add secondary photons as candidates
      42             :          @param photons                 if true, add secondary photons as candidates
      43             :          @param neutrinos               if true, add secondary neutrinos as candidates
      44             :          @param limit                   step size limit as fraction of mean free path
      45             :          */
      46             :         NuclearDecay(bool electrons = false, bool photons = false, bool neutrinos = false, double limit = 0.1);
      47             : 
      48             :         /** Limit the propagation step to a fraction of the mean free path
      49             :          * @param limit fraction of the mean free path
      50             :          */
      51             :         void setLimit(double limit);
      52             : 
      53             :         // decide if secondary electrons are added to the simulation    
      54             :         void setHaveElectrons(bool b);
      55             : 
      56             :         // decide if secondary photons are added to the simulation      
      57             :         void setHavePhotons(bool b);
      58             : 
      59             :         // decide if secondary neutrinos are added to the simulation    
      60             :         void setHaveNeutrinos(bool b);
      61             : 
      62             :         /** set a custom interaction tag to trace back this interaction
      63             :          * @param tag string that will be added to the candidate and output
      64             :          */
      65             :         void setInteractionTag(std::string tag);
      66             :         std::string getInteractionTag() const;
      67             : 
      68             :         void process(Candidate *candidate) const;
      69             :         void performInteraction(Candidate *candidate, int channel) const;
      70             :         void gammaEmission(Candidate *candidate, int channel) const;
      71             :         void betaDecay(Candidate *candidate, bool isBetaPlus) const;
      72             :         void nucleonEmission(Candidate *candidate, int dA, int dZ) const;
      73             : 
      74             :         /**
      75             :          Return the mean free path.
      76             :          This is not used in the simulation.
      77             :          @param id      PDG particle id
      78             :          @param gamma   Lorentz factor of particle
      79             :          @returns The mean free path [in meters]
      80             :          */
      81             :         double meanFreePath(int id, double gamma);
      82             : };
      83             : /** @}*/
      84             : 
      85             : } // namespace crpropa
      86             : 
      87             : #endif // CRPROPA_NUCLEARDECAY_H