G4HadronicInteraction.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 // Hadronic Interaction  base class
00029 // original by H.P. Wellisch
00030 // modified by J.L. Chuma, TRIUMF, 21-Mar-1997
00031 // Last modified: 04-Apr-1997
00032 // reimplemented 1.11.2003 JPW.
00033 // 23-Jan-2009 V.Ivanchenko move constructor and destructor to the body
00034 
00035 #include <iostream>
00036 
00037 #include "G4HadronicInteraction.hh"
00038 #include "G4SystemOfUnits.hh"
00039 #include "G4HadronicInteractionRegistry.hh"
00040 #include "G4HadronicException.hh"
00041 
00042 G4HadronicInteraction::G4HadronicInteraction(const G4String& modelName) :
00043   verboseLevel(0), theMinEnergy(0.0), theMaxEnergy(25.0*GeV), 
00044   isBlocked(false), recoilEnergyThreshold(0.0), theModelName(modelName),
00045   epCheckLevels(DBL_MAX, DBL_MAX)
00046 { 
00047   G4HadronicInteractionRegistry::Instance()->RegisterMe(this);
00048 }
00049 
00050 
00051 G4HadronicInteraction::~G4HadronicInteraction()
00052 {
00053   G4HadronicInteractionRegistry::Instance()->RemoveMe(this);
00054 }
00055 
00056 
00057 G4double 
00058 G4HadronicInteraction::SampleInvariantT(const G4ParticleDefinition*, 
00059                                         G4double, G4int, G4int)
00060 {
00061   return 0.0;
00062 }
00063  
00064 G4double G4HadronicInteraction::GetMinEnergy(
00065    const G4Material *aMaterial, const G4Element *anElement ) const
00066 {
00067   if( IsBlocked(aMaterial) ) { return 0.0; }
00068   if( IsBlocked(anElement) ) { return 0.0; }
00069   size_t length = theMinEnergyListElements.size();
00070   if(0 < length) {
00071     for(size_t i=0; i<length; ++i ) {
00072         if( anElement == theMinEnergyListElements[i].second )
00073           { return theMinEnergyListElements[i].first; }
00074     }
00075   }
00076   length = theMinEnergyList.size();
00077   if(0 < length) {
00078     for(size_t i=0; i<length; ++i ) {
00079       if( aMaterial == theMinEnergyList[i].second )
00080         { return theMinEnergyList[i].first; }
00081     }
00082   }
00083   if(IsBlocked()) { return 0.0; }
00084   if( verboseLevel > 1 ) {
00085     G4cout << "*** Warning from HadronicInteraction::GetMinEnergy" << G4endl
00086            << "    material " << aMaterial->GetName()
00087            << " not found in min energy List" << G4endl;
00088   } 
00089   return theMinEnergy;
00090 }
00091  
00092 void G4HadronicInteraction::SetMinEnergy(G4double anEnergy,
00093                                          const G4Element *anElement )
00094 {
00095   if( IsBlocked(anElement) ) {
00096     G4cout << "*** Warning from HadronicInteraction::SetMinEnergy" << G4endl
00097            << "    The model is not active for the Element  "
00098            << anElement->GetName() << "." << G4endl;
00099   }
00100   size_t length = theMinEnergyListElements.size();
00101   if(0 < length) {
00102     for(size_t i=0; i<length; ++i ) {
00103       if( anElement == theMinEnergyListElements[i].second )
00104         {
00105           theMinEnergyListElements[i].first = anEnergy;
00106           return;
00107         }
00108     }
00109   }
00110   theMinEnergyListElements.push_back(std::pair<G4double, const G4Element *>(anEnergy, anElement));
00111 }
00112  
00113 void G4HadronicInteraction::SetMinEnergy(G4double anEnergy,
00114                                          const G4Material *aMaterial )
00115 {
00116   if( IsBlocked(aMaterial) ) {
00117     G4cout << "*** Warning from HadronicInteraction::SetMinEnergy" << G4endl
00118            << "    The model is not active for the Material "
00119            << aMaterial->GetName() << "." << G4endl;
00120   }
00121   size_t length = theMinEnergyList.size();
00122   if(0 < length) {
00123     for(size_t i=0; i<length; ++i ) {
00124       if( aMaterial == theMinEnergyList[i].second )
00125         {
00126           theMinEnergyList[i].first = anEnergy;
00127           return;
00128         }
00129     }
00130   }
00131   theMinEnergyList.push_back(std::pair<G4double, const G4Material *>(anEnergy, aMaterial));
00132 }
00133  
00134 G4double G4HadronicInteraction::GetMaxEnergy(const G4Material *aMaterial, 
00135                                              const G4Element *anElement ) const
00136 {
00137   if( IsBlocked(aMaterial) ) { return 0.0; }
00138   if( IsBlocked(anElement) ) { return 0.0; }
00139   size_t length = theMaxEnergyListElements.size();
00140   if(0 < length) {
00141     for(size_t i=0; i<length; ++i ) {
00142         if( anElement == theMaxEnergyListElements[i].second )
00143           { return theMaxEnergyListElements[i].first; }
00144     }
00145   }
00146   length = theMaxEnergyList.size();
00147   if(0 < length) {
00148     for(size_t i=0; i<length; ++i ) {
00149       if( aMaterial == theMaxEnergyList[i].second )
00150         { return theMaxEnergyList[i].first; }
00151     }
00152   }
00153   if(IsBlocked()) { return 0.0; }
00154   if( verboseLevel > 1 ) {
00155     G4cout << "*** Warning from HadronicInteraction::GetMaxEnergy" << G4endl
00156            << "    material " << aMaterial->GetName()
00157            << " not found in min energy List" << G4endl;
00158   } 
00159   return theMaxEnergy;
00160 }
00161  
00162 void G4HadronicInteraction::SetMaxEnergy(G4double anEnergy,
00163                                          const G4Element *anElement ) 
00164 {
00165   if( IsBlocked(anElement) ) {
00166     G4cout << "*** Warning from HadronicInteraction::SetMaxEnergy" << G4endl
00167            << "Warning: The model is not active for the Element  "
00168            << anElement->GetName() << "." << G4endl;
00169   }
00170   size_t length = theMaxEnergyListElements.size();
00171   if(0 < length) {
00172     for(size_t i=0; i<length; ++i ) {
00173       if( anElement == theMaxEnergyListElements[i].second )
00174       {
00175         theMaxEnergyListElements[i].first = anEnergy;
00176         return;
00177       }
00178     }
00179   }
00180   theMaxEnergyListElements.push_back(std::pair<G4double, const G4Element *>(anEnergy, anElement));
00181 }
00182 
00183 void G4HadronicInteraction::SetMaxEnergy(G4double anEnergy,
00184                                          const G4Material *aMaterial )
00185 {
00186   if( IsBlocked(aMaterial) ) {
00187     G4cout << "*** Warning from HadronicInteraction::SetMaxEnergy" << G4endl
00188            << "Warning: The model is not active for the Material "
00189            << aMaterial->GetName() << "." << G4endl;
00190   }
00191   size_t length = theMaxEnergyList.size();
00192   if(0 < length) {
00193     for(size_t i=0; i<length; ++i ) {
00194       if( aMaterial == theMaxEnergyList[i].second )
00195         {
00196           theMaxEnergyList[i].first = anEnergy;
00197           return;
00198         }
00199     }
00200   }
00201   theMaxEnergyList.push_back(std::pair<G4double, const G4Material *>(anEnergy, aMaterial));
00202 }
00203 
00204 void G4HadronicInteraction::DeActivateFor( const G4Material *aMaterial )
00205 {
00206   theBlockedList.push_back(aMaterial);
00207 }
00208 
00209 void G4HadronicInteraction::DeActivateFor( const G4Element *anElement )
00210 {
00211   theBlockedListElements.push_back(anElement);
00212 }
00213 
00214 
00215 G4bool G4HadronicInteraction::IsBlocked(const G4Material* aMaterial) const
00216 {
00217   for (size_t i=0; i<theBlockedList.size(); ++i) {
00218     if (aMaterial == theBlockedList[i]) return true;
00219   }
00220   return false;
00221 }
00222 
00223 
00224 G4bool G4HadronicInteraction::IsBlocked(const G4Element* anElement) const
00225 {
00226   for (size_t i=0; i<theBlockedListElements.size(); ++i) {
00227     if (anElement == theBlockedListElements[i]) return true;
00228   }
00229   return false;
00230 }
00231 
00232 const std::pair<G4double, G4double> G4HadronicInteraction::GetFatalEnergyCheckLevels() const
00233 {
00234         // default level of Check
00235         return std::pair<G4double, G4double>(10.*perCent, 5 * GeV);
00236 }
00237 
00238 std::pair<G4double, G4double>
00239 G4HadronicInteraction::GetEnergyMomentumCheckLevels() const
00240 {
00241   return epCheckLevels;
00242 }
00243 
00244 
00245 void G4HadronicInteraction::ModelDescription(std::ostream& outFile) const
00246 {
00247   outFile << "The description for this model has not been written yet.\n";
00248 }
00249 
00250 /*
00251 G4HadronicInteraction::G4HadronicInteraction(const G4HadronicInteraction &right )
00252 { 
00253   *this = right; 
00254 }
00255     
00256 const G4HadronicInteraction& 
00257 G4HadronicInteraction::operator=(const G4HadronicInteraction &right )
00258 { 
00259   G4String text = "unintended use of G4HadronicInteraction::operator=";
00260   throw G4HadronicException(__FILE__, __LINE__, text); 
00261   return right;
00262 }
00263  */
00264 /* end of file */
00265  

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