Process the event and put new data products into it.
38 {
41
42 std::vector<ldmx::HcalHit> hcal_rec_hits;
43
45 input_pass_name_)};
46 std::unordered_map<unsigned int, std::vector<const ldmx::SimCalorimeterHit*>>
47 hits_by_id{};
48
49
50 for (const auto& hit : sim_hits) {
51 auto id{hit.getID()};
52 auto found{hits_by_id.find(id)};
53 if (found == hits_by_id.end()) {
54 hits_by_id[id] = std::vector<const ldmx::SimCalorimeterHit*>{&hit};
55 } else {
56 hits_by_id[id].push_back(&hit);
57 }
58 }
59 for (const auto& [barID, simhits_in_bar] : hits_by_id) {
61 double edep{};
62 double time{};
63 std::vector<double> pos{0, 0, 0};
64 for (auto hit : simhits_in_bar) {
65 edep += hit->getEdep();
66 double edep_hit = hit->getEdep();
67 time += hit->getTime() * edep_hit;
68 auto hit_pos{hit->getPosition()};
69 pos[0] += hit_pos[0] * edep_hit;
70 pos[1] += hit_pos[1] * edep_hit;
71 pos[2] += hit_pos[2] * edep_hit;
72 }
74
75
76 double mean_pe{(edep / mev_per_mip_) * pe_per_mip_};
77 double xpos{pos[0] / edep};
78 double ypos{pos[1] / edep};
79 double zpos{pos[2] / edep};
80 time /= edep;
81
82 auto orientation{hcal_geometry.getScintillatorOrientation(barID)};
83 double half_total_width{
84 hcal_geometry.getHalfTotalWidth(hit_id.section(), hit_id.layer())};
85 double scint_bar_length{hcal_geometry.getScintillatorLength(hit_id)};
86
87 auto strip_center{hcal_geometry.getStripCenterPosition(hit_id)};
88 if (hit_id.section() == ldmx::HcalID::HcalSection::BACK) {
89 double distance_along_bar =
90 (orientation ==
91 ldmx::HcalGeometry::ScintillatorOrientation::horizontal)
92 ? xpos
93 : ypos;
94 if (orientation ==
95 ldmx::HcalGeometry::ScintillatorOrientation::horizontal) {
96 ypos = strip_center.y();
97 xpos += (*position_resolution_smear_)(rng_);
98 } else {
99 xpos = strip_center.x();
100 ypos += (*position_resolution_smear_)(rng_);
101 }
102 zpos = strip_center.z();
103
104 mean_pe *= exp(1. / attenuation_length_);
105 double mean_pe_close =
106 mean_pe * exp(-1. *
107 ((half_total_width - distance_along_bar) /
108 (scint_bar_length * 0.5)) /
109 attenuation_length_);
110 double mean_pe_far =
111 mean_pe * exp(-1. *
112 ((half_total_width + distance_along_bar) /
113 (scint_bar_length * 0.5)) /
114 attenuation_length_);
115 int pe_close{
116 std::poisson_distribution<int>(mean_pe_close + mean_noise_)(rng_)};
117 int pe_far{
118 std::poisson_distribution<int>(mean_pe_far + mean_noise_)(rng_)};
119 rec_hit.
setPE(pe_close + pe_far);
120 rec_hit.
setMinPE(std::min(pe_close, pe_far));
121 } else {
122
123 int pe{std::poisson_distribution<int>(mean_pe + mean_noise_)(rng_)};
126
127
128
129 if (orientation ==
130 ldmx::HcalGeometry::ScintillatorOrientation::horizontal) {
131 xpos += (*position_resolution_smear_)(rng_);
132 ypos = strip_center.y();
133 zpos = strip_center.z();
134 } else if (orientation ==
135 ldmx::HcalGeometry::ScintillatorOrientation::vertical) {
136 xpos = strip_center.x();
137 ypos += (*position_resolution_smear_)(rng_);
138 zpos = strip_center.z();
139 } else if (orientation ==
140 ldmx::HcalGeometry::ScintillatorOrientation::depth) {
141 xpos = strip_center.x();
142 ypos = strip_center.y();
143 zpos += (*position_resolution_smear_)(rng_);
144 } else {
145 xpos = strip_center.x();
146 ypos = strip_center.y();
147 zpos = strip_center.z();
148 ldmx_log(warn) << "Bar orientation not found. Hit" << hit_id.raw()
149 << "positioned at bar center.";
150 }
151 }
152
153 rec_hit.
setID(hit_id.raw());
164 }
165 event.add(output_coll_name_, hcal_rec_hits);
166}
void setYPos(float ypos)
Set the Y position of the hit [mm].
void setID(int id)
Set the detector ID.
void setZPos(float zpos)
Set the Z position of the hit [mm].
void setXPos(float xpos)
Set the X position of the hit [mm].
void setTime(float time)
Set the time of the hit [ns].
void setEnergy(float energy)
Set the calorimetric energy of the hit, corrected for sampling factors [MeV].
void setNoise(bool yes)
Set if this hit is a noise hit.
static constexpr const char * CONDITIONS_OBJECT_NAME
Conditions object: The name of the python configuration calling this class (Hcal/python/HcalGeometry....
Stores reconstructed hit information from the HCAL.
void setSection(int section)
Set the section for this hit.
void setMinPE(float minpe)
Set the minimum number of photoelectrons estimated for this hit.
void setOrientation(int orientation)
Set if the bar is orientied in X / Y / Z meanig 0 / 1 / 2, respectively.
void setStrip(int strip)
Set the strip for this hit.
void setLayer(int layer)
Set the layer for this hit.
void setPE(float pe)
Set the number of photoelectrons estimated for this hit.
Implements detector ids for HCal subdetector.
Stores simulated calorimeter hit information.