00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020
00021
00022
00023
00024
00025
00026
00027
00028
00029
00030
00031
00032
00033
00034
00035
00036
00037
00038 #include "G4QDiscProcessMixer.hh"
00039 #include "G4PhysicalConstants.hh"
00040 #include "G4SystemOfUnits.hh"
00041 #include "G4HadronicDeprecate.hh"
00042
00043
00044
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
00059 G4QDiscProcessMixer::~G4QDiscProcessMixer()
00060 {
00061
00062
00063 }
00064
00065 void G4QDiscProcessMixer::AddDiscreteProcess(G4VDiscreteProcess* DP, G4double ME)
00066 {
00067 static const G4double maxEn = 1.E8*megaelectronvolt;
00068 if(!theDPVector.size())
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
00083 {
00084
00085
00086
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();
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
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();
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
00141
00142 return theDPVector[i]->first->PostStepDoIt(Track, Step);
00143 }
00144 return G4VDiscreteProcess::PostStepDoIt(Track, Step);
00145 }
00146
00147
00148
00149
00150
00151