Belle II Software development
NeuroTrigger Class Reference

Class to represent the CDC Neurotrigger. More...

#include <NeuroTrigger.h>

Classes

struct  Parameters
 Struct to keep neurotrigger parameters. More...
 

Public Member Functions

 NeuroTrigger ()
 Default constructor.
 
virtual ~NeuroTrigger ()
 Default destructor.
 
void initialize (const Parameters &p)
 Set parameters and get some network independent parameters.
 
void initialize (const NeuroTriggerParameters &p)
 
std::vector< unsigned > getRangeIndices (const NeuroTriggerParameters &p, unsigned isector)
 Get indices for sector ranges in parameter lists.
 
std::vector< unsigned > getRangeIndices (const Parameters &p, unsigned isector)
 
void save (const std::string &filename, const std::string &arrayname="MLPs")
 Save MLPs to file.
 
bool loadIDHist (const std::string &filename)
 function to load idhist from file
 
bool load (const std::string &filename, const std::string &arrayname="MLPs")
 Load MLPs from file.
 
void setConstants ()
 Loads parameters from the geometry and precalculates some constants that will be needed.
 
void setPrecision (const std::vector< unsigned > &precision)
 set fixed point precision
 
void initializeCollections (std::string hitCollectionName, std::string eventTimeName, const std::string &et_option)
 set the hit collection and event time to required and store the hit collection name
 
void initializeCollections (std::string hitCollectionName)
 
CDCTriggerMLPoperator[] (unsigned index)
 return reference to a neural network
 
const CDCTriggerMLPoperator[] (unsigned index) const
 return const reference to a neural network
 
unsigned nSectors () const
 return number of neural networks
 
void addMLP (const CDCTriggerMLP &newMLP)
 add an MLP to the list of networks
 
std::vector< int > selectMLPs (float phi0, float invpt, float theta)
 Select all matching expert MLPs based on the given track parameters.
 
std::vector< int > selectMLPsTrain (float phi0, float invpt, float theta)
 Select all matching expert MLPs based on the given track parameters.
 
int selectMLPbyPattern (std::vector< int > &MLPs, unsigned long pattern, const bool neurotrackinputmode)
 Select one MLP from a list of sector indices.
 
void updateTrack (const CDCTriggerTrack &track)
 Calculate 2D phi position and arclength for the given track and store them.
 
void updateTrackFix (const CDCTriggerTrack &track)
 Calculate 2D phi position and arclength for the given track and store them.
 
double getRelId (const CDCTriggerSegmentHit &hit)
 Calculate phi position of a hit relative to 2D track (scaled to number of wires).
 
int getLowestTime (unsigned isector, RelationVector< CDCTriggerSegmentHit > Hits, bool onlyAxials)
 helper function to get the fastest priority time of given ts array
 
void getEventTime (unsigned isector, const CDCTriggerTrack &track, std::string et_option, const bool)
 Read out the event time and store it.
 
void getEventTime (unsigned isector, const CDCTriggerTrack &track)
 DEPRECATED!! Read out the event time and store it.
 
std::string get_et_option ()
 Return value of m_et_option.
 
unsigned long getInputPattern (unsigned isector, const CDCTriggerTrack &track, const bool neurotrackinputmode)
 Calculate input pattern for MLP.
 
unsigned long getCompleteHitPattern (unsigned isector, const CDCTriggerTrack &track, const bool neurotrackinputmode)
 Get complete hit pattern of neurotrack.
 
unsigned long getPureDriftThreshold (unsigned isector, const CDCTriggerTrack &track, const bool neurotrackinputmode)
 Get the drift threshold bits, where the time of the TS was outside of the accepted time window and thus shifted to the allowed maximum within the borders.
 
std::vector< unsigned > selectHitsHWSim (unsigned isector, const CDCTriggerTrack &track)
 Select hits for each super layer from the ones related to input track.
 
std::vector< unsigned > selectHits (unsigned isector, const CDCTriggerTrack &track, bool returnAllRelevant=false)
 Select best hits for each super layer.
 
std::vector< float > getInputVector (unsigned isector, const std::vector< unsigned > &hitIds)
 Calculate input values for MLP.
 
std::vector< float > runMLP (unsigned isector, const std::vector< float > &input)
 Run an expert MLP.
 
std::vector< float > runMLPFix (unsigned isector, std::vector< float > input)
 Run an expert MLP with fixed point arithmetic.
 

Private Attributes

std::vector< CDCTriggerMLPm_MLPs = {}
 List of networks.
 
double m_radius [9][2] = {}
 Radius of the CDC layers with priority wires (2 per super layer)
 
unsigned m_TSoffset [10] = {}
 Number of track segments up to super layer.
 
double m_idRef [9][2] = {}
 2D phi position of current track scaled to number of wires
 
double m_alpha [9][2] = {}
 2D crossing angle of current track
 
int m_T0 = 0
 Event time of current event / track.
 
bool m_hasT0 = false
 Flag to show if stored event time is valid.
 
std::vector< unsigned > m_precision
 Fixed point precision in bit after radix point.
 
StoreArray< CDCTriggerSegmentHitm_segmentHits
 StoreArray containing the input track segment hits.
 
StoreObjPtr< BinnedEventT0m_eventTime
 StoreObjPtr containing the event time.
 
std::string m_hitCollectionName
 Name of the StoreArray containing the input track segment hits.
 
DBObjPtr< CDCTriggerNeuroConfigm_cdctriggerneuroconfig
 get NNT payload from database.
 

Detailed Description

Class to represent the CDC Neurotrigger.

The Neurotrigger consists of one or several Multi Layer Perceptrons. The input values are calculated from track segment hits and a 2D track estimate. The output is a scaled estimate of the z-vertex of the track. In case of several MLPs, each is an expert for a different track parameter region.

See also
Neurotrigger Modules:
NeuroTriggerTrainerModule for preparing training data and training,
NeuroTriggerModule for loading trained networks and using them.

Definition at line 40 of file NeuroTrigger.h.

Constructor & Destructor Documentation

◆ NeuroTrigger()

NeuroTrigger ( )
inline

Default constructor.

Definition at line 118 of file NeuroTrigger.h.

118{}

◆ ~NeuroTrigger()

virtual ~NeuroTrigger ( )
inlinevirtual

Default destructor.

Definition at line 121 of file NeuroTrigger.h.

121{}

Member Function Documentation

◆ addMLP()

void addMLP ( const CDCTriggerMLP newMLP)
inline

add an MLP to the list of networks

Definition at line 167 of file NeuroTrigger.h.

167{ m_MLPs.push_back(newMLP); }
std::vector< CDCTriggerMLP > m_MLPs
List of networks.
Definition: NeuroTrigger.h:293

◆ get_et_option()

std::string get_et_option ( )
inline

Return value of m_et_option.

Definition at line 229 of file NeuroTrigger.h.

230 {
231 std::string eto = m_MLPs[0].get_et_option();
232 for (unsigned int i = 0; i < m_MLPs.size(); ++i) {
233 if (m_MLPs[i].get_et_option() != eto) {
234 B2ERROR("Timing options in the expert networks in the CDC Neurotrigger differ!");
235 }
236 }
237 return eto;
238
239 }
std::string get_et_option()
Return value of m_et_option.
Definition: NeuroTrigger.h:229

◆ getCompleteHitPattern()

unsigned long getCompleteHitPattern ( unsigned  isector,
const CDCTriggerTrack track,
const bool  neurotrackinputmode 
)

Get complete hit pattern of neurotrack.

This does the same as the getInputPattern function, but also shows the axial hit bits. This function was made for the simulation of the hardware debug information "TSVector".

Definition at line 679 of file NeuroTrigger.cc.

680{
681 const CDCTriggerMLP& expert = m_MLPs[isector];
682 unsigned long chitPattern = 0;
683 vector<unsigned> nHits;
684 nHits.assign(9, 0);
685 // loop over axial hits related to input track
687 track.getRelationsTo<CDCTriggerSegmentHit>(m_hitCollectionName);
688 for (unsigned ihit = 0; ihit < axialHits.size(); ++ ihit) {
689 // skip hits with negative relation weight (not selected in finder)
690 if (axialHits.weight(ihit) < 0) {
691 continue;
692 }
693 unsigned short iSL = axialHits[ihit]->getISuperLayer();
694 // // skip stereo hits (should not be related to track, but check anyway)
695 if ((!neurotrackinputmode) && (iSL % 2 == 1)) continue;
696 double relId = getRelId(*axialHits[ihit]);
697 if (expert.isRelevant(relId, iSL)) {
698 if (nHits[iSL] < expert.getMaxHitsPerSL()) {
699 chitPattern |= 1 << (iSL + 9 * nHits[iSL]);
700 ++nHits[iSL];
701 }
702 }
703 }
704 if (!neurotrackinputmode) {
705 // loop over stereo hits
706 for (int ihit = 0; ihit < m_segmentHits.getEntries(); ++ ihit) {
707 unsigned short iSL = m_segmentHits[ihit]->getISuperLayer();
708 // skip axial hits
709 if (iSL % 2 == 0) continue;
710 // get priority time
711 double relId = getRelId(*m_segmentHits[ihit]);
712 if (expert.isRelevant(relId, iSL)) {
713 if (nHits[iSL] < expert.getMaxHitsPerSL()) {
714 chitPattern |= 1 << (iSL + 9 * nHits[iSL]);
715 ++nHits[iSL];
716 }
717 }
718 }
719 }
720 return chitPattern;
721}
Class to keep all parameters of an expert MLP for the neuro trigger.
Definition: CDCTriggerMLP.h:20
Combination of several CDCHits to a track segment hit for the trigger.
double getRelId(const CDCTriggerSegmentHit &hit)
Calculate phi position of a hit relative to 2D track (scaled to number of wires).
StoreArray< CDCTriggerSegmentHit > m_segmentHits
StoreArray containing the input track segment hits.
Definition: NeuroTrigger.h:316
std::string m_hitCollectionName
Name of the StoreArray containing the input track segment hits.
Definition: NeuroTrigger.h:320
Class for type safe access to objects that are referred to in relations.
size_t size() const
Get number of relations.
float weight(int index) const
Get weight with index.
int getEntries() const
Get the number of objects in the array.
Definition: StoreArray.h:216

◆ getEventTime() [1/2]

void getEventTime ( unsigned  isector,
const CDCTriggerTrack track 
)

DEPRECATED!! Read out the event time and store it.

If there is no valid event time, it can be determined from the shortest priority time of all hit candidates, if the option is enabled for the given sector.

◆ getEventTime() [2/2]

void getEventTime ( unsigned  isector,
const CDCTriggerTrack track,
std::string  et_option,
const bool  neuroinputmode = false 
)

Read out the event time and store it.

It can be given different options in the et_option ("EventTime option") parameter. The different options are: "etf_only" : only ETF info is used, otherwise an error is thrown. "fastestpriority" : event time is estimated by fastest priority time in selected track segments. if something fails, it is set to 0. "zero" : the event time is set to 0. "etf_or_fastestpriority" : the event time is obtained by the ETF, if not possible, the flag "fastestppriority" is used. "etf_or_zero" : the event time is obtained by the ETF, if

       m_T0 = 9999;

find shortest time of related and relevant axial hits RelationVector<CDCTriggerSegmentHit> Hits = track.getRelationsTo<CDCTriggerSegmentHit>(m_hitCollectionName); m_T0 = getLowestTime(isector, Hits, false); if (m_T0 < 9999) { m_hasT0 = true; } else { m_T0 = 0; m_hasT0 = false; }

  m_T0 = 9999;

find shortest time of related and relevant axial hits RelationVector<CDCTriggerSegmentHit> Hits = track.getRelationsTo<CDCTriggerSegmentHit>(m_hitCollectionName); m_T0 = getLowestTime(isector, Hits, true); if (m_T0 < 9999) { m_hasT0 = true; } else { m_T0 = 0; m_hasT0 = false; }

Definition at line 485 of file NeuroTrigger.cc.

486{
487
488 if (et_option != m_MLPs[isector].get_et_option()) {
489 B2DEBUG(20, "Used event time option is different to the one set in the MLP"
490 << LogVar("et_option", et_option) << LogVar("isector", isector)
491 << LogVar("et_option_mlp", m_MLPs[isector].get_et_option()));
492 }
493 if (et_option == "fastestpriority") {
494 B2DEBUG(200, "et_option is 'fastestpriority'");
495 m_T0 = 9999;
496 // find shortest time of related and relevant axial hits
498 track.getRelationsTo<CDCTriggerSegmentHit>(m_hitCollectionName);
499 m_T0 = getLowestTime(isector, Hits, false);
500 if (m_T0 < 9999) {
501 m_hasT0 = true;
502 } else {
503 m_T0 = 0;
504 m_hasT0 = false;
505 }
506
507 } else if (et_option == "fastest2d") {
508 m_T0 = 9999;
509 // find shortest time of related and relevant axial hits
511 track.getRelationsTo<CDCTriggerSegmentHit>(m_hitCollectionName);
512 m_T0 = getLowestTime(isector, Hits, true);
513 if (m_T0 < 9999) {
514 m_hasT0 = true;
515 } else {
516 m_T0 = 0;
517 m_hasT0 = false;
518 }
519 } else if (et_option == "zero") {
520 m_hasT0 = true;
521 m_T0 = 0;
522 } else if (et_option == "etf_only") {
523 bool hasT0 = (m_eventTime.isValid()) ? m_eventTime->hasBinnedEventT0(Const::CDC) : false;
524 if (hasT0) {
525 m_T0 = m_eventTime->getBinnedEventT0(Const::CDC);
526 m_hasT0 = true;
527 } else {
528 B2ERROR("No ETF output, but forced to use ETF!");
529 m_hasT0 = false;
530 m_T0 = 0;
531 }
532 } else if (et_option == "etf_or_fastestpriority") {
533 bool hasT0 = (m_eventTime.isValid()) ? m_eventTime->hasBinnedEventT0(Const::CDC) : false;
534 if (hasT0) {
535 m_T0 = m_eventTime->getBinnedEventT0(Const::CDC);
536 m_hasT0 = true;
537 } else {
538 getEventTime(isector, track, "fastestpriority", neuroinputmode);
551 }
552 } else if (et_option == "min_etf_fastestpriority") {
553 bool hasT0 = (m_eventTime.isValid()) ? m_eventTime->hasBinnedEventT0(Const::CDC) : false;
554 int T0_etf = 9999;
555 if (hasT0) {
556 T0_etf = m_eventTime->getBinnedEventT0(Const::CDC);
557 m_hasT0 = true;
558 }
559 getEventTime(isector, track, "fastestpriority", neuroinputmode);
560 if (m_T0 > T0_etf) {
561 m_T0 = T0_etf;
562 }
563 } else if (et_option == "etf_or_fastest2d") {
564 bool hasT0 = (m_eventTime.isValid()) ? m_eventTime->hasBinnedEventT0(Const::CDC) : false;
565 if (hasT0) {
566 m_T0 = m_eventTime->getBinnedEventT0(Const::CDC);
567 m_hasT0 = true;
568 } else {
569 getEventTime(isector, track, "fastest2d", neuroinputmode);
582 }
583 } else if (et_option == "etf_or_zero") {
584 bool hasT0 = (m_eventTime.isValid()) ? m_eventTime->hasBinnedEventT0(Const::CDC) : false;
585 if (hasT0) {
586 m_T0 = m_eventTime->getBinnedEventT0(Const::CDC);
587 m_hasT0 = true;
588 } else {
589 m_hasT0 = true;
590 m_T0 = 0;
591 }
592 } else if (et_option == "etfcc") {
593 if (track.getHasETFTime()) {
594 m_T0 = track.getETF_unpacked();
595 m_hasT0 = true;
596 } else {
597 m_T0 = 0;
598 m_hasT0 = false;
599 }
600
601 } else if (et_option == "etfcc_or_zero") {
602 if (track.getHasETFTime()) {
603 m_T0 = track.getETF_unpacked();
604 m_hasT0 = true;
605 } else {
606 m_T0 = 0;
607 m_hasT0 = true;
608 }
609
610 } else if (et_option == "etfcc_or_fastestpriority") {
611 if (track.getHasETFTime()) {
612 m_T0 = track.getETF_unpacked();
613 m_hasT0 = true;
614 } else {
615 getEventTime(isector, track, "fastestpriority", neuroinputmode);
616 }
617
618 } else if (et_option == "etfhwin") {
619 m_T0 = track.getETF_recalced();
620 m_hasT0 = true;
621
622 } else {
623 B2ERROR("No valid parameter for et_option (" << et_option << " )!");
624 }
625
626}
bool hasBinnedEventT0(const Const::EDetector detector) const
Check if one of the detectors in the given set has a binned t0 estimation.
int getBinnedEventT0(const Const::EDetector detector) const
Return the stored binned event t0 for the given detector or 0 if nothing stored.
int m_T0
Event time of current event / track.
Definition: NeuroTrigger.h:303
StoreObjPtr< BinnedEventT0 > m_eventTime
StoreObjPtr containing the event time.
Definition: NeuroTrigger.h:318
void getEventTime(unsigned isector, const CDCTriggerTrack &track, std::string et_option, const bool)
Read out the event time and store it.
int getLowestTime(unsigned isector, RelationVector< CDCTriggerSegmentHit > Hits, bool onlyAxials)
helper function to get the fastest priority time of given ts array
bool m_hasT0
Flag to show if stored event time is valid.
Definition: NeuroTrigger.h:305
bool isValid() const
Check whether the object was created.
Definition: StoreObjPtr.h:111
Class to store variables with their name which were sent to the logging service.

◆ getInputPattern()

unsigned long getInputPattern ( unsigned  isector,
const CDCTriggerTrack track,
const bool  neurotrackinputmode 
)

Calculate input pattern for MLP.

Parameters
isectorindex of the MLP that will use the input
trackaxial hit relations are taken from given track
neurotrackinputmodeinput mode
Returns
super layer pattern of hits in the current track

Definition at line 724 of file NeuroTrigger.cc.

725{
726 const CDCTriggerMLP& expert = m_MLPs[isector];
727 unsigned long hitPattern = 0;
728 vector<unsigned> nHits;
729 nHits.assign(9, 0);
730 // loop over axial hits related to input track
732 track.getRelationsTo<CDCTriggerSegmentHit>(m_hitCollectionName);
733 for (unsigned ihit = 0; ihit < axialHits.size(); ++ ihit) {
734 // skip hits with negative relation weight (not selected in finder)
735 if (axialHits.weight(ihit) < 0) {
736 continue;
737 }
738 unsigned short iSL = axialHits[ihit]->getISuperLayer();
739 // // skip stereo hits (should not be related to track, but check anyway)
740 if ((!neurotrackinputmode) && (iSL % 2 == 1)) continue;
741 // get priority time
742 int t = (m_hasT0) ? axialHits[ihit]->priorityTime() - m_T0 : 0;
743 if (t < 0) {
744 t = 0;
745 } else if (t > expert.getTMax()) {
746 t = expert.getTMax();
747 }
748 double relId = getRelId(*axialHits[ihit]);
749
750 if (expert.isRelevant(relId, iSL)) {
751 if (nHits[iSL] < expert.getMaxHitsPerSL()) {
752 hitPattern |= 1 << (iSL + 9 * nHits[iSL]);
753 ++nHits[iSL];
754 }
755 B2DEBUG(250, "hit in SL " << iSL);
756 } else {
757 B2DEBUG(250, "hit in SL " << iSL << " not relevant (relId = " << relId << ")");
758 }
759 }
760 if (!neurotrackinputmode) {
761 // loop over stereo hits
762 for (int ihit = 0; ihit < m_segmentHits.getEntries(); ++ ihit) {
763 unsigned short iSL = m_segmentHits[ihit]->getISuperLayer();
764 // skip axial hits
765 if (iSL % 2 == 0) continue;
766 // get priority time
767 int t = (m_hasT0) ? m_segmentHits[ihit]->priorityTime() - m_T0 : 0;
768 if (t < 0) {
769 t = 0;
770 } else if (t > expert.getTMax()) {
771 t = expert.getTMax();
772 }
773 double relId = getRelId(*m_segmentHits[ihit]);
774 if (expert.isRelevant(relId, iSL)) {
775 if (nHits[iSL] < expert.getMaxHitsPerSL()) {
776 hitPattern |= 1 << (iSL + 9 * nHits[iSL]);
777 ++nHits[iSL];
778 }
779 B2DEBUG(250, "hit in SL " << iSL);
780 } else {
781 B2DEBUG(250, "hit in SL " << iSL << " not relevant (relId = " << relId << ")");
782 }
783 }
784 }
785 B2DEBUG(250, "hitPattern " << hitPattern);
786 return hitPattern & expert.getSLpatternMask();
787}

◆ getInputVector()

vector< float > getInputVector ( unsigned  isector,
const std::vector< unsigned > &  hitIds 
)

Calculate input values for MLP.

Parameters
isectorindex of the MLP that will use the input
hitIdshit indices to be used for the input
Returns
scaled vector of input values (1 for each input node)

Definition at line 1007 of file NeuroTrigger.cc.

1008{
1009 const CDCTriggerMLP& expert = m_MLPs[isector];
1010 // prepare empty input vector and vectors to keep best drift times
1011 vector<float> inputVector;
1012 inputVector.assign(expert.nNodesLayer(0), 0.);
1013 // convert hits to input values
1014 vector<unsigned> nHits;
1015 nHits.assign(9, 0);
1016 for (unsigned ii = 0; ii < hitIds.size(); ++ii) {
1017 int ihit = hitIds[ii];
1018 unsigned short iSL = m_segmentHits[ihit]->getISuperLayer();
1019 unsigned short iRef = iSL + 9 * nHits[iSL];
1020 ++nHits[iSL];
1021 int priot = m_segmentHits[ihit]->priorityTime();
1022 int t = (m_hasT0) ? priot - m_T0 : 0;
1023 if (t < 0) {
1024 t = 0;
1025 } else if (t > expert.getTMax()) {
1026 t = expert.getTMax();
1027 }
1028 int LR = m_segmentHits[ihit]->getLeftRight();
1029 double relId = getRelId(*m_segmentHits[ihit]);
1030 int priority = m_segmentHits[ihit]->getPriorityPosition();
1031 // get scaled input values (scaling such that the absolute value of all inputs is < 1)
1032 inputVector[3 * iRef] = expert.scaleId(relId, iSL);
1033 float scaleT = pow(2, floor(log2(1. / expert.getTMax())));
1034 inputVector[3 * iRef + 1] = (((LR >> 1) & 1) - (LR & 1)) * t * scaleT;
1035 inputVector[3 * iRef + 2] = m_alpha[iSL][int(priority < 3)] * 0.5;
1036 }
1037 return inputVector;
1038}
double m_alpha[9][2]
2D crossing angle of current track
Definition: NeuroTrigger.h:301

◆ getLowestTime()

int getLowestTime ( unsigned  isector,
RelationVector< CDCTriggerSegmentHit Hits,
bool  onlyAxials = false 
)

helper function to get the fastest priority time of given ts array

Definition at line 461 of file NeuroTrigger.cc.

462{
463 int tlow = 9999;
464 B2DEBUG(200, "looping over axials:");
465 for (unsigned ihit = 0; ihit < Hits.size(); ++ihit) {
466 // skip hits with negative relation weight (not selected in finder)
467 if (Hits.weight(ihit) < 0) continue;
468 unsigned short iSL = Hits[ihit]->getISuperLayer();
469 if (iSL % 2 == 1 && onlyAxials) {continue;}
470 // get shortest time of relevant hits
471 B2DEBUG(200, " check drifttime: SL" + std::to_string(iSL) + ",ID = " + std::to_string(Hits[ihit]->getSegmentID()) + ", t = " +
472 std::to_string(Hits[ihit]->priorityTime()));
473 double relId = getRelId(*Hits[ihit]);
474 if (m_MLPs[isector].isRelevant(relId, iSL) &&
475 Hits[ihit]->priorityTime() < tlow) {
476 tlow = Hits[ihit]->priorityTime();
477 B2DEBUG(200, " new tlow: " << std::to_string(tlow));
478 }
479 }
480 return tlow;
481}

◆ getPureDriftThreshold()

unsigned long getPureDriftThreshold ( unsigned  isector,
const CDCTriggerTrack track,
const bool  neurotrackinputmode 
)

Get the drift threshold bits, where the time of the TS was outside of the accepted time window and thus shifted to the allowed maximum within the borders.

Note, that to get the same values as from the unpacker, this value has to be combined with the (complement of the) TSVector.

Definition at line 629 of file NeuroTrigger.cc.

630{
631 const CDCTriggerMLP& expert = m_MLPs[isector];
632 unsigned long driftth = 0;
633 vector<unsigned> nHits;
634 nHits.assign(9, 0);
635 // loop over axial hits related to input track
637 track.getRelationsTo<CDCTriggerSegmentHit>(m_hitCollectionName);
638 for (unsigned ihit = 0; ihit < axialHits.size(); ++ ihit) {
639 // skip hits with negative relation weight (not selected in finder)
640 if (axialHits.weight(ihit) < 0) {
641 continue;
642 }
643 unsigned short iSL = axialHits[ihit]->getISuperLayer();
644 // // skip stereo hits (should not be related to track, but check anyway)
645 if ((!neurotrackinputmode) && (iSL % 2 == 1)) continue;
646 // get priority time
647 int t = (m_hasT0) ? axialHits[ihit]->priorityTime() - m_T0 : 0;
648 double relId = getRelId(*axialHits[ihit]);
649 if (t < 0 || t > expert.getTMax()) {
650 if (expert.isRelevant(relId, iSL)) {
651 if (nHits[iSL] < expert.getMaxHitsPerSL()) {
652 driftth |= 1 << (8 - iSL + 9 * nHits[iSL]);
653 ++nHits[iSL];
654 }
655 }
656 }
657 }
658 if (!neurotrackinputmode) {
659 // loop over stereo hits
660 for (int ihit = 0; ihit < m_segmentHits.getEntries(); ++ ihit) {
661 unsigned short iSL = m_segmentHits[ihit]->getISuperLayer();
662 // skip axial hits
663 if (iSL % 2 == 0) continue;
664 // get priority time
665 int t = (m_hasT0) ? m_segmentHits[ihit]->priorityTime() - m_T0 : 0;
666 double relId = getRelId(*m_segmentHits[ihit]);
667 if (t < 0 || t > expert.getTMax()) {
668 if (expert.isRelevant(relId, iSL)) {
669 if (nHits[iSL] < expert.getMaxHitsPerSL()) {
670 driftth |= 1 << (8 - iSL + 9 * nHits[iSL]);
671 ++nHits[iSL];
672 }
673 }
674 }
675 }
676 }
677 return driftth;
678}

◆ getRangeIndices()

vector< unsigned > getRangeIndices ( const NeuroTriggerParameters p,
unsigned  isector 
)

Get indices for sector ranges in parameter lists.

Definition at line 260 of file NeuroTrigger.cc.

261{
262 // the indices can be used for both rangeuse and rangetrain, because the size of those arrays should be the same (it is checked in the initialize function).
263 std::vector<unsigned> indices = {0, 0, 0, 0};
264 if (p.phiRangeUse.size() * p.invptRangeUse.size() * p.thetaRangeUse.size() * p.SLpattern.size() == p.nMLP) {
265 indices[0] = isector % p.phiRangeUse.size();
266 indices[1] = (isector / p.phiRangeUse.size()) % p.invptRangeUse.size();
267 indices[2] = (isector / (p.phiRangeUse.size() * p.invptRangeUse.size())) % p.thetaRangeUse.size();
268 indices[3] = isector / (p.phiRangeUse.size() * p.invptRangeUse.size() * p.thetaRangeUse.size());
269 } else {
270 indices[0] = (p.phiRangeUse.size() == 1) ? 0 : isector;
271 indices[1] = (p.invptRangeUse.size() == 1) ? 0 : isector;
272 indices[2] = (p.thetaRangeUse.size() == 1) ? 0 : isector;
273 indices[3] = (p.SLpattern.size() == 1) ? 0 : isector;
274 }
275 return indices;
276}

◆ getRelId()

double getRelId ( const CDCTriggerSegmentHit hit)

Calculate phi position of a hit relative to 2D track (scaled to number of wires).

Definition at line 448 of file NeuroTrigger.cc.

449{
450 int iSL = hit.getISuperLayer();
451 int priority = hit.getPriorityPosition();
452 // (((priority >> 1) & 1) - (priority & 1)) is 0, -1, 1, 0 for priority = 0, 1, 2, 3
453 double relId = hit.getSegmentID() + 0.5 * (((priority >> 1) & 1) - (priority & 1))
454 - m_TSoffset[iSL] - m_idRef[iSL][int(priority < 3)];
455 relId = remainder(relId, (m_TSoffset[iSL + 1] - m_TSoffset[iSL]));
456 return relId;
457}
unsigned m_TSoffset[10]
Number of track segments up to super layer.
Definition: NeuroTrigger.h:297
double m_idRef[9][2]
2D phi position of current track scaled to number of wires
Definition: NeuroTrigger.h:299

◆ initialize()

void initialize ( const NeuroTriggerParameters p)

Definition at line 35 of file NeuroTrigger.cc.

36{
37 // check parameters
38 bool okay = true;
39 // ensure that length of lists matches number of sectors
40 if (p.nHidden.size() != 1 && p.nHidden.size() != p.nMLP) {
41 B2ERROR("Number of nHidden lists should be 1 or " << p.nMLP);
42 okay = false;
43 }
44 if (p.outputScale.size() != 1 && p.outputScale.size() != p.nMLP) {
45 B2ERROR("Number of outputScale lists should be 1 or " << p.nMLP);
46 okay = false;
47 }
48 bool rangeProduct = (p.phiRangeUse.size() * p.invptRangeUse.size() * p.thetaRangeUse.size() * p.SLpattern.size() == p.nMLP);
49 if (!rangeProduct) {
50 if (p.phiRangeUse.size() != 1 && p.phiRangeUse.size() != p.nMLP) {
51 B2ERROR("Number of phiRangeUse.lists should be 1 or " << p.nMLP);
52 okay = false;
53 }
54 if (p.invptRangeUse.size() != 1 && p.invptRangeUse.size() != p.nMLP) {
55 B2ERROR("Number of invptRangeUse.lists should be 1 or " << p.nMLP);
56 okay = false;
57 }
58 if (p.thetaRangeUse.size() != 1 && p.thetaRangeUse.size() != p.nMLP) {
59 B2ERROR("Number of thetaRangeUse.lists should be 1 or " << p.nMLP);
60 okay = false;
61 }
62 if (p.SLpattern.size() != 1 && p.SLpattern.size() != p.nMLP) {
63 B2ERROR("Number of SLpattern lists should be 1 or " << p.nMLP);
64 okay = false;
65 }
66 }
67 // ensure that length of maxHitsPerSL and SLpatternMask lists matches SLpattern
68 if (p.maxHitsPerSL.size() != 1 && p.maxHitsPerSL.size() != p.SLpattern.size()) {
69 B2ERROR("Number of maxHitsPerSL lists should be 1 or " << p.SLpattern.size());
70 okay = false;
71 }
72 if (p.SLpatternMask.size() != 1 && p.SLpatternMask.size() != p.SLpattern.size()) {
73 B2ERROR("Number of SLpatternMask lists should be 1 or " << p.SLpattern.size());
74 okay = false;
75 }
76 // ensure that number of target nodes is valid
77 if (p.targetZ.isSet() || p.targetTheta.isSet()) {
78 unsigned short nTarget = int(p.targetZ) + int(p.targetTheta);
79 if (nTarget < 1) {
80 B2ERROR("No outputs! Turn on either targetZ or targetTheta.");
81 okay = false;
82 }
83 }
84 // ensure that sector ranges are valid
85 for (unsigned iPhi = 0; iPhi < p.phiRangeUse.size(); ++iPhi) {
86 if (p.phiRangeUse[iPhi].size() != 2) {
87 B2ERROR("phiRangeUse should be exactly 2 values");
88 okay = false;
89 continue;
90 }
91 if (p.phiRangeUse[iPhi][0] >= p.phiRangeUse[iPhi][1]) {
92 B2ERROR("phiRangeUse should be smaller than phiRangeUse");
93 okay = false;
94 }
95 if (p.phiRangeUse[iPhi][0] < -360. || p.phiRangeUse[iPhi][1] > 360. ||
96 (p.phiRangeUse[iPhi][1] - p.phiRangeUse[iPhi][0]) > 360.) {
97 B2ERROR("phiRangeUse should be in [-360, 360], with maximal width of 360");
98 okay = false;
99 }
100 }
101 for (unsigned iPt = 0; iPt < p.invptRangeUse.size(); ++iPt) {
102 if (p.invptRangeUse[iPt].size() != 2) {
103 B2ERROR("invptRangeUse should be exactly 2 values");
104 okay = false;
105 }
106 if (p.invptRangeUse[iPt][0] >= p.invptRangeUse[iPt][1]) {
107 B2ERROR("invptRangeUse should be smaller than invptRangeUse");
108 okay = false;
109 }
110 }
111 for (unsigned iTheta = 0; iTheta < p.thetaRangeUse.size(); ++iTheta) {
112 if (p.thetaRangeUse[iTheta].size() != 2) {
113 B2ERROR("thetaRangeUse should be exactly 2 values");
114 okay = false;
115 continue;
116 }
117 if (p.thetaRangeUse[iTheta][0] >= p.thetaRangeUse[iTheta][1]) {
118 B2ERROR("thetaRangeUse should be smaller than thetaRangeUse");
119 okay = false;
120 }
121 if (p.thetaRangeUse[iTheta][0] < 0. || p.thetaRangeUse[iTheta][1] > 180.) {
122 B2ERROR("thetaRangeUse should be in [0, 180]");
123 okay = false;
124 }
125 }
126 int nTarget = (int) p.targetZ + (int) p.targetTheta;
127 if (p.outputScale.size() == p.nMLP || p.outputScale.size() == 1) {
128 for (unsigned iScale = 0; iScale < p.outputScale.size(); ++iScale) {
129 if (p.outputScale[iScale].size() != 2 * static_cast<unsigned int>(nTarget)) {
130 B2ERROR("outputScale should be exactly " << 2 * nTarget << " values");
131 okay = false;
132 }
133 }
134 } else {
135 B2ERROR("the size of outputscale should be 1 or match the number of experts");
136 }
137 // ensure that train sectors are valid
138 if (p.phiRangeUse.size() != p.phiRangeTrain.size()) {
139 B2ERROR("Number of phiRangeUse.lists and phiRangeTrain lists should be equal.");
140 okay = false;
141 } else {
142 for (unsigned iPhi = 0; iPhi < p.phiRangeUse.size(); ++iPhi) {
143 if (p.phiRangeTrain[iPhi].size() != 2) {
144 B2ERROR("phiRangeTrain should be exactly 2 values.");
145 okay = false;
146 } else if (p.phiRangeTrain[iPhi][0] > p.phiRangeUse[iPhi][0] ||
147 p.phiRangeTrain[iPhi][1] < p.phiRangeUse[iPhi][1]) {
148 B2ERROR("phiRangeTrain should be wider than phiRangeUse.or equal.");
149 okay = false;
150 }
151 }
152 }
153 if (p.invptRangeUse.size() != p.invptRangeTrain.size()) {
154 B2ERROR("Number of invptRangeUse.lists and invptRangeTrain lists should be equal.");
155 okay = false;
156 } else {
157 for (unsigned iPt = 0; iPt < p.invptRangeUse.size(); ++iPt) {
158 if (p.invptRangeTrain[iPt].size() != 2) {
159 B2ERROR("invptRangeTrain should be exactly 2 values.");
160 okay = false;
161 } else if (p.invptRangeTrain[iPt][0] > p.invptRangeUse[iPt][0] ||
162 p.invptRangeTrain[iPt][1] < p.invptRangeUse[iPt][1]) {
163 B2ERROR("invptRangeTrain should be wider than invptRangeUse.or equal.");
164 okay = false;
165 }
166 }
167 }
168 if (p.thetaRangeUse.size() != p.thetaRangeTrain.size()) {
169 B2ERROR("Number of thetaRangeUse.lists and thetaRangeTrain lists should be equal.");
170 okay = false;
171 } else {
172 for (unsigned iTheta = 0; iTheta < p.thetaRangeUse.size(); ++iTheta) {
173 if (p.thetaRangeTrain[iTheta].size() != 2) {
174 B2ERROR("thetaRangeTrain should be exactly 2 values.");
175 okay = false;
176 } else if (p.thetaRangeTrain[iTheta][0] > p.thetaRangeUse[iTheta][0] ||
177 p.thetaRangeTrain[iTheta][1] < p.thetaRangeUse[iTheta][1]) {
178 B2ERROR("thetaRangeTrain should be wider than thetaRangeUse.or equal.");
179 okay = false;
180 }
181 }
182 }
183
184 if (!okay) return;
185
186 // initialize MLPs
187 m_MLPs.clear();
188 for (unsigned iMLP = 0; iMLP < p.nMLP; ++iMLP) {
189 //get indices for sector parameters
190 //this is important for cases, where we have experts specialized on different geometrical sectors as well as the pattern mask. since they are all in one array, we need the specific index of the expert. E.g. p.maxhitspersl cloud look like: [<expert-trained-on-slpattern0+thetabigger90>,<expert-trained-on-slpattern1+thetabigger90>,<expert-trained-on-slpattern0+thetasmaller90>,<expert-trained-on-slpattern1+thetasmaller90>]
191 vector<unsigned> indices = getRangeIndices(p, iMLP);
192 //get number of nodes for each layer
193 unsigned short maxHits = (p.maxHitsPerSL.size() == 1) ? p.maxHitsPerSL[0] : p.maxHitsPerSL[indices[3]];
194 unsigned long SLpattern = p.SLpattern[indices[3]];
195 unsigned long SLpatternMask = (p.SLpatternMask.size() == 1) ? p.SLpatternMask[0] : p.SLpatternMask[indices[3]];
196 unsigned short nInput = 27 * maxHits;
197 vector<NNTParam<float>> nHidden = (p.nHidden.size() == 1) ? p.nHidden[0] : p.nHidden[iMLP];
198 vector<unsigned short> nNodes = {nInput};
199 for (unsigned iHid = 0; iHid < nHidden.size(); ++iHid) {
200 if (p.multiplyHidden) {
201 nNodes.push_back(nHidden[iHid] * nNodes[0]);
202 } else {
203 nNodes.push_back(nHidden[iHid]);
204 }
205 }
206 nNodes.push_back(nTarget);
207 unsigned short targetVars = int(p.targetZ) + (int(p.targetTheta) << 1);
208 // the parameters stored in the parameterset are not advanced enough to be vectors, they can only be single data types. the workaround was to make every variable contained in the (nested) vector an NNTParam. for the further use, those have to be converted to float vecors, which is done by the tcastvector function.
209 vector<float> phiRangeUse = p.tcastvector<float>(p.phiRangeUse)[indices[0]];
210 vector<float> invptRangeUse = p.tcastvector<float>(p.invptRangeUse)[indices[1]];
211 vector<float> thetaRangeUse = p.tcastvector<float>(p.thetaRangeUse)[indices[2]];
212 vector<float> phiRangeTrain = p.tcastvector<float>(p.phiRangeTrain)[indices[0]];
213 vector<float> invptRangeTrain = p.tcastvector<float>(p.invptRangeTrain)[indices[1]];
214 vector<float> thetaRangeTrain = p.tcastvector<float>(p.thetaRangeTrain)[indices[2]];
215 B2DEBUG(50, "Ranges for sector " << iMLP
216 << ": phiRange [" << phiRangeUse[0] << ", " << phiRangeUse[1]
217 << "], invptRange [" << invptRangeUse[0] << ", " << invptRangeUse[1]
218 << "], thetaRange [" << thetaRangeUse[0] << ", " << thetaRangeUse[1]
219 << "], SLpattern " << SLpattern);
220 //get scaling values
221 vector<float> outputScale = (p.outputScale.size() == 1) ? p.tcastvector<float>(p.outputScale)[0] : p.tcastvector<float>
222 (p.outputScale)[iMLP];
223 //convert phi and theta from degree to radian
224 phiRangeUse[0] *= Unit::deg;
225 phiRangeUse[1] *= Unit::deg;
226 thetaRangeUse[0] *= Unit::deg;
227 thetaRangeUse[1] *= Unit::deg;
228 phiRangeTrain[0] *= Unit::deg;
229 phiRangeTrain[1] *= Unit::deg;
230 thetaRangeTrain[0] *= Unit::deg;
231 thetaRangeTrain[1] *= Unit::deg;
232 if (p.targetTheta) {
233 outputScale[2 * int(p.targetZ)] *= Unit::deg;
234 outputScale[2 * int(p.targetZ) + 1] *= Unit::deg;
235 }
236 //create new MLP
237 m_MLPs.push_back(CDCTriggerMLP(nNodes, targetVars, outputScale,
238 phiRangeUse, invptRangeUse, thetaRangeUse,
239 phiRangeTrain, invptRangeTrain, thetaRangeTrain,
240 maxHits, SLpattern, SLpatternMask, p.tMax,
241 p.et_option()));
242 }
243
244 if (p.IDRanges.size() == p.nMLP) {
245 for (auto exp : p.IDRanges) {
246 // first entry is the expert number, after that follow the idranges for all the superlayers
247 std::vector<float> irange = {exp.begin() + 1, exp.end()};
248 m_MLPs[static_cast<int>(exp[0])].setRelID(irange);
249 }
250 } else if (p.IDRanges.size() == 0) {
251 B2WARNING("idranges have not been initialized yet, did you forget it?");
252 } else {
253 B2ERROR("number of idranges should match the number of experts!");
254 }
255 // load some values from the geometry that will be needed for the input
256 setConstants();
257}
Class to represent a complete set to describe a Neurotrigger.
std::vector< unsigned > getRangeIndices(const NeuroTriggerParameters &p, unsigned isector)
Get indices for sector ranges in parameter lists.
void setConstants()
Loads parameters from the geometry and precalculates some constants that will be needed.
static const double deg
degree to radians
Definition: Unit.h:109

◆ initializeCollections() [1/2]

void initializeCollections ( std::string  hitCollectionName)

Definition at line 307 of file NeuroTrigger.cc.

308{
309 m_segmentHits.isRequired(hitCollectionName);
310 m_hitCollectionName = hitCollectionName;
311}
bool isRequired(const std::string &name="")
Ensure this array/object has been registered previously.

◆ initializeCollections() [2/2]

void initializeCollections ( std::string  hitCollectionName,
std::string  eventTimeName,
const std::string &  et_option 
)

set the hit collection and event time to required and store the hit collection name

Definition at line 297 of file NeuroTrigger.cc.

298{
299 m_segmentHits.isRequired(hitCollectionName);
300 if (!((et_option == "fastestpriority") || (et_option == "etfhwin") || (et_option == "zero") || (et_option == "fastest2d"))) {
301 m_eventTime.isRequired(eventTimeName);
302 }
303 m_hitCollectionName = hitCollectionName;
304}

◆ load()

bool load ( const std::string &  filename,
const std::string &  arrayname = "MLPs" 
)

Load MLPs from file.

Parameters
filenamename of the TFile to read from
arraynamename of the TObjArray holding the MLPs in the file
Returns
true if the MLPs were loaded correctly

Definition at line 1162 of file NeuroTrigger.cc.

1163{
1164 if (filename.size() < 1) {
1165 m_MLPs.clear();
1166 m_MLPs = m_cdctriggerneuroconfig->getMLPs();
1167 if (m_MLPs.size() == 0) {
1168 B2ERROR("Could not load Neurotrigger weights from database!");
1169 return false;
1170 }
1171 B2INFO("Loaded Neurotrigger MLP weights from database: " + m_cdctriggerneuroconfig->getNNName());
1172 B2DEBUG(100, "loaded " << m_MLPs.size() << " networks from database");
1173 // load some values from the geometry that will be needed for the input
1174 setConstants();
1175 return true;
1176 } else {
1177 TFile datafile(filename.c_str(), "READ");
1178 if (!datafile.IsOpen()) {
1179 B2WARNING("Could not open file " << filename);
1180 return false;
1181 }
1182
1183 TObjArray* MLPs = (TObjArray*)datafile.Get(arrayname.c_str());
1184 if (!MLPs) {
1185 datafile.Close();
1186 B2WARNING("File " << filename << " does not contain key " << arrayname);
1187 return false;
1188 }
1189 m_MLPs.clear();
1190 for (int isector = 0; isector < MLPs->GetEntriesFast(); ++isector) {
1191 CDCTriggerMLP* expert = dynamic_cast<CDCTriggerMLP*>(MLPs->At(isector));
1192 if (expert) m_MLPs.push_back(*expert);
1193 else B2WARNING("Wrong type " << MLPs->At(isector)->ClassName() << ", ignoring this entry.");
1194 }
1195 MLPs->Clear();
1196 delete MLPs;
1197 datafile.Close();
1198 B2DEBUG(100, "loaded " << m_MLPs.size() << " networks");
1199
1200 B2INFO("Loaded Neurotrigger MLP weights from file: " + filename);
1201 // load some values from the geometry that will be needed for the input
1202 setConstants();
1203
1204 return true;
1205 }
1206}
DBObjPtr< CDCTriggerNeuroConfig > m_cdctriggerneuroconfig
get NNT payload from database.
Definition: NeuroTrigger.h:322

◆ loadIDHist()

bool loadIDHist ( const std::string &  filename)

function to load idhist from file

Definition at line 1143 of file NeuroTrigger.cc.

1144{
1145 std::ifstream gzipfile(filename, ios_base::in | ios_base::binary);
1146 boost::iostreams::filtering_istream arrayStream;
1147 arrayStream.push(boost::iostreams::gzip_decompressor());
1148 arrayStream.push(gzipfile);
1150 if (gzipfile.is_open()) {
1151 while (arrayStream >> hline) {
1152 for (unsigned i = 0; i < 18; ++i) {
1153 m_MLPs[hline.exPert].relevantID[i] = hline.relID[i];
1154 }
1155 }
1156 } else { return false;}
1157 return true;
1158}

◆ nSectors()

unsigned nSectors ( ) const
inline

return number of neural networks

Definition at line 164 of file NeuroTrigger.h.

164{ return m_MLPs.size(); }

◆ operator[]() [1/2]

CDCTriggerMLP & operator[] ( unsigned  index)
inline

return reference to a neural network

Definition at line 159 of file NeuroTrigger.h.

159{ return m_MLPs[index]; }

◆ operator[]() [2/2]

const CDCTriggerMLP & operator[] ( unsigned  index) const
inline

return const reference to a neural network

Definition at line 161 of file NeuroTrigger.h.

161{ return m_MLPs[index]; }

◆ runMLP()

vector< float > runMLP ( unsigned  isector,
const std::vector< float > &  input 
)

Run an expert MLP.

Parameters
isectorindex of the MLP
inputvector of input values
Returns
unscaled output values (z vertex in cm and/or theta in radian)

Definition at line 1041 of file NeuroTrigger.cc.

1042{
1043 const CDCTriggerMLP& expert = m_MLPs[isector];
1044 vector<float> weights = expert.getWeights();
1045 vector<float> layerinput = input;
1046 vector<float> layeroutput = {};
1047 unsigned iw = 0;
1048 for (unsigned il = 1; il < expert.nLayers(); ++il) {
1049 //add bias input
1050 layerinput.push_back(1.);
1051 //prepare output
1052 layeroutput.clear();
1053 layeroutput.assign(expert.nNodesLayer(il), 0.);
1054 //loop over outputs
1055 for (unsigned io = 0; io < layeroutput.size(); ++io) {
1056 //loop over inputs
1057 for (unsigned ii = 0; ii < layerinput.size(); ++ii) {
1058 layeroutput[io] += layerinput[ii] * weights[iw++];
1059 }
1060 //apply activation function
1061 layeroutput[io] = tanh(layeroutput[io] / 2.);
1062 }
1063 //output is new input
1064 layerinput = layeroutput;
1065 }
1066 return expert.unscaleTarget(layeroutput);
1067}

◆ runMLPFix()

vector< float > runMLPFix ( unsigned  isector,
std::vector< float >  input 
)

Run an expert MLP with fixed point arithmetic.

Definition at line 1070 of file NeuroTrigger.cc.

1071{
1072 unsigned precisionInput = m_precision[3];
1073 unsigned precisionWeights = m_precision[4];
1074 unsigned precisionLUT = m_precision[5];
1075 unsigned precisionTanh = m_precision[3];
1076 unsigned dp = precisionInput + precisionWeights - precisionLUT;
1077
1078 const CDCTriggerMLP& expert = m_MLPs[isector];
1079 // transform inputs to fixed point (cut off to input precision)
1080 vector<long> inputFix(input.size(), 0);
1081 for (unsigned ii = 0; ii < input.size(); ++ii) {
1082 inputFix[ii] = long(input[ii] * (1 << precisionInput));
1083 }
1084 // transform weights to fixed point (round to weight precision)
1085 vector<float> weights = expert.getWeights();
1086 vector<long> weightsFix(weights.size(), 0);
1087 for (unsigned iw = 0; iw < weights.size(); ++iw) {
1088 weightsFix[iw] = long(round(weights[iw] * (1 << precisionWeights)));
1089 }
1090 // maximum input value for the tanh LUT
1091 unsigned xMax = unsigned(ceil(atanh(1. - 1. / (1 << (precisionTanh + 1))) *
1092 (1 << (precisionLUT + 1))));
1093
1094 // run MLP
1095 vector<long> layerinput = inputFix;
1096 vector<long> layeroutput = {};
1097 unsigned iw = 0;
1098 for (unsigned il = 1; il < expert.nLayers(); ++il) {
1099 // add bias input
1100 layerinput.push_back(1 << precisionInput);
1101 // prepare output
1102 layeroutput.clear();
1103 layeroutput.assign(expert.nNodesLayer(il), 0);
1104 // loop over outputs
1105 for (unsigned io = 0; io < layeroutput.size(); ++io) {
1106 // loop over inputs
1107 for (unsigned ii = 0; ii < layerinput.size(); ++ii) {
1108 layeroutput[io] += layerinput[ii] * weightsFix[iw++];
1109 }
1110 // apply activation function -> LUT, calculated on the fly here
1111 unsigned long bin = abs(layeroutput[io]) >> dp;
1112 // correction to get symmetrical rounding errors
1113 float x = (bin + 0.5 - 1. / (1 << (dp + 1))) / (1 << precisionLUT);
1114 long tanhLUT = (bin < xMax) ? long(round(tanh(x / 2.) * (1 << precisionTanh))) : (1 << precisionTanh);
1115 layeroutput[io] = (layeroutput[io] < 0) ? -tanhLUT : tanhLUT;
1116 }
1117 // output is new input
1118 layerinput = layeroutput;
1119 }
1120
1121 // transform output back to float before unscaling
1122 vector<float> output(layeroutput.size(), 0.);
1123 for (unsigned io = 0; io < output.size(); ++io) {
1124 output[io] = layeroutput[io] / float(1 << precisionTanh);
1125 }
1126 return expert.unscaleTarget(output);
1127}
std::vector< unsigned > m_precision
Fixed point precision in bit after radix point.
Definition: NeuroTrigger.h:313

◆ save()

void save ( const std::string &  filename,
const std::string &  arrayname = "MLPs" 
)

Save MLPs to file.

Parameters
filenamename of the TFile to write to
arraynamename of the TObjArray holding the MLPs in the file

Definition at line 1130 of file NeuroTrigger.cc.

1131{
1132 B2INFO("Saving networks to file " << filename << ", array " << arrayname);
1133 TFile datafile(filename.c_str(), "UPDATE");
1134 TObjArray* MLPs = new TObjArray(m_MLPs.size());
1135 for (unsigned isector = 0; isector < m_MLPs.size(); ++isector) {
1136 MLPs->Add(&m_MLPs[isector]);
1137 }
1138 MLPs->Write(arrayname.c_str(), TObject::kSingleKey | TObject::kOverwrite);
1139 datafile.Close();
1140 MLPs->Clear();
1141 delete MLPs;
1142}

◆ selectHits()

vector< unsigned > selectHits ( unsigned  isector,
const CDCTriggerTrack track,
bool  returnAllRelevant = false 
)

Select best hits for each super layer.

Parameters
isectorindex of the MLP that will use the input
trackaxial hit relations are taken from given track
returnAllRelevantif true, return all relevant hits instead of selecting the best (for making relations)
Returns
list of selected hit indices

Definition at line 869 of file NeuroTrigger.cc.

871{
872 const CDCTriggerMLP& expert = m_MLPs[isector];
873 vector<unsigned> selectedHitIds = {};
874 // prepare vectors to keep best drift times, left/right and selected hit IDs
875 vector<int> tMin;
876 tMin.assign(expert.nNodesLayer(0), expert.getTMax());
877 vector<bool> LRknown;
878 LRknown.assign(expert.nNodesLayer(0), false);
879 vector<int> hitIds;
880 hitIds.assign(expert.nNodesLayer(0), -1);
881 vector<unsigned> nHits;
882 nHits.assign(9, 0);
883
884 // loop over axial hits related to input track
886 track.getRelationsTo<CDCTriggerSegmentHit>(m_hitCollectionName);
887 B2DEBUG(250, "start axial hit loop");
888 for (unsigned ihit = 0; ihit < axialHits.size(); ++ihit) {
889 // skip hits with negative relation weight (not selected in finder)
890 if (axialHits.weight(ihit) < 0) continue;
891 unsigned short iSL = axialHits[ihit]->getISuperLayer();
892 // skip stereo hits (should not be related to track, but check anyway)
893 if (iSL % 2 == 1) continue;
894 if (expert.getSLpatternUnmasked() != 0 &&
895 !((expert.getSLpatternUnmasked() >> iSL) & 1)) {
896 B2DEBUG(250, "skipping hit in SL " << iSL);
897 continue;
898 }
899 // get priority time and apply time window cut
900
901 int t = (m_hasT0) ? axialHits[ihit]->priorityTime() - m_T0 : 0;
902 if (t < 0) {
903 t = 0;
904 } else if (t > expert.getTMax()) {
905 t = expert.getTMax();
906 }
907 double relId = getRelId(*axialHits[ihit]);
908 if (expert.isRelevant(relId, iSL)) {
909 // get reference hit (worst of existing hits)
910 unsigned short iRef = iSL;
911 if (expert.getMaxHitsPerSL() > 1) {
912 if (nHits[iSL] < expert.getMaxHitsPerSL() &&
913 (expert.getSLpatternUnmasked() >> (iSL + 9 * nHits[iSL])) & 1) {
914 iRef += 9 * nHits[iSL];
915 ++nHits[iSL];
916 } else {
917 for (unsigned compare = iSL; compare < iSL + 9 * nHits[iSL]; compare += 9) {
918 if ((LRknown[iRef] && !LRknown[compare]) ||
919 (LRknown[iRef] == LRknown[compare] && tMin[iRef] < tMin[compare]))
920 iRef = compare;
921 }
922 }
923 }
924 // choose best hit (LR known before LR unknown, then shortest drift time)
925 bool useHit = false;
926 if (LRknown[iRef]) {
927 useHit = (axialHits[ihit]->LRknown() && t <= tMin[iRef]);
928 } else {
929 useHit = (axialHits[ihit]->LRknown() || t <= tMin[iRef]);
930 }
931 B2DEBUG(250, "relevant wire SL " << iSL << " LR " << axialHits[ihit]->getLeftRight()
932 << " t " << t << " iRef " << iRef << " useHit " << useHit);
933 if (useHit) {
934 // keep drift time and LR
935 LRknown[iRef] = axialHits[ihit]->LRknown();
936 tMin[iRef] = t;
937 hitIds[iRef] = axialHits[ihit]->getArrayIndex();
938 }
939 if (returnAllRelevant) selectedHitIds.push_back(axialHits[ihit]->getArrayIndex());
940 }
941 }
942
943 // loop over stereo hits, choosing up to MaxHitsPerSL per superlayer
944 B2DEBUG(250, "start stereo hit loop");
945 for (int ihit = 0; ihit < m_segmentHits.getEntries(); ++ ihit) {
946 unsigned short iSL = m_segmentHits[ihit]->getISuperLayer();
947 // skip axial hits
948 if (iSL % 2 == 0) continue;
949 if (expert.getSLpatternUnmasked() != 0 &&
950 !((expert.getSLpatternUnmasked() >> iSL) & 1)) {
951 B2DEBUG(250, "skipping hit in SL " << iSL);
952 continue;
953 }
954 // get priority time and apply time window cut
955 int t = (m_hasT0) ? m_segmentHits[ihit]->priorityTime() - m_T0 : 0;
956 if (t < 0) {
957 t = 0;
958 } else if (t > expert.getTMax()) {
959 t = expert.getTMax();
960 }
961 double relId = getRelId(*m_segmentHits[ihit]);
962 if (expert.isRelevant(relId, iSL)) {
963 // get reference hit (worst of existing hits)
964 unsigned short iRef = iSL;
965 if (expert.getMaxHitsPerSL() > 1) {
966 if (nHits[iSL] < expert.getMaxHitsPerSL() &&
967 (expert.getSLpatternUnmasked() >> (iSL + 9 * nHits[iSL])) & 1) {
968 iRef += 9 * nHits[iSL];
969 ++nHits[iSL];
970 } else {
971 for (unsigned compare = iSL; compare < iSL + 9 * nHits[iSL]; compare += 9) {
972 if ((LRknown[iRef] && !LRknown[compare]) ||
973 (LRknown[iRef] == LRknown[compare] && tMin[iRef] < tMin[compare]))
974 iRef = compare;
975 }
976 }
977 }
978 // choose best hit (LR known before LR unknown, then shortest drift time)
979 bool useHit = false;
980 if (LRknown[iRef]) {
981 useHit = (m_segmentHits[ihit]->LRknown() && t <= tMin[iRef]);
982 } else {
983 useHit = (m_segmentHits[ihit]->LRknown() || t <= tMin[iRef]);
984 }
985 B2DEBUG(250, "relevant wire SL " << iSL << " LR " << m_segmentHits[ihit]->getLeftRight()
986 << " t " << t << " iRef " << iRef << " useHit " << useHit);
987 if (useHit) {
988 // keep drift time and LR
989 LRknown[iRef] = m_segmentHits[ihit]->LRknown();
990 tMin[iRef] = t;
991 hitIds[iRef] = ihit;
992 }
993 if (returnAllRelevant) selectedHitIds.push_back(ihit);
994 }
995 }
996
997 // save selected hit Ids
998 if (!returnAllRelevant) {
999 for (unsigned iHit = 0; iHit < hitIds.size(); ++iHit) {
1000 if (hitIds[iHit] >= 0) selectedHitIds.push_back(hitIds[iHit]);
1001 }
1002 }
1003 return selectedHitIds;
1004}

◆ selectHitsHWSim()

vector< unsigned > selectHitsHWSim ( unsigned  isector,
const CDCTriggerTrack track 
)

Select hits for each super layer from the ones related to input track.

Parameters
isectorindex of the MLP that will use the input
trackall hit relations are taken from given track
Returns
list of selected hit indices

Definition at line 790 of file NeuroTrigger.cc.

791{
792 const CDCTriggerMLP& expert = m_MLPs[isector];
793 vector<unsigned> selectedHitIds = {};
794 // prepare vectors to keep best drift times, left/right and selected hit IDs
795 vector<int> tMin;
796 tMin.assign(expert.nNodesLayer(0), expert.getTMax());
797 vector<bool> LRknown;
798 LRknown.assign(expert.nNodesLayer(0), false);
799 vector<int> hitIds;
800 hitIds.assign(expert.nNodesLayer(0), -1);
801 vector<unsigned> nHits;
802 nHits.assign(9, 0);
803
804 // loop over all hits related to input track
806 track.getRelationsTo<CDCTriggerSegmentHit>(m_hitCollectionName);
807 B2DEBUG(250, "start hit loop over all related hits");
812 for (unsigned ihit = 0; ihit < allHits.size(); ++ihit) {
813 // skip hits with negative relation weight (not selected in finder)
814 if (allHits.weight(ihit) < 0) continue;
815 unsigned short iSL = allHits[ihit]->getISuperLayer();
816 if (expert.getSLpatternUnmasked() != 0 &&
817 !((expert.getSLpatternUnmasked() >> iSL) & 1)) {
818 B2DEBUG(250, "skipping hit in SL " << iSL);
819 continue;
820 }
821 // get priority time and apply time window cut
822 int t = (m_hasT0) ? allHits[ihit]->priorityTime() - m_T0 : 0;
823 if (t < 0) {
824 t = 0;
825 } else if (t > expert.getTMax()) {
826 t = expert.getTMax();
827 }
828 double relId = getRelId(*allHits[ihit]);
829 if (expert.isRelevant(relId, iSL)) {
830 // get reference hit (worst of existing hits)
831 unsigned short iRef = iSL;
832 if (expert.getMaxHitsPerSL() > 1) {
833 if (nHits[iSL] < expert.getMaxHitsPerSL() &&
834 (expert.getSLpatternUnmasked() >> (iSL + 9 * nHits[iSL])) & 1) {
835 iRef += 9 * nHits[iSL];
836 ++nHits[iSL];
837 } else {
838 for (unsigned compare = iSL; compare < iSL + 9 * nHits[iSL]; compare += 9) {
839 if ((LRknown[iRef] && !LRknown[compare]) ||
840 (LRknown[iRef] == LRknown[compare] && tMin[iRef] < tMin[compare]))
841 iRef = compare;
842 }
843 }
844 }
845 // choose best hit (LR known before LR unknown, then shortest drift time)
846 bool useHit = false;
847 if (LRknown[iRef]) {
848 useHit = (allHits[ihit]->LRknown() && t <= tMin[iRef]);
849 } else {
850 useHit = (allHits[ihit]->LRknown() || t <= tMin[iRef]);
851 }
852 B2DEBUG(250, "relevant wire SL " << iSL << " LR " << allHits[ihit]->getLeftRight()
853 << " t " << t << " iRef " << iRef << " useHit " << useHit);
854 if (useHit) {
855 // keep drift time and LR
856 LRknown[iRef] = allHits[ihit]->LRknown();
857 tMin[iRef] = t;
858 hitIds[iRef] = allHits[ihit]->getArrayIndex();
859 }
860 }
861 }
862 // save selected hit Ids
863 for (unsigned iHit = 0; iHit < hitIds.size(); ++iHit) {
864 if (hitIds[iHit] >= 0) selectedHitIds.push_back(hitIds[iHit]);
865 }
866 return selectedHitIds;
867}

◆ selectMLPbyPattern()

int selectMLPbyPattern ( std::vector< int > &  MLPs,
unsigned long  pattern,
const bool  neurotrackinputmode 
)

Select one MLP from a list of sector indices.

The selected expert either matches the given sector pattern, or has no pattern restriction. An unrestricted expert is returned only if there is no exactly matching expert.

Returns
index of the selected MLP, -1 if no matching MLP is found

Definition at line 365 of file NeuroTrigger.cc.

366{
367 if (MLPs.size() == 0) {
368 return -1;
369 }
370
371 int bestIndex = -1;
372 for (unsigned i = 0; i < MLPs.size(); ++i) {
373 int isector = MLPs[i];
374 unsigned long sectorPattern = m_MLPs[isector].getSLpattern();
375 B2DEBUG(250, "hitPattern " << pattern << " sectorPattern " << sectorPattern);
376 // no hit pattern restriction -> keep looking for exact match
377 if (sectorPattern == 0) {
378 B2DEBUG(250, "found match for general sector");
379 bestIndex = isector;
380 }
381 // exact match -> keep this sector and stop searching
382 if (pattern == sectorPattern) {
383 B2DEBUG(250, "found match for hit pattern " << pattern);
384 bestIndex = isector;
385 break;
386 }
387 }
388
389 if (bestIndex < 0) {
390 if (neurotrackinputmode) {
391 B2DEBUG(150, "No sector found to match pattern, using sector 0" << pattern << ".");
392 bestIndex = 0;
393 } else {
394 B2DEBUG(150, "No sector found to match pattern " << pattern << ".");
395 }
396 }
397 return bestIndex;
398}

◆ selectMLPs()

vector< int > selectMLPs ( float  phi0,
float  invpt,
float  theta 
)

Select all matching expert MLPs based on the given track parameters.

If the sectors are overlapping, there may be more than one matching expert. During training this is intended, afterwards sectors should be redefined to be unique. For unique geometrical sectors, this function can still find several experts with different sector patterns.

Returns
indices of the selected MLPs, empty if the track does not fit any sector

Definition at line 338 of file NeuroTrigger.cc.

339{
340 vector<int> indices = {};
341
342 if (m_MLPs.size() == 0) {
343 B2WARNING("Trying to select MLP before initializing MLPs.");
344 return indices;
345 }
346
347 // find all matching sectors
348 for (unsigned isector = 0; isector < m_MLPs.size(); ++isector) {
349 if (m_MLPs[isector].inPhiRangeUse(phi0) && m_MLPs[isector].inInvptRangeUse(invpt)
350 && m_MLPs[isector].inThetaRangeUse(theta)) {
351 indices.push_back(isector);
352 }
353 }
354
355 if (indices.size() == 0) {
356 B2DEBUG(150, "Track does not match any sector.");
357 B2DEBUG(150, "invpt=" << invpt << ", phi=" << phi0 * 180. / M_PI << ", theta=" << theta * 180. / M_PI);
358 }
359
360 return indices;
361}

◆ selectMLPsTrain()

vector< int > selectMLPsTrain ( float  phi0,
float  invpt,
float  theta 
)

Select all matching expert MLPs based on the given track parameters.

If the sectors are overlapping, there may be more than one matching expert. During training this is intended, afterwards sectors should be redefined to be unique. For unique geometrical sectors, this function can still find several experts with different sector patterns.

Returns
indices of the selected MLPs, empty if the track does not fit any sector

Definition at line 313 of file NeuroTrigger.cc.

314{
315 vector<int> indices = {};
316
317 if (m_MLPs.size() == 0) {
318 B2WARNING("Trying to select MLP before initializing MLPs.");
319 return indices;
320 }
321
322 // find all matching sectors
323 for (unsigned isector = 0; isector < m_MLPs.size(); ++isector) {
324 if (m_MLPs[isector].inPhiRangeTrain(phi0) && m_MLPs[isector].inInvptRangeTrain(invpt)
325 && m_MLPs[isector].inThetaRangeTrain(theta)) {
326 indices.push_back(isector);
327 }
328 }
329
330 if (indices.size() == 0) {
331 B2DEBUG(150, "Track does not match any sector.");
332 B2DEBUG(150, "invpt=" << invpt << ", phi=" << phi0 * 180. / M_PI << ", theta=" << theta * 180. / M_PI);
333 }
334
335 return indices;
336}

◆ setConstants()

void setConstants ( )

Loads parameters from the geometry and precalculates some constants that will be needed.

Definition at line 280 of file NeuroTrigger.cc.

281{
283 int layerId = 3;
284 int nTS = 0;
285 for (int iSL = 0; iSL < 9; ++iSL) {
286 m_TSoffset[iSL] = nTS;
287 nTS += cdc.nWiresInLayer(layerId);
288 m_TSoffset[iSL + 1] = nTS;
289 for (int priority = 0; priority < 2; ++ priority) {
290 m_radius[iSL][priority] = cdc.senseWireR(layerId + priority);
291 }
292 layerId += (iSL > 0 ? 6 : 7);
293 }
294}
The Class for CDC Geometry Parameters.
static CDCGeometryPar & Instance(const CDCGeometry *=nullptr)
Static method to get a reference to the CDCGeometryPar instance.
double m_radius[9][2]
Radius of the CDC layers with priority wires (2 per super layer)
Definition: NeuroTrigger.h:295

◆ setPrecision()

void setPrecision ( const std::vector< unsigned > &  precision)
inline

set fixed point precision

Definition at line 151 of file NeuroTrigger.h.

151{ m_precision = precision; }

◆ updateTrack()

void updateTrack ( const CDCTriggerTrack track)

Calculate 2D phi position and arclength for the given track and store them.

Definition at line 402 of file NeuroTrigger.cc.

403{
404 double omega = track.getOmega(); // signed track curvature
405 for (int iSL = 0; iSL < 9; ++iSL) {
406 for (int priority = 0; priority < 2; ++priority) {
407 m_alpha[iSL][priority] = (m_radius[iSL][priority] * abs(omega) < 2.) ?
408 asin(m_radius[iSL][priority] * omega / 2.) :
409 M_PI_2;
410 m_idRef[iSL][priority] = remainder(((track.getPhi0() - m_alpha[iSL][priority]) *
411 (m_TSoffset[iSL + 1] - m_TSoffset[iSL]) / 2. / M_PI),
412 (m_TSoffset[iSL + 1] - m_TSoffset[iSL]));
413 }
414 }
415}

◆ updateTrackFix()

void updateTrackFix ( const CDCTriggerTrack track)

Calculate 2D phi position and arclength for the given track and store them.

Definition at line 418 of file NeuroTrigger.cc.

419{
420 //unsigned precisionPhi = m_precision[0];
421 //unsigned precisionAlpha = m_precision[0];
422 unsigned precisionPhiAlpha = m_precision[0];
423 unsigned precisionScale = m_precision[1];
424 unsigned precisionId = m_precision[2];
425
426 double omega = track.getOmega();
427 long phi = round(track.getPhi0() * (1 << precisionPhiAlpha));
428
429 for (int iSL = 0; iSL < 9; ++iSL) {
430 for (int priority = 0; priority < 2; ++priority) {
431 // LUT, calculated on the fly here
432 m_alpha[iSL][priority] =
433 (m_radius[iSL][priority] * abs(omega) < 2) ?
434 round(asin(m_radius[iSL][priority] * omega / 2) * (1 << precisionPhiAlpha)) :
435 round(M_PI_2 * (1 << precisionPhiAlpha));
436 long dphi = phi - (long(m_alpha[iSL][priority]));
437 m_idRef[iSL][priority] =
438 long(dphi * round((m_TSoffset[iSL + 1] - m_TSoffset[iSL]) / 2. / M_PI * (1 << precisionScale)) /
439 (long(1) << (precisionPhiAlpha + precisionScale - precisionId)));
440 // unscale after rounding
441 m_alpha[iSL][priority] /= (1 << precisionPhiAlpha);
442 m_idRef[iSL][priority] /= (1 << precisionId);
443 }
444 }
445}

Member Data Documentation

◆ m_alpha

double m_alpha[9][2] = {}
private

2D crossing angle of current track

Definition at line 301 of file NeuroTrigger.h.

◆ m_cdctriggerneuroconfig

DBObjPtr<CDCTriggerNeuroConfig> m_cdctriggerneuroconfig
private

get NNT payload from database.

Definition at line 322 of file NeuroTrigger.h.

◆ m_eventTime

StoreObjPtr<BinnedEventT0> m_eventTime
private

StoreObjPtr containing the event time.

Definition at line 318 of file NeuroTrigger.h.

◆ m_hasT0

bool m_hasT0 = false
private

Flag to show if stored event time is valid.

Definition at line 305 of file NeuroTrigger.h.

◆ m_hitCollectionName

std::string m_hitCollectionName
private

Name of the StoreArray containing the input track segment hits.

Definition at line 320 of file NeuroTrigger.h.

◆ m_idRef

double m_idRef[9][2] = {}
private

2D phi position of current track scaled to number of wires

Definition at line 299 of file NeuroTrigger.h.

◆ m_MLPs

std::vector<CDCTriggerMLP> m_MLPs = {}
private

List of networks.

Definition at line 293 of file NeuroTrigger.h.

◆ m_precision

std::vector<unsigned> m_precision
private

Fixed point precision in bit after radix point.

8 values:

  • 2D track parameters: omega, phi
  • geometrical values derived from track: crossing angle, reference wire ID
  • scale factor: radian to wire ID
  • MLP values: nodes, weights, activation function LUT input (LUT output = nodes)

Definition at line 313 of file NeuroTrigger.h.

◆ m_radius

double m_radius[9][2] = {}
private

Radius of the CDC layers with priority wires (2 per super layer)

Definition at line 295 of file NeuroTrigger.h.

◆ m_segmentHits

StoreArray<CDCTriggerSegmentHit> m_segmentHits
private

StoreArray containing the input track segment hits.

Definition at line 316 of file NeuroTrigger.h.

◆ m_T0

int m_T0 = 0
private

Event time of current event / track.

Definition at line 303 of file NeuroTrigger.h.

◆ m_TSoffset

unsigned m_TSoffset[10] = {}
private

Number of track segments up to super layer.

Definition at line 297 of file NeuroTrigger.h.


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