#include <G4DNARuddIonisationExtendedModel.hh>
Inheritance diagram for G4DNARuddIonisationExtendedModel:
Public Member Functions | |
G4DNARuddIonisationExtendedModel (const G4ParticleDefinition *p=0, const G4String &nam="DNARuddIonisationExtendedModel") | |
virtual | ~G4DNARuddIonisationExtendedModel () |
virtual void | Initialise (const G4ParticleDefinition *, const G4DataVector &) |
virtual G4double | CrossSectionPerVolume (const G4Material *material, const G4ParticleDefinition *p, G4double ekin, G4double emin, G4double emax) |
virtual void | SampleSecondaries (std::vector< G4DynamicParticle * > *, const G4MaterialCutsCouple *, const G4DynamicParticle *, G4double tmin, G4double maxEnergy) |
Protected Attributes | |
G4ParticleChangeForGamma * | fParticleChangeForGamma |
Definition at line 46 of file G4DNARuddIonisationExtendedModel.hh.
G4DNARuddIonisationExtendedModel::G4DNARuddIonisationExtendedModel | ( | const G4ParticleDefinition * | p = 0 , |
|
const G4String & | nam = "DNARuddIonisationExtendedModel" | |||
) |
Definition at line 47 of file G4DNARuddIonisationExtendedModel.cc.
References fParticleChangeForGamma, G4cout, G4endl, and G4VEmModel::SetDeexcitationFlag().
00049 :G4VEmModel(nam),isInitialised(false) 00050 { 00051 // nistwater = G4NistManager::Instance()->FindOrBuildMaterial("G4_WATER"); 00052 fpWaterDensity = 0; 00053 00054 slaterEffectiveCharge[0]=0.; 00055 slaterEffectiveCharge[1]=0.; 00056 slaterEffectiveCharge[2]=0.; 00057 sCoefficient[0]=0.; 00058 sCoefficient[1]=0.; 00059 sCoefficient[2]=0.; 00060 00061 lowEnergyLimitForA[1] = 0 * eV; 00062 lowEnergyLimitForA[2] = 0 * eV; 00063 lowEnergyLimitForA[3] = 0 * eV; 00064 lowEnergyLimitOfModelForA[1] = 100 * eV; 00065 lowEnergyLimitOfModelForA[4] = 1 * keV; 00066 lowEnergyLimitOfModelForA[5] = 0.5 * MeV; // For A = 3 or above, limit is MeV/uma 00067 killBelowEnergyForA[1] = lowEnergyLimitOfModelForA[1]; 00068 killBelowEnergyForA[4] = lowEnergyLimitOfModelForA[4]; 00069 killBelowEnergyForA[5] = lowEnergyLimitOfModelForA[5]; 00070 00071 00072 verboseLevel= 0; 00073 // Verbosity scale: 00074 // 0 = nothing 00075 // 1 = warning for energy non-conservation 00076 // 2 = details of energy budget 00077 // 3 = calculation of cross sections, file openings, sampling of atoms 00078 // 4 = entering in methods 00079 00080 if( verboseLevel>0 ) 00081 { 00082 G4cout << "Rudd ionisation model is constructed " << G4endl; 00083 } 00084 00085 //Mark this model as "applicable" for atomic deexcitation 00086 SetDeexcitationFlag(true); 00087 fAtomDeexcitation = 0; 00088 fParticleChangeForGamma = 0; 00089 }
G4DNARuddIonisationExtendedModel::~G4DNARuddIonisationExtendedModel | ( | ) | [virtual] |
Definition at line 93 of file G4DNARuddIonisationExtendedModel.cc.
00094 { 00095 // Cross section 00096 00097 std::map< G4String,G4DNACrossSectionDataSet*,std::less<G4String> >::iterator pos; 00098 for (pos = tableData.begin(); pos != tableData.end(); ++pos) 00099 { 00100 G4DNACrossSectionDataSet* table = pos->second; 00101 delete table; 00102 } 00103 00104 // The following removal is forbidden G4VEnergyLossModel takes care of deletion 00105 // however coverity will signal this as an error 00106 //if (fAtomDeexcitation) {delete fAtomDeexcitation;} 00107 00108 }
G4double G4DNARuddIonisationExtendedModel::CrossSectionPerVolume | ( | const G4Material * | material, | |
const G4ParticleDefinition * | p, | |||
G4double | ekin, | |||
G4double | emin, | |||
G4double | emax | |||
) | [virtual] |
Reimplemented from G4VEmModel.
Definition at line 374 of file G4DNARuddIonisationExtendedModel.cc.
References FatalException, G4cout, G4endl, G4Exception(), G4Material::GetIndex(), G4DNAGenericIonsManager::GetIon(), G4ParticleDefinition::GetParticleName(), G4DNAGenericIonsManager::Instance(), and G4Proton::ProtonDefinition().
00379 { 00380 if (verboseLevel > 3) 00381 G4cout << "Calling CrossSectionPerVolume() of G4DNARuddIonisationExtendedModel" << G4endl; 00382 00383 // Calculate total cross section for model 00384 00385 G4DNAGenericIonsManager *instance; 00386 instance = G4DNAGenericIonsManager::Instance(); 00387 00388 if ( 00389 particleDefinition != G4Proton::ProtonDefinition() 00390 && 00391 particleDefinition != instance->GetIon("hydrogen") 00392 && 00393 particleDefinition != instance->GetIon("alpha++") 00394 && 00395 particleDefinition != instance->GetIon("alpha+") 00396 && 00397 particleDefinition != instance->GetIon("helium") 00398 && 00399 particleDefinition != instance->GetIon("carbon") 00400 && 00401 particleDefinition != instance->GetIon("nitrogen") 00402 && 00403 particleDefinition != instance->GetIon("oxygen") 00404 && 00405 particleDefinition != instance->GetIon("iron") 00406 ) 00407 00408 return 0; 00409 00410 G4double lowLim = 0; 00411 00412 if ( particleDefinition == G4Proton::ProtonDefinition() 00413 || particleDefinition == instance->GetIon("hydrogen") 00414 ) 00415 00416 lowLim = lowEnergyLimitOfModelForA[1]; 00417 00418 else if ( particleDefinition == instance->GetIon("alpha++") 00419 || particleDefinition == instance->GetIon("alpha+") 00420 || particleDefinition == instance->GetIon("helium") 00421 ) 00422 00423 lowLim = lowEnergyLimitOfModelForA[4]; 00424 00425 else lowLim = lowEnergyLimitOfModelForA[5]; 00426 00427 G4double highLim = 0; 00428 G4double sigma=0; 00429 00430 00431 G4double waterDensity = (*fpWaterDensity)[material->GetIndex()]; 00432 00433 if(waterDensity!= 0.0) 00434 // if (material == nistwater || material->GetBaseMaterial() == nistwater) 00435 { 00436 const G4String& particleName = particleDefinition->GetParticleName(); 00437 00438 std::map< G4String,G4double,std::less<G4String> >::iterator pos2; 00439 pos2 = highEnergyLimit.find(particleName); 00440 00441 if (pos2 != highEnergyLimit.end()) 00442 { 00443 highLim = pos2->second; 00444 } 00445 00446 if (k <= highLim) 00447 { 00448 00449 //SI : XS must not be zero otherwise sampling of secondaries method ignored 00450 00451 if (k < lowLim) k = lowLim; 00452 00453 // 00454 00455 std::map< G4String,G4DNACrossSectionDataSet*,std::less<G4String> >::iterator pos; 00456 pos = tableData.find(particleName); 00457 00458 if (pos != tableData.end()) 00459 { 00460 G4DNACrossSectionDataSet* table = pos->second; 00461 if (table != 0) 00462 { 00463 sigma = table->FindValue(k); 00464 } 00465 } 00466 else 00467 { 00468 G4Exception("G4DNARuddIonisationExtendedModel::CrossSectionPerVolume","em0002", 00469 FatalException,"Model not applicable to particle type."); 00470 } 00471 00472 } // if (k >= lowLim && k < highLim) 00473 00474 if (verboseLevel > 2) 00475 { 00476 G4cout << "__________________________________" << G4endl; 00477 G4cout << "°°° G4DNARuddIonisationExtendedModel - XS INFO START" << G4endl; 00478 G4cout << "°°° Kinetic energy(eV)=" << k/eV << " particle : " << particleDefinition->GetParticleName() << G4endl; 00479 G4cout << "°°° Cross section per water molecule (cm^2)=" << sigma/cm/cm << G4endl; 00480 G4cout << "°°° Cross section per water molecule (cm^-1)=" << sigma*waterDensity/(1./cm) << G4endl; 00481 // G4cout << " - Cross section per water molecule (cm^-1)=" << sigma*material->GetAtomicNumDensityVector()[1]/(1./cm) << G4endl; 00482 G4cout << "°°° G4DNARuddIonisationExtendedModel - XS INFO END" << G4endl; 00483 00484 } 00485 00486 } // if (waterMaterial) 00487 00488 return sigma*waterDensity; 00489 // return sigma*material->GetAtomicNumDensityVector()[1]; 00490 00491 }
void G4DNARuddIonisationExtendedModel::Initialise | ( | const G4ParticleDefinition * | , | |
const G4DataVector & | ||||
) | [virtual] |
Implements G4VEmModel.
Definition at line 112 of file G4DNARuddIonisationExtendedModel.cc.
References G4LossTableManager::AtomDeexcitation(), fParticleChangeForGamma, G4cout, G4endl, G4ParticleDefinition::GetAtomicMass(), G4DNAGenericIonsManager::GetIon(), G4Material::GetMaterial(), G4DNAMolecularMaterial::GetNumMolPerVolTableFor(), G4VEmModel::GetParticleChangeForGamma(), G4ParticleDefinition::GetParticleName(), G4VEmModel::HighEnergyLimit(), G4LossTableManager::Instance(), G4DNAMolecularMaterial::Instance(), G4DNAGenericIonsManager::Instance(), G4DNACrossSectionDataSet::LoadData(), G4VEmModel::LowEnergyLimit(), G4InuclParticleNames::proton, G4Proton::ProtonDefinition(), G4VEmModel::SetHighEnergyLimit(), and G4VEmModel::SetLowEnergyLimit().
00114 { 00115 00116 if (verboseLevel > 3) 00117 G4cout << "Calling G4DNARuddIonisationExtendedModel::Initialise()" << G4endl; 00118 00119 // Energy limits 00120 00121 G4String fileProton("dna/sigma_ionisation_p_rudd"); 00122 G4String fileHydrogen("dna/sigma_ionisation_h_rudd"); 00123 G4String fileAlphaPlusPlus("dna/sigma_ionisation_alphaplusplus_rudd"); 00124 G4String fileAlphaPlus("dna/sigma_ionisation_alphaplus_rudd"); 00125 G4String fileHelium("dna/sigma_ionisation_he_rudd"); 00126 G4String fileCarbon("dna/sigma_ionisation_c_rudd"); 00127 G4String fileNitrogen("dna/sigma_ionisation_n_rudd"); 00128 G4String fileOxygen("dna/sigma_ionisation_o_rudd"); 00129 G4String fileIron("dna/sigma_ionisation_fe_rudd"); 00130 00131 G4DNAGenericIonsManager *instance; 00132 instance = G4DNAGenericIonsManager::Instance(); 00133 G4ParticleDefinition* protonDef = G4Proton::ProtonDefinition(); 00134 G4ParticleDefinition* hydrogenDef = instance->GetIon("hydrogen"); 00135 G4ParticleDefinition* alphaPlusPlusDef = instance->GetIon("alpha++"); 00136 G4ParticleDefinition* alphaPlusDef = instance->GetIon("alpha+"); 00137 G4ParticleDefinition* heliumDef = instance->GetIon("helium"); 00138 G4ParticleDefinition* carbonDef = instance->GetIon("carbon"); 00139 G4ParticleDefinition* nitrogenDef = instance->GetIon("nitrogen"); 00140 G4ParticleDefinition* oxygenDef = instance->GetIon("oxygen"); 00141 G4ParticleDefinition* ironDef = instance->GetIon("iron"); 00142 00143 G4String proton; 00144 G4String hydrogen; 00145 G4String alphaPlusPlus; 00146 G4String alphaPlus; 00147 G4String helium; 00148 G4String carbon; 00149 G4String nitrogen; 00150 G4String oxygen; 00151 G4String iron; 00152 00153 G4double scaleFactor = 1 * m*m; 00154 00155 // LIMITS AND DATA 00156 00157 proton = protonDef->GetParticleName(); 00158 tableFile[proton] = fileProton; 00159 lowEnergyLimit[proton] = lowEnergyLimitForA[1]; 00160 highEnergyLimit[proton] = 500. * keV; 00161 00162 // Cross section 00163 00164 G4DNACrossSectionDataSet* tableProton = new G4DNACrossSectionDataSet(new G4LogLogInterpolation, 00165 eV, 00166 scaleFactor ); 00167 tableProton->LoadData(fileProton); 00168 tableData[proton] = tableProton; 00169 00170 // ********************************************************************************************** 00171 00172 hydrogen = hydrogenDef->GetParticleName(); 00173 tableFile[hydrogen] = fileHydrogen; 00174 00175 lowEnergyLimit[hydrogen] = lowEnergyLimitForA[1]; 00176 highEnergyLimit[hydrogen] = 100. * MeV; 00177 00178 // Cross section 00179 00180 G4DNACrossSectionDataSet* tableHydrogen = new G4DNACrossSectionDataSet(new G4LogLogInterpolation, 00181 eV, 00182 scaleFactor ); 00183 tableHydrogen->LoadData(fileHydrogen); 00184 00185 tableData[hydrogen] = tableHydrogen; 00186 00187 // ********************************************************************************************** 00188 00189 alphaPlusPlus = alphaPlusPlusDef->GetParticleName(); 00190 tableFile[alphaPlusPlus] = fileAlphaPlusPlus; 00191 00192 lowEnergyLimit[alphaPlusPlus] = lowEnergyLimitForA[4]; 00193 highEnergyLimit[alphaPlusPlus] = 400. * MeV; 00194 00195 // Cross section 00196 00197 G4DNACrossSectionDataSet* tableAlphaPlusPlus = new G4DNACrossSectionDataSet(new G4LogLogInterpolation, 00198 eV, 00199 scaleFactor ); 00200 tableAlphaPlusPlus->LoadData(fileAlphaPlusPlus); 00201 00202 tableData[alphaPlusPlus] = tableAlphaPlusPlus; 00203 00204 // ********************************************************************************************** 00205 00206 alphaPlus = alphaPlusDef->GetParticleName(); 00207 tableFile[alphaPlus] = fileAlphaPlus; 00208 00209 lowEnergyLimit[alphaPlus] = lowEnergyLimitForA[4]; 00210 highEnergyLimit[alphaPlus] = 400. * MeV; 00211 00212 // Cross section 00213 00214 G4DNACrossSectionDataSet* tableAlphaPlus = new G4DNACrossSectionDataSet(new G4LogLogInterpolation, 00215 eV, 00216 scaleFactor ); 00217 tableAlphaPlus->LoadData(fileAlphaPlus); 00218 tableData[alphaPlus] = tableAlphaPlus; 00219 00220 // ********************************************************************************************** 00221 00222 helium = heliumDef->GetParticleName(); 00223 tableFile[helium] = fileHelium; 00224 00225 lowEnergyLimit[helium] = lowEnergyLimitForA[4]; 00226 highEnergyLimit[helium] = 400. * MeV; 00227 00228 // Cross section 00229 00230 G4DNACrossSectionDataSet* tableHelium = new G4DNACrossSectionDataSet(new G4LogLogInterpolation, 00231 eV, 00232 scaleFactor ); 00233 tableHelium->LoadData(fileHelium); 00234 tableData[helium] = tableHelium; 00235 00236 // ********************************************************************************************** 00237 00238 carbon = carbonDef->GetParticleName(); 00239 tableFile[carbon] = fileCarbon; 00240 00241 lowEnergyLimit[carbon] = lowEnergyLimitForA[5] * particle->GetAtomicMass(); 00242 highEnergyLimit[carbon] = 1e6* particle->GetAtomicMass() * MeV; 00243 00244 // Cross section 00245 00246 G4DNACrossSectionDataSet* tableCarbon = new G4DNACrossSectionDataSet(new G4LogLogInterpolation, 00247 eV, 00248 scaleFactor ); 00249 tableCarbon->LoadData(fileCarbon); 00250 tableData[carbon] = tableCarbon; 00251 00252 // ********************************************************************************************** 00253 00254 oxygen = oxygenDef->GetParticleName(); 00255 tableFile[oxygen] = fileOxygen; 00256 00257 lowEnergyLimit[oxygen] = lowEnergyLimitForA[5]* particle->GetAtomicMass(); 00258 highEnergyLimit[oxygen] = 1e6* particle->GetAtomicMass()* MeV; 00259 00260 // Cross section 00261 00262 G4DNACrossSectionDataSet* tableOxygen = new G4DNACrossSectionDataSet(new G4LogLogInterpolation, 00263 eV, 00264 scaleFactor ); 00265 tableOxygen->LoadData(fileOxygen); 00266 tableData[oxygen] = tableOxygen; 00267 00268 // ********************************************************************************************** 00269 00270 nitrogen = nitrogenDef->GetParticleName(); 00271 tableFile[nitrogen] = fileNitrogen; 00272 00273 lowEnergyLimit[nitrogen] = lowEnergyLimitForA[5]* particle->GetAtomicMass(); 00274 highEnergyLimit[nitrogen] = 1e6* particle->GetAtomicMass()* MeV; 00275 00276 // Cross section 00277 00278 G4DNACrossSectionDataSet* tableNitrogen = new G4DNACrossSectionDataSet(new G4LogLogInterpolation, 00279 eV, 00280 scaleFactor ); 00281 tableNitrogen->LoadData(fileNitrogen); 00282 tableData[nitrogen] = tableNitrogen; 00283 00284 // ********************************************************************************************** 00285 00286 iron = ironDef->GetParticleName(); 00287 tableFile[iron] = fileIron; 00288 00289 lowEnergyLimit[iron] = lowEnergyLimitForA[5]* particle->GetAtomicMass(); 00290 highEnergyLimit[iron] = 1e6* particle->GetAtomicMass()* MeV; 00291 00292 // Cross section 00293 00294 G4DNACrossSectionDataSet* tableIron = new G4DNACrossSectionDataSet(new G4LogLogInterpolation, 00295 eV, 00296 scaleFactor ); 00297 tableIron->LoadData(fileIron); 00298 tableData[iron] = tableIron; 00299 00300 // ZF Following lines can be replaced by: 00301 SetLowEnergyLimit(lowEnergyLimit[particle->GetParticleName()]); 00302 SetHighEnergyLimit(highEnergyLimit[particle->GetParticleName()]); 00303 // at least for HZE 00304 00305 /* 00306 if (particle==protonDef) 00307 { 00308 SetLowEnergyLimit(lowEnergyLimit[proton]); 00309 SetHighEnergyLimit(highEnergyLimit[proton]); 00310 } 00311 00312 if (particle==hydrogenDef) 00313 { 00314 SetLowEnergyLimit(lowEnergyLimit[hydrogen]); 00315 SetHighEnergyLimit(highEnergyLimit[hydrogen]); 00316 } 00317 00318 if (particle==heliumDef) 00319 { 00320 SetLowEnergyLimit(lowEnergyLimit[helium]); 00321 SetHighEnergyLimit(highEnergyLimit[helium]); 00322 } 00323 00324 if (particle==alphaPlusDef) 00325 { 00326 SetLowEnergyLimit(lowEnergyLimit[alphaPlus]); 00327 SetHighEnergyLimit(highEnergyLimit[alphaPlus]); 00328 } 00329 00330 if (particle==alphaPlusPlusDef) 00331 { 00332 SetLowEnergyLimit(lowEnergyLimit[alphaPlusPlus]); 00333 SetHighEnergyLimit(highEnergyLimit[alphaPlusPlus]); 00334 } 00335 00336 if (particle==carbonDef) 00337 { 00338 SetLowEnergyLimit(lowEnergyLimit[carbon]); 00339 SetHighEnergyLimit(highEnergyLimit[carbon]); 00340 } 00341 00342 if (particle==oxygenDef) 00343 { 00344 SetLowEnergyLimit(lowEnergyLimit[oxygen]); 00345 SetHighEnergyLimit(highEnergyLimit[oxygen]); 00346 }*/ 00347 00348 //---------------------------------------------------------------------- 00349 00350 if( verboseLevel>0 ) 00351 { 00352 G4cout << "Rudd ionisation model is initialized " << G4endl 00353 << "Energy range: " 00354 << LowEnergyLimit() / eV << " eV - " 00355 << HighEnergyLimit() / keV << " keV for " 00356 << particle->GetParticleName() 00357 << G4endl; 00358 } 00359 00360 // Initialize water density pointer 00361 fpWaterDensity = G4DNAMolecularMaterial::Instance()->GetNumMolPerVolTableFor(G4Material::GetMaterial("G4_WATER")); 00362 00363 // 00364 00365 fAtomDeexcitation = G4LossTableManager::Instance()->AtomDeexcitation(); 00366 00367 if (isInitialised) { return; } 00368 fParticleChangeForGamma = GetParticleChangeForGamma(); 00369 isInitialised = true; 00370 }
void G4DNARuddIonisationExtendedModel::SampleSecondaries | ( | std::vector< G4DynamicParticle * > * | , | |
const G4MaterialCutsCouple * | , | |||
const G4DynamicParticle * | , | |||
G4double | tmin, | |||
G4double | maxEnergy | |||
) | [virtual] |
Implements G4VEmModel.
Definition at line 495 of file G4DNARuddIonisationExtendedModel.cc.
References G4InuclSpecialFunctions::bindingEnergy(), G4DNAChemistryManager::CreateWaterMolecule(), eIonizedMolecule, G4Electron::Electron(), fKShell, fParticleChangeForGamma, fStopAndKill, G4cout, G4endl, G4VAtomDeexcitation::GenerateParticles(), G4ParticleDefinition::GetAtomicMass(), G4VAtomDeexcitation::GetAtomicShell(), G4ParticleChangeForGamma::GetCurrentTrack(), G4DynamicParticle::GetDefinition(), G4DynamicParticle::GetKineticEnergy(), G4DynamicParticle::GetMomentumDirection(), G4ParticleDefinition::GetParticleName(), G4DNAChemistryManager::Instance(), G4DNAWaterIonisationStructure::IonisationEnergy(), G4VParticleChange::ProposeLocalEnergyDeposit(), G4ParticleChangeForGamma::ProposeMomentumDirection(), G4VParticleChange::ProposeTrackStatus(), and G4ParticleChangeForGamma::SetProposedKineticEnergy().
00500 { 00501 if (verboseLevel > 3) 00502 G4cout << "Calling SampleSecondaries() of G4DNARuddIonisationExtendedModel" << G4endl; 00503 00504 G4double lowLim = 0; 00505 G4double highLim = 0; 00506 00507 // ZF: the following line summarizes the commented part 00508 if(particle->GetDefinition()->GetAtomicMass() <= 4) lowLim = killBelowEnergyForA[particle->GetDefinition()->GetAtomicMass()]; 00509 else lowLim = killBelowEnergyForA[5]*particle->GetDefinition()->GetAtomicMass(); 00510 00511 /* 00512 if(particle->GetDefinition()->GetAtomicMass() >= 5) lowLim = killBelowEnergyForA[5]*particle->GetDefinition()->GetAtomicMass(); 00513 00514 00515 if ( particle->GetDefinition() == G4Proton::ProtonDefinition() 00516 || particle->GetDefinition() == instance->GetIon("hydrogen") 00517 ) 00518 00519 lowLim = killBelowEnergyForA[1]; 00520 00521 if ( particle->GetDefinition() == instance->GetIon("alpha++") 00522 || particle->GetDefinition() == instance->GetIon("alpha+") 00523 || particle->GetDefinition() == instance->GetIon("helium") 00524 ) 00525 00526 lowLim = killBelowEnergyForA[4];*/ 00527 00528 G4double k = particle->GetKineticEnergy(); 00529 00530 const G4String& particleName = particle->GetDefinition()->GetParticleName(); 00531 00532 // SI - the following is useless since lowLim is already defined 00533 /* 00534 std::map< G4String,G4double,std::less<G4String> >::iterator pos1; 00535 pos1 = lowEnergyLimit.find(particleName); 00536 00537 if (pos1 != lowEnergyLimit.end()) 00538 { 00539 lowLim = pos1->second; 00540 } 00541 */ 00542 00543 std::map< G4String,G4double,std::less<G4String> >::iterator pos2; 00544 pos2 = highEnergyLimit.find(particleName); 00545 00546 if (pos2 != highEnergyLimit.end())highLim = pos2->second; 00547 00548 if (k >= lowLim && k < highLim) 00549 { 00550 G4ParticleDefinition* definition = particle->GetDefinition(); 00551 G4ParticleMomentum primaryDirection = particle->GetMomentumDirection(); 00552 /* 00553 G4double particleMass = definition->GetPDGMass(); 00554 G4double totalEnergy = k + particleMass; 00555 G4double pSquare = k*(totalEnergy+particleMass); 00556 G4double totalMomentum = std::sqrt(pSquare); 00557 */ 00558 00559 G4int ionizationShell = RandomSelect(k,particleName); 00560 00561 // sample deexcitation 00562 // here we assume that H_{2}O electronic levels are the same of Oxigen. 00563 // this can be considered true with a rough 10% error in energy on K-shell, 00564 00565 G4int secNumberInit = 0; // need to know at a certain point the enrgy of secondaries 00566 G4int secNumberFinal = 0; // So I'll make the diference and then sum the energies 00567 G4double bindingEnergy = 0; 00568 bindingEnergy = waterStructure.IonisationEnergy(ionizationShell); 00569 00570 if(fAtomDeexcitation) { 00571 G4int Z = 8; 00572 G4AtomicShellEnumerator as = fKShell; 00573 00574 if (ionizationShell <5 && ionizationShell >1) 00575 { 00576 as = G4AtomicShellEnumerator(4-ionizationShell); 00577 } 00578 else if (ionizationShell <2) 00579 { 00580 as = G4AtomicShellEnumerator(3); 00581 } 00582 00583 00584 // DEBUG 00585 // if (ionizationShell == 4) { 00586 // 00587 // G4cout << "Z: " << Z << " as: " << as 00588 // << " ionizationShell: " << ionizationShell << " bindingEnergy: "<< bindingEnergy/eV << G4endl; 00589 // G4cout << "Press <Enter> key to continue..." << G4endl; 00590 // G4cin.ignore(); 00591 // } 00592 00593 00594 const G4AtomicShell* shell = fAtomDeexcitation->GetAtomicShell(Z, as); 00595 secNumberInit = fvect->size(); 00596 fAtomDeexcitation->GenerateParticles(fvect, shell, Z, 0, 0); 00597 secNumberFinal = fvect->size(); 00598 } 00599 00600 00601 G4double secondaryKinetic = RandomizeEjectedElectronEnergy(definition,k,ionizationShell); 00602 00603 G4double cosTheta = 0.; 00604 G4double phi = 0.; 00605 RandomizeEjectedElectronDirection(definition, k,secondaryKinetic, cosTheta, phi, ionizationShell); 00606 00607 G4double sinTheta = std::sqrt(1.-cosTheta*cosTheta); 00608 G4double dirX = sinTheta*std::cos(phi); 00609 G4double dirY = sinTheta*std::sin(phi); 00610 G4double dirZ = cosTheta; 00611 G4ThreeVector deltaDirection(dirX,dirY,dirZ); 00612 deltaDirection.rotateUz(primaryDirection); 00613 00614 // Ignored for ions on electrons 00615 /* 00616 G4double deltaTotalMomentum = std::sqrt(secondaryKinetic*(secondaryKinetic + 2.*electron_mass_c2 )); 00617 00618 G4double finalPx = totalMomentum*primaryDirection.x() - deltaTotalMomentum*deltaDirection.x(); 00619 G4double finalPy = totalMomentum*primaryDirection.y() - deltaTotalMomentum*deltaDirection.y(); 00620 G4double finalPz = totalMomentum*primaryDirection.z() - deltaTotalMomentum*deltaDirection.z(); 00621 G4double finalMomentum = std::sqrt(finalPx*finalPx+finalPy*finalPy+finalPz*finalPz); 00622 finalPx /= finalMomentum; 00623 finalPy /= finalMomentum; 00624 finalPz /= finalMomentum; 00625 00626 G4ThreeVector direction; 00627 direction.set(finalPx,finalPy,finalPz); 00628 00629 fParticleChangeForGamma->ProposeMomentumDirection(direction.unit()) ; 00630 */ 00631 00632 fParticleChangeForGamma->ProposeMomentumDirection(primaryDirection); 00633 G4double scatteredEnergy = k-bindingEnergy-secondaryKinetic; 00634 G4double deexSecEnergy = 0; 00635 for (G4int j=secNumberInit; j < secNumberFinal; j++) { 00636 00637 deexSecEnergy = deexSecEnergy + (*fvect)[j]->GetKineticEnergy(); 00638 00639 } 00640 00641 fParticleChangeForGamma->SetProposedKineticEnergy(scatteredEnergy); 00642 fParticleChangeForGamma->ProposeLocalEnergyDeposit(k-scatteredEnergy-secondaryKinetic-deexSecEnergy); 00643 00644 G4DynamicParticle* dp = new G4DynamicParticle (G4Electron::Electron(),deltaDirection,secondaryKinetic) ; 00645 fvect->push_back(dp); 00646 00647 const G4Track * theIncomingTrack = fParticleChangeForGamma->GetCurrentTrack(); 00648 G4DNAChemistryManager::Instance()->CreateWaterMolecule(eIonizedMolecule, 00649 ionizationShell, 00650 theIncomingTrack); 00651 } 00652 00653 // SI - not useful since low energy of model is 0 eV 00654 00655 if (k < lowLim) 00656 { 00657 fParticleChangeForGamma->SetProposedKineticEnergy(0.); 00658 fParticleChangeForGamma->ProposeTrackStatus(fStopAndKill); 00659 fParticleChangeForGamma->ProposeLocalEnergyDeposit(k); 00660 } 00661 00662 }
Definition at line 72 of file G4DNARuddIonisationExtendedModel.hh.
Referenced by G4DNARuddIonisationExtendedModel(), Initialise(), and SampleSecondaries().