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