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 #include "G4EmDNAPhysics.hh"
00029
00030 #include "G4SystemOfUnits.hh"
00031
00032 #include "G4DNAGenericIonsManager.hh"
00033
00034
00035
00036 #include "G4DNAElastic.hh"
00037 #include "G4DNAChampionElasticModel.hh"
00038 #include "G4DNAScreenedRutherfordElasticModel.hh"
00039
00040 #include "G4DNAExcitation.hh"
00041 #include "G4DNAAttachment.hh"
00042 #include "G4DNAVibExcitation.hh"
00043 #include "G4DNAIonisation.hh"
00044 #include "G4DNAChargeDecrease.hh"
00045 #include "G4DNAChargeIncrease.hh"
00046
00047
00048
00049 #include "G4Electron.hh"
00050 #include "G4Proton.hh"
00051 #include "G4GenericIon.hh"
00052
00053
00054
00055 #include "G4Positron.hh"
00056 #include "G4eMultipleScattering.hh"
00057 #include "G4UrbanMscModel95.hh"
00058 #include "G4eIonisation.hh"
00059 #include "G4eBremsstrahlung.hh"
00060 #include "G4eplusAnnihilation.hh"
00061
00062 #include "G4Gamma.hh"
00063 #include "G4PhotoElectricEffect.hh"
00064 #include "G4LivermorePhotoElectricModel.hh"
00065 #include "G4ComptonScattering.hh"
00066 #include "G4LivermoreComptonModel.hh"
00067 #include "G4GammaConversion.hh"
00068 #include "G4LivermoreGammaConversionModel.hh"
00069 #include "G4RayleighScattering.hh"
00070 #include "G4LivermoreRayleighModel.hh"
00071
00072
00073 #include "G4LossTableManager.hh"
00074 #include "G4UAtomicDeexcitation.hh"
00075 #include "G4PhysicsListHelper.hh"
00076 #include "G4BuilderType.hh"
00077
00078
00079 #include "G4PhysicsConstructorFactory.hh"
00080
00081 G4_DECLARE_PHYSCONSTR_FACTORY(G4EmDNAPhysics);
00082
00083
00084
00085
00086 G4EmDNAPhysics::G4EmDNAPhysics(G4int ver)
00087 : G4VPhysicsConstructor("G4EmDNAPhysics"), verbose(ver)
00088 {
00089 SetPhysicsType(bElectromagnetic);
00090 }
00091
00092
00093
00094 G4EmDNAPhysics::G4EmDNAPhysics(G4int ver, const G4String&)
00095 : G4VPhysicsConstructor("G4EmDNAPhysics"), verbose(ver)
00096 {
00097 SetPhysicsType(bElectromagnetic);
00098 }
00099
00100
00101
00102 G4EmDNAPhysics::~G4EmDNAPhysics()
00103 {}
00104
00105
00106
00107 void G4EmDNAPhysics::ConstructParticle()
00108 {
00109
00110 G4Gamma::Gamma();
00111
00112
00113 G4Electron::Electron();
00114 G4Positron::Positron();
00115
00116
00117 G4Proton::Proton();
00118
00119 G4GenericIon::GenericIonDefinition();
00120
00121 G4DNAGenericIonsManager * genericIonsManager;
00122 genericIonsManager=G4DNAGenericIonsManager::Instance();
00123 genericIonsManager->GetIon("alpha++");
00124 genericIonsManager->GetIon("alpha+");
00125 genericIonsManager->GetIon("helium");
00126 genericIonsManager->GetIon("hydrogen");
00127 genericIonsManager->GetIon("carbon");
00128 genericIonsManager->GetIon("nitrogen");
00129 genericIonsManager->GetIon("oxygen");
00130 genericIonsManager->GetIon("iron");
00131
00132 }
00133
00134
00135
00136 void G4EmDNAPhysics::ConstructProcess()
00137 {
00138 G4PhysicsListHelper* ph = G4PhysicsListHelper::GetPhysicsListHelper();
00139
00140 theParticleIterator->reset();
00141 while( (*theParticleIterator)() )
00142 {
00143 G4ParticleDefinition* particle = theParticleIterator->value();
00144 G4String particleName = particle->GetParticleName();
00145
00146 if (particleName == "e-") {
00147
00148
00149
00150 G4DNAElastic* theDNAElasticProcess = new G4DNAElastic("e-_G4DNAElastic");
00151 theDNAElasticProcess->SetEmModel(new G4DNAChampionElasticModel());
00152
00153
00154
00155
00156 ph->RegisterProcess(theDNAElasticProcess, particle);
00157
00158
00159 ph->RegisterProcess(new G4DNAExcitation("e-_G4DNAExcitation"), particle);
00160
00161
00162 ph->RegisterProcess(new G4DNAIonisation("e-_G4DNAIonisation"), particle);
00163
00164
00165 ph->RegisterProcess(new G4DNAVibExcitation("e-_G4DNAVibExcitation"), particle);
00166
00167
00168 ph->RegisterProcess(new G4DNAAttachment("e-_G4DNAAttachment"), particle);
00169
00170 } else if ( particleName == "proton" ) {
00171 ph->RegisterProcess(new G4DNAExcitation("proton_G4DNAExcitation"), particle);
00172 ph->RegisterProcess(new G4DNAIonisation("proton_G4DNAIonisation"), particle);
00173 ph->RegisterProcess(new G4DNAChargeDecrease("proton_G4DNAChargeDecrease"), particle);
00174
00175 } else if ( particleName == "hydrogen" ) {
00176 ph->RegisterProcess(new G4DNAExcitation("hydrogen_G4DNAExcitation"), particle);
00177 ph->RegisterProcess(new G4DNAIonisation("hydrogen_G4DNAIonisation"), particle);
00178 ph->RegisterProcess(new G4DNAChargeIncrease("hydrogen_G4DNAChargeIncrease"), particle);
00179
00180 } else if ( particleName == "alpha" ) {
00181 ph->RegisterProcess(new G4DNAExcitation("alpha_G4DNAExcitation"), particle);
00182 ph->RegisterProcess(new G4DNAIonisation("alpha_G4DNAIonisation"), particle);
00183 ph->RegisterProcess(new G4DNAChargeDecrease("alpha_G4DNAChargeDecrease"), particle);
00184
00185 } else if ( particleName == "alpha+" ) {
00186 ph->RegisterProcess(new G4DNAExcitation("alpha+_G4DNAExcitation"), particle);
00187 ph->RegisterProcess(new G4DNAIonisation("alpha+_G4DNAIonisation"), particle);
00188 ph->RegisterProcess(new G4DNAChargeDecrease("alpha+_G4DNAChargeDecrease"), particle);
00189 ph->RegisterProcess(new G4DNAChargeIncrease("alpha+_G4DNAChargeIncrease"), particle);
00190
00191 } else if ( particleName == "helium" ) {
00192 ph->RegisterProcess(new G4DNAExcitation("helium_G4DNAExcitation"), particle);
00193 ph->RegisterProcess(new G4DNAIonisation("helium_G4DNAIonisation"), particle);
00194 ph->RegisterProcess(new G4DNAChargeIncrease("helium_G4DNAChargeIncrease"), particle);
00195
00196
00197
00198 } else if ( particleName == "carbon" ) {
00199 ph->RegisterProcess(new G4DNAIonisation("carbon_G4DNAIonisation"), particle);
00200
00201 } else if ( particleName == "nitrogen" ) {
00202 ph->RegisterProcess(new G4DNAIonisation("nitrogen_G4DNAIonisation"), particle);
00203
00204 } else if ( particleName == "oxygen" ) {
00205 ph->RegisterProcess(new G4DNAIonisation("oxygen_G4DNAIonisation"), particle);
00206
00207 } else if ( particleName == "iron" ) {
00208 ph->RegisterProcess(new G4DNAIonisation("iron_G4DNAIonisation"), particle);
00209
00210 }
00211
00212
00213
00214
00215
00216
00217
00218 else if (particleName == "e+") {
00219
00220
00221
00222 G4eMultipleScattering* msc = new G4eMultipleScattering();
00223 msc->AddEmModel(0, new G4UrbanMscModel95());
00224 msc->SetStepLimitType(fUseDistanceToBoundary);
00225 G4eIonisation* eIoni = new G4eIonisation();
00226 eIoni->SetStepFunction(0.2, 100*um);
00227
00228 ph->RegisterProcess(msc, particle);
00229 ph->RegisterProcess(eIoni, particle);
00230 ph->RegisterProcess(new G4eBremsstrahlung(), particle);
00231 ph->RegisterProcess(new G4eplusAnnihilation(), particle);
00232
00233 } else if (particleName == "gamma") {
00234
00235 G4double LivermoreHighEnergyLimit = GeV;
00236
00237 G4PhotoElectricEffect* thePhotoElectricEffect = new G4PhotoElectricEffect();
00238 G4LivermorePhotoElectricModel* theLivermorePhotoElectricModel =
00239 new G4LivermorePhotoElectricModel();
00240 theLivermorePhotoElectricModel->SetHighEnergyLimit(LivermoreHighEnergyLimit);
00241 thePhotoElectricEffect->AddEmModel(0, theLivermorePhotoElectricModel);
00242 ph->RegisterProcess(thePhotoElectricEffect, particle);
00243
00244 G4ComptonScattering* theComptonScattering = new G4ComptonScattering();
00245 G4LivermoreComptonModel* theLivermoreComptonModel =
00246 new G4LivermoreComptonModel();
00247 theLivermoreComptonModel->SetHighEnergyLimit(LivermoreHighEnergyLimit);
00248 theComptonScattering->AddEmModel(0, theLivermoreComptonModel);
00249 ph->RegisterProcess(theComptonScattering, particle);
00250
00251 G4GammaConversion* theGammaConversion = new G4GammaConversion();
00252 G4LivermoreGammaConversionModel* theLivermoreGammaConversionModel =
00253 new G4LivermoreGammaConversionModel();
00254 theLivermoreGammaConversionModel->SetHighEnergyLimit(LivermoreHighEnergyLimit);
00255 theGammaConversion->AddEmModel(0, theLivermoreGammaConversionModel);
00256 ph->RegisterProcess(theGammaConversion, particle);
00257
00258 G4RayleighScattering* theRayleigh = new G4RayleighScattering();
00259 G4LivermoreRayleighModel* theRayleighModel = new G4LivermoreRayleighModel();
00260 theRayleighModel->SetHighEnergyLimit(LivermoreHighEnergyLimit);
00261 theRayleigh->AddEmModel(0, theRayleighModel);
00262 ph->RegisterProcess(theRayleigh, particle);
00263 }
00264
00265
00266
00267 }
00268
00269
00270
00271 G4VAtomDeexcitation* de = new G4UAtomicDeexcitation();
00272 G4LossTableManager::Instance()->SetAtomDeexcitation(de);
00273 de->SetFluo(true);
00274 }
00275
00276