Ignore:
Timestamp:
Sep 30, 2010, 2:47:17 PM (14 years ago)
Author:
garnier
Message:

tag geant4.9.4 beta 1 + modifs locales

File:
1 edited

Legend:

Unmodified
Added
Removed

  • trunk/examples/extended/medical/DICOM/src/DicomPhysicsList.cc

    r1230 r1337  
    2424// ********************************************************************
    2525//
    26 // The code was written by :
    27 //      *Louis Archambault louis.archambault@phy.ulaval.ca,
    28 //      *Luc Beaulieu beaulieu@phy.ulaval.ca
    29 //      +Vincent Hubert-Tremblay at tigre.2@sympatico.ca
    30 //
    31 //
    32 // *Centre Hospitalier Universitaire de Quebec (CHUQ),
    33 // Hotel-Dieu de Quebec, departement de Radio-oncologie
    34 // 11 cote du palais. Quebec, QC, Canada, G1R 2J6
    35 // tel (418) 525-4444 #6720
    36 // fax (418) 691 5268
    37 //
    38 // + Université Laval, Québec (QC) Canada
    39 //
    40 // History: 30.11.07  P.Arce default cut changed to 1 mm
    41 //*******************************************************
    42 
    43 #include "DicomPhysicsList.hh"
     26// -------------------------------------------------------------------
     27// $Id: DicomPhysicsList.cc,v 1.11 2009/10/26 11:20:49 chauvie Exp $
     28// -------------------------------------------------------------------
    4429
    4530#include "G4ParticleDefinition.hh"
    46 #include "G4ParticleWithCuts.hh"
    4731#include "G4ProcessManager.hh"
    4832#include "G4ParticleTypes.hh"
    49 #include "G4ParticleTable.hh"
    50 #include "G4Material.hh"
    51 #include "G4UnitsTable.hh"
    52 #include "G4ios.hh"
     33#include "G4StepLimiter.hh"
     34#include "G4BaryonConstructor.hh"                     
     35#include "G4IonConstructor.hh"   
     36#include "G4MesonConstructor.hh"         
     37
     38#include "DicomPhysicsList.hh"
     39
     40//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
    5341
    5442DicomPhysicsList::DicomPhysicsList():  G4VUserPhysicsList()
    5543{
    56   defaultCutValue = 1.*mm;
    57   cutForGamma     = 1.*mm;
     44  defaultCutValue = 0.01*micrometer;
     45  cutForGamma     = defaultCutValue;
    5846  cutForElectron  = defaultCutValue;
    5947  cutForPositron  = defaultCutValue;
    60 
    61   SetVerboseLevel(0);
    62 }
     48  SetVerboseLevel(1);
     49}
     50
     51//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
    6352
    6453DicomPhysicsList::~DicomPhysicsList()
    6554{}
    6655
     56//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
     57
    6758void DicomPhysicsList::ConstructParticle()
    6859{
    69   // In this method, static member functions should be called
    70   // for all particles which you want to use.
    71   // This ensures that objects of these particle types will be
    72   // created in the program.
    73 
    7460  ConstructBosons();
    7561  ConstructLeptons();
    76 }
     62  ConstructBaryons();
     63}
     64
     65//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
    7766
    7867void DicomPhysicsList::ConstructBosons()
    79 {
     68{ 
    8069  // gamma
    8170  G4Gamma::GammaDefinition();
    8271
    83 }
     72  // optical photon
     73  G4OpticalPhoton::OpticalPhotonDefinition();
     74}
     75 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
    8476
    8577void DicomPhysicsList::ConstructLeptons()
    … …  
    9082}
    9183
     84//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
     85
     86void DicomPhysicsList::ConstructBaryons()
     87{
     88  //  baryons
     89  G4BaryonConstructor bConstructor;
     90  bConstructor.ConstructParticle();
     91
     92  G4IonConstructor iConstructor;
     93  iConstructor.ConstructParticle();
     94
     95  G4MesonConstructor mConstructor;
     96  mConstructor.ConstructParticle();
     97}
     98
     99//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
     100
    92101void DicomPhysicsList::ConstructProcess()
    93102{
    94103  AddTransportation();
    95104  ConstructEM();
    96 }
    97 
    98 #include "G4MultipleScattering.hh"
     105  ConstructHad();
     106  ConstructGeneral();
     107}
     108
     109// *** Processes and models
     110
    99111// gamma
    100 #include "G4LowEnergyRayleigh.hh"
    101 #include "G4LowEnergyPhotoElectric.hh"
    102 #include "G4LowEnergyCompton.hh"
    103 #include "G4LowEnergyGammaConversion.hh"
     112
     113#include "G4PhotoElectricEffect.hh"
     114#include "G4LivermorePhotoElectricModel.hh"
     115
     116#include "G4ComptonScattering.hh"
     117#include "G4LivermoreComptonModel.hh"
     118
     119#include "G4GammaConversion.hh"
     120#include "G4LivermoreGammaConversionModel.hh"
     121
     122#include "G4RayleighScattering.hh"
     123#include "G4LivermoreRayleighModel.hh"
     124
    104125// e-
    105 #include "G4LowEnergyIonisation.hh"
    106 #include "G4LowEnergyBremsstrahlung.hh"
     126
     127#include "G4eMultipleScattering.hh"
     128#include "G4UniversalFluctuation.hh"
     129
     130#include "G4eIonisation.hh"
     131#include "G4LivermoreIonisationModel.hh"
     132
     133#include "G4eBremsstrahlung.hh"
     134#include "G4LivermoreBremsstrahlungModel.hh"
     135
    107136// e+
    108 #include "G4eIonisation.hh"
    109 #include "G4eBremsstrahlung.hh"
     137
    110138#include "G4eplusAnnihilation.hh"
    111139
     140// mu
     141
     142#include "G4MuIonisation.hh"
     143#include "G4MuBremsstrahlung.hh"
     144#include "G4MuPairProduction.hh"
     145
     146// hadrons
     147
     148#include "G4hMultipleScattering.hh"
     149#include "G4MscStepLimitType.hh"
     150
     151#include "G4hBremsstrahlung.hh"
     152#include "G4hPairProduction.hh"
     153
     154#include "G4hIonisation.hh"
     155#include "G4ionIonisation.hh"
     156#include "G4IonParametrisedLossModel.hh"
     157
     158//
     159
     160#include "G4LossTableManager.hh"
     161#include "G4EmProcessOptions.hh"
     162
     163
     164//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
     165
    112166void DicomPhysicsList::ConstructEM()
    113167{
     168  theParticleIterator->reset();
     169
     170  while( (*theParticleIterator)() ){
     171
     172    G4ParticleDefinition* particle = theParticleIterator->value();
     173    G4ProcessManager* pmanager = particle->GetProcessManager();
     174    G4String particleName = particle->GetParticleName();
     175
     176    if (particleName == "gamma") {
     177
     178
     179      G4PhotoElectricEffect* thePhotoElectricEffect = new G4PhotoElectricEffect();
     180      G4LivermorePhotoElectricModel* theLivermorePhotoElectricModel = new G4LivermorePhotoElectricModel();
     181      thePhotoElectricEffect->AddEmModel(0, theLivermorePhotoElectricModel);
     182      pmanager->AddDiscreteProcess(thePhotoElectricEffect);
     183
     184      G4ComptonScattering* theComptonScattering = new G4ComptonScattering();
     185      G4LivermoreComptonModel* theLivermoreComptonModel = new G4LivermoreComptonModel();
     186      theComptonScattering->AddEmModel(0, theLivermoreComptonModel);
     187      pmanager->AddDiscreteProcess(theComptonScattering);
     188
     189      G4GammaConversion* theGammaConversion = new G4GammaConversion();
     190      G4LivermoreGammaConversionModel* theLivermoreGammaConversionModel = new G4LivermoreGammaConversionModel();
     191      theGammaConversion->AddEmModel(0, theLivermoreGammaConversionModel);
     192      pmanager->AddDiscreteProcess(theGammaConversion);
     193
     194      G4RayleighScattering* theRayleigh = new G4RayleighScattering();
     195      G4LivermoreRayleighModel* theRayleighModel = new G4LivermoreRayleighModel();
     196      theRayleigh->AddEmModel(0, theRayleighModel);
     197      pmanager->AddDiscreteProcess(theRayleigh);
     198     
     199      pmanager->AddProcess(new G4StepLimiter(), -1, -1, 5);
     200     
     201    } else if (particleName == "e-") {
     202     
     203      G4eMultipleScattering* msc = new G4eMultipleScattering();
     204      msc->SetStepLimitType(fUseDistanceToBoundary);
     205      pmanager->AddProcess(msc,                   -1, 1, 1);
     206     
     207      // Ionisation
     208      G4eIonisation* eIoni = new G4eIonisation();
     209      eIoni->AddEmModel(0, new G4LivermoreIonisationModel(), new G4UniversalFluctuation() );
     210      eIoni->SetStepFunction(0.2, 100*um); //     
     211      pmanager->AddProcess(eIoni,                 -1, 2, 2);
     212     
     213      // Bremsstrahlung
     214      G4eBremsstrahlung* eBrem = new G4eBremsstrahlung();
     215      eBrem->AddEmModel(0, new G4LivermoreBremsstrahlungModel());
     216      pmanager->AddProcess(eBrem,                 -1,-3, 3);
     217
     218      pmanager->AddProcess(new G4StepLimiter(), -1, -1, 4);
     219     
     220    } else if (particleName == "e+") {
     221
     222      // Identical to G4EmStandardPhysics_option3
     223   
     224      G4eMultipleScattering* msc = new G4eMultipleScattering();
     225      msc->SetStepLimitType(fUseDistanceToBoundary);
     226      pmanager->AddProcess(msc,                   -1, 1, 1);
     227
     228      G4eIonisation* eIoni = new G4eIonisation();
     229      eIoni->SetStepFunction(0.2, 100*um);     
     230      pmanager->AddProcess(eIoni,                 -1, 2, 2);
     231     
     232      pmanager->AddProcess(new G4eBremsstrahlung, -1,-3, 3);
     233     
     234      pmanager->AddProcess(new G4eplusAnnihilation,0,-1, 4);
     235     
     236      pmanager->AddProcess(new G4StepLimiter(), -1, -1, 5);
     237
     238    } else if (particleName == "GenericIon") {
     239
     240      pmanager->AddProcess(new G4hMultipleScattering, -1, 1, 1);
     241
     242      G4ionIonisation* ionIoni = new G4ionIonisation();
     243      ionIoni->SetEmModel(new G4IonParametrisedLossModel());
     244      ionIoni->SetStepFunction(0.1, 20*um);
     245      pmanager->AddProcess(ionIoni,                   -1, 2, 2);
     246
     247      pmanager->AddProcess(new G4StepLimiter(), -1, -1, 3);
     248
     249    } else if (particleName == "alpha" ||
     250               particleName == "He3" ) {
     251
     252      // Identical to G4EmStandardPhysics_option3
     253     
     254      pmanager->AddProcess(new G4hMultipleScattering, -1, 1, 1);
     255
     256      G4ionIonisation* ionIoni = new G4ionIonisation();
     257      ionIoni->SetStepFunction(0.1, 20*um);
     258      pmanager->AddProcess(ionIoni,                   -1, 2, 2);
     259
     260      pmanager->AddProcess(new G4StepLimiter(), -1, -1, 3);
     261    }
     262
     263   //
     264
     265  }
     266}
     267
     268//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
     269
     270#include "G4HadronElasticProcess.hh"
     271#include "G4LElastic.hh"
     272
     273#include "G4AlphaInelasticProcess.hh"
     274#include "G4BinaryLightIonReaction.hh"
     275#include "G4TripathiCrossSection.hh"
     276#include "G4IonsShenCrossSection.hh"
     277#include "G4LEAlphaInelastic.hh"
     278
     279void DicomPhysicsList::ConstructHad()
     280{
     281
     282  G4HadronElasticProcess * theElasticProcess = new G4HadronElasticProcess;
     283  theElasticProcess->RegisterMe( new G4LElastic() );
     284
    114285  theParticleIterator->reset();
    115286  while( (*theParticleIterator)() )
    116     {
    117       G4ParticleDefinition* particle = theParticleIterator->value();
    118       G4ProcessManager* pmanager = particle->GetProcessManager();
    119       G4String particleName = particle->GetParticleName();
    120 
    121       //processes
    122       lowePhot = new  G4LowEnergyPhotoElectric("LowEnPhotoElec");
    123       loweIon  = new G4LowEnergyIonisation("LowEnergyIoni");
    124       loweBrem = new G4LowEnergyBremsstrahlung("LowEnBrem");
    125 
    126       if (particleName == "gamma")
    127         {
    128           //gamma
    129           pmanager->AddDiscreteProcess(new G4LowEnergyRayleigh);
    130           pmanager->AddDiscreteProcess(lowePhot);
    131           pmanager->AddDiscreteProcess(new G4LowEnergyCompton);
    132           pmanager->AddDiscreteProcess(new G4LowEnergyGammaConversion);
    133 
    134         }
    135       else if (particleName == "e-")
    136         {
    137           //electron
    138           pmanager->AddProcess(new G4MultipleScattering, -1, 1,1);
    139           pmanager->AddProcess(loweIon,     -1, 2,2);
    140           pmanager->AddProcess(loweBrem,    -1,-1,3);
    141 
    142         }
    143       else if (particleName == "e+")
    144         {
    145           //positron
    146           pmanager->AddProcess(new G4MultipleScattering, -1, 1,1);
    147           pmanager->AddProcess(new G4eIonisation,        -1, 2,2);
    148           pmanager->AddProcess(new G4eBremsstrahlung,    -1,-1,3);
    149           pmanager->AddProcess(new G4eplusAnnihilation,   0,-1,4);
    150 
    151         }
    152     }
    153 }
    154 #include "G4Decay.hh"
     287  {
     288    G4ParticleDefinition* particle = theParticleIterator->value();
     289    G4ProcessManager* pManager = particle->GetProcessManager();
     290
     291    if (particle->GetParticleName() == "alpha")
     292       {
     293
     294          // INELASTIC SCATTERING
     295          // Binary Cascade
     296          G4BinaryLightIonReaction* theBC = new G4BinaryLightIonReaction();
     297          theBC -> SetMinEnergy(80.*MeV);
     298          theBC -> SetMaxEnergy(40.*GeV);
     299 
     300          // TRIPATHI CROSS SECTION
     301          // Implementation of formulas in analogy to NASA technical paper 3621 by
     302          // Tripathi, et al. Cross-sections for ion ion scattering
     303          G4TripathiCrossSection* TripathiCrossSection = new G4TripathiCrossSection;
     304 
     305          // IONS SHEN CROSS SECTION
     306          // Implementation of formulas
     307          // Shen et al. Nuc. Phys. A 491 130 (1989)
     308          // Total Reaction Cross Section for Heavy-Ion Collisions
     309          G4IonsShenCrossSection* aShen = new G4IonsShenCrossSection;
     310 
     311          // Final state production model for Alpha inelastic scattering below 20 GeV
     312          G4LEAlphaInelastic* theAIModel = new G4LEAlphaInelastic;
     313          theAIModel -> SetMaxEnergy(100.*MeV);
     314
     315          G4AlphaInelasticProcess * theIPalpha = new G4AlphaInelasticProcess;             
     316          theIPalpha->AddDataSet(TripathiCrossSection);
     317          theIPalpha->AddDataSet(aShen);
     318
     319          // Register the Alpha Inelastic and Binary Cascade Model
     320          theIPalpha->RegisterMe(theAIModel);
     321          theIPalpha->RegisterMe(theBC);
     322         
     323          // Activate the alpha inelastic scattering using the alpha inelastic and binary cascade model
     324          pManager -> AddDiscreteProcess(theIPalpha);
     325         
     326          // Activate the Hadron Elastic Process
     327          pManager -> AddDiscreteProcess(theElasticProcess);
     328           
     329       }
     330  }
     331
     332}
     333
     334//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
     335
    155336void DicomPhysicsList::ConstructGeneral()
    156 {
    157   // Add Decay Process
    158   G4Decay* theDecayProcess = new G4Decay();
    159   theParticleIterator->reset();
    160   while( (*theParticleIterator)() )
    161     {
    162       G4ParticleDefinition* particle = theParticleIterator->value();
    163       G4ProcessManager* pmanager = particle->GetProcessManager();
    164       if (theDecayProcess->IsApplicable(*particle))
    165         {
    166           pmanager ->AddProcess(theDecayProcess);
    167           // set ordering for PostStepDoIt and AtRestDoIt
    168           pmanager ->SetProcessOrdering(theDecayProcess, idxPostStep);
    169           pmanager ->SetProcessOrdering(theDecayProcess, idxAtRest);
    170         }
    171     }
    172 }
     337{ }
     338
     339//....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
    173340
    174341void DicomPhysicsList::SetCuts()
    175342{
    176   if (verboseLevel >0)
    177     {
    178       G4cout << "DicomPhysicsList::SetCuts:";
    179       G4cout << "CutLength : " << G4BestUnit(defaultCutValue,"Length") << G4endl;
    180     }
    181 
     343  if (verboseLevel >0){
     344    G4cout << "DicomPhysicsList::SetCuts:";
     345    G4cout << "CutLength : " << G4BestUnit(defaultCutValue,"Length") << G4endl;
     346  } 
     347 
    182348  // set cut values for gamma at first and for e- second and next for e+,
    183   // because some processes for e+/e- need cut values for gamma
     349  // because some processes for e+/e- need cut values for gamma 
    184350  SetCutValue(cutForGamma, "gamma");
    185351  SetCutValue(cutForElectron, "e-");
    186352  SetCutValue(cutForPositron, "e+");
    187 
    188 
    189   if (verboseLevel>0)
    190     DumpCutValuesTable();
    191 }
    192 
    193 void DicomPhysicsList::SetGammaCut(G4double val)
    194 {
    195   ResetCuts();
    196   cutForGamma = val;
    197 }
    198 
    199 void DicomPhysicsList::SetElectronCut(G4double val)
    200 {
    201   //  ResetCuts();
    202   cutForElectron = val;
    203 }
    204 
    205 void DicomPhysicsList::SetPositronCut(G4double val)
    206 {
    207   //  ResetCuts();
    208   cutForPositron = val;
    209 }
    210 
     353 
     354  if (verboseLevel>0) DumpCutValuesTable();
     355
     356}
     357
Note: See TracChangeset for help on using the changeset viewer.