Renaming AliHLTReconstructorBase to AliHLTPluginBase to reflect the
[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 #ifndef HAVE_NOT_ESD_COPY
253   if (fName.IsNull()) {
254     // this is the first event, create the file name
255     fName="";
256     if (!fDirectory.IsNull()) {
257       fName+=fDirectory; fName+="/";
258     }
259     fName+="Ali"; fName+=GetPrefix();
260     if (fDt!=kAliHLTDataTypeESDObject &&
261         fDt!=kAliHLTDataTypeESDTree) {
262
263       HLTWarning("non-standard ESD type %s", AliHLTComponent::DataType2Text(fDt).c_str());
264       TString id;
265       id.Insert(0, fDt.fID, kAliHLTComponentDataTypefIDsize);
266       id.Remove(TString::kTrailing, ' ');
267       id.ToUpper();
268       fName+="_"; fName+=id; fName+=".root";
269     } else {
270       fName+="ESDs.root";
271     }
272
273     fpFile=new TFile(fName, "RECREATE");
274     fpTree=new TTree("esdTree", "Tree with HLT ESD objects");
275     fpTree->SetDirectory(0);
276     fpEsd=new AliESDEvent;
277     if (fpEsd) {
278       fpEsd->CreateStdContent();
279       *fpEsd=*pSrcESD;
280       if (fpTree) {
281         fpEsd->WriteToTree(fpTree);
282       }
283     }
284   }
285
286   if (fpFile && fpTree && fpEsd) {
287     // synchronize and add empty events
288     fpEsd->Reset();
289     int nofCurrentEvents=fpTree->GetEntries();
290     if (nofCurrentEvents<eventno) {
291       iResult=1; // indicate tree to be written
292       HLTDebug("adding %d empty events to file %s", eventno-nofCurrentEvents, fName.Data());
293       for (int i=nofCurrentEvents; i<eventno; i++) {
294         fpTree->Fill();
295       }
296     }
297
298     if (iResult>=0 && pSrcESD) {
299       int nofObjects=fpEsd->GetList()->GetEntries();
300       *fpEsd=*pSrcESD;
301       if (nofObjects!=fpEsd->GetList()->GetEntries()) {
302         // The source ESD contains object not present in the target ESD
303         // before. Those objects will not be written to the tree since
304         // the branch layout has been created earlier.
305         // Create new tree with the additional branches, copy the entries
306         // of the current tree into the new tree, and continue.
307         TTree* pNewTree=new TTree("esdTree", "Tree with HLT ESD objects");
308         pNewTree->SetDirectory(0);
309         AliESDEvent* readESD=new AliESDEvent;
310         readESD->CreateStdContent();
311         readESD->ReadFromTree(fpTree);
312         fpEsd->Reset();
313         fpEsd->WriteToTree(pNewTree);
314         HLTDebug("cloning tree with %d entries", fpTree->GetEntries());
315         for (int event=0; event<fpTree->GetEntries(); event++) {
316           fpTree->GetEntry(event);
317           *fpEsd=*readESD;
318           pNewTree->Fill();
319           fpEsd->Reset();
320         }
321         fpFile->Close();
322         delete fpFile;
323         delete readESD;
324         delete fpTree;
325         fpFile=new TFile(fName, "RECREATE");
326         fpTree=pNewTree;
327         *fpEsd=*pSrcESD;
328         HLTDebug("new ESD with %d objects", fpEsd->GetList()->GetEntries());
329       }
330       fpTree->Fill();
331       iResult=1; // indicate tree to be written
332     }
333
334     if (iResult>0) {
335       fpFile->cd();
336       fpTree->GetUserInfo()->Add(fpEsd);
337       fpTree->Write(fpTree->GetName(),TObject::kOverwrite);
338       fpTree->GetUserInfo()->Clear();
339     }
340   }
341 #else //HAVE_NOT_ESD_COPY
342   // this is the old workaround, necessary for older AliRoot versions
343   // version<=v4-12-Release
344
345   // we need to copy the ESD, I did not find an approptiate
346   // method, the workaround is to save the ESD in a temporary
347   // tree, read the content back into the ESD structure
348   // used for filling.
349   // Unfortunately the following code crashes at the second event.
350   // The expert on the ESD (Christian Klein Boesig) does not have
351   // a solution either. It seems to be a problem in ROOT.
352   //  TTree* dummy=new TTree("dummy","dummy");
353   //  dummy->SetDirectory(0);
354   //  pESD->WriteToTree(dummy);
355   //  dummy->Fill();
356   //  dummy->GetUserInfo()->Add(pESD);
357   //  fpEsd->ReadFromTree(dummy);
358   //  dummy->GetEvent(0);
359   //  fpEsd->WriteToTree(fpTree);
360   //  fpTree->Fill();
361   //  dummy->GetUserInfo()->Clear();
362   //  delete dummy;
363   //
364   // The only way is via TChain, which is working on files only at the
365   // time of writing.
366   // We use temporary files for the new event to be copied into the
367   // existing tree.
368   //
369   if (fName.IsNull()) {
370     // this is the first event, create the file on disk and write ESD
371     TString origin;
372     origin.Insert(0, fDt.fOrigin, kAliHLTComponentDataTypefOriginSize);
373     origin.Remove(TString::kTrailing, ' ');
374     origin.ToUpper();
375     fName="";
376     if (!fDirectory.IsNull()) {
377       fName+=fDirectory; fName+="/";
378     }
379     fName+="AliHLT"; fName+=origin;
380     if (fDt!=kAliHLTDataTypeESDObject &&
381         fDt!=kAliHLTDataTypeESDTree) {
382
383       HLTWarning("non-standard ESD type %s", AliHLTComponent::DataType2Text(fDt).c_str());
384       TString id;
385       id.Insert(0, fDt.fID, kAliHLTComponentDataTypefIDsize);
386       id.Remove(TString::kTrailing, ' ');
387       id.ToUpper();
388       fName+="_"; fName+=id; fName+=".root";
389     } else {
390       fName+="ESDs.root";
391     }
392
393     if (!gSystem->AccessPathName(fName)) {
394       // file exists, delete
395       TString shellcmd="rm -f ";
396       shellcmd+=fName;
397       gSystem->Exec(shellcmd);
398     }
399   }
400
401   TChain chain("esdTree");
402   TList cleanup;
403   cleanup.SetOwner();
404
405   int nofCurrentEvents=0;
406   if (iResult>=0) {
407     if (!gSystem->AccessPathName(fName)) {
408       // these are the other events, use the target file and temporary files to merge
409       // with TChain
410       chain.Add(fName);
411
412       if (eventno>=0) {
413         TFile file(fName);
414         if (!file.IsZombie()) {
415           TTree* pSrcTree;
416           file.GetObject("esdTree", pSrcTree);
417           if (pSrcTree) {
418             nofCurrentEvents=pSrcTree->GetEntries();
419           }
420           file.Close();
421         }
422       }
423     }
424   }
425
426   // synchronize and add empty events
427   if (nofCurrentEvents<eventno) {
428     iResult=1; // indicate files to merge
429     TTree* pTgtTree=new TTree("esdTree", "Tree with HLT ESD objects");
430     if (pTgtTree) {
431       pTgtTree->SetDirectory(0);
432       AliESDEvent* pTmpESD=new AliESDEvent;
433       if (pTmpESD) {
434         TString tmpfilename;
435         FILE* pTmpFile=gSystem->TempFileName(tmpfilename);
436         if (pTmpFile) {
437           fclose(pTmpFile);
438           pTmpFile=NULL;
439           cleanup.Add(new TObjString(tmpfilename));
440           TFile emptyevents(tmpfilename, "RECREATE");
441           if (!emptyevents.IsZombie()) {
442             pTmpESD->CreateStdContent();
443             pTmpESD->WriteToTree(pTgtTree);
444             HLTDebug("adding %d empty events to file %s", eventno-nofCurrentEvents, fName.Data());
445             for (int i=nofCurrentEvents; i<eventno; i++) {
446               pTgtTree->Fill();
447             }
448             pTgtTree->GetUserInfo()->Add(pTmpESD);
449             emptyevents.cd();
450             pTgtTree->Write();
451             emptyevents.Close();
452             chain.Add(tmpfilename);
453             pTgtTree->GetUserInfo()->Clear();
454           }
455         }
456         delete pTmpESD;
457       } else {
458         iResult=-ENOMEM;
459       }
460       delete pTgtTree;
461     } else {
462       iResult=-ENOMEM;
463     }
464   }
465
466   if (iResult>=0 && pSrcESD) {
467     // add the new event to the chain
468     iResult=1; // indicate files to merge
469     TString tmpfilename=WriteTempFile(pSrcESD);
470     if (!tmpfilename.IsNull()) {
471       chain.Add(tmpfilename);
472       cleanup.Add(new TObjString(tmpfilename));
473     }
474   }
475
476   if (iResult>0) {
477     // build temporary file name for chain output
478     TString tgtName;
479     FILE* pTmpFile=gSystem->TempFileName(tgtName);
480     if (pTmpFile) {
481       fclose(pTmpFile);
482       pTmpFile=NULL;
483
484       // there have been problems with the memory consumption when using
485       // TChain::Merge
486       // but using a separate loop soemtimes crashes in AliESDEvent::ReadFromTree
487       // since this is for backward compatiblity only, we take the TChain::Merge
488       chain.Merge(tgtName);
489 //       TFile tgtFile(tgtName, "RECREATE");
490 //       TTree* pTgtTree=new TTree("esdTree", "Tree with HLT ESD objects");
491 //       AliESDEvent* pTgtESD=new AliESDEvent;
492 //       if (pTgtTree && pTgtESD) {
493 //      pTgtESD->ReadFromTree(&chain);
494 //      pTgtESD->WriteToTree(pTgtTree);
495
496 //      int nofEvents=chain.GetEntries();
497 //      for (int event=0; event<nofEvents; event++) {
498 //        chain.GetEntry(event);
499 //        pTgtTree->Fill();
500 //      }
501
502 //      pTgtTree->GetUserInfo()->Add(pTgtESD);
503 //      tgtFile.cd();
504 //      pTgtTree->Write();
505 //      pTgtTree->GetUserInfo()->Clear();
506 //       } else {
507 //      iResult=-ENOMEM;
508 //       }
509
510 //       if (pTgtTree) delete pTgtTree;
511 //       if (pTgtESD) delete pTgtESD;
512 //       tgtFile.Close();
513
514       // rename the merged file to the original file
515       TString shellcmd="mv ";
516       shellcmd+=tgtName + " " + fName;
517       if (gSystem->Exec(shellcmd)==0) {
518         HLTDebug("renaming %s to %s", tgtName.Data(), fName.Data());
519       } else {
520         HLTError("can not rename temporary file %s to %s", tgtName.Data(), fName.Data());
521       }
522     } else {
523       HLTError("can not get temporary file name from system");
524       iResult=-EBADF;
525     }
526   }
527
528   // delete temporary files
529   // the list objects are cleaned up by the TList destructor as the
530   // list is owner
531   TIter entry(&cleanup);
532   while (TObject* pObj=entry.Next()) {
533     if (dynamic_cast<TObjString*>(pObj)) {
534       TString shellcmd="rm -f ";
535       shellcmd+=(dynamic_cast<TObjString*>(pObj))->GetString();
536       gSystem->Exec(shellcmd);
537     }
538   }
539 #endif //HAVE_NOT_ESD_COPY
540
541   return iResult;
542 }
543
544 TString AliHLTEsdManagerImplementation::AliHLTEsdListEntry::WriteTempFile(AliESDEvent* pESD) const
545 {
546   // see header file for class documentation
547   int iResult=0;
548   TString tmpfilename;
549   FILE* pTmpFile=gSystem->TempFileName(tmpfilename);
550   if (pTmpFile) {
551     fclose(pTmpFile);
552     pTmpFile=NULL;
553
554     TFile file(tmpfilename, "RECREATE");
555     if (!file.IsZombie()) {
556       TTree* pTree=AliHLTEsdManagerImplementation::EmbedIntoTree(pESD);
557       if (pTree) {
558         file.cd();
559         if (pTree->Write()>0) {
560         } else {
561           HLTError("can not write esd tree to temporary file %s", tmpfilename.Data());
562         }
563
564         pTree->GetUserInfo()->Clear();
565         delete pTree;
566       } else {
567         iResult=-ENOMEM;
568       }
569       file.Close();
570     } else {
571       HLTError("can not open file %s", tmpfilename.Data());
572       iResult=-EBADF;
573     }
574   } else {
575     HLTError("can not get temporary file name from system");
576     iResult=-EBADF;
577   }
578
579   if (iResult<0) {
580     if (gSystem->AccessPathName(tmpfilename)==0) {
581       TString shellcmd="rm -f ";
582       shellcmd+=tmpfilename;
583       gSystem->Exec(shellcmd);
584     }
585     tmpfilename="";
586   }
587   return tmpfilename;
588 }
589
590 void AliHLTEsdManagerImplementation::AliHLTEsdListEntry::SetDirectory(const char* directory)
591 {
592   // see header file for class documentation
593   if (!directory) return;
594   if (!fName.IsNull()) {
595     HLTWarning("ESD entry already in writing mode (%s), ignoring directory", fName.Data());
596     return;
597   }
598   fDirectory=directory;
599 }
600
601 void AliHLTEsdManagerImplementation::AliHLTEsdListEntry::Delete()
602 {
603   // see header file for class documentation
604   if (fName.IsNull()) return;
605   if (gSystem->AccessPathName(fName)!=0) return;
606
607   TString shellcmd="rm -f ";
608   shellcmd+=fName;
609   gSystem->Exec(shellcmd);
610   fName="";
611 }
612
613 const char* AliHLTEsdManagerImplementation::AliHLTEsdListEntry::GetFileName() const
614 {
615   // see header file for class documentation
616   return fName.Data();
617 }
618
619 const char* AliHLTEsdManagerImplementation::AliHLTEsdListEntry::GetPrefix()
620 {
621   // see header file for class documentation
622   if (fPrefix.IsNull()) {
623     fPrefix.Insert(0, fDt.fOrigin, kAliHLTComponentDataTypefOriginSize);
624     fPrefix.Remove(TString::kTrailing, ' ');
625     fPrefix.ToUpper();
626     if (!fPrefix.Contains("HLT")) {
627       fPrefix.Insert(0, "HLT");
628     }
629   }
630   return fPrefix.Data();
631 }
632
633 int AliHLTEsdManagerImplementation::Merge(AliESDEvent* pTgt, AliESDEvent* pSrc) const
634 {
635   // see header file for class documentation
636   int iResult=0;
637   if (!pTgt || !pSrc) return -EINVAL;
638
639   TIter next(pSrc->GetList());
640   TObject* pSrcObject=NULL;
641   TString name;
642   static int warningCount=0;
643   while ((pSrcObject=next())) {
644     if(!pSrcObject->InheritsFrom("TCollection")){
645       // simple objects
646     } else if(pSrcObject->InheritsFrom("TClonesArray")){
647       TClonesArray* pTClA=dynamic_cast<TClonesArray*>(pSrcObject);
648       if (pTClA!=NULL && pTClA->GetEntriesFast()>0) {
649         name=pTClA->GetName();
650         TObject* pTgtObject=pTgt->GetList()->FindObject(name);
651         TClonesArray* pTgtArray=NULL;
652         if (pTgtObject!=NULL && pTgtObject->InheritsFrom("TClonesArray")){
653           pTgtArray=dynamic_cast<TClonesArray*>(pTgtObject);
654           if (pTgtArray) {
655             TString classType=pTClA->Class()->GetName();
656             if (classType.CompareTo(pTgtArray->Class()->GetName())==0) {
657               if (pTgtArray->GetEntries()==0) {
658                 pTgtArray->ExpandCreate(pTClA->GetEntries());
659                 for(int i=0; i<pTClA->GetEntriesFast(); ++i){
660                   (*pTClA)[i]->Copy(*((*pTgtArray)[i]));
661                 }
662               } else {
663                 if (warningCount++<10) {
664                   HLTWarning("TClonesArray \"%s\"  in target ESD %p is already filled with %d entries",
665                              name.Data(), pTgt, pTgtArray->GetEntries());
666                 }
667                 iResult=-EBUSY;
668               }
669             } else {
670               if (warningCount++<10) {
671                 HLTWarning("TClonesArray \"%s\" exists in target ESD %p, but describes incompatible class type %s instead of %s",
672                            name.Data(), pTgt, pTgtArray->GetClass()->GetName(), pTClA->GetClass()->GetName());
673               }
674               iResult=-EBUSY;
675             }
676           } else {
677             if (warningCount++<10) {
678               HLTError("internal error: dynamic cast failed for object %s %p", pTgtObject->GetName(), pTgtObject);
679             }
680             iResult=-EBUSY;
681           }
682         } else if (pTgtObject) {
683           if (warningCount++<10) {
684             HLTWarning("object \"%s\" does already exist in target ESD %p, but is %s rather than TClonesArray"
685                        " skipping data",
686                        name.Data(), pTgt, pTgtObject->Class()->GetName());
687           }
688           iResult=-EBUSY;
689         } else {
690           if (warningCount++<10) {
691             HLTWarning("object \"%s\" does not exist in target ESD %p, data can not be copied because it will be lost when filling the tree",
692                        name.Data(), pTgt);
693           }
694           iResult=-ENOENT;
695         }
696       }
697     }
698   }
699   return iResult;
700 }