ATLAS Offline Software
Loading...
Searching...
No Matches
Trig::StaticBunchCrossingTool Class Reference

Stand-alone bunch group tool knowing some static configurations. More...

#include <StaticBunchCrossingTool.h>

Inheritance diagram for Trig::StaticBunchCrossingTool:
Collaboration diagram for Trig::StaticBunchCrossingTool:

Public Types

enum  BeamType { Beam1 = 0 , Beam2 = 1 , Crossing = 2 }
 Types of the return values of the bcIntensity function. More...
enum  BunchCrossingType {
  Empty = 0 , FirstEmpty = 1 , MiddleEmpty = 2 , Single = 100 ,
  Front = 200 , Middle = 201 , Tail = 202 , Unpaired = 300
}
 Simplified type for a given bunch crossing. More...
enum  BunchDistanceType { NanoSec , BunchCrossings , FilledBunches }
 Enumeration specifying the units in which to expect the bunch distance type. More...
enum  BunchFillType {
  CollidingBunch = 0 , UnpairedBunch = 1 , EmptyBunch = 2 , UnpairedBeam1 = 3 ,
  UnpairedBeam2 = 4
}
 Enumeration specifying what kind of bunch to use in the gap functions. More...
typedef unsigned int bcid_type
 Declare the interface that this class provides.
Definition of the StoreGate-like object's definition
typedef ServiceHandle< StoreGateSvcMetaStore_t
 Type of the metadata store object in Athena.
typedef const ServiceHandle< StoreGateSvc > & MetaStorePtr_t
 Type of the metadata store pointer in standalone mode.

Public Member Functions

 StaticBunchCrossingTool (const std::string &name="StaticBunchCrossingTool")
 Create a proper constructor for Athena.
virtual StatusCode initialize ()
 Function initialising the tool.
StatusCode loadConfig (int bgkey)
 Load a hard-coded bunch group key.
StatusCode loadConfig (const std::vector< int > &filledBunches, const std::vector< float > &filledIntensities=std::vector< float >(), const std::vector< int > &unpairedBunches=std::vector< int >(), const std::vector< float > &unpairedIntensities=std::vector< float >())
 Load a configuration specified by the user.
StatusCode loadConfig (const std::vector< float > &bunches)
 Load the configuration specified by the user.
virtual void print () const
 Print the state of the tool.
virtual StatusCode sysInitialize ()
 Function initialising the tool in the correct way in Athena.
ServiceHandle< StoreGateSvc > & evtStore ()
 The standard StoreGateSvc (event store) Returns (kind of) a pointer to the StoreGateSvc.
const ServiceHandle< StoreGateSvc > & detStore () const
 The standard StoreGateSvc/DetectorStore Returns (kind of) a pointer to the StoreGateSvc.
virtual StatusCode sysStart () override
 Handle START transition.
virtual std::vector< Gaudi::DataHandle * > inputHandles () const override
 Return this algorithm's input handles.
virtual std::vector< Gaudi::DataHandle * > outputHandles () const override
 Return this algorithm's output handles.
Gaudi::Details::PropertyBase & declareProperty (Gaudi::Property< T, V, H > &t)
void updateVHKA (Gaudi::Details::PropertyBase &)
MsgStream & msg () const
bool msgLvl (const MSG::Level lvl) const
Functions implementing the IBunchCrossingTool interface
virtual bool isFilled (bcid_type bcid) const
 The simplest query: Is the bunch crossing filled or not?
virtual bool isInTrain (bcid_type bcid) const
 Function deciding if a given bunch crossing is in a filled train.
virtual bool isUnpaired (bcid_type bcid) const
 Function deciding if a given bunch crossing has an unpaired bunch.
virtual bool isBeam1 (bcid_type bcid) const
 Function deciding if there was a bunch from "beam 1" in this bunch crossing.
virtual bool isBeam2 (bcid_type bcid) const
 Function deciding if there was a bunch from "beam 2" in this bunch crossing.
virtual float bcIntensity (bcid_type bcid, BeamType type=Crossing) const
 Function returning the "intensity" of a given bunch crossing.
virtual BunchCrossingType bcType (bcid_type bcid) const
 Get the type of the specific bunch crossing.
virtual int distanceFromFront (bcid_type bcid, BunchDistanceType type=NanoSec) const
 The distance of the specific bunch crossing from the front of the train.
virtual int distanceFromTail (bcid_type bcid, BunchDistanceType type=NanoSec) const
 The distance of the specific bunch crossing from the tail of the train.
virtual int gapBeforeTrain (bcid_type bcid, BunchDistanceType type=NanoSec) const
 Gap before the train this BCID is in.
virtual int gapAfterTrain (bcid_type bcid, BunchDistanceType type=NanoSec) const
 Gap after the train this BCID is in.
virtual int gapBeforeBunch (bcid_type bcid, BunchDistanceType dtype=NanoSec, BunchFillType ftype=CollidingBunch) const
 Gap before a particular bunch.
virtual int gapAfterBunch (bcid_type bcid, BunchDistanceType dtype=NanoSec, BunchFillType ftype=CollidingBunch) const
 Gap after a particular bunch.
virtual std::vector< bool > bunchesInFront (bcid_type bcid, int bunches=10) const
 Function returning whether the previous bunches were filled, and how.
virtual std::vector< bool > bunchesAfter (bcid_type bcid=0, int bunches=10) const
 Function returning whether the following bunches were filled, and how.
virtual std::vector< float > bunchIntInFront (bcid_type bcid, int bunches=10, BeamType type=Crossing) const
 Function returning the intensities of the bunch crossings before the reference.
virtual std::vector< float > bunchIntAfter (bcid_type bcid, int bunches=10, BeamType type=Crossing) const
 Function returning the intensities of the bunch crossings after the reference.
virtual unsigned int numberOfFilledBunches () const
 Get the number of filled bunches in the current configuration.
virtual unsigned int numberOfUnpairedBunches () const
 Get the number of unpaired bunches in the current configuration.
virtual unsigned int numberOfBunchTrains () const
 Get the number of the bunch trains in the current configuration.
virtual int bunchTrainSpacing (BunchDistanceType type=NanoSec) const
 Get the bunch spacing in the trains.
Functions providing access to the input/output metadata
MetaStorePtr_t inputMetaStore () const
 Accessor for the input metadata store.
MetaStorePtr_t outputMetaStore () const
 Accessor for the output metadata store.
Additional helper functions, not directly mimicking Athena
template<class T>
const T * getProperty (const std::string &name) const
 Get one of the tool's properties.
const std::string & msg_level_name () const __attribute__((deprecated))
 A deprecated function for getting the message level's name.
const std::string & getName (const void *ptr) const
 Get the name of an object that is / should be in the event store.
SG::sgkey_t getKey (const void *ptr) const
 Get the (hashed) key of an object that is in the event store.

Protected Member Functions

StatusCode loadSingleBunches (const std::vector< int > &bunches, const std::vector< float > &bunch_int1=std::vector< float >(), const std::vector< float > &bunch_int2=std::vector< float >())
 Interpret the configuration for single bunches.
StatusCode loadBunchTrains (const std::vector< int > &bunches, const std::vector< float > &bunch_int1=std::vector< float >(), const std::vector< float > &bunch_int2=std::vector< float >())
 Interpret the configuration for bunch trains.
StatusCode loadUnpairedBunches (const std::vector< int > &beam1, const std::vector< int > &beam2, const std::vector< float > &bunch_int1=std::vector< float >(), const std::vector< float > &bunch_int2=std::vector< float >())
 Interpret the configuration for unpaired bunches.
void printConfig () const
 Function printing the configuration of the tool.
void renounceArray (SG::VarHandleKeyArray &handlesArray)
 remove all handles from I/O resolution
std::enable_if_t< std::is_void_v< std::result_of_t< decltype(&T::renounce)(T)> > &&!std::is_base_of_v< SG::VarHandleKeyArray, T > &&std::is_base_of_v< Gaudi::DataHandle, T >, void > renounce (T &h)
void extraDeps_update_handler (Gaudi::Details::PropertyBase &ExtraDeps)
 Add StoreName to extra input/output deps as needed.
Callback functions helping in metadata reading/writing
void setUseIncidents (const bool flag)
virtual void handle (const Incident &inc)
 Function receiving incidents from IncidentSvc/TEvent.
virtual StatusCode beginInputFile ()
 Function called when a new input file is opened.
virtual StatusCode endInputFile ()
 Function called when the currently open input file got completely processed.
virtual StatusCode beginEvent ()
 Function called when a new events is loaded.
virtual StatusCode metaDataStop ()
 Function called when the tool should write out its metadata.

Protected Attributes

Variables holding the decoded bunch structure
std::set< Trig::BunchCrossingm_filledBunches
 List of colliding bunches.
std::set< Trig::BunchCrossingm_singleBunches
 Internal list of single bunches.
std::set< Trig::BunchCrossingm_unpairedBunches
 Internal list of unpaired bunches.
std::set< Trig::BunchTrainm_bunchTrains
 Internal list of bunch trains.
Configurable tool properties
int m_maxBunchSpacing
 The maximum bunch spacing that the tool should consider.
int m_frontLength
 Length of the "front" of a bunch train.
int m_tailLength
 Length of the "tail" of a bunch train.

Private Types

typedef ServiceHandle< StoreGateSvcStoreGateSvc_t

Private Member Functions

int bunchSpacing (const std::vector< int > &bunches) const
 Get the apparent bunch spacing in the current configuration.
Gaudi::Details::PropertyBase & declareGaudiProperty (Gaudi::Property< T, V, H > &hndl, const SG::VarHandleKeyType &)
 specialization for handling Gaudi::Property<SG::VarHandleKey>

Private Attributes

int m_bgkey
 Default key to be loaded.
std::map< int, std::vector< int > > m_knownBGKeys
 All the hard-coded configs.
MetaStore_t m_inputMetaStore
 Object accessing the input metadata store.
MetaStore_t m_outputMetaStore
 Object accessing the output metadata store.
bool m_beginInputFileCalled
 Flag helping to discover when the tool misses the opening of the first input file.
bool m_useIncidents
StoreGateSvc_t m_evtStore
 Pointer to StoreGate (event store by default)
StoreGateSvc_t m_detStore
 Pointer to StoreGate (detector store by default)
std::vector< SG::VarHandleKeyArray * > m_vhka
bool m_varHandleArraysDeclared

Detailed Description

Stand-alone bunch group tool knowing some static configurations.

   This tool can be used out of the box to get information about some
   of the bunch layouts which were used to take data. It takes its
   knowledge from the "StaticConfigs.h" file.

   It can be used in a standalone analysis code to analyse ntuples
   for instance.
Author
Attila Krasznahorkay Attil.nosp@m.a.Kr.nosp@m.aszna.nosp@m.hork.nosp@m.ay@ce.nosp@m.rn.c.nosp@m.h
Revision
749252
Date
2016-05-24 11:30:51 +0200 (Tue, 24 May 2016)

Definition at line 35 of file StaticBunchCrossingTool.h.

Member Typedef Documentation

◆ bcid_type

typedef unsigned int Trig::IBunchCrossingTool::bcid_type
inherited

Declare the interface that this class provides.

Convenience type definition

Definition at line 47 of file IBunchCrossingTool.h.

◆ MetaStore_t

Type of the metadata store object in Athena.

Definition at line 66 of file AsgMetadataTool.h.

◆ MetaStorePtr_t

Type of the metadata store pointer in standalone mode.

Definition at line 68 of file AsgMetadataTool.h.

◆ StoreGateSvc_t

typedef ServiceHandle<StoreGateSvc> AthCommonDataStore< AthCommonMsg< AlgTool > >::StoreGateSvc_t
privateinherited

Definition at line 388 of file AthCommonDataStore.h.

Member Enumeration Documentation

◆ BeamType

Types of the return values of the bcIntensity function.

The different information sources provide the "bunch intensity" information in quite different ways, and the information can mean different things actually. This enumeration is used to specify what exact type of information the user is looking for.

Enumerator
Beam1 

The returned intensity should be for "beam 1".

Beam2 

The returned intensity should be for "beam 2".

Crossing 

The returned intensity should describe the BC.

Definition at line 118 of file IBunchCrossingTool.h.

118 {
119 Beam1 = 0,
120 Beam2 = 1,
121 Crossing = 2
122 };
@ Beam1
The returned intensity should be for "beam 1".
@ Beam2
The returned intensity should be for "beam 2".
@ Crossing
The returned intensity should describe the BC.

◆ BunchCrossingType

Simplified type for a given bunch crossing.

This enumeration can specify what kind of bunch crossing one BCID belongs to. The types could easily be extended later on.

Enumerator
Empty 

An empty bunch far away from filled bunches.

FirstEmpty 

The first empty bunch after a train.

MiddleEmpty 

An empty BCID in the middle of a train.

Single 

This is a filled, single bunch (not in a train)

Front 

The BCID belongs to the first few bunches in a train.

Middle 

The BCID belongs to the middle bunches in a train.

Tail 

The BCID belongs to the last few bunces in a train.

Unpaired 

This is an unpaired bunch (either beam1 or beam2)

Definition at line 145 of file IBunchCrossingTool.h.

145 {
146 Empty = 0,
147 FirstEmpty = 1,
148 MiddleEmpty = 2,
149 Single = 100,
150 Front = 200,
151 Middle = 201,
152 Tail = 202,
153 Unpaired = 300
154 };
@ Unpaired
This is an unpaired bunch (either beam1 or beam2)
@ FirstEmpty
The first empty bunch after a train.
@ Tail
The BCID belongs to the last few bunces in a train.
@ Empty
An empty bunch far away from filled bunches.
@ Front
The BCID belongs to the first few bunches in a train.
@ Middle
The BCID belongs to the middle bunches in a train.
@ Single
This is a filled, single bunch (not in a train)
@ MiddleEmpty
An empty BCID in the middle of a train.

◆ BunchDistanceType

Enumeration specifying the units in which to expect the bunch distance type.

To make it clear for the following functions what units to interpret their return values in, it is possible to request their return values in different units.

Enumerator
NanoSec 

Distance in nanoseconds.

BunchCrossings 

Distance in units of 25 nanoseconds.

FilledBunches 

Distance in units of filled bunches (depends on filling scheme)

Definition at line 174 of file IBunchCrossingTool.h.

174 {
175 NanoSec,
179 };
@ BunchCrossings
Distance in units of 25 nanoseconds.
@ NanoSec
Distance in nanoseconds.
@ FilledBunches
Distance in units of filled bunches (depends on filling scheme)

◆ BunchFillType

Enumeration specifying what kind of bunch to use in the gap functions.

The following functions can be used to calculate the gap before and after a specific BCID to some other bunch type. The gap can actually be wrt. two different types of bunches. The user may be interested between the space of two filled bunches, the space between an unpaired bunch and the previous filled bunch, the space between two unpaired bunches, or the space between a filled bunch and the previous unpaired bunch.

The empty type is just put here for completeness. Maybe once we'll be using 25 ns spacing in the bunch trains, this will be a useful parameter as well.

This enumeration helps in answeing all of these questions.

Enumerator
CollidingBunch 

The gap should be calculated wrt. the closest colling bunch.

UnpairedBunch 

The gap should be calculated wrt. the closest unpaired bunch.

EmptyBunch 

The gap should be calculated wrt. the closest empty bunch.

UnpairedBeam1 

The gap should be calculated wrt.

the closest unpaired bunch from beam 1

UnpairedBeam2 

The gap should be calculated wrt.

the closest unpaired bunch from beam 2

Definition at line 274 of file IBunchCrossingTool.h.

274 {
276 CollidingBunch = 0,
278 UnpairedBunch = 1,
280 EmptyBunch = 2,
283 UnpairedBeam1 = 3,
286 UnpairedBeam2 = 4
287 };
@ CollidingBunch
The gap should be calculated wrt. the closest colling bunch.
@ UnpairedBunch
The gap should be calculated wrt. the closest unpaired bunch.
@ EmptyBunch
The gap should be calculated wrt. the closest empty bunch.
@ UnpairedBeam1
The gap should be calculated wrt.
@ UnpairedBeam2
The gap should be calculated wrt.

Constructor & Destructor Documentation

◆ StaticBunchCrossingTool()

Trig::StaticBunchCrossingTool::StaticBunchCrossingTool ( const std::string & name = "StaticBunchCrossingTool")

Create a proper constructor for Athena.

The constructor reads all the configurations out of "StaticConfigs.h" and sets up the tool using those.

Default constructor

Definition at line 17 of file StaticBunchCrossingTool.cxx.

18 : BunchCrossingToolBase( name ),
20
21 // Declare the property of the tool:
22 declareProperty( "BGKey", m_bgkey = 0 );
23
24 //
25 // Learn the static bunch group configurations:
26 //
27 for( int i = 0; i < BGK_CONF_N; ++i ) {
28 for( int j = 0; j < BGK_CONF_NUM[ i ]; ++j ) {
29 m_knownBGKeys[ BGK_CONF_KEY[ i ] ].push_back( BGK_CONF[ i ][ j ] );
30 }
31 }
32 }
static const int BGK_CONF_KEY[BGK_CONF_N]
static const int BGK_CONF_NUM[BGK_CONF_N]
static const int *const BGK_CONF[BGK_CONF_N]
static const int BGK_CONF_N
Gaudi::Details::PropertyBase & declareProperty(Gaudi::Property< T, V, H > &t)
BunchCrossingToolBase(const std::string &name="BunchCrossingToolBase")
Default constructor.
int m_bgkey
Default key to be loaded.
std::map< int, std::vector< int > > m_knownBGKeys
All the hard-coded configs.

Member Function Documentation

◆ bcIntensity()

float Trig::BunchCrossingToolBase::bcIntensity ( bcid_type bcid,
BeamType type = Crossing ) const
virtualinherited

Function returning the "intensity" of a given bunch crossing.

Implements Trig::IBunchCrossingTool.

Definition at line 106 of file BunchCrossingToolBase.cxx.

107 {
108
109 // Check if this is a colliding bunch:
110 std::set< BunchCrossing >::const_iterator itr;
111 if( ( itr = m_filledBunches.find( bcid ) ) != m_filledBunches.end() ) {
112 switch( type ) {
113 case Crossing:
114 if( itr->intensityBeam2() > 0.001 ) {
115 ATH_MSG_WARNING( "Crossing intensity not available, ask for "
116 << "separate beam intensities instead" );
117 return 0.0;
118 } else {
119 return itr->intensityBeam1();
120 }
121 case Beam1:
122 if( std::abs( itr->intensityBeam2() ) < 0.001 ) {
123 ATH_MSG_WARNING( "Separate beam intensities not available, ask "
124 << "for the crossing intensity instead" );
125 return 0.0;
126 } else {
127 return itr->intensityBeam1();
128 }
129 break;
130 case Beam2:
131 return itr->intensityBeam2();
132 break;
133 default:
134 ATH_MSG_ERROR( "Unknown intensity type requested (" << type
135 << ")" );
136 return -1.0;
137 }
138 }
139
140 // Check if this is an unpaired bunch:
141 if( ( itr = m_unpairedBunches.find( bcid ) ) !=
142 m_unpairedBunches.end() ) {
143 switch( type ) {
144 case Beam1:
145 return itr->intensityBeam1();
146 break;
147 case Beam2:
148 return itr->intensityBeam2();
149 break;
150 case Crossing:
151 ATH_MSG_WARNING( "Crossing intensity requested for unpaired bunch ("
152 << "bcid=" << bcid << ")" );
153 return 0.0;
154 default:
155 ATH_MSG_ERROR( "Unknown intensity type requested (" << type
156 << ")" );
157 return -1.0;
158 }
159 }
160
161 // If neither, then its intensity is 0.0 by definition:
162 return 0.0;
163 }
#define ATH_MSG_ERROR(x)
#define ATH_MSG_WARNING(x)
std::set< Trig::BunchCrossing > m_unpairedBunches
Internal list of unpaired bunches.
std::set< Trig::BunchCrossing > m_filledBunches
List of colliding bunches.

◆ bcType()

IBunchCrossingTool::BunchCrossingType Trig::BunchCrossingToolBase::bcType ( bcid_type bcid) const
virtualinherited

Get the type of the specific bunch crossing.

Implements Trig::IBunchCrossingTool.

Definition at line 166 of file BunchCrossingToolBase.cxx.

166 {
167
168 // First the obvious check:
169 if( ! isFilled( bcid ) ) {
170 // Check if it's an unpaired bunch:
171 if( isUnpaired( bcid ) ) {
172 return Unpaired;
173 }
174 // If the previous bunch crossing is the tail of a bunch train:
175 if( ! distanceFromTail( bcid - 1, BunchCrossings ) ) {
176 return FirstEmpty;
177 }
178 // Check if it's in the middle of a bunch train:
179 std::set< BunchTrain >::const_iterator itr = m_bunchTrains.begin();
180 std::set< BunchTrain >::const_iterator end = m_bunchTrains.end();
181 for( ; itr != end; ++itr ) {
182 if( itr->isInside( bcid ) ) {
183 return MiddleEmpty;
184 }
185 }
186 // If none of the above are true, it has to be a "simple" empty bunch:
187 return Empty;
188 }
189
190 // Now we know that the bunch has to be a filled one...
191
192 // If it's not in a train, it has to be a single filled bunch:
193 if( ! isInTrain( bcid ) ) return Single;
194
195 // Let's check if it is close to the front of a bunch train:
196 int distance = distanceFromFront( bcid, NanoSec );
197 if( ( distance >= 0 ) && ( distance <= m_frontLength ) ) {
198 return Front;
199 }
200 // Now let's check if it's close to the tail of a bunch train:
202 if( ( distance >= 0 ) && ( distance <= m_tailLength ) ) {
203 return Tail;
204 }
205
206 // If none of the above are true, it has to be in the middle of a train:
207 return Middle;
208 }
virtual bool isFilled(bcid_type bcid) const
The simplest query: Is the bunch crossing filled or not?
virtual int distanceFromFront(bcid_type bcid, BunchDistanceType type=NanoSec) const
The distance of the specific bunch crossing from the front of the train.
virtual bool isUnpaired(bcid_type bcid) const
Function deciding if a given bunch crossing has an unpaired bunch.
virtual int distanceFromTail(bcid_type bcid, BunchDistanceType type=NanoSec) const
The distance of the specific bunch crossing from the tail of the train.
int m_frontLength
Length of the "front" of a bunch train.
std::set< Trig::BunchTrain > m_bunchTrains
Internal list of bunch trains.
virtual bool isInTrain(bcid_type bcid) const
Function deciding if a given bunch crossing is in a filled train.
int m_tailLength
Length of the "tail" of a bunch train.
int distance(const BunchCrossing bc1, const BunchCrossing bc2)
I need this function only for technical reasons.

◆ beginEvent()

StatusCode asg::AsgMetadataTool::beginEvent ( )
protectedvirtualinherited

◆ beginInputFile()

StatusCode asg::AsgMetadataTool::beginInputFile ( )
protectedvirtualinherited

Function called when a new input file is opened.

Dummy implementation that can be overridden by the derived tool.

Reimplemented in AsgElectronEfficiencyCorrectionTool, BookkeeperDumperTool, BookkeeperTool, PMGTools::PMGTruthWeightTool, TauAnalysisTools::TauEfficiencyCorrectionsTool, TauAnalysisTools::TauSmearingTool, Trig::TrigDecisionTool, Trig::xAODBunchCrossingTool, TrigConf::xAODConfigTool, xAODMaker::TriggerMenuMetaDataTool, and xAODMaker::TruthMetaDataTool.

Definition at line 185 of file AsgMetadataTool.cxx.

185 {
186
187 // Return gracefully:
188 return StatusCode::SUCCESS;
189 }

◆ bunchesAfter()

std::vector< bool > Trig::BunchCrossingToolBase::bunchesAfter ( bcid_type bcid = 0,
int bunches = 10 ) const
virtualinherited

Function returning whether the following bunches were filled, and how.

Implements Trig::IBunchCrossingTool.

Definition at line 734 of file BunchCrossingToolBase.cxx.

735 {
736
737 // The only thing we have to be careful about is the bunches near the
738 // "turnover" region of the BCIDs. That's why I use the BunchCrossing
739 // class here:
740 std::vector< bool > result;
741 for( int i = 0; i < bunches; ++i ) {
742 result.push_back( isFilled( BunchCrossing( bcid ) +
743 BunchCrossing( i ) ) );
744 }
745
746 return result;
747 }

◆ bunchesInFront()

std::vector< bool > Trig::BunchCrossingToolBase::bunchesInFront ( bcid_type bcid,
int bunches = 10 ) const
virtualinherited

Function returning whether the previous bunches were filled, and how.

Implements Trig::IBunchCrossingTool.

Definition at line 718 of file BunchCrossingToolBase.cxx.

719 {
720
721 // The only thing we have to be careful about is the bunches near the
722 // "turnover" region of the BCIDs. That's why I use the BunchCrossing
723 // class here:
724 std::vector< bool > result;
725 for( int i = 0; i < bunches; ++i ) {
726 result.push_back( isFilled( BunchCrossing( bcid ) -
727 BunchCrossing( i ) ) );
728 }
729
730 return result;
731 }

◆ bunchIntAfter()

std::vector< float > Trig::BunchCrossingToolBase::bunchIntAfter ( bcid_type bcid,
int bunches = 10,
BeamType type = Crossing ) const
virtualinherited

Function returning the intensities of the bunch crossings after the reference.

Implements Trig::IBunchCrossingTool.

Definition at line 781 of file BunchCrossingToolBase.cxx.

783 {
784
785 std::vector< float > result;
786 for( int i = 0; i < bunches; ++i ) {
787 std::set< BunchCrossing >::const_iterator itr =
788 m_filledBunches.find( BunchCrossing( bcid ) + BunchCrossing( i ) );
789 if( itr != m_filledBunches.end() ) {
790 switch( type ) {
791 case Beam1:
792 case Crossing:
793 result.push_back( itr->intensityBeam1() );
794 break;
795 case Beam2:
796 result.push_back( itr->intensityBeam2() );
797 break;
798 default:
799 ATH_MSG_ERROR( "Unknown intensity type requested ("
800 << type << ")" );
801 return result;
802 }
803 } else {
804 result.push_back( 0.0 );
805 }
806 }
807
808 return result;
809 }

◆ bunchIntInFront()

std::vector< float > Trig::BunchCrossingToolBase::bunchIntInFront ( bcid_type bcid,
int bunches = 10,
BeamType type = Crossing ) const
virtualinherited

Function returning the intensities of the bunch crossings before the reference.

Implements Trig::IBunchCrossingTool.

Definition at line 750 of file BunchCrossingToolBase.cxx.

752 {
753
754 std::vector< float > result;
755 for( int i = 0; i < bunches; ++i ) {
756 std::set< BunchCrossing >::const_iterator itr =
757 m_filledBunches.find( BunchCrossing( bcid ) - BunchCrossing( i ) );
758 if( itr != m_filledBunches.end() ) {
759 switch( type ) {
760 case Beam1:
761 case Crossing:
762 result.push_back( itr->intensityBeam1() );
763 break;
764 case Beam2:
765 result.push_back( itr->intensityBeam2() );
766 break;
767 default:
768 ATH_MSG_ERROR( "Unknown intensity type requested ("
769 << type << ")" );
770 return result;
771 }
772 } else {
773 result.push_back( 0.0 );
774 }
775 }
776
777 return result;
778 }

◆ bunchSpacing()

int Trig::BunchCrossingToolBase::bunchSpacing ( const std::vector< int > & bunches) const
privateinherited

Get the apparent bunch spacing in the current configuration.

This function tries to figure out what's the right bunch spacing for the current configuration, before the code would try to properly interpret the configuration.

This allows us to set a "MaxBunchSpacing" property of a relatively large number, but still be able to interpret configurations with very small gaps between bunch trains. (Like the 8b4e 2017 configurations.)

Parameters
bunchesThe filled bunches in the current configuration
Returns
The apparent bunch spacing of the configuration in bunch crossings

Definition at line 1225 of file BunchCrossingToolBase.cxx.

1226 {
1227
1228 // The maximum BC spacing to start from.
1229 const int maxSpacing = m_maxBunchSpacing / BunchCrossing::BUNCH_SPACING;
1230
1231 // Iterate downwards from the maximum spacing, searching for the minimum
1232 // spacing with which bunches exist.
1233 int result = maxSpacing;
1234 for( ; result > 0; --result ) {
1235
1236 // Test how many bunches have neighbors inside of the current window.
1237 const int nbunches =
1238 std::count_if( bunches.begin(), bunches.end(),
1239 count_bunch_neighbors( bunches, result ) );
1240 ATH_MSG_VERBOSE( "Number of bunches with " << result
1241 << " BC neighbors: " << nbunches );
1242
1243 // If none of them do, then we're finished.
1244 if( ! nbunches ) {
1245 // If we're not at the max spacing anymore, then it means that in
1246 // the previous step there were still bunches with neighbors in this
1247 // window. That's the spacing we need then.
1248 //
1249 // But if we're still at the maximum spacing, then let's just keep
1250 // that as the final answer.
1251 if( result != maxSpacing ) {
1252 ++result;
1253 }
1254 break;
1255 }
1256 }
1257 // If we went down to 0, then the right answer is 1. This is just how this
1258 // algorithm works...
1259 if( result == 0 ) {
1260 result = 1;
1261 }
1262 ATH_MSG_DEBUG( "Bunch spacing: " << result << " BCs" );
1263
1264 // Return the smallest spacing found:
1265 return result;
1266 }
#define ATH_MSG_VERBOSE(x)
#define ATH_MSG_DEBUG(x)
int m_maxBunchSpacing
The maximum bunch spacing that the tool should consider.
static const int BUNCH_SPACING
Minimum spacing between the bunches, in nanoseconds.

◆ bunchTrainSpacing()

int Trig::BunchCrossingToolBase::bunchTrainSpacing ( BunchDistanceType type = NanoSec) const
virtualinherited

Get the bunch spacing in the trains.

Implements Trig::IBunchCrossingTool.

Definition at line 827 of file BunchCrossingToolBase.cxx.

827 {
828
829 // Check if there are bunch trains in the current configurations:
830 if( m_bunchTrains.size() ) {
831 switch( type ) {
832 case NanoSec:
833 return m_bunchTrains.begin()->spacing();
834 break;
835 case BunchCrossings:
836 return ( m_bunchTrains.begin()->spacing() /
838 break;
839 case FilledBunches:
840 ATH_MSG_ERROR( "Function should not be called with argument: "
841 "FilledBunches" );
842 return -1;
843 break;
844 default:
845 ATH_MSG_ERROR( "Function called with unknown argument: " << type );
846 return -1;
847 }
848 } else {
849 // Return -1 if there are no bunch trains in the configuration:
850 return -1;
851 }
852
853 ATH_MSG_FATAL( "The code should never reach this line. Check the code!" );
854 return -1;
855 }
#define ATH_MSG_FATAL(x)

◆ declareGaudiProperty()

Gaudi::Details::PropertyBase & AthCommonDataStore< AthCommonMsg< AlgTool > >::declareGaudiProperty ( Gaudi::Property< T, V, H > & hndl,
const SG::VarHandleKeyType &  )
inlineprivateinherited

specialization for handling Gaudi::Property<SG::VarHandleKey>

Definition at line 156 of file AthCommonDataStore.h.

158 {
160 hndl.value(),
161 hndl.documentation());
162
163 }

◆ declareProperty()

Gaudi::Details::PropertyBase & AthCommonDataStore< AthCommonMsg< AlgTool > >::declareProperty ( Gaudi::Property< T, V, H > & t)
inlineinherited

Definition at line 145 of file AthCommonDataStore.h.

145 {
146 typedef typename SG::HandleClassifier<T>::type htype;
148 }
Gaudi::Details::PropertyBase & declareGaudiProperty(Gaudi::Property< T, V, H > &hndl, const SG::VarHandleKeyType &)
specialization for handling Gaudi::Property<SG::VarHandleKey>

◆ detStore()

const ServiceHandle< StoreGateSvc > & AthCommonDataStore< AthCommonMsg< AlgTool > >::detStore ( ) const
inlineinherited

The standard StoreGateSvc/DetectorStore Returns (kind of) a pointer to the StoreGateSvc.

Definition at line 95 of file AthCommonDataStore.h.

◆ distanceFromFront()

int Trig::BunchCrossingToolBase::distanceFromFront ( bcid_type bcid,
BunchDistanceType type = NanoSec ) const
virtualinherited

The distance of the specific bunch crossing from the front of the train.

This is one of the most tricky functions of this class.

When the user wants to ask the distance of a bunch crossing from the front of its train in units of nano seconds or bunch crossings, the main part of the logic is propagated to the algebra defined for the BunchCrossing class.

But when the user wants to know the distance in terms of filled bunches, the code has to treat 3 different bunch train configurations:

  • First is when the bunch train doesn't spread across the "BCID turnover". The code in this case just has to see how many steps it is from bunchtrain->begin().
  • If the bunch train spreads across the "BCID turnover", and the bcid is after BCID==0, the code has to count the filled bunches from bunchtrain->train_front() to bunchtrain->end(), plus the filled bunches from bunchtrain->begin() to the bcid in question.
  • If the bunch train spreads across the "BCID turnover", and the bcid is "before" BCID==0, the code has to count the filled bunches from bunchtrain->train_front() to the bcid in question.
Parameters
bcidThe bcid that should be checked
typeThe type of the requested return value
Returns
The distance of the bcid in question from the front of its bunch train

Implements Trig::IBunchCrossingTool.

Definition at line 239 of file BunchCrossingToolBase.cxx.

240 {
241
242 // Look for this BCID in the list of bunch trains:
243 std::set< BunchTrain >::const_iterator itr = m_bunchTrains.begin();
244 std::set< BunchTrain >::const_iterator end = m_bunchTrains.end();
245 for( ; itr != end; ++itr ) {
246 BunchTrain::const_iterator element;
247 if( ( element = itr->find( bcid ) ) != itr->end() ) {
248 switch( type ) {
249
250 case NanoSec:
251 return BunchCrossing::BUNCH_SPACING * ( BunchCrossing( bcid ) -
252 *( itr->train_front() ) );
253 break;
254 case BunchCrossings:
255 return ( BunchCrossing( bcid ) - *( itr->train_front() ) );
256 break;
257 case FilledBunches:
258 if( *( itr->train_front() ) > *( itr->train_back() ) ) {
259 if( BunchCrossing( bcid ) <= *( itr->train_back() ) ) {
260 return ( std::count_if( itr->begin(), element,
261 std::bind( std::not_equal_to< BunchCrossing >(),
262 *element, std::placeholders::_1 ) ) +
263 std::count_if( itr->train_front(), itr->end(),
264 std::bind( std::not_equal_to< BunchCrossing >(),
265 *element, std::placeholders::_1 ) ) );
266 } else {
267 return std::count_if( itr->train_front(), element,
268 std::bind( std::not_equal_to< BunchCrossing >(),
269 *element, std::placeholders::_1 ) );
270 }
271 } else {
272 return std::count_if( itr->begin(), element,
273 std::bind( std::not_equal_to< BunchCrossing >(),
274 *element, std::placeholders::_1 ) );
275 }
276 break;
277 default:
278 ATH_MSG_ERROR( "BunchDistanceType not understood!" );
279 return -1;
280 }
281 }
282 }
283
284 // If the bunch crossing is not part of a train:
285 return -1;
286 }

◆ distanceFromTail()

int Trig::BunchCrossingToolBase::distanceFromTail ( bcid_type bcid,
BunchDistanceType type = NanoSec ) const
virtualinherited

The distance of the specific bunch crossing from the tail of the train.

This is one of the most tricky functions of this class.

When the user wants to ask the distance of a bunch crossing from the tail of its train in units of nano seconds or bunch crossings, the main part of the logic is propagated to the algebra defined for the BunchCrossing class.

But when the user wants to know the distance in terms of filled bunches, the code has to treat 3 different bunch train configurations:

  • First is when the bunch train doesn't spread across the "BCID turnover". The code in this case just has to see how many steps it is from bunchtrain->end().
  • If the bunch train spreads across the "BCID turnover", and the bcid is "before" BCID==0, the code has to count the filled bunches from the bcid in question to bunchtrain->end(), plus the filled bunches from bunchtrain->begin() to bunchtrain->train_back().
  • If the bunch train spreads across the "BCID turnover", and the bcid is after BCID==0, the code has to count the filled bunches from the bcid in question to bunchtrain->train_back().
Parameters
bcidThe bcid that should be checked
typeThe type of the requested return value
Returns
The distance of the bcid in question from the tail of its bunch train

Implements Trig::IBunchCrossingTool.

Definition at line 316 of file BunchCrossingToolBase.cxx.

317 {
318
319 // Look for this BCID in the list of bunch trains:
320 std::set< BunchTrain >::const_iterator itr = m_bunchTrains.begin();
321 std::set< BunchTrain >::const_iterator end = m_bunchTrains.end();
322 for( ; itr != end; ++itr ) {
323 BunchTrain::const_iterator element;
324 if( ( element = itr->find( bcid ) ) != itr->end() ) {
325 switch( type ) {
326
327 case NanoSec:
328 return BunchCrossing::BUNCH_SPACING * ( *( itr->train_back() ) -
329 BunchCrossing( bcid ) );
330 break;
331 case BunchCrossings:
332 return ( *( itr->train_back() ) - BunchCrossing( bcid ) );
333 break;
334 case FilledBunches:
335 if( *( itr->train_front() ) > *( itr->train_back() ) ) {
336 if( BunchCrossing( bcid ) > *( itr->train_back() ) ) {
337 return ( std::count_if( element, itr->end(),
338 std::bind( std::not_equal_to< BunchCrossing >(),
339 *element, std::placeholders::_1 ) ) +
340 std::count_if( itr->begin(), ++( itr->train_back() ),
341 std::bind( std::not_equal_to< BunchCrossing >(),
342 *element, std::placeholders::_1 ) ) );
343 } else {
344 return std::count_if( element, ++( itr->train_back() ),
345 std::bind( std::not_equal_to< BunchCrossing >(),
346 *element, std::placeholders::_1 ) );
347 }
348 } else {
349 return std::count_if( element, itr->end(),
350 std::bind( std::not_equal_to< BunchCrossing >(),
351 *element, std::placeholders::_1 ) );
352 }
353 break;
354 default:
355 ATH_MSG_ERROR( "BunchDistanceType not understood!" );
356 return -1;
357 }
358 }
359 }
360
361 // If the bunch crossing is not part of a train:
362 return -1;
363 }

◆ endInputFile()

StatusCode asg::AsgMetadataTool::endInputFile ( )
protectedvirtualinherited

Function called when the currently open input file got completely processed.

Dummy implementation that can be overridden by the derived tool.

Reimplemented in BookkeeperDumperTool, BookkeeperTool, xAODMaker::TriggerMenuMetaDataTool, and xAODMaker::TruthMetaDataTool.

Definition at line 193 of file AsgMetadataTool.cxx.

193 {
194
195 // Return gracefully:
196 return StatusCode::SUCCESS;
197 }

◆ evtStore()

ServiceHandle< StoreGateSvc > & AthCommonDataStore< AthCommonMsg< AlgTool > >::evtStore ( )
inlineinherited

The standard StoreGateSvc (event store) Returns (kind of) a pointer to the StoreGateSvc.

Definition at line 85 of file AthCommonDataStore.h.

◆ extraDeps_update_handler()

void AthCommonDataStore< AthCommonMsg< AlgTool > >::extraDeps_update_handler ( Gaudi::Details::PropertyBase & ExtraDeps)
protectedinherited

Add StoreName to extra input/output deps as needed.

use the logic of the VarHandleKey to parse the DataObjID keys supplied via the ExtraInputs and ExtraOuputs Properties to add the StoreName if it's not explicitly given

◆ gapAfterBunch()

int Trig::BunchCrossingToolBase::gapAfterBunch ( bcid_type bcid,
BunchDistanceType dtype = NanoSec,
BunchFillType ftype = CollidingBunch ) const
virtualinherited

Gap after a particular bunch.

The function creates a smart BunchCrossing out of the BCID provided, then it goes looking for the next bunch crossing of the specified type.

Finally it just calculates the distance between the two bunch crossings in the requested units. It can return the size of the gap either in nanoseconds or in BCIDs (25 ns steps).

Parameters
bcidThe bcid that should be investigated
dtypeThe type of the requested return value
ftypeThe type of the previous bunch to consider
Returns
The gap after the specified bcid

Implements Trig::IBunchCrossingTool.

Definition at line 608 of file BunchCrossingToolBase.cxx.

610 {
611
612 // Construct this "smart" BCID:
613 const BunchCrossing bunch( bcid );
614
615 // Search for the first next bunch that fulfills the requirement:
616 BunchCrossing next_bunch( bunch );
617 ++next_bunch;
618 int loop_counter = 0;
619 switch( ftype ) {
620
621 case CollidingBunch:
622 // There should always be filled bunches in the configuration, but
623 // let's make sure that we don't get into an endless loop:
624 while( ( ! isFilled( next_bunch ) ) &&
625 ( loop_counter < BunchCrossing::MAX_BCID ) ) {
626 ++next_bunch;
627 ++loop_counter;
628 }
629 if( loop_counter == BunchCrossing::MAX_BCID ) {
630 ATH_MSG_ERROR( "Failed to calculate gap after BCID "
631 << bcid << " to a filled bunch! This shouldn't have "
632 << "happened!" );
633 return -1;
634 }
635 break;
636 case UnpairedBunch:
637 // There are no unpaired bunches in every configuration, so make sure
638 // we don't get into an endless loop:
639 while( ( ! isUnpaired( next_bunch ) ) &&
640 ( loop_counter < BunchCrossing::MAX_BCID ) ) {
641 ++next_bunch;
642 ++loop_counter;
643 }
644 // Return "-1" if there are no unpaired bunches in the configuration:
645 if( loop_counter == BunchCrossing::MAX_BCID ) {
646 return -1;
647 }
648 break;
649 case UnpairedBeam1:
650 // There are no unpaired bunches from beam 1 in every configuration, so
651 // make sure we don't get into an endless loop:
652 while( ( ! ( isUnpaired( next_bunch ) && isBeam1( next_bunch ) ) ) &&
653 ( loop_counter < BunchCrossing::MAX_BCID ) ) {
654 ++next_bunch;
655 ++loop_counter;
656 }
657 // Return "-1" if there are no unpaired bunches from beam 1 in the
658 // configuration:
659 if( loop_counter == BunchCrossing::MAX_BCID ) {
660 return -1;
661 }
662 break;
663 case UnpairedBeam2:
664 // There are no unpaired bunches from beam 2 in every configuration, so
665 // make sure we don't get into an endless loop:
666 while( ( ! ( isUnpaired( next_bunch ) && isBeam2( next_bunch ) ) ) &&
667 ( loop_counter < BunchCrossing::MAX_BCID ) ) {
668 ++next_bunch;
669 ++loop_counter;
670 }
671 // Return "-1" if there are no unpaired bunches from beam 2 in the
672 // configuration:
673 if( loop_counter == BunchCrossing::MAX_BCID ) {
674 return -1;
675 }
676 break;
677 case EmptyBunch:
678 // There should always be empty bunches in the configuration, but let's
679 // make sure that we don't get into an endless loop:
680 while( ( isFilled( next_bunch ) || isUnpaired( next_bunch ) ) &&
681 ( loop_counter < BunchCrossing::MAX_BCID ) ) {
682 ++next_bunch;
683 ++loop_counter;
684 }
685 if( loop_counter == BunchCrossing::MAX_BCID ) {
686 ATH_MSG_ERROR( "Failed to calculate gap after BCID "
687 << bcid << " to an empty bunch! This shouldn't have "
688 << "happened!" );
689 return -1;
690 }
691 break;
692 default:
693 ATH_MSG_ERROR( "Unknown bunch fill type specified: "
694 << ftype );
695 return -1;
696 }
697
698 // Now return the results:
699 switch( dtype ) {
700
701 case NanoSec:
702 return BunchCrossing::BUNCH_SPACING * bunch.gapTo( next_bunch );
703 break;
704 case BunchCrossings:
705 return bunch.gapTo( next_bunch );
706 break;
707 default:
708 ATH_MSG_ERROR( "You can only use NanoSec or BunchCrossings for type "
709 "for gapBeforeBunch" );
710 return -1;
711 }
712
713 // This should actually never be reached:
714 return -1;
715 }
virtual bool isBeam1(bcid_type bcid) const
Function deciding if there was a bunch from "beam 1" in this bunch crossing.
virtual bool isBeam2(bcid_type bcid) const
Function deciding if there was a bunch from "beam 2" in this bunch crossing.
static const int MAX_BCID
The maximum number of bunches that can be in the LHC.

◆ gapAfterTrain()

int Trig::BunchCrossingToolBase::gapAfterTrain ( bcid_type bcid,
BunchDistanceType type = NanoSec ) const
virtualinherited

Gap after the train this BCID is in.

The function first finds the train that the specified BCID is in.

Then it calculates the size of the gap after this train. It can return the size of the gap either in nanoseconds or in BCIDs (25 ns steps).

Parameters
bcidThe bcid whose train should be investigated
typeThe type of the requested return value
Returns
The gap after the train of the specified bcid

Implements Trig::IBunchCrossingTool.

Definition at line 429 of file BunchCrossingToolBase.cxx.

430 {
431
432 // Find this BCID in the list of bunch trains:
433 std::set< BunchTrain >::const_iterator itr = m_bunchTrains.begin();
434 std::set< BunchTrain >::const_iterator end = m_bunchTrains.end();
435 std::set< BunchTrain >::const_iterator train = m_bunchTrains.end();
436 for( ; itr != end; ++itr ) {
437 if( itr->find( bcid ) != itr->end() ) {
438 train = itr;
439 break;
440 }
441 }
442
443 // If we didn't find this BCID in a train, let's return right away:
444 if( train == end ) {
445 return -1;
446 }
447
448 // Search for the first filled bunch before the front of this train:
449 BunchCrossing train_front = *( train->train_back() );
450 ++train_front;
451 while( ! isFilled( train_front ) ) {
452 ++train_front;
453 }
454
455 // Now return the results:
456 switch( type ) {
457
458 case NanoSec:
460 train_front.distance( *( train->train_back() ) );
461 break;
462 case BunchCrossings:
463 return train_front.distance( *( train->train_back() ) );
464 break;
465 default:
466 ATH_MSG_ERROR( "You can only use NanoSec or BunchCrossings for type "
467 "for gapAfterTrain" );
468 return -1;
469 }
470
471 // This should actually never be reached:
472 return -1;
473 }

◆ gapBeforeBunch()

int Trig::BunchCrossingToolBase::gapBeforeBunch ( bcid_type bcid,
BunchDistanceType dtype = NanoSec,
BunchFillType ftype = CollidingBunch ) const
virtualinherited

Gap before a particular bunch.

The function creates a smart BunchCrossing out of the BCID provided, then it goes looking for the previous bunch crossing of the specified type.

Finally it just calculates the distance between the two bunch crossings in the requested units. It can return the size of the gap either in nanoseconds or in BCIDs (25 ns steps).

Parameters
bcidThe bcid that should be investigated
dtypeThe type of the requested return value
ftypeThe type of the previous bunch to consider
Returns
The gap before the specified bcid

Implements Trig::IBunchCrossingTool.

Definition at line 487 of file BunchCrossingToolBase.cxx.

489 {
490
491 // Construct this "smart" BCID:
492 const BunchCrossing bunch( bcid );
493
494 // Search for the first previous bunch that fulfills the requirement:
495 BunchCrossing prev_bunch( bunch );
496 --prev_bunch;
497 int loop_counter = 0;
498 switch( ftype ) {
499
500 case CollidingBunch:
501 // There should always be filled bunches in the configuration, but
502 // let's make sure that we don't get into an endless loop:
503 while( ( ! isFilled( prev_bunch ) ) &&
504 ( loop_counter < BunchCrossing::MAX_BCID ) ) {
505 --prev_bunch;
506 ++loop_counter;
507 }
508 if( loop_counter == BunchCrossing::MAX_BCID ) {
509 ATH_MSG_ERROR( "Failed to calculate gap before BCID "
510 << bcid << " to a filled bunch! This shouldn't have "
511 << "happened!" );
512 return -1;
513 }
514 break;
515 case UnpairedBunch:
516 // There are no unpaired bunches in every configuration, so make sure
517 // we don't get into an endless loop:
518 while( ( ! isUnpaired( prev_bunch ) ) &&
519 ( loop_counter < BunchCrossing::MAX_BCID ) ) {
520 --prev_bunch;
521 ++loop_counter;
522 }
523 // Return "-1" if there are no unpaired bunches in the configuration:
524 if( loop_counter == BunchCrossing::MAX_BCID ) {
525 return -1;
526 }
527 break;
528 case UnpairedBeam1:
529 // There are no unpaired bunches from beam 1 in every configuration, so
530 // make sure we don't get into an endless loop:
531 while( ( ! ( isUnpaired( prev_bunch ) && isBeam1( prev_bunch ) ) ) &&
532 ( loop_counter < BunchCrossing::MAX_BCID ) ) {
533 --prev_bunch;
534 ++loop_counter;
535 }
536 // Return "-1" if there are no unpaired bunches from beam 1 in the
537 // configuration:
538 if( loop_counter == BunchCrossing::MAX_BCID ) {
539 return -1;
540 }
541 break;
542 case UnpairedBeam2:
543 // There are no unpaired bunches from beam 2 in every configuration, so
544 // make sure we don't get into an endless loop:
545 while( ( ! ( isUnpaired( prev_bunch ) && isBeam2( prev_bunch ) ) ) &&
546 ( loop_counter < BunchCrossing::MAX_BCID ) ) {
547 --prev_bunch;
548 ++loop_counter;
549 }
550 // Return "-1" if there are no unpaired bunches from beam 2 in the
551 // configuration:
552 if( loop_counter == BunchCrossing::MAX_BCID ) {
553 return -1;
554 }
555 break;
556 case EmptyBunch:
557 // There should always be empty bunches in the configuration, but let's
558 // make sure that we don't get into an endless loop:
559 while( ( isFilled( prev_bunch ) || isUnpaired( prev_bunch ) ) &&
560 ( loop_counter < BunchCrossing::MAX_BCID ) ) {
561 --prev_bunch;
562 ++loop_counter;
563 }
564 if( loop_counter == BunchCrossing::MAX_BCID ) {
565 ATH_MSG_ERROR( "Failed to calculate gap before BCID "
566 << bcid << " to an empty bunch! This shouldn't have "
567 << "happened!" );
568 return -1;
569 }
570 break;
571 default:
572 ATH_MSG_ERROR( "Unknown bunch fill type specified: "
573 << ftype );
574 return -1;
575 }
576
577 // Now return the results:
578 switch( dtype ) {
579
580 case NanoSec:
581 return BunchCrossing::BUNCH_SPACING * bunch.gapFrom( prev_bunch );
582 break;
583 case BunchCrossings:
584 return bunch.gapFrom( prev_bunch );
585 break;
586 default:
587 ATH_MSG_ERROR( "You can only use NanoSec or BunchCrossings for type "
588 "for gapBeforeBunch" );
589 return -1;
590 }
591
592 // This should actually never be reached:
593 return -1;
594 }

◆ gapBeforeTrain()

int Trig::BunchCrossingToolBase::gapBeforeTrain ( bcid_type bcid,
BunchDistanceType type = NanoSec ) const
virtualinherited

Gap before the train this BCID is in.

The function first finds the train that the specified BCID is in.

Then it calculates the size of the gap before this train. It can return the size of the gap either in nanoseconds or in BCIDs (25 ns steps).

Parameters
bcidThe bcid whose train should be investigated
typeThe type of the requested return value
Returns
The gap before the train of the specified bcid

Implements Trig::IBunchCrossingTool.

Definition at line 374 of file BunchCrossingToolBase.cxx.

375 {
376
377 // Find this BCID in the list of bunch trains:
378 std::set< BunchTrain >::const_iterator itr = m_bunchTrains.begin();
379 std::set< BunchTrain >::const_iterator end = m_bunchTrains.end();
380 std::set< BunchTrain >::const_iterator train = m_bunchTrains.end();
381 for( ; itr != end; ++itr ) {
382 if( itr->find( bcid ) != itr->end() ) {
383 train = itr;
384 break;
385 }
386 }
387
388 // If we didn't find this BCID in a train, let's return right away:
389 if( train == end ) {
390 return -1;
391 }
392
393 // Search for the first filled bunch before the front of this train:
394 BunchCrossing train_tail = *( train->train_front() );
395 --train_tail;
396 while( ! isFilled( train_tail ) ) {
397 --train_tail;
398 }
399
400 // Now return the results:
401 switch( type ) {
402
403 case NanoSec:
405 train->train_front()->distance( train_tail );
406 break;
407 case BunchCrossings:
408 return train->train_front()->distance( train_tail );
409 break;
410 default:
411 ATH_MSG_ERROR( "You can only use NanoSec or BunchCrossings for type "
412 "for gapBeforeTrain" );
413 return -1;
414 }
415
416 // This should actually never be reached:
417 return -1;
418 }

◆ getKey()

SG::sgkey_t asg::AsgTool::getKey ( const void * ptr) const
inherited

Get the (hashed) key of an object that is in the event store.

This is a bit of a special one. StoreGateSvc and xAOD::TEvent both provide ways for getting the SG::sgkey_t key for an object that is in the store, based on a bare pointer. But they provide different interfaces for doing so.

In order to allow tools to efficiently perform this operation, they can use this helper function.

See also
asg::AsgTool::getName
Parameters
ptrThe bare pointer to the object that the event store should know about
Returns
The hashed key of the object in the store. If not found, an invalid (zero) key.

Definition at line 119 of file AsgTool.cxx.

119 {
120
121#ifdef XAOD_STANDALONE
122 // In case we use @c xAOD::TEvent, we have a direct function call
123 // for this.
124 return evtStore()->event()->getKey( ptr );
125#else
126 const SG::DataProxy* proxy = evtStore()->proxy( ptr );
127 return ( proxy == nullptr ? 0 : proxy->sgkey() );
128#endif // XAOD_STANDALONE
129 }
ServiceHandle< StoreGateSvc > & evtStore()

◆ getName()

const std::string & asg::AsgTool::getName ( const void * ptr) const
inherited

Get the name of an object that is / should be in the event store.

This is a bit of a special one. StoreGateSvc and xAOD::TEvent both provide ways for getting the std::string name for an object that is in the store, based on a bare pointer. But they provide different interfaces for doing so.

In order to allow tools to efficiently perform this operation, they can use this helper function.

See also
asg::AsgTool::getKey
Parameters
ptrThe bare pointer to the object that the event store should know about
Returns
The string name of the object in the store. If not found, an empty string.

Definition at line 106 of file AsgTool.cxx.

106 {
107
108#ifdef XAOD_STANDALONE
109 // In case we use @c xAOD::TEvent, we have a direct function call
110 // for this.
111 return evtStore()->event()->getName( ptr );
112#else
113 const SG::DataProxy* proxy = evtStore()->proxy( ptr );
114 static const std::string dummy = "";
115 return ( proxy == nullptr ? dummy : proxy->name() );
116#endif // XAOD_STANDALONE
117 }

◆ getProperty()

template<class T>
const T * asg::AsgTool::getProperty ( const std::string & name) const
inherited

Get one of the tool's properties.

◆ handle()

void asg::AsgMetadataTool::handle ( const Incident & inc)
protectedvirtualinherited

Function receiving incidents from IncidentSvc/TEvent.

Reimplemented in Trig::TrigDecisionTool.

Definition at line 135 of file AsgMetadataTool.cxx.

135 {
136
137 // Tell the user what's happening:
138 ATH_MSG_VERBOSE( "Callback received with incident: " << inc.type() );
139
140 // Call the appropriate member function:
141 if( inc.type() == IncidentType::BeginInputFile ) {
143 if( beginInputFile().isFailure() ) {
144 ATH_MSG_FATAL( "Failed to call beginInputFile()" );
145 throw std::runtime_error( "Couldn't call beginInputFile()" );
146 }
147 } else if( inc.type() == IncidentType::EndInputFile ) {
148 if( endInputFile().isFailure() ) {
149 ATH_MSG_FATAL( "Failed to call endInputFile()" );
150 throw std::runtime_error( "Couldn't call endInputFile()" );
151 }
152 } else if( inc.type() == IncidentType::BeginEvent ) {
153 // If the tool didn't catch the begin input file incident for the
154 // first input file of the job, then call the appropriate function
155 // now.
156 if( ! m_beginInputFileCalled ) {
158 if( beginInputFile().isFailure() ) {
159 ATH_MSG_FATAL( "Failed to call beginInputFile()" );
160 throw std::runtime_error( "Couldn't call beginInputFile()" );
161 }
162 }
163 if( beginEvent().isFailure() ) {
164 ATH_MSG_FATAL( "Failed to call beginEvent()" );
165 throw std::runtime_error( "Couldn't call beginEvent()" );
166 }
167
168 #ifdef XAOD_STANDALONE
169 } else if( inc.type() == IncidentType::MetaDataStop ) {
170 if( metaDataStop().isFailure() ) {
171 ATH_MSG_FATAL( "Failed to call metaDataStop()" );
172 throw std::runtime_error( "Couldn't call metaDataStop()" );
173 }
174
175 #endif // XAOD_STANDALONE
176 } else {
177 ATH_MSG_WARNING( "Unknown incident type received in AsgMetaDataTool: " << inc.type() );
178 }
179
180 return;
181 }
virtual StatusCode beginInputFile()
Function called when a new input file is opened.
virtual StatusCode beginEvent()
Function called when a new events is loaded.
bool m_beginInputFileCalled
Flag helping to discover when the tool misses the opening of the first input file.
virtual StatusCode endInputFile()
Function called when the currently open input file got completely processed.
virtual StatusCode metaDataStop()
Function called when the tool should write out its metadata.

◆ initialize()

StatusCode Trig::StaticBunchCrossingTool::initialize ( void )
virtual

Function initialising the tool.

Reimplemented from asg::AsgTool.

Definition at line 34 of file StaticBunchCrossingTool.cxx.

34 {
35
36 // Let the user know what's happening:
37 ATH_MSG_INFO( "Initialising tool" );
38 ATH_MSG_INFO( " Bunch group config : " << m_bgkey );
39
40 // Load the configuration specified by the tool's property:
42
43 // Return gracefully:
44 return StatusCode::SUCCESS;
45 }
#define ATH_CHECK
Evaluate an expression and check for errors.
#define ATH_MSG_INFO(x)
StatusCode loadConfig(int bgkey)
Load a hard-coded bunch group key.

◆ inputHandles()

virtual std::vector< Gaudi::DataHandle * > AthCommonDataStore< AthCommonMsg< AlgTool > >::inputHandles ( ) const
overridevirtualinherited

Return this algorithm's input handles.

We override this to include handle instances from key arrays if they have not yet been declared. See comments on updateVHKA.

◆ inputMetaStore()

AsgMetadataTool::MetaStorePtr_t asg::AsgMetadataTool::inputMetaStore ( ) const
inherited

Accessor for the input metadata store.

Definition at line 93 of file AsgMetadataTool.cxx.

93 {
94
95#ifdef XAOD_STANDALONE
96 return &m_inputMetaStore;
97#else // XAOD_STANDALONE
98 return m_inputMetaStore;
99#endif // XAOD_STANDALONE
100 }
MetaStore_t m_inputMetaStore
Object accessing the input metadata store.

◆ isBeam1()

bool Trig::BunchCrossingToolBase::isBeam1 ( bcid_type bcid) const
virtualinherited

Function deciding if there was a bunch from "beam 1" in this bunch crossing.

Implements Trig::IBunchCrossingTool.

Definition at line 66 of file BunchCrossingToolBase.cxx.

66 {
67
68 // Check if this is a filled bunch:
69 std::set< BunchCrossing >::const_iterator itr_fill =
70 m_filledBunches.find( bcid );
71 if( itr_fill != m_filledBunches.end() ) {
72 return true;
73 }
74
75 // Check if this is an unpaired bunch with the bunch in beam 1:
76 std::set< BunchCrossing >::const_iterator itr_unp =
77 m_unpairedBunches.find( bcid );
78 if( ( itr_unp != m_unpairedBunches.end() ) &&
79 ( itr_unp->intensityBeam1() > 0.001 ) ) {
80 return true;
81 } else {
82 return false;
83 }
84 }

◆ isBeam2()

bool Trig::BunchCrossingToolBase::isBeam2 ( bcid_type bcid) const
virtualinherited

Function deciding if there was a bunch from "beam 2" in this bunch crossing.

Implements Trig::IBunchCrossingTool.

Definition at line 86 of file BunchCrossingToolBase.cxx.

86 {
87
88 // Check if this is a filled bunch:
89 std::set< BunchCrossing >::const_iterator itr_fill =
90 m_filledBunches.find( bcid );
91 if( itr_fill != m_filledBunches.end() ) {
92 return true;
93 }
94
95 // Check if this is an unpaired bunch with the bunch in beam 1:
96 std::set< BunchCrossing >::const_iterator itr_unp =
97 m_unpairedBunches.find( bcid );
98 if( ( itr_unp != m_unpairedBunches.end() ) &&
99 ( itr_unp->intensityBeam2() > 0.001 ) ) {
100 return true;
101 } else {
102 return false;
103 }
104 }

◆ isFilled()

bool Trig::BunchCrossingToolBase::isFilled ( bcid_type bcid) const
virtualinherited

The simplest query: Is the bunch crossing filled or not?

Implements Trig::IBunchCrossingTool.

Definition at line 32 of file BunchCrossingToolBase.cxx.

32 {
33
34 // Check if this is a filled bunch crossing:
35 if( m_filledBunches.find( bcid ) != m_filledBunches.end() ) {
36 return true;
37 } else {
38 return false;
39 }
40 }

◆ isInTrain()

bool Trig::BunchCrossingToolBase::isInTrain ( bcid_type bcid) const
virtualinherited

Function deciding if a given bunch crossing is in a filled train.

Implements Trig::IBunchCrossingTool.

Definition at line 42 of file BunchCrossingToolBase.cxx.

42 {
43
44 // Look for this BCID in the list of bunch trains:
45 std::set< BunchTrain >::const_iterator itr = m_bunchTrains.begin();
46 std::set< BunchTrain >::const_iterator end = m_bunchTrains.end();
47 for( ; itr != end; ++itr ) {
48 if( itr->find( bcid ) != itr->end() ) {
49 return true;
50 }
51 }
52
53 return false;
54 }

◆ isUnpaired()

bool Trig::BunchCrossingToolBase::isUnpaired ( bcid_type bcid) const
virtualinherited

Function deciding if a given bunch crossing has an unpaired bunch.

Implements Trig::IBunchCrossingTool.

Definition at line 56 of file BunchCrossingToolBase.cxx.

56 {
57
58 // Check if this is an unpaired bunch crossing:
59 if( m_unpairedBunches.find( bcid ) != m_unpairedBunches.end() ) {
60 return true;
61 } else {
62 return false;
63 }
64 }

◆ loadBunchTrains()

StatusCode Trig::BunchCrossingToolBase::loadBunchTrains ( const std::vector< int > & bunches,
const std::vector< float > & bunch_int1 = std::vector< float >(),
const std::vector< float > & bunch_int2 = std::vector< float >() )
protectedinherited

Interpret the configuration for bunch trains.

This function takes care of identifying the bunch trains in the configuration.

The algorithm is quite simple. It starts off with all the filled bunches that have not been identified as single bunches by the loadSingleBunches(...) function. It takes the first available bunch crossing, then loops over the rest of them. When it finds a BC that's "close enough" to the current bunch train, then the bunch train is extended. From there on the algorithm continues with this extended bunch train. When the loop reaches the end of the available bunches, the created bunch train is added to the cache, and the algorithms starts again with the first still available bunch crossing.

The bunch intensity parameter is optional. If it's not specified, the code will assign an intensity of "1.0" to all the bunch crossings.

Parameters
bunchesThe filled bunch crossings
bunch_intThe "intensities" of the paired bunches
Returns
StatusCode::SUCCESS if the interpretation was successful, StatusCode::FAILURE otherwise

Definition at line 986 of file BunchCrossingToolBase.cxx.

989 {
990
991 // Do a small sanity check:
992 if( ( ( bunches.size() != bunch_int1.size() ) && bunch_int1.size() ) ||
993 ( ( bunches.size() != bunch_int2.size() ) && bunch_int2.size() ) ) {
994 ATH_MSG_ERROR( "Received vectors of different sizes for the bunch "
995 "IDs and bunch intensities\n"
996 "Function can not work like this..." );
997 return StatusCode::FAILURE;
998 }
999 if( ! bunch_int1.size() ) {
1000 ATH_MSG_DEBUG( "Not using bunch intensity for the calculation" );
1001 }
1002 if( ( ! bunch_int2.size() ) && bunch_int1.size() ) {
1003 ATH_MSG_DEBUG( "Using 'bunch crossing intensity' for the "
1004 "calculation" );
1005 }
1006
1007 //
1008 // Calculate the allowed maximum BCID separation between single bunches:
1009 //
1010 const int maxBCSpacing = bunchSpacing( bunches );
1011
1012 // Reset the cache:
1013 m_bunchTrains.clear();
1014
1015 //
1016 // Create a cache of the bunches which have not been identified as a
1017 // single bunch:
1018 //
1019 std::set< BunchCrossing > cache;
1020 std::vector< int >::const_iterator b_itr = bunches.begin();
1021 std::vector< int >::const_iterator b_end = bunches.end();
1022 std::vector< float >::const_iterator i1_itr = bunch_int1.begin();
1023 std::vector< float >::const_iterator i2_itr = bunch_int2.begin();
1024 for( ; b_itr != b_end; ++b_itr ) {
1025 if( std::find( m_singleBunches.begin(), m_singleBunches.end(),
1026 BunchCrossing( *b_itr ) ) == m_singleBunches.end() ) {
1027 cache.insert( BunchCrossing( *b_itr,
1028 ( float )( bunch_int1.size() ? *i1_itr :
1029 1.0 ),
1030 ( float )( bunch_int2.size() ? *i2_itr :
1031 0.0 ) ) );
1032 }
1033 if( bunch_int1.size() ) ++i1_itr;
1034 if( bunch_int2.size() ) ++i2_itr;
1035 }
1036
1037 ATH_MSG_VERBOSE( "Bunches considered for trains: " << cache );
1038
1039 //
1040 // Continue the loop until we have unassigned bunches:
1041 //
1042 while( cache.size() ) {
1043
1044 // Create a new bunch train object:
1045 BunchTrain bt;
1046 bt.insert( *cache.begin() );
1047
1048 // Try finding attachable bunches until no other attachable bunch can
1049 // be found:
1050 size_t prev_size = 0;
1051 while( prev_size != cache.size() ) {
1052
1053 // Let's remember the size of the cache:
1054 prev_size = cache.size();
1055
1056 // Find all the bunches that should be a part of this train:
1057 std::set< BunchCrossing >::const_iterator c_itr = cache.begin();
1058 std::set< BunchCrossing >::const_iterator c_end = cache.end();
1059 for( ; c_itr != c_end; ++c_itr ) {
1060 if( bt.distance( *c_itr ) <= maxBCSpacing ) {
1061 ATH_MSG_VERBOSE( "Adding BunchCrossing " << *c_itr
1062 << " to Bunch Train " << bt );
1063 bt.insert( *c_itr );
1064 }
1065 }
1066
1067 // Now remove the selected bunches from the cache:
1068 BunchTrain::const_iterator itr = bt.begin();
1069 BunchTrain::const_iterator end = bt.end();
1070 for( ; itr != end; ++itr ) {
1071 cache.erase( *itr );
1072 }
1073 }
1074
1075 // Finally, remember this train:
1076 if( ! bt.validate() ) {
1077 ATH_MSG_ERROR( "Found a strange bunch train: " << bt );
1078 ATH_MSG_ERROR( "Keeping in it the list of trains!" );
1079 }
1080 m_bunchTrains.insert( bt );
1081 }
1082
1083 //
1084 // Check if the spacing in the bunch trains is the same. (It should be for
1085 // all real configurations.)
1086 //
1087 std::set< BunchTrain >::const_iterator train_itr = m_bunchTrains.begin();
1088 std::set< BunchTrain >::const_iterator train_end = m_bunchTrains.end();
1089 int spacing = -1;
1090 for( ; train_itr != train_end; ++train_itr ) {
1091 if( spacing < 0 ) {
1092 spacing = train_itr->spacing();
1093 continue;
1094 }
1095 if( train_itr->spacing() != spacing ) {
1096 ATH_MSG_WARNING( "The spacing seems to be different between the "
1097 "trains" );
1098 ATH_MSG_WARNING( "This should probably not happen" );
1099 }
1100 }
1101
1102 ATH_MSG_DEBUG( "Bunch trains found: " << m_bunchTrains );
1103
1104 return StatusCode::SUCCESS;
1105 }
std::set< Trig::BunchCrossing > m_singleBunches
Internal list of single bunches.
int bunchSpacing(const std::vector< int > &bunches) const
Get the apparent bunch spacing in the current configuration.

◆ loadConfig() [1/3]

StatusCode Trig::StaticBunchCrossingTool::loadConfig ( const std::vector< float > & bunches)

Load the configuration specified by the user.

This helper function is here to make it possible to extract the bunch configuration from an MC file in an Athena environment, and then use that as a hardcoded configuration to initialise this tool.

In this case the argument is the intensities of the bunches. From which the code uses the same logic as MCBunchCrossingTool to find the filled bunches.

Parameters
bunchesThe intensities of the BCIDs. Starting with BCID 0.
Returns
StatusCode::SUCCESS if loading the configuration was successful, StatusCode::FAILURE otherwise

Definition at line 128 of file StaticBunchCrossingTool.cxx.

129 {
130
131 // Check that the user specified something useful:
132 if( ! bunches.size() ) {
133 ATH_MSG_ERROR( "Empty container received" );
134 return StatusCode::FAILURE;
135 }
136
137 // Translate it into vectors of filled bunches and intensities:
138 std::vector< int > filled_bunches;
139 std::vector< float > filled_intensities;
140
141 // Minimum bunch intensity to consider a bunch filled:
142 static const float MIN_BUNCH_INTENSITY = 0.1;
143
144 // Check if the pattern "fits into" the LHC:
145 if( BunchCrossing::MAX_BCID % bunches.size() ) {
146
147 ATH_MSG_INFO( "Bunch pattern doesn't \"fit into\" "
149 // The loop doesn't go all the way up to MAX_BCID/2 in order not
150 // to produce "weird" patterns half way. This should be pretty safe
151 // to do, because the MC BCIDs will only be in the range defined by
152 // the pattern from the metadata.
153 for( int i = 0; i < ( BunchCrossing::MAX_BCID / 2 - 20 ); ++i ) {
154 const int pos1 = i % bunches.size();
155 const int pos2 = bunches.size() - 1 - ( i % bunches.size() );
156 if( bunches[ pos1 ] > MIN_BUNCH_INTENSITY ) {
157 filled_bunches.push_back( i );
158 filled_intensities.push_back( bunches[ pos1 ] );
159 }
160 if( bunches[ pos2 ] > MIN_BUNCH_INTENSITY ) {
161 filled_bunches.push_back( BunchCrossing::MAX_BCID - 1 - i );
162 filled_intensities.push_back( bunches[ pos2 ] );
163 }
164 }
165
166 } else {
167
168 // If the sample size fits into the number of available bunches,
169 // the algorithm is pretty simple:
170 ATH_MSG_INFO( "Bunch pattern \"fits into\" "
172 for( int i = 0; i < BunchCrossing::MAX_BCID; ++i ) {
173 const int pos = i % bunches.size();
174 if( bunches[ pos ] > MIN_BUNCH_INTENSITY ) {
175 filled_bunches.push_back( i );
176 filled_intensities.push_back( bunches[ pos ] );
177 }
178 }
179
180 }
181
182 //
183 // Now let the base class interpret the information:
184 //
185 ATH_CHECK( loadSingleBunches( filled_bunches, filled_intensities ) );
186 ATH_CHECK( loadBunchTrains( filled_bunches, filled_intensities ) );
187 ATH_CHECK( loadUnpairedBunches( std::vector< int >(),
188 std::vector< int >() ) );
189
190 // Print the configuration to give some feedback to the user:
191 printConfig();
192
193 // Return gracefully:
194 return StatusCode::SUCCESS;
195 }
StatusCode loadBunchTrains(const std::vector< int > &bunches, const std::vector< float > &bunch_int1=std::vector< float >(), const std::vector< float > &bunch_int2=std::vector< float >())
Interpret the configuration for bunch trains.
StatusCode loadSingleBunches(const std::vector< int > &bunches, const std::vector< float > &bunch_int1=std::vector< float >(), const std::vector< float > &bunch_int2=std::vector< float >())
Interpret the configuration for single bunches.
void printConfig() const
Function printing the configuration of the tool.
StatusCode loadUnpairedBunches(const std::vector< int > &beam1, const std::vector< int > &beam2, const std::vector< float > &bunch_int1=std::vector< float >(), const std::vector< float > &bunch_int2=std::vector< float >())
Interpret the configuration for unpaired bunches.

◆ loadConfig() [2/3]

StatusCode Trig::StaticBunchCrossingTool::loadConfig ( const std::vector< int > & filledBunches,
const std::vector< float > & filledIntensities = std::vector< float >(),
const std::vector< int > & unpairedBunches = std::vector< int >(),
const std::vector< float > & unpairedIntensities = std::vector< float >() )

Load a configuration specified by the user.

To make it possible to quickly try out new configurations, this public function configures the underlying tool from a list of bunches (and optionally bunch intensities) given by the user.

Only used for testing at the moment.

Parameters
bunchesThe filled bunches which should be interpreted
intensitiesThe intensities of these filled bunches (optional)
Returns
StatusCode::SUCCESS if loading the configuration was successful, StatusCode::FAILURE otherwise

Definition at line 94 of file StaticBunchCrossingTool.cxx.

98 {
99
100 //
101 // Let the base class interpret the configuration:
102 //
103 ATH_CHECK( loadSingleBunches( filledBunches, filledIntensities ) );
104 ATH_CHECK( loadBunchTrains( filledBunches, filledIntensities ) );
105 ATH_CHECK( loadUnpairedBunches( unpairedBunches, unpairedBunches,
106 unpairedIntensities,
107 unpairedIntensities ) );
108
109 // Print the configuration to give some feedback to the user:
110 printConfig();
111
112 return StatusCode::SUCCESS;
113 }

◆ loadConfig() [3/3]

StatusCode Trig::StaticBunchCrossingTool::loadConfig ( int bgkey)

Load a hard-coded bunch group key.

This function tries to load one of the pre-defined configurations.

The code knows about all the bunch-group settings up to BGKey=105 at the moment.

Parameters
bgkeyBunch group key of the requested configuration
Returns
StatusCode::SUCCESS if loading the configuration was successful, StatusCode::FAILURE otherwise

Definition at line 56 of file StaticBunchCrossingTool.cxx.

56 {
57
58 //
59 // Check if we have this configuration:
60 //
61 std::map< int, std::vector< int > >::const_iterator bgset;
62 if( ( bgset = m_knownBGKeys.find( bgkey ) ) == m_knownBGKeys.end() ) {
63 ATH_MSG_ERROR( "Couldn't find configuration for bunch group key: "
64 << bgkey );
65 return StatusCode::FAILURE;
66 }
67
68 //
69 // Now let the base class interpret the configuration:
70 //
71 ATH_CHECK( loadSingleBunches( bgset->second ) );
72 ATH_CHECK( loadBunchTrains( bgset->second ) );
73 ATH_CHECK( loadUnpairedBunches( std::vector< int >(),
74 std::vector< int >() ) );
75
76 // Print the configuration to give some feedback to the user:
78
79 return StatusCode::SUCCESS;
80 }

◆ loadSingleBunches()

StatusCode Trig::BunchCrossingToolBase::loadSingleBunches ( const std::vector< int > & bunches,
const std::vector< float > & bunch_int1 = std::vector< float >(),
const std::vector< float > & bunch_int2 = std::vector< float >() )
protectedinherited

Interpret the configuration for single bunches.

This function takes care of selecting the single bunches from the paired bunches.

These will then not be taken into account in the train finding algorithm.

The bunch intensity parameter is optional. If it's not specified, the code will assign an intensity of "1.0" to all the bunch crossings.

Parameters
bunchesThe paired bunches
bunch_intThe "intensities" of the paired bunches
Returns
StatusCode::SUCCESS if the interpretation was successful, StatusCode::FAILURE otherwise

Definition at line 870 of file BunchCrossingToolBase.cxx.

873 {
874
875 // Do a small sanity check:
876 if( ( ( bunches.size() != bunch_int1.size() ) && bunch_int1.size() ) ||
877 ( ( bunches.size() != bunch_int2.size() ) && bunch_int2.size() ) ) {
878 ATH_MSG_ERROR( "Received vectors of different sizes for the bunch "
879 "IDs and bunch intensities\n"
880 "Function can not work like this..." );
881 return StatusCode::FAILURE;
882 }
883 if( ! bunch_int1.size() ) {
884 ATH_MSG_DEBUG( "Not using bunch intensity for the calculation" );
885 }
886 if( ( ! bunch_int2.size() ) && bunch_int1.size() ) {
887 ATH_MSG_DEBUG( "Using 'bunch crossing intensity' for the "
888 "calculation" );
889 }
890
891 //
892 // Calculate the allowed maximum BCID separation between single bunches:
893 //
894 const int maxBCSpacing = bunchSpacing( bunches );
895
896 // Reset the cache:
897 m_filledBunches.clear();
898 m_singleBunches.clear();
899
900 //
901 // Loop over the paired bunch crossings:
902 //
903 std::vector< int >::const_iterator b_itr = bunches.begin();
904 std::vector< int >::const_iterator b_end = bunches.end();
905 std::vector< float >::const_iterator i1_itr = bunch_int1.begin();
906 std::vector< float >::const_iterator i2_itr = bunch_int2.begin();
907 for( ; b_itr != b_end; ++b_itr ) {
908
909 // Evaluate the intensity of this paired bunch:
910 const float intensity1 = bunch_int1.size() ? *i1_itr : 1.0;
911 const float intensity2 = bunch_int2.size() ? *i2_itr : 0.0;
912
913 // It's definitely a filled bunch, so let's remember it as such:
914 m_filledBunches.insert( BunchCrossing( *b_itr, intensity1,
915 intensity2 ) );
916
917 ATH_MSG_VERBOSE( "Evaluating bunch crossing: " << *b_itr );
918
919 //
920 // This expression counts how many of the paired bunches fulfill:
921 //
922 // distance( ref_bcid, bcid ) <= maxBCSpacing
923 //
924 // Since the calculation takes the reference bcid into account as well,
925 // the count is always >= 1. I have to use the specialised
926 // distance(...) function, because a regular
927 //
928 // bcid - maxBCSpacing <= ref_bcid <= bcid + maxBCSpacing
929 //
930 // expression wouldn't give the correct answer at the "turnover" of the
931 // bcid numbering. (When evaluating bcid 1 and
932 // BunchCrossing::MAX_BCID.)
933 //
934 const int neighbours =
935 std::count_if( bunches.begin(), bunches.end(),
936 [ maxBCSpacing, &b_itr ]( int bunch ) {
937 return ( Trig::distance( bunch, *b_itr ) <=
938 maxBCSpacing );
939 } );
940
941 //
942 // Now decide if we want to consider this bunch crossing as a single
943 // bunch or not:
944 //
945 ATH_MSG_VERBOSE( " Bunch neighbours: " << neighbours );
946 if( neighbours == 1 ) {
947 ATH_MSG_VERBOSE( " Bunch crossing " << *b_itr
948 << " seems to be a single bunch" );
949 m_singleBunches.insert( BunchCrossing( *b_itr, intensity1,
950 intensity2 ) );
951 }
952
953 // Only step through the intensity vector(s) if it has some elements:
954 if( bunch_int1.size() ) ++i1_itr;
955 if( bunch_int2.size() ) ++i2_itr;
956 }
957
958 //
959 // Finally some debugging message for the end:
960 //
961 ATH_MSG_DEBUG( "Single bunches found: " << m_singleBunches );
962
963 return StatusCode::SUCCESS;
964 }

◆ loadUnpairedBunches()

StatusCode Trig::BunchCrossingToolBase::loadUnpairedBunches ( const std::vector< int > & beam1,
const std::vector< int > & beam2,
const std::vector< float > & bunch_int1 = std::vector< float >(),
const std::vector< float > & bunch_int2 = std::vector< float >() )
protectedinherited

Interpret the configuration for unpaired bunches.

This function just caches the unpaired bunches internally.

It doesn't have to do anything as fancy as the other two load functions, it just takes the BCIDs as they are.

Parameters
bunchesThe unpaired bunch crossings
bunch_intThe "intensities" of the unpaired bunches
Returns
StatusCode::SUCCESS if the caching was successful, StatusCode::FAILURE otherwise

Definition at line 1117 of file BunchCrossingToolBase.cxx.

1121 {
1122
1123 // Do a small sanity check:
1124 if( ( ( beam1.size() != bunch_int1.size() ) && bunch_int1.size() ) ||
1125 ( ( beam2.size() != bunch_int2.size() ) && bunch_int2.size() ) ) {
1126 ATH_MSG_ERROR( "Received vectors of different sizes for the bunch "
1127 "IDs and bunch intensities\n"
1128 "Function can not work like this..." );
1129 return StatusCode::FAILURE;
1130 }
1131 if( ( ! bunch_int1.size() ) && ( ! bunch_int2.size() ) ) {
1132 ATH_MSG_DEBUG( "Not using bunch intensity for the calculation" );
1133 }
1134
1135 // Reset the cache:
1136 m_unpairedBunches.clear();
1137
1138 //
1139 // Add the unpaired bunches from beam 1:
1140 //
1141 std::vector< int >::const_iterator b_itr = beam1.begin();
1142 std::vector< int >::const_iterator b_end = beam1.end();
1143 std::vector< float >::const_iterator i_itr = bunch_int1.begin();
1144 for( ; b_itr != b_end; ++b_itr ) {
1145
1146 // Evaluate the intensity of this unpaired bunch:
1147 const float intensity = bunch_int1.size() ? *i_itr : 1.0;
1148
1149 // Nothing fancy to do, just put it in the cache:
1150 m_unpairedBunches.insert( BunchCrossing( *b_itr, intensity, 0.0 ) );
1151
1152 // Only step through the intensity vector if it has some elements:
1153 if( bunch_int1.size() ) ++i_itr;
1154 }
1155
1156 //
1157 // Add the unpaired bunches from beam 2:
1158 //
1159 b_itr = beam2.begin();
1160 b_end = beam2.end();
1161 i_itr = bunch_int2.begin();
1162 for( ; b_itr != b_end; ++b_itr ) {
1163
1164 // Evaluate the intensity of this unpaired bunch:
1165 const float intensity = bunch_int2.size() ? *i_itr : 1.0;
1166
1167 // Check if this BCID is already known as an unpaired BCID:
1168 std::set< BunchCrossing >::iterator itr =
1169 m_unpairedBunches.find( *b_itr );
1170 if( itr != m_unpairedBunches.end() ) {
1171 // Modify the BCID not to correspond to a particular beam. Most
1172 // implementations don't treat beam 1 and beam 2 separately.
1173 BunchCrossing bc( *itr );
1174 bc.setIntensityBeam2( intensity );
1175 m_unpairedBunches.erase( itr );
1176 m_unpairedBunches.insert( bc );
1177 } else {
1178 // Nothing fancy to do, just put it in the cache:
1179 m_unpairedBunches.insert( BunchCrossing( *b_itr, 0.0, intensity ) );
1180 }
1181
1182 // Only step through the intensity vector if it has some elements:
1183 if( bunch_int2.size() ) ++i_itr;
1184 }
1185
1186 //
1187 // Finally some debugging message for the end:
1188 //
1189 ATH_MSG_DEBUG( "Unpaired bunches found: " << m_unpairedBunches );
1190
1191 return StatusCode::SUCCESS;
1192 }

◆ metaDataStop()

StatusCode asg::AsgMetadataTool::metaDataStop ( )
protectedvirtualinherited

Function called when the tool should write out its metadata.

Dummy implementation that can be overridden by the derived tool.

Reimplemented in BookkeeperDumperTool, BookkeeperTool, xAODMaker::TriggerMenuMetaDataTool, and xAODMaker::TruthMetaDataTool.

Definition at line 209 of file AsgMetadataTool.cxx.

209 {
210
211 // Return gracefully:
212 return StatusCode::SUCCESS;
213 }

◆ msg()

MsgStream & AthCommonMsg< AlgTool >::msg ( ) const
inlineinherited

Definition at line 24 of file AthCommonMsg.h.

24 {
25 return this->msgStream();
26 }

◆ msg_level_name()

const std::string & asg::AsgTool::msg_level_name ( ) const
inherited

A deprecated function for getting the message level's name.

Instead of using this, weirdly named function, user code should get the string name of the current minimum message level (in case they really need it...), with:

MSG::name( msg().level() )

This function's name doesn't follow the ATLAS coding rules, and as such will be removed in the not too distant future.

Returns
The string name of the current minimum message level that's printed

Definition at line 101 of file AsgTool.cxx.

101 {
102
103 return MSG::name( msg().level() );
104 }
MsgStream & msg() const
const std::string & name(Level lvl)
Convenience function for translating message levels to strings.
Definition MsgLevel.cxx:19

◆ msgLvl()

bool AthCommonMsg< AlgTool >::msgLvl ( const MSG::Level lvl) const
inlineinherited

Definition at line 30 of file AthCommonMsg.h.

30 {
31 return this->msgLevel(lvl);
32 }

◆ numberOfBunchTrains()

unsigned int Trig::BunchCrossingToolBase::numberOfBunchTrains ( ) const
virtualinherited

Get the number of the bunch trains in the current configuration.

Implements Trig::IBunchCrossingTool.

Definition at line 821 of file BunchCrossingToolBase.cxx.

821 {
822
823 return m_bunchTrains.size();
824 }

◆ numberOfFilledBunches()

unsigned int Trig::BunchCrossingToolBase::numberOfFilledBunches ( ) const
virtualinherited

Get the number of filled bunches in the current configuration.

Implements Trig::IBunchCrossingTool.

Definition at line 811 of file BunchCrossingToolBase.cxx.

811 {
812
813 return m_filledBunches.size();
814 }

◆ numberOfUnpairedBunches()

unsigned int Trig::BunchCrossingToolBase::numberOfUnpairedBunches ( ) const
virtualinherited

Get the number of unpaired bunches in the current configuration.

Implements Trig::IBunchCrossingTool.

Definition at line 816 of file BunchCrossingToolBase.cxx.

816 {
817
818 return m_unpairedBunches.size();
819 }

◆ outputHandles()

virtual std::vector< Gaudi::DataHandle * > AthCommonDataStore< AthCommonMsg< AlgTool > >::outputHandles ( ) const
overridevirtualinherited

Return this algorithm's output handles.

We override this to include handle instances from key arrays if they have not yet been declared. See comments on updateVHKA.

◆ outputMetaStore()

AsgMetadataTool::MetaStorePtr_t asg::AsgMetadataTool::outputMetaStore ( ) const
inherited

Accessor for the output metadata store.

Definition at line 102 of file AsgMetadataTool.cxx.

102 {
103
104#ifdef XAOD_STANDALONE
105 return &m_outputMetaStore;
106#else // XAOD_STANDALONE
107 return m_outputMetaStore;
108#endif // XAOD_STANDALONE
109 }
MetaStore_t m_outputMetaStore
Object accessing the output metadata store.

◆ print()

◆ printConfig()

void Trig::BunchCrossingToolBase::printConfig ( ) const
protectedinherited

Function printing the configuration of the tool.

This function is used to print the overall configuration in a format very similar to how atlas-runqery.cern.ch started showing this information recently.

Definition at line 1199 of file BunchCrossingToolBase.cxx.

1199 {
1200
1201 ATH_MSG_INFO( "No. of coll. bunches : " << m_filledBunches.size() );
1202 ATH_MSG_INFO( "No. of unpaired bunches: " << m_unpairedBunches.size() );
1203 ATH_MSG_INFO( "No. of bunch trains : " << m_bunchTrains.size() );
1204 if( m_bunchTrains.size() ) {
1205 ATH_MSG_INFO( "Bunch spacing in trains: "
1206 << m_bunchTrains.begin()->spacing()
1207 << " ns" );
1208 }
1209
1210 return;
1211 }

◆ renounce()

std::enable_if_t< std::is_void_v< std::result_of_t< decltype(&T::renounce)(T)> > &&!std::is_base_of_v< SG::VarHandleKeyArray, T > &&std::is_base_of_v< Gaudi::DataHandle, T >, void > AthCommonDataStore< AthCommonMsg< AlgTool > >::renounce ( T & h)
inlineprotectedinherited

Definition at line 380 of file AthCommonDataStore.h.

381 {
382 h.renounce();
384 }
std::enable_if_t< std::is_void_v< std::result_of_t< decltype(&T::renounce)(T)> > &&!std::is_base_of_v< SG::VarHandleKeyArray, T > &&std::is_base_of_v< Gaudi::DataHandle, T >, void > renounce(T &h)

◆ renounceArray()

void AthCommonDataStore< AthCommonMsg< AlgTool > >::renounceArray ( SG::VarHandleKeyArray & handlesArray)
inlineprotectedinherited

remove all handles from I/O resolution

Definition at line 364 of file AthCommonDataStore.h.

364 {
366 }

◆ setUseIncidents()

void asg::AsgMetadataTool::setUseIncidents ( const bool flag)
inlineprotectedinherited

Definition at line 132 of file AsgMetadataTool.h.

133 {
135 }
bool flag
Definition master.py:29

◆ sysInitialize()

StatusCode asg::AsgMetadataTool::sysInitialize ( )
virtualinherited

Function initialising the tool in the correct way in Athena.

This function is used to set up the callbacks from IncidentSvc in Athena at the right time during initialisation, without the user having to do anything special in his/her code.

Reimplemented from AthCommonDataStore< AthCommonMsg< AlgTool > >.

Definition at line 115 of file AsgMetadataTool.cxx.

115 {
116
117#ifndef XAOD_STANDALONE
118 if (m_useIncidents) {
119 // Connect to the IncidentSvc:
120 ServiceHandle< IIncidentSvc > incSvc( "IncidentSvc", name() );
121 ATH_CHECK( incSvc.retrieve() );
122
123 // Set up the right callbacks: don't rethrow exceptions, any failure and we should end
124 incSvc->addListener( this, IncidentType::BeginEvent, 0, false );
125 }
126 // Let the base class do its thing:
127 ATH_CHECK( AlgTool::sysInitialize() );
128
129#endif // not XAOD_STANDALONE
130
131 // Return gracefully:
132 return StatusCode::SUCCESS;
133 }

◆ sysStart()

virtual StatusCode AthCommonDataStore< AthCommonMsg< AlgTool > >::sysStart ( )
overridevirtualinherited

Handle START transition.

We override this in order to make sure that conditions handle keys can cache a pointer to the conditions container.

◆ updateVHKA()

void AthCommonDataStore< AthCommonMsg< AlgTool > >::updateVHKA ( Gaudi::Details::PropertyBase & )
inlineinherited

Definition at line 308 of file AthCommonDataStore.h.

308 {
309 // debug() << "updateVHKA for property " << p.name() << " " << p.toString()
310 // << " size: " << m_vhka.size() << endmsg;
311 for (auto &a : m_vhka) {
313 for (auto k : keys) {
314 k->setOwner(this);
315 }
316 }
317 }
std::vector< SG::VarHandleKeyArray * > m_vhka

Member Data Documentation

◆ m_beginInputFileCalled

bool asg::AsgMetadataTool::m_beginInputFileCalled
privateinherited

Flag helping to discover when the tool misses the opening of the first input file.

Definition at line 126 of file AsgMetadataTool.h.

◆ m_bgkey

int Trig::StaticBunchCrossingTool::m_bgkey
private

Default key to be loaded.

Definition at line 63 of file StaticBunchCrossingTool.h.

◆ m_bunchTrains

std::set< Trig::BunchTrain > Trig::BunchCrossingToolBase::m_bunchTrains
protectedinherited

Internal list of bunch trains.

Definition at line 156 of file BunchCrossingToolBase.h.

◆ m_detStore

StoreGateSvc_t AthCommonDataStore< AthCommonMsg< AlgTool > >::m_detStore
privateinherited

Pointer to StoreGate (detector store by default)

Definition at line 393 of file AthCommonDataStore.h.

◆ m_evtStore

StoreGateSvc_t AthCommonDataStore< AthCommonMsg< AlgTool > >::m_evtStore
privateinherited

Pointer to StoreGate (event store by default)

Definition at line 390 of file AthCommonDataStore.h.

◆ m_filledBunches

std::set< Trig::BunchCrossing > Trig::BunchCrossingToolBase::m_filledBunches
protectedinherited

List of colliding bunches.

Definition at line 150 of file BunchCrossingToolBase.h.

◆ m_frontLength

int Trig::BunchCrossingToolBase::m_frontLength
protectedinherited

Length of the "front" of a bunch train.

Definition at line 166 of file BunchCrossingToolBase.h.

◆ m_inputMetaStore

MetaStore_t asg::AsgMetadataTool::m_inputMetaStore
privateinherited

Object accessing the input metadata store.

Definition at line 119 of file AsgMetadataTool.h.

◆ m_knownBGKeys

std::map< int, std::vector< int > > Trig::StaticBunchCrossingTool::m_knownBGKeys
private

All the hard-coded configs.

Definition at line 66 of file StaticBunchCrossingTool.h.

◆ m_maxBunchSpacing

int Trig::BunchCrossingToolBase::m_maxBunchSpacing
protectedinherited

The maximum bunch spacing that the tool should consider.

Definition at line 164 of file BunchCrossingToolBase.h.

◆ m_outputMetaStore

MetaStore_t asg::AsgMetadataTool::m_outputMetaStore
privateinherited

Object accessing the output metadata store.

Definition at line 121 of file AsgMetadataTool.h.

◆ m_singleBunches

std::set< Trig::BunchCrossing > Trig::BunchCrossingToolBase::m_singleBunches
protectedinherited

Internal list of single bunches.

Definition at line 152 of file BunchCrossingToolBase.h.

◆ m_tailLength

int Trig::BunchCrossingToolBase::m_tailLength
protectedinherited

Length of the "tail" of a bunch train.

Definition at line 168 of file BunchCrossingToolBase.h.

◆ m_unpairedBunches

std::set< Trig::BunchCrossing > Trig::BunchCrossingToolBase::m_unpairedBunches
protectedinherited

Internal list of unpaired bunches.

Definition at line 154 of file BunchCrossingToolBase.h.

◆ m_useIncidents

bool asg::AsgMetadataTool::m_useIncidents
privateinherited

Definition at line 128 of file AsgMetadataTool.h.

◆ m_varHandleArraysDeclared

bool AthCommonDataStore< AthCommonMsg< AlgTool > >::m_varHandleArraysDeclared
privateinherited

Definition at line 399 of file AthCommonDataStore.h.

◆ m_vhka

std::vector<SG::VarHandleKeyArray*> AthCommonDataStore< AthCommonMsg< AlgTool > >::m_vhka
privateinherited

Definition at line 398 of file AthCommonDataStore.h.


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