BDSIM
BDSIM is a Geant4 extension toolkit for simulation of particle transport in accelerator beamlines.
Loading...
Searching...
No Matches
BDSTrajectory.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 BDSTRAJECTORY_H
20#define BDSTRAJECTORY_H
21#include "BDSTrajectoryOptions.hh"
22#include "BDSTrajectoryPoint.hh"
23#include "G4Trajectory.hh"
24
25#include <ostream>
26#include <vector>
27
28class G4Step;
29class G4Track;
30class G4TrajectoryContainer;
31class G4VTrajectoryPoint;
32
33typedef std::vector<BDSTrajectoryPoint*> BDSTrajectoryPointsContainer;
34
43class BDSTrajectory: public G4Trajectory
44{
45public:
46 BDSTrajectory() = delete;
47 BDSTrajectory(const G4Track* aTrack,
48 G4bool interactiveIn,
49 const BDS::TrajectoryOptions& storageOptionsIn);
52
53 virtual ~BDSTrajectory();
54
55 inline void* operator new(size_t);
56 inline void operator delete(void*);
57 inline int operator == (const BDSTrajectory& right) const {return (this==&right);}
58
61 virtual void AppendStep(const G4Step* aStep);
62
65 void AppendStep(const BDSTrajectoryPoint* pointIn);
66
70 void CleanPoint(BDSTrajectoryPoint* point) const;
71
73 virtual void MergeTrajectory(G4VTrajectory* secondTrajectory);
74
76 virtual G4VTrajectoryPoint* GetPoint(G4int i) const {return (*fpBDSPointsContainer)[i];}
77
79 virtual int GetPointEntries() const {return (int)fpBDSPointsContainer->size();}
80
82 virtual G4bool IsPrimary() const {return false;}
83
86 inline void SetTrajIndex(G4int trajIndexIn) {trajIndex = trajIndexIn;}
87 inline G4int GetTrajIndex() const {return trajIndex;}
88
91 inline void SetParentIndex(G4int parentIndexIn) {parentIndex = parentIndexIn;}
92 inline G4int GetParentIndex() const {return parentIndex;}
93
95 inline void SetParentStepIndex(G4int parentStepIndexIn) {parentStepIndex = parentStepIndexIn;}
96 inline G4int GetParentStepIndex() const {return parentStepIndex;}
97
99 inline G4int GetDepth() const {return depth;}
100 inline void SetDepth(G4int depthIn) {depth = depthIn;}
101
103 inline void SetParent(BDSTrajectory* parentIn) {parent = parentIn;}
104 inline BDSTrajectory* GetParent() const {return parent;}
105 inline G4int GetCreatorProcessType() const {return creatorProcessType;}
106 inline G4int GetCreatorProcessSubType() const {return creatorProcessSubType;}
107
109 friend std::ostream& operator<< (std::ostream &out, BDSTrajectory const &t);
110
115 BDSTrajectoryPoint* LastInteraction() const;
116
117protected:
118 G4int creatorProcessType;
119 G4int creatorProcessSubType;
120 G4double weight;
121 G4bool interactive;
122 G4bool suppressTransportationAndNotInteractive;
123 const BDS::TrajectoryOptions storageOptions;
124 BDSTrajectory* parent;
125 G4int trajIndex;
126 G4int parentIndex;
127 G4int parentStepIndex;
128 G4int depth;
129
133 BDSTrajectoryPointsContainer* fpBDSPointsContainer;
134};
135
136extern G4Allocator<BDSTrajectory> bdsTrajectoryAllocator;
137
138inline void* BDSTrajectory::operator new(size_t)
139{
140 void* aTrajectory;
141 aTrajectory = (void*)bdsTrajectoryAllocator.MallocSingle();
142 return aTrajectory;
143}
144
145inline void BDSTrajectory::operator delete(void* aTrajectory)
146{bdsTrajectoryAllocator.FreeSingle((BDSTrajectory*)aTrajectory);}
147
148
149#endif
A Point in a trajectory with extra information.
Trajectory information from track including last scatter etc.
virtual int GetPointEntries() const
Get number of trajectory points in this trajectory.
BDSTrajectoryPoint * FirstInteraction() const
virtual void MergeTrajectory(G4VTrajectory *secondTrajectory)
Merge another trajectory into this one.
void SetTrajIndex(G4int trajIndexIn)
void SetParentStepIndex(G4int parentStepIndexIn)
The index of the step along the parent trajectory from which this one was created.
virtual G4VTrajectoryPoint * GetPoint(G4int i) const
Access a point - use this class's container.
BDSTrajectory(BDSTrajectory &)=delete
copy constructor is not needed
friend std::ostream & operator<<(std::ostream &out, BDSTrajectory const &t)
Output stream.
void SetParent(BDSTrajectory *parentIn)
Record the parent trajectory.
void CleanPoint(BDSTrajectoryPoint *point) const
virtual void AppendStep(const G4Step *aStep)
virtual G4bool IsPrimary() const
Method to identify which one is a primary. Overridden in derived class.
void SetParentIndex(G4int parentIndexIn)
BDSTrajectoryPointsContainer * fpBDSPointsContainer
G4int GetDepth() const
Depth in the tree. Will be filled later once all trajectories are created and sorted.