Geant4 9.6.0
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
G4ProtonAntiProtonAtRestChips.hh
Go to the documentation of this file.
1//
2// ********************************************************************
3// * License and Disclaimer *
4// * *
5// * The Geant4 software is copyright of the Copyright Holders of *
6// * the Geant4 Collaboration. It is provided under the terms and *
7// * conditions of the Geant4 Software License, included in the file *
8// * LICENSE and available at http://cern.ch/geant4/license . These *
9// * include a list of copyright holders. *
10// * *
11// * Neither the authors of this software system, nor their employing *
12// * institutes,nor the agencies providing financial support for this *
13// * work make any representation or warranty, express or implied, *
14// * regarding this software system or assume any liability for its *
15// * use. Please see the license in the file LICENSE and URL above *
16// * for the full disclaimer and the limitation of liability. *
17// * *
18// * This code implementation is the result of the scientific and *
19// * technical work of the GEANT4 collaboration. *
20// * By using, copying, modifying or distributing the software (or *
21// * any work based on the software) you agree to acknowledge its *
22// * use in resulting scientific publications, and indicate your *
23// * acceptance of all terms of the Geant4 Software license. *
24// ********************************************************************
25//
26// (Why only antiproton-proton, when the antiproton-nucleus is made? - M.K.)
27// 17.02.2009 M.Kossov, now it is recommended to use the G4QCaptureAtRest process
28#ifndef G4ProtonAntiProtonAtRestChips_h
29#define G4ProtonAntiProtonAtRestChips_h
30
32
33#include "globals.hh"
34#include "G4VRestProcess.hh"
35#include "G4ParticleTable.hh"
36#include "G4Quasmon.hh"
37#include "G4QHadronVector.hh"
38#include "G4ParticleChange.hh"
39#include "G4LorentzVector.hh"
40#include "G4DynamicParticle.hh"
41#include "G4IonTable.hh"
42#include "G4Neutron.hh"
47
48
50{
51 private:
52 // hide assignment operator as private
55
56 public:
57
59 "AntiProtonAnnihilationAtRest")
60 : G4VRestProcess (processName, fHadronic)
61 {
62 G4HadronicDeprecate("G4ProtonAntiProtonAtRestChips");
64 }
65
67
69 {
70 return ( &aParticle == G4AntiProton::AntiProtonDefinition() );
71 }
72
73 // null physics table
75
78
79 // zero mean lifetime
81 G4ForceCondition* condition) {return 0.0;}
82
84
85 private:
87 G4StopElementSelector theSelector; // Assume identical laws as for muons
88};
89
90inline
92AtRestDoIt(const G4Track& aTrack, const G4Step&aStep)
93{
94 // Create target
95 G4Element * theTarget = theSelector.GetElement(aTrack.GetMaterial());
96 G4Nucleus aTargetNucleus(theTarget->GetA_asInt(), theTarget->GetZ_asInt());
97
98 // Check model validity - note this will be a sub-branch in the ordinary stopping @@@@@@
99 // in the long haul. @@@@@@
101 {
102 throw G4HadronicException(__FILE__, __LINE__,
103 "Calling G4ProtonAntiProtonAtRestChips with particle other than p-bar!!!");
104 }
105 if(aTargetNucleus.GetZ_asInt() != 1)
106 {
107 throw G4HadronicException(__FILE__, __LINE__,
108 "Calling G4ProtonAntiProtonAtRestChips for target other than Hydrogen!!!");
109 }
110
111 // Call chips
112 return theModel.ApplyYourself(aTrack, aTargetNucleus);
113}
114
118{
122#ifdef CHIPSdebug
123 if ((currentInteractionLength <0.0) || (verboseLevel>2))
124 {
125 G4cout << "G4ProtonAntiProtonAtRestChips::AtRestGetPhysicalInteractionLength ";
126 G4cout << "[ " << GetProcessName() << "]" <<G4endl;
127 track.GetDynamicParticle()->DumpInfo();
128 G4cout << " in Material " << track.GetMaterial()->GetName() <<G4endl;
129 G4cout << "MeanLifeTime = " << currentInteractionLength/CLHEP::ns << "[ns]" <<G4endl;
130 }
131#endif
133}
134#endif
G4double condition(const G4ErrorSymMatrix &m)
G4ForceCondition
@ NotForced
#define G4HadronicDeprecate(name)
@ fHadronAtRest
@ fHadronic
double G4double
Definition: G4Types.hh:64
bool G4bool
Definition: G4Types.hh:67
#define G4endl
Definition: G4ios.hh:52
G4DLLIMPORT std::ostream G4cout
static G4AntiProton * AntiProton()
Definition: G4AntiProton.cc:93
static G4AntiProton * AntiProtonDefinition()
Definition: G4AntiProton.cc:88
G4HadFinalState * ApplyYourself(const G4HadProjectile &aTrack, G4Nucleus &aTargetNucleus, G4HadFinalState *aChange=0)
void DumpInfo(G4int mode=0) const
G4ParticleDefinition * GetDefinition() const
const G4String & GetName() const
Definition: G4Material.hh:177
G4int GetZ_asInt() const
Definition: G4Nucleus.hh:115
void BuildPhysicsTable(const G4ParticleDefinition &)
G4VParticleChange * AtRestDoIt(const G4Track &, const G4Step &)
G4ProtonAntiProtonAtRestChips(const G4String &processName="AntiProtonAnnihilationAtRest")
G4double AtRestGetPhysicalInteractionLength(const G4Track &track, G4ForceCondition *condition)
G4double GetMeanLifeTime(const G4Track &aTrack, G4ForceCondition *condition)
G4bool IsApplicable(const G4ParticleDefinition &aParticle)
Definition: G4Step.hh:78
G4Element * GetElement(const G4Material *aMaterial)
G4Material * GetMaterial() const
const G4DynamicParticle * GetDynamicParticle() const
G4double currentInteractionLength
Definition: G4VProcess.hh:297
virtual void ResetNumberOfInteractionLengthLeft()
Definition: G4VProcess.cc:92
G4int verboseLevel
Definition: G4VProcess.hh:368
G4double theNumberOfInteractionLengthLeft
Definition: G4VProcess.hh:293
void SetProcessSubType(G4int)
Definition: G4VProcess.hh:403
const G4String & GetProcessName() const
Definition: G4VProcess.hh:379