]> git.uio.no Git - u/mrichter/AliRoot.git/blob - HLT/rec/AliHLTEsdManagerImplementation.cxx
several bugfixes in order to get the HLTGlobalTrigger decision safely into the ESD
[u/mrichter/AliRoot.git] / HLT / rec / AliHLTEsdManagerImplementation.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   AliHLTEsdManagerImplementation.cxx
20     @author Matthias Richter
21     @date   
22     @brief  Manager for merging and writing of HLT ESDs
23 */
24
25 #include "AliHLTEsdManagerImplementation.h"
26 #include "AliHLTComponent.h"
27 #include "AliESDEvent.h"
28 #include "AliHLTMessage.h"
29 #include "AliESDEvent.h"
30 #include "AliESDtrack.h"
31 #include "TFile.h"
32 #include "TTree.h"
33 #include "TClass.h"
34 #include "TObject.h"
35 #include "TObjectTable.h"
36 #include "TSystem.h"
37 #include "TChain.h"
38 #include "TList.h"
39
40 /** ROOT macro for the implementation of ROOT specific class methods */
41 ClassImp(AliHLTEsdManagerImplementation)
42
43 AliHLTEsdManagerImplementation::AliHLTEsdManagerImplementation()
44   :
45   fESDs(),
46   fDirectory()
47 {
48   // see header file for class documentation
49   // or
50   // refer to README to build package
51   // or
52   // visit http://web.ift.uib.no/~kjeks/doc/alice-hlt
53 }
54
55 AliHLTEsdManagerImplementation::~AliHLTEsdManagerImplementation()
56 {
57   // see header file for class documentation
58   for (unsigned int i=0; i<fESDs.size(); i++) {
59     if (fESDs[i]) {
60       delete fESDs[i];
61     }
62     fESDs[i]=NULL;
63   }
64 }
65
66 AliHLTEsdManagerImplementation::AliHLTEsdListEntry* AliHLTEsdManagerImplementation::Find(AliHLTComponentDataType dt) const
67 {
68   // see header file for class documentation
69   AliHLTEsdListEntry* pEntry=NULL;
70   for (unsigned int i=0; i<fESDs.size(); i++) {
71     if (fESDs[i] && *(fESDs[i])==dt) {
72       pEntry=const_cast<AliHLTEsdListEntry*>(fESDs[i]);
73     }
74   }
75   return pEntry;
76 }
77
78 int AliHLTEsdManagerImplementation::WriteESD(const AliHLTUInt8_t* pBuffer, AliHLTUInt32_t size,
79                                AliHLTComponentDataType dt, AliESDEvent* tgtesd, int eventno)
80 {
81   // see header file for class documentation
82   if (!pBuffer && size<=0) return -EINVAL;
83   int iResult=0;
84   AliHLTUInt32_t firstWord=*((AliHLTUInt32_t*)pBuffer);
85   if (firstWord==size-sizeof(AliHLTUInt32_t)) {
86     HLTDebug("create object from block %s size %d", AliHLTComponent::DataType2Text(dt).c_str(), size);
87     AliHLTMessage msg(const_cast<AliHLTUInt8_t*>(pBuffer), size);
88     TClass* objclass=msg.GetClass();
89     TObject* pObj=msg.ReadObject(objclass);
90     if (pObj && objclass) {
91       HLTDebug("object %p type %s created", pObj, objclass->GetName());
92       AliESDEvent* pESD=dynamic_cast<AliESDEvent*>(pObj);
93       TTree* pTree=NULL;
94       if (!pESD) {
95         pTree=dynamic_cast<TTree*>(pObj);
96         if (pTree) {
97           pESD=new AliESDEvent;
98           pESD->CreateStdContent();
99           if (pTree->GetEntries()>0) {
100             if (pTree->GetEntries()>1) {
101               HLTWarning("only one entry allowed for ESD embedded into tree, data block %s contains tree with %d entires, taking first entry",
102                          AliHLTComponent::DataType2Text(dt).c_str(), pTree->GetEntries());
103             }
104             pESD->ReadFromTree(pTree);
105             pTree->GetEvent(0);
106           }
107         } else {
108           HLTWarning("tree of data block %s has no events, skipping data block", AliHLTComponent::DataType2Text(dt).c_str());
109         }
110       }
111       if (pESD) {
112         AliHLTEsdListEntry* entry=Find(dt);
113         if (!entry) {
114           if ((entry=new AliHLTEsdListEntry(dt))!=NULL) {
115             if (!fDirectory.IsNull()) {
116               entry->SetDirectory(fDirectory);
117             }
118             fESDs.push_back(entry);
119           }
120         }
121         if (entry) {
122           if (tgtesd) {
123 #if !defined(HAVE_NOT_ESD_COPY)
124             Merge(tgtesd, pESD);
125 #else //HAVE_NOT_ESD_COPY
126             static bool warningPrinted=false;
127             if (!warningPrinted) {
128               HLTWarning("old version of AliESDEvent does not provide assignment operator, skip merging to global hltEsd");
129             }
130             warningPrinted=true;
131 #endif //HAVE_NOT_ESD_COPY
132           }
133           entry->WriteESD(pESD, eventno);
134         } else {
135           HLTError("internal mismatch, can not create list entry");
136           iResult=-ENOMEM;
137         }
138       } else {
139         HLTWarning("data block %s is not of class type AliESDEvent, ignoring ...", AliHLTComponent::DataType2Text(dt).c_str());
140       }
141       if (pTree) {
142         // ESD has been created and must be cleaned up
143         pESD->Reset();
144         delete pESD;
145         pESD=NULL;
146       }
147       delete pObj;
148       pObj=NULL;
149     } else {
150     }
151   }
152   return iResult;
153 }
154
155 int AliHLTEsdManagerImplementation::PadESDs(int eventno)
156 {
157   // see header file for class documentation
158   int iResult=0;
159   for (unsigned int i=0; i<fESDs.size(); i++) {
160     if (fESDs[i]) {
161       int res=fESDs[i]->WriteESD(NULL, eventno);
162       if (res<0 && iResult>=0) iResult=res;
163     }
164   }
165   return iResult;
166 }
167
168 void AliHLTEsdManagerImplementation::SetDirectory(const char* directory)
169 {
170   // see header file for class documentation
171   if (!directory) return;
172   fDirectory=directory;
173   for (unsigned int i=0; i<fESDs.size(); i++) {
174     if (fESDs[i]) {
175       fESDs[i]->SetDirectory(directory);
176     }
177   }
178 }
179
180 TString AliHLTEsdManagerImplementation::GetFileNames(AliHLTComponentDataType dt) const
181 {
182   TString result;
183   for (unsigned int i=0; i<fESDs.size(); i++) {
184     if (fESDs[i] && *(fESDs[i])==dt) {
185       if (!result.IsNull()) result+=" ";
186       result+=fESDs[i]->GetFileName();
187     }
188   }
189   return result;
190 }
191
192 TTree* AliHLTEsdManagerImplementation::EmbedIntoTree(AliESDEvent* pESD, const char* name, const char* title)
193 {
194   // see header file for class documentation
195   int iResult=0;
196   TTree* pTree=new TTree(name, title);
197   if (pTree) {
198     pESD->WriteToTree(pTree);
199     pTree->Fill();
200     pTree->GetUserInfo()->Add(pESD);
201   } else {
202     iResult=-ENOMEM;
203   }
204
205   if (iResult<0) {
206     pTree->GetUserInfo()->Clear();
207     delete pTree;
208   }
209
210   return pTree;
211 }
212
213 AliHLTEsdManagerImplementation::AliHLTEsdListEntry::AliHLTEsdListEntry(AliHLTComponentDataType dt)
214   :
215   fName(),
216   fDirectory(),
217   fDt(dt),
218   fpFile(NULL),
219   fpTree(NULL),
220   fpEsd(NULL),
221   fPrefix()
222 {
223   // see header file for class documentation
224 }
225
226 AliHLTEsdManagerImplementation::AliHLTEsdListEntry::~AliHLTEsdListEntry()
227 {
228   // see header file for class documentation
229   if (fpEsd) delete fpEsd;
230   fpEsd=NULL;
231
232   if (fpTree) delete fpTree;
233   fpTree=NULL;
234
235   if (fpFile) {
236     fpFile->Close();
237     delete fpFile;
238   }
239   fpFile=NULL;
240 }
241
242 bool AliHLTEsdManagerImplementation::AliHLTEsdListEntry::operator==(AliHLTComponentDataType dt) const
243 {
244   // see header file for class documentation
245   return fDt==dt;
246 }
247
248 int AliHLTEsdManagerImplementation::AliHLTEsdListEntry::WriteESD(AliESDEvent* pSrcESD, int eventno)
249 {
250   // see header file for class documentation
251   int iResult=0;
252
253   // Matthias 2009-06-06: writing of individual ESD files for the different origins was a
254   // first attempt when functionality was missing in the AliRoot framework and remained as
255   // debugging feature. ESD merging is now implemented and data written to the hltEsd, so
256   // the feature is now disabled because it causes increasing memory consumption. Presumably
257   // not because of a memory leak but the way the internal TTree is used and kept in memory. 
258   return 0;
259
260 #ifndef HAVE_NOT_ESD_COPY
261   if (fName.IsNull()) {
262     // this is the first event, create the file name
263     fName="";
264     if (!fDirectory.IsNull()) {
265       fName+=fDirectory; fName+="/";
266     }
267     fName+="Ali"; fName+=GetPrefix();
268     if (fDt!=kAliHLTDataTypeESDObject &&
269         fDt!=kAliHLTDataTypeESDTree) {
270
271       HLTWarning("non-standard ESD type %s", AliHLTComponent::DataType2Text(fDt).c_str());
272       TString id;
273       id.Insert(0, fDt.fID, kAliHLTComponentDataTypefIDsize);
274       id.Remove(TString::kTrailing, ' ');
275       id.ToUpper();
276       fName+="_"; fName+=id; fName+=".root";
277     } else {
278       fName+="ESDs.root";
279     }
280
281     fpFile=new TFile(fName, "RECREATE");
282     fpTree=new TTree("esdTree", "Tree with HLT ESD objects");
283     fpTree->SetDirectory(0);
284     fpEsd=new AliESDEvent;
285     if (fpEsd) {
286       fpEsd->CreateStdContent();
287       *fpEsd=*pSrcESD;
288       if (fpTree) {
289         fpEsd->WriteToTree(fpTree);
290       }
291     }
292   }
293
294   if (fpFile && fpTree && fpEsd) {
295     // synchronize and add empty events
296     fpEsd->Reset();
297     int nofCurrentEvents=fpTree->GetEntries();
298     if (nofCurrentEvents<eventno) {
299       iResult=1; // indicate tree to be written
300       HLTDebug("adding %d empty events to file %s", eventno-nofCurrentEvents, fName.Data());
301       for (int i=nofCurrentEvents; i<eventno; i++) {
302         fpTree->Fill();
303       }
304     }
305
306     if (iResult>=0 && pSrcESD) {
307       int nofObjects=fpEsd->GetList()->GetEntries();
308       *fpEsd=*pSrcESD;
309       if (nofObjects!=fpEsd->GetList()->GetEntries()) {
310         // The source ESD contains object not present in the target ESD
311         // before. Those objects will not be written to the tree since
312         // the branch layout has been created earlier.
313         // Create new tree with the additional branches, copy the entries
314         // of the current tree into the new tree, and continue.
315         TTree* pNewTree=new TTree("esdTree", "Tree with HLT ESD objects");
316         pNewTree->SetDirectory(0);
317         AliESDEvent* readESD=new AliESDEvent;
318         readESD->CreateStdContent();
319         readESD->ReadFromTree(fpTree);
320         fpEsd->Reset();
321         fpEsd->WriteToTree(pNewTree);
322         HLTDebug("cloning tree with %d entries", fpTree->GetEntries());
323         for (int event=0; event<fpTree->GetEntries(); event++) {
324           fpTree->GetEntry(event);
325           *fpEsd=*readESD;
326           pNewTree->Fill();
327           fpEsd->Reset();
328         }
329         fpFile->Close();
330         delete fpFile;
331         delete readESD;
332         delete fpTree;
333         fpFile=new TFile(fName, "RECREATE");
334         fpTree=pNewTree;
335         *fpEsd=*pSrcESD;
336         HLTDebug("new ESD with %d objects", fpEsd->GetList()->GetEntries());
337       }
338       fpTree->Fill();
339       iResult=1; // indicate tree to be written
340     }
341
342     if (iResult>0) {
343       fpFile->cd();
344       fpTree->GetUserInfo()->Add(fpEsd);
345       fpTree->Write(fpTree->GetName(),TObject::kOverwrite);
346       fpTree->GetUserInfo()->Clear();
347     }
348   }
349 #else //HAVE_NOT_ESD_COPY
350   // this is the old workaround, necessary for older AliRoot versions
351   // version<=v4-12-Release
352
353   // we need to copy the ESD, I did not find an approptiate
354   // method, the workaround is to save the ESD in a temporary
355   // tree, read the content back into the ESD structure
356   // used for filling.
357   // Unfortunately the following code crashes at the second event.
358   // The expert on the ESD (Christian Klein Boesig) does not have
359   // a solution either. It seems to be a problem in ROOT.
360   //  TTree* dummy=new TTree("dummy","dummy");
361   //  dummy->SetDirectory(0);
362   //  pESD->WriteToTree(dummy);
363   //  dummy->Fill();
364   //  dummy->GetUserInfo()->Add(pESD);
365   //  fpEsd->ReadFromTree(dummy);
366   //  dummy->GetEvent(0);
367   //  fpEsd->WriteToTree(fpTree);
368   //  fpTree->Fill();
369   //  dummy->GetUserInfo()->Clear();
370   //  delete dummy;
371   //
372   // The only way is via TChain, which is working on files only at the
373   // time of writing.
374   // We use temporary files for the new event to be copied into the
375   // existing tree.
376   //
377   if (fName.IsNull()) {
378     // this is the first event, create the file on disk and write ESD
379     TString origin;
380     origin.Insert(0, fDt.fOrigin, kAliHLTComponentDataTypefOriginSize);
381     origin.Remove(TString::kTrailing, ' ');
382     origin.ToUpper();
383     fName="";
384     if (!fDirectory.IsNull()) {
385       fName+=fDirectory; fName+="/";
386     }
387     fName+="AliHLT"; fName+=origin;
388     if (fDt!=kAliHLTDataTypeESDObject &&
389         fDt!=kAliHLTDataTypeESDTree) {
390
391       HLTWarning("non-standard ESD type %s", AliHLTComponent::DataType2Text(fDt).c_str());
392       TString id;
393       id.Insert(0, fDt.fID, kAliHLTComponentDataTypefIDsize);
394       id.Remove(TString::kTrailing, ' ');
395       id.ToUpper();
396       fName+="_"; fName+=id; fName+=".root";
397     } else {
398       fName+="ESDs.root";
399     }
400
401     if (!gSystem->AccessPathName(fName)) {
402       // file exists, delete
403       TString shellcmd="rm -f ";
404       shellcmd+=fName;
405       gSystem->Exec(shellcmd);
406     }
407   }
408
409   TChain chain("esdTree");
410   TList cleanup;
411   cleanup.SetOwner();
412
413   int nofCurrentEvents=0;
414   if (iResult>=0) {
415     if (!gSystem->AccessPathName(fName)) {
416       // these are the other events, use the target file and temporary files to merge
417       // with TChain
418       chain.Add(fName);
419
420       if (eventno>=0) {
421         TFile file(fName);
422         if (!file.IsZombie()) {
423           TTree* pSrcTree;
424           file.GetObject("esdTree", pSrcTree);
425           if (pSrcTree) {
426             nofCurrentEvents=pSrcTree->GetEntries();
427           }
428           file.Close();
429         }
430       }
431     }
432   }
433
434   // synchronize and add empty events
435   if (nofCurrentEvents<eventno) {
436     iResult=1; // indicate files to merge
437     TTree* pTgtTree=new TTree("esdTree", "Tree with HLT ESD objects");
438     if (pTgtTree) {
439       pTgtTree->SetDirectory(0);
440       AliESDEvent* pTmpESD=new AliESDEvent;
441       if (pTmpESD) {
442         TString tmpfilename;
443         FILE* pTmpFile=gSystem->TempFileName(tmpfilename);
444         if (pTmpFile) {
445           fclose(pTmpFile);
446           pTmpFile=NULL;
447           cleanup.Add(new TObjString(tmpfilename));
448           TFile emptyevents(tmpfilename, "RECREATE");
449           if (!emptyevents.IsZombie()) {
450             pTmpESD->CreateStdContent();
451             pTmpESD->WriteToTree(pTgtTree);
452             HLTDebug("adding %d empty events to file %s", eventno-nofCurrentEvents, fName.Data());
453             for (int i=nofCurrentEvents; i<eventno; i++) {
454               pTgtTree->Fill();
455             }
456             pTgtTree->GetUserInfo()->Add(pTmpESD);
457             emptyevents.cd();
458             pTgtTree->Write();
459             emptyevents.Close();
460             chain.Add(tmpfilename);
461             pTgtTree->GetUserInfo()->Clear();
462           }
463         }
464         delete pTmpESD;
465       } else {
466         iResult=-ENOMEM;
467       }
468       delete pTgtTree;
469     } else {
470       iResult=-ENOMEM;
471     }
472   }
473
474   if (iResult>=0 && pSrcESD) {
475     // add the new event to the chain
476     iResult=1; // indicate files to merge
477     TString tmpfilename=WriteTempFile(pSrcESD);
478     if (!tmpfilename.IsNull()) {
479       chain.Add(tmpfilename);
480       cleanup.Add(new TObjString(tmpfilename));
481     }
482   }
483
484   if (iResult>0) {
485     // build temporary file name for chain output
486     TString tgtName;
487     FILE* pTmpFile=gSystem->TempFileName(tgtName);
488     if (pTmpFile) {
489       fclose(pTmpFile);
490       pTmpFile=NULL;
491
492       // there have been problems with the memory consumption when using
493       // TChain::Merge
494       // but using a separate loop soemtimes crashes in AliESDEvent::ReadFromTree
495       // since this is for backward compatiblity only, we take the TChain::Merge
496       chain.Merge(tgtName);
497 //       TFile tgtFile(tgtName, "RECREATE");
498 //       TTree* pTgtTree=new TTree("esdTree", "Tree with HLT ESD objects");
499 //       AliESDEvent* pTgtESD=new AliESDEvent;
500 //       if (pTgtTree && pTgtESD) {
501 //      pTgtESD->ReadFromTree(&chain);
502 //      pTgtESD->WriteToTree(pTgtTree);
503
504 //      int nofEvents=chain.GetEntries();
505 //      for (int event=0; event<nofEvents; event++) {
506 //        chain.GetEntry(event);
507 //        pTgtTree->Fill();
508 //      }
509
510 //      pTgtTree->GetUserInfo()->Add(pTgtESD);
511 //      tgtFile.cd();
512 //      pTgtTree->Write();
513 //      pTgtTree->GetUserInfo()->Clear();
514 //       } else {
515 //      iResult=-ENOMEM;
516 //       }
517
518 //       if (pTgtTree) delete pTgtTree;
519 //       if (pTgtESD) delete pTgtESD;
520 //       tgtFile.Close();
521
522       // rename the merged file to the original file
523       TString shellcmd="mv ";
524       shellcmd+=tgtName + " " + fName;
525       if (gSystem->Exec(shellcmd)==0) {
526         HLTDebug("renaming %s to %s", tgtName.Data(), fName.Data());
527       } else {
528         HLTError("can not rename temporary file %s to %s", tgtName.Data(), fName.Data());
529       }
530     } else {
531       HLTError("can not get temporary file name from system");
532       iResult=-EBADF;
533     }
534   }
535
536   // delete temporary files
537   // the list objects are cleaned up by the TList destructor as the
538   // list is owner
539   TIter entry(&cleanup);
540   while (TObject* pObj=entry.Next()) {
541     if (dynamic_cast<TObjString*>(pObj)) {
542       TString shellcmd="rm -f ";
543       shellcmd+=(dynamic_cast<TObjString*>(pObj))->GetString();
544       gSystem->Exec(shellcmd);
545     }
546   }
547 #endif //HAVE_NOT_ESD_COPY
548
549   return iResult;
550 }
551
552 TString AliHLTEsdManagerImplementation::AliHLTEsdListEntry::WriteTempFile(AliESDEvent* pESD) const
553 {
554   // see header file for class documentation
555   int iResult=0;
556   TString tmpfilename;
557   FILE* pTmpFile=gSystem->TempFileName(tmpfilename);
558   if (pTmpFile) {
559     fclose(pTmpFile);
560     pTmpFile=NULL;
561
562     TFile file(tmpfilename, "RECREATE");
563     if (!file.IsZombie()) {
564       TTree* pTree=AliHLTEsdManagerImplementation::EmbedIntoTree(pESD);
565       if (pTree) {
566         file.cd();
567         if (pTree->Write()>0) {
568         } else {
569           HLTError("can not write esd tree to temporary file %s", tmpfilename.Data());
570         }
571
572         pTree->GetUserInfo()->Clear();
573         delete pTree;
574       } else {
575         iResult=-ENOMEM;
576       }
577       file.Close();
578     } else {
579       HLTError("can not open file %s", tmpfilename.Data());
580       iResult=-EBADF;
581     }
582   } else {
583     HLTError("can not get temporary file name from system");
584     iResult=-EBADF;
585   }
586
587   if (iResult<0) {
588     if (gSystem->AccessPathName(tmpfilename)==0) {
589       TString shellcmd="rm -f ";
590       shellcmd+=tmpfilename;
591       gSystem->Exec(shellcmd);
592     }
593     tmpfilename="";
594   }
595   return tmpfilename;
596 }
597
598 void AliHLTEsdManagerImplementation::AliHLTEsdListEntry::SetDirectory(const char* directory)
599 {
600   // see header file for class documentation
601   if (!directory) return;
602   if (!fName.IsNull()) {
603     HLTWarning("ESD entry already in writing mode (%s), ignoring directory", fName.Data());
604     return;
605   }
606   fDirectory=directory;
607 }
608
609 void AliHLTEsdManagerImplementation::AliHLTEsdListEntry::Delete()
610 {
611   // see header file for class documentation
612   if (fName.IsNull()) return;
613   if (gSystem->AccessPathName(fName)!=0) return;
614
615   TString shellcmd="rm -f ";
616   shellcmd+=fName;
617   gSystem->Exec(shellcmd);
618   fName="";
619 }
620
621 const char* AliHLTEsdManagerImplementation::AliHLTEsdListEntry::GetFileName() const
622 {
623   // see header file for class documentation
624   return fName.Data();
625 }
626
627 const char* AliHLTEsdManagerImplementation::AliHLTEsdListEntry::GetPrefix()
628 {
629   // see header file for class documentation
630   if (fPrefix.IsNull()) {
631     fPrefix.Insert(0, fDt.fOrigin, kAliHLTComponentDataTypefOriginSize);
632     fPrefix.Remove(TString::kTrailing, ' ');
633     fPrefix.ToUpper();
634     if (!fPrefix.Contains("HLT")) {
635       fPrefix.Insert(0, "HLT");
636     }
637   }
638   return fPrefix.Data();
639 }
640
641 int AliHLTEsdManagerImplementation::Merge(AliESDEvent* pTgt, AliESDEvent* pSrc) const
642 {
643   // see header file for class documentation
644   int iResult=0;
645   if (!pTgt || !pSrc) return -EINVAL;
646
647   TIter next(pSrc->GetList());
648   TObject* pSrcObject=NULL;
649   static int warningCount=0;
650   while ((pSrcObject=next())) {
651     if(!pSrcObject->InheritsFrom("TCollection")){
652       // simple objects
653       TString name=pSrcObject->GetName();
654       if(pSrcObject->InheritsFrom("AliHLTTriggerDecision")){
655         //pSrcObject->Print();
656         TObject* pTgtObject=pTgt->FindListObject(name);
657         if (pTgtObject) {
658           pSrcObject->Copy(*pTgtObject);
659         } else {
660           pTgt->AddObject(pSrcObject->Clone());
661         }
662       } else {
663         // TODO: implement the handling of other objects, some kind of mapping
664       }
665     } else if(pSrcObject->InheritsFrom("TClonesArray")){
666       TClonesArray* pTClA=dynamic_cast<TClonesArray*>(pSrcObject);
667       if (pTClA!=NULL && pTClA->GetEntriesFast()>0) {
668         TString name=pTClA->GetName();
669         TObject* pTgtObject=pTgt->GetList()->FindObject(name);
670         TClonesArray* pTgtArray=NULL;
671         if (pTgtObject!=NULL && pTgtObject->InheritsFrom("TClonesArray")){
672           pTgtArray=dynamic_cast<TClonesArray*>(pTgtObject);
673           if (pTgtArray) {
674             TString classType=pTClA->Class()->GetName();
675             if (classType.CompareTo(pTgtArray->Class()->GetName())==0) {
676               if (pTgtArray->GetEntries()==0) {
677                 pTgtArray->ExpandCreate(pTClA->GetEntries());
678                 for(int i=0; i<pTClA->GetEntriesFast(); ++i){
679                   (*pTClA)[i]->Copy(*((*pTgtArray)[i]));
680                 }
681               } else {
682                 if (warningCount++<10) {
683                   HLTWarning("TClonesArray \"%s\"  in target ESD %p is already filled with %d entries",
684                              name.Data(), pTgt, pTgtArray->GetEntries());
685                 }
686                 iResult=-EBUSY;
687               }
688             } else {
689               if (warningCount++<10) {
690                 HLTWarning("TClonesArray \"%s\" exists in target ESD %p, but describes incompatible class type %s instead of %s",
691                            name.Data(), pTgt, pTgtArray->GetClass()->GetName(), pTClA->GetClass()->GetName());
692               }
693               iResult=-EBUSY;
694             }
695           } else {
696             if (warningCount++<10) {
697               HLTError("internal error: dynamic cast failed for object %s %p", pTgtObject->GetName(), pTgtObject);
698             }
699             iResult=-EBUSY;
700           }
701         } else if (pTgtObject) {
702           if (warningCount++<10) {
703             HLTWarning("object \"%s\" does already exist in target ESD %p, but is %s rather than TClonesArray"
704                        " skipping data",
705                        name.Data(), pTgt, pTgtObject->Class()->GetName());
706           }
707           iResult=-EBUSY;
708         } else {
709           if (warningCount++<10) {
710             HLTWarning("object \"%s\" does not exist in target ESD %p, data can not be copied because it will be lost when filling the tree",
711                        name.Data(), pTgt);
712           }
713           iResult=-ENOENT;
714         }
715       }
716     }
717   }
718   return iResult;
719 }