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 // CPA100 elastic model class for electrons 27 // 28 // Based on the work of M. Terrissol and M. C. Bordage 29 // 30 // Users are requested to cite the following papers: 31 // - M. Terrissol, A. Baudre, Radiat. Prot. Dosim. 31 (1990) 175-177 32 // - M.C. Bordage, J. Bordes, S. Edel, M. Terrissol, X. Franceries, 33 // M. Bardies, N. Lampe, S. Incerti, Phys. Med. 32 (2016) 1833-1840 34 // 35 // Authors of this class: 36 // M.C. Bordage, M. Terrissol, S. Edel, J. Bordes, S. Incerti 37 // 38 // 15.01.2014: creation 39 // Based on the study by S. Zein et. al. Nucl. Inst. Meth. B 488 (2021) 70-82 40 // 1/2/2023 : Hoang added modification 41 42 #ifndef G4DNACPA100ElasticModel_h 43 #define G4DNACPA100ElasticModel_h 1 44 45 #include "G4DNACrossSectionDataSet.hh" 46 #include "G4Electron.hh" 47 #include "G4LogLogInterpolation.hh" 48 #include "G4NistManager.hh" 49 #include "G4ParticleChangeForGamma.hh" 50 #include "G4ProductionCutsTable.hh" 51 #include "G4VDNAModel.hh" 52 53 #include <map> 54 55 class G4DNACPA100ElasticModel : public G4VDNAModel 56 { 57 using TriDimensionMap = 58 std::map<std::size_t, std::map<const G4ParticleDefinition*, 59 std::map<G4double, std::map<G4double, G4double>>>>; 60 using VecMap = 61 std::map<std::size_t, 62 std::map<const G4ParticleDefinition*, std::map<G4double, std::vector<G4double>>>>; 63 64 public: 65 explicit G4DNACPA100ElasticModel(const G4ParticleDefinition* p = nullptr, 66 const G4String& nam = "DNACPA100ElasticModel"); 67 68 ~G4DNACPA100ElasticModel() override = default; 69 70 void Initialise(const G4ParticleDefinition*, const G4DataVector&) override; 71 72 G4double CrossSectionPerVolume(const G4Material* material, const G4ParticleDefinition* p, 73 G4double ekin, G4double emin, G4double emax) override; 74 75 void SampleSecondaries(std::vector<G4DynamicParticle*>*, const G4MaterialCutsCouple*, 76 const G4DynamicParticle*, G4double tmin, G4double maxEnergy) override; 77 78 inline void SelectStationary(G4bool input); 79 80 G4DNACPA100ElasticModel& operator=(const G4DNACPA100ElasticModel& right) = delete; 81 82 G4DNACPA100ElasticModel(const G4DNACPA100ElasticModel&) = delete; 83 84 inline G4double GetElasticLevel(const std::size_t& l) 85 { 86 return fLevels[l]; 87 } 88 89 private: 90 G4ParticleChangeForGamma* fParticleChangeForGamma = nullptr; 91 G4bool statCode = false; 92 G4bool isInitialised = false; 93 G4int verboseLevel = 0; 94 95 G4double Theta(const G4ParticleDefinition* p, G4double k, G4double integrDiff, 96 const std::size_t& materialID); 97 98 G4double LinLinInterpolate(G4double e1, G4double e2, G4double e, G4double xs1, G4double xs2); 99 100 G4double LinLogInterpolate(G4double e1, G4double e2, G4double e, G4double xs1, G4double xs2); 101 102 G4double LogLogInterpolate(G4double e1, G4double e2, G4double e, G4double xs1, G4double xs2); 103 104 G4double QuadInterpolator(G4double e11, G4double e12, G4double e21, G4double e22, G4double x11, 105 G4double x12, G4double x21, G4double x22, G4double t1, G4double t2, 106 G4double t, G4double e); 107 108 G4double fKillBelowEnergy = 0.; 109 TriDimensionMap diffCrossSectionData; 110 VecMap eValuesVect; 111 std::map<std::size_t, std::map<const G4ParticleDefinition*, std::vector<G4double>>> tValuesVec; 112 113 G4double RandomizeCosTheta(G4double k, const std::size_t& materialID); 114 115 void ReadDiffCSFile(const std::size_t& id, const G4ParticleDefinition* p, const G4String& file, 116 const G4double&) override; 117 118 const G4Material* fpGuanine = nullptr; 119 const G4Material* fpG4_WATER = nullptr; 120 const G4Material* fpDeoxyribose = nullptr; 121 const G4Material* fpCytosine = nullptr; 122 const G4Material* fpThymine = nullptr; 123 const G4Material* fpAdenine = nullptr; 124 const G4Material* fpPhosphate = nullptr; 125 const G4ParticleDefinition* fpParticle = nullptr; 126 127 G4DNACPA100ElasticModel* fpModelData = nullptr; 128 129 std::map<std::size_t /*MatIndex*/, G4double> fLevels; 130 }; 131 132 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 133 134 inline void G4DNACPA100ElasticModel::SelectStationary(G4bool input) 135 { 136 statCode = input; 137 } 138 139 #endif 140