KaliVeda  1.13/01
Heavy-Ion Analysis Toolkit
KVTemplateEvent.h
Go to the documentation of this file.
1 /***************************************************************************
2  kvevent.h - description
3  -------------------
4  begin : Sun May 19 2002
5  copyright : (C) 2002 by J.D. Frankland
6  email : frankland@ganil.fr
7 
8 $Id: KVEvent.h,v 1.29 2008/12/17 11:23:12 ebonnet Exp $
9  ***************************************************************************/
10 
11 /***************************************************************************
12  * *
13  * This program is free software; you can redistribute it and/or modify *
14  * it under the terms of the GNU General Public License as published by *
15  * the Free Software Foundation; either version 2 of the License, or *
16  * (at your option) any later version. *
17  * *
18  ***************************************************************************/
19 
20 #ifndef KVEVENT_H
21 #define KVEVENT_H
22 
23 #define KVEVENT_PART_INDEX_OOB "Particle index %d out of bounds [1,%d]"
24 
25 #include "TTree.h"
26 #include "TVector3.h"
27 #include "TClonesArray.h"
28 #include "KVEvent.h"
29 #include "KVConfig.h"
30 #include "TRotation.h"
31 #include "TLorentzRotation.h"
33 #include "KVNameValueList.h"
34 #include "TMethodCall.h"
35 
36 #include <KVFrameTransform.h>
37 #include <KVIntegerList.h>
38 #include <TH1.h>
39 #include <TObjString.h>
40 #include <iterator>
41 
42 class KVIntegerList;
43 
89 template <typename Particle>
90 class KVTemplateEvent: public KVEvent {
91 
92 public:
108 
117  class Iterator {
118  public:
119  typedef std::forward_iterator_tag iterator_category;
120  typedef Particle value_type;
121  typedef std::ptrdiff_t difference_type;
122  typedef Particle* pointer;
123  typedef Particle& reference;
124 
125  enum Type { // type of iterator
126  Null, // null value
127  All, // include all particles
128  OK, // include particles which are "OK"
129  Group, // include particles belonging to group
130  Bad // type mismatch: iterator goes straight to end
131  };
132 
133  private:
134  KVTemplateParticleCondition<Particle> fSelection;//condition for selecting particles
135  TIter fIter;//iterator over TClonesArray
136  Type fType;//iterator type
137  mutable Bool_t fIterating;//=kTRUE when iteration in progress
139  {
144 
145  if (fType == Type::Bad)
146  return kFALSE;
147  return fSelection.Test(current());
148  }
149  Particle* current() const
150  {
152  return static_cast<Particle*>(*fIter);
153  }
154  public:
156  : fSelection(),
157  fIter(static_cast<TIterator*>(nullptr)),
158  fType(Type::Null),
160  {}
161  Iterator(const Iterator& i)
162  : fSelection(i.fSelection),
163  fIter(i.fIter),
164  fType(i.fType),
166  {}
167 
169  : fSelection(selection), fIter(e->GetParticleArray()), fIterating(kTRUE)
170  {
173 
174  if (!e->GetParticleArray()->GetClass()->InheritsFrom(Particle::Class())) {
175  ::Warning("KVTemplateEvent::Iterator", "KVTemplateEvent<%s>::Iterator for %s particles requested for event containing %s particles. Iteration is aborted.",
176  Particle::Class()->GetName(), Particle::Class()->GetName(), e->GetParticleArray()->GetClass()->GetName());
177  fType = Bad;
178  }
179 
180  fIter.Begin();// set iterator to first particle of event corresponding to selection
181  while ((current() != nullptr) && !AcceptableIteration()) ++fIter;
182  if (current() == nullptr) fIterating = kFALSE;
183  }
184 
187  {
190 
191  if (!e.GetParticleArray()->GetClass()->InheritsFrom(Particle::Class())) {
192  ::Warning("KVTemplateEvent::Iterator", "KVTemplateEvent<%s>::Iterator for %s particles requested for event containing %s particles. Iteration is aborted.",
193  Particle::Class()->GetName(), Particle::Class()->GetName(), e.GetParticleArray()->GetClass()->GetName());
194  fType = Bad;
195  }
196 
197  fIter.Begin();// set iterator to first particle of event corresponding to selection
198  while ((current() != nullptr) && !AcceptableIteration()) ++fIter;
199  if (current() == nullptr) fIterating = kFALSE;
200  }
201 
202  Iterator(const KVEvent* e, Type t = Type::All, TString grp = "")
204  {
206 
207  if (!e->GetParticleArray()->GetClass()->InheritsFrom(Particle::Class())) {
208  ::Warning("KVTemplateEvent::Iterator", "KVTemplateEvent<%s>::Iterator for %s particles requested for event containing %s particles. Iteration is aborted.",
209  Particle::Class()->GetName(), Particle::Class()->GetName(), e->GetParticleArray()->GetClass()->GetName());
210  fType = Bad;
211  }
212  if (fType == Type::OK) {
213  fSelection.Set("ok", [](const Particle * n) {
214  return n->IsOK();
215  });
216  }
217  else if (fType == Type::Group) {
218  fSelection.Set("group", [grp](const Particle * n) {
219  return n->BelongsToGroup(grp);
220  });
221  }
222 
223  fIter.Begin();// set iterator to first particle of event corresponding to selection
224  while ((current() != nullptr) && !AcceptableIteration()) ++fIter;
225  if (current() == nullptr) fIterating = kFALSE;
226  }
227 
228  Iterator(const KVEvent& e, Type t = Type::All, TString grp = "")
230  {
232 
233  if (!e.GetParticleArray()->GetClass()->InheritsFrom(Particle::Class())) {
234  ::Warning("KVTemplateEvent::Iterator", "KVTemplateEvent<%s>::Iterator for %s particles requested for event containing %s particles. Iteration is aborted.",
235  Particle::Class()->GetName(), Particle::Class()->GetName(), e.GetParticleArray()->GetClass()->GetName());
236  fType = Bad;
237  }
238  if (fType == Type::OK) {
239  fSelection.Set("ok", [](const Particle * n) {
240  return n->IsOK();
241  });
242  }
243  else if (fType == Type::Group) {
244  fSelection.Set("group", [grp](const Particle * n) {
245  return n->BelongsToGroup(grp);
246  });
247  }
248 
249  fIter.Begin();// set iterator to first particle of event corresponding to selection
250  while ((current() != nullptr) && !AcceptableIteration()) ++fIter;
251  if (current() == nullptr) fIterating = kFALSE;
252  }
253 
254  Particle& operator* () const
255  {
257 
258  return *(current());
259  }
260  template<typename PointerType = Particle>
261  PointerType * get_pointer() const
262  {
263  return dynamic_cast<PointerType*>(current());
264  }
265  template<typename ReferenceType = Particle>
266  ReferenceType & get_reference() const
267  {
268  return dynamic_cast<ReferenceType&>(*current());
269  }
270  template<typename PointerType = Particle>
271  const PointerType * get_const_pointer() const
272  {
273  return dynamic_cast<const PointerType*>(current());
274  }
275  template<typename ReferenceType = Particle>
276  const ReferenceType & get_const_reference() const
277  {
278  return dynamic_cast<const ReferenceType&>(*current());
279  }
280  Bool_t operator!= (const Iterator& it) const
281  {
283  return current() != it.current();
284  }
285  Bool_t operator== (const Iterator& it) const
286  {
288  return current() == it.current();
289  }
291  {
294  if (current() == nullptr) fIterating = kFALSE;
295  ++fIter;
296  while (current() != nullptr && !AcceptableIteration()) ++fIter;
297  return *this;
298  }
300  {
303  Iterator tmp(*this);
304  operator++();
305  return tmp;
306  }
308  {
310  if (this != &rhs) { // check self-assignment based on address of object
311  fSelection = rhs.fSelection;
312  fIter = rhs.fIter;
313  fType = rhs.fType;
314  fIterating = rhs.fIterating;
315  }
316  return *this;
317  }
318 
319  static Iterator End()
320  {
321  return Iterator();
322  }
323 
324  void Reset(Type t = Type::Null, TString grp = "")
325  {
331  if (t != Type::Null) {
332  fType = t;
333  if (fType == Type::OK) {
334  fSelection.Set("ok", [](const Particle * n) {
335  return n->IsOK();
336  });
337  }
338  else if (fType == Type::Group) {
339  fSelection.Set("group", [grp](const Particle * n) {
340  return n->BelongsToGroup(grp);
341  });
342  }
343  }
344  fIter.Begin();
345  fIterating = kTRUE;
346  while ((current() != nullptr) && !AcceptableIteration()) ++fIter;
347  }
349  {
351  return fIterating;
352  }
353  void SetIsIterating(Bool_t on = kTRUE)
354  {
356  fIterating = on;
357  }
358  };
359 protected:
360  mutable Iterator fIter;
361 
362 public:
364  : KVEvent(Particle::Class(), mult)
365  {}
366 
367  Particle* AddParticle()
368  {
373 
374  Int_t mult = GetMult();
375 #ifdef __WITHOUT_TCA_CONSTRUCTED_AT
376  Particle* tmp = (Particle*) ConstructedAt(mult, "C");
377 #else
378  Particle* tmp = (Particle*) fParticles->ConstructedAt(mult, "C");
379 #endif
380  if (!tmp) {
381  Error("AddParticle", "Allocation failure, Mult=%d", mult);
382  return 0;
383  }
384  return tmp;
385  }
386  Particle* GetParticle(Int_t npart) const
387  {
389 
390  if (npart < 1 || npart > fParticles->GetEntriesFast()) {
391  Error("GetParticle", KVEVENT_PART_INDEX_OOB, npart,
393  return 0;
394  }
395 
396  return (Particle*)((*fParticles)[npart - 1]);
397  }
398  virtual Int_t GetMult(Option_t* opt = "") const
399  {
408 
409  if (TString(opt) == "") return KVEvent::GetMult();
410  Int_t fMultOK = 0;
411  for (Iterator it = GetNextParticleIterator(opt); it != end(); ++it) ++fMultOK;
412  return fMultOK;
413  }
414  Int_t GetMultiplicity(Int_t Z, Int_t A = 0, Option_t* opt = "")
415  {
421 
422  if (A > 0) return (Int_t)GetSum("IsIsotope", "int,int", Form("%d,%d", Z, A), opt);
423  return (Int_t)GetSum("IsElement", "int", Form("%d", Z), opt);
424  }
425  void GetMultiplicities(Int_t mult[], const TString& species, Option_t* opt = "")
426  {
438 
439  std::unique_ptr<TObjArray> spec(species.Tokenize(", "));// remove any spaces
440  Int_t nspec = spec->GetEntries();
441  memset(mult, 0, nspec * sizeof(Int_t)); // set multiplicities to zero
442  for (Iterator it = GetNextParticleIterator(opt); it != end(); ++it) {
443  for (int i = 0; i < nspec; i++) {
444  if (((TObjString*)(*spec)[i])->String() == (*it).GetSymbol()) mult[i] += 1;
445  }
446  }
447  }
448  Double_t GetSum(const Char_t* Nucleus_method, Option_t* opt = "")
449  {
457 
458  Double_t fSum = 0;
459  TMethodCall mt;
460  mt.InitWithPrototype(Particle::Class(), Nucleus_method, "");
461 
462  if (mt.IsValid()) {
464  if (mt.ReturnType() == TMethodCall::kLong) {
465  Long_t ret;
466  for (; it != end(); ++it) {
467  Particle* tmp = it.get_pointer();
468  mt.Execute(tmp, "", ret);
469  fSum += ret;
470  }
471  }
472  else if (mt.ReturnType() == TMethodCall::kDouble) {
473  Double_t ret;
474  for (; it != end(); ++it) {
475  Particle* tmp = it.get_pointer();
476  mt.Execute(tmp, "", ret);
477  fSum += ret;
478  }
479  }
480  }
481 
482  return fSum;
483  }
484  Double_t GetSum(const Char_t* Nucleus_method, const Char_t* method_prototype, const Char_t* args, Option_t* opt = "")
485  {
492 
493  Double_t fSum = 0;
494  TMethodCall mt;
495  mt.InitWithPrototype(Particle::Class(), Nucleus_method, method_prototype);
496 
497  if (mt.IsValid()) {
499  if (mt.ReturnType() == TMethodCall::kLong) {
500  Long_t ret;
501  for (; it != end(); ++it) {
502  Particle* tmp = it.get_pointer();
503  mt.Execute(tmp, args, ret);
504  fSum += ret;
505  }
506  }
507  else if (mt.ReturnType() == TMethodCall::kDouble) {
508  Double_t ret;
509  for (; it != end(); ++it) {
510  Particle* tmp = it.get_pointer();
511  mt.Execute(tmp, args, ret);
512  fSum += ret;
513  }
514  }
515  }
516 
517  return fSum;
518  }
519  void FillHisto(TH1* h, const Char_t* Nucleus_method, Option_t* opt = "")
520  {
527 
528  TMethodCall mt;
529  mt.InitWithPrototype(Particle::Class(), Nucleus_method, "");
530 
531  if (mt.IsValid()) {
533  if (mt.ReturnType() == TMethodCall::kLong) {
534  Long_t ret;
535  for (; it != end(); ++it) {
536  Particle* tmp = it.get_pointer();
537  mt.Execute(tmp, "", ret);
538  h->Fill((Double_t)ret);
539  }
540  }
541  else if (mt.ReturnType() == TMethodCall::kDouble) {
542  Double_t ret;
543  for (; it != end(); ++it) {
544  Particle* tmp = it.get_pointer();
545  mt.Execute(tmp, "", ret);
546  h->Fill(ret);
547  }
548  }
549  }
550  }
551  void FillHisto(TH1* h, const Char_t* Nucleus_method, const Char_t* method_prototype, const Char_t* args, Option_t* opt = "")
552  {
558 
559  TMethodCall mt;
560  mt.InitWithPrototype(Particle::Class(), Nucleus_method, method_prototype);
561 
562  if (mt.IsValid()) {
564  if (mt.ReturnType() == TMethodCall::kLong) {
565  Long_t ret;
566  for (; it != end(); ++it) {
567  Particle* tmp = it.get_pointer();
568  mt.Execute(tmp, args, ret);
569  h->Fill((Double_t)ret);
570  }
571  }
572  else if (mt.ReturnType() == TMethodCall::kDouble) {
573  Double_t ret;
574  for (; it != end(); ++it) {
575  Particle* tmp = it.get_pointer();
576  mt.Execute(tmp, args, ret);
577  h->Fill(ret);
578  }
579  }
580  }
581  }
582  virtual void Print(Option_t* t = "") const
583  {
586 
587  std::cout << "\nKVEvent with " << GetMult(t) << " particles :" << std::endl;
588  std::cout << "------------------------------------" << std::endl;
589  fParameters.Print();
590  for (Iterator it = GetNextParticleIterator(t); it != end(); ++it)(*it).Print();
591  }
592  virtual void ls(Option_t* t = "") const
593  {
594  Print(t);
595  }
596  Particle* GetParticleWithName(const Char_t* name) const
597  {
600 
601  Particle* tmp = (Particle*)fParticles->FindObject(name);
602  return tmp;
603  }
604  Particle* GetParticle(const Char_t* group_name) const
605  {
607 
608  Iterator it = GetNextParticleIterator(group_name);
609  Particle* tmp = it.get_pointer();
610  if (!tmp) Warning("GetParticle", "Particle not found: %s", group_name);
611  return tmp;
612  }
613 
614  Iterator begin() const
615  {
617  return Iterator(this);
618  }
619  Iterator end() const
620  {
622  return Iterator::End();
623  }
624 
625  Particle* GetNextParticle(Option_t* opt = "") const
626  {
649 
650  TString Opt(opt);
651  Opt.ToUpper();
652 
654  if (fIter.IsIterating()) return &(*(fIter++));
655 
657  Bool_t ok_iter = (Opt == "OK");
658  Bool_t grp_iter = (!ok_iter && Opt.Length());
659 
660  if (ok_iter) fIter = Iterator(this, Iterator::Type::OK);
661  else if (grp_iter) fIter = Iterator(this, Iterator::Type::Group, Opt);
662  else fIter = Iterator(this);
663  return &(*(fIter++));
664  }
665  void ResetGetNextParticle() const
666  {
670 
672  }
673 
675  {
683 
684  for (Iterator it = begin(); it != end(); ++it) {
685  (*it).ResetEnergy();
686  }
687  }
688 
689  void DefineGroup(const Char_t* groupname, const Char_t* from = "")
690  {
695 
696  for (Iterator it = GetNextParticleIterator(from); it != end(); ++it) {
697  (*it).AddGroup(groupname);
698  }
699  }
700  void DefineGroup(const Char_t* groupname, KVTemplateParticleCondition<Particle>* cond, const Char_t* from = "")
701  {
705 
706  for (Iterator it = GetNextParticleIterator(from); it != end(); ++it) {
707  if (cond->Test(it.get_const_pointer()))(*it).AddGroup(groupname, from);
708  }
709  }
710 
711  void SetFrame(const Char_t* frame, const KVFrameTransform& ft)
712  {
718 
719  for (auto& part : *this) {
720  part.SetFrame(frame, ft);
721  }
722  }
723  void SetFrame(const Char_t* newframe, const Char_t* oldframe, const KVFrameTransform& ft)
724  {
734 
735  for (auto& part : *this) {
736  part.SetFrame(newframe, oldframe, ft);
737  }
738  }
739  void ChangeFrame(const KVFrameTransform& ft, const KVString& name = "")
740  {
746 
747 
748  for (auto& part : *this) {
749  part.ChangeFrame(ft, name);
750  }
751  if (name != "") SetParameter("defaultFrame", name);
752  }
753  void ChangeDefaultFrame(const Char_t* newdef, const Char_t* defname = "")
754  {
760 
761  for (auto& part : *this) {
762  part.ChangeDefaultFrame(newdef, defname);
763  }
764  SetParameter("defaultFrame", newdef);
765  }
767  {
772 
773  for (auto& part : *this) {
774  part.UpdateAllFrames();
775  }
776  }
777 
778  template<typename U = Particle>
779  typename std::enable_if<std::is_base_of<KVNucleus, U>::value>::type
781  {
787 
788  IL->Clear();
789  for (Iterator it = GetNextParticleIterator(opt); it != end(); ++it) IL->Add((*it).GetZ());
790  IL->SetPopulation(1);
791  IL->CheckForUpdate();
792  }
793 
794  void GetMasses(std::vector<Double_t>& mass)
795  {
797 
798  mass.clear();
799  mass.reserve(GetMult());
800  int i = 0;
801  for (Iterator it = begin(); it != end(); ++it) mass[i++] = (*it).GetMass();
802  }
803 
804  template<typename U = Particle>
805  typename std::enable_if<std::is_base_of<KVNucleus, U>::value>::type
806  GetGSMasses(std::vector<Double_t>& mass)
807  {
809 
810  mass.clear();
811  mass.reserve(GetMult());
812  int i = 0;
813  for (Iterator it = begin(); it != end(); ++it) mass[i++] = (*it).GetMassGS();
814  }
815 
816  template<typename U = Particle>
817  typename std::enable_if<std::is_base_of<KVNucleus, U>::value, Double_t>::type
819  {
834 
835  Double_t sumM = 0;
836  Particle CN;
837  Int_t M = GetMult();
838  for (int i = 1; i <= M; i++) {
839  sumM += GetParticle(i)->GetMass();
840  CN += *(GetParticle(i));
841  }
842  return CN.GetMassGS() - sumM;
843  }
844  template<typename U = Particle>
845  typename std::enable_if < !std::is_base_of<KVNucleus, U>::value, Double_t >::type
847  {
849  return 0;
850  }
852  {
853  return get_channel_qvalue();
854  }
855  template<typename U = Particle>
856  typename std::enable_if<std::is_base_of<KVNucleus, U>::value, Double_t>::type
858  {
871 
872  Double_t sumM = 0;
873  Particle CN;
874  Int_t M = GetMult();
875  for (int i = 1; i <= M; i++) {
876  sumM += GetParticle(i)->GetMassGS();
877  CN += *(GetParticle(i));
878  }
879  return CN.GetMassGS() - sumM;
880  }
881  template<typename U = Particle>
882  typename std::enable_if<std::is_base_of<KVNucleus, U>::value, KVString>::type
884  {
892  fParticles->Sort();
893  KVString partition;
894 
895  KVNameValueList nvl;
896  partition = "";
897  for (Iterator it = begin(); it != end(); ++it) {
898  TString st = (*it).GetSymbol();
899  Int_t pop = TMath::Max(nvl.GetIntValue(st.Data()), 0);
900  pop += 1;
901  nvl.SetValue(st.Data(), pop);
902  }
903  for (Int_t ii = 0; ii < nvl.GetEntries(); ii += 1) {
904  Int_t pop = nvl.GetIntValue(ii);
905  if (pop == 1) partition += nvl.GetNameAt(ii);
906  else partition += Form("%s(%d)", nvl.GetNameAt(ii), pop);
907  if (ii < nvl.GetEntries() - 1) partition += " ";
908  }
909  return partition;
910  }
911 
912  template<typename U = Particle>
913  typename std::enable_if < !std::is_base_of<KVNucleus, U>::value, KVString >::type
915  {
917  return "";
918  }
919 
921  {
922  return get_partition_name();
923  }
924 
925  void SetFrameName(const KVString& name)
926  {
933 
934  for (Iterator it = begin(); it != end(); ++it) {
935  (*it).SetFrameName(name);
936  }
937  SetParameter("defaultFrame", name);
938  }
939 
940  struct EventIterator {
942  EventIterator(const KVEvent& event, typename Iterator::Type t = Iterator::Type::All, const TString& grp = "")
943  : it(event, t, grp)
944  {}
945  EventIterator(const KVEvent* event, typename Iterator::Type t = Iterator::Type::All, const TString& grp = "")
946  : it(event, t, grp)
947  {}
949  : it(event, selection)
950  {}
952  : it(event, selection)
953  {}
954  Iterator begin() const
955  {
956  return it;
957  }
958  Iterator end() const
959  {
960  return Iterator::End();
961  }
962  };
963 
964  struct EventOKIterator : public EventIterator {
967  {}
970  {}
971  };
972 
974  EventGroupIterator(const KVEvent& event, const TString& grp) :
975  EventIterator(event, Iterator::Type::Group, grp)
976  {}
977  EventGroupIterator(const KVEvent* event, const TString& grp) :
978  EventIterator(event, Iterator::Type::Group, grp)
979  {}
980  };
981 
983  {
989 
990  TString Opt(opt);
991  Opt.ToUpper();
992  Bool_t ok_iter = (Opt == "OK");
993  Bool_t grp_iter = (!ok_iter && Opt.Length());
994 
995  if (ok_iter) return EventOKIterator(this).begin();
996  else if (grp_iter) return EventGroupIterator(this, Opt).begin();
997  return EventIterator(this).begin();
998  }
999 
1001  {
1009  return EventIterator(this, c);
1010  }
1011 
1012  ClassDef(KVTemplateEvent, 0) //Base class for all types of multiparticle event
1013 };
1014 
1015 #endif
1016 
int Int_t
long Long_t
KVIonRangeTableMaterial * M
KVTemplateEvent< KVNucleus >::EventOKIterator EventOKIterator
KVTemplateEvent< KVNucleus >::EventGroupIterator EventGroupIterator
KVTemplateEvent< KVNucleus >::EventIterator EventIterator
#define KVEVENT_PART_INDEX_OOB
#define c(i)
#define e(i)
char Char_t
const Bool_t kFALSE
bool Bool_t
double Double_t
const Bool_t kTRUE
const char Option_t
#define ClassDef(name, id)
char * Form(const char *fmt,...)
Abstract base class container for multi-particle events.
Definition: KVEvent.h:66
KVNameValueList fParameters
general-purpose list of parameters
Definition: KVEvent.h:71
TClonesArray * fParticles
array of particles in event
Definition: KVEvent.h:70
const TClonesArray * GetParticleArray() const
Definition: KVEvent.h:137
virtual Int_t GetMult(Option_t *opt="") const
Definition: KVEvent.h:141
void SetParameter(const Char_t *name, ValType value) const
Definition: KVEvent.h:221
Utility class for kinematical transformations of KVParticle class.
Handle a list of positive integers (partition)
Definition: KVIntegerList.h:68
void Add(TArrayI *tab)
void Fill(Double_t* tab,Int_t mult);
void Clear(Option_t *option="")
Classe dérivée de TNamed, Reinitialisation de l'object.
void SetPopulation(Int_t pop)
Initialise la population à "pop".
Handles lists of named parameters with different types, a list of KVNamedParameter objects.
virtual void Print(Option_t *opt="") const
Int_t GetIntValue(const Char_t *name) const
void SetValue(const Char_t *name, value_type value)
const Char_t * GetNameAt(Int_t idx) const
Int_t GetEntries() const
Extension of ROOT TString class which allows backwards compatibility with ROOT v3....
Definition: KVString.h:72
Class used for iterating over particles in events.
const PointerType * get_const_pointer() const
Iterator(const KVEvent *e, const KVTemplateParticleCondition< Particle > &selection)
Particle & operator*() const
const Iterator & operator++()
std::forward_iterator_tag iterator_category
PointerType * get_pointer() const
Iterator(const Iterator &i)
ReferenceType & get_reference() const
Iterator(const KVEvent &e, Type t=Type::All, TString grp="")
std::ptrdiff_t difference_type
Iterator(const KVEvent &e, const KVTemplateParticleCondition< Particle > &selection)
KVTemplateParticleCondition< Particle > fSelection
Bool_t operator==(const Iterator &it) const
Iterator(const KVEvent *e, Type t=Type::All, TString grp="")
Bool_t operator!=(const Iterator &it) const
Particle * current() const
Iterator & operator=(const Iterator &rhs)
void Reset(Type t=Type::Null, TString grp="")
const ReferenceType & get_const_reference() const
void SetIsIterating(Bool_t on=kTRUE)
Base class for event classes (containers for different types of particle objects)
std::enable_if< !std::is_base_of< KVNucleus, U >::value, KVString >::type get_partition_name()
EventIterator ConditionalIterator(const KVTemplateParticleCondition< Particle > &c)
void ResetGetNextParticle() const
void ChangeFrame(const KVFrameTransform &ft, const KVString &name="")
Particle * GetNextParticle(Option_t *opt="") const
std::enable_if< std::is_base_of< KVNucleus, U >::value >::type GetGSMasses(std::vector< Double_t > &mass)
void GetMasses(std::vector< Double_t > &mass)
virtual void ls(Option_t *t="") const
void DefineGroup(const Char_t *groupname, const Char_t *from="")
Int_t GetMultiplicity(Int_t Z, Int_t A=0, Option_t *opt="")
void SetFrame(const Char_t *newframe, const Char_t *oldframe, const KVFrameTransform &ft)
Particle * GetParticle(Int_t npart) const
Iterator GetNextParticleIterator(Option_t *opt) const
KVString GetPartitionName()
KVTemplateEvent(Int_t mult=50)
internal iterator used by GetNextParticle()
Particle * GetParticle(const Char_t *group_name) const
void FillHisto(TH1 *h, const Char_t *Nucleus_method, const Char_t *method_prototype, const Char_t *args, Option_t *opt="")
Double_t GetSum(const Char_t *Nucleus_method, Option_t *opt="")
Iterator begin() const
std::enable_if< !std::is_base_of< KVNucleus, U >::value, Double_t >::type get_channel_qvalue() const
std::enable_if< std::is_base_of< KVNucleus, U >::value, KVString >::type get_partition_name()
std::enable_if< std::is_base_of< KVNucleus, U >::value, Double_t >::type get_channel_qvalue() const
void SetFrame(const Char_t *frame, const KVFrameTransform &ft)
Particle * GetParticleWithName(const Char_t *name) const
void DefineGroup(const Char_t *groupname, KVTemplateParticleCondition< Particle > *cond, const Char_t *from="")
void SetFrameName(const KVString &name)
std::enable_if< std::is_base_of< KVNucleus, U >::value, Double_t >::type GetGSChannelQValue() const
Double_t GetChannelQValue() const
void GetMultiplicities(Int_t mult[], const TString &species, Option_t *opt="")
std::enable_if< std::is_base_of< KVNucleus, U >::value >::type FillIntegerList(KVIntegerList *IL, Option_t *opt)
Iterator end() const
void FillHisto(TH1 *h, const Char_t *Nucleus_method, Option_t *opt="")
virtual Int_t GetMult(Option_t *opt="") const
void ChangeDefaultFrame(const Char_t *newdef, const Char_t *defname="")
Double_t GetSum(const Char_t *Nucleus_method, const Char_t *method_prototype, const Char_t *args, Option_t *opt="")
virtual void Print(Option_t *t="") const
Particle * AddParticle()
Bool_t Test(const ParticleType *nuc) const
void Set(const KVString &name, const LambdaFunc &F)
TObject * ConstructedAt(Int_t idx)
virtual void Sort(Int_t upto=kMaxInt)
virtual Int_t Fill(const char *name, Double_t w)
TIter & Begin()
EReturnType ReturnType()
static const EReturnType kLong
void InitWithPrototype(const char *function, const char *proto, ROOT::EFunctionMatchMode mode=ROOT::kConversionMatch)
Bool_t IsValid() const
static const EReturnType kDouble
void Execute()
virtual const char * GetName() const
Int_t GetEntriesFast() const
virtual TObject * FindObject(const char *name) const
virtual void Warning(const char *method, const char *msgfmt,...) const
virtual void Error(const char *method, const char *msgfmt,...) const
Ssiz_t Length() const
void ToUpper()
TObjArray * Tokenize(const TString &delim) const
const char * Data() const
Type
auto All(const RVec< T > &v) -> decltype(v[0]==false)
const Int_t n
TH1 * h
Double_t Max(Double_t a, Double_t b)
const char * Class
EventGroupIterator(const KVEvent &event, const TString &grp)
EventGroupIterator(const KVEvent *event, const TString &grp)
EventIterator(const KVEvent &event, typename Iterator::Type t=Iterator::Type::All, const TString &grp="")
EventIterator(const KVEvent *event, typename Iterator::Type t=Iterator::Type::All, const TString &grp="")
EventIterator(const KVEvent *event, const KVTemplateParticleCondition< Particle > &selection)
EventIterator(const KVEvent &event, const KVTemplateParticleCondition< Particle > &selection)
EventOKIterator(const KVEvent &event)
EventOKIterator(const KVEvent *event)