Geant4 Cross Reference

Cross-Referencing   Geant4
Geant4/processes/electromagnetic/dna/models/include/G4DNABornIonisationModel1.hh

Version: [ ReleaseNotes ] [ 1.0 ] [ 1.1 ] [ 2.0 ] [ 3.0 ] [ 3.1 ] [ 3.2 ] [ 4.0 ] [ 4.0.p1 ] [ 4.0.p2 ] [ 4.1 ] [ 4.1.p1 ] [ 5.0 ] [ 5.0.p1 ] [ 5.1 ] [ 5.1.p1 ] [ 5.2 ] [ 5.2.p1 ] [ 5.2.p2 ] [ 6.0 ] [ 6.0.p1 ] [ 6.1 ] [ 6.2 ] [ 6.2.p1 ] [ 6.2.p2 ] [ 7.0 ] [ 7.0.p1 ] [ 7.1 ] [ 7.1.p1 ] [ 8.0 ] [ 8.0.p1 ] [ 8.1 ] [ 8.1.p1 ] [ 8.1.p2 ] [ 8.2 ] [ 8.2.p1 ] [ 8.3 ] [ 8.3.p1 ] [ 8.3.p2 ] [ 9.0 ] [ 9.0.p1 ] [ 9.0.p2 ] [ 9.1 ] [ 9.1.p1 ] [ 9.1.p2 ] [ 9.1.p3 ] [ 9.2 ] [ 9.2.p1 ] [ 9.2.p2 ] [ 9.2.p3 ] [ 9.2.p4 ] [ 9.3 ] [ 9.3.p1 ] [ 9.3.p2 ] [ 9.4 ] [ 9.4.p1 ] [ 9.4.p2 ] [ 9.4.p3 ] [ 9.4.p4 ] [ 9.5 ] [ 9.5.p1 ] [ 9.5.p2 ] [ 9.6 ] [ 9.6.p1 ] [ 9.6.p2 ] [ 9.6.p3 ] [ 9.6.p4 ] [ 10.0 ] [ 10.0.p1 ] [ 10.0.p2 ] [ 10.0.p3 ] [ 10.0.p4 ] [ 10.1 ] [ 10.1.p1 ] [ 10.1.p2 ] [ 10.1.p3 ] [ 10.2 ] [ 10.2.p1 ] [ 10.2.p2 ] [ 10.2.p3 ] [ 10.3 ] [ 10.3.p1 ] [ 10.3.p2 ] [ 10.3.p3 ] [ 10.4 ] [ 10.4.p1 ] [ 10.4.p2 ] [ 10.4.p3 ] [ 10.5 ] [ 10.5.p1 ] [ 10.6 ] [ 10.6.p1 ] [ 10.6.p2 ] [ 10.6.p3 ] [ 10.7 ] [ 10.7.p1 ] [ 10.7.p2 ] [ 10.7.p3 ] [ 10.7.p4 ] [ 11.0 ] [ 11.0.p1 ] [ 11.0.p2 ] [ 11.0.p3, ] [ 11.0.p4 ] [ 11.1 ] [ 11.1.1 ] [ 11.1.2 ] [ 11.1.3 ] [ 11.2 ] [ 11.2.1 ] [ 11.2.2 ] [ 11.3.0 ]

  1 //
  2 // ********************************************************************
  3 // * License and Disclaimer                                           *
  4 // *                                                                  *
  5 // * The  Geant4 software  is  copyright of the Copyright Holders  of *
  6 // * the Geant4 Collaboration.  It is provided  under  the terms  and *
  7 // * conditions of the Geant4 Software License,  included in the file *
  8 // * LICENSE and available at  http://cern.ch/geant4/license .  These *
  9 // * include a list of copyright holders.                             *
 10 // *                                                                  *
 11 // * Neither the authors of this software system, nor their employing *
 12 // * institutes,nor the agencies providing financial support for this *
 13 // * work  make  any representation or  warranty, express or implied, *
 14 // * regarding  this  software system or assume any liability for its *
 15 // * use.  Please see the license in the file  LICENSE  and URL above *
 16 // * for the full disclaimer and the limitation of liability.         *
 17 // *                                                                  *
 18 // * This  code  implementation is the result of  the  scientific and *
 19 // * technical work of the GEANT4 collaboration.                      *
 20 // * By using,  copying,  modifying or  distributing the software (or *
 21 // * any work based  on the software)  you  agree  to acknowledge its *
 22 // * use  in  resulting  scientific  publications,  and indicate your *
 23 // * acceptance of all terms of the Geant4 Software license.          *
 24 // ********************************************************************
 25 //
 26 //
 27 
 28 #ifndef G4DNABornIonisationModel1_h
 29 #define G4DNABornIonisationModel1_h 1
 30 
 31 #include "G4VEmModel.hh"
 32 #include "G4ParticleChangeForGamma.hh"
 33 #include "G4ProductionCutsTable.hh"
 34 
 35 #include "G4DNACrossSectionDataSet.hh"
 36 #include "G4Electron.hh"
 37 #include "G4Proton.hh"
 38 #include "G4DNAGenericIonsManager.hh"
 39 
 40 #include "G4LogLogInterpolation.hh"
 41 
 42 #include "G4DNAWaterIonisationStructure.hh"
 43 #include "G4VAtomDeexcitation.hh"
 44 #include "G4NistManager.hh"
 45 
 46 
 47 class G4DNABornIonisationModel1 : public G4VEmModel
 48 {
 49 
 50 public:
 51 
 52   G4DNABornIonisationModel1(const G4ParticleDefinition* p = nullptr,
 53                const G4String& nam = "DNABornIonisationModel");
 54 
 55   ~G4DNABornIonisationModel1() override;
 56    
 57   G4DNABornIonisationModel1 & operator=(const  G4DNABornIonisationModel1 &right) = delete;
 58   G4DNABornIonisationModel1(const  G4DNABornIonisationModel1&) = delete;
 59 
 60   void Initialise(const G4ParticleDefinition*, const G4DataVector& = *(new G4DataVector())) override;
 61 
 62   G4double CrossSectionPerVolume(  const G4Material* material,
 63              const G4ParticleDefinition* p,
 64              G4double ekin,
 65              G4double emin,
 66              G4double emax) override;
 67 
 68   void SampleSecondaries(std::vector<G4DynamicParticle*>*,
 69          const G4MaterialCutsCouple*,
 70          const G4DynamicParticle*,
 71          G4double tmin,
 72          G4double maxEnergy) override;
 73 
 74   G4double GetPartialCrossSection(const G4Material*,
 75                                           G4int /*level*/,
 76                                           const G4ParticleDefinition*,
 77                                           G4double /*kineticEnergy*/) override;
 78         
 79   G4double DifferentialCrossSection(G4ParticleDefinition * aParticleDefinition, G4double k, G4double energyTransfer, G4int shell);
 80 
 81   G4double TransferedEnergy(G4ParticleDefinition * aParticleDefinition,
 82                             G4double incomingParticleEnergy, G4int shell, G4double random) ;
 83 
 84   inline void SelectFasterComputation(G4bool input); 
 85 
 86   inline void SelectStationary(G4bool input); 
 87 
 88   inline void SelectSPScaling(G4bool input); 
 89 
 90 protected:
 91 
 92   G4ParticleChangeForGamma* fParticleChangeForGamma;
 93 
 94 private:
 95 
 96   G4bool fasterCode;
 97   G4bool statCode;
 98   G4bool spScaling;
 99 
100   // Water density table
101   const std::vector<G4double>* fpMolWaterDensity;
102 
103   // Deexcitation manager to produce fluo photons and e-
104   G4VAtomDeexcitation*      fAtomDeexcitation;
105 
106   std::map<G4String,G4double,std::less<G4String> > lowEnergyLimit;
107   std::map<G4String,G4double,std::less<G4String> > highEnergyLimit;
108 
109   // TODO :
110 //  std::map<const G4ParticleDefinition*,std::pair<G4double,G4double> > fEnergyLimits;
111 
112 
113   G4bool isInitialised{false};
114   G4int verboseLevel;
115   
116   // Cross section
117 
118   using MapFile = std::map<G4String, G4String, std::less<G4String>>;
119   MapFile tableFile; // useful ?
120 
121   using MapData = std::map<G4String, G4DNACrossSectionDataSet *, std::less<G4String>>;
122   MapData tableData;
123   
124   // Final state
125   
126   G4DNAWaterIonisationStructure waterStructure;
127 
128   G4double RandomizeEjectedElectronEnergy(G4ParticleDefinition * aParticleDefinition, G4double incomingParticleEnergy, G4int shell) ;
129 
130   G4double RandomizeEjectedElectronEnergyFromCumulatedDcs(G4ParticleDefinition * aParticleDefinition, G4double incomingParticleEnergy, G4int shell) ;
131    
132   G4double Interpolate(G4double e1, G4double e2, G4double e, G4double xs1, G4double xs2);
133    
134   G4double QuadInterpolator( G4double e11, 
135            G4double e12, 
136            G4double e21, 
137            G4double e22, 
138            G4double x11,
139            G4double x12, 
140            G4double x21, 
141            G4double x22, 
142            G4double t1, 
143            G4double t2, 
144            G4double t, 
145            G4double e);
146 
147   using TriDimensionMap = std::map<G4double, std::map<G4double, G4double>>;
148   
149   TriDimensionMap eDiffCrossSectionData[6];
150   TriDimensionMap eNrjTransfData[6]; // for cumulated dcs
151   
152   TriDimensionMap pDiffCrossSectionData[6];
153   TriDimensionMap pNrjTransfData[6]; // for cumulated dcs
154   
155   std::vector<G4double> eTdummyVec;
156   std::vector<G4double> pTdummyVec;
157 
158   using VecMap = std::map<G4double, std::vector<G4double>>;
159   
160   VecMap eVecm;
161   VecMap pVecm;
162   
163   VecMap eProbaShellMap[6]; // for cumulated dcs
164   VecMap pProbaShellMap[6]; // for cumulated dcs
165   
166   // Partial cross section
167   
168   G4int RandomSelect(G4double energy,const G4String& particle );
169 
170 };
171 
172 inline void G4DNABornIonisationModel1::SelectFasterComputation (G4bool input)
173 { 
174     fasterCode = input; 
175 }    
176 
177 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
178 
179 inline void G4DNABornIonisationModel1::SelectStationary (G4bool input)
180 { 
181     statCode = input; 
182 }    
183 
184 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
185 
186 inline void G4DNABornIonisationModel1::SelectSPScaling (G4bool input)
187 { 
188     spScaling = input; 
189 }    
190 
191 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo....
192 
193 #endif
194