Logo ROOT  
Reference Guide
TTreeCache.cxx
Go to the documentation of this file.
1 // @(#)root/tree:$Id$
2 // Author: Rene Brun 04/06/2006
3 
4 /*************************************************************************
5  * Copyright (C) 1995-2018, 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 TTreeCache
13 \ingroup tree
14 \brief A cache to speed-up the reading of ROOT datasets
15 
16 # A cache to speed-up the reading of ROOT datasets
17 
18 ## Table of Contents
19 - [Motivation](\ref motivation)
20 - [General Description](\ref description)
21 - [Changes in behaviour](\ref changesbehaviour)
22 - [Self-optimization](\ref cachemisses)
23 - [Examples of usage](\ref examples)
24 - [Check performance and stats](\ref checkPerf)
25 
26 \anchor motivation
27 ## Motivation: why having a cache is needed?
28 
29 When writing a TTree, the branch buffers are kept in memory.
30 A typical branch buffersize (before compression) is typically 32 KBytes.
31 After compression, the zipped buffer may be just a few Kbytes.
32 The branch buffers cannot be much larger in case of TTrees with several
33 hundred or thousand branches.
34 
35 When writing, this does not generate a performance problem because branch
36 buffers are always written sequentially and, thanks to OS optimisations,
37 content is flushed to the output file when a few MBytes of data are available.
38 On the other hand, when reading, one may hit performance problems because of
39 latencies e.g imposed by network.
40 For example in a WAN with 10ms latency, reading 1000 buffers of 10 KBytes each
41 with no cache will imply 10s penalty where a local read of the 10 MBytes would
42 take about 1 second.
43 
44 The TreeCache tries to prefetch all the buffers for the selected branches
45 in order to transfer a few multi-Megabytes large buffers instead of many
46 multi-kilobytes small buffers. In addition, TTreeCache can sort the blocks to
47 be read in increasing order such that the file is read sequentially.
48 
49 Systems like xrootd, dCache or httpd take advantage of the TTreeCache in
50 reading ahead as much data as they can and return to the application
51 the maximum data specified in the cache and have the next chunk of data ready
52 when the next request comes.
53 
54 ### Are there cases for which the usage of TTreeCache is detrimental for performance?
55 Yes, some corner cases. For example, when reading only a small fraction of all
56 entries such that not all branch buffers are read.
57 
58 \anchor description
59 ## General Description
60 This class acts as a file cache, registering automatically the baskets from
61 the branches being processed via direct manipulation of TTrees or with tools
62 such as TTree::Draw, TTree::Process, TSelector, TTreeReader and RDataFrame
63 when in the learning phase. The learning phase is by default 100 entries.
64 It can be changed via TTreeCache::SetLearnEntries.
65 
66 The usage of a TTreeCache can considerably improve the runtime performance at
67 the price of a modest investment in memory, in particular when the TTree is
68 accessed remotely, e.g. via a high latency network.
69 
70 For each TTree being processed a TTreeCache object is created.
71 This object is automatically deleted when the Tree is deleted or
72 when the file is deleted.
73 The user can change the size of the cache with the TTree::SetCacheSize method
74 (by default the size is 30 Megabytes). This feature can be controlled with the
75 environment variable `ROOT_TTREECACHE_SIZE` or the TTreeCache.Size option.
76 The entry range for which the cache is active can also be set with the
77 SetEntryRange method.
78 
79 \anchor changesbehaviour
80 ## Changes of behavior when using TChain and TEventList
81 
82 The usage of TChain or TEventList have influence on the behaviour of the cache:
83 
84 - Special case of a TChain
85  Once the training is done on the first Tree, the list of branches
86  in the cache is kept for the following files.
87 
88 - Special case of a TEventlist
89  if the Tree or TChain has a TEventlist, only the buffers
90  referenced by the list are put in the cache.
91 
92 The learning phase is started or restarted when:
93  - TTree automatically creates a cache.
94  - TTree::SetCacheSize is called with a non-zero size and a cache
95  did not previously exist
96  - TTreeCache::StartLearningPhase is called.
97  - TTreeCache::SetEntryRange is called
98  * and the learning is not yet finished
99  * and has not been set to manual
100  * and the new minimun entry is different.
101 
102 The learning period is stopped (and prefetching is started) when:
103  - TTreeCache::StopLearningPhase is called.
104  - An entry outside the 'learning' range is requested
105  The 'learning range is from fEntryMin (default to 0) to
106  fEntryMin + fgLearnEntries.
107  - A 'cached' TChain switches over to a new file.
108 
109 
110 \anchor cachemisses
111 ## Self-optimization in presence of cache misses
112 
113 The TTreeCache can optimize its behavior on a cache miss. When
114 miss optimization is enabled (see the SetOptimizeMisses method),
115 it tracks all branches utilized after the learning phase which caused a cache
116 miss.
117 When one cache miss occurs, all the utilized branches are be prefetched
118 for that event. This optimization utilizes the observation that infrequently
119 accessed branches are often accessed together.
120 An example scenario where such behavior is desirable, is an analysis where
121 a set of collections are read only for a few events in which a certain
122 condition is respected, e.g. a trigger fired.
123 
124 ### Additional memory and CPU usage when optimizing for cache misses
125 When this mode is enabled, the memory dedicated to the cache can increase
126 by at most a factor two in the case of cache miss.
127 Additionally, on the first miss of an event, we must iterate through all the
128 "active branches" for the miss cache and find the correct basket.
129 This can be potentially a CPU-expensive operation compared to, e.g., the
130 latency of a SSD. This is why the miss cache is currently disabled by default.
131 
132 \anchor examples
133 ## Example usages of TTreeCache
134 
135 A few use cases are discussed below. A cache may be created with automatic
136 sizing when a TTree is used:
137 
138 In some applications, e.g. central processing workflows of experiments, the list
139 of branches to read is known a priori. For these cases, the TTreeCache can be
140 instructed about the branches which will be read via explicit calls to the TTree
141 or TTreeCache interfaces.
142 In less streamlined applications such as analysis, predicting the branches which
143 will be read can be difficult. In such cases, ROOT I/O flags used branches
144 automatically when a branch buffer is read during the learning phase.
145 
146 In the examples below, portions of analysis code are shown.
147 The few statements involving the TreeCache are marked with `//<<<`
148 
149 ### ROOT::RDataFrame and TTreeReader Examples
150 
151 If you use RDataFrame or TTreeReader, the system will automatically cache the
152 best set of branches: no action is required by the user.
153 
154 ### TTree::Draw Example
155 
156 The TreeCache is automatically used by TTree::Draw. The method knows
157 which branches are used in the query and it puts automatically these branches
158 in the cache. The entry range is also inferred automatically.
159 
160 ### TTree::Process and TSelectors Examples
161 
162 The user must enable the cache and tell the system which branches to cache
163 and also specify the entry range. It is important to specify the entry range
164 in case only a subset of the events is processed to avoid wasteful caching.
165 
166 #### Reading all branches
167 
168 ~~~ {.cpp}
169  TTree *T;
170  f->GetObject(T, "mytree");
171  auto nentries = T->GetEntries();
172  auto cachesize = 10000000U; // 10 MBytes
173  T->SetCacheSize(cachesize); //<<<
174  T->AddBranchToCache("*", true); //<<< add all branches to the cache
175  T->Process("myselector.C+");
176  // In the TSelector::Process function we read all branches
177  T->GetEntry(i);
178  // ... Here the entry is processed
179 ~~~
180 
181 #### Reading a subset of all branches
182 
183 In the Process function we read a subset of the branches.
184 Only the branches used in the first entry will be put in the cache
185 ~~~ {.cpp}
186  TTree *T;
187  f->GetObject(T, "mytree");
188  // We want to process only the 200 first entries
189  auto nentries=200UL;
190  auto efirst = 0;
191  auto elast = efirst+nentries;
192  auto cachesize = 10000000U; // 10 MBytes
193  TTreeCache::SetLearnEntries(1); //<<< we can take the decision after 1 entry
194  T->SetCacheSize(cachesize); //<<<
195  T->SetCacheEntryRange(efirst,elast); //<<<
196  T->Process("myselector.C+","",nentries,efirst);
197  // In the TSelector::Process we read only 2 branches
198  auto b1 = T->GetBranch("branch1");
199  b1->GetEntry(i);
200  if (somecondition) return;
201  auto b2 = T->GetBranch("branch2");
202  b2->GetEntry(i);
203  ... Here the entry is processed
204 ~~~
205 ### Custom event loop
206 
207 #### Always using the same two branches
208 
209 In this example, exactly two branches are always used: those need to be
210 prefetched.
211 ~~~ {.cpp}
212  TTree *T;
213  f->GetObject(T, "mytree");
214  auto b1 = T->GetBranch("branch1");
215  auto b2 = T->GetBranch("branch2");
216  auto nentries = T->GetEntries();
217  auto cachesize = 10000000U; //10 MBytes
218  T->SetCacheSize(cachesize); //<<<
219  T->AddBranchToCache(b1, true); //<<< add branch1 and branch2 to the cache
220  T->AddBranchToCache(b2, true); //<<<
221  T->StopCacheLearningPhase(); //<<< we do not need the system to guess anything
222  for (auto i : TSeqL(nentries)) {
223  T->LoadTree(i); //<<< important call when calling TBranch::GetEntry after
224  b1->GetEntry(i);
225  if (some condition not met) continue;
226  b2->GetEntry(i);
227  if (some condition not met) continue;
228  // Here we read the full event only in some rare cases.
229  // There is no point in caching the other branches as it might be
230  // more economical to read only the branch buffers really used.
231  T->GetEntry(i);
232  ... Here the entry is processed
233  }
234 ~~~
235 #### Always using at least the same two branches
236 
237 In this example, two branches are always used: in addition, some analysis
238 functions are invoked and those may trigger the reading of other branches which
239 are a priori not known.
240 There is no point in prefetching branches that will be used very rarely: we can
241 rely on the system to cache the right branches.
242 ~~~ {.cpp}
243  TTree *T;
244  f->GetObject(T, "mytree");
245  auto nentries = T->GetEntries();
246  auto cachesize = 10000000; //10 MBytes
247  T->SetCacheSize(cachesize); //<<<
248  T->SetCacheLearnEntries(5); //<<< we can take the decision after 5 entries
249  auto b1 = T->GetBranch("branch1");
250  auto b2 = T->GetBranch("branch2");
251  for (auto i : TSeqL(nentries)) {
252  T->LoadTree(i);
253  b1->GetEntry(i);
254  if (some condition not met) continue;
255  b2->GetEntry(i);
256  // At this point we may call a user function where a few more branches
257  // will be read conditionally. These branches will be put in the cache
258  // if they have been used in the first 10 entries
259  if (some condition not met) continue;
260  // Here we read the full event only in some rare cases.
261  // There is no point in caching the other branches as it might be
262  // more economical to read only the branch buffers really used.
263  T->GetEntry(i);
264  .. process the rare but interesting cases.
265  ... Here the entry is processed
266  }
267 ~~~
268 
269 \anchor checkPerf
270 ## How can the usage and performance of TTreeCache be verified?
271 
272 Once the event loop terminated, the number of effective system reads for a
273 given file can be checked with a code like the following:
274 ~~~ {.cpp}
275  printf("Reading %lld bytes in %d transactions\n",myTFilePtr->GetBytesRead(), f->GetReadCalls());
276 ~~~
277 
278 Another handy command is:
279 ~~~ {.cpp}
280 myTreeOrChain.GetTree()->PrintCacheStats();
281 ~~~
282 
283 */
284 
285 #include "TSystem.h"
286 #include "TEnv.h"
287 #include "TTreeCache.h"
288 #include "TChain.h"
289 #include "TList.h"
290 #include "TBranch.h"
291 #include "TBranchElement.h"
292 #include "TEventList.h"
293 #include "TObjArray.h"
294 #include "TObjString.h"
295 #include "TRegexp.h"
296 #include "TLeaf.h"
297 #include "TFriendElement.h"
298 #include "TFile.h"
299 #include "TMath.h"
300 #include "TBranchCacheInfo.h"
301 #include "TVirtualPerfStats.h"
302 #include <limits.h>
303 
305 
307 
308 ////////////////////////////////////////////////////////////////////////////////
309 /// Default Constructor.
310 
311 TTreeCache::TTreeCache() : TFileCacheRead(), fPrefillType(GetConfiguredPrefillType())
312 {
313 }
314 
315 ////////////////////////////////////////////////////////////////////////////////
316 /// Constructor.
317 
319  : TFileCacheRead(tree->GetCurrentFile(), buffersize, tree), fEntryMax(tree->GetEntriesFast()), fEntryNext(0),
320  fBrNames(new TList), fTree(tree), fPrefillType(GetConfiguredPrefillType())
321 {
323  Int_t nleaves = tree->GetListOfLeaves()->GetEntriesFast();
324  fBranches = new TObjArray(nleaves);
325 }
326 
327 ////////////////////////////////////////////////////////////////////////////////
328 /// Destructor. (in general called by the TFile destructor)
329 
331 {
332  // Informe the TFile that we have been deleted (in case
333  // we are deleted explicitly by legacy user code).
334  if (fFile) fFile->SetCacheRead(0, fTree);
335 
336  delete fBranches;
337  if (fBrNames) {fBrNames->Delete(); delete fBrNames; fBrNames=0;}
338 }
339 
340 ////////////////////////////////////////////////////////////////////////////////
341 /// Add a branch discovered by actual usage to the list of branches to be stored
342 /// in the cache this function is called by TBranch::GetBasket
343 /// If we are not longer in the training phase this is an error.
344 /// Returns:
345 /// - 0 branch added or already included
346 /// - -1 on error
347 
348 Int_t TTreeCache::LearnBranch(TBranch *b, Bool_t subbranches /*= kFALSE*/)
349 {
350  if (!fIsLearning) {
351  return -1;
352  }
353 
354  // Reject branch that are not from the cached tree.
355  if (!b || fTree->GetTree() != b->GetTree()) return -1;
356 
357  // Is this the first addition of a branch (and we are learning and we are in
358  // the expected TTree), then prefill the cache. (We expect that in future
359  // release the Prefill-ing will be the default so we test for that inside the
360  // LearnPrefill call).
361  if (!fLearnPrefilling && fNbranches == 0) LearnPrefill();
362 
363  return AddBranch(b, subbranches);
364 }
365 
366 ////////////////////////////////////////////////////////////////////////////////
367 /// Add a branch to the list of branches to be stored in the cache
368 /// this function is called by the user via TTree::AddBranchToCache.
369 /// The branch is added even if we are outside of the training phase.
370 /// Returns:
371 /// - 0 branch added or already included
372 /// - -1 on error
373 
374 Int_t TTreeCache::AddBranch(TBranch *b, Bool_t subbranches /*= kFALSE*/)
375 {
376  // Reject branch that are not from the cached tree.
377  if (!b || fTree->GetTree() != b->GetTree()) return -1;
378 
379  //Is branch already in the cache?
380  Bool_t isNew = kTRUE;
381  for (int i=0;i<fNbranches;i++) {
382  if (fBranches->UncheckedAt(i) == b) {isNew = kFALSE; break;}
383  }
384  if (isNew) {
385  fTree = b->GetTree();
387  const char *bname = b->GetName();
388  if (fTree->IsA() == TChain::Class()) {
389  // If we have a TChain, we will need to use the branch name
390  // and we better disambiguate them (see atlasFlushed.root for example)
391  // in order to cache all the requested branches.
392  // We do not do this all the time as GetMother is slow (it contains
393  // a linear search from list of top level branch).
394  TString build;
395  const char *mothername = b->GetMother()->GetName();
396  if (b != b->GetMother() && mothername[strlen(mothername)-1] != '.') {
397  // Maybe we ought to prefix the name to avoid ambiguity.
398  auto bem = dynamic_cast<TBranchElement*>(b->GetMother());
399  if (bem->GetType() < 3) {
400  // Not a collection.
401  build = mothername;
402  build.Append(".");
403  if (strncmp(bname,build.Data(),build.Length()) != 0) {
404  build.Append(bname);
405  bname = build.Data();
406  }
407  }
408  }
409  }
410  fBrNames->Add(new TObjString(bname));
411  fNbranches++;
412  if (gDebug > 0) printf("Entry: %lld, registering branch: %s\n",b->GetTree()->GetReadEntry(),b->GetName());
413  }
414 
415  // process subbranches
416  Int_t res = 0;
417  if (subbranches) {
418  TObjArray *lb = b->GetListOfBranches();
419  Int_t nb = lb->GetEntriesFast();
420  for (Int_t j = 0; j < nb; j++) {
421  TBranch* branch = (TBranch*) lb->UncheckedAt(j);
422  if (!branch) continue;
423  if (AddBranch(branch, subbranches)<0) {
424  res = -1;
425  }
426  }
427  }
428  return res;
429 }
430 
431 ////////////////////////////////////////////////////////////////////////////////
432 /// Add a branch to the list of branches to be stored in the cache
433 /// this is to be used by user (thats why we pass the name of the branch).
434 /// It works in exactly the same way as TTree::SetBranchStatus so you
435 /// probably want to look over there for details about the use of bname
436 /// with regular expressions.
437 /// The branches are taken with respect to the Owner of this TTreeCache
438 /// (i.e. the original Tree)
439 /// NB: if bname="*" all branches are put in the cache and the learning phase stopped
440 /// Returns:
441 /// - 0 branch added or already included
442 /// - -1 on error
443 
444 Int_t TTreeCache::AddBranch(const char *bname, Bool_t subbranches /*= kFALSE*/)
445 {
446  TBranch *branch, *bcount;
447  TLeaf *leaf, *leafcount;
448 
449  Int_t i;
450  Int_t nleaves = (fTree->GetListOfLeaves())->GetEntriesFast();
451  TRegexp re(bname,kTRUE);
452  Int_t nb = 0;
453  Int_t res = 0;
454 
455  // first pass, loop on all branches
456  // for leafcount branches activate/deactivate in function of status
457  Bool_t all = kFALSE;
458  if (!strcmp(bname,"*")) all = kTRUE;
459  for (i=0;i<nleaves;i++) {
460  leaf = (TLeaf*)(fTree->GetListOfLeaves())->UncheckedAt(i);
461  branch = (TBranch*)leaf->GetBranch();
462  TString s = branch->GetName();
463  if (!all) { //Regexp gives wrong result for [] in name
464  TString longname;
465  longname.Form("%s.%s",fTree->GetName(),branch->GetName());
466  if (strcmp(bname,branch->GetName())
467  && longname != bname
468  && s.Index(re) == kNPOS) continue;
469  }
470  nb++;
471  if (AddBranch(branch, subbranches)<0) {
472  res = -1;
473  }
474  leafcount = leaf->GetLeafCount();
475  if (leafcount && !all) {
476  bcount = leafcount->GetBranch();
477  if (AddBranch(bcount, subbranches)<0) {
478  res = -1;
479  }
480  }
481  }
482  if (nb==0 && strchr(bname,'*')==0) {
483  branch = fTree->GetBranch(bname);
484  if (branch) {
485  if (AddBranch(branch, subbranches)<0) {
486  res = -1;
487  }
488  ++nb;
489  }
490  }
491 
492  //search in list of friends
493  UInt_t foundInFriend = 0;
494  if (fTree->GetListOfFriends()) {
495  TIter nextf(fTree->GetListOfFriends());
496  TFriendElement *fe;
497  TString name;
498  while ((fe = (TFriendElement*)nextf())) {
499  TTree *t = fe->GetTree();
500  if (t==0) continue;
501 
502  // If the alias is present replace it with the real name.
503  char *subbranch = (char*)strstr(bname,fe->GetName());
504  if (subbranch!=bname) subbranch = 0;
505  if (subbranch) {
506  subbranch += strlen(fe->GetName());
507  if ( *subbranch != '.' ) subbranch = 0;
508  else subbranch ++;
509  }
510  if (subbranch) {
511  name.Form("%s.%s",t->GetName(),subbranch);
512  if (name != bname && AddBranch(name, subbranches)<0) {
513  res = -1;
514  }
515  ++foundInFriend;
516  }
517  }
518  }
519  if (!nb && !foundInFriend) {
520  if (gDebug > 0) printf("AddBranch: unknown branch -> %s \n", bname);
521  Error("AddBranch", "unknown branch -> %s", bname);
522  return -1;
523  }
524  //if all branches are selected stop the learning phase
525  if (*bname == '*') {
526  fEntryNext = -1; // We are likely to have change the set of branches, so for the [re-]reading of the cluster.
528  }
529  return res;
530 }
531 
532 ////////////////////////////////////////////////////////////////////////////////
533 /// Remove a branch to the list of branches to be stored in the cache
534 /// this function is called by TBranch::GetBasket.
535 /// Returns:
536 /// - 0 branch dropped or not in cache
537 /// - -1 on error
538 
539 Int_t TTreeCache::DropBranch(TBranch *b, Bool_t subbranches /*= kFALSE*/)
540 {
541  if (!fIsLearning) {
542  return -1;
543  }
544 
545  // Reject branch that are not from the cached tree.
546  if (!b || fTree->GetTree() != b->GetTree()) return -1;
547 
548  //Is branch already in the cache?
549  if (fBranches->Remove(b)) {
550  --fNbranches;
551  if (gDebug > 0) printf("Entry: %lld, un-registering branch: %s\n",b->GetTree()->GetReadEntry(),b->GetName());
552  }
553  delete fBrNames->Remove(fBrNames->FindObject(b->GetName()));
554 
555  // process subbranches
556  Int_t res = 0;
557  if (subbranches) {
558  TObjArray *lb = b->GetListOfBranches();
559  Int_t nb = lb->GetEntriesFast();
560  for (Int_t j = 0; j < nb; j++) {
561  TBranch* branch = (TBranch*) lb->UncheckedAt(j);
562  if (!branch) continue;
563  if (DropBranch(branch, subbranches)<0) {
564  res = -1;
565  }
566  }
567  }
568  return res;
569 }
570 
571 ////////////////////////////////////////////////////////////////////////////////
572 /// Remove a branch to the list of branches to be stored in the cache
573 /// this is to be used by user (thats why we pass the name of the branch).
574 /// It works in exactly the same way as TTree::SetBranchStatus so you
575 /// probably want to look over there for details about the use of bname
576 /// with regular expressions.
577 /// The branches are taken with respect to the Owner of this TTreeCache
578 /// (i.e. the original Tree)
579 /// NB: if bname="*" all branches are put in the cache and the learning phase stopped
580 /// Returns:
581 /// - 0 branch dropped or not in cache
582 /// - -1 on error
583 
584 Int_t TTreeCache::DropBranch(const char *bname, Bool_t subbranches /*= kFALSE*/)
585 {
586  TBranch *branch, *bcount;
587  TLeaf *leaf, *leafcount;
588 
589  Int_t i;
590  Int_t nleaves = (fTree->GetListOfLeaves())->GetEntriesFast();
591  TRegexp re(bname,kTRUE);
592  Int_t nb = 0;
593  Int_t res = 0;
594 
595  // first pass, loop on all branches
596  // for leafcount branches activate/deactivate in function of status
597  Bool_t all = kFALSE;
598  if (!strcmp(bname,"*")) all = kTRUE;
599  for (i=0;i<nleaves;i++) {
600  leaf = (TLeaf*)(fTree->GetListOfLeaves())->UncheckedAt(i);
601  branch = (TBranch*)leaf->GetBranch();
602  TString s = branch->GetName();
603  if (!all) { //Regexp gives wrong result for [] in name
604  TString longname;
605  longname.Form("%s.%s",fTree->GetName(),branch->GetName());
606  if (strcmp(bname,branch->GetName())
607  && longname != bname
608  && s.Index(re) == kNPOS) continue;
609  }
610  nb++;
611  if (DropBranch(branch, subbranches)<0) {
612  res = -1;
613  }
614  leafcount = leaf->GetLeafCount();
615  if (leafcount && !all) {
616  bcount = leafcount->GetBranch();
617  if (DropBranch(bcount, subbranches)<0) {
618  res = -1;
619  }
620  }
621  }
622  if (nb==0 && strchr(bname,'*')==0) {
623  branch = fTree->GetBranch(bname);
624  if (branch) {
625  if (DropBranch(branch, subbranches)<0) {
626  res = -1;
627  }
628  ++nb;
629  }
630  }
631 
632  //search in list of friends
633  UInt_t foundInFriend = 0;
634  if (fTree->GetListOfFriends()) {
635  TIter nextf(fTree->GetListOfFriends());
636  TFriendElement *fe;
637  TString name;
638  while ((fe = (TFriendElement*)nextf())) {
639  TTree *t = fe->GetTree();
640  if (t==0) continue;
641 
642  // If the alias is present replace it with the real name.
643  char *subbranch = (char*)strstr(bname,fe->GetName());
644  if (subbranch!=bname) subbranch = 0;
645  if (subbranch) {
646  subbranch += strlen(fe->GetName());
647  if ( *subbranch != '.' ) subbranch = 0;
648  else subbranch ++;
649  }
650  if (subbranch) {
651  name.Form("%s.%s",t->GetName(),subbranch);
652  if (DropBranch(name, subbranches)<0) {
653  res = -1;
654  }
655  ++foundInFriend;
656  }
657  }
658  }
659  if (!nb && !foundInFriend) {
660  if (gDebug > 0) printf("DropBranch: unknown branch -> %s \n", bname);
661  Error("DropBranch", "unknown branch -> %s", bname);
662  return -1;
663  }
664  //if all branches are selected stop the learning phase
665  if (*bname == '*') {
666  fEntryNext = -1; // We are likely to have change the set of branches, so for the [re-]reading of the cluster.
667  }
668  return res;
669 }
670 
671 ////////////////////////////////////////////////////////////////////////////////
672 /// Start of methods for the miss cache.
673 ////////////////////////////////////////////////////////////////////////////////
674 
675 ////////////////////////////////////////////////////////////////////////////////
676 /// Enable / disable the miss cache.
677 ///
678 /// The first time this is called on a TTreeCache object, the corresponding
679 /// data structures will be allocated. Subsequent enable / disables will
680 /// simply turn the functionality on/off.
682 {
683 
684  if (opt && !fMissCache) {
685  ResetMissCache();
686  }
687  fOptimizeMisses = opt;
688 }
689 
690 ////////////////////////////////////////////////////////////////////////////////
691 /// Reset all the miss cache training.
692 ///
693 /// The contents of the miss cache will be emptied as well as the list of
694 /// branches used.
696 {
697 
698  fLastMiss = -1;
699  fFirstMiss = -1;
700 
701  if (!fMissCache) {
702  fMissCache.reset(new MissCache());
703  }
704  fMissCache->clear();
705 }
706 
707 ////////////////////////////////////////////////////////////////////////////////
708 /// For the event currently being fetched into the miss cache, find the IO
709 /// (offset / length tuple) to pull in the current basket for a given branch.
710 ///
711 /// Returns:
712 /// - IOPos describing the IO operation necessary for the basket on this branch
713 /// - On failure, IOPos.length will be set to 0.
715 {
716  if (R__unlikely(b.GetDirectory() == 0)) {
717  // printf("Branch at %p has no valid directory.\n", &b);
718  return IOPos{0, 0};
719  }
720  if (R__unlikely(b.GetDirectory()->GetFile() != fFile)) {
721  // printf("Branch at %p is in wrong file (branch file %p, my file %p).\n", &b, b.GetDirectory()->GetFile(),
722  // fFile);
723  return IOPos{0, 0};
724  }
725 
726  // printf("Trying to find a basket for branch %p\n", &b);
727  // Pull in metadata about branch; make sure it is valid
728  Int_t *lbaskets = b.GetBasketBytes();
729  Long64_t *entries = b.GetBasketEntry();
730  if (R__unlikely(!lbaskets || !entries)) {
731  // printf("No baskets or entries.\n");
732  return IOPos{0, 0};
733  }
734  // Int_t blistsize = b.GetListOfBaskets()->GetSize();
735  Int_t blistsize = b.GetWriteBasket();
736  if (R__unlikely(blistsize <= 0)) {
737  // printf("Basket list is size 0.\n");
738  return IOPos{0, 0};
739  }
740 
741  // Search for the basket that contains the event of interest. Unlike the primary cache, we
742  // are only interested in a single basket per branch - we don't try to fill the cache.
743  Long64_t basketOffset = TMath::BinarySearch(blistsize, entries, entry);
744  if (basketOffset < 0) { // No entry found.
745  // printf("No entry offset found for entry %ld\n", fTree->GetReadEntry());
746  return IOPos{0, 0};
747  }
748 
749  // Check to see if there's already a copy of this basket in memory. If so, don't fetch it
750  if ((basketOffset < blistsize) && b.GetListOfBaskets()->UncheckedAt(basketOffset)) {
751 
752  // printf("Basket is already in memory.\n");
753  return IOPos{0, 0};
754  }
755 
756  Long64_t pos = b.GetBasketSeek(basketOffset);
757  Int_t len = lbaskets[basketOffset];
758  if (R__unlikely(pos <= 0 || len <= 0)) {
759  /*printf("Basket returned was invalid (basketOffset=%ld, pos=%ld, len=%d).\n", basketOffset, pos, len);
760  for (int idx=0; idx<blistsize; idx++) {
761  printf("Basket entry %d, first event %d, pos %ld\n", idx, entries[idx], b.GetBasketSeek(idx));
762  }*/
763  return IOPos{0, 0};
764  } // Sanity check
765  // Do not cache a basket if it is bigger than the cache size!
766  if (R__unlikely(len > fBufferSizeMin)) {
767  // printf("Basket size is greater than the cache size.\n");
768  return IOPos{0, 0};
769  }
770 
771  return {pos, len};
772 }
773 
774 ////////////////////////////////////////////////////////////////////////////////
775 /// Given a particular IO description (offset / length) representing a 'miss' of
776 /// the TTreeCache's primary cache, calculate all the corresponding IO that
777 /// should be performed.
778 ///
779 /// `all` indicates that this function should search the set of _all_ branches
780 /// in this TTree. When set to false, we only search through branches that
781 /// have previously incurred a miss.
782 ///
783 /// Returns:
784 /// - TBranch pointer corresponding to the basket that will be retrieved by
785 /// this IO operation.
786 /// - If no corresponding branch could be found (or an error occurs), this
787 /// returns nullptr.
789 {
790  if (R__unlikely((pos < 0) || (len < 0))) {
791  return nullptr;
792  }
793 
794  int count = all ? (fTree->GetListOfLeaves())->GetEntriesFast() : fMissCache->fBranches.size();
795  fMissCache->fEntries.reserve(count);
796  fMissCache->fEntries.clear();
797  Bool_t found_request = kFALSE;
798  TBranch *resultBranch = nullptr;
799  Long64_t entry = fTree->GetReadEntry();
800 
801  std::vector<std::pair<size_t, Int_t>> basketsInfo;
802  auto perfStats = GetTree()->GetPerfStats();
803 
804  // printf("Will search %d branches for basket at %ld.\n", count, pos);
805  for (int i = 0; i < count; i++) {
806  TBranch *b =
807  all ? static_cast<TBranch *>(static_cast<TLeaf *>((fTree->GetListOfLeaves())->UncheckedAt(i))->GetBranch())
808  : fMissCache->fBranches[i];
809  IOPos iopos = FindBranchBasketPos(*b, entry);
810  if (iopos.fLen == 0) { // Error indicator
811  continue;
812  }
813  if (iopos.fPos == pos && iopos.fLen == len) {
814  found_request = kTRUE;
815  resultBranch = b;
816  // Note that we continue to iterate; fills up the rest of the entries in the cache.
817  }
818  // At this point, we are ready to push back a new offset
819  fMissCache->fEntries.emplace_back(std::move(iopos));
820 
821  if (R__unlikely(perfStats)) {
822  Int_t blistsize = b->GetWriteBasket();
823  Int_t basketNumber = -1;
824  for (Int_t bn = 0; bn < blistsize; ++bn) {
825  if (iopos.fPos == b->GetBasketSeek(bn)) {
826  basketNumber = bn;
827  break;
828  }
829  }
830  if (basketNumber >= 0)
831  basketsInfo.emplace_back((size_t)i, basketNumber);
832  }
833  }
834  if (R__unlikely(!found_request)) {
835  // We have gone through all the branches in this file and the requested basket
836  // doesn't appear to be in any of them. Likely a logic error / bug.
837  fMissCache->fEntries.clear();
838  }
839  if (R__unlikely(perfStats)) {
840  for (auto &info : basketsInfo) {
841  perfStats->SetLoadedMiss(info.first, info.second);
842  }
843  }
844  return resultBranch;
845 }
846 
847 ////////////////////////////////////////////////////////////////////////////////
848 ///
849 /// Process a cache miss; (pos, len) isn't in the buffer.
850 ///
851 /// The first time we have a miss, we buffer as many baskets we can (up to the
852 /// maximum size of the TTreeCache) in memory from all branches that are not in
853 /// the prefetch list.
854 ///
855 /// Subsequent times, we fetch all the buffers corresponding to branches that
856 /// had previously seen misses. If it turns out the (pos, len) isn't in the
857 /// list of branches, we treat this as if it was the first miss.
858 ///
859 /// Returns true if we were able to pull the data into the miss cache.
860 ///
862 {
863 
864  Bool_t firstMiss = kFALSE;
865  if (fFirstMiss == -1) {
867  firstMiss = kTRUE;
868  }
870  // The first time this is executed, we try to pull in as much data as we can.
871  TBranch *b = CalculateMissEntries(pos, len, firstMiss);
872  if (!b) {
873  if (!firstMiss) {
874  // TODO: this recalculates for *all* branches, throwing away the above work.
875  b = CalculateMissEntries(pos, len, kTRUE);
876  }
877  if (!b) {
878  // printf("ProcessMiss: pos %ld does not appear to correspond to a buffer in this file.\n", pos);
879  // We have gone through all the branches in this file and the requested basket
880  // doesn't appear to be in any of them. Likely a logic error / bug.
881  fMissCache->fEntries.clear();
882  return kFALSE;
883  }
884  }
885  // TODO: this should be a set.
886  fMissCache->fBranches.push_back(b);
887 
888  // OK, sort the entries
889  std::sort(fMissCache->fEntries.begin(), fMissCache->fEntries.end());
890 
891  // Now, fetch the buffer.
892  std::vector<Long64_t> positions;
893  positions.reserve(fMissCache->fEntries.size());
894  std::vector<Int_t> lengths;
895  lengths.reserve(fMissCache->fEntries.size());
896  ULong64_t cumulative = 0;
897  for (auto &mcentry : fMissCache->fEntries) {
898  positions.push_back(mcentry.fIO.fPos);
899  lengths.push_back(mcentry.fIO.fLen);
900  mcentry.fIndex = cumulative;
901  cumulative += mcentry.fIO.fLen;
902  }
903  fMissCache->fData.reserve(cumulative);
904  // printf("Reading %lu bytes into miss cache for %lu entries.\n", cumulative, fEntries->size());
905  fNMissReadPref += fMissCache->fEntries.size();
906  fFile->ReadBuffers(&(fMissCache->fData[0]), &(positions[0]), &(lengths[0]), fMissCache->fEntries.size());
908 
909  return kTRUE;
910 }
911 
912 ////////////////////////////////////////////////////////////////////////////////
913 /// Given an IO operation (pos, len) that was a cache miss in the primary TTC,
914 /// try the operation again with the miss cache.
915 ///
916 /// Returns true if the IO operation was successful and the contents of buf
917 /// were populated with the requested data.
918 ///
920 {
921 
922  if (!fOptimizeMisses) {
923  return kFALSE;
924  }
925  if (R__unlikely((pos < 0) || (len < 0))) {
926  return kFALSE;
927  }
928 
929  // printf("Checking the miss cache for offset=%ld, length=%d\n", pos, len);
930 
931  // First, binary search to see if the desired basket is already cached.
932  MissCache::Entry mcentry{IOPos{pos, len}};
933  auto iter = std::lower_bound(fMissCache->fEntries.begin(), fMissCache->fEntries.end(), mcentry);
934 
935  if (iter != fMissCache->fEntries.end()) {
936  if (len > iter->fIO.fLen) {
937  ++fNMissReadMiss;
938  return kFALSE;
939  }
940  auto offset = iter->fIndex;
941  memcpy(buf, &(fMissCache->fData[offset]), len);
942  // printf("Returning data from pos=%ld in miss cache.\n", offset);
943  ++fNMissReadOk;
944  return kTRUE;
945  }
946 
947  // printf("Data not in miss cache.\n");
948 
949  // Update the cache, looking for this (pos, len).
950  if (!ProcessMiss(pos, len)) {
951  // printf("Unable to pull data into miss cache.\n");
952  ++fNMissReadMiss;
953  return kFALSE;
954  }
955 
956  // OK, we updated the cache with as much information as possible. Search again for
957  // the entry we want.
958  iter = std::lower_bound(fMissCache->fEntries.begin(), fMissCache->fEntries.end(), mcentry);
959 
960  if (iter != fMissCache->fEntries.end()) {
961  auto offset = iter->fIndex;
962  // printf("Expecting data at offset %ld in miss cache.\n", offset);
963  memcpy(buf, &(fMissCache->fData[offset]), len);
964  ++fNMissReadOk;
965  return kTRUE;
966  }
967 
968  // This must be a logic bug. ProcessMiss should return false if (pos, len)
969  // wasn't put into fEntries.
970  ++fNMissReadMiss;
971  return kFALSE;
972 }
973 
974 ////////////////////////////////////////////////////////////////////////////////
975 /// End of methods for miss cache.
976 ////////////////////////////////////////////////////////////////////////////////
977 
978 namespace {
979 struct BasketRanges {
980  struct Range {
981  Long64_t fMin; ///< Inclusive minimum
982  Long64_t fMax; ///< Inclusive maximum
983 
984  Range() : fMin(-1), fMax(-1) {}
985 
986  void UpdateMin(Long64_t min)
987  {
988  if (fMin == -1 || min < fMin)
989  fMin = min;
990  }
991 
992  void UpdateMax(Long64_t max)
993  {
994  if (fMax == -1 || fMax < max)
995  fMax = max;
996  }
997 
998  Bool_t Contains(Long64_t entry) { return (fMin <= entry && entry <= fMax); }
999  };
1000 
1001  std::vector<Range> fRanges;
1002  std::map<Long64_t,size_t> fMinimums;
1003  std::map<Long64_t,size_t> fMaximums;
1004 
1005  BasketRanges(size_t nBranches) { fRanges.resize(nBranches); }
1006 
1007  void Update(size_t branchNumber, Long64_t min, Long64_t max)
1008  {
1009  Range &range = fRanges.at(branchNumber);
1010  auto old(range);
1011 
1012  range.UpdateMin(min);
1013  range.UpdateMax(max);
1014 
1015  if (old.fMax != range.fMax) {
1016  if (old.fMax != -1) {
1017  auto maxIter = fMaximums.find(old.fMax);
1018  if (maxIter != fMaximums.end()) {
1019  if (maxIter->second == 1) {
1020  fMaximums.erase(maxIter);
1021  } else {
1022  --(maxIter->second);
1023  }
1024  }
1025  }
1026  ++(fMaximums[max]);
1027  }
1028  }
1029 
1030  void Update(size_t branchNumber, size_t basketNumber, Long64_t *entries, size_t nb, size_t max)
1031  {
1032  Update(branchNumber, entries[basketNumber],
1033  (basketNumber < (nb - 1)) ? (entries[basketNumber + 1] - 1) : max - 1);
1034  }
1035 
1036  // Check that fMaximums and fMinimums are properly set
1037  bool CheckAllIncludeRange()
1038  {
1039  Range result;
1040  for (const auto &r : fRanges) {
1041  if (result.fMin == -1 || result.fMin < r.fMin) {
1042  if (r.fMin != -1)
1043  result.fMin = r.fMin;
1044  }
1045  if (result.fMax == -1 || r.fMax < result.fMax) {
1046  if (r.fMax != -1)
1047  result.fMax = r.fMax;
1048  }
1049  }
1050  // if (result.fMax < result.fMin) {
1051  // // No overlapping range.
1052  // }
1053 
1054  Range allIncludedRange(AllIncludedRange());
1055 
1056  return (result.fMin == allIncludedRange.fMin && result.fMax == allIncludedRange.fMax);
1057  }
1058 
1059  // This returns a Range object where fMin is the maximum of all the minimun entry
1060  // number loaded for each branch and fMax is the minimum of all the maximum entry
1061  // number loaded for each branch.
1062  // As such it is valid to have fMin > fMax, this is the case where there
1063  // are no overlap between the branch's range. For example for 2 branches
1064  // where we have for one the entry [50,99] and for the other [0,49] then
1065  // we will have fMin = max(50,0) = 50 and fMax = min(99,49) = 49
1066  Range AllIncludedRange()
1067  {
1068  Range result;
1069  if (!fMinimums.empty())
1070  result.fMin = fMinimums.rbegin()->first;
1071  if (!fMaximums.empty())
1072  result.fMax = fMaximums.begin()->first;
1073  return result;
1074  }
1075 
1076  // Returns the number of branches with at least one baskets registered.
1077  UInt_t BranchesRegistered()
1078  {
1079  UInt_t result = 0;
1080  for (const auto &r : fRanges) {
1081  if (r.fMin != -1 && r.fMax != -1)
1082  ++result;
1083  }
1084  return result;
1085  }
1086 
1087  // Returns true if at least one of the branch's range contains
1088  // the entry.
1089  Bool_t Contains(Long64_t entry)
1090  {
1091  for (const auto &r : fRanges) {
1092  if (r.fMin != -1 && r.fMax != -1)
1093  if (r.fMin <= entry && entry <= r.fMax)
1094  return kTRUE;
1095  }
1096  return kFALSE;
1097  }
1098 
1099  void Print()
1100  {
1101  for (size_t i = 0; i < fRanges.size(); ++i) {
1102  if (fRanges[i].fMin != -1 || fRanges[i].fMax != -1)
1103  Printf("Range #%zu : %lld to %lld", i, fRanges[i].fMin, fRanges[i].fMax);
1104  }
1105  }
1106 };
1107 } // Anonymous namespace.
1108 
1109 ////////////////////////////////////////////////////////////////////////////////
1110 /// Fill the cache buffer with the branches in the cache.
1111 
1113 {
1114 
1115  if (fNbranches <= 0) return kFALSE;
1117  Long64_t entry = tree->GetReadEntry();
1118  Long64_t fEntryCurrentMax = 0;
1119 
1120  if (entry != -1 && (entry < fEntryMin || fEntryMax < entry))
1121  return kFALSE;
1122 
1123  if (fEnablePrefetching) { // Prefetching mode
1124  if (fIsLearning) { // Learning mode
1125  if (fEntryNext >= 0 && entry >= fEntryNext) {
1126  // entry is outside the learn range, need to stop the learning
1127  // phase. Doing so may trigger a recursive call to FillBuffer in
1128  // the process of filling both prefetching buffers
1130  fIsManual = kFALSE;
1131  }
1132  }
1133  if (fIsLearning) { // Learning mode
1134  // The learning phase should start from the minimum entry in the cache
1135  entry = fEntryMin;
1136  }
1137  if (fFirstTime) {
1138  //try to detect if it is normal or reverse read
1139  fFirstEntry = entry;
1140  }
1141  else {
1142  if (fFirstEntry == entry) return kFALSE;
1143  // Set the read direction
1144  if (!fReadDirectionSet) {
1145  if (entry < fFirstEntry) {
1146  fReverseRead = kTRUE;
1148  }
1149  else if (entry > fFirstEntry) {
1152  }
1153  }
1154 
1155  if (fReverseRead) {
1156  // Reverse reading with prefetching
1157  if (fEntryCurrent >0 && entry < fEntryNext) {
1158  // We can prefetch the next buffer
1159  if (entry >= fEntryCurrent) {
1160  entry = fEntryCurrent - tree->GetAutoFlush() * fFillTimes;
1161  }
1162  if (entry < 0) entry = 0;
1163  }
1164  else if (fEntryCurrent >= 0) {
1165  // We are still reading from the oldest buffer, no need to prefetch a new one
1166  return kFALSE;
1167  }
1168  if (entry < 0) return kFALSE;
1170  }
1171  else {
1172  // Normal reading with prefetching
1173  if (fEnablePrefetching) {
1174  if (entry < 0 && fEntryNext > 0) {
1175  entry = fEntryCurrent;
1176  } else if (entry >= fEntryCurrent) {
1177  if (entry < fEntryNext) {
1178  entry = fEntryNext;
1179  }
1180  }
1181  else {
1182  // We are still reading from the oldest buffer,
1183  // no need to prefetch a new one
1184  return kFALSE;
1185  }
1187  }
1188  }
1189  }
1190  }
1191 
1192  // Set to true to enable all debug output without having to set gDebug
1193  // Replace this once we have a per module and/or per class debugging level/setting.
1194  static constexpr bool showMore = kFALSE;
1195 
1196  static const auto PrintAllCacheInfo = [](TObjArray *branches) {
1197  for (Int_t i = 0; i < branches->GetEntries(); i++) {
1198  TBranch *b = (TBranch *)branches->UncheckedAt(i);
1199  b->PrintCacheInfo();
1200  }
1201  };
1202 
1203  if (showMore || gDebug > 6)
1204  Info("FillBuffer", "***** Called for entry %lld", entry);
1205 
1206  if (!fIsLearning && fEntryCurrent <= entry && entry < fEntryNext) {
1207  // Check if all the basket in the cache have already be used and
1208  // thus we can reuse the cache.
1209  Bool_t allUsed = kTRUE;
1210  for (Int_t i = 0; i < fNbranches; ++i) {
1212  if (!b->fCacheInfo.AllUsed()) {
1213  allUsed = kFALSE;
1214  break;
1215  }
1216  }
1217  if (allUsed) {
1218  fEntryNext = entry;
1219  if (showMore || gDebug > 5)
1220  Info("FillBuffer", "All baskets used already, so refresh the cache early at entry %lld", entry);
1221  }
1222  if (gDebug > 8)
1223  PrintAllCacheInfo(fBranches);
1224  }
1225 
1226  // If the entry is in the range we previously prefetched, there is
1227  // no point in retrying. Note that this will also return false
1228  // during the training phase (fEntryNext is then set intentional to
1229  // the end of the training phase).
1230  if (fEntryCurrent <= entry && entry < fEntryNext) return kFALSE;
1231 
1232  // Triggered by the user, not the learning phase
1233  if (entry == -1)
1234  entry = 0;
1235 
1236  Bool_t resetBranchInfo = kFALSE;
1237  if (entry < fCurrentClusterStart || fNextClusterStart <= entry) {
1238  // We are moving on to another set of clusters.
1239  resetBranchInfo = kTRUE;
1240  if (showMore || gDebug > 6)
1241  Info("FillBuffer", "*** Will reset the branch information about baskets");
1242  } else if (showMore || gDebug > 6) {
1243  Info("FillBuffer", "*** Info we have on the set of baskets");
1244  PrintAllCacheInfo(fBranches);
1245  }
1246 
1247  fEntryCurrentMax = fEntryCurrent;
1248  TTree::TClusterIterator clusterIter = tree->GetClusterIterator(entry);
1249 
1250  auto entryCurrent = clusterIter();
1251  auto entryNext = clusterIter.GetNextEntry();
1252 
1253  if (entryNext < fEntryMin || fEntryMax < entryCurrent) {
1254  // There is no overlap between the cluster we found [entryCurrent, entryNext[
1255  // and the authorized range [fEntryMin, fEntryMax]
1256  // so we have nothing to do
1257  return kFALSE;
1258  }
1259 
1260  // If there is overlap between the found cluster and the authorized range
1261  // update the cache data members with the information about the current cluster.
1262  fEntryCurrent = entryCurrent;
1263  fEntryNext = entryNext;
1264 
1265  auto firstClusterEnd = fEntryNext;
1266  if (showMore || gDebug > 6)
1267  Info("FillBuffer", "Looking at cluster spanning from %lld to %lld", fEntryCurrent, fEntryNext);
1268 
1270  if (fEntryMax <= 0) fEntryMax = tree->GetEntries();
1272 
1273  if ( fEnablePrefetching ) {
1274  if ( entry == fEntryMax ) {
1275  // We are at the end, no need to do anything else
1276  return kFALSE;
1277  }
1278  }
1279 
1280  if (resetBranchInfo) {
1281  // We earlier thought we were onto the next set of clusters.
1282  if (fCurrentClusterStart != -1 || fNextClusterStart != -1) {
1283  if (!(fEntryCurrent < fCurrentClusterStart || fEntryCurrent >= fNextClusterStart)) {
1284  Error("FillBuffer", "Inconsistency: fCurrentClusterStart=%lld fEntryCurrent=%lld fNextClusterStart=%lld "
1285  "but fEntryCurrent should not be in between the two",
1287  }
1288  }
1289 
1290  // Start the next cluster set.
1292  fNextClusterStart = firstClusterEnd;
1293  }
1294 
1295  // Check if owner has a TEventList set. If yes we optimize for this
1296  // Special case reading only the baskets containing entries in the
1297  // list.
1298  TEventList *elist = fTree->GetEventList();
1299  Long64_t chainOffset = 0;
1300  if (elist) {
1301  if (fTree->IsA() ==TChain::Class()) {
1302  TChain *chain = (TChain*)fTree;
1303  Int_t t = chain->GetTreeNumber();
1304  chainOffset = chain->GetTreeOffset()[t];
1305  }
1306  }
1307 
1308  //clear cache buffer
1309  Int_t ntotCurrentBuf = 0;
1310  if (fEnablePrefetching){ //prefetching mode
1311  if (fFirstBuffer) {
1313  ntotCurrentBuf = fNtot;
1314  }
1315  else {
1317  ntotCurrentBuf = fBNtot;
1318  }
1319  }
1320  else {
1322  ntotCurrentBuf = fNtot;
1323  }
1324 
1325  //store baskets
1326  BasketRanges ranges((showMore || gDebug > 6) ? fNbranches : 0);
1327  BasketRanges reqRanges(fNbranches);
1328  BasketRanges memRanges((showMore || gDebug > 6) ? fNbranches : 0);
1329  Int_t clusterIterations = 0;
1330  Long64_t minEntry = fEntryCurrent;
1331  Int_t prevNtot;
1332  Long64_t maxReadEntry = minEntry; // If we are stopped before the end of the 2nd pass, this marker will where we need to start next time.
1333  Int_t nReadPrefRequest = 0;
1334  auto perfStats = GetTree()->GetPerfStats();
1335  do {
1336  prevNtot = ntotCurrentBuf;
1337  Long64_t lowestMaxEntry = fEntryMax; // The lowest maximum entry in the TTreeCache for each branch for each pass.
1338 
1339  struct collectionInfo {
1340  Int_t fClusterStart{-1}; // First basket belonging to the current cluster
1341  Int_t fCurrent{0}; // Currently visited basket
1342  Bool_t fLoadedOnce{kFALSE};
1343 
1344  void Rewind() { fCurrent = (fClusterStart >= 0) ? fClusterStart : 0; }
1345  };
1346  std::vector<collectionInfo> cursor(fNbranches);
1347  Bool_t reachedEnd = kFALSE;
1348  Bool_t skippedFirst = kFALSE;
1349  Bool_t oncePerBranch = kFALSE;
1350  Int_t nDistinctLoad = 0;
1351  Bool_t progress = kTRUE;
1352  enum ENarrow {
1353  kFull = 0,
1354  kNarrow = 1
1355  };
1356  enum EPass {
1357  kStart = 1,
1358  kRegular = 2,
1359  kRewind = 3
1360  };
1361 
1362  auto CollectBaskets = [this, elist, chainOffset, entry, clusterIterations, resetBranchInfo, perfStats,
1363  &cursor, &lowestMaxEntry, &maxReadEntry, &minEntry,
1364  &reachedEnd, &skippedFirst, &oncePerBranch, &nDistinctLoad, &progress,
1365  &ranges, &memRanges, &reqRanges,
1366  &ntotCurrentBuf, &nReadPrefRequest](EPass pass, ENarrow narrow, Long64_t maxCollectEntry) {
1367  // The first pass we add one basket per branches around the requested entry
1368  // then in the second pass we add the other baskets of the cluster.
1369  // This is to support the case where the cache is too small to hold a full cluster.
1370  Int_t nReachedEnd = 0;
1371  Int_t nSkipped = 0;
1372  auto oldnReadPrefRequest = nReadPrefRequest;
1373  std::vector<Int_t> potentialVetoes;
1374 
1375  if (showMore || gDebug > 7)
1376  Info("CollectBaskets", "Called with pass=%d narrow=%d maxCollectEntry=%lld", pass, narrow, maxCollectEntry);
1377 
1378  Bool_t filled = kFALSE;
1379  for (Int_t i = 0; i < fNbranches; ++i) {
1381  if (b->GetDirectory()==0 || b->TestBit(TBranch::kDoNotProcess))
1382  continue;
1383  if (b->GetDirectory()->GetFile() != fFile)
1384  continue;
1385  potentialVetoes.clear();
1386  if (pass == kStart && !cursor[i].fLoadedOnce && resetBranchInfo) {
1387  // First check if we have any cluster that is currently in the
1388  // cache but was not used and would be reloaded in the next
1389  // cluster.
1390  b->fCacheInfo.GetUnused(potentialVetoes);
1391  if (showMore || gDebug > 7) {
1392  TString vetolist;
1393  for(auto v : potentialVetoes) {
1394  vetolist += v;
1395  vetolist.Append(' ');
1396  }
1397  if (!potentialVetoes.empty())
1398  Info("FillBuffer", "*** Potential Vetos for branch #%d: %s", i, vetolist.Data());
1399  }
1400  b->fCacheInfo.Reset();
1401  }
1402  Int_t nb = b->GetMaxBaskets();
1403  Int_t *lbaskets = b->GetBasketBytes();
1404  Long64_t *entries = b->GetBasketEntry();
1405  if (!lbaskets || !entries)
1406  continue;
1407  //we have found the branch. We now register all its baskets
1408  // from the requested offset to the basket below fEntryMax
1409  Int_t blistsize = b->GetListOfBaskets()->GetSize();
1410 
1411  auto maxOfBasket = [this, nb, entries](int j) {
1412  return ((j < (nb - 1)) ? (entries[j + 1] - 1) : fEntryMax - 1);
1413  };
1414 
1415  if (pass == kRewind)
1416  cursor[i].Rewind();
1417  for (auto &j = cursor[i].fCurrent; j < nb; j++) {
1418  // This basket has already been read, skip it
1419 
1420  if (j < blistsize && b->GetListOfBaskets()->UncheckedAt(j)) {
1421 
1422  if (showMore || gDebug > 6) {
1423  ranges.Update(i, entries[j], maxOfBasket(j));
1424  memRanges.Update(i, entries[j], maxOfBasket(j));
1425  }
1426  if (entries[j] <= entry && entry <= maxOfBasket(j)) {
1427  b->fCacheInfo.SetIsInCache(j);
1428  b->fCacheInfo.SetUsed(j);
1429  if (narrow) {
1430  // In narrow mode, we would select 'only' this basket,
1431  // so we are done for this round, let's 'consume' this
1432  // basket and go.
1433  ++nReachedEnd;
1434  ++j;
1435  break;
1436  }
1437  }
1438  continue;
1439  }
1440 
1441  // Important: do not try to read maxCollectEntry, otherwise we might jump to the next autoflush
1442  if (entries[j] >= maxCollectEntry) {
1443  ++nReachedEnd;
1444  break; // break out of the for each branch loop.
1445  }
1446 
1447  Long64_t pos = b->GetBasketSeek(j);
1448  Int_t len = lbaskets[j];
1449  if (pos <= 0 || len <= 0)
1450  continue;
1451  if (len > fBufferSizeMin) {
1452  // Do not cache a basket if it is bigger than the cache size!
1453  if ((showMore || gDebug > 7) &&
1454  (!(entries[j] < minEntry && (j < nb - 1 && entries[j + 1] <= minEntry))))
1455  Info("FillBuffer", "Skipping branch %s basket %d is too large for the cache: %d > %d",
1456  b->GetName(), j, len, fBufferSizeMin);
1457  continue;
1458  }
1459 
1460  if (nReadPrefRequest && entries[j] > (reqRanges.AllIncludedRange().fMax + 1)) {
1461  // There is a gap between this basket and the max of the 'lowest' already loaded basket
1462  // If we are tight in memory, reading this basket may prevent reading the basket (for the other branches)
1463  // that covers this gap, forcing those baskets to be read uncached (because the cache wont be reloaded
1464  // until we use this basket).
1465  // eg. We could end up with the cache contain
1466  // b1: [428, 514[ // 'this' basket and we can assume [321 to 428[ is already in memory
1467  // b2: [400, 424[
1468  // and when reading entry 425 we will read b2's basket uncached.
1469 
1470  if (showMore || gDebug > 8)
1471  Info("FillBuffer", "Skipping for now due to gap %d/%d with %lld > %lld", i, j, entries[j],
1472  (reqRanges.AllIncludedRange().fMax + 1));
1473  break; // Without consuming the basket.
1474  }
1475 
1476  if (entries[j] < minEntry && (j<nb-1 && entries[j+1] <= minEntry))
1477  continue;
1478 
1479  // We are within the range
1480  if (cursor[i].fClusterStart == -1)
1481  cursor[i].fClusterStart = j;
1482 
1483  if (elist) {
1484  Long64_t emax = fEntryMax;
1485  if (j<nb-1)
1486  emax = entries[j + 1] - 1;
1487  if (!elist->ContainsRange(entries[j]+chainOffset,emax+chainOffset))
1488  continue;
1489  }
1490 
1491  if (b->fCacheInfo.HasBeenUsed(j) || b->fCacheInfo.IsInCache(j) || b->fCacheInfo.IsVetoed(j)) {
1492  // We already cached and used this basket during this cluster range,
1493  // let's not redo it
1494  if (showMore || gDebug > 7)
1495  Info("FillBuffer", "Skipping basket to avoid redo: %d/%d veto: %d", i, j, b->fCacheInfo.IsVetoed(j));
1496  continue;
1497  }
1498 
1499  if (std::find(std::begin(potentialVetoes), std::end(potentialVetoes), j) != std::end(potentialVetoes)) {
1500  // This basket was in the previous cache/cluster and was not used,
1501  // let's not read it again. I.e. we bet that it will continue to not
1502  // be used. At worst it will be used and thus read by itself.
1503  // Usually in this situation the basket is large so the penalty for
1504  // (re)reading it uselessly is high and the penalty to read it by
1505  // itself is 'small' (i.e. size bigger than latency).
1506  b->fCacheInfo.Veto(j);
1507  if (showMore || gDebug > 7)
1508  Info("FillBuffer", "Veto-ing cluster %d [%lld,%lld[ in branch %s #%d", j, entries[j],
1509  maxOfBasket(j) + 1, b->GetName(), i);
1510  continue;
1511  }
1512 
1513  if (narrow) {
1514  if ((((entries[j] > entry)) || (j < nb - 1 && entries[j + 1] <= entry))) {
1515  // Keep only the basket that contains the entry
1516  if (j == cursor[i].fClusterStart && entry > entries[j])
1517  ++nSkipped;
1518  if (entries[j] > entry)
1519  break;
1520  else
1521  continue;
1522  }
1523  }
1524 
1525  if ((ntotCurrentBuf + len) > fBufferSizeMin) {
1526  // Humm ... we are going to go over the requested size.
1527  if (clusterIterations > 0 && cursor[i].fLoadedOnce) {
1528  // We already have a full cluster and now we would go over the requested
1529  // size, let's stop caching (and make sure we start next time from the
1530  // end of the previous cluster).
1531  if (showMore || gDebug > 5) {
1532  Info(
1533  "FillBuffer",
1534  "Breaking early because %d is greater than %d at cluster iteration %d will restart at %lld",
1535  (ntotCurrentBuf + len), fBufferSizeMin, clusterIterations, minEntry);
1536  }
1537  fEntryNext = minEntry;
1538  filled = kTRUE;
1539  break;
1540  } else {
1541  if (pass == kStart || !cursor[i].fLoadedOnce) {
1542  if ((ntotCurrentBuf + len) > 4 * fBufferSizeMin) {
1543  // Okay, so we have not even made one pass and we already have
1544  // accumulated request for more than twice the memory size ...
1545  // So stop for now, and will restart at the same point, hoping
1546  // that the basket will still be in memory and not asked again ..
1547  fEntryNext = maxReadEntry;
1548 
1549  if (showMore || gDebug > 5) {
1550  Info("FillBuffer", "Breaking early because %d is greater than 4*%d at cluster iteration "
1551  "%d pass %d will restart at %lld",
1552  (ntotCurrentBuf + len), fBufferSizeMin, clusterIterations, pass, fEntryNext);
1553  }
1554  filled = kTRUE;
1555  break;
1556  }
1557  } else {
1558  // We have made one pass through the branches and thus already
1559  // requested one basket per branch, let's stop prefetching
1560  // now.
1561  if ((ntotCurrentBuf + len) > 2 * fBufferSizeMin) {
1562  fEntryNext = maxReadEntry;
1563  if (showMore || gDebug > 5) {
1564  Info("FillBuffer", "Breaking early because %d is greater than 2*%d at cluster iteration "
1565  "%d pass %d will restart at %lld",
1566  (ntotCurrentBuf + len), fBufferSizeMin, clusterIterations, pass, fEntryNext);
1567  }
1568  filled = kTRUE;
1569  break;
1570  }
1571  }
1572  }
1573  }
1574 
1575  ++nReadPrefRequest;
1576 
1577  reqRanges.Update(i, j, entries, nb, fEntryMax);
1578  if (showMore || gDebug > 6)
1579  ranges.Update(i, j, entries, nb, fEntryMax);
1580 
1581  b->fCacheInfo.SetIsInCache(j);
1582 
1583  if (showMore || gDebug > 6)
1584  Info("FillBuffer", "*** Registering branch %d basket %d %s", i, j, b->GetName());
1585 
1586  if (!cursor[i].fLoadedOnce) {
1587  cursor[i].fLoadedOnce = kTRUE;
1588  ++nDistinctLoad;
1589  }
1590  if (R__unlikely(perfStats)) {
1591  perfStats->SetLoaded(i, j);
1592  }
1593 
1594  // Actual registering the basket for loading from the file.
1595  if (fEnablePrefetching){
1596  if (fFirstBuffer) {
1597  TFileCacheRead::Prefetch(pos,len);
1598  ntotCurrentBuf = fNtot;
1599  }
1600  else {
1602  ntotCurrentBuf = fBNtot;
1603  }
1604  }
1605  else {
1606  TFileCacheRead::Prefetch(pos,len);
1607  ntotCurrentBuf = fNtot;
1608  }
1609 
1610  if ( ( j < (nb-1) ) && entries[j+1] > maxReadEntry ) {
1611  // Info("FillBuffer","maxCollectEntry incremented from %lld to %lld", maxReadEntry, entries[j+1]);
1612  maxReadEntry = entries[j+1];
1613  }
1614  if (ntotCurrentBuf > 4 * fBufferSizeMin) {
1615  // Humm something wrong happened.
1616  Warning("FillBuffer", "There is more data in this cluster (starting at entry %lld to %lld, "
1617  "current=%lld) than usual ... with %d %.3f%% of the branches we already have "
1618  "%d bytes (instead of %d)",
1619  fEntryCurrent, fEntryNext, entries[j], i, (100.0 * i) / ((float)fNbranches), ntotCurrentBuf,
1620  fBufferSizeMin);
1621  }
1622  if (pass == kStart) {
1623  // In the first pass, we record one basket per branch and move on to the next branch.
1624  auto high = maxOfBasket(j);
1625  if (high < lowestMaxEntry)
1626  lowestMaxEntry = high;
1627  // 'Consume' the baskets (i.e. avoid looking at it during a subsequent pass)
1628  ++j;
1629  break;
1630  } else if ((j + 1) == nb || entries[j + 1] >= maxReadEntry || entries[j + 1] >= lowestMaxEntry) {
1631  // In the other pass, load the baskets until we get to the maximum loaded so far.
1632  auto high = maxOfBasket(j);
1633  if (high < lowestMaxEntry)
1634  lowestMaxEntry = high;
1635  // 'Consume' the baskets (i.e. avoid looking at it during a subsequent pass)
1636  ++j;
1637  break;
1638  }
1639  }
1640 
1641  if (cursor[i].fCurrent == nb) {
1642  ++nReachedEnd;
1643  }
1644 
1645  if (gDebug > 0)
1646  Info("CollectBaskets",
1647  "Entry: %lld, registering baskets branch %s, fEntryNext=%lld, fNseek=%d, ntotCurrentBuf=%d",
1648  minEntry, ((TBranch *)fBranches->UncheckedAt(i))->GetName(), fEntryNext, fNseek, ntotCurrentBuf);
1649  }
1650  reachedEnd = (nReachedEnd == fNbranches);
1651  skippedFirst = (nSkipped > 0);
1652  oncePerBranch = (nDistinctLoad == fNbranches);
1653  progress = nReadPrefRequest - oldnReadPrefRequest;
1654  return filled;
1655  };
1656 
1657  // First collect all the basket containing the request entry.
1658  bool full = kFALSE;
1659 
1660  full = CollectBaskets(kStart, kNarrow, fEntryNext);
1661 
1662  // Then fill out from all but the 'largest' branch to even out
1663  // the range across branches;
1664  while (!full && !reachedEnd && progress) { // used to be restricted to !oncePerBranch
1665  full = CollectBaskets(kStart, kFull, std::min(maxReadEntry, fEntryNext));
1666  }
1667 
1668  resetBranchInfo = kFALSE; // Make sure the 2nd cluster iteration does not erase the info.
1669 
1670  // Then fill out to the end of the cluster.
1671  if (!full && !fReverseRead) {
1672  do {
1673  full = CollectBaskets(kRegular, kFull, fEntryNext);
1674  } while (!full && !reachedEnd && progress);
1675  }
1676 
1677  // The restart from the start of the cluster.
1678  if (!full && skippedFirst) {
1679  full = CollectBaskets(kRewind, kFull, fEntryNext);
1680  while (!full && !reachedEnd && progress) {
1681  full = CollectBaskets(kRegular, kFull, fEntryNext);
1682  }
1683  }
1684 
1685  clusterIterations++;
1686 
1687  minEntry = clusterIter.Next();
1688  if (fIsLearning) {
1689  fFillTimes++;
1690  }
1691 
1692  // Continue as long as we still make progress (prevNtot < ntotCurrentBuf), that the next entry range to be looked
1693  // at,
1694  // which start at 'minEntry', is not past the end of the requested range (minEntry < fEntryMax)
1695  // and we guess that we not going to go over the requested amount of memory by asking for another set
1696  // of entries (fBufferSizeMin > ((Long64_t)ntotCurrentBuf*(clusterIterations+1))/clusterIterations).
1697  // ntotCurrentBuf / clusterIterations is the average size we are accumulated so far at each loop.
1698  // and thus (ntotCurrentBuf / clusterIterations) * (clusterIterations+1) is a good guess at what the next total
1699  // size
1700  // would be if we run the loop one more time. ntotCurrentBuf and clusterIterations are Int_t but can sometimes
1701  // be 'large' (i.e. 30Mb * 300 intervals) and can overflow the numerical limit of Int_t (i.e. become
1702  // artificially negative). To avoid this issue we promote ntotCurrentBuf to a long long (64 bits rather than 32
1703  // bits)
1704  if (!((fBufferSizeMin > ((Long64_t)ntotCurrentBuf * (clusterIterations + 1)) / clusterIterations) &&
1705  (prevNtot < ntotCurrentBuf) && (minEntry < fEntryMax))) {
1706  if (showMore || gDebug > 6)
1707  Info("FillBuffer", "Breaking because %d <= %lld || (%d >= %d) || %lld >= %lld", fBufferSizeMin,
1708  ((Long64_t)ntotCurrentBuf * (clusterIterations + 1)) / clusterIterations, prevNtot, ntotCurrentBuf,
1709  minEntry, fEntryMax);
1710  break;
1711  }
1712 
1713  //for the reverse reading case
1714  if (!fIsLearning && fReverseRead) {
1715  if (clusterIterations >= fFillTimes)
1716  break;
1717  if (minEntry >= fEntryCurrentMax && fEntryCurrentMax > 0)
1718  break;
1719  }
1720  fEntryNext = clusterIter.GetNextEntry();
1723  } while (kTRUE);
1724 
1725  if (showMore || gDebug > 6) {
1726  Info("FillBuffer", "Mem ranges");
1727  memRanges.Print();
1728  Info("FillBuffer", "Combined ranges");
1729  ranges.Print();
1730  Info("FillBuffer", "Requested ranges");
1731  reqRanges.Print();
1732  PrintAllCacheInfo(fBranches);
1733  }
1734 
1735  if (nReadPrefRequest == 0) {
1736  // Nothing was added in the cache. This usually indicates that the baskets
1737  // contains the requested entry are either already in memory or are too large
1738  // on their own to fit in the cache.
1739  if (showMore || gDebug > 5) {
1740  Info("FillBuffer", "For entry %lld, nothing was added to the cache.", entry);
1741  }
1742  } else if (fEntryNext < firstClusterEnd && !reqRanges.Contains(entry)) {
1743  // Something went very wrong and even-though we searched for the baskets
1744  // holding 'entry' we somehow ended up with a range of entries that does
1745  // validate. So we must have been unable to find or fit the needed basket.
1746  // And thus even-though, we know the corresponding baskets wont be in the cache,
1747  // Let's make it official that 'entry' is within the range of this TTreeCache ('s search.)
1748 
1749  // Without this, the next read will be flagged as 'out-of-range' and then we start at
1750  // the exact same point as this FillBuffer execution resulting in both the requested
1751  // entry still not being part of the cache **and** the beginning of the cluster being
1752  // read **again**.
1753 
1754  if (showMore || gDebug > 5) {
1755  Error("FillBuffer", "Reset the next entry because the currently loaded range does not contains the request "
1756  "entry: %lld. fEntryNext updated from %lld to %lld. %d",
1757  entry, fEntryNext, firstClusterEnd, nReadPrefRequest);
1758  reqRanges.Print();
1759  }
1760 
1761  fEntryNext = firstClusterEnd;
1762  } else {
1763  if (showMore || gDebug > 5) {
1764  Info("FillBuffer", "Complete adding %d baskets from %d branches taking in memory %d out of %d",
1765  nReadPrefRequest, reqRanges.BranchesRegistered(), ntotCurrentBuf, fBufferSizeMin);
1766  }
1767  }
1768 
1769  fNReadPref += nReadPrefRequest;
1770  if (fEnablePrefetching) {
1771  if (fIsLearning) {
1773  }
1774  if (!fIsLearning && fFirstTime){
1775  // First time we add autoFlush entries , after fFillTimes * autoFlush
1776  // only in reverse prefetching mode
1777  fFirstTime = kFALSE;
1778  }
1779  }
1780  fIsLearning = kFALSE;
1781  return kTRUE;
1782 }
1783 
1784 ////////////////////////////////////////////////////////////////////////////////
1785 /// Return the desired prefill type from the environment or resource variable
1786 /// - 0 - No prefill
1787 /// - 1 - All branches
1788 
1790 {
1791  const char *stcp;
1792  Int_t s = 0;
1793 
1794  if (!(stcp = gSystem->Getenv("ROOT_TTREECACHE_PREFILL")) || !*stcp) {
1795  s = gEnv->GetValue("TTreeCache.Prefill", 1);
1796  } else {
1797  s = TString(stcp).Atoi();
1798  }
1799 
1800  return static_cast<TTreeCache::EPrefillType>(s);
1801 }
1802 
1803 ////////////////////////////////////////////////////////////////////////////////
1804 /// Give the total efficiency of the primary cache... defined as the ratio
1805 /// of blocks found in the cache vs. the number of blocks prefetched
1806 /// ( it could be more than 1 if we read the same block from the cache more
1807 /// than once )
1808 ///
1809 /// Note: This should eb used at the end of the processing or we will
1810 /// get incomplete stats
1811 
1813 {
1814  if ( !fNReadPref )
1815  return 0;
1816 
1817  return ((Double_t)fNReadOk / (Double_t)fNReadPref);
1818 }
1819 
1820 ////////////////////////////////////////////////////////////////////////////////
1821 /// The total efficiency of the 'miss cache' - defined as the ratio
1822 /// of blocks found in the cache versus the number of blocks prefetched
1823 
1825 {
1826  if (!fNMissReadPref) {
1827  return 0;
1828  }
1829  return static_cast<double>(fNMissReadOk) / static_cast<double>(fNMissReadPref);
1830 }
1831 
1832 ////////////////////////////////////////////////////////////////////////////////
1833 /// This will indicate a sort of relative efficiency... a ratio of the
1834 /// reads found in the cache to the number of reads so far
1835 
1837 {
1838  if ( !fNReadOk && !fNReadMiss )
1839  return 0;
1840 
1841  return ((Double_t)fNReadOk / (Double_t)(fNReadOk + fNReadMiss));
1842 }
1843 
1844 ////////////////////////////////////////////////////////////////////////////////
1845 /// Relative efficiency of the 'miss cache' - ratio of the reads found in cache
1846 /// to the number of reads so far.
1847 
1849 {
1850  if (!fNMissReadOk && !fNMissReadMiss) {
1851  return 0;
1852  }
1853 
1854  return static_cast<double>(fNMissReadOk) / static_cast<double>(fNMissReadOk + fNMissReadMiss);
1855 }
1856 
1857 ////////////////////////////////////////////////////////////////////////////////
1858 /// Static function returning the number of entries used to train the cache
1859 /// see SetLearnEntries
1860 
1862 {
1863  return fgLearnEntries;
1864 }
1865 
1866 ////////////////////////////////////////////////////////////////////////////////
1867 /// Print cache statistics. Like:
1868 ///
1869 /// ~~~ {.cpp}
1870 /// ******TreeCache statistics for file: cms2.root ******
1871 /// Number of branches in the cache ...: 1093
1872 /// Cache Efficiency ..................: 0.997372
1873 /// Cache Efficiency Rel...............: 1.000000
1874 /// Learn entries......................: 100
1875 /// Reading............................: 72761843 bytes in 7 transactions
1876 /// Readahead..........................: 256000 bytes with overhead = 0 bytes
1877 /// Average transaction................: 10394.549000 Kbytes
1878 /// Number of blocks in current cache..: 210, total size: 6280352
1879 /// ~~~
1880 ///
1881 /// - if option = "a" the list of blocks in the cache is printed
1882 /// see also class TTreePerfStats.
1883 /// - if option contains 'cachedbranches', the list of branches being
1884 /// cached is printed.
1885 
1886 void TTreeCache::Print(Option_t *option) const
1887 {
1888  TString opt = option;
1889  opt.ToLower();
1890  printf("******TreeCache statistics for tree: %s in file: %s ******\n",fTree ? fTree->GetName() : "no tree set",fFile ? fFile->GetName() : "no file set");
1891  if (fNbranches <= 0) return;
1892  printf("Number of branches in the cache ...: %d\n",fNbranches);
1893  printf("Cache Efficiency ..................: %f\n",GetEfficiency());
1894  printf("Cache Efficiency Rel...............: %f\n",GetEfficiencyRel());
1895  printf("Secondary Efficiency ..............: %f\n", GetMissEfficiency());
1896  printf("Secondary Efficiency Rel ..........: %f\n", GetMissEfficiencyRel());
1897  printf("Learn entries......................: %d\n",TTreeCache::GetLearnEntries());
1898  if ( opt.Contains("cachedbranches") ) {
1899  opt.ReplaceAll("cachedbranches","");
1900  printf("Cached branches....................:\n");
1901  const TObjArray *cachedBranches = this->GetCachedBranches();
1902  Int_t nbranches = cachedBranches->GetEntriesFast();
1903  for (Int_t i = 0; i < nbranches; ++i) {
1904  TBranch* branch = (TBranch*) cachedBranches->UncheckedAt(i);
1905  printf("Branch name........................: %s\n",branch->GetName());
1906  }
1907  }
1908  TFileCacheRead::Print(opt);
1909 }
1910 
1911 ////////////////////////////////////////////////////////////////////////////////
1912 /// Old method ReadBuffer before the addition of the prefetch mechanism.
1913 
1915  //Is request already in the cache?
1916  if (TFileCacheRead::ReadBuffer(buf,pos,len) == 1){
1917  fNReadOk++;
1918  return 1;
1919  }
1920 
1921  static const auto recordMiss = [](TVirtualPerfStats *perfStats, TObjArray *branches, Bool_t bufferFilled,
1922  Long64_t basketpos) {
1923  if (gDebug > 6)
1924  ::Info("TTreeCache::ReadBufferNormal", "Cache miss after an %s FillBuffer: pos=%lld",
1925  bufferFilled ? "active" : "inactive", basketpos);
1926  for (Int_t i = 0; i < branches->GetEntries(); ++i) {
1927  TBranch *b = (TBranch *)branches->UncheckedAt(i);
1928  Int_t blistsize = b->GetListOfBaskets()->GetSize();
1929  for (Int_t j = 0; j < blistsize; ++j) {
1930  if (basketpos == b->GetBasketSeek(j)) {
1931  if (gDebug > 6)
1932  ::Info("TTreeCache::ReadBufferNormal", " Missing basket: %d for %s", j, b->GetName());
1933  perfStats->SetMissed(i, j);
1934  }
1935  }
1936  }
1937  };
1938 
1939  //not found in cache. Do we need to fill the cache?
1940  Bool_t bufferFilled = FillBuffer();
1941  if (bufferFilled) {
1942  Int_t res = TFileCacheRead::ReadBuffer(buf,pos,len);
1943 
1944  if (res == 1)
1945  fNReadOk++;
1946  else if (res == 0) {
1947  fNReadMiss++;
1948  auto perfStats = GetTree()->GetPerfStats();
1949  if (perfStats)
1950  recordMiss(perfStats, fBranches, bufferFilled, pos);
1951  }
1952 
1953  return res;
1954  }
1955 
1956  if (CheckMissCache(buf, pos, len)) {
1957  return 1;
1958  }
1959 
1960  fNReadMiss++;
1961  auto perfStats = GetTree()->GetPerfStats();
1962  if (perfStats)
1963  recordMiss(perfStats, fBranches, bufferFilled, pos);
1964 
1965  return 0;
1966 }
1967 
1968 ////////////////////////////////////////////////////////////////////////////////
1969 /// Used to read a chunk from a block previously fetched. It will call FillBuffer
1970 /// even if the cache lookup succeeds, because it will try to prefetch the next block
1971 /// as soon as we start reading from the current block.
1972 
1974 {
1975  if (TFileCacheRead::ReadBuffer(buf, pos, len) == 1){
1976  //call FillBuffer to prefetch next block if necessary
1977  //(if we are currently reading from the last block available)
1978  FillBuffer();
1979  fNReadOk++;
1980  return 1;
1981  }
1982 
1983  //keep on prefetching until request is satisfied
1984  // try to prefetch a couple of times and if request is still not satisfied then
1985  // fall back to normal reading without prefetching for the current request
1986  Int_t counter = 0;
1987  while (1) {
1988  if(TFileCacheRead::ReadBuffer(buf, pos, len)) {
1989  break;
1990  }
1991  FillBuffer();
1992  fNReadMiss++;
1993  counter++;
1994  if (counter>1) {
1995  return 0;
1996  }
1997  }
1998 
1999  fNReadOk++;
2000  return 1;
2001 }
2002 
2003 ////////////////////////////////////////////////////////////////////////////////
2004 /// Read buffer at position pos if the request is in the list of
2005 /// prefetched blocks read from fBuffer.
2006 /// Otherwise try to fill the cache from the list of selected branches,
2007 /// and recheck if pos is now in the list.
2008 /// Returns:
2009 /// - -1 in case of read failure,
2010 /// - 0 in case not in cache,
2011 /// - 1 in case read from cache.
2012 /// This function overloads TFileCacheRead::ReadBuffer.
2013 
2015 {
2016  if (!fEnabled) return 0;
2017 
2018  if (fEnablePrefetching)
2019  return TTreeCache::ReadBufferPrefetch(buf, pos, len);
2020  else
2021  return TTreeCache::ReadBufferNormal(buf, pos, len);
2022 }
2023 
2024 ////////////////////////////////////////////////////////////////////////////////
2025 /// This will simply clear the cache
2026 
2028 {
2029  for (Int_t i = 0; i < fNbranches; ++i) {
2031  if (b->GetDirectory()==0 || b->TestBit(TBranch::kDoNotProcess))
2032  continue;
2033  if (b->GetDirectory()->GetFile() != fFile)
2034  continue;
2035  b->fCacheInfo.Reset();
2036  }
2037  fEntryCurrent = -1;
2038  fEntryNext = -1;
2039  fCurrentClusterStart = -1;
2040  fNextClusterStart = -1;
2041 
2043 
2044  if (fEnablePrefetching) {
2045  fFirstTime = kTRUE;
2047  }
2048 }
2049 
2050 ////////////////////////////////////////////////////////////////////////////////
2051 /// Change the underlying buffer size of the cache.
2052 /// If the change of size means some cache content is lost, or if the buffer
2053 /// is now larger, setup for a cache refill the next time there is a read
2054 /// Returns:
2055 /// - 0 if the buffer content is still available
2056 /// - 1 if some or all of the buffer content has been made unavailable
2057 /// - -1 on error
2058 
2060 {
2061  Int_t prevsize = GetBufferSize();
2062  Int_t res = TFileCacheRead::SetBufferSize(buffersize);
2063  if (res < 0) {
2064  return res;
2065  }
2066 
2067  if (res == 0 && buffersize <= prevsize) {
2068  return res;
2069  }
2070 
2071  // if content was removed from the buffer, or the buffer was enlarged then
2072  // empty the prefetch lists and prime to fill the cache again
2073 
2075  if (fEnablePrefetching) {
2077  }
2078 
2079  fEntryCurrent = -1;
2080  if (!fIsLearning) {
2081  fEntryNext = -1;
2082  }
2083 
2084  return 1;
2085 }
2086 
2087 ////////////////////////////////////////////////////////////////////////////////
2088 /// Set the minimum and maximum entry number to be processed
2089 /// this information helps to optimize the number of baskets to read
2090 /// when prefetching the branch buffers.
2091 
2093 {
2094  // This is called by TTreePlayer::Process in an automatic way...
2095  // don't restart it if the user has specified the branches.
2096  Bool_t needLearningStart = (fEntryMin != emin) && fIsLearning && !fIsManual;
2097 
2098  fEntryMin = emin;
2099  fEntryMax = emax;
2101  if (gDebug > 0)
2102  Info("SetEntryRange", "fEntryMin=%lld, fEntryMax=%lld, fEntryNext=%lld",
2104 
2105  if (needLearningStart) {
2106  // Restart learning
2108  }
2109 }
2110 
2111 ////////////////////////////////////////////////////////////////////////////////
2112 /// Overload to make sure that the object specific
2113 
2115 {
2116  // The infinite recursion is 'broken' by the fact that
2117  // TFile::SetCacheRead remove the entry from fCacheReadMap _before_
2118  // calling SetFile (and also by setting fFile to zero before the calling).
2119  if (fFile) {
2120  TFile *prevFile = fFile;
2121  fFile = 0;
2122  prevFile->SetCacheRead(0, fTree, action);
2123  }
2124  TFileCacheRead::SetFile(file, action);
2125 }
2126 
2127 ////////////////////////////////////////////////////////////////////////////////
2128 /// Static function to set the number of entries to be used in learning mode
2129 /// The default value for n is 10. n must be >= 1
2130 
2132 {
2133  if (n < 1) n = 1;
2134  fgLearnEntries = n;
2135 }
2136 
2137 ////////////////////////////////////////////////////////////////////////////////
2138 /// Set whether the learning period is started with a prefilling of the
2139 /// cache and which type of prefilling is used.
2140 /// The two value currently supported are:
2141 /// - TTreeCache::kNoPrefill disable the prefilling
2142 /// - TTreeCache::kAllBranches fill the cache with baskets from all branches.
2143 /// The default prefilling behavior can be controlled by setting
2144 /// TTreeCache.Prefill or the environment variable ROOT_TTREECACHE_PREFILL.
2145 
2147 {
2148  fPrefillType = type;
2149 }
2150 
2151 ////////////////////////////////////////////////////////////////////////////////
2152 /// The name should be enough to explain the method.
2153 /// The only additional comments is that the cache is cleaned before
2154 /// the new learning phase.
2155 
2157 {
2158  fIsLearning = kTRUE;
2159  fIsManual = kFALSE;
2160  fNbranches = 0;
2161  if (fBrNames) fBrNames->Delete();
2163  fEntryCurrent = -1;
2164 }
2165 
2166 ////////////////////////////////////////////////////////////////////////////////
2167 /// This is the counterpart of StartLearningPhase() and can be used to stop
2168 /// the learning phase. It's useful when the user knows exactly what branches
2169 /// they are going to use.
2170 /// For the moment it's just a call to FillBuffer() since that method
2171 /// will create the buffer lists from the specified branches.
2172 
2174 {
2175  if (fIsLearning) {
2176  // This will force FillBuffer to read the buffers.
2177  fEntryNext = -1;
2178  fIsLearning = kFALSE;
2179  }
2180  fIsManual = kTRUE;
2181 
2182  auto perfStats = GetTree()->GetPerfStats();
2183  if (perfStats)
2184  perfStats->UpdateBranchIndices(fBranches);
2185 
2186  //fill the buffers only once during learning
2187  if (fEnablePrefetching && !fOneTime) {
2188  fIsLearning = kTRUE;
2189  FillBuffer();
2190  fOneTime = kTRUE;
2191  }
2192 }
2193 
2194 ////////////////////////////////////////////////////////////////////////////////
2195 /// Update pointer to current Tree and recompute pointers to the branches in the cache.
2196 
2198 {
2199 
2200  fTree = tree;
2201 
2202  fEntryMin = 0;
2203  fEntryMax = fTree->GetEntries();
2204 
2205  fEntryCurrent = -1;
2206 
2207  if (fBrNames->GetEntries() == 0 && fIsLearning) {
2208  // We still need to learn.
2210  } else {
2211  // We learnt from a previous file.
2212  fIsLearning = kFALSE;
2213  fEntryNext = -1;
2214  }
2215  fNbranches = 0;
2216 
2217  TIter next(fBrNames);
2218  TObjString *os;
2219  while ((os = (TObjString*)next())) {
2220  TBranch *b = fTree->GetBranch(os->GetName());
2221  if (!b) {
2222  continue;
2223  }
2225  fNbranches++;
2226  }
2227 
2228  auto perfStats = GetTree()->GetPerfStats();
2229  if (perfStats)
2230  perfStats->UpdateBranchIndices(fBranches);
2231 }
2232 
2233 ////////////////////////////////////////////////////////////////////////////////
2234 /// Perform an initial prefetch, attempting to read as much of the learning
2235 /// phase baskets for all branches at once
2236 
2238 {
2239  // This is meant for the learning phase
2240  if (!fIsLearning) return;
2241 
2242  // This should be called before reading entries, otherwise we'll
2243  // always exit here, since TBranch adds itself before reading
2244  if (fNbranches > 0) return;
2245 
2246  // Is the LearnPrefill enabled (using an Int_t here to allow for future
2247  // extension to alternative Prefilling).
2248  if (fPrefillType == kNoPrefill) return;
2249 
2250  Long64_t entry = fTree ? fTree->GetReadEntry() : 0;
2251 
2252  // Return early if we are out of the requested range.
2253  if (entry < fEntryMin || entry > fEntryMax) return;
2254 
2256 
2257 
2258  // Force only the learn entries to be cached by temporarily setting min/max
2259  // to the learning phase entry range
2260  // But save all the old values, so we can restore everything to how it was
2261  Long64_t eminOld = fEntryMin;
2262  Long64_t emaxOld = fEntryMax;
2263  Long64_t ecurrentOld = fEntryCurrent;
2264  Long64_t enextOld = fEntryNext;
2265  auto currentClusterStartOld = fCurrentClusterStart;
2266  auto nextClusterStartOld = fNextClusterStart;
2267 
2268  fEntryMin = std::max(fEntryMin, fEntryCurrent);
2269  fEntryMax = std::min(fEntryMax, fEntryNext);
2270 
2271  // We check earlier that we are within the authorized range, but
2272  // we might still be out of the (default) learning range and since
2273  // this is called before any branch is added to the cache, this means
2274  // that the user's first GetEntry is this one which is outside of the
2275  // learning range ... so let's do something sensible-ish.
2276  // Note: we probably should also fix the learning range but we may
2277  // or may not have enough information to know if we can move it
2278  // (for example fEntryMin (eminOld right now) might be the default or user provided)
2279  if (entry < fEntryMin) fEntryMin = entry;
2280  if (entry > fEntryMax) fEntryMax = entry;
2281 
2282  // Add all branches to be cached. This also sets fIsManual, stops learning,
2283  // and makes fEntryNext = -1 (which forces a cache fill, which is good)
2284  AddBranch("*");
2285  fIsManual = kFALSE; // AddBranch sets fIsManual, so we reset it
2286 
2287  // Now, fill the buffer with the learning phase entry range
2288  FillBuffer();
2289 
2290  // Leave everything the way we found it
2291  fIsLearning = kTRUE;
2292  DropBranch("*"); // This doesn't work unless we're already learning
2293 
2294  // Restore entry values
2295  fEntryMin = eminOld;
2296  fEntryMax = emaxOld;
2297  fEntryCurrent = ecurrentOld;
2298  fEntryNext = enextOld;
2299  fCurrentClusterStart = currentClusterStartOld;
2300  fNextClusterStart = nextClusterStartOld;
2301 
2303 }
TTreeCache::StartLearningPhase
void StartLearningPhase()
The name should be enough to explain the method.
Definition: TTreeCache.cxx:2156
TFileCacheRead::Prefetch
virtual void Prefetch(Long64_t pos, Int_t len)
Add block of length len at position pos in the list of blocks to be prefetched.
Definition: TFileCacheRead.cxx:201
TTreeCache::SetFile
virtual void SetFile(TFile *file, TFile::ECacheAction action=TFile::kDisconnect)
Overload to make sure that the object specific.
Definition: TTreeCache.cxx:2114
TTreeCache::EPrefillType
EPrefillType
Definition: TTreeCache.h:35
TEventList::ContainsRange
virtual Bool_t ContainsRange(Long64_t entrymin, Long64_t entrymax)
Return TRUE if list contains entries from entrymin to entrymax included.
Definition: TEventList.cxx:174
TTreeCache::SetBufferSize
virtual Int_t SetBufferSize(Int_t buffersize)
Change the underlying buffer size of the cache.
Definition: TTreeCache.cxx:2059
n
const Int_t n
Definition: legend1.C:16
TVirtualPerfStats.h
ROOT::Math::IntegOptionsUtil::Print
void Print(std::ostream &os, const OptionType &opt)
Definition: IntegratorOptions.cxx:91
TTreeCache::GetConfiguredPrefillType
EPrefillType GetConfiguredPrefillType() const
Return the desired prefill type from the environment or resource variable.
Definition: TTreeCache.cxx:1789
TTreeCache::fEnabled
Bool_t fEnabled
! cache enabled for cached reading
Definition: TTreeCache.h:63
kTRUE
const Bool_t kTRUE
Definition: RtypesCore.h:100
Range
Ta Range(0, 0, 1, 1)
TEventList
A TEventList object is a list of selected events (entries) in a TTree.
Definition: TEventList.h:31
TObjArray
An array of TObjects.
Definition: TObjArray.h:37
TTreeCache::LearnPrefill
virtual void LearnPrefill()
Perform an initial prefetch, attempting to read as much of the learning phase baskets for all branche...
Definition: TTreeCache.cxx:2237
TTreeCache::FillBuffer
virtual Bool_t FillBuffer()
Fill the cache buffer with the branches in the cache.
Definition: TTreeCache.cxx:1112
Option_t
const char Option_t
Definition: RtypesCore.h:66
TTreeCache::fNReadPref
Int_t fNReadPref
Number of blocks that were prefetched.
Definition: TTreeCache.h:49
TString::Atoi
Int_t Atoi() const
Return integer value of string.
Definition: TString.cxx:1943
kNPOS
const Ssiz_t kNPOS
Definition: RtypesCore.h:124
TCollection::GetEntries
virtual Int_t GetEntries() const
Definition: TCollection.h:177
TBranchElement
A Branch for the case of an object.
Definition: TBranchElement.h:39
TTree::GetListOfFriends
virtual TList * GetListOfFriends() const
Definition: TTree.h:485
TChain::GetTreeOffset
Long64_t * GetTreeOffset() const
Definition: TChain.h:117
gEnv
R__EXTERN TEnv * gEnv
Definition: TEnv.h:171
TList::FindObject
virtual TObject * FindObject(const char *name) const
Find an object in this list using its name.
Definition: TList.cxx:578
TTreeCache::fNMissReadPref
Int_t fNMissReadPref
Number of blocks read into the secondary ("miss") cache.
Definition: TTreeCache.h:50
TList::Delete
virtual void Delete(Option_t *option="")
Remove all objects from the list AND delete all heap based objects.
Definition: TList.cxx:470
TString::Data
const char * Data() const
Definition: TString.h:369
TTreeCache::ReadBuffer
virtual Int_t ReadBuffer(char *buf, Long64_t pos, Int_t len)
Read buffer at position pos if the request is in the list of prefetched blocks read from fBuffer.
Definition: TTreeCache.cxx:2014
TFileCacheRead::GetBufferSize
virtual Int_t GetBufferSize() const
Definition: TFileCacheRead.h:88
TFile::ECacheAction
ECacheAction
TTreeCache flushing semantics.
Definition: TFile.h:71
tree
Definition: tree.py:1
ClassImp
#define ClassImp(name)
Definition: Rtypes.h:364
TObjArray::Remove
virtual TObject * Remove(TObject *obj)
Remove object from array.
Definition: TObjArray.cxx:719
TTreeCache::fNReadMiss
Int_t fNReadMiss
Number of blocks read and not found in the cache.
Definition: TTreeCache.h:47
TObjString.h
TTreeCache::TTreeCache
TTreeCache()
Default Constructor.
Definition: TTreeCache.cxx:311
r
ROOT::R::TRInterface & r
Definition: Object.C:4
TBranch.h
TObject::Info
virtual void Info(const char *method, const char *msgfmt,...) const
Issue info message.
Definition: TObject.cxx:867
Long64_t
long long Long64_t
Definition: RtypesCore.h:80
TObject::Error
virtual void Error(const char *method, const char *msgfmt,...) const
Issue error message.
Definition: TObject.cxx:893
TTreeCache::fOptimizeMisses
Bool_t fOptimizeMisses
! true if we should optimize cache misses.
Definition: TTreeCache.h:72
TTreeCache::MissCache
Definition: TTreeCache.h:85
TTree
A TTree represents a columnar dataset.
Definition: TTree.h:79
TTreeCache::fNMissReadOk
Int_t fNMissReadOk
Number of blocks read, not found in the primary cache, and found in the secondary cache.
Definition: TTreeCache.h:46
TTreeCache::fReverseRead
Bool_t fReverseRead
! reading in reverse mode
Definition: TTreeCache.h:58
TTreeCache::fTree
TTree * fTree
! pointer to the current Tree
Definition: TTreeCache.h:53
TGeant4Unit::s
static constexpr double s
Definition: TGeant4SystemOfUnits.h:162
Int_t
int Int_t
Definition: RtypesCore.h:45
TVirtualPerfStats
Provides the interface for the PROOF internal performance measurement and event tracing.
Definition: TVirtualPerfStats.h:32
TTreeCache::CalculateMissEntries
TBranch * CalculateMissEntries(Long64_t, int, bool)
Given an file read, try to determine the corresponding branch.
Definition: TTreeCache.cxx:788
TTreeCache::fFirstTime
Bool_t fFirstTime
! save the fact that we processes the first entry
Definition: TTreeCache.h:60
TFileCacheRead::SecondPrefetch
virtual void SecondPrefetch(Long64_t, Int_t)
Definition: TFileCacheRead.cxx:258
TString::Contains
Bool_t Contains(const char *pat, ECaseCompare cmp=kExact) const
Definition: TString.h:624
TString::Length
Ssiz_t Length() const
Definition: TString.h:410
TList.h
TEnv::GetValue
virtual Int_t GetValue(const char *name, Int_t dflt) const
Returns the integer value for a resource.
Definition: TEnv.cxx:491
TTreeCache::fIsManual
Bool_t fIsManual
! true if cache is StopLearningPhase was used
Definition: TTreeCache.h:55
TFriendElement
A TFriendElement TF describes a TTree object TF in a file.
Definition: TFriendElement.h:33
TTreeCache::UpdateBranches
virtual void UpdateBranches(TTree *tree)
Update pointer to current Tree and recompute pointers to the branches in the cache.
Definition: TTreeCache.cxx:2197
TTreeCache::kNoPrefill
@ kNoPrefill
Definition: TTreeCache.h:35
TTreeCache::fFillTimes
Int_t fFillTimes
! how many times we can fill the current buffer
Definition: TTreeCache.h:59
TTreeCache::fFirstEntry
Long64_t fFirstEntry
! save the value of the first entry
Definition: TTreeCache.h:61
TObjArray::UncheckedAt
TObject * UncheckedAt(Int_t i) const
Definition: TObjArray.h:90
TFriendElement.h
TEnv.h
TString
Basic string class.
Definition: TString.h:136
Bool_t
bool Bool_t
Definition: RtypesCore.h:63
v
@ v
Definition: rootcling_impl.cxx:3664
b
#define b(i)
Definition: RSha256.hxx:100
TTreeCache::Print
virtual void Print(Option_t *option="") const
Print cache statistics.
Definition: TTreeCache.cxx:1886
TFile.h
TObjString::GetName
const char * GetName() const
Returns name of object.
Definition: TObjString.h:38
TBranchCacheInfo.h
bool
TString::ReplaceAll
TString & ReplaceAll(const TString &s1, const TString &s2)
Definition: TString.h:692
TTreeCache::fIsLearning
Bool_t fIsLearning
! true if cache is in learning mode
Definition: TTreeCache.h:54
TTreeCache::GetTree
TTree * GetTree() const
Definition: TTreeCache.h:149
TTreeCache::fBrNames
TList * fBrNames
! list of branch names in the cache
Definition: TTreeCache.h:52
TTreeCache::fPrefillType
EPrefillType fPrefillType
Whether a pre-filling is enabled (and if applicable which type)
Definition: TTreeCache.h:64
TTreeCache::fFirstBuffer
Bool_t fFirstBuffer
! true if first buffer is used for prefetching
Definition: TTreeCache.h:56
TObjString
Collectable string class.
Definition: TObjString.h:28
TTreeCache::SetOptimizeMisses
void SetOptimizeMisses(Bool_t opt)
Start of methods for the miss cache.
Definition: TTreeCache.cxx:681
TTree::GetTree
virtual TTree * GetTree() const
Definition: TTree.h:512
TBranch
A TTree is a list of TBranches.
Definition: TBranch.h:89
TTree::GetBranch
virtual TBranch * GetBranch(const char *name)
Return pointer to the branch with the given name in this tree or its friends.
Definition: TTree.cxx:5241
TRegexp.h
TString::Form
void Form(const char *fmt,...)
Formats a string using a printf style format descriptor.
Definition: TString.cxx:2311
TTreeCache::fReadDirectionSet
Bool_t fReadDirectionSet
! read direction established
Definition: TTreeCache.h:62
TTreeCache::fMissCache
std::unique_ptr< MissCache > fMissCache
! Cache contents for misses
Definition: TTreeCache.h:105
TTreeCache::fNextClusterStart
Long64_t fNextClusterStart
! End+1 of the cluster(s) where the current content was picked out
Definition: TTreeCache.h:43
TTreeCache.h
TChain.h
TBranchElement.h
TLeaf.h
TTreeCache::fFirstMiss
Long64_t fFirstMiss
! set to the event # of the first miss.
Definition: TTreeCache.h:73
TTreeCache::StopLearningPhase
virtual void StopLearningPhase()
This is the counterpart of StartLearningPhase() and can be used to stop the learning phase.
Definition: TTreeCache.cxx:2173
TSystem.h
TTreeCache::fNReadOk
Int_t fNReadOk
Number of blocks read and found in the cache.
Definition: TTreeCache.h:45
TFileCacheRead::ReadBuffer
virtual Int_t ReadBuffer(char *buf, Long64_t pos, Int_t len)
Read buffer at position pos.
Definition: TFileCacheRead.cxx:363
TLeaf
A TLeaf describes individual elements of a TBranch See TBranch structure in TTree.
Definition: TLeaf.h:57
TTree::GetReadEntry
virtual Long64_t GetReadEntry() const
Definition: TTree.h:504
TTree::TClusterIterator::GetNextEntry
Long64_t GetNextEntry()
Definition: TTree.h:303
TLeaf::GetLeafCount
virtual TLeaf * GetLeafCount() const
If this leaf stores a variable-sized array or a multi-dimensional array whose last dimension has vari...
Definition: TLeaf.h:120
TFriendElement::GetTree
virtual TTree * GetTree()
Return pointer to friend TTree.
Definition: TFriendElement.cxx:209
TObjArray::GetEntriesFast
Int_t GetEntriesFast() const
Definition: TObjArray.h:64
TTree::GetPerfStats
virtual TVirtualPerfStats * GetPerfStats() const
Definition: TTree.h:501
TTreeCache::AddBranch
virtual Int_t AddBranch(TBranch *b, Bool_t subgbranches=kFALSE)
Add a branch to the list of branches to be stored in the cache this function is called by the user vi...
Definition: TTreeCache.cxx:374
TTreeCache::MissCache::Entry
Definition: TTreeCache.h:86
TTreeCache::GetCachedBranches
const TObjArray * GetCachedBranches() const
Definition: TTreeCache.h:139
kFALSE
const Bool_t kFALSE
Definition: RtypesCore.h:101
TString::Append
TString & Append(const char *cs)
Definition: TString.h:564
TTreeCache::fEntryMax
Long64_t fEntryMax
! last entry in the cache
Definition: TTreeCache.h:39
TTreeCache::fNMissReadMiss
Int_t fNMissReadMiss
Number of blocks read and not found in either cache.
Definition: TTreeCache.h:48
TChain::GetTreeNumber
virtual Int_t GetTreeNumber() const
Definition: TChain.h:116
TFileCacheRead::SetBufferSize
virtual Int_t SetBufferSize(Int_t buffersize)
Sets the buffer size.
Definition: TFileCacheRead.cxx:709
R__unlikely
#define R__unlikely(expr)
Definition: RConfig.hxx:597
TObjArray::AddAtAndExpand
virtual void AddAtAndExpand(TObject *obj, Int_t idx)
Add object at position idx.
Definition: TObjArray.cxx:235
gDebug
Int_t gDebug
Definition: TROOT.cxx:590
TTreeCache::GetMissEfficiencyRel
Double_t GetMissEfficiencyRel() const
Relative efficiency of the 'miss cache' - ratio of the reads found in cache to the number of reads so...
Definition: TTreeCache.cxx:1848
UInt_t
unsigned int UInt_t
Definition: RtypesCore.h:46
TTreeCache::GetMissEfficiency
Double_t GetMissEfficiency() const
The total efficiency of the 'miss cache' - defined as the ratio of blocks found in the cache versus t...
Definition: TTreeCache.cxx:1824
TTreeCache::fLearnPrefilling
Bool_t fLearnPrefilling
! true if we are in the process of executing LearnPrefill
Definition: TTreeCache.h:68
TTreeCache::fNbranches
Int_t fNbranches
! Number of branches in the cache
Definition: TTreeCache.h:44
TTree::TClusterIterator::Next
Long64_t Next()
Move on to the next cluster and return the starting entry of this next cluster.
Definition: TTree.cxx:648
TObject::Warning
virtual void Warning(const char *method, const char *msgfmt,...) const
Issue warning message.
Definition: TObject.cxx:879
TTreeCache::LearnBranch
virtual Int_t LearnBranch(TBranch *b, Bool_t subgbranches=kFALSE)
Add a branch discovered by actual usage to the list of branches to be stored in the cache this functi...
Definition: TTreeCache.cxx:348
TTreeCache::fEntryMin
Long64_t fEntryMin
! first entry in the cache
Definition: TTreeCache.h:38
TFile
A ROOT file is a suite of consecutive data records (TKey instances) with a well defined format.
Definition: TFile.h:54
TTreeCache::fEntryNext
Long64_t fEntryNext
! next entry number where cache must be filled
Definition: TTreeCache.h:41
TTreeCache::ProcessMiss
Bool_t ProcessMiss(Long64_t pos, int len)
! Given a file read not in the miss cache, handle (possibly) loading the data.
Definition: TTreeCache.cxx:861
TTreeCache::fLastMiss
Long64_t fLastMiss
! set to the event # of the last miss.
Definition: TTreeCache.h:74
TFileCacheRead::fBNtot
Int_t fBNtot
Definition: TFileCacheRead.h:59
TTreeCache::ResetCache
virtual void ResetCache()
This will simply clear the cache.
Definition: TTreeCache.cxx:2027
unsigned int
TFileCacheRead
A cache when reading files over the network.
Definition: TFileCacheRead.h:22
TTreeCache::SetEntryRange
virtual void SetEntryRange(Long64_t emin, Long64_t emax)
Set the minimum and maximum entry number to be processed this information helps to optimize the numbe...
Definition: TTreeCache.cxx:2092
TRegexp
Regular expression class.
Definition: TRegexp.h:31
TMath::BinarySearch
Long64_t BinarySearch(Long64_t n, const T *array, T value)
Definition: TMathBase.h:278
TFileCacheRead::Print
virtual void Print(Option_t *option="") const
Print cache statistics.
Definition: TFileCacheRead.cxx:325
TFileCacheRead::fNtot
Int_t fNtot
Total size of prefetched blocks.
Definition: TFileCacheRead.h:40
Printf
void Printf(const char *fmt,...)
gSystem
R__EXTERN TSystem * gSystem
Definition: TSystem.h:559
TFile::ReadBuffers
virtual Bool_t ReadBuffers(char *buf, Long64_t *pos, Int_t *len, Int_t nbuf)
Read the nbuf blocks described in arrays pos and len.
Definition: TFile.cxx:1738
TTree::TClusterIterator
Helper class to iterate over cluster of baskets.
Definition: TTree.h:266
TTreeCache::GetEfficiency
Double_t GetEfficiency() const
Give the total efficiency of the primary cache...
Definition: TTreeCache.cxx:1812
TObjArray::AddAt
virtual void AddAt(TObject *obj, Int_t idx)
Add object at position ids.
Definition: TObjArray.cxx:254
TSystem::Getenv
virtual const char * Getenv(const char *env)
Get environment variable.
Definition: TSystem.cxx:1662
TTreeCache::IOPos
Definition: TTreeCache.h:78
TTreeCache::fgLearnEntries
static Int_t fgLearnEntries
number of entries used for learning mode
Definition: TTreeCache.h:65
ULong64_t
unsigned long long ULong64_t
Definition: RtypesCore.h:81
TLeaf::GetBranch
TBranch * GetBranch() const
Definition: TLeaf.h:115
TTreeCache
A cache to speed-up the reading of ROOT datasets.
Definition: TTreeCache.h:32
TVirtualPerfStats::UpdateBranchIndices
virtual void UpdateBranchIndices(TObjArray *branches)=0
Double_t
double Double_t
Definition: RtypesCore.h:59
TObjArray.h
TTreeCache::fOneTime
Bool_t fOneTime
! used in the learning phase
Definition: TTreeCache.h:57
TList::Remove
virtual TObject * Remove(TObject *obj)
Remove object from the list.
Definition: TList.cxx:822
TTreeCache::fBranches
TObjArray * fBranches
! List of branches to be stored in the cache
Definition: TTreeCache.h:51
file
Definition: file.py:1
TList::Add
virtual void Add(TObject *obj)
Definition: TList.h:87
TTreeCache::SetLearnEntries
static void SetLearnEntries(Int_t n=10)
Static function to set the number of entries to be used in learning mode The default value for n is 1...
Definition: TTreeCache.cxx:2131
TTreeCache::GetLearnEntries
static Int_t GetLearnEntries()
Static function returning the number of entries used to train the cache see SetLearnEntries.
Definition: TTreeCache.cxx:1861
name
char name[80]
Definition: TGX11.cxx:110
TFileCacheRead::fEnablePrefetching
Bool_t fEnablePrefetching
reading by prefetching asynchronously
Definition: TFileCacheRead.h:37
TIter
Definition: TCollection.h:233
TTree::GetEventList
TEventList * GetEventList() const
Definition: TTree.h:468
TTreeCache::fEntryCurrent
Long64_t fEntryCurrent
! current lowest entry number in the cache
Definition: TTreeCache.h:40
TFileCacheRead::SetFile
virtual void SetFile(TFile *file, TFile::ECacheAction action=TFile::kDisconnect)
Set the file using this cache and reset the current blocks (if any).
Definition: TFileCacheRead.cxx:544
TTreeCache::CheckMissCache
Bool_t CheckMissCache(char *buf, Long64_t pos, int len)
Check the miss cache for a particular buffer, fetching if deemed necessary.
Definition: TTreeCache.cxx:919
TTreeCache::ReadBufferNormal
virtual Int_t ReadBufferNormal(char *buf, Long64_t pos, Int_t len)
Old method ReadBuffer before the addition of the prefetch mechanism.
Definition: TTreeCache.cxx:1914
TFileCacheRead::fNseek
Int_t fNseek
Number of blocks to be prefetched.
Definition: TFileCacheRead.h:39
TChain
A chain is a collection of files containing TTree objects.
Definition: TChain.h:33
TTreeCache::SetLearnPrefill
virtual void SetLearnPrefill(EPrefillType type=kNoPrefill)
Set whether the learning period is started with a prefilling of the cache and which type of prefillin...
Definition: TTreeCache.cxx:2146
TTreeCache::DropBranch
virtual Int_t DropBranch(TBranch *b, Bool_t subbranches=kFALSE)
Remove a branch to the list of branches to be stored in the cache this function is called by TBranch:...
Definition: TTreeCache.cxx:539
TNamed::GetName
virtual const char * GetName() const
Returns name of object.
Definition: TNamed.h:47
type
int type
Definition: TGX11.cxx:121
TBranch::kDoNotProcess
@ kDoNotProcess
Active bit for branches.
Definition: TBranch.h:101
TTreeCache::fCurrentClusterStart
Long64_t fCurrentClusterStart
! Start of the cluster(s) where the current content was picked out
Definition: TTreeCache.h:42
TFileCacheRead::fFile
TFile * fFile
Pointer to file.
Definition: TFileCacheRead.h:51
Class
void Class()
Definition: Class.C:29
TString::ToLower
void ToLower()
Change string to lower-case.
Definition: TString.cxx:1147
TTreeCache::ResetMissCache
void ResetMissCache()
Reset all the miss cache training.
Definition: TTreeCache.cxx:695
TTreeCache::GetEfficiencyRel
Double_t GetEfficiencyRel() const
This will indicate a sort of relative efficiency...
Definition: TTreeCache.cxx:1836
TTreeCache::~TTreeCache
virtual ~TTreeCache()
Destructor. (in general called by the TFile destructor)
Definition: TTreeCache.cxx:330
TTreeCache::ReadBufferPrefetch
virtual Int_t ReadBufferPrefetch(char *buf, Long64_t pos, Int_t len)
Used to read a chunk from a block previously fetched.
Definition: TTreeCache.cxx:1973
TVirtualPerfStats::SetMissed
virtual void SetMissed(TBranch *b, size_t basketNumber)=0
TTree::GetListOfLeaves
virtual TObjArray * GetListOfLeaves()
Definition: TTree.h:484
TTree::GetEntries
virtual Long64_t GetEntries() const
Definition: TTree.h:458
TFileCacheRead::fBufferSizeMin
Int_t fBufferSizeMin
Original size of fBuffer.
Definition: TFileCacheRead.h:26
TList
A doubly linked list.
Definition: TList.h:44
TTreeCache::FindBranchBasketPos
IOPos FindBranchBasketPos(TBranch &, Long64_t entry)
Given a branch and an entry, determine the file location (offset / size) of the corresponding basket.
Definition: TTreeCache.cxx:714
TMath.h
TFile::SetCacheRead
virtual void SetCacheRead(TFileCacheRead *cache, TObject *tree=0, ECacheAction action=kDisconnect)
Set a pointer to the read cache.
Definition: TFile.cxx:2282
TFileCacheRead::fIsTransferred
Bool_t fIsTransferred
True when fBuffer contains something valid.
Definition: TFileCacheRead.h:54
int
TEventList.h