diff --git a/Event/DummyProducers/CMakeLists.txt b/Event/DummyProducers/CMakeLists.txt new file mode 100644 index 0000000000000000000000000000000000000000..341cc88e16d07f2a2c9ea10e29578151bd9f1016 --- /dev/null +++ b/Event/DummyProducers/CMakeLists.txt @@ -0,0 +1,24 @@ +############################################################################### +# (c) Copyright 2019 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. # +############################################################################### +################################################################################ +# Package: DummyProducers +################################################################################ +gaudi_subdir(DummyProducers) + +gaudi_depends_on_subdirs(Event/TrackEvent) + +find_package(Boost) +find_package(ROOT) +include_directories(SYSTEM ${Boost_INCLUDE_DIRS} ${ROOT_INCLUDE_DIRS}) + +gaudi_add_module(DummyProducers + src/*.cpp + LINK_LIBRARIES TrackEvent) \ No newline at end of file diff --git a/Event/DummyProducers/src/ProducePrFittedForwardTracks.cpp b/Event/DummyProducers/src/ProducePrFittedForwardTracks.cpp new file mode 100644 index 0000000000000000000000000000000000000000..9ab650a74b24d0a32c5512dd1c32d57a2df8cd34 --- /dev/null +++ b/Event/DummyProducers/src/ProducePrFittedForwardTracks.cpp @@ -0,0 +1,28 @@ +/*****************************************************************************\ +* (c) Copyright 2019 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/GeneratePrFittedForwardTracks.h" +#include "GaudiAlg/Producer.h" + +using namespace LHCb::Pr::Fitted::Forward; +using base_t = Gaudi::Functional::Producer<Tracks()>; + +struct ProducePrFittedForwardTracks final : public base_t { + ProducePrFittedForwardTracks( std::string const& name, ISvcLocator* svcLoc ) + : base_t( name, svcLoc, KeyValue{"Output", ""} ) {} + + Tracks operator()() const override { return generate_tracks( m_nTracks, ( m_eventCount++ ).value() ); } + +private: + mutable Gaudi::Accumulators::Counter<unsigned int> m_eventCount{this, "Event"}; + Gaudi::Property<std::size_t> m_nTracks{this, "NumberToGenerate", 100, "Number of objects to generate"}; +}; + +DECLARE_COMPONENT( ProducePrFittedForwardTracks ) \ No newline at end of file diff --git a/Event/RecEvent/CMakeLists.txt b/Event/RecEvent/CMakeLists.txt index 62d17e52465de5915ff5891656b1bf1274427977..a89ee4953752d5e03c5843b158352e90beb4e754 100644 --- a/Event/RecEvent/CMakeLists.txt +++ b/Event/RecEvent/CMakeLists.txt @@ -36,5 +36,5 @@ gaudi_add_dictionary(RecEvent dict/dictionary.h dict/selection.xml LINK_LIBRARIE gaudi_add_unit_test(TestMuonZipInfrastructure tests/src/TestMuonZipInfrastructure.cpp - LINK_LIBRARIES LHCbMathLib + LINK_LIBRARIES LHCbMathLib TrackEvent TYPE Boost) diff --git a/Event/RecEvent/tests/src/TestMuonZipInfrastructure.cpp b/Event/RecEvent/tests/src/TestMuonZipInfrastructure.cpp index 8d3fbadf51b234bf4dfce7eafc5e1e7c787cc6aa..d2c3cf2caf287855d3045b141128574e4cda0c48 100644 --- a/Event/RecEvent/tests/src/TestMuonZipInfrastructure.cpp +++ b/Event/RecEvent/tests/src/TestMuonZipInfrastructure.cpp @@ -8,6 +8,7 @@ * granted to it by virtue of its status as an Intergovernmental Organization * * or submit itself to any jurisdiction. * \*****************************************************************************/ +#include "Event/GeneratePrFittedForwardTracks.h" #include "Event/MuonPID_v2.h" #include "Event/PrFittedForwardTracks.h" #include "Event/PrIterableFittedForwardTracks.h" @@ -34,31 +35,6 @@ using Tracks = LHCb::Pr::Fitted::Forward::Tracks; using MuonPIDs = LHCb::Pr::Muon::PIDs; -Tracks make_tracks( size_t ntracks ) { - using dType = SIMDWrapper::scalar::types; - using F = dType::float_v; - using I = dType::int_v; - std::mt19937 gen( 42 ); // Random engine with fixed seed - std::normal_distribution<float> xy_dist{0.f, 0.1f}, z_dist{0.f, 1.f}, txy_dist{0.f, 5e-3f}; - Tracks tracks{nullptr}; - for ( unsigned int i = 0; i < ntracks; ++i ) { - F x{xy_dist( gen )}, y{xy_dist( gen )}, z{z_dist( gen )}, tx{txy_dist( gen )}, ty{txy_dist( gen )}; - Vec3<F> beam_pos{x, y, z}, beam_dir{tx, ty, 1.f}, - covX{1e-5f /* x_err**2 */, 0.f /* cov(x, tx) */, 1e-8 /* tx_err**2 */}, - covY{1e-5f /* y_err**2 */, 0.f /* cov(y, ty) */, 1e-8 /* ty_err**2 */}; - tracks.store_trackFT( i, I{0} ); // Ancestor index. - tracks.store_QoP( i, F{1.f / 5e4} ); // q/p for q = +1, p = 50 GeV - tracks.store_beamStatePos( i, beam_pos ); - tracks.store_beamStateDir( i, beam_dir ); - tracks.store_chi2( i, F{1.f} ); // chi2/dof - tracks.store_chi2nDof( i, I{42} ); // chi2 dof - tracks.store_covX( i, covX ); - tracks.store_covY( i, covY ); - ++tracks.size(); - } - return tracks; -} - BOOST_AUTO_TEST_CASE( test_emplace_back ) { using dType_scalar = SIMDWrapper::scalar::types; LHCb::Pr::Muon::PIDs pids{Zipping::generateZipIdentifier()}; @@ -111,7 +87,7 @@ std::tuple<size_t, MuonPIDs> create_muonids( Tracks const& tracks ) { BOOST_AUTO_TEST_CASE( test_transform_from_tracks ) { unsigned int ntracks = 999; - auto tracks = make_tracks( ntracks ); + auto tracks = LHCb::Pr::Fitted::Forward::generate_tracks( ntracks ); auto [nmuons, muonPIDs] = create_muonids( tracks ); // Check if counted correctly using vectorized version, compiler flag dependent. @@ -123,7 +99,7 @@ BOOST_AUTO_TEST_CASE( test_transform_from_tracks ) { BOOST_AUTO_TEST_CASE( test_constructor_from_tracks_and_lambda ) { unsigned int ntracks = 999; - auto tracks = make_tracks( ntracks ); + auto tracks = LHCb::Pr::Fitted::Forward::generate_tracks( ntracks ); auto iterable_tracks_scalar = LHCb::Pr::make_zip<SIMDWrapper::InstructionSet::Scalar, true>( tracks ); auto ismuon = []( int t ) { return ( t % 3 == 0 ); }; @@ -144,7 +120,7 @@ BOOST_AUTO_TEST_CASE( test_constructor_from_tracks_and_lambda ) { BOOST_AUTO_TEST_CASE( test_zipping_of_tracks_and_muonids ) { unsigned int ntracks = 999; - auto tracks = make_tracks( ntracks ); + auto tracks = LHCb::Pr::Fitted::Forward::generate_tracks( ntracks ); auto [nmuons, muonPIDs] = create_muonids( tracks ); // Test zipping @@ -182,7 +158,7 @@ BOOST_AUTO_TEST_CASE( test_timing_of_ismuon ) { std::cout << "Test timing of scalar and vectorised loops" << std::endl; unsigned int ntracks = 999; - auto tracks = make_tracks( ntracks ); + auto tracks = LHCb::Pr::Fitted::Forward::generate_tracks( ntracks ); auto [nmuons, muonPIDs] = create_muonids( tracks ); auto iterable_muonPIDs = LHCb::Pr::make_zip( muonPIDs ); diff --git a/Event/TrackEvent/Event/GeneratePrFittedForwardTracks.h b/Event/TrackEvent/Event/GeneratePrFittedForwardTracks.h new file mode 100644 index 0000000000000000000000000000000000000000..c46acf94099de8ab5b8df48c38a7b29d010d49f3 --- /dev/null +++ b/Event/TrackEvent/Event/GeneratePrFittedForwardTracks.h @@ -0,0 +1,20 @@ +/*****************************************************************************\ +* (c) Copyright 2019 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. * +\*****************************************************************************/ +#pragma once +#include "Event/PrFittedForwardTracks.h" + +namespace LHCb::Pr::Fitted::Forward { + /** Helper for tests and benchmarks that generates a + * LHCb::Pr::Fitted::Forward::Tracks container that is populated with + * fake-but-sane data. + */ + Tracks generate_tracks( std::size_t nTracks, unsigned int seed = 0 ); +} // namespace LHCb::Pr::Fitted::Forward \ No newline at end of file diff --git a/Event/TrackEvent/Event/PrIterableFittedForwardTracks.h b/Event/TrackEvent/Event/PrIterableFittedForwardTracks.h index a2d7e06a19114fc7c6f7f989ec1249f4890b32a7..4b2ea11fcfe6bac03d7db9fad45bd93271ea541d 100644 --- a/Event/TrackEvent/Event/PrIterableFittedForwardTracks.h +++ b/Event/TrackEvent/Event/PrIterableFittedForwardTracks.h @@ -120,6 +120,7 @@ namespace LHCb::Pr::Fitted::Forward { } } + auto nDoF() const { return cast( this->m_tracks->template chi2nDof<IType>( this->offset() ) ); } auto chi2PerDoF() const { return cast( this->m_tracks->template chi2<FType>( this->offset() ) ); } auto lhcbIDs( LHCb::Pr::Velo::Hits const& velo_hits ) const { diff --git a/Event/TrackEvent/Event/PrProxyHelpers.h b/Event/TrackEvent/Event/PrProxyHelpers.h index d963d42df7c750af7a4c5aa11f0fa5559bdb4893..e841b23da9d9c445e197924a98ea2330a1014985 100644 --- a/Event/TrackEvent/Event/PrProxyHelpers.h +++ b/Event/TrackEvent/Event/PrProxyHelpers.h @@ -47,10 +47,10 @@ namespace LHCb::Pr::detail { : m_position{std::move( position )}, m_slopes{std::move( slopes )} {} PosType const& position() const { return m_position; } SlopeType const& slopes() const { return m_slopes; } - auto x() const { return position().x; } - auto y() const { return position().y; } - auto z() const { return position().z; } - auto tx() const { return slopes().x; } - auto ty() const { return slopes().y; } + auto x() const { return position().X(); } + auto y() const { return position().Y(); } + auto z() const { return position().Z(); } + auto tx() const { return slopes().X(); } + auto ty() const { return slopes().Y(); } }; } // namespace LHCb::Pr::detail \ No newline at end of file diff --git a/Event/TrackEvent/src/GeneratePrFittedForwardTracks.cpp b/Event/TrackEvent/src/GeneratePrFittedForwardTracks.cpp new file mode 100644 index 0000000000000000000000000000000000000000..5028ca44d9152431ddee69fd097d26968cd1cb65 --- /dev/null +++ b/Event/TrackEvent/src/GeneratePrFittedForwardTracks.cpp @@ -0,0 +1,128 @@ +/*****************************************************************************\ +* (c) Copyright 2019 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/GeneratePrFittedForwardTracks.h" + +#include <TMath.h> + +#include <random> + +LHCb::Pr::Fitted::Forward::Tracks LHCb::Pr::Fitted::Forward::generate_tracks( std::size_t nTracks, unsigned int seed ) { + using dType = SIMDWrapper::scalar::types; + using F = dType::float_v; + using I = dType::int_v; + // Try to ensure that: + // - increasing nTracks from N -> M does not affect the first N tracks + // - adding a new field to the output would not change the existing fields + // Achieve this by using a bunch of different random engines with fixed seeds + seed *= 30; // avoid reusing a seed in event n-1 in event n, giving some spare slots if this needs to be made more + // sophisticated in future + std::mt19937 gen_x{seed + 1}, gen_y{seed + 2}, gen_z{seed + 3}, gen_phi{seed + 4}, gen_eta{seed + 5}, + gen_cov_x_x{seed + 6}, gen_corr_x_tx{seed + 7}, gen_cov_tx_tx{seed + 8}, gen_qoverp{seed + 9}, + gen_ndof{seed + 10}, gen_chi2{seed + 11}; + // The various parametrisations were fitted to the output of the VeloKalman + // algorithm in the default Moore configuration in early November 2019. The + // parametrisations were chosen based on the available distributions in the + // standard library. They aim to reproduce the gross features of the + // distributions without worrying too much about capturing the details. + + // Helpers for sampling (x, y) + auto make_xy_maker = [part_choice = std::uniform_real_distribution<float>{0.f, 1.f}, + core_dist = std::normal_distribution<float>{0.f, 4.51e-2f}, + tail_dist = std::cauchy_distribution<float>{0.f, 9.78e-1f}]() mutable { + return [&]( auto&& gen ) { + for ( ;; ) { + auto x = part_choice( gen ) > 5e-2 ? core_dist( gen ) : tail_dist( gen ); + if ( std::abs( x ) < 10.f ) return x; + } + }; + }; + auto get_x = make_xy_maker(); + auto get_y = make_xy_maker(); + // Helper for sampling z + auto get_z = [part_choice = std::uniform_real_distribution<float>{0.f, 1.f}, + core = std::normal_distribution<float>{0.f, 64.f}, + tail = std::normal_distribution<float>{313.f, 200.f}]( auto&& gen ) mutable { + for ( ;; ) { + auto z = part_choice( gen ) > 1.3e-2 ? core( gen ) : tail( gen ); + if ( z > -400.f && z < 900.f ) return z; + } + }; + auto get_phi = [sign_gen = std::uniform_real_distribution<float>{0.f, 1.f}, + gaus = std::normal_distribution<float>{TMath::PiOver2(), 1.5f}]( auto&& gen ) mutable { + auto sign = sign_gen( gen ) > 0.5f ? 1.f : -1.f; + return sign * gaus( gen ); + }; + auto get_eta = [dist = std::gamma_distribution<float>{5.f, 0.4f}]( auto&& gen ) mutable { + for ( ;; ) { + auto eta = 1.4f + dist( gen ); + if ( eta < 5.f ) return eta; + } + }; + auto get_qoverp = [sign_choice = std::uniform_real_distribution<float>{0.f, 1.f}, + dist = std::lognormal_distribution<float>{-2.1e-1f, 0.986f}]( auto&& gen ) mutable { + for ( ;; ) { + auto abs_qoverp = 1e-4f * dist( gen ); + if ( abs_qoverp > 1e-6 && abs_qoverp < 4.5e-4 ) return ( sign_choice( gen ) > 0.5 ? 1.f : -1.f ) * abs_qoverp; + } + }; + auto get_ndof = [dist = std::poisson_distribution<int>{36}]( auto&& gen ) mutable { + for ( ;; ) { + auto ndof = dist( gen ) - 22; + if ( ndof >= 2 && ndof % 2 == 0 ) return ndof; + } + }; + auto get_cov_tx_tx = [dist = std::gamma_distribution<float>{1.3f, 1.9e-2f}]( auto&& gen ) mutable { + return 1e-6f * ( 1e-3f + dist( gen ) ); + }; + auto get_cov_x_tx = [dist = std::lognormal_distribution<float>{-3.3, 0.555f}]( auto&& gen, float cov_x_x, + float cov_tx_tx ) mutable { + for ( ;; ) { + auto corr_x_tx = -0.999f + dist( gen ); + if ( corr_x_tx < +0.999f ) return corr_x_tx * std::sqrt( cov_x_x * cov_tx_tx ); + } + }; + auto get_cov_x_x = [dist = std::gamma_distribution<float>{2.5f, 4e-4f}]( auto&& gen ) mutable { return dist( gen ); }; + auto get_chi2 = [dist = std::chi_squared_distribution<float>{}]( auto&& gen, int ndof ) mutable { + dist.param( std::chi_squared_distribution<float>::param_type( ndof ) ); + return dist( gen ); + }; + + Tracks tracks{/* ancestor relations */ nullptr}; + for ( int i = 0; i < int( nTracks ); ++i ) { + // Generate (x, y, z) directly + F x{get_x( gen_x )}, y{get_y( gen_y )}, z{get_z( gen_z )}; + // Generate (eta, phi) then derive (tx, ty) + float phi{get_phi( gen_phi )}, eta{get_eta( gen_eta )}; + float tan_theta{( 2.f * std::exp( eta ) ) / ( std::exp( 2.f * eta ) - 1.f )}; + F tx{tan_theta * std::cos( phi )}, ty{tan_theta * std::sin( phi )}; + // Generate cov(x, x) and cov(tx, tx) + float cov_x_x{get_cov_x_x( gen_cov_x_x )}, cov_tx_tx{get_cov_tx_tx( gen_cov_tx_tx )}; + // Generate corr(x, tx) and calculate cov(x, tx) + F cov_x_tx{get_cov_x_tx( gen_corr_x_tx, cov_x_x, cov_tx_tx )}; + // Generate the number of degrees of freedom + int ndof{get_ndof( gen_ndof )}; + // Generate the chi2 itself + float chi2{get_chi2( gen_chi2, ndof )}; + Vec3<F> beam_pos{x, y, z}, beam_dir{tx, ty, 1.f}, + covX{F{cov_x_x} /* x_err**2 */, cov_x_tx /* cov(x, tx) */, F{cov_tx_tx} /* tx_err**2 */}, + covY{F{cov_x_x} /* y_err**2 */, cov_x_tx /* cov(y, ty) */, F{cov_tx_tx} /* ty_err**2 */}; + tracks.store_trackFT( i, I{i} ); // Ancestor index. + tracks.store_QoP( i, F{get_qoverp( gen_qoverp )} ); + tracks.store_beamStatePos( i, beam_pos ); + tracks.store_beamStateDir( i, beam_dir ); + tracks.store_chi2( i, F{chi2 / ndof} ); // chi2/dof + tracks.store_chi2nDof( i, I{ndof} ); // chi2 dof + tracks.store_covX( i, covX ); + tracks.store_covY( i, covY ); + ++tracks.size(); + } + return tracks; +} \ No newline at end of file diff --git a/Kernel/SOAExtensions/CMakeLists.txt b/Kernel/SOAExtensions/CMakeLists.txt index 10015c24d59338c0b0ebcb7c56c7af7c8cc656f7..e9007167f1a608a43f8660679c82f790d43574cc 100644 --- a/Kernel/SOAExtensions/CMakeLists.txt +++ b/Kernel/SOAExtensions/CMakeLists.txt @@ -16,8 +16,7 @@ gaudi_subdir(SOAExtensions) gaudi_depends_on_subdirs( Kernel/SOAContainer Kernel/LHCbMath - GaudiAlg - ) + GaudiAlg) gaudi_install_headers(SOAExtensions) # as transitive dependency, projects that include SOAExtensions must also include Kernel/SOAContainer