Commit ce356a1a authored by Rosen Matev's avatar Rosen Matev 🌲
Browse files

Merge branch 'pkoppenb-DecayTreeFitterAlg' into 'master'

DecayTreeFitterAlg to access refitted candidates

See merge request !970
parents 4601a20a 5d2974ef
Pipeline #3276552 passed with stage
in 32 seconds
......@@ -43,7 +43,7 @@ gaudi_add_module(DaVinciFilters
Phys::DaVinciInterfacesLib
Phys::DaVinciKernelLib
Phys::DaVinciTypesLib
Phys::DecayTreeFitter
Phys::DecayTreeFitterLib
Phys::LoKiArrayFunctorsLib
Phys::LoKiPhysLib
Rec::TrackInterfacesLib
......
......@@ -13,35 +13,36 @@ Phys/DecayTreeFitter
--------------------
#]=======================================================================]
gaudi_add_library(DecayTreeFitter
gaudi_add_library(DecayTreeFitterLib
SOURCES
src/Constraint.cpp
src/ConvertedPhoton.cpp
src/DecayChain.cpp
src/FitParams.cpp
src/Fitter.cpp
src/HelixUtils.cpp
src/InteractionPoint.cpp
src/InternalParticle.cpp
src/InternalRecoTrack.cpp
src/JetMomentum.cpp
src/KalmanCalculator.cpp
src/MergedConstraint.cpp
src/MissingParticle.cpp
src/ParticleBase.cpp
src/RecoComposite.cpp
src/RecoMergedPi0.cpp
src/RecoParticle.cpp
src/RecoPhoton.cpp
src/RecoResonance.cpp
src/RecoTrack.cpp
src/Resonance.cpp
src/VtxErrCode.cpp
src/lib/Constraint.cpp
src/lib/ConvertedPhoton.cpp
src/lib/DecayChain.cpp
src/lib/FitParams.cpp
src/lib/Fitter.cpp
src/lib/HelixUtils.cpp
src/lib/InteractionPoint.cpp
src/lib/InternalParticle.cpp
src/lib/InternalRecoTrack.cpp
src/lib/JetMomentum.cpp
src/lib/KalmanCalculator.cpp
src/lib/MergedConstraint.cpp
src/lib/MissingParticle.cpp
src/lib/ParticleBase.cpp
src/lib/RecoComposite.cpp
src/lib/RecoMergedPi0.cpp
src/lib/RecoParticle.cpp
src/lib/RecoPhoton.cpp
src/lib/RecoResonance.cpp
src/lib/RecoTrack.cpp
src/lib/Resonance.cpp
src/lib/VtxErrCode.cpp
LINK
PUBLIC
Gaudi::GaudiKernel
LHCb::LHCbMathLib
Phys::DaVinciTypesLib
Phys::DaVinciInterfacesLib
Rec::TrackInterfacesLib
PRIVATE
Boost::headers
......@@ -55,8 +56,19 @@ gaudi_add_library(DecayTreeFitter
Rec::TrackKernel
)
gaudi_add_module(DecayTreeFitter
SOURCES
src/component/DecayTreeFitterAlg.cpp
LINK
Gaudi::GaudiAlgLib
LHCb::LHCbKernel
LHCb::PhysEvent
DecayTreeFitterLib
LHCb::LinkerEvent
)
gaudi_add_dictionary(DecayTreeFitterDict
HEADERFILES dict/DecayTreeFitterDict.h
SELECTION dict/DecayTreeFitterDict.xml
LINK DecayTreeFitter
LINK DecayTreeFitterLib
)
/*****************************************************************************\
* (c) Copyright 2021 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 files
#include "DecayTreeFitter/Fitter.h"
#include "DetDesc/IDetectorElement.h"
#include "Event/Particle.h"
#include "GaudiAlg/Transformer.h"
#include "Kernel/DecayTree.h"
#include "Kernel/IParticleDescendants.h"
#include "Kernel/IPrintDecay.h"
#include "Kernel/ParticleID.h"
#include "Relations/Relation1D.h"
#include "TrackInterfaces/ITrackStateProvider.h"
#include <cassert>
namespace {
using Table = LHCb::Relation1D<LHCb::Particle, LHCb::Particle>;
using out_t = std::tuple<LHCb::Particles, Table>;
} // namespace
// ============================================================================
/** @class DecayTreeFitterAlg
* Algorithm that runs DecayTreeFitter on an input and stores the output.
* Can be used stand-alone as part of a selection, or in conjuction with MAP_INPUT functor
*
* @param Input Location of particles
* @returns Output Location of DecayTreeFit-ed particles
* @returns OutputRelations location of relation table between in and out particles
*
* Properties:
* - "PrintTree" = false: print decay tree (useful for debugging)
*
* Use:
* @code
* from PyConf.Algorithms import DecayTreeFitterAlg
* import Functors as F
* DTF = DecayTreeFitterAlg(Input=dimuons, MassConstraints=["J/psi(1S)"])
* DTFRelations = DTF.OutputRelations # Relations
* functor = F.MAP_INPUT(Functor=F.MASS, Relations=DTFRelations)
* @endcode
*
* @see DecayTreeFitter::Fitter
* @see FunTuple
* @see ParticleMapper_Rel
* @author Patrick Koppenburg
* @date 2021-07-19
*/
class DecayTreeFitterAlg : public Gaudi::Functional::MultiTransformer<out_t( LHCb::Particles const& )> {
public:
// ==========================================================================
/** the standard constructor
* @param name algorithm instance name
* @param pSvc service locator
*/
DecayTreeFitterAlg( const std::string& name, ISvcLocator* pSvc )
: MultiTransformer( name, pSvc, {KeyValue{"Input", ""}},
{KeyValue{"Output", ""}, KeyValue{"OutputRelations", ""}} ) {}
// ==========================================================================
/// the standard execution of the algorithm
// ==========================================================================
out_t operator()( const LHCb::Particles& parts ) const override {
auto lhcb = getDet<IDetectorElement>( m_standardGeometry_address );
if ( !lhcb ) { throw GaudiException( "Could not load geometry", name(), StatusCode::FAILURE ); }
auto& geometry = *lhcb->geometry();
// output
LHCb::Particles parts_out;
Table p2p_out;
for ( auto* p : parts ) {
DecayTreeFitter::Fitter fitter( *p, m_stateprovider ); // @todo also foresee the case with a PV given
if ( !m_massConstraints.empty() ) {
for ( const auto& C : m_massConstraints ) { fitter.setMassConstraint( LHCb::ParticleID( C ) ); }
}
fitter.fit( geometry );
LHCb::DecayTree fittedTree = fitter.getFittedTree();
if ( fittedTree.valid() ) {
m_fittedCount++;
if ( m_print ) m_printDecay->printTree( fittedTree.head(), -1 );
// fill relations. First get decay tree of input
LHCb::Particle::ConstVector parts_in; // input tree, flattened
parts_in.push_back( p ); // input mother
LHCb::Particle::ConstVector dauts_in = m_descendants->descendants( p );
parts_in.insert( parts_in.end(), dauts_in.begin(), dauts_in.end() ); // insert descendants
// search for fitted particles in flattened tree
for ( const auto inP : parts_in ) {
// 1. Find refitted partner in tree
auto* clone = fittedTree.findClone( *inP );
if ( !clone ) {
err() << "No refitted particle found for a " << inP->particleID().pid() << endmsg;
throw GaudiException( "No refitted particle found for", name(), StatusCode::FAILURE );
}
LHCb::Particle* refittedP = clone->clone(); // get fitted particle
// 2. save partner in output
LHCb::Particles::key_type key =
parts_out.insert( refittedP ); // insert, but not yet release as this would delete the clone map
const LHCb::Particle* savedP = static_cast<const LHCb::Particle*>(
parts_out.containedObject( key ) ); // get pointer to newly instered particle
// 3. relate
p2p_out.relate( inP, savedP ).ignore();
}
fittedTree.release(); // release now
} else {
m_failedCount++;
Warning( "Failed fit", StatusCode::SUCCESS, 0 ).ignore();
}
}
if ( !parts.empty() ) {
m_eventCount++;
m_candidateCount += parts.size();
m_fittedCount += parts_out.size();
}
return std::make_tuple( std::move( parts_out ), p2p_out );
}
private:
// ==========================================================================
// tools
ToolHandle<ITrackStateProvider> m_stateprovider{this, "StateProvider", "TrackStateProvider"};
ToolHandle<IParticleDescendants> m_descendants{this, "Descendants", "ParticleDescendants"};
ToolHandle<IPrintDecay> m_printDecay{this, "PrintDecayTreeTool", "PrintDecayTreeTool"};
// properties
Gaudi::Property<std::vector<int>> m_massConstraints{this, "MassConstraints"};
Gaudi::Property<std::string> m_standardGeometry_address{this, "StandardGeometryTop", "/dd/Structure/LHCb"};
// Gaudi::Property<bool> m_originVertex{this, "ConstrainToOriginVertex", false};
Gaudi::Property<bool> m_print{this, "PrintTree", false};
// counters
mutable Gaudi::Accumulators::Counter<> m_eventCount{this, "Events"};
mutable Gaudi::Accumulators::StatCounter<unsigned int> m_candidateCount{this, "Input Particles"};
mutable Gaudi::Accumulators::StatCounter<unsigned int> m_fittedCount{this, "Fitted Particles"};
mutable Gaudi::Accumulators::StatCounter<unsigned int> m_failedCount{this, "Failed fits"};
mutable Gaudi::Accumulators::StatCounter<unsigned int> m_savedCount{this, "saved Particles"};
// ==========================================================================
};
// ============================================================================
/// declare the factory (needed for instantiation)
DECLARE_COMPONENT( DecayTreeFitterAlg )
// ============================================================================
// The END
// ============================================================================
Markdown is supported
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment