Belle II Software  release-06-02-00
ROEVariables.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 include
10 #include <analysis/variables/ROEVariables.h>
11 #include <analysis/variables/Variables.h>
12 #include <analysis/utility/PCmsLabTransform.h>
13 
14 #include <analysis/VariableManager/Manager.h>
15 
16 // framework - DataStore
17 #include <framework/datastore/StoreArray.h>
18 #include <framework/datastore/StoreObjPtr.h>
19 
20 // dataobjects
21 #include <analysis/dataobjects/Particle.h>
22 #include <analysis/dataobjects/ParticleList.h>
23 
24 #include <mdst/dataobjects/ECLCluster.h>
25 
26 // framework aux
27 #include <framework/logging/Logger.h>
28 #include <framework/gearbox/Const.h>
29 
30 // utility
31 #include <analysis/utility/ReferenceFrame.h>
32 
33 #include <TRandom.h>
34 #include <TMath.h>
35 
36 #include <iostream>
37 
38 using namespace std;
39 
40 namespace Belle2 {
45  namespace Variable {
46 
47  double isInRestOfEvent(const Particle* particle)
48  {
49 
50  StoreObjPtr<RestOfEvent> roeobjptr;
51  if (not roeobjptr.isValid())
52  return 0;
53 
54  const RestOfEvent* roe = &(*roeobjptr);
55 
56  return isInThisRestOfEvent(particle, roe);
57  }
58 
59  double isCloneOfSignalSide(const Particle* particle)
60  {
61 
62  StoreObjPtr<RestOfEvent> roe;
63  if (not roe.isValid()) {
64  B2WARNING("Please use isCloneOfSignalSide variable in for_each ROE loop!");
65  return std::numeric_limits<float>::quiet_NaN();
66  }
67  auto* particleMC = particle->getMCParticle();
68  if (!particleMC) {
69  return 0.0;
70  }
71  auto* signal = roe->getRelatedFrom<Particle>();
72  auto signalFSPs = signal->getFinalStateDaughters();
73  for (auto* daughter : signalFSPs) {
74  auto* daughterMC = daughter->getMCParticle();
75  if (daughterMC == particleMC) {
76  return 1.0;
77  }
78  }
79  return 0.0;
80  }
81  double hasAncestorFromSignalSide(const Particle* particle)
82  {
83  StoreObjPtr<RestOfEvent> roe;
84  if (!roe.isValid()) {
85  B2WARNING("Please use hasAncestorFromSignalSide variable in for_each ROE loop!");
86  return std::numeric_limits<float>::quiet_NaN();
87  }
88  auto* particleMC = particle->getMCParticle();
89  if (!particleMC) {
90  return 0.0;
91  }
92  auto* signalReco = roe->getRelatedFrom<Particle>();
93  auto* signalMC = signalReco->getMCParticle();
94  MCParticle* ancestorMC = particleMC->getMother();
95  while (ancestorMC) {
96  if (ancestorMC == signalMC) {
97  return 1.0;
98  }
99  ancestorMC = ancestorMC->getMother();
100  }
101  return 0.0;
102  }
103 
104 
105 
106  Manager::FunctionPtr currentROEIsInList(const std::vector<std::string>& arguments)
107  {
108  if (arguments.size() != 1)
109  B2FATAL("Wrong number of arguments (1 required) for meta function currentROEIsInList");
110 
111  std::string listName = arguments[0];
112 
113  auto func = [listName](const Particle*) -> double {
114 
115  StoreObjPtr<ParticleList> particleList(listName);
116  if (!(particleList.isValid()))
117  {
118  B2FATAL("Invalid Listname " << listName << " given to currentROEIsInList!");
119  }
120  StoreObjPtr<RestOfEvent> roe("RestOfEvent");
121 
122  if (not roe.isValid())
123  return 0;
124 
125  auto* particle = roe->getRelatedFrom<Particle>();
126  if (particle == nullptr)
127  {
128  B2ERROR("Relation between particle and ROE doesn't exist! currentROEIsInList() variable has to be called from ROE loop");
129  return std::numeric_limits<float>::quiet_NaN();
130  }
131  return particleList->contains(particle) ? 1 : 0;
132 
133  };
134  return func;
135  }
136 
137  Manager::FunctionPtr particleRelatedToCurrentROE(const std::vector<std::string>& arguments)
138  {
139  if (arguments.size() != 1)
140  B2FATAL("Wrong number of arguments (1 required) for meta function particleRelatedToCurrentROE");
141 
142  const Variable::Manager::Var* var = Manager::Instance().getVariable(arguments[0]);
143  auto func = [var](const Particle*) -> double {
144 
145  StoreObjPtr<RestOfEvent> roe("RestOfEvent");
146 
147  if (not roe.isValid())
148  return std::numeric_limits<float>::quiet_NaN();
149 
150  auto* particle = roe->getRelatedFrom<Particle>();
151  if (particle == nullptr)
152  {
153  B2ERROR("Relation between particle and ROE doesn't exist! particleRelatedToCurrentROE() variable has to be called from ROE loop");
154  return std::numeric_limits<float>::quiet_NaN();
155  }
156  return var->function(particle);
157 
158  };
159  return func;
160  }
161 
162  Manager::FunctionPtr useROERecoilFrame(const std::vector<std::string>& arguments)
163  {
164  if (arguments.size() == 1) {
165  const Variable::Manager::Var* var = Manager::Instance().getVariable(arguments[0]);
166  auto func = [var](const Particle * particle) -> double {
167  // Here we prioritize old variable behaviour first:
168  const RestOfEvent* roe = particle->getRelatedTo<RestOfEvent>();
169  // if related ROE not found, get the StoreArray pointer
170  if (roe == nullptr)
171  {
172  StoreObjPtr<RestOfEvent> roeObjPtr("RestOfEvent");
173  if (roeObjPtr.isValid()) {
174  roe = &*roeObjPtr;
175  }
176  }
177  if (roe == nullptr)
178  {
179  B2ERROR("Neither relation between particle and ROE doesn't exist nor ROE object has not been found!");
180  return std::numeric_limits<float>::quiet_NaN();
181  }
182  PCmsLabTransform T;
183  TLorentzVector pRecoil = T.getBeamFourMomentum() - roe->get4Vector();
184  Particle tmp(pRecoil, 0);
185  UseReferenceFrame<RestFrame> frame(&tmp);
186  double result = var->function(particle);
187  return result;
188  };
189  return func;
190  } else {
191  B2WARNING("Wrong number of arguments for meta function useROERecoilFrame");
192  return nullptr;
193  }
194  }
195 
196  // only the helper function
197  double nRemainingTracksInROE(const Particle* particle, const std::string& maskName)
198  {
199  StoreObjPtr<RestOfEvent> roe("RestOfEvent");
200  if (not roe.isValid())
201  return 0.0;
202  int n_roe_tracks = roe->getNTracks(maskName);
203  int n_par_tracks = 0;
204  const auto& daughters = particle->getFinalStateDaughters();
205  for (const auto& daughter : daughters) {
206  if (daughter->getParticleSource() == Particle::EParticleSourceObject::c_Track && roe->hasParticle(daughter, maskName)) {
207  n_par_tracks++;
208  }
209  }
210  return n_roe_tracks - n_par_tracks;
211  }
212 
213  Manager::FunctionPtr nROE_RemainingTracksWithMask(const std::vector<std::string>& arguments)
214  {
215  std::string maskName;
216 
217  if (arguments.size() == 0)
218  maskName = "";
219  else if (arguments.size() == 1)
220  maskName = arguments[0];
221  else
222  B2FATAL("At most 1 argument (name of mask) accepted for meta function nROETracks");
223 
224  auto func = [maskName](const Particle * particle) -> double {
225  return nRemainingTracksInROE(particle, maskName);
226  };
227  return func;
228  }
229 
230  double nROE_RemainingTracks(const Particle* particle)
231  {
232  return nRemainingTracksInROE(particle);
233  }
234 
235  double nROE_KLMClusters(const Particle* particle)
236  {
237  // Get related ROE object
238  const RestOfEvent* roe = getRelatedROEObject(particle);
239 
240  if (!roe) {
241  B2ERROR("Relation between particle and ROE doesn't exist!");
242  return std::numeric_limits<float>::quiet_NaN();
243  }
244 
245  return roe->getNKLMClusters();
246  }
247 
248  double ROE_MC_E(const Particle* particle)
249  {
250  const MCParticle* mcp = particle->getMCParticle();
251 
252  if (!mcp)
253  return std::numeric_limits<float>::quiet_NaN();
254 
255  PCmsLabTransform T;
256  TLorentzVector boostvec = T.getBeamFourMomentum();
257  auto mcroe4vector = boostvec - mcp->get4Vector();
258  const auto& frame = ReferenceFrame::GetCurrent();
259  auto frameMCRoe4Vector = frame.getMomentum(mcroe4vector);
260  return frameMCRoe4Vector.Energy();
261  }
262 
263  double ROE_MC_P(const Particle* particle)
264  {
265  const MCParticle* mcp = particle->getMCParticle();
266 
267  if (!mcp)
268  return std::numeric_limits<float>::quiet_NaN();
269 
270  PCmsLabTransform T;
271  TLorentzVector boostvec = T.getBeamFourMomentum();
272  auto mcroe4vector = boostvec - mcp->get4Vector();
273  const auto& frame = ReferenceFrame::GetCurrent();
274  auto frameMCRoe4Vector = frame.getMomentum(mcroe4vector);
275  return frameMCRoe4Vector.Vect().Mag();
276  }
277 
278  double ROE_MC_Px(const Particle* particle)
279  {
280  const MCParticle* mcp = particle->getMCParticle();
281 
282  if (!mcp)
283  return std::numeric_limits<float>::quiet_NaN();
284 
285  PCmsLabTransform T;
286  TLorentzVector boostvec = T.getBeamFourMomentum();
287  auto mcroe4vector = boostvec - mcp->get4Vector();
288  const auto& frame = ReferenceFrame::GetCurrent();
289  auto frameMCRoe4Vector = frame.getMomentum(mcroe4vector);
290 
291  return frameMCRoe4Vector.Vect().X();
292  }
293 
294  double ROE_MC_Py(const Particle* particle)
295  {
296  const MCParticle* mcp = particle->getMCParticle();
297 
298  if (!mcp)
299  return std::numeric_limits<float>::quiet_NaN();
300 
301  PCmsLabTransform T;
302  TLorentzVector boostvec = T.getBeamFourMomentum();
303  auto mcroe4vector = boostvec - mcp->get4Vector();
304  const auto& frame = ReferenceFrame::GetCurrent();
305  auto frameMCRoe4Vector = frame.getMomentum(mcroe4vector);
306 
307  return frameMCRoe4Vector.Vect().Y();
308  }
309 
310  double ROE_MC_Pz(const Particle* particle)
311  {
312  const MCParticle* mcp = particle->getMCParticle();
313 
314  if (!mcp)
315  return std::numeric_limits<float>::quiet_NaN();
316 
317  PCmsLabTransform T;
318  TLorentzVector boostvec = T.getBeamFourMomentum();
319  auto mcroe4vector = boostvec - mcp->get4Vector();
320  const auto& frame = ReferenceFrame::GetCurrent();
321  auto frameMCRoe4Vector = frame.getMomentum(mcroe4vector);
322 
323  return frameMCRoe4Vector.Vect().Z();
324  }
325 
326  double ROE_MC_Pt(const Particle* particle)
327  {
328  const MCParticle* mcp = particle->getMCParticle();
329 
330  if (!mcp)
331  return std::numeric_limits<float>::quiet_NaN();
332 
333  PCmsLabTransform T;
334  TLorentzVector boostvec = T.getBeamFourMomentum();
335  auto mcroe4vector = boostvec - mcp->get4Vector();
336  const auto& frame = ReferenceFrame::GetCurrent();
337  auto frameMCRoe4Vector = frame.getMomentum(mcroe4vector);
338 
339  return frameMCRoe4Vector.Vect().Perp();
340  }
341 
342  double ROE_MC_PTheta(const Particle* particle)
343  {
344  const MCParticle* mcp = particle->getMCParticle();
345 
346  if (!mcp)
347  return std::numeric_limits<float>::quiet_NaN();
348 
349  PCmsLabTransform T;
350  TLorentzVector boostvec = T.getBeamFourMomentum();
351  auto mcroe4vector = boostvec - mcp->get4Vector();
352  const auto& frame = ReferenceFrame::GetCurrent();
353  auto frameMCRoe4Vector = frame.getMomentum(mcroe4vector);
354 
355  return frameMCRoe4Vector.Theta();
356  }
357 
358  double ROE_MC_M(const Particle* particle)
359  {
360  const MCParticle* mcp = particle->getMCParticle();
361 
362  if (!mcp)
363  return std::numeric_limits<float>::quiet_NaN();
364 
365  PCmsLabTransform T;
366  TLorentzVector boostvec = T.getBeamFourMomentum();
367 
368  return (boostvec - mcp->get4Vector()).Mag();
369  }
370 
371  Manager::FunctionPtr ROE_MC_MissingFlags(const std::vector<std::string>& arguments)
372  {
373  std::string maskName;
374 
375  if (arguments.size() == 0)
376  maskName = "";
377  else if (arguments.size() == 1)
378  maskName = arguments[0];
379  else
380  B2FATAL("At most 1 argument (name of mask) accepted for meta function nROETracks");
381 
382  auto func = [maskName](const Particle * particle) -> double {
383 
384  StoreArray<Particle> particles;
385 
386  //Get MC Particle of the B meson
387  const MCParticle* mcParticle = particle->getMCParticle();
388 
389  if (!mcParticle)
390  return std::numeric_limits<float>::quiet_NaN();
391 
392  // Get Mother
393  const MCParticle* mcMother = mcParticle->getMother();
394 
395  if (!mcMother)
396  return std::numeric_limits<float>::quiet_NaN();
397 
398  // Get daughters
399  std::vector<MCParticle*> mcDaughters = mcMother->getDaughters();
400 
401  if (mcDaughters.size() != 2)
402  return std::numeric_limits<float>::quiet_NaN();
403 
404  // Get the companion B meson
405  MCParticle* mcROE = nullptr;
406  if (mcDaughters[0]->getArrayIndex() == mcParticle->getArrayIndex())
407  mcROE = mcDaughters[1];
408  else
409  mcROE = mcDaughters[0];
410 
411  // Get related ROE object
412  const RestOfEvent* roe = getRelatedROEObject(particle);
413 
414  std::set<const MCParticle*> mcROEObjects;
415 
416  auto roeParticles = roe->getParticles(maskName);
417  for (auto* roeParticle : roeParticles)
418  {
419  auto* mcroeParticle = roeParticle->getMCParticle();
420  if (mcroeParticle != nullptr) {
421  mcROEObjects.insert(mcroeParticle);
422  }
423  }
424  int flags = 0;
425  checkMCParticleMissingFlags(mcROE, mcROEObjects, flags);
426 
427  return flags;
428  };
429  return func;
430  }
431 
432  Manager::FunctionPtr nROE_Tracks(const std::vector<std::string>& arguments)
433  {
434  std::string maskName;
435 
436  if (arguments.size() == 0)
437  maskName = "";
438  else if (arguments.size() == 1)
439  maskName = arguments[0];
440  else
441  B2FATAL("At most 1 argument (name of mask) accepted for meta function nROETracks");
442 
443  auto func = [maskName](const Particle * particle) -> double {
444 
445  // Get related ROE object
446  const RestOfEvent* roe = getRelatedROEObject(particle);
447 
448  if (!roe)
449  {
450  B2ERROR("Relation between particle and ROE doesn't exist!");
451  return std::numeric_limits<float>::quiet_NaN();
452  }
453 
454  return roe->getNTracks(maskName);
455  };
456  return func;
457  }
458 
459  Manager::FunctionPtr nROE_ECLClusters(const std::vector<std::string>& arguments)
460  {
461  std::string maskName;
462 
463  if (arguments.size() == 0)
464  maskName = "";
465  else if (arguments.size() == 1)
466  maskName = arguments[0];
467  else
468  B2FATAL("At most 1 argument (name of mask) accepted for meta function nROEECLClusters");
469 
470  auto func = [maskName](const Particle * particle) -> double {
471 
472  // Get related ROE object
473  const RestOfEvent* roe = getRelatedROEObject(particle);
474 
475  if (!roe)
476  {
477  B2ERROR("Relation between particle and ROE doesn't exist!");
478  return std::numeric_limits<float>::quiet_NaN();
479  }
480 
481  return roe->getNECLClusters(maskName);
482  };
483  return func;
484  }
485 
486  Manager::FunctionPtr nROE_NeutralECLClusters(const std::vector<std::string>& arguments)
487  {
488  std::string maskName;
489 
490  if (arguments.size() == 0)
491  maskName = "";
492  else if (arguments.size() == 1)
493  maskName = arguments[0];
494  else
495  B2FATAL("At most 1 argument (name of mask) accepted for meta function nROENeutralECLClusters");
496 
497  auto func = [maskName](const Particle * particle) -> double {
498 
499  // Get related ROE object
500  const RestOfEvent* roe = getRelatedROEObject(particle);
501 
502  if (!roe)
503  {
504  B2ERROR("Relation between particle and ROE doesn't exist!");
505  return std::numeric_limits<float>::quiet_NaN();
506  }
507 
508  return roe->getPhotons(maskName).size();
509  };
510  return func;
511  }
512 
513  Manager::FunctionPtr nROE_Photons(const std::vector<std::string>& arguments)
514  {
515  std::string maskName = "";
516 
517  if (arguments.size() == 1) {
518  maskName = arguments[0];
519  }
520  if (arguments.size() > 1) {
521  B2FATAL("At most 1 argument (name of mask) accepted for meta function nROE_Photons");
522  }
523  auto func = [maskName](const Particle * particle) -> double {
524 
525  // Get related ROE object
526  const RestOfEvent* roe = getRelatedROEObject(particle);
527 
528  if (!roe)
529  {
530  B2ERROR("Relation between particle and ROE doesn't exist!");
531  return std::numeric_limits<float>::quiet_NaN();
532  }
533 
534  // Get unused ECLClusters in ROE
535  auto roeClusters = roe->getPhotons(maskName);
536  int nPhotons = 0;
537 
538  // Select ECLClusters with photon hypothesis
539  for (auto& roeCluster : roeClusters)
540  if (roeCluster->getECLClusterEHypothesisBit() == ECLCluster::EHypothesisBit::c_nPhotons)
541  nPhotons++;
542  return nPhotons;
543  };
544  return func;
545  }
546 
547  Manager::FunctionPtr nROE_NeutralHadrons(const std::vector<std::string>& arguments)
548  {
549  std::string maskName = "";
550 
551  if (arguments.size() == 1) {
552  maskName = arguments[0];
553  }
554  if (arguments.size() > 1) {
555  B2FATAL("At most 1 argument (name of mask) accepted for meta function nROE_NeutralHadrons");
556  }
557  auto func = [maskName](const Particle * particle) -> double {
558 
559  // Get related ROE object
560  const RestOfEvent* roe = getRelatedROEObject(particle);
561 
562  if (!roe)
563  {
564  B2ERROR("Relation between particle and ROE doesn't exist!");
565  return std::numeric_limits<float>::quiet_NaN();
566  }
567 
568  return roe->getHadrons(maskName).size();
569  };
570  return func;
571  }
572 
573  Manager::FunctionPtr nROE_ChargedParticles(const std::vector<std::string>& arguments)
574  {
575  std::string maskName = "";
576  int pdgCode = 0;
577  if (arguments.size() == 1) {
578  maskName = arguments[0];
579  }
580  if (arguments.size() == 2) {
581  maskName = arguments[0];
582  try {
583  pdgCode = Belle2::convertString<int>(arguments[1]);
584  } catch (std::invalid_argument&) {
585  B2ERROR("First argument of nROE_ChargedParticles must be a PDG code");
586  return nullptr;
587  }
588  }
589  if (arguments.size() > 2) {
590  B2FATAL("Wrong number of arguments (2 optional) for meta function nROE_ChargedParticles");
591  }
592  auto func = [maskName, pdgCode](const Particle * particle) -> double {
593 
594  // Get related ROE object
595  const RestOfEvent* roe = getRelatedROEObject(particle);
596 
597  if (!roe)
598  {
599  B2ERROR("Relation between particle and ROE doesn't exist!");
600  return std::numeric_limits<float>::quiet_NaN();
601  }
602 
603  return roe->getChargedParticles(maskName, abs(pdgCode)).size();
604  };
605  return func;
606  }
607 
608  Manager::FunctionPtr nROE_Composites(const std::vector<std::string>& arguments)
609  {
610  std::string maskName = "";
611 
612  if (arguments.size() == 1) {
613  maskName = arguments[0];
614  }
615  if (arguments.size() > 1) {
616  B2FATAL("At most 1 argument (name of mask) accepted for meta function nROE_Composites");
617  }
618  auto func = [maskName](const Particle * particle) -> double {
619 
620  // Get related ROE object
621  const RestOfEvent* roe = getRelatedROEObject(particle);
622 
623  if (!roe)
624  {
625  B2ERROR("Relation between particle and ROE doesn't exist!");
626  return std::numeric_limits<float>::quiet_NaN();
627  }
628  int result = 0;
629  auto particles = roe->getParticles(maskName, false);
630 
631  for (auto roeParticle : particles)
632  {
633  if (roeParticle->getParticleSource() == Particle::c_Composite or
634  roeParticle->getParticleSource() == Particle::c_V0) {
635  result++;
636  }
637  }
638  return result;
639  };
640  return func;
641  }
642 
643 
644  Manager::FunctionPtr nROE_ParticlesInList(const std::vector<std::string>& arguments)
645  {
646  std::string pListName;
647  std::string maskName;
648 
649  if (arguments.size() == 1) {
650  pListName = arguments[0];
651  maskName = "";
652  } else if (arguments.size() == 2) {
653  pListName = arguments[0];
654  maskName = arguments[1];
655  } else
656  B2FATAL("Wrong number of arguments (1 or 2 required) for meta function nROE_ParticlesInList");
657 
658  auto func = [pListName, maskName](const Particle * particle) -> double {
659 
660  // Get related ROE object
661  const RestOfEvent* roe = getRelatedROEObject(particle);
662 
663  if (!roe)
664  {
665  B2ERROR("Relation between particle and ROE doesn't exist!");
666  return std::numeric_limits<float>::quiet_NaN();
667  }
668 
669  int nPart = 0;
670 
671  // Get particle list
672  StoreObjPtr<ParticleList> pList(pListName);
673  if (!pList.isValid())
674  B2FATAL("ParticleList " << pListName << " could not be found or is not valid!");
675 
676  for (unsigned int i = 0; i < pList->getListSize(); i++)
677  {
678  const Particle* part = pList->getParticle(i);
679  if (isInThisRestOfEvent(part, roe, maskName))
680  ++nPart;
681  }
682 
683  return nPart;
684  };
685  return func;
686  }
687 
688  Manager::FunctionPtr ROE_Charge(const std::vector<std::string>& arguments)
689  {
690  std::string maskName;
691 
692  if (arguments.size() == 0)
693  maskName = "";
694  else if (arguments.size() == 1)
695  maskName = arguments[0];
696  else
697  B2FATAL("At most 1 argument (name of mask) accepted for meta function ROECharge");
698 
699  auto func = [maskName](const Particle * particle) -> double {
700 
701  // Get related ROE object
702  const RestOfEvent* roe = getRelatedROEObject(particle);
703 
704  if (!roe)
705  {
706  B2ERROR("Relation between particle and ROE doesn't exist!");
707  return std::numeric_limits<float>::quiet_NaN();
708  }
709 
710  // Get tracks in ROE
711  auto roeParticles = roe->getParticles(maskName);
712  int roeCharge = 0;
713 
714  for (auto* roeParticle : roeParticles)
715  {
716  roeCharge += roeParticle->getCharge();
717  }
718 
719  return roeCharge;
720  };
721  return func;
722  }
723 
724  Manager::FunctionPtr ROE_ExtraEnergy(const std::vector<std::string>& arguments)
725  {
726  std::string maskName;
727 
728  if (arguments.size() == 0)
729  maskName = "";
730  else if (arguments.size() == 1)
731  maskName = arguments[0];
732  else
733  B2FATAL("At most 1 argument (name of mask) accepted for meta function extraEnergy");
734 
735  auto func = [maskName](const Particle * particle) -> double {
736 
737  // Get related ROE object
738  const RestOfEvent* roe = getRelatedROEObject(particle);
739 
740  if (!roe)
741  {
742  B2ERROR("Relation between particle and ROE doesn't exist!");
743  return std::numeric_limits<float>::quiet_NaN();
744  }
745 
746  double extraE = 0.0;
747 
748  auto roeClusters = roe->getPhotons(maskName);
749 
750  for (auto& roeCluster : roeClusters)
751  if (roeCluster->getECLClusterEHypothesisBit() == ECLCluster::EHypothesisBit::c_nPhotons)
752  extraE += roeCluster->getECLClusterEnergy();
753 
754  auto roeChargedParticles = roe->getChargedParticles(maskName);
755 
756  for (auto& roeChargedParticle : roeChargedParticles)
757  {
758  if (roeChargedParticle->getECLCluster())
759  extraE += roeChargedParticle->getECLClusterEnergy();
760  }
761 
762  return extraE;
763  };
764  return func;
765  }
766 
767  Manager::FunctionPtr ROE_NeutralExtraEnergy(const std::vector<std::string>& arguments)
768  {
769  std::string maskName;
770 
771  if (arguments.size() == 0)
772  maskName = "";
773  else if (arguments.size() == 1)
774  maskName = arguments[0];
775  else
776  B2FATAL("At most 1 argument (name of mask) accepted for meta function extraEnergy");
777 
778  auto func = [maskName](const Particle * particle) -> double {
779 
780  // Get related ROE object
781  const RestOfEvent* roe = getRelatedROEObject(particle);
782 
783  if (!roe)
784  {
785  B2ERROR("Relation between particle and ROE doesn't exist!");
786  return std::numeric_limits<float>::quiet_NaN();
787  }
788  auto roephotons = roe->getPhotons(maskName);
789  TLorentzVector total4vector;
790  for (auto* photon : roephotons)
791  {
792  total4vector += photon->get4Vector();
793  }
794  const auto& frame = ReferenceFrame::GetCurrent();
795  auto frameRoe4Vector = frame.getMomentum(total4vector);
796  return frameRoe4Vector.Energy();
797  };
798  return func;
799  }
800 
801  Manager::FunctionPtr ROE_E(const std::vector<std::string>& arguments)
802  {
803  std::string maskName;
804  if (arguments.size() == 0)
805  maskName = "";
806  else if (arguments.size() == 1)
807  maskName = arguments[0];
808  else
809  B2FATAL("At most 1 argument (name of mask) accepted for meta function ROE_E");
810  auto func = [maskName](const Particle * particle) -> double {
811  const RestOfEvent* roe = particle->getRelatedTo<RestOfEvent>();
812  if (!roe)
813  {
814  B2ERROR("Relation between particle and ROE doesn't exist!");
815  return std::numeric_limits<float>::quiet_NaN();
816  }
817  const auto& frame = ReferenceFrame::GetCurrent();
818  auto frameRoe4Vector = frame.getMomentum(roe->get4Vector(maskName));
819  return frameRoe4Vector.Energy();
820  };
821  return func;
822  }
823 
824  Manager::FunctionPtr ROE_M(const std::vector<std::string>& arguments)
825  {
826  std::string maskName;
827 
828  if (arguments.size() == 0)
829  maskName = "";
830  else if (arguments.size() == 1)
831  maskName = arguments[0];
832  else
833  B2FATAL("At most 1 argument (name of mask) accepted for meta function ROE_M");
834 
835  auto func = [maskName](const Particle * particle) -> double {
836 
837  // Get related ROE object
838  const RestOfEvent* roe = getRelatedROEObject(particle);
839 
840  if (!roe)
841  {
842  B2ERROR("Relation between particle and ROE doesn't exist!");
843  return std::numeric_limits<float>::quiet_NaN();
844  }
845 
846  return roe->get4Vector(maskName).Mag();
847  };
848  return func;
849  }
850 
851  Manager::FunctionPtr ROE_P(const std::vector<std::string>& arguments)
852  {
853  std::string maskName;
854 
855  if (arguments.size() == 0)
856  maskName = "";
857  else if (arguments.size() == 1)
858  maskName = arguments[0];
859  else
860  B2FATAL("At most 1 argument (name of mask) accepted for meta function ROE_P");
861 
862  auto func = [maskName](const Particle * particle) -> double {
863 
864  // Get related ROE object
865  const RestOfEvent* roe = getRelatedROEObject(particle);
866 
867  if (!roe)
868  {
869  B2ERROR("Relation between particle and ROE doesn't exist!");
870  return std::numeric_limits<float>::quiet_NaN();
871  }
872 
873  const auto& frame = ReferenceFrame::GetCurrent();
874  auto frameRoe4Vector = frame.getMomentum(roe->get4Vector(maskName));
875  return frameRoe4Vector.Vect().Mag();
876  };
877  return func;
878  }
879 
880  Manager::FunctionPtr ROE_Px(const std::vector<std::string>& arguments)
881  {
882  std::string maskName;
883 
884  if (arguments.size() == 0)
885  maskName = "";
886  else if (arguments.size() == 1)
887  maskName = arguments[0];
888  else
889  B2FATAL("At most 1 argument (name of mask) accepted for meta function ROE_Px");
890 
891  auto func = [maskName](const Particle * particle) -> double {
892 
893  // Get related ROE object
894  const RestOfEvent* roe = getRelatedROEObject(particle);
895 
896  if (!roe)
897  {
898  B2ERROR("Relation between particle and ROE doesn't exist!");
899  return std::numeric_limits<float>::quiet_NaN();
900  }
901 
902  const auto& frame = ReferenceFrame::GetCurrent();
903  auto frameRoe4Vector = frame.getMomentum(roe->get4Vector(maskName));
904  return frameRoe4Vector.Vect().X();
905  };
906  return func;
907  }
908 
909  Manager::FunctionPtr ROE_Py(const std::vector<std::string>& arguments)
910  {
911  std::string maskName;
912 
913  if (arguments.size() == 0)
914  maskName = "";
915  else if (arguments.size() == 1)
916  maskName = arguments[0];
917  else
918  B2FATAL("At most 1 argument (name of mask) accepted for meta function ROE_Py");
919 
920  auto func = [maskName](const Particle * particle) -> double {
921 
922  // Get related ROE object
923  const RestOfEvent* roe = getRelatedROEObject(particle);
924 
925  if (!roe)
926  {
927  B2ERROR("Relation between particle and ROE doesn't exist!");
928  return std::numeric_limits<float>::quiet_NaN();
929  }
930 
931  const auto& frame = ReferenceFrame::GetCurrent();
932  auto frameRoe4Vector = frame.getMomentum(roe->get4Vector(maskName));
933  return frameRoe4Vector.Vect().Y();
934  };
935  return func;
936  }
937 
938  Manager::FunctionPtr ROE_Pt(const std::vector<std::string>& arguments)
939  {
940  std::string maskName;
941 
942  if (arguments.size() == 0)
943  maskName = "";
944  else if (arguments.size() == 1)
945  maskName = arguments[0];
946  else
947  B2FATAL("At most 1 argument (name of mask) accepted for meta function ROE_Pt");
948 
949  auto func = [maskName](const Particle * particle) -> double {
950 
951  // Get related ROE object
952  const RestOfEvent* roe = getRelatedROEObject(particle);
953 
954  if (!roe)
955  {
956  B2ERROR("Relation between particle and ROE doesn't exist!");
957  return std::numeric_limits<float>::quiet_NaN();
958  }
959 
960  const auto& frame = ReferenceFrame::GetCurrent();
961  auto frameRoe4Vector = frame.getMomentum(roe->get4Vector(maskName));
962  return frameRoe4Vector.Vect().Perp();
963  };
964  return func;
965  }
966 
967  Manager::FunctionPtr ROE_Pz(const std::vector<std::string>& arguments)
968  {
969  std::string maskName;
970 
971  if (arguments.size() == 0)
972  maskName = "";
973  else if (arguments.size() == 1)
974  maskName = arguments[0];
975  else
976  B2FATAL("At most 1 argument (name of mask) accepted for meta function ROE_Pz");
977 
978  auto func = [maskName](const Particle * particle) -> double {
979 
980  // Get related ROE object
981  const RestOfEvent* roe = getRelatedROEObject(particle);
982 
983  if (!roe)
984  {
985  B2ERROR("Relation between particle and ROE doesn't exist!");
986  return std::numeric_limits<float>::quiet_NaN();
987  }
988 
989  const auto& frame = ReferenceFrame::GetCurrent();
990  auto frameRoe4Vector = frame.getMomentum(roe->get4Vector(maskName));
991  return frameRoe4Vector.Vect().Z();
992  };
993  return func;
994  }
995 
996  Manager::FunctionPtr ROE_PTheta(const std::vector<std::string>& arguments)
997  {
998  std::string maskName;
999 
1000  if (arguments.size() == 0)
1001  maskName = "";
1002  else if (arguments.size() == 1)
1003  maskName = arguments[0];
1004  else
1005  B2FATAL("At most 1 argument (name of mask) accepted for meta function ROE_PTheta");
1006 
1007  auto func = [maskName](const Particle * particle) -> double {
1008 
1009  // Get related ROE object
1010  const RestOfEvent* roe = getRelatedROEObject(particle);
1011 
1012  if (!roe)
1013  {
1014  B2ERROR("Relation between particle and ROE doesn't exist!");
1015  return std::numeric_limits<float>::quiet_NaN();
1016  }
1017 
1018  const auto& frame = ReferenceFrame::GetCurrent();
1019  auto frameRoe4Vector = frame.getMomentum(roe->get4Vector(maskName));
1020  return frameRoe4Vector.Theta();
1021  };
1022  return func;
1023  }
1024 
1025  Manager::FunctionPtr ROE_DeltaE(const std::vector<std::string>& arguments)
1026  {
1027  std::string maskName;
1028 
1029  if (arguments.size() == 0)
1030  maskName = "";
1031  else if (arguments.size() == 1)
1032  maskName = arguments[0];
1033  else
1034  B2FATAL("At most 1 argument (name of mask) accepted for meta function ROE_deltae");
1035 
1036  auto func = [maskName](const Particle * particle) -> double {
1037 
1038  // Get related ROE object
1039  const RestOfEvent* roe = getRelatedROEObject(particle);
1040 
1041  if (!roe)
1042  {
1043  B2ERROR("Relation between particle and ROE doesn't exist!");
1044  return std::numeric_limits<float>::quiet_NaN();
1045  }
1046 
1047  PCmsLabTransform T;
1048  TLorentzVector vec = T.rotateLabToCms() * roe->get4Vector(maskName);
1049  return vec.E() - T.getCMSEnergy() / 2;
1050  };
1051  return func;
1052  }
1053 
1054  Manager::FunctionPtr ROE_Mbc(const std::vector<std::string>& arguments)
1055  {
1056  std::string maskName;
1057 
1058  if (arguments.size() == 0)
1059  maskName = "";
1060  else if (arguments.size() == 1)
1061  maskName = arguments[0];
1062  else
1063  B2FATAL("At most 1 argument (name of mask) accepted for meta function ROE_mbc");
1064 
1065  auto func = [maskName](const Particle * particle) -> double {
1066 
1067  // Get related ROE object
1068  const RestOfEvent* roe = getRelatedROEObject(particle);
1069 
1070  if (!roe)
1071  {
1072  B2ERROR("Relation between particle and ROE doesn't exist!");
1073  return std::numeric_limits<float>::quiet_NaN();
1074  }
1075 
1076  PCmsLabTransform T;
1077  TLorentzVector vec = T.rotateLabToCms() * roe->get4Vector(maskName);
1078 
1079  double E = T.getCMSEnergy() / 2;
1080  double m2 = E * E - vec.Vect().Mag2();
1081  double mbc = m2 > 0 ? sqrt(m2) : 0;
1082 
1083  return mbc;
1084  };
1085  return func;
1086  }
1087 
1088  Manager::FunctionPtr bssMassDifference(const std::vector<std::string>& arguments)
1089  {
1090  std::string maskName;
1091 
1092  if (arguments.size() == 0)
1093  maskName = "";
1094  else if (arguments.size() == 1)
1095  maskName = arguments[0];
1096  else
1097  B2FATAL("At most 1 argument (name of mask) accepted for meta function bssMassDifference");
1098 
1099  auto func = [maskName](const Particle * particle) -> double {
1100 
1101  // Get related ROE object
1102  TLorentzVector neutrino4vec = missing4Vector(particle->getDaughter(0), maskName, "6");
1103  TLorentzVector sig4vec = particle->getDaughter(0)->get4Vector();
1104 
1105  TLorentzVector bsMom = neutrino4vec + sig4vec;
1106  TLorentzVector bssMom = bsMom + particle->getDaughter(1)->get4Vector();
1107 
1108  return bssMom.M() - bsMom.M();
1109  };
1110  return func;
1111  }
1112 
1113  Manager::FunctionPtr WE_DeltaE(const std::vector<std::string>& arguments)
1114  {
1115  std::string maskName;
1116  std::string opt;
1117 
1118  if (arguments.size() == 1) {
1119  maskName = "";
1120  opt = arguments[0];
1121  } else if (arguments.size() == 2) {
1122  maskName = arguments[0];
1123  opt = arguments[1];
1124  } else
1125  B2FATAL("Wrong number of arguments (2 required) for meta function correctedB_deltae");
1126 
1127  auto func = [maskName, opt](const Particle * particle) -> double {
1128 
1129  PCmsLabTransform T;
1130  TLorentzVector boostvec = T.getBeamFourMomentum();
1131  TLorentzVector sig4vec = T.rotateLabToCms() * particle->get4Vector();
1132  TLorentzVector sig4vecLAB = particle->get4Vector();
1133  TLorentzVector neutrino4vec = missing4Vector(particle, maskName, "1");
1134  TLorentzVector neutrino4vecLAB = missing4Vector(particle, maskName, "6");
1135 
1136  double deltaE = std::numeric_limits<float>::quiet_NaN();
1137 
1138  // Definition 0: CMS
1139  if (opt == "0")
1140  {
1141  double totalSigEnergy = (sig4vec + neutrino4vec).Energy();
1142  double E = T.getCMSEnergy() / 2;
1143  deltaE = totalSigEnergy - E;
1144  }
1145 
1146  // Definition 1: LAB
1147  else if (opt == "1")
1148  {
1149  double Ecms = T.getCMSEnergy();
1150  double s = Ecms * Ecms;
1151  deltaE = ((sig4vecLAB + neutrino4vecLAB) * boostvec - s / 2.0) / sqrt(s);
1152  }
1153 
1154  else
1155  B2FATAL("Option for correctedB_deltae variable should only be 0/1 (CMS/LAB)");
1156 
1157  return deltaE;
1158  };
1159  return func;
1160  }
1161 
1162  Manager::FunctionPtr WE_Mbc(const std::vector<std::string>& arguments)
1163  {
1164  std::string maskName;
1165  std::string opt;
1166 
1167  if (arguments.size() == 1) {
1168  maskName = "";
1169  opt = arguments[0];
1170  } else if (arguments.size() == 2) {
1171  maskName = arguments[0];
1172  opt = arguments[1];
1173  } else
1174  B2FATAL("Wrong number of arguments (2 required) for meta function correctedB_mbc");
1175 
1176  auto func = [maskName, opt](const Particle * particle) -> double {
1177 
1178  PCmsLabTransform T;
1179  TLorentzVector boostvec = T.getBeamFourMomentum();
1180  TLorentzVector sig4vec = T.rotateLabToCms() * particle->get4Vector();
1181  TLorentzVector sig4vecLAB = particle->get4Vector();
1182  TLorentzVector neutrino4vec;
1183 
1184  double mbc = std::numeric_limits<float>::quiet_NaN();
1185 
1186  // Definition 0: CMS
1187  if (opt == "0")
1188  {
1189  neutrino4vec = missing4Vector(particle, maskName, "1");
1190  TVector3 bmom = (sig4vec + neutrino4vec).Vect();
1191  double E = T.getCMSEnergy() / 2;
1192  double m2 = E * E - bmom.Mag2();
1193  mbc = m2 > 0 ? sqrt(m2) : 0;
1194  }
1195 
1196  // Definition 1: LAB
1197  else if (opt == "1")
1198  {
1199  neutrino4vec = missing4Vector(particle, maskName, "6");
1200  TVector3 bmom = (sig4vecLAB + neutrino4vec).Vect();
1201  double Ecms = T.getCMSEnergy();
1202  double s = Ecms * Ecms;
1203  double m2 = pow((s / 2.0 + bmom * boostvec.Vect()) / boostvec.Energy(), 2.0) - bmom.Mag2();
1204  mbc = m2 > 0 ? sqrt(m2) : 0;
1205  }
1206 
1207  // Definition 2: CMS with factor alpha (so that dE == 0)
1208  else if (opt == "2")
1209  {
1210  neutrino4vec = missing4Vector(particle, maskName, "7");
1211  TVector3 bmom = (sig4vec + neutrino4vec).Vect();
1212  double E = T.getCMSEnergy() / 2;
1213  double m2 = E * E - bmom.Mag2();
1214  mbc = m2 > 0 ? sqrt(m2) : 0;
1215  }
1216 
1217  else
1218  B2FATAL("Option for correctedB_mbc variable should only be 0/1/2 (CMS/LAB/CMS with factor)");
1219 
1220  return mbc;
1221  };
1222  return func;
1223  }
1224 
1225  Manager::FunctionPtr WE_MissM2(const std::vector<std::string>& arguments)
1226  {
1227  std::string maskName;
1228  std::string opt;
1229 
1230  if (arguments.size() == 1) {
1231  maskName = "";
1232  opt = arguments[0];
1233  } else if (arguments.size() == 2) {
1234  maskName = arguments[0];
1235  opt = arguments[1];
1236  } else
1237  B2FATAL("Wrong number of arguments (2 required) for meta function missM2");
1238 
1239  auto func = [maskName, opt](const Particle * particle) -> double {
1240 
1241  return missing4Vector(particle, maskName, opt).Mag2();
1242  };
1243  return func;
1244  }
1245 
1246  Manager::FunctionPtr WE_MissPTheta(const std::vector<std::string>& arguments)
1247  {
1248  std::string maskName;
1249  std::string opt;
1250 
1251  if (arguments.size() == 1) {
1252  maskName = "";
1253  opt = arguments[0];
1254  } else if (arguments.size() == 2) {
1255  maskName = arguments[0];
1256  opt = arguments[1];
1257  } else
1258  B2FATAL("Wrong number of arguments (2 required) for meta function WE_MissPTheta");
1259 
1260  auto func = [maskName, opt](const Particle * particle) -> double {
1261 
1262  // Get related ROE object
1263  const RestOfEvent* roe = getRelatedROEObject(particle);
1264 
1265  if (!roe)
1266  {
1267  B2ERROR("Relation between particle and ROE doesn't exist!");
1268  return std::numeric_limits<float>::quiet_NaN();
1269  }
1270 
1271  return missing4Vector(particle, maskName, opt).Theta();
1272  };
1273  return func;
1274  }
1275 
1276  Manager::FunctionPtr WE_MissP(const std::vector<std::string>& arguments)
1277  {
1278  std::string maskName;
1279  std::string opt;
1280 
1281  if (arguments.size() == 1) {
1282  maskName = "";
1283  opt = arguments[0];
1284  } else if (arguments.size() == 2) {
1285  maskName = arguments[0];
1286  opt = arguments[1];
1287  } else
1288  B2FATAL("Wrong number of arguments (2 required) for meta function WE_MissP");
1289 
1290  auto func = [maskName, opt](const Particle * particle) -> double {
1291 
1292  // Get related ROE object
1293  const RestOfEvent* roe = getRelatedROEObject(particle);
1294 
1295  if (!roe)
1296  {
1297  B2ERROR("Relation between particle and ROE doesn't exist!");
1298  return std::numeric_limits<float>::quiet_NaN();
1299  }
1300 
1301  return missing4Vector(particle, maskName, opt).Vect().Mag();
1302  };
1303  return func;
1304  }
1305 
1306  Manager::FunctionPtr WE_MissPx(const std::vector<std::string>& arguments)
1307  {
1308  std::string maskName;
1309  std::string opt;
1310 
1311  if (arguments.size() == 1) {
1312  maskName = "";
1313  opt = arguments[0];
1314  } else if (arguments.size() == 2) {
1315  maskName = arguments[0];
1316  opt = arguments[1];
1317  } else
1318  B2FATAL("Wrong number of arguments (2 required) for meta function WE_MissPx");
1319 
1320  auto func = [maskName, opt](const Particle * particle) -> double {
1321 
1322  // Get related ROE object
1323  const RestOfEvent* roe = getRelatedROEObject(particle);
1324 
1325  if (!roe)
1326  {
1327  B2ERROR("Relation between particle and ROE doesn't exist!");
1328  return std::numeric_limits<float>::quiet_NaN();
1329  }
1330 
1331  return missing4Vector(particle, maskName, opt).Vect().Px();
1332  };
1333  return func;
1334  }
1335 
1336  Manager::FunctionPtr WE_MissPy(const std::vector<std::string>& arguments)
1337  {
1338  std::string maskName;
1339  std::string opt;
1340 
1341  if (arguments.size() == 1) {
1342  maskName = "";
1343  opt = arguments[0];
1344  } else if (arguments.size() == 2) {
1345  maskName = arguments[0];
1346  opt = arguments[1];
1347  } else
1348  B2FATAL("Wrong number of arguments (2 required) for meta function WE_MissPy");
1349 
1350  auto func = [maskName, opt](const Particle * particle) -> double {
1351 
1352  // Get related ROE object
1353  const RestOfEvent* roe = getRelatedROEObject(particle);
1354 
1355  if (!roe)
1356  {
1357  B2ERROR("Relation between particle and ROE doesn't exist!");
1358  return std::numeric_limits<float>::quiet_NaN();
1359  }
1360 
1361  return missing4Vector(particle, maskName, opt).Vect().Py();
1362  };
1363  return func;
1364  }
1365 
1366  Manager::FunctionPtr WE_MissPz(const std::vector<std::string>& arguments)
1367  {
1368  std::string maskName;
1369  std::string opt;
1370 
1371  if (arguments.size() == 1) {
1372  maskName = "";
1373  opt = arguments[0];
1374  } else if (arguments.size() == 2) {
1375  maskName = arguments[0];
1376  opt = arguments[1];
1377  } else
1378  B2FATAL("Wrong number of arguments (2 required) for meta function WE_MissPz");
1379 
1380  auto func = [maskName, opt](const Particle * particle) -> double {
1381 
1382  // Get related ROE object
1383  const RestOfEvent* roe = getRelatedROEObject(particle);
1384 
1385  if (!roe)
1386  {
1387  B2ERROR("Relation between particle and ROE doesn't exist!");
1388  return std::numeric_limits<float>::quiet_NaN();
1389  }
1390 
1391  return missing4Vector(particle, maskName, opt).Vect().Pz();
1392  };
1393  return func;
1394  }
1395 
1396  Manager::FunctionPtr WE_MissE(const std::vector<std::string>& arguments)
1397  {
1398  std::string maskName;
1399  std::string opt;
1400 
1401  if (arguments.size() == 1) {
1402  maskName = "";
1403  opt = arguments[0];
1404  } else if (arguments.size() == 2) {
1405  maskName = arguments[0];
1406  opt = arguments[1];
1407  } else
1408  B2FATAL("Wrong number of arguments (2 required) for meta function WE_MissE");
1409 
1410  auto func = [maskName, opt](const Particle * particle) -> double {
1411 
1412  // Get related ROE object
1413  const RestOfEvent* roe = getRelatedROEObject(particle);
1414 
1415  if (!roe)
1416  {
1417  B2ERROR("Relation between particle and ROE doesn't exist!");
1418  return std::numeric_limits<float>::quiet_NaN();
1419  }
1420 
1421  return missing4Vector(particle, maskName, opt).Energy();
1422  };
1423  return func;
1424  }
1425 
1426  Manager::FunctionPtr WE_xiZ(const std::vector<std::string>& arguments)
1427  {
1428  std::string maskName;
1429 
1430  if (arguments.size() == 0)
1431  maskName = "";
1432  else if (arguments.size() == 1)
1433  maskName = arguments[0];
1434  else
1435  B2FATAL("At most 1 argument (name of mask) accepted for meta function xiZ");
1436 
1437  auto func = [maskName](const Particle * particle) -> double {
1438 
1439  // Get related ROE object
1440  const RestOfEvent* roe = getRelatedROEObject(particle);
1441 
1442  if (!roe)
1443  {
1444  B2ERROR("Relation between particle and ROE doesn't exist!");
1445  return std::numeric_limits<float>::quiet_NaN();
1446  }
1447 
1448  double pz = 0;
1449  double energy = 0;
1450 
1451  // Get all Tracks on reconstructed side
1452  std::vector<const Particle*> recTrackParticles = particle->getFinalStateDaughters();
1453 
1454  // Loop the reconstructed side
1455  for (auto& recTrackParticle : recTrackParticles)
1456  {
1457  pz += recTrackParticle->getPz();
1458  energy += recTrackParticle->getEnergy();
1459  }
1460 
1461  // Loop the ROE side
1462  auto roeParticles = roe->getChargedParticles(maskName);
1463  for (auto* roeParticle : roeParticles)
1464  {
1465  pz += roeParticle->getPz();
1466  energy += roeParticle->getEnergy();
1467  }
1468 
1469  return pz / energy;
1470  };
1471  return func;
1472  }
1473 
1474  Manager::FunctionPtr WE_MissM2OverMissE(const std::vector<std::string>& arguments)
1475  {
1476  std::string maskName;
1477 
1478  if (arguments.size() == 0)
1479  maskName = "";
1480  else if (arguments.size() == 1)
1481  maskName = arguments[0];
1482  else
1483  B2FATAL("At most 1 argument (name of mask) accepted for meta function WE_MissM2OverMissE");
1484 
1485  auto func = [maskName](const Particle * particle) -> double {
1486 
1487  // Get related ROE object
1488  const RestOfEvent* roe = getRelatedROEObject(particle);
1489 
1490  if (!roe)
1491  {
1492  B2ERROR("Relation between particle and ROE doesn't exist!");
1493  return std::numeric_limits<float>::quiet_NaN();
1494  }
1495 
1496  PCmsLabTransform T;
1497  TLorentzVector missing4Momentum;
1498  TLorentzVector boostvec = T.getBeamFourMomentum();
1499 
1500  return missing4Vector(particle, maskName, "5").Mag2() / (2.0 * missing4Vector(particle, maskName, "5").Energy());
1501  };
1502  return func;
1503  }
1504 
1505  Manager::FunctionPtr WE_q2lnuSimple(const std::vector<std::string>& arguments)
1506  {
1507  std::string maskName("");
1508  std::string option("1");
1509 
1510  if (arguments.size() == 1) {
1511  maskName = arguments[0];
1512  } else if (arguments.size() == 2) {
1513  maskName = arguments[0];
1514  option = arguments[1];
1515  } else if (arguments.size() > 2) {
1516  B2FATAL("Too many arguments. At most two arguments are allowed for meta function q2lnuSimple(maskname,option)");
1517  }
1518 
1519  auto func = [maskName, option](const Particle * particle) -> double {
1520 
1521  // Get related ROE object
1522  const RestOfEvent* roe = getRelatedROEObject(particle);
1523 
1524  if (!roe)
1525  {
1526  B2ERROR("Relation between particle and ROE doesn't exist!");
1527  return std::numeric_limits<float>::quiet_NaN();
1528  }
1529 
1530  int n = particle->getNDaughters();
1531 
1532  if (n < 1)
1533  return std::numeric_limits<float>::quiet_NaN();
1534 
1535  // Assumes lepton is the last particle in the reco decay chain!
1536  PCmsLabTransform T;
1537  const Particle* lep = particle->getDaughter(n - 1);
1538  TLorentzVector lep4vec = T.rotateLabToCms() * lep->get4Vector();
1539  TLorentzVector nu4vec = missing4Vector(particle, maskName, option);
1540 
1541  return (lep4vec + nu4vec).Mag2();
1542  };
1543  return func;
1544  }
1545 
1546  Manager::FunctionPtr WE_q2lnu(const std::vector<std::string>& arguments)
1547  {
1548  std::string maskName("");
1549  std::string option("7");
1550 
1551  if (arguments.size() == 1) {
1552  maskName = arguments[0];
1553  } else if (arguments.size() == 2) {
1554  maskName = arguments[0];
1555  option = arguments[1];
1556  } else if (arguments.size() > 2) {
1557  B2FATAL("Too many arguments. At most two arguments are allowed for meta function q2lnu(maskname, option)");
1558  }
1559 
1560  auto func = [maskName, option](const Particle * particle) -> double {
1561 
1562  // Get related ROE object
1563  const RestOfEvent* roe = getRelatedROEObject(particle);
1564 
1565  if (!roe)
1566  {
1567  B2ERROR("Relation between particle and ROE doesn't exist!");
1568  return std::numeric_limits<float>::quiet_NaN();
1569  }
1570 
1571  int n = particle->getNDaughters();
1572 
1573  if (n < 1)
1574  return std::numeric_limits<float>::quiet_NaN();
1575 
1576  PCmsLabTransform T;
1577  const Particle* lep = particle->getDaughter(n - 1);
1578  TLorentzVector lep_cm = T.rotateLabToCms() * lep->get4Vector();
1579 
1580  TLorentzVector Y_cm = T.rotateLabToCms() * particle->get4Vector();
1581  TLorentzVector neu_cm = missing4Vector(particle, maskName, option);
1582 
1583  double e_beam = T.getCMSEnergy() / 2.0;
1584 
1585  //just to make the formula simpler
1586  double bmass = particle->getPDGMass();
1587  double pB2 = e_beam * e_beam - bmass * bmass;
1588 
1589  //angle between the Y and the neutrino, from the Mbc=M_B constraint
1590  double cos_angle_nu = (pB2 - Y_cm.Vect().Mag2() - neu_cm.Vect().Mag2()) / (2.0 * Y_cm.Vect().Mag() * neu_cm.Vect().Mag());
1591  if (abs(cos_angle_nu) > 1)
1592  {
1593  return (lep_cm + neu_cm).Mag2();
1594  }
1595 
1596  double angle_nu = TMath::ACos(cos_angle_nu);
1597  //first get one random neutrino, on the allowed cone for the constraint
1598  TLorentzVector rotated_neu(-1 * Y_cm.Vect(), Y_cm.E()); //first get reverse Y
1599 
1600  double nu_theta = rotated_neu.Theta() + (TMath::Pi() - angle_nu);
1601  double nu_phi = rotated_neu.Phi();
1602  rotated_neu.SetTheta(nu_theta);
1603  rotated_neu.SetPhi(nu_phi);
1604  rotated_neu.SetRho(neu_cm.Rho());
1605  rotated_neu.SetE(neu_cm.E());
1606 
1607  TVector3 Yneu_norm = Y_cm.Vect().Cross(neu_cm.Vect());
1608  TVector3 Yrot_norm = Y_cm.Vect().Cross(rotated_neu.Vect());
1609  //angle between the two crossproducts -> angle between the two vectors perpendicular to the Y-p_inc and Y-B planes -> angle between the planes
1610  //this angle needs to come out as zero
1611 
1612  double rot_angle = Yneu_norm.Angle(Yrot_norm);
1613 
1614  TLorentzVector rotated_neu2(rotated_neu);
1615  //unfortunately don't -and probably can't- know in which direction to rotate without trying
1616  //so create a copy of the vector, and later choose the correct one
1617  //However, rotation by 180 degrees is never needed, direction of the cross-product vector assures that when after rotation
1618  //the B-vector is in the plane, it always is on the side closer to pcm_lv_inc.
1619  //rotate around Y into the Y-neutrino-plane
1620  rotated_neu.Rotate(rot_angle, Y_cm.Vect());
1621  rotated_neu2.Rotate(-rot_angle, Y_cm.Vect());
1622 
1623  double dot1 = rotated_neu.Vect().Dot(Yneu_norm);
1624  double dot2 = rotated_neu2.Vect().Dot(Yneu_norm);
1625 
1626  if (abs(dot2) < abs(dot1)) rotated_neu = rotated_neu2;
1627 
1628  return (lep_cm + rotated_neu).Mag2();
1629  };
1630  return func;
1631  }
1632 
1633  Manager::FunctionPtr WE_cosThetaEll(const std::vector<std::string>& arguments)
1634  {
1635  std::string maskName;
1636 
1637  if (arguments.size() == 0)
1638  maskName = "";
1639  else if (arguments.size() == 1)
1640  maskName = arguments[0];
1641  else
1642  B2FATAL("At most 1 argument (name of mask) accepted for meta function cosThetaEll");
1643 
1644  auto func = [maskName](const Particle * particle) -> double {
1645 
1646  PCmsLabTransform T;
1647  TLorentzVector boostvec = T.getBeamFourMomentum();
1648  TLorentzVector pNu = missing4Vector(particle, maskName, "6");
1649 
1650  TLorentzVector pLep;
1651  // TODO: avoid hardocoded values
1652  for (unsigned i = 0; i < particle->getNDaughters(); i++)
1653  {
1654  int absPDG = abs(particle->getDaughter(i)->getPDGCode());
1655  if (absPDG == Const::electron.getPDGCode() || absPDG == Const::muon.getPDGCode() || absPDG == 15) {
1656  pLep = particle->getDaughter(i)->get4Vector();
1657  break;
1658  }
1659  }
1660 
1661  TLorentzVector pW = pNu + pLep;
1662  TLorentzVector pB = particle->get4Vector() + pNu;
1663 
1664  // boost lepton and B momentum to W frame
1665  TVector3 boost2W = -(pW.BoostVector());
1666  pLep.Boost(boost2W);
1667  pB.Boost(boost2W);
1668 
1669  TVector3 lep3Vector = pLep.Vect();
1670  TVector3 B3Vector = pB.Vect();
1671  double numerator = lep3Vector.Dot(B3Vector);
1672  double denominator = (lep3Vector.Mag()) * (B3Vector.Mag());
1673 
1674  return numerator / denominator;
1675  };
1676  return func;
1677  }
1678 
1679  Manager::FunctionPtr passesROEMask(const std::vector<std::string>& arguments)
1680  {
1681  std::string maskName;
1682 
1683  if (arguments.size() == 0)
1684  maskName = "";
1685  else if (arguments.size() == 1)
1686  maskName = arguments[0];
1687  else
1688  B2FATAL("At most 1 argument (name of mask) accepted for meta function passesROEMask");
1689 
1690  auto func = [maskName](const Particle * particle) -> double {
1691 
1692  double result = 0;
1693 
1694  StoreObjPtr<RestOfEvent> roe("RestOfEvent");
1695  if (not roe.isValid())
1696  return std::numeric_limits<float>::quiet_NaN();
1697 
1698  if (roe->hasParticle(particle, maskName))
1699  {
1700  return 1.0;
1701  }
1702 
1703  return result;
1704  };
1705  return func;
1706  }
1707 
1708  double printROE(const Particle* particle)
1709  {
1710  const RestOfEvent* roe = getRelatedROEObject(particle);
1711 
1712  if (!roe) {
1713  B2ERROR("Relation between particle and ROE doesn't exist!");
1714  } else roe->print();
1715  return 0.0;
1716  }
1717 
1718 
1719  Manager::FunctionPtr pi0Prob(const std::vector<std::string>& arguments)
1720  {
1721  if (arguments.size() != 1)
1722  B2ERROR("Wrong number of arguments (1 required) for pi0Prob");
1723 
1724  std::string mode;
1725  mode = arguments[0];
1726 
1727  if (mode != "standard" and mode != "tight" and mode != "cluster" and mode != "both")
1728  B2ERROR("the given argument is not supported in pi0Prob!");
1729 
1730  auto func = [mode](const Particle * particle) -> double {
1731  if (mode == "standard")
1732  {
1733  if (particle->hasExtraInfo("Pi0ProbOrigin")) {
1734  return particle->getExtraInfo("Pi0ProbOrigin");
1735  } else {
1736  B2WARNING("Pi0ProbOrigin is not registerted in extraInfo! \n"
1737  "the function writePi0EtaVeto has to be executed to register this extraInfo.");
1738  return std::numeric_limits<float>::quiet_NaN();
1739  }
1740  }
1741 
1742  else if (mode == "tight")
1743  {
1744  if (particle->hasExtraInfo("Pi0ProbTightEnergyThreshold")) {
1745  return particle->getExtraInfo("Pi0ProbTightEnergyThreshold");
1746  } else {
1747  B2WARNING("Pi0ProbTightEnergyThreshold is not registerted in extraInfo! \n"
1748  "the function writePi0EtaVeto has to be executed to register this extraInfo.");
1749  return std::numeric_limits<float>::quiet_NaN();
1750  }
1751  }
1752 
1753  else if (mode == "cluster")
1754  {
1755  if (particle->hasExtraInfo("Pi0ProbLargeClusterSize")) {
1756  return particle->getExtraInfo("Pi0ProbLargeClusterSize");
1757  } else {
1758  B2WARNING("Pi0ProbLargeClusterSize is not registerted in extraInfo! \n"
1759  "the function writePi0EtaVeto has to be executed to register this extraInfo.");
1760  return std::numeric_limits<float>::quiet_NaN();
1761  }
1762  }
1763 
1764  else if (mode == "both")
1765  {
1766  if (particle->hasExtraInfo("Pi0ProbTightEnergyThresholdAndLargeClusterSize")) {
1767  return particle->getExtraInfo("Pi0ProbTightEnergyThresholdAndLargeClusterSize");
1768  } else {
1769  B2WARNING("Pi0ProbTightEnergyThresholdAndLargeClusterSize is not registerted in extraInfo! \n"
1770  "the function writePi0EtaVeto has to be executed to register this extraInfo.");
1771  return std::numeric_limits<float>::quiet_NaN();
1772  }
1773  }
1774 
1775  else
1776  {
1777  return std::numeric_limits<float>::quiet_NaN();
1778  }
1779  };
1780  return func;
1781  }
1782 
1783  Manager::FunctionPtr etaProb(const std::vector<std::string>& arguments)
1784  {
1785  if (arguments.size() != 1)
1786  B2ERROR("Wrong number of arguments (1 required) for etaProb");
1787 
1788  std::string mode;
1789  mode = arguments[0];
1790 
1791  if (mode != "standard" and mode != "tight" and mode != "cluster" and mode != "both")
1792  B2ERROR("the given argument is not supported in etaProb!");
1793 
1794  auto func = [mode](const Particle * particle) -> double {
1795  if (mode == "standard")
1796  {
1797  if (particle->hasExtraInfo("EtaProbOrigin")) {
1798  return particle->getExtraInfo("EtaProbOrigin");
1799  } else {
1800  B2WARNING("EtaProbOrigin is not registerted in extraInfo! \n"
1801  "the function writePi0EtaVeto has to be executed to register this extraInfo.");
1802  return std::numeric_limits<float>::quiet_NaN();
1803  }
1804  }
1805 
1806  else if (mode == "tight")
1807  {
1808  if (particle->hasExtraInfo("EtaProbTightEnergyThreshold")) {
1809  return particle->getExtraInfo("EtaProbTightEnergyThreshold");
1810  } else {
1811  B2WARNING("EtaProbTightEnergyThreshold is not registerted in extraInfo! \n"
1812  "the function writePi0EtaVeto has to be executed to register this extraInfo.");
1813  return std::numeric_limits<float>::quiet_NaN();
1814  }
1815  }
1816 
1817  else if (mode == "cluster")
1818  {
1819  if (particle->hasExtraInfo("EtaProbLargeClusterSize")) {
1820  return particle->getExtraInfo("EtaProbLargeClusterSize");
1821  } else {
1822  B2WARNING("EtaProbLargeClusterSize is not registerted in extraInfo! \n"
1823  "the function writePi0EtaVeto has to be executed to register this extraInfo.");
1824  return std::numeric_limits<float>::quiet_NaN();
1825  }
1826  }
1827 
1828  else if (mode == "both")
1829  {
1830  if (particle->hasExtraInfo("EtaProbTightEnergyThresholdAndLargeClusterSize")) {
1831  return particle->getExtraInfo("EtaProbTightEnergyThresholdAndLargeClusterSize");
1832  } else {
1833  B2WARNING("EtaProbTightEnergyThresholdAndLargeClusterSize is not registerted in extraInfo! \n"
1834  "the function writePi0EtaVeto has to be executed to register this extraInfo.");
1835  return std::numeric_limits<float>::quiet_NaN();
1836  }
1837  }
1838 
1839  else
1840  {
1841  return std::numeric_limits<float>::quiet_NaN();
1842  }
1843  };
1844  return func;
1845  }
1846 
1847  // ------------------------------------------------------------------------------
1848  // Below are some functions for ease of usage, they are not a part of variables
1849  // ------------------------------------------------------------------------------
1850 
1851  TLorentzVector missing4Vector(const Particle* particle, const std::string& maskName, const std::string& opt)
1852  {
1853  // Get related ROE object
1854  const RestOfEvent* roe = getRelatedROEObject(particle);
1855 
1856  if (!roe) {
1857  B2ERROR("Relation between particle and ROE doesn't exist!");
1858  TLorentzVector empty;
1859  return empty;
1860  }
1861 
1862  PCmsLabTransform T;
1863  TLorentzVector boostvec = T.getBeamFourMomentum();
1864 
1865  TLorentzVector rec4vecLAB = particle->get4Vector();
1866  TLorentzVector roe4vecLAB = roe->get4Vector(maskName);
1867 
1868  TLorentzVector rec4vec = T.rotateLabToCms() * rec4vecLAB;
1869  TLorentzVector roe4vec = T.rotateLabToCms() * roe4vecLAB;
1870 
1871  TLorentzVector miss4vec;
1872  double E_beam_cms = T.getCMSEnergy() / 2.0;
1873 
1874  // Definition 0: CMS, use energy and momentum of tracks and clusters
1875  if (opt == "0") {
1876  miss4vec.SetVect(- (rec4vec.Vect() + roe4vec.Vect()));
1877  miss4vec.SetE(2 * E_beam_cms - (rec4vec.E() + roe4vec.E()));
1878  }
1879 
1880  // Definition 1: CMS, same as 0, fix Emiss = pmiss
1881  else if (opt == "1") {
1882  miss4vec.SetVect(- (rec4vec.Vect() + roe4vec.Vect()));
1883  miss4vec.SetE(miss4vec.Vect().Mag());
1884  }
1885 
1886  // Definition 2: CMS, same as 0, fix Eroe = Ecms/2
1887  else if (opt == "2") {
1888  miss4vec.SetVect(- (rec4vec.Vect() + roe4vec.Vect()));
1889  miss4vec.SetE(E_beam_cms - rec4vec.E());
1890  }
1891 
1892  // Definition 3: CMS, use only energy and momentum of signal side
1893  else if (opt == "3") {
1894  miss4vec.SetVect(- rec4vec.Vect());
1895  miss4vec.SetE(E_beam_cms - rec4vec.E());
1896  }
1897 
1898  // Definition 4: CMS, same as 3, update with direction of ROE momentum
1899  else if (opt == "4") {
1900  TVector3 pB = - roe4vec.Vect();
1901  pB.SetMag(0.340);
1902  miss4vec.SetVect(pB - rec4vec.Vect());
1903  miss4vec.SetE(E_beam_cms - rec4vec.E());
1904  }
1905 
1906  // Definition 5: LAB, use energy and momentum of tracks and clusters from whole event
1907  else if (opt == "5") {
1908  miss4vec.SetVect(boostvec.Vect() - (rec4vecLAB.Vect() + roe4vecLAB.Vect()));
1909  miss4vec.SetE(boostvec.E() - (rec4vecLAB.E() + roe4vecLAB.E()));
1910  }
1911 
1912  // Definition 6: LAB, same as 5, fix Emiss = pmiss
1913  else if (opt == "6") {
1914  miss4vec.SetVect(boostvec.Vect() - (rec4vecLAB.Vect() + roe4vecLAB.Vect()));
1915  miss4vec.SetE(miss4vec.Vect().Mag());
1916  }
1917 
1918  // Definition 7: CMS, correct pmiss 3-momentum vector with factor alpha so that dE = 0 (used for Mbc calculation)
1919  else if (opt == "7") {
1920  miss4vec.SetVect(-(rec4vec.Vect() + roe4vec.Vect()));
1921  miss4vec.SetE(miss4vec.Vect().Mag());
1922  double factorAlpha = (E_beam_cms - rec4vec.E()) / miss4vec.E();
1923  miss4vec.SetRho(factorAlpha * miss4vec.Rho());
1924  miss4vec.SetE(miss4vec.Rho());
1925  }
1926 
1927  return miss4vec;
1928  }
1929 
1930  void checkMCParticleMissingFlags(const MCParticle* mcp, std::set<const MCParticle*> mcROEObjects, int& missingFlags)
1931  {
1932  std::vector<MCParticle*> daughters = mcp->getDaughters();
1933  for (auto& daughter : daughters) {
1934 
1935  if (!daughter->hasStatus(MCParticle::c_PrimaryParticle))
1936  continue;
1937 
1938  if (mcROEObjects.find(daughter) == mcROEObjects.end()) {
1939 
1940  int pdg = abs(daughter->getPDG());
1941 
1942  // photon
1943  if (pdg == Const::photon.getPDGCode() and (missingFlags & 1) == 0)
1944  missingFlags += 1;
1945 
1946  // electrons
1947  else if (pdg == Const::electron.getPDGCode() and (missingFlags & 2) == 0)
1948  missingFlags += 2;
1949 
1950  // muons
1951  else if (pdg == Const::muon.getPDGCode() and (missingFlags & 4) == 0)
1952  missingFlags += 4;
1953 
1954  // pions
1955  else if (pdg == Const::pion.getPDGCode() and (missingFlags & 8) == 0)
1956  missingFlags += 8;
1957 
1958  // kaons
1959  else if (pdg == Const::kaon.getPDGCode() and (missingFlags & 16) == 0)
1960  missingFlags += 16;
1961 
1962  // protons
1963  else if (pdg == Const::proton.getPDGCode() and (missingFlags & 32) == 0)
1964  missingFlags += 32;
1965 
1966  // neutrons
1967  else if (pdg == Const::neutron.getPDGCode() and (missingFlags & 64) == 0)
1968  missingFlags += 64;
1969 
1970  // kshort
1971  else if (pdg == Const::Kshort.getPDGCode() and ((missingFlags & 128) == 0 or (missingFlags & 256) == 0)) {
1972  std::vector<MCParticle*> ksDaug = daughter->getDaughters();
1973  if (ksDaug.size() == 2) {
1974  // K_S0 -> pi+ pi-
1975  if (abs(ksDaug[0]->getPDG()) == Const::pion.getPDGCode() and abs(ksDaug[1]->getPDG()) == Const::pion.getPDGCode()
1976  and (missingFlags & 128) == 0) {
1977  if (mcROEObjects.find(ksDaug[0]) == mcROEObjects.end() or mcROEObjects.find(ksDaug[1]) == mcROEObjects.end())
1978  missingFlags += 128;
1979  }
1980  // K_S0 -> pi0 pi0
1981  else if (abs(ksDaug[0]->getPDG()) == Const::pi0.getPDGCode() and abs(ksDaug[1]->getPDG()) == Const::pi0.getPDGCode()
1982  and (missingFlags & 256) == 0) {
1983  std::vector<MCParticle*> pi0Daug0 = ksDaug[0]->getDaughters();
1984  std::vector<MCParticle*> pi0Daug1 = ksDaug[1]->getDaughters();
1985  if (mcROEObjects.find(pi0Daug0[0]) == mcROEObjects.end() or
1986  mcROEObjects.find(pi0Daug0[1]) == mcROEObjects.end() or
1987  mcROEObjects.find(pi0Daug1[0]) == mcROEObjects.end() or
1988  mcROEObjects.find(pi0Daug1[1]) == mcROEObjects.end())
1989  missingFlags += 256;
1990  }
1991  }
1992  }
1993 
1994  // klong
1995  else if (pdg == Const::Klong.getPDGCode() and (missingFlags & 512) == 0)
1996  missingFlags += 512;
1997 
1998  // neutrinos, which are not in the Const::
1999  else if ((pdg == 12 or pdg == 14 or pdg == 16) and (missingFlags & 1024) == 0)
2000  missingFlags += 1024;
2001  }
2002  checkMCParticleMissingFlags(daughter, mcROEObjects, missingFlags);
2003  }
2004  }
2005 
2006  double isInThisRestOfEvent(const Particle* particle, const RestOfEvent* roe, const std::string& maskName)
2007  {
2008  if (particle->getParticleSource() == Particle::c_Composite or
2009  particle->getParticleSource() == Particle::c_V0) {
2010  std::vector<const Particle*> fspDaug = particle->getFinalStateDaughters();
2011  for (auto& i : fspDaug) {
2012  if (isInThisRestOfEvent(i, roe, maskName) == 0)
2013  return 0;
2014  }
2015  return 1.0;
2016  }
2017  return roe->hasParticle(particle, maskName);
2018  }
2019 
2020  const RestOfEvent* getRelatedROEObject(const Particle* particle, bool returnHostOnly)
2021  {
2022  // Get related ROE object
2023  const RestOfEvent* roe = particle->getRelatedTo<RestOfEvent>();
2024  if (!roe && !returnHostOnly) {
2025  roe = particle->getRelatedTo<RestOfEvent>("NestedRestOfEvents");
2026 
2027  }
2028  return roe;
2029  }
2030 
2031  VARIABLE_GROUP("Rest Of Event");
2032 
2033  REGISTER_VARIABLE("useROERecoilFrame(variable)", useROERecoilFrame,
2034  "Returns the value of the variable using the rest frame of the ROE recoil as current reference frame.\n"
2035  "Can be used inside for_each loop or outside of it if the particle has associated Rest of Event.\n"
2036  "E.g. ``useROERecoilFrame(E)`` returns the energy of a particle in the ROE recoil frame.");
2037 
2038  REGISTER_VARIABLE("isInRestOfEvent", isInRestOfEvent,
2039  "Returns 1 if a track, ecl or klmCluster associated to particle is in the current RestOfEvent object, 0 otherwise."
2040  "One can use this variable only in a for_each loop over the RestOfEvent StoreArray.");
2041 
2042  REGISTER_VARIABLE("isCloneOfSignalSide", isCloneOfSignalSide,
2043  "Returns 1 if a particle is a clone of signal side final state particles, 0 otherwise. "
2044  "Requires generator information and truth-matching. "
2045  "One can use this variable only in a ``for_each`` loop over the RestOfEvent StoreArray.");
2046 
2047  REGISTER_VARIABLE("hasAncestorFromSignalSide", hasAncestorFromSignalSide,
2048  "Returns 1 if a particle has ancestor from signal side, 0 otherwise. "
2049  "Requires generator information and truth-matching. "
2050  "One can use this variable only in a ``for_each`` loop over the RestOfEvent StoreArray.");
2051 
2052  REGISTER_VARIABLE("currentROEIsInList(particleList)", currentROEIsInList,
2053  "[Eventbased] Returns 1 the associated particle of the current ROE is contained in the given list or its charge-conjugated."
2054  "Useful to restrict the for_each loop over ROEs to ROEs of a certain ParticleList.");
2055 
2056  REGISTER_VARIABLE("nROE_RemainingTracks", nROE_RemainingTracks,
2057  "Returns number of tracks in ROE - number of tracks of given particle"
2058  "One can use this variable only in a for_each loop over the RestOfEvent StoreArray.");
2059 
2060  REGISTER_VARIABLE("nROE_RemainingTracks(maskName)", nROE_RemainingTracksWithMask,
2061  "Returns number of remaining tracks between the ROE (specified via a mask) and the given particle. For the given particle only tracks are counted which are in the RoE."
2062  "One can use this variable only in a for_each loop over the RestOfEvent StoreArray."
2063  "Is required for the specific FEI. :noindex:");
2064  // nROE_RemainingTracks is overloaded (two C++ functions sharing one
2065  // variable name) so one of the two needs to be made the indexed
2066  // variable in sphinx
2067 
2068  REGISTER_VARIABLE("nROE_KLMClusters", nROE_KLMClusters,
2069  "Returns number of all remaining KLM clusters in the related RestOfEvent object.");
2070 
2071  REGISTER_VARIABLE("nROE_Charged(maskName, PDGcode = 0)", nROE_ChargedParticles,
2072  "Returns number of all charged particles in the related RestOfEvent object. First optional argument is ROE mask name. "
2073  "Second argument is a PDG code to count only one charged particle species, independently of charge. "
2074  "For example: ``nROE_Charged(cleanMask, 321)`` will output number of kaons in Rest Of Event with ``cleanMask``. "
2075  "PDG code 0 is used to count all charged particles");
2076 
2077  REGISTER_VARIABLE("nROE_Photons(maskName)", nROE_Photons,
2078  "Returns number of all photons in the related RestOfEvent object, accepts 1 optional argument of ROE mask name. ");
2079 
2080  REGISTER_VARIABLE("nROE_NeutralHadrons(maskName)", nROE_NeutralHadrons,
2081  "Returns number of all neutral hadrons in the related RestOfEvent object, accepts 1 optional argument of ROE mask name. ");
2082 
2083  REGISTER_VARIABLE("particleRelatedToCurrentROE(var)", particleRelatedToCurrentROE,
2084  "[Eventbased] Returns variable applied to the particle which is related to the current RestOfEvent object"
2085  "One can use this variable only in a for_each loop over the RestOfEvent StoreArray.");
2086 
2087  REGISTER_VARIABLE("roeMC_E", ROE_MC_E,
2088  "Returns true energy of unused tracks and clusters in ROE, can be used with ``use***Frame()`` function.");
2089 
2090  REGISTER_VARIABLE("roeMC_M", ROE_MC_M,
2091  "Returns true invariant mass of unused tracks and clusters in ROE");
2092 
2093  REGISTER_VARIABLE("roeMC_P", ROE_MC_P,
2094  "Returns true momentum of unused tracks and clusters in ROE, can be used with ``use***Frame()`` function.");
2095 
2096  REGISTER_VARIABLE("roeMC_Px", ROE_MC_Px,
2097  "Returns x component of true momentum of unused tracks and clusters in ROE, can be used with ``use***Frame()`` function.");
2098 
2099  REGISTER_VARIABLE("roeMC_Py", ROE_MC_Py,
2100  "Returns y component of true momentum of unused tracks and clusters in ROE, can be used with ``use***Frame()`` function.");
2101 
2102  REGISTER_VARIABLE("roeMC_Pz", ROE_MC_Pz,
2103  "Returns z component of true momentum of unused tracks and clusters in ROE, can be used with ``use***Frame()`` function.");
2104 
2105  REGISTER_VARIABLE("roeMC_Pt", ROE_MC_Pt,
2106  "Returns transverse component of true momentum of unused tracks and clusters in ROE, can be used with ``use***Frame()`` function.");
2107 
2108  REGISTER_VARIABLE("roeMC_PTheta", ROE_MC_PTheta,
2109  "Returns polar angle of true momentum of unused tracks and clusters in ROE, can be used with ``use***Frame()`` function.");
2110 
2111  REGISTER_VARIABLE("roeMC_MissFlags(maskName)", ROE_MC_MissingFlags,
2112  "Returns flags corresponding to missing particles on ROE side.");
2113 
2114  REGISTER_VARIABLE("nROE_Tracks(maskName)", nROE_Tracks,
2115  "Returns number of tracks in the related RestOfEvent object that pass the selection criteria.");
2116 
2117  REGISTER_VARIABLE("nROE_ECLClusters(maskName)", nROE_ECLClusters,
2118  "Returns number of ECL clusters in the related RestOfEvent object that pass the selection criteria.");
2119 
2120  REGISTER_VARIABLE("nROE_NeutralECLClusters(maskName)", nROE_NeutralECLClusters,
2121  "Returns number of neutral ECL clusters in the related RestOfEvent object that pass the selection criteria.");
2122 
2123  REGISTER_VARIABLE("nROE_Composites(maskName)", nROE_Composites,
2124  "Returns number of composite particles or V0s in the related RestOfEvent object that pass the selection criteria.");
2125 
2126  REGISTER_VARIABLE("nROE_ParticlesInList(pListName)", nROE_ParticlesInList,
2127  "Returns the number of particles in ROE from the given particle list.\n"
2128  "Use of variable aliases is advised.");
2129 
2130  REGISTER_VARIABLE("roeCharge(maskName)", ROE_Charge,
2131  "Returns total charge of the related RestOfEvent object.");
2132 
2133  REGISTER_VARIABLE("roeEextra(maskName)", ROE_ExtraEnergy,
2134  "Returns extra energy from ECLClusters in the calorimeter that is not associated to the given Particle");
2135 
2136  REGISTER_VARIABLE("roeNeextra(maskName)", ROE_NeutralExtraEnergy,
2137  "Returns extra energy from neutral ECLClusters in the calorimeter that is not associated to the given Particle, can be used with ``use***Frame()`` function.");
2138 
2139  REGISTER_VARIABLE("roeE(maskName)", ROE_E,
2140  "Returns energy of unused tracks and clusters in ROE, can be used with ``use***Frame()`` function.");
2141 
2142  REGISTER_VARIABLE("roeM(maskName)", ROE_M,
2143  "Returns invariant mass of unused tracks and clusters in ROE");
2144 
2145  REGISTER_VARIABLE("roeP(maskName)", ROE_P,
2146  "Returns momentum of unused tracks and clusters in ROE, can be used with ``use***Frame()`` function.");
2147 
2148  REGISTER_VARIABLE("roePt(maskName)", ROE_Pt,
2149  "Returns transverse component of momentum of unused tracks and clusters in ROE, can be used with ``use***Frame()`` function.");
2150 
2151  REGISTER_VARIABLE("roePx(maskName)", ROE_Px,
2152  "Returns x component of momentum of unused tracks and clusters in ROE, can be used with ``use***Frame()`` function.");
2153 
2154  REGISTER_VARIABLE("roePy(maskName)", ROE_Py,
2155  "Returns y component of momentum of unused tracks and clusters in ROE, can be used with ``use***Frame()`` function.");
2156 
2157  REGISTER_VARIABLE("roePz(maskName)", ROE_Pz,
2158  "Returns z component of momentum of unused tracks and clusters in ROE, can be used with ``use***Frame()`` function.");
2159 
2160  REGISTER_VARIABLE("roePTheta(maskName)", ROE_PTheta,
2161  "Returns theta angle of momentum of unused tracks and clusters in ROE, can be used with ``use***Frame()`` function.");
2162 
2163  REGISTER_VARIABLE("roeDeltae(maskName)", ROE_DeltaE,
2164  "Returns energy difference of the related RestOfEvent object with respect to :math:`E_\\mathrm{cms}/2`.");
2165 
2166  REGISTER_VARIABLE("roeMbc(maskName)", ROE_Mbc,
2167  "Returns beam constrained mass of the related RestOfEvent object with respect to :math:`E_\\mathrm{cms}/2`.");
2168 
2169  REGISTER_VARIABLE("weDeltae(maskName, opt)", WE_DeltaE,
2170  "Returns the energy difference of the B meson, corrected with the missing neutrino momentum (reconstructed side + neutrino) with respect to :math:`E_\\mathrm{cms}/2`.");
2171 
2172  REGISTER_VARIABLE("weMbc(maskName, opt)", WE_Mbc,
2173  "Returns beam constrained mass of B meson, corrected with the missing neutrino momentum (reconstructed side + neutrino) with respect to :math:`E_\\mathrm{cms}/2`.");
2174 
2175  REGISTER_VARIABLE("weMissM2(maskName, opt)", WE_MissM2,
2176  "Returns the invariant mass squared of the missing momentum (see :b2:var:`weMissE` possible options)");
2177 
2178  REGISTER_VARIABLE("weMissPTheta(maskName, opt)", WE_MissPTheta,
2179  "Returns the polar angle of the missing momentum (see possible :b2:var:`weMissE` options)");
2180 
2181  REGISTER_VARIABLE("weMissP(maskName, opt)", WE_MissP,
2182  "Returns the magnitude of the missing momentum (see possible :b2:var:`weMissE` options)");
2183 
2184  REGISTER_VARIABLE("weMissPx(maskName, opt)", WE_MissPx,
2185  "Returns the x component of the missing momentum (see :b2:var:`weMissE` possible options)");
2186 
2187  REGISTER_VARIABLE("weMissPy(maskName, opt)", WE_MissPy,
2188  "Returns the y component of the missing momentum (see :b2:var:`weMissE` possible options)");
2189 
2190  REGISTER_VARIABLE("weMissPz(maskName, opt)", WE_MissPz,
2191  "Returns the z component of the missing momentum (see :b2:var:`weMissE` possible options)");
2192 
2193  REGISTER_VARIABLE("weMissE(maskName, opt)", WE_MissE,
2194  R"DOC(Returns the energy of the missing momentum, possible options ``opt`` are the following:
2195 
2196 - ``0``: CMS, use energy and momentum of charged particles and photons
2197 - ``1``: CMS, same as ``0``, fix :math:`E_\mathrm{miss} = p_\mathrm{miss}`
2198 - ``2``: CMS, same as ``0``, fix :math:`E_\mathrm{roe} = E_\mathrm{cms}/2`
2199 - ``3``: CMS, use only energy and momentum of signal side
2200 - ``4``: CMS, same as ``3``, update with direction of ROE momentum
2201 - ``5``: LAB, use energy and momentum of charged particles and photons from whole event
2202 - ``6``: LAB, same as ``5``, fix :math:`E_\mathrm{miss} = p_\mathrm{miss}``
2203 - ``7``: CMS, correct pmiss 3-momentum vector with factor alpha so that :math:`d_E = 0`` (used for :math:`M_\mathrm{bc}` calculation).)DOC");
2204 
2205  REGISTER_VARIABLE("weXiZ(maskName)", WE_xiZ,
2206  "Returns Xi_z in event (for Bhabha suppression and two-photon scattering)");
2207 
2208  REGISTER_VARIABLE("bssMassDifference(maskName)", bssMassDifference,
2209  "Bs* - Bs mass difference");
2210 
2211  REGISTER_VARIABLE("weCosThetaEll(maskName)", WE_cosThetaEll, R"DOC(
2212 
2213 Returns the angle between :math:`M` and lepton in :math:`W` rest frame in the decays of the type:
2214 :math:`M \to h_1 ... h_n \ell`, where W 4-momentum is given as
2215 
2216 .. math::
2217  p_W = p_\ell + p_\nu.
2218 
2219 The neutrino momentum is calculated from ROE taking into account the specified mask, and setting
2220 
2221 .. math::
2222  E_{\nu} = |p_{miss}|.
2223 
2224 )DOC");
2225 
2226  REGISTER_VARIABLE("weQ2lnuSimple(maskName,option)", WE_q2lnuSimple,
2227  "Returns the momentum transfer squared, :math:`q^2`, calculated in CMS as :math:`q^2 = (p_l + p_\\nu)^2`, \n"
2228  "where :math:`B \\to H_1\\dots H_n \\ell \\nu_\\ell`. Lepton is assumed to be the last reconstructed daughter. \n"
2229  "By default, option is set to ``1`` (see :b2:var:`weMissE`). Unless you know what you are doing, keep this default value.");
2230 
2231  REGISTER_VARIABLE("weQ2lnu(maskName,option)", WE_q2lnu,
2232  "Returns the momentum transfer squared, :math:`q^2`, calculated in CMS as :math:`q^2 = (p_l + p_\\nu)^2`, \n"
2233  "where :math:`B \\to H_1\\dots H_n \\ell \\nu_\\ell`. Lepton is assumed to be the last reconstructed daughter. \n"
2234  "This calculation uses constraints from dE = 0 and Mbc = Mb to correct the neutrino direction. \n"
2235  "By default, option is set to ``7`` (see :b2:var:`weMissE`). Unless you know what you are doing, keep this default value.");
2236 
2237  REGISTER_VARIABLE("weMissM2OverMissE(maskName)", WE_MissM2OverMissE,
2238  "Returns missing mass squared over missing energy");
2239 
2240  REGISTER_VARIABLE("passesROEMask(maskName)", passesROEMask,
2241  "Returns boolean value if a particle passes a certain mask or not. Only to be used in for_each path, otherwise returns quiet NaN.");
2242 
2243  REGISTER_VARIABLE("printROE", printROE,
2244  "For debugging, prints indices of all particles in the ROE and all masks. Returns 0.");
2245 
2246  REGISTER_VARIABLE("pi0Prob(mode)", pi0Prob,
2247  "Returns pi0 probability, where mode is used to specify the selection criteria for soft photon. \n"
2248  "The following strings are available. \n\n"
2249  "- ``standard``: loose energy cut and no clusterNHits cut are applied to soft photon \n"
2250  "- ``tight``: tight energy cut and no clusterNHits cut are applied to soft photon \n"
2251  "- ``cluster``: loose energy cut and clusterNHits cut are applied to soft photon \n"
2252  "- ``both``: tight energy cut and clusterNHits cut are applied to soft photon \n\n"
2253  "You can find more details in `writePi0EtaVeto` function in modularAnalysis.py.");
2254 
2255  REGISTER_VARIABLE("etaProb(mode)", etaProb,
2256  "Returns eta probability, where mode is used to specify the selection criteria for soft photon. \n"
2257  "The following strings are available. \n\n"
2258  "- ``standard``: loose energy cut and no clusterNHits cut are applied to soft photon \n"
2259  "- ``tight``: tight energy cut and no clusterNHits cut are applied to soft photon \n"
2260  "- ``cluster``: loose energy cut and clusterNHits cut are applied to soft photon \n"
2261  "- ``both``: tight energy cut and clusterNHits cut are applied to soft photon \n\n"
2262  "You can find more details in `writePi0EtaVeto` function in modularAnalysis.py.");
2263 
2264  }
2266 }
Belle2
Abstract base class for different kinds of events.
Definition: MillepedeAlgorithm.h:17