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 // 01 Oct 2011 A.M., S.I. - 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 A. Taborda et al. 34 // X-Ray Spectrom. 2011, 40, 127-134 35 // --------------------------------------------------------------------------------------- 36 37 #include <fstream> 38 #include <iomanip> 39 40 #include "G4ecpssrFormFactorLixsModel.hh" 41 #include "globals.hh" 42 #include "G4SystemOfUnits.hh" 43 #include "G4ios.hh" 44 #include "G4EMDataSet.hh" 45 #include "G4LinInterpolation.hh" 46 #include "G4Proton.hh" 47 #include "G4Alpha.hh" 48 49 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 50 51 G4ecpssrFormFactorLixsModel::G4ecpssrFormFactorLixsModel() 52 { 53 interpolation = new G4LinInterpolation(); 54 55 for (G4int i=11; i<93; i++) 56 { 57 protonL1DataSetMap[i] = new G4EMDataSet(i,interpolation); 58 protonL1DataSetMap[i]->LoadData("pixe/ecpssr/proton/l1-i01m001c01-"); 59 60 protonL2DataSetMap[i] = new G4EMDataSet(i,interpolation); 61 protonL2DataSetMap[i]->LoadData("pixe/ecpssr/proton/l2-i01m001c01-"); 62 63 protonL3DataSetMap[i] = new G4EMDataSet(i,interpolation); 64 protonL3DataSetMap[i]->LoadData("pixe/ecpssr/proton/l3-i01m001c01-"); 65 } 66 67 for (G4int i=11; i<93; i++) 68 { 69 alphaL1DataSetMap[i] = new G4EMDataSet(i,interpolation); 70 alphaL1DataSetMap[i]->LoadData("pixe/ecpssr/alpha/l1-i02m004c02-"); 71 72 alphaL2DataSetMap[i] = new G4EMDataSet(i,interpolation); 73 alphaL2DataSetMap[i]->LoadData("pixe/ecpssr/alpha/l2-i02m004c02-"); 74 75 alphaL3DataSetMap[i] = new G4EMDataSet(i,interpolation); 76 alphaL3DataSetMap[i]->LoadData("pixe/ecpssr/alpha/l3-i02m004c02-"); 77 } 78 79 } 80 81 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 82 83 G4ecpssrFormFactorLixsModel::~G4ecpssrFormFactorLixsModel() 84 { 85 protonL1DataSetMap.clear(); 86 alphaL1DataSetMap.clear(); 87 88 protonL2DataSetMap.clear(); 89 alphaL2DataSetMap.clear(); 90 91 protonL3DataSetMap.clear(); 92 alphaL3DataSetMap.clear(); 93 94 delete interpolation; 95 } 96 97 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 98 99 G4double G4ecpssrFormFactorLixsModel::CalculateL1CrossSection(G4int zTarget,G4double massIncident, G4double energyIncident) 100 { 101 G4Proton* aProton = G4Proton::Proton(); 102 G4Alpha* aAlpha = G4Alpha::Alpha(); 103 G4double sigma = 0; 104 105 if (energyIncident > 0.1*MeV && energyIncident < 100.*MeV && zTarget < 93 && zTarget > 10) { 106 107 if (massIncident == aProton->GetPDGMass()) 108 { 109 sigma = protonL1DataSetMap[zTarget]->FindValue(energyIncident/MeV); 110 if (sigma !=0 && energyIncident > protonL1DataSetMap[zTarget]->GetEnergies(0).back()*MeV) return 0.; 111 } 112 else if (massIncident == aAlpha->GetPDGMass()) 113 { 114 sigma = alphaL1DataSetMap[zTarget]->FindValue(energyIncident/MeV); 115 if (sigma !=0 && energyIncident > alphaL1DataSetMap[zTarget]->GetEnergies(0).back()*MeV) return 0.; 116 } 117 else 118 { 119 sigma = 0.; 120 } 121 } 122 123 // sigma is in internal units: it has been converted from 124 // the input file in barns bt the EmDataset 125 return sigma; 126 } 127 128 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 129 130 G4double G4ecpssrFormFactorLixsModel::CalculateL2CrossSection(G4int zTarget,G4double massIncident, G4double energyIncident) 131 { 132 G4Proton* aProton = G4Proton::Proton(); 133 G4Alpha* aAlpha = G4Alpha::Alpha(); 134 G4double sigma = 0; 135 136 if (energyIncident > 0.1*MeV && energyIncident < 100.*MeV && zTarget < 93 && zTarget > 10) { 137 138 if (massIncident == aProton->GetPDGMass()) 139 { 140 sigma = protonL2DataSetMap[zTarget]->FindValue(energyIncident/MeV); 141 if (sigma !=0 && energyIncident > protonL2DataSetMap[zTarget]->GetEnergies(0).back()*MeV) return 0.; 142 } 143 else if (massIncident == aAlpha->GetPDGMass()) 144 { 145 sigma = alphaL2DataSetMap[zTarget]->FindValue(energyIncident/MeV); 146 if (sigma !=0 && energyIncident > alphaL2DataSetMap[zTarget]->GetEnergies(0).back()*MeV) return 0.; 147 } 148 else 149 { 150 sigma = 0.; 151 } 152 } 153 154 // sigma is in internal units: it has been converted from 155 // the input file in barns bt the EmDataset 156 return sigma; 157 } 158 159 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 160 161 G4double G4ecpssrFormFactorLixsModel::CalculateL3CrossSection(G4int zTarget,G4double massIncident, G4double energyIncident) 162 { 163 G4Proton* aProton = G4Proton::Proton(); 164 G4Alpha* aAlpha = G4Alpha::Alpha(); 165 G4double sigma = 0; 166 167 if (energyIncident > 0.1*MeV && energyIncident < 100.*MeV && zTarget < 93 && zTarget > 10) { 168 169 if (massIncident == aProton->GetPDGMass()) 170 { 171 sigma = protonL3DataSetMap[zTarget]->FindValue(energyIncident/MeV); 172 if (sigma !=0 && energyIncident > protonL3DataSetMap[zTarget]->GetEnergies(0).back()*MeV) return 0.; 173 } 174 else if (massIncident == aAlpha->GetPDGMass()) 175 { 176 sigma = alphaL3DataSetMap[zTarget]->FindValue(energyIncident/MeV); 177 if (sigma !=0 && energyIncident > alphaL3DataSetMap[zTarget]->GetEnergies(0).back()*MeV) return 0.; 178 } 179 else 180 { 181 sigma = 0.; 182 } 183 } 184 185 // sigma is in internal units: it has been converted from 186 // the input file in barns bt the EmDataset 187 return sigma; 188 } 189