LDMX Software
dqm::VisiblesFeatureProducer Class Reference

Public Member Functions

 VisiblesFeatureProducer (const std::string &name, framework::Process &process)
 
void configure (framework::config::Parameters &parameters) override
 Callback for the EventProcessor to configure itself from the given set of parameters.
 
void analyze (const framework::Event &event) override
 Process the event and make histograms or summaries.
 
- Public Member Functions inherited from framework::Analyzer
 Analyzer (const std::string &name, Process &process)
 Class constructor.
 
virtual void process (Event &event) final
 Processing an event for an Analyzer is calling analyze.
 
virtual void beforeNewRun (ldmx::RunHeader &run_header) final
 Don't allow Analyzers to add parameters to the run header.
 
- Public Member Functions inherited from framework::EventProcessor
 DECLARE_FACTORY (EventProcessor, EventProcessor *, const std::string &, Process &)
 declare that we have a factory for this class
 
 EventProcessor (const std::string &name, Process &process)
 Class constructor.
 
virtual ~EventProcessor ()=default
 Class destructor.
 
virtual void onNewRun (const ldmx::RunHeader &run_header)
 Callback for the EventProcessor to take any necessary action when the run being processed changes.
 
virtual void onFileOpen (EventFile &event_file)
 Callback for the EventProcessor to take any necessary action when a new event input ROOT file is opened.
 
virtual void onFileClose (EventFile &event_file)
 Callback for the EventProcessor to take any necessary action when a event input ROOT file is closed.
 
virtual void onProcessStart ()
 Callback for the EventProcessor to take any necessary action when the processing of events starts, such as creating histograms.
 
virtual void onProcessEnd ()
 Callback for the EventProcessor to take any necessary action when the processing of events finishes, such as calculating job-summary quantities.
 
template<class T >
const T & getCondition (const std::string &condition_name)
 Access a conditions object for the current event.
 
TDirectory * getHistoDirectory ()
 Access/create a directory in the histogram file for this event processor to create histograms and analysis tuples.
 
void setStorageHint (framework::StorageControl::Hint hint)
 Mark the current event as having the given storage control hint from this module_.
 
void setStorageHint (framework::StorageControl::Hint hint, const std::string &purposeString)
 Mark the current event as having the given storage control hint from this module and the given purpose string.
 
int getLogFrequency () const
 Get the current logging frequency from the process.
 
int getRunNumber () const
 Get the run number from the process.
 
std::string getName () const
 Get the processor name.
 
void createHistograms (const std::vector< framework::config::Parameters > &histos)
 Internal function which is used to create histograms passed from the python configuration @parma histos vector of Parameters that configure histograms to create.
 

Private Member Functions

bool inList (std::vector< int > parents, int track_id)
 

Private Attributes

bool training_ {false}
 
std::string training_file_
 
double beam_energy_mev_ {0.}
 
std::string hcal_rec_collection_
 
std::string hcal_rec_pass_name_
 
std::string ecal_rec_collection_
 
std::string ecal_rec_pass_name_
 
bool recoil_from_tracking_ {false}
 
std::string track_collection_
 
std::string track_pass_name_
 
std::string sp_collection_
 
std::string sp_pass_name_
 
std::string sim_particles_pass_name_
 

Additional Inherited Members

- Protected Member Functions inherited from framework::EventProcessor
void abortEvent ()
 Abort the event immediately.
 
- Protected Attributes inherited from framework::EventProcessor
HistogramPool histograms_
 helper object for making and filling histograms
 
NtupleManagerntuple_ {NtupleManager::getInstance()}
 Manager for any ntuples.
 
logging::logger the_log_
 The logger for this EventProcessor.
 

Detailed Description

Definition at line 17 of file VisiblesFeatureProducer.h.

Constructor & Destructor Documentation

◆ VisiblesFeatureProducer()

dqm::VisiblesFeatureProducer::VisiblesFeatureProducer ( const std::string & name,
framework::Process & process )
inline

Definition at line 19 of file VisiblesFeatureProducer.h.

20 : Analyzer(name, process) {}
virtual void process(Event &event) final
Processing an event for an Analyzer is calling analyze.
Analyzer(const std::string &name, Process &process)
Class constructor.

Member Function Documentation

◆ analyze()

void dqm::VisiblesFeatureProducer::analyze ( const framework::Event & event)
overridevirtual

Process the event and make histograms or summaries.

Parameters
eventThe Event to analyze

Implements framework::Analyzer.

Definition at line 50 of file VisiblesFeatureProducer.cxx.

50 {
51 std::vector<double> bdt_features;
52
53 const auto &particle_map{event.getMap<int, ldmx::SimParticle>(
54 "SimParticles", sim_particles_pass_name_)};
55
56 // Get target scoring plane hits for recoil electron
57 // Use this to calculate the projected photon line vector
58 // This currently uses truth-level information, but it should be replaced
59 // by reconstructed tracker information, when available
60 std::vector<double> gamma_p(3);
61 std::vector<double> gamma_x0(3);
62
63 std::vector<double> recoil_p(3);
64 bool found_recoil_e = false;
65
66 if (recoil_from_tracking_) {
67 const auto &recoil_tracks{
68 event.getCollection<ldmx::Track>(track_collection_, track_pass_name_)};
69 // Fill this in later when you know how to use it
70 for (auto &track : recoil_tracks) {
71 // need to figure out how to best isolate candidate electron track
72 if (track.q() == 1 && track.getNhits() == 5) {
73 gamma_x0 = track.getPosition();
74 gamma_p[0] = -1. * track.getMomentum()[0];
75 gamma_p[1] = -1. * track.getMomentum()[1];
76 gamma_p[2] = beam_energy_mev_ - track.getMomentum()[2];
77 }
78 }
79 } else {
80 if (event.exists(sp_collection_, sp_pass_name_)) {
81 const auto &target_sp_hits = event.getCollection<ldmx::SimTrackerHit>(
82 sp_collection_, sp_pass_name_);
83 bool found_rec = false;
84 for (auto const &it : particle_map) {
85 for (auto const &sphit : target_sp_hits) {
86 if (sphit.getPosition()[2] > 0) {
87 if (it.first == sphit.getTrackID()) {
88 if (it.second.getPdgID() == 622) {
89 std::vector<float> x0_float = sphit.getPosition();
90 std::vector<double> x0_double(x0_float.begin(), x0_float.end());
91 gamma_x0 = x0_double;
92 gamma_p = sphit.getMomentum();
93 found_rec = true;
94 }
95 if (it.second.getPdgID() == 11 &&
96 inList(it.second.getParents(), 0)) {
97 if (!found_rec) {
98 std::vector<float> x0_float = sphit.getPosition();
99 std::vector<double> x0_double(x0_float.begin(),
100 x0_float.end());
101 gamma_x0 = x0_double;
102 gamma_p[0] = -1. * sphit.getMomentum()[0];
103 gamma_p[1] = -1. * sphit.getMomentum()[1];
104 gamma_p[2] = beam_energy_mev_ - sphit.getMomentum()[2];
105 found_rec = true;
106 }
107 recoil_p = sphit.getMomentum();
108 found_recoil_e = true;
109 }
110 }
111 }
112 }
113 }
114 }
115 }
116
117 double p_mag = 0.;
118 if (found_recoil_e) {
119 p_mag = std::sqrt(recoil_p[0] * recoil_p[0] + recoil_p[1] * recoil_p[1] +
120 recoil_p[2] * recoil_p[2]);
121 }
122
123 const auto &ecal_rec_hits = event.getCollection<ldmx::EcalHit>(
124 ecal_rec_collection_, ecal_rec_pass_name_);
125 const auto &hcal_rec_hits = event.getCollection<ldmx::HcalHit>(
126 hcal_rec_collection_, hcal_rec_pass_name_);
127
128 double ecal_energy = 0.;
129 double hcal_energy = 0.;
130 bool hcal_containment = true;
131
132 for (const ldmx::EcalHit &hit : ecal_rec_hits) {
133 if (hit.getEnergy() > 0.) {
134 ecal_energy += hit.getEnergy();
135 }
136 }
137 for (const ldmx::HcalHit &hit : hcal_rec_hits) {
138 if (hit.getEnergy() > 0.) {
139 ldmx::HcalID det_id(hit.getID());
140 if (det_id.getSection() != 0) {
141 continue;
142 }
143 if (det_id.getLayerID() == 1 && hit.getPE() > 5) {
144 hcal_containment = false;
145 }
146 hcal_energy += 12. * hit.getEnergy();
147 }
148 }
149
150 if (ecal_energy < 3160 && hcal_energy > 4840 && hcal_containment &&
151 p_mag < 2400) {
152 // initialize all of the features
153 int n_layers_hit = 0;
154 double x_std = 0.;
155 double y_std = 0.;
156 double z_std = 0.;
157 double x_mean = 0.;
158 double y_mean = 0.;
159 double r_mean = 0.;
160 int iso_hits = 0;
161 double iso_energy = 0.;
162 int n_readout_hits = 0;
163 double summed_det = 0.;
164 double r_mean_from_photon_track = 0.;
165
166 double z_mean = 0.; // need this when calculating z_std
167 std::vector<int> layers_hit;
168
169 for (const ldmx::HcalHit &hit : hcal_rec_hits) {
170 if (hit.getEnergy() > 0.) {
171 ldmx::HcalID det_id(hit.getID());
172 if (det_id.getSection() != 0) { // skip hits that aren't in main Hcal
173 continue;
174 }
175 if (fabs(hit.getXPos()) > 1000 || fabs(hit.getYPos()) > 1000) {
176 continue;
177 }
178
179 n_readout_hits += 1;
180 double hit_x = hit.getXPos();
181 double hit_y = hit.getYPos();
182 double hit_z = hit.getZPos();
183 double hit_r = sqrt(hit_x * hit_x + hit_y * hit_y);
184
185 summed_det += hit.getEnergy();
186
187 x_mean += hit_x * hit.getEnergy();
188 y_mean += hit_y * hit.getEnergy();
189 z_mean += hit_z * hit.getEnergy();
190 r_mean += hit_r * hit.getEnergy();
191
192 // check if this is a new layer in the collection
193 if (!(std::find(layers_hit.begin(), layers_hit.end(),
194 det_id.getLayerID()) != layers_hit.end())) {
195 layers_hit.push_back(det_id.getLayerID());
196 }
197
198 double proj_x =
199 gamma_x0[0] + (hit_z - gamma_x0[2]) * gamma_p[0] / gamma_p[2];
200 double proj_y =
201 gamma_x0[1] + (hit_z - gamma_x0[2]) * gamma_p[1] / gamma_p[2];
202
203 r_mean_from_photon_track +=
204 hit.getEnergy() * sqrt((hit_x - proj_x) * (hit_x - proj_x) +
205 (hit_y - proj_y) * (hit_y - proj_y));
206
207 // Calculate isolated hits
208 double closest_point = 9999.;
209 for (const ldmx::HcalHit &hit2 : hcal_rec_hits) {
210 if (hit2.getEnergy() > 0.) {
211 ldmx::HcalID det_i_d2(hit2.getID());
212 if (fabs(hit2.getXPos()) > 1000 || fabs(hit2.getYPos()) > 1000) {
213 continue;
214 }
215 if (det_i_d2.getLayerID() == det_id.getLayerID()) {
216 // Determine if a bar is vertical (along y-axis) or horizontal
217 // (along x-axis) Odd layers have horizontal strips Even layers
218 // have vertical strips
219 if (hit.isOrientationX()) {
220 if (fabs(hit2.getYPos() - hit_y) > 0) {
221 if (fabs(hit2.getYPos() - hit_y) < closest_point) {
222 closest_point = fabs(hit2.getYPos() - hit_y);
223 }
224 }
225 }
226 if (hit.isOrientationY()) {
227 if (fabs(hit2.getXPos() - hit_x) > 0) {
228 if (fabs(hit2.getXPos() - hit_x) < closest_point) {
229 closest_point = fabs(hit2.getXPos() - hit_x);
230 }
231 }
232 }
233 }
234 }
235 }
236 if (closest_point > 50.) {
237 iso_hits += 1;
238 iso_energy += hit.getEnergy();
239 }
240 }
241 }
242
243 n_layers_hit = layers_hit.size();
244
245 if (summed_det > 0.) {
246 x_mean /= summed_det;
247 y_mean /= summed_det;
248 z_mean /= summed_det;
249 r_mean /= summed_det;
250
251 r_mean_from_photon_track /= summed_det;
252 }
253
254 for (const ldmx::HcalHit &hit : hcal_rec_hits) {
255 if (hit.getEnergy() > 0.) {
256 if (fabs(hit.getXPos()) > 1000 || fabs(hit.getYPos()) > 1000) {
257 continue;
258 }
259 ldmx::HcalID det_id(hit.getID());
260 if (det_id.getSection() == 0) {
261 x_std += hit.getEnergy() * (hit.getXPos() - x_mean) *
262 (hit.getXPos() - x_mean);
263 y_std += hit.getEnergy() * (hit.getYPos() - y_mean) *
264 (hit.getYPos() - y_mean);
265 z_std += hit.getEnergy() * (hit.getZPos() - z_mean) *
266 (hit.getZPos() - z_mean);
267 }
268 }
269 }
270
271 if (summed_det > 0.) {
272 x_std = sqrt(x_std / summed_det);
273 y_std = sqrt(y_std / summed_det);
274 z_std = sqrt(z_std / summed_det);
275 }
276
277 // Fill histograms
278 histograms_.fill("layers_hit", n_layers_hit);
279 histograms_.fill("x_std", x_std);
280 histograms_.fill("y_std", y_std);
281 histograms_.fill("z_std", z_std);
282 histograms_.fill("x_mean", x_mean);
283 histograms_.fill("y_mean", y_mean);
284 histograms_.fill("r_mean", r_mean);
285 histograms_.fill("iso_hits", iso_hits);
286 histograms_.fill("iso_energy", iso_energy);
287 histograms_.fill("n_hits", n_readout_hits);
288 histograms_.fill("total_energy", hcal_energy);
289 histograms_.fill("photon_track", r_mean_from_photon_track);
290
291 bdt_features.push_back(n_layers_hit);
292 bdt_features.push_back(x_std);
293 bdt_features.push_back(y_std);
294 bdt_features.push_back(z_std);
295 bdt_features.push_back(x_mean);
296 bdt_features.push_back(y_mean);
297 bdt_features.push_back(r_mean);
298 bdt_features.push_back(iso_hits);
299 bdt_features.push_back(iso_energy);
300 bdt_features.push_back(n_readout_hits);
301 bdt_features.push_back(hcal_energy);
302 bdt_features.push_back(r_mean_from_photon_track);
303
304 if (training_) {
305 std::ofstream file(training_file_, std::ios::app);
306 if (!file.is_open()) {
307 ldmx_log(fatal) << "Error: Could not open file " << training_file_;
308 return;
309 }
310 for (int i = 0; i < bdt_features.size(); ++i) {
311 file << bdt_features[i] << (i + 1 == bdt_features.size() ? "\n" : ", ");
312 }
313 }
314 }
315
316 return;
317}
HistogramPool histograms_
helper object for making and filling histograms
bool exists(const std::string &name, const std::string &passName, bool unique=true) const
Check for the existence of an object or collection with the given name and pass name in the event.
Definition Event.cxx:92
void fill(const std::string &name, const T &val)
Fill a 1D histogram.
Stores reconstructed hit information from the ECAL.
Definition EcalHit.h:19
Stores reconstructed hit information from the HCAL.
Definition HcalHit.h:24
Implements detector ids for HCal subdetector.
Definition HcalID.h:19
Class representing a simulated particle.
Definition SimParticle.h:23
Represents a simulated tracker hit in the simulation.
Implementation of a track object.
Definition Track.h:53

References framework::Event::exists(), framework::HistogramPool::fill(), ldmx::HcalID::getLayerID(), and framework::EventProcessor::histograms_.

◆ configure()

void dqm::VisiblesFeatureProducer::configure ( framework::config::Parameters & parameters)
overridevirtual

Callback for the EventProcessor to configure itself from the given set of parameters.

The parameters a processor has access to are the member variables of the python class in the sequence that has className equal to the EventProcessor class name.

For an example, look at MyProcessor.

Parameters
parametersParameters for configuration.

Reimplemented from framework::EventProcessor.

Definition at line 22 of file VisiblesFeatureProducer.cxx.

23 {
24 training_ = parameters.get<bool>("training");
25 training_file_ = parameters.get<std::string>("training_file");
26
27 beam_energy_mev_ = parameters.get<double>("beam_energy");
28
29 // collection names
30 hcal_rec_collection_ = parameters.get<std::string>("hcal_rec_coll_name");
31 hcal_rec_pass_name_ = parameters.get<std::string>("hcal_rec_pass_name");
32
33 ecal_rec_collection_ = parameters.get<std::string>("ecal_rec_coll_name");
34 ecal_rec_pass_name_ = parameters.get<std::string>("ecal_rec_pass_name");
35
36 recoil_from_tracking_ = parameters.get<bool>("recoil_from_tracking");
37 track_collection_ = parameters.get<std::string>("track_collection");
38 track_pass_name_ = parameters.get<std::string>("track_pass_name");
39
40 sp_collection_ = parameters.get<std::string>("sp_coll_name");
41 sp_pass_name_ = parameters.get<std::string>("sp_pass_name");
42 sim_particles_pass_name_ =
43 parameters.get<std::string>("sim_particles_pass_name");
44}
const T & get(const std::string &name) const
Retrieve the parameter of the given name.
Definition Parameters.h:78

References framework::config::Parameters::get().

◆ inList()

bool dqm::VisiblesFeatureProducer::inList ( std::vector< int > parents,
int track_id )
private

Definition at line 46 of file VisiblesFeatureProducer.cxx.

46 {
47 return std::find(parents.begin(), parents.end(), track_id) != parents.end();
48}

Member Data Documentation

◆ beam_energy_mev_

double dqm::VisiblesFeatureProducer::beam_energy_mev_ {0.}
private

Definition at line 32 of file VisiblesFeatureProducer.h.

32{0.};

◆ ecal_rec_collection_

std::string dqm::VisiblesFeatureProducer::ecal_rec_collection_
private

Definition at line 36 of file VisiblesFeatureProducer.h.

◆ ecal_rec_pass_name_

std::string dqm::VisiblesFeatureProducer::ecal_rec_pass_name_
private

Definition at line 37 of file VisiblesFeatureProducer.h.

◆ hcal_rec_collection_

std::string dqm::VisiblesFeatureProducer::hcal_rec_collection_
private

Definition at line 34 of file VisiblesFeatureProducer.h.

◆ hcal_rec_pass_name_

std::string dqm::VisiblesFeatureProducer::hcal_rec_pass_name_
private

Definition at line 35 of file VisiblesFeatureProducer.h.

◆ recoil_from_tracking_

bool dqm::VisiblesFeatureProducer::recoil_from_tracking_ {false}
private

Definition at line 38 of file VisiblesFeatureProducer.h.

38{false};

◆ sim_particles_pass_name_

std::string dqm::VisiblesFeatureProducer::sim_particles_pass_name_
private

Definition at line 43 of file VisiblesFeatureProducer.h.

◆ sp_collection_

std::string dqm::VisiblesFeatureProducer::sp_collection_
private

Definition at line 41 of file VisiblesFeatureProducer.h.

◆ sp_pass_name_

std::string dqm::VisiblesFeatureProducer::sp_pass_name_
private

Definition at line 42 of file VisiblesFeatureProducer.h.

◆ track_collection_

std::string dqm::VisiblesFeatureProducer::track_collection_
private

Definition at line 39 of file VisiblesFeatureProducer.h.

◆ track_pass_name_

std::string dqm::VisiblesFeatureProducer::track_pass_name_
private

Definition at line 40 of file VisiblesFeatureProducer.h.

◆ training_

bool dqm::VisiblesFeatureProducer::training_ {false}
private

Definition at line 29 of file VisiblesFeatureProducer.h.

29{false};

◆ training_file_

std::string dqm::VisiblesFeatureProducer::training_file_
private

Definition at line 30 of file VisiblesFeatureProducer.h.


The documentation for this class was generated from the following files: