REST-for-Physics  v2.3
Rare Event Searches ToolKit for Physics
TRestGeant4PrimaryGeneratorInfo.cxx
1//
2// Created by lobis on 7/14/2022.
3//
4
5#include "TRestGeant4PrimaryGeneratorInfo.h"
6
7#include <TAxis.h>
8#include <TF1.h>
9#include <TF2.h>
10#include <TMath.h>
11#include <TRestStringOutput.h>
12
13#include <iostream>
14
15using namespace std;
16using namespace TRestGeant4PrimaryGeneratorTypes;
17
18string TRestGeant4PrimaryGeneratorTypes::SpatialGeneratorTypesToString(const SpatialGeneratorTypes& type) {
19 switch (type) {
20 case SpatialGeneratorTypes::CUSTOM:
21 return "Custom";
22 case SpatialGeneratorTypes::VOLUME:
23 return "Volume";
24 case SpatialGeneratorTypes::SURFACE:
25 return "Surface";
26 case SpatialGeneratorTypes::POINT:
27 return "Point";
28 case SpatialGeneratorTypes::COSMIC:
29 return "Cosmic";
30 case SpatialGeneratorTypes::SOURCE:
31 return "Source";
32 }
33 cout << "TRestGeant4PrimaryGeneratorTypes::SpatialGeneratorTypesToString - Error - Unknown "
34 "SpatialGeneratorTypes"
35 << endl;
36 exit(1);
37}
38
39SpatialGeneratorTypes TRestGeant4PrimaryGeneratorTypes::StringToSpatialGeneratorTypes(const string& type) {
40 if (TString(type).EqualTo(SpatialGeneratorTypesToString(SpatialGeneratorTypes::CUSTOM),
41 TString::ECaseCompare::kIgnoreCase)) {
42 return SpatialGeneratorTypes::CUSTOM;
43 } else if (TString(type).EqualTo(SpatialGeneratorTypesToString(SpatialGeneratorTypes::VOLUME),
44 TString::ECaseCompare::kIgnoreCase)) {
45 return SpatialGeneratorTypes::VOLUME;
46 } else if (TString(type).EqualTo(SpatialGeneratorTypesToString(SpatialGeneratorTypes::SURFACE),
47 TString::ECaseCompare::kIgnoreCase)) {
48 return SpatialGeneratorTypes::SURFACE;
49 } else if (TString(type).EqualTo(SpatialGeneratorTypesToString(SpatialGeneratorTypes::POINT),
50 TString::ECaseCompare::kIgnoreCase)) {
51 return SpatialGeneratorTypes::POINT;
52 } else if (TString(type).EqualTo(SpatialGeneratorTypesToString(SpatialGeneratorTypes::COSMIC),
53 TString::ECaseCompare::kIgnoreCase)) {
54 return SpatialGeneratorTypes::COSMIC;
55 } else if (TString(type).EqualTo(SpatialGeneratorTypesToString(SpatialGeneratorTypes::SOURCE),
56 TString::ECaseCompare::kIgnoreCase)) {
57 return SpatialGeneratorTypes::SOURCE;
58 } else {
59 cout << "TRestGeant4PrimaryGeneratorTypes::StringToSpatialGeneratorTypes - Error - Unknown "
60 "SpatialGeneratorTypes: "
61 << type << endl;
62 exit(1);
63 }
64}
65
66string TRestGeant4PrimaryGeneratorTypes::SpatialGeneratorShapesToString(const SpatialGeneratorShapes& type) {
67 switch (type) {
68 case SpatialGeneratorShapes::GDML:
69 return "GDML";
70 case SpatialGeneratorShapes::WALL:
71 return "Wall";
72 case SpatialGeneratorShapes::CIRCLE:
73 return "Circle";
74 case SpatialGeneratorShapes::BOX:
75 return "Box";
76 case SpatialGeneratorShapes::SPHERE:
77 return "Sphere";
78 case SpatialGeneratorShapes::CYLINDER:
79 return "Cylinder";
80 }
81 cout << "TRestGeant4PrimaryGeneratorTypes::SpatialGeneratorShapesToString - Error - Unknown "
82 "SpatialGeneratorShapes"
83 << endl;
84 exit(1);
85}
86
87SpatialGeneratorShapes TRestGeant4PrimaryGeneratorTypes::StringToSpatialGeneratorShapes(const string& type) {
88 if (TString(type).EqualTo(SpatialGeneratorShapesToString(SpatialGeneratorShapes::GDML),
89 TString::ECaseCompare::kIgnoreCase)) {
90 return SpatialGeneratorShapes::GDML;
91 } else if (TString(type).EqualTo(SpatialGeneratorShapesToString(SpatialGeneratorShapes::WALL),
92 TString::ECaseCompare::kIgnoreCase)) {
93 return SpatialGeneratorShapes::WALL;
94 } else if (TString(type).EqualTo(SpatialGeneratorShapesToString(SpatialGeneratorShapes::CIRCLE),
95 TString::ECaseCompare::kIgnoreCase)) {
96 return SpatialGeneratorShapes::CIRCLE;
97 } else if (TString(type).EqualTo(SpatialGeneratorShapesToString(SpatialGeneratorShapes::BOX),
98 TString::ECaseCompare::kIgnoreCase)) {
99 return SpatialGeneratorShapes::BOX;
100 } else if (TString(type).EqualTo(SpatialGeneratorShapesToString(SpatialGeneratorShapes::SPHERE),
101 TString::ECaseCompare::kIgnoreCase)) {
102 return SpatialGeneratorShapes::SPHERE;
103 } else if (TString(type).EqualTo(SpatialGeneratorShapesToString(SpatialGeneratorShapes::CYLINDER),
104 TString::ECaseCompare::kIgnoreCase)) {
105 return SpatialGeneratorShapes::CYLINDER;
106
107 } else {
108 cout << "TRestGeant4PrimaryGeneratorTypes::StringToSpatialGeneratorShapes - Error - Unknown "
109 "SpatialGeneratorShapes: "
110 << type << endl;
111 exit(1);
112 }
113}
114
115string TRestGeant4PrimaryGeneratorTypes::EnergyDistributionTypesToString(
116 const EnergyDistributionTypes& type) {
117 switch (type) {
118 case EnergyDistributionTypes::TH1D:
119 return "TH1D";
120 case EnergyDistributionTypes::TH2D:
121 return "TH2D";
122 case EnergyDistributionTypes::FORMULA:
123 return "Formula";
124 case EnergyDistributionTypes::FORMULA2:
125 return "Formula2";
126 case EnergyDistributionTypes::MONO:
127 return "Mono";
128 case EnergyDistributionTypes::FLAT:
129 return "Flat";
130 case EnergyDistributionTypes::LOG:
131 return "Log";
132 }
133 cout << "TRestGeant4PrimaryGeneratorTypes::EnergyDistributionTypesToString - Error - Unknown energy "
134 "distribution type"
135 << endl;
136 exit(1);
137}
138
139EnergyDistributionTypes TRestGeant4PrimaryGeneratorTypes::StringToEnergyDistributionTypes(
140 const string& type) {
141 if (TString(type).EqualTo(EnergyDistributionTypesToString(EnergyDistributionTypes::TH1D),
142 TString::ECaseCompare::kIgnoreCase)) {
143 return EnergyDistributionTypes::TH1D;
144 } else if (TString(type).EqualTo(EnergyDistributionTypesToString(EnergyDistributionTypes::TH2D),
145 TString::ECaseCompare::kIgnoreCase)) {
146 return EnergyDistributionTypes::TH2D;
147 } else if (TString(type).EqualTo(EnergyDistributionTypesToString(EnergyDistributionTypes::FORMULA),
148 TString::ECaseCompare::kIgnoreCase)) {
149 return EnergyDistributionTypes::FORMULA;
150 } else if (TString(type).EqualTo(EnergyDistributionTypesToString(EnergyDistributionTypes::FORMULA2),
151 TString::ECaseCompare::kIgnoreCase)) {
152 return EnergyDistributionTypes::FORMULA2;
153 } else if (TString(type).EqualTo(EnergyDistributionTypesToString(EnergyDistributionTypes::MONO),
154 TString::ECaseCompare::kIgnoreCase)) {
155 return EnergyDistributionTypes::MONO;
156 } else if (TString(type).EqualTo(EnergyDistributionTypesToString(EnergyDistributionTypes::FLAT),
157 TString::ECaseCompare::kIgnoreCase)) {
158 return EnergyDistributionTypes::FLAT;
159 } else if (TString(type).EqualTo(EnergyDistributionTypesToString(EnergyDistributionTypes::LOG),
160 TString::ECaseCompare::kIgnoreCase)) {
161 return EnergyDistributionTypes::LOG;
162 } else {
163 cout << "TRestGeant4PrimaryGeneratorTypes::StringToEnergyDistributionTypes - Error - Unknown "
164 "EnergyDistributionTypes: "
165 << type << endl;
166 exit(1);
167 }
168}
169
170string TRestGeant4PrimaryGeneratorTypes::EnergyDistributionFormulasToString(
171 const EnergyDistributionFormulas& type) {
172 switch (type) {
173 case EnergyDistributionFormulas::COSMIC_NEUTRONS:
174 return "CosmicNeutrons";
175 case EnergyDistributionFormulas::COSMIC_GAMMAS:
176 return "CosmicGammas";
177 case EnergyDistributionFormulas::FISSION_NEUTRONS_U238:
178 return "FissionNeutronsU238";
179 case EnergyDistributionFormulas::ENVIRONMENTAL_GAMMAS:
180 return "EnvironmentalGammas";
181 case EnergyDistributionFormulas::ENVIRONMENTAL_NEUTRONS:
182 return "EnvironmentalNeutrons";
183 }
184 cout << "TRestGeant4PrimaryGeneratorTypes::EnergyDistributionFormulasToString - Error - Unknown energy "
185 "distribution formula"
186 << endl;
187 exit(1);
188}
189
190EnergyDistributionFormulas TRestGeant4PrimaryGeneratorTypes::StringToEnergyDistributionFormulas(
191 const string& type) {
192 if (TString(type).EqualTo(EnergyDistributionFormulasToString(EnergyDistributionFormulas::COSMIC_NEUTRONS),
193 TString::ECaseCompare::kIgnoreCase)) {
194 return EnergyDistributionFormulas::COSMIC_NEUTRONS;
195 } else if (TString(type).EqualTo(
196 EnergyDistributionFormulasToString(EnergyDistributionFormulas::COSMIC_GAMMAS),
197 TString::ECaseCompare::kIgnoreCase)) {
198 return EnergyDistributionFormulas::COSMIC_GAMMAS;
199 } else if (TString(type).EqualTo(
200 EnergyDistributionFormulasToString(EnergyDistributionFormulas::FISSION_NEUTRONS_U238),
201 TString::ECaseCompare::kIgnoreCase)) {
202 return EnergyDistributionFormulas::FISSION_NEUTRONS_U238;
203 } else if (TString(type).EqualTo(
204 EnergyDistributionFormulasToString(EnergyDistributionFormulas::ENVIRONMENTAL_GAMMAS),
205 TString::ECaseCompare::kIgnoreCase)) {
206 return EnergyDistributionFormulas::ENVIRONMENTAL_GAMMAS;
207 } else if (TString(type).EqualTo(
208 EnergyDistributionFormulasToString(EnergyDistributionFormulas::ENVIRONMENTAL_NEUTRONS),
209 TString::ECaseCompare::kIgnoreCase)) {
210 return EnergyDistributionFormulas::ENVIRONMENTAL_NEUTRONS;
211 } else {
212 cout << "TRestGeant4PrimaryGeneratorTypes::StringToEnergyDistributionFormulas - Error - Unknown "
213 "energyDistributionFormulas: "
214 << type << endl;
215 exit(1);
216 }
217}
218
219TF1 TRestGeant4PrimaryGeneratorTypes::EnergyDistributionFormulasToRootFormula(
220 const EnergyDistributionFormulas& type) {
221 switch (type) {
222 case EnergyDistributionFormulas::COSMIC_NEUTRONS: {
223 // Formula from https://ieeexplore.ieee.org/document/1369506
224 const char* title = "Cosmic Neutrons at Sea Level";
225 auto distribution =
226 TF1(title,
227 "1.006E-6 * TMath::Exp(-0.3500 * TMath::Power(TMath::Log(x * 1E-3), 2) + 2.1451 * "
228 "TMath::Log(x * 1E-3)) + "
229 "1.011E-3 * TMath::Exp(-0.4106 * TMath::Power(TMath::Log(x * 1E-3), 2) - 0.6670 * "
230 "TMath::Log(x * 1E-3))",
231 1E2, 1E7);
232 distribution.SetNormalized(true); // Normalized, not really necessary
233 distribution.SetTitle(title);
234 distribution.GetXaxis()->SetTitle("Energy (keV)");
235 return distribution;
236 }
237 case EnergyDistributionFormulas::COSMIC_GAMMAS:
238 exit(1);
239 case EnergyDistributionFormulas::FISSION_NEUTRONS_U238: {
240 // Watt formula for U238. The spectrum parameters are from this reference
241 // https://arxiv.org/pdf/hep-ex/0312050
242 const char* title = "Watt formula: Neutrons from U238 fission";
243 auto distribution =
244 TF1(title, "TMath::Exp(-x/712.4) * TMath::SinH(TMath::Sqrt(0.0056405*x))", 0, 10000); // keV
245 distribution.SetNormalized(true);
246 distribution.SetTitle(title);
247 distribution.GetXaxis()->SetTitle("Energy (keV)");
248 return distribution;
249 }
250 case EnergyDistributionFormulas::ENVIRONMENTAL_GAMMAS: {
251 // Environmental gamma radiation approximation
252 const char* title = "Environmental Gammas";
253 auto distribution = TF1(title, "(1 / TMath::Power(x, 1.5)) * TMath::Exp(-x / 1500.)", 250, 10000);
254 distribution.SetNormalized(true);
255 distribution.SetTitle(title);
256 distribution.GetXaxis()->SetTitle("Energy (keV)");
257 return distribution;
258 }
259 case EnergyDistributionFormulas::ENVIRONMENTAL_NEUTRONS: {
260 // Environmental neutron radiation approximation
261 const char* title = "Environmental Neutrons";
262 auto distribution = TF1(title, "TMath::Exp(-x * 0.0005)", 0, 10000); // keV
263 distribution.SetNormalized(true);
264 distribution.SetTitle(title);
265 distribution.GetXaxis()->SetTitle("Energy (keV)");
266 return distribution;
267 }
268 }
269 cout << "TRestGeant4PrimaryGeneratorTypes::EnergyDistributionFormulasToRootFormula - Error - Unknown "
270 "energy distribution formula"
271 << endl;
272 exit(1);
273}
274
275string TRestGeant4PrimaryGeneratorTypes::AngularDistributionTypesToString(
276 const AngularDistributionTypes& type) {
277 switch (type) {
278 case AngularDistributionTypes::TH1D:
279 return "TH1D";
280 case AngularDistributionTypes::TH2D:
281 return "TH2D";
282 case AngularDistributionTypes::FORMULA:
283 return "Formula";
284 case AngularDistributionTypes::FORMULA2:
285 return "Formula2";
286 case AngularDistributionTypes::ISOTROPIC:
287 return "Isotropic";
288 case AngularDistributionTypes::FLUX:
289 return "Flux";
290 case AngularDistributionTypes::BACK_TO_BACK:
291 return "Back to back";
292 }
293 cout << "TRestGeant4PrimaryGeneratorTypes::AngularDistributionTypesToString - Error - Unknown angular "
294 "distribution type"
295 << endl;
296 exit(1);
297}
298
299AngularDistributionTypes TRestGeant4PrimaryGeneratorTypes::StringToAngularDistributionTypes(
300 const string& type) {
301 if (TString(type).EqualTo(AngularDistributionTypesToString(AngularDistributionTypes::TH1D),
302 TString::ECaseCompare::kIgnoreCase)) {
303 return AngularDistributionTypes::TH1D;
304 } else if (TString(type).EqualTo(AngularDistributionTypesToString(AngularDistributionTypes::TH2D),
305 TString::ECaseCompare::kIgnoreCase)) {
306 return AngularDistributionTypes::TH2D;
307 } else if (TString(type).EqualTo(AngularDistributionTypesToString(AngularDistributionTypes::FORMULA),
308 TString::ECaseCompare::kIgnoreCase)) {
309 return AngularDistributionTypes::FORMULA;
310 } else if (TString(type).EqualTo(AngularDistributionTypesToString(AngularDistributionTypes::FORMULA2),
311 TString::ECaseCompare::kIgnoreCase)) {
312 return AngularDistributionTypes::FORMULA2;
313 } else if (TString(type).EqualTo(AngularDistributionTypesToString(AngularDistributionTypes::ISOTROPIC),
314 TString::ECaseCompare::kIgnoreCase)) {
315 return AngularDistributionTypes::ISOTROPIC;
316 } else if (TString(type).EqualTo(AngularDistributionTypesToString(AngularDistributionTypes::FLUX),
317 TString::ECaseCompare::kIgnoreCase)) {
318 return AngularDistributionTypes::FLUX;
319 } else if (TString(type).EqualTo(AngularDistributionTypesToString(AngularDistributionTypes::BACK_TO_BACK),
320 TString::ECaseCompare::kIgnoreCase)) {
321 return AngularDistributionTypes::BACK_TO_BACK;
322 } else {
323 cout << "TRestGeant4PrimaryGeneratorTypes::StringToAngularDistributionTypes - Error - Unknown "
324 "AngularDistributionTypes: "
325 << type << endl;
326 exit(1);
327 }
328}
329
330string TRestGeant4PrimaryGeneratorTypes::AngularDistributionFormulasToString(
331 const AngularDistributionFormulas& type) {
332 switch (type) {
333 case AngularDistributionFormulas::COS2:
334 return "Cos2";
335 case AngularDistributionFormulas::COS3:
336 return "Cos3";
337 case AngularDistributionFormulas::SIN_COS2:
338 return "SinCos2";
339 case AngularDistributionFormulas::SIN_2THETA:
340 return "Sin2theta";
341 }
342 cout << "TRestGeant4PrimaryGeneratorTypes::AngularDistributionFormulasToString - Error - Unknown angular "
343 "distribution formula"
344 << endl;
345 exit(1);
346}
347
348AngularDistributionFormulas TRestGeant4PrimaryGeneratorTypes::StringToAngularDistributionFormulas(
349 const string& type) {
350 if (TString(type).EqualTo(AngularDistributionFormulasToString(AngularDistributionFormulas::COS2),
351 TString::ECaseCompare::kIgnoreCase)) {
352 return AngularDistributionFormulas::COS2;
353 } else if (TString(type).EqualTo(AngularDistributionFormulasToString(AngularDistributionFormulas::COS3),
354 TString::ECaseCompare::kIgnoreCase)) {
355 return AngularDistributionFormulas::COS3;
356 } else if (TString(type).EqualTo(
357 AngularDistributionFormulasToString(AngularDistributionFormulas::SIN_COS2),
358 TString::ECaseCompare::kIgnoreCase)) {
359 return AngularDistributionFormulas::SIN_COS2;
360 } else if (TString(type).EqualTo(
361 AngularDistributionFormulasToString(AngularDistributionFormulas::SIN_2THETA),
362 TString::ECaseCompare::kIgnoreCase)) {
363 return AngularDistributionFormulas::SIN_2THETA;
364 } else {
365 cout << "TRestGeant4PrimaryGeneratorTypes::StringToAngularDistributionFormulas - Error - Unknown "
366 "AngularDistributionFormulas: "
367 << type << endl;
368 exit(1);
369 }
370}
371
372TF1 TRestGeant4PrimaryGeneratorTypes::AngularDistributionFormulasToRootFormula(
373 const AngularDistributionFormulas& type) {
374 switch (type) {
375 case AngularDistributionFormulas::COS2: {
376 auto cos2 = [](double* xs, double* ps) {
377 if (xs[0] >= 0 && xs[0] <= TMath::Pi() / 2) {
378 return TMath::Power(TMath::Cos(xs[0]), 2);
379 }
380 return 0.0;
381 };
382 const char* title = "AngularDistribution: Cos2";
383 auto f = TF1(title, cos2, 0.0, TMath::Pi());
384 f.SetTitle(title);
385 return f;
386 }
387 case AngularDistributionFormulas::COS3: {
388 auto cos3 = [](double* xs, double* ps) {
389 if (xs[0] >= 0 && xs[0] <= TMath::Pi() / 2) {
390 return TMath::Power(TMath::Cos(xs[0]), 3);
391 }
392 return 0.0;
393 };
394 const char* title = "AngularDistribution: Cos3";
395 auto f = TF1(title, cos3, 0.0, TMath::Pi());
396 f.SetTitle(title);
397 return f;
398 }
399 case AngularDistributionFormulas::SIN_COS2: {
400 auto sinCos2 = [](double* xs, double* ps) {
401 if (xs[0] >= 0 && xs[0] <= TMath::Pi() / 2) {
402 return TMath::Sin(xs[0]) * TMath::Power(TMath::Cos(xs[0]), 2);
403 }
404 return 0.0;
405 };
406 const char* title = "AngularDistribution: SinCos2";
407 auto f = TF1(title, sinCos2, 0.0, TMath::Pi());
408 f.SetTitle(title);
409 return f;
410 }
411 case AngularDistributionFormulas::SIN_2THETA: {
412 auto sin2theta = [](double* xs, double* ps) {
413 if (xs[0] >= 0 && xs[0] <= TMath::Pi() / 2) {
414 return TMath::Sin(2 * xs[0]);
415 }
416 return 0.0;
417 };
418 const char* title = "AngularDistribution: Sin2theta";
419 auto f = TF1(title, sin2theta, 0.0, TMath::Pi());
420 f.SetTitle(title);
421 return f;
422 }
423 }
424 cout << "TRestGeant4PrimaryGeneratorTypes::AngularDistributionFormulasToRootFormula - Error - Unknown "
425 "angular distribution formula"
426 << endl;
427 exit(1);
428}
429
430string TRestGeant4PrimaryGeneratorTypes::EnergyAndAngularDistributionFormulasToString(
431 const EnergyAndAngularDistributionFormulas& type) {
432 switch (type) {
433 case EnergyAndAngularDistributionFormulas::COSMIC_MUONS:
434 return "CosmicMuons";
435 }
436 cout << "TRestGeant4PrimaryGeneratorTypes::EnergyAndAngularDistributionFormulasToString - Error - "
437 "Unknown energy/angular distribution formula"
438 << endl;
439 exit(1);
440}
441
442EnergyAndAngularDistributionFormulas
443TRestGeant4PrimaryGeneratorTypes::StringToEnergyAndAngularDistributionFormulas(const string& type) {
444 if (TString(type).EqualTo(
445 EnergyAndAngularDistributionFormulasToString(EnergyAndAngularDistributionFormulas::COSMIC_MUONS),
446 TString::ECaseCompare::kIgnoreCase)) {
447 return EnergyAndAngularDistributionFormulas::COSMIC_MUONS;
448 } else {
449 cout << "TRestGeant4PrimaryGeneratorTypes::StringToEnergyAndAngularDistributionFormulas - Error - "
450 "Unknown AngularDistributionFormulas: "
451 << type << endl;
452 exit(1);
453 }
454}
455
456TF2 TRestGeant4PrimaryGeneratorTypes::EnergyAndAngularDistributionFormulasToRootFormula(
457 const EnergyAndAngularDistributionFormulas& type) {
458 switch (type) {
459 case EnergyAndAngularDistributionFormulas::COSMIC_MUONS: {
460 // Guan formula from https://arxiv.org/pdf/1509.06176.pdf
461 // muon rest mass is 105.7 MeV
462 // energy in keV
463 // Original Guan formula is differential in solid angle. The final factor integrates over
464 // azimuth and converts dOmega to 2*pi*sin(theta)*dtheta, so do not add another 2*pi here.
465 const char* title = "Cosmic Muons Energy and Angular";
466 auto f =
467 TF2(title,
468 "1E-6*0.14*TMath::Power(x*1E-6*(1.+3.64/"
469 "(x*1E-6*TMath::Power(TMath::Power((TMath::Power(TMath::Cos(y),2)+0.0105212-0.068287*"
470 "TMath::Power(TMath::Cos(y),0.958633)+0.0407253*TMath::Power(TMath::Cos(y),0.817285)"
471 ")/(0.982960),0.5),1.29))),-2.7)*(1./"
472 "(1.+(1.1*x*1E-6*TMath::Power((TMath::Power(TMath::Cos(y),2)+0.0105212-0.068287*TMath::"
473 "Power(TMath::Cos(y),0.958633)+0.0407253*TMath::Power(TMath::Cos(y),0.817285))/"
474 "(0.982960),0.5))/115.)+0.054/"
475 "(1.+(1.1*x*1E-6*TMath::Power((TMath::Power(TMath::Cos(y),2)+0.0105212-0.068287*TMath::"
476 "Power(TMath::Cos(y),0.958633)+0.0407253*TMath::Power(TMath::Cos(y),0.817285))/"
477 "(0.982960),0.5))/850.))*(2.*TMath::Sin(y)*TMath::Pi())",
478 2.0E5, 5.0E9, 0, TMath::Pi() / 2.);
479 f.SetTitle(title);
480 /*
481 * we need to increase the number of bins to get a smooth distribution
482 * this makes the first random generation slow, but it is only a constant factor which does not
483 * matter in long simulations
484 * If we don't include this, the distribution will be very discrete. Max value is 10000
485 */
486 return f;
487 }
488 }
489 cout << "TRestGeant4PrimaryGeneratorTypes::EnergyAndAngularDistributionFormulasToRootFormula - Error - "
490 "Unknown energy and angular distribution formula"
491 << endl;
492 exit(1);
493}
494
495void TRestGeant4PrimaryGeneratorInfo::Print() const {
496 const auto typeEnum = StringToSpatialGeneratorTypes(fSpatialGeneratorType.Data());
497 RESTMetadata << "Generator type: " << SpatialGeneratorTypesToString(typeEnum) << RESTendl;
498
499 if (typeEnum != SpatialGeneratorTypes::POINT) {
500 const auto shapeEnum = StringToSpatialGeneratorShapes(fSpatialGeneratorShape.Data());
501 RESTMetadata << "Generator shape: " << SpatialGeneratorShapesToString(shapeEnum);
502 if (shapeEnum == SpatialGeneratorShapes::GDML) {
503 RESTMetadata << "::" << fSpatialGeneratorFrom << RESTendl;
504 } else {
505 if (shapeEnum == SpatialGeneratorShapes::BOX) {
506 RESTMetadata << ", (length, width, height): ";
507 } else if (shapeEnum == SpatialGeneratorShapes::SPHERE) {
508 RESTMetadata << ", (radius, , ): ";
509 } else if (shapeEnum == SpatialGeneratorShapes::WALL) {
510 RESTMetadata << ", (length, width, ): ";
511 } else if (shapeEnum == SpatialGeneratorShapes::CIRCLE) {
512 RESTMetadata << ", (radius, , ): ";
513 } else if (shapeEnum == SpatialGeneratorShapes::CYLINDER) {
514 RESTMetadata << ", (radius, length, ): ";
515 }
516
517 RESTMetadata << fSpatialGeneratorSize.X() << ", " << fSpatialGeneratorSize.Y() << ", "
518 << fSpatialGeneratorSize.Z() << RESTendl;
519 }
520 }
521 RESTMetadata << "Generator center : (" << fSpatialGeneratorPosition.X() << ","
522 << fSpatialGeneratorPosition.Y() << "," << fSpatialGeneratorPosition.Z() << ") mm"
523 << RESTendl;
524 RESTMetadata << "Generator rotation : (" << fSpatialGeneratorRotationAxis.X() << ","
526 << "), angle: " << fSpatialGeneratorRotationValue * TMath::RadToDeg() << "ยบ" << RESTendl;
527}
TVector3 fSpatialGeneratorSize
The size of the generator. Stores up to three dimensions according to the shape box: length,...
TString fSpatialGeneratorFrom
The volume name where the events are generated, in case of volume or surface generator types.
Double_t fSpatialGeneratorRotationValue
degrees of rotation for generator virtual shape along the axis
TVector3 fSpatialGeneratorPosition
The position of the generator for virtual, and point, generator types.
TVector3 fSpatialGeneratorRotationAxis
A 3d-std::vector with the angles, measured in degrees, of a XYZ rotation applied to the virtual gener...
TString fSpatialGeneratorShape
Shape of spatial generator (wall, GDML, sphere, etc)
TString fSpatialGeneratorType
Type of spatial generator (point, surface, volume, custom)