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 // Hadronic Process: Very Low Energy Neutron X-Sections 27 // original by H.P. Wellisch, TRIUMF, 14-Feb-97 28 // Builds and has the Cross-section data for one element and channel. 29 // 30 // Bug fixes and workarounds in the destructor, F.W.Jones 06-Jul-1999 31 // 070612 Fix memory leaking by T. Koi 32 // 33 // 080520 Delete unnecessary dependencies by T. Koi 34 35 // P. Arce, June-2014 Conversion neutron_hp to particle_hp 36 // V. Ivanchenko, July-2023 Basic revision of particle HP classes 37 // 38 39 #ifndef G4ParticleHPChannel_h 40 #define G4ParticleHPChannel_h 1 41 42 #include "G4Element.hh" 43 #include "G4HadProjectile.hh" 44 #include "G4Material.hh" 45 #include "G4ParticleHPCaptureFS.hh" 46 #include "G4ParticleHPFinalState.hh" 47 #include "G4ParticleHPIsoData.hh" 48 #include "G4ParticleHPManager.hh" 49 #include "G4ParticleHPVector.hh" 50 #include "G4StableIsotopes.hh" 51 #include "G4WendtFissionFragmentGenerator.hh" 52 #include "globals.hh" 53 54 class G4ParticleDefinition; 55 56 class G4ParticleHPChannel 57 { 58 public: 59 60 G4ParticleHPChannel(G4ParticleDefinition* projectile = nullptr); 61 62 ~G4ParticleHPChannel(); 63 64 G4double GetXsec(G4double energy) const; 65 66 G4double GetWeightedXsec(G4double energy, G4int isoNumber) const; 67 68 G4double GetFSCrossSection(G4double energy, G4int isoNumber) const; 69 70 G4bool IsActive(G4int isoNumber) const 71 { 72 return active[isoNumber]; 73 } 74 75 G4bool HasFSData(G4int isoNumber) const 76 { 77 return theFinalStates[isoNumber]->HasFSData(); 78 } 79 80 G4bool HasAnyData(G4int isoNumber) const 81 { 82 return theFinalStates[isoNumber]->HasAnyData(); 83 } 84 85 G4bool Register(G4ParticleHPFinalState* theFS); 86 87 void Init(G4Element* theElement, const G4String& dirName); 88 89 void Init(G4Element* theElement, const G4String& dirName, 90 const G4String& fsType); 91 92 void UpdateData(G4int A, G4int Z, G4int index, G4double abundance, 93 G4ParticleDefinition* projectile) 94 { 95 UpdateData(A, Z, 0, index, abundance, projectile); 96 } 97 98 void UpdateData(G4int A, G4int Z, G4int M, G4int index, G4double abundance, 99 G4ParticleDefinition* projectile); 100 101 void Harmonise(G4ParticleHPVector*& theStore, G4ParticleHPVector* theNew); 102 103 G4HadFinalState* ApplyYourself(const G4HadProjectile& theTrack, 104 G4int isoNumber = -1, 105 G4bool isElastic = false); 106 107 G4int GetNiso() const { return niso; } 108 109 G4double GetN(G4int i) const { return theFinalStates[i]->GetN(); } 110 G4double GetZ(G4int i) const { return theFinalStates[i]->GetZ(); } 111 G4double GetM(G4int i) const { return theFinalStates[i]->GetM(); } 112 113 G4bool HasDataInAnyFinalState() const 114 { 115 G4bool result = false; 116 for (G4int i = 0; i < niso; ++i) { 117 if (theFinalStates[i]->HasAnyData()) { 118 result = true; 119 break; 120 } 121 } 122 return result; 123 } 124 125 void DumpInfo() const; 126 127 const G4String& GetFSType() { return theFSType; } 128 129 G4ParticleHPFinalState** GetFinalStates() const { return theFinalStates; } 130 131 // method added by M.Zmeskal 02/2024 - to be used in G4ParticleHPFissionURR 132 G4WendtFissionFragmentGenerator* GetWendtFissionGenerator() const; 133 134 G4ParticleHPChannel(G4ParticleHPChannel &) = delete; 135 G4ParticleHPChannel & operator=(const G4ParticleHPChannel &right) = delete; 136 137 protected: 138 139 G4ParticleHPManager* fManager; 140 141 private: 142 143 G4ParticleDefinition* theProjectile; 144 145 G4ParticleHPVector* theChannelData; 146 G4Element* theElement{nullptr}; 147 148 // total (element) cross-section for this channel 149 G4ParticleHPVector* theBuffer{nullptr}; 150 151 G4ParticleHPIsoData* theIsotopeWiseData{nullptr}; 152 // these are the isotope-wise cross-sections for each final state. 153 G4ParticleHPFinalState** theFinalStates{nullptr}; 154 // also these are isotope-wise pionters, parallel to the above. 155 156 G4WendtFissionFragmentGenerator* wendtFissionGenerator{nullptr}; 157 G4bool* active{nullptr}; 158 G4int niso{-1}; 159 G4int registerCount{-1}; 160 161 G4String theDir{""}; 162 G4String theFSType{""}; 163 }; 164 165 #endif 166