G4CHIPSElastic.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 // Geant4 Class : G4CHIPSElastic
00031 //
00032 // Author : V.Ivanchenko 29 June 2009 
00033 //  
00034 // Modified:
00035 // 13.01.10: M.Kosov: Use G4Q(Pr/Neut)ElasticCS instead of G4QElasticCS
00036 //
00037 //---------------------------------------------------------------------
00038 // CHIPS model of hadron elastic scattering
00039 //
00040 
00041 #include "G4CHIPSElastic.hh"
00042 #include "G4VQCrossSection.hh"
00043 #include "G4ParticleDefinition.hh"
00044 #include "G4QProtonElasticCrossSection.hh"
00045 #include "G4QNeutronElasticCrossSection.hh"
00046 
00047 #include "G4QAntiBaryonElasticCrossSection.hh"        // Uzhi
00048 #include "G4QPionPlusElasticCrossSection.hh"          // Uzhi
00049 #include "G4QPionMinusElasticCrossSection.hh"         // Uzhi
00050 #include "G4QKaonPlusElasticCrossSection.hh"          // Uzhi
00051 #include "G4QKaonMinusElasticCrossSection.hh"         // Uzhi
00052 #include <iostream>
00053 
00054 
00055 G4VQCrossSection* G4CHIPSElastic::pxsManager = 0;
00056 G4VQCrossSection* G4CHIPSElastic::nxsManager = 0;
00057 
00058 G4VQCrossSection* G4CHIPSElastic::PBARxsManager = 0;   // Uzhi
00059 G4VQCrossSection* G4CHIPSElastic::PIPxsManager = 0;
00060 G4VQCrossSection* G4CHIPSElastic::PIMxsManager = 0;
00061 G4VQCrossSection* G4CHIPSElastic::KPxsManager = 0;
00062 G4VQCrossSection* G4CHIPSElastic::KMxsManager = 0;
00063 
00064 G4CHIPSElastic::G4CHIPSElastic() : G4HadronElastic("hElasticCHIPS")
00065 {
00066   if(!pxsManager)
00067   {
00068     pxsManager    = G4QProtonElasticCrossSection::GetPointer();
00069     nxsManager    = G4QNeutronElasticCrossSection::GetPointer();
00070 
00071     PBARxsManager = G4QAntiBaryonElasticCrossSection::GetPointer(); // Uzhi
00072     PIPxsManager  = G4QPionPlusElasticCrossSection::GetPointer();   // Uzhi
00073     PIMxsManager  = G4QPionMinusElasticCrossSection::GetPointer();  // Uzhi
00074     KPxsManager   = G4QKaonPlusElasticCrossSection::GetPointer();   // Uzhi
00075     KMxsManager   = G4QKaonMinusElasticCrossSection::GetPointer();  // Uzhi
00076   }
00077   //Description();
00078 }
00079 
00080 G4CHIPSElastic::~G4CHIPSElastic()
00081 {}
00082 
00083 void G4CHIPSElastic::Description() const
00084 {
00085   char* dirName = getenv("G4PhysListDocDir");
00086   if (dirName) {
00087     std::ofstream outFile;
00088     G4String outFileName = GetModelName() + ".html";
00089     G4String pathName = G4String(dirName) + "/" + outFileName;
00090     outFile.open(pathName);
00091     outFile << "<html>\n";
00092     outFile << "<head>\n";
00093 
00094     outFile << "<title>Description of G4CHIPSElastic</title>\n";
00095     outFile << "</head>\n";
00096     outFile << "<body>\n";
00097 
00098     outFile << "The G4CHIPSElastic model performs hadron-nucleus elastic\n"
00099             << "scattering using the parameterized elastic cross sections\n"
00100             << "of M. Kossov\n";
00101 
00102     outFile << "</body>\n";
00103     outFile << "</html>\n";
00104     outFile.close();
00105   }
00106 }
00107 
00108 
00109 G4double 
00110 G4CHIPSElastic::SampleInvariantT(const G4ParticleDefinition* p, 
00111                                  G4double plab, G4int Z, G4int A)
00112 {
00113   G4int N = A - Z;
00114   if(Z == 1 && N == 2)      { N = 1; }
00115   else if(Z == 2 && N == 1) { N = 2; }
00116   G4int projPDG = p->GetPDGEncoding();
00117   G4double cs = 0.;
00118   if     (projPDG==2212) { cs = pxsManager->GetCrossSection(false,plab,Z,N,projPDG); }
00119   else if(projPDG==2112) { cs = nxsManager->GetCrossSection(false,plab,Z,N,projPDG); }
00120   else if(projPDG==-2212){ cs = PBARxsManager->GetCrossSection(false,plab,Z,N,projPDG); } //Pbar
00121   else if(projPDG== 211) { cs = PIPxsManager->GetCrossSection(false,plab,Z,N,projPDG); } // Pi+
00122   else if(projPDG==-211) { cs = PIMxsManager->GetCrossSection(false,plab,Z,N,projPDG); } // Pi-
00123   else if(projPDG== 321) { cs = KPxsManager->GetCrossSection(false,plab,Z,N,projPDG); } // K+
00124   else if(projPDG==-321) { cs = KMxsManager->GetCrossSection(false,plab,Z,N,projPDG); } // K-
00125 
00126   G4double t = 0.0;
00127   if(cs > 0.0)
00128   {
00129     if     (projPDG== 2212) { t = pxsManager->GetExchangeT(Z,N,projPDG); }
00130     else if(projPDG== 2112) { t = nxsManager->GetExchangeT(Z,N,projPDG); }
00131     else if(projPDG==-2212) { t = PBARxsManager->GetExchangeT(Z,N,projPDG); } // Pbar
00132     else if(projPDG==  211) { t = PIPxsManager->GetExchangeT(Z,N,projPDG); }  // Pi+
00133     else if(projPDG== -211) { t = PIMxsManager->GetExchangeT(Z,N,projPDG); }  // Pi-
00134     else if(projPDG==  321) { t = KPxsManager->GetExchangeT(Z,N,projPDG); }   // K+
00135     else if(projPDG== -321) { t = KMxsManager->GetExchangeT(Z,N,projPDG); }   // K-
00136   }
00137   else  { t = G4HadronElastic::SampleInvariantT(p, plab, Z, A); }
00138   return t;
00139 }
00140 

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