Geant4 11.2.2
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
G4ChannelingFastSimModel.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#ifndef G4ChannelingFastSimModel_h
28#define G4ChannelingFastSimModel_h 1
29
31#include "G4Step.hh"
32#include "G4TouchableHandle.hh"
33#include <vector>
36
38#include <unordered_map>
39#include "G4BaierKatkov.hh"
40#include "G4LogicalVolume.hh"
41#include "G4ParticleTable.hh"
42
43/** \file G4ChannelingFastSimModel.hh
44* \brief Definition of the G4ChannelingFastSimModel class
45* FastSimulation Channeling model: calculates charge particle trajectories
46* in oriented crystals in the field of crystal planes/axes either straight or bent.
47* It is also possible to simulate radiation using Baier-Katkov method.
48*/
49
51{
52public:
53 // Constructor, destructor
57
58 /// -- IsApplicable
60 /// -- ModelTrigger
61 G4bool ModelTrigger(const G4FastTrack &) override;
62 /// -- User method DoIt
63 void DoIt(const G4FastTrack&, G4FastStep&) override;
64
65 ///special functions
66 void Input(const G4Material* crystal, const G4String &lattice);
67
69
71
72 G4BaierKatkov* GetRadiationModel() {return fBaierKatkov;}
73
75
76 ///set cuts
77 void SetLowKineticEnergyLimit(G4double ekinetic, const G4String& particleName)
78 {fLowEnergyLimit[particleTable->FindParticle(particleName)->
79 GetParticleDefinitionID()] = ekinetic;}
80 void SetLindhardAngleNumberHighLimit(G4double angleNumber, const G4String& particleName)
81 {fLindhardAngleNumberHighLimit[particleTable->FindParticle(particleName)->
82 GetParticleDefinitionID()]=angleNumber;}
83
85 {fDefaultLowEnergyLimit=ekinetic;}
87 {fDefaultLindhardAngleNumberHighLimit=angleNumber;}
88
89
90 /// get the maximal number of photons that can be produced per fastStep
91 /// Caution: is redundant, if the radiation model is not activated
93 {fMaxPhotonsProducedPerStep=nPhotons;}
94
95 ///get cuts
97 {return GetLowKineticEnergyLimit(particleTable->
98 FindParticle(particleName)->
99 GetParticleDefinitionID());}
101 {return GetLindhardAngleNumberHighLimit(particleTable->
102 FindParticle(particleName)->
103 GetParticleDefinitionID());}
104 //the same functions but using particleDefinitionID (needed for faster model execution)
106 {return (fLowEnergyLimit.count(particleDefinitionID) == 1)
107 ? fLowEnergyLimit[particleDefinitionID]
108 : fDefaultLowEnergyLimit;}
110 {return (fLindhardAngleNumberHighLimit.count(particleDefinitionID) == 1)
111 ? fLindhardAngleNumberHighLimit[particleDefinitionID]
112 : fDefaultLindhardAngleNumberHighLimit;}
113
114 /// get the maximal number of photons that can be produced per fastStep
115 G4int GetMaxPhotonsProducedPerStep(){return fMaxPhotonsProducedPerStep;}
116
117private:
118
119 G4ChannelingFastSimCrystalData* fCrystalData{nullptr};
120 G4BaierKatkov* fBaierKatkov{nullptr};
121
123
124 ///flag of radiation model
125 G4bool fRad = false;
126
127 /// maps of cuts
128 std::unordered_map<G4int, G4double> fLowEnergyLimit;
129 std::unordered_map<G4int, G4double> fLindhardAngleNumberHighLimit;
130
131 G4double fDefaultLowEnergyLimit = 200*CLHEP::MeV;
132 G4double fDefaultLindhardAngleNumberHighLimit = 100.;
133
134 /// the maximal number of photons that can be produced per fastStep
135 G4int fMaxPhotonsProducedPerStep=1000.;
136
137};
138#endif
139
140
141
142
Definition of the G4BaierKatkov class This class is designed for the calculation of radiation probabi...
Definition of the G4ChannelingFastSimCrystalData class The class inherits G4VChannelingFastSimCrystal...
double G4double
Definition G4Types.hh:83
bool G4bool
Definition G4Types.hh:86
int G4int
Definition G4Types.hh:85
void SetDefaultLindhardAngleNumberHighLimit(G4double angleNumber)
G4ChannelingFastSimCrystalData * GetCrystalData()
G4ChannelingFastSimModel(const G4String &, G4Region *)
void SetDefaultLowKineticEnergyLimit(G4double ekinetic)
void DoIt(const G4FastTrack &, G4FastStep &) override
– User method DoIt
void Input(const G4Material *crystal, const G4String &lattice)
special functions
G4bool IsApplicable(const G4ParticleDefinition &) override
– IsApplicable
void SetLindhardAngleNumberHighLimit(G4double angleNumber, const G4String &particleName)
G4double GetLindhardAngleNumberHighLimit(G4int particleDefinitionID)
void SetMaxPhotonsProducedPerStep(G4double nPhotons)
G4bool ModelTrigger(const G4FastTrack &) override
– ModelTrigger
void SetLowKineticEnergyLimit(G4double ekinetic, const G4String &particleName)
set cuts
G4double GetLowKineticEnergyLimit(const G4String &particleName)
get cuts
G4int GetMaxPhotonsProducedPerStep()
get the maximal number of photons that can be produced per fastStep
G4double GetLindhardAngleNumberHighLimit(const G4String &particleName)
G4double GetLowKineticEnergyLimit(G4int particleDefinitionID)
G4ParticleDefinition * FindParticle(G4int PDGEncoding)
static G4ParticleTable * GetParticleTable()