BDSIM
BDSIM is a Geant4 extension toolkit for simulation of particle transport in accelerator beamlines.
Loading...
Searching...
No Matches
BDSLaserCompton.cc
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#include "BDSComptonEngine.hh"
20#include "BDSGlobalConstants.hh"
21#include "BDSLaserCompton.hh"
22
23#include "globals.hh"
24#include "G4Gamma.hh"
25
26#include "CLHEP/Units/PhysicalConstants.h"
27
28BDSLaserCompton::BDSLaserCompton(const G4String& processName):
29 G4VDiscreteProcess(processName),
30 laserEnergy(0.0)
31{
33 laserWavelength = globals->GetLaserwireWavelength();
34 laserDirection = globals->GetLaserwireDir();
35 comptonEngine = new BDSComptonEngine();
36}
37
38BDSLaserCompton::~BDSLaserCompton()
39{
40 delete comptonEngine;
41}
42
43G4VParticleChange* BDSLaserCompton::PostStepDoIt(const G4Track& trackData,
44 const G4Step& stepData)
45{
46 aParticleChange.Initialize(trackData);
47
48 // ensure that Laserwire can only occur once in an event
49 //G4cout << "FireLaserCompton == " << FireLaserCompton << G4endl;
50
51 if(!FireLaserCompton)
52 {return G4VDiscreteProcess::PostStepDoIt(trackData,stepData);}
53
54 G4Material* aMaterial = trackData.GetMaterial();
55
56 if (aMaterial == BDSMaterials::Instance()->GetMaterial("LaserVac"))
57 {
58 laserWavelength = globals->GetLaserwireWavelength();
59 laserDirection = globals->GetLaserwireDir();
60
61 laserEnergy = CLHEP::twopi*CLHEP::hbarc/laserWavelength;
62
63 // point laserwire in x: P_x Py Pz E
64 G4LorentzVector Laser4mom(laserEnergy*laserDirection.unit(),laserEnergy);
65
66 const G4DynamicParticle* aDynamicParticle = trackData.GetDynamicParticle();
67
68 comptonEngine->SetIncomingElectron4Vec(aDynamicParticle->Get4Momentum());
69 comptonEngine->SetIncomingPhoton4Vec(Laser4mom);
70 comptonEngine->PerformCompton();
71
72 if(globals->GetLaserwireTrackPhotons())
73 {
74 // create G4DynamicParticle object for the Gamma
75 G4LorentzVector ScatGam = comptonEngine->GetScatteredGamma();
76 G4DynamicParticle* aGamma = new G4DynamicParticle (G4Gamma::Gamma(),
77 ScatGam.vect().unit(),// direction
78 ScatGam.e());
79
80 aParticleChange.SetNumberOfSecondaries(1);
81 aParticleChange.AddSecondary(aGamma);
82 if(!globals->GetLaserwireTrackElectrons())
83 {
84 aParticleChange.ProposeEnergy(0.);
85 aParticleChange.ProposeLocalEnergyDeposit (0.);
86 aParticleChange.ProposeTrackStatus(fStopAndKill);
87 }
88 }
89 else
90 {
91 aParticleChange.SetNumberOfSecondaries(0);
92 aParticleChange.ProposeLocalEnergyDeposit (0.);
93 }
94
95 G4double NewKinEnergy = comptonEngine->GetScatteredElectron().e()-CLHEP::electron_mass_c2;
96
97 G4LorentzVector ScatEl = comptonEngine->GetScatteredElectron();
98
99 if (NewKinEnergy > 0.)
100 {
101 aParticleChange.ProposeMomentumDirection(ScatEl.vect().unit());
102 aParticleChange.ProposeEnergy(NewKinEnergy);
103 aParticleChange.ProposeLocalEnergyDeposit (0.);
104 }
105 else
106 {
107 aParticleChange.ProposeEnergy( 0. );
108 aParticleChange.ProposeLocalEnergyDeposit (0.);
109 G4double charge= aDynamicParticle->GetCharge();
110 if (charge<0.) aParticleChange.ProposeTrackStatus(fStopAndKill);
111 else aParticleChange.ProposeTrackStatus(fStopButAlive);
112 }
113 }
114
115 FireLaserCompton = false;
116
117 return G4VDiscreteProcess::PostStepDoIt(trackData,stepData);
118}
Engine to calculate product of compton scattering process.
static BDSGlobalConstants * Instance()
Access method.
static BDSMaterials * Instance()
Singleton pattern access.
Definition: BDSMaterials.cc:39