Geant4  v4-10.4-release
 모두 클래스 네임스페이스들 파일들 함수 변수 타입정의 열거형 타입 열거형 멤버 Friends 매크로 그룹들 페이지들
G4HadronicProcess.hh
이 파일의 문서화 페이지로 가기
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 // $Id: G4HadronicProcess.hh 110586 2018-05-31 12:01:42Z gcosmo $
27 //
28 // -------------------------------------------------------------------
29 //
30 // GEANT4 Class header file
31 //
32 // G4HadronicProcess
33 //
34 // This is the top level Hadronic Process class
35 // The inelastic, elastic, capture, and fission processes
36 // should derive from this class
37 //
38 // original by H.P.Wellisch
39 // J.L. Chuma, TRIUMF, 10-Mar-1997
40 // Last modified: 04-Apr-1997
41 // 19-May-2008 V.Ivanchenko cleanup and added comments
42 // 05-Jul-2010 V.Ivanchenko cleanup commented lines
43 // 28-Jul-2012 M.Maire add function GetTargetDefinition()
44 // 14-Sep-2012 Inherit from RestDiscrete, use subtype code (now in ctor) to
45 // configure base-class
46 // 28-Sep-2012 M. Kelsey -- Undo inheritance change, keep new ctor
47 
48 #ifndef G4HadronicProcess_h
49 #define G4HadronicProcess_h 1
50 
51 #include "globals.hh"
52 #include "G4VDiscreteProcess.hh"
53 #include "G4EnergyRangeManager.hh"
54 #include "G4Nucleus.hh"
55 #include "G4ReactionProduct.hh"
56 #include "G4HadronicProcessType.hh"
58 #include <vector>
59 
60 class G4Track;
61 class G4Step;
62 class G4Element;
63 class G4ParticleChange;
68 
70 {
71 public:
72  G4HadronicProcess(const G4String& processName="Hadronic",
73  G4ProcessType procType=fHadronic);
74 
75  // Preferred signature for subclasses, specifying their subtype here
76  G4HadronicProcess(const G4String& processName,
77  G4HadronicProcessType subType);
78 
79  ~G4HadronicProcess() override;
80 
81  // register generator of secondaries
83 
84  // get cross section per element
86  const G4Element * elm,
87  const G4Material* mat = nullptr);
88 
89  // obsolete method to get cross section per element
90  inline
92  const G4Element * elm,
93  const G4Material* mat = nullptr)
94  { return GetElementCrossSection(part, elm, mat); }
95 
96  // generic PostStepDoIt recommended for all derived classes
97  G4VParticleChange* PostStepDoIt(const G4Track& aTrack,
98  const G4Step& aStep) override;
99 
100  // initialisation of physics tables and G4HadronicProcessStore
101  void PreparePhysicsTable(const G4ParticleDefinition&) override;
102 
103  // build physics tables and print out the configuration of the process
104  void BuildPhysicsTable(const G4ParticleDefinition&) override;
105 
106  // dump physics tables
108 
109  // add cross section data set
110  void AddDataSet(G4VCrossSectionDataSet * aDataSet);
111 
112  // access to the list of hadronic interactions
113  std::vector<G4HadronicInteraction*>& GetHadronicInteractionList();
114 
115  // get inverse cross section per volume
116  G4double GetMeanFreePath(const G4Track &aTrack, G4double,
117  G4ForceCondition *) override;
118 
119  // access to the target nucleus
120  inline const G4Nucleus* GetTargetNucleus() const
121  { return &targetNucleus; }
122 
123  // G4ParticleDefinition* GetTargetDefinition();
124  inline const G4Isotope* GetTargetIsotope()
125  { return targetNucleus.GetIsotope(); }
126 
127  void ProcessDescription(std::ostream& outFile) const override;
128 
129 protected:
130 
131  // generic method to choose secondary generator
132  // recommended for all derived classes
134  const G4HadProjectile & aHadProjectile, G4Nucleus& aTargetNucleus,
135  const G4Material* aMaterial, const G4Element* anElement)
136  { return theEnergyRangeManager.GetHadronicInteraction(aHadProjectile,
137  aTargetNucleus,
138  aMaterial,anElement);
139  }
140 
141  // access to the target nucleus
143  { return &targetNucleus; }
144 
145 public:
146 
147  // scale cross section
148  void BiasCrossSectionByFactor(G4double aScale);
149  void MultiplyCrossSectionBy(G4double factor);
151  { return aScaleFactor; }
152 
153  // Integral option
154  inline void SetIntegral(G4bool val)
155  { useIntegralXS = val; }
156 
157  // Energy-momentum non-conservation limits and reporting
158  inline void SetEpReportLevel(G4int level)
159  { epReportLevel = level; }
160 
161  inline void SetEnergyMomentumCheckLevels(G4double relativeLevel, G4double absoluteLevel)
162  { epCheckLevels.first = relativeLevel;
163  epCheckLevels.second = absoluteLevel;
164  levelsSetByProcess = true;
165  }
166 
167  inline std::pair<G4double, G4double> GetEnergyMomentumCheckLevels() const
168  { return epCheckLevels; }
169 
170  // access to the cross section data store
172  {return theCrossSectionDataStore;}
173 
174 protected:
175 
176  void DumpState(const G4Track&, const G4String&, G4ExceptionDescription&);
177 
178  // access to the chosen generator
180  { return theInteraction; }
181 
182  // access to the cross section data set
184  { return theLastCrossSection; }
185 
186  // fill result
187  void FillResult(G4HadFinalState* aR, const G4Track& aT);
188 
189  // Check the result for catastrophic energy non-conservation
191  const G4Nucleus& targetNucleus,
193 
194  // Check 4-momentum balance
196 
197 private:
198 
199  void InitialiseLocal();
200 
203 
204  // hide assignment operator as private
206  G4HadronicProcess(const G4HadronicProcess&) = delete;
207 
208  // Set E/p conservation check levels from environment variables
210 
211 protected:
212 
214 
216 
218 
220 
221 private:
222 
224 
226 
228 
230 
232 
234 
236 
238 
240 
242 
243  // Energy-momentum checking
244  std::pair<G4double, G4double> epCheckLevels;
246 
247  std::vector<G4VLeadingParticleBiasing*> theBias;
248 
250 
254 };
255 
256 #endif
257 
G4ParticleChange * theTotalResult
G4Nucleus * GetTargetNucleusPointer()
G4HadronicProcess & operator=(const G4HadronicProcess &right)=delete
G4HadProjectile thePro
std::ostringstream G4ExceptionDescription
Definition: G4Exception.hh:45
std::vector< ExP01TrackerHit * > a
Definition: ExP01Classes.hh:33
G4CrossSectionDataStore * theCrossSectionDataStore
G4HadronicInteraction * GetHadronicInteraction() const
void RegisterMe(G4HadronicInteraction *a)
const char * p
Definition: xmltok.h:285
std::pair< G4double, G4double > GetEnergyMomentumCheckLevels() const
G4double GetLastCrossSection()
G4EnergyRangeManager theEnergyRangeManager
G4VParticleChange * PostStepDoIt(const G4Track &aTrack, const G4Step &aStep) override
std::vector< G4VLeadingParticleBiasing * > theBias
G4HadronicProcessType
void PreparePhysicsTable(const G4ParticleDefinition &) override
void ProcessDescription(std::ostream &outFile) const override
G4double CrossSectionFactor() const
G4double XBiasSurvivalProbability()
void SetIntegral(G4bool val)
G4double GetElementCrossSection(const G4DynamicParticle *part, const G4Element *elm, const G4Material *mat=nullptr)
G4HadronicProcessStore * theProcessStore
G4HadronicInteraction * GetHadronicInteraction(const G4HadProjectile &aHadProjectile, G4Nucleus &aTargetNucleus, const G4Material *aMaterial, const G4Element *anElement) const
void BuildPhysicsTable(const G4ParticleDefinition &) override
double G4double
Definition: G4Types.hh:76
bool G4bool
Definition: G4Types.hh:79
TString part[npart]
Definition: Style.C:32
G4double theInitialNumberOfInteractionLength
G4double GetMeanFreePath(const G4Track &aTrack, G4double, G4ForceCondition *) override
void GetEnergyMomentumCheckEnvvars()
void DumpPhysicsTable(const G4ParticleDefinition &p)
void SetEpReportLevel(G4int level)
G4ProcessType
const G4Isotope * GetIsotope()
Definition: G4Nucleus.hh:119
std::pair< G4double, G4double > epCheckLevels
G4HadFinalState * CheckResult(const G4HadProjectile &thePro, const G4Nucleus &targetNucleus, G4HadFinalState *result)
void BiasCrossSectionByFactor(G4double aScale)
Definition: G4Step.hh:76
Float_t mat
void MultiplyCrossSectionBy(G4double factor)
void FillResult(G4HadFinalState *aR, const G4Track &aT)
~G4HadronicProcess() override
G4double G4ParticleHPJENDLHEData::G4double result
void DumpState(const G4Track &, const G4String &, G4ExceptionDescription &)
void CheckEnergyMomentumConservation(const G4Track &, const G4Nucleus &)
void AddDataSet(G4VCrossSectionDataSet *aDataSet)
G4CrossSectionDataStore * GetCrossSectionDataStore()
int G4int
Definition: G4Types.hh:78
G4double XBiasSecondaryWeight()
const G4Nucleus * GetTargetNucleus() const
const G4Isotope * GetTargetIsotope()
G4ForceCondition
G4HadronicProcess(const G4String &processName="Hadronic", G4ProcessType procType=fHadronic)
void SetEnergyMomentumCheckLevels(G4double relativeLevel, G4double absoluteLevel)
G4double GetMicroscopicCrossSection(const G4DynamicParticle *part, const G4Element *elm, const G4Material *mat=nullptr)
G4HadronicInteraction * theInteraction
G4HadronicInteraction * ChooseHadronicInteraction(const G4HadProjectile &aHadProjectile, G4Nucleus &aTargetNucleus, const G4Material *aMaterial, const G4Element *anElement)
std::vector< G4HadronicInteraction * > & GetHadronicInteractionList()