AliceO2Group · sawenzel · Jan 18, 2023 · Dec 21, 2022 · Jan 18, 2023
@@ -24,7 +24,8 @@ namespace conf
 enum class EG4Physics {
  kFTFP_BERT_optical = 0,         /* just ordinary */
  kFTFP_BERT_optical_biasing = 1, /* with biasing enabled */
-  kFTFP_INCLXX_optical = 2        /* special INCL++ version */
+  kFTFP_INCLXX_optical = 2,       /* special INCL++ version */
+  kFTFP_BERT_HP_optical = 3       /* enable low energy neutron transport */
 };

 // parameters to influence the G4 engine

@@ -37,7 +37,7 @@ struct SimCutParams : public o2::conf::ConfigurableParamHelper<SimCutParams> {
  float tunnelZ = 1900;       // Z-value from where we apply maxRTrackingZDC (default value taken from standard "hall" dimensions)

  float globalDensityFactor = 1.f; // global factor that scales all material densities for systematic studies
-
+  bool lowneut = false;
  O2ParamDef(SimCutParams, "SimCutParams");
 };


@@ -19,7 +19,8 @@ namespace conf

 namespace
 {
-static const std::string confstrings[3] = {"FTFP_BERT_EMV+optical", "FTFP_BERT_EMV+optical+biasing", "FTFP_INCLXX_EMV+optical"};
+static const std::string confstrings[4] = {"FTFP_BERT_EMV+optical", "FTFP_BERT_EMV+optical+biasing", "FTFP_INCLXX_EMV+optical",
+                                           "FTFP_BERT_HP_EMV+optical"};
 }

 std::string const& G4Params::getPhysicsConfigString() const

@@ -131,6 +131,9 @@ class MaterialManager
  /// Custom setting of process or cut given parameter name and value
  void SpecialCut(const char* modname, int localindex, ECut parID, Float_t val);

+  /// Set flag fro low energy neutron transport
+  void SetLowEnergyNeutronTransport(bool flag) { mLowNeut = flag; }
+
  /// load cuts and process flags from a data file (like AliRoot did)
  void loadCutsAndProcessesFromFile(const char* modname, const char* filename);
  void loadCutsAndProcessesFromJSON(ESpecial special = ESpecial::kFALSE, std::string const& filename = "");
@@ -264,6 +267,8 @@ class MaterialManager
  /// Decide whether special process and cut settings should be applied
  bool mApplySpecialProcesses = true;
  bool mApplySpecialCuts = true;
+  /// Flag for low energy neutron transport
+  bool mLowNeut = false;

 public:
  ClassDefNV(MaterialManager, 0);

@@ -254,6 +254,12 @@ void MaterialManager::Cut(ESpecial special, int globalindex, ECut cut, Float_t v
  if (val < 0.) {
    return;
  }
+  // if low energy neutron transport is requested setting kCUTNEU will set to 0.005eV
+  if (mLowNeut && cut == ECut::kCUTNEU) {
+    LOG(info) << "Due to low energy neutrons, neutron cut value " << val << " discarded and reset to 5e-12";
+    val = 5.e-12;
+  }
+
  auto it = mCutIDToName.find(cut);
  if (it == mCutIDToName.end()) {
    return;

@@ -13,6 +13,7 @@ TVirtualMC* FlukaRuntimeConfig()
  auto inpFile = params.scoringFile;
  auto userStepping = params.userStepping;
  auto hadronCut = params.activationHadronCut;
+  fluka->SetLowEnergyNeutronTransport(lowNeutron);
  if (!inpFile.empty()) {
    fluka->SetActivationSimulation(isAct, hadronCut);
    fluka->SetUserScoringFileName(inpFile.c_str());

@@ -133,6 +133,7 @@ void O2MCApplicationBase::InitGeometry()
  // load special cuts which might be given from the outside first.
  auto& matMgr = o2::base::MaterialManager::Instance();
  matMgr.loadCutsAndProcessesFromJSON(o2::base::MaterialManager::ESpecial::kTRUE);
+  matMgr.SetLowEnergyNeutronTransport(mCutParams.lowneut);
  // During the following, FairModule::SetSpecialPhysicsCuts will be called for each module
  FairMCApplication::InitGeometry();
  matMgr.writeCutsAndProcessesToJSON();