G4PhotoElectricAngularGeneratorSauterGavrila.cc

Go to the documentation of this file.
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 //
00027 // -------------------------------------------------------------------
00028 //
00029 // GEANT4 Class file
00030 //
00031 //
00032 // File name:     G4PhotoElectricAngularGeneratorSauterGavrila
00033 //
00034 // Creation date: 10 May 2004
00035 //
00036 // Modifications: 
00037 // 10 May 2003     P. Rodrigues    First implementation acording with new design
00038 //
00039 // Class Description: 
00040 //
00041 // Concrete class for PhotoElectric Electron Angular Distribution Generation 
00042 // This model is a re-implementation of the Photolectric angular distribution
00043 // developed my M. Maire for the Standard EM Physics G4PhotoElectricEffect 
00044 //
00045 // Class Description: End 
00046 //
00047 // -------------------------------------------------------------------
00048 //
00049 
00050 #include "G4PhotoElectricAngularGeneratorSauterGavrila.hh"
00051 #include "G4PhysicalConstants.hh"
00052 #include "Randomize.hh"
00053 
00054 G4PhotoElectricAngularGeneratorSauterGavrila::G4PhotoElectricAngularGeneratorSauterGavrila():
00055   G4VEmAngularDistribution("AngularGenSauterGavrilaLowE")
00056 {}
00057 
00058 G4PhotoElectricAngularGeneratorSauterGavrila::~G4PhotoElectricAngularGeneratorSauterGavrila() 
00059 {}
00060 
00061 G4ThreeVector& 
00062 G4PhotoElectricAngularGeneratorSauterGavrila::SampleDirection(
00063                          const G4DynamicParticle* dp,
00064                          G4double, G4int, const G4Material*)
00065 {
00066 
00067   // Compute Theta distribution of the emitted electron, with respect to the
00068   // incident Gamma.
00069   // The Sauter-Gavrila distribution for the K-shell is used. 
00070 
00071   G4double costeta = 1.;
00072   G4double Phi     = twopi * G4UniformRand();
00073   G4double cosphi = std::cos(Phi);
00074   G4double sinphi = std::sin(Phi);
00075   G4double sinteta = 0;
00076   G4double gamma   = 1. + dp->GetKineticEnergy()/electron_mass_c2;
00077 
00078   if (gamma > 5.) {
00079     fLocalDirection = dp->GetMomentumDirection(); 
00080     return fLocalDirection; 
00081     // Bugzilla 1120
00082     // SI on 05/09/2010 as suggested by JG 04/09/10 
00083   }
00084 
00085   G4double beta  = std::sqrt((gamma - 1)*(gamma + 1))/gamma;
00086   G4double b     = 0.5*gamma*(gamma - 1)*(gamma - 2);
00087     
00088   G4double rndm,term,greject,grejsup;
00089   if (gamma < 2.) grejsup = gamma*gamma*(1.+b-beta*b);
00090   else            grejsup = gamma*gamma*(1.+b+beta*b);
00091   
00092   do { rndm = 1.-2*G4UniformRand();
00093        costeta = (rndm+beta)/(rndm*beta+1.);
00094        term = 1.-beta*costeta;
00095        greject = (1.-costeta*costeta)*(1.+b*term)/(term*term);
00096   } while(greject < G4UniformRand()*grejsup);
00097        
00098   sinteta = std::sqrt((1 - costeta)*(1 + costeta));
00099   fLocalDirection.set(sinteta*cosphi, sinteta*sinphi, costeta);
00100   fLocalDirection.rotateUz(dp->GetMomentumDirection());
00101   return fLocalDirection;
00102 }
00103 
00104 void G4PhotoElectricAngularGeneratorSauterGavrila::PrintGeneratorInformation() const
00105 {
00106   G4cout << "\n" << G4endl;
00107   G4cout << "" << G4endl;
00108   G4cout << "Re-implementation of the photolectric angular distribution" << G4endl;
00109   G4cout << "developed my M. Maire for the Standard EM Physics G4PhotoElectricEffect" << G4endl;
00110   G4cout << "It computes the theta distribution of the emitted electron, with respect to the" << G4endl;
00111   G4cout << "incident Gamma, using the Sauter-Gavrila distribution for the K-shell\n" << G4endl;
00112 } 

Generated on Mon May 27 17:49:18 2013 for Geant4 by  doxygen 1.4.7