G4QDiscProcessMixer.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 //      ---------------- G4QDiscProcessMixer class -------------------
00029 //                 by Mikhail Kossov, Aug 2007.
00030 // G4QDiscProcessMixer class of the CHIPS Simulation Branch in GEANT4
00031 // ------------------------------------------------------------------------
00032 // Short description: universal mixer of processes (NOT models as in GHAD!)
00033 // depending on the application energy region (defined by users).
00034 // ------------------------------------------------------------------------
00035 
00036 //#define debug
00037 
00038 #include "G4QDiscProcessMixer.hh"
00039 #include "G4PhysicalConstants.hh"
00040 #include "G4SystemOfUnits.hh"
00041 #include "G4HadronicDeprecate.hh"
00042 
00043 
00044 // Constructor
00045 G4QDiscProcessMixer::G4QDiscProcessMixer(const G4String& name,
00046                                          const G4ParticleDefinition* particle,
00047                                          G4ProcessType pType):
00048   G4VDiscreteProcess(name, pType), DPParticle(particle)
00049 {
00050   G4HadronicDeprecate("G4QDiscProcessMixer");
00051 
00052 #ifdef debug
00053   G4cout<<"G4QDiscProcessMixer::Constructor is called processName="<<name<<G4endl;
00054 #endif
00055   if (verboseLevel>0) G4cout<<GetProcessName()<<" process is created "<<G4endl;
00056 }
00057 
00058 // Destructor
00059 G4QDiscProcessMixer::~G4QDiscProcessMixer()
00060 {
00061   // Now the responsibility of deleting is deligated to the user, who created them
00062   //for_each(theDPVector.begin(), theDPVector.end(), DeleteDiscreteProcess());
00063 }
00064 
00065 void G4QDiscProcessMixer::AddDiscreteProcess(G4VDiscreteProcess* DP, G4double ME)
00066 {
00067   static const G4double maxEn = 1.E8*megaelectronvolt; // Conditional INF
00068   if(!theDPVector.size()) // The first process in the DiscreteProcessVector (MaxE=INF)
00069   {
00070     std::pair<G4VDiscreteProcess*, G4double>* QDiscProc =
00071                                    new std::pair<G4VDiscreteProcess*, G4double>(DP, maxEn);
00072     theDPVector.push_back( QDiscProc );
00073   }
00074   else
00075   {
00076     if(ME < theDPVector[theDPVector.size()-1]->second)
00077     {
00078       std::pair<G4VDiscreteProcess*, G4double>* QDiscProc =
00079                                       new std::pair<G4VDiscreteProcess*, G4double>(DP, ME);
00080       theDPVector.push_back( QDiscProc );
00081     }
00082     else // Wrong Max Energy Order for the new process in the sequence of processes
00083     {
00084       // G4cerr<<"G4QDiscProcessMixer::AddDiscreteProcess:LastMaxE("<<theDPVector.size()-1
00085       //       <<")="<<theDPVector[theDPVector.size()-1]->second<<" <= MaxE="<<ME<<G4endl;
00086       // G4Exception("G4QDiscProcessMixer::AddDiscreteProcess: Wrong Max Energy Order");
00087       G4ExceptionDescription ed;
00088       ed << " LastMaxE(" << theDPVector.size()-1 << ")="
00089          << theDPVector[theDPVector.size()-1]->second << " <= MaxE=" << ME << G4endl;
00090       G4Exception("G4QDiscProcessMixer::AddDiscreteProcess()", "HAD_CHPS_0000",
00091                   FatalException, ed);
00092     }
00093   }
00094 }
00095 
00096 G4bool G4QDiscProcessMixer::IsApplicable(const G4ParticleDefinition& particle) 
00097 {
00098 #ifdef debug
00099   G4cout<<"G4QDiscProcessMixer::IsApplicable: part="<<particle.GetParticleName()<<" = "
00100         <<DPParticle->GetParticleName()<<G4endl;
00101 #endif
00102   if(particle == *DPParticle) return true;
00103   return false;
00104 }
00105 
00106 G4double G4QDiscProcessMixer::PostStepGetPhysicalInteractionLength(const G4Track& Track,
00107                                                                    G4double  PrevStSize,
00108                                                                    G4ForceCondition* F)
00109 {
00110   G4double kEn=Track.GetDynamicParticle()->GetKineticEnergy(); // Projectile kinetic energy
00111   G4int maxDP=theDPVector.size();
00112   if(maxDP) for(G4int i=maxDP-1; i>-1; i--) if(kEn < theDPVector[i]->second)
00113   {
00114 #ifdef debug
00115     G4cout<<"G4QDPMix::PSGetPIL:"<<i<<",kEn="<<kEn<<" < "<<theDPVector[i]->second<<", PIL="
00116         <<theDPVector[i]->first->PostStepGetPhysicalInteractionLength(Track,PrevStSize,F)
00117         <<G4endl;
00118 #endif
00119     return theDPVector[i]->first->PostStepGetPhysicalInteractionLength(Track,PrevStSize,F);
00120   }
00121   return DBL_MAX;
00122 }
00123 
00124 // compilation fake class (length is calculated in PostStepGetPhysicalInteractionLength)
00125 G4double G4QDiscProcessMixer::GetMeanFreePath(const G4Track&, G4double, G4ForceCondition*)
00126 {
00127   return DBL_MAX;
00128 }
00129 
00130 G4VParticleChange* G4QDiscProcessMixer::PostStepDoIt(const G4Track& Track,
00131                                                      const G4Step& Step)
00132 {
00133   G4double kEn=Track.GetDynamicParticle()->GetKineticEnergy(); // Projectile kinetic energy
00134   G4int maxDP=theDPVector.size();
00135   if(maxDP) for(G4int i=maxDP-1; i>-1; i--) if(kEn < theDPVector[i]->second)
00136   {
00137 #ifdef debug
00138     G4cout<<"G4QDPMix::PSDoIt: i="<<i<<",kEn="<<kEn<<" < "<<theDPVector[i]->second<<G4endl;
00139 #endif
00140     //EnMomConservation= theDPVector[i]->first->GetEnegryMomentumConservation();
00141     //nOfNeutrons      = theDPVector[i]->first->GetNumberOfNeutronsInTarget();
00142     return theDPVector[i]->first->PostStepDoIt(Track, Step);
00143   }
00144   return G4VDiscreteProcess::PostStepDoIt(Track, Step);
00145 }
00146 
00147 //G4LorentzVector G4QDiscProcessMixer::GetEnegryMomentumConservation()
00148 //                                                              {return EnMomConservation;}
00149 
00150 //G4int G4QDiscProcessMixer::GetNumberOfNeutronsInTarget()
00151 //                                                                    {return nOfNeutrons;}

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