ConfigurablePhysicsList.icc

Go to the documentation of this file.

// Configurable Physics List Class
//
// Ben Jones, MIT, 24/06/09
//
// Based on the QGSP_BERT physics list supplied with geant 4, but with
// options to switch on and off different physics processed from the config
// file.
//
// To include a new physics process, create a new physics builder registering
// the relevant particles and interactions.  Then add the builder to the
// GetPhysicsBuilders and GetDefaultSettings functions in this class,
// using the name by which the builder will be referenced in the config
// file.  Physics builders to be enabled are specified in the LArG4 config.
//

#include "art/Framework/Services/Registry/ServiceHandle.h"

#include "Geant4/G4ChargeExchangePhysics.hh"
#include "Geant4/G4EmExtraPhysics.hh"
#include "Geant4/G4ParticleDefinition.hh"
#include "Geant4/G4ParticleTable.hh"
#include "Geant4/G4ParticleTypes.hh"
#include "Geant4/G4ParticleWithCuts.hh"
#include "Geant4/G4ProcessManager.hh"
#include "Geant4/G4ProcessVector.hh"
#include "Geant4/globals.hh"

#include "Geant4/G4Material.hh"
#include "Geant4/G4MaterialTable.hh"
#include "Geant4/G4MuonNuclearProcess.hh"
#include "Geant4/G4ios.hh"
#include <iomanip>

#include "Geant4/G4Version.hh"
#if G4VERSION_NUMBER < 1060
#include "Geant4/G4DataQuestionaire.hh"
#else
// The G4DataQuestionaire.hh interface was removed in the 10.6, as described in the release notes:
// http://geant4-data.web.cern.ch/geant4-data/ReleaseNotes/ReleaseNotes4.10.6.html 8
// The class was header-only and only used internally by packaged physics lists to check dataset
// environment variables. These are checked anyway by processes that require them,
// so it should be possible to remove it from your application without an issue.
#endif

#include "larsim/LegacyLArG4/CustomPhysicsTable.hh"
#include "larsim/LegacyLArG4/MuNuclearSplittingProcess.h"
#include "larsim/LegacyLArG4/MuNuclearSplittingProcessXSecBias.h"
#include "larsim/Simulation/LArG4Parameters.h"

#include "cetlib_except/exception.h"
#include "messagefacility/MessageLogger/MessageLogger.h"

namespace larg4 {

  template <class T>
  TConfigurablePhysicsList<T>::TConfigurablePhysicsList(G4int const ver) : T()
  {
    mf::LogWarning logmsg{"ConfigurablePhysics"};
    logmsg << "Setting up Configurable Physics List.\n";
    // note - not sure exact purpose of these commands, but left in from QGSP_BERT list
    this->defaultCutValue = 0.7 * CLHEP::mm;
    this->SetVerboseLevel(ver);
#if G4VERSION_NUMBER < 1060
    G4DataQuestionaire it(photon);
#endif
    art::ServiceHandle<sim::LArG4Parameters> lgp;
    if (lgp->UseCustomPhysics()) {
      EnabledPhysics = lgp->EnabledPhysics();
      logmsg << "Custom physics list enabled, contents:\n";
      for (auto const& val : EnabledPhysics) {
        logmsg << "  " << val << '\n';
        assert(TheCustomPhysicsTable);
        if (!TheCustomPhysicsTable->IsPhysicsAvailable(val)) {
          logmsg << "\nERROR: physics not available: " << val << '\n';
          std::vector<std::string> availablePhysics(
            TheCustomPhysicsTable->GetAvailablePhysicsList());
          logmsg << "\n-- Available physics:\n";
          for (auto const& nm : availablePhysics) {
            if (TheCustomPhysicsTable->IsPhysicsAvailable(nm)) { logmsg << "  " << nm << '\n'; }
          }
          logmsg << "--\n";
          logmsg << "\nERROR: Throwing exception!\n";
          throw cet::exception("ConfigurablePhysics") << "Physics not available: " << val;
        }
      }
    }
    else {
      logmsg << "Custom physics list disabled, using default QGSP_BERT configuration.\n";
      EnabledPhysics = GetDefaultSettings();
    }
    for (auto const& PhysicsName : EnabledPhysics) {
      logmsg << "Registering physics: " << PhysicsName << '\n';
      G4VPhysicsConstructor* g4v = TheCustomPhysicsTable->GetPhysicsConstructor(PhysicsName);
      if (PhysicsName == "SynchrotronAndGN") {
        logmsg << PhysicsName << ": Turning on GammaNuclear, Synchrotron.\n";
        G4String on("on");
        ((G4EmExtraPhysics*)g4v)->GammaNuclear(on);
        ((G4EmExtraPhysics*)g4v)->Synch(on);
        // We want MuonNuclear off, since we're gonna activate it
        // within a Wrapped Process. This would double count these evts.
        if (!lgp->K0Bias()) {
          logmsg << PhysicsName << ": Turning on MuNuclear.\n";
          ((G4EmExtraPhysics*)g4v)->MuonNuclear(on);
        }
      }
      this->RegisterPhysics(g4v);
    }
  }

  template <class T>
  TConfigurablePhysicsList<T>::~TConfigurablePhysicsList()
  {}

  template <class T>
  void TConfigurablePhysicsList<T>::SetCuts()
  {

    art::ServiceHandle<sim::LArG4Parameters> lg4p;

    if (this->verboseLevel > 1) G4cout << "ConfigurablePhysicsList::SetCuts:";
    this->SetCutsWithDefault();

    // Set Proton Cut to  0, particularly important when using High Precision Packages
    // D.R. 02/25/19
    bool ModifyProtonCut = lg4p->ModifyProtonCut();
    if (ModifyProtonCut) {
      double theProtonCut = (double)(lg4p->NewProtonCut());
      this->SetCutValue(theProtonCut, "proton");
      mf::LogInfo("ConfigurablePhysicsList::SetCuts:") << "Setting Proton cut to: " << theProtonCut;
    }
  }

  template <class T>
  std::vector<std::string> TConfigurablePhysicsList<T>::GetDefaultSettings()
  {
    // Set default enabled physics processes (equivalent to QGSP_BERT physics list)
    std::vector<std::string> TheVector;
    TheVector.push_back("Em");
    TheVector.push_back("SynchrotronAndGN");
    TheVector.push_back("Decay");
    TheVector.push_back("Hadron");
    TheVector.push_back("HadronElastic");
    TheVector.push_back("Stopping");
    TheVector.push_back("Ion");
    TheVector.push_back("NeutronTrackingCut");
    return TheVector;
  }

}

// Sept 2009 - Ben Jones, MIT