8ce85f9834b9c47def0f8be33e4ad57c9cd66947
[u/mrichter/AliRoot.git] / HLT / rec / AliHLTEsdManager.cxx
1 // $Id$
2
3 //**************************************************************************
4 //* This file is property of and copyright by the ALICE HLT Project        * 
5 //* ALICE Experiment at CERN, All rights reserved.                         *
6 //*                                                                        *
7 //* Primary Authors: Matthias Richter <Matthias.Richter@ift.uib.no>        *
8 //*                  for The ALICE HLT Project.                            *
9 //*                                                                        *
10 //* Permission to use, copy, modify and distribute this software and its   *
11 //* documentation strictly for non-commercial purposes is hereby granted   *
12 //* without fee, provided that the above copyright notice appears in all   *
13 //* copies and that both the copyright notice and this permission notice   *
14 //* appear in the supporting documentation. The authors make no claims     *
15 //* about the suitability of this software for any purpose. It is          *
16 //* provided "as is" without express or implied warranty.                  *
17 //**************************************************************************
18
19 /** @file   AliHLTEsdManager.cxx
20     @author Matthias Richter
21     @date   
22     @brief  Manager for merging and writing of HLT ESDs
23 */
24
25 #include "AliHLTEsdManager.h"
26 #include "AliHLTComponent.h"
27 #include "AliESDEvent.h"
28 #include "AliHLTMessage.h"
29 #include "AliESDEvent.h"
30 #include "TFile.h"
31 #include "TTree.h"
32 #include "TClass.h"
33 #include "TObject.h"
34 #include "TObjectTable.h"
35 #include "TSystem.h"
36 #include "TChain.h"
37 #include "TList.h"
38
39 /** ROOT macro for the implementation of ROOT specific class methods */
40 ClassImp(AliHLTEsdManager)
41
42 AliHLTEsdManager::AliHLTEsdManager()
43   :
44   fESDs(),
45   fDirectory()
46 {
47   // see header file for class documentation
48   // or
49   // refer to README to build package
50   // or
51   // visit http://web.ift.uib.no/~kjeks/doc/alice-hlt
52 }
53
54 AliHLTEsdManager::~AliHLTEsdManager()
55 {
56   // see header file for class documentation
57   for (unsigned int i=0; i<fESDs.size(); i++) {
58     if (fESDs[i]) {
59       delete fESDs[i];
60     }
61     fESDs[i]=NULL;
62   }
63 }
64
65 AliHLTEsdManager::AliHLTEsdListEntry* AliHLTEsdManager::Find(AliHLTComponentDataType dt) const
66 {
67   // see header file for class documentation
68   AliHLTEsdListEntry* pEntry=NULL;
69   for (unsigned int i=0; i<fESDs.size(); i++) {
70     if (fESDs[i] && *(fESDs[i])==dt) {
71       pEntry=const_cast<AliHLTEsdListEntry*>(fESDs[i]);
72     }
73   }
74   return pEntry;
75 }
76
77 int AliHLTEsdManager::WriteESD(const AliHLTUInt8_t* pBuffer, AliHLTUInt32_t size,
78                                AliHLTComponentDataType dt, AliESDEvent* tgtesd, int eventno)
79 {
80   // see header file for class documentation
81   if (!pBuffer && size<=0) return -EINVAL;
82   int iResult=0;
83   AliHLTUInt32_t firstWord=*((AliHLTUInt32_t*)pBuffer);
84   if (firstWord==size-sizeof(AliHLTUInt32_t)) {
85     HLTDebug("create object from block %s size %d", AliHLTComponent::DataType2Text(dt).c_str(), size);
86     AliHLTMessage msg(const_cast<AliHLTUInt8_t*>(pBuffer), size);
87     TClass* objclass=msg.GetClass();
88     TObject* pObj=msg.ReadObject(objclass);
89     if (pObj && objclass) {
90       HLTDebug("object %p type %s created", pObj, objclass->GetName());
91       AliESDEvent* pESD=dynamic_cast<AliESDEvent*>(pObj);
92       TTree* pTree=NULL;
93       if (!pESD) {
94         pTree=dynamic_cast<TTree*>(pObj);
95         if (pTree) {
96           pESD=new AliESDEvent;
97           pESD->CreateStdContent();
98           if (pTree->GetEntries()>0) {
99             if (pTree->GetEntries()>1) {
100               HLTWarning("only one entry allowed for ESD embedded into tree, data block %s contains tree with %d entires, taking first entry",
101                          AliHLTComponent::DataType2Text(dt).c_str(), pTree->GetEntries());
102             }
103             pESD->ReadFromTree(pTree);
104             pTree->GetEvent(0);
105           }
106         } else {
107           HLTWarning("tree of data block %s has no events, skipping data block", AliHLTComponent::DataType2Text(dt).c_str());
108         }
109       }
110       if (pESD) {
111         AliHLTEsdListEntry* entry=Find(dt);
112         if (!entry) {
113           AliHLTEsdListEntry* newEntry=new AliHLTEsdListEntry(dt);
114           if (!fDirectory.IsNull()) {
115             newEntry->SetDirectory(fDirectory);
116           }
117           fESDs.push_back(newEntry);
118         }
119         if (tgtesd) {
120 #ifdef HAVE_ESD_COPY
121           *tgtesd=*pESD;
122 #else
123           TTree* pTmpTree=AliHLTEsdManager::EmbedIntoTree(pESD);
124           if (pTmpTree) {
125             tgtesd->ReadFromTree(pTmpTree);
126             pTmpTree->GetEvent(0);
127             pTmpTree->GetUserInfo()->Clear();
128             delete pTmpTree;
129             HLTDebug("data block %s written to target ESD", AliHLTComponent::DataType2Text(dt).c_str());
130           } else {
131             iResult=-ENOMEM;
132           }
133 #endif
134         } else {
135         entry=Find(dt);
136         if (entry) {
137           entry->WriteESD(pESD, eventno);
138         } else {
139           HLTError("internal mismatch, can not create list entry");
140           iResult=-ENOMEM;
141         }
142         }
143       } else {
144         HLTWarning("data block %s is not of class type AliESDEvent, ignoring ...", AliHLTComponent::DataType2Text(dt).c_str());
145       }
146       if (pTree) {
147         // ESD has been created and must be cleaned up
148         pESD->Reset();
149         delete pESD;
150         pESD=NULL;
151       }
152       delete pObj;
153       pObj=NULL;
154     } else {
155     }
156   }
157   return iResult;
158 }
159
160 int AliHLTEsdManager::PadESDs(int eventno)
161 {
162   // see header file for class documentation
163   int iResult=0;
164   for (unsigned int i=0; i<fESDs.size(); i++) {
165     if (fESDs[i]) {
166       int res=fESDs[i]->WriteESD(NULL, eventno);
167       if (res<0 && iResult>=0) iResult=res;
168     }
169   }
170   return iResult;
171 }
172
173 void AliHLTEsdManager::SetDirectory(const char* directory)
174 {
175   // see header file for class documentation
176   if (!directory) return;
177   fDirectory=directory;
178   for (unsigned int i=0; i<fESDs.size(); i++) {
179     if (fESDs[i]) {
180       fESDs[i]->SetDirectory(directory);
181     }
182   }
183 }
184
185 TString AliHLTEsdManager::GetFileNames(AliHLTComponentDataType dt) const
186 {
187   TString result;
188   for (unsigned int i=0; i<fESDs.size(); i++) {
189     if (fESDs[i] && *(fESDs[i])==dt) {
190       if (!result.IsNull()) result+=" ";
191       result+=fESDs[i]->GetFileName();
192     }
193   }
194   return result;
195 }
196
197 TTree* AliHLTEsdManager::EmbedIntoTree(AliESDEvent* pESD, const char* name, const char* title)
198 {
199   // see header file for class documentation
200   int iResult=0;
201   TTree* pTree=new TTree(name, title);
202   if (pTree) {
203     pESD->WriteToTree(pTree);
204     pTree->Fill();
205     pTree->GetUserInfo()->Add(pESD);
206   } else {
207     iResult=-ENOMEM;
208   }
209
210   if (iResult<0) {
211     pTree->GetUserInfo()->Clear();
212     delete pTree;
213   }
214
215   return pTree;
216 }
217
218 AliHLTEsdManager::AliHLTEsdListEntry::AliHLTEsdListEntry(AliHLTComponentDataType dt)
219   :
220   fName(),
221   fDirectory(),
222   fDt(dt),
223   fpFile(NULL),
224   fpTree(NULL),
225   fpEsd(NULL)
226 {
227   // see header file for class documentation
228 }
229
230 AliHLTEsdManager::AliHLTEsdListEntry::~AliHLTEsdListEntry()
231 {
232   // see header file for class documentation
233   if (fpEsd) delete fpEsd;
234   fpEsd=NULL;
235
236   if (fpTree) delete fpTree;
237   fpTree=NULL;
238
239   if (fpFile) delete fpFile;
240   fpFile=NULL;
241 }
242
243 bool AliHLTEsdManager::AliHLTEsdListEntry::operator==(AliHLTComponentDataType dt) const
244 {
245   // see header file for class documentation
246   return fDt==dt;
247 }
248
249 int AliHLTEsdManager::AliHLTEsdListEntry::WriteESD(AliESDEvent* pSrcESD, int eventno)
250 {
251   // see header file for class documentation
252   int iResult=0;
253 #ifdef HAVE_ESD_COPY
254   if (fName.IsNull()) {
255     // this is the first event, create the file name
256     TString origin;
257     origin.Insert(0, fDt.fOrigin, kAliHLTComponentDataTypefOriginSize);
258     origin.Remove(TString::kTrailing, ' ');
259     origin.ToUpper();
260     fName="";
261     if (!fDirectory.IsNull()) {
262       fName+=fDirectory; fName+="/";
263     }
264     fName+="AliHLT"; fName+=origin;
265     if (fDt!=kAliHLTDataTypeESDObject &&
266         fDt!=kAliHLTDataTypeESDTree) {
267
268       HLTWarning("non-standard ESD type %s", AliHLTComponent::DataType2Text(fDt).c_str());
269       TString id;
270       id.Insert(0, fDt.fID, kAliHLTComponentDataTypefIDsize);
271       id.Remove(TString::kTrailing, ' ');
272       id.ToUpper();
273       fName+="_"; fName+=id; fName+=".root";
274     } else {
275       fName+="ESDs.root";
276     }
277
278     fpFile=new TFile(fName, "RECREATE");
279     fpTree=new TTree("esdTree", "Tree with HLT ESD objects");
280     fpEsd=new AliESDEvent;
281     if (fpEsd) {
282       fpEsd->CreateStdContent();
283       if (fpTree) {
284         fpEsd->WriteToTree(fpTree);
285       }
286     }
287   }
288
289   if (fpFile && fpTree && fpEsd) {
290     // synchronize and add empty events
291     fpEsd->Reset();
292     int nofCurrentEvents=fpTree->GetEntries();
293     if (nofCurrentEvents<eventno) {
294       iResult=1; // indicate tree to be written
295       HLTDebug("adding %d empty events to file %s", eventno-nofCurrentEvents, fName.Data());
296       for (int i=nofCurrentEvents; i<eventno; i++) {
297         fpTree->Fill();
298       }
299     }
300
301     if (iResult>=0 && pSrcESD) {
302       *fpEsd=*pSrcESD;
303       fpTree->Fill();
304       iResult=1; // indicate tree to be written
305     }
306
307     if (iResult>0) {
308       fpFile->cd();
309       fpTree->GetUserInfo()->Add(fpEsd);
310       fpTree->Write(fpTree->GetName(),TObject::kOverwrite);
311       fpTree->GetUserInfo()->Clear();
312     }
313   }
314 #else //HAVE_ESD_COPY
315   // this is the old workaround, necessary for older AliRoot versions
316   // version<=v4-12-Release
317
318   // we need to copy the ESD, I did not find an approptiate
319   // method, the workaround is to save the ESD in a temporary
320   // tree, read the content back into the ESD structure
321   // used for filling.
322   // Unfortunately the following code crashes at the second event.
323   // The expert on the ESD (Christian Klein Boesig) does not have
324   // a solution either. It seems to be a problem in ROOT.
325   //  TTree* dummy=new TTree("dummy","dummy");
326   //  dummy->SetDirectory(0);
327   //  pESD->WriteToTree(dummy);
328   //  dummy->Fill();
329   //  dummy->GetUserInfo()->Add(pESD);
330   //  fpEsd->ReadFromTree(dummy);
331   //  dummy->GetEvent(0);
332   //  fpEsd->WriteToTree(fpTree);
333   //  fpTree->Fill();
334   //  dummy->GetUserInfo()->Clear();
335   //  delete dummy;
336   //
337   // The only way is via TChain, which is working on files only at the
338   // time of writing.
339   // We use temporary files for the new event to be copied into the
340   // existing tree.
341   //
342   if (fName.IsNull()) {
343     // this is the first event, create the file on disk and write ESD
344     TString origin;
345     origin.Insert(0, fDt.fOrigin, kAliHLTComponentDataTypefOriginSize);
346     origin.Remove(TString::kTrailing, ' ');
347     origin.ToUpper();
348     fName="";
349     if (!fDirectory.IsNull()) {
350       fName+=fDirectory; fName+="/";
351     }
352     fName+="AliHLT"; fName+=origin;
353     if (fDt!=kAliHLTDataTypeESDObject &&
354         fDt!=kAliHLTDataTypeESDTree) {
355
356       HLTWarning("non-standard ESD type %s", AliHLTComponent::DataType2Text(fDt).c_str());
357       TString id;
358       id.Insert(0, fDt.fID, kAliHLTComponentDataTypefIDsize);
359       id.Remove(TString::kTrailing, ' ');
360       id.ToUpper();
361       fName+="_"; fName+=id; fName+=".root";
362     } else {
363       fName+="ESDs.root";
364     }
365
366     if (!gSystem->AccessPathName(fName)) {
367       // file exists, delete
368       TString shellcmd="rm -f ";
369       shellcmd+=fName;
370       gSystem->Exec(shellcmd);
371     }
372   }
373
374   TChain chain("esdTree");
375   TList cleanup;
376   cleanup.SetOwner();
377
378   int nofCurrentEvents=0;
379   if (iResult>=0) {
380     if (!gSystem->AccessPathName(fName)) {
381       // these are the other events, use the target file and temporary files to merge
382       // with TChain
383       chain.Add(fName);
384
385       if (eventno>=0) {
386         TFile file(fName);
387         if (!file.IsZombie()) {
388           TTree* pSrcTree;
389           file.GetObject("esdTree", pSrcTree);
390           if (pSrcTree) {
391             nofCurrentEvents=pSrcTree->GetEntries();
392           }
393           file.Close();
394         }
395       }
396     }
397   }
398
399   // synchronize and add empty events
400   if (nofCurrentEvents<eventno) {
401     iResult=1; // indicate files to merge
402     TTree* pTgtTree=new TTree("esdTree", "Tree with HLT ESD objects");
403     if (pTgtTree) {
404       pTgtTree->SetDirectory(0);
405       AliESDEvent* pTmpESD=new AliESDEvent;
406       if (pTmpESD) {
407         TString tmpfilename;
408         FILE* pTmpFile=gSystem->TempFileName(tmpfilename);
409         if (pTmpFile) {
410           fclose(pTmpFile);
411           pTmpFile=NULL;
412           cleanup.Add(new TObjString(tmpfilename));
413           TFile emptyevents(tmpfilename, "RECREATE");
414           if (!emptyevents.IsZombie()) {
415             pTmpESD->CreateStdContent();
416             pTmpESD->WriteToTree(pTgtTree);
417             HLTDebug("adding %d empty events to file %s", eventno-nofCurrentEvents, fName.Data());
418             for (int i=nofCurrentEvents; i<eventno; i++) {
419               pTgtTree->Fill();
420             }
421             pTgtTree->GetUserInfo()->Add(pTmpESD);
422             emptyevents.cd();
423             pTgtTree->Write();
424             emptyevents.Close();
425             chain.Add(tmpfilename);
426             pTgtTree->GetUserInfo()->Clear();
427           }
428         }
429         delete pTmpESD;
430       } else {
431         iResult=-ENOMEM;
432       }
433       delete pTgtTree;
434     } else {
435       iResult=-ENOMEM;
436     }
437   }
438
439   if (iResult>=0 && pSrcESD) {
440     // add the new event to the chain
441     iResult=1; // indicate files to merge
442     TString tmpfilename=WriteTempFile(pSrcESD);
443     if (!tmpfilename.IsNull()) {
444       chain.Add(tmpfilename);
445       cleanup.Add(new TObjString(tmpfilename));
446     }
447   }
448
449   if (iResult>0) {
450     // build temporary file name for chain output
451     TString tgtName;
452     FILE* pTmpFile=gSystem->TempFileName(tgtName);
453     if (pTmpFile) {
454       fclose(pTmpFile);
455       pTmpFile=NULL;
456
457       // there have been problems with the memory consumption when using
458       // TChain::Merge
459       // but using a separate loop soemtimes crashes in AliESDEvent::ReadFromTree
460       // since this is for backward compatiblity only, we take the TChain::Merge
461       chain.Merge(tgtName);
462 //       TFile tgtFile(tgtName, "RECREATE");
463 //       TTree* pTgtTree=new TTree("esdTree", "Tree with HLT ESD objects");
464 //       AliESDEvent* pTgtESD=new AliESDEvent;
465 //       if (pTgtTree && pTgtESD) {
466 //      pTgtESD->ReadFromTree(&chain);
467 //      pTgtESD->WriteToTree(pTgtTree);
468
469 //      int nofEvents=chain.GetEntries();
470 //      for (int event=0; event<nofEvents; event++) {
471 //        chain.GetEntry(event);
472 //        pTgtTree->Fill();
473 //      }
474
475 //      pTgtTree->GetUserInfo()->Add(pTgtESD);
476 //      tgtFile.cd();
477 //      pTgtTree->Write();
478 //      pTgtTree->GetUserInfo()->Clear();
479 //       } else {
480 //      iResult=-ENOMEM;
481 //       }
482
483 //       if (pTgtTree) delete pTgtTree;
484 //       if (pTgtESD) delete pTgtESD;
485 //       tgtFile.Close();
486
487       // rename the merged file to the original file
488       TString shellcmd="mv ";
489       shellcmd+=tgtName + " " + fName;
490       if (gSystem->Exec(shellcmd)==0) {
491         HLTDebug("renaming %s to %s", tgtName.Data(), fName.Data());
492       } else {
493         HLTError("can not rename temporary file %s to %s", tgtName.Data(), fName.Data());
494       }
495     } else {
496       HLTError("can not get temporary file name from system");
497       iResult=-EBADF;
498     }
499   }
500
501   // delete temporary files
502   // the list objects are cleaned up by the TList destructor as the
503   // list is owner
504   TIter entry(&cleanup);
505   while (TObject* pObj=entry.Next()) {
506     if (dynamic_cast<TObjString*>(pObj)) {
507       TString shellcmd="rm -f ";
508       shellcmd+=(dynamic_cast<TObjString*>(pObj))->GetString();
509       gSystem->Exec(shellcmd);
510     }
511   }
512 #endif //HAVE_ESD_COPY
513
514   return iResult;
515 }
516
517 TString AliHLTEsdManager::AliHLTEsdListEntry::WriteTempFile(AliESDEvent* pESD) const
518 {
519   // see header file for class documentation
520   int iResult;
521   TString tmpfilename;
522   FILE* pTmpFile=gSystem->TempFileName(tmpfilename);
523   if (pTmpFile) {
524     fclose(pTmpFile);
525     pTmpFile=NULL;
526
527     TFile file(tmpfilename, "RECREATE");
528     if (!file.IsZombie()) {
529       TTree* pTree=AliHLTEsdManager::EmbedIntoTree(pESD);
530       if (pTree) {
531         file.cd();
532         if (pTree->Write()>0) {
533         } else {
534           HLTError("can not write esd tree to temporary file %s", tmpfilename.Data());
535         }
536
537         pTree->GetUserInfo()->Clear();
538         delete pTree;
539       } else {
540         iResult=-ENOMEM;
541       }
542       file.Close();
543     } else {
544       HLTError("can not open file %s", tmpfilename.Data());
545       iResult=-EBADF;
546     }
547   } else {
548     HLTError("can not get temporary file name from system");
549     iResult=-EBADF;
550   }
551
552   if (iResult<0) {
553     if (gSystem->AccessPathName(tmpfilename)==0) {
554       TString shellcmd="rm -f ";
555       shellcmd+=tmpfilename;
556       gSystem->Exec(shellcmd);
557     }
558     tmpfilename="";
559   }
560   return tmpfilename;
561 }
562
563 void AliHLTEsdManager::AliHLTEsdListEntry::SetDirectory(const char* directory)
564 {
565   // see header file for class documentation
566   if (!directory) return;
567   if (!fName.IsNull()) {
568     HLTWarning("ESD entry already in writing mode (%s), ignoring directory", fName.Data());
569     return;
570   }
571   fDirectory=directory;
572 }
573
574 void AliHLTEsdManager::AliHLTEsdListEntry::Delete()
575 {
576   // see header file for class documentation
577   if (fName.IsNull()) return;
578   if (gSystem->AccessPathName(fName)!=0) return;
579
580   TString shellcmd="rm -f ";
581   shellcmd+=fName;
582   gSystem->Exec(shellcmd);
583   fName="";
584 }
585
586 const char* AliHLTEsdManager::AliHLTEsdListEntry::GetFileName() const
587 {
588   // see header file for class documentation
589   return fName.Data();
590 }