BDSIM
BDSIM is a Geant4 extension toolkit for simulation of particle transport in accelerator beamlines.
Loading...
Searching...
No Matches
BDSBH4D.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 "BDSBH4D.hh"
20
21#include "TFile.h"
22
23#ifdef USE_BOOST
24#include "BDSBH4DTypeDefs.hh"
25#include <boost/format.hpp>
26#include <boost/histogram.hpp>
27#endif
28
29#include <iostream>
30#include <sstream>
31#include <string>
32#include <vector>
33
34#ifdef USE_BOOST
35templateClassImp(BDSBH4D)
36
37template <>
39 BDSBH4DBase(3,3,3,3, 0,1, 0,1, 0,1, 0.1,0.230, "BDSBH4D","BDSBH4D","linear")
40{
41 h = boost::histogram::make_histogram_with(std::vector<double>(),
42 boost::histogram::axis::regular<double> {3, 0.0, 1.0, "x"},
43 boost::histogram::axis::regular<double> {3, 0.0, 1.0, "y"},
44 boost::histogram::axis::regular<double> {3, 0.0, 1.0, "z"},
45 boost::histogram::axis::regular<double> {3, 1.0, 230.0, "energy"});
46
47 h_err = boost::histogram::make_histogram_with(std::vector<double>(),
48 boost::histogram::axis::regular<double> {3, 0.0, 1.0, "x"},
49 boost::histogram::axis::regular<double> {3, 0.0, 1.0, "y"},
50 boost::histogram::axis::regular<double> {3, 0.0, 1.0, "z"},
51 boost::histogram::axis::regular<double> {3, 1.0, 230.0, "energy"});
52}
53
54template <>
56 BDSBH4DBase(3,3,3,3, 0,1, 0,1, 0,1, 0.1,0.230, "BDSBH4D","BDSBH4D","log")
57{
58 h = boost::histogram::make_histogram_with(std::vector<double>(),
59 boost::histogram::axis::regular<double> {3, 0.0, 1.0, "x"},
60 boost::histogram::axis::regular<double> {3, 0.0, 1.0, "y"},
61 boost::histogram::axis::regular<double> {3, 0.0, 1.0, "z"},
62 boost::histogram::axis::regular<double, boost::histogram::axis::transform::log> {3, 1.0, 230.0, "energy"});
63
64 h_err = boost::histogram::make_histogram_with(std::vector<double>(),
65 boost::histogram::axis::regular<double> {3, 0.0, 1.0, "x"},
66 boost::histogram::axis::regular<double> {3, 0.0, 1.0, "y"},
67 boost::histogram::axis::regular<double> {3, 0.0, 1.0, "z"},
68 boost::histogram::axis::regular<double, boost::histogram::axis::transform::log> {3, 1.0, 230.0, "energy"});
69}
70
71template <>
73 BDSBH4DBase(3,3,3,0,1, 0,1, 0,1, "BDSBH4D","BDSBH4D","user",std::vector<double>{0.001,0.1,0.230})
74{
75 h = boost::histogram::make_histogram_with(std::vector<double>(),
76 boost::histogram::axis::regular<double> {3, 0.0, 1.0, "x"},
77 boost::histogram::axis::regular<double> {3, 0.0, 1.0, "y"},
78 boost::histogram::axis::regular<double> {3, 0.0, 1.0, "z"},
79 boost::histogram::axis::variable<double> {std::vector<double>{0.001,0.230}, "energy"});
80
81 h_err = boost::histogram::make_histogram_with(std::vector<double>(),
82 boost::histogram::axis::regular<double> {3, 0.0, 1.0, "x"},
83 boost::histogram::axis::regular<double> {3, 0.0, 1.0, "y"},
84 boost::histogram::axis::regular<double> {3, 0.0, 1.0, "z"},
85 boost::histogram::axis::variable<double> {std::vector<double>{0.001,0.1,0.230}, "energy"});
86}
87
88template <>
89BDSBH4D<boost_histogram_linear>::BDSBH4D(std::string& name, std::string& title, const std::string& eScale,
90 unsigned int nxbins, double xmin, double xmax,
91 unsigned int nybins, double ymin, double ymax,
92 unsigned int nzbins, double zmin, double zmax,
93 unsigned int nebins, double emin, double emax):
94 BDSBH4DBase(nxbins, nybins, nzbins, nebins,
95 xmin, xmax, ymin, ymax, zmin, zmax, emin, emax,
96 name, title, eScale)
97{
98 h = boost::histogram::make_histogram_with(std::vector<double>(),
99 boost::histogram::axis::regular<double> {nxbins, xmin, xmax, "x"},
100 boost::histogram::axis::regular<double> {nybins, ymin, ymax, "y"},
101 boost::histogram::axis::regular<double> {nzbins, zmin, zmax, "z"},
102 boost::histogram::axis::regular<double> {nebins, emin, emax, "energy"});
103
104 h_err = boost::histogram::make_histogram_with(std::vector<double>(),
105 boost::histogram::axis::regular<double> {nxbins, xmin, xmax, "x"},
106 boost::histogram::axis::regular<double> {nybins, ymin, ymax, "y"},
107 boost::histogram::axis::regular<double> {nzbins, zmin, zmax, "z"},
108 boost::histogram::axis::regular<double> {nebins, emin, emax, "energy"});
109}
110
111template <>
112BDSBH4D<boost_histogram_log>::BDSBH4D(std::string& name, std::string& title, const std::string& eScale,
113 unsigned int nxbins, double xmin, double xmax,
114 unsigned int nybins, double ymin, double ymax,
115 unsigned int nzbins, double zmin, double zmax,
116 unsigned int nebins, double emin, double emax):
117 BDSBH4DBase(nxbins, nybins, nzbins, nebins,
118 xmin, xmax, ymin, ymax, zmin, zmax, emin, emax,
119 name, title, eScale)
120{
121 h = boost::histogram::make_histogram_with(std::vector<double>(),
122 boost::histogram::axis::regular<double> {nxbins, xmin, xmax, "x"},
123 boost::histogram::axis::regular<double> {nybins, ymin, ymax, "y"},
124 boost::histogram::axis::regular<double> {nzbins, zmin, zmax, "z"},
125 boost::histogram::axis::regular<double, boost::histogram::axis::transform::log> {nebins, emin, emax, "energy"});
126
127 h_err = boost::histogram::make_histogram_with(std::vector<double>(),
128 boost::histogram::axis::regular<double> {nxbins, xmin, xmax, "x"},
129 boost::histogram::axis::regular<double> {nybins, ymin, ymax, "y"},
130 boost::histogram::axis::regular<double> {nzbins, zmin, zmax, "z"},
131 boost::histogram::axis::regular<double, boost::histogram::axis::transform::log> {nebins, emin, emax, "energy"});
132}
133
134template <>
135BDSBH4D<boost_histogram_variable>::BDSBH4D(std::string& name, std::string& title, const std::string& eScale,
136 const std::vector<double>& eBinEdgesIn,
137 unsigned int nxbins, double xmin, double xmax,
138 unsigned int nybins, double ymin, double ymax,
139 unsigned int nzbins, double zmin, double zmax):
140 BDSBH4DBase(nxbins, nybins, nzbins, xmin, xmax, ymin, ymax, zmin, zmax, name, title, eScale, eBinEdgesIn)
141{
142 h = boost::histogram::make_histogram_with(std::vector<double>(),
143 boost::histogram::axis::regular<double> {nxbins, xmin, xmax, "x"},
144 boost::histogram::axis::regular<double> {nybins, ymin, ymax, "y"},
145 boost::histogram::axis::regular<double> {nzbins, zmin, zmax, "z"},
146 boost::histogram::axis::variable<double> {eBinEdgesIn, "energy"});
147
148 h_err = boost::histogram::make_histogram_with(std::vector<double>(),
149 boost::histogram::axis::regular<double> {nxbins, xmin, xmax, "x"},
150 boost::histogram::axis::regular<double> {nybins, ymin, ymax, "y"},
151 boost::histogram::axis::regular<double> {nzbins, zmin, zmax, "z"},
152 boost::histogram::axis::variable<double> {eBinEdgesIn, "energy"});
153}
154
155template<class T>
157{;}
158
159template <class T>
161{
162 h.reset();
163 h_entries = 0;
164}
165
166template <class T>
167BDSBH4D<T>* BDSBH4D<T>::Clone(const char* newname) const
168{
169 auto clone = new BDSBH4D<T>(*this);
170 clone->SetName(newname);
171 return clone;
172}
173
174template <class T>
175void BDSBH4D<T>::Fill_BDSBH4D(double xValue,
176 double yValue,
177 double zValue,
178 double eValue)
179{
180 h(xValue, yValue, zValue, eValue);
181}
182
183template <class T>
184void BDSBH4D<T>::Set_BDSBH4D(int x,
185 int y,
186 int z,
187 int e,
188 double value)
189{
190 h.at(x, y, z, e) = value;
191}
192
193template <class T>
195 int y,
196 int z,
197 int e,
198 double value)
199{
200 h_err.at(x, y, z, e) = value;
201}
202
203template <class T>
204void BDSBH4D<T>::Add_BDSBH4D(BDSBH4DBase* otherHistogram)
205{
206 auto tmp = dynamic_cast<BDSBH4D<T>*>(otherHistogram);
207 h += tmp->h;
208}
209
210template <class T>
211double BDSBH4D<T>::At(int x, int y, int z, int e)
212{
213 return h.at(x, y, z, e);
214}
215
216template <class T>
217double BDSBH4D<T>::LowBinEdgeAt(int x, int y, int z, int e)
218{
219 std::ostringstream os4;
220 for (auto&& i : indexed(h))
221 {
222 if (i.index(0) == x and i.index(1) == y and i.index(2) == z and i.index(3) == e)
223 {return i.bin(3).lower();}
224 }
225 return 0;
226}
227
228template <class T>
229double BDSBH4D<T>::HighBinEdgeAt(int x, int y, int z, int e)
230{
231 std::ostringstream os4;
232 for (auto&& i : indexed(h))
233 {
234 if (i.index(0) == x and i.index(1) == y and i.index(2) == z and i.index(3) == e)
235 {return i.bin(3).upper();}
236 }
237 return 0;
238}
239
240template <class T>
241double BDSBH4D<T>::AtError(int x, int y, int z, int e)
242{
243 return h_err.at(x, y, z, e);
244}
245
246template <class T>
247void BDSBH4D<T>::Print_BDSBH4D(bool with_zero_values)
248{
249 std::ostringstream os4;
250 for (auto&& x : indexed(this->h))
251 {
252 if (!with_zero_values and *x != 0)
253 {
254 os4 << boost::format("(%i, %i, %i, %i) [%.3f, %.3f) [%.3f, %.3f) [%.3f, %.3f) [%.3f, %.3f) %i\n")
255 % x.index(0) % x.index(1) % x.index(2) % x.index(3)
256 % x.bin(0).lower() % x.bin(0).upper()
257 % x.bin(1).lower() % x.bin(1).upper()
258 % x.bin(2).lower() % x.bin(2).upper()
259 % x.bin(3).lower() % x.bin(3).upper()
260 % *x;
261 }
262 else if (with_zero_values)
263 {
264 os4 << boost::format("(%i, %i, %i, %i) [%.3f, %.3f) [%.3f, %.3f) [%.3f, %.3f) [%.3f, %.3f) %i\n")
265 % x.index(0) % x.index(1) % x.index(2) % x.index(3)
266 % x.bin(0).lower() % x.bin(0).upper()
267 % x.bin(1).lower() % x.bin(1).upper()
268 % x.bin(2).lower() % x.bin(2).upper()
269 % x.bin(3).lower() % x.bin(3).upper()
270 % *x;
271 }
272 }
273 std::cout << os4.str() << std::flush;
274}
275
276template <class T>
277void BDSBH4D<T>::Print_BDSBH4D(int x, int y, int z, int e)
278{
279 std::ostringstream os4;
280 for (auto&& i : indexed(this->h))
281 {
282 if (i.index(0) == x and i.index(1) == y and i.index(2) ==z and i.index(3) ==e)
283 {
284 os4 << boost::format("(%i, %i, %i, %i) [%.3f, %.3f) [%.3f, %.3f) [%.3f, %.3f) [%.3f, %.3f) %i\n")
285 % i.index(0) % i.index(1) % i.index(2) % i.index(3)
286 % i.bin(0).lower() % i.bin(0).upper()
287 % i.bin(1).lower() % i.bin(1).upper()
288 % i.bin(2).lower() % i.bin(2).upper()
289 % i.bin(3).lower() % i.bin(3).upper()
290 % *i;
291 }
292 }
293 std::cout << os4.str() << std::flush;
294}
295
297template class BDSBH4D<boost_histogram_log>;
299
300#endif
Base class for the 4D histogram classes.
Definition: BDSBH4DBase.hh:37
4D histogram classes with linear, logarithmic and user-defined energy binning.
Definition: BDSBH4D.hh:41