include/Geant4/G4UrbanAdjointMscModel.hh

0001 //
0002 // ********************************************************************
0003 // * License and Disclaimer                                           *
0004 // *                                                                  *
0005 // * The  Geant4 software  is  copyright of the Copyright Holders  of *
0006 // * the Geant4 Collaboration.  It is provided  under  the terms  and *
0007 // * conditions of the Geant4 Software License,  included in the file *
0008 // * LICENSE and available at  http://cern.ch/geant4/license .  These *
0009 // * include a list of copyright holders.                             *
0010 // *                                                                  *
0011 // * Neither the authors of this software system, nor their employing *
0012 // * institutes,nor the agencies providing financial support for this *
0013 // * work  make  any representation or  warranty, express or implied, *
0014 // * regarding  this  software system or assume any liability for its *
0015 // * use.  Please see the license in the file  LICENSE  and URL above *
0016 // * for the full disclaimer and the limitation of liability.         *
0017 // *                                                                  *
0018 // * This  code  implementation is the result of  the  scientific and *
0019 // * technical work of the GEANT4 collaboration.                      *
0020 // * By using,  copying,  modifying or  distributing the software (or *
0021 // * any work based  on the software)  you  agree  to acknowledge its *
0022 // * use  in  resulting  scientific  publications,  and indicate your *
0023 // * acceptance of all terms of the Geant4 Software license.          *
0024 // ********************************************************************
0025 //
0026 // File name:     G4UrbanAdjointMscModel
0027 //
0028 // Author:        Laszlo Urban
0029 //
0030 // Class Description:
0031 //   Implementation of the model of multiple scattering based on
0032 //   H.W.Lewis Phys Rev 78 (1950) 526 and L.Urban model
0033 // -------------------------------------------------------------------
0034
0035 #ifndef G4UrbanAdjointMscModel_h
0036 #define G4UrbanAdjointMscModel_h 1
0037
0038 #include "G4Electron.hh"
0039 #include "G4Exp.hh"
0040 #include "G4Log.hh"
0041 #include "G4MscStepLimitType.hh"
0042 #include "G4VMscModel.hh"
0043
0044 class G4LossTableManager;
0045 class G4MaterialCutsCouple;
0046 class G4ParticleChangeForMSC;
0047 class G4ParticleDefinition;
0048 class G4SafetyHelper;
0049
0050 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0051
0052 class G4UrbanAdjointMscModel : public G4VMscModel
0053 {
0054  public:
0055   explicit G4UrbanAdjointMscModel(const G4String& nam = "UrbanMsc");
0056
0057   ~G4UrbanAdjointMscModel() override;
0058
0059   void Initialise(const G4ParticleDefinition*, const G4DataVector&) override;
0060
0061   void StartTracking(G4Track*) override;
0062
0063   G4double ComputeCrossSectionPerAtom(const G4ParticleDefinition* particle,
0064                                       G4double KineticEnergy,
0065                                       G4double AtomicNumber,
0066                                       G4double AtomicWeight = 0.,
0067                                       G4double cut          = 0.,
0068                                       G4double emax         = DBL_MAX) override;
0069
0070   G4ThreeVector& SampleScattering(const G4ThreeVector&,
0071                                   G4double safety) override;
0072
0073   G4double ComputeTruePathLengthLimit(const G4Track& track,
0074                                       G4double& currentMinimalStep) override;
0075
0076   G4double ComputeGeomPathLength(G4double truePathLength) override;
0077
0078   G4double ComputeTrueStepLength(G4double geomStepLength) override;
0079
0080   G4double ComputeTheta0(G4double truePathLength, G4double KineticEnergy);
0081
0082   inline void SetNewDisplacementFlag(G4bool);
0083
0084   G4UrbanAdjointMscModel& operator=(const G4UrbanAdjointMscModel& right) =
0085     delete;
0086   G4UrbanAdjointMscModel(const G4UrbanAdjointMscModel&) = delete;
0087
0088  private:
0089   G4double SampleCosineTheta(G4double trueStepLength, G4double KineticEnergy);
0090
0091   void SampleDisplacement(G4double sinTheta, G4double phi);
0092
0093   void SampleDisplacementNew(G4double sinTheta, G4double phi);
0094
0095   inline void SetParticle(const G4ParticleDefinition*);
0096
0097   inline void UpdateCache();
0098
0099   inline G4double Randomizetlimit();
0100
0101   inline G4double SimpleScattering(G4double xmeanth, G4double x2meanth);
0102
0103   CLHEP::HepRandomEngine* rndmEngineMod;
0104
0105   const G4ParticleDefinition* particle;
0106   const G4ParticleDefinition* positron;
0107   G4ParticleChangeForMSC* fParticleChange;
0108
0109   const G4MaterialCutsCouple* couple;
0110   G4LossTableManager* theManager;
0111
0112   G4double mass;
0113   G4double charge, ChargeSquare;
0114   G4double masslimite, lambdalimit, fr;
0115
0116   G4double taubig;
0117   G4double tausmall;
0118   G4double taulim;
0119   G4double currentTau;
0120   G4double tlimit;
0121   G4double tlimitmin;
0122   G4double tlimitminfix, tlimitminfix2;
0123   G4double tgeom;
0124
0125   G4double geombig;
0126   G4double geommin;
0127   G4double geomlimit;
0128   G4double skindepth;
0129   G4double smallstep;
0130
0131   G4double presafety;
0132
0133   G4double lambda0;
0134   G4double lambdaeff;
0135   G4double tPathLength;
0136   G4double zPathLength;
0137   G4double par1, par2, par3;
0138
0139   G4double stepmin;
0140
0141   G4double currentKinEnergy;
0142   G4double currentRange;
0143   G4double rangeinit;
0144   G4double currentRadLength;
0145
0146   G4double Zold;
0147   G4double Zeff, Z2, Z23, lnZ;
0148   G4double coeffth1, coeffth2;
0149   G4double coeffc1, coeffc2, coeffc3, coeffc4;
0150
0151   G4double rangecut;
0152   G4double drr, finalr;
0153
0154   G4int currentMaterialIndex;
0155
0156   G4bool firstStep;
0157   G4bool insideskin;
0158
0159   G4bool latDisplasmentbackup;
0160   G4bool displacementFlag;
0161 };
0162
0163 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0164 inline void G4UrbanAdjointMscModel::SetNewDisplacementFlag(G4bool val)
0165 {
0166   displacementFlag = val;
0167 }
0168
0169 inline void G4UrbanAdjointMscModel::SetParticle(const G4ParticleDefinition* p)
0170 {
0171   const G4ParticleDefinition* p1 = p;
0172
0173   if(p->GetParticleName() == "adj_e-")
0174     p1 = G4Electron::Electron();
0175
0176   if(p1 != particle)
0177   {
0178     particle     = p1;
0179     mass         = p1->GetPDGMass();
0180     charge       = p1->GetPDGCharge() / CLHEP::eplus;
0181     ChargeSquare = charge * charge;
0182   }
0183 }
0184
0185 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0186 inline G4double G4UrbanAdjointMscModel::Randomizetlimit()
0187 {
0188   G4double temptlimit = tlimit;
0189   if(tlimit > tlimitmin)
0190   {
0191     G4double delta = tlimit - tlimitmin;
0192     do
0193     {
0194       temptlimit = G4RandGauss::shoot(rndmEngineMod, tlimit, 0.1 * delta);
0195       // Loop checking, 10-Apr-2016, Laszlo Urban
0196     } while((temptlimit < tlimit - delta) || (temptlimit > tlimit + delta));
0197   }
0198   else
0199   {
0200     temptlimit = tlimitmin;
0201   }
0202
0203   return temptlimit;
0204 }
0205
0206 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0207 inline void G4UrbanAdjointMscModel::UpdateCache()
0208 {
0209   lnZ = G4Log(Zeff);
0210   // correction in theta0 formula
0211   G4double w    = G4Exp(lnZ / 6.);
0212   G4double facz = 0.990395 + w * (-0.168386 + w * 0.093286);
0213   coeffth1      = facz * (1. - 8.7780e-2 / Zeff);
0214   coeffth2      = facz * (4.0780e-2 + 1.7315e-4 * Zeff);
0215
0216   // tail parameters
0217   G4double Z13 = w * w;
0218   coeffc1      = 2.3785 - Z13 * (4.1981e-1 - Z13 * 6.3100e-2);
0219   coeffc2      = 4.7526e-1 + Z13 * (1.7694 - Z13 * 3.3885e-1);
0220   coeffc3      = 2.3683e-1 - Z13 * (1.8111 - Z13 * 3.2774e-1);
0221   coeffc4      = 1.7888e-2 + Z13 * (1.9659e-2 - Z13 * 2.6664e-3);
0222
0223   Z2  = Zeff * Zeff;
0224   Z23 = Z13 * Z13;
0225
0226   Zold = Zeff;
0227 }
0228
0229 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
0230 inline G4double G4UrbanAdjointMscModel::SimpleScattering(G4double xmeanth,
0231                                                          G4double x2meanth)
0232 {
0233   // 'large angle scattering'
0234   // 2 model functions with correct xmean and x2mean
0235   G4double a =
0236     (2. * xmeanth + 9. * x2meanth - 3.) / (2. * xmeanth - 3. * x2meanth + 1.);
0237   G4double prob = (a + 2.) * xmeanth / a;
0238
0239   // sampling
0240   G4double cth = 1.;
0241   if(rndmEngineMod->flat() < prob)
0242   {
0243     cth = -1. + 2. * G4Exp(G4Log(rndmEngineMod->flat()) / (a + 1.));
0244   }
0245   else
0246   {
0247     cth = -1. + 2. * rndmEngineMod->flat();
0248   }
0249   return cth;
0250 }
0251
0252 #endif