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

Strip RichSmartID time info when used for MC associations

parent a079afd5
No related branches found
No related tags found
1 merge request!4281Strip RichSmartID time info when used for MC associations
......@@ -28,8 +28,10 @@ namespace Rich::Future::Rec::MC {
public:
/// Contructor from the underlying MC relations
Helper( const TkToMCPRels& rels, const LHCb::MCRichDigitSummarys& histories )
: SmartIDUtils( histories ), TrackToMCParticle( rels ) {}
Helper( const TkToMCPRels& rels, //
const LHCb::MCRichDigitSummarys& histories, //
const bool EnableDetDescMCDD4HepFix = true )
: SmartIDUtils( histories, EnableDetDescMCDD4HepFix ), TrackToMCParticle( rels ) {}
public:
using SmartIDUtils::mcParticles;
......
......@@ -35,7 +35,7 @@ namespace Rich::Future::MC::Relations {
public:
/// Constructor from a summaries object
SmartIDUtils( const LHCb::MCRichDigitSummarys& histories );
SmartIDUtils( const LHCb::MCRichDigitSummarys& histories, const bool EnableDetDescMCDD4HepFix = true );
public:
/// Returns a vector of the MCParticles associated to a given pixel cluster
......
......@@ -70,23 +70,9 @@ namespace Rich::Future::MC {
// Collect all the hits in the same pixel
std::map<LHCb::RichSmartID, std::vector<LHCb::MCRichHit*>> hitsPerPix;
for ( const auto mchit : mchits ) {
auto id = mchit->sensDetID();
const auto id = mchit->sensDetID();
if ( mchit && id.isValid() ) {
#ifdef USE_DD4HEP
if ( m_detdescMCinput ) {
// If built for DD4HEP apply correction to PMT module numbers to account
// for different numbering scheme between DD4HEP and DetDesc.
// Option needs to be explicitly activated only when input is known to
// be DetDesc based MC.
// ***** To eventually be removed when DetDesc finally dies completely *****
const Rich::DetectorArray<Rich::PanelArray<LHCb::RichSmartID::DataType>> mod_corr{{{0, 0}, {6, 18}}};
const auto pdMod = id.pdMod() + mod_corr[id.rich()][id.panel()];
const auto pdNumInMod = id.pdNumInMod();
id.setPD( pdMod, pdNumInMod );
}
#endif
// Apply all sorts MC based filtering here, e.g. signal only, timing windows etc.
// To be extended as needed
if ( m_rejectBackground && mchit->isBackground() ) {
......@@ -109,6 +95,20 @@ namespace Rich::Future::MC {
// temporary copy of SmartID
auto newID = id;
#ifdef USE_DD4HEP
if ( m_detdescMCinput ) {
// If built for DD4HEP apply correction to PMT module numbers to account
// for different numbering scheme between DD4HEP and DetDesc.
// Option needs to be explicitly activated only when input is known to
// be DetDesc based MC.
// ***** To eventually be removed when DetDesc finally dies completely *****
const Rich::DetectorArray<Rich::PanelArray<LHCb::RichSmartID::DataType>> mod_corr{{{0, 0}, {6, 18}}};
const auto pdMod = id.pdMod() + mod_corr[id.rich()][id.panel()];
const auto pdNumInMod = id.pdNumInMod();
newID.setPD( pdMod, pdNumInMod );
}
#endif
// Add best time for this ID
if ( m_includeTimeInfo ) {
// attempt to select the 'best' time for this ID ...
......@@ -124,15 +124,18 @@ namespace Rich::Future::MC {
}
}
// Add ID to decoded list
// Add (dd4hep corrected) ID to decoded list
ids.emplace_back( newID );
// Form MC summaries for each MCHit
for ( const auto hit : hits ) {
auto summary = std::make_unique<LHCb::MCRichDigitSummary>();
summary->setRichSmartID( id ); // Use version without time
// do not want dd4hep corrected ID here as utility does its own correction.
summary->setRichSmartID( id );
summary->setMCParticle( hit->mcParticle() );
summary->setHistory( hit->mcRichDigitHistoryCode() );
auto history = hit->mcRichDigitHistoryCode();
history.setSignalEvent( true ); // normally done by Boole for signal event only
summary->setHistory( history );
summaries.add( summary.release() );
}
}
......
......@@ -22,7 +22,7 @@ MCHitUtils::MCHitUtils( const LHCb::MCRichHits& mchits ) {
// Build the mappings to MCHits
// Get the pixel level ID
const auto pixelID = hit->sensDetID().pixelID();
const auto pixelID = hit->sensDetID().pixelID().channelDataOnly();
// append hit to list for this ID
m_smartIDsToHits[pixelID].emplace_back( hit );
......
......@@ -57,11 +57,11 @@ const LHCb::MCParticle* Helper::trueRecPhoton( const LHCb::MCParticle& mcPart, /
const LHCb::RichSmartID id ) const {
// Get MCParticles for the channel
const auto mcParts = mcParticles( id );
const auto mcParts = mcParticles( id.channelDataOnly() );
// Loop over all MCParticles associated to the pixel
for ( const auto* MCP : mcParts ) {
if ( &mcPart == MCP ) return MCP;
if ( &mcPart == MCP ) { return MCP; }
}
// if get here no association
......@@ -98,11 +98,11 @@ LHCb::MCParticle::ConstVector Helper::trueCherenkovPhoton( const LHCb::Track&
// loop over the track MCPs
for ( const auto tkMCP : tkMCPs ) {
if ( !tkMCP ) continue;
if ( !tkMCP ) { continue; }
// save primary ID
saveMCP( tkMCP, cluster.primaryID() );
saveMCP( tkMCP, cluster.primaryID().channelDataOnly() );
// loop over secondary IDs
for ( const auto S : cluster.secondaryIDs() ) { saveMCP( tkMCP, S ); }
for ( const auto S : cluster.secondaryIDs() ) { saveMCP( tkMCP, S.channelDataOnly() ); }
}
return trueMCPs;
......
......@@ -17,7 +17,7 @@
using namespace Rich::Future::MC::Relations;
SmartIDUtils::SmartIDUtils( const LHCb::MCRichDigitSummarys& histories ) {
SmartIDUtils::SmartIDUtils( const LHCb::MCRichDigitSummarys& histories, const bool EnableDetDescMCDD4HepFix ) {
// Fill the summary map
// If this class is used a lot might want to write an algorithm
// to create the map and save it in the TES, to avoid it being created
......@@ -25,16 +25,18 @@ SmartIDUtils::SmartIDUtils( const LHCb::MCRichDigitSummarys& histories ) {
for ( const auto* sum : histories ) {
if ( sum ) {
auto id = sum->richSmartID();
if ( EnableDetDescMCDD4HepFix ) {
#ifdef USE_DD4HEP
// If built for DD4HEP apply correction to PMT module numbers to account
// for different numbering scheme between DD4HEP and DetDesc.
// This is only required when processing DetDesc MC samples with dd4hep
// geometry and should be removed once dd4hep MC samples are available.
const Rich::DetectorArray<Rich::PanelArray<LHCb::RichSmartID::DataType>> mod_corr{{{0, 0}, {6, 18}}};
const auto pdMod = id.pdMod() + mod_corr[id.rich()][id.panel()];
const auto pdNumInMod = id.pdNumInMod();
id.setPD( pdMod, pdNumInMod );
// If built for DD4HEP apply correction to PMT module numbers to account
// for different numbering scheme between DD4HEP and DetDesc.
// This is only required when processing DetDesc MC samples with dd4hep
// geometry and should be removed once dd4hep MC samples are available.
const Rich::DetectorArray<Rich::PanelArray<LHCb::RichSmartID::DataType>> mod_corr{{{0, 0}, {6, 18}}};
const auto pdMod = id.pdMod() + mod_corr[id.rich()][id.panel()];
const auto pdNumInMod = id.pdNumInMod();
id.setPD( pdMod, pdNumInMod );
#endif
}
m_sumMap[id].push_back( sum );
}
}
......@@ -58,7 +60,7 @@ SmartIDUtils::mcDigitHistoryCodes( const Rich::PDPixelCluster& cluster ) const {
std::vector<LHCb::MCRichDigitHistoryCode> SmartIDUtils::mcDigitHistoryCodes( const LHCb::RichSmartID id ) const {
std::vector<LHCb::MCRichDigitHistoryCode> codes;
// Does this hit have an entry in the pixel summary map
const auto iEn = m_sumMap.find( id );
const auto iEn = m_sumMap.find( id.channelDataOnly() );
if ( iEn != m_sumMap.end() ) {
// reserve size
codes.reserve( ( *iEn ).second.size() );
......@@ -73,12 +75,12 @@ std::vector<LHCb::MCRichDigitHistoryCode> SmartIDUtils::mcDigitHistoryCodes( con
LHCb::MCParticle::ConstVector SmartIDUtils::mcParticles( const Rich::PDPixelCluster& cluster ) const {
// Get the MCps for the primary ID
auto mcParts = mcParticles( cluster.primaryID() );
auto mcParts = mcParticles( cluster.primaryID().channelDataOnly() );
// now loop over the cluster secondary smart IDs
for ( const auto& S : cluster.secondaryIDs() ) {
// Get the MCPs for this channel
const auto tmp_vect = mcParticles( S );
const auto tmp_vect = mcParticles( S.channelDataOnly() );
// reserve space
mcParts.reserve( mcParts.size() + tmp_vect.size() );
// save if not already present
......@@ -95,7 +97,7 @@ LHCb::MCParticle::ConstVector SmartIDUtils::mcParticles( const LHCb::RichSmartID
LHCb::MCParticle::ConstVector mcParts;
// Does this hit have an entry in the pixel summary map
const auto iEn = m_sumMap.find( id );
const auto iEn = m_sumMap.find( id.channelDataOnly() );
if ( iEn != m_sumMap.end() ) {
mcParts.reserve( ( *iEn ).second.size() );
// Loop over the associated summaries
......@@ -103,10 +105,10 @@ LHCb::MCParticle::ConstVector SmartIDUtils::mcParticles( const LHCb::RichSmartID
// Get the MCP
const auto* mcP = sum->mcParticle();
// protect against null references
if ( !mcP ) continue;
if ( !mcP ) { continue; }
// Add to vector, once per MCParticle
const auto iFind = std::find( mcParts.begin(), mcParts.end(), mcP );
if ( mcParts.end() == iFind ) mcParts.push_back( mcP );
if ( mcParts.end() == iFind ) { mcParts.push_back( mcP ); }
}
}
......@@ -115,14 +117,17 @@ LHCb::MCParticle::ConstVector SmartIDUtils::mcParticles( const LHCb::RichSmartID
const LHCb::MCParticle* SmartIDUtils::trueCherenkovRadiation( const LHCb::RichSmartID id,
const Rich::RadiatorType rad ) const {
// Strip time info from ID
const auto id_noT = id.channelDataOnly();
// Test if hit is background
if ( isBackground( id ) ) return nullptr;
if ( isBackground( id_noT ) ) { return nullptr; }
// Test if hit is from correct radiator
if ( !isCherenkovRadiation( id, rad ) ) return nullptr;
if ( !isCherenkovRadiation( id_noT, rad ) ) { return nullptr; }
// All OK so find correct MCParticle
const auto mcParts = mcParticles( id );
const auto mcParts = mcParticles( id_noT );
// Should do something better than just return first, but what ??
return ( mcParts.empty() ? nullptr : mcParts.front() );
......@@ -136,7 +141,7 @@ const LHCb::MCParticle* SmartIDUtils::trueCherenkovRadiation( const Rich::PDPixe
// now try the secondary IDs
for ( const auto S : cluster.secondaryIDs() ) {
mcP = trueCherenkovRadiation( S, rad );
if ( mcP ) break;
if ( mcP ) { break; }
}
}
// return
......@@ -147,7 +152,7 @@ bool SmartIDUtils::isBackground( const LHCb::RichSmartID id ) const {
bool isback = true;
// try via summary objects
const auto iEn = m_sumMap.find( id );
const auto iEn = m_sumMap.find( id.channelDataOnly() );
if ( iEn != m_sumMap.end() ) {
// returns true if hit is only background (no signal contribution)
isback = std::none_of( ( *iEn ).second.begin(), ( *iEn ).second.end(),
......@@ -161,7 +166,7 @@ bool SmartIDUtils::isCherenkovRadiation( const LHCb::RichSmartID id, const Rich:
bool isCKRad = false;
// first, try via summary objects
const auto iEn = m_sumMap.find( id );
const auto iEn = m_sumMap.find( id.channelDataOnly() );
if ( iEn != m_sumMap.end() ) {
isCKRad = std::any_of( ( *iEn ).second.begin(), ( *iEn ).second.end(), [&rad]( const auto* sum ) {
const auto& code = sum->history();
......
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