3 //**************************************************************************
4 //* This file is property of and copyright by the ALICE HLT Project *
5 //* ALICE Experiment at CERN, All rights reserved. *
7 //* Primary Authors: Matthias Richter <Matthias.Richter@ift.uib.no> *
8 //* for The ALICE HLT Project. *
10 //* Permission to use, copy, modify and distribute this software and its *
11 //* documentation strictly for non-commercial purposes is hereby granted *
12 //* without fee, provided that the above copyright notice appears in all *
13 //* copies and that both the copyright notice and this permission notice *
14 //* appear in the supporting documentation. The authors make no claims *
15 //* about the suitability of this software for any purpose. It is *
16 //* provided "as is" without express or implied warranty. *
17 //**************************************************************************
19 // @file AliHLTOUTComponent.cxx
20 // @author Matthias Richter
22 // @brief The HLTOUT data sink component similar to HLTOUT nodes
23 // @note Used in the AliRoot environment only.
31 #include "AliHLTOUTComponent.h"
32 #include "AliHLTOUT.h"
33 #include "AliHLTHOMERLibManager.h"
34 #include "AliHLTHOMERWriter.h"
35 #include "AliHLTErrorGuard.h"
36 #include "AliDAQ.h" // equipment Ids
37 #include "AliRawDataHeader.h" // Common Data Header
38 #include <TDatime.h> // seed for TRandom
39 #include <TRandom.h> // random int generation for DDL no
45 /** ROOT macro for the implementation of ROOT specific class methods */
46 ClassImp(AliHLTOUTComponent)
48 AliHLTOUTComponent::AliHLTOUTComponent(EType type)
49 : AliHLTOfflineDataSink()
52 , fIdFirstDDL(7680) // 0x1e<<8
56 , fDigitFileName("HLT.Digits.root")
59 , fppDigitArrays(NULL)
64 // see header file for class documentation
66 // refer to README to build package
68 // visit http://web.ift.uib.no/~kjeks/doc/alice-hlt
70 // I guess DDL definitions should never change any more
71 assert(fNofDDLs==AliDAQ::NumberOfDdls("HLT"));
72 fNofDDLs=AliDAQ::NumberOfDdls("HLT");
74 /* AliDAQ::DdlIDOffset returns wrong offset for HLT links
75 assert(fIdFirstDDL==AliDAQ::DdlIDOffset("HLT"));
76 fIdFirstDDL=AliDAQ::DdlIDOffset("HLT");
79 if (fType!=kGlobal && fType!=kDigits && fType!=kRaw) {
80 ALIHLTERRORGUARD(1, "invalid component type %d", fType);
84 int AliHLTOUTComponent::fgOptions=kWriteRawFiles|kWriteDigits;
86 AliHLTOUTComponent::~AliHLTOUTComponent()
88 // see header file for class documentation
89 if (fpLibManager) delete fpLibManager;
93 const char* AliHLTOUTComponent::GetComponentID()
95 // see header file for class documentation
97 case kDigits: return "HLTOUTdigits";
98 case kRaw: return "HLTOUTraw";
106 void AliHLTOUTComponent::GetInputDataTypes( vector<AliHLTComponentDataType>& list)
108 // see header file for class documentation
110 list.push_back(kAliHLTAnyDataType);
113 AliHLTComponent* AliHLTOUTComponent::Spawn()
115 // see header file for class documentation
116 return new AliHLTOUTComponent(fType);
119 int AliHLTOUTComponent::DoInit( int argc, const char** argv )
121 // see header file for class documentation
126 fOptions|=kWriteDigits; fOptions&=~kWriteRawFiles;
127 HLTInfo("initializing HLTOUT component for digits generation");
130 fOptions|=kWriteRawFiles; fOptions&=~kWriteDigits;
131 HLTInfo("initializing HLTOUT component for raw data generation");
138 if ((iResult=ConfigureFromArgumentString(argc, argv))<0) return iResult;
140 // Make sure there is no library manager before we try and create a new one.
146 // Create a new library manager and allocate the appropriate number of
147 // HOMER writers for the HLTOUT component.
148 fpLibManager=new AliHLTHOMERLibManager;
151 for (writerNo=0; writerNo<fNofDDLs; writerNo++) {
152 AliHLTMonitoringWriter* pWriter=fpLibManager->OpenWriter();
154 HLTDebug("HOMER writer %p added", pWriter);
155 fWriters.push_back(pWriter);
157 HLTError("can not open HOMER writer");
169 int AliHLTOUTComponent::ScanConfigurationArgument(int argc, const char** argv)
171 // see header file for class documentation
172 if (argc<=0) return 0;
174 TString argument=argv[i];
178 // specify number of ddl links
179 if (argument.CompareTo("-links")==0) {
180 if (++i>=argc) return -EPROTO;
181 TString parameter(argv[i]);
182 parameter.Remove(TString::kLeading, ' '); // remove all blanks
183 if (parameter.IsDigit()) {
184 fNofDDLs=parameter.Atoi();
186 HLTError("wrong parameter for argument %s, number expected", argument.Data());
194 if (argument.CompareTo("-digitfile")==0) {
195 if (++i>=argc) return -EPROTO;
196 fDigitFileName=argv[i];
203 if (argument.Contains(key)) {
204 argument.ReplaceAll(key, "");
205 if (argument.IsNull()) {
206 fOptions|=kWriteRawFiles;
207 } else if (argument.CompareTo("=off")==0) {
208 fOptions&=~kWriteRawFiles;
209 } else if (argument.CompareTo("=on")==0) {
210 fOptions|=kWriteRawFiles;
212 HLTError("invalid parameter for argument %s: possible %s=off/%s=on", key, key, key);
221 if (argument.Contains(key)) {
222 argument.ReplaceAll(key, "");
223 if (argument.IsNull()) {
224 fOptions|=kWriteDigits;
225 } else if (argument.CompareTo("=off")==0) {
226 fOptions&=~kWriteDigits;
227 } else if (argument.CompareTo("=on")==0) {
228 fOptions|=kWriteDigits;
230 HLTError("invalid parameter for argument %s: possible %s=off/%s=on", key, key, key);
241 int AliHLTOUTComponent::DoDeinit()
243 // see header file for class documentation
247 AliHLTMonitoringWriterPVector::iterator element=fWriters.begin();
248 while (element!= fWriters.end()) {
250 // wanted to have a dynamic_cast<AliHLTHOMERWriter*> here, but this results into
251 // undefined symbol when loading the library
252 if (*element!=NULL) {
254 fpLibManager->DeleteWriter((AliHLTHOMERWriter*)(*element));
256 HLTError("writer instance is NULL");
258 element=fWriters.erase(element);
272 fpDigitFile->Close();
277 if (fppDigitArrays) {
278 for (int i=0; i<fNofDDLs; i++) {
279 if (fppDigitArrays[i]) delete fppDigitArrays[i];
281 delete[] fppDigitArrays;
288 int AliHLTOUTComponent::DumpEvent( const AliHLTComponentEventData& evtData,
289 const AliHLTComponentBlockData* blocks,
290 AliHLTComponentTriggerData& /*trigData*/ )
292 // see header file for class documentation
294 HLTInfo("write %d output block(s)", evtData.fBlockCnt);
297 AliHLTUInt32_t eventType=gkAliEventTypeUnknown;
298 bool bIsDataEvent=IsDataEvent(&eventType);
300 homer_uint64 homerHeader[kCount_64b_Words];
301 HOMERBlockDescriptor homerDescriptor(homerHeader);
302 for (int n=0; n<(int)evtData.fBlockCnt; n++ ) {
303 if (blocks[n].fDataType==kAliHLTDataTypeEvent ||
304 blocks[n].fDataType==kAliHLTDataTypeSOR ||
305 blocks[n].fDataType==kAliHLTDataTypeEOR ||
306 blocks[n].fDataType==kAliHLTDataTypeComConf ||
307 blocks[n].fDataType==kAliHLTDataTypeUpdtDCS)
309 // the special events have to be ignored.
313 (blocks[n].fDataType!=kAliHLTDataTypeComponentTable))
315 // In simulation, there are no SOR and EOR events created. Thats
316 // why all data blocks of those events are currently ignored.
317 // Strictly speaking, components should not create output blocks
318 // on the SOR/EOR event
320 // Exeptions: some blocks are added, the buffer must be prepared and
321 // kept since the pointers will be invalid
322 // - kAliHLTDataTypeComponentTable component table entries
325 memset( homerHeader, 0, sizeof(homer_uint64)*kCount_64b_Words );
326 homerDescriptor.Initialize();
327 // for some traditional reason the TCPDumpSubscriber swaps the bytes
328 // of the data type id and data type origin. Actually I do not understand
329 // the corresponding code line
330 // homerBlock.SetType( blocks[n].fDataType.fID );
331 // this compiles in the PubSub framework and in addition does a byte swap
333 homer_uint64 origin=0;
334 memcpy(&id, blocks[n].fDataType.fID, sizeof(homer_uint64));
335 memcpy(((AliHLTUInt8_t*)&origin)+sizeof(homer_uint32), blocks[n].fDataType.fOrigin, sizeof(homer_uint32));
336 homerDescriptor.SetType(AliHLTOUT::ByteSwap64(id));
337 homerDescriptor.SetSubType1(AliHLTOUT::ByteSwap64(origin));
338 homerDescriptor.SetSubType2(blocks[n].fSpecification);
339 homerDescriptor.SetBlockSize(blocks[n].fSize);
341 writerNo=ShuffleWriters(fWriters, blocks[n].fSize);
343 assert(writerNo>=0 && writerNo<(int)fWriters.size());
344 // I'm puzzled by the different headers, buffers etc. used in the
345 // HOMER writer/data. In additional, there is no type check as there
346 // are void pointers used and names mixed.
347 // It seems that HOMERBlockDescriptor is just a tool to set the
348 // different fields in the homer header, which is an array of 64 bit
350 fWriters[writerNo]->AddBlock(homerHeader, blocks[n].fPtr);
355 if (iResult>=0 && !bIsDataEvent && fNofDDLs>=2) {
356 // data blocks from a special event are kept to be added to the
357 // following event. In the current implementation at least 2 DDLs
358 // are required to allow to keep the blocks of the SOR event and
359 // include it in the first event. If only one writer is available
360 // the blocks are ignored. For the moment this is not expexted to
361 // be a problem since components should not gererate anything on
362 // SOR/EOR. The only case is the list of AliHLTComponentTableEntry
363 // transmitted for component statistics in debug mode.
364 if (fReservedWriter>=0) {
365 HLTWarning("overriding previous buffer of non-data event data blocks");
367 const AliHLTUInt8_t* pBuffer=NULL;
369 // TODO: not yet clear whether it is smart to send the event id of
370 // this special event or if it should be set from the id of the
371 // following event where the data will be added
372 if (blockCount>0 && (bufferSize=FillOutputBuffer(evtData.fEventID, fWriters[writerNo], pBuffer))>0) {
373 fReservedWriter=writerNo;
374 fReservedData=bufferSize;
376 fWriters[writerNo]->Clear();
377 } else if (iResult>=0 && !bIsDataEvent && fNofDDLs<2 && blockCount>0) {
378 HLTWarning("ignoring %d block(s) for special event of type %d: at least 2 DDLs are required", blockCount, eventType);
381 if (iResult>=0 && bIsDataEvent) {
382 iResult=Write(GetEventCount(), GetRunLoader());
389 int AliHLTOUTComponent::FillESD(int /*eventNo*/, AliRunLoader* /*runLoader*/, AliESDEvent* /*esd*/)
391 // see header file for class documentation
392 // 2010-04-14 nothing to do any more. The data is written at the end of
397 int AliHLTOUTComponent::Write(int eventNo, AliRunLoader* runLoader)
399 // see header file for class documentation
402 if (fWriters.size()==0) return 0;
404 if (fReservedWriter>=0) {
405 if (fOptions&kWriteDigits) WriteDigitArray(fReservedWriter, &fBuffer[0], fReservedData);
406 if (fOptions&kWriteRawFiles) WriteRawFile(eventNo, runLoader, fReservedWriter, &fBuffer[0], fReservedData);
410 // search for the writer with the biggest data volume in order to allocate the
411 // output buffer of sufficient size
413 for (size_t i=0; i<fWriters.size(); i++) {
414 if ((int)i==fReservedWriter) continue;
417 if (sorted.size()==0 || fWriters[i]->GetTotalMemorySize()<=fWriters[sorted[0]]->GetTotalMemorySize()) {
420 sorted.insert(sorted.begin(), i);
426 vector<int>::iterator ddlno=sorted.begin();
427 while (ddlno!=sorted.end()) {
428 const AliHLTUInt8_t* pBuffer=NULL;
431 if ((bufferSize=FillOutputBuffer(eventNo, fWriters[*ddlno], pBuffer))>0) {
432 if (fOptions&kWriteDigits) WriteDigitArray(*ddlno, pBuffer, bufferSize);
433 if (fOptions&kWriteRawFiles) WriteRawFile(eventNo, runLoader, *ddlno, pBuffer, bufferSize);
435 fWriters[*ddlno]->Clear();
438 if (fOptions&kWriteDigits) WriteDigits(eventNo, runLoader);
442 int AliHLTOUTComponent::ShuffleWriters(AliHLTMonitoringWriterPVector &list, AliHLTUInt32_t /*size*/)
444 // see header file for class documentation
446 assert(list.size()>0);
447 if (list.size()==0) return iResult;
450 for (i=0; i<list.size(); i++) {
451 if ((int)i==fReservedWriter) continue;
452 if (list[i]->GetTotalMemorySize()==0)
453 writers.push_back(i);
454 else if (iResult<0 ||
455 list[i]->GetTotalMemorySize()<list[iResult]->GetTotalMemorySize())
459 if (writers.size()>0) {
461 if (writers.size()>0) {
462 // shuffle among the empty writers
465 rand.SetSeed(dt.Get()*(iResult+1));
466 i=rand.Integer(writers.size()-1);
467 assert(i>0 && i<writers.size()-1);
471 // take the writer with the least data volume
477 int AliHLTOUTComponent::FillOutputBuffer(int eventNo, AliHLTMonitoringWriter* pWriter, const AliHLTUInt8_t* &pBuffer)
479 // see header file for class documentation
481 unsigned int bufferSize=0;
483 // space for common data header
484 bufferSize+=sizeof(AliRawDataHeader);
485 assert(sizeof(AliRawDataHeader)==32);
487 // space for HLT event header
488 bufferSize+=sizeof(AliHLTOUT::AliHLTOUTEventHeader);
490 // space for payload from the writer
491 if (pWriter) bufferSize+=pWriter->GetTotalMemorySize();
493 // payload data must be aligned to 32bit
494 bufferSize=(bufferSize+3)/4;
497 if (bufferSize>fBuffer.size())
498 fBuffer.resize(bufferSize);
500 // reset the last 32bit word, rest will be overwritten
501 memset(&fBuffer[bufferSize-4], 0, 4);
503 if (bufferSize<=fBuffer.size()) {
504 AliRawDataHeader* pCDH=reinterpret_cast<AliRawDataHeader*>(&fBuffer[0]);
505 AliHLTOUT::AliHLTOUTEventHeader* pHLTH=reinterpret_cast<AliHLTOUT::AliHLTOUTEventHeader*>(&fBuffer[sizeof(AliRawDataHeader)]);
506 *pCDH = AliRawDataHeader(); // Fill with default values.
507 memset(pHLTH, 0, sizeof(AliHLTOUT::AliHLTOUTEventHeader));
511 pWriter->Copy(&fBuffer[sizeof(AliRawDataHeader)+sizeof(AliHLTOUT::AliHLTOUTEventHeader)], 0, 0, 0, 0);
512 pHLTH->fLength=pWriter->GetTotalMemorySize();
513 // set status bit to indicate HLT payload
514 pCDH->fStatusMiniEventID|=0x1<<(AliHLTOUT::kCDHStatusFlagsOffset+AliHLTOUT::kCDHFlagsHLTPayload);
516 pHLTH->fLength+=sizeof(AliHLTOUT::AliHLTOUTEventHeader);
517 // pHLTH->fEventIDLow is already set to zero in memset above.
518 pHLTH->fEventIDLow = eventNo;
519 // version does not really matter since we do not add decision data
520 pHLTH->fVersion=AliHLTOUT::kVersion1;
522 pCDH->fSize=bufferSize;
523 pCDH->fStatusMiniEventID|=0x1<<(AliHLTOUT::kCDHStatusFlagsOffset + AliHLTOUT::kCDHFlagsHLTPayload);
526 iResult=(int)bufferSize;
535 int AliHLTOUTComponent::WriteDigitArray(int hltddl, const AliHLTUInt8_t* pBuffer, unsigned int bufferSize)
537 // see header file for class documentation
539 assert(hltddl<fNofDDLs);
540 if (hltddl>=fNofDDLs) return -ERANGE;
542 if (!fppDigitArrays) {
543 fppDigitArrays=new TArrayC*[fNofDDLs];
544 if (fppDigitArrays) {
545 for (int i=0; i<fNofDDLs; i++) {
546 fppDigitArrays[i]=new TArrayC(0);
550 if (fppDigitArrays && fppDigitArrays[hltddl]) {
551 fppDigitArrays[hltddl]->Set(bufferSize, reinterpret_cast<const Char_t*>(pBuffer));
558 int AliHLTOUTComponent::WriteDigits(int /*eventNo*/, AliRunLoader* /*runLoader*/)
560 // see header file for class documentation
563 fpDigitFile=new TFile(fDigitFileName, "RECREATE");
565 if (fpDigitFile && !fpDigitFile->IsZombie()) {
567 fpDigitTree=new TTree("rawhltout","HLTOUT raw data");
568 if (fpDigitTree && fppDigitArrays) {
569 for (int i=0; i<fNofDDLs; i++) {
570 const char* branchName=AliDAQ::DdlFileName("HLT", i);
571 if (fppDigitArrays[i]) fpDigitTree->Branch(branchName, "TArrayC", &fppDigitArrays[i], 32000/*just as the default*/, 0);
577 int res=fpDigitTree->Fill();
578 HLTDebug("writing digit tree: %d", res);
580 res=fpDigitTree->Write("",TObject::kOverwrite);
581 HLTDebug("writing digit tree: %d", res);
585 fpDigitTree->Write("",TObject::kOverwrite);
587 if (fppDigitArrays) for (int i=0; i<fNofDDLs; i++) {
588 if (fppDigitArrays[i]) fppDigitArrays[i]->Set(0);
592 const char* errorMsg="";
593 if (GetEventCount()==5) {
594 errorMsg=" (suppressing further error messages)";
596 if (GetEventCount()<5) {
597 HLTError("can not open HLT digit file %s%s", fDigitFileName.Data(), errorMsg);
604 int AliHLTOUTComponent::WriteRawFile(int eventNo, AliRunLoader* /*runLoader*/, int hltddl, const AliHLTUInt8_t* pBuffer, unsigned int bufferSize)
606 // see header file for class documentation
608 const char* fileName=AliDAQ::DdlFileName("HLT", hltddl);
609 assert(fileName!=NULL);
611 filePath.Form("raw%d/", eventNo);
612 if (gSystem->AccessPathName(filePath)) {
613 TString command="mkdir "; command+=filePath;
614 gSystem->Exec(command);
618 ios::openmode filemode=(ios::openmode)0;
619 ofstream rawfile(filePath.Data(), filemode);
620 if (rawfile.good()) {
621 if (pBuffer && bufferSize>0) {
622 rawfile.write(reinterpret_cast<const char*>(pBuffer), bufferSize);
624 HLTWarning("writing zero length raw data file %s");
626 HLTDebug("wrote %d byte(s) to file %s", bufferSize, filePath.Data());
628 HLTError("can not open file %s for writing", filePath.Data());
636 void AliHLTOUTComponent::SetGlobalOption(unsigned int options)
638 // see header file for class documentation
642 void AliHLTOUTComponent::ClearGlobalOption(unsigned int options)
644 // see header file for class documentation
648 bool AliHLTOUTComponent::TestGlobalOption(unsigned int option)
650 // see header file for class documentation
651 return (fgOptions&option)!=0;