Geant4  10.03.p02
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4FTFModel.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 //
27 // $Id: G4FTFModel.hh 100828 2016-11-02 15:25:59Z gcosmo $
28 // GEANT4 tag $Name: $
29 //
30 // Class Description
31 // Final state production code for hadron inelastic scattering above 3 GeV
32 // based on the modeling ansatz used in FRITIOF.
33 // To be used in your physics list in case you need this physics.
34 // In this case you want to register an object of this class with an object
35 // of G4TheoFSGenerator.
36 // Class Description - End
37 
38 #ifndef G4FTFModel_h
39 #define G4FTFModel_h 1
40 
41 // ------------------------------------------------------------
42 // GEANT 4 class header file
43 //
44 // ---------------- G4FTFModel ----------------
45 // by Gunter Folger, May 1998.
46 // class implementing the excitation in the FTF Parton String Model
47 // ------------------------------------------------------------
48 
49 #include "G4VPartonStringModel.hh"
50 #include "G4FTFParameters.hh"
51 #include "G4FTFParticipants.hh"
52 #include "G4ExcitedStringVector.hh"
54 #include "G4ElasticHNScattering.hh"
55 #include "G4FTFAnnihilation.hh"
56 #include "G4Proton.hh"
57 #include "G4Neutron.hh"
58 
59 class G4VSplitableHadron;
60 class G4ExcitedString;
61 
62 
64  public:
65  G4FTFModel( const G4String& modelName = "FTF" );
66  ~G4FTFModel();
67 
68  void Init( const G4Nucleus& aNucleus, const G4DynamicParticle& aProjectile );
73 
74  virtual void ModelDescription( std::ostream& ) const;
75 
76  private:
77  G4FTFModel( const G4FTFModel& right );
78  const G4FTFModel& operator=( const G4FTFModel& right );
79  int operator==( const G4FTFModel& right ) const;
80  int operator!=( const G4FTFModel& right ) const;
81 
82  void StoreInvolvedNucleon();
83  void ReggeonCascade();
84  G4bool PutOnMassShell();
85  G4bool ExciteParticipants();
86  G4ExcitedStringVector* BuildStrings();
87  void GetResiduals();
88  G4bool AdjustNucleons( G4VSplitableHadron* SelectedAntiBaryon,
89  G4Nucleon* ProjectileNucleon,
90  G4VSplitableHadron* SelectedTargetNucleon,
91  G4Nucleon* TargetNucleon,
92  G4bool Annihilation );
93  G4ThreeVector GaussianPt( G4double AveragePt2, G4double maxPtSquare ) const;
94 
95  G4bool ComputeNucleusProperties( G4V3DNucleus* nucleus, G4LorentzVector& nucleusMomentum,
96  G4LorentzVector& residualMomentum, G4double& sumMasses,
97  G4double& residualExcitationEnergy, G4double& residualMass,
98  G4int& residualMassNumber, G4int& residualCharge );
99  // Utility method used by PutOnMassShell.
100 
101  G4bool GenerateDeltaIsobar( const G4double sqrtS, const G4int numberOfInvolvedNucleons,
102  G4Nucleon* involvedNucleons[], G4double& sumMasses );
103  // Utility method used by PutOnMassShell.
104 
105  G4bool SamplingNucleonKinematics( G4double averagePt2, const G4double maxPt2,
106  G4double dCor, G4V3DNucleus* nucleus,
107  const G4LorentzVector& pResidual,
108  const G4double residualMass, const G4int residualMassNumber,
109  const G4int numberOfInvolvedNucleons,
110  G4Nucleon* involvedNucleons[], G4double& mass2 );
111 
112  // Utility method used by PutOnMassShell.
113 
114  G4bool CheckKinematics( const G4double sValue, const G4double sqrtS,
115  const G4double projectileMass2, const G4double targetMass2,
116  const G4double nucleusY, const G4bool isProjectileNucleus,
117  const G4int numberOfInvolvedNucleons, G4Nucleon* involvedNucleons[],
118  G4double& targetWminus, G4double& projectileWplus, G4bool& success );
119  // Utility method used by PutOnMassShell.
120 
121  G4bool FinalizeKinematics( const G4double w, const G4bool isProjectileNucleus,
122  const G4LorentzRotation& boostFromCmsToLab,
123  const G4double residualMass, const G4int residualMassNumber,
124  const G4int numberOfInvolvedNucleons,
125  G4Nucleon* involvedNucleons[],
126  G4LorentzVector& residual4Momentum );
127  // Utility method used by PutOnMassShell.
128 
129  G4ReactionProduct theProjectile;
130  G4FTFParticipants theParticipants;
131 
132  G4Nucleon* TheInvolvedNucleonsOfTarget[250];
133  G4int NumberOfInvolvedNucleonsOfTarget;
134 
135  G4Nucleon* TheInvolvedNucleonsOfProjectile[250];
136  G4int NumberOfInvolvedNucleonsOfProjectile;
137 
138  G4FTFParameters* theParameters;
139  G4DiffractiveExcitation* theExcitation;
140  G4ElasticHNScattering* theElastic;
141  G4FTFAnnihilation* theAnnihilation;
142 
143  std::vector< G4VSplitableHadron* > theAdditionalString;
144 
145  G4double LowEnergyLimit;
146  G4bool HighEnergyInter;
147 
148  G4LorentzVector ProjectileResidual4Momentum;
149  G4int ProjectileResidualMassNumber;
150  G4int ProjectileResidualCharge;
151  G4double ProjectileResidualExcitationEnergy;
152 
153  G4LorentzVector TargetResidual4Momentum;
154  G4int TargetResidualMassNumber;
155  G4int TargetResidualCharge;
156  G4double TargetResidualExcitationEnergy;
157 };
158 
159 
161  return theParticipants.GetWoundedNucleus();
162 }
163 
164 
166  return theParticipants.GetWoundedNucleus();
167 }
168 
169 
171  return theParticipants.GetProjectileNucleus();
172 }
173 
174 #endif
175 
std::vector< G4ExcitedString * > G4ExcitedStringVector
G4V3DNucleus * GetWoundedNucleus() const
Definition: G4FTFModel.hh:160
G4FTFModel(const G4String &modelName="FTF")
Definition: G4FTFModel.cc:72
G4V3DNucleus * GetTargetNucleus() const
Definition: G4FTFModel.hh:165
int G4int
Definition: G4Types.hh:78
G4V3DNucleus * GetProjectileNucleus() const
Definition: G4FTFModel.hh:170
bool G4bool
Definition: G4Types.hh:79
virtual G4V3DNucleus * GetProjectileNucleus() const
virtual G4V3DNucleus * GetWoundedNucleus() const
void Init(const G4Nucleus &aNucleus, const G4DynamicParticle &aProjectile)
Definition: G4FTFModel.cc:150
virtual void ModelDescription(std::ostream &) const
Definition: G4FTFModel.cc:3366
double G4double
Definition: G4Types.hh:76
G4ExcitedStringVector * GetStrings()
Definition: G4FTFModel.cc:272