Belle II Software development
ParticleVertexFitterModule.cc
1/**************************************************************************
2 * basf2 (Belle II Analysis Software Framework) *
3 * Author: The Belle II Collaboration *
4 * *
5 * See git log for contributors and copyright holders. *
6 * This file is licensed under LGPL-3.0, see LICENSE.md. *
7 **************************************************************************/
8
9// Own header.
10#include <analysis/modules/ParticleVertexFitter/ParticleVertexFitterModule.h>
11
12// framework aux
13#include <framework/gearbox/Unit.h>
14#include <framework/gearbox/Const.h>
15#include <framework/logging/Logger.h>
16#include <framework/particledb/EvtGenDatabasePDG.h>
17
18// dataobjects
19#include <analysis/dataobjects/Particle.h>
20#include <analysis/dataobjects/Btube.h>
21#include <mdst/dataobjects/V0.h>
22
23// utilities
24#include <analysis/utility/CLHEPToROOT.h>
25#include <analysis/utility/PCmsLabTransform.h>
26#include <analysis/utility/ParticleCopy.h>
27#include <analysis/utility/ROOTToCLHEP.h>
28
29// Magnetic field
30#include <framework/geometry/BFieldManager.h>
31
32// KFit
33#include <analysis/VertexFitting/KFit/KFitConst.h>
34#include <analysis/VertexFitting/KFit/MassFitKFit.h>
35#include <analysis/VertexFitting/KFit/FourCFitKFit.h>
36#include <analysis/VertexFitting/KFit/MassFourCFitKFit.h>
37#include <analysis/VertexFitting/KFit/MassPointingVertexFitKFit.h>
38#include <analysis/VertexFitting/KFit/MassVertexFitKFit.h>
39#include <analysis/VertexFitting/KFit/VertexFitKFit.h>
40#include <analysis/VertexFitting/KFit/MakeMotherKFit.h>
41#include <analysis/VertexFitting/KFit/RecoilMassKFit.h>
42
43// Rave
44#include <analysis/VertexFitting/RaveInterface/RaveSetup.h>
45#include <analysis/VertexFitting/RaveInterface/RaveVertexFitter.h>
46#include <analysis/VertexFitting/RaveInterface/RaveKinematicVertexFitter.h>
47
48#include <TVector.h>
49#include <TRotation.h>
50#include <TMath.h>
51
52using namespace std;
53using namespace Belle2;
54
55
56//-----------------------------------------------------------------
57// Register module
58//-----------------------------------------------------------------
59
60REG_MODULE(ParticleVertexFitter);
61
62//-----------------------------------------------------------------
63// Implementation
64//-----------------------------------------------------------------
65
67 m_Bfield(0)
68{
69 // set module description (e.g. insert text)
70 setDescription("Vertex fitter for modular analysis");
72
73 // Add parameters
74 addParam("listName", m_listName, "name of particle list", string(""));
75 addParam("confidenceLevel", m_confidenceLevel,
76 "Confidence level to accept the fit. Particle candidates with "
77 "p-value less than confidenceLevel are removed from the particle "
78 "list. If set to -1, all candidates are kept; if set to 0, "
79 "the candidates failing the fit are removed.",
80 0.001);
81 addParam("vertexFitter", m_vertexFitter, "KFit or Rave", string("KFit"));
82 addParam("fitType", m_fitType, "type of the kinematic fit (vertex, massvertex, mass)", string("vertex"));
83 addParam("withConstraint", m_withConstraint,
84 "additional constraint on vertex: ipprofile, iptube, mother, iptubecut, pointing, btube",
85 string(""));
86 addParam("decayString", m_decayString, "specifies which daughter particles are included in the kinematic fit", string(""));
87 addParam("updateDaughters", m_updateDaughters, "true: update the daughters after the vertex fit", false);
88 addParam("smearing", m_smearing, "smear IP tube width by given length", 0.002);
89 addParam("recoilMass", m_recoilMass, "recoil invariant mass (GeV)", 0.);
90 addParam("massConstraintList", m_massConstraintList,
91 "Type::[int]. List of daughter particles to mass constrain with int = pdg code. (only for MassFourCKFit)", {});
92}
93
95{
96 // Particle list with name m_listName has to exist
97 m_plist.isRequired(m_listName);
98
99 // magnetic field
100 m_Bfield = BFieldManager::getFieldInTesla(ROOT::Math::XYZVector(0, 0, 0)).Z();
101
102 // RAVE setup
103 if (m_vertexFitter == "Rave")
105
106 B2DEBUG(1, "ParticleVertexFitterModule : magnetic field = " << m_Bfield);
107
108
109 if (m_decayString != "")
111
112 B2INFO("ParticleVertexFitter: Performing " << m_fitType << " fit on " << m_listName << " using " << m_vertexFitter);
113 if (m_decayString != "")
114 B2INFO("ParticleVertexFitter: Using specified decay string: " << m_decayString);
115 if (m_withConstraint != "")
116 B2INFO("ParticleVertexFitter: Additional " << m_withConstraint << " will be applied");
117
118}
119
121{
122 // TODO: set magnetic field for each run
123 // m_Bfield = BFieldManager::getFieldInTesla(ROOT::Math::XYZVector(0, 0, 0)).Z();
124 // TODO: set IP spot size for each run
125}
126
128{
129 if (m_vertexFitter == "Rave")
131
132 m_BeamSpotCenter = m_beamSpotDB->getIPPosition();
133 m_beamSpotCov.ResizeTo(3, 3);
134 if (m_withConstraint == "ipprofile") m_beamSpotCov = m_beamSpotDB->getCovVertex();
135 if (m_withConstraint == "iptube") {
136 if (m_smearing > 0 && m_vertexFitter == "KFit") {
138 } else {
140 }
141 }
142 if (m_withConstraint == "iptubecut") { // for development purpose only
143 m_BeamSpotCenter = ROOT::Math::XYZVector(0.001, 0., .013);
145 }
146 if ((m_vertexFitter == "Rave") && (m_withConstraint == "ipprofile" || m_withConstraint == "iptube"
147 || m_withConstraint == "mother" || m_withConstraint == "iptubecut" || m_withConstraint == "btube"))
149
150 std::vector<unsigned int> toRemove;
151 unsigned int nParticles = m_plist->getListSize();
152
153 for (unsigned iPart = 0; iPart < nParticles; iPart++) {
154 Particle* particle = m_plist->getParticle(iPart);
155 m_hasCovMatrix = false;
156 if (m_updateDaughters == true) {
157 if (m_decayString.empty() || m_vertexFitter == "KFit")
159 else B2ERROR("Daughters update works only when all daughters are selected. Daughters will not be updated");
160 }
161
162 if (m_withConstraint == "mother") {
163 m_BeamSpotCenter = particle->getVertex();
164 m_beamSpotCov = particle->getVertexErrorMatrix();
165 }
166
167 TMatrixFSym mother_errMatrix(7);
168 mother_errMatrix = particle->getMomentumVertexErrorMatrix();
169 for (int k = 0; k < 7; k++) {
170 for (int j = 0; j < 7; j++) {
171 if (mother_errMatrix[k][j] > 0) {
172 m_hasCovMatrix = true;
173 }
174 }
175 }
176
177 bool hasTube = true;
178 if (m_withConstraint == "btube") {
179 Btube* Ver = particle->getRelatedTo<Btube>();
180 if (!Ver) {
181 hasTube = false;
182 toRemove.push_back(particle->getArrayIndex());
183 } else {
184 m_BeamSpotCenter.SetXYZ(Ver->getTubeCenter()(0, 0), Ver->getTubeCenter()(1, 0), Ver->getTubeCenter()(2, 0));
186 }
187 }
188 bool ok = false;
189 if (hasTube) {
190 ok = doVertexFit(particle);
191 }
192 if (!ok)
193 particle->setPValue(-1);
194 if (particle->getPValue() < m_confidenceLevel)
195 toRemove.push_back(particle->getArrayIndex());
196
197 }
198 m_plist->removeParticles(toRemove);
199
200 //free memory allocated by rave. initialize() would be enough, except that we must clean things up before program end...
201 if (m_vertexFitter == "Rave")
203}
204
206{
207 // steering starts here
208
209 if (m_Bfield == 0) {
210 B2FATAL("ParticleVertexFitter: No magnetic field");
211 }
212
213 if (m_withConstraint != "ipprofile" &&
214 m_withConstraint != "iptube" &&
215 m_withConstraint != "mother" &&
216 m_withConstraint != "iptubecut" &&
217 m_withConstraint != "pointing" &&
218 m_withConstraint != "btube" &&
219 m_withConstraint != "")
220 B2FATAL("ParticleVertexFitter: " << m_withConstraint << " ***invalid Constraint ");
221
222 bool ok = false;
223 // fits with KFit
224 if (m_vertexFitter == "KFit") {
225
226 if (m_decayString != "" and m_fitType != "vertex")
227 B2FATAL("ParticleVertexFitter: KFit does not support yet selection of daughters via decay string except for vertex fit!");
228
229 // vertex fit
230 if (m_fitType == "vertex") {
231 if (m_withConstraint == "ipprofile") {
232 ok = doKVertexFit(mother, true, false);
233 } else if (m_withConstraint == "iptube") {
234 ok = doKVertexFit(mother, false, true);
235 } else {
236 ok = doKVertexFit(mother, false, false);
237 }
238 }
239
240 // mass-constrained vertex fit
241 if (m_fitType == "massvertex") {
242 if (m_withConstraint == "ipprofile" || m_withConstraint == "iptube" || m_withConstraint == "iptubecut") {
243 B2FATAL("ParticleVertexFitter: Invalid options - mass-constrained fit using KFit does not work with iptube or ipprofile constraint.");
244 } else if (m_withConstraint == "pointing") {
245 ok = doKMassPointingVertexFit(mother);
246 } else {
247 ok = doKMassVertexFit(mother);
248 }
249 }
250
251 // mass fit
252 if (m_fitType == "mass") {
253 if (m_withConstraint == "ipprofile" || m_withConstraint == "iptube" || m_withConstraint == "iptubecut") {
254 B2FATAL("ParticleVertexFitter: Invalid options - mass fit using KFit does not work with iptube or ipprofile constraint.");
255 } else {
256 ok = doKMassFit(mother);
257 }
258 }
259
260 // four C fit
261 if (m_fitType == "fourC") {
262 if (m_withConstraint == "ipprofile" || m_withConstraint == "iptube" || m_withConstraint == "iptubecut") {
263 B2FATAL("ParticleVertexFitter: Invalid options - four C fit using KFit does not work with iptube or ipprofile constraint.");
264 } else {
265 ok = doKFourCFit(mother);
266 }
267 }
268
269 // four mass C fit
270 if (m_fitType == "massfourC") {
271 if (m_withConstraint == "ipprofile" || m_withConstraint == "iptube" || m_withConstraint == "iptubecut") {
272 B2FATAL("ParticleVertexFitter: Invalid options - four C fit using KFit does not work with iptube or ipprofile constraint.");
273 } else {
274 ok = doKMassFourCFit(mother);
275 }
276 }
277
278 // recoil mass C fit
279 if (m_fitType == "recoilmass") {
280 if (m_withConstraint == "ipprofile" || m_withConstraint == "iptube" || m_withConstraint == "iptubecut") {
281 B2FATAL("ParticleVertexFitter: Invalid options - recoil mass fit using KFit does not work with iptube or ipprofile constraint.");
282 } else {
283 ok = doKRecoilMassFit(mother);
284 }
285 }
286
287 // invalid KFit fit type
288 if (m_fitType != "vertex"
289 && m_fitType != "massvertex"
290 && m_fitType != "mass"
291 && m_fitType != "fourC"
292 && m_fitType != "massfourC"
293 && m_fitType != "recoilmass")
294 B2FATAL("ParticleVertexFitter: " << m_fitType << " ***invalid fit type for the vertex fitter ");
295 }
296
297 // fits using Rave
298 if (m_vertexFitter == "Rave") {
299 try {
300 ok = doRaveFit(mother);
301 } catch (const rave::CheckedFloatException&) {
302 B2ERROR("Invalid inputs (nan/inf)?");
303 ok = false;
304 }
305 }
306
307 // invalid fitter
308 if (m_vertexFitter != "KFit" && m_vertexFitter != "Rave")
309 B2FATAL("ParticleVertexFitter: " << m_vertexFitter << " ***invalid vertex fitter ");
310
311 if (!ok) return false;
312
313 // steering ends here
314
315 //if (mother->getPValue() < m_confidenceLevel) return false;
316 return true;
317
318}
319
320bool ParticleVertexFitterModule::fillFitParticles(const Particle* mother, std::vector<const Particle*>& fitChildren,
321 std::vector<const Particle*>& twoPhotonChildren)
322{
323 if (m_decayString.empty()) {
324 // if decayString is empty, just use all primary daughters
325 for (unsigned ichild = 0; ichild < mother->getNDaughters(); ichild++) {
326 const Particle* child = mother->getDaughter(ichild);
327 // This if allows to skip the daughters, which cannot be used in the fits, particularly K_L0 from KLM.
328 // Useful for fully-inclusive particles.
329 if (mother->getProperty() == Particle::PropertyFlags::c_IsUnspecified and child->getPValue() < 0) {
330 continue;
331 }
332 fitChildren.push_back(child);
333 }
334 } else {
335 fitChildren = m_decaydescriptor.getSelectionParticles(mother);
336 }
337
338 auto itr = fitChildren.begin();
339 while (itr != fitChildren.end()) {
340 const Particle* child = *itr;
341
342 if (child->getPValue() < 0) {
343 B2WARNING("Daughter with PDG code " << child->getPDGCode() << " does not have a valid error matrix.");
344 return false; // error matrix not valid
345 }
346 bool isTwoPhotonParticle = false;
347 if (m_hasCovMatrix == false) {
348 if (child->getPDGCode() == Const::pi0.getPDGCode() or child->getPDGCode() == 221) { // pi0 or eta
349 if (child->getNDaughters() == 2) {
350 if (child->getDaughter(0)->getPDGCode() == Const::photon.getPDGCode()
351 && child->getDaughter(1)->getPDGCode() == Const::photon.getPDGCode()) {
352 isTwoPhotonParticle = true;
353 }
354 }
355 }
356 }
357 if (isTwoPhotonParticle) {
358 // move children from fitChildren to twoPhotonChildren
359 twoPhotonChildren.push_back(child);
360 itr = fitChildren.erase(itr);
361 } else {
362 itr++;
363 }
364 }
365
366 return true;
367}
368
369bool ParticleVertexFitterModule::fillNotFitParticles(const Particle* mother, std::vector<const Particle*>& notFitChildren,
370 const std::vector<const Particle*>& fitChildren)
371{
372 if (fitChildren.empty())
373 B2WARNING("[ParticleVertexFitterModule::fillNotFitParticles] fitChildren is empty! Please call fillFitParticles firstly");
374 if (!notFitChildren.empty())
375 B2WARNING("[ParticleVertexFitterModule::fillNotFitParticles] notFitChildren is NOT empty!"
376 << " The function should be called only once");
377
378 if (m_decayString.empty())
379 // if decayString is empty, just use all primary daughters
380 return true;
381
382 std::function<bool(const Particle*)> funcCheckInFit =
383 [&funcCheckInFit, &notFitChildren, fitChildren](const Particle * part) {
384
385 // check if the given particle in fitChildren
386 // if it is included, return true
387 if (std::find(fitChildren.begin(), fitChildren.end(), part) != fitChildren.end())
388 return true;
389
390 // if not, firstly check if particle has children
391 if (part->getNDaughters() == 0)
392 // if it has no children (=final-state-particle), return false
393 return false;
394
395 // here, the given particle is not in fitChildren and has children
396 bool isAnyChildrenInFit = false;
397 vector<const Particle*> notFitChildren_tmp;
398 for (unsigned ichild = 0; ichild < part->getNDaughters(); ichild++) {
399 // call funcCheckInFit recursively for all children
400 const Particle* child = part->getDaughter(ichild);
401 bool isChildrenInFit = funcCheckInFit(child);
402 isAnyChildrenInFit = isChildrenInFit or isAnyChildrenInFit;
403
404 // if the child is not in fitChildren, fill the child in a temporary vector
405 if (!isChildrenInFit)
406 notFitChildren_tmp.push_back(child);
407 }
408
409 // if there are a sister in fitChildren, the children in the temporary vector will be filled in notFitChildren
410 if (isAnyChildrenInFit)
411 notFitChildren.insert(notFitChildren.end(), notFitChildren_tmp.begin(), notFitChildren_tmp.end());
412
413 // if no children in fitChildren, the given particle should be filled instead of all children.
414
415 return isAnyChildrenInFit;
416 };
417
418
419 // call funcCheckInFit for all primary children
420 for (unsigned ichild = 0; ichild < mother->getNDaughters(); ichild++) {
421 const Particle* child = mother->getDaughter(ichild);
422 bool isGivenParticleOrAnyChildrenInFit = funcCheckInFit(child);
423 if (!isGivenParticleOrAnyChildrenInFit)
424 notFitChildren.push_back(child);
425 }
426
427 return true;
428}
429
431 const analysis::VertexFitKFit& kv)
432{
433 // TODO: something like setGammaError is necessary
434 // this is just workaround for the moment
435
436 const Particle* g1Orig = preFit->getDaughter(0);
437 const Particle* g2Orig = preFit->getDaughter(1);
438 Particle g1Temp(g1Orig->get4Vector(), 22);
439 Particle g2Temp(g2Orig->get4Vector(), 22);
440
441 TMatrixFSym g1ErrMatrix = g1Orig->getMomentumVertexErrorMatrix();
442 TMatrixFSym g2ErrMatrix = g2Orig->getMomentumVertexErrorMatrix();
443
444 ROOT::Math::XYZVector pos(kv.getVertex().x(), kv.getVertex().y(), kv.getVertex().z());
445 CLHEP::HepSymMatrix posErrorMatrix = kv.getVertexError();
446
447 TMatrixFSym errMatrix(3);
448 for (int i = 0; i < 3; i++)
449 for (int j = 0; j < 3; j++)
450 errMatrix(i, j) = posErrorMatrix[i][j];
451
452 g1ErrMatrix.SetSub(4, errMatrix);
453 g2ErrMatrix.SetSub(4, errMatrix);
454
455 g1Temp.updateMomentum(g1Orig->get4Vector(), pos, g1ErrMatrix, 1.0);
456 g2Temp.updateMomentum(g2Orig->get4Vector(), pos, g2ErrMatrix, 1.0);
457
458 // perform the mass fit for the two-photon particle
461
462 km.addParticle(&g1Temp);
463 km.addParticle(&g2Temp);
464
465 km.setVertex(kv.getVertex());
467 km.setInvariantMass(preFit->getPDGMass());
468
469 int err = km.doFit();
470 if (err != 0) {
471 return false;
472 }
473
474 // The update of the daughters is disabled for this mass fit.
475 bool updateDaughters = m_updateDaughters;
476 m_updateDaughters = false;
477 bool ok = makeKMassMother(km, postFit);
478 m_updateDaughters = updateDaughters;
479
480 return ok;
481}
482
483bool ParticleVertexFitterModule::doKVertexFit(Particle* mother, bool ipProfileConstraint, bool ipTubeConstraint)
484{
485 if ((mother->getNDaughters() < 2 && !ipTubeConstraint) || mother->getNDaughters() < 1) return false;
486
487 std::vector<const Particle*> fitChildren;
488 std::vector<const Particle*> twoPhotonChildren;
489 bool validChildren = fillFitParticles(mother, fitChildren, twoPhotonChildren);
490
491 if (!validChildren)
492 return false;
493
494 std::vector<const Particle*> notFitChildren;
495 fillNotFitParticles(mother, notFitChildren, fitChildren);
496
497
498 if (twoPhotonChildren.size() > 1) {
499 B2FATAL("[ParticleVertexFitterModule::doKVertexFit] Vertex fit using KFit does not support fit with multiple particles decaying to two photons like pi0 (yet).");
500 }
501
502 if ((fitChildren.size() < 2 && !ipTubeConstraint) || fitChildren.size() < 1) {
503 B2WARNING("[ParticleVertexFitterModule::doKVertexFit] Number of particles with valid error matrix entering the vertex fit using KFit is too low.");
504 return false;
505 }
506
507 // Initialise the Fitter
510
511 if (mother->getV0()) {
513 }
514
515 for (auto& child : fitChildren)
516 kv.addParticle(child);
517
518 if (ipProfileConstraint)
520
521 if (ipTubeConstraint)
522 addIPTubeToKFit(kv);
523
524 // Perform vertex fit using only the particles with valid error matrices
525 int err = kv.doFit();
526 if (err != 0)
527 return false;
528
529 double chi2_track = getChi2TracksLBoost(kv);
530 unsigned track_count = kv.getTrackCount();
531 mother->writeExtraInfo("chiSquared_trackL", chi2_track);
532 mother->writeExtraInfo("kFit_nTracks", track_count);
533
534 bool ok = false;
535 if (twoPhotonChildren.size() == 0)
536 // in the case daughters do not include pi0 - this is it (fit done)
537 ok = makeKVertexMother(kv, mother);
538 else if (twoPhotonChildren.size() == 1) {
539 // there is a daughter reconstructed from two photons so without position information
540 // 1. determine vertex based on all other valid daughters
541 // 2. set position and error matrix of two-photon daughter to previously determined vertex
542 // 3. redo the fit using all particles (including two-photon particle this time)
543
544 const Particle* twoPhotonDaughter = twoPhotonChildren[0];
545 Particle fixedTwoPhotonDaughter(twoPhotonDaughter->get4Vector(), twoPhotonDaughter->getPDGCode());
546 ok = redoTwoPhotonDaughterMassFit(&fixedTwoPhotonDaughter, twoPhotonDaughter, kv);
547 if (!ok)
548 return false;
549
550 // finally perform the fit using all daughter particles
553
554 for (auto& child : fitChildren)
555 kv2.addParticle(child);
556
557 kv2.addParticle(&fixedTwoPhotonDaughter);
558
559 if (ipProfileConstraint)
561
562 err = kv2.doFit();
563
564 if (err != 0)
565 return false;
566
567 ok = makeKVertexMother(kv2, mother);
568 }
569
570 // update 4-vector using not-fit-particles
571 ROOT::Math::PxPyPzEVector total4Vector(mother->get4Vector());
572 for (auto& child : notFitChildren)
573 total4Vector += child->get4Vector();
574 mother->set4Vector(total4Vector);
575
576 return ok;
577}
578
580{
581 if (mother->getNDaughters() < 2) return false;
582
583 std::vector<const Particle*> fitChildren;
584 std::vector<const Particle*> twoPhotonChildren;
585 bool validChildren = fillFitParticles(mother, fitChildren, twoPhotonChildren);
586
587 if (!validChildren)
588 return false;
589
590 if (twoPhotonChildren.size() > 1) {
591 B2FATAL("[ParticleVertexFitterModule::doKVertexFit] MassVertex fit using KFit does not support fit with multiple particles decaying to two photons like pi0 (yet).");
592 }
593
594 if (fitChildren.size() < 2) {
595 B2WARNING("[ParticleVertexFitterModule::doKVertexFit] Number of particles with valid error matrix entering the vertex fit using KFit is less than 2.");
596 return false;
597 }
598
599 bool ok = false;
600 if (twoPhotonChildren.size() == 0) {
601 // Initialise the Fitter
604
605 if (mother->getV0()) {
606 HepPoint3D V0vertex_heppoint(mother->getV0()->getFittedVertexX(),
607 mother->getV0()->getFittedVertexY(),
608 mother->getV0()->getFittedVertexZ());
609 kmv.setInitialVertex(V0vertex_heppoint);
610 }
611
612 for (auto child : fitChildren)
613 kmv.addParticle(child);
614
615 kmv.setInvariantMass(mother->getPDGMass());
616 int err = kmv.doFit();
617 if (err != 0)
618 return false;
619
620 // in the case daughters do not include particles with two photon daughters like pi0 - this is it (fit done)
621 ok = makeKMassVertexMother(kmv, mother);
622 } else if (twoPhotonChildren.size() == 1) {
623 // there is a daughter reconstructed from two photons so without position information
624 // 1. determine vertex based on all other valid daughters
625 // 2. set position and error matrix of two-photon daughter to previously determined vertex
626 // 3. redo the fit using all particles (including two-photon particle this time)
627
630
631 for (auto child : fitChildren)
632 kv.addParticle(child);
633
634 // Perform vertex fit using only the particles with valid error matrices
635 int err = kv.doFit();
636 if (err != 0)
637 return false;
638
639 const Particle* twoPhotonDaughter = twoPhotonChildren[0];
640 Particle fixedTwoPhotonDaughter(twoPhotonDaughter->get4Vector(), twoPhotonDaughter->getPDGCode());
641 ok = redoTwoPhotonDaughterMassFit(&fixedTwoPhotonDaughter, twoPhotonDaughter, kv);
642 if (!ok)
643 return false;
644
645 // finally perform the fit using all daughter particles
648
649 for (auto child : fitChildren)
650 kmv2.addParticle(child);
651 kmv2.addParticle(&fixedTwoPhotonDaughter);
652
653 kmv2.setInvariantMass(mother->getPDGMass());
654 err = kmv2.doFit();
655
656 if (err != 0)
657 return false;
658
659 ok = makeKMassVertexMother(kmv2, mother);
660 }
661
662 return ok;
663
664}
665
667{
668 if (!(mother->hasExtraInfo("prodVertX") && mother->hasExtraInfo("prodVertY") && mother->hasExtraInfo("prodVertZ"))) {
669 return false;
670 }
671
672 if (mother->getNDaughters() < 2) return false;
673
674 std::vector<const Particle*> fitChildren;
675 std::vector<const Particle*> twoPhotonChildren;
676 bool validChildren = fillFitParticles(mother, fitChildren, twoPhotonChildren);
677
678 if (!validChildren)
679 return false;
680
681 if (twoPhotonChildren.size() > 0) {
682 B2FATAL("[ParticleVertexFitterModule::doKMassPointingVertexFit] MassPointingVertex fit using KFit does not support fit with two-photon daughters (yet).");
683 }
684
685 if (fitChildren.size() < 2) {
686 B2WARNING("[ParticleVertexFitterModule::doKMassPointingVertexFit] Number of particles with valid error matrix entering the vertex fit using KFit is less than 2.");
687 return false;
688 }
689
690 bool ok = false;
691 // Initialise the Fitter
694
695 for (auto child : fitChildren)
696 kmpv.addParticle(child);
697
698 kmpv.setInvariantMass(mother->getPDGMass());
699 HepPoint3D productionVertex(mother->getExtraInfo("prodVertX"),
700 mother->getExtraInfo("prodVertY"),
701 mother->getExtraInfo("prodVertZ"));
702 kmpv.setProductionVertex(productionVertex);
703 int err = kmpv.doFit();
704 if (err != 0) return false;
705
706 ok = makeKMassPointingVertexMother(kmpv, mother);
707
708 return ok;
709}
710
712{
713 if (mother->getNDaughters() < 2) return false;
714
717
718 for (unsigned ichild = 0; ichild < mother->getNDaughters(); ichild++) {
719 const Particle* child = mother->getDaughter(ichild);
720
721 if (child->getPValue() < 0) return false; // error matrix not valid
722
723 km.addParticle(child);
724 }
725
726 // apply mass constraint
727 km.setInvariantMass(mother->getPDGMass());
728
729 int err = km.doFit();
730
731 if (err != 0) return false;
732
733 bool ok = makeKMassMother(km, mother);
734
735 return ok;
736}
737
739{
740 if (mother->getNDaughters() < 2) return false;
741
744
745 for (unsigned ichild = 0; ichild < mother->getNDaughters(); ichild++) {
746 const Particle* child = mother->getDaughter(ichild);
747
748 if (child->getNDaughters() > 0) {
749 bool err = addChildofParticletoKFit(kf, child);
750 if (!err) return false;
751 } else {
752 if (child->getPValue() < 0) return false; // error matrix not valid
753
754 kf.addParticle(child);
755 }
756 }
757
758 // apply four momentum constraint
760 kf.setFourMomentum(T.getBeamFourMomentum());
761
762 int err = kf.doFit();
763
764 if (err != 0) return false;
765
766 bool ok = makeKFourCMother(kf, mother);
767
768 return ok;
769}
770
772{
773 if (mother->getNDaughters() < 2) return false;
774
777
778 for (unsigned ichild = 0; ichild < mother->getNDaughters(); ichild++) {
779 const Particle* child = mother->getDaughter(ichild);
780
781 if (child->getNDaughters() > 0) {
782 bool massconstraint = std::find(m_massConstraintList.begin(), m_massConstraintList.end(),
783 std::abs(child->getPDGCode())) != m_massConstraintList.end();
784 std::vector<unsigned> childId;
785 bool err = addChildofParticletoMassKFit(kf, child, childId);
786 if (massconstraint) kf.addMassConstraint(child->getPDGMass(), childId);
787 if (!err) return false;
788 } else {
789 if (child->getPValue() < 0) return false; // error matrix not valid
790 kf.addParticle(child);
791 }
792 }
793
794 // apply four momentum constraint
796 kf.setFourMomentum(T.getBeamFourMomentum());
797
798 int err = kf.doFit();
799
800 if (err != 0) return false;
801
802 bool ok = makeMassKFourCMother(kf, mother);
803
804 return ok;
805}
806
808{
811
812 for (unsigned ichild = 0; ichild < mother->getNDaughters(); ichild++) {
813 const Particle* child = mother->getDaughter(ichild);
814
815 if (child->getPValue() < 0) return false; // error matrix not valid
816
817 kf.addParticle(child);
818 }
819
820 // apply four momentum constraint
822 kf.setFourMomentum(T.getBeamFourMomentum());
823
824 // apply recoil mass constraint
826
827 int err = kf.doFit();
828
829 if (err != 0) return false;
830
831 bool ok = makeKRecoilMassMother(kf, mother);
832
833 return ok;
834}
835
837 Particle* mother)
838{
839 enum analysis::KFitError::ECode fitError;
840 fitError = kv.updateMother(mother);
841 if (fitError != analysis::KFitError::kNoError)
842 return false;
843 if (m_decayString.empty() && m_updateDaughters == true) {
844 // update daughter momenta as well
845 // the order of daughters in the *fitter is the same as in the mother Particle
846
847 std::vector<Particle*> daughters = mother->getDaughters();
848
849 unsigned track_count = kv.getTrackCount();
850 if (daughters.size() != track_count)
851 return false;
852
853 for (unsigned iChild = 0; iChild < track_count; iChild++) {
854 double a = -1 * Const::speedOfLight * 1e-4 * m_Bfield * daughters[iChild]->getCharge();
855 double dx = kv.getVertex().x() - kv.getTrackPosition(iChild).x();
856 double dy = kv.getVertex().y() - kv.getTrackPosition(iChild).y();
857
858 ROOT::Math::PxPyPzEVector i4Vector(kv.getTrackMomentum(iChild).x() - a * dy,
859 kv.getTrackMomentum(iChild).y() + a * dx,
860 kv.getTrackMomentum(iChild).z(),
861 kv.getTrackMomentum(iChild).t());
862 daughters[iChild]->set4VectorDividingByMomentumScaling(i4Vector);
863
864 daughters[iChild]->setVertex(
865 CLHEPToROOT::getXYZVector(kv.getTrackPosition(iChild)));
866 daughters[iChild]->setMomentumVertexErrorMatrix(
867 CLHEPToROOT::getTMatrixFSym(kv.getTrackError(iChild)));
868 }
869
870 } else if (m_updateDaughters == true) { // if decayString is not empty
871 // first, update only the fit children
872 std::vector<const Particle*> fitChildren = m_decaydescriptor.getSelectionParticles(mother);
873
874 unsigned track_count = kv.getTrackCount();
875 if (fitChildren.size() != track_count)
876 return false;
877
878 for (unsigned iChild = 0; iChild < track_count; iChild++) {
879 auto daughter = const_cast<Particle*>(fitChildren[iChild]);
880
881 double a = -1 * Const::speedOfLight * 1e-4 * m_Bfield * daughter->getCharge();
882 double dx = kv.getVertex().x() - kv.getTrackPosition(iChild).x();
883 double dy = kv.getVertex().y() - kv.getTrackPosition(iChild).y();
884
885 ROOT::Math::PxPyPzEVector i4Vector(kv.getTrackMomentum(iChild).x() - a * dy,
886 kv.getTrackMomentum(iChild).y() + a * dx,
887 kv.getTrackMomentum(iChild).z(),
888 kv.getTrackMomentum(iChild).t());
889 daughter->set4VectorDividingByMomentumScaling(i4Vector);
890
891 daughter->setVertex(CLHEPToROOT::getXYZVector(kv.getTrackPosition(iChild)));
892 daughter->setMomentumVertexErrorMatrix(CLHEPToROOT::getTMatrixFSym(kv.getTrackError(iChild)));
893 }
894
895 // then, update other particles that have a fit-child in decay
896 std::function<bool(Particle*)> funcUpdateMomentum =
897 [&funcUpdateMomentum, fitChildren](Particle * part) {
898
899 if (part->getNDaughters() == 0) {
900 // check if part is included in fitChildren
901 if (std::find(fitChildren.begin(), fitChildren.end(), part) != fitChildren.end())
902 return true;
903 else
904 return false;
905 }
906
907 bool includeFitChildren = false;
908
909 // Update daughters' momentum
910 for (auto daughter : part->getDaughters())
911 includeFitChildren = funcUpdateMomentum(daughter) || includeFitChildren;
912
913 if (includeFitChildren) {
914 // Using updated daughters, update part's momentum
915 ROOT::Math::PxPyPzEVector sum4Vector;
916 for (auto daughter : part->getDaughters())
917 sum4Vector += daughter->get4Vector();
918
919 part->set4VectorDividingByMomentumScaling(sum4Vector);
920 }
921
922 return includeFitChildren;
923 };
924
925 // Update all daughters
926 for (auto daughter : mother->getDaughters())
927 funcUpdateMomentum(daughter);
928
929 }
930
931
932 return true;
933}
934
936 Particle* mother)
937{
938 enum analysis::KFitError::ECode fitError;
939 fitError = kmv.updateMother(mother);
940 if (fitError != analysis::KFitError::kNoError)
941 return false;
942 if (m_decayString.empty() && m_updateDaughters == true) {
943 // update daughter momenta as well
944 // the order of daughters in the *fitter is the same as in the mother Particle
945
946 std::vector<Particle*> daughters = mother->getDaughters();
947
948 unsigned track_count = kmv.getTrackCount();
949 if (daughters.size() != track_count)
950 return false;
951
952 for (unsigned iChild = 0; iChild < track_count; iChild++) {
953 double a = -1 * Const::speedOfLight * 1e-4 * m_Bfield * daughters[iChild]->getCharge();
954 double dx = kmv.getVertex().x() - kmv.getTrackPosition(iChild).x();
955 double dy = kmv.getVertex().y() - kmv.getTrackPosition(iChild).y();
956
957 ROOT::Math::PxPyPzEVector i4Vector(kmv.getTrackMomentum(iChild).x() - a * dy,
958 kmv.getTrackMomentum(iChild).y() + a * dx,
959 kmv.getTrackMomentum(iChild).z(),
960 kmv.getTrackMomentum(iChild).t());
961 daughters[iChild]->set4VectorDividingByMomentumScaling(i4Vector);
962
963 daughters[iChild]->setVertex(
964 CLHEPToROOT::getXYZVector(kmv.getTrackPosition(iChild)));
965 daughters[iChild]->setMomentumVertexErrorMatrix(
966 CLHEPToROOT::getTMatrixFSym(kmv.getTrackError(iChild)));
967 }
968 }
969
970 return true;
971}
972
974 Particle* mother)
975{
976 enum analysis::KFitError::ECode fitError;
977 fitError = kmpv.updateMother(mother);
978 if (fitError != analysis::KFitError::kNoError) {
979 return false;
980 }
981
982 if (m_decayString.empty() && m_updateDaughters == true) {
983 // update daughter momenta as well
984 // the order of daughters in the *fitter is the same as in the mother Particle
985
986 std::vector<Particle*> daughters = mother->getDaughters();
987
988 unsigned track_count = kmpv.getTrackCount();
989 if (daughters.size() != track_count)
990 return false;
991
992 for (unsigned iChild = 0; iChild < track_count; iChild++) {
993 double a = -1 * Const::speedOfLight * 1e-4 * m_Bfield * daughters[iChild]->getCharge();
994 double dx = kmpv.getVertex().x() - kmpv.getTrackPosition(iChild).x();
995 double dy = kmpv.getVertex().y() - kmpv.getTrackPosition(iChild).y();
996
997 ROOT::Math::PxPyPzEVector i4Vector(kmpv.getTrackMomentum(iChild).x() - a * dy,
998 kmpv.getTrackMomentum(iChild).y() + a * dx,
999 kmpv.getTrackMomentum(iChild).z(),
1000 kmpv.getTrackMomentum(iChild).t());
1001 daughters[iChild]->set4VectorDividingByMomentumScaling(i4Vector);
1002
1003 daughters[iChild]->setVertex(
1004 CLHEPToROOT::getXYZVector(kmpv.getTrackPosition(iChild)));
1005 daughters[iChild]->setMomentumVertexErrorMatrix(
1006 CLHEPToROOT::getTMatrixFSym(kmpv.getTrackError(iChild)));
1007 }
1008 }
1009
1010 return true;
1011}
1012
1013
1015 Particle* mother)
1016{
1017 enum analysis::KFitError::ECode fitError;
1018 fitError = km.updateMother(mother);
1019 if (fitError != analysis::KFitError::kNoError)
1020 return false;
1021 if (m_decayString.empty() && m_updateDaughters == true) {
1022 // update daughter momenta as well
1023 // the order of daughters in the *fitter is the same as in the mother Particle
1024
1025 std::vector<Particle*> daughters = mother->getDaughters();
1026
1027 unsigned track_count = km.getTrackCount();
1028 if (daughters.size() != track_count)
1029 return false;
1030
1031 for (unsigned iChild = 0; iChild < track_count; iChild++) {
1032 double a = -1 * Const::speedOfLight * 1e-4 * m_Bfield * daughters[iChild]->getCharge();
1033 double dx = km.getVertex().x() - km.getTrackPosition(iChild).x();
1034 double dy = km.getVertex().y() - km.getTrackPosition(iChild).y();
1035
1036 ROOT::Math::PxPyPzEVector i4Vector(km.getTrackMomentum(iChild).x() - a * dy,
1037 km.getTrackMomentum(iChild).y() + a * dx,
1038 km.getTrackMomentum(iChild).z(),
1039 km.getTrackMomentum(iChild).t());
1040 daughters[iChild]->set4VectorDividingByMomentumScaling(i4Vector);
1041
1042 daughters[iChild]->setVertex(
1043 CLHEPToROOT::getXYZVector(km.getTrackPosition(iChild)));
1044 daughters[iChild]->setMomentumVertexErrorMatrix(
1045 CLHEPToROOT::getTMatrixFSym(km.getTrackError(iChild)));
1046 }
1047 }
1048
1049 return true;
1050}
1051
1052
1053
1055{
1056 enum analysis::KFitError::ECode fitError;
1057 fitError = kf.updateMother(mother);
1058 if (fitError != analysis::KFitError::kNoError)
1059 return false;
1060 mother->addExtraInfo("FourCFitProb", kf.getCHIsq());
1061 mother->addExtraInfo("FourCFitChi2", kf.getNDF());
1062 if (m_decayString.empty() && m_updateDaughters == true) {
1063 // update daughter momenta as well
1064 // the order of daughters in the *fitter is the same as in the mother Particle
1065
1066 std::vector<Particle*> daughters = mother->getDaughters();
1067
1068 const unsigned nd = daughters.size();
1069 unsigned l = 0;
1070 std::vector<std::vector<unsigned>> pars;
1071 std::vector<Particle*> allparticles;
1072 for (unsigned ichild = 0; ichild < nd; ichild++) {
1073 const Particle* daughter = mother->getDaughter(ichild);
1074 std::vector<unsigned> pard;
1075 if (daughter->getNDaughters() > 0) {
1076 updateMapOfTrackAndDaughter(l, pars, pard, allparticles, daughter);
1077 pars.push_back(pard);
1078 allparticles.push_back(daughters[ichild]);
1079 } else {
1080 pard.push_back(l);
1081 pars.push_back(pard);
1082 allparticles.push_back(daughters[ichild]);
1083 l++;
1084 }
1085 }
1086
1087 unsigned track_count = kf.getTrackCount();
1088 if (l != track_count)
1089 return false;
1090
1091 for (unsigned iDaug = 0; iDaug < allparticles.size(); iDaug++) {
1092 ROOT::Math::PxPyPzEVector childMoms;
1093 ROOT::Math::XYZVector childPoss;
1094 TMatrixFSym childErrMatrixs(7);
1095 for (unsigned int iChild : pars[iDaug]) {
1096 childMoms = childMoms +
1097 CLHEPToROOT::getLorentzVector(
1098 kf.getTrackMomentum(iChild));
1099 childPoss = childPoss +
1100 CLHEPToROOT::getXYZVector(
1101 kf.getTrackPosition(iChild));
1102 TMatrixFSym childErrMatrix =
1103 CLHEPToROOT::getTMatrixFSym(kf.getTrackError(iChild));
1104 childErrMatrixs = childErrMatrixs + childErrMatrix;
1105 }
1106 allparticles[iDaug]->set4Vector(childMoms);
1107 allparticles[iDaug]->setVertex(childPoss);
1108 allparticles[iDaug]->setMomentumVertexErrorMatrix(childErrMatrixs);
1109 }
1110 }
1111
1112 return true;
1113}
1114
1116{
1117 enum analysis::KFitError::ECode fitError;
1118 fitError = kf.updateMother(mother);
1119 if (fitError != analysis::KFitError::kNoError)
1120 return false;
1121 mother->addExtraInfo("MassFourCFitProb", TMath::Prob(kf.getCHIsq(), kf.getNDF()));
1122 mother->addExtraInfo("MassFourCFitChi2", kf.getCHIsq());
1123 mother->addExtraInfo("MassFourCFitNDF", kf.getNDF());
1124 if (m_decayString.empty() && m_updateDaughters == true) {
1125 // update daughter momenta as well
1126 // the order of daughters in the *fitter is the same as in the mother Particle
1127
1128 std::vector<Particle*> daughters = mother->getDaughters();
1129
1130 const unsigned nd = daughters.size();
1131 unsigned l = 0;
1132 std::vector<std::vector<unsigned>> pars;
1133 std::vector<Particle*> allparticles;
1134 for (unsigned ichild = 0; ichild < nd; ichild++) {
1135 const Particle* daughter = mother->getDaughter(ichild);
1136 std::vector<unsigned> pard;
1137 if (daughter->getNDaughters() > 0) {
1138 updateMapOfTrackAndDaughter(l, pars, pard, allparticles, daughter);
1139 pars.push_back(pard);
1140 allparticles.push_back(daughters[ichild]);
1141 } else {
1142 pard.push_back(l);
1143 pars.push_back(pard);
1144 allparticles.push_back(daughters[ichild]);
1145 l++;
1146 }
1147 }
1148
1149 unsigned track_count = kf.getTrackCount();
1150 if (l != track_count)
1151 return false;
1152
1153 for (unsigned iDaug = 0; iDaug < allparticles.size(); iDaug++) {
1154 ROOT::Math::PxPyPzEVector childMoms;
1155 ROOT::Math::XYZVector childPoss;
1156 TMatrixFSym childErrMatrixs(7);
1157 for (unsigned int iChild : pars[iDaug]) {
1158 childMoms = childMoms +
1159 CLHEPToROOT::getLorentzVector(
1160 kf.getTrackMomentum(iChild));
1161 childPoss = childPoss +
1162 CLHEPToROOT::getXYZVector(
1163 kf.getTrackPosition(iChild));
1164 TMatrixFSym childErrMatrix =
1165 CLHEPToROOT::getTMatrixFSym(kf.getTrackError(iChild));
1166 childErrMatrixs = childErrMatrixs + childErrMatrix;
1167 }
1168 allparticles[iDaug]->set4Vector(childMoms);
1169 allparticles[iDaug]->setVertex(childPoss);
1170 allparticles[iDaug]->setMomentumVertexErrorMatrix(childErrMatrixs);
1171 }
1172 }
1173
1174 return true;
1175}
1176
1178{
1179 enum analysis::KFitError::ECode fitError;
1180 fitError = kf.updateMother(mother);
1181 if (fitError != analysis::KFitError::kNoError)
1182 return false;
1183 mother->addExtraInfo("RecoilMassFitProb", TMath::Prob(kf.getCHIsq(), kf.getNDF()));
1184 mother->addExtraInfo("RecoilMassFitChi2", kf.getCHIsq());
1185 mother->addExtraInfo("RecoilMassFitNDF", kf.getNDF());
1186 if (m_decayString.empty() && m_updateDaughters == true) {
1187 // update daughter momenta as well
1188 // the order of daughters in the *fitter is the same as in the mother Particle
1189
1190 std::vector<Particle*> daughters = mother->getDaughters();
1191
1192 unsigned track_count = kf.getTrackCount();
1193 if (daughters.size() != track_count)
1194 return false;
1195
1196 for (unsigned iChild = 0; iChild < track_count; iChild++) {
1197 double a = -1 * Const::speedOfLight * 1e-4 * m_Bfield * daughters[iChild]->getCharge();
1198 double dx = kf.getVertex().x() - kf.getTrackPosition(iChild).x();
1199 double dy = kf.getVertex().y() - kf.getTrackPosition(iChild).y();
1200
1201 ROOT::Math::PxPyPzEVector i4Vector(kf.getTrackMomentum(iChild).x() - a * dy,
1202 kf.getTrackMomentum(iChild).y() + a * dx,
1203 kf.getTrackMomentum(iChild).z(),
1204 kf.getTrackMomentum(iChild).t());
1205 daughters[iChild]->set4VectorDividingByMomentumScaling(i4Vector);
1206
1207 daughters[iChild]->setVertex(
1208 CLHEPToROOT::getXYZVector(kf.getTrackPosition(iChild)));
1209 daughters[iChild]->setMomentumVertexErrorMatrix(
1210 CLHEPToROOT::getTMatrixFSym(kf.getTrackError(iChild)));
1211 }
1212 }
1213
1214 return true;
1215}
1216
1217
1218void ParticleVertexFitterModule::updateMapOfTrackAndDaughter(unsigned& l, std::vector<std::vector<unsigned>>& pars,
1219 std::vector<unsigned>& parm, std::vector<Particle*>& allparticles, const Particle* daughter)
1220{
1221 std::vector <Particle*> daughters = daughter->getDaughters();
1222 for (unsigned ichild = 0; ichild < daughter->getNDaughters(); ichild++) {
1223 const Particle* child = daughter->getDaughter(ichild);
1224 std::vector<unsigned> pard;
1225 if (child->getNDaughters() > 0) {
1226 updateMapOfTrackAndDaughter(l, pars, pard, allparticles, child);
1227 parm.insert(parm.end(), pard.begin(), pard.end());
1228 pars.push_back(pard);
1229 allparticles.push_back(daughters[ichild]);
1230 } else {
1231 pard.push_back(l);
1232 parm.push_back(l);
1233 pars.push_back(pard);
1234 allparticles.push_back(daughters[ichild]);
1235 l++;
1236 }
1237 }
1238}
1239
1240
1242{
1243 if ((m_decayString.empty() ||
1244 (m_withConstraint == "" && m_fitType != "mass")) && mother->getNDaughters() < 2) return false;
1246 if (m_withConstraint == "ipprofile" || m_withConstraint == "iptube" || m_withConstraint == "mother"
1247 || m_withConstraint == "iptubecut" || m_withConstraint == "btube")
1249
1251 if (m_fitType == "mass") rf.setVertFit(false);
1252
1253 if (m_decayString.empty()) {
1254 rf.addMother(mother);
1255 } else {
1256 std::vector<const Particle*> tracksVertex = m_decaydescriptor.getSelectionParticles(mother);
1257 std::vector<std::string> tracksName = m_decaydescriptor.getSelectionNames();
1258
1259 if (allSelectedDaughters(mother, tracksVertex)) {
1260 for (auto& itrack : tracksVertex) {
1261 if (itrack != mother) rf.addTrack(itrack);
1262 }
1263 rf.setMother(mother);
1264 } else {
1265
1267 bool mothSel = false;
1268 int nTrk = 0;
1269 for (unsigned itrack = 0; itrack < tracksVertex.size(); itrack++) {
1270 if (tracksVertex[itrack] != mother) {
1271 rsf.addTrack(tracksVertex[itrack]);
1272 B2DEBUG(1, "ParticleVertexFitterModule: Adding particle " << tracksName[itrack] << " to vertex fit ");
1273 nTrk++;
1274 }
1275 if (tracksVertex[itrack] == mother) mothSel = true;
1276 }
1277
1278
1279 // Fit one particle constrained to originate from the beam spot
1280 bool mothIPfit = false;
1281 if (tracksVertex.size() == 1 && mothSel == true && m_withConstraint != "" && nTrk == 0) {
1282 rsf.addTrack(tracksVertex[0]);
1283 if (tracksVertex[0] != mother)
1284 B2FATAL("ParticleVertexFitterModule: FATAL Error in IP constrained mother fit");
1285 nTrk++;
1286 mothIPfit = true;
1287 }
1288
1289
1290 ROOT::Math::XYZVector pos;
1291 TMatrixDSym RerrMatrix(3);
1292 int nvert = 0;
1293
1294 // one track fit is not kinematic
1295 if (nTrk == 1) {
1297 for (auto& itrack : tracksVertex) {
1298 rsg.addTrack(itrack);
1299 nvert = rsg.fit("kalman");
1300 if (nvert > 0) {
1301 pos = rsg.getPos(0);
1302 RerrMatrix = rsg.getCov(0);
1303 double prob = rsg.getPValue(0);
1304 ROOT::Math::PxPyPzEVector mom(mother->get4Vector());
1305 TMatrixDSym errMatrix(7);
1306 for (int i = 0; i < 7; i++) {
1307 for (int j = 0; j < 7; j++) {
1308 if (i > 3 && j > 3) {errMatrix[i][j] = RerrMatrix[i - 4][j - 4];}
1309 else {errMatrix[i][j] = 0;}
1310 }
1311 }
1312 if (mothIPfit) {
1313 mother->writeExtraInfo("prodVertX", pos.X());
1314 mother->writeExtraInfo("prodVertY", pos.Y());
1315 mother->writeExtraInfo("prodVertZ", pos.Z());
1316 mother->writeExtraInfo("prodVertSxx", RerrMatrix[0][0]);
1317 mother->writeExtraInfo("prodVertSxy", RerrMatrix[0][1]);
1318 mother->writeExtraInfo("prodVertSxz", RerrMatrix[0][2]);
1319 mother->writeExtraInfo("prodVertSyx", RerrMatrix[1][0]);
1320 mother->writeExtraInfo("prodVertSyy", RerrMatrix[1][1]);
1321 mother->writeExtraInfo("prodVertSyz", RerrMatrix[1][2]);
1322 mother->writeExtraInfo("prodVertSzx", RerrMatrix[2][0]);
1323 mother->writeExtraInfo("prodVertSzy", RerrMatrix[2][1]);
1324 mother->writeExtraInfo("prodVertSzz", RerrMatrix[2][2]);
1325 } else {
1326 mother->updateMomentum(mom, pos, errMatrix, prob);
1327 }
1328 return true;
1329 } else {return false;}
1330 }
1331 } else {
1332 nvert = rsf.fit();
1333 }
1334
1335 if (nvert > 0) {
1336 pos = rsf.getPos();
1337 RerrMatrix = rsf.getVertexErrorMatrix();
1338 double prob = rsf.getPValue();
1339 ROOT::Math::PxPyPzEVector mom(mother->get4Vector());
1340 TMatrixDSym errMatrix(7);
1341 for (int i = 0; i < 7; i++) {
1342 for (int j = 0; j < 7; j++) {
1343 if (i > 3 && j > 3) {errMatrix[i][j] = RerrMatrix[i - 4][j - 4];}
1344 else {errMatrix[i][j] = 0;}
1345 }
1346 }
1347 mother->updateMomentum(mom, pos, errMatrix, prob);
1348 } else {return false;}
1349
1350
1351 if (mothSel && nTrk > 1) {
1353 rf.addMother(mother);
1354 int nKfit = rf.fit();
1355 rf.updateMother();
1357
1358 if (nKfit > 0) {return true;}
1359 else return false;
1360 } else return true;
1361 }
1362 }
1363
1364 bool okFT = false;
1365 if (m_fitType == "vertex") {
1366 okFT = true;
1367 int nVert = rf.fit();
1368 rf.updateMother();
1369 if (m_decayString.empty() && m_updateDaughters == true) rf.updateDaughters();
1370 if (nVert != 1) return false;
1371 }
1372 if (m_fitType == "mass") {
1373 // add protection
1374 okFT = true;
1375 rf.setMassConstFit(true);
1376 rf.setVertFit(false);
1377 int nVert = rf.fit();
1378 rf.updateMother();
1379 if (nVert != 1) return false;
1380 };
1381 if (m_fitType == "massvertex") {
1382 okFT = true;
1383 rf.setMassConstFit(true);
1384 int nVert = rf.fit();
1385 rf.updateMother();
1386 if (m_decayString.empty() && m_updateDaughters == true) rf.updateDaughters();
1387 if (nVert != 1) return false;
1388 };
1389 if (!okFT) {
1390 B2FATAL("fitType : " << m_fitType << " ***invalid fit type ");
1391 }
1392
1393 return true;
1394}
1395
1397 const std::vector<const Particle*>& tracksVertex)
1398{
1399 bool isAll = false;
1400 if (mother->getNDaughters() == 0) return false;
1401
1402 int nNotIncluded = mother->getNDaughters();
1403
1404 for (unsigned i = 0; i < mother->getNDaughters(); i++) {
1405 bool dauOk = false;
1406 for (auto& vi : tracksVertex) {
1407 if (vi == mother->getDaughter(i)) {
1408 nNotIncluded = nNotIncluded - 1;
1409 dauOk = true;
1410 }
1411 }
1412 if (!dauOk) {
1413 if (allSelectedDaughters(mother->getDaughter(i), tracksVertex)) nNotIncluded--;
1414 }
1415 }
1416 if (nNotIncluded == 0) isAll = true;
1417 return isAll;
1418}
1419
1421{
1422 for (unsigned ichild = 0; ichild < particle->getNDaughters(); ichild++) {
1423 const Particle* child = particle->getDaughter(ichild);
1424 if (child->getNDaughters() > 0) addChildofParticletoKFit(kf, child);
1425 else {
1426 if (child->getPValue() < 0) return false; // error matrix not valid
1427
1428 kf.addParticle(child);
1429 }
1430 }
1431 return true;
1432}
1433
1435 std::vector<unsigned>& particleId)
1436{
1437 for (unsigned ichild = 0; ichild < particle->getNDaughters(); ichild++) {
1438 const Particle* child = particle->getDaughter(ichild);
1439 if (child->getNDaughters() > 0) {
1440 bool massconstraint = std::find(m_massConstraintList.begin(), m_massConstraintList.end(),
1441 std::abs(child->getPDGCode())) != m_massConstraintList.end();
1442 std::vector<unsigned> childId;
1443 addChildofParticletoMassKFit(kf, child, childId);
1444 if (massconstraint) kf.addMassConstraint(child->getPDGMass(), childId);
1445 particleId.insert(particleId.end(), childId.begin(), childId.end());
1446 } else {
1447 if (child->getPValue() < 0) return false; // error matrix not valid
1448 kf.addParticle(child);
1449 particleId.push_back(kf.getTrackCount() - 1);
1450 }
1451 }
1452 return true;
1453}
1454
1456{
1457 HepPoint3D pos = ROOTToCLHEP::getPoint3D(m_BeamSpotCenter);
1458 CLHEP::HepSymMatrix covMatrix = ROOTToCLHEP::getHepSymMatrix(m_beamSpotCov);
1459
1460 kv.setIpProfile(pos, covMatrix);
1461}
1462
1464{
1465 CLHEP::HepSymMatrix err(7, 0);
1466
1467 for (int i = 0; i < 3; i++) {
1468 for (int j = 0; j < 3; j++) {
1469 err[i + 4][j + 4] = m_beamSpotCov(i, j);
1470 }
1471 }
1472
1474 ROOT::Math::PxPyPzEVector iptube_mom = T.getBeamFourMomentum();
1475
1477 ROOTToCLHEP::getHepLorentzVector(iptube_mom),
1478 ROOTToCLHEP::getPoint3D(m_BeamSpotCenter),
1479 err,
1480 0.);
1481}
1482
1484{
1486
1487 ROOT::Math::XYZVector boostDir = T.getBoostVector().Unit();
1488
1489 TMatrixDSym beamSpotCov = m_beamSpotDB->getCovVertex();
1490 beamSpotCov(2, 2) = cut * cut;
1491 double thetab = boostDir.Theta();
1492 double phib = boostDir.Phi();
1493
1494 double stb = TMath::Sin(thetab);
1495 double ctb = TMath::Cos(thetab);
1496 double spb = TMath::Sin(phib);
1497 double cpb = TMath::Cos(phib);
1498
1499
1500 TMatrix rz(3, 3); rz(2, 2) = 1;
1501 rz(0, 0) = cpb; rz(0, 1) = spb;
1502 rz(1, 0) = -1 * spb; rz(1, 1) = cpb;
1503
1504 TMatrix ry(3, 3); ry(1, 1) = 1;
1505 ry(0, 0) = ctb; ry(0, 2) = -1 * stb;
1506 ry(2, 0) = stb; ry(2, 2) = ctb;
1507
1508 TMatrix r(3, 3); r.Mult(rz, ry);
1509 TMatrix rt(3, 3); rt.Transpose(r);
1510
1511 TMatrix TubePart(3, 3); TubePart.Mult(rt, beamSpotCov);
1512 TMatrix Tube(3, 3); Tube.Mult(TubePart, r);
1513
1514 m_beamSpotCov(0, 0) = Tube(0, 0); m_beamSpotCov(0, 1) = Tube(0, 1); m_beamSpotCov(0, 2) = Tube(0, 2);
1515 m_beamSpotCov(1, 0) = Tube(1, 0); m_beamSpotCov(1, 1) = Tube(1, 1); m_beamSpotCov(1, 2) = Tube(1, 2);
1516 m_beamSpotCov(2, 0) = Tube(2, 0); m_beamSpotCov(2, 1) = Tube(2, 1); m_beamSpotCov(2, 2) = Tube(2, 2);
1517}
1518
1520{
1521 TMatrixDSym beamSpotCov = m_beamSpotDB->getCovVertex();
1522 for (int i = 0; i < 3; i++)
1523 beamSpotCov(i, i) += width * width;
1524
1525 m_beamSpotCov = beamSpotCov;
1526}
1527
1529{
1530 double chi2TrackL = 0;
1531
1532 for (int iTrack = 0; iTrack < kv.getTrackCount(); iTrack++) {
1533
1534 analysis::KFitTrack trk_i = kv.getTrack(iTrack); // KFitTrack contains parameters before/after fit.
1535
1536 TMatrixFSym err = CLHEPToROOT::getTMatrixFSym(trk_i.getError(analysis::KFitConst::kBeforeFit)); // px, py, pz, E, x, y, z
1537
1538 ROOT::Math::XYZVector x_before = CLHEPToROOT::getXYZVector(trk_i.getPosition(analysis::KFitConst::kBeforeFit));
1539 ROOT::Math::XYZVector x_after = CLHEPToROOT::getXYZVector(trk_i.getPosition());
1540 ROOT::Math::XYZVector dPos = x_after - x_before;
1541
1543 ROOT::Math::XYZVector boost3 = T.getBoostVector().Unit();
1544 TVectorD boostD(0, 6, 0., 0., 0., 0., boost3.X(), boost3.Y(), boost3.Z(), "END");
1545
1546 double dLBoost = dPos.Dot(boost3);
1547
1548 chi2TrackL += TMath::Power(dLBoost, 2) / err.Similarity(boostD);
1549 }
1550 return chi2TrackL;
1551}
static ROOT::Math::XYZVector getFieldInTesla(const ROOT::Math::XYZVector &pos)
return the magnetic field at a given position in Tesla.
For each MCParticle with hits in the CDC, this class stores some summarising information on those hit...
Definition Btube.h:22
TMatrixFSym getTubeMatrix() const
Returns Btube matrix.
Definition Btube.h:73
Eigen::Matrix< double, 3, 1 > getTubeCenter() const
Returns Btube center.
Definition Btube.h:59
static const ParticleType pi0
neutral pion particle
Definition Const.h:674
static const double speedOfLight
[cm/ns]
Definition Const.h:695
static const ParticleType photon
photon particle
Definition Const.h:673
void setDescription(const std::string &description)
Sets the description of the module.
Definition Module.cc:214
void setPropertyFlags(unsigned int propertyFlags)
Sets the flags for the module properties.
Definition Module.cc:208
Module()
Constructor.
Definition Module.cc:30
@ c_ParallelProcessingCertified
This module can be run in parallel processing mode safely (All I/O must be done through the data stor...
Definition Module.h:80
Class to hold Lorentz transformations from/to CMS and boost vector.
ROOT::Math::PxPyPzEVector getBeamFourMomentum() const
Returns LAB four-momentum of e+e-, i.e.
ROOT::Math::XYZVector getBoostVector() const
Returns boost vector (beta=p/E)
TMatrixDSym m_beamSpotCov
Beam spot covariance matrix.
bool makeKMassMother(analysis::MassFitKFit &kv, Particle *p)
Update mother particle after mass fit using KFit.
std::string m_withConstraint
additional constraint on vertex
bool doKMassPointingVertexFit(Particle *p)
Mass-constrained vertex fit with additional pointing constraint using KFit.
void smearBeamSpot(double width)
smear beam spot covariance
bool makeMassKFourCMother(analysis::MassFourCFitKFit &kv, Particle *p)
Update mother particle after MassFourC fit using KFit.
void updateMapOfTrackAndDaughter(unsigned &l, std::vector< std::vector< unsigned > > &pars, std::vector< unsigned > &pard, std::vector< Particle * > &allparticles, const Particle *daughter)
update the map of daughter and tracks, find out which tracks belong to each daughter.
bool doKVertexFit(Particle *p, bool ipProfileConstraint, bool ipTubeConstraint)
Unconstrained vertex fit using KFit.
bool doKMassVertexFit(Particle *p)
Mass-constrained vertex fit using KFit.
virtual void initialize() override
Initialize the Module.
void addIPTubeToKFit(analysis::VertexFitKFit &kv)
Adds IPTube constraint to the vertex fit using KFit.
bool addChildofParticletoMassKFit(analysis::MassFourCFitKFit &kf, const Particle *particle, std::vector< unsigned > &particleId)
Adds given particle's child to the MassFourCFitKFit.
virtual void event() override
Event processor.
bool makeKMassPointingVertexMother(analysis::MassPointingVertexFitKFit &kv, Particle *p)
Update mother particle after mass-constrained vertex fit with additional pointing constraint using KF...
bool m_updateDaughters
flag for daughters update
std::string m_decayString
daughter particles selection
std::vector< int > m_massConstraintList
PDG codes of the particles to be mass constraint (massfourC)
bool m_hasCovMatrix
flag for mother covariance matrix (PseudoFitter)
bool makeKVertexMother(analysis::VertexFitKFit &kv, Particle *p)
Update mother particle after unconstrained vertex fit using KFit.
bool redoTwoPhotonDaughterMassFit(Particle *postFit, const Particle *preFit, const analysis::VertexFitKFit &kv)
Combines preFit particle and vertex information from vertex fit kv to create new postFit particle.
bool makeKMassVertexMother(analysis::MassVertexFitKFit &kv, Particle *p)
Update mother particle after mass-constrained vertex fit using KFit.
bool makeKFourCMother(analysis::FourCFitKFit &kv, Particle *p)
Update mother particle after FourC fit using KFit.
bool fillFitParticles(const Particle *mother, std::vector< const Particle * > &fitChildren, std::vector< const Particle * > &twoPhotonChildren)
Fills valid particle's children (with valid error matrix) in the vector of Particles that will enter ...
bool doKFourCFit(Particle *p)
FourC fit using KFit.
ROOT::Math::XYZVector m_BeamSpotCenter
Beam spot position.
virtual void beginRun() override
Called when entering a new run.
bool doKRecoilMassFit(Particle *p)
RecoilMass fit using KFit.
DBObjPtr< BeamSpot > m_beamSpotDB
Beam spot database object.
bool doVertexFit(Particle *p)
Main steering routine.
std::string m_vertexFitter
Vertex Fitter name.
bool doRaveFit(Particle *mother)
Fit using Rave.
bool makeKRecoilMassMother(analysis::RecoilMassKFit &kf, Particle *p)
Update mother particle after RecoilMass fit using KFit.
bool doKMassFit(Particle *p)
Mass fit using KFit.
double m_recoilMass
recoil mass for constraint
double m_confidenceLevel
required fit confidence level
bool doKMassFourCFit(Particle *p)
MassFourC fit using KFit.
DecayDescriptor m_decaydescriptor
Decay descriptor of decays to look for.
bool fillNotFitParticles(const Particle *mother, std::vector< const Particle * > &notFitChildren, const std::vector< const Particle * > &fitChildren)
Fills valid particle's children (with valid error matrix) in the vector of Particles that will not en...
void findConstraintBoost(double cut)
calculate iptube constraint (quasi cylinder along boost direction) for RAVE fit
double m_Bfield
magnetic field from data base
bool addChildofParticletoKFit(analysis::FourCFitKFit &kv, const Particle *particle)
Adds given particle's child to the FourCFitKFit.
double getChi2TracksLBoost(const analysis::VertexFitKFit &kv)
calculate the chi2 using only lboost information of tracks
double m_smearing
smearing width applied to IP tube
std::string m_fitType
type of the kinematic fit
void addIPProfileToKFit(analysis::VertexFitKFit &kv)
Adds IPProfile constraint to the vertex fit using KFit.
StoreObjPtr< ParticleList > m_plist
particle list
bool allSelectedDaughters(const Particle *mother, const std::vector< const Particle * > &tracksVertex)
check if all the Daughters (o grand-daughters) are selected for the vertex fit
Class to store reconstructed particles.
Definition Particle.h:76
void writeExtraInfo(const std::string &name, const double value)
Sets the user defined extraInfo.
Definition Particle.cc:1393
bool hasExtraInfo(const std::string &name) const
Return whether the extra info with the given name is set.
Definition Particle.cc:1351
const V0 * getV0() const
Returns the pointer to the V0 object that was used to create this Particle (if ParticleType == c_V0).
Definition Particle.cc:956
int getPDGCode(void) const
Returns PDG code.
Definition Particle.h:465
int getProperty() const
Returns particle property as a bit pattern The values are defined in the PropertyFlags enum and descr...
Definition Particle.h:518
unsigned getNDaughters(void) const
Returns number of daughter particles.
Definition Particle.h:747
double getPDGMass(void) const
Returns uncertainty on the invariant mass (requires valid momentum error matrix)
Definition Particle.cc:635
ROOT::Math::PxPyPzEVector get4Vector() const
Returns Lorentz vector.
Definition Particle.h:567
void addExtraInfo(const std::string &name, double value)
Sets the user-defined data of given name to the given value.
Definition Particle.cc:1421
void set4Vector(const ROOT::Math::PxPyPzEVector &p4)
Sets Lorentz vector.
Definition Particle.h:282
std::vector< Particle * > getDaughters() const
Returns a vector of pointers to daughter particles.
Definition Particle.cc:668
void updateMomentum(const ROOT::Math::PxPyPzEVector &p4, const ROOT::Math::XYZVector &vertex, const TMatrixFSym &errMatrix, double pValue)
Sets Lorentz vector, position, 7x7 error matrix and p-value.
Definition Particle.h:397
TMatrixFSym getMomentumVertexErrorMatrix() const
Returns 7x7 error matrix.
Definition Particle.cc:451
@ c_IsUnspecified
Ordinary particles.
Definition Particle.h:120
const Particle * getDaughter(unsigned i) const
Returns a pointer to the i-th daughter particle.
Definition Particle.cc:662
double getExtraInfo(const std::string &name) const
Return given value if set.
Definition Particle.cc:1374
double getFittedVertexX() const
Get the X coordinate of the fitted vertex.
Definition V0.h:75
double getFittedVertexY() const
Get the Y coordinate of the fitted vertex.
Definition V0.h:78
double getFittedVertexZ() const
Get the Z coordinate of the fitted vertex.
Definition V0.h:81
ROOT::Math::XYZVector getFittedVertexPosition() const
Get the fitted vertex position.
Definition V0.h:84
FourCFitKFit is a derived class from KFitBase to perform 4 momentum-constraint kinematical fit.
enum KFitError::ECode setFourMomentum(const ROOT::Math::PxPyPzEVector &m)
Set an 4 Momentum for the FourC-constraint fit.
double getCHIsq(void) const override
Get a chi-square of the fit.
enum KFitError::ECode updateMother(Particle *mother)
Update mother particle.
enum KFitError::ECode doFit(void)
Perform a four momentum-constraint fit.
const CLHEP::HepSymMatrix getTrackError(const int id) const
Get an error matrix of the track.
Definition KFitBase.cc:168
const CLHEP::HepLorentzVector getTrackMomentum(const int id) const
Get a Lorentz vector of the track.
Definition KFitBase.cc:154
const HepPoint3D getTrackPosition(const int id) const
Get a position of the track.
Definition KFitBase.cc:161
virtual int getNDF(void) const
Get an NDF of the fit.
Definition KFitBase.cc:114
enum KFitError::ECode setMagneticField(const double mf)
Change a magnetic field from the default value KFitConst::kDefaultMagneticField.
Definition KFitBase.cc:93
const KFitTrack getTrack(const int id) const
Get a specified track object.
Definition KFitBase.cc:175
enum KFitError::ECode addParticle(const Particle *particle)
Add a particle to the fitter.
Definition KFitBase.cc:59
int getTrackCount(void) const
Get the number of added tracks.
Definition KFitBase.cc:107
ECode
ECode is a error code enumerate.
Definition KFitError.h:33
KFitTrack is a container of the track information (Lorentz vector, position, and error matrix),...
Definition KFitTrack.h:38
const CLHEP::HepSymMatrix getError(const int flag=KFitConst::kAfterFit) const
Get an error matrix of the track.
Definition KFitTrack.cc:171
const HepPoint3D getPosition(const int flag=KFitConst::kAfterFit) const
Get a position of the track.
Definition KFitTrack.cc:163
MassFitKFit is a derived class from KFitBase to perform mass-constraint kinematical fit.
Definition MassFitKFit.h:33
enum KFitError::ECode setVertex(const HepPoint3D &v)
Set an initial vertex position for the mass-constraint fit.
enum KFitError::ECode updateMother(Particle *mother)
Update mother particle.
enum KFitError::ECode doFit(void)
Perform a mass-constraint fit.
enum KFitError::ECode setInvariantMass(const double m)
Set an invariant mass for the mass-constraint fit.
enum KFitError::ECode setVertexError(const CLHEP::HepSymMatrix &e)
Set an initial vertex error matrix for the mass-constraint fit.
const HepPoint3D getVertex(const int flag=KFitConst::kAfterFit) const
Get a vertex position.
MassFourCFitKFit is a derived class from KFitBase to perform mass and 4 momentum-constraint kinematic...
enum KFitError::ECode setFourMomentum(const ROOT::Math::PxPyPzEVector &m)
Set an 4 Momentum for the mass-four-constraint fit.
enum KFitError::ECode addMassConstraint(const double m, std::vector< unsigned > &childTrackId)
Set an invariant mass of daughter particle for the mass-four-momentum-constraint fit.
double getCHIsq(void) const override
Get a chi-square of the fit.
enum KFitError::ECode updateMother(Particle *mother)
Update mother particle.
enum KFitError::ECode doFit(void)
Perform a mass-four-momentum-constraint fit.
MassPointingVertexFitKFit is a derived class from KFitBase It performs a kinematical fit with three c...
enum KFitError::ECode updateMother(Particle *mother)
Update mother particle.
enum KFitError::ECode doFit(void)
Perform a mass-vertex-pointing constraint fit.
enum KFitError::ECode setInvariantMass(const double m)
Set an invariant mass for the mass-vertex-pointing constraint fit.
enum KFitError::ECode setProductionVertex(const HepPoint3D &v)
Set the production vertex of the particle.
const HepPoint3D getVertex(const int flag=KFitConst::kAfterFit) const
Get a vertex position.
MassVertexFitKFit is a derived class from KFitBase to perform mass-vertex-constraint kinematical fit.
enum KFitError::ECode setInitialVertex(const HepPoint3D &v)
Set an initial vertex point for the mass-vertex constraint fit.
enum KFitError::ECode updateMother(Particle *mother)
Update mother particle.
enum KFitError::ECode doFit(void)
Perform a mass-vertex-constraint fit.
enum KFitError::ECode setInvariantMass(const double m)
Set an invariant mass for the mass-vertex constraint fit.
const HepPoint3D getVertex(const int flag=KFitConst::kAfterFit) const
Get a vertex position.
The RaveKinematicVertexFitter class is part of the RaveInterface together with RaveSetup.
void addMother(const Particle *aMotherParticlePtr)
All daughters of the argument of this function will be used as input for the vertex fit.
ROOT::Math::XYZVector getPos()
get the position of the fitted vertex.
void addTrack(const Particle *aParticlePtr)
add a track (in the format of a Particle) to set of tracks that should be fitted to a vertex
void setMother(const Particle *aMotherParticlePtr)
Set Mother particle for Vertex/momentum update.
int fit()
do the kinematic vertex fit with all tracks previously added with the addTrack or addMother function.
void setVertFit(bool isVertFit=true)
Set vertex fit: set false in case of mass fit only.
void setMassConstFit(bool isConstFit=true)
Set mass constrained fit.
double getPValue()
get the p value of the fitted vertex.
TMatrixDSym getVertexErrorMatrix()
get the covariance matrix (3x3) of the of the fitted vertex position.
void updateDaughters()
update the Daughters particles
void unsetBeamSpot()
unset beam spot constraint
Definition RaveSetup.cc:82
static void initialize(int verbosity=1, double MagneticField=1.5)
Set everything up so everything needed for vertex fitting is there.
Definition RaveSetup.cc:35
static RaveSetup * getInstance()
get the pointer to the instance to get/set any of options stored in RaveSetup
Definition RaveSetup.h:40
void setBeamSpot(const ROOT::Math::XYZVector &beamSpot, const TMatrixDSym &beamSpotCov)
The beam spot position and covariance is known you can set it here so that and a vertex in the beam s...
Definition RaveSetup.cc:74
void reset()
frees memory allocated by initialize().
Definition RaveSetup.cc:60
The RaveVertexFitter class is part of the RaveInterface together with RaveSetup.
TMatrixDSym getCov(VecSize vertexId=0) const
get the covariance matrix (3x3) of the of the fitted vertex position.
int fit(std::string options="default")
do the vertex fit with all tracks previously added with the addTrack or addMother function.
ROOT::Math::XYZVector getPos(VecSize vertexId=0) const
get the position of the fitted vertex.
void addTrack(const Particle *const aParticlePtr)
add a track (in the format of a Particle) to set of tracks that should be fitted to a vertex
double getPValue(VecSize vertexId=0) const
get the p value of the fitted vertex.
RecoilMassKFit is a derived class from KFitBase to perform a kinematical fit with a recoil mass const...
enum KFitError::ECode setFourMomentum(const ROOT::Math::PxPyPzEVector &m)
Set a recoil mass .
double getCHIsq(void) const override
Get a chi-square of the fit.
enum KFitError::ECode updateMother(Particle *mother)
Update mother particle.
enum KFitError::ECode setRecoilMass(const double m)
Set an invariant mass for the four momentum-constraint fit.
enum KFitError::ECode doFit(void)
Perform a recoil-mass constraint fit.
const HepPoint3D getVertex(const int flag=KFitConst::kAfterFit) const
Get a vertex position.
VertexFitKFit is a derived class from KFitBase to perform vertex-constraint kinematical fit.
enum KFitError::ECode setIpTubeProfile(const CLHEP::HepLorentzVector &p, const HepPoint3D &x, const CLHEP::HepSymMatrix &e, const double q)
Set a virtual IP-tube track for the vertex constraint fit.
enum KFitError::ECode setInitialVertex(const HepPoint3D &v)
Set an initial vertex point for the vertex-vertex constraint fit.
enum KFitError::ECode updateMother(Particle *mother)
Update mother particle.
const CLHEP::HepSymMatrix getVertexError(void) const
Get a fitted vertex error matrix.
enum KFitError::ECode doFit(void)
Perform a vertex-constraint fit.
enum KFitError::ECode setIpProfile(const HepPoint3D &ip, const CLHEP::HepSymMatrix &ipe)
Set an IP-ellipsoid shape for the vertex constraint fit.
const HepPoint3D getVertex(const int flag=KFitConst::kAfterFit) const
Get a vertex position.
void addParam(const std::string &name, T &paramVariable, const std::string &description, const T &defaultValue)
Adds a new parameter to the module.
Definition Module.h:559
#define REG_MODULE(moduleName)
Register the given module (without 'Module' suffix) with the framework.
Definition Module.h:649
void copyDaughters(Particle *mother)
Function copies all (grand-)^n-daughter particles of the argument mother Particle.
Abstract base class for different kinds of events.
STL namespace.
static const int kBeforeFit
Input parameter to specify before-fit when setting/getting a track attribute.
Definition KFitConst.h:33