ProposalProcessBase.hpp

/*
 * (c) Copyright 2020 CORSIKA Project, corsika-project@lists.kit.edu
 *
 * This software is distributed under the terms of the GNU General Public
 * Licence version 3 (GPL Version 3). See file LICENSE for a full version of
 * the license.
 */

#pragma once

#include <PROPOSAL/PROPOSAL.h>

#include <corsika/framework/core/ParticleProperties.hpp>
#include <corsika/framework/random/RNGManager.hpp>

#include <array>

namespace corsika::proposal {

  static constexpr std::array<Code, 7> tracked{
      Code::Photon, Code::Electron, Code::Positron, Code::MuMinus,
      Code::MuPlus, Code::TauPlus,  Code::TauMinus,
  };

  static std::map<Code, PROPOSAL::ParticleDef> particle = {
      {Code::Photon, PROPOSAL::GammaDef()},   {Code::Electron, PROPOSAL::EMinusDef()},
      {Code::Positron, PROPOSAL::EPlusDef()}, {Code::MuMinus, PROPOSAL::MuMinusDef()},
      {Code::MuPlus, PROPOSAL::MuPlusDef()},  {Code::TauMinus, PROPOSAL::TauMinusDef()},
      {Code::TauPlus, PROPOSAL::TauPlusDef()}};

  template <typename T>
  static auto cross_builder =
      [](PROPOSAL::Medium& m,
         corsika::units::si::HEPEnergyType
             emCut) { 
        auto p_cut = std::make_shared<const PROPOSAL::EnergyCutSettings>(
            0.5 * emCut / 1_MeV, 1, false);
        return PROPOSAL::GetStdCrossSections(T(), m, p_cut, true);
      };

  static std::map<Code, std::function<PROPOSAL::crosssection_list_t(
                            PROPOSAL::Medium&, corsika::units::si::HEPEnergyType)>>
      cross = {{Code::Photon, cross_builder<PROPOSAL::GammaDef>},
               {Code::Electron, cross_builder<PROPOSAL::EMinusDef>},
               {Code::Positron, cross_builder<PROPOSAL::EPlusDef>},
               {Code::MuMinus, cross_builder<PROPOSAL::MuMinusDef>},
               {Code::MuPlus, cross_builder<PROPOSAL::MuPlusDef>},
               {Code::TauMinus, cross_builder<PROPOSAL::TauMinusDef>},
               {Code::TauPlus, cross_builder<PROPOSAL::TauPlusDef>}};

  class ProposalProcessBase {
  protected:
    default_prng_type& RNG_ = RNGManager<>::getInstance().getRandomStream("proposal");

    std::unordered_map<std::size_t, PROPOSAL::Medium>
        media; 

    template <typename TEnvironment>
    ProposalProcessBase(TEnvironment const& _env);

    bool canInteract(Code pcode) const;

    using calc_key_t = std::pair<std::size_t, Code>;

    struct hash {
      size_t operator()(const calc_key_t& p) const noexcept;
    };

    virtual void buildCalculator(Code, NuclearComposition const&) = 0;

    template <typename Particle, typename Calculators>
    auto getCalculator(Particle& vP, Calculators& calc) {
      const auto& comp = vP.getNode()->getModelProperties().getNuclearComposition();
      auto calc_it = calc.find(std::make_pair(comp.getHash(), vP.getPID()));
      if (calc_it != calc.end()) return calc_it;
      buildCalculator(vP.getPID(), comp);
      return getCalculator(vP, calc);
    }
  };
} // namespace corsika::proposal

#include <corsika/detail/modules/proposal/ProposalProcessBase.inl>