CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
FastTimerSD.cc
Go to the documentation of this file.
2 
4 
11 
15 
18 
19 #include <vector>
20 
21 #include "G4Track.hh"
22 #include "G4Step.hh"
23 #include "G4StepPoint.hh"
24 #include "G4VPhysicalVolume.hh"
25 
26 #include <iostream>
27 
28 //#define EDM_ML_DEBUG
29 //-------------------------------------------------------------------
31  const edm::EventSetup& es,
32  const SensitiveDetectorCatalog& clg,
33  edm::ParameterSet const& p,
34  const SimTrackManager* manager)
35  : TimingSD(name, es, clg, p, manager), ftcons(nullptr) {
36  //Parameters
37  edm::ParameterSet m_p = p.getParameter<edm::ParameterSet>("FastTimerSD");
38  int verbn = m_p.getUntrackedParameter<int>("Verbosity");
39 
40  SetVerboseLevel(verbn);
41 
43  es.get<IdealGeometryRecord>().get(cpv);
44 
45  std::string attribute = "ReadOutName";
46  DDSpecificsMatchesValueFilter filter{DDValue(attribute, name, 0)};
47  DDFilteredView fv(*cpv, filter);
48  fv.firstChild();
50  std::vector<int> temp = dbl_to_int(getDDDArray("Type", sv));
51  type_ = temp[0];
52 
53  setTimeFactor(100.);
54 
55  edm::LogInfo("FastTimerSim") << "FastTimerSD: Instantiation completed for " << name << " of type " << type_;
56 }
57 
59 
60 uint32_t FastTimerSD::setDetUnitId(const G4Step* aStep) {
61  //Find the depth segment
62  const G4ThreeVector& global = getGlobalEntryPoint();
63  const G4ThreeVector& local = getLocalEntryPoint();
64  int iz = (global.z() > 0) ? 1 : -1;
65  std::pair<int, int> izphi = ((ftcons) ? ((type_ == 1) ? (ftcons->getZPhi(std::abs(local.z()), local.phi()))
66  : (ftcons->getEtaPhi(local.perp(), local.phi())))
67  : (std::pair<int, int>(0, 0)));
68  uint32_t id = FastTimeDetId(type_, izphi.first, izphi.second, iz).rawId();
69 #ifdef EDM_ML_DEBUG
70  edm::LogVerbatim("FastTimerSD") << "Volume " << aStep->GetPreStepPoint()->GetPhysicalVolume()->GetName() << ": "
71  << global.z() << " Iz(eta)phi " << izphi.first << ":" << izphi.second << ":" << iz
72  << " id " << std::hex << id << std::dec;
73 #endif
74  return id;
75 }
76 
77 void FastTimerSD::update(const BeginOfJob* job) {
78  const edm::EventSetup* es = (*job)();
80  es->get<IdealGeometryRecord>().get(fdc);
81  if (fdc.isValid()) {
82  ftcons = &(*fdc);
83  } else {
84  edm::LogError("FastTimerSim") << "FastTimerSD : Cannot find FastTimeDDDConstants";
85  throw cms::Exception("Unknown", "FastTimerSD") << "Cannot find FastTimeDDDConstants\n";
86  }
87 #ifdef EDM_ML_DEBUG
88  edm::LogVerbatim("FastTimerSD") << "FastTimerSD::Initialized with FastTimeDDDConstants\n";
89 #endif
90 }
91 
92 std::vector<double> FastTimerSD::getDDDArray(const std::string& str, const DDsvalues_type& sv) {
93  DDValue value(str);
94  if (DDfetch(&sv, value)) {
95  const std::vector<double>& fvec = value.doubles();
96  int nval = fvec.size();
97  if (nval < 1) {
98  edm::LogError("FastTimerSim") << "FastTimerSD : # of " << str << " bins " << nval << " < 1 ==> illegal";
99  throw cms::Exception("DDException") << "FastTimerSD: cannot get array " << str;
100  }
101  return fvec;
102  } else {
103  edm::LogError("FastTimerSim") << "FastTimerSD: cannot get array " << str;
104  throw cms::Exception("DDException") << "FastTimerSD: cannot get array " << str;
105  }
106 }
T getParameter(std::string const &) const
T getUntrackedParameter(std::string const &, T const &) const
const std::vector< double > & doubles() const
a reference to the double-valued values stored in the given instance of DDValue
Definition: DDValue.cc:111
void setTimeFactor(double)
Definition: TimingSD.cc:86
std::pair< int, int > getEtaPhi(double r, double phi) const
const G4ThreeVector & getLocalEntryPoint() const
Definition: TimingSD.h:66
#define nullptr
constexpr uint32_t rawId() const
get the raw id
Definition: DetId.h:57
bool DDfetch(const DDsvalues_type *, DDValue &)
helper for retrieving DDValues from DDsvalues_type *.
Definition: DDsvalues.cc:79
void update(const BeginOfJob *) override
This routine will be called when the appropriate signal arrives.
Definition: FastTimerSD.cc:77
const FastTimeDDDConstants * ftcons
Definition: FastTimerSD.h:39
std::vector< int > dbl_to_int(const std::vector< double > &vecdbl)
Converts a std::vector of doubles to a std::vector of int.
Definition: DDutils.h:7
std::pair< int, int > getZPhi(double z, double phi) const
std::vector< std::pair< unsigned int, DDValue > > DDsvalues_type
Definition: DDsvalues.h:12
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
uint32_t setDetUnitId(const G4Step *) override
Definition: FastTimerSD.cc:60
const G4ThreeVector & getGlobalEntryPoint() const
Definition: TimingSD.h:67
~FastTimerSD() override
Definition: FastTimerSD.cc:58
std::vector< double > getDDDArray(const std::string &, const DDsvalues_type &)
Definition: FastTimerSD.cc:92
DDsvalues_type mergedSpecifics() const
FastTimerSD(const std::string &, const edm::EventSetup &, const SensitiveDetectorCatalog &, edm::ParameterSet const &, const SimTrackManager *)
Definition: FastTimerSD.cc:30
T get() const
Definition: EventSetup.h:73
bool firstChild()
set the current node to the first child ...
#define str(s)