Skip to content
Snippets Groups Projects
Commit 5c4e33ef authored by Christopher Rob Jones's avatar Christopher Rob Jones
Browse files

RichPID: Small clean up and modernisation

parent d02217c8
No related branches found
No related tags found
1 merge request!4288RichPID: Small clean up and modernisation
Pipeline #6738100 passed
......@@ -35,13 +35,8 @@ namespace LHCb {
// Namespace for locations in TDS
namespace RichPIDLocation {
inline const std::string Default = "Rec/Rich/PIDs";
inline const std::string Offline = "Rec/Rich/PIDs";
inline const std::string OfflineGlobal = "Rec/Rich/GlobalPIDs";
inline const std::string OfflineLocal = "Rec/Rich/LocalPIDs";
inline const std::string HLT = "Rec/Rich/HltPIDs";
inline const std::string HLTGlobal = "Rec/Rich/HltGlobalPIDs";
inline const std::string HLTLocal = "Rec/Rich/HltLocalPIDs";
inline static const std::string Default = "Rec/Rich/PIDs";
inline static const std::string Offline = Default;
} // namespace RichPIDLocation
/** @class RichPID RichPID.h
......@@ -54,181 +49,197 @@ namespace LHCb {
class RichPID : public KeyedObject<int> {
public:
/// typedef for std::vector of RichPID
// Types for containers of PIDs
using Vector = std::vector<RichPID*>;
using ConstVector = std::vector<const RichPID*>;
using Container = KeyedContainer<RichPID, Containers::HashMap>;
using Selection = SharedObjectsContainer<RichPID>;
using Range = Gaudi::Range_<ConstVector>;
// Data types
using FltType = float;
using DLL = FltType;
using DLLs = std::vector<DLL>;
/// typedef for KeyedContainer of RichPID
typedef KeyedContainer<RichPID, Containers::HashMap> Container;
using Selection = SharedObjectsContainer<RichPID>;
using Range = Gaudi::Range_<ConstVector>;
public:
/// Constructor from data fields
RichPID( Rich::ParticleIDType pid, const LHCb::Track* track, std::vector<float> dllValues );
RichPID( const Rich::ParticleIDType pid, const LHCb::Track* track, DLLs&& dllValues )
: m_particleLLValues( std::forward<DLLs>( dllValues ) ), m_track( track ) {
setBestParticleID( pid );
}
/// Copy constructor. Creates a new RichPID object with the same pid information
/// Copy constructor
RichPID( const LHCb::RichPID& pid )
: KeyedObject<int>()
, m_pidResultCode( pid.m_pidResultCode )
, m_particleLLValues( pid.m_particleLLValues )
, m_track( pid.m_track ) {}
/// Move constructor
RichPID( LHCb::RichPID&& pid )
: KeyedObject<int>()
, m_pidResultCode( pid.m_pidResultCode )
, m_particleLLValues( std::move( pid.m_particleLLValues ) )
, m_track( pid.m_track ) {}
/// Default Constructor
RichPID() : m_particleLLValues( Rich::NParticleTypes, 0 ) {}
RichPID() : m_particleLLValues( Rich::NParticleTypes, 0.0 ) {}
public:
// Retrieve pointer to class definition structure
const CLID& clID() const override { return LHCb::RichPID::classID(); }
static const CLID& classID() { return CLID_RichPID; }
public:
/// Returns the raw gaussian probability value for a given particle species
float particleRawProb( Rich::ParticleIDType type ) const;
FltType particleRawProb( const Rich::ParticleIDType type ) const {
const auto dll = particleDeltaLL( type );
return ( dll > 0 ? 1.0F - (FltType)gsl_sf_erf( std::sqrt( dll ) ) : 1.0F );
}
/// Returns the normalised gaussian probability value for a given particle species
float particleNormProb( Rich::ParticleIDType type ) const;
FltType particleNormProb( const Rich::ParticleIDType type ) const {
FltType norm = 0;
for ( const auto pid : Rich::Particles() ) { norm += particleRawProb( pid ); }
return ( norm > 0 ? particleRawProb( type ) / norm : 0 );
}
/// Returns the delta LL value for a given hypothesis
float particleDeltaLL( Rich::ParticleIDType type ) const { return m_particleLLValues[type]; }
DLL particleDeltaLL( const Rich::ParticleIDType type ) const noexcept { return m_particleLLValues[type]; }
/// Sets the particle delta LL value for the given hypothesis
void setParticleDeltaLL( Rich::ParticleIDType type, float deltaLL ) { m_particleLLValues[type] = deltaLL; }
void setParticleDeltaLL( const Rich::ParticleIDType type, const DLL deltaLL ) noexcept {
m_particleLLValues[type] = deltaLL;
}
/// Sets the particle delta LL value for the given hypothesis
void setParticleLLValues( std::vector<float>&& value ) { m_particleLLValues = std::move( value ); }
void setParticleLLValues( DLLs&& value ) {
assert( value.size() == m_particleLLValues.size() );
m_particleLLValues = std::move( value );
}
/// Boolean method to check if the given radiator was used to create this PID result (i.e. if the associated track
/// was found to traverse the radiator in a manner that would have produced detectable Cherenkov photons
bool traversedRadiator( Rich::RadiatorType radiator ) const;
bool traversedRadiator( const Rich::RadiatorType radiator ) const noexcept {
return ( Rich::Aerogel == radiator ? this->usedAerogel() : //
Rich::Rich1Gas == radiator ? this->usedRich1Gas() : //
Rich::Rich2Gas == radiator ? this->usedRich2Gas() : false );
}
/// Verify if a given hypothesis was above threshold and the associated track present in any radiator
bool isAboveThreshold( Rich::ParticleIDType type ) const;
bool isAboveThreshold( const Rich::ParticleIDType type ) const;
/// Set a given hypothesis above threshold and the associated track present in any radiator
void setAboveThreshold( Rich::ParticleIDType type, bool flag );
void setAboveThreshold( const Rich::ParticleIDType type, const bool flag );
/// Returns the signed sigma separation beween 2 particle hypotheses (first relative to last)
float nSigmaSeparation( Rich::ParticleIDType firstPart, Rich::ParticleIDType lastPart ) const;
FltType nSigmaSeparation( const Rich::ParticleIDType firstPart, const Rich::ParticleIDType lastPart ) const {
const auto dLL = m_particleLLValues[lastPart] - m_particleLLValues[firstPart];
return std::sqrt( std::abs( 2.0F * dLL ) ) * ( dLL > 0 ? 1.0F : -1.0F );
}
/// Returns true if the given mass hypothesis is within the given number of sigmas separation from the best PID type
bool isConsistentNSigma( Rich::ParticleIDType type, float nsigma ) const;
/// Textual representation of PID type
std::string pidType() const;
bool isConsistentNSigma( const Rich::ParticleIDType type, const FltType nsigma ) const {
return ( nSigmaSeparation( bestParticleID(), type ) > nsigma );
}
/// The best Particle ID
Rich::ParticleIDType bestParticleID() const;
Rich::ParticleIDType bestParticleID() const noexcept {
// Shift by -1 to convert packed representation to Rich::ParticleIDType
return ( Rich::ParticleIDType )( ( ( m_pidResultCode & packedBestParticleIDMask ) >> packedBestParticleIDBits ) -
1 );
}
/// set the best particle ID
void setBestParticleID( Rich::ParticleIDType type );
void setBestParticleID( const Rich::ParticleIDType type ) noexcept {
// Shift bit-packed representation by +1 to start numbering from 0
const auto val = (unsigned int)type + 1;
m_pidResultCode &= ~packedBestParticleIDMask;
m_pidResultCode |= ( ( ( (unsigned int)val ) << packedBestParticleIDBits ) & packedBestParticleIDMask );
}
/// Print this RichPID data object in a human readable way
std::ostream& fillStream( std::ostream& s ) const override;
/// Retrieve const Bit-packed information (Best particle ID and History) for the RichPID result
unsigned int pidResultCode() const { return m_pidResultCode; }
unsigned int pidResultCode() const noexcept { return m_pidResultCode; }
/// Update Bit-packed information (Best particle ID and History) for the RichPID result
void setPidResultCode( unsigned int value ) { m_pidResultCode = value; }
void setPidResultCode( const unsigned int value ) noexcept { m_pidResultCode = value; }
/// Retrieve Information from aerogel was used to form this PID result
bool usedAerogel() const { return test<usedAerogelMask>(); }
bool usedAerogel() const noexcept { return test<usedAerogelMask>(); }
/// Update Information from aerogel was used to form this PID result
void setUsedAerogel( bool value ) { set<usedAerogelMask>( value ); }
void setUsedAerogel( const bool value ) noexcept { set<usedAerogelMask>( value ); }
/// Retrieve Information from Rich1 gas was used to form this PID result
bool usedRich1Gas() const { return test<usedRich1GasMask>(); }
bool usedRich1Gas() const noexcept { return test<usedRich1GasMask>(); }
/// Update Information from Rich1 gas was used to form this PID result
void setUsedRich1Gas( bool value ) { set<usedRich1GasMask>( value ); }
void setUsedRich1Gas( const bool value ) noexcept { set<usedRich1GasMask>( value ); }
/// Retrieve Information from Rich2 gas was used to form this PID result
bool usedRich2Gas() const { return test<usedRich2GasMask>(); }
bool usedRich2Gas() const noexcept { return test<usedRich2GasMask>(); }
/// Update Information from Rich2 gas was used to form this PID result
void setUsedRich2Gas( bool value ) { set<usedRich2GasMask>( value ); }
void setUsedRich2Gas( const bool value ) noexcept { set<usedRich2GasMask>( value ); }
/// Retrieve The electron hypothesis is above threshold in at least one active radiator
bool electronHypoAboveThres() const { return test<electronHypoAboveThresMask>(); }
bool electronHypoAboveThres() const noexcept { return test<electronHypoAboveThresMask>(); }
/// Update The electron hypothesis is above threshold in at least one active radiator
void setElectronHypoAboveThres( bool value ) { set<electronHypoAboveThresMask>( value ); }
void setElectronHypoAboveThres( const bool value ) noexcept { set<electronHypoAboveThresMask>( value ); }
/// Retrieve The muon hypothesis is above threshold in at least one active radiator
bool muonHypoAboveThres() const { return test<muonHypoAboveThresMask>(); }
bool muonHypoAboveThres() const noexcept { return test<muonHypoAboveThresMask>(); }
/// Update The muon hypothesis is above threshold in at least one active radiator
void setMuonHypoAboveThres( bool value ) { set<muonHypoAboveThresMask>( value ); }
void setMuonHypoAboveThres( const bool value ) noexcept { set<muonHypoAboveThresMask>( value ); }
/// Retrieve The pion hypothesis is above threshold in at least one active radiator
bool pionHypoAboveThres() const { return test<pionHypoAboveThresMask>(); }
bool pionHypoAboveThres() const noexcept { return test<pionHypoAboveThresMask>(); }
/// Update The pion hypothesis is above threshold in at least one active radiator
void setPionHypoAboveThres( bool value ) { set<pionHypoAboveThresMask>( value ); }
void setPionHypoAboveThres( const bool value ) noexcept { set<pionHypoAboveThresMask>( value ); }
/// Retrieve The kaon hypothesis is above threshold in at least one active radiator
bool kaonHypoAboveThres() const { return test<kaonHypoAboveThresMask>(); }
bool kaonHypoAboveThres() const noexcept { return test<kaonHypoAboveThresMask>(); }
/// Update The kaon hypothesis is above threshold in at least one active radiator
void setKaonHypoAboveThres( bool value ) { set<kaonHypoAboveThresMask>( value ); }
void setKaonHypoAboveThres( const bool value ) noexcept { set<kaonHypoAboveThresMask>( value ); }
/// Retrieve The proton hypothesis is above threshold in at least one active radiator
bool protonHypoAboveThres() const { return test<protonHypoAboveThresMask>(); }
bool protonHypoAboveThres() const noexcept { return test<protonHypoAboveThresMask>(); }
/// Update The proton hypothesis is above threshold in at least one active radiator
void setProtonHypoAboveThres( bool value ) { set<protonHypoAboveThresMask>( value ); }
/// Retrieve RICH Offline Global PID result
bool offlineGlobal() const { return test<offlineGlobalMask>(); }
/// Update RICH Offline Global PID result
void setOfflineGlobal( bool value ) { return set<offlineGlobalMask>( value ); }
/// Retrieve RICH Offline Local PID result
bool offlineLocal() const { return test<offlineLocalMask>(); }
/// Update RICH Offline Local PID result
void setOfflineLocal( bool value ) { set<offlineLocalMask>( value ); }
/// Retrieve RICH Offline Ring Refit PID result
bool ringRefit() const { return test<ringRefitMask>(); }
/// Update RICH Offline Ring Refit PID result
void setRingRefit( bool value ) { set<ringRefitMask>( value ); }
/// Retrieve RICH HLT Local PID result
bool hltLocal() const { return test<hltLocalMask>(); }
/// Update RICH HLT Local PID result
void setHltLocal( bool value ) { set<hltLocalMask>( value ); }
/// Retrieve RICH HLT Global PID result
bool hltGlobal() const { return test<hltGlobalMask>(); }
/// Update RICH HLT Global PID result
void setHltGlobal( bool value ) { set<hltGlobalMask>( value ); }
void setProtonHypoAboveThres( const bool value ) noexcept { set<protonHypoAboveThresMask>( value ); }
/// Retrieve The deuteron hypothesis is above threshold in at least one active radiator
bool deuteronHypoAboveThres() const { return test<deuteronHypoAboveThresMask>(); }
bool deuteronHypoAboveThres() const noexcept { return test<deuteronHypoAboveThresMask>(); }
/// Update The deuteron hypothesis is above threshold in at least one active radiator
void setDeuteronHypoAboveThres( bool value ) { set<deuteronHypoAboveThresMask>( value ); }
void setDeuteronHypoAboveThres( const bool value ) noexcept { set<deuteronHypoAboveThresMask>( value ); }
/// Retrieve non-const Vector of particle hypothesis log likelihood values
DLLs& particleLLValues() noexcept { return m_particleLLValues; }
/// Retrieve const Vector of particle hypothesis log likelihood values
const std::vector<float>& particleLLValues() const { return m_particleLLValues; }
const DLLs& particleLLValues() const noexcept { return m_particleLLValues; }
/// Update Vector of particle hypothesis log likelihood values
void setParticleLLValues( const std::vector<float>& value ) { m_particleLLValues = value; }
void setParticleLLValues( const DLLs& value ) noexcept { m_particleLLValues = value; }
/// Retrieve (const) Associated reconstructed Track
const LHCb::Track* track() const { return m_track; }
const LHCb::Track* track() const noexcept { return m_track; }
/// Update Associated reconstructed Track
void setTrack( const SmartRef<LHCb::Track>& value ) { m_track = value; }
void setTrack( const SmartRef<LHCb::Track>& value ) noexcept { m_track = value; }
/// Update (pointer) Associated reconstructed Track
void setTrack( const LHCb::Track* value ) { m_track = value; }
void setTrack( const LHCb::Track* value ) noexcept { m_track = value; }
public:
/// overload ostream output
friend std::ostream& operator<<( std::ostream& str, const RichPID& obj ) { return obj.fillStream( str ); }
public:
......@@ -254,11 +265,11 @@ namespace LHCb {
private:
template <unsigned int mask>
void set( bool value ) {
void set( const bool value ) noexcept {
m_pidResultCode = ( m_pidResultCode & ~mask ) | ( -value & mask );
}
template <unsigned int mask>
bool test() const {
bool test() const noexcept {
return ( m_pidResultCode & mask ) != 0;
}
......@@ -278,17 +289,18 @@ namespace LHCb {
pionHypoAboveThresMask = 0x200U,
kaonHypoAboveThresMask = 0x400U,
protonHypoAboveThresMask = 0x800U,
offlineGlobalMask = 0x1000U,
offlineLocalMask = 0x2000U,
ringRefitMask = 0x4000U,
hltLocalMask = 0x8000U,
hltGlobalMask = 0x10000U,
// Redundant types removed here ...
deuteronHypoAboveThresMask = 0x20000U
};
unsigned int m_pidResultCode{0}; ///< Bit-packed information (Best particle ID and History) for the RichPID result
std::vector<float> m_particleLLValues; ///< Vector of particle hypothesis log likelihood values
SmartRef<LHCb::Track> m_track; ///< Associated reconstructed Track
/// Bit-packed information (Best particle ID and History) for the RichPID result
unsigned int m_pidResultCode{0};
/// Vector of particle hypothesis log likelihood values
DLLs m_particleLLValues;
/// Associated reconstructed Track
SmartRef<LHCb::Track> m_track;
}; // class RichPID
......@@ -296,54 +308,3 @@ namespace LHCb {
using RichPIDs = RichPID::Container;
} // namespace LHCb
// -----------------------------------------------------------------------------
// end of class
// -----------------------------------------------------------------------------
// Including forward declarations
inline LHCb::RichPID::RichPID( const Rich::ParticleIDType pid, const LHCb::Track* track, std::vector<float> dllValues )
: m_pidResultCode( 0 ), m_particleLLValues( std::move( dllValues ) ), m_track( track ) {
setBestParticleID( pid );
}
inline float LHCb::RichPID::particleRawProb( const Rich::ParticleIDType type ) const {
const auto dll = particleDeltaLL( type );
return dll > 0 ? 1.0F - (float)gsl_sf_erf( std::sqrt( dll ) ) : 1.0F;
}
inline float LHCb::RichPID::particleNormProb( const Rich::ParticleIDType type ) const {
float norm = 0;
for ( const auto pid : Rich::Particles() ) { norm += particleRawProb( pid ); }
return ( norm > 0 ? particleRawProb( type ) / norm : 0 );
}
inline float LHCb::RichPID::nSigmaSeparation( const Rich::ParticleIDType firstPart,
const Rich::ParticleIDType lastPart ) const {
const auto dLL = m_particleLLValues[lastPart] - m_particleLLValues[firstPart];
return std::sqrt( std::abs( 2.0F * dLL ) ) * ( dLL > 0 ? 1.0F : -1.0F );
}
inline bool LHCb::RichPID::isConsistentNSigma( const Rich::ParticleIDType type, const float nsigma ) const {
return ( nSigmaSeparation( bestParticleID(), type ) > nsigma );
}
inline Rich::ParticleIDType LHCb::RichPID::bestParticleID() const {
// Shift by -1 to convert packed representation to Rich::ParticleIDType
return ( Rich::ParticleIDType )( ( ( m_pidResultCode & packedBestParticleIDMask ) >> packedBestParticleIDBits ) - 1 );
}
inline void LHCb::RichPID::setBestParticleID( const Rich::ParticleIDType type ) {
// Shift bit-packed representation by +1 to start numbering from 0
const auto val = (unsigned int)type + 1;
m_pidResultCode &= ~packedBestParticleIDMask;
m_pidResultCode |= ( ( ( (unsigned int)val ) << packedBestParticleIDBits ) & packedBestParticleIDMask );
}
......@@ -30,48 +30,23 @@
// Boost
#include "boost/format.hpp"
std::string LHCb::RichPID::pidType() const {
std::string hist;
int cnt = 0;
if ( this->offlineGlobal() ) {
hist += "OfflineGlobal";
++cnt;
}
if ( this->offlineLocal() ) {
hist += ( cnt > 0 ? "+" : "" );
hist += "OfflineLocal";
++cnt;
}
if ( this->ringRefit() ) {
hist += ( cnt > 0 ? "+" : "" );
hist += "RingRefit";
++cnt;
}
if ( this->hltGlobal() ) {
hist += ( cnt > 0 ? "+" : "" );
hist += "HltGlobal";
++cnt;
}
if ( this->hltLocal() ) {
hist += ( cnt > 0 ? "+" : "" );
hist += "HltLocal";
++cnt;
}
return hist;
}
bool LHCb::RichPID::isAboveThreshold( const Rich::ParticleIDType type ) const {
return ( Rich::Electron == type
? this->electronHypoAboveThres()
: Rich::Muon == type
? this->muonHypoAboveThres()
: Rich::Pion == type
? this->pionHypoAboveThres()
: Rich::Kaon == type
? this->kaonHypoAboveThres()
: Rich::Proton == type
? this->protonHypoAboveThres()
: Rich::Deuteron == type ? this->deuteronHypoAboveThres() : false );
switch ( type ) {
case Rich::Electron:
return this->electronHypoAboveThres();
case Rich::Muon:
return this->muonHypoAboveThres();
case Rich::Pion:
return this->pionHypoAboveThres();
case Rich::Kaon:
return this->kaonHypoAboveThres();
case Rich::Proton:
return this->protonHypoAboveThres();
case Rich::Deuteron:
return this->deuteronHypoAboveThres();
default:
return false;
}
}
void LHCb::RichPID::setAboveThreshold( const Rich::ParticleIDType type, const bool flag ) {
......@@ -90,13 +65,6 @@ void LHCb::RichPID::setAboveThreshold( const Rich::ParticleIDType type, const bo
}
}
bool LHCb::RichPID::traversedRadiator( const Rich::RadiatorType radiator ) const {
return ( Rich::Aerogel == radiator
? this->usedAerogel()
: Rich::Rich1Gas == radiator ? this->usedRich1Gas()
: Rich::Rich2Gas == radiator ? this->usedRich2Gas() : false );
}
std::ostream& LHCb::RichPID::fillStream( std::ostream& s ) const {
s << "[ ";
......@@ -105,7 +73,7 @@ std::ostream& LHCb::RichPID::fillStream( std::ostream& s ) const {
const std::string iF = "%5i"; // ints
// PID type
s << "Key" << boost::format( iF ) % key() << " " << pidType();
s << "Key" << boost::format( iF ) % key();
// Track info
if ( track() ) {
......
0% Loading or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment