]>
Commit | Line | Data |
---|---|---|
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 | } |