BDSIM
BDSIM is a Geant4 extension toolkit for simulation of particle transport in accelerator beamlines.
EventAnalysisOrbit.cc
1/*
2Beam Delivery Simulation (BDSIM) Copyright (C) Royal Holloway,
3University of London 2001 - 2022.
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 "EventAnalysisOrbit.hh"
20#include "RBDSException.hh"
21#include "SamplerAnalysis.hh"
22
23#include <iomanip>
24#include <iostream>
25#include <string>
26#include <vector>
27
28#include "TChain.h"
29#include "TFile.h"
30#include "TString.h"
31#include "TTree.h"
32
33ClassImp(EventAnalysisOrbit)
34
37{;}
38
40 TChain* chainIn,
41 bool perEntryAnalysis,
42 bool processSamplersIn,
43 bool debugIn,
44 double printModuloFraction,
45 bool emittanceOnTheFlyIn,
46 long int eventStartIn,
47 long int eventEndIn):
48 EventAnalysis(eventIn, chainIn, perEntryAnalysis, processSamplersIn,
49 debugIn, printModuloFraction, emittanceOnTheFlyIn,
50 eventStartIn, eventEndIn)
51{;}
52
54{
55 ss.clear();
56 x.clear();
57 xp.clear();
58 y.clear();
59 yp.clear();
60 elementName.clear();
61}
62
64{
65 if (index > entries-1)
66 {
67 std::string errString = "Orbit index: " + std::to_string(index) +
68 " greater than number of events: " + std::to_string(entries) + " in file minus one";
69 throw RBDSException(errString);
70 }
71
72 int nSamplers = (int)samplerAnalyses.size();
73
74 std::cout << "Getting orbit " << index << std::endl;
75 chain->GetEntry(index);
76 std::cout << "Loaded" << std::endl;
77
78 int counter = 0;
79 for (auto s : samplerAnalyses)
80 {
81 std::cout << "\rSampler #" << std::setw(6) << counter << " of " << nSamplers;
82 std::cout.flush();
83 if (s->s->n > 0)
84 {// valid entry on that sampler
85 ss.push_back(s->s->S);
86 x.push_back(s->s->x[0]);
87 xp.push_back(s->s->xp[0]);
88 y.push_back(s->s->y[0]);
89 yp.push_back(s->s->yp[0]);
90 elementName.push_back(s->s->samplerName);
91 }
92 counter++;
93 }
94 std::cout << std::endl;
95}
96
98{
99 f->cd("/");
100 TTree* orbitTree = new TTree("Orbit", "Orbit");
101
102 double dss;
103 double dx;
104 double dxp;
105 double dy;
106 double dyp;
107 std::string delementName;
108
109 orbitTree->Branch("s", &dss, "s/D");
110 orbitTree->Branch("x", &dx, "x/D");
111 orbitTree->Branch("xp", &dxp, "xp/D");
112 orbitTree->Branch("y", &dy, "y/D");
113 orbitTree->Branch("yp", &dyp, "yp/D");
114 orbitTree->Branch("elementName", &delementName);
115
116 for (int i = 0; i < (int)ss.size(); ++i)
117 {
118 dss = ss.at(i);
119 dx = x.at(i);
120 dxp = xp.at(i);
121 dy = y.at(i);
122 dyp = yp.at(i);
123 delementName = elementName.at(i);
124 orbitTree->Fill();
125 }
126
127 orbitTree->Write();
128 f->Close();
129}
long int entries
Number of entries in the chain.
Definition: Analysis.hh:99
Simple analysis to pull out first hit in each sampler.
std::vector< std::string > elementName
Temporary storage for orbit.
std::vector< double > y
Temporary storage for orbit.
std::vector< double > ss
Temporary storage for orbit.
std::vector< double > xp
Temporary storage for orbit.
void Clear()
Empty the member vectors of their data.
void ExtractOrbit(int index)
Extract an orbit from the data.
std::vector< double > x
Temporary storage for orbit.
std::vector< double > yp
Temporary storage for orbit.
void WriteOrbit(TFile *f)
Write orbit to a ROOT file.
Event level analysis.
std::vector< SamplerAnalysis * > samplerAnalyses
Holder for sampler analysis objects.
Event loader.
Definition: Event.hh:50
General exception with possible name of object and message.