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        18846 :         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