#include <G4DNASancheSolvatationModel.hh>
Inheritance diagram for G4DNASancheSolvatationModel:
Public Member Functions | |
G4DNASancheSolvatationModel (const G4ParticleDefinition *p=0, const G4String &nam="DNASancheSolvatationModel") | |
virtual | ~G4DNASancheSolvatationModel () |
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) |
void | SetVerbose (int) |
Protected Member Functions | |
G4ThreeVector | RadialDistributionOfProducts (G4double Rrms) const |
Protected Attributes | |
const std::vector< G4double > * | fpWaterDensity |
G4ParticleChangeForGamma * | fParticleChangeForGamma |
G4bool | fIsInitialised |
G4int | fVerboseLevel |
Definition at line 51 of file G4DNASancheSolvatationModel.hh.
G4DNASancheSolvatationModel::G4DNASancheSolvatationModel | ( | const G4ParticleDefinition * | p = 0 , |
|
const G4String & | nam = "DNASancheSolvatationModel" | |||
) |
Definition at line 49 of file G4DNASancheSolvatationModel.cc.
References G4DNAWaterExcitationStructure::ExcitationEnergy(), fParticleChangeForGamma, fpWaterDensity, fVerboseLevel, G4VEmModel::SetHighEnergyLimit(), and G4VEmModel::SetLowEnergyLimit().
00050 : 00051 G4VEmModel(nam),fIsInitialised(false) 00052 { 00053 fVerboseLevel = 0 ; 00054 SetLowEnergyLimit(0.); 00055 G4DNAWaterExcitationStructure exStructure ; 00056 SetHighEnergyLimit(exStructure.ExcitationEnergy(0)); 00057 fParticleChangeForGamma = 0; 00058 // fNistWater = G4NistManager::Instance()->FindOrBuildMaterial("G4_WATER"); 00059 fpWaterDensity = 0; 00060 }
G4DNASancheSolvatationModel::~G4DNASancheSolvatationModel | ( | ) | [virtual] |
G4double G4DNASancheSolvatationModel::CrossSectionPerVolume | ( | const G4Material * | material, | |
const G4ParticleDefinition * | p, | |||
G4double | ekin, | |||
G4double | emin, | |||
G4double | emax | |||
) | [virtual] |
Reimplemented from G4VEmModel.
Definition at line 94 of file G4DNASancheSolvatationModel.cc.
References DBL_MAX, fVerboseLevel, G4cout, G4endl, G4Material::GetIndex(), and G4VEmModel::HighEnergyLimit().
00099 { 00100 #ifdef G4VERBOSE 00101 if (fVerboseLevel > 1) 00102 G4cout << "Calling CrossSectionPerVolume() of G4SancheSolvatationModel" << G4endl; 00103 #endif 00104 00105 if(ekin > HighEnergyLimit()) 00106 { 00107 return 0.0; 00108 } 00109 00110 G4double waterDensity = (*fpWaterDensity)[material->GetIndex()]; 00111 00112 if(waterDensity!= 0.0) 00113 // if (material == nistwater || material->GetBaseMaterial() == nistwater) 00114 { 00115 if (ekin <= HighEnergyLimit()) 00116 { 00117 return DBL_MAX; 00118 } 00119 } 00120 return 0. ; 00121 }
void G4DNASancheSolvatationModel::Initialise | ( | const G4ParticleDefinition * | , | |
const G4DataVector & | ||||
) | [virtual] |
Implements G4VEmModel.
Definition at line 67 of file G4DNASancheSolvatationModel.cc.
References G4Electron::ElectronDefinition(), FatalErrorInArgument, fIsInitialised, fParticleChangeForGamma, fpWaterDensity, fVerboseLevel, G4cout, G4endl, G4Exception(), G4Material::GetMaterial(), G4DNAMolecularMaterial::GetNumMolPerVolTableFor(), G4VEmModel::GetParticleChangeForGamma(), and G4DNAMolecularMaterial::Instance().
00069 { 00070 #ifdef G4VERBOSE 00071 if (fVerboseLevel) 00072 G4cout << "Calling G4SancheSolvatationModel::Initialise()" << G4endl; 00073 #endif 00074 if (particleDefinition != G4Electron::ElectronDefinition()) 00075 { 00076 G4ExceptionDescription exceptionDescription ; 00077 exceptionDescription << "Attempting to calculate cross section for wrong particle"; 00078 G4Exception("G4DNASancheSolvatationModel::CrossSectionPerVolume","G4DNASancheSolvatationModel001", 00079 FatalErrorInArgument,exceptionDescription); 00080 return; 00081 } 00082 00083 if(!fIsInitialised) 00084 { 00085 fIsInitialised = true; 00086 fParticleChangeForGamma = GetParticleChangeForGamma(); 00087 } 00088 00089 // Initialize water density pointer 00090 fpWaterDensity = G4DNAMolecularMaterial::Instance()->GetNumMolPerVolTableFor(G4Material::GetMaterial("G4_WATER")); 00091 }
G4ThreeVector G4DNASancheSolvatationModel::RadialDistributionOfProducts | ( | G4double | Rrms | ) | const [protected] |
Definition at line 124 of file G4DNASancheSolvatationModel.cc.
References G4UniformRand, G4INCL::Math::pi, and G4INCL::Math::sign().
Referenced by SampleSecondaries().
00125 { 00126 G4double sigma = std::sqrt(1.57)/2*expectationValue; 00127 00128 G4double XValueForfMax = std::sqrt(2.*sigma*sigma); 00129 G4double fMaxValue = std::sqrt(2./3.14) * 1./(sigma*sigma*sigma) * 00130 (XValueForfMax*XValueForfMax)* 00131 std::exp(-1./2. * (XValueForfMax*XValueForfMax)/(sigma*sigma)); 00132 00133 G4double R; 00134 00135 do 00136 { 00137 G4double aRandomfValue = fMaxValue*G4UniformRand(); 00138 00139 G4double sign; 00140 if(G4UniformRand() > 0.5) 00141 { 00142 sign = +1.; 00143 } 00144 else 00145 { 00146 sign = -1; 00147 } 00148 00149 R = expectationValue + sign*3.*sigma* G4UniformRand(); 00150 G4double f = std::sqrt(2./3.14) * 1/std::pow(sigma, 3) * R*R * std::exp(-1./2. * R*R/(sigma*sigma)); 00151 00152 if(aRandomfValue < f) 00153 { 00154 break; 00155 } 00156 } 00157 while(1); 00158 00159 G4double costheta = (2.*G4UniformRand()-1.); 00160 G4double theta = std::acos (costheta); 00161 G4double phi = 2.*pi*G4UniformRand(); 00162 00163 G4double xDirection = R*std::cos(phi)* std::sin(theta); 00164 G4double yDirection = R*std::sin(theta)*std::sin(phi); 00165 G4double zDirection = R*costheta; 00166 G4ThreeVector RandDirection = G4ThreeVector(xDirection, yDirection, zDirection); 00167 00168 return RandDirection; 00169 }
void G4DNASancheSolvatationModel::SampleSecondaries | ( | std::vector< G4DynamicParticle * > * | , | |
const G4MaterialCutsCouple * | , | |||
const G4DynamicParticle * | , | |||
G4double | tmin, | |||
G4double | maxEnergy | |||
) | [virtual] |
Implements G4VEmModel.
Definition at line 172 of file G4DNASancheSolvatationModel.cc.
References G4DNAChemistryManager::CreateSolvatedElectron(), fParticleChangeForGamma, fStopAndKill, fVerboseLevel, G4cout, G4endl, G4ParticleChangeForGamma::GetCurrentTrack(), G4DynamicParticle::GetKineticEnergy(), G4Track::GetPosition(), G4VEmModel::HighEnergyLimit(), G4DNAChemistryManager::Instance(), G4VParticleChange::ProposeLocalEnergyDeposit(), G4VParticleChange::ProposeTrackStatus(), RadialDistributionOfProducts(), and G4ParticleChangeForGamma::SetProposedKineticEnergy().
00177 { 00178 #ifdef G4VERBOSE 00179 if (fVerboseLevel) 00180 G4cout << "Calling SampleSecondaries() of G4SancheSolvatationModel" << G4endl; 00181 #endif 00182 G4double k = particle->GetKineticEnergy(); 00183 00184 if (k <= HighEnergyLimit()) 00185 { 00186 G4double r_mean = 00187 (-0.003*std::pow(k/eV,6) + 0.0749*std::pow(k/eV,5) - 0.7197*std::pow(k/eV,4) 00188 + 3.1384*std::pow(k/eV,3) - 5.6926*std::pow(k/eV,2) + 5.6237*k/eV - 0.7883)*nanometer; 00189 00190 G4ThreeVector displacement = RadialDistributionOfProducts (r_mean); 00191 00192 //______________________________________________________________ 00193 const G4Track * theIncomingTrack = fParticleChangeForGamma->GetCurrentTrack(); 00194 G4ThreeVector finalPosition(theIncomingTrack->GetPosition()+displacement); 00195 00196 G4DNAChemistryManager::Instance()->CreateSolvatedElectron(theIncomingTrack,&finalPosition); 00197 00198 fParticleChangeForGamma->SetProposedKineticEnergy(25.e-3*eV); 00199 fParticleChangeForGamma->ProposeTrackStatus(fStopAndKill); 00200 fParticleChangeForGamma->ProposeLocalEnergyDeposit(k); 00201 } 00202 }
void G4DNASancheSolvatationModel::SetVerbose | ( | int | ) | [inline] |
Definition at line 89 of file G4DNASancheSolvatationModel.hh.
References fVerboseLevel.
00090 { 00091 fVerboseLevel = flag; 00092 }
G4bool G4DNASancheSolvatationModel::fIsInitialised [protected] |
Definition at line 79 of file G4DNASancheSolvatationModel.hh.
Referenced by G4DNASancheSolvatationModel(), Initialise(), and SampleSecondaries().
const std::vector<G4double>* G4DNASancheSolvatationModel::fpWaterDensity [protected] |
Definition at line 76 of file G4DNASancheSolvatationModel.hh.
Referenced by G4DNASancheSolvatationModel(), and Initialise().
G4int G4DNASancheSolvatationModel::fVerboseLevel [protected] |
Definition at line 82 of file G4DNASancheSolvatationModel.hh.
Referenced by CrossSectionPerVolume(), G4DNASancheSolvatationModel(), Initialise(), SampleSecondaries(), and SetVerbose().