Newer
Older
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
//==========================================================================
// AIDA Detector description implementation
//--------------------------------------------------------------------------
// Copyright (C) Organisation europeenne pour la Recherche nucleaire (CERN)
// All rights reserved.
//
// For the licensing terms see $DD4hepINSTALL/LICENSE.
// For the list of contributors see $DD4hepINSTALL/doc/CREDITS.
//
// Author : M.Frank
//
//==========================================================================
// Framework include files
#include "DDG4/Geant4UserLimits.h"
#include "DDG4/Geant4Particle.h"
#include "DD4hep/InstanceCount.h"
#include "DD4hep/DD4hepUnits.h"
#include "DD4hep/Primitives.h"
// Geant 4 include files
#include "G4Track.hh"
#include "CLHEP/Units/SystemOfUnits.h"
// C/C++ include files
#include <stdexcept>
using namespace std;
using namespace dd4hep::sim;
/// Access value according to track
double Geant4UserLimits::Handler::value(const G4Track& track) const {
if ( !particleLimits.empty() ) {
auto i = particleLimits.find(track.GetDefinition());
if ( i != particleLimits.end() ) {
return (*i).second;
}
}
return defaultValue;
}
/// Set the handler value(s)
void Geant4UserLimits::Handler::set(const string& particles, double val) {
if ( particles == "*" ) {
defaultValue = val;
return;
}
auto defs = Geant4ParticleHandle::g4DefinitionsRegEx(particles);
for(auto* d : defs)
particleLimits[d] = val;
}
/// Initializing Constructor
Geant4UserLimits::Geant4UserLimits(LimitSet ls)
: G4UserLimits(ls.name()), limits(ls)
{
const auto& lim = limits.limits();
InstanceCount::increment(this);
/// Set defaults
maxStepLength.defaultValue = fMaxStep;
maxTrackLength.defaultValue = fMaxTrack;
maxTime.defaultValue = fMaxTime;
minEKine.defaultValue = fMinEkine;
minRange.defaultValue = fMinRange;
/// Overwrite with values if present:
for(const Limit& l : lim) {
if (l.name == "step_length_max")
maxStepLength.set(l.particles, l.value*CLHEP::mm/dd4hep::mm);
else if (l.name == "track_length_max")
maxTrackLength.set(l.particles, l.value*CLHEP::mm/dd4hep::mm);
else if (l.name == "time_max")
maxTime.set(l.particles, l.value*CLHEP::ns/dd4hep::ns);
else if (l.name == "ekin_min")
minEKine.set(l.particles, l.value*CLHEP::MeV/dd4hep::MeV);
else if (l.name == "range_min")
minRange.set(l.particles, l.value);
else
throw runtime_error("Unknown Geant4 user limit: " + l.toString());
}
}
/// Standard destructor
Geant4UserLimits::~Geant4UserLimits() {
InstanceCount::decrement(this);
}
/// Setters may not be called!
void Geant4UserLimits::SetMaxAllowedStep(G4double /* ustepMax */) {
dd4hep::notImplemented(string(__PRETTY_FUNCTION__)+" May not be called!");
}
void Geant4UserLimits::SetUserMaxTrackLength(G4double /* utrakMax */) {
dd4hep::notImplemented(string(__PRETTY_FUNCTION__)+" May not be called!");
}
void Geant4UserLimits::SetUserMaxTime(G4double /* utimeMax */) {
dd4hep::notImplemented(string(__PRETTY_FUNCTION__)+" May not be called!");
}
void Geant4UserLimits::SetUserMinEkine(G4double /* uekinMin */) {
dd4hep::notImplemented(string(__PRETTY_FUNCTION__)+" May not be called!");
}
void Geant4UserLimits::SetUserMinRange(G4double /* urangMin */) {
dd4hep::notImplemented(string(__PRETTY_FUNCTION__)+" May not be called!");
}