Geant4  10.01.p02
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 90572 2015-06-04 09:32:32Z 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 
65  public:
66  G4FTFModel( const G4String& modelName = "FTF" );
67  ~G4FTFModel();
68 
69  void Init( const G4Nucleus& aNucleus, const G4DynamicParticle& aProjectile );
74 
75  virtual void ModelDescription( std::ostream& ) const;
76 
77  private:
78  G4FTFModel( const G4FTFModel& right );
79  const G4FTFModel& operator=( const G4FTFModel& right );
80  int operator==( const G4FTFModel& right ) const;
81  int operator!=( const G4FTFModel& right ) const;
82 
83  void StoreInvolvedNucleon();
84  void ReggeonCascade();
88  void GetResiduals();
89  G4bool AdjustNucleons( G4VSplitableHadron* SelectedAntiBaryon,
90  G4Nucleon* ProjectileNucleon,
91  G4VSplitableHadron* SelectedTargetNucleon,
92  G4Nucleon* TargetNucleon,
93  G4bool Annihilation );
94  G4ThreeVector GaussianPt( G4double AveragePt2, G4double maxPtSquare ) const;
95 
96  G4bool ComputeNucleusProperties( G4V3DNucleus* nucleus, G4LorentzVector& nucleusMomentum,
97  G4LorentzVector& residualMomentum, G4double& sumMasses,
98  G4double& residualExcitationEnergy, G4double& residualMass,
99  G4int& residualMassNumber, G4int& residualCharge );
100  // Utility method used by PutOnMassShell.
101 
102  G4bool GenerateDeltaIsobar( const G4double sqrtS, const G4int numberOfInvolvedNucleons,
103  G4Nucleon* involvedNucleons[], G4double& sumMasses );
104  // Utility method used by PutOnMassShell.
105 
106  G4bool SamplingNucleonKinematics( G4double averagePt2, const G4double maxPt2,
107  G4double dCor, G4V3DNucleus* nucleus,
108  const G4LorentzVector& pResidual,
109  const G4double residualMass, const G4int residualMassNumber,
110  const G4int numberOfInvolvedNucleons,
111  G4Nucleon* involvedNucleons[], G4double& mass2 );
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 
131 
134 
137 
142 
143  std::vector< G4VSplitableHadron* > theAdditionalString;
144 
147 
152 
157 
158 };
159 
160 
163 }
164 
165 
168 }
169 
170 
173 }
174 
175 #endif
G4ElasticHNScattering * theElastic
Definition: G4FTFModel.hh:140
std::vector< G4ExcitedString * > G4ExcitedStringVector
G4V3DNucleus * GetWoundedNucleus() const
Definition: G4FTFModel.hh:161
G4bool PutOnMassShell()
Definition: G4FTFModel.cc:518
G4FTFModel(const G4String &modelName="FTF")
Definition: G4FTFModel.cc:69
CLHEP::Hep3Vector G4ThreeVector
CLHEP::HepLorentzRotation G4LorentzRotation
G4bool AdjustNucleons(G4VSplitableHadron *SelectedAntiBaryon, G4Nucleon *ProjectileNucleon, G4VSplitableHadron *SelectedTargetNucleon, G4Nucleon *TargetNucleon, G4bool Annihilation)
Definition: G4FTFModel.cc:981
void ReggeonCascade()
Definition: G4FTFModel.cc:413
G4ReactionProduct theProjectile
Definition: G4FTFModel.hh:129
int operator!=(const G4FTFModel &right) const
G4V3DNucleus * GetTargetNucleus() const
Definition: G4FTFModel.hh:166
G4int TargetResidualMassNumber
Definition: G4FTFModel.hh:154
int G4int
Definition: G4Types.hh:78
G4FTFParameters * theParameters
Definition: G4FTFModel.hh:138
G4FTFAnnihilation * theAnnihilation
Definition: G4FTFModel.hh:141
G4bool ExciteParticipants()
Definition: G4FTFModel.cc:808
G4Nucleon * TheInvolvedNucleonsOfTarget[250]
Definition: G4FTFModel.hh:132
const G4FTFModel & operator=(const G4FTFModel &right)
std::vector< G4VSplitableHadron * > theAdditionalString
Definition: G4FTFModel.hh:143
void GetResiduals()
Definition: G4FTFModel.cc:2426
G4FTFParticipants theParticipants
Definition: G4FTFModel.hh:130
G4V3DNucleus * GetProjectileNucleus() const
Definition: G4FTFModel.hh:171
bool G4bool
Definition: G4Types.hh:79
G4ThreeVector GaussianPt(G4double AveragePt2, G4double maxPtSquare) const
Definition: G4FTFModel.cc:2597
G4int ProjectileResidualMassNumber
Definition: G4FTFModel.hh:149
G4DiffractiveExcitation * theExcitation
Definition: G4FTFModel.hh:139
virtual G4V3DNucleus * GetProjectileNucleus() const
virtual G4V3DNucleus * GetWoundedNucleus() const
G4double LowEnergyLimit
Definition: G4FTFModel.hh:145
G4ExcitedStringVector * BuildStrings()
Definition: G4FTFModel.cc:2103
G4bool ComputeNucleusProperties(G4V3DNucleus *nucleus, G4LorentzVector &nucleusMomentum, G4LorentzVector &residualMomentum, G4double &sumMasses, G4double &residualExcitationEnergy, G4double &residualMass, G4int &residualMassNumber, G4int &residualCharge)
Definition: G4FTFModel.cc:2617
G4Nucleon * TheInvolvedNucleonsOfProjectile[250]
Definition: G4FTFModel.hh:135
void Init(const G4Nucleus &aNucleus, const G4DynamicParticle &aProjectile)
Definition: G4FTFModel.cc:141
virtual void ModelDescription(std::ostream &) const
Definition: G4FTFModel.cc:3025
G4int ProjectileResidualCharge
Definition: G4FTFModel.hh:150
G4bool CheckKinematics(const G4double sValue, const G4double sqrtS, const G4double projectileMass2, const G4double targetMass2, const G4double nucleusY, const G4bool isProjectileNucleus, const G4int numberOfInvolvedNucleons, G4Nucleon *involvedNucleons[], G4double &targetWminus, G4double &projectileWplus, G4bool &success)
Definition: G4FTFModel.cc:2879
G4int NumberOfInvolvedNucleonsOfProjectile
Definition: G4FTFModel.hh:136
void StoreInvolvedNucleon()
Definition: G4FTFModel.cc:362
G4bool SamplingNucleonKinematics(G4double averagePt2, const G4double maxPt2, G4double dCor, G4V3DNucleus *nucleus, const G4LorentzVector &pResidual, const G4double residualMass, const G4int residualMassNumber, const G4int numberOfInvolvedNucleons, G4Nucleon *involvedNucleons[], G4double &mass2)
Definition: G4FTFModel.cc:2762
G4bool GenerateDeltaIsobar(const G4double sqrtS, const G4int numberOfInvolvedNucleons, G4Nucleon *involvedNucleons[], G4double &sumMasses)
Definition: G4FTFModel.cc:2696
G4bool FinalizeKinematics(const G4double w, const G4bool isProjectileNucleus, const G4LorentzRotation &boostFromCmsToLab, const G4double residualMass, const G4int residualMassNumber, const G4int numberOfInvolvedNucleons, G4Nucleon *involvedNucleons[], G4LorentzVector &residual4Momentum)
Definition: G4FTFModel.cc:2953
G4int TargetResidualCharge
Definition: G4FTFModel.hh:155
G4LorentzVector TargetResidual4Momentum
Definition: G4FTFModel.hh:153
G4bool HighEnergyInter
Definition: G4FTFModel.hh:146
double G4double
Definition: G4Types.hh:76
G4LorentzVector ProjectileResidual4Momentum
Definition: G4FTFModel.hh:148
G4double ProjectileResidualExcitationEnergy
Definition: G4FTFModel.hh:151
G4ExcitedStringVector * GetStrings()
Definition: G4FTFModel.cc:259
G4int NumberOfInvolvedNucleonsOfTarget
Definition: G4FTFModel.hh:133
G4double TargetResidualExcitationEnergy
Definition: G4FTFModel.hh:156
CLHEP::HepLorentzVector G4LorentzVector
int operator==(const G4FTFModel &right) const