ROOT  6.06/09
Reference Guide
TEventList.cxx
Go to the documentation of this file.
1 // @(#)root/tree:$Id$
2 // Author: Rene Brun 11/02/97
3 
4 /*************************************************************************
5  * Copyright (C) 1995-2000, Rene Brun and Fons Rademakers. *
6  * All rights reserved. *
7  * *
8  * For the licensing terms see $ROOTSYS/LICENSE. *
9  * For the list of contributors see $ROOTSYS/README/CREDITS. *
10  *************************************************************************/
11 
12 /** \class TEventList
13 A TEventList object is a list of selected events (entries) in a TTree.
14 
15 A TEventList is automatically generated by TTree::Draw: example
16 ~~~ {.cpp}
17  tree->Draw(">>elist1","x<0 && y> 0");
18 ~~~
19 In this example, a TEventList object named "elist1" will be
20 generated. (Previous contents are overwritten).
21 ~~~ {.cpp}
22  tree->Draw(">>+elist1","x<0 && y> 0");
23 ~~~
24 In this example, selected entries are added to the list.
25 
26 The TEventList object is added to the list of objects in the current
27 directory.
28 
29 Use TTree:SetEventList(TEventList *list) to inform TTree that you
30 want to use the list as input. The following code gets a pointer to
31 the TEventList object created in the above commands:
32 ~~~ {.cpp}
33  TEventList *list = (TEventList*)gDirectory->Get("elist1");
34 ~~~
35 - Use function Enter to enter an element in the list
36 - The function Add may be used to merge two lists.
37 - The function Subtract may be used to subtract two lists.
38 - The function Reset may be used to reset a list.
39 - Use TEventList::Print(option) to print the contents.
40  (option "all" prints all the list entries).
41 - Operators + and - correspond to functions Add and Subtract.
42 - A TEventList object can be saved on a file via the Write function.
43 */
44 
45 #include "TEventList.h"
46 #include "TCut.h"
47 #include "TClass.h"
48 #include "TFile.h"
49 #include "TMath.h"
50 
52 
53 ////////////////////////////////////////////////////////////////////////////////
54 /// Default constructor for a EventList.
55 
57 {
58  fN = 0;
59  fSize = 100;
60  fDelta = 100;
61  fList = 0;
62  fDirectory = 0;
63  fReapply = kFALSE;
64 }
65 
66 ////////////////////////////////////////////////////////////////////////////////
67 /// Create a EventList.
68 ///
69 /// This Eventlist is added to the list of objects in current directory.
70 
71 TEventList::TEventList(const char *name, const char *title, Int_t initsize, Int_t delta)
72  :TNamed(name,title), fReapply(kFALSE)
73 {
74  fN = 0;
75  if (initsize > 100) fSize = initsize;
76  else fSize = 100;
77  if (delta > 100) fDelta = delta;
78  else fDelta = 100;
79  fList = 0;
81  if (fDirectory) fDirectory->Append(this);
82 }
83 
84 ////////////////////////////////////////////////////////////////////////////////
85 /// Copy constructor.
86 
88 {
89  fN = list.fN;
90  fSize = list.fSize;
91  fDelta = list.fDelta;
92  fList = new Long64_t[fSize];
93  for (Int_t i=0; i<fN; i++)
94  fList[i] = list.fList[i];
95  fReapply = list.fReapply;
96  fDirectory = 0;
97 }
98 
99 ////////////////////////////////////////////////////////////////////////////////
100 /// Default destructor for a EventList.
101 
103 {
104  delete [] fList; fList = 0;
105  if (fDirectory) fDirectory->Remove(this);
106  fDirectory = 0;
107 }
108 
109 ////////////////////////////////////////////////////////////////////////////////
110 /// Merge contents of alist with this list.
111 ///
112 /// Both alist and this list are assumed to be sorted prior to this call
113 
114 void TEventList::Add(const TEventList *alist)
115 {
116  Int_t i;
117  Int_t an = alist->GetN();
118  if (!an) return;
119  Long64_t *alst = alist->GetList();
120  if (!fList) {
121  fList = new Long64_t[an];
122  for (i=0;i<an;i++) fList[i] = alst[i];
123  fN = an;
124  fSize = an;
125  return;
126  }
127  Int_t newsize = fN + an;
128  Long64_t *newlist = new Long64_t[newsize];
129  Int_t newpos, alpos;
130  newpos = alpos = 0;
131  for (i=0;i<fN;i++) {
132  while (alpos < an && fList[i] > alst[alpos]) {
133  newlist[newpos] = alst[alpos];
134  newpos++;
135  alpos++;
136  }
137  if (alpos < an && fList[i] == alst[alpos]) alpos++;
138  newlist[newpos] = fList[i];
139  newpos++;
140  }
141  while (alpos < an) {
142  newlist[newpos] = alst[alpos];
143  newpos++;
144  alpos++;
145  }
146  delete [] fList;
147  fN = newpos;
148  fSize = newsize;
149  fList = newlist;
150 
151  TCut orig = GetTitle();
152  TCut added = alist->GetTitle();
153  TCut updated = orig || added;
154  SetTitle(updated.GetTitle());
155 }
156 
157 ////////////////////////////////////////////////////////////////////////////////
158 /// Return TRUE if list contains entry.
159 
161 {
162  if (GetIndex(entry) < 0) return kFALSE;
163  return kTRUE;
164 }
165 
166 ////////////////////////////////////////////////////////////////////////////////
167 /// Return TRUE if list contains entries from entrymin to entrymax included.
168 
170 {
171  Long64_t imax = TMath::BinarySearch(fN,fList,entrymax);
172  //printf("ContainsRange: entrymin=%lld, entrymax=%lld,imax=%lld, fList[imax]=%lld\n",entrymin,entrymax,imax,fList[imax]);
173 
174  if (fList[imax] < entrymin) return kFALSE;
175  return kTRUE;
176 }
177 
178 ////////////////////////////////////////////////////////////////////////////////
179 /// Called by TKey and others to automatically add us to a directory when we are read from a file.
180 
182 {
183  SetDirectory(dir);
184 }
185 
186 ////////////////////////////////////////////////////////////////////////////////
187 /// Enter element entry into the list.
188 
190 {
191  if (!fList) {
192  fList = new Long64_t[fSize];
193  fList[0] = entry;
194  fN = 1;
195  return;
196  }
197  if (fN>0 && entry==fList[fN-1]) return;
198  if (fN >= fSize) {
199  Int_t newsize = TMath::Max(2*fSize,fN+fDelta);
200  Resize(newsize-fSize);
201  }
202  if(fN==0 || entry>fList[fN-1]) {
203  fList[fN] = entry;
204  ++fN;
205  } else {
206  Int_t pos = TMath::BinarySearch(fN, fList, entry);
207  if(pos>=0 && entry==fList[pos])
208  return;
209  ++pos;
210  memmove( &(fList[pos+1]), &(fList[pos]), 8*(fN-pos));
211  fList[pos] = entry;
212  ++fN;
213  }
214 }
215 
216 ////////////////////////////////////////////////////////////////////////////////
217 /// Return value of entry at index in the list.
218 /// Return -1 if index is not in the list range.
219 
221 {
222  if (!fList) return -1;
223  if (index < 0 || index >= fN) return -1;
224  return fList[index];
225 }
226 
227 ////////////////////////////////////////////////////////////////////////////////
228 /// Return index in the list of element with value entry
229 /// array is supposed to be sorted prior to this call.
230 /// If match is found, function returns position of element.
231 /// If no match found, function returns -1.
232 
234 {
235  Long64_t nabove, nbelow, middle;
236  nabove = fN+1;
237  nbelow = 0;
238  while(nabove-nbelow > 1) {
239  middle = (nabove+nbelow)/2;
240  if (entry == fList[middle-1]) return middle-1;
241  if (entry < fList[middle-1]) nabove = middle;
242  else nbelow = middle;
243  }
244  return -1;
245 }
246 
247 ////////////////////////////////////////////////////////////////////////////////
248 /// Remove elements from this list that are NOT present in alist.
249 
251 {
252  if (!alist) return;
253  if (!fList) return;
254 
255  Long64_t *newlist = new Long64_t[fN];
256  Int_t newpos = 0;
257  Int_t i;
258  for (i=0;i<fN;i++) {
259  if (alist->GetIndex(fList[i]) >= 0) {
260  newlist[newpos] = fList[i];
261  newpos++;
262  }
263  }
264  delete [] fList;
265  fN = newpos;
266  fList = newlist;
267 
268  TCut orig = GetTitle();
269  TCut removed = alist->GetTitle();
270  TCut updated = orig && removed;
271  SetTitle(updated.GetTitle());
272 }
273 
274 ////////////////////////////////////////////////////////////////////////////////
275 /// Merge entries in all the TEventList in the collection in this event list.
276 
278 {
279  if (!list) return -1;
280  TIter next(list);
281 
282  //first loop to count the number of entries
283  TEventList *el;
284  Int_t nevents = 0;
285  while ((el = (TEventList*)next())) {
286  if (!el->InheritsFrom(TEventList::Class())) {
287  Error("Add","Attempt to add object of class: %s to a %s",el->ClassName(),this->ClassName());
288  return -1;
289  }
290  Add(el);
291  nevents += el->GetN();
292  }
293 
294  return nevents;
295 }
296 
297 ////////////////////////////////////////////////////////////////////////////////
298 /// Print contents of this list.
299 
300 void TEventList::Print(Option_t *option) const
301 {
302  printf("EventList:%s/%s, number of entries =%d, size=%d\n",GetName(),GetTitle(),fN,fSize);
303  if (!strstr(option,"all")) return;
304  Int_t i;
305  Int_t nbuf = 0;
306  char element[10];
307  char *line = new char[100];
308  snprintf(line,100,"%5d : ",0);
309  for (i=0;i<fN;i++) {
310  nbuf++;
311  if (nbuf > 10) {
312  printf("%s\n",line);
313  snprintf(line,100,"%5d : ",i);
314  nbuf = 1;
315  }
316  snprintf(element,10,"%7lld ",fList[i]);
317  strlcat(line,element,100);
318  }
319  if (nbuf) printf("%s\n",line);
320  delete [] line;
321 }
322 
323 ////////////////////////////////////////////////////////////////////////////////
324 /// Reset number of entries in event list.
325 
327 {
328  fN = 0;
329 }
330 
331 ////////////////////////////////////////////////////////////////////////////////
332 /// Resize list by delta entries.
333 
335 {
336  if (!delta) delta = fDelta;
337  fSize += delta;
338  Long64_t *newlist = new Long64_t[fSize];
339  for (Int_t i=0;i<fN;i++) newlist[i] = fList[i];
340  delete [] fList;
341  fList = newlist;
342 }
343 
344 ////////////////////////////////////////////////////////////////////////////////
345 /// Remove reference to this EventList from current directory and add
346 /// reference to new directory dir. dir can be 0 in which case the list
347 /// does not belong to any directory.
348 
350 {
351  if (fDirectory == dir) return;
352  if (fDirectory) fDirectory->Remove(this);
353  fDirectory = dir;
354  if (fDirectory) fDirectory->Append(this);
355 }
356 
357 ////////////////////////////////////////////////////////////////////////////////
358 /// Change the name of this TEventList.
359 
360 void TEventList::SetName(const char *name)
361 {
362  // TEventLists are named objects in a THashList.
363  // We must update the hashlist if we change the name
364  if (fDirectory) fDirectory->Remove(this);
365  fName = name;
366  if (fDirectory) fDirectory->Append(this);
367 }
368 
369 ////////////////////////////////////////////////////////////////////////////////
370 /// Sort list entries in increasing order
371 
373 {
374  Int_t *index = new Int_t[fN];
375  Long64_t *newlist = new Long64_t[fSize];
376  Int_t i,ind;
377  TMath::Sort(fN,fList,index); //sort in decreasing order
378  for (i=0;i<fN;i++) {
379  ind = index[fN-i-1];
380  newlist[i] = fList[ind];
381  }
382  for (i=fN;i<fSize;i++) {
383  newlist[i] = 0;
384  }
385  delete [] index;
386  delete [] fList;
387  fList = newlist;
388 }
389 
390 ////////////////////////////////////////////////////////////////////////////////
391 /// Stream an object of class TEventList.
392 
393 void TEventList::Streamer(TBuffer &b)
394 {
395  if (b.IsReading()) {
396  UInt_t R__s, R__c;
397  Version_t R__v = b.ReadVersion(&R__s, &R__c);
398  fDirectory = 0;
399  if (R__v > 1) {
400  b.ReadClassBuffer(TEventList::Class(), this, R__v, R__s, R__c);
402  return;
403  }
404  //====process old versions before automatic schema evolution
405  TNamed::Streamer(b);
406  b >> fN;
407  b >> fSize;
408  b >> fDelta;
409  if (fN) {
410  Int_t *tlist = new Int_t[fSize];
411  b.ReadFastArray(tlist,fN);
412  fList = new Long64_t[fSize];
413  for (Int_t i=0;i<fN;i++) fList[i] = tlist[i];
414  delete [] tlist;
415  }
417  b.CheckByteCount(R__s, R__c, TEventList::IsA());
418  //====end of old versions
419 
420  } else {
422  }
423 }
424 
425 ////////////////////////////////////////////////////////////////////////////////
426 /// Remove elements from this list that are present in alist.
427 
429 {
430  if (!alist) return;
431  if (!fList) return;
432 
433  Long64_t *newlist = new Long64_t[fN];
434  Int_t newpos = 0;
435  Int_t i;
436  for (i=0;i<fN;i++) {
437  if (alist->GetIndex(fList[i]) < 0) {
438  newlist[newpos] = fList[i];
439  newpos++;
440  }
441  }
442  delete [] fList;
443  fN = newpos;
444  fList = newlist;
445 
446  TCut orig = GetTitle();
447  TCut removed = alist->GetTitle();
448  TCut updated = orig && !removed;
449  SetTitle(updated.GetTitle());
450 }
451 
452 ////////////////////////////////////////////////////////////////////////////////
453 /// Assingment.
454 
456 {
457  if (this != &list) {
458  TNamed::operator=(list);
459  if (fSize < list.fSize) {
460  delete [] fList;
461  fList = new Long64_t[list.fSize];
462  }
463  fN = list.fN;
464  fSize = list.fSize;
465  fDelta = list.fDelta;
466  for (Int_t i=0; i<fN; i++)
467  fList[i] = list.fList[i];
468  }
469  return *this;
470 }
471 
472 ////////////////////////////////////////////////////////////////////////////////
473 /// Addition.
474 
475 TEventList operator+(const TEventList &list1, const TEventList &list2)
476 {
477  TEventList newlist = list1;
478  newlist.Add(&list2);
479  return newlist;
480 }
481 
482 ////////////////////////////////////////////////////////////////////////////////
483 /// Substraction
484 
485 TEventList operator-(const TEventList &list1, const TEventList &list2)
486 {
487  TEventList newlist = list1;
488  newlist.Subtract(&list2);
489  return newlist;
490 }
491 
492 ////////////////////////////////////////////////////////////////////////////////
493 /// Intersection.
494 
495 TEventList operator*(const TEventList &list1, const TEventList &list2)
496 {
497  TEventList newlist = list1;
498  newlist.Intersect(&list2);
499  return newlist;
500 }
501 
virtual void Enter(Long64_t entry)
Enter element entry into the list.
Definition: TEventList.cxx:189
virtual const char * GetTitle() const
Returns title of object.
Definition: TNamed.h:52
virtual Bool_t Contains(Long64_t entry)
Return TRUE if list contains entry.
Definition: TEventList.cxx:160
virtual Int_t WriteClassBuffer(const TClass *cl, void *pointer)=0
long long Long64_t
Definition: RtypesCore.h:69
Bool_t IsReading() const
Definition: TBuffer.h:81
virtual Bool_t InheritsFrom(const char *classname) const
Returns kTRUE if object inherits from class "classname".
Definition: TObject.cxx:487
short Version_t
Definition: RtypesCore.h:61
ClassImp(TSeqCollection) Int_t TSeqCollection TIter next(this)
Return index of object in collection.
TLine * line
const char Option_t
Definition: RtypesCore.h:62
#define gDirectory
Definition: TDirectory.h:218
virtual void Print(Option_t *option="") const
Print contents of this list.
Definition: TEventList.cxx:300
virtual void Intersect(const TEventList *list)
Remove elements from this list that are NOT present in alist.
Definition: TEventList.cxx:250
Buffer base class used for serializing objects.
Definition: TBuffer.h:40
virtual Int_t CheckByteCount(UInt_t startpos, UInt_t bcnt, const TClass *clss)=0
TDirectory * fDirectory
Definition: TEventList.h:41
ClassImp(TEventList) TEventList
Default constructor for a EventList.
Definition: TEventList.cxx:51
int Int_t
Definition: RtypesCore.h:41
bool Bool_t
Definition: RtypesCore.h:59
Int_t fDelta
Definition: TEventList.h:38
virtual void DirectoryAutoAdd(TDirectory *)
Called by TKey and others to automatically add us to a directory when we are read from a file...
Definition: TEventList.cxx:181
const Bool_t kFALSE
Definition: Rtypes.h:92
TEventList & operator=(const TEventList &list)
Assingment.
Definition: TEventList.cxx:455
virtual void Sort()
Sort list entries in increasing order.
Definition: TEventList.cxx:372
Int_t fSize
Definition: TEventList.h:37
virtual void Resize(Int_t delta=0)
Resize list by delta entries.
Definition: TEventList.cxx:334
virtual Int_t GetN() const
Definition: TEventList.h:58
void Class()
Definition: Class.C:29
The TNamed class is the base class for all named ROOT classes.
Definition: TNamed.h:33
virtual Bool_t ContainsRange(Long64_t entrymin, Long64_t entrymax)
Return TRUE if list contains entries from entrymin to entrymax included.
Definition: TEventList.cxx:169
void Sort(Index n, const Element *a, Index *index, Bool_t down=kTRUE)
Definition: TMath.h:1002
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
Definition: TObject.cxx:918
TEventList()
Pointer to directory holding this tree.
virtual void SetDirectory(TDirectory *dir)
Remove reference to this EventList from current directory and add reference to new directory dir...
Definition: TEventList.cxx:349
A specialized string object used for TTree selections.
Definition: TCut.h:27
virtual Long64_t GetEntry(Int_t index) const
Return value of entry at index in the list.
Definition: TEventList.cxx:220
TEventList operator*(const TEventList &list1, const TEventList &list2)
Intersection.
Definition: TEventList.cxx:495
TNamed & operator=(const TNamed &rhs)
TNamed assignment operator.
Definition: TNamed.cxx:40
virtual const char * ClassName() const
Returns name of class to which the object belongs.
Definition: TObject.cxx:187
Int_t fN
Definition: TEventList.h:36
Collection abstract base class.
Definition: TCollection.h:48
TClass * IsA() const
unsigned int UInt_t
Definition: RtypesCore.h:42
virtual void Append(TObject *obj, Bool_t replace=kFALSE)
Append object to this directory.
Definition: TDirectory.cxx:151
A TEventList object is a list of selected events (entries) in a TTree.
Definition: TEventList.h:33
virtual ~TEventList()
Default destructor for a EventList.
Definition: TEventList.cxx:102
virtual const char * GetName() const
Returns name of object.
Definition: TNamed.h:51
TEventList operator+(const TEventList &list1, const TEventList &list2)
Addition.
Definition: TEventList.cxx:475
Long64_t entry
TString fName
Definition: TNamed.h:36
virtual void ReadFastArray(Bool_t *b, Int_t n)=0
virtual Int_t ReadClassBuffer(const TClass *cl, void *pointer, const TClass *onfile_class=0)=0
virtual Int_t Merge(TCollection *list)
Merge entries in all the TEventList in the collection in this event list.
Definition: TEventList.cxx:277
virtual TObject * Remove(TObject *)
Remove an object from the in-memory list.
Describe directory structure in memory.
Definition: TDirectory.h:41
ClassImp(TMCParticle) void TMCParticle printf(": p=(%7.3f,%7.3f,%9.3f) ;", fPx, fPy, fPz)
virtual void SetName(const char *name)
Change the name of this TEventList.
Definition: TEventList.cxx:360
#define name(a, b)
Definition: linkTestLib0.cpp:5
virtual void Add(const TEventList *list)
Merge contents of alist with this list.
Definition: TEventList.cxx:114
virtual Int_t GetIndex(Long64_t entry) const
Return index in the list of element with value entry array is supposed to be sorted prior to this cal...
Definition: TEventList.cxx:233
Bool_t fReapply
Definition: TEventList.h:39
virtual Long64_t * GetList() const
Definition: TEventList.h:57
virtual void Subtract(const TEventList *list)
Remove elements from this list that are present in alist.
Definition: TEventList.cxx:428
Short_t Max(Short_t a, Short_t b)
Definition: TMathBase.h:202
Long64_t * fList
Definition: TEventList.h:40
void ResetBit(UInt_t f)
Definition: TObject.h:172
const Bool_t kTRUE
Definition: Rtypes.h:91
virtual void SetTitle(const char *title="")
Change (i.e. set) the title of the TNamed.
Definition: TNamed.cxx:152
TEventList operator-(const TEventList &list1, const TEventList &list2)
Substraction.
Definition: TEventList.cxx:485
Long64_t BinarySearch(Long64_t n, const T *array, T value)
Definition: TMath.h:944
virtual Version_t ReadVersion(UInt_t *start=0, UInt_t *bcnt=0, const TClass *cl=0)=0
virtual void Reset(Option_t *option="")
Reset number of entries in event list.
Definition: TEventList.cxx:326