BDSIM
BDSIM is a Geant4 extension toolkit for simulation of particle transport in accelerator beamlines.
Loading...
Searching...
No Matches
BDSLaserCompton.hh
1/*
2Beam Delivery Simulation (BDSIM) Copyright (C) Royal Holloway,
3University of London 2001 - 2023.
4
5This file is part of BDSIM.
6
7BDSIM is free software: you can redistribute it and/or modify
8it under the terms of the GNU General Public License as published
9by the Free Software Foundation version 3 of the License.
10
11BDSIM is distributed in the hope that it will be useful, but
12WITHOUT ANY WARRANTY; without even the implied warranty of
13MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14GNU General Public License for more details.
15
16You should have received a copy of the GNU General Public License
17along with BDSIM. If not, see <http://www.gnu.org/licenses/>.
18*/
19#ifndef BDSLASERCOMPTON_H
20#define BDSLASERCOMPTON_H
21
22#include "globals.hh"
23#include "G4VDiscreteProcess.hh"
24#include "G4Electron.hh"
25#include "G4Positron.hh"
26#include "BDSMaterials.hh"
27
28#include <limits>
29
32
33class G4Step;
34class G4Track;
35
36// flag initiated in BDSEventAction
37extern G4bool FireLaserCompton;
38
53class BDSLaserCompton: public G4VDiscreteProcess
54{
55public:
56
57 explicit BDSLaserCompton(const G4String& processName = "eLaser");
58
59 virtual ~BDSLaserCompton();
60
62 virtual G4bool IsApplicable(const G4ParticleDefinition&);
63
65 virtual G4double GetMeanFreePath(const G4Track& track,
66 G4double previousStepSize,
67 G4ForceCondition* condition );
68
69 virtual G4VParticleChange *PostStepDoIt(const G4Track& track,
70 const G4Step& step);
71
72 inline void SetLaserDirection(G4ThreeVector aDirection);
73 inline G4ThreeVector GetLaserDirection() const;
74
75 inline void SetLaserWavelength(G4double aWavelength);
76 inline G4double GetLaserWavelength() const;
77
78protected:
79 G4double ComputeMeanFreePath( const G4ParticleDefinition* ParticleType,
80 G4double KineticEnergy,
81 const G4Material* aMaterial);
82
83private:
84 // assignment and copy constructor not implemented nor used
85 BDSLaserCompton & operator=(const BDSLaserCompton &right);
87
88 BDSGlobalConstants* globals;
89
90 G4double laserWavelength;
91 G4ThreeVector laserDirection;
92 G4double laserEnergy;
93 BDSComptonEngine* comptonEngine;
94};
95
96inline G4bool BDSLaserCompton::IsApplicable(const G4ParticleDefinition& particle)
97{
98 return ( (&particle == G4Electron::Electron()) || (&particle == G4Positron::Positron()) );
99}
100
101inline G4double BDSLaserCompton::GetMeanFreePath(const G4Track& track,
102 G4double /*PreviousStepSize*/,
103 G4ForceCondition* ForceCondition)
104{
105 if ( track.GetMaterial() == BDSMaterials::Instance()->GetMaterial("LaserVac") &&
106 FireLaserCompton )
107 {*ForceCondition = Forced;}
108 return std::numeric_limits<double>::max();
109}
110
111inline void BDSLaserCompton::SetLaserDirection(G4ThreeVector laserDirectionIn)
112{laserDirection = laserDirectionIn;}
113
114inline G4ThreeVector BDSLaserCompton::GetLaserDirection() const
115{return laserDirection;}
116
117inline void BDSLaserCompton::SetLaserWavelength(G4double wavelengthIn)
118{laserWavelength = wavelengthIn;}
119
120inline G4double BDSLaserCompton::GetLaserWavelength() const
121{return laserWavelength;}
122
123#endif
Engine to calculate product of compton scattering process.
A class that holds global options and constants.
Laser compton scattering process to achieve a laserwire.
virtual G4bool IsApplicable(const G4ParticleDefinition &)
Overridden from G4VProcess - whether applicable to a particular particle.
virtual G4double GetMeanFreePath(const G4Track &track, G4double previousStepSize, G4ForceCondition *condition)
Overridden from G4VDiscreteProcess.
static BDSMaterials * Instance()
Singleton pattern access.
Definition: BDSMaterials.cc:39
G4Material * GetMaterial(G4String material) const
Get material by name.