Hall-D Software  alpha
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
DTreeInterface.cc
Go to the documentation of this file.
1 #include "DTreeInterface.h"
2 
3 /************************************************* STATIC-VARIABLE-ACCESSING PRIVATE MEMBER FUNCTIONS *************************************************/
4 
5 map<string, int>& DTreeInterface::Get_NumWritersByFileMap(void) const
6 {
7  // string is file name
8  // must be read/used entirely in global root lock: when 0, close the file (modifying files needs global root lock)
9  static map<string, int> locNumWritersByFileMap;
10  return locNumWritersByFileMap;
11 }
12 
13 map<string, size_t>& DTreeInterface::Get_FundamentalArraySizeMap(TTree* locTree) const
14 {
15  //A "global" lock (shared across all threads & DTreeInterface objects) is necessary to read the outer map
16  //This lock is acquired locally: DO NOT CALL THIS FUNCTION WHILE WITHIN A LOCK!!
17  //However, only a file-lock is necessary to read the inner map once it is acquired:
18  //Only the given tree is viewing/modifying it
19  japp->WriteLock("DTreeInterface_SizeMap"); //LOCK MAP
20 
21  static map<TTree*, map<string, size_t> > locFundamentalArraySizeMap;
22  map<string, size_t>& locTreeSpecificMap = locFundamentalArraySizeMap[locTree];
23 
24  japp->Unlock("DTreeInterface_SizeMap"); //UNLOCK MAP
25  return locTreeSpecificMap;
26 }
27 
28 /********************************************************************* INITIALIZE *********************************************************************/
29 
30 DTreeInterface* DTreeInterface::Create_DTreeInterface(string locTreeName, string locFileName)
31 {
32  if(locFileName == "hd_root.root")
33  {
34  cout << "WARNING: SAVING TREES TO hd_root.root IS NOT SUPPORTED (or wise). ";
35  cout << "RETURNING NULL FROM DTreeInterface::Create_DTreeInterface()" << endl;
36  return NULL;
37  }
38  return new DTreeInterface(locTreeName, locFileName);
39 }
40 
41 //Constructor
42 DTreeInterface::DTreeInterface(string locTreeName, string locFileName) : dFileName(locFileName),
43 dTreeIndex_MajorBranchName("0"), dTreeIndex_MinorBranchName("0")
44 {
45  japp->RootWriteLock();
46  {
47  map<string, int>& locNumWritersByFileMap = Get_NumWritersByFileMap();
48  if(locNumWritersByFileMap.find(dFileName) == locNumWritersByFileMap.end())
49  locNumWritersByFileMap[dFileName] = 1;
50  else
51  ++locNumWritersByFileMap[dFileName];
52  }
53  japp->RootUnLock();
54 
55  GetOrCreate_FileAndTree(locTreeName);
56 }
57 
58 //Destructor
60 {
61  japp->RootWriteLock();
62  {
63  map<string, int>& locNumWritersByFileMap = Get_NumWritersByFileMap();
64 
65  --locNumWritersByFileMap[dFileName];
66  if(locNumWritersByFileMap[dFileName] != 0)
67  {
68  japp->RootUnLock();
69  return;
70  }
71 
72  //Build index if requested
74  dTree->BuildIndex(dTreeIndex_MajorBranchName.c_str(), dTreeIndex_MinorBranchName.c_str());
75 
76  //Save to output
77  TFile* locOutputFile = (TFile*)gROOT->GetListOfFiles()->FindObject(dFileName.c_str());
78  locOutputFile->Write(0, TObject::kOverwrite);
79  locOutputFile->Close();
80  delete locOutputFile;
81  }
82  japp->RootUnLock();
83 }
84 
86 {
87  japp->RootWriteLock();
88  {
89  TDirectory* locCurrentDir = gDirectory;
90 
91  //see if root file exists already
92  TFile* locOutputFile = (TFile*)gROOT->GetListOfFiles()->FindObject(dFileName.c_str());
93  if(locOutputFile == nullptr)
94  locOutputFile = new TFile(dFileName.c_str(), "RECREATE");
95  locOutputFile->cd();
96 
97  //see if ttree exists already. if not, create it
98  dTree = (TTree*)gDirectory->Get(locTreeName.c_str());
99  if(dTree == nullptr)
100  {
101  dTree = new TTree(locTreeName.c_str(), locTreeName.c_str());
102  dTree->SetAutoFlush(dAutoFlush);
103  }
104 
105  locCurrentDir->cd();
106  }
107  japp->RootUnLock();
108 }
109 
110 /****************************************************************** CREATE BRANCHES *******************************************************************/
111 
113 {
114  japp->WriteLock(dFileName); //LOCK FILE
115  bool locBranchesCreatedFlag = (dTree->GetNbranches() > 0);
116  japp->Unlock(dFileName); //UNLOCK FILE
117  return locBranchesCreatedFlag;
118 }
119 
120 const TList* DTreeInterface::Get_UserInfo(void) const
121 {
122  return dTree->GetUserInfo();
123 }
124 
125 bool DTreeInterface::Create_Branches(const DTreeBranchRegister& locTreeBranchRegister)
126 {
127  //return value: true if created, false if previously created (nothing done here)
128  //if false, user can safely delete memory in the branch register UserInfo TList
129 
130  //MUST CARRY AROUND A REFERENCE TO THIS. ONLY READ/MODIFY THE MAP WITHIN A FILE LOCK.
131  map<string, size_t>& locFundamentalArraySizeMap = Get_FundamentalArraySizeMap(dTree);
132 
133  japp->WriteLock(dFileName); //LOCK FILE
134  {
135  //if there are branches already, don't do anything
136  if(dTree->GetNbranches() > 0)
137  {
138  japp->Unlock(dFileName); //UNLOCK FILE
139  return false;
140  }
141 
142  //set user info
143  TList* locInputUserInfo = locTreeBranchRegister.Get_UserInfo();
144  TList* locTreeUserInfo = dTree->GetUserInfo();
145  for(Int_t loc_i = 0; loc_i < locInputUserInfo->GetSize(); ++loc_i)
146  locTreeUserInfo->Add(locInputUserInfo->At(loc_i));
147 
148  //loop over branches
149  //branches that are arrays, for which the branch with the size was not created it yet
150  //save them, and create them at the end
151  vector<string> locPostponedBranchNames;
152  for(const auto& locBranchName : locTreeBranchRegister.dBranchNames)
153  {
154  //check if need to postpone this branch
155  auto locSizeNameIterator = locTreeBranchRegister.dArraySizeNameMap.find(locBranchName);
156  if(locSizeNameIterator != locTreeBranchRegister.dArraySizeNameMap.end())
157  {
158  //fundamental array: make sure array size branch has been created first. if not, postpone it
159  if(dTree->GetBranch(locSizeNameIterator->second.c_str()) == nullptr)
160  {
161  locPostponedBranchNames.push_back(locBranchName);
162  continue;
163  }
164  }
165 
166  //it's ok: create it
167  Create_Branch(locTreeBranchRegister, locBranchName, locFundamentalArraySizeMap);
168  }
169 
170  //create postponed branches
171  for(string& locBranchName : locPostponedBranchNames)
172  Create_Branch(locTreeBranchRegister, locBranchName, locFundamentalArraySizeMap);
173  }
174  japp->Unlock(dFileName); //UNLOCK FILE
175 
176  return true;
177 }
178 
179 void DTreeInterface::Create_Branch(const DTreeBranchRegister& locTreeBranchRegister, string locBranchName, map<string, size_t>& locFundamentalArraySizeMap)
180 {
181  const type_index& locTypeIndex = locTreeBranchRegister.dBranchTypeMap.find(locBranchName)->second;
182 
183  //array size
184  auto locSizeIterator = locTreeBranchRegister.dInitialArraySizeMap.find(locBranchName);
185  if(locSizeIterator == locTreeBranchRegister.dInitialArraySizeMap.end())
186  {
187  //not an array
188  Create_Branch(locBranchName, locTypeIndex, 0, "");
189  return;
190  }
191  size_t locArraySize = locSizeIterator->second;
192 
193  //array size name
194  auto locSizeNameIterator = locTreeBranchRegister.dArraySizeNameMap.find(locBranchName);
195  if(locSizeNameIterator == locTreeBranchRegister.dArraySizeNameMap.end())
196  Create_Branch(locBranchName, locTypeIndex, locArraySize, ""); //clones array
197  else //fundamental array
198  {
199  Create_Branch(locBranchName, locTypeIndex, locArraySize, locSizeNameIterator->second);
200  locFundamentalArraySizeMap[locBranchName] = locArraySize;
201  }
202 }
203 
204 void DTreeInterface::Create_Branch(string locBranchName, type_index locTypeIndex, size_t locArraySize, string locArraySizeName)
205 {
206  //Fundamental types
207  if(locTypeIndex == type_index(typeid(Char_t)))
208  Create_Branch<Char_t>(locBranchName, locArraySize, locArraySizeName);
209  else if(locTypeIndex == type_index(typeid(UChar_t)))
210  Create_Branch<UChar_t>(locBranchName, locArraySize, locArraySizeName);
211  else if(locTypeIndex == type_index(typeid(Short_t)))
212  Create_Branch<Short_t>(locBranchName, locArraySize, locArraySizeName);
213  else if(locTypeIndex == type_index(typeid(UShort_t)))
214  Create_Branch<UShort_t>(locBranchName, locArraySize, locArraySizeName);
215  else if(locTypeIndex == type_index(typeid(Int_t)))
216  Create_Branch<Int_t>(locBranchName, locArraySize, locArraySizeName);
217  else if(locTypeIndex == type_index(typeid(UInt_t)))
218  Create_Branch<UInt_t>(locBranchName, locArraySize, locArraySizeName);
219  else if(locTypeIndex == type_index(typeid(Float_t)))
220  Create_Branch<Float_t>(locBranchName, locArraySize, locArraySizeName);
221  else if(locTypeIndex == type_index(typeid(Double_t)))
222  Create_Branch<Double_t>(locBranchName, locArraySize, locArraySizeName);
223  else if(locTypeIndex == type_index(typeid(Long64_t)))
224  Create_Branch<Long64_t>(locBranchName, locArraySize, locArraySizeName);
225  else if(locTypeIndex == type_index(typeid(ULong64_t)))
226  Create_Branch<ULong64_t>(locBranchName, locArraySize, locArraySizeName);
227  else if(locTypeIndex == type_index(typeid(Bool_t)))
228  Create_Branch<Bool_t>(locBranchName, locArraySize, locArraySizeName);
229 
230  //TObject
231  else if(locTypeIndex == type_index(typeid(TVector2)))
232  Create_Branch<TVector2>(locBranchName, locArraySize, locArraySizeName);
233  else if(locTypeIndex == type_index(typeid(TVector3)))
234  Create_Branch<TVector3>(locBranchName, locArraySize, locArraySizeName);
235  else if(locTypeIndex == type_index(typeid(TLorentzVector)))
236  Create_Branch<TLorentzVector>(locBranchName, locArraySize, locArraySizeName);
237 }
238 
239 /**************************************************************** FILL BRANCHES & TREE ****************************************************************/
240 
241 void DTreeInterface::Fill(DTreeFillData& locTreeFillData)
242 {
243  //MUST CARRY AROUND A REFERENCE TO THIS. ONLY READ/MODIFY THE MAP WITHIN A FILE LOCK.
244  map<string, size_t>& locFundamentalArraySizeMap = Get_FundamentalArraySizeMap(dTree);
245  japp->WriteLock(dFileName); //LOCK FILE
246  {
247  //loop over branches
248  for(auto locBranchPair : locTreeFillData.dFillData)
249  {
250  //type
251  auto locBranchName = locBranchPair.first;
252  auto locTypeIndex = locBranchPair.second.first;
253  auto locFillBaseClass = locBranchPair.second.second;
254 
255  if(dTree->GetBranch(locBranchName.c_str()) == NULL)
256  {
257  cout << "WARNING, CANNOT FILL DATA, BRANCH " << locBranchName << " DOES NOT EXIST." << endl;
258  continue;
259  }
260 
261  //check if is array. if not, fill
262  auto locLargestIndexFilledIterator = locTreeFillData.dArrayLargestIndexFilledMap.find(locBranchName);
263  bool locIsArrayFlag = (locLargestIndexFilledIterator != locTreeFillData.dArrayLargestIndexFilledMap.end());
264  if(!locIsArrayFlag)
265  {
266  Fill(locBranchName, locTypeIndex, locFillBaseClass->Get(0), false);
267  continue;
268  }
269 
270  //is array, get how many to fill
271  auto& locLargestIndexFilled = locLargestIndexFilledIterator->second;
272 
273  //increase array size if necessary, or decrease if too large from last event
274  auto locFundamentalArraySizeIterator = locFundamentalArraySizeMap.find(locBranchName);
275  if(locFundamentalArraySizeIterator != locFundamentalArraySizeMap.end())
276  {
277  size_t locCurrentArraySize = locFundamentalArraySizeMap[locBranchName]; //may not be in map! (tobj)
278  if((locLargestIndexFilled + 1) > int(locCurrentArraySize))
279  {
280  Change_ArraySize(locBranchName, locTypeIndex, locLargestIndexFilled + 1);
281  locFundamentalArraySizeMap[locBranchName] = locLargestIndexFilled + 1;
282  }
283 /* else if(locCurrentArraySize > dMaxArraySize)
284  {
285  Change_ArraySize(locBranchName, locTypeIndex, dMaxArraySize);
286  locFundamentalArraySizeMap[locBranchName] = dMaxArraySize;
287  } */
288  }
289  else //is clones array: clear it
290  {
291  auto locClonesArray = Get_Pointer_TClonesArray(locBranchName);
292  locClonesArray->Clear(); //empties array
293 // if(size_t(locClonesArray->GetEntriesFast()) > dMaxArraySize)
294 // locClonesArray->Expand(dMaxArraySize);
295  }
296 
297  //fill array
298  for(int locArrayIndex = 0; locArrayIndex <= locLargestIndexFilled; ++locArrayIndex)
299  Fill(locBranchName, locTypeIndex, locFillBaseClass->Get(locArrayIndex), true, locArrayIndex);
300 
301  //reset DTreeFillData for next event!
302  locLargestIndexFilled = -1;
303  }
304 
305  //fill tree
306  dTree->Fill();
307  }
308  japp->Unlock(dFileName); //UNLOCK FILE
309 
310  //Reset fill vectors if too large!!
311  for(auto locBranchPair : locTreeFillData.dFillData)
312  locBranchPair.second.second->Check_Capacity();
313 }
314 
315 void DTreeInterface::Change_ArraySize(string locBranchName, type_index locTypeIndex, size_t locNewArraySize)
316 {
317  //Fundamental types
318  if(locTypeIndex == type_index(typeid(Char_t)))
319  Change_ArraySize<Char_t>(locBranchName, locNewArraySize);
320  else if(locTypeIndex == type_index(typeid(UChar_t)))
321  Change_ArraySize<UChar_t>(locBranchName, locNewArraySize);
322  else if(locTypeIndex == type_index(typeid(Short_t)))
323  Change_ArraySize<Short_t>(locBranchName, locNewArraySize);
324  else if(locTypeIndex == type_index(typeid(UShort_t)))
325  Change_ArraySize<UShort_t>(locBranchName, locNewArraySize);
326  else if(locTypeIndex == type_index(typeid(Int_t)))
327  Change_ArraySize<Int_t>(locBranchName, locNewArraySize);
328  else if(locTypeIndex == type_index(typeid(UInt_t)))
329  Change_ArraySize<UInt_t>(locBranchName, locNewArraySize);
330  else if(locTypeIndex == type_index(typeid(Float_t)))
331  Change_ArraySize<Float_t>(locBranchName, locNewArraySize);
332  else if(locTypeIndex == type_index(typeid(Double_t)))
333  Change_ArraySize<Double_t>(locBranchName, locNewArraySize);
334  else if(locTypeIndex == type_index(typeid(Long64_t)))
335  Change_ArraySize<Long64_t>(locBranchName, locNewArraySize);
336  else if(locTypeIndex == type_index(typeid(ULong64_t)))
337  Change_ArraySize<ULong64_t>(locBranchName, locNewArraySize);
338  else if(locTypeIndex == type_index(typeid(Bool_t)))
339  Change_ArraySize<Bool_t>(locBranchName, locNewArraySize);
340 }
341 
342 void DTreeInterface::Fill(string locBranchName, type_index locTypeIndex, void* locVoidPointer, bool locIsArrayFlag, size_t locArrayIndex)
343 {
344  //Fundamental types
345  if(locTypeIndex == type_index(typeid(Char_t)))
346  Get_Pointer_Fundamental<Char_t>(locBranchName)[locArrayIndex] = *(static_cast<Char_t*>(locVoidPointer));
347  else if(locTypeIndex == type_index(typeid(UChar_t)))
348  Get_Pointer_Fundamental<UChar_t>(locBranchName)[locArrayIndex] = *(static_cast<UChar_t*>(locVoidPointer));
349  else if(locTypeIndex == type_index(typeid(Short_t)))
350  Get_Pointer_Fundamental<Short_t>(locBranchName)[locArrayIndex] = *(static_cast<Short_t*>(locVoidPointer));
351  else if(locTypeIndex == type_index(typeid(UShort_t)))
352  Get_Pointer_Fundamental<UShort_t>(locBranchName)[locArrayIndex] = *(static_cast<UShort_t*>(locVoidPointer));
353  else if(locTypeIndex == type_index(typeid(Int_t)))
354  Get_Pointer_Fundamental<Int_t>(locBranchName)[locArrayIndex] = *(static_cast<Int_t*>(locVoidPointer));
355  else if(locTypeIndex == type_index(typeid(UInt_t)))
356  Get_Pointer_Fundamental<UInt_t>(locBranchName)[locArrayIndex] = *(static_cast<UInt_t*>(locVoidPointer));
357  else if(locTypeIndex == type_index(typeid(Float_t)))
358  Get_Pointer_Fundamental<Float_t>(locBranchName)[locArrayIndex] = *(static_cast<Float_t*>(locVoidPointer));
359  else if(locTypeIndex == type_index(typeid(Double_t)))
360  Get_Pointer_Fundamental<Double_t>(locBranchName)[locArrayIndex] = *(static_cast<Double_t*>(locVoidPointer));
361  else if(locTypeIndex == type_index(typeid(Long64_t)))
362  Get_Pointer_Fundamental<Long64_t>(locBranchName)[locArrayIndex] = *(static_cast<Long64_t*>(locVoidPointer));
363  else if(locTypeIndex == type_index(typeid(ULong64_t)))
364  Get_Pointer_Fundamental<ULong64_t>(locBranchName)[locArrayIndex] = *(static_cast<ULong64_t*>(locVoidPointer));
365  else if(locTypeIndex == type_index(typeid(Bool_t)))
366  Get_Pointer_Fundamental<Bool_t>(locBranchName)[locArrayIndex] = *(static_cast<Bool_t*>(locVoidPointer));
367 
368  //TObject
369  else if(locTypeIndex == type_index(typeid(TVector3)))
370  {
371  TVector3& locObject = *(static_cast<TVector3*>(locVoidPointer));
372  Fill_TObject<TVector3>(locBranchName, locObject, locIsArrayFlag, locArrayIndex);
373  }
374  else if(locTypeIndex == type_index(typeid(TVector2)))
375  {
376  TVector2& locObject = *(static_cast<TVector2*>(locVoidPointer));
377  Fill_TObject<TVector2>(locBranchName, locObject, locIsArrayFlag, locArrayIndex);
378  }
379  else if(locTypeIndex == type_index(typeid(TLorentzVector)))
380  {
381  TLorentzVector& locObject = *(static_cast<TLorentzVector*>(locVoidPointer));
382  Fill_TObject<TLorentzVector>(locBranchName, locObject, locIsArrayFlag, locArrayIndex);
383  }
384 }
map< string, size_t > & Get_FundamentalArraySizeMap(TTree *locTree) const
map< string, size_t > dInitialArraySizeMap
void Change_ArraySize(string locBranchName, type_index locTypeIndex, size_t locNewArraySize)
void Create_Branch(const DTreeBranchRegister &locTreeBranchRegister, string locBranchName, map< string, size_t > &locFundamentalArraySizeMap)
const TList * Get_UserInfo(void) const
void GetOrCreate_FileAndTree(string locTreeName)
string dTreeIndex_MinorBranchName
bool Create_Branches(const DTreeBranchRegister &locTreeBranchRegister)
map< string, string > dArraySizeNameMap
map< string, pair< type_index, DFillBaseClass * > > dFillData
JApplication * japp
vector< string > dBranchNames
void Fill(DTreeFillData &locTreeFillData)
static DTreeInterface * Create_DTreeInterface(string locTreeName, string locFileName)
map< string, int > & Get_NumWritersByFileMap(void) const
DTreeInterface(void)
TList * Get_UserInfo(void) const
bool Get_BranchesCreatedFlag(void) const
map< string, int > dArrayLargestIndexFilledMap
TClonesArray * Get_Pointer_TClonesArray(string locBranchName)
Long64_t dAutoFlush
map< string, type_index > dBranchTypeMap
string dTreeIndex_MajorBranchName