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)
CDCTriggerMLP &	operator[] (unsigned index)
	return reference to a neural network
const CDCTriggerMLP &	operator[] (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< CDCTriggerMLP >	m_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< CDCTriggerSegmentHit >	m_segmentHits
	StoreArray containing the input track segment hits.
StoreObjPtr< BinnedEventT0 >	m_eventTime
	StoreObjPtr containing the event time.
std::string	m_hitCollectionName
	Name of the StoreArray containing the input track segment hits.
DBObjPtr< CDCTriggerNeuroConfig >	m_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.

{}

~NeuroTrigger()

virtual ~NeuroTrigger ( )

inlinevirtual

Default destructor.

Definition at line 121 of file NeuroTrigger.h.

{}

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.

{ m_MLPs.push_back(newMLP); }

get_et_option()

std::string get_et_option ( )

inline

Return value of m_et_option.

Definition at line 229 of file NeuroTrigger.h.

    {
      std::string eto = m_MLPs[0].get_et_option();
      for (unsigned int i = 0; i < m_MLPs.size(); ++i) {
        if (m_MLPs[i].get_et_option() != eto) {
          B2ERROR("Timing options in the expert networks in the CDC Neurotrigger differ!");
        }
      }
      return eto;
 
    }

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.

{
  const CDCTriggerMLP& expert = m_MLPs[isector];
  unsigned long chitPattern = 0;
  vector<unsigned> nHits;
  nHits.assign(9, 0);
  // loop over axial hits related to input track
  RelationVector<CDCTriggerSegmentHit> axialHits =
    track.getRelationsTo<CDCTriggerSegmentHit>(m_hitCollectionName);
  for (unsigned ihit = 0; ihit < axialHits.size(); ++ ihit) {
    // skip hits with negative relation weight (not selected in finder)
    if (axialHits.weight(ihit) < 0) {
      continue;
    }
    unsigned short iSL = axialHits[ihit]->getISuperLayer();
    // // skip stereo hits (should not be related to track, but check anyway)
    if ((!neurotrackinputmode) && (iSL % 2 == 1)) continue;
    double relId = getRelId(*axialHits[ihit]);
    if (expert.isRelevant(relId, iSL)) {
      if (nHits[iSL] < expert.getMaxHitsPerSL()) {
        chitPattern |= 1 << (iSL + 9 * nHits[iSL]);
        ++nHits[iSL];
      }
    }
  }
  if (!neurotrackinputmode) {
    // loop over stereo hits
    for (int ihit = 0; ihit < m_segmentHits.getEntries(); ++ ihit) {
      unsigned short iSL = m_segmentHits[ihit]->getISuperLayer();
      // skip axial hits
      if (iSL % 2 == 0) continue;
      // get priority time
      double relId = getRelId(*m_segmentHits[ihit]);
      if (expert.isRelevant(relId, iSL)) {
        if (nHits[iSL] < expert.getMaxHitsPerSL()) {
          chitPattern |= 1 << (iSL + 9 * nHits[iSL]);
          ++nHits[iSL];
        }
      }
    }
  }
  return chitPattern;
}

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.

{
 
  if (et_option != m_MLPs[isector].get_et_option()) {
    B2DEBUG(20, "Used event time option is different to the one set in the MLP"
            << LogVar("et_option", et_option) << LogVar("isector", isector)
            << LogVar("et_option_mlp", m_MLPs[isector].get_et_option()));
  }
  if (et_option == "fastestpriority") {
    B2DEBUG(200, "et_option is 'fastestpriority'");
    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;
    }
 
  } else if (et_option == "fastest2d") {
    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;
    }
  } else if (et_option == "zero") {
    m_hasT0 = true;
    m_T0 = 0;
  } else if (et_option == "etf_only") {
    bool hasT0 = (m_eventTime.isValid()) ? m_eventTime->hasBinnedEventT0(Const::CDC) : false;
    if (hasT0) {
      m_T0 = m_eventTime->getBinnedEventT0(Const::CDC);
      m_hasT0 = true;
    } else {
      B2ERROR("No ETF output, but forced to use ETF!");
      m_hasT0 = false;
      m_T0 = 0;
    }
  } else if (et_option == "etf_or_fastestpriority") {
    bool hasT0 = (m_eventTime.isValid()) ? m_eventTime->hasBinnedEventT0(Const::CDC) : false;
    if (hasT0) {
      m_T0 = m_eventTime->getBinnedEventT0(Const::CDC);
      m_hasT0 = true;
    } else {
      getEventTime(isector, track, "fastestpriority", neuroinputmode);
    }
  } else if (et_option == "min_etf_fastestpriority") {
    bool hasT0 = (m_eventTime.isValid()) ? m_eventTime->hasBinnedEventT0(Const::CDC) : false;
    int T0_etf = 9999;
    if (hasT0) {
      T0_etf = m_eventTime->getBinnedEventT0(Const::CDC);
      m_hasT0 = true;
    }
    getEventTime(isector, track, "fastestpriority", neuroinputmode);
    if (m_T0 > T0_etf) {
      m_T0 = T0_etf;
    }
  } else if (et_option == "etf_or_fastest2d") {
    bool hasT0 = (m_eventTime.isValid()) ? m_eventTime->hasBinnedEventT0(Const::CDC) : false;
    if (hasT0) {
      m_T0 = m_eventTime->getBinnedEventT0(Const::CDC);
      m_hasT0 = true;
    } else {
      getEventTime(isector, track, "fastest2d", neuroinputmode);
    }
  } else if (et_option == "etf_or_zero") {
    bool hasT0 = (m_eventTime.isValid()) ? m_eventTime->hasBinnedEventT0(Const::CDC) : false;
    if (hasT0) {
      m_T0 = m_eventTime->getBinnedEventT0(Const::CDC);
      m_hasT0 = true;
    } else {
      m_hasT0 = true;
      m_T0 = 0;
    }
  } else if (et_option == "etfcc") {
    if (track.getHasETFTime()) {
      m_T0 = track.getETF_unpacked();
      m_hasT0 = true;
    } else {
      m_T0 = 0;
      m_hasT0 = false;
    }
 
  } else if (et_option == "etfcc_or_zero") {
    if (track.getHasETFTime()) {
      m_T0 = track.getETF_unpacked();
      m_hasT0 = true;
    } else {
      m_T0 = 0;
      m_hasT0 = true;
    }
 
  } else if (et_option == "etfcc_or_fastestpriority") {
    if (track.getHasETFTime()) {
      m_T0 = track.getETF_unpacked();
      m_hasT0 = true;
    } else {
      getEventTime(isector, track, "fastestpriority", neuroinputmode);
    }
 
  } else if (et_option == "etfhwin") {
    m_T0 = track.getETF_recalced();
    m_hasT0 = true;
 
  } else {
    B2ERROR("No valid parameter for et_option (" << et_option << " )!");
  }
 
}

getInputPattern()

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

Calculate input pattern for MLP.

Parameters

isector	index of the MLP that will use the input
track	axial hit relations are taken from given track
neurotrackinputmode	input mode

Returns

super layer pattern of hits in the current track

Definition at line 724 of file NeuroTrigger.cc.

{
  const CDCTriggerMLP& expert = m_MLPs[isector];
  unsigned long hitPattern = 0;
  vector<unsigned> nHits;
  nHits.assign(9, 0);
  // loop over axial hits related to input track
  RelationVector<CDCTriggerSegmentHit> axialHits =
    track.getRelationsTo<CDCTriggerSegmentHit>(m_hitCollectionName);
  for (unsigned ihit = 0; ihit < axialHits.size(); ++ ihit) {
    // skip hits with negative relation weight (not selected in finder)
    if (axialHits.weight(ihit) < 0) {
      continue;
    }
    unsigned short iSL = axialHits[ihit]->getISuperLayer();
    // // skip stereo hits (should not be related to track, but check anyway)
    if ((!neurotrackinputmode) && (iSL % 2 == 1)) continue;
    double relId = getRelId(*axialHits[ihit]);
 
    if (expert.isRelevant(relId, iSL)) {
      if (nHits[iSL] < expert.getMaxHitsPerSL()) {
        hitPattern |= 1 << (iSL + 9 * nHits[iSL]);
        ++nHits[iSL];
      }
      B2DEBUG(250, "hit in SL " << iSL);
    } else {
      B2DEBUG(250, "hit in SL " << iSL << " not relevant (relId = " << relId << ")");
    }
  }
  if (!neurotrackinputmode) {
    // loop over stereo hits
    for (int ihit = 0; ihit < m_segmentHits.getEntries(); ++ ihit) {
      unsigned short iSL = m_segmentHits[ihit]->getISuperLayer();
      // skip axial hits
      if (iSL % 2 == 0) continue;
      double relId = getRelId(*m_segmentHits[ihit]);
      if (expert.isRelevant(relId, iSL)) {
        if (nHits[iSL] < expert.getMaxHitsPerSL()) {
          hitPattern |= 1 << (iSL + 9 * nHits[iSL]);
          ++nHits[iSL];
        }
        B2DEBUG(250, "hit in SL " << iSL);
      } else {
        B2DEBUG(250, "hit in SL " << iSL << " not relevant (relId = " << relId << ")");
      }
    }
  }
  B2DEBUG(250, "hitPattern " << hitPattern);
  return hitPattern & expert.getSLpatternMask();
}

getInputVector()

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

Calculate input values for MLP.

Parameters

isector	index of the MLP that will use the input
hitIds	hit indices to be used for the input

Returns

scaled vector of input values (1 for each input node)

Definition at line 993 of file NeuroTrigger.cc.

{
  const CDCTriggerMLP& expert = m_MLPs[isector];
  // prepare empty input vector and vectors to keep best drift times
  vector<float> inputVector;
  inputVector.assign(expert.nNodesLayer(0), 0.);
  // convert hits to input values
  vector<unsigned> nHits;
  nHits.assign(9, 0);
  for (unsigned ii = 0; ii < hitIds.size(); ++ii) {
    int ihit = hitIds[ii];
    unsigned short iSL = m_segmentHits[ihit]->getISuperLayer();
    unsigned short iRef = iSL + 9 * nHits[iSL];
    ++nHits[iSL];
    int priot = m_segmentHits[ihit]->priorityTime();
    int t = (m_hasT0) ? priot - m_T0 : 0;
    if (t < 0) {
      t = 0;
    } else if (t > expert.getTMax()) {
      t = expert.getTMax();
    }
    int LR = m_segmentHits[ihit]->getLeftRight();
    double relId = getRelId(*m_segmentHits[ihit]);
    int priority = m_segmentHits[ihit]->getPriorityPosition();
    // get scaled input values (scaling such that the absolute value of all inputs is < 1)
    inputVector[3 * iRef] = expert.scaleId(relId, iSL);
    float scaleT = pow(2, floor(log2(1. / expert.getTMax())));
    inputVector[3 * iRef + 1] = (((LR >> 1) & 1) - (LR & 1)) * t * scaleT;
    inputVector[3 * iRef + 2] = m_alpha[iSL][int(priority < 3)] * 0.5;
  }
  return inputVector;
}

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.

{
  int tlow = 9999;
  B2DEBUG(200, "looping over axials:");
  for (unsigned ihit = 0; ihit < Hits.size(); ++ihit) {
    // skip hits with negative relation weight (not selected in finder)
    if (Hits.weight(ihit) < 0) continue;
    unsigned short iSL = Hits[ihit]->getISuperLayer();
    if (iSL % 2 == 1 && onlyAxials) {continue;}
    // get shortest time of relevant hits
    B2DEBUG(200, "  check drifttime: SL" + std::to_string(iSL) + ",ID = " + std::to_string(Hits[ihit]->getSegmentID()) + ", t = " +
            std::to_string(Hits[ihit]->priorityTime()));
    double relId = getRelId(*Hits[ihit]);
    if (m_MLPs[isector].isRelevant(relId, iSL) &&
        Hits[ihit]->priorityTime() < tlow) {
      tlow = Hits[ihit]->priorityTime();
      B2DEBUG(200, "    new tlow: " << std::to_string(tlow));
    }
  }
  return tlow;
}

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.

{
  const CDCTriggerMLP& expert = m_MLPs[isector];
  unsigned long driftth = 0;
  vector<unsigned> nHits;
  nHits.assign(9, 0);
  // loop over axial hits related to input track
  RelationVector<CDCTriggerSegmentHit> axialHits =
    track.getRelationsTo<CDCTriggerSegmentHit>(m_hitCollectionName);
  for (unsigned ihit = 0; ihit < axialHits.size(); ++ ihit) {
    // skip hits with negative relation weight (not selected in finder)
    if (axialHits.weight(ihit) < 0) {
      continue;
    }
    unsigned short iSL = axialHits[ihit]->getISuperLayer();
    // // skip stereo hits (should not be related to track, but check anyway)
    if ((!neurotrackinputmode) && (iSL % 2 == 1)) continue;
    // get priority time
    int t = (m_hasT0) ? axialHits[ihit]->priorityTime() - m_T0 : 0;
    double relId = getRelId(*axialHits[ihit]);
    if (t < 0 || t > expert.getTMax()) {
      if (expert.isRelevant(relId, iSL)) {
        if (nHits[iSL] < expert.getMaxHitsPerSL()) {
          driftth |= 1 << (8 - iSL + 9 * nHits[iSL]);
          ++nHits[iSL];
        }
      }
    }
  }
  if (!neurotrackinputmode) {
    // loop over stereo hits
    for (int ihit = 0; ihit < m_segmentHits.getEntries(); ++ ihit) {
      unsigned short iSL = m_segmentHits[ihit]->getISuperLayer();
      // skip axial hits
      if (iSL % 2 == 0) continue;
      // get priority time
      int t = (m_hasT0) ? m_segmentHits[ihit]->priorityTime() - m_T0 : 0;
      double relId = getRelId(*m_segmentHits[ihit]);
      if (t < 0 || t > expert.getTMax()) {
        if (expert.isRelevant(relId, iSL)) {
          if (nHits[iSL] < expert.getMaxHitsPerSL()) {
            driftth |= 1 << (8 - iSL + 9 * nHits[iSL]);
            ++nHits[iSL];
          }
        }
      }
    }
  }
  return driftth;
}

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.

{
  // 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).
  std::vector<unsigned> indices = {0, 0, 0, 0};
  if (p.phiRangeUse.size() * p.invptRangeUse.size() * p.thetaRangeUse.size() * p.SLpattern.size() == p.nMLP) {
    indices[0] = isector % p.phiRangeUse.size();
    indices[1] = (isector / p.phiRangeUse.size()) % p.invptRangeUse.size();
    indices[2] = (isector / (p.phiRangeUse.size() * p.invptRangeUse.size())) % p.thetaRangeUse.size();
    indices[3] = isector / (p.phiRangeUse.size() * p.invptRangeUse.size() * p.thetaRangeUse.size());
  } else {
    indices[0] = (p.phiRangeUse.size() == 1) ? 0 : isector;
    indices[1] = (p.invptRangeUse.size() == 1) ? 0 : isector;
    indices[2] = (p.thetaRangeUse.size() == 1) ? 0 : isector;
    indices[3] = (p.SLpattern.size() == 1) ? 0 : isector;
  }
  return indices;
}

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.

{
  int iSL = hit.getISuperLayer();
  int priority = hit.getPriorityPosition();
  // (((priority >> 1) & 1) - (priority & 1)) is 0, -1, 1, 0 for priority = 0, 1, 2, 3
  double relId = hit.getSegmentID() + 0.5 * (((priority >> 1) & 1) - (priority & 1))
                 - m_TSoffset[iSL] - m_idRef[iSL][int(priority < 3)];
  relId = remainder(relId, (m_TSoffset[iSL + 1] - m_TSoffset[iSL]));
  return relId;
}

initialize()

void initialize

(

const NeuroTriggerParameters &

p

)

Definition at line 35 of file NeuroTrigger.cc.

{
  // check parameters
  bool okay = true;
  // ensure that length of lists matches number of sectors
  if (p.nHidden.size() != 1 && p.nHidden.size() != p.nMLP) {
    B2ERROR("Number of nHidden lists should be 1 or " << p.nMLP);
    okay = false;
  }
  if (p.outputScale.size() != 1 && p.outputScale.size() != p.nMLP) {
    B2ERROR("Number of outputScale lists should be 1 or " << p.nMLP);
    okay = false;
  }
  bool rangeProduct = (p.phiRangeUse.size() * p.invptRangeUse.size() * p.thetaRangeUse.size() * p.SLpattern.size() == p.nMLP);
  if (!rangeProduct) {
    if (p.phiRangeUse.size() != 1 && p.phiRangeUse.size() != p.nMLP) {
      B2ERROR("Number of phiRangeUse.lists should be 1 or " << p.nMLP);
      okay = false;
    }
    if (p.invptRangeUse.size() != 1 && p.invptRangeUse.size() != p.nMLP) {
      B2ERROR("Number of invptRangeUse.lists should be 1 or " << p.nMLP);
      okay = false;
    }
    if (p.thetaRangeUse.size() != 1 && p.thetaRangeUse.size() != p.nMLP) {
      B2ERROR("Number of thetaRangeUse.lists should be 1 or " << p.nMLP);
      okay = false;
    }
    if (p.SLpattern.size() != 1 && p.SLpattern.size() != p.nMLP) {
      B2ERROR("Number of SLpattern lists should be 1 or " << p.nMLP);
      okay = false;
    }
  }
  // ensure that length of maxHitsPerSL and SLpatternMask lists matches SLpattern
  if (p.maxHitsPerSL.size() != 1 && p.maxHitsPerSL.size() != p.SLpattern.size()) {
    B2ERROR("Number of maxHitsPerSL lists should be 1 or " << p.SLpattern.size());
    okay = false;
  }
  if (p.SLpatternMask.size() != 1 && p.SLpatternMask.size() != p.SLpattern.size()) {
    B2ERROR("Number of SLpatternMask lists should be 1 or " << p.SLpattern.size());
    okay = false;
  }
  // ensure that number of target nodes is valid
  if (p.targetZ.isSet() || p.targetTheta.isSet()) {
    unsigned short nTarget = int(p.targetZ) + int(p.targetTheta);
    if (nTarget < 1) {
      B2ERROR("No outputs! Turn on either targetZ or targetTheta.");
      okay = false;
    }
  }
  // ensure that sector ranges are valid
  for (unsigned iPhi = 0; iPhi < p.phiRangeUse.size(); ++iPhi) {
    if (p.phiRangeUse[iPhi].size() != 2) {
      B2ERROR("phiRangeUse should be exactly 2 values");
      okay = false;
      continue;
    }
    if (p.phiRangeUse[iPhi][0] >= p.phiRangeUse[iPhi][1]) {
      B2ERROR("phiRangeUse should be smaller than phiRangeUse");
      okay = false;
    }
    if (p.phiRangeUse[iPhi][0] < -360. || p.phiRangeUse[iPhi][1] > 360. ||
        (p.phiRangeUse[iPhi][1] - p.phiRangeUse[iPhi][0]) > 360.) {
      B2ERROR("phiRangeUse should be in [-360, 360], with maximal width of 360");
      okay = false;
    }
  }
  for (unsigned iPt = 0; iPt < p.invptRangeUse.size(); ++iPt) {
    if (p.invptRangeUse[iPt].size() != 2) {
      B2ERROR("invptRangeUse should be exactly 2 values");
      okay = false;
    }
    if (p.invptRangeUse[iPt][0] >= p.invptRangeUse[iPt][1]) {
      B2ERROR("invptRangeUse should be smaller than invptRangeUse");
      okay = false;
    }
  }
  for (unsigned iTheta = 0; iTheta < p.thetaRangeUse.size(); ++iTheta) {
    if (p.thetaRangeUse[iTheta].size() != 2) {
      B2ERROR("thetaRangeUse should be exactly 2 values");
      okay = false;
      continue;
    }
    if (p.thetaRangeUse[iTheta][0] >= p.thetaRangeUse[iTheta][1]) {
      B2ERROR("thetaRangeUse should be smaller than thetaRangeUse");
      okay = false;
    }
    if (p.thetaRangeUse[iTheta][0] < 0. || p.thetaRangeUse[iTheta][1] > 180.) {
      B2ERROR("thetaRangeUse should be in [0, 180]");
      okay = false;
    }
  }
  int nTarget = (int) p.targetZ + (int) p.targetTheta;
  if (p.outputScale.size() == p.nMLP || p.outputScale.size() == 1) {
    for (unsigned iScale = 0; iScale < p.outputScale.size(); ++iScale) {
      if (p.outputScale[iScale].size() != 2 * static_cast<unsigned int>(nTarget)) {
        B2ERROR("outputScale should be exactly " << 2 * nTarget << " values");
        okay = false;
      }
    }
  } else {
    B2ERROR("the size of outputscale should be 1 or match the number of experts");
  }
  // ensure that train sectors are valid
  if (p.phiRangeUse.size() != p.phiRangeTrain.size()) {
    B2ERROR("Number of phiRangeUse.lists and phiRangeTrain lists should be equal.");
    okay = false;
  } else {
    for (unsigned iPhi = 0; iPhi < p.phiRangeUse.size(); ++iPhi) {
      if (p.phiRangeTrain[iPhi].size() != 2) {
        B2ERROR("phiRangeTrain should be exactly 2 values.");
        okay = false;
      } else if (p.phiRangeTrain[iPhi][0] > p.phiRangeUse[iPhi][0] ||
                 p.phiRangeTrain[iPhi][1] < p.phiRangeUse[iPhi][1]) {
        B2ERROR("phiRangeTrain should be wider than phiRangeUse.or equal.");
        okay = false;
      }
    }
  }
  if (p.invptRangeUse.size() != p.invptRangeTrain.size()) {
    B2ERROR("Number of invptRangeUse.lists and invptRangeTrain lists should be equal.");
    okay = false;
  } else {
    for (unsigned iPt = 0; iPt < p.invptRangeUse.size(); ++iPt) {
      if (p.invptRangeTrain[iPt].size() != 2) {
        B2ERROR("invptRangeTrain should be exactly 2 values.");
        okay = false;
      } else if (p.invptRangeTrain[iPt][0] > p.invptRangeUse[iPt][0] ||
                 p.invptRangeTrain[iPt][1] < p.invptRangeUse[iPt][1]) {
        B2ERROR("invptRangeTrain should be wider than invptRangeUse.or equal.");
        okay = false;
      }
    }
  }
  if (p.thetaRangeUse.size() != p.thetaRangeTrain.size()) {
    B2ERROR("Number of thetaRangeUse.lists and thetaRangeTrain lists should be equal.");
    okay = false;
  } else {
    for (unsigned iTheta = 0; iTheta < p.thetaRangeUse.size(); ++iTheta) {
      if (p.thetaRangeTrain[iTheta].size() != 2) {
        B2ERROR("thetaRangeTrain should be exactly 2 values.");
        okay = false;
      } else if (p.thetaRangeTrain[iTheta][0] > p.thetaRangeUse[iTheta][0] ||
                 p.thetaRangeTrain[iTheta][1] < p.thetaRangeUse[iTheta][1]) {
        B2ERROR("thetaRangeTrain should be wider than thetaRangeUse.or equal.");
        okay = false;
      }
    }
  }
 
  if (!okay) return;
 
  // initialize MLPs
  m_MLPs.clear();
  for (unsigned iMLP = 0; iMLP < p.nMLP; ++iMLP) {
    //get indices for sector parameters
    //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>]
    vector<unsigned> indices = getRangeIndices(p, iMLP);
    //get number of nodes for each layer
    unsigned short maxHits = (p.maxHitsPerSL.size() == 1) ? p.maxHitsPerSL[0] : p.maxHitsPerSL[indices[3]];
    unsigned long SLpattern = p.SLpattern[indices[3]];
    unsigned long SLpatternMask = (p.SLpatternMask.size() == 1) ? p.SLpatternMask[0] : p.SLpatternMask[indices[3]];
    unsigned short nInput = 27 * maxHits;
    vector<NNTParam<float>> nHidden = (p.nHidden.size() == 1) ? p.nHidden[0] : p.nHidden[iMLP];
    vector<unsigned short> nNodes = {nInput};
    for (unsigned iHid = 0; iHid < nHidden.size(); ++iHid) {
      if (p.multiplyHidden) {
        nNodes.push_back(nHidden[iHid] * nNodes[0]);
      } else {
        nNodes.push_back(nHidden[iHid]);
      }
    }
    nNodes.push_back(nTarget);
    unsigned short targetVars = int(p.targetZ) + (int(p.targetTheta) << 1);
    // 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.
    vector<float> phiRangeUse = p.tcastvector<float>(p.phiRangeUse)[indices[0]];
    vector<float> invptRangeUse = p.tcastvector<float>(p.invptRangeUse)[indices[1]];
    vector<float> thetaRangeUse = p.tcastvector<float>(p.thetaRangeUse)[indices[2]];
    vector<float> phiRangeTrain = p.tcastvector<float>(p.phiRangeTrain)[indices[0]];
    vector<float> invptRangeTrain = p.tcastvector<float>(p.invptRangeTrain)[indices[1]];
    vector<float> thetaRangeTrain = p.tcastvector<float>(p.thetaRangeTrain)[indices[2]];
    B2DEBUG(50, "Ranges for sector " << iMLP
            << ": phiRange [" << phiRangeUse[0] << ", " << phiRangeUse[1]
            << "], invptRange [" << invptRangeUse[0] << ", " << invptRangeUse[1]
            << "], thetaRange [" << thetaRangeUse[0] << ", " << thetaRangeUse[1]
            << "], SLpattern " << SLpattern);
    //get scaling values
    vector<float> outputScale = (p.outputScale.size() == 1) ? p.tcastvector<float>(p.outputScale)[0] : p.tcastvector<float>
                                (p.outputScale)[iMLP];
    //convert phi and theta from degree to radian
    phiRangeUse[0] *= Unit::deg;
    phiRangeUse[1] *= Unit::deg;
    thetaRangeUse[0] *= Unit::deg;
    thetaRangeUse[1] *= Unit::deg;
    phiRangeTrain[0] *= Unit::deg;
    phiRangeTrain[1] *= Unit::deg;
    thetaRangeTrain[0] *= Unit::deg;
    thetaRangeTrain[1] *= Unit::deg;
    if (p.targetTheta) {
      outputScale[2 * int(p.targetZ)] *= Unit::deg;
      outputScale[2 * int(p.targetZ) + 1] *= Unit::deg;
    }
    //create new MLP
    m_MLPs.push_back(CDCTriggerMLP(nNodes, targetVars, outputScale,
                                   phiRangeUse, invptRangeUse, thetaRangeUse,
                                   phiRangeTrain, invptRangeTrain, thetaRangeTrain,
                                   maxHits, SLpattern, SLpatternMask, p.tMax,
                                   p.et_option()));
  }
 
  if (p.IDRanges.size() == p.nMLP) {
    for (auto exp : p.IDRanges) {
      // first entry is the expert number, after that follow the idranges for all the superlayers
      std::vector<float> irange = {exp.begin() + 1, exp.end()};
      m_MLPs[static_cast<int>(exp[0])].setRelID(irange);
    }
  } else if (p.IDRanges.size() == 0) {
    B2WARNING("idranges have not been initialized yet, did you forget it?");
  } else {
    B2ERROR("number of idranges should match the number of experts!");
  }
  // load some values from the geometry that will be needed for the input
  setConstants();
}

initializeCollections() [1/2]

void initializeCollections

(

std::string

hitCollectionName

)

Definition at line 307 of file NeuroTrigger.cc.

{
  m_segmentHits.isRequired(hitCollectionName);
  m_hitCollectionName = hitCollectionName;
}

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.

{
  m_segmentHits.isRequired(hitCollectionName);
  if (!((et_option == "fastestpriority") || (et_option == "etfhwin") || (et_option == "zero") || (et_option == "fastest2d"))) {
    m_eventTime.isRequired(eventTimeName);
  }
  m_hitCollectionName = hitCollectionName;
}

load()

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

Load MLPs from file.

Parameters

filename	name of the TFile to read from
arrayname	name of the TObjArray holding the MLPs in the file

Returns

true if the MLPs were loaded correctly

Definition at line 1148 of file NeuroTrigger.cc.

{
  if (filename.size() < 1) {
    m_MLPs.clear();
    m_MLPs = m_cdctriggerneuroconfig->getMLPs();
    if (m_MLPs.size() == 0) {
      B2ERROR("Could not load Neurotrigger weights from database!");
      return false;
    }
    B2INFO("Loaded Neurotrigger MLP weights from database: " +  m_cdctriggerneuroconfig->getNNName());
    B2DEBUG(100, "loaded " << m_MLPs.size() << " networks from database");
    // load some values from the geometry that will be needed for the input
    setConstants();
    return true;
  } else {
    TFile datafile(filename.c_str(), "READ");
    if (!datafile.IsOpen()) {
      B2WARNING("Could not open file " << filename);
      return false;
    }
 
    TObjArray* MLPs = (TObjArray*)datafile.Get(arrayname.c_str());
    if (!MLPs) {
      datafile.Close();
      B2WARNING("File " << filename << " does not contain key " << arrayname);
      return false;
    }
    m_MLPs.clear();
    for (int isector = 0; isector < MLPs->GetEntriesFast(); ++isector) {
      CDCTriggerMLP* expert = dynamic_cast<CDCTriggerMLP*>(MLPs->At(isector));
      if (expert) m_MLPs.push_back(*expert);
      else B2WARNING("Wrong type " << MLPs->At(isector)->ClassName() << ", ignoring this entry.");
    }
    MLPs->Clear();
    delete MLPs;
    datafile.Close();
    B2DEBUG(100, "loaded " << m_MLPs.size() << " networks");
 
    B2INFO("Loaded Neurotrigger MLP weights from file: " +  filename);
    // load some values from the geometry that will be needed for the input
    setConstants();
 
    return true;
  }
}

loadIDHist()

bool loadIDHist

(

const std::string &

filename

)

function to load idhist from file

Definition at line 1129 of file NeuroTrigger.cc.

{
  std::ifstream gzipfile(filename, ios_base::in | ios_base::binary);
  boost::iostreams::filtering_istream arrayStream;
  arrayStream.push(boost::iostreams::gzip_decompressor());
  arrayStream.push(gzipfile);
  CDCTriggerMLPData::HeaderSet hline;
  if (gzipfile.is_open()) {
    while (arrayStream >> hline) {
      for (unsigned i = 0; i < 18; ++i) {
        m_MLPs[hline.exPert].relevantID[i] = hline.relID[i];
      }
    }
  } else { return false;}
  return true;
}

nSectors()

unsigned nSectors ( ) const

inline

return number of neural networks

Definition at line 164 of file NeuroTrigger.h.

{ 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.

{ 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.

{ return m_MLPs[index]; }

runMLP()

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

Run an expert MLP.

Parameters

isector	index of the MLP
input	vector of input values

Returns

unscaled output values (z vertex in cm and/or theta in radian)

Definition at line 1027 of file NeuroTrigger.cc.

{
  const CDCTriggerMLP& expert = m_MLPs[isector];
  vector<float> weights = expert.getWeights();
  vector<float> layerinput = input;
  vector<float> layeroutput = {};
  unsigned iw = 0;
  for (unsigned il = 1; il < expert.nLayers(); ++il) {
    //add bias input
    layerinput.push_back(1.);
    //prepare output
    layeroutput.clear();
    layeroutput.assign(expert.nNodesLayer(il), 0.);
    //loop over outputs
    for (unsigned io = 0; io < layeroutput.size(); ++io) {
      //loop over inputs
      for (unsigned ii = 0; ii < layerinput.size(); ++ii) {
        layeroutput[io] += layerinput[ii] * weights[iw++];
      }
      //apply activation function
      layeroutput[io] = tanh(layeroutput[io] / 2.);
    }
    //output is new input
    layerinput = layeroutput;
  }
  return expert.unscaleTarget(layeroutput);
}

runMLPFix()

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

Run an expert MLP with fixed point arithmetic.

Definition at line 1056 of file NeuroTrigger.cc.

{
  unsigned precisionInput = m_precision[3];
  unsigned precisionWeights = m_precision[4];
  unsigned precisionLUT = m_precision[5];
  unsigned precisionTanh = m_precision[3];
  unsigned dp = precisionInput + precisionWeights - precisionLUT;
 
  const CDCTriggerMLP& expert = m_MLPs[isector];
  // transform inputs to fixed point (cut off to input precision)
  vector<long> inputFix(input.size(), 0);
  for (unsigned ii = 0; ii < input.size(); ++ii) {
    inputFix[ii] = long(input[ii] * (1 << precisionInput));
  }
  // transform weights to fixed point (round to weight precision)
  vector<float> weights = expert.getWeights();
  vector<long> weightsFix(weights.size(), 0);
  for (unsigned iw = 0; iw < weights.size(); ++iw) {
    weightsFix[iw] = long(round(weights[iw] * (1 << precisionWeights)));
  }
  // maximum input value for the tanh LUT
  unsigned xMax = unsigned(ceil(atanh(1. - 1. / (1 << (precisionTanh + 1))) *
                                (1 << (precisionLUT + 1))));
 
  // run MLP
  vector<long> layerinput = inputFix;
  vector<long> layeroutput = {};
  unsigned iw = 0;
  for (unsigned il = 1; il < expert.nLayers(); ++il) {
    // add bias input
    layerinput.push_back(1 << precisionInput);
    // prepare output
    layeroutput.clear();
    layeroutput.assign(expert.nNodesLayer(il), 0);
    // loop over outputs
    for (unsigned io = 0; io < layeroutput.size(); ++io) {
      // loop over inputs
      for (unsigned ii = 0; ii < layerinput.size(); ++ii) {
        layeroutput[io] += layerinput[ii] * weightsFix[iw++];
      }
      // apply activation function -> LUT, calculated on the fly here
      unsigned long bin = abs(layeroutput[io]) >> dp;
      // correction to get symmetrical rounding errors
      float x = (bin + 0.5 - 1. / (1 << (dp + 1))) / (1 << precisionLUT);
      long tanhLUT = (bin < xMax) ? long(round(tanh(x / 2.) * (1 << precisionTanh))) : (1 << precisionTanh);
      layeroutput[io] = (layeroutput[io] < 0) ? -tanhLUT : tanhLUT;
    }
    // output is new input
    layerinput = layeroutput;
  }
 
  // transform output back to float before unscaling
  vector<float> output(layeroutput.size(), 0.);
  for (unsigned io = 0; io < output.size(); ++io) {
    output[io] = layeroutput[io] / float(1 << precisionTanh);
  }
  return expert.unscaleTarget(output);
}

save()

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

Save MLPs to file.

Parameters

filename	name of the TFile to write to
arrayname	name of the TObjArray holding the MLPs in the file

Definition at line 1116 of file NeuroTrigger.cc.

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

selectHits()

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

Select best hits for each super layer.

Parameters

isector	index of the MLP that will use the input
track	axial hit relations are taken from given track
returnAllRelevant	if true, return all relevant hits instead of selecting the best (for making relations)

Returns

list of selected hit indices

Definition at line 855 of file NeuroTrigger.cc.

{
  const CDCTriggerMLP& expert = m_MLPs[isector];
  vector<unsigned> selectedHitIds = {};
  // prepare vectors to keep best drift times, left/right and selected hit IDs
  vector<int> tMin;
  tMin.assign(expert.nNodesLayer(0), expert.getTMax());
  vector<bool> LRknown;
  LRknown.assign(expert.nNodesLayer(0), false);
  vector<int> hitIds;
  hitIds.assign(expert.nNodesLayer(0), -1);
  vector<unsigned> nHits;
  nHits.assign(9, 0);
 
  // loop over axial hits related to input track
  RelationVector<CDCTriggerSegmentHit> axialHits =
    track.getRelationsTo<CDCTriggerSegmentHit>(m_hitCollectionName);
  B2DEBUG(250, "start axial hit loop");
  for (unsigned ihit = 0; ihit < axialHits.size(); ++ihit) {
    // skip hits with negative relation weight (not selected in finder)
    if (axialHits.weight(ihit) < 0) continue;
    unsigned short iSL = axialHits[ihit]->getISuperLayer();
    // skip stereo hits (should not be related to track, but check anyway)
    if (iSL % 2 == 1) continue;
    if (expert.getSLpatternUnmasked() != 0 &&
        !((expert.getSLpatternUnmasked() >> iSL) & 1)) {
      B2DEBUG(250, "skipping hit in SL " << iSL);
      continue;
    }
    // get priority time and apply time window cut
 
    int t = (m_hasT0) ? axialHits[ihit]->priorityTime() - m_T0 : 0;
    if (t < 0) {
      t = 0;
    } else if (t > expert.getTMax()) {
      t = expert.getTMax();
    }
    double relId = getRelId(*axialHits[ihit]);
    if (expert.isRelevant(relId, iSL)) {
      // get reference hit (worst of existing hits)
      unsigned short iRef = iSL;
      if (expert.getMaxHitsPerSL() > 1) {
        if (nHits[iSL] < expert.getMaxHitsPerSL() &&
            (expert.getSLpatternUnmasked() >> (iSL + 9 * nHits[iSL])) & 1) {
          iRef += 9 * nHits[iSL];
          ++nHits[iSL];
        } else {
          for (unsigned compare = iSL; compare < iSL + 9 * nHits[iSL]; compare += 9) {
            if ((LRknown[iRef] && !LRknown[compare]) ||
                (LRknown[iRef] == LRknown[compare] && tMin[iRef] < tMin[compare]))
              iRef = compare;
          }
        }
      }
      // choose best hit (LR known before LR unknown, then shortest drift time)
      bool useHit = false;
      if (LRknown[iRef]) {
        useHit = (axialHits[ihit]->LRknown() && t <= tMin[iRef]);
      } else {
        useHit = (axialHits[ihit]->LRknown() || t <= tMin[iRef]);
      }
      B2DEBUG(250, "relevant wire SL " << iSL << " LR " << axialHits[ihit]->getLeftRight()
              << " t " << t << " iRef " << iRef << " useHit " << useHit);
      if (useHit) {
        // keep drift time and LR
        LRknown[iRef] = axialHits[ihit]->LRknown();
        tMin[iRef] = t;
        hitIds[iRef] = axialHits[ihit]->getArrayIndex();
      }
      if (returnAllRelevant) selectedHitIds.push_back(axialHits[ihit]->getArrayIndex());
    }
  }
 
  // loop over stereo hits, choosing up to MaxHitsPerSL per superlayer
  B2DEBUG(250, "start stereo hit loop");
  for (int ihit = 0; ihit < m_segmentHits.getEntries(); ++ ihit) {
    unsigned short iSL = m_segmentHits[ihit]->getISuperLayer();
    // skip axial hits
    if (iSL % 2 == 0) continue;
    if (expert.getSLpatternUnmasked() != 0 &&
        !((expert.getSLpatternUnmasked() >> iSL) & 1)) {
      B2DEBUG(250, "skipping hit in SL " << iSL);
      continue;
    }
    // get priority time and apply time window cut
    int t = (m_hasT0) ? m_segmentHits[ihit]->priorityTime() - m_T0 : 0;
    if (t < 0) {
      t = 0;
    } else if (t > expert.getTMax()) {
      t = expert.getTMax();
    }
    double relId = getRelId(*m_segmentHits[ihit]);
    if (expert.isRelevant(relId, iSL)) {
      // get reference hit (worst of existing hits)
      unsigned short iRef = iSL;
      if (expert.getMaxHitsPerSL() > 1) {
        if (nHits[iSL] < expert.getMaxHitsPerSL() &&
            (expert.getSLpatternUnmasked() >> (iSL + 9 * nHits[iSL])) & 1) {
          iRef += 9 * nHits[iSL];
          ++nHits[iSL];
        } else {
          for (unsigned compare = iSL; compare < iSL + 9 * nHits[iSL]; compare += 9) {
            if ((LRknown[iRef] && !LRknown[compare]) ||
                (LRknown[iRef] == LRknown[compare] && tMin[iRef] < tMin[compare]))
              iRef = compare;
          }
        }
      }
      // choose best hit (LR known before LR unknown, then shortest drift time)
      bool useHit = false;
      if (LRknown[iRef]) {
        useHit = (m_segmentHits[ihit]->LRknown() && t <= tMin[iRef]);
      } else {
        useHit = (m_segmentHits[ihit]->LRknown() || t <= tMin[iRef]);
      }
      B2DEBUG(250, "relevant wire SL " << iSL << " LR " << m_segmentHits[ihit]->getLeftRight()
              << " t " << t << " iRef " << iRef << " useHit " << useHit);
      if (useHit) {
        // keep drift time and LR
        LRknown[iRef] = m_segmentHits[ihit]->LRknown();
        tMin[iRef] = t;
        hitIds[iRef] = ihit;
      }
      if (returnAllRelevant) selectedHitIds.push_back(ihit);
    }
  }
 
  // save selected hit Ids
  if (!returnAllRelevant) {
    for (unsigned iHit = 0; iHit < hitIds.size(); ++iHit) {
      if (hitIds[iHit] >= 0) selectedHitIds.push_back(hitIds[iHit]);
    }
  }
  return selectedHitIds;
}

selectHitsHWSim()

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

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

Parameters

isector	index of the MLP that will use the input
track	all hit relations are taken from given track

Returns

list of selected hit indices

Definition at line 776 of file NeuroTrigger.cc.

{
  const CDCTriggerMLP& expert = m_MLPs[isector];
  vector<unsigned> selectedHitIds = {};
  // prepare vectors to keep best drift times, left/right and selected hit IDs
  vector<int> tMin;
  tMin.assign(expert.nNodesLayer(0), expert.getTMax());
  vector<bool> LRknown;
  LRknown.assign(expert.nNodesLayer(0), false);
  vector<int> hitIds;
  hitIds.assign(expert.nNodesLayer(0), -1);
  vector<unsigned> nHits;
  nHits.assign(9, 0);
 
  // loop over all hits related to input track
  RelationVector<CDCTriggerSegmentHit> allHits =
    track.getRelationsTo<CDCTriggerSegmentHit>(m_hitCollectionName);
  B2DEBUG(250, "start hit loop over all related hits");
  for (unsigned ihit = 0; ihit < allHits.size(); ++ihit) {
    // skip hits with negative relation weight (not selected in finder)
    if (allHits.weight(ihit) < 0) continue;
    unsigned short iSL = allHits[ihit]->getISuperLayer();
    if (expert.getSLpatternUnmasked() != 0 &&
        !((expert.getSLpatternUnmasked() >> iSL) & 1)) {
      B2DEBUG(250, "skipping hit in SL " << iSL);
      continue;
    }
    // get priority time and apply time window cut
    int t = (m_hasT0) ? allHits[ihit]->priorityTime() - m_T0 : 0;
    if (t < 0) {
      t = 0;
    } else if (t > expert.getTMax()) {
      t = expert.getTMax();
    }
    double relId = getRelId(*allHits[ihit]);
    if (expert.isRelevant(relId, iSL)) {
      // get reference hit (worst of existing hits)
      unsigned short iRef = iSL;
      if (expert.getMaxHitsPerSL() > 1) {
        if (nHits[iSL] < expert.getMaxHitsPerSL() &&
            (expert.getSLpatternUnmasked() >> (iSL + 9 * nHits[iSL])) & 1) {
          iRef += 9 * nHits[iSL];
          ++nHits[iSL];
        } else {
          for (unsigned compare = iSL; compare < iSL + 9 * nHits[iSL]; compare += 9) {
            if ((LRknown[iRef] && !LRknown[compare]) ||
                (LRknown[iRef] == LRknown[compare] && tMin[iRef] < tMin[compare]))
              iRef = compare;
          }
        }
      }
      // choose best hit (LR known before LR unknown, then shortest drift time)
      bool useHit = false;
      if (LRknown[iRef]) {
        useHit = (allHits[ihit]->LRknown() && t <= tMin[iRef]);
      } else {
        useHit = (allHits[ihit]->LRknown() || t <= tMin[iRef]);
      }
      B2DEBUG(250, "relevant wire SL " << iSL << " LR " << allHits[ihit]->getLeftRight()
              << " t " << t << " iRef " << iRef << " useHit " << useHit);
      if (useHit) {
        // keep drift time and LR
        LRknown[iRef] = allHits[ihit]->LRknown();
        tMin[iRef] = t;
        hitIds[iRef] = allHits[ihit]->getArrayIndex();
      }
    }
  }
  // save selected hit Ids
  for (unsigned iHit = 0; iHit < hitIds.size(); ++iHit) {
    if (hitIds[iHit] >= 0) selectedHitIds.push_back(hitIds[iHit]);
  }
  return selectedHitIds;
}

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.

{
  if (MLPs.size() == 0) {
    return -1;
  }
 
  int bestIndex = -1;
  for (unsigned i = 0; i < MLPs.size(); ++i) {
    int isector = MLPs[i];
    unsigned long sectorPattern = m_MLPs[isector].getSLpattern();
    B2DEBUG(250, "hitPattern " << pattern << " sectorPattern " << sectorPattern);
    // no hit pattern restriction -> keep looking for exact match
    if (sectorPattern == 0) {
      B2DEBUG(250, "found match for general sector");
      bestIndex = isector;
    }
    // exact match -> keep this sector and stop searching
    if (pattern == sectorPattern) {
      B2DEBUG(250, "found match for hit pattern " << pattern);
      bestIndex = isector;
      break;
    }
  }
 
  if (bestIndex < 0) {
    if (neurotrackinputmode) {
      B2DEBUG(150, "No sector found to match pattern, using sector 0" << pattern << ".");
      bestIndex = 0;
    } else {
      B2DEBUG(150, "No sector found to match pattern " << pattern << ".");
    }
  }
  return bestIndex;
}

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.

{
  vector<int> indices = {};
 
  if (m_MLPs.size() == 0) {
    B2WARNING("Trying to select MLP before initializing MLPs.");
    return indices;
  }
 
  // find all matching sectors
  for (unsigned isector = 0; isector < m_MLPs.size(); ++isector) {
    if (m_MLPs[isector].inPhiRangeUse(phi0) && m_MLPs[isector].inInvptRangeUse(invpt)
        && m_MLPs[isector].inThetaRangeUse(theta)) {
      indices.push_back(isector);
    }
  }
 
  if (indices.size() == 0) {
    B2DEBUG(150, "Track does not match any sector.");
    B2DEBUG(150, "invpt=" << invpt << ", phi=" << phi0 * 180. / M_PI << ", theta=" << theta * 180. / M_PI);
  }
 
  return indices;
}

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.

{
  vector<int> indices = {};
 
  if (m_MLPs.size() == 0) {
    B2WARNING("Trying to select MLP before initializing MLPs.");
    return indices;
  }
 
  // find all matching sectors
  for (unsigned isector = 0; isector < m_MLPs.size(); ++isector) {
    if (m_MLPs[isector].inPhiRangeTrain(phi0) && m_MLPs[isector].inInvptRangeTrain(invpt)
        && m_MLPs[isector].inThetaRangeTrain(theta)) {
      indices.push_back(isector);
    }
  }
 
  if (indices.size() == 0) {
    B2DEBUG(150, "Track does not match any sector.");
    B2DEBUG(150, "invpt=" << invpt << ", phi=" << phi0 * 180. / M_PI << ", theta=" << theta * 180. / M_PI);
  }
 
  return indices;
}

setConstants()

void setConstants

(

)

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

Definition at line 280 of file NeuroTrigger.cc.

{
  CDCGeometryPar& cdc = CDCGeometryPar::Instance();
  int layerId = 3;
  int nTS = 0;
  for (int iSL = 0; iSL < 9; ++iSL) {
    m_TSoffset[iSL] = nTS;
    nTS += cdc.nWiresInLayer(layerId);
    m_TSoffset[iSL + 1] = nTS;
    for (int priority = 0; priority < 2; ++ priority) {
      m_radius[iSL][priority] = cdc.senseWireR(layerId + priority);
    }
    layerId += (iSL > 0 ? 6 : 7);
  }
}

setPrecision()

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

inline

set fixed point precision

Definition at line 151 of file NeuroTrigger.h.

{ 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.

{
  double omega = track.getOmega(); // signed track curvature
  for (int iSL = 0; iSL < 9; ++iSL) {
    for (int priority = 0; priority < 2; ++priority) {
      m_alpha[iSL][priority] = (m_radius[iSL][priority] * abs(omega) < 2.) ?
                               asin(m_radius[iSL][priority] * omega / 2.) :
                               M_PI_2;
      m_idRef[iSL][priority] = remainder(((track.getPhi0() - m_alpha[iSL][priority]) *
                                          (m_TSoffset[iSL + 1] - m_TSoffset[iSL]) / 2. / M_PI),
                                         (m_TSoffset[iSL + 1] - m_TSoffset[iSL]));
    }
  }
}

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.

{
  //unsigned precisionPhi = m_precision[0];
  //unsigned precisionAlpha = m_precision[0];
  unsigned precisionPhiAlpha = m_precision[0];
  unsigned precisionScale = m_precision[1];
  unsigned precisionId = m_precision[2];
 
  double omega = track.getOmega();
  long phi = round(track.getPhi0() * (1 << precisionPhiAlpha));
 
  for (int iSL = 0; iSL < 9; ++iSL) {
    for (int priority = 0; priority < 2; ++priority) {
      // LUT, calculated on the fly here
      m_alpha[iSL][priority] =
        (m_radius[iSL][priority] * abs(omega) < 2) ?
        round(asin(m_radius[iSL][priority] * omega / 2) * (1 << precisionPhiAlpha)) :
        round(M_PI_2 * (1 << precisionPhiAlpha));
      long dphi = phi - (long(m_alpha[iSL][priority]));
      m_idRef[iSL][priority] =
        long(dphi * round((m_TSoffset[iSL + 1] - m_TSoffset[iSL]) / 2. / M_PI * (1 << precisionScale)) /
             (long(1) << (precisionPhiAlpha + precisionScale - precisionId)));
      // unscale after rounding
      m_alpha[iSL][priority] /= (1 << precisionPhiAlpha);
      m_idRef[iSL][priority] /= (1 << precisionId);
    }
  }
}

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:

• trg/cdc/include/NeuroTrigger.h
• trg/cdc/src/NeuroTrigger.cc