Geant4 Cross Reference |
1 // 2 // ******************************************************************** 3 // * License and Disclaimer * 4 // * * 5 // * The Geant4 software is copyright of the Copyright Holders of * 6 // * the Geant4 Collaboration. It is provided under the terms and * 7 // * conditions of the Geant4 Software License, included in the file * 8 // * LICENSE and available at http://cern.ch/geant4/license . These * 9 // * include a list of copyright holders. * 10 // * * 11 // * Neither the authors of this software system, nor their employing * 12 // * institutes,nor the agencies providing financial support for this * 13 // * work make any representation or warranty, express or implied, * 14 // * regarding this software system or assume any liability for its * 15 // * use. Please see the license in the file LICENSE and URL above * 16 // * for the full disclaimer and the limitation of liability. * 17 // * * 18 // * This code implementation is the result of the scientific and * 19 // * technical work of the GEANT4 collaboration. * 20 // * By using, copying, modifying or distributing the software (or * 21 // * any work based on the software) you agree to acknowledge its * 22 // * use in resulting scientific publications, and indicate your * 23 // * acceptance of all terms of the Geant4 Software license. * 24 // ******************************************************************** 25 // 26 // 27 // 080721 Create adjust_final_state method by T. Koi 28 // 080801 Introduce theNDLDataA,Z which has A and Z of NDL data by T. Koi 29 // 30 // P. Arce, June-2014 Conversion neutron_hp to particle_hp 31 // V. Ivanchenko, July-2023 Basic revision of particle HP classes 32 // 33 34 #ifndef G4ParticleHPFinalState_h 35 #define G4ParticleHPFinalState_h 36 37 #include "G4Cache.hh" 38 #include "G4HadFinalState.hh" 39 #include "G4HadProjectile.hh" 40 #include "G4Material.hh" 41 #include "G4Neutron.hh" 42 #include "G4IonTable.hh" 43 #include "G4ParticleHPManager.hh" 44 #include "G4ParticleHPNames.hh" 45 #include "G4ParticleHPVector.hh" 46 47 class G4ParticleDefinition; 48 49 class G4ParticleHPFinalState 50 { 51 public: 52 53 G4ParticleHPFinalState(); 54 virtual ~G4ParticleHPFinalState(); 55 56 inline void Init(G4double A, G4double Z, const G4String& dirName, 57 const G4String& aFSType, G4ParticleDefinition* p) 58 { 59 theProjectile = p; 60 Init(A, Z, 0, dirName, aFSType, p); 61 } 62 virtual void Init(G4double A, G4double Z, G4int M, const G4String& dirName, 63 const G4String& aFSType, G4ParticleDefinition*) = 0; 64 65 virtual G4HadFinalState* ApplyYourself(const G4HadProjectile&) 66 { 67 throw G4HadronicException( 68 __FILE__, __LINE__, 69 "G4ParticleHPFinalState::ApplyYourself(..) needs implementation"); 70 return nullptr; 71 } 72 73 // this would better be done templating G4ParticleHPChannel..., 74 // 75 virtual G4ParticleHPFinalState* New() = 0; 76 77 inline G4bool HasXsec() const { return hasXsec; } 78 inline G4bool HasFSData() const { return hasFSData; } 79 inline G4bool HasAnyData() const { return hasAnyData; } 80 81 virtual G4double GetXsec(G4double) const { return 0.; } 82 virtual G4ParticleHPVector* GetXsec() const { return nullptr; } 83 84 void SetA_Z(G4double anA, G4double aZ, G4int aM = 0) 85 { 86 theBaseA = G4lrint(anA); 87 theBaseZ = G4lrint(aZ); 88 theBaseM = aM; 89 } 90 G4double GetZ() const { return theBaseZ; } 91 G4double GetN() const { return theBaseA; } 92 G4double GetA() const { return theBaseA; } 93 G4int GetM() const { return theBaseM; } 94 95 inline void SetAZMs(const G4ParticleHPDataUsed& used) 96 { 97 theNDLDataA = G4lrint(used.GetA()); 98 theNDLDataZ = G4lrint(used.GetZ()); 99 theNDLDataM = used.GetM(); 100 } 101 102 inline void SetAZMs(G4double anA, G4double aZ, G4int aM, const G4ParticleHPDataUsed& used) 103 { 104 theBaseA = G4lrint(anA); 105 theBaseZ = G4lrint(aZ); 106 theBaseM = aM; 107 theNDLDataA = G4lrint(used.GetA()); 108 theNDLDataZ = G4lrint(used.GetZ()); 109 theNDLDataM = used.GetM(); 110 } 111 112 inline void SetProjectile(G4ParticleDefinition* projectile) 113 { 114 theProjectile = projectile; 115 } 116 117 G4ParticleHPFinalState& operator=(const G4ParticleHPFinalState& right) = delete; 118 G4ParticleHPFinalState(const G4ParticleHPFinalState&) = delete; 119 120 protected: 121 122 void adjust_final_state(G4LorentzVector); 123 124 G4ParticleDefinition* theProjectile{nullptr}; 125 G4ParticleHPManager* fManager; 126 G4IonTable* ionTable; 127 128 G4int theBaseA{0}; 129 G4int theBaseZ{0}; 130 G4int theBaseM{0}; 131 G4int theNDLDataZ{0}; 132 G4int theNDLDataA{0}; 133 G4int theNDLDataM{0}; 134 135 G4int secID{-1}; 136 // Creator model ID for the secondaries created by this class or derived ones 137 138 G4bool hasXsec{true}; 139 G4bool hasFSData{true}; 140 G4bool hasAnyData{true}; 141 G4ParticleHPNames theNames; 142 143 G4Cache<G4HadFinalState*> theResult; 144 145 }; 146 147 #endif 148