LDMX Software
trigscint::TrigScintTrackProducer Class Reference

making tracks from trigger scintillator clusters More...

#include <TrigScintTrackProducer.h>

Public Member Functions

 TrigScintTrackProducer (const std::string &name, framework::Process &process)
 
void configure (framework::config::Parameters &ps) override
 Callback for the EventProcessor to configure itself from the given set of parameters.
 
void produce (framework::Event &event) override
 Process the event and put new data products into it.
 
void onProcessStart () override
 Callback for the EventProcessor to take any necessary action when the processing of events starts, such as creating histograms.
 
void onProcessEnd () override
 Callback for the EventProcessor to take any necessary action when the processing of events finishes, such as calculating job-summary quantities.
 
- Public Member Functions inherited from framework::Producer
 Producer (const std::string &name, Process &process)
 Class constructor.
 
virtual void process (Event &event) final
 Processing an event for a Producer is calling produce.
 
- 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 beforeNewRun (ldmx::RunHeader &run_header)
 Callback for Producers to add parameters to the run header before conditions are initialized.
 
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.
 
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

ldmx::TrigScintTrack makeTrack (std::vector< ldmx::TrigScintCluster > clusters)
 
void matchXYTracks (std::vector< ldmx::TrigScintTrack > &tracks)
 

Private Attributes

std::vector< ldmx::TrigScintTracktracks_
 
double max_delta_ {0.}
 
int verbose_ {0}
 
std::string seeding_collection_
 
std::vector< std::string > input_collections_
 
std::string output_collection_
 
std::string pass_name_ {""}
 
bool skip_last_ {false}
 
int vert_bar_start_idx_ {52}
 
int n_bars_y_ {16}
 
int n_bars_x_ {8}
 
float bar_width_y_ {3.}
 
float bar_gap_y_ {2.1}
 
float bar_width_x_ {3.}
 
float bar_gap_x_ {0.1}
 
float x_conv_factor_
 
float x_start_
 
float y_conv_factor_
 
float y_start_
 

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

making tracks from trigger scintillator clusters

Definition at line 19 of file TrigScintTrackProducer.h.

Constructor & Destructor Documentation

◆ TrigScintTrackProducer()

trigscint::TrigScintTrackProducer::TrigScintTrackProducer ( const std::string & name,
framework::Process & process )
inline

Definition at line 21 of file TrigScintTrackProducer.h.

22 : Producer(name, process) {}
Producer(const std::string &name, Process &process)
Class constructor.
virtual void process(Event &event) final
Processing an event for a Producer is calling produce.

Member Function Documentation

◆ configure()

void trigscint::TrigScintTrackProducer::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 8 of file TrigScintTrackProducer.cxx.

8 {
9 max_delta_ = ps.get<double>(
10 "delta_max"); // max distance to consider adding in a cluster to track
11 seeding_collection_ = ps.get<std::string>(
12 "seeding_collection"); // probably tagger pad, "TriggerPadTagClusters"
13 input_collections_ = ps.get<std::vector<std::string>>(
14 "further_input_collections"); // {"TriggerPadUpClusters" ,
15 // "TriggerPadDnClusters" }
16 output_collection_ = ps.get<std::string>("output_collection");
17 pass_name_ = ps.get<std::string>("input_pass_name");
18 verbose_ = ps.get<int>("verbosity");
19 vert_bar_start_idx_ = ps.get<int>("vertical_bar_start_index");
20 n_bars_y_ = ps.get<int>("number_horizontal_bars");
21 bar_width_y_ = ps.get<double>("horizontal_bar_width");
22 bar_gap_y_ = ps.get<double>("horizontal_bar_gap");
23 n_bars_x_ = ps.get<int>("number_vertical_bars");
24 bar_width_x_ = ps.get<double>("vertical_bar_width");
25 bar_gap_x_ = ps.get<double>("vertical_bar_gap");
26 skip_last_ = ps.get<bool>("allow_skip_last_collection");
27
28 // TO DO: allow any number of input collections
29
30 if (verbose_) {
31 ldmx_log(info) << "In TrigScintTrackProducer: configure done!" << std::endl;
32 ldmx_log(info) << "Got parameters: \nSeeding: " << seeding_collection_
33 << "\nTolerance: " << max_delta_
34 << "\nInput: " << input_collections_.at(0) << " and "
35 << input_collections_.at(1)
36 << "\nInput pass name: " << pass_name_
37 << "\nAllow tracks with no hit in last collection: "
38 << skip_last_
39 << "\nVertical bar start index: " << vert_bar_start_idx_
40 << "\nNumber of horizontal bars: " << n_bars_y_
41 << "\nHorizontal bar width: " << bar_width_y_
42 << "\nHorizontal bar gap: " << bar_gap_y_
43 << "\nNumber of vertical bars: " << n_bars_x_
44 << "\nVertical bar width: " << bar_width_x_
45 << "\nVertical bar gap: " << bar_gap_x_
46 << "\nOutput: " << output_collection_
47 << "\nVerbosity: " << verbose_;
48 }
49 // each bar only goes half this distance up (overlap/zig-zag)
50 y_conv_factor_ = (bar_width_y_ + bar_gap_y_) / 2.;
51 // half height of pad
52 y_start_ = -(n_bars_y_ * (bar_width_y_ + bar_gap_y_) - bar_gap_y_) / 2.;
53 // each bar goes entire distance sideways (no overlap)
54 x_conv_factor_ = bar_width_x_ + bar_gap_x_;
55 // half width of pad
56 x_start_ = -(n_bars_x_ * (bar_width_x_ + bar_gap_x_) - bar_gap_x_) / 2.;
57
58 return;
59}

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

◆ makeTrack()

ldmx::TrigScintTrack trigscint::TrigScintTrackProducer::makeTrack ( std::vector< ldmx::TrigScintCluster > clusters)
private

Definition at line 474 of file TrigScintTrackProducer.cxx.

475 {
476 // for now let's keep a straight, unweighted centroid
477 // consider the possibility that at least one cluster has a centroid
478 // identically == 0. then we need to shift them by 1 if we want to do energy
479 // weighted track centroid later. but no need now
481 float centroid = 0;
482 float centroid_x = 0;
483 float centroid_y = 0;
484 float beam_efrac = 0;
485 float pe = 0;
486 for (uint i = 0; i < clusters.size(); i++) {
487 centroid += (clusters.at(i)).getCentroid();
488 centroid_x += (clusters.at(i)).getCentroidX();
489 centroid_y += (clusters.at(i)).getCentroidY();
490 tr.addConstituent(clusters.at(i));
491 beam_efrac += (clusters.at(i)).getBeamEfrac();
492 pe += (clusters.at(i)).getPE();
493 }
494 centroid /= clusters.size();
495 centroid_x /= clusters.size();
496 if (centroid >= vert_bar_start_idx_) {
497 if (verbose_) {
498 ldmx_log(debug)
499 << " -- In makeTrack made vertical bar track with centroid "
500 << centroid << " and y flag sum " << centroid_y;
501 // try commenting this to check if that helps with an out-of-bounds
502 // problem
503 // << " from clusters with y centroids";
504 // for (uint i = 0; i < clusters.size(); i++)
505 // ldmx_log(debug) << "\tpad " << i << ": centroidY "
506 // << (clusters.at(i)).getCentroidY();
507 }
508 // then the sum of centroid y is 0, 2, 4 or 6
509 // we have 4 divisions, so, the center of it should be divNb/8
510 // (or rather, that's where channel nBars/8 begins)
511 // and then a factor 2 for the zig-zag pattern
512 centroid_y = (centroid_y + 1) * 2 * n_bars_y_ / 8.;
513 // TODO: here we could instead just use quadrant indices 0-3 by dividing by
514 // 2 but that would mean that in the raw, x and y track centroidY would mean
515 // different things
516 if (verbose_) ldmx_log(debug) << " -- new centroidY = " << centroid_y;
517 } else
518 centroid_y /= clusters.size();
519
520 beam_efrac /= clusters.size();
521 pe /= clusters.size();
522
523 float residual = 0;
524 for (uint i = 0; i < clusters.size(); i++)
525 residual += ((clusters.at(i)).getCentroid() - centroid) *
526 ((clusters.at(i)).getCentroid() - centroid);
527 residual = sqrt(residual / clusters.size());
528
529 tr.setCentroid(centroid);
530 tr.setCentroidX(centroid_x);
531 tr.setCentroidY(centroid_y);
532 tr.setResidual(residual);
533 tr.setBeamEfrac(beam_efrac);
534 tr.setPE(pe);
535
536 if (verbose_) {
537 ldmx_log(debug) << " -- In makeTrack made track with centroid "
538 << centroid << " and residual " << residual << " and pe "
539 << pe << " from clusters with centroids";
540 for (uint i = 0; i < clusters.size(); i++)
541 ldmx_log(debug) << "\tpad " << i << ": centroid "
542 << (clusters.at(i)).getCentroid();
543 }
544
545 return tr;
546}
Represents a track of trigger scintillator clusters.
void setCentroidX(float centroid)
Set the x centroid of the track.
void setResidual(float resid)
Set the detector ID residual of the track.
void setCentroidY(float centroid)
Set the y centroid of the track.
void setPE(float pe)
Set the average cluster pe of the track.
void addConstituent(TrigScintCluster cl)
Add a cluster to the list of track constituents.
void setCentroid(float centroid)
Set the detector ID centroid of the track.
void setBeamEfrac(float e)
Set beam energy fraction of hit.

◆ matchXYTracks()

void trigscint::TrigScintTrackProducer::matchXYTracks ( std::vector< ldmx::TrigScintTrack > & tracks)
private

Definition at line 549 of file TrigScintTrackProducer.cxx.

550 {
551 // map quadrant nb to track (can be multiple per quadrant)
552 std::multimap<int, int>
553 y_idx_quad_map; // key = quad, val = track index in collection
554 std::multimap<int, int> x_idx_quad_map;
555
556 std::multimap<int, ldmx::TrigScintTrack> y_quad_map;
557 std::multimap<int, ldmx::TrigScintTrack> x_quad_map;
558 // map track in quadrant back to index in entire track collection
559 // used for updating collection track variables
560 std::map<ldmx::TrigScintTrack, int> y_track_map;
561 std::map<ldmx::TrigScintTrack, int> x_track_map;
562
563 uint trk_idx = -1;
564 for (auto trk : tracks) {
565 trk_idx++;
566 // 1. get the y bar tracks with centroidX = -1
567 if (trk.getCentroidX() == -1) {
568 if (verbose_)
569 ldmx_log(debug) << " -- In matchXYTracks found y track at "
570 << trk.getCentroidY() << "; mapping to quad "
571 << (int)trk.getCentroidY() / 8 << " with trk index "
572 << trk_idx;
573 // 2. order them... or map them to quadrants. note that there are 2 layers
574 // so 2*n_bars_y_/4 channels per quadrant
575 y_quad_map.insert(std::make_pair((int)(trk.getCentroidY() / 8), trk));
576 y_track_map[trk] = trk_idx;
577 y_idx_quad_map.insert(
578 std::make_pair((int)(trk.getCentroidY() / 8), trk_idx));
579
580 } else { // 3. get the remaining tracks (from vertical bars) and map them
581 // (back) to (middle of) quadrants
582 x_quad_map.insert(std::make_pair((int)(trk.getCentroidY() / 8), trk));
583 x_track_map[trk] = trk_idx;
584 x_idx_quad_map.insert(
585 std::make_pair((int)(trk.getCentroidY() / 8), trk_idx));
586 if (verbose_)
587 ldmx_log(debug) << " -- In matchXYTracks found x track at (x,y) = ("
588 << trk.getCentroidX() << ", " << trk.getCentroidY()
589 << "); mapping to quad " << (int)trk.getCentroidY() / 8
590 << " with trk index " << trk_idx;
591 }
592 }
593
594 // 4a
595 //
596 // 1) here use the geometry? if we can assume perfect alignment we can take
597 // width and nBars and take nBars/2 as origin
598 // --- now do the matching ---
599
600 // if there is no useful matching to be done: these are the pad width wide
601 // numbers
602 float x0 = 0;
603 // this should be half the pad... could also set
604 // it to full beam spot width
605 float sx0 = fabs(x_start_);
606
607 // y_start_ is half the pad, so this should be half a quadrant
608 float sy0 = fabs(y_start_) / 4.;
609
610 // assume at least one y track. will have to figure out if there is ever a
611 // reason to use an isolated x track in its place.
612 for (auto yitr = y_quad_map.begin(); yitr != y_quad_map.end(); ++yitr) {
613 int n_yin_quad = y_quad_map.count((*yitr).first);
614 int n_xin_quad = x_quad_map.count((*yitr).first);
615 float y{-9999.}, sy{-9999.}, x{-9999.}, x1{-9999.}, x2{-9999.}, sx1{-9999.},
616 sx2{-9999.}, y1{-9999.}, y2{-9999.}, sy1{-9999.}, sy2{-9999.};
617 // quad midpoint:
618 float y0 = (*yitr).first * 8 + sy0;
619 float sx = 1. / 2 *
620 x_conv_factor_; // rely on x precision being one single bar
621 // width; always used unless x is undeterminable
622
623 // check all x first
624 // do the easiest first:
625 if (n_xin_quad == 0) { // then there's no hope of setting a better x here
626 // just use the beam spot width... and center of pad
627 x = x0;
628 sx = sx0;
629 if (verbose_)
630 ldmx_log(debug) << "\t\t\t no x info in quad " << (*yitr).first
631 << "; will set x to middle of pad, pad half-width as "
632 "precision: set (x, sx)=("
633 << x << ", " << sx << ")";
634 } // 0 x tracks in quadrant
635 else if (n_xin_quad ==
636 1) { // slightly harder: 1 x track -- might be easy if
637 // it's just one y track; if several, need to
638 // think about overlaps. but in overlap case, just
639 // revert to setting x0 and sx0, when we know
640 auto xitr = x_quad_map.find((*yitr).first);
641 x = ((*xitr).second).getCentroidX() * x_conv_factor_ + x_start_;
642
643 if (verbose_)
644 ldmx_log(debug) << "\t\t\t 1 x in quad " << (*yitr).first
645 << ", getting (x, sx)=(" << x << ", " << sx << ")";
646 } // 1 x track in quadrant
647 else if (n_xin_quad == 2) { // finally if we have two tracks, get x1 and x2
648 // and decide later how to use them
649 // don't think we want to experiment with discerning three overlapping
650 // tracks, so not >= 2
651 // continue; //debugging: skip for now -- didn't help
652 auto xitr1 = x_quad_map.lower_bound((*yitr).first);
653 auto xitr2 = x_quad_map.upper_bound((*yitr).first);
654 xitr2--; // upper_bound points to next element
655
656 if (xitr1 != xitr2) { // should be true already but...
657 x1 = ((*xitr1).second).getCentroidX() * x_conv_factor_ + x_start_;
658 x2 = ((*xitr2).second).getCentroidX() * x_conv_factor_ + x_start_;
659 sx1 = x_conv_factor_ / 2.; // 1 bar width
660 sx2 = sx1;
661 x = (x1 + x2) / 2.;
662 sx = fabs(x1 - x2) / 2 *
663 x_conv_factor_; // rely on x precision being one single pad width
664 if (verbose_)
665 ldmx_log(debug) << "\t\t -- 2 x in quad: setting y track x "
666 "coordinate to midpoint";
667 }
668 } // if 2 x tracks in quad
669
670 // ok! over y:
671 // can skip 0 y case by construction
672 if (n_yin_quad == 1) { // we can already now tell what the y coordinate and
673 // its precision is
674 y = ((*yitr).second).getCentroidY() * y_conv_factor_ + y_start_;
675 sy = ((*yitr).second).getResidual() * y_conv_factor_;
676 // if all clusters lined up, assign
677 // precision of 1 bar width
678 if (sy == 0) sy = 1. / 2 * y_conv_factor_;
679
680 if (n_xin_quad <= 1) {
681 // 4. every quadrant which just has one of each --> done ;
682 // b) set the sx, sy of the x track now, using the residuals from the
683 // other b1) special case: no x tracks; then x, sx have been set above
684 auto xidx = x_idx_quad_map.find((*yitr).first);
685 auto yidx = y_idx_quad_map.find((*yitr).first);
686 tracks.at((*xidx).second).setPosition(x, y);
687 tracks.at((*xidx).second).setSigmaXY(sx, sy);
688 tracks.at((*yidx).second).setPosition(x, y);
689 tracks.at((*yidx).second).setSigmaXY(sx, sy);
690 if (verbose_)
691 ldmx_log(debug) << "\t\t\t in quad " << (*yitr).first
692 << ", set (x, y) = (" << x << ", " << y
693 << ") and (sx, sy) = " << sx << ", " << sy << ")";
694 continue;
695 }
696 } // 1 y, 0 or 1 x track in quadrant
697
698 if (verbose_)
699 ldmx_log(debug) << "\t\t in quad " << (*yitr).first
700 << ", not single x,y tracks: " << n_xin_quad
701 << " of x and " << n_yin_quad << " of y";
702
703 if (n_yin_quad == 2) { // let's start here and see if we can do >= 2 later
704 // here one could do sth to avoid checking the other y track again in the
705 // outermost loop over y
706 auto yitr1 = y_quad_map.lower_bound((*yitr).first);
707 auto yitr2 = y_quad_map.upper_bound((*yitr).first);
708 yitr2--; // back up once
709 y1 = ((*yitr1).second).getCentroidY() * y_conv_factor_ + y_start_;
710 y2 = ((*yitr2).second).getCentroidY() * y_conv_factor_ + y_start_;
711 sy1 = ((*yitr1).second).getResidual() * y_conv_factor_;
712 sy2 = ((*yitr2).second).getResidual() * y_conv_factor_;
713 y = (y1 + y2) / 2.;
714 sy = fabs(y1 - y2) / 2 * y_conv_factor_;
715 if (verbose_)
716 ldmx_log(debug)
717 << "\t\t -- 2 y in quad: setting x track y coordinate to midpoint";
718 } // 2y in quad
719
720 if (n_yin_quad == 1 &&
721 n_xin_quad == 2) { // don't think we want to experiment with discerning
722 // three overlapping tracks, so not >= 2
723
724 // first: set the y track coordinates to x = the mid of x tracks, y = y
725 // of y track
726 auto yidx = y_idx_quad_map.find((*yitr).first);
727 tracks.at((*yidx).second).setPosition(x, y);
728 tracks.at((*yidx).second).setSigmaXY(sx, sy);
729
730 int min_overlap_pe = 250;
731 if (((*yitr).second).getPE() < min_overlap_pe) {
732 // can't tell, really, that either of these belong to the y track. so.
733 // let them keep their own x coordinate but set y to quadrant midpoint,
734 // with uncertainty +/- half quadrant width (1/8 of pad height)
735 y = y0;
736 sy = sy0;
737 if (verbose_)
738 ldmx_log(debug) << "\t\t -- Can't tell which x track should be "
739 "matched to single y track. Setting both x track "
740 "coordinates to y quadrant value:";
741 } // if can't assume overlap
742 else if (verbose_)
743 ldmx_log(debug) << "\t\t -- Found large PE count ("
744 << ((*yitr).second).getPE() << " > " << min_overlap_pe
745 << "), suggesting overlap! Setting both x track "
746 "coordinates to y track value:";
747
748 // consider making two x tracks out if this one, and, anyway have to set
749 // their average as the y track x cocordinate
750 // EXPERIMENTAL : apply only to x tracks, which can be disregarded for
751 // electron counting
752 if (verbose_)
753 ldmx_log(debug) << "\t\t -- (x1, x2, y) = (" << x1 << ", " << x2
754 << ", " << y << ") and (sx1, sx2, sy) = " << sx1 << ", "
755 << sx2 << ", " << sy << ")";
756
757 // now set x track coordinates according to overlap check result
758 auto xidx1 = x_idx_quad_map.lower_bound((*yitr).first);
759 auto xidx2 = x_idx_quad_map.upper_bound((*yitr).first);
760 xidx2--; // upper_bound points to (last+1) element
761 tracks.at((*xidx1).second).setPosition(x1, y);
762 tracks.at((*xidx1).second).setSigmaXY(sx1, sy);
763 tracks.at((*xidx2).second).setPosition(x2, y);
764 tracks.at((*xidx2).second).setSigmaXY(sx2, sy);
765
766 } // 1 y, 2 x tracks in the quadrant
767 else if (n_yin_quad == 2 && n_xin_quad == 1) {
768 // 5b) if there are more y than x: could be an overlap
769
770 // first: set the x track coordinates to x = x of x track, y = the mid of
771 // y tracks
772 auto xidx = x_idx_quad_map.find((*yitr).first);
773 tracks.at((*xidx).second).setPosition(x, y);
774 tracks.at((*xidx).second).setSigmaXY(sx, sy);
775
776 auto xitr = x_quad_map.lower_bound((*yitr).first);
777 int min_overlap_pe = 300;
778 if (((*xitr).second).getPE() < min_overlap_pe) {
779 if (verbose_)
780 ldmx_log(debug)
781 << "\t\t just 1 x track with not-unusual PE in the quad -- can't "
782 "match; setting mid-point values for x ";
783 x = x0;
784 sx = sx0;
785 } // if can't assume overlap
786 else {
787 // consider making two x tracks out if this one, and, anyway have to set
788 // their average as the y track x cocordinate
789 // EXPERIMENTAL : apply only to x tracks, which can be disregarded for
790 // electron counting
791 if (verbose_)
792 ldmx_log(debug) << "\t\t -- Found large PE count ("
793 << ((*xitr).second).getPE() << " > " << min_overlap_pe
794 << ") in x track, suggesting overlap! Setting both y "
795 "track coordinates to x track value:";
796 } // if can assume overlap
797 if (verbose_)
798 ldmx_log(debug) << "\t\t -- (x, y1, y2) = (" << x << ", " << y1 << ", "
799 << y2 << ") and (sx, sy1, sy2) = " << sx << ", " << sy1
800 << ", " << sy2 << ")";
801
802 auto yidx1 = y_idx_quad_map.lower_bound((*yitr).first);
803 auto yidx2 = y_idx_quad_map.upper_bound((*yitr).first);
804 yidx2--; // upper_bound points to next element
805 tracks.at((*yidx1).second).setPosition(x, y1);
806 tracks.at((*yidx1).second).setSigmaXY(sx, sy1);
807 tracks.at((*yidx2).second).setPosition(x, y2);
808 tracks.at((*yidx2).second).setSigmaXY(sx, sy2);
809
810 } // 2 y and 1 x track in quad
811 else if (n_yin_quad == 2 && n_xin_quad == 2) {
812 // MIDPONTS ALL OVER!
813 auto xidx1 = x_idx_quad_map.lower_bound((*yitr).first);
814 auto xidx2 = x_idx_quad_map.upper_bound((*yitr).first);
815 xidx2--;
816 auto yidx1 = y_idx_quad_map.lower_bound((*yitr).first);
817 auto yidx2 = y_idx_quad_map.upper_bound((*yitr).first);
818 yidx2--;
819
820 if (y_idx_quad_map.find((*yitr).first) == y_idx_quad_map.end())
821 ldmx_log(error) << "The two y tracks in the same quadrant at "
822 << (*yitr).first
823 << " appear to not be found in the y track map! "
824 "investigate. Note that yidx1.first = "
825 << (*yidx1).first
826 << " and yidx2.first = " << (*yidx2).first;
827 else {
828 tracks.at((*xidx1).second).setPosition(x1, y);
829 tracks.at((*xidx1).second).setSigmaXY(sx1, sy);
830 tracks.at((*xidx2).second).setPosition(x2, y);
831 tracks.at((*xidx2).second).setSigmaXY(sx2, sy);
832
833 tracks.at((*yidx1).second).setPosition(x, y1);
834 tracks.at((*yidx1).second).setSigmaXY(sx, sy1);
835 tracks.at((*yidx2).second).setPosition(x, y2);
836 tracks.at((*yidx2).second).setSigmaXY(sx, sy2);
837
838 if (verbose_)
839 ldmx_log(debug) << "\t\t -- in a 2 x 2 situaiton; midpoint y: " << y
840 << " for both x tracks, midpoint x: " << x
841 << " for both y tracks";
842 }
843 } // if 2 y, 2 x tracks
844
845 if (n_xin_quad > 2) {
846 if (verbose_)
847 ldmx_log(debug) << "\t\t -*-*-*- more than 2 x tracks in the same quad "
848 "-- nothing done about the x,y coordinates in this "
849 "situation -- implement if needed!!";
850 }
851 if (n_yin_quad > 2) {
852 if (verbose_)
853 ldmx_log(debug) << "\t\t -*-*-*- more than 2 y tracks in the same quad "
854 "-- nothing done about the x,y coordinates in this "
855 "situation -- implement if needed!!";
856 }
857
858 } // over y tracks
859
860 y_quad_map.clear();
861 x_quad_map.clear();
862
863 // return tracks;
864}

◆ onProcessEnd()

void trigscint::TrigScintTrackProducer::onProcessEnd ( )
overridevirtual

Callback for the EventProcessor to take any necessary action when the processing of events finishes, such as calculating job-summary quantities.

Reimplemented from framework::EventProcessor.

Definition at line 872 of file TrigScintTrackProducer.cxx.

872 {
873 ldmx_log(debug) << "Process ends!";
874
875 return;
876}

◆ onProcessStart()

void trigscint::TrigScintTrackProducer::onProcessStart ( )
overridevirtual

Callback for the EventProcessor to take any necessary action when the processing of events starts, such as creating histograms.

Reimplemented from framework::EventProcessor.

Definition at line 866 of file TrigScintTrackProducer.cxx.

866 {
867 ldmx_log(debug) << "Process starts!";
868
869 return;
870}

◆ produce()

void trigscint::TrigScintTrackProducer::produce ( framework::Event & event)
overridevirtual

Process the event and put new data products into it.

Parameters
eventThe Event to process.

Implements framework::Producer.

Definition at line 61 of file TrigScintTrackProducer.cxx.

61 {
62 // parameters.
63 // one pad cluster collection to use as seed
64 // a vector with the other two
65 // a maximum distance between seed centroid and other pad clusters
66 // allowed to belong to the same track
67 // an output collection name a verbosity controller
68
69 if (verbose_) {
70 ldmx_log(debug)
71 << "TrigScintTrackProducer: produce() starts! Event number: "
72 << event.getEventHeader().getEventNumber();
73 }
74 if (!event.exists(seeding_collection_, pass_name_)) {
75 ldmx_log(info) << "No collection called " << seeding_collection_
76 << "; skipping event";
77 // << "; still, not skipping event";
78 return;
79 }
80
81 if (!event.exists(seeding_collection_, pass_name_)) {
82 ldmx_log(info) << "No collection called " << seeding_collection_
83 << "; skipping event";
84 return;
85 }
86 const auto seeds{event.getCollection<ldmx::TrigScintCluster>(
87 seeding_collection_, pass_name_)};
88 uint num_seeds = seeds.size();
89
90 if (verbose_) {
91 ldmx_log(debug) << "Got track seeding cluster collection "
92 << seeding_collection_ << " with " << num_seeds
93 << " entries ";
94 }
95
96 if (!event.exists(input_collections_.at(0), pass_name_)) {
97 ldmx_log(info) << "No collection called " << input_collections_.at(0)
98 << "; skipping event";
99 // << "; still, not skipping event";
100
101 return;
102 }
103 const auto clusters_pad1{event.getCollection<ldmx::TrigScintCluster>(
104 input_collections_.at(0), pass_name_)};
105
106 if (!event.exists(input_collections_.at(1), pass_name_)) {
107 ldmx_log(info) << "No collection called "
108 << input_collections_.at(1)
109 // << "; still, not skipping event";
110 << "; skipping event";
111 std::vector<ldmx::TrigScintTrack> empty{};
112 event.add(output_collection_, empty);
113 return;
114 }
115
116 const auto clusters_pad2{event.getCollection<ldmx::TrigScintCluster>(
117 input_collections_.at(1), pass_name_)};
118
119 if (verbose_) {
120 ldmx_log(debug) << "Got the other two pad collections:"
121 << input_collections_.at(0) << " with "
122 << clusters_pad1.size() << " entries, and "
123 << input_collections_.at(1) << " with "
124 << clusters_pad2.size() << " entries.";
125 }
126
127 std::vector<ldmx::TrigScintTrack> cleaned_tracks;
128 std::vector<ldmx::TrigScintTrack> cleaned_tracks_y;
129 std::vector<ldmx::TrigScintTrack> cleaned_tracks_x;
130
131 // loop over the clusters in the seeding pad collection, if there are clusters
132 // in all pads
133 // bool skipDn = false;
134 if (num_seeds && clusters_pad1.size()) {
135 // could check this explicitly here: and then just get out of all checks on
136 // the dn pad immediately
137 // if (! clusters_pad2.size())
138 // skipDn = true ;
139 for (const auto &seed : seeds) {
140 // for each seed, search through the other two pads to match all clusters
141 // with centroids within tolerance to tracks
142 float centroid = seed.getCentroid();
143
144 std::vector<ldmx::TrigScintTrack> track_candidates;
145
146 if (verbose_ > 1) {
147 ldmx_log(debug) << "Got seed with centroid " << centroid;
148 }
149
150 // reset for each seed
151 // bool madeTrack = false;
152
153 for (const auto &cluster1 : clusters_pad1) {
154 if (verbose_ > 1) {
155 ldmx_log(debug) << "\tGot pad1 cluster with centroid "
156 << cluster1.getCentroid();
157 }
158 if (fabs(cluster1.getCentroid() - centroid) < max_delta_ ||
159 (centroid >= vert_bar_start_idx_ // then in vertical bars
160 && seed.getCentroidX() == cluster1.getCentroidX())) {
161 // use geometry y overlap scheme to see if this is really a match in x
162 // should be done in a map
163
164 if (centroid >= vert_bar_start_idx_ &&
165 seed.getCentroidY() < cluster1.getCentroidY()) {
166 // impossible combination
167 if (verbose_ > 1) {
168 ldmx_log(debug) << "\tSkipping impossible x cluster combination "
169 "with y flags (tag up) ("
170 << seed.getCentroidY() << " "
171 << cluster1.getCentroidY() << ")";
172 }
173 continue;
174 }
175
176 // else: first (possible) match! loop through next pad too
177
178 if (verbose_ > 1) {
179 ldmx_log(debug) << "\t\tIt is close enough!. Check pad2";
180 }
181
182 // try making third pad clusters an optional part of track
183
184 std::vector<ldmx::TrigScintCluster> cluster_vec = {seed, cluster1};
185
186 bool has_match_dn = false;
187
188 for (const auto &cluster2 : clusters_pad2) {
189 if (verbose_ > 1) {
190 ldmx_log(debug) << "\tGot pad2 cluster with centroid "
191 << cluster2.getCentroid();
192 }
193
194 if (fabs(cluster2.getCentroid() - centroid) < max_delta_ ||
195 (centroid >= vert_bar_start_idx_ // then in vertical bars
196 && seed.getCentroidX() == cluster2.getCentroidX())) {
197 // use geometry y overlap scheme to see if this is really a match
198 // in x
199
200 if (centroid >= vert_bar_start_idx_ &&
201 (seed.getCentroidY() < cluster2.getCentroidY() ||
202 cluster1.getCentroidY() >
203 cluster2.getCentroidY())) { // impossible
204 if (verbose_ > 1) {
205 ldmx_log(debug)
206 << "\tSkipping impossible x cluster combination with y "
207 "flags (tag up dn) ("
208 << seed.getCentroidY() << " " << cluster1.getCentroidY()
209 << " " << cluster2.getCentroidY() << ")";
210 }
211 continue;
212 }
213
214 // first match! loop through next pad too
215
216 if (verbose_ > 1) {
217 ldmx_log(debug) << "\t\tIt is close enough!. Make a track";
218 }
219
220 // only make this vector now! this ensures against hanging
221 // clusters with indices from earlier in the loop
222 std::vector<ldmx::TrigScintCluster> three_cluster_vec = {
223 seed, cluster1, cluster2};
224
225 /*
226 // here we could break if we didn't want to allow all possible
227 combinations madeTrack=true; break; //we're done with this
228 iteration once there's a track made
229 */
230 // make a track
231 ldmx::TrigScintTrack track = makeTrack(three_cluster_vec);
232 track_candidates.push_back(track);
233 has_match_dn = true;
234 } // if match in pad2
235 } // over clusters in pad2
236 // if there was no match to this in pad 2, make a track with just
237 // these two clusters
238 if (!has_match_dn && skip_last_) {
239 // we allow skipping last pad if needed
240 ldmx::TrigScintTrack track = makeTrack(cluster_vec);
241 track_candidates.push_back(track);
242 }
243
244 } // if possible (x,)y match in pad1
245 /*
246//same here
247if (madeTrack)
248break;
249*/
250
251 } // over clusters in pad1
252
253 // continue to next seed if 0 track candidates
254 if (track_candidates.size() == 0) continue;
255
256 int keep_idx = 0;
257 float min_residual = 1000; // some large number
258
259 // no need to choose between only one candidate track
260 if (track_candidates.size() > 1) {
261 // now for each seed, pick only the track with the smallest residual.
262
263 if (verbose_) {
264 ldmx_log(debug) << "Got " << track_candidates.size()
265 << " tracks to check.";
266 }
267
268 for (uint idx = 0; idx < track_candidates.size(); idx++) {
269 if ((track_candidates.at(idx)).getResidual() < min_residual) {
270 keep_idx = (int)idx;
271 min_residual =
272 (track_candidates.at(idx)).getResidual(); // update minimum
273
274 if (verbose_ > 1) {
275 ldmx_log(debug)
276 << "Track at index " << idx
277 << " has smallest residual so far: " << min_residual;
278 }
279
280 } // finding min residual
281 } // over track candidates
282 } // if more than 1 to choose from
283
284 // store the track at keepIdx, if there was one we made it this far and
285 // keepIdx is 0 or has been updated to the smallest residual track idx
286 // if (keepIdx >= 0) {
287 tracks_.push_back(track_candidates.at(keep_idx));
288 if (verbose_) {
289 ldmx_log(debug) << "Kept track at index " << keep_idx;
290 ldmx_log(trace) << track_candidates.at(keep_idx);
291 }
292 //}
293 } // over seeds
294
295 // done here if there were no tracks found
296 if (tracks_.size() == 0) {
297 if (verbose_) {
298 ldmx_log(debug) << "No tracks found!";
299 }
300 std::vector<ldmx::TrigScintTrack> empty{};
301 event.add(output_collection_, empty);
302 return;
303 }
304 // now, if there are multiple seeds sharing the same downstream hits, this
305 // should also be remedied with a selection on min residual.
306
307 // The logic of this loop kind of assumes I can remove tracks immediately --
308 // that way I can do pairwise checks between more tracks within a single
309 // loop. But for now I haven't figured out how to erase elements in a fool
310 // proof way. So I iterate over a vector...
311
312 std::vector keep_indices(tracks_.size(), 1);
313 if (verbose_ > 1)
314 ldmx_log(debug) << "vector of indices to keep has size "
315 << keep_indices.size();
316
317 for (uint idx = tracks_.size() - 1; idx > 0; idx--) {
318 // since we start in one end, we only have to check matches in one
319 // direction
320 ldmx::TrigScintTrack track = tracks_.at(idx);
321 for (int idx_comp = idx - 1; idx_comp >= 0; idx_comp--) {
322 if (verbose_ > 1)
323 ldmx_log(debug) << "In track disambiguation loop, idx points at "
324 << idx << " and prev idx points at " << idx_comp;
325
326 ldmx::TrigScintTrack next_track = tracks_.at(idx_comp);
327
328 // no need to start pulling constituents from tracks that are
329 // ridiculously far apart
330 if (fabs(track.getCentroid() - next_track.getCentroid()) <
331 3 * max_delta_) {
332 std::vector<ldmx::TrigScintCluster> consts_1 =
333 track.getConstituents();
334 std::vector<ldmx::TrigScintCluster> consts_2 =
335 next_track.getConstituents();
336 if (verbose_ > 1)
337 ldmx_log(debug)
338 << "In track disambiguation loop, got the two tracks, "
339 "with nConstituents "
340 << consts_1.size() << " and " << consts_2.size()
341 << ", respectively. ";
342 // let's do "if either cluster is shared" right now... but could also
343 // have it settable to use a stricter cut: an AND
344 if (consts_1[1].getCentroid() == consts_2[1].getCentroid() ||
345 (consts_1.size() > 2 && consts_2.size() > 2 &&
346 consts_1[2].getCentroid() ==
347 consts_2[2]
348 .getCentroid())) { // we have overlap downstream of the
349 // seeding pad. probably, one cluster
350 // in seeding pad is noise
351
352 if (verbose_ > 1) {
353 ldmx_log(debug) << "Found overlap! Tracks at index " << idx
354 << " and " << idx_comp;
355 ldmx_log(trace) << tracks_.at(idx);
356 ldmx_log(trace) << tracks_.at(idx_comp);
357 }
358
359 if ((tracks_.at(idx)).getResidual() <
360 (tracks_.at(idx_comp)).getResidual()) {
361 // next track (lower index) is a worse choice, remove its flag for
362 // keeping
363 keep_indices.at(idx_comp) = 0;
364 } else // prefer next track over current. remove current track's
365 // keep
366 // flag
367 keep_indices.at(idx) = 0;
368 /*}
369 else {
370 tracks_.erase(itNext);
371 // removeIdx.push_back(idx+1);
372 // we might see the same index two times in the loop in this
373 case, if there are three seeds sharing the same clusters
374 downstream.
375 // then the third only gets removed if it's even worse than
376 the second.
377 // one could deal with this with an extra overlap check. not
378 sure we will be in this situation any time soon though.
379 }*/
380 } // over matching/overlapping tracks
381 } // over tracks close enough to share constituents
382 } // over constructed tracks at other indices, to match
383 } // over constructed tracks
384
385 for (uint idx = 0; idx < tracks_.size(); idx++) {
386 if (verbose_ > 1) {
387 ldmx_log(debug) << "keep flag for idx " << idx << " is "
388 << keep_indices.at(idx);
389 }
390 if (keep_indices.at(idx)) { // this hasn't been flagged for removal
391
392 cleaned_tracks.push_back(tracks_.at(idx));
393
394 if (verbose_) {
395 ldmx_log(debug) << "After cleaning, keeping track at index " << idx
396 << ": Centroid = " << (tracks_.at(idx)).getCentroid()
397 << "; CentroidX = "
398 << (tracks_.at(idx)).getCentroidX()
399 << "; CentroidY = "
400 << (tracks_.at(idx)).getCentroidY()
401 << "; track PE = " << (tracks_.at(idx)).getPE()
402 << tracks_.at(idx);
403 }
404 } // if index flagged for keeping
405 } // over all (uniquely seeded) tracks in the event
406
407 if (verbose_) {
408 for (uint idx = 0; idx < tracks_.size(); idx++) {
409 ldmx_log(debug) << "Keeping track at index " << idx << ":"
410 << tracks_.at(idx);
411 }
412 }
413
414 if (verbose_) {
415 ldmx_log(debug) << "Running track x,y matching ";
416 }
417
418 if (cleaned_tracks.size() > 0) {
419 matchXYTracks(cleaned_tracks);
420 std::vector<ldmx::TrigScintTrack> matched_tracks =
421 cleaned_tracks; // don't know why this copying needs to happen but it
422 // does
423 // std::vector<ldmx::TrigScintTrack> matchXYTracks( cleanedTracks
424 //); std::vector<ldmx::TrigScintTrack> matchedTracks =
425 // matchXYTracks( cleanedTracks );
426 for (auto trk : matched_tracks) {
427 /* for (uint idx = 0; idx < tracks_.size(); idx++) {
428 if (verbose_ > 1) {
429 ldmx_log(debug) << "keep flag for idx " << idx << " is "
430 << keepIndices.at(idx);
431 }
432 if (keepIndices.at(idx)) { // this hasn't
433 been flagged for removal
434 //check if channel nb is above that of horizontal bars
435 if (tracks_.at(idx).getCentroid() >= vert_bar_start_idx_)
436 */
437 if (trk.getCentroid() >= vert_bar_start_idx_)
438 cleaned_tracks_x.push_back(trk); // acks_.at(idx));
439 else
440 cleaned_tracks_y.push_back(trk); // acks_.at(idx));
441 // cleanedTracksY.push_back(trk);
442 if (verbose_ > 1) {
443 float centr = trk.getCentroid(); // tracks_.at(idx).getCentroid(); //
444 std::string coll_str = centr >= vert_bar_start_idx_ ? "X" : "Y";
445 coll_str = output_collection_ + coll_str;
446 ldmx_log(debug) << "saving track with centroid " << centr
447 << " to output track collection " << coll_str;
448 }
449 // }
450 }
451 }
452
453 } // if there are clusters in all pads
454 else if (verbose_) {
455 ldmx_log(info)
456 << "Not all pads had clusters; (maybe) skipping tracking attempt";
457 }
458
459 if (verbose_) {
460 ldmx_log(debug) << "Done with tracking step. ";
461 }
462
463 event.add(output_collection_, cleaned_tracks);
464 // event.add(output_collection_, matchedTracks);
465
466 event.add(output_collection_ + "Y", cleaned_tracks_y);
467 event.add(output_collection_ + "X", cleaned_tracks_x);
468
469 tracks_.resize(0);
470
471 return;
472}
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
Stores cluster information from the trigger scintillator pads.
float getCentroid() const
Get the detector ID centroid of the track.
std::vector< ldmx::TrigScintCluster > getConstituents() const
Get the cluster constituents of the track.

References framework::Event::exists(), ldmx::TrigScintTrack::getCentroid(), and ldmx::TrigScintTrack::getConstituents().

Member Data Documentation

◆ bar_gap_x_

float trigscint::TrigScintTrackProducer::bar_gap_x_ {0.1}
private

Definition at line 91 of file TrigScintTrackProducer.h.

91{0.1}; // mm

◆ bar_gap_y_

float trigscint::TrigScintTrackProducer::bar_gap_y_ {2.1}
private

Definition at line 89 of file TrigScintTrackProducer.h.

89{2.1}; // mm

◆ bar_width_x_

float trigscint::TrigScintTrackProducer::bar_width_x_ {3.}
private

Definition at line 90 of file TrigScintTrackProducer.h.

90{3.}; // mm

◆ bar_width_y_

float trigscint::TrigScintTrackProducer::bar_width_y_ {3.}
private

Definition at line 88 of file TrigScintTrackProducer.h.

88{3.}; // mm

◆ input_collections_

std::vector<std::string> trigscint::TrigScintTrackProducer::input_collections_
private

Definition at line 54 of file TrigScintTrackProducer.h.

◆ max_delta_

double trigscint::TrigScintTrackProducer::max_delta_ {0.}
private

Definition at line 45 of file TrigScintTrackProducer.h.

45{0.};

◆ n_bars_x_

int trigscint::TrigScintTrackProducer::n_bars_x_ {8}
private

Definition at line 72 of file TrigScintTrackProducer.h.

72{8};

◆ n_bars_y_

int trigscint::TrigScintTrackProducer::n_bars_y_ {16}
private

Definition at line 69 of file TrigScintTrackProducer.h.

69{16};

◆ output_collection_

std::string trigscint::TrigScintTrackProducer::output_collection_
private

Definition at line 57 of file TrigScintTrackProducer.h.

◆ pass_name_

std::string trigscint::TrigScintTrackProducer::pass_name_ {""}
private

Definition at line 60 of file TrigScintTrackProducer.h.

60{""};

◆ seeding_collection_

std::string trigscint::TrigScintTrackProducer::seeding_collection_
private

Definition at line 51 of file TrigScintTrackProducer.h.

◆ skip_last_

bool trigscint::TrigScintTrackProducer::skip_last_ {false}
private

Definition at line 63 of file TrigScintTrackProducer.h.

63{false};

◆ tracks_

std::vector<ldmx::TrigScintTrack> trigscint::TrigScintTrackProducer::tracks_
private

Definition at line 33 of file TrigScintTrackProducer.h.

◆ verbose_

int trigscint::TrigScintTrackProducer::verbose_ {0}
private

Definition at line 48 of file TrigScintTrackProducer.h.

48{0};

◆ vert_bar_start_idx_

int trigscint::TrigScintTrackProducer::vert_bar_start_idx_ {52}
private

Definition at line 66 of file TrigScintTrackProducer.h.

66{52};

◆ x_conv_factor_

float trigscint::TrigScintTrackProducer::x_conv_factor_
private

Definition at line 93 of file TrigScintTrackProducer.h.

◆ x_start_

float trigscint::TrigScintTrackProducer::x_start_
private

Definition at line 94 of file TrigScintTrackProducer.h.

◆ y_conv_factor_

float trigscint::TrigScintTrackProducer::y_conv_factor_
private

Definition at line 95 of file TrigScintTrackProducer.h.

◆ y_start_

float trigscint::TrigScintTrackProducer::y_start_
private

Definition at line 96 of file TrigScintTrackProducer.h.


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