3 //**************************************************************************
4 //* This file is property of and copyright by the ALICE HLT Project *
5 //* ALICE Experiment at CERN, All rights reserved. *
7 //* Primary Authors: Matthias Richter <Matthias.Richter@ift.uib.no> *
8 //* for The ALICE HLT Project. *
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 //**************************************************************************
19 /** @file AliHLTEsdManager.cxx
20 @author Matthias Richter
22 @brief Manager for merging and writing of HLT ESDs
25 #include "AliHLTEsdManager.h"
26 #include "AliHLTComponent.h"
27 #include "AliESDEvent.h"
28 #include "AliHLTMessage.h"
29 #include "AliESDEvent.h"
34 #include "TObjectTable.h"
39 /** ROOT macro for the implementation of ROOT specific class methods */
40 ClassImp(AliHLTEsdManager)
42 AliHLTEsdManager::AliHLTEsdManager()
47 // see header file for class documentation
49 // refer to README to build package
51 // visit http://web.ift.uib.no/~kjeks/doc/alice-hlt
54 AliHLTEsdManager::~AliHLTEsdManager()
56 // see header file for class documentation
57 for (unsigned int i=0; i<fESDs.size(); i++) {
65 AliHLTEsdManager::AliHLTEsdListEntry* AliHLTEsdManager::Find(AliHLTComponentDataType dt) const
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]);
77 int AliHLTEsdManager::WriteESD(const AliHLTUInt8_t* pBuffer, AliHLTUInt32_t size,
78 AliHLTComponentDataType dt, AliESDEvent* tgtesd, int eventno)
80 // see header file for class documentation
81 if (!pBuffer && size<=0) return -EINVAL;
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);
94 pTree=dynamic_cast<TTree*>(pObj);
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());
103 pESD->ReadFromTree(pTree);
107 HLTWarning("tree of data block %s has no events, skipping data block", AliHLTComponent::DataType2Text(dt).c_str());
111 AliHLTEsdListEntry* entry=Find(dt);
113 AliHLTEsdListEntry* newEntry=new AliHLTEsdListEntry(dt);
114 if (!fDirectory.IsNull()) {
115 newEntry->SetDirectory(fDirectory);
117 fESDs.push_back(newEntry);
123 TTree* pTmpTree=AliHLTEsdManager::EmbedIntoTree(pESD);
125 tgtesd->ReadFromTree(pTmpTree);
126 pTmpTree->GetEvent(0);
127 pTmpTree->GetUserInfo()->Clear();
129 HLTDebug("data block %s written to target ESD", AliHLTComponent::DataType2Text(dt).c_str());
137 entry->WriteESD(pESD, eventno);
139 HLTError("internal mismatch, can not create list entry");
144 HLTWarning("data block %s is not of class type AliESDEvent, ignoring ...", AliHLTComponent::DataType2Text(dt).c_str());
147 // ESD has been created and must be cleaned up
159 int AliHLTEsdManager::PadESDs(int eventno)
161 // see header file for class documentation
163 for (unsigned int i=0; i<fESDs.size(); i++) {
165 int res=fESDs[i]->WriteESD(NULL, eventno);
166 if (res<0 && iResult>=0) iResult=res;
172 void AliHLTEsdManager::SetDirectory(const char* directory)
174 // see header file for class documentation
175 if (!directory) return;
176 fDirectory=directory;
177 for (unsigned int i=0; i<fESDs.size(); i++) {
179 fESDs[i]->SetDirectory(directory);
184 TString AliHLTEsdManager::GetFileNames(AliHLTComponentDataType dt) const
187 for (unsigned int i=0; i<fESDs.size(); i++) {
188 if (fESDs[i] && *(fESDs[i])==dt) {
189 if (!result.IsNull()) result+=" ";
190 result+=fESDs[i]->GetFileName();
196 TTree* AliHLTEsdManager::EmbedIntoTree(AliESDEvent* pESD, const char* name, const char* title)
198 // see header file for class documentation
200 TTree* pTree=new TTree(name, title);
202 pESD->WriteToTree(pTree);
204 pTree->GetUserInfo()->Add(pESD);
210 pTree->GetUserInfo()->Clear();
217 AliHLTEsdManager::AliHLTEsdListEntry::AliHLTEsdListEntry(AliHLTComponentDataType dt)
226 // see header file for class documentation
229 AliHLTEsdManager::AliHLTEsdListEntry::~AliHLTEsdListEntry()
231 // see header file for class documentation
234 bool AliHLTEsdManager::AliHLTEsdListEntry::operator==(AliHLTComponentDataType dt) const
236 // see header file for class documentation
240 int AliHLTEsdManager::AliHLTEsdListEntry::WriteESD(AliESDEvent* pSrcESD, int eventno)
242 // see header file for class documentation
245 if (fName.IsNull()) {
246 // this is the first event, create the file name
248 origin.Insert(0, fDt.fOrigin, kAliHLTComponentDataTypefOriginSize);
249 origin.Remove(TString::kTrailing, ' ');
252 if (!fDirectory.IsNull()) {
253 fName+=fDirectory; fName+="/";
255 fName+="AliHLT"; fName+=origin;
256 if (fDt!=kAliHLTDataTypeESDObject &&
257 fDt!=kAliHLTDataTypeESDTree) {
259 HLTWarning("non-standard ESD type %s", AliHLTComponent::DataType2Text(fDt).c_str());
261 id.Insert(0, fDt.fID, kAliHLTComponentDataTypefIDsize);
262 id.Remove(TString::kTrailing, ' ');
264 fName+="_"; fName+=id; fName+=".root";
269 fpFile=new TFile(fName, "RECREATE");
270 fpTree=new TTree("esdTree", "Tree with HLT ESD objects");
271 fpEsd=new AliESDEvent;
273 fpEsd->CreateStdContent();
275 fpEsd->WriteToTree(fpTree);
280 if (fpFile && fpTree && fpEsd) {
281 // synchronize and add empty events
283 int nofCurrentEvents=fpTree->GetEntries();
284 if (nofCurrentEvents<eventno) {
285 iResult=1; // indicate tree to be written
286 HLTDebug("adding %d empty events to file %s", eventno-nofCurrentEvents, fName.Data());
287 for (int i=nofCurrentEvents; i<eventno; i++) {
292 if (iResult>=0 && pSrcESD) {
295 iResult=1; // indicate tree to be written
300 fpTree->GetUserInfo()->Add(fpEsd);
301 fpTree->Write(fpTree->GetName(),TObject::kOverwrite);
302 fpTree->GetUserInfo()->Clear();
305 #else //HAVE_ESD_COPY
306 // this is the old workaround, necessary for older AliRoot versions
307 // version<=v4-12-Release
309 // we need to copy the ESD, I did not find an approptiate
310 // method, the workaround is to save the ESD in a temporary
311 // tree, read the content back into the ESD structure
313 // Unfortunately the following code crashes at the second event.
314 // The expert on the ESD (Christian Klein Boesig) does not have
315 // a solution either. It seems to be a problem in ROOT.
316 // TTree* dummy=new TTree("dummy","dummy");
317 // dummy->SetDirectory(0);
318 // pESD->WriteToTree(dummy);
320 // dummy->GetUserInfo()->Add(pESD);
321 // fpEsd->ReadFromTree(dummy);
322 // dummy->GetEvent(0);
323 // fpEsd->WriteToTree(fpTree);
325 // dummy->GetUserInfo()->Clear();
328 // The only way is via TChain, which is working on files only at the
330 // We use temporary files for the new event to be copied into the
333 if (fName.IsNull()) {
334 // this is the first event, create the file on disk and write ESD
336 origin.Insert(0, fDt.fOrigin, kAliHLTComponentDataTypefOriginSize);
337 origin.Remove(TString::kTrailing, ' ');
340 if (!fDirectory.IsNull()) {
341 fName+=fDirectory; fName+="/";
343 fName+="AliHLT"; fName+=origin;
344 if (fDt!=kAliHLTDataTypeESDObject &&
345 fDt!=kAliHLTDataTypeESDTree) {
347 HLTWarning("non-standard ESD type %s", AliHLTComponent::DataType2Text(fDt).c_str());
349 id.Insert(0, fDt.fID, kAliHLTComponentDataTypefIDsize);
350 id.Remove(TString::kTrailing, ' ');
352 fName+="_"; fName+=id; fName+=".root";
357 if (!gSystem->AccessPathName(fName)) {
358 // file exists, delete
359 TString shellcmd="rm -f ";
361 gSystem->Exec(shellcmd);
365 TChain chain("esdTree");
369 int nofCurrentEvents=0;
371 if (!gSystem->AccessPathName(fName)) {
372 // these are the other events, use the target file and temporary files to merge
378 if (!file.IsZombie()) {
380 file.GetObject("esdTree", pSrcTree);
382 nofCurrentEvents=pSrcTree->GetEntries();
390 // synchronize and add empty events
391 if (nofCurrentEvents<eventno) {
392 iResult=1; // indicate files to merge
393 TTree* pTgtTree=new TTree("esdTree", "Tree with HLT ESD objects");
395 pTgtTree->SetDirectory(0);
396 AliESDEvent* pTmpESD=new AliESDEvent;
399 FILE* pTmpFile=gSystem->TempFileName(tmpfilename);
403 cleanup.Add(new TObjString(tmpfilename));
404 TFile emptyevents(tmpfilename, "RECREATE");
405 if (!emptyevents.IsZombie()) {
406 pTmpESD->CreateStdContent();
407 pTmpESD->WriteToTree(pTgtTree);
408 HLTDebug("adding %d empty events to file %s", eventno-nofCurrentEvents, fName.Data());
409 for (int i=nofCurrentEvents; i<eventno; i++) {
412 pTgtTree->GetUserInfo()->Add(pTmpESD);
416 chain.Add(tmpfilename);
417 pTgtTree->GetUserInfo()->Clear();
430 if (iResult>=0 && pSrcESD) {
431 // add the new event to the chain
432 iResult=1; // indicate files to merge
433 TString tmpfilename=WriteTempFile(pSrcESD);
434 if (!tmpfilename.IsNull()) {
435 chain.Add(tmpfilename);
436 cleanup.Add(new TObjString(tmpfilename));
441 // build temporary file name for chain output
443 FILE* pTmpFile=gSystem->TempFileName(tgtName);
448 // there have been problems with the memory consumption when using
450 // but using a separate loop soemtimes crashes in AliESDEvent::ReadFromTree
451 // since this is for backward compatiblity only, we take the TChain::Merge
452 chain.Merge(tgtName);
453 // TFile tgtFile(tgtName, "RECREATE");
454 // TTree* pTgtTree=new TTree("esdTree", "Tree with HLT ESD objects");
455 // AliESDEvent* pTgtESD=new AliESDEvent;
456 // if (pTgtTree && pTgtESD) {
457 // pTgtESD->ReadFromTree(&chain);
458 // pTgtESD->WriteToTree(pTgtTree);
460 // int nofEvents=chain.GetEntries();
461 // for (int event=0; event<nofEvents; event++) {
462 // chain.GetEntry(event);
466 // pTgtTree->GetUserInfo()->Add(pTgtESD);
468 // pTgtTree->Write();
469 // pTgtTree->GetUserInfo()->Clear();
474 // if (pTgtTree) delete pTgtTree;
475 // if (pTgtESD) delete pTgtESD;
478 // rename the merged file to the original file
479 TString shellcmd="mv ";
480 shellcmd+=tgtName + " " + fName;
481 if (gSystem->Exec(shellcmd)==0) {
482 HLTDebug("renaming %s to %s", tgtName.Data(), fName.Data());
484 HLTError("can not rename temporary file %s to %s", tgtName.Data(), fName.Data());
487 HLTError("can not get temporary file name from system");
492 // delete temporary files
493 // the list objects are cleaned up by the TList destructor as the
495 TIter entry(&cleanup);
496 while (TObject* pObj=entry.Next()) {
497 if (dynamic_cast<TObjString*>(pObj)) {
498 TString shellcmd="rm -f ";
499 shellcmd+=(dynamic_cast<TObjString*>(pObj))->GetString();
500 gSystem->Exec(shellcmd);
503 #endif //HAVE_ESD_COPY
508 TString AliHLTEsdManager::AliHLTEsdListEntry::WriteTempFile(AliESDEvent* pESD) const
510 // see header file for class documentation
513 FILE* pTmpFile=gSystem->TempFileName(tmpfilename);
518 TFile file(tmpfilename, "RECREATE");
519 if (!file.IsZombie()) {
520 TTree* pTree=AliHLTEsdManager::EmbedIntoTree(pESD);
523 if (pTree->Write()>0) {
525 HLTError("can not write esd tree to temporary file %s", tmpfilename.Data());
528 pTree->GetUserInfo()->Clear();
535 HLTError("can not open file %s", tmpfilename.Data());
539 HLTError("can not get temporary file name from system");
544 if (gSystem->AccessPathName(tmpfilename)==0) {
545 TString shellcmd="rm -f ";
546 shellcmd+=tmpfilename;
547 gSystem->Exec(shellcmd);
554 void AliHLTEsdManager::AliHLTEsdListEntry::SetDirectory(const char* directory)
556 // see header file for class documentation
557 if (!directory) return;
558 if (!fName.IsNull()) {
559 HLTWarning("ESD entry already in writing mode (%s), ignoring directory", fName.Data());
562 fDirectory=directory;
565 void AliHLTEsdManager::AliHLTEsdListEntry::Delete()
567 // see header file for class documentation
568 if (fName.IsNull()) return;
569 if (gSystem->AccessPathName(fName)!=0) return;
571 TString shellcmd="rm -f ";
573 gSystem->Exec(shellcmd);
577 const char* AliHLTEsdManager::AliHLTEsdListEntry::GetFileName() const
579 // see header file for class documentation