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 #ifndef G4ProtonAntiProtonAtRestChips_h
00029 #define G4ProtonAntiProtonAtRestChips_h
00030
00031 #include <CLHEP/Units/SystemOfUnits.h>
00032
00033 #include "globals.hh"
00034 #include "G4VRestProcess.hh"
00035 #include "G4ParticleTable.hh"
00036 #include "G4Quasmon.hh"
00037 #include "G4QHadronVector.hh"
00038 #include "G4ParticleChange.hh"
00039 #include "G4LorentzVector.hh"
00040 #include "G4DynamicParticle.hh"
00041 #include "G4IonTable.hh"
00042 #include "G4Neutron.hh"
00043 #include "G4StopElementSelector.hh"
00044 #include "G4ChiralInvariantPhaseSpace.hh"
00045 #include "G4HadronicProcessType.hh"
00046 #include "G4HadronicDeprecate.hh"
00047
00048
00049 class G4ProtonAntiProtonAtRestChips : public G4VRestProcess
00050 {
00051 private:
00052
00053 G4ProtonAntiProtonAtRestChips& operator=(const G4ProtonAntiProtonAtRestChips &right);
00054 G4ProtonAntiProtonAtRestChips(const G4ProtonAntiProtonAtRestChips& );
00055
00056 public:
00057
00058 G4ProtonAntiProtonAtRestChips(const G4String& processName=
00059 "AntiProtonAnnihilationAtRest")
00060 : G4VRestProcess (processName, fHadronic)
00061 {
00062 G4HadronicDeprecate("G4ProtonAntiProtonAtRestChips");
00063 SetProcessSubType(fHadronAtRest);
00064 }
00065
00066 ~G4ProtonAntiProtonAtRestChips() {}
00067
00068 G4bool IsApplicable(const G4ParticleDefinition& aParticle)
00069 {
00070 return ( &aParticle == G4AntiProton::AntiProtonDefinition() );
00071 }
00072
00073
00074 void BuildPhysicsTable(const G4ParticleDefinition&){}
00075
00076 G4double AtRestGetPhysicalInteractionLength(const G4Track&track,
00077 G4ForceCondition*condition);
00078
00079
00080 G4double GetMeanLifeTime(const G4Track& aTrack,
00081 G4ForceCondition* condition) {return 0.0;}
00082
00083 G4VParticleChange* AtRestDoIt(const G4Track&, const G4Step&);
00084
00085 private:
00086 G4ChiralInvariantPhaseSpace theModel;
00087 G4StopElementSelector theSelector;
00088 };
00089
00090 inline
00091 G4VParticleChange * G4ProtonAntiProtonAtRestChips::
00092 AtRestDoIt(const G4Track& aTrack, const G4Step&aStep)
00093 {
00094
00095 G4Element * theTarget = theSelector.GetElement(aTrack.GetMaterial());
00096 G4Nucleus aTargetNucleus(theTarget->GetA_asInt(), theTarget->GetZ_asInt());
00097
00098
00099
00100 if(aTrack.GetDynamicParticle()->GetDefinition() != G4AntiProton::AntiProton())
00101 {
00102 throw G4HadronicException(__FILE__, __LINE__,
00103 "Calling G4ProtonAntiProtonAtRestChips with particle other than p-bar!!!");
00104 }
00105 if(aTargetNucleus.GetZ_asInt() != 1)
00106 {
00107 throw G4HadronicException(__FILE__, __LINE__,
00108 "Calling G4ProtonAntiProtonAtRestChips for target other than Hydrogen!!!");
00109 }
00110
00111
00112 return theModel.ApplyYourself(aTrack, aTargetNucleus);
00113 }
00114
00115 G4double G4ProtonAntiProtonAtRestChips::
00116 AtRestGetPhysicalInteractionLength(const G4Track&track,
00117 G4ForceCondition*condition)
00118 {
00119 ResetNumberOfInteractionLengthLeft();
00120 *condition = NotForced;
00121 currentInteractionLength = GetMeanLifeTime(track, condition);
00122 #ifdef CHIPSdebug
00123 if ((currentInteractionLength <0.0) || (verboseLevel>2))
00124 {
00125 G4cout << "G4ProtonAntiProtonAtRestChips::AtRestGetPhysicalInteractionLength ";
00126 G4cout << "[ " << GetProcessName() << "]" <<G4endl;
00127 track.GetDynamicParticle()->DumpInfo();
00128 G4cout << " in Material " << track.GetMaterial()->GetName() <<G4endl;
00129 G4cout << "MeanLifeTime = " << currentInteractionLength/CLHEP::ns << "[ns]" <<G4endl;
00130 }
00131 #endif
00132 return theNumberOfInteractionLengthLeft * currentInteractionLength;
00133 }
00134 #endif