From 6cd9bed62670c4c26ba9d7c6a56604fb05d4b692 Mon Sep 17 00:00:00 2001
From: nmchugh <niall.thomas.mchugh@cern.ch>
Date: Thu, 13 Mar 2025 00:32:44 +0000
Subject: [PATCH 1/2] add initial version of pv splitting
---
Rec/RecAlgs/CMakeLists.txt | 1 +
Rec/RecAlgs/src/LumiPVs_nobeamline.cpp | 112 ++++++++-----
Rec/RecAlgs/src/PVSplitter.cpp | 222 +++++++++++++++++++++++++
3 files changed, 295 insertions(+), 40 deletions(-)
create mode 100644 Rec/RecAlgs/src/PVSplitter.cpp
diff --git a/Rec/RecAlgs/CMakeLists.txt b/Rec/RecAlgs/CMakeLists.txt
index 93fdeb5ec75..e1582852bac 100644
--- a/Rec/RecAlgs/CMakeLists.txt
+++ b/Rec/RecAlgs/CMakeLists.txt
@@ -21,6 +21,7 @@ gaudi_add_module(RecAlgs
src/UniqueIDGeneratorAlg.cpp
src/LumiSummaryTuple.cpp
src/LumiPVs_nobeamline.cpp
+ src/PVSplitter.cpp
LINK
AIDA::aida
Boost::headers
diff --git a/Rec/RecAlgs/src/LumiPVs_nobeamline.cpp b/Rec/RecAlgs/src/LumiPVs_nobeamline.cpp
index 3235b435a54..f1c3367de7b 100644
--- a/Rec/RecAlgs/src/LumiPVs_nobeamline.cpp
+++ b/Rec/RecAlgs/src/LumiPVs_nobeamline.cpp
@@ -21,64 +21,96 @@
#include "DetDesc/DetectorElement.h"
#include "TrackInterfaces/IPVOfflineTool.h"
-using InputVeloTracks = LHCb::Event::v1::Tracks;
+using InputPVs = LHCb::RecVertices;
using InputHlt = LHCb::HltDecReports;
using RecPV = LHCb::RecVertex;
+
+namespace {
+ int count_backwards_tracks( RecPV* pv ) {
+ int nbw = 0;
+ for ( auto tr: pv->tracks()) {
+ if( tr->isVeloBackward() ) nbw++;
+ }
+ return nbw;
+ }
+
+ // convert from LHCb::RecVertices -> std::vector<LHCb::RecVertex*>
+ void convert ( const InputPVs& oldvec, std::vector<RecPV*>& newvec ) {
+ std::copy( oldvec.begin(), oldvec.end(), std::back_inserter( newvec ) );
+ }
+}
+
class LumiPVs_nobeamline
: public virtual LHCb::Algorithm::Consumer<
- void( const LHCb::ODIN&, const InputVeloTracks&, const InputHlt&, const DetectorElement& ),
- LHCb::Algorithm::Traits::usesBaseAndConditions<GaudiTupleAlg, DetectorElement>> {
+ void( const LHCb::ODIN&, const InputPVs&, const InputPVs&, const InputPVs&, const InputHlt& ),
+ LHCb::Algorithm::Traits::usesBaseAndConditions<GaudiTupleAlg>> {
+
private:
- Gaudi::Property<std::string> m_tupleName{ this, "BeamGasTuple", "BeamGasTuple" };
- ToolHandle<IPVOfflineTool> m_toolpv{ this, "PVOfflineTool", "PVOfflineTool" };
- LHCb::RecVertices reconstructMultiPVFromTracks( std::vector<const LHCb::Event::v1::Track*>& tracks2use,
- Gaudi::XYZPoint const& beamSpot ) const;
+ Gaudi::Property<std::string> m_tupleName{this, "BeamGasTuple", "BeamGasTuple"};
public:
- virtual void operator()( const LHCb::ODIN&, const InputVeloTracks&, const InputHlt&,
- const DetectorElement& ) const override;
+ virtual void operator()( const LHCb::ODIN&, const InputPVs&, const InputPVs&, const InputPVs&, const InputHlt& ) const override;
LumiPVs_nobeamline( const std::string& name, ISvcLocator* pSvcLocator );
};
DECLARE_COMPONENT_WITH_ID( LumiPVs_nobeamline, "LumiPVs_nobeamline" );
LumiPVs_nobeamline::LumiPVs_nobeamline( const std::string& name, ISvcLocator* pSvcLocator )
: Consumer( name, pSvcLocator,
- { KeyValue{ "ODIN", "" }, KeyValue{ "InputVeloTracks", "" }, KeyValue{ "InputHlt", "" },
- KeyValue{ "StandardGeometryTop", LHCb::standard_geometry_top } } ) {}
-
-void LumiPVs_nobeamline::operator()( const LHCb::ODIN& odin, const InputVeloTracks& tracks, const InputHlt& decreport,
- const DetectorElement& lhcb ) const {
+ {KeyValue{"ODIN", ""},
+ KeyValue{"InputPVs", ""},
+ KeyValue{"InputSplitPVs1", ""},
+ KeyValue{"InputSplitPVs2", ""},
+ KeyValue{"InputHlt", ""}} ) {}
+void LumiPVs_nobeamline::operator()( const LHCb::ODIN& odin, const InputPVs& _pvs, const InputPVs& _split1,
+ const InputPVs& _split2, const InputHlt& decreport ) const {
auto tuple = nTuple( m_tupleName );
- std::vector<const LHCb::Track*> selectedtracks;
- std::vector<RecPV> vtxlist;
- std::copy_if( tracks.begin(), tracks.end(), std::back_inserter( selectedtracks ),
- []( const LHCb::Track* tr ) { return tr->hasVelo(); } );
+ // convert to vectors, as I can't find a way to simultaneously iterate over KeyedContainers
+ std::vector<RecPV*> pvs, split1, split2;
+ convert( _pvs, pvs );
+ convert( _split1, split1 );
+ convert( _split2, split2 );
+
+ StatusCode sc;
+ std::string suffix;
+ RecPV* pv;
+ for ( int i = 0; i < static_cast<int>( pvs.size() ); i++ ) {
+ sc = tuple->column( "runNumber", odin.runNumber() );
+ sc &= tuple->column( "gpsTime", static_cast<unsigned long long>( odin.gpsTime() ) );
+ sc &= tuple->column( "eventNumber", static_cast<unsigned long long>( odin.eventNumber() ) );
+ sc &= tuple->column( "eventType", odin.eventType() );
+ sc &= tuple->column( "triggerType", odin.triggerType() );
+ sc &= tuple->column( "orbitNumber", odin.orbitNumber() );
+ sc &= tuple->column( "bxType", static_cast<std::uint16_t>( odin.bunchCrossingType() ) );
+ sc &= tuple->column( "bxId", odin.bunchId() );
+ sc &= tuple->column( "step", odin.calibrationStep() );
+ for ( const auto& kvp : decreport ) { sc &= tuple->column( kvp.first, kvp.second.decision() ); }
- m_toolpv->reconstructMultiPVFromTracks( selectedtracks, vtxlist, *lhcb.geometry() ).ignore();
+ for ( int j = 0; j < 3; j++) {
+ if ( j == 0 ) {
+ suffix = "";
+ pv = pvs.at(i);
+ } else if ( j == 1 ) {
+ suffix = "_1";
+ pv = split1.at(i);
+ } else {
+ suffix = "_2";
+ pv = split2.at(i);
+ }
- auto sc = tuple->column( "runNumber", odin.runNumber() );
- sc = tuple->column( "gpsTime", static_cast<unsigned long long>( odin.gpsTime() ) );
- sc = tuple->column( "eventNumber", static_cast<unsigned long long>( odin.eventNumber() ) );
- sc = tuple->column( "eventType", odin.eventType() );
- sc = tuple->column( "triggerType", odin.triggerType() );
- sc = tuple->column( "orbitNumber", odin.orbitNumber() );
- sc = tuple->column( "bxType", static_cast<std::uint16_t>( odin.bunchCrossingType() ) );
- sc = tuple->column( "bxId", odin.bunchId() );
- sc = tuple->column( "step", odin.calibrationStep() );
- sc = tuple->column( "nVeloTracks", static_cast<unsigned long long>( tracks.size() ) );
- for ( const auto& kvp : decreport ) { sc = tuple->column( kvp.first, kvp.second.decision() ); }
+ sc &= tuple->column( "PV_ntracks"+suffix, static_cast<unsigned long long>( pv->tracks().size() ) );
+ sc &= tuple->column( "PV_ntracksbw"+suffix, count_backwards_tracks( pv ) );
+ sc &= tuple->column( "PV_chi2"+suffix, pv->chi2() );
+ sc &= tuple->column( "PV_ndof"+suffix, pv->nDoF() );
+ sc &= tuple->column( "PVX"+suffix, pv->position().x() );
+ sc &= tuple->column( "PVY"+suffix, pv->position().y() );
+ sc &= tuple->column( "PVZ"+suffix, pv->position().z() );
+ }
- sc = tuple->farray(
- { { "PV_nTracks", +[]( const RecPV& pv ) -> double { return pv.tracks().size(); } },
- { "PV_chi2ndof", +[]( const RecPV& pv ) { return pv.chi2() / pv.nDoF(); } },
- { "PVX", +[]( const RecPV& pv ) { return pv.position().x(); } },
- { "PVY", +[]( const RecPV& pv ) { return pv.position().y(); } },
- { "PVZ", +[]( const RecPV& pv ) { return pv.position().z(); } },
- { "PVR",
- +[]( const RecPV& pv ) { return sqrt( pow( pv.position().x(), 2 ) + pow( pv.position().y(), 2 ) ); } } },
- vtxlist.begin(), vtxlist.end(), "nPVs", 1000 );
+ if ( sc.isFailure() ) continue;
- sc.andThen( [&] { return tuple->write(); } ).orThrow( "Failed to fill ntuple", "LumiPV_nobeamline Alg" );
+ sc = tuple->write();
+ if ( sc.isFailure() ) throw GaudiException( "Cannot fill ntuple", "LumiPVs_nobeamline", StatusCode::FAILURE );
+ }
}
diff --git a/Rec/RecAlgs/src/PVSplitter.cpp b/Rec/RecAlgs/src/PVSplitter.cpp
new file mode 100644
index 00000000000..d3e47be6350
--- /dev/null
+++ b/Rec/RecAlgs/src/PVSplitter.cpp
@@ -0,0 +1,222 @@
+/*****************************************************************************\
+* (c) Copyright 2000-2018 CERN for the benefit of the LHCb Collaboration *
+* *
+* This software is distributed under the terms of the GNU General Public *
+* Licence version 3 (GPL Version 3), copied verbatim in the file "COPYING". *
+* *
+* In applying this licence, CERN does not waive the privileges and immunities *
+* granted to it by virtue of its status as an Intergovernmental Organization *
+* or submit itself to any jurisdiction. *
+\*****************************************************************************/
+#include "Event/RecVertex.h"
+#include "Event/Track.h"
+#include "LHCbAlgs/Transformer.h"
+#include "TrackInterfaces/IPVOfflineTool.h"
+#include "Kernel/HitPattern.h"
+#include "Detector/LHCb/DeLHCb.h"
+#include "DetDesc/DetectorElement.h"
+
+#include <random>
+#include <cmath>
+
+using Track = LHCb::Track;
+using TrackVec = std::vector<const LHCb::Track*>;
+using Tracks = LHCb::Tracks;
+using RecPV = LHCb::RecVertex;
+using RecPVs = LHCb::RecVertices;
+using OutputPVs = std::tuple<RecPVs, RecPVs, RecPVs>;
+
+namespace {
+ using base_t = LHCb::Algorithm::MultiTransformer< OutputPVs( const Tracks&, LHCb::Detector::DeLHCb const& ),
+ LHCb::Algorithm::Traits::usesConditions<LHCb::Detector::DeLHCb> >;
+
+ int random_offset( const Track* track ) {
+ // bit hacky but this should be randomish, right ...?
+ const int rnd = static_cast<int>( track->lhcbIDs()[0].lhcbID() );
+ return rnd % 2;
+ }
+
+ bool is_left_track( const Track* tr ) {
+ LHCb::HitPattern hits{tr->lhcbIDs()};
+ if ( hits.numVeloA() == hits.numVeloC() ) {
+ // pseudo-randomly allocate "mixed" tracks
+ return static_cast<bool>( random_offset( tr ) );
+ }
+ return hits.numVeloA() > hits.numVeloC();
+ }
+
+ bool is_top_track( const Track* tr ) {
+ // the old version checked for sin(phi) > 0 of each VELO sensor/hit -
+ // make sure this naive port does the same thing...
+ return ( std::sin( tr->phi() ) > 0 );
+ }
+
+ TrackVec unpack_tracks( const RecPV* pv, const bool shuffle ) {
+ const SmartRefVector<Track>& _tracks = pv->tracks();
+ TrackVec tracks;
+ for ( auto tr: _tracks ) {
+ tracks.push_back( tr );
+ }
+
+ if ( shuffle ) {
+ const int seed = static_cast<int>( tracks[0]->lhcbIDs()[0].lhcbID() );
+ std::mt19937 rng(seed);
+ std::shuffle(tracks.begin(), tracks.end(), rng);
+ }
+ return tracks;
+ }
+
+ void random_split_subsets( TrackVec& tracks1, TrackVec& tracks2, TrackVec* splittracks ) {
+ if ( tracks1.size() < 5 || tracks2.size() < 5 ) return;
+
+ // randomly divide two sets of tracks between two containers
+ int offset = random_offset( tracks1[0] );
+ for ( auto tr: tracks1 ) {
+ splittracks[offset % 2].push_back( tr );
+ offset++;
+ }
+ offset = random_offset( tracks2[0] );
+ for ( auto tr: tracks2 ) {
+ splittracks[offset % 2].push_back( tr );
+ offset++;
+ }
+ }
+}
+
+class PVSplitter : public base_t {
+
+public:
+ PVSplitter( const std::string& name, ISvcLocator* pSvcLocator )
+ : base_t ( name, pSvcLocator,
+ { KeyValue{"InputVeloTracks", ""},
+ KeyValue{"StandardGeometryTop", LHCb::standard_geometry_top} },
+ { KeyValue{"OutputPVs", "Rec/Vertex/SplitPVs/PVs"},
+ KeyValue{"OutputSplitPVs1", "Rec/Vertex/SplitPVs/SplitPVs1"},
+ KeyValue{"OutputSplitPVs2", "Rec/Vertex/SplitPVs/SplitPVs2"} }) {}
+
+ OutputPVs operator()( const Tracks&, LHCb::Detector::DeLHCb const& ) const override;
+
+private:
+ ToolHandle<IPVOfflineTool> m_toolpv{this, "PVOfflineTool", "PVOfflineTool"};
+ Gaudi::Property<std::string> m_splitMethod{this, "SplitMethod", "random"};
+
+ StatusCode split_pv( const RecPV, RecPV*, LHCb::Detector::DeLHCb const& ) const;
+
+ typedef void ( PVSplitter::*SplitTracksFunc )( const TrackVec, TrackVec* ) const;
+ SplitTracksFunc m_splitFunc;
+ void splitTracksRandomly( const TrackVec, TrackVec* ) const;
+ void splitTracksByVeloHalf( const TrackVec, TrackVec* ) const;
+ void splitTracksEqualVeloHalf( const TrackVec, TrackVec* ) const;
+ void splitTracksByTopBottom( const TrackVec, TrackVec* ) const;
+ void splitTracksEqualTopBottom( const TrackVec, TrackVec* ) const;
+
+ StatusCode initialize() override;
+};
+
+DECLARE_COMPONENT( PVSplitter )
+
+StatusCode PVSplitter::initialize() {
+ return MultiTransformer::initialize().andThen( [&] {
+ if ( m_splitMethod == "random" )
+ m_splitFunc = &PVSplitter::splitTracksRandomly;
+ else if ( m_splitMethod == "byvelohalf" )
+ m_splitFunc = &PVSplitter::splitTracksByVeloHalf;
+ else if ( m_splitMethod == "equalvelohalf" )
+ m_splitFunc = &PVSplitter::splitTracksEqualVeloHalf;
+ else if ( m_splitMethod == "bytopbottom" )
+ m_splitFunc = &PVSplitter::splitTracksByTopBottom;
+ else if ( m_splitMethod == "equaltopbottom" )
+ m_splitFunc = &PVSplitter::splitTracksEqualTopBottom;
+ else
+ throw GaudiException( "Unrecognised splitting algorithm", "PVSplitter", StatusCode::FAILURE );
+ } );
+}
+
+OutputPVs PVSplitter::operator()( const Tracks& _tracks, LHCb::Detector::DeLHCb const& lhcb ) const {
+ std::vector<RecPV> pvs;
+ std::vector<const Track*> tracks;
+ std::copy( _tracks.begin(), _tracks.end(), std::back_inserter( tracks ) );
+ m_toolpv->reconstructMultiPVFromTracks( tracks, pvs, *lhcb.geometry() ).ignore();
+
+ RecPVs retpvs;
+ RecPVs splitpvs1;
+ RecPVs splitpvs2;
+ StatusCode sc;
+ for ( auto pv: pvs ) {
+ // require 5 tracks / split pv => 10 total
+ if (pv.tracks().size() < 10) continue;
+
+ RecPV splitpvs[2];
+ sc = split_pv( pv, splitpvs, lhcb );
+ if ( sc.isFailure() ) continue;
+
+ retpvs.insert( new RecPV{pv} );
+ splitpvs1.insert( new RecPV{splitpvs[0]} );
+ splitpvs2.insert( new RecPV{splitpvs[1]} );
+ }
+
+ return OutputPVs{ std::move( retpvs ), std::move( splitpvs1 ), std::move( splitpvs2 ) };
+}
+
+StatusCode PVSplitter::split_pv( const RecPV pv, RecPV* splitpvs, LHCb::Detector::DeLHCb const& lhcb ) const {
+ TrackVec splittracks[2];
+ TrackVec tracks = unpack_tracks( &pv, true );
+ ( this->*m_splitFunc )( tracks, splittracks );
+
+ StatusCode sc;
+ for ( int i = 0; i < 2; i++ ){
+ sc = m_toolpv->reconstructSinglePVFromTracks( pv.position(), splittracks[i], splitpvs[i], *lhcb.geometry() );
+ if ( sc.isFailure() ) return sc;
+ }
+ return sc;
+}
+
+void PVSplitter::splitTracksRandomly( const TrackVec tracks, TrackVec* splittracks ) const {
+ // randomly 0 or 1, to avoid first split pv always having >= ntracks in second split pv
+ int offset = random_offset( tracks[0] );
+ for ( auto tr: tracks ) {
+ splittracks[offset % 2].push_back( tr );
+ offset++;
+ }
+}
+
+void PVSplitter::splitTracksByVeloHalf( const TrackVec tracks, TrackVec* splittracks ) const {
+ for (auto tr: tracks ) {
+ if ( is_left_track( tr ) )
+ splittracks[0].push_back( tr );
+ else
+ splittracks[1].push_back( tr );
+ }
+}
+
+void PVSplitter::splitTracksEqualVeloHalf( const TrackVec tracks, TrackVec* splittracks ) const {
+ TrackVec left, right;
+ for( auto tr: tracks ) {
+ if ( is_left_track( tr ) )
+ left.push_back( tr );
+ else
+ right.push_back( tr );
+ }
+ random_split_subsets( left, right, splittracks );
+}
+
+void PVSplitter::splitTracksByTopBottom( const TrackVec tracks, TrackVec* splittracks ) const {
+ for (auto tr: tracks ) {
+ if ( is_top_track( tr ) )
+ splittracks[0].push_back( tr );
+ else
+ splittracks[1].push_back( tr );
+ }
+}
+
+void PVSplitter::splitTracksEqualTopBottom( const TrackVec tracks, TrackVec* splittracks ) const {
+ TrackVec top, bottom;
+ for( auto tr: tracks ) {
+ if ( is_top_track( tr ) )
+ top.push_back( tr );
+ else
+ bottom.push_back( tr );
+ }
+ random_split_subsets( top, bottom, splittracks );
+}
+
--
GitLab
From b79072c7ae319067ef379edc5d555f0dceb18733 Mon Sep 17 00:00:00 2001
From: Gitlab CI <noreply@cern.ch>
Date: Thu, 13 Mar 2025 00:33:43 +0000
Subject: [PATCH 2/2] pre-commit fixes
patch generated by https://gitlab.cern.ch/lhcb/Rec/-/jobs/52659815
---
Rec/RecAlgs/src/LumiPVs_nobeamline.cpp | 58 +++++++-------
Rec/RecAlgs/src/PVSplitter.cpp | 100 ++++++++++++-------------
2 files changed, 76 insertions(+), 82 deletions(-)
diff --git a/Rec/RecAlgs/src/LumiPVs_nobeamline.cpp b/Rec/RecAlgs/src/LumiPVs_nobeamline.cpp
index f1c3367de7b..69decb073ee 100644
--- a/Rec/RecAlgs/src/LumiPVs_nobeamline.cpp
+++ b/Rec/RecAlgs/src/LumiPVs_nobeamline.cpp
@@ -21,46 +21,44 @@
#include "DetDesc/DetectorElement.h"
#include "TrackInterfaces/IPVOfflineTool.h"
-using InputPVs = LHCb::RecVertices;
-using InputHlt = LHCb::HltDecReports;
-using RecPV = LHCb::RecVertex;
+using InputPVs = LHCb::RecVertices;
+using InputHlt = LHCb::HltDecReports;
+using RecPV = LHCb::RecVertex;
namespace {
int count_backwards_tracks( RecPV* pv ) {
int nbw = 0;
- for ( auto tr: pv->tracks()) {
- if( tr->isVeloBackward() ) nbw++;
+ for ( auto tr : pv->tracks() ) {
+ if ( tr->isVeloBackward() ) nbw++;
}
return nbw;
}
// convert from LHCb::RecVertices -> std::vector<LHCb::RecVertex*>
- void convert ( const InputPVs& oldvec, std::vector<RecPV*>& newvec ) {
+ void convert( const InputPVs& oldvec, std::vector<RecPV*>& newvec ) {
std::copy( oldvec.begin(), oldvec.end(), std::back_inserter( newvec ) );
}
-}
+} // namespace
class LumiPVs_nobeamline
- : public virtual LHCb::Algorithm::Consumer<
- void( const LHCb::ODIN&, const InputPVs&, const InputPVs&, const InputPVs&, const InputHlt& ),
- LHCb::Algorithm::Traits::usesBaseAndConditions<GaudiTupleAlg>> {
+ : public virtual LHCb::Algorithm::Consumer<void( const LHCb::ODIN&, const InputPVs&, const InputPVs&,
+ const InputPVs&, const InputHlt& ),
+ LHCb::Algorithm::Traits::usesBaseAndConditions<GaudiTupleAlg>> {
private:
- Gaudi::Property<std::string> m_tupleName{this, "BeamGasTuple", "BeamGasTuple"};
+ Gaudi::Property<std::string> m_tupleName{ this, "BeamGasTuple", "BeamGasTuple" };
public:
- virtual void operator()( const LHCb::ODIN&, const InputPVs&, const InputPVs&, const InputPVs&, const InputHlt& ) const override;
+ virtual void operator()( const LHCb::ODIN&, const InputPVs&, const InputPVs&, const InputPVs&,
+ const InputHlt& ) const override;
LumiPVs_nobeamline( const std::string& name, ISvcLocator* pSvcLocator );
};
DECLARE_COMPONENT_WITH_ID( LumiPVs_nobeamline, "LumiPVs_nobeamline" );
LumiPVs_nobeamline::LumiPVs_nobeamline( const std::string& name, ISvcLocator* pSvcLocator )
: Consumer( name, pSvcLocator,
- {KeyValue{"ODIN", ""},
- KeyValue{"InputPVs", ""},
- KeyValue{"InputSplitPVs1", ""},
- KeyValue{"InputSplitPVs2", ""},
- KeyValue{"InputHlt", ""}} ) {}
+ { KeyValue{ "ODIN", "" }, KeyValue{ "InputPVs", "" }, KeyValue{ "InputSplitPVs1", "" },
+ KeyValue{ "InputSplitPVs2", "" }, KeyValue{ "InputHlt", "" } } ) {}
void LumiPVs_nobeamline::operator()( const LHCb::ODIN& odin, const InputPVs& _pvs, const InputPVs& _split1,
const InputPVs& _split2, const InputHlt& decreport ) const {
@@ -72,9 +70,9 @@ void LumiPVs_nobeamline::operator()( const LHCb::ODIN& odin, const InputPVs& _pv
convert( _split1, split1 );
convert( _split2, split2 );
- StatusCode sc;
+ StatusCode sc;
std::string suffix;
- RecPV* pv;
+ RecPV* pv;
for ( int i = 0; i < static_cast<int>( pvs.size() ); i++ ) {
sc = tuple->column( "runNumber", odin.runNumber() );
sc &= tuple->column( "gpsTime", static_cast<unsigned long long>( odin.gpsTime() ) );
@@ -87,25 +85,25 @@ void LumiPVs_nobeamline::operator()( const LHCb::ODIN& odin, const InputPVs& _pv
sc &= tuple->column( "step", odin.calibrationStep() );
for ( const auto& kvp : decreport ) { sc &= tuple->column( kvp.first, kvp.second.decision() ); }
- for ( int j = 0; j < 3; j++) {
+ for ( int j = 0; j < 3; j++ ) {
if ( j == 0 ) {
suffix = "";
- pv = pvs.at(i);
+ pv = pvs.at( i );
} else if ( j == 1 ) {
suffix = "_1";
- pv = split1.at(i);
+ pv = split1.at( i );
} else {
suffix = "_2";
- pv = split2.at(i);
+ pv = split2.at( i );
}
- sc &= tuple->column( "PV_ntracks"+suffix, static_cast<unsigned long long>( pv->tracks().size() ) );
- sc &= tuple->column( "PV_ntracksbw"+suffix, count_backwards_tracks( pv ) );
- sc &= tuple->column( "PV_chi2"+suffix, pv->chi2() );
- sc &= tuple->column( "PV_ndof"+suffix, pv->nDoF() );
- sc &= tuple->column( "PVX"+suffix, pv->position().x() );
- sc &= tuple->column( "PVY"+suffix, pv->position().y() );
- sc &= tuple->column( "PVZ"+suffix, pv->position().z() );
+ sc &= tuple->column( "PV_ntracks" + suffix, static_cast<unsigned long long>( pv->tracks().size() ) );
+ sc &= tuple->column( "PV_ntracksbw" + suffix, count_backwards_tracks( pv ) );
+ sc &= tuple->column( "PV_chi2" + suffix, pv->chi2() );
+ sc &= tuple->column( "PV_ndof" + suffix, pv->nDoF() );
+ sc &= tuple->column( "PVX" + suffix, pv->position().x() );
+ sc &= tuple->column( "PVY" + suffix, pv->position().y() );
+ sc &= tuple->column( "PVZ" + suffix, pv->position().z() );
}
if ( sc.isFailure() ) continue;
diff --git a/Rec/RecAlgs/src/PVSplitter.cpp b/Rec/RecAlgs/src/PVSplitter.cpp
index d3e47be6350..67305d0432d 100644
--- a/Rec/RecAlgs/src/PVSplitter.cpp
+++ b/Rec/RecAlgs/src/PVSplitter.cpp
@@ -8,16 +8,16 @@
* granted to it by virtue of its status as an Intergovernmental Organization *
* or submit itself to any jurisdiction. *
\*****************************************************************************/
+#include "DetDesc/DetectorElement.h"
+#include "Detector/LHCb/DeLHCb.h"
#include "Event/RecVertex.h"
#include "Event/Track.h"
+#include "Kernel/HitPattern.h"
#include "LHCbAlgs/Transformer.h"
#include "TrackInterfaces/IPVOfflineTool.h"
-#include "Kernel/HitPattern.h"
-#include "Detector/LHCb/DeLHCb.h"
-#include "DetDesc/DetectorElement.h"
-#include <random>
#include <cmath>
+#include <random>
using Track = LHCb::Track;
using TrackVec = std::vector<const LHCb::Track*>;
@@ -27,8 +27,8 @@ using RecPVs = LHCb::RecVertices;
using OutputPVs = std::tuple<RecPVs, RecPVs, RecPVs>;
namespace {
- using base_t = LHCb::Algorithm::MultiTransformer< OutputPVs( const Tracks&, LHCb::Detector::DeLHCb const& ),
- LHCb::Algorithm::Traits::usesConditions<LHCb::Detector::DeLHCb> >;
+ using base_t = LHCb::Algorithm::MultiTransformer<OutputPVs( const Tracks&, LHCb::Detector::DeLHCb const& ),
+ LHCb::Algorithm::Traits::usesConditions<LHCb::Detector::DeLHCb>>;
int random_offset( const Track* track ) {
// bit hacky but this should be randomish, right ...?
@@ -37,7 +37,7 @@ namespace {
}
bool is_left_track( const Track* tr ) {
- LHCb::HitPattern hits{tr->lhcbIDs()};
+ LHCb::HitPattern hits{ tr->lhcbIDs() };
if ( hits.numVeloA() == hits.numVeloC() ) {
// pseudo-randomly allocate "mixed" tracks
return static_cast<bool>( random_offset( tr ) );
@@ -53,15 +53,13 @@ namespace {
TrackVec unpack_tracks( const RecPV* pv, const bool shuffle ) {
const SmartRefVector<Track>& _tracks = pv->tracks();
- TrackVec tracks;
- for ( auto tr: _tracks ) {
- tracks.push_back( tr );
- }
+ TrackVec tracks;
+ for ( auto tr : _tracks ) { tracks.push_back( tr ); }
if ( shuffle ) {
- const int seed = static_cast<int>( tracks[0]->lhcbIDs()[0].lhcbID() );
- std::mt19937 rng(seed);
- std::shuffle(tracks.begin(), tracks.end(), rng);
+ const int seed = static_cast<int>( tracks[0]->lhcbIDs()[0].lhcbID() );
+ std::mt19937 rng( seed );
+ std::shuffle( tracks.begin(), tracks.end(), rng );
}
return tracks;
}
@@ -71,44 +69,43 @@ namespace {
// randomly divide two sets of tracks between two containers
int offset = random_offset( tracks1[0] );
- for ( auto tr: tracks1 ) {
+ for ( auto tr : tracks1 ) {
splittracks[offset % 2].push_back( tr );
offset++;
}
offset = random_offset( tracks2[0] );
- for ( auto tr: tracks2 ) {
+ for ( auto tr : tracks2 ) {
splittracks[offset % 2].push_back( tr );
offset++;
}
}
-}
+} // namespace
class PVSplitter : public base_t {
public:
PVSplitter( const std::string& name, ISvcLocator* pSvcLocator )
- : base_t ( name, pSvcLocator,
- { KeyValue{"InputVeloTracks", ""},
- KeyValue{"StandardGeometryTop", LHCb::standard_geometry_top} },
- { KeyValue{"OutputPVs", "Rec/Vertex/SplitPVs/PVs"},
- KeyValue{"OutputSplitPVs1", "Rec/Vertex/SplitPVs/SplitPVs1"},
- KeyValue{"OutputSplitPVs2", "Rec/Vertex/SplitPVs/SplitPVs2"} }) {}
+ : base_t( name, pSvcLocator,
+ { KeyValue{ "InputVeloTracks", "" }, KeyValue{ "StandardGeometryTop", LHCb::standard_geometry_top } },
+ { KeyValue{ "OutputPVs", "Rec/Vertex/SplitPVs/PVs" },
+ KeyValue{ "OutputSplitPVs1", "Rec/Vertex/SplitPVs/SplitPVs1" },
+ KeyValue{ "OutputSplitPVs2", "Rec/Vertex/SplitPVs/SplitPVs2" } } ) {}
OutputPVs operator()( const Tracks&, LHCb::Detector::DeLHCb const& ) const override;
private:
- ToolHandle<IPVOfflineTool> m_toolpv{this, "PVOfflineTool", "PVOfflineTool"};
- Gaudi::Property<std::string> m_splitMethod{this, "SplitMethod", "random"};
+ ToolHandle<IPVOfflineTool> m_toolpv{ this, "PVOfflineTool", "PVOfflineTool" };
+ Gaudi::Property<std::string> m_splitMethod{ this, "SplitMethod", "random" };
StatusCode split_pv( const RecPV, RecPV*, LHCb::Detector::DeLHCb const& ) const;
typedef void ( PVSplitter::*SplitTracksFunc )( const TrackVec, TrackVec* ) const;
SplitTracksFunc m_splitFunc;
- void splitTracksRandomly( const TrackVec, TrackVec* ) const;
- void splitTracksByVeloHalf( const TrackVec, TrackVec* ) const;
- void splitTracksEqualVeloHalf( const TrackVec, TrackVec* ) const;
- void splitTracksByTopBottom( const TrackVec, TrackVec* ) const;
- void splitTracksEqualTopBottom( const TrackVec, TrackVec* ) const;
+ void splitTracksRandomly( const TrackVec, TrackVec* ) const;
+ void splitTracksByVeloHalf( const TrackVec, TrackVec* ) const;
+ void splitTracksEqualVeloHalf( const TrackVec, TrackVec* ) const;
+ void splitTracksByTopBottom( const TrackVec, TrackVec* ) const;
+ void splitTracksEqualTopBottom( const TrackVec, TrackVec* ) const;
StatusCode initialize() override;
};
@@ -133,26 +130,26 @@ StatusCode PVSplitter::initialize() {
}
OutputPVs PVSplitter::operator()( const Tracks& _tracks, LHCb::Detector::DeLHCb const& lhcb ) const {
- std::vector<RecPV> pvs;
+ std::vector<RecPV> pvs;
std::vector<const Track*> tracks;
std::copy( _tracks.begin(), _tracks.end(), std::back_inserter( tracks ) );
m_toolpv->reconstructMultiPVFromTracks( tracks, pvs, *lhcb.geometry() ).ignore();
- RecPVs retpvs;
- RecPVs splitpvs1;
- RecPVs splitpvs2;
+ RecPVs retpvs;
+ RecPVs splitpvs1;
+ RecPVs splitpvs2;
StatusCode sc;
- for ( auto pv: pvs ) {
+ for ( auto pv : pvs ) {
// require 5 tracks / split pv => 10 total
- if (pv.tracks().size() < 10) continue;
+ if ( pv.tracks().size() < 10 ) continue;
RecPV splitpvs[2];
sc = split_pv( pv, splitpvs, lhcb );
if ( sc.isFailure() ) continue;
- retpvs.insert( new RecPV{pv} );
- splitpvs1.insert( new RecPV{splitpvs[0]} );
- splitpvs2.insert( new RecPV{splitpvs[1]} );
+ retpvs.insert( new RecPV{ pv } );
+ splitpvs1.insert( new RecPV{ splitpvs[0] } );
+ splitpvs2.insert( new RecPV{ splitpvs[1] } );
}
return OutputPVs{ std::move( retpvs ), std::move( splitpvs1 ), std::move( splitpvs2 ) };
@@ -164,7 +161,7 @@ StatusCode PVSplitter::split_pv( const RecPV pv, RecPV* splitpvs, LHCb::Detector
( this->*m_splitFunc )( tracks, splittracks );
StatusCode sc;
- for ( int i = 0; i < 2; i++ ){
+ for ( int i = 0; i < 2; i++ ) {
sc = m_toolpv->reconstructSinglePVFromTracks( pv.position(), splittracks[i], splitpvs[i], *lhcb.geometry() );
if ( sc.isFailure() ) return sc;
}
@@ -174,14 +171,14 @@ StatusCode PVSplitter::split_pv( const RecPV pv, RecPV* splitpvs, LHCb::Detector
void PVSplitter::splitTracksRandomly( const TrackVec tracks, TrackVec* splittracks ) const {
// randomly 0 or 1, to avoid first split pv always having >= ntracks in second split pv
int offset = random_offset( tracks[0] );
- for ( auto tr: tracks ) {
+ for ( auto tr : tracks ) {
splittracks[offset % 2].push_back( tr );
offset++;
}
}
void PVSplitter::splitTracksByVeloHalf( const TrackVec tracks, TrackVec* splittracks ) const {
- for (auto tr: tracks ) {
+ for ( auto tr : tracks ) {
if ( is_left_track( tr ) )
splittracks[0].push_back( tr );
else
@@ -191,7 +188,7 @@ void PVSplitter::splitTracksByVeloHalf( const TrackVec tracks, TrackVec* splittr
void PVSplitter::splitTracksEqualVeloHalf( const TrackVec tracks, TrackVec* splittracks ) const {
TrackVec left, right;
- for( auto tr: tracks ) {
+ for ( auto tr : tracks ) {
if ( is_left_track( tr ) )
left.push_back( tr );
else
@@ -201,7 +198,7 @@ void PVSplitter::splitTracksEqualVeloHalf( const TrackVec tracks, TrackVec* spli
}
void PVSplitter::splitTracksByTopBottom( const TrackVec tracks, TrackVec* splittracks ) const {
- for (auto tr: tracks ) {
+ for ( auto tr : tracks ) {
if ( is_top_track( tr ) )
splittracks[0].push_back( tr );
else
@@ -210,13 +207,12 @@ void PVSplitter::splitTracksByTopBottom( const TrackVec tracks, TrackVec* splitt
}
void PVSplitter::splitTracksEqualTopBottom( const TrackVec tracks, TrackVec* splittracks ) const {
- TrackVec top, bottom;
- for( auto tr: tracks ) {
- if ( is_top_track( tr ) )
- top.push_back( tr );
- else
- bottom.push_back( tr );
- }
+ TrackVec top, bottom;
+ for ( auto tr : tracks ) {
+ if ( is_top_track( tr ) )
+ top.push_back( tr );
+ else
+ bottom.push_back( tr );
+ }
random_split_subsets( top, bottom, splittracks );
}
-
--
GitLab