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 // History: 27 // ----------- 28 // 10 Nov 2021 S. Guatelli & S. Bakr, Wollongong University - 1st implementation 29 // 30 // Class description 31 // ---------------- 32 // Computation of K, L & M shell ECPSSR ionisation cross sections for protons and alphas 33 // Based on the work of 34 // - S. Bakr et al. (2021) NIM B, 507:11-19. 35 // - S. Bakr et al (2018), NIMB B, 436: 285-291. 36 // --------------------------------------------------------------------------------------- 37 38 #include <fstream> 39 #include <iomanip> 40 41 #include "G4ANSTOecpssrLixsModel.hh" 42 43 #include "globals.hh" 44 #include "G4SystemOfUnits.hh" 45 #include "G4ios.hh" 46 47 #include "G4EMDataSet.hh" 48 #include "G4LinInterpolation.hh" 49 #include "G4Proton.hh" 50 #include "G4Alpha.hh" 51 52 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 53 54 G4ANSTOecpssrLixsModel::G4ANSTOecpssrLixsModel() 55 { 56 G4cout << "Using ANSTO L Cross Sections! "<< G4endl; 57 58 interpolation = new G4LinInterpolation(); 59 60 for (G4int i=26; i<93; i++) 61 { 62 protonL1DataSetMap[i] = new G4EMDataSet(i,interpolation); 63 protonL1DataSetMap[i]->LoadData("pixe_ANSTO/proton/l1-"); 64 65 protonL2DataSetMap[i] = new G4EMDataSet(i,interpolation); 66 protonL2DataSetMap[i]->LoadData("pixe_ANSTO/proton/l2-"); 67 68 protonL3DataSetMap[i] = new G4EMDataSet(i,interpolation); 69 protonL3DataSetMap[i]->LoadData("pixe_ANSTO/proton/l3-"); 70 } 71 72 for (G4int i=26; i<93; i++) 73 { 74 alphaL1DataSetMap[i] = new G4EMDataSet(i,interpolation); 75 alphaL1DataSetMap[i]->LoadData("pixe_ANSTO/alpha/l1-"); 76 77 alphaL2DataSetMap[i] = new G4EMDataSet(i,interpolation); 78 alphaL2DataSetMap[i]->LoadData("pixe_ANSTO/alpha/l2-"); 79 80 alphaL3DataSetMap[i] = new G4EMDataSet(i,interpolation); 81 alphaL3DataSetMap[i]->LoadData("pixe_ANSTO/alpha/l3-"); 82 } 83 84 } 85 86 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 87 88 G4ANSTOecpssrLixsModel::~G4ANSTOecpssrLixsModel() 89 { 90 protonL1DataSetMap.clear(); 91 alphaL1DataSetMap.clear(); 92 93 protonL2DataSetMap.clear(); 94 alphaL2DataSetMap.clear(); 95 96 protonL3DataSetMap.clear(); 97 alphaL3DataSetMap.clear(); 98 99 delete interpolation; 100 } 101 102 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 103 104 G4double G4ANSTOecpssrLixsModel::CalculateL1CrossSection(G4int zTarget,G4double massIncident, G4double energyIncident) 105 { 106 G4Proton* aProton = G4Proton::Proton(); 107 G4Alpha* aAlpha = G4Alpha::Alpha(); 108 G4double sigma = 0; 109 110 if (massIncident == aProton->GetPDGMass()) 111 { 112 113 if (energyIncident > 0.2*MeV && energyIncident < 5.*MeV && zTarget < 93 && zTarget > 25) { 114 115 sigma = protonL1DataSetMap[zTarget]->FindValue(energyIncident/MeV); 116 if (sigma !=0 && energyIncident > protonL1DataSetMap[zTarget]->GetEnergies(0).back()*MeV) return 0.; 117 } 118 119 else if (massIncident == aAlpha->GetPDGMass()) 120 { 121 122 if (energyIncident > 0.2*MeV && energyIncident < 40.*MeV && zTarget < 93 && zTarget > 25) { 123 124 sigma = alphaL1DataSetMap[zTarget]->FindValue(energyIncident/MeV); 125 if (sigma !=0 && energyIncident > alphaL1DataSetMap[zTarget]->GetEnergies(0).back()*MeV) return 0.; 126 } 127 else 128 { 129 sigma = 0.; 130 } 131 } 132 } 133 // sigma is in internal units: it has been converted from 134 // the input file in barns bt the EmDataset 135 return sigma; 136 } 137 138 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 139 140 G4double G4ANSTOecpssrLixsModel::CalculateL2CrossSection(G4int zTarget,G4double massIncident, G4double energyIncident) 141 { 142 G4Proton* aProton = G4Proton::Proton(); 143 G4Alpha* aAlpha = G4Alpha::Alpha(); 144 G4double sigma = 0; 145 146 if (massIncident == aProton->GetPDGMass()) 147 { 148 if (energyIncident > 0.2*MeV && energyIncident < 5.*MeV && zTarget < 93 && zTarget > 25) { 149 150 sigma = protonL2DataSetMap[zTarget]->FindValue(energyIncident/MeV); 151 if (sigma !=0 && energyIncident > protonL2DataSetMap[zTarget]->GetEnergies(0).back()*MeV) return 0.; 152 } 153 154 else if (massIncident == aAlpha->GetPDGMass()) 155 { 156 if (energyIncident > 0.2*MeV && energyIncident < 40.*MeV && zTarget < 93 && zTarget > 25) { 157 158 sigma = alphaL2DataSetMap[zTarget]->FindValue(energyIncident/MeV); 159 if (sigma !=0 && energyIncident > alphaL2DataSetMap[zTarget]->GetEnergies(0).back()*MeV) return 0.; 160 } 161 else 162 { 163 sigma = 0.; 164 } 165 } 166 } 167 // sigma is in internal units: it has been converted from 168 // the input file in barns bt the EmDataset 169 return sigma; 170 } 171 172 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 173 174 G4double G4ANSTOecpssrLixsModel::CalculateL3CrossSection(G4int zTarget,G4double massIncident, G4double energyIncident) 175 { 176 G4Proton* aProton = G4Proton::Proton(); 177 G4Alpha* aAlpha = G4Alpha::Alpha(); 178 G4double sigma = 0; 179 180 if (massIncident == aProton->GetPDGMass()) 181 { 182 if (energyIncident > 0.2*MeV && energyIncident < 5.*MeV && zTarget < 93 && zTarget > 25) { 183 184 sigma = protonL3DataSetMap[zTarget]->FindValue(energyIncident/MeV); 185 if (sigma !=0 && energyIncident > protonL3DataSetMap[zTarget]->GetEnergies(0).back()*MeV) return 0.; 186 } 187 } 188 189 else if (massIncident == aAlpha->GetPDGMass()) 190 { 191 if (energyIncident > 0.2*MeV && energyIncident < 40.*MeV && zTarget < 93 && zTarget > 25) { 192 193 sigma = alphaL3DataSetMap[zTarget]->FindValue(energyIncident/MeV); 194 if (sigma !=0 && energyIncident > alphaL3DataSetMap[zTarget]->GetEnergies(0).back()*MeV) return 0.; 195 } 196 } 197 198 else 199 { 200 sigma = 0.; 201 } 202 203 204 // sigma is in internal units: it has been converted from 205 // the input file in barns bt the EmDataset 206 return sigma; 207 } 208