12 #include <analysis/variables/V0DaughterTrackVariables.h>
13 #include <analysis/VariableManager/Manager.h>
15 #include <analysis/variables/TrackVariables.h>
18 #include <framework/datastore/StoreObjPtr.h>
19 #include <framework/datastore/StoreArray.h>
20 #include <framework/dataobjects/Helix.h>
23 #include <mdst/dataobjects/Track.h>
24 #include <mdst/dataobjects/MCParticle.h>
25 #include <mdst/dataobjects/TrackFitResult.h>
26 #include <mdst/dataobjects/V0.h>
27 #include <mdst/dataobjects/EventLevelTrackingInfo.h>
28 #include <mdst/dataobjects/HitPatternCDC.h>
29 #include <mdst/dataobjects/HitPatternVXD.h>
30 #include <mdst/dataobjects/ECLCluster.h>
33 #include <framework/logging/Logger.h>
50 const TrackFitResult* getTrackFitResultFromV0DaughterParticle(
const Particle* particle,
const double daughterID)
52 const int dID = int(std::lround(daughterID));
53 if (not(dID == 0 || dID == 1))
56 const Particle* daughter = particle->getDaughter(dID);
60 const TrackFitResult* trackFit = daughter->getRelated<TrackFitResult>();
67 double getV0DaughterTrackDetNHits(
const Particle* particle,
const double daughterID,
const Const::EDetector& det)
69 auto trackFit = getTrackFitResultFromV0DaughterParticle(particle, daughterID);
71 return std::numeric_limits<double>::quiet_NaN();
74 if (det == Const::EDetector::CDC) {
75 return trackFit->getHitPatternCDC().getNHits();
76 }
else if (det == Const::EDetector::SVD) {
77 return trackFit->getHitPatternVXD().getNSVDHits();
78 }
else if (det == Const::EDetector::PXD) {
79 return trackFit->getHitPatternVXD().getNPXDHits();
81 return std::numeric_limits<double>::quiet_NaN();
85 double v0DaughterTrackNCDCHits(
const Particle* part,
const std::vector<double>& daughterID)
87 return getV0DaughterTrackDetNHits(part, daughterID[0], Const::EDetector::CDC);
90 double v0DaughterTrackNSVDHits(
const Particle* part,
const std::vector<double>& daughterID)
92 return getV0DaughterTrackDetNHits(part, daughterID[0], Const::EDetector::SVD);
95 double v0DaughterTrackNPXDHits(
const Particle* part,
const std::vector<double>& daughterID)
97 return getV0DaughterTrackDetNHits(part, daughterID[0], Const::EDetector::PXD);
100 double v0DaughterTrackNVXDHits(
const Particle* part,
const std::vector<double>& daughterID)
102 return v0DaughterTrackNPXDHits(part, daughterID) + v0DaughterTrackNSVDHits(part, daughterID);
105 double v0DaughterTrackFirstSVDLayer(
const Particle* part,
const std::vector<double>& daughterID)
107 auto trackFit = getTrackFitResultFromV0DaughterParticle(part, daughterID[0]);
109 return std::numeric_limits<double>::quiet_NaN();
111 return trackFit->getHitPatternVXD().getFirstSVDLayer();
114 double v0DaughterTrackFirstPXDLayer(
const Particle* part,
const std::vector<double>& daughterID)
116 auto trackFit = getTrackFitResultFromV0DaughterParticle(part, daughterID[0]);
118 return std::numeric_limits<double>::quiet_NaN();
120 return trackFit->getHitPatternVXD().getFirstPXDLayer(HitPatternVXD::PXDMode::normal);
123 double v0DaughterTrackFirstCDCLayer(
const Particle* part,
const std::vector<double>& daughterID)
125 auto trackFit = getTrackFitResultFromV0DaughterParticle(part, daughterID[0]);
127 return std::numeric_limits<double>::quiet_NaN();
129 return trackFit->getHitPatternCDC().getFirstLayer();
132 double v0DaughterTrackLastCDCLayer(
const Particle* part,
const std::vector<double>& daughterID)
134 auto trackFit = getTrackFitResultFromV0DaughterParticle(part, daughterID[0]);
136 return std::numeric_limits<double>::quiet_NaN();
138 return trackFit->getHitPatternCDC().getLastLayer();
141 double v0DaughterTrackPValue(
const Particle* part,
const std::vector<double>& daughterID)
143 auto trackFit = getTrackFitResultFromV0DaughterParticle(part, daughterID[0]);
145 return std::numeric_limits<float>::quiet_NaN();
147 return trackFit->getPValue();
152 double getV0DaughterTrackParamAtIndex(
const Particle* particle,
const double daughterID,
const int tauIndex)
154 auto trackFit = getTrackFitResultFromV0DaughterParticle(particle, daughterID);
156 return std::numeric_limits<double>::quiet_NaN();
158 return trackFit->getTau().at(tauIndex);
161 double v0DaughterTrackD0(
const Particle* part,
const std::vector<double>& daughterID)
163 return getV0DaughterTrackParamAtIndex(part, daughterID[0], 0);
166 double v0DaughterTrackPhi0(
const Particle* part,
const std::vector<double>& daughterID)
168 return getV0DaughterTrackParamAtIndex(part, daughterID[0], 1);
171 double v0DaughterTrackOmega(
const Particle* part,
const std::vector<double>& daughterID)
173 return getV0DaughterTrackParamAtIndex(part, daughterID[0], 2);
176 double v0DaughterTrackZ0(
const Particle* part,
const std::vector<double>& daughterID)
178 return getV0DaughterTrackParamAtIndex(part, daughterID[0], 3);
181 double v0DaughterTrackTanLambda(
const Particle* part,
const std::vector<double>& daughterID)
183 return getV0DaughterTrackParamAtIndex(part, daughterID[0], 4);
188 double getV0DaughterTrackParamErrorAtIndex(
const Particle* particle,
const double daughterID,
const int tauIndex)
190 auto trackFit = getTrackFitResultFromV0DaughterParticle(particle, daughterID);
192 return std::numeric_limits<double>::quiet_NaN();
194 double errorSquared = trackFit->getCovariance5()[tauIndex][tauIndex];
195 if (errorSquared > 0)
196 return sqrt(errorSquared);
198 return std::numeric_limits<double>::quiet_NaN();
201 double v0DaughterTrackD0Error(
const Particle* part,
const std::vector<double>& daughterID)
203 return getV0DaughterTrackParamErrorAtIndex(part, daughterID[0], 0);
206 double v0DaughterTrackPhi0Error(
const Particle* part,
const std::vector<double>& daughterID)
208 return getV0DaughterTrackParamErrorAtIndex(part, daughterID[0], 1);
211 double v0DaughterTrackOmegaError(
const Particle* part,
const std::vector<double>& daughterID)
213 return getV0DaughterTrackParamErrorAtIndex(part, daughterID[0], 2);
216 double v0DaughterTrackZ0Error(
const Particle* part,
const std::vector<double>& daughterID)
218 return getV0DaughterTrackParamErrorAtIndex(part, daughterID[0], 3);
221 double v0DaughterTrackTanLambdaError(
const Particle* part,
const std::vector<double>& daughterID)
223 return getV0DaughterTrackParamErrorAtIndex(part, daughterID[0], 4);
228 double getHelixParameterPullOfV0DaughterWithTrueVertexAsPivotAtIndex(
const Particle* particle,
const double daughterID,
231 if (!particle) {
return std::numeric_limits<double>::quiet_NaN(); }
233 const int dID = int(std::lround(daughterID));
234 if (not(dID == 0 || dID == 1)) {
return std::numeric_limits<double>::quiet_NaN(); }
236 const MCParticle* mcparticle_v0 = particle->getMCParticle();
237 if (!mcparticle_v0) {
return std::numeric_limits<double>::quiet_NaN(); }
239 if (!(particle->getDaughter(dID))) {
return std::numeric_limits<double>::quiet_NaN(); }
241 const MCParticle* mcparticle = particle->getDaughter(dID)->getMCParticle();
242 if (!mcparticle) {
return std::numeric_limits<double>::quiet_NaN(); }
244 const TrackFitResult* trackFit = getTrackFitResultFromV0DaughterParticle(particle, daughterID);
245 if (!trackFit) {
return std::numeric_limits<double>::quiet_NaN(); }
248 const TVector3 mcProdVertex = mcparticle->getVertex();
249 const TVector3 mcMomentum = mcparticle->getMomentum();
250 const double mcParticleCharge = mcparticle->getCharge();
251 const double BzAtProdVertex = BFieldManager::getField(TVector3(mcProdVertex.X(), mcProdVertex.Y(),
253 Helix mcHelix = Helix(mcProdVertex, mcMomentum, mcParticleCharge, BzAtProdVertex);
254 mcHelix.passiveMoveBy(mcProdVertex);
255 const std::vector<double> mcHelixPars = { mcHelix.getD0(), mcHelix.getPhi0(), mcHelix.getOmega(),
256 mcHelix.getZ0(), mcHelix.getTanLambda()
260 UncertainHelix measHelix = trackFit->getUncertainHelix();
261 measHelix.passiveMoveBy(mcProdVertex);
262 const TMatrixDSym measCovariance = measHelix.getCovariance();
263 const std::vector<double> measHelixPars = {measHelix.getD0(), measHelix.getPhi0(), measHelix.getOmega(),
264 measHelix.getZ0(), measHelix.getTanLambda()
266 const std::vector<double> measErrSquare = {measCovariance[0][0], measCovariance[1][1], measCovariance[2][2],
267 measCovariance[3][3], measCovariance[4][4]
270 if (measErrSquare.at(tauIndex) > 0)
271 return (mcHelixPars.at(tauIndex) - measHelixPars.at(tauIndex)) / std::sqrt(measErrSquare.at(tauIndex));
273 return std::numeric_limits<double>::quiet_NaN();
276 double v0DaughterHelixWithTrueVertexAsPivotD0Pull(
const Particle* part,
const std::vector<double>& daughterID)
278 return getHelixParameterPullOfV0DaughterWithTrueVertexAsPivotAtIndex(part, daughterID[0], 0);
281 double v0DaughterHelixWithTrueVertexAsPivotPhi0Pull(
const Particle* part,
const std::vector<double>& daughterID)
283 return getHelixParameterPullOfV0DaughterWithTrueVertexAsPivotAtIndex(part, daughterID[0], 1);
286 double v0DaughterHelixWithTrueVertexAsPivotOmegaPull(
const Particle* part,
const std::vector<double>& daughterID)
288 return getHelixParameterPullOfV0DaughterWithTrueVertexAsPivotAtIndex(part, daughterID[0], 2);
291 double v0DaughterHelixWithTrueVertexAsPivotZ0Pull(
const Particle* part,
const std::vector<double>& daughterID)
293 return getHelixParameterPullOfV0DaughterWithTrueVertexAsPivotAtIndex(part, daughterID[0], 3);
296 double v0DaughterHelixWithTrueVertexAsPivotTanLambdaPull(
const Particle* part,
const std::vector<double>& daughterID)
298 return getHelixParameterPullOfV0DaughterWithTrueVertexAsPivotAtIndex(part, daughterID[0], 4);
303 double getHelixParameterPullOfV0DaughterWithOriginAsPivotAtIndex(
const Particle* particle,
const double daughterID,
306 if (!particle) {
return std::numeric_limits<double>::quiet_NaN(); }
308 const int dID = int(std::lround(daughterID));
309 if (not(dID == 0 || dID == 1)) {
return std::numeric_limits<double>::quiet_NaN(); }
311 const MCParticle* mcparticle_v0 = particle->getMCParticle();
312 if (!mcparticle_v0) {
return std::numeric_limits<double>::quiet_NaN(); }
314 if (!(particle->getDaughter(dID))) {
return std::numeric_limits<double>::quiet_NaN(); }
316 const MCParticle* mcparticle = particle->getDaughter(dID)->getMCParticle();
317 if (!mcparticle) {
return std::numeric_limits<double>::quiet_NaN(); }
319 const TrackFitResult* trackFit = getTrackFitResultFromV0DaughterParticle(particle, daughterID);
320 if (!trackFit) {
return std::numeric_limits<double>::quiet_NaN(); }
323 const TVector3 mcProdVertex = mcparticle->getVertex();
324 const TVector3 mcMomentum = mcparticle->getMomentum();
325 const double mcParticleCharge = mcparticle->getCharge();
326 const double BzAtOrigin = BFieldManager::getField(TVector3(0, 0, 0)).Z() /
Belle2::Unit::T;
328 Helix mcHelix = Helix(mcProdVertex, mcMomentum, mcParticleCharge, BzAtOrigin);
329 const std::vector<double> mcHelixPars = { mcHelix.getD0(), mcHelix.getPhi0(), mcHelix.getOmega(),
330 mcHelix.getZ0(), mcHelix.getTanLambda()
334 UncertainHelix measHelix = trackFit->getUncertainHelix();
335 const TMatrixDSym measCovariance = measHelix.getCovariance();
336 const std::vector<double> measHelixPars = {measHelix.getD0(), measHelix.getPhi0(), measHelix.getOmega(),
337 measHelix.getZ0(), measHelix.getTanLambda()
339 const std::vector<double> measErrSquare = {measCovariance[0][0], measCovariance[1][1], measCovariance[2][2],
340 measCovariance[3][3], measCovariance[4][4]
343 if (measErrSquare.at(tauIndex) > 0)
344 return (mcHelixPars.at(tauIndex) - measHelixPars.at(tauIndex)) / std::sqrt(measErrSquare.at(tauIndex));
346 return std::numeric_limits<double>::quiet_NaN();
349 double v0DaughterHelixWithOriginAsPivotD0Pull(
const Particle* part,
const std::vector<double>& daughterID)
351 return getHelixParameterPullOfV0DaughterWithOriginAsPivotAtIndex(part, daughterID[0], 0);
354 double v0DaughterHelixWithOriginAsPivotPhi0Pull(
const Particle* part,
const std::vector<double>& daughterID)
356 return getHelixParameterPullOfV0DaughterWithOriginAsPivotAtIndex(part, daughterID[0], 1);
359 double v0DaughterHelixWithOriginAsPivotOmegaPull(
const Particle* part,
const std::vector<double>& daughterID)
361 return getHelixParameterPullOfV0DaughterWithOriginAsPivotAtIndex(part, daughterID[0], 2);
364 double v0DaughterHelixWithOriginAsPivotZ0Pull(
const Particle* part,
const std::vector<double>& daughterID)
366 return getHelixParameterPullOfV0DaughterWithOriginAsPivotAtIndex(part, daughterID[0], 3);
369 double v0DaughterHelixWithOriginAsPivotTanLambdaPull(
const Particle* part,
const std::vector<double>& daughterID)
371 return getHelixParameterPullOfV0DaughterWithOriginAsPivotAtIndex(part, daughterID[0], 4);
374 double v0DaughterTrackParam5AtIPPerigee(
const Particle* part,
const std::vector<double>& params)
376 auto trackFit = getTrackFitResultFromV0DaughterParticle(part, params[0]);
378 return std::numeric_limits<double>::quiet_NaN();
381 const int paramID = int(std::lround(params[1]));
382 if (not(0 <= paramID && paramID < 5))
383 return std::numeric_limits<double>::quiet_NaN();
385 std::vector<float> tau = trackFit->getTau();
389 double v0DaughterTrackParamCov5x5AtIPPerigee(
const Particle* part,
const std::vector<double>& params)
391 auto trackFit = getTrackFitResultFromV0DaughterParticle(part, params[0]);
393 return std::numeric_limits<double>::quiet_NaN();
396 const int paramID = int(std::lround(params[1]));
397 if (not(0 <= paramID && paramID < 15))
398 return std::numeric_limits<double>::quiet_NaN();
400 std::vector<float> cov = trackFit->getCov();
405 VARIABLE_GROUP(
"V0Daughter");
407 REGISTER_VARIABLE(
"v0DaughterNCDCHits(i)", v0DaughterTrackNCDCHits,
"Number of CDC hits associated to the i-th daughter track");
408 REGISTER_VARIABLE(
"v0DaughterNSVDHits(i)", v0DaughterTrackNSVDHits,
"Number of SVD hits associated to the i-th daughter track");
409 REGISTER_VARIABLE(
"v0DaughterNPXDHits(i)", v0DaughterTrackNPXDHits,
"Number of PXD hits associated to the i-th daughter track");
410 REGISTER_VARIABLE(
"v0DaughterNVXDHits(i)", v0DaughterTrackNVXDHits,
"Number of PXD+SVD hits associated to the i-th daughter track");
411 REGISTER_VARIABLE(
"v0DaughterFirstSVDLayer(i)", v0DaughterTrackFirstSVDLayer,
412 "First activated SVD layer associated to the i-th daughter track");
413 REGISTER_VARIABLE(
"v0DaughterFirstPXDLayer(i)", v0DaughterTrackFirstPXDLayer,
414 "First activated PXD layer associated to the i-th daughter track");
415 REGISTER_VARIABLE(
"v0DaughterFirstCDCLayer(i)", v0DaughterTrackFirstCDCLayer,
416 "First activated CDC layer associated to the i-th daughter track");
417 REGISTER_VARIABLE(
"v0DaughterLastCDCLayer(i)", v0DaughterTrackLastCDCLayer,
418 "Last CDC layer associated to the i-th daughter track");
419 REGISTER_VARIABLE(
"v0DaughterPValue(i)", v0DaughterTrackPValue,
420 "chi2 probalility of the i-th daughter track fit");
422 REGISTER_VARIABLE(
"v0DaughterD0(i)", v0DaughterTrackD0,
"d0 of the i-th daughter track fit");
423 REGISTER_VARIABLE(
"v0DaughterPhi0(i)", v0DaughterTrackPhi0,
"phi0 of the i-th daughter track fit");
424 REGISTER_VARIABLE(
"v0DaughterOmega(i)", v0DaughterTrackOmega,
"omega of the i-th daughter track fit");
425 REGISTER_VARIABLE(
"v0DaughterZ0(i)", v0DaughterTrackZ0,
"z0 of the i-th daughter track fit");
426 REGISTER_VARIABLE(
"v0DaughterTanLambda(i)", v0DaughterTrackTanLambda,
"tan(lambda) of the i-th daughter track fit");
428 REGISTER_VARIABLE(
"v0DaughterD0Error(i)", v0DaughterTrackD0Error,
"d0 error of the i-th daughter track fit");
429 REGISTER_VARIABLE(
"v0DaughterPhi0Error(i)", v0DaughterTrackPhi0Error,
"phi0 error of the i-th daughter track fit");
430 REGISTER_VARIABLE(
"v0DaughterOmegaError(i)", v0DaughterTrackOmegaError,
"omega error of the i-th daughter track fit");
431 REGISTER_VARIABLE(
"v0DaughterZ0Error(i)", v0DaughterTrackZ0Error,
"z0 error of the i-th daughter track fit");
432 REGISTER_VARIABLE(
"v0DaughterTanLambdaError(i)", v0DaughterTrackTanLambdaError,
"tan(lambda) error of the i-th daughter track fit");
434 REGISTER_VARIABLE(
"v0DaughterD0PullWithTrueVertexAsPivot(i)", v0DaughterHelixWithTrueVertexAsPivotD0Pull,
435 "d0 pull of the i-th daughter track with the true V0 vertex as the track pivot");
436 REGISTER_VARIABLE(
"v0DaughterPhi0PullWithTrueVertexAsPivot(i)", v0DaughterHelixWithTrueVertexAsPivotPhi0Pull,
437 "phi0 pull of the i-th daughter track with the true V0 vertex as the track pivot");
438 REGISTER_VARIABLE(
"v0DaughterOmegaPullWithTrueVertexAsPivot(i)", v0DaughterHelixWithTrueVertexAsPivotOmegaPull,
439 "omega pull of the i-th daughter track with the true V0 vertex as the track pivot");
440 REGISTER_VARIABLE(
"v0DaughterZ0PullWithTrueVertexAsPivot(i)", v0DaughterHelixWithTrueVertexAsPivotZ0Pull,
441 "z0 pull of the i-th daughter track with the true V0 vertex as the track pivot");
442 REGISTER_VARIABLE(
"v0DaughterTanLambdaPullWithTrueVertexAsPivot(i)", v0DaughterHelixWithTrueVertexAsPivotTanLambdaPull,
443 "tan(lambda) pull of the i-th daughter track with the true V0 vertex as the track pivot");
445 REGISTER_VARIABLE(
"v0DaughterD0PullWithOriginAsPivot(i)", v0DaughterHelixWithOriginAsPivotD0Pull,
446 "d0 pull of the i-th daughter track with the origin as the track pivot");
447 REGISTER_VARIABLE(
"v0DaughterPhi0PullWithOriginAsPivot(i)", v0DaughterHelixWithOriginAsPivotPhi0Pull,
448 "phi0 pull of the i-th daughter track with the origin as the track pivot");
449 REGISTER_VARIABLE(
"v0DaughterOmegaPullWithOriginAsPivot(i)", v0DaughterHelixWithOriginAsPivotOmegaPull,
450 "omega pull of the i-th daughter track with the origin as the track pivot");
451 REGISTER_VARIABLE(
"v0DaughterZ0PullWithOriginAsPivot(i)", v0DaughterHelixWithOriginAsPivotZ0Pull,
452 "z0 pull of the i-th daughter track with the origin as the track pivot");
453 REGISTER_VARIABLE(
"v0DaughterTanLambdaPullWithOriginAsPivot(i)", v0DaughterHelixWithOriginAsPivotTanLambdaPull,
454 "tan(lambda) pull of the i-th daughter track with the origin as the track pivot");
456 REGISTER_VARIABLE(
"v0DaughterTau(i,j)", v0DaughterTrackParam5AtIPPerigee,
457 "j-th track parameter (at IP perigee) of the i-th daughter track. "
458 "j: 0:d0, 1:phi0, 2:omega, 3:z0, 4:tanLambda");
459 REGISTER_VARIABLE(
"v0DaughterCov(i,j)", v0DaughterTrackParamCov5x5AtIPPerigee,
460 "j-th element of the 15 covariance matrix elements (at IP perigee) of the i-th daughter track. "
461 "(0,0), (0,1) ... (1,1), (1,2) ... (2,2) ...");