Back to home page

EIC code displayed by LXR

 
 

    


Warning, file /geant4/examples/extended/hadronic/FlukaCern/ProcessLevel/FinalState/include/HadronicGenerator.hh was not indexed or was modified since last indexation (in which case cross-reference links may be missing, inaccurate or erroneous).

0001 //
0002 // ********************************************************************
0003 // * License and Disclaimer                                           *
0004 // *                                                                  *
0005 // * The  Geant4 software  is  copyright of the Copyright Holders  of *
0006 // * the Geant4 Collaboration.  It is provided  under  the terms  and *
0007 // * conditions of the Geant4 Software License,  included in the file *
0008 // * LICENSE and available at  http://cern.ch/geant4/license .  These *
0009 // * include a list of copyright holders.                             *
0010 // *                                                                  *
0011 // * Neither the authors of this software system, nor their employing *
0012 // * institutes,nor the agencies providing financial support for this *
0013 // * work  make  any representation or  warranty, express or implied, *
0014 // * regarding  this  software system or assume any liability for its *
0015 // * use.  Please see the license in the file  LICENSE  and URL above *
0016 // * for the full disclaimer and the limitation of liability.         *
0017 // *                                                                  *
0018 // * This  code  implementation is the result of  the  scientific and *
0019 // * technical work of the GEANT4 collaboration.                      *
0020 // * By using,  copying,  modifying or  distributing the software (or *
0021 // * any work based  on the software)  you  agree  to acknowledge its *
0022 // * use  in  resulting  scientific  publications,  and indicate your *
0023 // * acceptance of all terms of the Geant4 Software license.          *
0024 // ********************************************************************
0025 //
0026 /// \file HadronicGenerator.hh
0027 /// \brief Definition of the HadronicGenerator class
0028 //
0029 //------------------------------------------------------------------------
0030 // Class: HadronicGenerator
0031 // Author: Alberto Ribon (CERN EP/SFT), May 2020
0032 // Modified: G. Hugo, 8 December 2022
0033 //
0034 /// This class shows how to use Geant4 as a generator for simulating
0035 /// inelastic hadron-nuclear interactions.
0036 /// Some of the most used hadronic models are currently supported in
0037 /// this class:
0038 /// - the hadronic string models Fritiof (FTF) and Quark-Gluon-String (QGS)
0039 ///   coupled with Precompound/de-excitation
0040 /// - the intranuclear cascade models: Bertini (BERT), Binary Cascade (BIC),
0041 ///                                    and Liege (INCL)
0042 /// Combinations of two models - in a transition energy interval, with a
0043 /// linear probability as a function of the energy - are also available to
0044 /// "mimic" the transition between hadronic models as in the most common
0045 /// Geant4 reference physics lists.
0046 ///
0047 /// The current version of this class does NOT support:
0048 /// -  hadron elastic interactions
0049 /// -  neutron capture and fission
0050 /// -  precise low-energy inelastic interactions of neutrons and
0051 ///    charged particles (i.e. ParticleHP)
0052 /// -  gamma/lepton-nuclear inelastic interactions
0053 ///
0054 /// This class does NOT use the Geant4 run-manager, and therefore should
0055 /// be usable in a multi-threaded application, with one instance of this
0056 /// class in each thread.
0057 ///
0058 /// This class has been inspired by test30 (whose author is Vladimir
0059 /// Ivanchenko), with various simplifications and restricted to hadronic
0060 /// inelastic interactions.
0061 //------------------------------------------------------------------------
0062 
0063 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0064 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0065 
0066 #ifndef HadronicGenerator_h
0067 #define HadronicGenerator_h 1
0068 
0069 #include "G4HadronicProcess.hh"
0070 #include "G4ThreeVector.hh"
0071 #include "G4ios.hh"
0072 #include "globals.hh"
0073 
0074 #include <iomanip>
0075 #include <map>
0076 
0077 class G4ParticleDefinition;
0078 class G4VParticleChange;
0079 class G4ParticleTable;
0080 class G4Material;
0081 class G4HadronicInteraction;
0082 
0083 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0084 
0085 class HadronicGenerator
0086 {
0087     // This class provides the functionality of a "hadronic generator"
0088     // for Geant4 final-state inelastic hadronic collisions.
0089     // Only a few of the available Geant4 final-state hadronic inelastic
0090     // "physics cases" are currently available in this class - but it can
0091     // be extended to other cases if needed.
0092     // It is important to notice that this class does NOT use the Geant4
0093     // run-manager, so it should work fine in a multi-threaded environment,
0094     // with a separate instance of this class in each thread.
0095   public:
0096     explicit HadronicGenerator(const G4String physicsCase = "FTFP_BERT_ATL");
0097     // Currently supported final-state hadronic inelastic "physics cases":
0098     // -   Hadronic models :       CernFLUKAHadronInelastic,
0099     //                             BERT, BIC, IonBIC, INCL, FTFP, QGSP
0100     // -  "Physics-list proxies" : FTFP_BERT_ATL (default), FTFP_BERT,
0101     //                             QGSP_BERT, QGSP_BIC, FTFP_INCLXX
0102     //    (i.e. they are not real, complete physics lists - for instance
0103     //     they do not have: transportation, electromagnetic physics,
0104     //     hadron elastic scattering, neutron fission and capture, etc. -
0105     //     however, they cover all hadron types and all energies by
0106     //     combining different hadronic models, i.e. there are transitions
0107     //     between two hadronic models in well-defined energy intervals,
0108     //     e.g. "FTFP_BERT" has the transition between BERT and FTFP
0109     //     hadronic models; moreover, the transition intervals used in
0110     //     our "physics cases"might not be the same as in the corresponding
0111     //     physics lists).
0112 
0113     ~HadronicGenerator();
0114 
0115     inline G4bool IsPhysicsCaseSupported() const;
0116     // Returns "true" if the physicsCase is supported; "false" otherwise.
0117 
0118     G4bool IsApplicable(const G4String& nameProjectile, const G4double projectileEnergy) const;
0119     G4bool IsApplicable(G4ParticleDefinition* projectileDefinition,
0120                         const G4double projectileEnergy) const;
0121     // Returns "true" if the specified projectile (either by name or particle definition)
0122     // of given energy is applicable, "false" otherwise.
0123 
0124     G4VParticleChange* GenerateInteraction(const G4String& nameProjectile,
0125                                            const G4double projectileEnergy,
0126                                            const G4ThreeVector& projectileDirection,
0127                                            G4Material* targetMaterial);
0128     G4VParticleChange* GenerateInteraction(G4ParticleDefinition* projectileDefinition,
0129                                            const G4double projectileEnergy,
0130                                            const G4ThreeVector& projectileDirection,
0131                                            G4Material* targetMaterial);
0132     // This is the main method provided by the class:
0133     // in input it receives the projectile (either by name or particle definition),
0134     // its energy, its direction and the target material, and it returns one sampled
0135     // final-state of the inelastic hadron-nuclear collision as modelled by the
0136     // final-state hadronic inelastic "physics case" specified in the constructor.
0137     // If the required hadronic collision is not possible, then the method returns
0138     // immediately an empty "G4VParticleChange", i.e. without secondaries produced.
0139 
0140     inline G4HadronicProcess* GetHadronicProcess() const;
0141     inline G4HadronicInteraction* GetHadronicInteraction() const;
0142     // Returns the hadronic process and the hadronic interaction, respectively,
0143     // that handled the last call of "GenerateInteraction".
0144 
0145     G4double GetImpactParameter() const;
0146     G4int GetNumberOfTargetSpectatorNucleons() const;
0147     G4int GetNumberOfProjectileSpectatorNucleons() const;
0148     G4int GetNumberOfNNcollisions() const;
0149     // In the case of hadronic interactions handled by the FTF model, returns,
0150     // respectively, the impact parameter, the number of target/projectile
0151     // spectator nucleons, and the number of nucleon-nucleon collisions,
0152     // else, returns a negative value (-999).
0153 
0154   private:
0155     G4String fPhysicsCase;
0156     G4bool fPhysicsCaseIsSupported = false;
0157     G4HadronicProcess* fLastHadronicProcess = nullptr;
0158     G4ParticleTable* fPartTable = nullptr;
0159     std::map<G4ParticleDefinition*, G4HadronicProcess*> fProcessMap;
0160 };
0161 
0162 inline G4bool HadronicGenerator::IsPhysicsCaseSupported() const
0163 {
0164   return fPhysicsCaseIsSupported;
0165 }
0166 
0167 inline G4HadronicProcess* HadronicGenerator::GetHadronicProcess() const
0168 {
0169   return fLastHadronicProcess;
0170 }
0171 
0172 inline G4HadronicInteraction* HadronicGenerator::GetHadronicInteraction() const
0173 {
0174   return fLastHadronicProcess == nullptr ? nullptr : fLastHadronicProcess->GetHadronicInteraction();
0175 }
0176 
0177 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0178 
0179 #endif