Changeset 1315 for trunk/source/processes/hadronic/models/cascade/cascade/include/G4InuclParticle.hh
- Timestamp:
- Jun 18, 2010, 11:42:07 AM (14 years ago)
- File:
-
- 1 edited
Legend:
- Unmodified
- Added
- Removed
-
trunk/source/processes/hadronic/models/cascade/cascade/include/G4InuclParticle.hh
r962 r1315 1 #ifndef G4INUCL_PARTICLE_HH 2 #define G4INUCL_PARTICLE_HH 1 3 // 2 4 // ******************************************************************** … … 23 25 // * acceptance of all terms of the Geant4 Software license. * 24 26 // ******************************************************************** 27 // $Id: G4InuclParticle.hh,v 1.19 2010/05/21 18:07:30 mkelsey Exp $ 28 // Geant4 tag: $Name: geant4-09-04-beta-cand-01 $ 25 29 // 26 #ifndef G4INUCL_PARTICLE_HH 27 #define G4INUCL_PARTICLE_HH 30 // 20100112 M. Kelsey -- Remove G4CascadeMomentum, use G4LorentzVector directly 31 // 20100409 M. Kelsey -- Drop unused string argument from ctors. 32 // 20100519 M. Kelsey -- Add public access to G4DynamicParticle content 28 33 29 #ifndef GLOB 34 #include "G4DynamicParticle.hh" 35 #include "G4LorentzVector.hh" 30 36 #include "globals.hh" 31 #endif32 37 33 #include <iostream>34 #include <vector>35 #include "G4CascadeMomentum.hh"36 37 // Notice: no cc-file for G4InuclParticle38 38 39 39 class G4InuclParticle { 40 public: 41 G4InuclParticle() : modelId(0) {} 40 42 41 public: 42 G4InuclParticle() { 43 setModel(0); // default model 44 }; 43 explicit G4InuclParticle(const G4LorentzVector& mom) : modelId(0) { 44 pDP.Set4Momentum(mom*GeV/MeV); // From Bertini to G4 units 45 } 45 46 46 virtual ~G4InuclParticle() { }; 47 48 G4InuclParticle(const G4CascadeMomentum& mom) { 49 setMomentum(mom); 50 setModel(0); 51 }; 47 virtual ~G4InuclParticle() {} 52 48 53 void setMomentum(const G4CascadeMomentum& mom) {54 momentum = mom;55 };49 // Copy and assignment constructors for use with std::vector<> 50 G4InuclParticle(const G4InuclParticle& right) 51 : pDP(right.pDP), modelId(right.modelId) {} 56 52 53 G4InuclParticle& operator=(const G4InuclParticle& right); 57 54 58 const G4CascadeMomentum& getMomentum() const { 59 return momentum; 60 }; 55 // This is no longer required, as setMomentum() handles mass adjustment 56 void setEnergy() { ; } 61 57 62 G4double getMomModule() const { 63 return std::sqrt(momentum[1] * momentum[1] + 64 momentum[2] * momentum[2] + 65 momentum[3] * momentum[3]); 66 }; 67 68 virtual void printParticle() const { 69 G4cout << " px " << momentum[1] << " py " << momentum[2] << 70 " pz " << momentum[3] << 71 " pmod " << std::sqrt(momentum[1] * momentum[1] + 72 momentum[2] * momentum[2] + 73 momentum[3] * momentum[3]) 74 << " E " << momentum[0] 75 << " creator model " << modelId << G4endl; 76 }; 58 // These are call-throughs to G4DynamicParticle 59 void setMomentum(const G4LorentzVector& mom); 77 60 78 void setModel(G4int model) { 79 modelId = model; 80 }; 61 void setMass(G4double mass) { pDP.SetMass(mass*GeV/MeV); } 81 62 82 G4int getModel() { 83 return modelId; 84 }; 63 G4double getMass() const { 64 return pDP.GetMass()*MeV/GeV; // From G4 to Bertini units 65 } 66 67 G4double getCharge() const { 68 return pDP.GetCharge(); 69 } 70 71 G4double getKineticEnergy() const { 72 return pDP.GetKineticEnergy()*MeV/GeV; // From G4 to Bertini units 73 } 74 75 G4double getEnergy() const { 76 return pDP.GetTotalEnergy()*MeV/GeV; // From G4 to Bertini units 77 } 78 79 G4double getMomModule() const { 80 return pDP.GetTotalMomentum()*MeV/GeV; // From G4 to Bertini units 81 } 82 83 G4LorentzVector getMomentum() const { 84 return pDP.Get4Momentum()*MeV/GeV; // From G4 to Bertini units 85 } 86 87 virtual void printParticle() const; 88 89 void setModel(G4int model) { modelId = model; } 90 91 G4int getModel() const { return modelId; } 92 93 G4ParticleDefinition* getDefinition() const { 94 return pDP.GetDefinition(); 95 } 96 97 const G4DynamicParticle& getDynamicParticle() const { 98 return pDP; 99 } 85 100 86 101 protected: 87 G4CascadeMomentum momentum; 102 // Special constructors for subclasses to set particle type correctly 103 explicit G4InuclParticle(G4ParticleDefinition* pd) : modelId(0) { 104 setDefinition(pd); 105 } 106 107 // FIXME: Bertini code doesn't pass valid 4-vectors, so force mass value 108 // from supplied PartDefn, with required unit conversions 109 G4InuclParticle(G4ParticleDefinition* pd, const G4LorentzVector& mom); 110 111 // NOTE: Momentum forced along Z direction 112 G4InuclParticle(G4ParticleDefinition* pd, G4double ekin) 113 : pDP(pd,G4ThreeVector(0.,0.,1.),ekin*GeV/MeV), modelId(0) {} 114 115 void setDefinition(G4ParticleDefinition* pd) { pDP.SetDefinition(pd); } 88 116 89 117 private: 90 G4 int modelId; // used to indicate model that created instance of G4InuclParticle118 G4DynamicParticle pDP; // Carries all the kinematics and info 91 119 120 G4int modelId; 121 // used to indicate model that created instance of G4InuclParticle 92 122 // 0 default 93 123 // 1 bullet
Note: See TracChangeset
for help on using the changeset viewer.