00001 // 00002 // ******************************************************************** 00003 // * License and Disclaimer * 00004 // * * 00005 // * The Geant4 software is copyright of the Copyright Holders of * 00006 // * the Geant4 Collaboration. It is provided under the terms and * 00007 // * conditions of the Geant4 Software License, included in the file * 00008 // * LICENSE and available at http://cern.ch/geant4/license . These * 00009 // * include a list of copyright holders. * 00010 // * * 00011 // * Neither the authors of this software system, nor their employing * 00012 // * institutes,nor the agencies providing financial support for this * 00013 // * work make any representation or warranty, express or implied, * 00014 // * regarding this software system or assume any liability for its * 00015 // * use. Please see the license in the file LICENSE and URL above * 00016 // * for the full disclaimer and the limitation of liability. * 00017 // * * 00018 // * This code implementation is the result of the scientific and * 00019 // * technical work of the GEANT4 collaboration. * 00020 // * By using, copying, modifying or distributing the software (or * 00021 // * any work based on the software) you agree to acknowledge its * 00022 // * use in resulting scientific publications, and indicate your * 00023 // * acceptance of all terms of the Geant4 Software license. * 00024 // ******************************************************************** 00025 // 00026 // $Id: G4WentzelVIRelModel.hh 66592 2012-12-23 09:34:55Z vnivanch $ 00027 // 00028 // ------------------------------------------------------------------- 00029 // 00030 // 00031 // GEANT4 Class header file 00032 // 00033 // 00034 // File name: G4WentzelVIRelModel 00035 // 00036 // Author: V.Ivanchenko 00037 // 00038 // Creation date: 08.06.2012 from G4WentzelVIModel 00039 // 00040 // Modifications: 00041 // 00042 // Class Description: 00043 // 00044 // Implementation of the model of multiple scattering based on 00045 // G.Wentzel, Z. Phys. 40 (1927) 590. 00046 // H.W.Lewis, Phys Rev 78 (1950) 526. 00047 // J.M. Fernandez-Varea et al., NIM B73 (1993) 447. 00048 // L.Urban, CERN-OPEN-2006-077. 00049 00050 // ------------------------------------------------------------------- 00051 // 00052 00053 #ifndef G4WentzelVIRelModel_h 00054 #define G4WentzelVIRelModel_h 1 00055 00056 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 00057 00058 #include "G4VMscModel.hh" 00059 #include "G4MaterialCutsCouple.hh" 00060 #include "G4WentzelVIRelXSection.hh" 00061 00062 class G4ParticleDefinition; 00063 class G4LossTableManager; 00064 class G4NistManager; 00065 class G4Pow; 00066 00067 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 00068 00069 class G4WentzelVIRelModel : public G4VMscModel 00070 { 00071 00072 public: 00073 00074 G4WentzelVIRelModel(const G4String& nam = "WentzelVIUni"); 00075 00076 virtual ~G4WentzelVIRelModel(); 00077 00078 virtual void Initialise(const G4ParticleDefinition*, const G4DataVector&); 00079 00080 void StartTracking(G4Track*); 00081 00082 virtual G4double ComputeCrossSectionPerAtom(const G4ParticleDefinition*, 00083 G4double KineticEnergy, 00084 G4double AtomicNumber, 00085 G4double AtomicWeight=0., 00086 G4double cut = DBL_MAX, 00087 G4double emax= DBL_MAX); 00088 00089 virtual G4ThreeVector& SampleScattering(const G4ThreeVector&, 00090 G4double safety); 00091 00092 virtual G4double ComputeTruePathLengthLimit(const G4Track& track, 00093 G4double& currentMinimalStep); 00094 00095 virtual G4double ComputeGeomPathLength(G4double truePathLength); 00096 00097 virtual G4double ComputeTrueStepLength(G4double geomStepLength); 00098 00099 private: 00100 00101 G4double ComputeXSectionPerVolume(); 00102 00103 inline void SetupParticle(const G4ParticleDefinition*); 00104 00105 inline void DefineMaterial(const G4MaterialCutsCouple*); 00106 00107 // hide assignment operator 00108 G4WentzelVIRelModel & operator=(const G4WentzelVIRelModel &right); 00109 G4WentzelVIRelModel(const G4WentzelVIRelModel&); 00110 00111 G4LossTableManager* theManager; 00112 G4NistManager* fNistManager; 00113 G4ParticleChangeForMSC* fParticleChange; 00114 G4WentzelVIRelXSection* wokvi; 00115 G4Pow* fG4pow; 00116 00117 const G4DataVector* currentCuts; 00118 00119 G4double tlimitminfix; 00120 G4double invsqrt12; 00121 00122 // cache kinematics 00123 G4double preKinEnergy; 00124 G4double tPathLength; 00125 G4double zPathLength; 00126 G4double lambdaeff; 00127 G4double currentRange; 00128 00129 // data for single scattering mode 00130 G4double xtsec; 00131 std::vector<G4double> xsecn; 00132 std::vector<G4double> prob; 00133 G4int nelments; 00134 00135 G4double numlimit; 00136 00137 // cache material 00138 G4int currentMaterialIndex; 00139 const G4MaterialCutsCouple* currentCouple; 00140 const G4Material* currentMaterial; 00141 00142 // single scattering parameters 00143 G4double cosThetaMin; 00144 G4double cosThetaMax; 00145 G4double cosTetMaxNuc; 00146 00147 // projectile 00148 const G4ParticleDefinition* particle; 00149 G4double lowEnergyLimit; 00150 00151 // flags 00152 G4bool inside; 00153 G4bool singleScatteringMode; 00154 }; 00155 00156 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 00157 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 00158 00159 inline 00160 void G4WentzelVIRelModel::DefineMaterial(const G4MaterialCutsCouple* cup) 00161 { 00162 if(cup != currentCouple) { 00163 currentCouple = cup; 00164 SetCurrentCouple(cup); 00165 currentMaterial = cup->GetMaterial(); 00166 currentMaterialIndex = currentCouple->GetIndex(); 00167 } 00168 } 00169 00170 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo.... 00171 00172 inline void G4WentzelVIRelModel::SetupParticle(const G4ParticleDefinition* p) 00173 { 00174 // Initialise mass and charge 00175 if(p != particle) { 00176 particle = p; 00177 wokvi->SetupParticle(p); 00178 } 00179 } 00180 00181 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo...... 00182 00183 #endif 00184