10 #include <analysis/variables/V0DaughterTrackVariables.h>
13 #include <analysis/VariableManager/Manager.h>
15 #include <analysis/dataobjects/Particle.h>
16 #include <analysis/variables/TrackVariables.h>
19 #include <framework/dataobjects/Helix.h>
22 #include <mdst/dataobjects/Track.h>
23 #include <mdst/dataobjects/MCParticle.h>
24 #include <mdst/dataobjects/TrackFitResult.h>
25 #include <mdst/dataobjects/HitPatternCDC.h>
26 #include <mdst/dataobjects/HitPatternVXD.h>
29 #include <framework/logging/Logger.h>
44 static const double realNaN = std::numeric_limits<double>::quiet_NaN();
48 double getV0DaughterTrackDetNHits(
const Particle* particle,
const double daughterID,
const Const::EDetector& det)
50 auto daughter = particle->getDaughter(daughterID);
51 return trackNHits(daughter, det);
54 double v0DaughterTrackNCDCHits(
const Particle* part,
const std::vector<double>& daughterID)
56 return getV0DaughterTrackDetNHits(part, daughterID[0], Const::EDetector::CDC);
59 double v0DaughterTrackNSVDHits(
const Particle* part,
const std::vector<double>& daughterID)
61 return getV0DaughterTrackDetNHits(part, daughterID[0], Const::EDetector::SVD);
64 double v0DaughterTrackNPXDHits(
const Particle* part,
const std::vector<double>& daughterID)
66 return getV0DaughterTrackDetNHits(part, daughterID[0], Const::EDetector::PXD);
69 double v0DaughterTrackNVXDHits(
const Particle* part,
const std::vector<double>& daughterID)
71 return v0DaughterTrackNPXDHits(part, daughterID) + v0DaughterTrackNSVDHits(part, daughterID);
74 double v0DaughterTrackNRemovedHits(
const Particle* part,
const std::vector<double>& daughterID)
77 if (part->getParticleSource() != Particle::EParticleSourceObject::c_V0)
79 auto daughter = part->getDaughter(daughterID[0]);
80 const Track* track = daughter->getTrack();
81 if (!track)
return realNaN;
82 const TrackFitResult* trackFit = track->getTrackFitResultWithClosestMass(Const::ChargedStable(abs(
83 daughter->getPDGCode())));
84 if (!trackFit)
return realNaN;
85 double nHitsBeforeRemoval = trackFit->getHitPatternCDC().getNHits()
86 + trackFit->getHitPatternVXD().getNSVDHits()
87 + trackFit->getHitPatternVXD().getNPXDHits();
88 double nHitsAfterRemoval = trackNVXDHits(daughter) + trackNCDCHits(daughter);
89 return nHitsBeforeRemoval - nHitsAfterRemoval;
92 double v0DaughterTrackFirstSVDLayer(
const Particle* part,
const std::vector<double>& daughterID)
94 auto daughter = part->getDaughter(daughterID[0]);
95 return trackFirstSVDLayer(daughter);
98 double v0DaughterTrackFirstPXDLayer(
const Particle* part,
const std::vector<double>& daughterID)
100 auto daughter = part->getDaughter(daughterID[0]);
101 return trackFirstPXDLayer(daughter);
104 double v0DaughterTrackFirstCDCLayer(
const Particle* part,
const std::vector<double>& daughterID)
106 auto daughter = part->getDaughter(daughterID[0]);
107 return trackFirstCDCLayer(daughter);
110 double v0DaughterTrackLastCDCLayer(
const Particle* part,
const std::vector<double>& daughterID)
112 auto daughter = part->getDaughter(daughterID[0]);
113 return trackLastCDCLayer(daughter);
116 double v0DaughterTrackPValue(
const Particle* part,
const std::vector<double>& daughterID)
118 auto daughter = part->getDaughter(daughterID[0]);
119 return trackPValue(daughter);
122 double v0DaughterTrackD0(
const Particle* part,
const std::vector<double>& daughterID)
124 auto daughter = part->getDaughter(daughterID[0]);
125 return trackD0(daughter);
128 double v0DaughterTrackPhi0(
const Particle* part,
const std::vector<double>& daughterID)
130 auto daughter = part->getDaughter(daughterID[0]);
131 return trackPhi0(daughter);
134 double v0DaughterTrackOmega(
const Particle* part,
const std::vector<double>& daughterID)
136 auto daughter = part->getDaughter(daughterID[0]);
137 return trackOmega(daughter);
140 double v0DaughterTrackZ0(
const Particle* part,
const std::vector<double>& daughterID)
142 auto daughter = part->getDaughter(daughterID[0]);
143 return trackZ0(daughter);
146 double v0DaughterTrackTanLambda(
const Particle* part,
const std::vector<double>& daughterID)
148 auto daughter = part->getDaughter(daughterID[0]);
149 return trackTanLambda(daughter);
152 double v0DaughterTrackD0Error(
const Particle* part,
const std::vector<double>& daughterID)
154 auto daughter = part->getDaughter(daughterID[0]);
155 return trackD0Error(daughter);
158 double v0DaughterTrackPhi0Error(
const Particle* part,
const std::vector<double>& daughterID)
160 auto daughter = part->getDaughter(daughterID[0]);
161 return trackPhi0Error(daughter);
164 double v0DaughterTrackOmegaError(
const Particle* part,
const std::vector<double>& daughterID)
166 auto daughter = part->getDaughter(daughterID[0]);
167 return trackOmegaError(daughter);
170 double v0DaughterTrackZ0Error(
const Particle* part,
const std::vector<double>& daughterID)
172 auto daughter = part->getDaughter(daughterID[0]);
173 return trackZ0Error(daughter);
176 double v0DaughterTrackTanLambdaError(
const Particle* part,
const std::vector<double>& daughterID)
178 auto daughter = part->getDaughter(daughterID[0]);
179 return trackTanLambdaError(daughter);
182 double v0DaughterD0(
const Particle* particle,
const std::vector<double>& daughterID)
185 return std::numeric_limits<float>::quiet_NaN();
187 ROOT::Math::XYZVector v0Vertex = particle->getVertex();
189 const Particle* daug = particle->getDaughter(daughterID[0]);
191 const TrackFitResult* trackFit = daug->getTrackFitResult();
192 if (!trackFit)
return std::numeric_limits<float>::quiet_NaN();
194 UncertainHelix helix = trackFit->getUncertainHelix();
195 helix.passiveMoveBy(v0Vertex);
197 return helix.getD0();
200 double v0DaughterD0Diff(
const Particle* particle)
202 return v0DaughterD0(particle, {0}) - v0DaughterD0(particle, {1});
205 double v0DaughterZ0(
const Particle* particle,
const std::vector<double>& daughterID)
208 return std::numeric_limits<float>::quiet_NaN();
210 ROOT::Math::XYZVector v0Vertex = particle->getVertex();
212 const Particle* daug = particle->getDaughter(daughterID[0]);
214 const TrackFitResult* trackFit = daug->getTrackFitResult();
215 if (!trackFit)
return std::numeric_limits<float>::quiet_NaN();
217 UncertainHelix helix = trackFit->getUncertainHelix();
218 helix.passiveMoveBy(v0Vertex);
220 return helix.getZ0();
223 double v0DaughterZ0Diff(
const Particle* particle)
225 return v0DaughterZ0(particle, {0}) - v0DaughterZ0(particle, {1});
230 double getHelixParameterPullOfV0DaughterWithTrueVertexAsPivotAtIndex(
const Particle* particle,
const double daughterID,
233 if (!particle) {
return std::numeric_limits<double>::quiet_NaN(); }
235 const int dID = int(std::lround(daughterID));
236 if (not(dID == 0 || dID == 1)) {
return std::numeric_limits<double>::quiet_NaN(); }
239 if (!mcparticle_v0) {
return std::numeric_limits<double>::quiet_NaN(); }
241 if (!(particle->getDaughter(dID))) {
return std::numeric_limits<double>::quiet_NaN(); }
244 if (!mcparticle) {
return std::numeric_limits<double>::quiet_NaN(); }
246 const TrackFitResult* trackFit = particle->getDaughter(dID)->getTrackFitResult();
247 if (!trackFit) {
return std::numeric_limits<double>::quiet_NaN(); }
250 const ROOT::Math::XYZVector mcProdVertex = mcparticle->getVertex();
251 const ROOT::Math::XYZVector mcMomentum = mcparticle->getMomentum();
252 const double mcParticleCharge = mcparticle->getCharge();
253 const double BzAtProdVertex = BFieldManager::getFieldInTesla(mcProdVertex).Z();
254 Helix mcHelix = Helix(mcProdVertex, mcMomentum, mcParticleCharge, BzAtProdVertex);
255 mcHelix.passiveMoveBy(mcProdVertex);
256 const std::vector<double> mcHelixPars = { mcHelix.getD0(), mcHelix.getPhi0(), mcHelix.getOmega(),
257 mcHelix.getZ0(), mcHelix.getTanLambda()
261 UncertainHelix measHelix = trackFit->getUncertainHelix();
262 measHelix.passiveMoveBy(mcProdVertex);
263 const TMatrixDSym measCovariance = measHelix.getCovariance();
264 const std::vector<double> measHelixPars = {measHelix.getD0(), measHelix.getPhi0(), measHelix.getOmega(),
265 measHelix.getZ0(), measHelix.getTanLambda()
267 const std::vector<double> measErrSquare = {measCovariance[0][0], measCovariance[1][1], measCovariance[2][2],
268 measCovariance[3][3], measCovariance[4][4]
271 if (measErrSquare.at(tauIndex) > 0)
272 return (mcHelixPars.at(tauIndex) - measHelixPars.at(tauIndex)) / std::sqrt(measErrSquare.at(tauIndex));
274 return std::numeric_limits<double>::quiet_NaN();
277 double v0DaughterHelixWithTrueVertexAsPivotD0Pull(
const Particle* part,
const std::vector<double>& daughterID)
279 return getHelixParameterPullOfV0DaughterWithTrueVertexAsPivotAtIndex(part, daughterID[0], 0);
282 double v0DaughterHelixWithTrueVertexAsPivotPhi0Pull(
const Particle* part,
const std::vector<double>& daughterID)
284 return getHelixParameterPullOfV0DaughterWithTrueVertexAsPivotAtIndex(part, daughterID[0], 1);
287 double v0DaughterHelixWithTrueVertexAsPivotOmegaPull(
const Particle* part,
const std::vector<double>& daughterID)
289 return getHelixParameterPullOfV0DaughterWithTrueVertexAsPivotAtIndex(part, daughterID[0], 2);
292 double v0DaughterHelixWithTrueVertexAsPivotZ0Pull(
const Particle* part,
const std::vector<double>& daughterID)
294 return getHelixParameterPullOfV0DaughterWithTrueVertexAsPivotAtIndex(part, daughterID[0], 3);
297 double v0DaughterHelixWithTrueVertexAsPivotTanLambdaPull(
const Particle* part,
const std::vector<double>& daughterID)
299 return getHelixParameterPullOfV0DaughterWithTrueVertexAsPivotAtIndex(part, daughterID[0], 4);
302 double v0DaughterHelixWithOriginAsPivotD0Pull(
const Particle* part,
const std::vector<double>& daughterID)
304 auto daughter = part->getDaughter(daughterID[0]);
305 return getHelixD0Pull(daughter);
308 double v0DaughterHelixWithOriginAsPivotPhi0Pull(
const Particle* part,
const std::vector<double>& daughterID)
310 auto daughter = part->getDaughter(daughterID[0]);
311 return getHelixPhi0Pull(daughter);
314 double v0DaughterHelixWithOriginAsPivotOmegaPull(
const Particle* part,
const std::vector<double>& daughterID)
316 auto daughter = part->getDaughter(daughterID[0]);
317 return getHelixOmegaPull(daughter);
320 double v0DaughterHelixWithOriginAsPivotZ0Pull(
const Particle* part,
const std::vector<double>& daughterID)
322 auto daughter = part->getDaughter(daughterID[0]);
323 return getHelixZ0Pull(daughter);
326 double v0DaughterHelixWithOriginAsPivotTanLambdaPull(
const Particle* part,
const std::vector<double>& daughterID)
328 auto daughter = part->getDaughter(daughterID[0]);
329 return getHelixTanLambdaPull(daughter);
332 double v0DaughterTrackParam5AtIPPerigee(
const Particle* part,
const std::vector<double>& params)
334 auto daughter = part->getDaughter(params[0]);
336 return std::numeric_limits<double>::quiet_NaN();
338 auto trackFit = daughter->getTrackFitResult();
340 return std::numeric_limits<double>::quiet_NaN();
343 const int paramID = int(std::lround(params[1]));
344 if (not(0 <= paramID && paramID < 5))
345 return std::numeric_limits<double>::quiet_NaN();
347 std::vector<float> tau = trackFit->getTau();
351 double v0DaughterTrackParamCov5x5AtIPPerigee(
const Particle* part,
const std::vector<double>& params)
353 auto daughter = part->getDaughter(params[0]);
355 return std::numeric_limits<double>::quiet_NaN();
357 auto trackFit = daughter->getTrackFitResult();
359 return std::numeric_limits<double>::quiet_NaN();
362 const int paramID = int(std::lround(params[1]));
363 if (not(0 <= paramID && paramID < 15))
364 return std::numeric_limits<double>::quiet_NaN();
366 std::vector<float> cov = trackFit->getCov();
370 int convertedPhotonErrorChecks(
const Particle* gamma,
const std::vector<double>& daughterIndices)
373 if (daughterIndices.size() != 2) {
374 B2ERROR(
"Invalid number of daughter indices. Please specify exactly two valid daughter indices.");
379 int daughterIndex1 = int(daughterIndices[0]);
380 int daughterIndex2 = int(daughterIndices[1]);
381 if (
int(gamma->getNDaughters()) <= std::max(daughterIndex1, daughterIndex2)) {
382 B2ERROR(
"Invalid daughter indices provided. Particle does not have that many daughters.");
387 if (!gamma->getDaughter(daughterIndex1)->getTrack()) {
388 B2ERROR(
"There is no track associated with daughter index " << daughterIndex1);
391 if (!gamma->getDaughter(daughterIndex2)->getTrack()) {
392 B2ERROR(
"There is no track associated with daughter index " << daughterIndex2);
397 if (fabs(gamma->getDaughter(daughterIndex1)->getPDGCode()) != 11) {
398 B2INFO(
"The first track provided has not been reconstructed as an electron/positron. It has PDG code " << gamma->getDaughter(
399 daughterIndex1)->getPDGCode() <<
". However, this is still fully admissible.");
401 if (fabs(gamma->getDaughter(daughterIndex2)->getPDGCode()) != 11) {
402 B2INFO(
"The second track provided has not been reconstructed as an electron/positron. It has PDG code " << gamma->getDaughter(
403 daughterIndex1)->getPDGCode() <<
".However, this is still fully admissible.");
410 int convertedPhotonLoadHelixParams(
const Particle* gamma,
int daughterIndex1,
int daughterIndex2,
double& Phi01,
double& D01,
411 double& Omega1,
double& Z01,
double& TanLambda1,
double& Phi02,
double& D02,
double& Omega2,
double& Z02,
416 Helix e1Helix = gamma->getDaughter(daughterIndex1)->getTrackFitResult()->getHelix();
418 Phi01 = e1Helix.getPhi0();
419 D01 = e1Helix.getD0() ;
420 Omega1 = e1Helix.getOmega();
421 Z01 = e1Helix.getZ0();
422 TanLambda1 = e1Helix.getTanLambda();
425 Helix e2Helix = gamma->getDaughter(daughterIndex2)->getTrackFitResult()->getHelix();
427 Phi02 = e2Helix.getPhi0();
428 D02 = e2Helix.getD0() ;
429 Omega2 = e2Helix.getOmega();
430 Z02 = e2Helix.getZ0();
431 TanLambda2 = e2Helix.getTanLambda();
434 if (Omega1 == 0) {
return -1;}
435 else if (Omega2 == 0) {
return -2;}
440 double convertedPhotonInvariantMass(
const Particle* gamma,
const std::vector<double>& daughterIndices)
443 int errFlag = convertedPhotonErrorChecks(gamma, daughterIndices);
444 if (errFlag == -1) {
return std::numeric_limits<double>::quiet_NaN();}
447 double Phi01, D01, Omega1, Z01, TanLambda1, Phi02, D02, Omega2, Z02, TanLambda2;
448 int daughterIndex1 = int(daughterIndices[0]);
449 int daughterIndex2 = int(daughterIndices[1]);
450 convertedPhotonLoadHelixParams(gamma, daughterIndex1, daughterIndex2, Phi01, D01, Omega1, Z01, TanLambda1, Phi02, D02,
451 Omega2, Z02, TanLambda2);
456 double sinlam1 = TanLambda1 / sqrt(1 + (TanLambda1 * TanLambda1));
457 double coslam1 = 1 / sqrt(1 + (TanLambda1 * TanLambda1));
458 double sinlam2 = TanLambda2 / sqrt(1 + (TanLambda2 * TanLambda2));
459 double coslam2 = 1 / sqrt(1 + (TanLambda2 * TanLambda2));
463 double p1 = gamma->getDaughter(daughterIndex1)->getMomentumMagnitude();
464 double pt1 = p1 * coslam1, pz1 = p1 * sinlam1;
465 double e1 = sqrt((p1 * p1) + (Const::electronMass * Const::electronMass));
467 double p2 = gamma->getDaughter(daughterIndex2)->getMomentumMagnitude();
468 double pt2 = p2 * coslam2, pz2 = p2 * sinlam2;
469 double e2 = sqrt((p2 * p2) + (Const::electronMass * Const::electronMass));
472 double vtxMass = sqrt(pow(e1 + e2, 2.0) - pow(pt1 + pt2, 2.0) - pow(pz1 + pz2, 2.0));
476 double convertedPhotonDelTanLambda(
const Particle* gamma,
const std::vector<double>& daughterIndices)
479 int errFlag = convertedPhotonErrorChecks(gamma, daughterIndices);
480 if (errFlag == -1) {
return std::numeric_limits<double>::quiet_NaN();}
483 double Phi01, D01, Omega1, Z01, TanLambda1, Phi02, D02, Omega2, Z02, TanLambda2;
484 int daughterIndex1 = int(daughterIndices[0]);
485 int daughterIndex2 = int(daughterIndices[1]);
486 convertedPhotonLoadHelixParams(gamma, daughterIndex1, daughterIndex2, Phi01, D01, Omega1, Z01, TanLambda1, Phi02, D02,
487 Omega2, Z02, TanLambda2);
492 return (TanLambda2 - TanLambda1);
495 double convertedPhotonDelR(
const Particle* gamma,
const std::vector<double>& daughterIndices)
498 int errFlag = convertedPhotonErrorChecks(gamma, daughterIndices);
499 if (errFlag == -1) {
return std::numeric_limits<double>::quiet_NaN();}
502 double Phi01, D01, Omega1, Z01, TanLambda1, Phi02, D02, Omega2, Z02, TanLambda2;
503 int daughterIndex1 = int(daughterIndices[0]);
504 int daughterIndex2 = int(daughterIndices[1]);
505 errFlag = convertedPhotonLoadHelixParams(gamma, daughterIndex1, daughterIndex2, Phi01, D01, Omega1, Z01, TanLambda1, Phi02, D02,
506 Omega2, Z02, TanLambda2);
508 B2ERROR(
"First track provided has curvature zero. Calculation of convertedPhotonDelR failed.");
509 return std::numeric_limits<double>::quiet_NaN();
512 B2ERROR(
"Second track provided has curvature zero. Calculation of convertedPhotonDelR failed.");
513 return std::numeric_limits<double>::quiet_NaN();
517 double radius1 = 1 / Omega1;
518 double radius2 = 1 / Omega2;
520 TVector2 center1((radius1 + D01) * sin(Phi01), -1 * (radius1 + D01) * cos(Phi01));
521 TVector2 center2((radius2 + D02) * sin(Phi02), -1 * (radius2 + D02) * cos(Phi02));
522 TVector2 cenDiff = center1 - center2;
524 double delR = fabs(radius1) + fabs(radius2) - cenDiff.Mod();
528 std::pair<double, double> convertedPhotonZ1Z2(
const Particle* gamma,
const std::vector<double>& daughterIndices)
531 int errFlag = convertedPhotonErrorChecks(gamma, daughterIndices);
532 if (errFlag == -1) {
return std::pair<double, double>(std::numeric_limits<double>::quiet_NaN(), std::numeric_limits<double>::quiet_NaN());}
535 double Phi01, D01, Omega1, Z01, TanLambda1, Phi02, D02, Omega2, Z02, TanLambda2;
536 int daughterIndex1 = int(daughterIndices[0]);
537 int daughterIndex2 = int(daughterIndices[1]);
538 errFlag = convertedPhotonLoadHelixParams(gamma, daughterIndex1, daughterIndex2, Phi01, D01, Omega1, Z01, TanLambda1, Phi02, D02,
539 Omega2, Z02, TanLambda2);
541 B2ERROR(
"First track provided has curvature zero. Calculation of convertedPhotonZ1Z2 failed.");
542 return std::pair<double, double>(std::numeric_limits<double>::quiet_NaN(), std::numeric_limits<double>::quiet_NaN());
545 B2ERROR(
"Second track provided has curvature zero. Calculation of convertedPhotonZ1Z2 failed.");
546 return std::pair<double, double>(std::numeric_limits<double>::quiet_NaN(), std::numeric_limits<double>::quiet_NaN());
551 double radius1 = 1 / Omega1;
552 double radius2 = 1 / Omega2;
554 TVector2 center1((radius1 + D01) * sin(Phi01), -1 * (radius1 + D01) * cos(Phi01));
555 TVector2 center2((radius2 + D02) * sin(Phi02), -1 * (radius2 + D02) * cos(Phi02));
557 TVector2 n1 = center1 - center2; n1 = n1.Unit();
558 TVector2 n2 = -1 * n1;
559 n1 = copysign(1.0, Omega1) * n1;
560 n2 = copysign(1.0, Omega2) * n2;
563 double phiN1 = atan2(n1.X(), -n1.Y());
564 double phiN2 = atan2(n2.X(), -n2.Y());
565 double Phi01Intersect = phiN1 - Phi01;
566 double Phi02Intersect = phiN2 - Phi02;
568 double z1 = Z01 - (radius1 * TanLambda1 * Phi01Intersect);
569 double z2 = Z02 - (radius2 * TanLambda2 * Phi02Intersect);
570 std::pair<double, double> z1z2(z1, z2);
574 double convertedPhotonDelZ(
const Particle* gamma,
const std::vector<double>& daughterIndices)
576 std::pair<double, double> z1z2 = convertedPhotonZ1Z2(gamma, daughterIndices);
577 double z1 = z1z2.first;
double z2 = z1z2.second;
581 double convertedPhotonZ(
const Particle* gamma,
const std::vector<double>& daughterIndices)
583 std::pair<double, double> z1z2 = convertedPhotonZ1Z2(gamma, daughterIndices);
584 double z1 = z1z2.first;
double z2 = z1z2.second;
585 return (z1 + z2) * 0.5;
588 TVector2 convertedPhotonXY(
const Particle* gamma,
const std::vector<double>& daughterIndices)
591 int errFlag = convertedPhotonErrorChecks(gamma, daughterIndices);
592 if (errFlag == -1) {
return TVector2(std::numeric_limits<double>::quiet_NaN(), std::numeric_limits<double>::quiet_NaN());}
595 double Phi01, D01, Omega1, Z01, TanLambda1, Phi02, D02, Omega2, Z02, TanLambda2;
596 int daughterIndex1 = int(daughterIndices[0]);
597 int daughterIndex2 = int(daughterIndices[1]);
598 errFlag = convertedPhotonLoadHelixParams(gamma, daughterIndex1, daughterIndex2, Phi01, D01, Omega1, Z01, TanLambda1, Phi02, D02,
599 Omega2, Z02, TanLambda2);
601 B2ERROR(
"First track provided has curvature zero. Calculation of convertedPhotonXY failed.");
602 return TVector2(std::numeric_limits<double>::quiet_NaN(), std::numeric_limits<double>::quiet_NaN());
605 B2ERROR(
"Second track provided has curvature zero. Calculation of convertedPhotonXY failed.");
606 return TVector2(std::numeric_limits<double>::quiet_NaN(), std::numeric_limits<double>::quiet_NaN());
610 double radius1 = 1 / Omega1;
611 double radius2 = 1 / Omega2;
613 TVector2 center1((radius1 + D01) * sin(Phi01), -1 * (radius1 + D01) * cos(Phi01));
614 TVector2 center2((radius2 + D02) * sin(Phi02), -1 * (radius2 + D02) * cos(Phi02));
615 TVector2 cenDiff = center2 - center1;
616 double delR = fabs(radius1) + fabs(radius2) - cenDiff.Mod();
619 TVector2 n1 = cenDiff.Unit();
620 TVector2 vtxXY = center1 + ((fabs(radius1) - (delR / 2)) * n1);
624 double convertedPhotonX(
const Particle* gamma,
const std::vector<double>& daughterIndices)
626 auto vtxXY = convertedPhotonXY(gamma, daughterIndices);
630 double convertedPhotonY(
const Particle* gamma,
const std::vector<double>& daughterIndices)
632 auto vtxXY = convertedPhotonXY(gamma, daughterIndices);
636 double convertedPhotonRho(
const Particle* gamma,
const std::vector<double>& daughterIndices)
638 auto vtxXY = convertedPhotonXY(gamma, daughterIndices);
642 B2Vector3D convertedPhoton3Momentum(
const Particle* gamma,
const std::vector<double>& daughterIndices)
645 int errFlag = convertedPhotonErrorChecks(gamma, daughterIndices);
646 if (errFlag == -1) {
return B2Vector3D(std::numeric_limits<double>::quiet_NaN(), std::numeric_limits<double>::quiet_NaN(), std::numeric_limits<double>::quiet_NaN());}
649 double Phi01, D01, Omega1, Z01, TanLambda1, Phi02, D02, Omega2, Z02, TanLambda2;
650 int daughterIndex1 = int(daughterIndices[0]);
651 int daughterIndex2 = int(daughterIndices[1]);
652 errFlag = convertedPhotonLoadHelixParams(gamma, daughterIndex1, daughterIndex2, Phi01, D01, Omega1, Z01, TanLambda1, Phi02, D02,
656 B2ERROR(
"First track provided has curvature zero. Calculation of convertedPhoton3Momentum failed.");
657 return B2Vector3D(std::numeric_limits<double>::quiet_NaN(), std::numeric_limits<double>::quiet_NaN(),
658 std::numeric_limits<double>::quiet_NaN());
661 B2ERROR(
"Second track provided has curvature zero. Calculation of convertedPhoton3Momentum failed.");
662 return B2Vector3D(std::numeric_limits<double>::quiet_NaN(), std::numeric_limits<double>::quiet_NaN(),
663 std::numeric_limits<double>::quiet_NaN());
668 double radius1 = 1 / Omega1;
669 double radius2 = 1 / Omega2;
671 TVector2 center1((radius1 + D01) * sin(Phi01), -1 * (radius1 + D01) * cos(Phi01));
672 TVector2 center2((radius2 + D02) * sin(Phi02), -1 * (radius2 + D02) * cos(Phi02));
673 TVector2 n1 = center1 - center2; n1 = n1.Unit();
674 TVector2 n2 = -1 * n1;
675 n1 = copysign(1.0, Omega1) * n1;
676 n2 = copysign(1.0, Omega2) * n2;
679 double phiN1 = atan2(n1.X(), -n1.Y());
680 double phiN2 = atan2(n2.X(), -n2.Y());
683 double sinlam1 = TanLambda1 / sqrt(1 + (TanLambda1 * TanLambda1));
684 double coslam1 = 1 / sqrt(1 + (TanLambda1 * TanLambda1));
685 double sinlam2 = TanLambda2 / sqrt(1 + (TanLambda2 * TanLambda2));
686 double coslam2 = 1 / sqrt(1 + (TanLambda2 * TanLambda2));
689 double p1 = gamma->getDaughter(daughterIndex1)->getMomentumMagnitude();
690 B2Vector3D e1Momentum(coslam1 * cos(phiN1), coslam1 * sin(phiN1), sinlam1);
691 double p2 = gamma->getDaughter(daughterIndex2)->getMomentumMagnitude();
692 B2Vector3D e2Momentum(coslam2 * cos(phiN2), coslam2 * sin(phiN2), sinlam2);
693 B2Vector3D gammaMomentum = (e1Momentum * p1) + (e2Momentum * p2);
695 return gammaMomentum;
698 double convertedPhotonPx(
const Particle* gamma,
const std::vector<double>& daughterIndices)
700 auto gammaMomentum = convertedPhoton3Momentum(gamma, daughterIndices);
701 return gammaMomentum.Px();
704 double convertedPhotonPy(
const Particle* gamma,
const std::vector<double>& daughterIndices)
706 auto gammaMomentum = convertedPhoton3Momentum(gamma, daughterIndices);
707 return gammaMomentum.Py();
710 double convertedPhotonPz(
const Particle* gamma,
const std::vector<double>& daughterIndices)
712 auto gammaMomentum = convertedPhoton3Momentum(gamma, daughterIndices);
713 return gammaMomentum.Pz();
716 int v0DaughtersShareInnermostHit(
const Particle* part)
720 auto daughterPlus = part->getDaughter(0);
721 auto daughterMinus = part->getDaughter(1);
722 if (!daughterPlus || !daughterMinus)
724 auto trackFitPlus = daughterPlus->getTrackFitResult();
725 auto trackFitMinus = daughterMinus->getTrackFitResult();
726 if (!trackFitPlus || !trackFitMinus)
728 int flagPlus = trackFitPlus->getHitPatternVXD().getInnermostHitShareStatus();
729 int flagMinus = trackFitMinus->getHitPatternVXD().getInnermostHitShareStatus();
730 if (flagPlus != flagMinus)
735 bool v0DaughtersShareInnermostUHit(
const Particle* part)
737 return ((v0DaughtersShareInnermostHit(part) / 2) == 1);
740 bool v0DaughtersShareInnermostVHit(
const Particle* part)
742 return ((v0DaughtersShareInnermostHit(part) % 2) == 1);
745 VARIABLE_GROUP(
"V0Daughter");
747 REGISTER_VARIABLE(
"v0DaughterNCDCHits(i)", v0DaughterTrackNCDCHits,
"Number of CDC hits associated to the i-th daughter track");
748 MAKE_DEPRECATED(
"v0DaughterNCDCHits",
false,
"light-2104-poseidon", R
"DOC(
749 The same value can be calculated with the more generic variable `nCDCHits`,
750 so replace the current call with ``daughter(i, nCDCHits)``.)DOC");
751 REGISTER_VARIABLE("v0DaughterNSVDHits(i)", v0DaughterTrackNSVDHits,
"Number of SVD hits associated to the i-th daughter track");
752 MAKE_DEPRECATED(
"v0DaughterNSVDHits",
false,
"light-2104-poseidon", R
"DOC(
753 The same value can be calculated with the more generic variable `nSVDHits`,
754 so replace the current call with ``daughter(i, nSVDHits)``.)DOC");
755 REGISTER_VARIABLE("v0DaughterNPXDHits(i)", v0DaughterTrackNPXDHits,
"Number of PXD hits associated to the i-th daughter track");
756 MAKE_DEPRECATED(
"v0DaughterNPXDHits",
false,
"light-2104-poseidon", R
"DOC(
757 The same value can be calculated with the more generic variable `nPXDHits`,
758 so replace the current call with ``daughter(i, nPXDHits)``.)DOC");
759 REGISTER_VARIABLE("v0DaughterNVXDHits(i)", v0DaughterTrackNVXDHits,
"Number of PXD+SVD hits associated to the i-th daughter track");
760 MAKE_DEPRECATED(
"v0DaughterNVXDHits",
false,
"light-2104-poseidon", R
"DOC(
761 The same value can be calculated with the more generic variable `nVXDHits`,
762 so replace the current call with ``daughter(i, nVXDHits)``.)DOC");
763 REGISTER_VARIABLE("v0DaughterNRemovedHits(i)", v0DaughterTrackNRemovedHits,
764 "The number of the i-th daughter track hits removed in V0Finder. Returns 0 if called for something other than V0 daughters.");
765 REGISTER_VARIABLE(
"v0DaughterFirstSVDLayer(i)", v0DaughterTrackFirstSVDLayer,
766 "First activated SVD layer associated to the i-th daughter track");
767 MAKE_DEPRECATED(
"v0DaughterFirstSVDLayer",
false,
"light-2104-poseidon", R
"DOC(
768 The same value can be calculated with the more generic variable `firstSVDLayer`,
769 so replace the current call with ``daughter(i, firstSVDLayer)``.)DOC");
770 REGISTER_VARIABLE("v0DaughterFirstPXDLayer(i)", v0DaughterTrackFirstPXDLayer,
771 "First activated PXD layer associated to the i-th daughter track");
772 MAKE_DEPRECATED(
"v0DaughterFirstPXDLayer",
false,
"light-2104-poseidon", R
"DOC(
773 The same value can be calculated with the more generic variable `firstPXDLayer`,
774 so replace the current call with ``daughter(i, firstPXDLayer)``.)DOC");
775 REGISTER_VARIABLE("v0DaughterFirstCDCLayer(i)", v0DaughterTrackFirstCDCLayer,
776 "First activated CDC layer associated to the i-th daughter track");
777 MAKE_DEPRECATED(
"v0DaughterFirstCDCLayer",
false,
"light-2104-poseidon", R
"DOC(
778 The same value can be calculated with the more generic variable `firstCDCLayer`,
779 so replace the current call with ``daughter(i, firstCDCLayer)``.)DOC");
780 REGISTER_VARIABLE("v0DaughterLastCDCLayer(i)", v0DaughterTrackLastCDCLayer,
781 "Last CDC layer associated to the i-th daughter track");
782 MAKE_DEPRECATED(
"v0DaughterLastCDCLayer",
false,
"light-2104-poseidon", R
"DOC(
783 The same value can be calculated with the more generic variable `lastCDCLayer`,
784 so replace the current call with ``daughter(i, lastCDCLayer)``.)DOC");
785 REGISTER_VARIABLE("v0DaughterPValue(i)", v0DaughterTrackPValue,
786 "chi2 probalility of the i-th daughter track fit");
787 MAKE_DEPRECATED(
"v0DaughterPValue",
false,
"light-2104-poseidon", R
"DOC(
788 The same value can be calculated with the more generic variable `pValue`,
789 so replace the current call with ``daughter(i, pValue)``.)DOC");
791 REGISTER_VARIABLE("v0DaughterD0(i)", v0DaughterTrackD0,
"d0 of the i-th daughter track fit",
"cm");
793 The same value can be calculated with the more generic variable `d0`,
794 so replace the current call with ``daughter(i, d0)``.)DOC");
795 REGISTER_VARIABLE("v0DaughterPhi0(i)", v0DaughterTrackPhi0,
"phi0 of the i-th daughter track fit",
"rad");
797 The same value can be calculated with the more generic variable `phi0`,
798 so replace the current call with ``daughter(i, phi0)``.)DOC");
799 REGISTER_VARIABLE("v0DaughterOmega(i)", v0DaughterTrackOmega,
"omega of the i-th daughter track fit",
800 ":math:`\\text{cm}^{-1}`");
801 MAKE_DEPRECATED(
"v0DaughterOmega",
false,
"light-2104-poseidon", R
"DOC(
802 The same value can be calculated with the more generic variable `omega`,
803 so replace the current call with ``daughter(i, omega)``.)DOC");
804 REGISTER_VARIABLE("v0DaughterZ0(i)", v0DaughterTrackZ0,
"z0 of the i-th daughter track fit",
"cm");
806 The same value can be calculated with the more generic variable `z0`,
807 so replace the current call with ``daughter(i, z0)``.)DOC");
808 REGISTER_VARIABLE("v0DaughterTanLambda(i)", v0DaughterTrackTanLambda,
"tan(lambda) of the i-th daughter track fit");
809 MAKE_DEPRECATED(
"v0DaughterTanLambda",
false,
"light-2104-poseidon", R
"DOC(
810 The same value can be calculated with the more generic variable `tanLambda`,
811 so replace the current call with ``daughter(i, tanLambda)``.)DOC");
813 REGISTER_VARIABLE("v0DaughterD0Error(i)", v0DaughterTrackD0Error,
"d0 error of the i-th daughter track fit",
"cm");
814 MAKE_DEPRECATED(
"v0DaughterD0Error",
false,
"light-2104-poseidon", R
"DOC(
815 The same value can be calculated with the more generic variable `d0Err`,
816 so replace the current call with ``daughter(i, d0Err)``.)DOC");
817 REGISTER_VARIABLE("v0DaughterPhi0Error(i)", v0DaughterTrackPhi0Error,
"phi0 error of the i-th daughter track fit",
"rad");
818 MAKE_DEPRECATED(
"v0DaughterPhi0Error",
false,
"light-2104-poseidon", R
"DOC(
819 The same value can be calculated with the more generic variable `phi0Err`,
820 so replace the current call with ``daughter(i, phi0Err)``.)DOC");
821 REGISTER_VARIABLE("v0DaughterOmegaError(i)", v0DaughterTrackOmegaError,
"omega error of the i-th daughter track fit",
822 ":math:`\\text{cm}^{-1}`");
823 MAKE_DEPRECATED(
"v0DaughterOmegaError",
false,
"light-2104-poseidon", R
"DOC(
824 The same value can be calculated with the more generic variable `omegaErr`,
825 so replace the current call with ``daughter(i, omegaErr)``.)DOC");
826 REGISTER_VARIABLE("v0DaughterZ0Error(i)", v0DaughterTrackZ0Error,
"z0 error of the i-th daughter track fit",
"cm");
827 MAKE_DEPRECATED(
"v0DaughterZ0Error",
false,
"light-2104-poseidon", R
"DOC(
828 The same value can be calculated with the more generic variable `z0Err`,
829 so replace the current call with ``daughter(i, z0Err)``.)DOC");
830 REGISTER_VARIABLE("v0DaughterTanLambdaError(i)", v0DaughterTrackTanLambdaError,
"tan(lambda) error of the i-th daughter track fit");
831 MAKE_DEPRECATED(
"v0DaughterTanLambdaError",
false,
"light-2104-poseidon", R
"DOC(
832 The same value can be calculated with the more generic variable `tanLambdaErr`,
833 so replace the current call with ``daughter(i, tanLambdaErr)``.)DOC");
836 REGISTER_VARIABLE("V0d0(id)", v0DaughterD0,
837 "Return the d0 impact parameter of a V0's daughter with daughterID index with the V0 vertex point as a pivot for the track.",
"cm");
838 REGISTER_VARIABLE(
"V0Deltad0", v0DaughterD0Diff,
839 "Return the difference between d0 impact parameters of V0's daughters with the V0 vertex point as a pivot for the track.",
"cm");
840 REGISTER_VARIABLE(
"V0z0(id)", v0DaughterZ0,
841 "Return the z0 impact parameter of a V0's daughter with daughterID index with the V0 vertex point as a pivot for the track.",
"cm");
842 REGISTER_VARIABLE(
"V0Deltaz0", v0DaughterZ0Diff,
843 "Return the difference between z0 impact parameters of V0's daughters with the V0 vertex point as a pivot for the track.",
"cm");
846 REGISTER_VARIABLE(
"v0DaughterD0PullWithTrueVertexAsPivot(i)", v0DaughterHelixWithTrueVertexAsPivotD0Pull,
847 "d0 pull of the i-th daughter track with the true V0 vertex as the track pivot");
848 REGISTER_VARIABLE(
"v0DaughterPhi0PullWithTrueVertexAsPivot(i)", v0DaughterHelixWithTrueVertexAsPivotPhi0Pull,
849 "phi0 pull of the i-th daughter track with the true V0 vertex as the track pivot");
850 REGISTER_VARIABLE(
"v0DaughterOmegaPullWithTrueVertexAsPivot(i)", v0DaughterHelixWithTrueVertexAsPivotOmegaPull,
851 "omega pull of the i-th daughter track with the true V0 vertex as the track pivot");
852 REGISTER_VARIABLE(
"v0DaughterZ0PullWithTrueVertexAsPivot(i)", v0DaughterHelixWithTrueVertexAsPivotZ0Pull,
853 "z0 pull of the i-th daughter track with the true V0 vertex as the track pivot");
854 REGISTER_VARIABLE(
"v0DaughterTanLambdaPullWithTrueVertexAsPivot(i)", v0DaughterHelixWithTrueVertexAsPivotTanLambdaPull,
855 "tan(lambda) pull of the i-th daughter track with the true V0 vertex as the track pivot");
857 REGISTER_VARIABLE(
"v0DaughterD0PullWithOriginAsPivot(i)", v0DaughterHelixWithOriginAsPivotD0Pull,
858 "d0 pull of the i-th daughter track with the origin as the track pivot");
859 MAKE_DEPRECATED(
"v0DaughterD0PullWithOriginAsPivot",
false,
"light-2104-poseidon", R
"DOC(
860 The same value can be calculated with the more generic variable `d0Pull`,
861 so replace the current call with ``daughter(i, d0Pull)``.)DOC");
862 REGISTER_VARIABLE("v0DaughterPhi0PullWithOriginAsPivot(i)", v0DaughterHelixWithOriginAsPivotPhi0Pull,
863 "phi0 pull of the i-th daughter track with the origin as the track pivot");
864 MAKE_DEPRECATED(
"v0DaughterPhi0PullWithOriginAsPivot",
false,
"light-2104-poseidon", R
"DOC(
865 The same value can be calculated with the more generic variable `phi0Pull`,
866 so replace the current call with ``daughter(i, phi0Pull)``.)DOC");
867 REGISTER_VARIABLE("v0DaughterOmegaPullWithOriginAsPivot(i)", v0DaughterHelixWithOriginAsPivotOmegaPull,
868 "omega pull of the i-th daughter track with the origin as the track pivot");
869 MAKE_DEPRECATED(
"v0DaughterOmegaPullWithOriginAsPivot",
false,
"light-2104-poseidon", R
"DOC(
870 The same value can be calculated with the more generic variable `omegaPull`,
871 so replace the current call with ``daughter(i, omegaPull)``.)DOC");
872 REGISTER_VARIABLE("v0DaughterZ0PullWithOriginAsPivot(i)", v0DaughterHelixWithOriginAsPivotZ0Pull,
873 "z0 pull of the i-th daughter track with the origin as the track pivot");
874 MAKE_DEPRECATED(
"v0DaughterZ0PullWithOriginAsPivot",
false,
"light-2104-poseidon", R
"DOC(
875 The same value can be calculated with the more generic variable `z0Pull`,
876 so replace the current call with ``daughter(i, z0Pull)``.)DOC");
877 REGISTER_VARIABLE("v0DaughterTanLambdaPullWithOriginAsPivot(i)", v0DaughterHelixWithOriginAsPivotTanLambdaPull,
878 "tan(lambda) pull of the i-th daughter track with the origin as the track pivot");
879 MAKE_DEPRECATED(
"v0DaughterTanLambdaPullWithOriginAsPivot",
false,
"light-2104-poseidon", R
"DOC(
880 The same value can be calculated with the more generic variable `tanLambdaPull`,
881 so replace the current call with ``daughter(i, tanLambdaPull)``.)DOC");
883 REGISTER_VARIABLE("v0DaughterTau(i,j)", v0DaughterTrackParam5AtIPPerigee,
884 "j-th track parameter (at IP perigee) of the i-th daughter track. "
885 "j: 0:d0, 1:phi0, 2:omega, 3:z0, 4:tanLambda",
"cm, rad, :math:`\\text{cm}^{-1}`, cm, unitless");
886 REGISTER_VARIABLE(
"v0DaughterCov(i,j)", v0DaughterTrackParamCov5x5AtIPPerigee,
887 "j-th element of the 15 covariance matrix elements (at IP perigee) of the i-th daughter track. "
888 "(0,0), (0,1) ... (1,1), (1,2) ... (2,2) ..."
889 "index order is: 0:d0, 1:phi0, 2:omega, 3:z0, 4:tanLambda",
"cm, rad, :math:`\\text{cm}^{-1}`, cm, unitless");
891 REGISTER_VARIABLE(
"convertedPhotonInvariantMass(i,j)", convertedPhotonInvariantMass,
892 "Invariant mass of the i-j daughter system as defined in https://indico.belle2.org/event/3644/contributions/18622/attachments/9401/14443/Photon_vertexin_B2GM.pdf, assuming it's a converted photon",
893 "GeV/:math:`\\text{c}^2`");
894 REGISTER_VARIABLE(
"convertedPhotonDelTanLambda(i,j)", convertedPhotonDelTanLambda,
895 "Discriminating variable Delta-TanLambda calculated for daughters (i,j) as defined in https://indico.belle2.org/event/3644/contributions/18622/attachments/9401/14443/Photon_vertexin_B2GM.pdf, assuming it's a converted photon");
896 REGISTER_VARIABLE(
"convertedPhotonDelR(i,j)", convertedPhotonDelR,
897 "Discriminating variable Delta-R calculated for daughters (i,j) as defined in https://indico.belle2.org/event/3644/contributions/18622/attachments/9401/14443/Photon_vertexin_B2GM.pdf, assuming it's a converted photon",
899 REGISTER_VARIABLE(
"convertedPhotonDelZ(i,j)", convertedPhotonDelZ,
900 "Discriminating variable Delta-Z calculated for daughters (i,j) as defined in https://indico.belle2.org/event/3644/contributions/18622/attachments/9401/14443/Photon_vertexin_B2GM.pdf, assuming it's a converted photon",
902 REGISTER_VARIABLE(
"convertedPhotonX(i,j)", convertedPhotonX,
903 "Estimate of vertex X coordinate calculated for daughters (i,j) as defined in https://indico.belle2.org/event/3644/contributions/18622/attachments/9401/14443/Photon_vertexin_B2GM.pdf, assuming it's a converted photon",
905 REGISTER_VARIABLE(
"convertedPhotonY(i,j)", convertedPhotonY,
906 "Estimate of vertex Y coordinate calculated for daughters (i,j) as defined in https://indico.belle2.org/event/3644/contributions/18622/attachments/9401/14443/Photon_vertexin_B2GM.pdf, assuming it's a converted photon",
908 REGISTER_VARIABLE(
"convertedPhotonZ(i,j)", convertedPhotonZ,
909 "Estimate of vertex Z coordinate calculated for daughters (i,j) as defined in https://indico.belle2.org/event/3644/contributions/18622/attachments/9401/14443/Photon_vertexin_B2GM.pdf, assuming it's a converted photon",
911 REGISTER_VARIABLE(
"convertedPhotonRho(i,j)", convertedPhotonRho,
912 "Estimate of vertex Rho calculated for daughters (i,j) as defined in https://indico.belle2.org/event/3644/contributions/18622/attachments/9401/14443/Photon_vertexin_B2GM.pdf, assuming it's a converted photon",
914 REGISTER_VARIABLE(
"convertedPhotonPx(i,j)", convertedPhotonPx,
915 "Estimate of x-component of photon momentum calculated for daughters (i,j) as defined in https://indico.belle2.org/event/3644/contributions/18622/attachments/9401/14443/Photon_vertexin_B2GM.pdf, assuming it's a converted photon",
917 REGISTER_VARIABLE(
"convertedPhotonPy(i,j)", convertedPhotonPy,
918 "Estimate of y-component of photon momentum calculated for daughters (i,j) as defined in https://indico.belle2.org/event/3644/contributions/18622/attachments/9401/14443/Photon_vertexin_B2GM.pdf, assuming it's a converted photon",
920 REGISTER_VARIABLE(
"convertedPhotonPz(i,j)", convertedPhotonPz,
921 "Estimate of z-component of photon momentum calculated for daughters (i,j) as defined in https://indico.belle2.org/event/3644/contributions/18622/attachments/9401/14443/Photon_vertexin_B2GM.pdf, assuming it's a converted photon",
924 REGISTER_VARIABLE(
"v0DaughtersShare1stHit", v0DaughtersShareInnermostHit,
925 "flag for V0 daughters sharing the first(innermost) VXD hit. 0x1(0x2) bit represents V/z(U/r-phi)-hit share.");
926 REGISTER_VARIABLE(
"v0DaughtersShare1stUHit", v0DaughtersShareInnermostUHit,
927 "flag for V0 daughters sharing the first(innermost) VXD U-side hit.");
928 REGISTER_VARIABLE(
"v0DaughtersShare1stVHit", v0DaughtersShareInnermostVHit,
929 "flag for V0 daughters sharing the first(innermost) VXD V-side hit.");
const MCParticle * getMCParticle() const
Returns the pointer to the MCParticle object that was used to create this Particle (ParticleType == c...
const Particle * getDaughter(unsigned i) const
Returns a pointer to the i-th daughter particle.
#define MAKE_DEPRECATED(name, make_fatal, version, description)
Registers a variable as deprecated.
B2Vector3< double > B2Vector3D
typedef for common usage with double
Abstract base class for different kinds of events.