G4VLowEnergyDiscretePhotonProcess.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 // $Id$
00027 //
00028 // --------------------------------------------------------------
00029 //
00030 // File name:     G4VLowEnergyDiscretePhotonProcess.cc
00031 //
00032 // Author:        Capra Riccardo
00033 //
00034 // Creation date: May 2005
00035 //
00036 // History:
00037 // -----------
00038 // 02 May 2005  R. Capra         1st implementation
00039 //
00040 //----------------------------------------------------------------
00041 
00042   
00043 
00044 #include "G4VLowEnergyDiscretePhotonProcess.hh"
00045 
00046 #include "G4PhysicalConstants.hh"
00047 #include "G4SystemOfUnits.hh"
00048 #include "G4String.hh"
00049 #include "G4CrossSectionHandler.hh"
00050 #include "G4CompositeEMDataSet.hh"
00051 #include "G4Gamma.hh"
00052 #include "G4Track.hh"
00053 #include "G4DynamicParticle.hh"
00054 #include "G4ThreeVector.hh"
00055 #include "Randomize.hh" // G4UniformRand
00056 
00057 G4VLowEnergyDiscretePhotonProcess :: G4VLowEnergyDiscretePhotonProcess(const G4String& processName, 
00058                                                                        const G4String& aCrossSectionFileName, 
00059                                                                        const G4String& aScatterFileName, 
00060                                                                        G4VDataSetAlgorithm* aScatterInterpolation, 
00061                                                                        G4double aLowEnergyLimit, 
00062                                                                        G4double aHighEnergyLimit)
00063   :
00064   G4VLowEnergyTestableDiscreteProcess(processName),
00065   lowEnergyLimit(aLowEnergyLimit),
00066   highEnergyLimit(aHighEnergyLimit),
00067   crossSectionFileName(aCrossSectionFileName),
00068   meanFreePathTable(0)
00069 {
00070   crossSectionHandler = new G4CrossSectionHandler();
00071   scatterFunctionData = new G4CompositeEMDataSet(aScatterInterpolation, 1., 1.);
00072   scatterFunctionData->LoadData(aScatterFileName);
00073  
00074   if (verboseLevel > 0)
00075     {
00076       G4cout << GetProcessName() << " is created " << G4endl
00077              << "Energy range: "
00078              << lowEnergyLimit / keV << " keV - "
00079              << highEnergyLimit / GeV << " GeV"
00080              << G4endl;
00081     }
00082 }
00083 
00084 
00085 
00086 G4VLowEnergyDiscretePhotonProcess::~G4VLowEnergyDiscretePhotonProcess(void)
00087 {
00088   if (meanFreePathTable)
00089     delete meanFreePathTable;
00090  
00091   delete crossSectionHandler;
00092   delete scatterFunctionData;
00093 }
00094 
00095 
00096 
00097 
00098 
00099 G4bool G4VLowEnergyDiscretePhotonProcess::IsApplicable(const G4ParticleDefinition& particleDefinition)
00100 {
00101   return (&particleDefinition)==G4Gamma::Gamma(); 
00102 }
00103 
00104 
00105 
00106 void G4VLowEnergyDiscretePhotonProcess::BuildPhysicsTable(const G4ParticleDefinition&  /*photon*/)
00107 {
00108   crossSectionHandler->Clear();
00109   crossSectionHandler->LoadData(crossSectionFileName);
00110 
00111   if (meanFreePathTable)
00112     delete meanFreePathTable; 
00113   meanFreePathTable=crossSectionHandler->BuildMeanFreePathForMaterials();
00114 }
00115 
00116 
00117 
00118 
00119 
00120 G4double G4VLowEnergyDiscretePhotonProcess::GetMeanFreePath(const G4Track& aTrack, G4double /*previousStepSize*/, G4ForceCondition*  /*condition*/)
00121 {
00122   G4double photonEnergy;
00123   photonEnergy = aTrack.GetDynamicParticle()->GetKineticEnergy();
00124 
00125   if (photonEnergy < lowEnergyLimit)
00126     return DBL_MAX;
00127  
00128   if (photonEnergy > highEnergyLimit)
00129     photonEnergy=highEnergyLimit;
00130  
00131   size_t materialIndex;
00132   materialIndex = aTrack.GetMaterialCutsCouple()->GetIndex(); 
00133 
00134   return meanFreePathTable->FindValue(photonEnergy, materialIndex);
00135 }
00136 
00137 
00138 
00139 
00140 
00141 G4ThreeVector G4VLowEnergyDiscretePhotonProcess::GetPhotonPolarization(const G4DynamicParticle&  photon)
00142 {
00143   G4ThreeVector photonMomentumDirection;
00144   G4ThreeVector photonPolarization;
00145 
00146   photonPolarization = photon.GetPolarization(); 
00147   photonMomentumDirection = photon.GetMomentumDirection();
00148 
00149   if ((!photonPolarization.isOrthogonal(photonMomentumDirection, 1e-6)) || photonPolarization.mag()==0.)
00150     {
00151       // if |photonPolarization|==0. or |photonPolarization * photonDirection0| > 1e-6 * |photonPolarization ^ photonDirection0|
00152       // then polarization is choosen randomly.
00153   
00154       G4ThreeVector e1(photonMomentumDirection.orthogonal().unit());
00155       G4ThreeVector e2(photonMomentumDirection.cross(e1).unit());
00156   
00157       G4double angle(G4UniformRand() * twopi);
00158   
00159       e1*=std::cos(angle);
00160       e2*=std::sin(angle);
00161   
00162       photonPolarization=e1+e2;
00163     }
00164   else if (photonPolarization.howOrthogonal(photonMomentumDirection) != 0.)
00165     {
00166       // if |photonPolarization * photonDirection0| != 0.
00167       // then polarization is made orthonormal;
00168   
00169       photonPolarization=photonPolarization.perpPart(photonMomentumDirection);
00170     }
00171  
00172   return photonPolarization.unit();
00173 }

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