BDSIM
BDSIM is a Geant4 extension toolkit for simulation of particle transport in accelerator beamlines.
Loading...
Searching...
No Matches
BDSTrajectoryPrimary.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 BDSTRAJECTORYPRIMARY_H
20#define BDSTRAJECTORYPRIMARY_H
21
22#include "BDSTrajectory.hh"
23
24#include "G4Allocator.hh"
25#include "G4Types.hh"
26
27#include <ostream>
28
30class G4Step;
31class G4Track;
32
33namespace BDS
34{
35 struct TrajectoryOptions;
36}
37
49{
50public:
51 BDSTrajectoryPrimary(const G4Track* aTrack,
52 G4bool interactive,
53 const BDS::TrajectoryOptions& storageOptionsIn,
54 G4bool storeTrajectoryPointsIn);
55
58
59 virtual ~BDSTrajectoryPrimary();
60
62 virtual G4bool IsPrimary() const {return true;}
63
64 inline void* operator new(size_t);
65 inline void operator delete(void*);
66 inline int operator == (const BDSTrajectoryPrimary& right) const {return (this==&right);}
67
70 virtual void AppendStep(const G4Step* aStep);
71
75 virtual void MergeTrajectory(G4VTrajectory* secondTrajectory);
76
78 friend std::ostream& operator<< (std::ostream &out, BDSTrajectoryPrimary const &t);
79
81 const BDSTrajectoryPoint* FirstHit() const {return firstHit;}
82 const BDSTrajectoryPoint* LastPoint() const {return lastPoint;}
84
87 G4bool HasHitSomething() const {return firstHit;}
88
91 static G4bool hasScatteredThisTurn;
92
93protected:
96
97private:
98 const G4bool storeTrajectoryPoints;
99
101};
102
103extern G4Allocator<BDSTrajectoryPrimary> bdsTrajectoryPrimaryAllocator;
104
105inline void* BDSTrajectoryPrimary::operator new(size_t)
106{
107 void* aTrajectory;
108 aTrajectory = (void*)bdsTrajectoryPrimaryAllocator.MallocSingle();
109 return aTrajectory;
110}
111
112inline void BDSTrajectoryPrimary::operator delete(void* aTrajectory)
113{bdsTrajectoryPrimaryAllocator.FreeSingle((BDSTrajectoryPrimary*)aTrajectory);}
114
115#endif
A Point in a trajectory with extra information.
Trajectory information for only the primary.
friend std::ostream & operator<<(std::ostream &out, BDSTrajectoryPrimary const &t)
Output stream.
const BDSTrajectoryPoint * LastPoint() const
Accessor.
virtual void AppendStep(const G4Step *aStep)
G4bool HasHitSomething() const
BDSTrajectoryPoint * firstHit
Point owned by this class for the first scattering point.
const BDSTrajectoryPoint * FirstHit() const
Accessor.
static G4bool hasScatteredThisTurn
BDSTrajectoryPoint * lastPoint
Point owned by this class for the last-most point.
BDSTrajectoryPrimary(BDSTrajectoryPrimary &)=delete
copy constructor is not needed
const G4bool storeTrajectoryPoints
Whether to use base class to store all points.
BDSTrajectoryPrimary()=delete
No default constructor required.
virtual G4bool IsPrimary() const
Override and return true always.
virtual void MergeTrajectory(G4VTrajectory *secondTrajectory)
Trajectory information from track including last scatter etc.
Return either G4Tubs or G4CutTubs depending on flat face.