G4ChargeExchangeProcess.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 // $Id$
00028 //
00029 //
00030 // Geant4 Hadron Charge Exchange Process -- source file
00031 //
00032 // Created 21 April 2006 V.Ivanchenko
00033 //
00034 // Modified:
00035 // 24-Apr-06 V.Ivanchenko add neutron scattering on hydrogen from CHIPS
00036 // 07-Jun-06 V.Ivanchenko fix problem of rotation of final state
00037 // 25-Jul-06 V.Ivanchenko add 19 MeV low energy for CHIPS
00038 // 23-Jan-07 V.Ivanchenko add cross section interfaces with Z and A
00039 //                        and do not use CHIPS for cross sections
00040 // 14-Sep-12 M.Kelsey -- Pass subType code to base ctor
00041 
00042 #include "G4ChargeExchangeProcess.hh"
00043 #include "globals.hh"
00044 #include "G4SystemOfUnits.hh"
00045 #include "G4CrossSectionDataStore.hh"
00046 #include "G4HadronElasticDataSet.hh"
00047 #include "G4Element.hh"
00048 #include "G4ElementVector.hh"
00049 #include "G4IsotopeVector.hh"
00050 #include "G4Neutron.hh"
00051 #include "G4Proton.hh"
00052 #include "G4PhysicsLinearVector.hh"
00053 
00054 G4ChargeExchangeProcess::G4ChargeExchangeProcess(const G4String& procName)
00055   : G4HadronicProcess(procName,fChargeExchange), first(true)
00056 {
00057   thEnergy = 20.*MeV;
00058   pPDG = 0;
00059   verboseLevel= 1;
00060   AddDataSet(new G4HadronElasticDataSet);
00061   theProton   = G4Proton::Proton();
00062   theNeutron  = G4Neutron::Neutron();
00063   theAProton  = G4AntiProton::AntiProton();
00064   theANeutron = G4AntiNeutron::AntiNeutron();
00065   thePiPlus   = G4PionPlus::PionPlus();
00066   thePiMinus  = G4PionMinus::PionMinus();
00067   thePiZero   = G4PionZero::PionZero();
00068   theKPlus    = G4KaonPlus::KaonPlus();
00069   theKMinus   = G4KaonMinus::KaonMinus();
00070   theK0S      = G4KaonZeroShort::KaonZeroShort();
00071   theK0L      = G4KaonZeroLong::KaonZeroLong();
00072   theL        = G4Lambda::Lambda();
00073   theAntiL    = G4AntiLambda::AntiLambda();
00074   theSPlus    = G4SigmaPlus::SigmaPlus();
00075   theASPlus   = G4AntiSigmaPlus::AntiSigmaPlus();
00076   theSMinus   = G4SigmaMinus::SigmaMinus();
00077   theASMinus  = G4AntiSigmaMinus::AntiSigmaMinus();
00078   theS0       = G4SigmaZero::SigmaZero();
00079   theAS0      = G4AntiSigmaZero::AntiSigmaZero();
00080   theXiMinus  = G4XiMinus::XiMinus();
00081   theXi0      = G4XiZero::XiZero();
00082   theAXiMinus = G4AntiXiMinus::AntiXiMinus();
00083   theAXi0     = G4AntiXiZero::AntiXiZero();
00084   theOmega    = G4OmegaMinus::OmegaMinus();
00085   theAOmega   = G4AntiOmegaMinus::AntiOmegaMinus();
00086   theD        = G4Deuteron::Deuteron();
00087   theT        = G4Triton::Triton();
00088   theA        = G4Alpha::Alpha();
00089   theHe3      = G4He3::He3();
00090 }
00091 
00092 G4ChargeExchangeProcess::~G4ChargeExchangeProcess()
00093 {
00094   delete factors;
00095 }
00096 
00097 void G4ChargeExchangeProcess::
00098 BuildPhysicsTable(const G4ParticleDefinition& aParticleType)
00099 {
00100   if(first) {
00101     first = false;
00102     theParticle = &aParticleType;
00103     pPDG = theParticle->GetPDGEncoding();
00104 
00105     store = G4HadronicProcess::GetCrossSectionDataStore();
00106 
00107     const size_t n = 10;
00108     if(theParticle == thePiPlus || theParticle == thePiMinus ||
00109        theParticle == theKPlus  || theParticle == theKMinus ||
00110        theParticle == theK0S    || theParticle == theK0L) {
00111 
00112       G4double F[n] = {0.33,0.27,0.29,0.31,0.27,0.18,0.13,0.1,0.09,0.07};
00113       factors = new G4PhysicsLinearVector(0.0,2.0*GeV,n);
00114       for(size_t i=0; i<n; i++) {factors->PutValue(i,F[i]);}
00115 
00116     } else {
00117 
00118       G4double F[n] = {0.50,0.45,0.40,0.35,0.30,0.25,0.06,0.04,0.005,0.0};
00119       factors = new G4PhysicsLinearVector(0.0,4.0*GeV,n);
00120       for(size_t i=0; i<n; i++) {factors->PutValue(i,F[i]);}
00121     }
00122     //factors->SetSpline(true);
00123 
00124     if(verboseLevel>1)
00125       G4cout << "G4ChargeExchangeProcess for "
00126              << theParticle->GetParticleName()
00127              << G4endl;
00128   }
00129   G4HadronicProcess::BuildPhysicsTable(aParticleType);
00130 }
00131 
00132 G4double G4ChargeExchangeProcess::GetElementCrossSection(
00133                                   const G4DynamicParticle* dp,
00134                                   const G4Element* elm,
00135                                   const G4Material* mat)
00136 {
00137   // gives the microscopic cross section in GEANT4 internal units
00138   G4double Z = elm->GetZ();
00139   G4int iz = G4int(Z);
00140   G4double x = 0.0;
00141 
00142   // The process is effective only above the threshold
00143   if(iz == 1 || dp->GetKineticEnergy() < thEnergy) return x;
00144 
00145   if(verboseLevel>1)
00146     G4cout << "G4ChargeExchangeProcess compute GHAD CS for element "
00147            << elm->GetName()
00148            << G4endl;
00149   x = store->GetCrossSection(dp, elm, mat);
00150 
00151   if(verboseLevel>1)
00152     G4cout << "G4ChargeExchangeProcess cross(mb)= " << x/millibarn
00153            << "  E(MeV)= " << dp->GetKineticEnergy()
00154            << "  " << theParticle->GetParticleName()
00155            << "  in Z= " << iz
00156            << G4endl;
00157   G4bool b;
00158   G4double A = elm->GetN();
00159   G4double ptot = dp->GetTotalMomentum();
00160   x *= factors->GetValue(ptot, b)/std::pow(A, 0.42);
00161   if(theParticle == thePiPlus || theParticle == theProton ||
00162      theParticle == theKPlus  || theParticle == theANeutron)
00163     { x *= (1.0 - Z/A); }
00164 
00165   else if(theParticle == thePiMinus || theParticle == theNeutron ||
00166           theParticle == theKMinus  || theParticle == theAProton)
00167     { x *= Z/A; }
00168 
00169   if(theParticle->GetPDGMass() < GeV) {
00170     if(ptot > 2.*GeV) x *= 4.0*GeV*GeV/(ptot*ptot);
00171   }
00172 
00173   if(verboseLevel>1) 
00174     G4cout << "Corrected cross(mb)= " << x/millibarn << G4endl;
00175 
00176   return x;
00177 }
00178 
00179 G4bool G4ChargeExchangeProcess::
00180 IsApplicable(const G4ParticleDefinition& aParticleType)
00181 {
00182   const G4ParticleDefinition* p = &aParticleType;
00183   return (p == thePiPlus || p == thePiMinus ||
00184           p == theProton || p == theNeutron ||
00185           p == theAProton|| p == theANeutron||
00186           p == theKPlus  || p == theKMinus  ||
00187           p == theK0S    || p == theK0L     ||
00188           p == theL);
00189 }
00190 
00191 void G4ChargeExchangeProcess::
00192 DumpPhysicsTable(const G4ParticleDefinition& aParticleType)
00193 {
00194   store->DumpPhysicsTable(aParticleType);
00195 }

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