Belle II Software development
TrackVariables.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/variables/TrackVariables.h>
11
12// include VariableManager
13#include <analysis/VariableManager/Manager.h>
14
15#include <analysis/dataobjects/Particle.h>
16#include <analysis/utility/DetectorSurface.h>
17
18// framework - DataStore
19#include <framework/datastore/StoreObjPtr.h>
20#include <framework/dataobjects/Helix.h>
21
22// dataobjects from the MDST
23#include <mdst/dbobjects/BeamSpot.h>
24#include <mdst/dataobjects/Track.h>
25#include <mdst/dataobjects/MCParticle.h>
26#include <mdst/dataobjects/TrackFitResult.h>
27#include <mdst/dataobjects/EventLevelTrackingInfo.h>
28#include <mdst/dataobjects/HitPatternVXD.h>
29#include <mdst/dataobjects/ECLCluster.h>
30
31// framework aux
32#include <framework/logging/Logger.h>
33
34#include <cmath>
35
36namespace Belle2 {
41 namespace Variable {
42
43 static const ROOT::Math::XYZVector vecNaN(Const::doubleNaN, Const::doubleNaN, Const::doubleNaN);
44
45 double trackNHits(const Particle* part, const Const::EDetector& det)
46 {
47 auto trackFit = part->getTrackFitResult();
48 if (!trackFit) return Const::doubleNaN;
49
50 if (det == Const::EDetector::CDC) {
51 return trackFit->getHitPatternCDC().getNHits();
52 } else if (det == Const::EDetector::SVD) {
53 return trackFit->getHitPatternVXD().getNSVDHits();
54 } else if (det == Const::EDetector::PXD) {
55 return trackFit->getHitPatternVXD().getNPXDHits();
56 } else {
57 return Const::doubleNaN;
58 }
59 }
60
61 double trackNCDCHits(const Particle* part)
62 {
63 return trackNHits(part, Const::EDetector::CDC);
64 }
65
66 double trackNSVDHits(const Particle* part)
67 {
68 return trackNHits(part, Const::EDetector::SVD);
69 }
70
71 double trackNPXDHits(const Particle* part)
72 {
73 return trackNHits(part, Const::EDetector::PXD);
74 }
75
76 double trackNVXDHits(const Particle* part)
77 {
78 return trackNPXDHits(part) + trackNSVDHits(part);
79 }
80
81 double trackNDF(const Particle* part)
82 {
83 auto trackFit = part->getTrackFitResult();
84 if (!trackFit) return Const::doubleNaN;
85 return trackFit->getNDF();
86 }
87
88 double trackChi2(const Particle* part)
89 {
90 auto trackFit = part->getTrackFitResult();
91 if (!trackFit) return Const::doubleNaN;
92 return trackFit->getChi2();
93 }
94
95 double trackFirstSVDLayer(const Particle* part)
96 {
97 auto trackFit = part->getTrackFitResult();
98 if (!trackFit) return Const::doubleNaN;
99 return trackFit->getHitPatternVXD().getFirstSVDLayer();
100 }
101
102 double trackFirstPXDLayer(const Particle* part)
103 {
104 auto trackFit = part->getTrackFitResult();
105 if (!trackFit) return Const::doubleNaN;
106 return trackFit->getHitPatternVXD().getFirstPXDLayer(HitPatternVXD::PXDMode::normal);
107 }
108
109 double trackFirstCDCLayer(const Particle* part)
110 {
111 auto trackFit = part->getTrackFitResult();
112 if (!trackFit) return Const::doubleNaN;
113 return trackFit->getHitPatternCDC().getFirstLayer();
114 }
115
116 double trackLastCDCLayer(const Particle* part)
117 {
118 auto trackFit = part->getTrackFitResult();
119 if (!trackFit) return Const::doubleNaN;
120 return trackFit->getHitPatternCDC().getLastLayer();
121 }
122
123 double trackD0(const Particle* part)
124 {
125 auto trackFit = part->getTrackFitResult();
126 if (!trackFit) return Const::doubleNaN;
127 return trackFit->getD0();
128 }
129
130 double trackPhi0(const Particle* part)
131 {
132 auto trackFit = part->getTrackFitResult();
133 if (!trackFit) return Const::doubleNaN;
134 return trackFit->getPhi0();
135 }
136
137 double trackOmega(const Particle* part)
138 {
139 auto trackFit = part->getTrackFitResult();
140 if (!trackFit) return Const::doubleNaN;
141 return trackFit->getOmega();
142 }
143
144 double trackZ0(const Particle* part)
145 {
146 auto trackFit = part->getTrackFitResult();
147 if (!trackFit) return Const::doubleNaN;
148 return trackFit->getZ0();
149 }
150
151 double trackTanLambda(const Particle* part)
152 {
153 auto trackFit = part->getTrackFitResult();
154 if (!trackFit) return Const::doubleNaN;
155 return trackFit->getTanLambda();
156 }
157
158 double trackD0FromIP(const Particle* part)
159 {
160 auto trackFit = part->getTrackFitResult();
161 if (!trackFit) return Const::doubleNaN;
162 static DBObjPtr<BeamSpot> beamSpotDB;
163 if (!beamSpotDB.isValid())
164 return Const::doubleNaN;
165 auto helix = trackFit->getHelix();
166 helix.passiveMoveBy(beamSpotDB->getIPPosition());
167 return helix.getD0();
168 }
169
170 double trackZ0FromIP(const Particle* part)
171 {
172 auto trackFit = part->getTrackFitResult();
173 if (!trackFit) return Const::doubleNaN;
174 static DBObjPtr<BeamSpot> beamSpotDB;
175 if (!beamSpotDB.isValid())
176 return Const::doubleNaN;
177 auto helix = trackFit->getHelix();
178 helix.passiveMoveBy(beamSpotDB->getIPPosition());
179 return helix.getZ0();
180 }
181
182 double trackPhi0FromIP(const Particle* part)
183 {
184 auto trackFit = part->getTrackFitResult();
185 if (!trackFit) return Const::doubleNaN;
186 static DBObjPtr<BeamSpot> beamSpotDB;
187 if (!beamSpotDB.isValid())
188 return Const::doubleNaN;
189 auto helix = trackFit->getHelix();
190 helix.passiveMoveBy(beamSpotDB->getIPPosition());
191 return helix.getPhi0();
192 }
193
194 double trackD0Error(const Particle* part)
195 {
196 auto trackFit = part->getTrackFitResult();
197 if (!trackFit) return Const::doubleNaN;
198 double errorSquared = trackFit->getCovariance5()[0][0];
199 if (errorSquared <= 0) return Const::doubleNaN;
200 return sqrt(errorSquared);
201 }
202
203 double trackPhi0Error(const Particle* part)
204 {
205 auto trackFit = part->getTrackFitResult();
206 if (!trackFit) return Const::doubleNaN;
207
208 double errorSquared = trackFit->getCovariance5()[1][1];
209 if (errorSquared <= 0) return Const::doubleNaN;
210 return sqrt(errorSquared);
211 }
212
213 double trackOmegaError(const Particle* part)
214 {
215 auto trackFit = part->getTrackFitResult();
216 if (!trackFit) return Const::doubleNaN;
217
218 double errorSquared = trackFit->getCovariance5()[2][2];
219 if (errorSquared <= 0) return Const::doubleNaN;
220 return sqrt(errorSquared);
221 }
222
223 double trackZ0Error(const Particle* part)
224 {
225 auto trackFit = part->getTrackFitResult();
226 if (!trackFit) return Const::doubleNaN;
227
228 double errorSquared = trackFit->getCovariance5()[3][3];
229 if (errorSquared <= 0) return Const::doubleNaN;
230 return sqrt(errorSquared);
231 }
232
233 double trackTanLambdaError(const Particle* part)
234 {
235 auto trackFit = part->getTrackFitResult();
236 if (!trackFit) return Const::doubleNaN;
237
238 double errorSquared = trackFit->getCovariance5()[4][4];
239 if (errorSquared <= 0) return Const::doubleNaN;
240 return sqrt(errorSquared);
241 }
242
243 double trackFitCovariance(const Particle* particle, const std::vector<double>& indices)
244 {
245 if (indices.size() != 2) {
246 B2FATAL("Exactly two indices must be provided to the variable trackFitCovariance!");
247 }
248 if (*(std::min_element(indices.begin(), indices.end())) < 0 or *(std::max_element(indices.begin(), indices.end())) > 4) {
249 B2FATAL("The indices provided to the variable trackFitCovariance must be in the range 0 - 4!");
250 }
251 auto trackFit = particle->getTrackFitResult();
252 if (!trackFit) return Const::doubleNaN;
253 return trackFit->getCovariance5()[indices[0]][indices[1]];
254 }
255
256 double trackPValue(const Particle* part)
257 {
258 auto trackFit = part->getTrackFitResult();
259 if (!trackFit) return Const::doubleNaN;
260 return trackFit->getPValue();
261 }
262
263 double trackFitHypothesisPDG(const Particle* part)
264 {
265 auto trackFit = part->getTrackFitResult();
266 if (!trackFit) return Const::doubleNaN;
267 return trackFit->getParticleType().getPDGCode();
268 }
269
270 double trackNECLClusters(const Particle* part)
271 {
272 const Track* track = part->getTrack();
273 if (!track) return Const::doubleNaN;
274
275 // count the number of nPhotons hypothesis ecl clusters
276 int count = 0;
277 for (const ECLCluster& cluster : track->getRelationsTo<ECLCluster>())
278 if (cluster.hasHypothesis(ECLCluster::EHypothesisBit::c_nPhotons))
279 ++count;
280 return count;
281 }
282
283 // used in trackHelixExtTheta and trackHelixExtPhi
284 ROOT::Math::XYZVector getPositionOnHelix(const TrackFitResult* trackFit, const std::vector<double>& pars)
285 {
286 const double r = pars[0];
287 const double zfwd = pars[1];
288 const double zbwd = pars[2];
289
290 // get helix and parameters
291 const double z0 = trackFit->getZ0();
292 const double tanlambda = trackFit->getTanLambda();
293 const Helix h = trackFit->getHelix();
294
295 // extrapolate to radius
296 const double arcLength = h.getArcLength2DAtCylindricalR(r);
297 const double lHelixRadius = arcLength > 0 ? arcLength : std::numeric_limits<double>::max();
298
299 // extrapolate to FWD z
300 const double lFWD = (zfwd - z0) / tanlambda > 0 ? (zfwd - z0) / tanlambda : std::numeric_limits<double>::max();
301
302 // extrapolate to BWD z
303 const double lBWD = (zbwd - z0) / tanlambda > 0 ? (zbwd - z0) / tanlambda : std::numeric_limits<double>::max();
304
305 // pick smallest arclength
306 const double l = std::min({lHelixRadius, lFWD, lBWD});
307
308 return h.getPositionAtArcLength2D(l);
309 }
310
311 ROOT::Math::XYZVector getPositionOnHelix(const Particle* part, const std::vector<double>& pars)
312 {
313 if (pars.size() == 4 and pars[3]) {
314 const Track* track = part->getTrack();
315 if (!track)
316 return vecNaN;
317
318 auto highestProbMass = part->getMostLikelyTrackFitResult().first;
319 const TrackFitResult* trackFit = track->getTrackFitResultWithClosestMass(highestProbMass);
320 if (!trackFit) return vecNaN;
321 else return getPositionOnHelix(trackFit, pars);
322 } else {
323 const TrackFitResult* trackFit = part->getTrackFitResult();
324 if (!trackFit) return vecNaN;
325 else return getPositionOnHelix(trackFit, pars);
326 }
327 }
328
329 // returns extrapolated theta position based on helix parameters
330 double trackHelixExtTheta(const Particle* part, const std::vector<double>& pars)
331 {
332 const auto nParams = pars.size();
333 if (nParams != 3 && nParams != 4) {
334 B2FATAL("Exactly three (+1 optional) parameters (r, zfwd, zbwd, [useHighestProbMass]) required for helixExtTheta.");
335 }
336
337 ROOT::Math::XYZVector position = getPositionOnHelix(part, pars);
338 if (position == vecNaN) return Const::doubleNaN;
339 return position.Theta();
340 }
341
342 // returns extrapolated phi position based on helix parameters
343 double trackHelixExtPhi(const Particle* part, const std::vector<double>& pars)
344 {
345 const auto nParams = pars.size();
346 if (nParams != 3 && nParams != 4) {
347 B2FATAL("Exactly three (+1 optional) parameters (r, zfwd, zbwd, [useHighestProbMass]) required for helixExtPhi.");
348 }
349
350 ROOT::Math::XYZVector position = getPositionOnHelix(part, pars);
351 if (position == vecNaN) return Const::doubleNaN;
352 return position.Phi();
353 }
354
355 Manager::FunctionPtr trackHelixExtThetaOnDet(const std::vector<std::string>& arguments)
356 {
357 if (arguments.size() != 1 && arguments.size() != 2)
358 B2FATAL("Exactly one (+1 optional) parameter (detector_surface_name, [useHighestProbMass]) is required for helixExtThetaOnDet.");
359
360 std::vector<double> parameters(3);
361 const std::string det = arguments[0];
363 parameters[0] = DetectorSurface::detToSurfBoundaries.at(det).m_rho;
364 parameters[1] = DetectorSurface::detToSurfBoundaries.at(det).m_zfwd;
365 parameters[2] = DetectorSurface::detToSurfBoundaries.at(det).m_zbwd;
367 parameters[0] = DetectorSurface::detLayerToSurfBoundaries.at(det).m_rho;
368 parameters[1] = DetectorSurface::detLayerToSurfBoundaries.at(det).m_zfwd;
369 parameters[2] = DetectorSurface::detLayerToSurfBoundaries.at(det).m_zbwd;
370 } else
371 B2FATAL("Given detector surface name is not supported.");
372
373 if (arguments.size() == 2)
374 parameters.push_back(std::stod(arguments[1]));
375
376 auto func = [parameters](const Particle * part) -> double {
377
378 ROOT::Math::XYZVector position = getPositionOnHelix(part, parameters);
379 if (position == vecNaN) return Const::doubleNaN;
380 return position.Theta();
381 };
382 return func;
383 }
384
385 Manager::FunctionPtr trackHelixExtPhiOnDet(const std::vector<std::string>& arguments)
386 {
387 if (arguments.size() != 1 && arguments.size() != 2)
388 B2FATAL("Exactly one (+1 optional) parameter (detector_surface_name, [useHighestProbMass]) is required for helixExtPhiOnDet.");
389
390 std::vector<double> parameters(3);
391 const std::string det = arguments[0];
393 parameters[0] = DetectorSurface::detToSurfBoundaries.at(det).m_rho;
394 parameters[1] = DetectorSurface::detToSurfBoundaries.at(det).m_zfwd;
395 parameters[2] = DetectorSurface::detToSurfBoundaries.at(det).m_zbwd;
397 parameters[0] = DetectorSurface::detLayerToSurfBoundaries.at(det).m_rho;
398 parameters[1] = DetectorSurface::detLayerToSurfBoundaries.at(det).m_zfwd;
399 parameters[2] = DetectorSurface::detLayerToSurfBoundaries.at(det).m_zbwd;
400 } else
401 B2FATAL("Given detector surface name is not supported.");
402
403 if (arguments.size() == 2)
404 parameters.push_back(std::stod(arguments[1]));
405
406 auto func = [parameters](const Particle * part) -> double {
407
408 ROOT::Math::XYZVector position = getPositionOnHelix(part, parameters);
409 if (position == vecNaN) return Const::doubleNaN;
410 return position.Phi();
411 };
412 return func;
413 }
414
415
416 /***************************************************
417 * Event level tracking quantities
418 */
419
420 // The number of CDC hits in the event not assigned to any track
421 double nExtraCDCHits(const Particle*)
422 {
423 StoreObjPtr<EventLevelTrackingInfo> elti;
424 if (!elti) return Const::doubleNaN;
425 return elti->getNCDCHitsNotAssigned();
426 }
427
428 // The number of CDC hits in the event not assigned to any track nor very
429 // likely beam background (i.e. hits that survive a cleanup selection)
430 double nExtraCDCHitsPostCleaning(const Particle*)
431 {
432 StoreObjPtr<EventLevelTrackingInfo> elti;
433 if (!elti) return Const::doubleNaN;
434 return elti->getNCDCHitsNotAssignedPostCleaning();
435 }
436
437 // Check for the presence of a non-assigned hit in the specified CDC layer
438 double hasExtraCDCHitsInLayer(const Particle*, const std::vector<double>& layer)
439 {
440 StoreObjPtr<EventLevelTrackingInfo> elti;
441 if (!elti) return Const::doubleNaN;
442 int ilayer = std::lround(layer[0]);
443 return elti->hasCDCLayer(ilayer);
444 }
445
446 // Check for the presence of a non-assigned hit in the specified CDC SuperLayer
447 double hasExtraCDCHitsInSuperLayer(const Particle*, const std::vector<double>& layer)
448 {
449 StoreObjPtr<EventLevelTrackingInfo> elti;
450 if (!elti) return Const::doubleNaN;
451 int ilayer = std::lround(layer[0]);
452 return elti->hasCDCSLayer(ilayer);
453 }
454
455 // The number of segments that couldn't be assigned to any track
456 double nExtraCDCSegments(const Particle*)
457 {
458 StoreObjPtr<EventLevelTrackingInfo> elti;
459 if (!elti) return Const::doubleNaN;
460 return elti->getNCDCSegments();
461 }
462
463 // The number of VXD hits not assigned to any track in the specified layer
464 double nExtraVXDHitsInLayer(const Particle*, const std::vector<double>& layer)
465 {
466 StoreObjPtr<EventLevelTrackingInfo> elti;
467 if (!elti) return Const::doubleNaN;
468 int ilayer = std::lround(layer[0]);
469 return elti->getNVXDClustersInLayer(ilayer);
470 }
471
472 // The number of VXD hits not assigned to any track
473 double nExtraVXDHits(const Particle*)
474 {
475 StoreObjPtr<EventLevelTrackingInfo> elti;
476 if (!elti) return Const::doubleNaN;
477 double out = 0.0;
478 for (uint16_t ilayer = 1; ilayer < 7; ++ilayer)
479 out += elti->getNVXDClustersInLayer(ilayer);
480 return out;
481 }
482
483 // time of first SVD sample relative to event T0
484 double svdFirstSampleTime(const Particle*)
485 {
486 StoreObjPtr<EventLevelTrackingInfo> elti;
487 if (!elti) return Const::doubleNaN;
488 return elti->getSVDFirstSampleTime();
489 }
490
491 // A flag set by the tracking if there is reason to assume there was a track
492 // in the event missed by the tracking or the track finding was (partly) aborted
493 // for this event. Further information about this can be obtained from the flagBlock
494 // of the EventLevelTrackingInfo object.
495 double trackFindingFailureFlag(const Particle*)
496 {
497 StoreObjPtr<EventLevelTrackingInfo> elti;
498 if (!elti) return Const::doubleNaN;
499 return elti->hasAnErrorFlag();
500 }
501
502 // Build a Helix from a MCParticle's kinematics. Does NOT check for null pointer!!!
503 Helix getMCHelix(const MCParticle* mcparticle)
504 {
505 const ROOT::Math::XYZVector mcProdVertex = mcparticle->getVertex();
506 const ROOT::Math::XYZVector mcMomentum = mcparticle->getMomentum();
507 const double BzAtProdVertex = BFieldManager::getFieldInTesla(mcProdVertex).Z();
508 const double mcParticleCharge = mcparticle->getCharge();
509 return Helix(mcProdVertex, mcMomentum, mcParticleCharge, BzAtProdVertex);
510 }
511
512 double getMCHelixParameterAtIndex(const Particle* particle, const int index)
513 {
514 if (!particle) return Const::doubleNaN;
515
516 const MCParticle* mcparticle = particle->getMCParticle();
517 if (!mcparticle) return Const::doubleNaN;
518
519 const Helix mcHelix(getMCHelix(mcparticle));
520 const std::vector<double> mcHelixPars{mcHelix.getD0(), mcHelix.getPhi0(), mcHelix.getOmega(), mcHelix.getZ0(), mcHelix.getTanLambda()};
521 return mcHelixPars.at(index);
522 }
523
524 double getHelixParameterPullAtIndex(const Particle* particle, const int index)
525 {
526 if (!particle) return Const::doubleNaN;
527
528 const MCParticle* mcparticle = particle->getMCParticle();
529 if (!mcparticle) return Const::doubleNaN;
530
531 const TrackFitResult* trackfit = particle->getTrackFitResult();
532 if (!trackfit) return Const::doubleNaN;
533
534 const UncertainHelix measHelix = trackfit->getUncertainHelix();
535 const TMatrixDSym& measCovariance = measHelix.getCovariance();
536 const Helix mcHelix(getMCHelix(mcparticle));
537
538 const std::vector<double> mcHelixPars = {mcHelix.getD0(), mcHelix.getPhi0(), mcHelix.getOmega(), mcHelix.getZ0(), mcHelix.getTanLambda()};
539 const std::vector<double> measHelixPars = {measHelix.getD0(), measHelix.getPhi0(), measHelix.getOmega(), measHelix.getZ0(), measHelix.getTanLambda()};
540
541 return (mcHelixPars.at(index) - measHelixPars.at(index)) / std::sqrt(measCovariance(index, index));
542 }
543
544 double getHelixMCD0(const Particle* part) { return getMCHelixParameterAtIndex(part, 0); }
545 double getHelixMCPhi0(const Particle* part) { return getMCHelixParameterAtIndex(part, 1); }
546 double getHelixMCOmega(const Particle* part) { return getMCHelixParameterAtIndex(part, 2); }
547 double getHelixMCZ0(const Particle* part) { return getMCHelixParameterAtIndex(part, 3); }
548 double getHelixMCTanLambda(const Particle* part) { return getMCHelixParameterAtIndex(part, 4); }
549
550 double getHelixD0Pull(const Particle* part)
551 {
552 return getHelixParameterPullAtIndex(part, 0);
553 }
554
555 double getHelixPhi0Pull(const Particle* part)
556 {
557 return getHelixParameterPullAtIndex(part, 1);
558 }
559
560 double getHelixOmegaPull(const Particle* part)
561 {
562 return getHelixParameterPullAtIndex(part, 2);
563 }
564
565 double getHelixZ0Pull(const Particle* part)
566 {
567 return getHelixParameterPullAtIndex(part, 3);
568 }
569
570 double getHelixTanLambdaPull(const Particle* part)
571 {
572 return getHelixParameterPullAtIndex(part, 4);
573 }
574
575 double getTrackTime(const Particle* part)
576 {
577 const Track* track = part->getTrack();
578 if (!track) return Const::doubleNaN;
579 return track->getTrackTime();
580 }
581
582 double isTrackFlippedAndRefitted(const Particle* part)
583 {
584 auto track = part->getTrack();
585 if (!track) return Const::doubleNaN;
586 return track->isFlippedAndRefitted() ? 1 : 0;
587 }
588
589 double getTrackLength(const Particle* part)
590 {
591 auto trackFit = part->getTrackFitResult();
592 if (!trackFit) return Const::doubleNaN;
593
594 const double lastCDCLayer = trackLastCDCLayer(part);
595 if (std::isnan(lastCDCLayer) or lastCDCLayer < 0)
596 return Const::doubleNaN;
597
598 const double r = DetectorSurface::cdcWireRadiuses.at((int)lastCDCLayer);
599
600 return trackFit->getHelix().getArcLength2DAtCylindricalR(r);
601 }
602
603
604 VARIABLE_GROUP("Tracking");
605
606 REGISTER_VARIABLE("mcD0", getHelixMCD0, R"DOC(
607Returns the MC value of :math:`d_0`, the signed distance to the
608point-of-closest-approach (POCA) in the :math:`r-\phi` plane.
609
610.. seealso:: :b2:var:`d0`
611
612Returns NaN if the particle is not related to any MCParticle.
613
614)DOC", "cm");
615 REGISTER_VARIABLE("mcPhi0", getHelixMCPhi0, R"DOC(
616Returns the MC value of :math:`\phi_0`, the angle of the transverse momentum
617in the :math:`r-\phi` plane.
618
619.. seealso:: :b2:var:`phi0`
620
621Returns NaN if the particle is not related to any MCParticle.
622
623)DOC", "rad");
624 REGISTER_VARIABLE("mcOmega", getHelixMCOmega, R"DOC(
625Returns the MC value of :math:`\omega`, the curvature of the track.
626
627.. seealso:: :b2:var:`omega`
628
629Returns NaN if the particle is not related to any MCParticle.
630
631)DOC", ":math:`\\text{cm}^{-1}`");
632 REGISTER_VARIABLE("mcZ0", getHelixMCZ0, R"DOC(
633Returns the MC value of :math:`z_0`, the z-coordinate of the
634point-of-closest-approach (POCA).
635
636.. seealso:: :b2:var:`z0`
637
638Returns NaN if the particle is not related to any MCParticle.
639
640)DOC", "cm");
641 REGISTER_VARIABLE("mcTanLambda", getHelixMCTanLambda, R"DOC(
642Returns the MC value of :math:`\tan\lambda`, the slope of the track in the
643:math:`r-z` plane.
644
645.. seealso:: :b2:var:`tanLambda`
646
647Returns NaN if the particle is not related to any MCParticle.
648)DOC");
649
650 REGISTER_VARIABLE("d0Pull", getHelixD0Pull, R"DOC(
651The pull of the tracking parameter :math:`d_0` for the reconstructed
652pattern-recognition track, with respect to the MC track. That is:
653
654.. math::
655
656 \frac{d_0^\textrm{MC} - d_0^\textrm{PR}}{\sigma_{d_0; \textrm{PR}}}
657
658.. seealso:: :b2:var:`d0`, :b2:var:`d0Err`
659
660Returns NaN if no MC particle is related or if called on something other than a
661track-based particle.
662 )DOC");
663 REGISTER_VARIABLE("phi0Pull", getHelixPhi0Pull, R"DOC(
664The pull of the tracking parameter :math:`\phi_0` for the reconstructed
665pattern-recognition track, with respect to the MC track. That is:
666
667.. math::
668
669 \frac{\phi_0^\textrm{MC} - \phi_0^\textrm{PR}}{\sigma_{\phi_0; \textrm{PR}}}
670
671.. seealso:: :b2:var:`phi0`, :b2:var:`phi0Err`
672
673Returns NaN if no MC particle is related or if called on something other than a
674track-based particle.
675 )DOC");
676 REGISTER_VARIABLE("omegaPull", getHelixOmegaPull, R"DOC(
677The pull of the tracking parameter :math:`\omega` for the reconstructed
678pattern-recognition track, with respect to the MC track. That is:
679
680.. math::
681
682 \frac{\omega^\textrm{MC} - \omega^\textrm{PR}}{\sigma_{\omega; \textrm{PR}}}
683
684.. seealso:: :b2:var:`omega`, :b2:var:`omegaErr`
685
686Returns NaN if no MC particle is related or if called on something other than a
687track-based particle.
688 )DOC");
689 REGISTER_VARIABLE("z0Pull", getHelixZ0Pull, R"DOC(
690The pull of the tracking parameter :math:`z_0` for the reconstructed
691pattern-recognition track, with respect to the MC track. That is:
692
693.. math::
694
695 \frac{z_0^\textrm{MC} - z_0^\textrm{PR}}{\sigma_{z_0; \textrm{PR}}}
696
697.. seealso:: :b2:var:`z0`, :b2:var:`z0Err`
698
699Returns NaN if no MC particle is related or if called on something other than a
700track-based particle.
701 )DOC");
702 REGISTER_VARIABLE("tanLambdaPull", getHelixTanLambdaPull, R"DOC(
703The pull of the tracking parameter :math:`\tan\lambda` for the reconstructed
704pattern-recognition track, with respect to the MC track. That is:
705
706.. math::
707
708 \frac{(\tan\lambda)^\textrm{MC} - (\tan\lambda)^\textrm{PR}}{\sigma_{\tan\lambda; \textrm{PR}}}
709
710.. seealso:: :b2:var:`tanLambda`, :b2:var:`tanLambdaErr`
711
712Returns NaN if no MC particle is related or if called on something other than a
713track-based particle.
714 )DOC");
715 REGISTER_VARIABLE("nCDCHits", trackNCDCHits,
716 "The number of CDC hits associated to the track. Returns NaN if called for something other than a track-based particle.");
717 REGISTER_VARIABLE("nSVDHits", trackNSVDHits,
718 "The number of SVD hits associated to the track. Returns NaN if called for something other than a track-based particle.");
719 REGISTER_VARIABLE("nPXDHits", trackNPXDHits,
720 "The number of PXD hits associated to the track. Returns NaN if called for something other than a track-based particle.");
721 REGISTER_VARIABLE("nVXDHits", trackNVXDHits,
722 "The number of PXD and SVD hits associated to the track. Returns NaN if called for something other than a track-based particle.");
723 REGISTER_VARIABLE("ndf", trackNDF, R"DOC(
724Returns the number of degrees of freedom of the track fit.
725
726.. note::
727
728 Note that this is not simply the number of hits -5 due to outlier hit
729 rejection.
730
731Returns NaN if called for something other than a track-based particle, or for
732mdst files processed with basf2 versions older than ``release-05-01``.
733 )DOC");
734 REGISTER_VARIABLE("chi2", trackChi2, R"DOC(
735Returns the :math:`\chi^2` of the track fit. This is actually computed based on
736:b2:var:`pValue` and :b2:var:`ndf`.
737
738.. note:: Note that for :b2:var:`pValue` exactly equal to 0 it returns infinity.
739
740Returns NaN if called for something other than a track-based particle, or for
741mdst files processed with basf2 versions older than ``release-05-01``.
742 )DOC");
743 REGISTER_VARIABLE("firstSVDLayer", trackFirstSVDLayer,
744 "The first activated SVD layer associated to the track. Returns NaN if called for something other than a track-based particle.");
745 REGISTER_VARIABLE("firstPXDLayer", trackFirstPXDLayer,
746 "The first activated PXD layer associated to the track. Returns NaN if called for something other than a track-based particle.");
747 REGISTER_VARIABLE("firstCDCLayer", trackFirstCDCLayer,
748 "The first activated CDC layer associated to the track. Returns NaN if called for something other than a track-based particle.");
749 REGISTER_VARIABLE("lastCDCLayer", trackLastCDCLayer,
750 "The last CDC layer associated to the track. Returns NaN if called for something other than a track-based particle.");
751 REGISTER_VARIABLE("d0", trackD0, R"DOC(
752Returns the tracking parameter :math:`d_0`, the signed distance to the
753point-of-closest-approach (POCA) in the :math:`r-\phi` plane.
754
755.. note::
756
757 Tracking parameters are with respect to the origin (0,0,0). For the
758 POCA with respect to the measured beam interaction point, see
759 :b2:var:`dr` (you probably want this unless you're doing a tracking
760 study or some debugging) and :b2:var:`d0FromIP`.
761
762Returns NaN if called for something other than a track-based particle.
763
764)DOC", "cm");
765 REGISTER_VARIABLE("phi0", trackPhi0, R"DOC(
766Returns the tracking parameter :math:`\phi_0`, the angle of the transverse
767momentum in the :math:`r-\phi` plane.
768
769.. note::
770
771 Tracking parameters are with respect to the origin (0,0,0). For the
772 POCA with respect to the measured beam interaction point, see
773 :b2:var:`phi0FromIP`.
774
775Returns NaN if called for something other than a track-based particle.
776
777)DOC", "rad");
778 REGISTER_VARIABLE("omega", trackOmega, R"DOC(
779Returns the tracking parameter :math:`\omega`, the curvature of the track.
780
781Returns NaN if called for something other than a track-based particle.
782
783)DOC", ":math:`\\text{cm}^{-1}`");
784 REGISTER_VARIABLE("z0", trackZ0, R"DOC(
785Returns the tracking parameter :math:`z_0`, the z-coordinate of the
786point-of-closest-approach (POCA).
787
788.. note::
789
790 Tracking parameters are with respect to the origin (0,0,0). For the
791 POCA with respect to the measured beam interaction point, see
792 :b2:var:`dz` (you probably want this unless you're doing a tracking
793 study or some debugging) and :b2:var:`z0FromIP`.
794
795Returns NaN if called for something other than a track-based particle.
796
797)DOC", "cm");
798 REGISTER_VARIABLE("tanLambda", trackTanLambda, R"DOC(
799Returns :math:`\tan\lambda`, the slope of the track in the :math:`r-z` plane.
800
801Returns NaN if called for something other than a track-based particle.
802 )DOC");
803 REGISTER_VARIABLE("d0FromIP", trackD0FromIP, R"DOC(
804Returns the tracking parameter :math:`d_0`, the signed distance to the
805point-of-closest-approach (POCA) in the :math:`r-\phi` plane, with respect to the measured beam interaction point.
806
807Returns NaN if called for something other than a track-based particle.
808
809)DOC", "cm");
810 REGISTER_VARIABLE("z0FromIP", trackZ0FromIP, R"DOC(
811Returns the tracking parameter :math:`z_0`, the z-coordinate of the
812point-of-closest-approach (POCA), with respect to the measured beam interaction point.
813
814Returns NaN if called for something other than a track-based particle.
815
816)DOC", "cm");
817 REGISTER_VARIABLE("phi0FromIP", trackPhi0FromIP, R"DOC(
818Returns the tracking parameter :math:`\phi_0`, the angle of the transverse
819momentum in the :math:`r-\phi` plane, with respect to the measured beam interaction point.
820
821Returns NaN if called for something other than a track-based particle.
822
823)DOC", "rad");
824 REGISTER_VARIABLE("d0Err", trackD0Error, R"DOC(
825Returns the uncertainty on :math:`d_0`, the signed distance to the
826point-of-closest-approach (POCA) in the :math:`r-\phi` plane.
827
828.. seealso:: :b2:var:`d0`, :b2:var:`d0Pull`
829
830Returns NaN if called for something other than a track-based particle.
831
832)DOC", "cm");
833 REGISTER_VARIABLE("phi0Err", trackPhi0Error, R"DOC(
834Returns the uncertainty on :math:`\phi_0`, the angle of the transverse momentum
835in the :math:`r-\phi` plane.
836
837.. seealso:: :b2:var:`phi0`, :b2:var:`phi0Pull`
838
839Returns NaN if called for something other than a track-based particle.
840
841)DOC", "rad");
842 REGISTER_VARIABLE("omegaErr", trackOmegaError, R"DOC(
843Returns the uncertainty on :math:`\omega`, the curvature of the track.
844
845.. seealso:: :b2:var:`omega`, :b2:var:`omegaPull`
846
847Returns NaN if called for something other than a track-based particle.
848
849)DOC", ":math:`\\text{cm}^{-1}`");
850 REGISTER_VARIABLE("z0Err", trackZ0Error, R"DOC(
851Returns the uncertainty on :math:`z_0`, the z-coordinate of the
852point-of-closest-approach (POCA).
853
854.. seealso:: :b2:var:`z0`, :b2:var:`z0Pull`
855
856Returns NaN if called for something other than a track-based particle."
857
858)DOC", "cm");
859 REGISTER_VARIABLE("tanLambdaErr", trackTanLambdaError, R"DOC(
860Returns the uncertainty on :math:`\tan\lambda`, the slope of the track in the
861:math:`r-z` plane.
862
863.. seealso:: :b2:var:`tanLambda`, :b2:var:`tanLambdaPull`
864
865Returns NaN if called for something other than a track-based particle.
866 )DOC");
867 REGISTER_VARIABLE("trackFitCovariance(i, j)", trackFitCovariance, R"DOC(
868 The track fit covariance matrix element corresponding to the two indices is returned.
869 This is the association between integers and parameters:
870
871 * 0: :math:`d_0`
872 * 1: :math:`\phi_0`
873 * 2: :math:`\omega`
874 * 3: :math:`z_0`
875 * 4: :math:`\tan\lambda`
876
877 .. note::
878
879 The covariance is returned. This means that the return value can be negative.
880 Furthermore, it's the squared value of the track fit error variables :b2:var:`d0Err`, etc.
881 when selecting the diagonal entries.
882
883 )DOC");
884 REGISTER_VARIABLE("pValue", trackPValue, R"DOC(
885The :math:`\chi^2` probability of the **track** fit.
886
887.. note::
888
889 This is the p-value of the track-fit. It does not get updated after
890 vertex fitting or kinematic fitting, and is meaningless for composite
891 particles.
892
893 See :b2:var:`chiProb` (you probably want this for high-level analysis).
894
895Returns NaN if called for something other than a track-based particle.
896 )DOC");
897 REGISTER_VARIABLE("trackFitHypothesisPDG", trackFitHypothesisPDG, R"DOC(
898Returns the PDG code of the track hypothesis actually used for the fit.
899Returns NaN if called for something other than a track-based particle.
900 )DOC");
901 REGISTER_VARIABLE("trackNECLClusters", trackNECLClusters, R"DOC(
902Returns a count of the number of ECLClusters matched to the track. This is
903always 0 or 1 with newer versions of ECL reconstruction.
904
905.. note::
906
907 For high-level analysis it is recommended to require the presence of a
908 matched ECL cluster along with a minimum energy requirement. A
909 track-based particle will have a clusterE if it is matched (NaN if
910 there is no cluster match for the track.
911
912 .. code-block:: python
913
914 import modularAnalysis as ma
915 # minimum energy of 200 MeV
916 ma.fillParticleList("e+:clusters", "clusterE > 0.2", path)
917
918 # these two are equivalent
919 ma.fillParticleList("e+:unmatched", "isNAN(clusterE) == 1", path)
920 ma.fillParticleList("e+:unmatched2", "trackNECLClusters == 0", path)
921
922Returns NaN if called for something other than a track-based particle.
923 )DOC");
924 REGISTER_VARIABLE("helixExtTheta(radius [cm], z fwd [cm], z bwd [cm], useHighestProbMass=0)", trackHelixExtTheta,
925 R"DOC(Returns theta of extrapolated helix parameters. If ``useHighestProbMass=1`` is set, the extrapolation will
926 use the track fit result for the mass hypothesis with the highest pValue.
927
928 )DOC", "rad");
929 REGISTER_VARIABLE("helixExtPhi(radius, z fwd, z bwd, useHighestProbMass=0)", trackHelixExtPhi,
930 "Returns phi of extrapolated helix parameters. If ``useHighestProbMass=1`` is set, the extrapolation will use the track fit result for the mass hypothesis with the highest pValue.\n\n",
931 "rad");
932
933 REGISTER_METAVARIABLE("helixExtThetaOnDet(detector_surface_name, useHighestProbMass=0)", trackHelixExtThetaOnDet,
934 R"DOC(Returns theta of extrapolated helix parameters on the given detector surface. The unit of angle is ``rad``.
935 If ``useHighestProbMass=1`` is set, the extrapolation will use the track fit result for the mass hypothesis with the highest pValue.
936 The supported detector surface names are ``{'CDC', 'TOP', 'ARICH', 'ECL', 'KLM'}``.
937 Also, the detector name with number of meaningful-layer is supported, e.g. ``'CDC8'``: last superlayer of CDC, ``'ECL1'``: mid-point of ECL.
938
939 ..note:: You can find more information in `modularAnalysis.calculateTrackIsolation`.
940 )DOC", Manager::VariableDataType::c_double);
941 REGISTER_METAVARIABLE("helixExtPhiOnDet(detector_surface_name, useHighestProbMass=0)", trackHelixExtPhiOnDet,
942 R"DOC(Returns phi of extrapolated helix parameters on the given detector surface. The unit of angle is ``rad``.
943 If ``useHighestProbMass=1`` is set, the extrapolation will use the track fit result for the mass hypothesis with the highest pValue.
944 The supported detector surface names are ``{'CDC', 'TOP', 'ARICH', 'ECL', 'KLM'}``.
945 Also, the detector name with number of meaningful-layer is supported, e.g. ``'CDC8'``: last superlayer of CDC, ``'ECL1'``: mid-point of ECL.
946
947 ..note:: You can find more information in `modularAnalysis.calculateTrackIsolation`.
948 )DOC", Manager::VariableDataType::c_double);
949
950
951 REGISTER_VARIABLE("nExtraCDCHits", nExtraCDCHits, R"DOC(
952[Eventbased] The number of CDC hits in the event not assigned to any track.
953
954Returns NaN if there is no event-level tracking information available.
955 )DOC");
956 REGISTER_VARIABLE("nExtraCDCHitsPostCleaning", nExtraCDCHitsPostCleaning, R"DOC(
957[Eventbased] Returns a count of the number of CDC hits in the event not assigned
958to any track nor very likely beam background (i.e. hits that survive a cleanup
959selection).
960
961Returns NaN if there is no event-level tracking information available.
962 )DOC");
963 REGISTER_VARIABLE("hasExtraCDCHitsInLayer(i)", hasExtraCDCHitsInLayer, R"DOC(
964[Eventbased] Returns 1 if a non-assigned hit exists in the specified CDC layer,
9650 otherwise.
966
967Returns NaN if there is no event-level tracking information available.
968 )DOC");
969 REGISTER_VARIABLE("hasExtraCDCHitsInSuperLayer(i)", hasExtraCDCHitsInSuperLayer, R"DOC(
970[Eventbased] Returns 1 if a non-assigned hit exists in the specified CDC
971SuperLayer, 0 otherwise.
972
973Returns NaN if there is no event-level tracking information available.
974 )DOC");
975 REGISTER_VARIABLE("nExtraCDCSegments", nExtraCDCSegments, R"DOC(
976[Eventbased] Returns the number of CDC segments not assigned to any track.
977
978Returns NaN if there is no event-level tracking information available.
979 )DOC");
980 // TODO: once the Tracking group fill the dataobject these can be
981 // uncommented - at the moment they are not filled, so leave out
982 //REGISTER_VARIABLE("nExtraVXDHitsInLayer(i)", nExtraVXDHitsInLayer,
983 //"[Eventbased] The number VXD hits not assigned in the specified VXD layer");
984 //REGISTER_VARIABLE("nExtraVXDHits", nExtraVXDHits, "[Eventbased] The number of VXD hits not assigned to any track");
985 //REGISTER_VARIABLE("svdFirstSampleTime", svdFirstSampleTime, "[Eventbased] The time of first SVD sample relatvie to event T0");
986 REGISTER_VARIABLE("trackFindingFailureFlag", trackFindingFailureFlag, R"DOC(
987[Eventbased] Returns a flag set by the tracking if there is reason to assume
988there was a track in the event missed by the tracking, or the track finding was
989(partly) aborted for this event.
990
991Returns NaN if there is no event-level tracking information available.
992 )DOC");
993
994 REGISTER_VARIABLE("isTrackFlippedAndRefitted", isTrackFlippedAndRefitted, R"DOC(
995Returns 1 if the charged final state particle comes from a track that has been flipped and refitted
996at the end of the tracking chain, in particular before the outer detector reconstruction.
997 )DOC");
998
999 REGISTER_VARIABLE("trackTime", getTrackTime, R"DOC(
1000Returns the time at which the track is produced relative to the time of the collision (given by SVD EventT0).
1001Both the time of the collision and the track time are computed using only SVD hits.
1002Returns NaN if SVD EventT0 is NaN, or if no SVD Hits are attached to the track.
1003For more details, see :ref:`Time Extraction <tracking_eventTimeExtraction>` page.
1004
1005)DOC", "ns");
1006
1007 REGISTER_VARIABLE("trackLength", getTrackLength, R"DOC(
1008Returns the arc length of the helix for the TrackFitResult associated with the particle.
1009The arc length is measured from the track origin to the radius of the CDC layer in which the Track has a hit.
1010Returns NaN if the particle has no CDC Hits.
1011
1012)DOC", "cm");
1013
1014
1015 }
1017}
static ROOT::Math::XYZVector getFieldInTesla(const ROOT::Math::XYZVector &pos)
return the magnetic field at a given position in Tesla.
EDetector
Enum for identifying the detector components (detector and subdetector).
Definition Const.h:42
static const double doubleNaN
quiet_NaN
Definition Const.h:703
@ c_nPhotons
CR is split into n photons (N1)
Definition ECLCluster.h:41
std::function< VarVariant(const Particle *)> FunctionPtr
functions stored take a const Particle* and return VarVariant.
Definition Manager.h:112
double sqrt(double a)
sqrt for double
Definition beamHelpers.h:28
Abstract base class for different kinds of events.
static const std::unordered_map< int, double > cdcWireRadiuses
CDC sense wire radiuses Values are take from cdc/data/CDC.xml.
static const std::unordered_map< std::string, DetSurfCylBoundaries > detToSurfBoundaries
Map that associates to each detector its valid cylindrical surface's boundaries.
static const std::unordered_map< std::string, DetSurfCylBoundaries > detLayerToSurfBoundaries
Map that associates to each detector layer its valid cylindrical surface's boundaries.