#include <G4LivermoreBremsstrahlungModel.hh>
Inheritance diagram for G4LivermoreBremsstrahlungModel:
Public Member Functions | |
G4LivermoreBremsstrahlungModel (const G4ParticleDefinition *p=0, const G4String &processName="LowEnBrem") | |
virtual | ~G4LivermoreBremsstrahlungModel () |
virtual void | Initialise (const G4ParticleDefinition *, const G4DataVector &) |
virtual G4double | ComputeCrossSectionPerAtom (const G4ParticleDefinition *, G4double kinEnergy, G4double Z, G4double A=0, G4double cut=0, G4double emax=DBL_MAX) |
virtual void | SampleSecondaries (std::vector< G4DynamicParticle * > *, const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double tmin, G4double maxEnergy) |
virtual G4double | ComputeDEDXPerVolume (const G4Material *, const G4ParticleDefinition *, G4double kineticEnergy, G4double cutEnergy) |
virtual G4double | MinEnergyCut (const G4ParticleDefinition *, const G4MaterialCutsCouple *) |
void | SetVerboseLevel (G4int vl) |
Protected Attributes | |
G4ParticleChangeForLoss * | fParticleChange |
Definition at line 57 of file G4LivermoreBremsstrahlungModel.hh.
G4LivermoreBremsstrahlungModel::G4LivermoreBremsstrahlungModel | ( | const G4ParticleDefinition * | p = 0 , |
|
const G4String & | processName = "LowEnBrem" | |||
) |
Definition at line 72 of file G4LivermoreBremsstrahlungModel.cc.
References G4VEmModel::SetAngularDistribution(), and G4VEmModel::SetHighEnergyLimit().
00074 :G4VEmModel(nam),fParticleChange(0),isInitialised(false), 00075 crossSectionHandler(0),energySpectrum(0) 00076 { 00077 fIntrinsicLowEnergyLimit = 10.0*eV; 00078 fIntrinsicHighEnergyLimit = 100.0*GeV; 00079 fNBinEnergyLoss = 360; 00080 // SetLowEnergyLimit(fIntrinsicLowEnergyLimit); 00081 SetHighEnergyLimit(fIntrinsicHighEnergyLimit); 00082 // 00083 verboseLevel = 0; 00084 SetAngularDistribution(new G4Generator2BS()); 00085 // 00086 //generatorName = "tsai"; 00087 //angularDistribution = new G4ModifiedTsai("TsaiGenerator"); //default generator 00088 // 00089 //TsaiAngularDistribution = new G4ModifiedTsai("TsaiGenerator"); 00090 // 00091 }
G4LivermoreBremsstrahlungModel::~G4LivermoreBremsstrahlungModel | ( | ) | [virtual] |
Definition at line 95 of file G4LivermoreBremsstrahlungModel.cc.
00096 { 00097 if (crossSectionHandler) delete crossSectionHandler; 00098 if (energySpectrum) delete energySpectrum; 00099 energyBins.clear(); 00100 //delete angularDistribution; 00101 //delete TsaiAngularDistribution; 00102 }
G4double G4LivermoreBremsstrahlungModel::ComputeCrossSectionPerAtom | ( | const G4ParticleDefinition * | , | |
G4double | kinEnergy, | |||
G4double | Z, | |||
G4double | A = 0 , |
|||
G4double | cut = 0 , |
|||
G4double | emax = DBL_MAX | |||
) | [virtual] |
Reimplemented from G4VEmModel.
Definition at line 189 of file G4LivermoreBremsstrahlungModel.cc.
References FatalException, G4cout, G4endl, G4Exception(), and G4BremsstrahlungCrossSectionHandler::GetCrossSectionAboveThresholdForElement().
00194 { 00195 G4int iZ = (G4int) Z; 00196 if (!crossSectionHandler) 00197 { 00198 G4Exception("G4LivermoreBremsstrahlungModel::ComputeCrossSectionPerAtom", 00199 "em1007",FatalException,"The cross section handler is not correctly initialized"); 00200 return 0; 00201 } 00202 00203 //The cut is already included in the crossSectionHandler 00204 G4double cs = 00205 crossSectionHandler->GetCrossSectionAboveThresholdForElement(energy,cutEnergy,iZ); 00206 00207 if (verboseLevel > 1) 00208 { 00209 G4cout << "G4LivermoreBremsstrahlungModel " << G4endl; 00210 G4cout << "Cross section for gamma emission > " << cutEnergy/keV << " keV at " << 00211 energy/keV << " keV and Z = " << iZ << " --> " << cs/barn << " barn" << G4endl; 00212 } 00213 return cs; 00214 }
G4double G4LivermoreBremsstrahlungModel::ComputeDEDXPerVolume | ( | const G4Material * | , | |
const G4ParticleDefinition * | , | |||
G4double | kineticEnergy, | |||
G4double | cutEnergy | |||
) | [virtual] |
Reimplemented from G4VEmModel.
Definition at line 219 of file G4LivermoreBremsstrahlungModel.cc.
References G4VEnergySpectrum::AverageEnergy(), G4VCrossSectionHandler::FindValue(), G4cout, G4endl, G4Material::GetAtomicNumDensityVector(), G4Material::GetElementVector(), and G4Material::GetNumberOfElements().
00223 { 00224 G4double sPower = 0.0; 00225 00226 const G4ElementVector* theElementVector = material->GetElementVector(); 00227 size_t NumberOfElements = material->GetNumberOfElements() ; 00228 const G4double* theAtomicNumDensityVector = 00229 material->GetAtomicNumDensityVector(); 00230 00231 // loop for elements in the material 00232 for (size_t iel=0; iel<NumberOfElements; iel++ ) 00233 { 00234 G4int iZ = (G4int)((*theElementVector)[iel]->GetZ()); 00235 G4double e = energySpectrum->AverageEnergy(iZ, 0.0,cutEnergy, 00236 kineticEnergy); 00237 G4double cs= crossSectionHandler->FindValue(iZ,kineticEnergy); 00238 sPower += e * cs * theAtomicNumDensityVector[iel]; 00239 } 00240 00241 if (verboseLevel > 2) 00242 { 00243 G4cout << "G4LivermoreBremsstrahlungModel " << G4endl; 00244 G4cout << "Stopping power < " << cutEnergy/keV << " keV at " << 00245 kineticEnergy/keV << " keV = " << sPower/(keV/mm) << " keV/mm" << G4endl; 00246 } 00247 00248 return sPower; 00249 }
void G4LivermoreBremsstrahlungModel::Initialise | ( | const G4ParticleDefinition * | , | |
const G4DataVector & | ||||
) | [virtual] |
Implements G4VEmModel.
Definition at line 106 of file G4LivermoreBremsstrahlungModel.cc.
References G4VCrossSectionHandler::BuildMeanFreePathForMaterials(), G4VCrossSectionHandler::Clear(), G4Electron::Electron(), FatalException, fParticleChange, G4cout, G4endl, G4Exception(), G4VEmModel::GetParticleChangeForLoss(), G4VEmModel::HighEnergyLimit(), G4VCrossSectionHandler::Initialise(), G4VCrossSectionHandler::LoadShellData(), G4VEmModel::LowEnergyLimit(), G4VEnergySpectrum::PrintData(), and G4VCrossSectionHandler::PrintData().
00108 { 00109 //Check that the Livermore Bremsstrahlung is NOT attached to e+ 00110 if (particle != G4Electron::Electron()) 00111 { 00112 G4Exception("G4LivermoreBremsstrahlungModel::Initialise", 00113 "em0002",FatalException,"Livermore Bremsstrahlung Model is applicable only to electrons"); 00114 } 00115 //Prepare energy spectrum 00116 if (energySpectrum) 00117 { 00118 delete energySpectrum; 00119 energySpectrum = 0; 00120 } 00121 00122 energyBins.clear(); 00123 for(size_t i=0; i<15; i++) 00124 { 00125 G4double x = 0.1*((G4double)i); 00126 if(i == 0) x = 0.01; 00127 if(i == 10) x = 0.95; 00128 if(i == 11) x = 0.97; 00129 if(i == 12) x = 0.99; 00130 if(i == 13) x = 0.995; 00131 if(i == 14) x = 1.0; 00132 energyBins.push_back(x); 00133 } 00134 const G4String dataName("/brem/br-sp.dat"); 00135 energySpectrum = new G4eBremsstrahlungSpectrum(energyBins,dataName); 00136 00137 if (verboseLevel > 0) 00138 G4cout << "G4eBremsstrahlungSpectrum is initialized" << G4endl; 00139 00140 //Initialize cross section handler 00141 if (crossSectionHandler) 00142 { 00143 delete crossSectionHandler; 00144 crossSectionHandler = 0; 00145 } 00146 G4VDataSetAlgorithm* interpolation = 0;//new G4SemiLogInterpolation(); 00147 crossSectionHandler = new G4BremsstrahlungCrossSectionHandler(energySpectrum,interpolation); 00148 crossSectionHandler->Initialise(0,LowEnergyLimit(),HighEnergyLimit(), 00149 fNBinEnergyLoss); 00150 crossSectionHandler->Clear(); 00151 crossSectionHandler->LoadShellData("brem/br-cs-"); 00152 //This is used to retrieve cross section values later on 00153 G4VEMDataSet* p = crossSectionHandler->BuildMeanFreePathForMaterials(&cuts); 00154 delete p; 00155 00156 if (verboseLevel > 0) 00157 { 00158 G4cout << "Livermore Bremsstrahlung model is initialized " << G4endl 00159 << "Energy range: " 00160 << LowEnergyLimit() / keV << " keV - " 00161 << HighEnergyLimit() / GeV << " GeV" 00162 << G4endl; 00163 } 00164 00165 if (verboseLevel > 1) 00166 { 00167 G4cout << "Cross section data: " << G4endl; 00168 crossSectionHandler->PrintData(); 00169 G4cout << "Parameters: " << G4endl; 00170 energySpectrum->PrintData(); 00171 } 00172 00173 if(isInitialised) return; 00174 fParticleChange = GetParticleChangeForLoss(); 00175 isInitialised = true; 00176 }
G4double G4LivermoreBremsstrahlungModel::MinEnergyCut | ( | const G4ParticleDefinition * | , | |
const G4MaterialCutsCouple * | ||||
) | [virtual] |
void G4LivermoreBremsstrahlungModel::SampleSecondaries | ( | std::vector< G4DynamicParticle * > * | , | |
const G4MaterialCutsCouple * | , | |||
const G4DynamicParticle * | , | |||
G4double | tmin, | |||
G4double | maxEnergy | |||
) | [virtual] |
Implements G4VEmModel.
Definition at line 253 of file G4LivermoreBremsstrahlungModel.cc.
References fParticleChange, G4cout, G4endl, G4Gamma::Gamma(), G4VEmModel::GetAngularDistribution(), G4DynamicParticle::GetKineticEnergy(), G4MaterialCutsCouple::GetMaterial(), G4DynamicParticle::GetMomentumDirection(), G4VParticleChange::ProposeLocalEnergyDeposit(), G4ParticleChangeForLoss::ProposeMomentumDirection(), G4VEmAngularDistribution::SampleDirection(), G4VEnergySpectrum::SampleEnergy(), G4VCrossSectionHandler::SelectRandomAtom(), and G4ParticleChangeForLoss::SetProposedKineticEnergy().
00258 { 00259 00260 G4double kineticEnergy = aDynamicParticle->GetKineticEnergy(); 00261 00262 // this is neede for pathalogical cases of no ionisation 00263 if (kineticEnergy <= fIntrinsicLowEnergyLimit) 00264 { 00265 fParticleChange->SetProposedKineticEnergy(0.); 00266 fParticleChange->ProposeLocalEnergyDeposit(kineticEnergy); 00267 return; 00268 } 00269 00270 //Sample material 00271 G4int Z = crossSectionHandler->SelectRandomAtom(couple, kineticEnergy); 00272 00273 //Sample gamma energy 00274 G4double tGamma = energySpectrum->SampleEnergy(Z, energyCut, kineticEnergy, kineticEnergy); 00275 //nothing happens 00276 if (tGamma == 0.) { return; } 00277 00278 G4double totalEnergy = kineticEnergy + electron_mass_c2; 00279 G4double finalEnergy = kineticEnergy - tGamma; // electron final energy 00280 00281 //Sample gamma direction 00282 G4ThreeVector gammaDirection = 00283 GetAngularDistribution()->SampleDirection(aDynamicParticle, 00284 totalEnergy-tGamma, 00285 Z, 00286 couple->GetMaterial()); 00287 00288 G4ThreeVector electronDirection = aDynamicParticle->GetMomentumDirection(); 00289 00290 //Update the incident particle 00291 if (finalEnergy < 0.) 00292 { 00293 // Kinematic problem 00294 tGamma = kineticEnergy; 00295 fParticleChange->SetProposedKineticEnergy(0.); 00296 } 00297 else 00298 { 00299 G4double momentum = std::sqrt((totalEnergy + electron_mass_c2)*kineticEnergy); 00300 G4double finalX = momentum*electronDirection.x() - tGamma*gammaDirection.x(); 00301 G4double finalY = momentum*electronDirection.y() - tGamma*gammaDirection.y(); 00302 G4double finalZ = momentum*electronDirection.z() - tGamma*gammaDirection.z(); 00303 G4double norm = 1./std::sqrt(finalX*finalX + finalY*finalY + finalZ*finalZ); 00304 00305 fParticleChange->ProposeMomentumDirection(finalX*norm, finalY*norm, finalZ*norm); 00306 fParticleChange->SetProposedKineticEnergy(finalEnergy); 00307 } 00308 00309 //Generate the bremsstrahlung gamma 00310 G4DynamicParticle* aGamma= new G4DynamicParticle (G4Gamma::Gamma(), 00311 gammaDirection, tGamma); 00312 fvect->push_back(aGamma); 00313 00314 if (verboseLevel > 1) 00315 { 00316 G4cout << "-----------------------------------------------------------" << G4endl; 00317 G4cout << "Energy balance from G4LivermoreBremsstrahlung" << G4endl; 00318 G4cout << "Incoming primary energy: " << kineticEnergy/keV << " keV" << G4endl; 00319 G4cout << "-----------------------------------------------------------" << G4endl; 00320 G4cout << "Outgoing primary energy: " << finalEnergy/keV << " keV" << G4endl; 00321 G4cout << "Gamma ray " << tGamma/keV << " keV" << G4endl; 00322 G4cout << "Total final state: " << (finalEnergy+tGamma)/keV << " keV" << G4endl; 00323 G4cout << "-----------------------------------------------------------" << G4endl; 00324 } 00325 if (verboseLevel > 0) 00326 { 00327 G4double energyDiff = std::fabs(finalEnergy+tGamma-kineticEnergy); 00328 if (energyDiff > 0.05*keV) 00329 G4cout << "G4LivermoreBremsstrahlung WARNING: problem with energy conservation: " 00330 << (finalEnergy+tGamma)/keV << " keV (final) vs. " 00331 << kineticEnergy/keV << " keV (initial)" << G4endl; 00332 } 00333 return; 00334 }
void G4LivermoreBremsstrahlungModel::SetVerboseLevel | ( | G4int | vl | ) | [inline] |
Definition at line 92 of file G4LivermoreBremsstrahlungModel.hh.
Referenced by Initialise(), and SampleSecondaries().