Belle II Software development
EventShapeVariables.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/EventShapeVariables.h>
11
12#include <analysis/dataobjects/Particle.h>
13#include <analysis/dataobjects/EventShapeContainer.h>
14
15#include <analysis/utility/ReferenceFrame.h>
16
17#include <framework/logging/Logger.h>
18#include <framework/datastore/StoreObjPtr.h>
19#include <framework/utilities/Conversion.h>
20
21#include <TVectorF.h>
22
23#include <boost/algorithm/string.hpp>
24
25namespace Belle2 {
30 namespace Variable {
31
32 double foxWolframR(const Particle*, const std::vector<double>& index)
33 {
34 if (index.size() != 1)
35 B2FATAL("foxWolframR cannot be called without providing the moment order");
36
37 int order = std::lround(index[0]);
38
39 if (order < 0 || order > 8)
40 B2FATAL("The Fox-Wolfram moment order must be within 0 and 8.");
41
42 StoreObjPtr<EventShapeContainer> evtShapeCont;
43 if (!evtShapeCont) {
44 B2ERROR("No EventShapeContainer object has been found in the datastore");
45 return Const::doubleNaN;
46 }
47 if (evtShapeCont->getFWMoment(0) == 0) {
48 B2INFO("The 0th-order FoxWolfram moment is zero");
49 return Const::doubleNaN;
50 }
51 return evtShapeCont->getFWMoment(order) / evtShapeCont->getFWMoment(0);
52 }
53
54 double foxWolframH(const Particle*, const std::vector<double>& index)
55 {
56 if (index.size() != 1)
57 B2FATAL("foxWolframH cannot be called without providing the moment order");
58
59 int order = std::lround(index[0]);
60
61 if (order < 0 || order > 8)
62 B2FATAL("The Fox-Wolfram moment order must be within 0 and 8.");
63
64 StoreObjPtr<EventShapeContainer> evtShapeCont;
65 if (!evtShapeCont) {
66 B2ERROR("No EventShapeContainer object has been found in the datastore");
67 return Const::doubleNaN;
68 }
69 return evtShapeCont->getFWMoment(order);
70 }
71
72 Manager::FunctionPtr harmonicMoment(const std::vector<std::string>& arguments)
73 {
74 if (arguments.size() != 2)
75 B2FATAL("harmonicMoment requires two arguments: the harmonic order (0-8) and the reference axis name (thrust or collision)");
76
77 int order = -1;
78 try {
79 order = Belle2::convertString<int>(arguments[0]);
80 } catch (std::invalid_argument&) {
81 B2FATAL("First argument of harmonicMoment must be an integer");
82 }
83 std::string axisName = arguments[1];
84 boost::to_lower(axisName);
85
86 if (order < 0 || order > 8)
87 B2FATAL("The Fox-Wolfram moment order must be within 0 and 8.");
88
89 if (axisName != "thrust" && axisName != "collision")
90 B2FATAL("Invalid axis name " << arguments[1] << ". The valid options are thrust and collision");
91
92 auto func = [order, axisName](const Particle*) -> double{
93 StoreObjPtr<EventShapeContainer> evtShapeCont;
94 if (!evtShapeCont)
95 {
96 B2ERROR("No EventShapeContainer object has been found in the datastore");
97 return Const::doubleNaN;
98 }
99 if (axisName == "thrust")
100 return evtShapeCont->getHarmonicMomentThrust(order);
101 else
102 return evtShapeCont->getHarmonicMomentCollision(order);
103 };
104 return func;
105 }
106
107 Manager::FunctionPtr cleoCone(const std::vector<std::string>& arguments)
108 {
109 if (arguments.size() != 2)
110 B2FATAL("cleoCone requires two arguments: the cone order (0-9) and the reference axis name (thrust or collision)");
111
112 int order = -1;
113 try {
114 order = Belle2::convertString<int>(arguments[0]);
115 } catch (std::invalid_argument&) {
116 B2FATAL("Argument of cleoCone must be an integer");
117 }
118 std::string axisName = arguments[1];
119 boost::to_lower(axisName);
120
121 if (order < 0 || order > 8)
122 B2FATAL("The CLEO cone order must be within 0 and 8.");
123
124 if (axisName != "thrust" && axisName != "collision")
125 B2FATAL("Invalid axis name " << arguments[1] << ". The valid options are thrust and collision");
126
127 auto func = [order, axisName](const Particle*) -> double{
128 StoreObjPtr<EventShapeContainer> evtShapeCont;
129 if (!evtShapeCont)
130 {
131 B2ERROR("No EventShapeContainer object has been found in the datastore");
132 return Const::doubleNaN;
133 }
134 if (axisName == "thrust")
135 return evtShapeCont->getCleoConeThrust(order);
136 else
137 return evtShapeCont->getCleoConeCollision(order);
138 };
139 return func;
140 }
141
142 double foxWolframR1(const Particle*)
143 {
144 StoreObjPtr<EventShapeContainer> evtShapeCont;
145 if (!evtShapeCont) {
146 B2ERROR("No EventShapeContainer object has been found in the datastore");
147 return Const::doubleNaN;
148 }
149 if (evtShapeCont->getFWMoment(0) == 0) {
150 B2INFO("The 0th-order FoxWolfram moment is zero");
151 return Const::doubleNaN;
152 }
153 return evtShapeCont->getFWMoment(1) / evtShapeCont->getFWMoment(0);
154 }
155
156 double foxWolframR2(const Particle*)
157 {
158 StoreObjPtr<EventShapeContainer> evtShapeCont;
159 if (!evtShapeCont) {
160 B2ERROR("No EventShapeContainer object has been found in the datastore");
161 return Const::doubleNaN;
162 }
163 if (evtShapeCont->getFWMoment(0) == 0) {
164 B2INFO("The 0th-order FoxWolfram moment is zero");
165 return Const::doubleNaN;
166 }
167 return evtShapeCont->getFWMoment(2) / evtShapeCont->getFWMoment(0);
168 }
169
170 double foxWolframR3(const Particle*)
171 {
172 StoreObjPtr<EventShapeContainer> evtShapeCont;
173 if (!evtShapeCont) {
174 B2ERROR("No EventShapeContainer object has been found in the datastore");
175 return Const::doubleNaN;
176 }
177 if (evtShapeCont->getFWMoment(0) == 0) {
178 B2INFO("The 0th-order FoxWolfram moment is zero");
179 return Const::doubleNaN;
180 }
181 return evtShapeCont->getFWMoment(3) / evtShapeCont->getFWMoment(0);
182 }
183
184 double foxWolframR4(const Particle*)
185 {
186 StoreObjPtr<EventShapeContainer> evtShapeCont;
187 if (!evtShapeCont) {
188 B2ERROR("No EventShapeContainer object has been found in the datastore");
189 return Const::doubleNaN;
190 }
191 if (evtShapeCont->getFWMoment(0) == 0) {
192 B2INFO("The 0th-order FoxWolfram moment is zero");
193 return Const::doubleNaN;
194 }
195 return evtShapeCont->getFWMoment(4) / evtShapeCont->getFWMoment(0);
196 }
197
198 double harmonicMomentThrust0(const Particle*)
199 {
200 StoreObjPtr<EventShapeContainer> evtShapeCont;
201 if (!evtShapeCont) {
202 B2ERROR("No EventShapeContainer object has been found in the datastore");
203 return Const::doubleNaN;
204 }
205 return evtShapeCont->getHarmonicMomentThrust(0);
206 }
207
208 double harmonicMomentThrust1(const Particle*)
209 {
210 StoreObjPtr<EventShapeContainer> evtShapeCont;
211 if (!evtShapeCont) {
212 B2ERROR("No EventShapeContainer object has been found in the datastore");
213 return Const::doubleNaN;
214 }
215 return evtShapeCont->getHarmonicMomentThrust(1);
216 }
217
218 double harmonicMomentThrust2(const Particle*)
219 {
220 StoreObjPtr<EventShapeContainer> evtShapeCont;
221 if (!evtShapeCont) {
222 B2ERROR("No EventShapeContainer object has been found in the datastore");
223 return Const::doubleNaN;
224 }
225 return evtShapeCont->getHarmonicMomentThrust(2);
226 }
227
228 double harmonicMomentThrust3(const Particle*)
229 {
230 StoreObjPtr<EventShapeContainer> evtShapeCont;
231 if (!evtShapeCont) {
232 B2ERROR("No EventShapeContainer object has been found in the datastore");
233 return Const::doubleNaN;
234 }
235 return evtShapeCont->getHarmonicMomentThrust(3);
236 }
237
238 double harmonicMomentThrust4(const Particle*)
239 {
240 StoreObjPtr<EventShapeContainer> evtShapeCont;
241 if (!evtShapeCont) {
242 B2ERROR("No EventShapeContainer object has been found in the datastore");
243 return Const::doubleNaN;
244 }
245 return evtShapeCont->getHarmonicMomentThrust(4);
246 }
247
248 double cleoConeThrust0(const Particle*)
249 {
250 StoreObjPtr<EventShapeContainer> evtShapeCont;
251 if (!evtShapeCont) {
252 B2ERROR("No EventShapeContainer object has been found in the datastore");
253 return Const::doubleNaN;
254 }
255 return evtShapeCont->getCleoConeThrust(0);
256 }
257
258 double cleoConeThrust1(const Particle*)
259 {
260 StoreObjPtr<EventShapeContainer> evtShapeCont;
261 if (!evtShapeCont) {
262 B2ERROR("No EventShapeContainer object has been found in the datastore");
263 return Const::doubleNaN;
264 }
265 return evtShapeCont->getCleoConeThrust(1);
266 }
267
268 double cleoConeThrust2(const Particle*)
269 {
270 StoreObjPtr<EventShapeContainer> evtShapeCont;
271 if (!evtShapeCont) {
272 B2ERROR("No EventShapeContainer object has been found in the datastore");
273 return Const::doubleNaN;
274 }
275 return evtShapeCont->getCleoConeThrust(2);
276 }
277
278 double cleoConeThrust3(const Particle*)
279 {
280 StoreObjPtr<EventShapeContainer> evtShapeCont;
281 if (!evtShapeCont) {
282 B2ERROR("No EventShapeContainer object has been found in the datastore");
283 return Const::doubleNaN;
284 }
285 return evtShapeCont->getCleoConeThrust(3);
286 }
287
288 double cleoConeThrust4(const Particle*)
289 {
290 StoreObjPtr<EventShapeContainer> evtShapeCont;
291 if (!evtShapeCont) {
292 B2ERROR("No EventShapeContainer object has been found in the datastore");
293 return Const::doubleNaN;
294 }
295 return evtShapeCont->getCleoConeThrust(4);
296 }
297
298 double cleoConeThrust5(const Particle*)
299 {
300 StoreObjPtr<EventShapeContainer> evtShapeCont;
301 if (!evtShapeCont) {
302 B2ERROR("No EventShapeContainer object has been found in the datastore");
303 return Const::doubleNaN;
304 }
305 return evtShapeCont->getCleoConeThrust(5);
306 }
307
308 double cleoConeThrust6(const Particle*)
309 {
310 StoreObjPtr<EventShapeContainer> evtShapeCont;
311 if (!evtShapeCont) {
312 B2ERROR("No EventShapeContainer object has been found in the datastore");
313 return Const::doubleNaN;
314 }
315 return evtShapeCont->getCleoConeThrust(6);
316 }
317
318 double cleoConeThrust7(const Particle*)
319 {
320 StoreObjPtr<EventShapeContainer> evtShapeCont;
321 if (!evtShapeCont) {
322 B2ERROR("No EventShapeContainer object has been found in the datastore");
323 return Const::doubleNaN;
324 }
325 return evtShapeCont->getCleoConeThrust(7);
326 }
327
328 double cleoConeThrust8(const Particle*)
329 {
330 StoreObjPtr<EventShapeContainer> evtShapeCont;
331 if (!evtShapeCont) {
332 B2ERROR("No EventShapeContainer object has been found in the datastore");
333 return Const::doubleNaN;
334 }
335 return evtShapeCont->getCleoConeThrust(8);
336 }
337
338 double sphericity(const Particle*)
339 {
340 StoreObjPtr<EventShapeContainer> evtShapeCont;
341 if (!evtShapeCont) {
342 B2ERROR("No EventShapeContainer object has been found in the datastore");
343 return Const::doubleNaN;
344 }
345 // (3/2)(lamda_2 + lambda_3)
346 return 1.5 * (evtShapeCont->getSphericityEigenvalue(1) + evtShapeCont->getSphericityEigenvalue(2)) ;
347 }
348
349 double aplanarity(const Particle*)
350 {
351 StoreObjPtr<EventShapeContainer> evtShapeCont;
352 if (!evtShapeCont) {
353 B2ERROR("No EventShapeContainer object has been found in the datastore");
354 return Const::doubleNaN;
355 }
356 // (3/2)(lambda_3)
357 return 1.5 * evtShapeCont->getSphericityEigenvalue(2);
358 }
359
360 double forwardHemisphereMass(const Particle*)
361 {
362 StoreObjPtr<EventShapeContainer> evtShapeCont;
363 if (!evtShapeCont) {
364 B2ERROR("No EventShapeContainer object has been found in the datastore");
365 return Const::doubleNaN;
366 }
367 return evtShapeCont->getForwardHemisphere4Momentum().M();
368 }
369
370 double forwardHemisphereX(const Particle*)
371 {
372 StoreObjPtr<EventShapeContainer> evtShapeCont;
373 if (!evtShapeCont) {
374 B2ERROR("No EventShapeContainer object has been found in the datastore");
375 return Const::doubleNaN;
376 }
377 return evtShapeCont->getForwardHemisphere4Momentum().px();
378 }
379
380 double forwardHemisphereY(const Particle*)
381 {
382 StoreObjPtr<EventShapeContainer> evtShapeCont;
383 if (!evtShapeCont) {
384 B2ERROR("No EventShapeContainer object has been found in the datastore");
385 return Const::doubleNaN;
386 }
387 return evtShapeCont->getForwardHemisphere4Momentum().py();
388 }
389
390 double forwardHemisphereZ(const Particle*)
391 {
392 StoreObjPtr<EventShapeContainer> evtShapeCont;
393 if (!evtShapeCont) {
394 B2ERROR("No EventShapeContainer object has been found in the datastore");
395 return Const::doubleNaN;
396 }
397 return evtShapeCont->getForwardHemisphere4Momentum().pz();
398 }
399
400 double forwardHemisphereMomentum(const Particle*)
401 {
402 StoreObjPtr<EventShapeContainer> evtShapeCont;
403 if (!evtShapeCont) {
404 B2ERROR("No EventShapeContainer object has been found in the datastore");
405 return Const::doubleNaN;
406 }
407 return evtShapeCont->getForwardHemisphere4Momentum().P();
408 }
409
410 double forwardHemisphereEnergy(const Particle*)
411 {
412 StoreObjPtr<EventShapeContainer> evtShapeCont;
413 if (!evtShapeCont) {
414 B2ERROR("No EventShapeContainer object has been found in the datastore");
415 return Const::doubleNaN;
416 }
417 return evtShapeCont->getForwardHemisphere4Momentum().E();
418 }
419
420 double backwardHemisphereMass(const Particle*)
421 {
422 StoreObjPtr<EventShapeContainer> evtShapeCont;
423 if (!evtShapeCont) {
424 B2ERROR("No EventShapeContainer object has been found in the datastore");
425 return Const::doubleNaN;
426 }
427 return evtShapeCont->getBackwardHemisphere4Momentum().M();
428 }
429
430 double backwardHemisphereX(const Particle*)
431 {
432 StoreObjPtr<EventShapeContainer> evtShapeCont;
433 if (!evtShapeCont) {
434 B2ERROR("No EventShapeContainer object has been found in the datastore");
435 return Const::doubleNaN;
436 }
437 return evtShapeCont->getBackwardHemisphere4Momentum().px();
438 }
439
440 double backwardHemisphereY(const Particle*)
441 {
442 StoreObjPtr<EventShapeContainer> evtShapeCont;
443 if (!evtShapeCont) {
444 B2ERROR("No EventShapeContainer object has been found in the datastore");
445 return Const::doubleNaN;
446 }
447 return evtShapeCont->getBackwardHemisphere4Momentum().py();
448 }
449
450 double backwardHemisphereZ(const Particle*)
451 {
452 StoreObjPtr<EventShapeContainer> evtShapeCont;
453 if (!evtShapeCont) {
454 B2ERROR("No EventShapeContainer object has been found in the datastore");
455 return Const::doubleNaN;
456 }
457 return evtShapeCont->getBackwardHemisphere4Momentum().pz();
458 }
459
460 double backwardHemisphereMomentum(const Particle*)
461 {
462 StoreObjPtr<EventShapeContainer> evtShapeCont;
463 if (!evtShapeCont) {
464 B2ERROR("No EventShapeContainer object has been found in the datastore");
465 return Const::doubleNaN;
466 }
467 return evtShapeCont->getBackwardHemisphere4Momentum().P();
468 }
469
470 double backwardHemisphereEnergy(const Particle*)
471 {
472 StoreObjPtr<EventShapeContainer> evtShapeCont;
473 if (!evtShapeCont) {
474 B2ERROR("No EventShapeContainer object has been found in the datastore");
475 return Const::doubleNaN;
476 }
477 return evtShapeCont->getBackwardHemisphere4Momentum().E();
478 }
479
480 double thrust(const Particle*)
481 {
482 StoreObjPtr<EventShapeContainer> evtShapeCont;
483 if (!evtShapeCont) {
484 B2ERROR("No EventShapeContainer object has been found in the datastore");
485 return Const::doubleNaN;
486 }
487 return evtShapeCont->getThrust();
488 }
489
490 double thrustAxisX(const Particle*)
491 {
492 StoreObjPtr<EventShapeContainer> evtShapeCont;
493 if (!evtShapeCont) {
494 B2ERROR("No EventShapeContainer object has been found in the datastore");
495 return Const::doubleNaN;
496 }
497 return evtShapeCont->getThrustAxis().X();
498 }
499
500 double thrustAxisY(const Particle*)
501 {
502 StoreObjPtr<EventShapeContainer> evtShapeCont;
503 if (!evtShapeCont) {
504 B2ERROR("No EventShapeContainer object has been found in the datastore");
505 return Const::doubleNaN;
506 }
507 return evtShapeCont->getThrustAxis().Y();
508 }
509
510 double thrustAxisZ(const Particle*)
511 {
512 StoreObjPtr<EventShapeContainer> evtShapeCont;
513 if (!evtShapeCont) {
514 B2ERROR("No EventShapeContainer object has been found in the datastore");
515 return Const::doubleNaN;
516 }
517 return evtShapeCont->getThrustAxis().Z();
518 }
519
520 double thrustAxisCosTheta(const Particle*)
521 {
522 StoreObjPtr<EventShapeContainer> evtShapeCont;
523 if (!evtShapeCont) {
524 B2ERROR("No EventShapeContainer object has been found in the datastore");
525 return Const::doubleNaN;
526 }
527 return cos(evtShapeCont->getThrustAxis().Theta());
528 }
529
530 Manager::FunctionPtr useThrustFrame(const std::vector<std::string>& arguments)
531 {
532 if (arguments.size() == 1) {
533 auto variableName = arguments[0];
534
535 const Variable::Manager::Var* var = Manager::Instance().getVariable(variableName);
536
537 auto func = [var](const Particle * particle) -> double {
538 StoreObjPtr<EventShapeContainer> evtShapeCont;
539 if (!evtShapeCont)
540 {
541 B2ERROR("No EventShapeContainer object has been found in the datastore");
542 return Const::doubleNaN;
543 }
544
545 ROOT::Math::XYZVector newZ = evtShapeCont->getThrustAxis();
546 ROOT::Math::XYZVector newY(0, 0, 0);
547 if (newZ.Z() == 0 and newZ.Y() == 0)
548 newY.SetX(1);
549 else
550 {
551 newY.SetY(newZ.Z());
552 newY.SetZ(-newZ.Y());
553 }
554 ROOT::Math::XYZVector newX = newY.Cross(newZ);
555
556 UseReferenceFrame<CMSRotationFrame> signalframe(newX, newY, newZ);
557
558 return std::get<double>(var->function(particle));
559 };
560 return func;
561 } else {
562 B2FATAL("Wrong number of arguments for meta function useThrustFrame. It only takes one argument, the variable name.");
563 }
564 }
565
566 VARIABLE_GROUP("EventShape");
567
568 REGISTER_VARIABLE("foxWolframR(i)", foxWolframR, R"DOC(
569[Eventbased] Ratio of the i-th to the 0-th order Fox Wolfram moments. The order ``i`` can go from 0 up to 8th.
570
571.. warning:: You have to run the Event Shape builder module for this variable to be meaningful.
572.. seealso:: :ref:`analysis_eventshape` and `modularAnalysis.buildEventShape`.
573)DOC");
574 REGISTER_VARIABLE("foxWolframH(i)", foxWolframH, R"DOC(
575[Eventbased] Returns i-th order Fox Wolfram moment. The order ``i`` can go from 0 up to 8th."
576
577.. warning:: You have to run the Event Shape builder module for this variable to be meaningful.
578.. seealso:: :ref:`analysis_eventshape` and `modularAnalysis.buildEventShape`.
579
580)DOC",":math:`\\text{GeV}^2`/:math:`\\text{c}^2`");
581 REGISTER_METAVARIABLE("harmonicMoment(i, axisName)", harmonicMoment, R"DOC(
582[Eventbased] Returns i-th order harmonic moment, calculated with respect to the axis ``axisName``.
583The order ``i`` can go from 0 up to 8th, the ``axisName`` can be either 'thrust' or 'collision'.
584The unit of the harmonic moment is ``GeV/c``.
585
586.. warning:: You have to run the Event Shape builder module for this variable to be meaningful.
587.. seealso:: :ref:`analysis_eventshape` and `modularAnalysis.buildEventShape`.
588)DOC", Manager::VariableDataType::c_double);
589 REGISTER_METAVARIABLE("cleoCone(i, axisName)", cleoCone, R"DOC(
590[Eventbased] Returns i-th order Cleo cone, calculated with respect to the axis ``axisName``. The order ``i`` can go from 0 up to 8th, the ``axisName`` can be either 'thrust' or 'collision'.
591
592.. warning:: You have to run the Event Shape builder module for this variable to be meaningful.
593.. seealso:: :ref:`analysis_eventshape` and `modularAnalysis.buildEventShape`.
594)DOC", Manager::VariableDataType::c_double);
595 REGISTER_METAVARIABLE("useThrustFrame(variable)", useThrustFrame, R"DOC(
596Evaluates a variable value in the thrust reference frame.
597
598.. warning:: You have to run the Event Shape builder module for this variable to be meaningful.
599.. seealso:: :ref:`analysis_eventshape` and `modularAnalysis.buildEventShape`.
600)DOC", Manager::VariableDataType::c_double);
601
602
603 REGISTER_VARIABLE("foxWolframR1", foxWolframR1, R"DOC(
604[Eventbased] Ratio of the 1-st to the 0-th order Fox Wolfram moments. This is just an alias of foxWolframR(1) defined for the user's convenience.
605
606.. warning:: You have to run the Event Shape builder module for this variable to be meaningful.
607.. seealso:: :ref:`analysis_eventshape` and `modularAnalysis.buildEventShape`.
608)DOC");
609 REGISTER_VARIABLE("foxWolframR2", foxWolframR2, R"DOC(
610[Eventbased] Ratio of the 2-nd to the 0-th order Fox Wolfram moments. This is just an alias of foxWolframR(2) defined for the user's convenience.
611
612.. warning:: You have to run the Event Shape builder module for this variable to be meaningful.
613.. seealso:: :ref:`analysis_eventshape` and `modularAnalysis.buildEventShape`.
614)DOC");
615 REGISTER_VARIABLE("foxWolframR3", foxWolframR3, R"DOC(
616[Eventbased] Ratio of the 3-rd to the 0-th order Fox Wolfram moments. This is just an alias of foxWolframR(3) defined for the user's convenience.
617
618.. warning:: You have to run the Event Shape builder module for this variable to be meaningful.
619.. seealso:: :ref:`analysis_eventshape` and `modularAnalysis.buildEventShape`.
620)DOC");
621 REGISTER_VARIABLE("foxWolframR4", foxWolframR4, R"DOC(
622[Eventbased] Ratio of the 4-th to the 0-th order Fox Wolfram moments. This is just an alias of foxWolframR(4) defined for the user's convenience.
623
624.. warning:: You have to run the Event Shape builder module for this variable to be meaningful.
625.. seealso:: :ref:`analysis_eventshape` and `modularAnalysis.buildEventShape`.
626)DOC");
627
628 REGISTER_VARIABLE("harmonicMomentThrust0", harmonicMomentThrust0, R"DOC(
629[Eventbased] Harmonic moment of the 0th order calculated with respect to the thrust axis.
630
631.. warning:: You have to run the Event Shape builder module for this variable to be meaningful.
632.. seealso:: :ref:`analysis_eventshape` and `modularAnalysis.buildEventShape`.
633
634)DOC","GeV/c");
635 REGISTER_VARIABLE("harmonicMomentThrust1", harmonicMomentThrust1, R"DOC(
636[Eventbased] Harmonic moment of the 1st order calculated with respect to the thrust axis.
637
638.. warning:: You have to run the Event Shape builder module for this variable to be meaningful.
639.. seealso:: :ref:`analysis_eventshape` and `modularAnalysis.buildEventShape`.
640
641)DOC","GeV/c");
642 REGISTER_VARIABLE("harmonicMomentThrust2", harmonicMomentThrust2, R"DOC(
643[Eventbased] Harmonic moment of the 2nd order calculated with respect to the thrust axis.
644
645.. warning:: You have to run the Event Shape builder module for this variable to be meaningful.
646.. seealso:: :ref:`analysis_eventshape` and `modularAnalysis.buildEventShape`.
647
648)DOC","GeV/c");
649 REGISTER_VARIABLE("harmonicMomentThrust3", harmonicMomentThrust3, R"DOC(
650[Eventbased] Harmonic moment of the 3rd order calculated with respect to the thrust axis.
651
652.. warning:: You have to run the Event Shape builder module for this variable to be meaningful.
653.. seealso:: :ref:`analysis_eventshape` and `modularAnalysis.buildEventShape`.
654
655)DOC","GeV/c");
656 REGISTER_VARIABLE("harmonicMomentThrust4", harmonicMomentThrust4, R"DOC(
657[Eventbased] Harmonic moment of the 4th order calculated with respect to the thrust axis.
658
659.. warning:: You have to run the Event Shape builder module for this variable to be meaningful.
660.. seealso:: :ref:`analysis_eventshape` and `modularAnalysis.buildEventShape`.
661
662)DOC","GeV/c");
663
664 REGISTER_VARIABLE("cleoConeThrust0", cleoConeThrust0, R"DOC(
665[Eventbased] 0th Cleo cone calculated with respect to the thrust axis.
666
667.. warning:: You have to run the Event Shape builder module for this variable to be meaningful.
668.. seealso:: :ref:`analysis_eventshape` and `modularAnalysis.buildEventShape`.
669)DOC");
670 REGISTER_VARIABLE("cleoConeThrust1", cleoConeThrust1, R"DOC(
671[Eventbased] 1st Cleo cone calculated with respect to the thrust axis.
672
673.. warning:: You have to run the Event Shape builder module for this variable to be meaningful.
674.. seealso:: :ref:`analysis_eventshape` and `modularAnalysis.buildEventShape`.
675)DOC");
676 REGISTER_VARIABLE("cleoConeThrust2", cleoConeThrust2, R"DOC(
677[Eventbased] 2nd Cleo cone calculated with respect to the thrust axis.
678
679.. warning:: You have to run the Event Shape builder module for this variable to be meaningful.
680.. seealso:: :ref:`analysis_eventshape` and `modularAnalysis.buildEventShape`.
681)DOC");
682 REGISTER_VARIABLE("cleoConeThrust3", cleoConeThrust3, R"DOC(
683[Eventbased] 3rd Cleo cone calculated with respect to the thrust axis.
684
685.. warning:: You have to run the Event Shape builder module for this variable to be meaningful.
686.. seealso:: :ref:`analysis_eventshape` and `modularAnalysis.buildEventShape`.
687)DOC");
688 REGISTER_VARIABLE("cleoConeThrust4", cleoConeThrust4, R"DOC(
689[Eventbased] 4th Cleo cone calculated with respect to the thrust axis.
690
691.. warning:: You have to run the Event Shape builder module for this variable to be meaningful.
692.. seealso:: :ref:`analysis_eventshape` and `modularAnalysis.buildEventShape`.
693)DOC");
694 REGISTER_VARIABLE("cleoConeThrust5", cleoConeThrust5, R"DOC(
695[Eventbased] 5th Cleo cone calculated with respect to the thrust axis.
696
697.. warning:: You have to run the Event Shape builder module for this variable to be meaningful.
698.. seealso:: :ref:`analysis_eventshape` and `modularAnalysis.buildEventShape`.
699)DOC");
700 REGISTER_VARIABLE("cleoConeThrust6", cleoConeThrust6, R"DOC(
701[Eventbased] 6th Cleo cone calculated with respect to the thrust axis.
702
703.. warning:: You have to run the Event Shape builder module for this variable to be meaningful.
704.. seealso:: :ref:`analysis_eventshape` and `modularAnalysis.buildEventShape`.
705)DOC");
706 REGISTER_VARIABLE("cleoConeThrust7", cleoConeThrust7, R"DOC(
707[Eventbased] 7th Cleo cone calculated with respect to the thrust axis.
708
709.. warning:: You have to run the Event Shape builder module for this variable to be meaningful.
710.. seealso:: :ref:`analysis_eventshape` and `modularAnalysis.buildEventShape`.
711)DOC");
712 REGISTER_VARIABLE("cleoConeThrust8", cleoConeThrust8, R"DOC(
713[Eventbased] 8th Cleo cone calculated with respect to the thrust axis.
714
715.. warning:: You have to run the Event Shape builder module for this variable to be meaningful.
716.. seealso:: :ref:`analysis_eventshape` and `modularAnalysis.buildEventShape`.
717)DOC");
718
719
720 REGISTER_VARIABLE("sphericity", sphericity, R"DOC(
721[Eventbased] Event sphericity, defined as the linear combination of the sphericity eigenvalues :math:`\lambda_i`: :math:`S = (3/2)(\lambda_2+\lambda_3)`.
722
723.. warning:: You have to run the Event Shape builder module for this variable to be meaningful.
724.. seealso:: :ref:`analysis_eventshape` and `modularAnalysis.buildEventShape`.
725)DOC");
726 REGISTER_VARIABLE("aplanarity", aplanarity, R"DOC(
727[Eventbased] Event aplanarity, defined as the 3/2 of the third sphericity eigenvalue.
728
729.. warning:: You have to run the Event Shape builder module for this variable to be meaningful.
730.. seealso:: :ref:`analysis_eventshape` and `modularAnalysis.buildEventShape`.
731)DOC");
732
733 REGISTER_VARIABLE("thrust", thrust, R"DOC(
734[Eventbased] Event thrust.
735
736.. warning:: You have to run the Event Shape builder module for this variable to be meaningful.
737.. seealso:: :ref:`analysis_eventshape` and `modularAnalysis.buildEventShape`.
738)DOC");
739 REGISTER_VARIABLE("thrustAxisX", thrustAxisX, R"DOC(
740[Eventbased] X component of the thrust axis.
741
742.. warning:: You have to run the Event Shape builder module for this variable to be meaningful.
743.. seealso:: :ref:`analysis_eventshape` and `modularAnalysis.buildEventShape`.
744)DOC");
745 REGISTER_VARIABLE("thrustAxisY", thrustAxisY, R"DOC(
746[Eventbased] Y component of the thrust axis.
747
748.. warning:: You have to run the Event Shape builder module for this variable to be meaningful.
749.. seealso:: :ref:`analysis_eventshape` and `modularAnalysis.buildEventShape`.
750)DOC");
751 REGISTER_VARIABLE("thrustAxisZ", thrustAxisZ, R"DOC(
752[Eventbased] Z component of the thrust axis.
753
754.. warning:: You have to run the Event Shape builder module for this variable to be meaningful.
755.. seealso:: :ref:`analysis_eventshape` and `modularAnalysis.buildEventShape`.
756)DOC");
757 REGISTER_VARIABLE("thrustAxisCosTheta", thrustAxisCosTheta, R"DOC(
758[Eventbased] Cosine of the polar angle component of the thrust axis.
759
760.. warning:: You have to run the Event Shape builder module for this variable to be meaningful.
761.. seealso:: :ref:`analysis_eventshape` and `modularAnalysis.buildEventShape`.
762)DOC");
763
764 REGISTER_VARIABLE("forwardHemisphereMass", forwardHemisphereMass, R"DOC(
765[Eventbased] Invariant mass of the particles flying in the same direction as the thrust axis.
766
767.. warning:: You have to run the Event Shape builder module for this variable to be meaningful.
768.. seealso:: :ref:`analysis_eventshape` and `modularAnalysis.buildEventShape`.
769
770)DOC","GeV/:math:`\\text{c}^2`");
771 REGISTER_VARIABLE("forwardHemisphereX", forwardHemisphereX, R"DOC(
772[Eventbased] X component of the total momentum of the particles flying in the same direction as the thrust axis.
773
774.. warning:: You have to run the Event Shape builder module for this variable to be meaningful.
775.. seealso:: :ref:`analysis_eventshape` and `modularAnalysis.buildEventShape`.
776
777)DOC","GeV/c");
778 REGISTER_VARIABLE("forwardHemisphereY", forwardHemisphereY, R"DOC(
779[Eventbased] Y component of the total momentum of the particles flying in the same direction as the thrust axis.
780
781.. warning:: You have to run the Event Shape builder module for this variable to be meaningful.
782.. seealso:: :ref:`analysis_eventshape` and `modularAnalysis.buildEventShape`.
783
784)DOC","GeV/c");
785 REGISTER_VARIABLE("forwardHemisphereZ", forwardHemisphereZ, R"DOC(
786[Eventbased] Z component of the total momentum of the particles flying in the same direction of the thrust axis.
787
788.. warning:: You have to run the Event Shape builder module for this variable to be meaningful.
789.. seealso:: :ref:`analysis_eventshape` and `modularAnalysis.buildEventShape`.
790
791)DOC","GeV/c");
792 REGISTER_VARIABLE("forwardHemisphereMomentum", forwardHemisphereMomentum, R"DOC(
793[Eventbased] Total momentum of the particles flying in the same direction as the thrust axis.
794
795.. warning:: You have to run the Event Shape builder module for this variable to be meaningful.
796.. seealso:: :ref:`analysis_eventshape` and `modularAnalysis.buildEventShape`.
797
798)DOC","GeV/c");
799 REGISTER_VARIABLE("forwardHemisphereEnergy", forwardHemisphereEnergy, R"DOC(
800[Eventbased] Total energy of the particles flying in the same direction as the thrust axis.
801
802.. warning:: You have to run the Event Shape builder module for this variable to be meaningful.
803.. seealso:: :ref:`analysis_eventshape` and `modularAnalysis.buildEventShape`.
804
805)DOC","GeV");
806 REGISTER_VARIABLE("backwardHemisphereMass", backwardHemisphereMass, R"DOC(
807[Eventbased] Invariant mass of the particles flying in the direction opposite to the thrust axis.
808
809.. warning:: You have to run the Event Shape builder module for this variable to be meaningful.
810.. seealso:: :ref:`analysis_eventshape` and `modularAnalysis.buildEventShape`.
811
812)DOC","GeV/:math:`\\text{c}^2`");
813 REGISTER_VARIABLE("backwardHemisphereX", backwardHemisphereX, R"DOC(
814[Eventbased] X component of the total momentum of the particles flying in the direction opposite to the thrust axis.
815
816.. warning:: You have to run the Event Shape builder module for this variable to be meaningful.
817.. seealso:: :ref:`analysis_eventshape` and `modularAnalysis.buildEventShape`.
818
819)DOC","GeV/c");
820 REGISTER_VARIABLE("backwardHemisphereY", backwardHemisphereY, R"DOC(
821[Eventbased] Y component of the total momentum of the particles flying in the direction opposite to the thrust axis.
822
823.. warning:: You have to run the Event Shape builder module for this variable to be meaningful.
824.. seealso:: :ref:`analysis_eventshape` and `modularAnalysis.buildEventShape`.
825
826)DOC","GeV/c");
827 REGISTER_VARIABLE("backwardHemisphereZ", backwardHemisphereZ, R"DOC(
828[Eventbased] Z component of the total momentum of the particles flying in the direction opposite to the thrust axis.
829
830.. warning:: You have to run the Event Shape builder module for this variable to be meaningful.
831.. seealso:: :ref:`analysis_eventshape` and `modularAnalysis.buildEventShape`.
832
833)DOC","GeV/c");
834 REGISTER_VARIABLE("backwardHemisphereMomentum", backwardHemisphereMomentum, R"DOC(
835[Eventbased] Total momentum of the particles flying in the direction opposite to the thrust axis.
836
837.. warning:: You have to run the Event Shape builder module for this variable to be meaningful.
838.. seealso:: :ref:`analysis_eventshape` and `modularAnalysis.buildEventShape`.
839
840)DOC","GeV/c");
841 REGISTER_VARIABLE("backwardHemisphereEnergy", backwardHemisphereEnergy, R"DOC(
842[Eventbased] Total energy of the particles flying in the direction opposite to the thrust axis.
843
844.. warning:: You have to run the Event Shape builder module for this variable to be meaningful.
845.. seealso:: :ref:`analysis_eventshape` and `modularAnalysis.buildEventShape`.
846
847)DOC","GeV");
848
849 }
851}
static const double doubleNaN
quiet_NaN
Definition: Const.h:703
std::function< VarVariant(const Particle *)> FunctionPtr
functions stored take a const Particle* and return VarVariant.
Definition: Manager.h:113
const Var * getVariable(std::string name)
Get the variable belonging to the given key.
Definition: Manager.cc:57
static Manager & Instance()
get singleton instance.
Definition: Manager.cc:25
Abstract base class for different kinds of events.