Belle II Software development
RecoTrack Class Reference

representation of all charged final states FIXME rename since this name is taken in tracking More...

#include <RecoTrack.h>

Inheritance diagram for RecoTrack:
RecoParticle ParticleBase

Public Types

enum  TFParticleType {
  kInteractionPoint ,
  kOrigin ,
  kComposite ,
  kRecoResonance ,
  kInternalParticle ,
  kRecoTrack ,
  kResonance ,
  kRecoNeutral ,
  kMissingParticle
}
 particle types More...
 
typedef std::vector< Constraintconstraintlist
 alias
 
typedef std::vector< std::pair< const ParticleBase *, int > > indexmap
 alias
 

Public Member Functions

 RecoTrack (Belle2::Particle *bc, const ParticleBase *mother)
 constructor
 
virtual ~RecoTrack ()
 destructor
 
virtual ErrCode initParticleWithMother (FitParams &fitparams) override
 init with mother particle (replacing initPar2)
 
virtual ErrCode initMotherlessParticle (FitParams &fitparams) override
 init without mother particle
 
ErrCode initCovariance (FitParams &fitparams) const override
 init covariance matrix of this particle constraint
 
ErrCode updFltToMotherCopy (const FitParams *fitparams)
 update m_flt
 
virtual ErrCode projectRecoConstraint (const FitParams &, Projection &) const override
 project the constraint (calculate residuals)
 
ErrCode updateParams (double flt)
 updated the cached parameters
 
virtual int dimM () const override
 dimension (5)
 
virtual int type () const override
 type of the constraint
 
virtual int nFinalChargedCandidates () const override
 number of final charged candidates
 
virtual void addToConstraintList (constraintlist &alist, int depth) const override
 add to the list of constraints
 
ErrCode updFltToMother (const FitParams &fitparams)
 update flight length to mother
 
void setFlightLength (double flt)
 setter for the flight length
 
virtual std::string parname (int index) const override
 name
 
virtual int dim () const override
 this here sets the size in the state vector we are only interested in the momenta of photons, neutral hadrons, and tracks as the position were the track ends or the cluster is, is not relevant for physics
 
virtual int momIndex () const override
 get momentum index
 
virtual bool hasEnergy () const override
 has an energy in the statevector?
 
virtual ErrCode projectConstraint (Constraint::Type, const FitParams &, Projection &) const override
 abstract abstract projection
 
virtual void updateIndex (int &offset)
 this sets the index for momentum, position, etc.
 
const ParticleBaselocate (Belle2::Particle *particle) const
 get particle base from basf2 particle
 
Belle2::Particleparticle () const
 get basf2 particle
 
int index () const
 get index
 
const ParticleBasemother () const
 getMother() / hasMother()
 
virtual ErrCode projectGeoConstraint (const FitParams &, Projection &) const
 project geometrical constraint
 
virtual ErrCode projectMassConstraintParticle (const FitParams &, Projection &) const
 project mass constraint using the particles parameters
 
virtual ErrCode projectMassConstraintDaughters (const FitParams &, Projection &) const
 project mass constraint using the parameters of the daughters
 
virtual ErrCode projectMassConstraint (const FitParams &, Projection &) const
 project mass constraint abstract
 
virtual void forceP4Sum (FitParams &) const
 force p4 sum conservation all along the tree
 
virtual int posIndex () const
 get vertex index (in statevector!)
 
virtual int tauIndex () const
 get tau index
 
virtual bool hasPosition () const
 get false
 
int eneIndex () const
 get energy index
 
virtual double chiSquare (const FitParams &) const
 get chi2
 
int charge () const
 get charge
 
virtual ParticleBaseaddDaughter (Belle2::Particle *, const ConstraintConfiguration &config, bool forceFitAll=false)
 add daughter
 
virtual void removeDaughter (const ParticleBase *pb)
 remove daughter
 
virtual void retrieveIndexMap (indexmap &anindexmap) const
 get index map
 
void setMother (const ParticleBase *m)
 set mother
 
void collectVertexDaughters (std::vector< ParticleBase * > &particles, int posindex)
 get vertex daughters
 

Static Public Member Functions

static ParticleBasecreateParticle (Belle2::Particle *particle, const ParticleBase *mother, const ConstraintConfiguration &config, bool forceFitAll=false)
 create the according treeFitter particle obj for a basf2 particle type
 
static ParticleBasecreateOrigin (Belle2::Particle *daughter, const ConstraintConfiguration &config, bool forceFitAll)
 create a custom origin particle or a beamspot
 

Protected Types

typedef std::vector< ParticleBase * > ParticleContainer
 just an alias
 

Protected Member Functions

ErrCode initTau (FitParams &par) const
 initialises tau as a length
 
void setIndex (int i)
 set Index (in statevector)
 

Static Protected Member Functions

static bool isAResonance (Belle2::Particle *particle)
 controls if a particle is treated as a resonance(lifetime=0) or a particle that has a finite lifetime.
 

Protected Attributes

Belle2::Particlem_particle
 pointer to framework type
 
const ParticleBasem_mother
 motherparticle
 
std::vector< ParticleBase * > m_daughters
 daughter container
 
bool m_isStronglyDecayingResonance
 decay length less than 1 micron
 
const ConstraintConfigurationm_config
 has all the constraint config
 

Private Attributes

double m_bfield
 B field along z.
 
const Belle2::TrackFitResultm_trackfit
 trackfit result from reconstruction
 
bool m_cached
 flag to mark the particle as initialised
 
double m_flt
 helix arc length at vertex
 
Eigen::Matrix< double, 1, 5 > m_params
 column vector to store the measurement
 
Eigen::Matrix< double, 5, 5 > m_covariance
 only lower triangle filled!
 
const float m_momentumScalingFactor
 scale the momenta by this correction factor
 
int m_index
 index
 
std::string m_name
 name
 

Detailed Description

representation of all charged final states FIXME rename since this name is taken in tracking

Definition at line 18 of file RecoTrack.h.

Member Typedef Documentation

◆ constraintlist

typedef std::vector<Constraint> constraintlist
inherited

alias

Definition at line 48 of file ParticleBase.h.

◆ indexmap

typedef std::vector< std::pair<const ParticleBase*, int> > indexmap
inherited

alias

Definition at line 51 of file ParticleBase.h.

◆ ParticleContainer

typedef std::vector<ParticleBase*> ParticleContainer
protectedinherited

just an alias

Definition at line 174 of file ParticleBase.h.

Member Enumeration Documentation

◆ TFParticleType

enum TFParticleType
inherited

particle types

Definition at line 27 of file ParticleBase.h.

27 {kInteractionPoint,
28 kOrigin,
29 kComposite,
30 kRecoResonance,
31 kInternalParticle,
32 kRecoTrack,
33 kResonance,
34 kRecoNeutral,
35 kMissingParticle
36 };

Constructor & Destructor Documentation

◆ RecoTrack()

RecoTrack ( Belle2::Particle * bc,
const ParticleBase * mother )

constructor

Definition at line 23 of file RecoTrack.cc.

23 :
24 RecoParticle(particle, mother),
25 m_bfield(0),
26 m_trackfit(particle->getTrackFitResult()),
27 m_cached(false),
28 m_flt(0),
29 m_params(5),
30 m_covariance(5, 5),
31 m_momentumScalingFactor(particle->getEffectiveMomentumScale())
32 {
33 m_bfield = Belle2::BFieldManager::getFieldInTesla(ROOT::Math::XYZVector(0, 0, 0)).Z(); //Bz in Tesla
34 m_covariance = Eigen::Matrix<double, 5, 5>::Zero(5, 5);
35 }
static ROOT::Math::XYZVector getFieldInTesla(const ROOT::Math::XYZVector &pos)
return the magnetic field at a given position in Tesla.

◆ ~RecoTrack()

virtual ~RecoTrack ( )
inlinevirtual

destructor

Definition at line 25 of file RecoTrack.h.

25{};

Member Function Documentation

◆ addDaughter()

ParticleBase * addDaughter ( Belle2::Particle * cand,
const ConstraintConfiguration & config,
bool forceFitAll = false )
virtualinherited

add daughter

Definition at line 63 of file ParticleBase.cc.

64 {
65 auto newDaughter = ParticleBase::createParticle(cand, this, config, forceFitAll);
66 m_daughters.push_back(newDaughter);
67 return m_daughters.back();
68 }

◆ addToConstraintList()

virtual void addToConstraintList ( constraintlist & alist,
int depth ) const
inlineoverridevirtual

add to the list of constraints

Implements ParticleBase.

Definition at line 55 of file RecoTrack.h.

56 {
57 alist.push_back(Constraint(this, Constraint::track, depth, dimM(), 1)) ;
58 }

◆ charge()

int charge ( ) const
inlineinherited

get charge

Definition at line 140 of file ParticleBase.h.

141 {
142 if (m_particle->getPDGCode()) {
143 double fltcharge = m_particle->getCharge();
144 return fltcharge < 0 ? int(fltcharge - 0.5) : int(fltcharge + 0.5);
145 } else {
146 return m_particle->getCharge() > 0 ? 1 : (m_particle->getCharge() < 0 ? -1 : 0);
147 }
148 }

◆ chiSquare()

double chiSquare ( const FitParams & fitparams) const
virtualinherited

get chi2

Definition at line 229 of file ParticleBase.cc.

230 {
231 double rc = 0;
232 for (auto* daughter : m_daughters) {
233 rc += daughter->chiSquare(fitparams);
234 }
235 return rc;
236 }

◆ collectVertexDaughters()

void collectVertexDaughters ( std::vector< ParticleBase * > & particles,
int posindex )
inherited

get vertex daughters

Definition at line 154 of file ParticleBase.cc.

155 {
156 if (mother() && mother()->posIndex() == posindex) {
157 particles.push_back(this);
158 }
159
160 for (auto* daughter : m_daughters) {
161 daughter->collectVertexDaughters(particles, posindex);
162 }
163 }

◆ createOrigin()

ParticleBase * createOrigin ( Belle2::Particle * daughter,
const ConstraintConfiguration & config,
bool forceFitAll )
staticinherited

create a custom origin particle or a beamspot

Definition at line 91 of file ParticleBase.cc.

96 {
97 return new Origin(daughter, config, forceFitAll);
98 }

◆ createParticle()

ParticleBase * createParticle ( Belle2::Particle * particle,
const ParticleBase * mother,
const ConstraintConfiguration & config,
bool forceFitAll = false )
staticinherited

create the according treeFitter particle obj for a basf2 particle type

Definition at line 100 of file ParticleBase.cc.

102 {
103 ParticleBase* rc = nullptr;
104
105 if (!mother) { // If there is no mother, this is the 'head of tree' particle (is never a resonance)
106 rc = new InternalParticle(particle, nullptr, config, forceFitAll);
107 } else if (particle->hasExtraInfo("bremsCorrected") // Has Bremsstrahlungs-recovery
108 && particle->getExtraInfo("bremsCorrected") != 0) { // and gammas are attached
109 rc = new Composite(particle, mother, config, true);
110 // if no gamma is attached, it is treated as a RecoTrack
111 } else if (particle->hasExtraInfo("treeFitterTreatMeAsInvisible")
112 && particle->getExtraInfo("treeFitterTreatMeAsInvisible") == 1) { // dummy particles with invisible flag
113 rc = new RecoResonance(particle, mother, config);
114 } else if (particle->getTrack()) { // external reconstructed track
115 rc = new RecoTrack(particle, mother);
116 } else if (particle->getParticleSource() == Belle2::Particle::EParticleSourceObject::c_ECLCluster
117 or particle->getParticleSource() == Belle2::Particle::EParticleSourceObject::c_KLMCluster) {
118 // external reconstructed neutral final-state particle
119 rc = new RecoNeutral(particle, mother);
120 } else if (particle->getMdstArrayIndex()) { // external composite e.g. V0
121 rc = new InternalParticle(particle, mother, config, forceFitAll);
122 } else { // 'internal' particles
123 if (isAResonance(particle)) {
124 rc = new Resonance(particle, mother, config, forceFitAll);
125 } else {
126 rc = new InternalParticle(particle, mother, config, forceFitAll);
127 }
128 }
129 return rc;
130 }

◆ dim()

virtual int dim ( ) const
inlineoverridevirtualinherited

this here sets the size in the state vector we are only interested in the momenta of photons, neutral hadrons, and tracks as the position were the track ends or the cluster is, is not relevant for physics

Implements ParticleBase.

Reimplemented in RecoNeutral.

Definition at line 39 of file RecoParticle.h.

39{ return 3; }

◆ dimM()

virtual int dimM ( ) const
inlineoverridevirtual

dimension (5)

Implements RecoParticle.

Definition at line 46 of file RecoTrack.h.

46{ return 5 ; }

◆ eneIndex()

int eneIndex ( ) const
inlineinherited

get energy index

Definition at line 134 of file ParticleBase.h.

134{ return hasEnergy() ? momIndex() + 3 : -1 ; }

◆ forceP4Sum()

virtual void forceP4Sum ( FitParams & ) const
inlinevirtualinherited

force p4 sum conservation all along the tree

Reimplemented in InternalParticle.

Definition at line 112 of file ParticleBase.h.

112{} ;

◆ hasEnergy()

virtual bool hasEnergy ( ) const
inlineoverridevirtualinherited

has an energy in the statevector?

Reimplemented from ParticleBase.

Reimplemented in RecoNeutral.

Definition at line 45 of file RecoParticle.h.

45{ return false; }

◆ hasPosition()

virtual bool hasPosition ( ) const
inlinevirtualinherited

get false

Reimplemented in Composite, InternalParticle, and Resonance.

Definition at line 131 of file ParticleBase.h.

131{ return false ; }

◆ index()

int index ( ) const
inlineinherited

get index

Definition at line 91 of file ParticleBase.h.

91{ return m_index ; }

◆ initCovariance()

ErrCode initCovariance ( FitParams & fitparams) const
overridevirtual

init covariance matrix of this particle constraint

Reimplemented from ParticleBase.

Definition at line 56 of file RecoTrack.cc.

57 {
58 // we only need a rough estimate of the covariance
59 TMatrixFSym p4Err = particle()->getMomentumErrorMatrix();
60 const int momindex = momIndex();
61
62 for (int row = 0; row < 3; ++row) {
63 fitparams.getCovariance()(momindex + row, momindex + row) = 1000 * p4Err[row][row];
64 }
65
66 return ErrCode(ErrCode::Status::success);
67 }

◆ initMotherlessParticle()

ErrCode initMotherlessParticle ( FitParams & fitparams)
overridevirtual

init without mother particle

Reimplemented from RecoParticle.

Definition at line 51 of file RecoTrack.cc.

52 {
53 return ErrCode(ErrCode::Status::success);
54 }

◆ initParticleWithMother()

ErrCode initParticleWithMother ( FitParams & fitparams)
overridevirtual

init with mother particle (replacing initPar2)

Implements ParticleBase.

Definition at line 37 of file RecoTrack.cc.

38 {
39 //initPar2
40 if (m_flt == 0) {
41 const_cast<RecoTrack*>(this)->updFltToMother(fitparams);
42 }
43 ROOT::Math::XYZVector recoP = m_trackfit->getHelix(m_momentumScalingFactor).getMomentumAtArcLength2D(m_flt, m_bfield);
44 const int momindex = momIndex();
45 fitparams.getStateVector()(momindex) = recoP.X();
46 fitparams.getStateVector()(momindex + 1) = recoP.Y();
47 fitparams.getStateVector()(momindex + 2) = recoP.Z();
48 return ErrCode(ErrCode::Status::success);
49 }

◆ initTau()

ErrCode initTau ( FitParams & par) const
protectedinherited

initialises tau as a length

Definition at line 428 of file ParticleBase.cc.

429 {
430 const int tauindex = tauIndex();
431 if (tauindex >= 0 && hasPosition()) {
432
433 const int posindex = posIndex();
434 const int mother_ps_index = mother()->posIndex();
435 const int dim = m_config->m_originDimension; // TODO can we configure this to be particle specific?
436
437 // tau has different meaning depending on the dimension of the constraint
438 // 2-> use x-y projection
439 const Eigen::Matrix < double, 1, -1, 1, 1, 3 > vertex_dist =
440 fitparams.getStateVector().segment(posindex, dim) - fitparams.getStateVector().segment(mother_ps_index, dim);
441 const Eigen::Matrix < double, 1, -1, 1, 1, 3 >
442 mom = fitparams.getStateVector().segment(posindex, dim);
443
444 // if an intermediate vertex is not well defined by a track or so it will be initialised with 0
445 // same for the momentum of for example B0, it might be initialised with 0
446 // in those cases use pdg value
447 const double mom_norm = mom.norm();
448 const double dot = std::abs(vertex_dist.dot(mom));
449 const double tau = dot / mom_norm;
450 if (0 == mom_norm || 0 == dot) {
451 const double mass = m_particle->getPDGMass();
452 if (mass > 0)
453 fitparams.getStateVector()(tauindex) = m_particle->getPDGLifetime() * 1e9 * Belle2::Const::speedOfLight / mass;
454 else
455 fitparams.getStateVector()(tauindex) = 0;
456 } else {
457 fitparams.getStateVector()(tauindex) = tau;
458 }
459 }
460
461 return ErrCode(ErrCode::Status::success);
462 }
static const double speedOfLight
[cm/ns]
Definition Const.h:695
T dot(GeneralVector< T > a, GeneralVector< T > b)
dot product of two general vectors

◆ isAResonance()

bool isAResonance ( Belle2::Particle * particle)
staticprotectedinherited

controls if a particle is treated as a resonance(lifetime=0) or a particle that has a finite lifetime.

A finite life time means it will register a geo constraint for this particle

Definition at line 132 of file ParticleBase.cc.

133 {
134 bool rc = false ;
135 const int pdgcode = std::abs(particle->getPDGCode());
136
137 if (pdgcode && !(particle->getMdstArrayIndex())) {
138 switch (pdgcode) {
139 case 22: //photon conversion
140 rc = false;
141 break ;
142
143 case -11: //bremsstrahlung
144 case 11:
145 rc = true ;
146 break ;
147 default: //everything with boosted flight length less than 1 micrometer
148 rc = (pdgcode && particle->getPDGLifetime() < 1e-14);
149 }
150 }
151 return rc ;
152 }

◆ locate()

const ParticleBase * locate ( Belle2::Particle * particle) const
inherited

get particle base from basf2 particle

Definition at line 209 of file ParticleBase.cc.

210 {
211 const ParticleBase* rc = (m_particle == particle) ? this : nullptr;
212 if (!rc) {
213 for (auto* daughter : m_daughters) {
214 rc = daughter->locate(particle);
215 if (rc) {break;}
216 }
217 }
218 return rc;
219 }

◆ momIndex()

virtual int momIndex ( ) const
inlineoverridevirtualinherited

get momentum index

Reimplemented from ParticleBase.

Definition at line 42 of file RecoParticle.h.

42{ return index(); }

◆ mother()

const ParticleBase * mother ( ) const
inlineinherited

getMother() / hasMother()

Definition at line 94 of file ParticleBase.h.

94{ return m_mother; };

◆ nFinalChargedCandidates()

virtual int nFinalChargedCandidates ( ) const
inlineoverridevirtual

number of final charged candidates

Reimplemented from ParticleBase.

Definition at line 52 of file RecoTrack.h.

52{ return 1 ; }

◆ parname()

std::string parname ( int index) const
overridevirtualinherited

name

Reimplemented from ParticleBase.

Definition at line 23 of file RecoParticle.cc.

24 {
25 return ParticleBase::parname(index + 4);
26 }

◆ particle()

Belle2::Particle * particle ( ) const
inlineinherited

get basf2 particle

Definition at line 88 of file ParticleBase.h.

88{ return m_particle ; }

◆ posIndex()

virtual int posIndex ( ) const
inlinevirtualinherited

get vertex index (in statevector!)

Reimplemented in Composite, InternalParticle, Origin, RecoResonance, and Resonance.

Definition at line 118 of file ParticleBase.h.

118{ return -1 ; }

◆ projectConstraint()

ErrCode projectConstraint ( Constraint::Type type,
const FitParams & fitparams,
Projection & p ) const
overridevirtualinherited

abstract abstract projection

Reimplemented from ParticleBase.

Definition at line 28 of file RecoParticle.cc.

29 {
30 ErrCode status ;
31 switch (type) {
32 case Constraint::track :
33 case Constraint::neutralHadron :
34 status |= projectRecoConstraint(fitparams, p);
35 break ;
36 default:
37 status |= ParticleBase::projectConstraint(type, fitparams, p);
38 }
39 return status;
40 }

◆ projectGeoConstraint()

ErrCode projectGeoConstraint ( const FitParams & fitparams,
Projection & p ) const
virtualinherited

project geometrical constraint

the direction of the momentum is very well known from the kinematic constraints that is why we do not extract the distance as a vector here

Definition at line 247 of file ParticleBase.cc.

248 {
249 assert(m_config);
250 // only allow 2d for head of tree particles that are beam constrained
251 const int dim = m_config->m_originDimension == 2 && std::abs(m_particle->getPDGCode()) == m_config->m_headOfTreePDG ? 2 : 3;
252 const int posindexmother = mother()->posIndex();
253 const int posindex = posIndex();
254 const int tauindex = tauIndex();
255 const int momindex = momIndex();
256
257 const double tau = fitparams.getStateVector()(tauindex);
258 Eigen::VectorXd x_vec = fitparams.getStateVector().segment(posindex, dim);
259 Eigen::VectorXd x_m = fitparams.getStateVector().segment(posindexmother, dim);
260 Eigen::VectorXd p_vec = fitparams.getStateVector().segment(momindex, dim);
261 const double mom = p_vec.norm();
262 const double mom3 = mom * mom * mom;
263
264 if (3 == dim) {
265 // we can already set these
266 //diagonal momentum
267 p.getH()(0, momindex) = tau * (p_vec(1) * p_vec(1) + p_vec(2) * p_vec(2)) / mom3 ;
268 p.getH()(1, momindex + 1) = tau * (p_vec(0) * p_vec(0) + p_vec(2) * p_vec(2)) / mom3 ;
269 p.getH()(2, momindex + 2) = tau * (p_vec(0) * p_vec(0) + p_vec(1) * p_vec(1)) / mom3 ;
270
271 //offdiagonal momentum
272 p.getH()(0, momindex + 1) = - tau * p_vec(0) * p_vec(1) / mom3 ;
273 p.getH()(0, momindex + 2) = - tau * p_vec(0) * p_vec(2) / mom3 ;
274
275 p.getH()(1, momindex + 0) = - tau * p_vec(1) * p_vec(0) / mom3 ;
276 p.getH()(1, momindex + 2) = - tau * p_vec(1) * p_vec(2) / mom3 ;
277
278 p.getH()(2, momindex + 0) = - tau * p_vec(2) * p_vec(0) / mom3 ;
279 p.getH()(2, momindex + 1) = - tau * p_vec(2) * p_vec(1) / mom3 ;
280
281 } else if (2 == dim) {
282
283 // NOTE THAT THESE ARE DIFFERENT IN 2d
284 p.getH()(0, momindex) = tau * (p_vec(1) * p_vec(1)) / mom3 ;
285 p.getH()(1, momindex + 1) = tau * (p_vec(0) * p_vec(0)) / mom3 ;
286
287 //offdiagonal momentum
288 p.getH()(0, momindex + 1) = - tau * p_vec(0) * p_vec(1) / mom3 ;
289 p.getH()(1, momindex + 0) = - tau * p_vec(1) * p_vec(0) / mom3 ;
290 } else {
291 B2FATAL("Dimension of Geometric constraint is not 2 or 3. This will crash many things. You should feel bad.");
292 }
293
294 // --- Fill residuals and remaining Jacobian entries ---
298 Eigen::VectorXd residuals = x_m + tau * p_vec / mom - x_vec;
299 p.getResiduals().segment(0, dim) = residuals;
300
301 for (int row = 0; row < dim; ++row) {
302 p.getH()(row, posindexmother + row) = 1;
303 p.getH()(row, posindex + row) = -1;
304 p.getH()(row, tauindex) = p_vec(row) / mom;
305 }
306
307 return ErrCode(ErrCode::Status::success);
308 }

◆ projectMassConstraint()

ErrCode projectMassConstraint ( const FitParams & fitparams,
Projection & p ) const
virtualinherited

project mass constraint abstract

Definition at line 407 of file ParticleBase.cc.

409 {
410 assert(m_config);
411 if (m_config->m_massConstraintType == 0) {
412 return projectMassConstraintParticle(fitparams, p);
413 } else {
414 return projectMassConstraintDaughters(fitparams, p);
415 }
416 }

◆ projectMassConstraintDaughters()

ErrCode projectMassConstraintDaughters ( const FitParams & fitparams,
Projection & p ) const
virtualinherited

project mass constraint using the parameters of the daughters

be aware that the signs here are important E-|p|-m extracts a negative mass and messes with the momentum !

Definition at line 310 of file ParticleBase.cc.

312 {
313 double mass = 0;
314 if (particle()->hasExtraInfo("treeFitterMassConstraintValue")) {
315 mass = particle()->getExtraInfo("treeFitterMassConstraintValue");
316 } else mass = particle()->getPDGMass();
317 const double mass2 = mass * mass;
318 double px = 0;
319 double py = 0;
320 double pz = 0;
321 double E = 0;
322
323 // the parameters of the daughters must be used otherwise the mass constraint does not have an effect on the extracted daughter momenta
324 for (const auto* daughter : m_daughters) {
325 const int momindex = daughter->momIndex();
326 // in most cases the daughters will be final states so we cache the value to use it in the energy column
327 const double px_daughter = fitparams.getStateVector()(momindex);
328 const double py_daughter = fitparams.getStateVector()(momindex + 1);
329 const double pz_daughter = fitparams.getStateVector()(momindex + 2);
330
331 px += px_daughter;
332 py += py_daughter;
333 pz += pz_daughter;
334 if (daughter->hasEnergy()) {
335 E += fitparams.getStateVector()(momindex + 3);
336 } else {
337 // final states dont have an energy index
338 double m = 0;
339 if (daughter->particle()->hasExtraInfo("treeFitterMassConstraintValue")) {
340 m = daughter->particle()->getExtraInfo("treeFitterMassConstraintValue");
341 } else m = daughter->particle()->getPDGMass();
342 E += std::sqrt(m * m + px_daughter * px_daughter + py_daughter * py_daughter + pz_daughter * pz_daughter);
343 }
344 }
345
349 p.getResiduals()(0) = mass2 - E * E + px * px + py * py + pz * pz;
350
351 for (const auto* daughter : m_daughters) {
352 //dr/dx = d/dx m2-{E1+E2+...}^2+{p1+p2+...}^2 = 2*x (x= E or p)
353 const int momindex = daughter->momIndex();
354 Eigen::Vector3d mom(px, py, pz);
355 p.getH().block<1, 3>(0, momindex) = 2.0 * mom;
356
357 if (daughter->hasEnergy()) {
358 p.getH()(0, momindex + 3) = -2.0 * E;
359 } else {
360 Eigen::Vector3d mom_daughter = fitparams.getStateVector().segment<3>(momindex);
361 double m = 0;
362 if (daughter->particle()->hasExtraInfo("treeFitterMassConstraintValue")) {
363 m = daughter->particle()->getExtraInfo("treeFitterMassConstraintValue");
364 } else m = daughter->particle()->getPDGMass();
365
366 const double E_daughter = std::sqrt(m * m + mom_daughter.squaredNorm());
367 const double E_by_E_daughter = E / E_daughter;
368 p.getH().block<1, 3>(0, momindex) -= 2.0 * E_by_E_daughter * mom_daughter;
369 }
370
371 }
372 return ErrCode(ErrCode::Status::success);
373 }
R E
internal precision of FFTW codelets

◆ projectMassConstraintParticle()

ErrCode projectMassConstraintParticle ( const FitParams & fitparams,
Projection & p ) const
virtualinherited

project mass constraint using the particles parameters

be aware that the signs here are important E-|p|-m extracts a negative mass and messes with the momentum !

Definition at line 375 of file ParticleBase.cc.

377 {
378 double mass = 0;
379 if (particle()->hasExtraInfo("treeFitterMassConstraintValue")) mass = particle()->getExtraInfo("treeFitterMassConstraintValue");
380 else mass = particle()->getPDGMass();
381 const double mass2 = mass * mass;
382 const int momindex = momIndex();
383 Eigen::Vector4d state = fitparams.getStateVector().segment<4>(momindex);
384 Eigen::Vector3d p3 = state.head<3>();
385 const double px = state(0);
386 const double py = state(1);
387 const double pz = state(2);
388 const double E = state(3);
389
393 p.getResiduals()(0) = mass2 - E * E + px * px + py * py + pz * pz;
394
395 p.getH().block<1, 3>(0, momindex) = 2.0 * p3;
396 p.getH()(0, momindex + 3) = -2.0 * E;
397
398 // TODO 0 in most cases -> needs special treatment if width=0 to not crash chi2 calculation
399 // const double width = TDatabasePDG::Instance()->GetParticle(particle()->getPDGCode())->Width();
400 // transport measurement uncertainty into residual system
401 // f' = sigma_x^2 * (df/dx)^2
402 // p.getV()(0) = width * width * 4 * mass2;
403
404 return ErrCode(ErrCode::Status::success);
405 }

◆ projectRecoConstraint()

ErrCode projectRecoConstraint ( const FitParams & fitparams,
Projection & p ) const
overridevirtual

project the constraint (calculate residuals)

Implements RecoParticle.

Definition at line 98 of file RecoTrack.cc.

99 {
100 ErrCode status;
101 const int posindexmother = mother()->posIndex();
102 const int momindex = momIndex();
103
104 Eigen::Matrix<double, 1, 6> positionAndMom = Eigen::Matrix<double, 1, 6>::Zero(1, 6);
105 positionAndMom.segment(0, 3) = fitparams.getStateVector().segment(posindexmother, 3);
106 positionAndMom.segment(3, 3) = fitparams.getStateVector().segment(momindex, 3);
107 Eigen::Matrix<double, 5, 6> jacobian = Eigen::Matrix<double, 5, 6>::Zero(5, 6);
108
109 const Belle2::Helix helix = Belle2::Helix(ROOT::Math::XYZVector(positionAndMom(0), positionAndMom(1), positionAndMom(2)),
110 ROOT::Math::XYZVector(positionAndMom(3), positionAndMom(4), positionAndMom(5)),
111 charge(), m_bfield);
112
113 HelixUtils::getJacobianToCartesianFrameworkHelix(jacobian,
114 positionAndMom(0),
115 positionAndMom(1),
116 positionAndMom(2),
117 positionAndMom(3),
118 positionAndMom(4),
119 positionAndMom(5),
120 m_bfield,
121 charge()
122 );
123 if (!m_cached) {
124 auto* nonconst = const_cast<RecoTrack*>(this);
125 if (m_flt == 0) { nonconst->updFltToMother(fitparams); }
126 nonconst->updateParams(m_flt);
127 }
128
129 Eigen::Matrix<double, 1, 5> helixpars(5);
130 helixpars(0) = helix.getD0();
131 helixpars(1) = helix.getPhi0();
132 helixpars(2) = helix.getOmega();
133 helixpars(3) = helix.getZ0();
134 helixpars(4) = helix.getTanLambda();
135
136 p.getResiduals().segment(0, 5) = m_params - helixpars;
137
138 //account for periodic boundary in phi residual
139 double phiResidual = p.getResiduals().segment(0, 5)(1);
140 phiResidual = std::fmod(phiResidual + pi, twoPi);
141 if (phiResidual < 0) phiResidual += twoPi;
142 phiResidual -= pi;
143 p.getResiduals().segment(0, 5)(1) = phiResidual;
144
145 p.getV().triangularView<Eigen::Lower>() = m_covariance.triangularView<Eigen::Lower>();
146
147 p.getH().block<5, 3>(0, posindexmother) = -1.0 * jacobian.block<5, 3>(0, 0);
148 p.getH().block<5, 3>(0, momindex) = -1.0 * jacobian.block<5, 3>(0, 3);
149
150 return status;
151 }
double charge(int pdgCode)
Returns electric charge of a particle with given pdg code.
Definition EvtPDLUtil.cc:44

◆ removeDaughter()

void removeDaughter ( const ParticleBase * pb)
virtualinherited

remove daughter

Definition at line 71 of file ParticleBase.cc.

72 {
73 auto iter = std::find(m_daughters.begin(), m_daughters.end(), pb);
74 if (iter != m_daughters.end()) {
75 delete *iter;
76 m_daughters.erase(iter);
77 } else {
78 B2ERROR("Cannot remove particle, because not found ...");
79 }
80 }

◆ retrieveIndexMap()

void retrieveIndexMap ( indexmap & anindexmap) const
virtualinherited

get index map

Definition at line 221 of file ParticleBase.cc.

222 {
223 map.push_back(std::pair<const ParticleBase*, int>(this, index()));
224 for (auto* daughter : m_daughters) {
225 daughter->retrieveIndexMap(map);
226 }
227 }

◆ setFlightLength()

void setFlightLength ( double flt)
inline

setter for the flight length

Definition at line 64 of file RecoTrack.h.

64{ m_flt = flt ; }

◆ setIndex()

void setIndex ( int i)
inlineprotectedinherited

set Index (in statevector)

Definition at line 185 of file ParticleBase.h.

185{ m_index = i ; }

◆ setMother()

void setMother ( const ParticleBase * m)
inlineinherited

set mother

Definition at line 160 of file ParticleBase.h.

160{ m_mother = m ; }

◆ tauIndex()

virtual int tauIndex ( ) const
inlinevirtualinherited

get tau index

Reimplemented in Composite, InternalParticle, Origin, RecoResonance, and Resonance.

Definition at line 121 of file ParticleBase.h.

121{ return -1 ; }

◆ type()

virtual int type ( ) const
inlineoverridevirtual

type of the constraint

Implements ParticleBase.

Definition at line 49 of file RecoTrack.h.

49{ return kRecoTrack ; }

◆ updateIndex()

void updateIndex ( int & offset)
virtualinherited

this sets the index for momentum, position, etc.

in the statevector

Definition at line 82 of file ParticleBase.cc.

83 {
84 for (auto* daughter : m_daughters) {
85 daughter->updateIndex(offset);
86 }
87 m_index = offset;
88 offset += dim();
89 }

◆ updateParams()

ErrCode updateParams ( double flt)

updated the cached parameters

Definition at line 79 of file RecoTrack.cc.

80 {
81 m_flt = flt;
82 std::vector<float> trackParameter = m_trackfit->getTau();
83 for (unsigned int i = 0; i < trackParameter.size(); ++i) {
84 m_params(i) = trackParameter[i];
85 if (i == 2) m_params(2) /= m_momentumScalingFactor; // index 2 is the curvature of the track
86 }
87 TMatrixDSym cov = m_trackfit->getCovariance5();
88 for (int row = 0; row < 5; ++row) {
89 for (int col = 0; col <= row; ++col) {
90 m_covariance(row, col) = cov[row][col];
91 }
92 }
93
94 m_cached = true;
95 return ErrCode(ErrCode::Status::success);
96 }

◆ updFltToMother()

ErrCode updFltToMother ( const FitParams & fitparams)

update flight length to mother

Definition at line 70 of file RecoTrack.cc.

71 {
72 const int posindexmother = mother()->posIndex();
73 m_flt = m_trackfit->getHelix(m_momentumScalingFactor).getArcLength2DAtXY(
74 fitparams.getStateVector()(posindexmother),
75 fitparams.getStateVector()(posindexmother + 1));
76 return ErrCode(ErrCode::Status::success);
77 } ;

Member Data Documentation

◆ m_bfield

double m_bfield
private

B field along z.

Definition at line 69 of file RecoTrack.h.

◆ m_cached

bool m_cached
private

flag to mark the particle as initialised

Definition at line 75 of file RecoTrack.h.

◆ m_config

const ConstraintConfiguration* m_config
protectedinherited

has all the constraint config

Definition at line 200 of file ParticleBase.h.

◆ m_covariance

Eigen::Matrix<double, 5, 5> m_covariance
private

only lower triangle filled!

Definition at line 84 of file RecoTrack.h.

◆ m_daughters

std::vector<ParticleBase*> m_daughters
protectedinherited

daughter container

Definition at line 194 of file ParticleBase.h.

◆ m_flt

double m_flt
private

helix arc length at vertex

Definition at line 78 of file RecoTrack.h.

◆ m_index

int m_index
privateinherited

index

Definition at line 204 of file ParticleBase.h.

◆ m_isStronglyDecayingResonance

bool m_isStronglyDecayingResonance
protectedinherited

decay length less than 1 micron

Definition at line 197 of file ParticleBase.h.

◆ m_momentumScalingFactor

const float m_momentumScalingFactor
private

scale the momenta by this correction factor

Definition at line 87 of file RecoTrack.h.

◆ m_mother

const ParticleBase* m_mother
protectedinherited

motherparticle

Definition at line 191 of file ParticleBase.h.

◆ m_name

std::string m_name
privateinherited

name

Definition at line 207 of file ParticleBase.h.

◆ m_params

Eigen::Matrix<double, 1, 5> m_params
private

column vector to store the measurement

Definition at line 81 of file RecoTrack.h.

◆ m_particle

Belle2::Particle* m_particle
protectedinherited

pointer to framework type

Definition at line 188 of file ParticleBase.h.

◆ m_trackfit

const Belle2::TrackFitResult* m_trackfit
private

trackfit result from reconstruction

Definition at line 72 of file RecoTrack.h.


The documentation for this class was generated from the following files: