00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031 #include "G4DNAMeltonAttachmentModel.hh"
00032 #include "G4SystemOfUnits.hh"
00033 #include "G4DNAChemistryManager.hh"
00034 #include "G4DNAMolecularMaterial.hh"
00035
00036
00037
00038 using namespace std;
00039
00040
00041
00042 G4DNAMeltonAttachmentModel::G4DNAMeltonAttachmentModel(const G4ParticleDefinition*,
00043 const G4String& nam)
00044 :G4VEmModel(nam),isInitialised(false)
00045 {
00046
00047 fpWaterDensity = 0;
00048
00049 lowEnergyLimit = 4 * eV;
00050 lowEnergyLimitOfModel = 4 * eV;
00051 highEnergyLimit = 13 * eV;
00052 SetLowEnergyLimit(lowEnergyLimit);
00053 SetHighEnergyLimit(highEnergyLimit);
00054
00055 verboseLevel= 0;
00056
00057
00058
00059
00060
00061
00062
00063 if( verboseLevel>0 )
00064 {
00065 G4cout << "Melton Attachment model is constructed " << G4endl
00066 << "Energy range: "
00067 << lowEnergyLimit / eV << " eV - "
00068 << highEnergyLimit / eV << " eV"
00069 << G4endl;
00070 }
00071 fParticleChangeForGamma = 0;
00072 fDissociationFlag = true;
00073 }
00074
00075
00076
00077 G4DNAMeltonAttachmentModel::~G4DNAMeltonAttachmentModel()
00078 {
00079
00080
00081 std::map< G4String,G4DNACrossSectionDataSet*,std::less<G4String> >::iterator pos;
00082
00083 for (pos = tableData.begin(); pos != tableData.end(); ++pos)
00084 {
00085 G4DNACrossSectionDataSet* table = pos->second;
00086 delete table;
00087 }
00088
00089
00090
00091 }
00092
00093
00094
00095 void G4DNAMeltonAttachmentModel::Initialise(const G4ParticleDefinition* ,
00096 const G4DataVector& )
00097 {
00098
00099 if (verboseLevel > 3)
00100 G4cout << "Calling G4DNAMeltonAttachmentModel::Initialise()" << G4endl;
00101
00102
00103
00104 if (LowEnergyLimit() < lowEnergyLimit)
00105 {
00106 G4cout << "G4DNAMeltonAttachmentModel: low energy limit increased from " <<
00107 LowEnergyLimit()/eV << " eV to " << lowEnergyLimit/eV << " eV" << G4endl;
00108 SetLowEnergyLimit(lowEnergyLimit);
00109 }
00110
00111 if (HighEnergyLimit() > highEnergyLimit)
00112 {
00113 G4cout << "G4DNAMeltonAttachmentModel: high energy limit decreased from " <<
00114 HighEnergyLimit()/eV << " eV to " << highEnergyLimit/eV << " eV" << G4endl;
00115 SetHighEnergyLimit(highEnergyLimit);
00116 }
00117
00118
00119
00120 G4double scaleFactor = 1e-18*cm*cm;
00121
00122 G4String fileElectron("dna/sigma_attachment_e_melton");
00123
00124 G4ParticleDefinition* electronDef = G4Electron::ElectronDefinition();
00125 G4String electron;
00126
00127
00128
00129
00130
00131 electron = electronDef->GetParticleName();
00132
00133 tableFile[electron] = fileElectron;
00134
00135 G4DNACrossSectionDataSet* tableE = new G4DNACrossSectionDataSet(new G4LogLogInterpolation, eV,scaleFactor );
00136 tableE->LoadData(fileElectron);
00137 tableData[electron] = tableE;
00138
00139
00140
00141 if (verboseLevel > 2)
00142 G4cout << "Loaded cross section data for Melton Attachment model" << G4endl;
00143
00144 if( verboseLevel>0 )
00145 {
00146 G4cout << "Melton Attachment model is initialized " << G4endl
00147 << "Energy range: "
00148 << LowEnergyLimit() / eV << " eV - "
00149 << HighEnergyLimit() / eV << " eV"
00150 << G4endl;
00151 }
00152
00153 fpWaterDensity = G4DNAMolecularMaterial::Instance()->GetNumMolPerVolTableFor(G4Material::GetMaterial("G4_WATER"));
00154
00155 if (isInitialised) { return; }
00156 fParticleChangeForGamma = GetParticleChangeForGamma();
00157 isInitialised = true;
00158
00159 }
00160
00161
00162
00163 G4double G4DNAMeltonAttachmentModel::CrossSectionPerVolume(const G4Material* material,
00164 const G4ParticleDefinition* particleDefinition,
00165 G4double ekin,
00166 G4double,
00167 G4double)
00168 {
00169 if (verboseLevel > 3)
00170 G4cout << "Calling CrossSectionPerVolume() of G4DNAMeltonAttachmentModel" << G4endl;
00171
00172
00173
00174 G4double sigma=0;
00175
00176 G4double waterDensity = (*fpWaterDensity)[material->GetIndex()];
00177
00178 if(waterDensity!= 0.0)
00179
00180 {
00181 const G4String& particleName = particleDefinition->GetParticleName();
00182
00183 if (ekin >= lowEnergyLimit && ekin < highEnergyLimit)
00184 {
00185
00186 std::map< G4String,G4DNACrossSectionDataSet*,std::less<G4String> >::iterator pos;
00187 pos = tableData.find(particleName);
00188
00189 if (pos != tableData.end())
00190 {
00191 G4DNACrossSectionDataSet* table = pos->second;
00192 if (table != 0)
00193 {
00194 sigma = table->FindValue(ekin);
00195 }
00196 }
00197 else
00198 {
00199 G4Exception("G4DNAMeltonAttachmentModel::ComputeCrossSectionPerVolume","em0002",
00200 FatalException,"Model not applicable to particle type.");
00201 }
00202 }
00203
00204 if (verboseLevel > 2)
00205 {
00206 G4cout << "__________________________________" << G4endl;
00207 G4cout << "°°° G4DNAMeltonAttachmentModel - XS INFO START" << G4endl;
00208 G4cout << "°°° Kinetic energy(eV)=" << ekin/eV << " particle : " << particleDefinition->GetParticleName() << G4endl;
00209 G4cout << "°°° Cross section per water molecule (cm^2)=" << sigma/cm/cm << G4endl;
00210 G4cout << "°°° Cross section per water molecule (cm^-1)=" << sigma*waterDensity/(1./cm) << G4endl;
00211
00212 G4cout << "°°° G4DNAMeltonAttachmentModel - XS INFO END" << G4endl;
00213 }
00214
00215 }
00216
00217 return sigma*waterDensity;
00218
00219 }
00220
00221
00222
00223 void G4DNAMeltonAttachmentModel::SampleSecondaries(std::vector<G4DynamicParticle*>* ,
00224 const G4MaterialCutsCouple* ,
00225 const G4DynamicParticle* aDynamicElectron,
00226 G4double,
00227 G4double)
00228 {
00229
00230 if (verboseLevel > 3)
00231 G4cout << "Calling SampleSecondaries() of G4DNAMeltonAttachmentModel" << G4endl;
00232
00233
00234
00235 G4double electronEnergy0 = aDynamicElectron->GetKineticEnergy();
00236 fParticleChangeForGamma->SetProposedKineticEnergy(0.);
00237 fParticleChangeForGamma->ProposeTrackStatus(fStopAndKill);
00238 fParticleChangeForGamma->ProposeLocalEnergyDeposit(electronEnergy0);
00239
00240 if(fDissociationFlag)
00241 {
00242 G4DNAChemistryManager::Instance()->CreateWaterMolecule(eDissociativeAttachment,-1,
00243 fParticleChangeForGamma->GetCurrentTrack());
00244 }
00245 return ;
00246 }