removing deprecated code for backward compatibility, all supported AliRoot version...
[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 ///         This is an implementation of the abstract interface AliHLTEsdManager
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 "AliESDFMD.h"
32 #include "AliESDVZERO.h"
33 #include "AliESDTZERO.h"
34 #include "AliESDCaloCells.h"
35 #include "AliMultiplicity.h"
36 #include "AliESDACORDE.h"
37 #include "TFile.h"
38 #include "TTree.h"
39 #include "TClass.h"
40 #include "TObject.h"
41 #include "TList.h"
42
43 const float kAliESDVZERODefaultTime = -1024.;
44 const float kAliESDVZERODefaultTimeGlitch = 1e-6;
45 const float kAliESDZDCDefaultEMEnergy = 0.;
46 const float kAliESDZDCDefaultEMEnergyGlitch = 1e-6;
47
48 /** ROOT macro for the implementation of ROOT specific class methods */
49 ClassImp(AliHLTEsdManagerImplementation)
50
51 AliHLTEsdManagerImplementation::AliHLTEsdManagerImplementation()
52   :
53   fESDs()
54   , fDirectory()
55   , fWriteLocal(false)
56 {
57   // see header file for class documentation
58   // or
59   // refer to README to build package
60   // or
61   // visit http://web.ift.uib.no/~kjeks/doc/alice-hlt
62
63   CheckClassConditions();
64 }
65
66 AliHLTEsdManagerImplementation::~AliHLTEsdManagerImplementation()
67 {
68   // see header file for class documentation
69   for (unsigned int i=0; i<fESDs.size(); i++) {
70     if (fESDs[i]) {
71       delete fESDs[i];
72     }
73     fESDs[i]=NULL;
74   }
75 }
76
77 int AliHLTEsdManagerImplementation::SetOption(const char* option)
78 {
79   // see header file for class documentation
80   int iResult=0;
81   TString strOptions=option;
82   TObjArray* pTokens=strOptions.Tokenize(" ");
83   if (pTokens) {
84     if (pTokens->GetEntriesFast()>0) {
85       for (int n=0; n<pTokens->GetEntriesFast(); n++) {
86         TString data=((TObjString*)pTokens->At(n))->GetString();
87         if (data.IsNull()) continue;
88
89         if (data.CompareTo("-writelocal")==0) {
90           fWriteLocal=true;
91         } else if (data.Contains("-directory=")) {
92           data.ReplaceAll("-directory=", "");
93           SetDirectory(data.Data());
94         } else {
95           HLTError("unknown argument %s", data.Data());
96           iResult=-EINVAL;
97           break;
98         }
99       }
100     }
101     delete pTokens;
102   }
103   return iResult;
104 }
105
106 AliHLTEsdManagerImplementation::AliHLTEsdListEntry* AliHLTEsdManagerImplementation::Find(AliHLTComponentDataType dt) const
107 {
108   // see header file for class documentation
109   AliHLTEsdListEntry* pEntry=NULL;
110   for (unsigned int i=0; i<fESDs.size(); i++) {
111     if (fESDs[i] && *(fESDs[i])==dt) {
112       pEntry=const_cast<AliHLTEsdListEntry*>(fESDs[i]);
113     }
114   }
115   return pEntry;
116 }
117
118 int AliHLTEsdManagerImplementation::WriteESD(const AliHLTUInt8_t* pBuffer, AliHLTUInt32_t size,
119                                AliHLTComponentDataType dt, AliESDEvent* tgtesd, int eventno)
120 {
121   // see header file for class documentation
122   if (!pBuffer && size<=0) return -EINVAL;
123   int iResult=0;
124   AliHLTUInt32_t firstWord=*((AliHLTUInt32_t*)pBuffer);
125   if (firstWord==size-sizeof(AliHLTUInt32_t)) {
126     HLTDebug("create object from block %s size %d", AliHLTComponent::DataType2Text(dt).c_str(), size);
127     AliHLTMessage msg(const_cast<AliHLTUInt8_t*>(pBuffer), size);
128     TClass* objclass=msg.GetClass();
129     TObject* pObj=msg.ReadObject(objclass);
130     if (pObj && objclass) {
131       HLTDebug("object %p type %s created", pObj, objclass->GetName());
132       AliESDEvent* pESD=dynamic_cast<AliESDEvent*>(pObj);
133       TTree* pTree=NULL;
134       if (!pESD) {
135         pTree=dynamic_cast<TTree*>(pObj);
136         if (pTree) {
137           pESD=new AliESDEvent;
138           pESD->CreateStdContent();
139           if (pTree->GetEntries()>0) {
140             if (pTree->GetEntries()>1) {
141               HLTWarning("only one entry allowed for ESD embedded into tree, data block %s contains tree with %d entires, taking first entry",
142                          AliHLTComponent::DataType2Text(dt).c_str(), pTree->GetEntries());
143             }
144             pESD->ReadFromTree(pTree);
145             pTree->GetEvent(0);
146           }
147         } else {
148           HLTWarning("tree of data block %s has no events, skipping data block", AliHLTComponent::DataType2Text(dt).c_str());
149         }
150       }
151       if (pESD) {
152         AliHLTEsdListEntry* entry=Find(dt);
153         if (!entry) {
154           if ((entry=new AliHLTEsdListEntry(dt))!=NULL) {
155             if (!fDirectory.IsNull()) {
156               entry->SetDirectory(fDirectory);
157             }
158             fESDs.push_back(entry);
159           }
160         }
161         if (entry) {
162           if (tgtesd) {
163             Merge(tgtesd, pESD);
164           }
165
166           // Matthias 2009-06-06: writing of individual ESD files for the different origins was a
167           // first attempt when functionality was missing in the AliRoot framework and remained as
168           // debugging feature. ESD merging is now implemented and data written to the hltEsd, so
169           // the feature is now disabled by default because it causes increasing memory consumption.
170           // Presumably not because of a memory leak but the way the internal TTree is used and kept
171           // in memory.
172           // Writing of local files can be optionally switched on as e.g. by the EsdCollector component.
173           if (fWriteLocal) entry->WriteESD(pESD, eventno);
174         } else {
175           HLTError("internal mismatch, can not create list entry");
176           iResult=-ENOMEM;
177         }
178       } else {
179         HLTWarning("data block %s is not of class type AliESDEvent, ignoring ...", AliHLTComponent::DataType2Text(dt).c_str());
180       }
181       if (pTree) {
182         // ESD has been created and must be cleaned up
183         pESD->Reset();
184         delete pESD;
185         pESD=NULL;
186       }
187       delete pObj;
188       pObj=NULL;
189     } else {
190     }
191   }
192   return iResult;
193 }
194
195 int AliHLTEsdManagerImplementation::PadESDs(int eventno)
196 {
197   // see header file for class documentation
198   int iResult=0;
199   for (unsigned int i=0; i<fESDs.size(); i++) {
200     if (fESDs[i]) {
201       int res=fESDs[i]->WriteESD(NULL, eventno);
202       if (res<0 && iResult>=0) iResult=res;
203     }
204   }
205   return iResult;
206 }
207
208 void AliHLTEsdManagerImplementation::SetDirectory(const char* directory)
209 {
210   // see header file for class documentation
211   if (!directory) return;
212   fDirectory=directory;
213   for (unsigned int i=0; i<fESDs.size(); i++) {
214     if (fESDs[i]) {
215       fESDs[i]->SetDirectory(directory);
216     }
217   }
218 }
219
220 TString AliHLTEsdManagerImplementation::GetFileNames(AliHLTComponentDataType dt) const
221 {
222   TString result;
223   for (unsigned int i=0; i<fESDs.size(); i++) {
224     if (fESDs[i] && *(fESDs[i])==dt) {
225       if (!result.IsNull()) result+=" ";
226       result+=fESDs[i]->GetFileName();
227     }
228   }
229   return result;
230 }
231
232 TTree* AliHLTEsdManagerImplementation::EmbedIntoTree(AliESDEvent* pESD, const char* name, const char* title)
233 {
234   // see header file for class documentation
235   int iResult=0;
236   TTree* pTree=new TTree(name, title);
237   if (pTree) {
238     pESD->WriteToTree(pTree);
239     pTree->Fill();
240     pTree->GetUserInfo()->Add(pESD);
241   } else {
242     iResult=-ENOMEM;
243   }
244
245   if (iResult<0) {
246     pTree->GetUserInfo()->Clear();
247     delete pTree;
248   }
249
250   return pTree;
251 }
252
253 AliHLTEsdManagerImplementation::AliHLTEsdListEntry::AliHLTEsdListEntry(AliHLTComponentDataType dt)
254   :
255   fName(),
256   fDirectory(),
257   fDt(dt),
258   fpFile(NULL),
259   fpTree(NULL),
260   fpEsd(NULL),
261   fPrefix()
262 {
263   // see header file for class documentation
264 }
265
266 AliHLTEsdManagerImplementation::AliHLTEsdListEntry::~AliHLTEsdListEntry()
267 {
268   // see header file for class documentation
269   if (fpEsd) delete fpEsd;
270   fpEsd=NULL;
271
272   if (fpTree) delete fpTree;
273   fpTree=NULL;
274
275   if (fpFile) {
276     fpFile->Close();
277     delete fpFile;
278   }
279   fpFile=NULL;
280 }
281
282 bool AliHLTEsdManagerImplementation::AliHLTEsdListEntry::operator==(AliHLTComponentDataType dt) const
283 {
284   // see header file for class documentation
285   return fDt==dt;
286 }
287
288 int AliHLTEsdManagerImplementation::AliHLTEsdListEntry::WriteESD(AliESDEvent* pSrcESD, int eventno)
289 {
290   // see header file for class documentation
291   int iResult=0;
292
293   if (fName.IsNull()) {
294     // this is the first event, create the file name
295     fName="";
296     if (!fDirectory.IsNull()) {
297       fName+=fDirectory; fName+="/";
298     }
299     fName+="Ali"; fName+=GetPrefix();
300     if (fDt!=kAliHLTDataTypeESDObject &&
301         fDt!=kAliHLTDataTypeESDTree) {
302
303       HLTWarning("non-standard ESD type %s", AliHLTComponent::DataType2Text(fDt).c_str());
304       TString id;
305       id.Insert(0, fDt.fID, kAliHLTComponentDataTypefIDsize);
306       id.Remove(TString::kTrailing, ' ');
307       id.ToUpper();
308       fName+="_"; fName+=id; fName+=".root";
309     } else {
310       fName+="ESDs.root";
311     }
312
313     fpFile=new TFile(fName, "RECREATE");
314     fpTree=new TTree("esdTree", "Tree with HLT ESD objects");
315     fpTree->SetDirectory(0);
316     fpEsd=new AliESDEvent;
317     if (fpEsd) {
318       fpEsd->CreateStdContent();
319       *fpEsd=*pSrcESD;
320       if (fpTree) {
321         fpEsd->WriteToTree(fpTree);
322       }
323     }
324   }
325
326   if (fpFile && fpTree && fpEsd) {
327     // synchronize and add empty events
328     fpEsd->Reset();
329     int nofCurrentEvents=fpTree->GetEntries();
330     if (nofCurrentEvents<eventno) {
331       iResult=1; // indicate tree to be written
332       HLTDebug("adding %d empty events to file %s", eventno-nofCurrentEvents, fName.Data());
333       for (int i=nofCurrentEvents; i<eventno; i++) {
334         fpTree->Fill();
335       }
336     }
337
338     if (iResult>=0 && pSrcESD) {
339       int nofObjects=fpEsd->GetList()->GetEntries();
340       *fpEsd=*pSrcESD;
341       if (nofObjects!=fpEsd->GetList()->GetEntries()) {
342         // The source ESD contains object not present in the target ESD
343         // before. Those objects will not be written to the tree since
344         // the branch layout has been created earlier.
345         // Create new tree with the additional branches, copy the entries
346         // of the current tree into the new tree, and continue.
347         TTree* pNewTree=new TTree("esdTree", "Tree with HLT ESD objects");
348         pNewTree->SetDirectory(0);
349         AliESDEvent* readESD=new AliESDEvent;
350         readESD->CreateStdContent();
351         readESD->ReadFromTree(fpTree);
352         fpEsd->Reset();
353         fpEsd->WriteToTree(pNewTree);
354         HLTDebug("cloning tree with %d entries", fpTree->GetEntries());
355         for (int event=0; event<fpTree->GetEntries(); event++) {
356           fpTree->GetEntry(event);
357           *fpEsd=*readESD;
358           pNewTree->Fill();
359           fpEsd->Reset();
360         }
361         fpFile->Close();
362         delete fpFile;
363         delete readESD;
364         delete fpTree;
365         fpFile=new TFile(fName, "RECREATE");
366         fpTree=pNewTree;
367         *fpEsd=*pSrcESD;
368         HLTDebug("new ESD with %d objects", fpEsd->GetList()->GetEntries());
369       }
370       fpTree->Fill();
371       iResult=1; // indicate tree to be written
372     }
373
374     if (iResult>0) {
375       fpFile->cd();
376       fpTree->GetUserInfo()->Add(fpEsd);
377       fpTree->Write(fpTree->GetName(),TObject::kOverwrite);
378       fpTree->GetUserInfo()->Clear();
379     }
380   }
381
382   return iResult;
383 }
384
385 void AliHLTEsdManagerImplementation::AliHLTEsdListEntry::SetDirectory(const char* directory)
386 {
387   // see header file for class documentation
388   if (!directory) return;
389   if (!fName.IsNull()) {
390     HLTWarning("ESD entry already in writing mode (%s), ignoring directory", fName.Data());
391     return;
392   }
393   fDirectory=directory;
394 }
395
396 const char* AliHLTEsdManagerImplementation::AliHLTEsdListEntry::GetFileName() const
397 {
398   // see header file for class documentation
399   return fName.Data();
400 }
401
402 const char* AliHLTEsdManagerImplementation::AliHLTEsdListEntry::GetPrefix()
403 {
404   // see header file for class documentation
405   if (fPrefix.IsNull()) {
406     fPrefix.Insert(0, fDt.fOrigin, kAliHLTComponentDataTypefOriginSize);
407     fPrefix.Remove(TString::kTrailing, ' ');
408     fPrefix.ToUpper();
409     if (!fPrefix.Contains("HLT")) {
410       fPrefix.Insert(0, "HLT");
411     }
412   }
413   return fPrefix.Data();
414 }
415
416 int AliHLTEsdManagerImplementation::Merge(AliESDEvent* pTgt, AliESDEvent* pSrc) const
417 {
418   // see header file for class documentation
419   int iResult=0;
420   if (!pTgt || !pSrc) return -EINVAL;
421
422   TIter next(pSrc->GetList());
423   TObject* pSrcObject=NULL;
424   static int warningCount=0;
425   while ((pSrcObject=next())) {
426     if(!pSrcObject->InheritsFrom("TCollection")){
427       // simple objects
428       // for every type of object we have to find out whether it is empty or not
429       // objects are only copied if non-empty, otherwise valid content would be
430       // overridden by empty objects during the merging
431       bool copy=false;
432       TString name=pSrcObject->GetName();
433       if(pSrcObject->InheritsFrom("AliHLTTriggerDecision")){
434         copy=true;
435       } else if (pSrcObject->IsA()==AliESDRun::Class()) {
436         AliESDRun* pESDRun=dynamic_cast<AliESDRun*>(pSrcObject);
437         // zero might be a valid run no in simulation, so we check also whether the CTP trigger classes are set
438         copy=(pESDRun && (pESDRun->GetRunNumber()>0 || !pESDRun->GetActiveTriggerClasses().IsNull()));
439       } else if (pSrcObject->IsA()==AliESDHeader::Class()) {
440         AliESDHeader* pESDHeader=dynamic_cast<AliESDHeader*>(pSrcObject);
441         copy=(pESDHeader && pESDHeader->GetTriggerMask()!=0);
442       } else if (pSrcObject->IsA()==AliESDVertex::Class()) {
443         AliESDVertex* pESDVertex=dynamic_cast<AliESDVertex*>(pSrcObject);
444         copy=(pESDVertex && pESDVertex->GetNContributors()>0);
445       } else if (pSrcObject->IsA()==AliESDTZERO::Class()) {
446         AliESDTZERO* pESDTZero=dynamic_cast<AliESDTZERO*>(pSrcObject);
447         copy=(pESDTZero && (pESDTZero->GetT0zVertex()!=0.0 || pESDTZero->GetT0()!=0.0));
448       } else if (pSrcObject->IsA()==AliESDVZERO::Class()) {
449         AliESDVZERO* pESDVZERO=dynamic_cast<AliESDVZERO*>(pSrcObject);
450         // object contains data if one of the times exceeds the default value
451         copy=pESDVZERO &&
452           (pESDVZERO->GetV0ATime() > kAliESDVZERODefaultTime+kAliESDVZERODefaultTimeGlitch ||
453            pESDVZERO->GetV0CTime() > kAliESDVZERODefaultTime+kAliESDVZERODefaultTimeGlitch);
454       } else if (pSrcObject->IsA()==AliESDFMD::Class()) {
455         AliESDFMD* pESDFMD=dynamic_cast<AliESDFMD*>(pSrcObject);
456         copy=(pESDFMD && false); // have to find an easy valid condition
457       } else if (pSrcObject->IsA()==AliESDZDC::Class()) {
458         AliESDZDC* pESDZDC=dynamic_cast<AliESDZDC*>(pSrcObject);
459         // object contains data if the EM energies are different from the default value
460         copy=pESDZDC &&
461           (TMath::Abs(pESDZDC->GetZDCEMEnergy(0)-kAliESDZDCDefaultEMEnergy) > kAliESDZDCDefaultEMEnergyGlitch ||
462            TMath::Abs(pESDZDC->GetZDCEMEnergy(1)-kAliESDZDCDefaultEMEnergy) > kAliESDZDCDefaultEMEnergyGlitch);
463       } else if (pSrcObject->IsA()==AliMultiplicity::Class()) {
464         AliMultiplicity* pMultiplicity=dynamic_cast<AliMultiplicity*>(pSrcObject);
465         copy=(pMultiplicity && pMultiplicity->GetNumberOfTracklets()>0);
466       } else if (pSrcObject->IsA()==AliESDCaloTrigger::Class()) {
467         AliESDCaloTrigger* pESDCaloTrigger=dynamic_cast<AliESDCaloTrigger*>(pSrcObject);
468         copy=(pESDCaloTrigger && false); // have to find an easy valid condition
469       } else if (pSrcObject->IsA()==AliESDCaloCells::Class()) {
470         AliESDCaloCells* pESDCaloCells=dynamic_cast<AliESDCaloCells*>(pSrcObject);
471         copy=(pESDCaloCells && false); // have to find an easy valid condition
472       } else if (pSrcObject->IsA()==AliESDACORDE::Class()) {
473         AliESDACORDE* pESDACORDE=dynamic_cast<AliESDACORDE*>(pSrcObject);
474         copy=(pESDACORDE && false); // have to find an easy valid condition
475       } else if (!AliHLTESDEventHelper::IsStdContent(name)) {
476         // this is likely to be ok as long as it is not any object of the std content.
477         copy=true;
478       } else {
479         HLTError("no merging implemented for object %s, omitting", name.Data());
480       }
481       if (copy) {
482         //pSrcObject->Print();
483         TObject* pTgtObject=pTgt->FindListObject(name);
484         if (pTgtObject) {
485           pSrcObject->Copy(*pTgtObject);
486         } else {
487           pTgt->AddObject(pSrcObject->Clone());
488         }
489       }
490     } else if(pSrcObject->InheritsFrom("TClonesArray")){
491       TClonesArray* pTClA=dynamic_cast<TClonesArray*>(pSrcObject);
492       if (pTClA!=NULL && pTClA->GetEntriesFast()>0) {
493         TString name=pTClA->GetName();
494         TObject* pTgtObject=pTgt->GetList()->FindObject(name);
495         TClonesArray* pTgtArray=NULL;
496         if (pTgtObject!=NULL && pTgtObject->InheritsFrom("TClonesArray")){
497           pTgtArray=dynamic_cast<TClonesArray*>(pTgtObject);
498           if (pTgtArray) {
499             TString classType=pTClA->Class()->GetName();
500             if (classType.CompareTo(pTgtArray->Class()->GetName())==0) {
501               if (pTgtArray->GetEntries()==0) {
502                 pTgtArray->ExpandCreate(pTClA->GetEntries());
503                 for(int i=0; i<pTClA->GetEntriesFast(); ++i){
504                   (*pTClA)[i]->Copy(*((*pTgtArray)[i]));
505                 }
506               } else {
507                 if (warningCount++<10) {
508                   HLTWarning("TClonesArray \"%s\"  in target ESD %p is already filled with %d entries",
509                              name.Data(), pTgt, pTgtArray->GetEntries());
510                 }
511                 iResult=-EBUSY;
512               }
513             } else {
514               if (warningCount++<10) {
515                 HLTWarning("TClonesArray \"%s\" exists in target ESD %p, but describes incompatible class type %s instead of %s",
516                            name.Data(), pTgt, pTgtArray->GetClass()->GetName(), pTClA->GetClass()->GetName());
517               }
518               iResult=-EBUSY;
519             }
520           } else {
521             if (warningCount++<10) {
522               HLTError("internal error: dynamic cast failed for object %s %p", pTgtObject->GetName(), pTgtObject);
523             }
524             iResult=-EBUSY;
525           }
526         } else if (pTgtObject) {
527           if (warningCount++<10) {
528             HLTWarning("object \"%s\" does already exist in target ESD %p, but is %s rather than TClonesArray"
529                        " skipping data",
530                        name.Data(), pTgt, pTgtObject->Class()->GetName());
531           }
532           iResult=-EBUSY;
533         } else {
534           if (warningCount++<10) {
535             HLTWarning("object \"%s\" does not exist in target ESD %p, data can not be copied because it will be lost when filling the tree",
536                        name.Data(), pTgt);
537           }
538           iResult=-ENOENT;
539         }
540       }
541     }
542   }
543   return iResult;
544 }
545
546 bool AliHLTEsdManagerImplementation::AliHLTESDEventHelper::IsStdContent(const char* key)
547 {
548   // check if the key denotes a std object
549   TString needle=key;
550   for (int i=0; i<kESDListN; i++) {
551     if (needle.CompareTo(fgkESDListName[i])==0) return true;
552   }
553   return false;
554 }
555
556 int AliHLTEsdManagerImplementation::CheckClassConditions() const
557 {
558   // this is a helper method which checks if some thing in the default
559   // object initialization changes during the evolution of the source code
560   // which is necessary for checking the validity of data in the object
561
562   // check AliESDVZERO
563   AliESDVZERO vzerodummy;
564   if (TMath::Abs(vzerodummy.GetV0ATime()-kAliESDVZERODefaultTime) > kAliESDVZERODefaultTimeGlitch ||
565       TMath::Abs(vzerodummy.GetV0CTime()-kAliESDVZERODefaultTime) > kAliESDVZERODefaultTimeGlitch) {
566     HLTWarning("initialization of AliESDVZERO has changed, default time values now %f/%f, "
567                "revision of condition might be needed",
568                vzerodummy.GetV0ATime(), vzerodummy.GetV0CTime());
569   }
570
571   // check AliESDZDC
572   AliESDZDC zdcdummy;
573   if (TMath::Abs(zdcdummy.GetZDCEMEnergy(0)-kAliESDZDCDefaultEMEnergy) > kAliESDZDCDefaultEMEnergyGlitch ||
574       TMath::Abs(zdcdummy.GetZDCEMEnergy(1)-kAliESDZDCDefaultEMEnergy) > kAliESDZDCDefaultEMEnergyGlitch) {
575     HLTWarning("initialization of AliESDZDC has changed, default em energy values now %f/%f, "
576                "revision of condition might be needed",
577                zdcdummy.GetZDCEMEnergy(0), zdcdummy.GetZDCEMEnergy(1));
578   }
579
580   return 0;
581 }