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)
63 , fRoundRobinCounter(0)
65 // see header file for class documentation
67 // refer to README to build package
69 // visit http://web.ift.uib.no/~kjeks/doc/alice-hlt
71 fIdFirstDDL=AliDAQ::DdlIDOffset("HLT");
72 fNofDDLs=AliDAQ::NumberOfDdls("HLT");
74 if (fType!=kGlobal && fType!=kDigits && fType!=kRaw) {
75 ALIHLTERRORGUARD(1, "invalid component type %d", fType);
79 int AliHLTOUTComponent::fgOptions=kWriteRawFiles|kWriteDigits;
81 AliHLTOUTComponent::~AliHLTOUTComponent()
84 if (fpLibManager) delete fpLibManager;
88 const char* AliHLTOUTComponent::GetComponentID()
90 // overloaded from AliHLTComponent: get component id
92 case kDigits: return "HLTOUTdigits";
93 case kRaw: return "HLTOUTraw";
101 void AliHLTOUTComponent::GetInputDataTypes( vector<AliHLTComponentDataType>& list)
103 // overloaded from AliHLTComponent: indicate input data types
105 list.push_back(kAliHLTAnyDataType);
108 AliHLTComponent* AliHLTOUTComponent::Spawn()
110 // overloaded from AliHLTComponent: create instance
111 return new AliHLTOUTComponent(fType);
114 int AliHLTOUTComponent::DoInit( int argc, const char** argv )
116 // overloaded from AliHLTComponent: initialization
121 fOptions|=kWriteDigits; fOptions&=~kWriteRawFiles;
122 HLTInfo("initializing HLTOUT component for digits generation");
125 fOptions|=kWriteRawFiles; fOptions&=~kWriteDigits;
126 HLTInfo("initializing HLTOUT component for raw data generation");
133 if ((iResult=ConfigureFromArgumentString(argc, argv))<0) return iResult;
135 // Create a new library manager and allocate the appropriate number of
136 // HOMER writers for the HLTOUT component.
137 if (!fpLibManager) fpLibManager=new AliHLTHOMERLibManager;
140 for (writerNo=0; writerNo<fNofDDLs; writerNo++) {
141 AliHLTMonitoringWriter* pWriter=fpLibManager->OpenWriter();
143 HLTDebug("HOMER writer %p added", pWriter);
144 fWriters.push_back(pWriter);
146 HLTError("can not open HOMER writer");
158 int AliHLTOUTComponent::ScanConfigurationArgument(int argc, const char** argv)
160 // overloaded from AliHLTComponent: argument scan
161 if (argc<=0) return 0;
163 TString argument=argv[i];
167 // specify number of ddl links
168 if (argument.CompareTo("-links")==0) {
169 if (++i>=argc) return -EPROTO;
170 TString parameter(argv[i]);
171 parameter.Remove(TString::kLeading, ' '); // remove all blanks
172 if (parameter.IsDigit()) {
173 fNofDDLs=parameter.Atoi();
175 HLTError("wrong parameter for argument %s, number expected", argument.Data());
183 if (argument.CompareTo("-digitfile")==0) {
184 if (++i>=argc) return -EPROTO;
185 fDigitFileName=argv[i];
192 if (argument.Contains(key)) {
193 argument.ReplaceAll(key, "");
194 if (argument.IsNull()) {
195 fOptions|=kWriteRawFiles;
196 } else if (argument.CompareTo("=off")==0) {
197 fOptions&=~kWriteRawFiles;
198 } else if (argument.CompareTo("=on")==0) {
199 fOptions|=kWriteRawFiles;
201 HLTError("invalid parameter for argument %s: possible %s=off/%s=on", key, key, key);
210 if (argument.Contains(key)) {
211 argument.ReplaceAll(key, "");
212 if (argument.IsNull()) {
213 fOptions|=kWriteDigits;
214 } else if (argument.CompareTo("=off")==0) {
215 fOptions&=~kWriteDigits;
216 } else if (argument.CompareTo("=on")==0) {
217 fOptions|=kWriteDigits;
219 HLTError("invalid parameter for argument %s: possible %s=off/%s=on", key, key, key);
226 // -distribute-blocks
227 key="-distribute-blocks";
228 if (argument.CompareTo(key)==0) {
229 fRoundRobinCounter=-1;
238 int AliHLTOUTComponent::DoDeinit()
240 // overloaded from AliHLTComponent: cleanup
244 AliHLTMonitoringWriterPVector::iterator element=fWriters.begin();
245 while (element!= fWriters.end()) {
247 // wanted to have a dynamic_cast<AliHLTHOMERWriter*> here, but this results into
248 // undefined symbol when loading the library
249 if (*element!=NULL) {
251 fpLibManager->DeleteWriter((AliHLTHOMERWriter*)(*element));
253 HLTError("writer instance is NULL");
255 element=fWriters.erase(element);
269 fpDigitFile->Close();
274 if (fppDigitArrays) {
275 for (int i=0; i<fNofDDLs; i++) {
276 if (fppDigitArrays[i]) delete fppDigitArrays[i];
278 delete[] fppDigitArrays;
285 int AliHLTOUTComponent::DumpEvent( const AliHLTComponentEventData& evtData,
286 const AliHLTComponentBlockData* blocks,
287 AliHLTComponentTriggerData& /*trigData*/ )
289 // overloaded from AliHLTDataSink: event processing
291 HLTInfo("write %d output block(s)", evtData.fBlockCnt);
294 AliHLTUInt32_t eventType=gkAliEventTypeUnknown;
295 bool bIsDataEvent=IsDataEvent(&eventType);
297 homer_uint64 homerHeader[kCount_64b_Words];
298 HOMERBlockDescriptor homerDescriptor(homerHeader);
299 for (int n=0; n<(int)evtData.fBlockCnt; n++ ) {
300 if (blocks[n].fDataType==kAliHLTDataTypeEvent ||
301 blocks[n].fDataType==kAliHLTDataTypeSOR ||
302 blocks[n].fDataType==kAliHLTDataTypeEOR ||
303 blocks[n].fDataType==kAliHLTDataTypeComConf ||
304 blocks[n].fDataType==kAliHLTDataTypeUpdtDCS)
306 // the special events have to be ignored.
310 (blocks[n].fDataType!=kAliHLTDataTypeComponentTable))
312 // In simulation, there are no SOR and EOR events created. Thats
313 // why all data blocks of those events are currently ignored.
314 // Strictly speaking, components should not create output blocks
315 // on the SOR/EOR event
317 // Exeptions: some blocks are added, the buffer must be prepared and
318 // kept since the pointers will be invalid
319 // - kAliHLTDataTypeComponentTable component table entries
322 memset( homerHeader, 0, sizeof(homer_uint64)*kCount_64b_Words );
323 homerDescriptor.Initialize();
324 // for some traditional reason the TCPDumpSubscriber swaps the bytes
325 // of the data type id and data type origin. Actually I do not understand
326 // the corresponding code line
327 // homerBlock.SetType( blocks[n].fDataType.fID );
328 // this compiles in the PubSub framework and in addition does a byte swap
330 homer_uint64 origin=0;
331 memcpy(&id, blocks[n].fDataType.fID, sizeof(homer_uint64));
332 memcpy(((AliHLTUInt8_t*)&origin)+sizeof(homer_uint32), blocks[n].fDataType.fOrigin, sizeof(homer_uint32));
333 homerDescriptor.SetType(AliHLTOUT::ByteSwap64(id));
334 homerDescriptor.SetSubType1(AliHLTOUT::ByteSwap64(origin));
335 homerDescriptor.SetSubType2(blocks[n].fSpecification);
336 homerDescriptor.SetBlockSize(blocks[n].fSize);
338 writerNo=ShuffleWriters(fWriters, blocks[n].fSize);
340 assert(writerNo>=0 && writerNo<(int)fWriters.size());
341 // I'm puzzled by the different headers, buffers etc. used in the
342 // HOMER writer/data. In additional, there is no type check as there
343 // are void pointers used and names mixed.
344 // It seems that HOMERBlockDescriptor is just a tool to set the
345 // different fields in the homer header, which is an array of 64 bit
347 fWriters[writerNo]->AddBlock(homerHeader, blocks[n].fPtr);
352 if (iResult>=0 && !bIsDataEvent && fNofDDLs>=2) {
353 // data blocks from a special event are kept to be added to the
354 // following event. In the current implementation at least 2 DDLs
355 // are required to allow to keep the blocks of the SOR event and
356 // include it in the first event. If only one writer is available
357 // the blocks are ignored. For the moment this is not expexted to
358 // be a problem since components should not gererate anything on
359 // SOR/EOR. The only case is the list of AliHLTComponentTableEntry
360 // transmitted for component statistics in debug mode.
361 if (fReservedWriter>=0) {
362 HLTWarning("overriding previous buffer of non-data event data blocks");
364 const AliHLTUInt8_t* pBuffer=NULL;
366 // TODO: not yet clear whether it is smart to send the event id of
367 // this special event or if it should be set from the id of the
368 // following event where the data will be added
369 if (blockCount>0 && (bufferSize=FillOutputBuffer(evtData.fEventID, fWriters[writerNo], pBuffer))>0) {
370 fReservedWriter=writerNo;
371 fReservedData=bufferSize;
373 fWriters[writerNo]->Clear();
374 } else if (iResult>=0 && !bIsDataEvent && fNofDDLs<2 && blockCount>0) {
375 HLTWarning("ignoring %d block(s) for special event of type %d: at least 2 DDLs are required", blockCount, eventType);
378 if (iResult>=0 && bIsDataEvent) {
379 iResult=Write(GetEventCount(), GetRunLoader());
382 if (fRoundRobinCounter>=0) {
383 if (++fRoundRobinCounter>=fNofDDLs) fRoundRobinCounter=0;
390 int AliHLTOUTComponent::FillESD(int /*eventNo*/, AliRunLoader* /*runLoader*/, AliESDEvent* /*esd*/)
392 // Nop. The data is written at the end of DumpEvent
396 int AliHLTOUTComponent::Write(int eventNo, AliRunLoader* runLoader)
398 // write digits and raw files for the current event
401 if (fWriters.size()==0) return 0;
403 if (fReservedWriter>=0) {
404 if (fOptions&kWriteDigits) WriteDigitArray(fReservedWriter, &fBuffer[0], fReservedData);
405 if (fOptions&kWriteRawFiles) WriteRawFile(eventNo, runLoader, fReservedWriter, &fBuffer[0], fReservedData);
409 // search for the writer with the biggest data volume in order to allocate the
410 // output buffer of sufficient size
412 for (size_t i=0; i<fWriters.size(); i++) {
413 if ((int)i==fReservedWriter) continue;
416 if (sorted.size()==0 || fWriters[i]->GetTotalMemorySize()<=fWriters[sorted[0]]->GetTotalMemorySize()) {
419 sorted.insert(sorted.begin(), i);
425 vector<int>::iterator ddlno=sorted.begin();
426 while (ddlno!=sorted.end()) {
427 const AliHLTUInt8_t* pBuffer=NULL;
430 if ((bufferSize=FillOutputBuffer(eventNo, fWriters[*ddlno], pBuffer))>0) {
431 if (fOptions&kWriteDigits) WriteDigitArray(*ddlno, pBuffer, bufferSize);
432 if (fOptions&kWriteRawFiles &&
433 (fRoundRobinCounter<0 || fRoundRobinCounter==*ddlno))
434 WriteRawFile(eventNo, runLoader, *ddlno, pBuffer, bufferSize);
436 fWriters[*ddlno]->Clear();
439 if (fOptions&kWriteDigits) WriteDigits(eventNo, runLoader);
443 int AliHLTOUTComponent::ShuffleWriters(AliHLTMonitoringWriterPVector &list, AliHLTUInt32_t /*size*/)
445 /// get a writer for the next block
446 /// in round robin mode (like the online HLTOUT) all blocks of one event go to the same link
447 /// this is now also the default behavior of the HLTOUTComponent and indicated by
448 /// fRoundRobinCounter>=0
449 /// Writers are selected randomly otherwise.
450 if (fRoundRobinCounter>=0) {
451 if (fRoundRobinCounter==fReservedWriter) {
452 if (++fRoundRobinCounter>=fNofDDLs) fRoundRobinCounter=0;
453 if (fRoundRobinCounter==fReservedWriter) {
454 HLTWarning("there are not enough links to use a reserved writer, discarding data in reserved writer %d (total %d)",
455 fReservedWriter, fNofDDLs);
459 return fRoundRobinCounter;
463 assert(list.size()>0);
464 if (list.size()==0) return iResult;
467 for (i=0; i<list.size(); i++) {
468 if ((int)i==fReservedWriter) continue;
469 if (list[i]->GetTotalMemorySize()==0)
470 writers.push_back(i);
471 else if (iResult<0 ||
472 list[i]->GetTotalMemorySize()<list[iResult]->GetTotalMemorySize())
476 if (writers.size()>0) {
478 if (writers.size()>0) {
479 // shuffle among the empty writers
482 rand.SetSeed(dt.Get()*(iResult+1));
483 i=rand.Integer(writers.size()-1);
484 assert(i>0 && i<writers.size()-1);
488 // take the writer with the least data volume
494 int AliHLTOUTComponent::FillOutputBuffer(int eventNo, AliHLTMonitoringWriter* pWriter, const AliHLTUInt8_t* &pBuffer)
496 // prepare the output buffer for writing, consists of
500 // buffer is allocated internally and data is valid until next call
502 unsigned int bufferSize=0;
504 // space for common data header
505 bufferSize+=sizeof(AliRawDataHeader);
506 assert(sizeof(AliRawDataHeader)==32);
508 // space for HLT event header
509 bufferSize+=sizeof(AliHLTOUT::AliHLTOUTEventHeader);
511 // space for payload from the writer
512 if (pWriter) bufferSize+=pWriter->GetTotalMemorySize();
514 // payload data must be aligned to 32bit
515 bufferSize=(bufferSize+3)/4;
518 if (bufferSize>fBuffer.size())
519 fBuffer.resize(bufferSize);
521 // reset the last 32bit word, rest will be overwritten
522 memset(&fBuffer[bufferSize-4], 0, 4);
524 if (bufferSize<=fBuffer.size()) {
525 AliRawDataHeader* pCDH=reinterpret_cast<AliRawDataHeader*>(&fBuffer[0]);
526 AliHLTOUT::AliHLTOUTEventHeader* pHLTH=reinterpret_cast<AliHLTOUT::AliHLTOUTEventHeader*>(&fBuffer[sizeof(AliRawDataHeader)]);
527 *pCDH = AliRawDataHeader(); // Fill with default values.
528 memset(pHLTH, 0, sizeof(AliHLTOUT::AliHLTOUTEventHeader));
532 pWriter->Copy(&fBuffer[sizeof(AliRawDataHeader)+sizeof(AliHLTOUT::AliHLTOUTEventHeader)], 0, 0, 0, 0);
533 pHLTH->fLength=pWriter->GetTotalMemorySize();
534 // set status bit to indicate HLT payload
535 pCDH->fStatusMiniEventID|=0x1<<(AliHLTOUT::kCDHStatusFlagsOffset+AliHLTOUT::kCDHFlagsHLTPayload);
537 pHLTH->fLength+=sizeof(AliHLTOUT::AliHLTOUTEventHeader);
538 // pHLTH->fEventIDLow is already set to zero in memset above.
539 pHLTH->fEventIDLow = eventNo;
540 // version does not really matter since we do not add decision data
541 pHLTH->fVersion=AliHLTOUT::kVersion1;
543 pCDH->fSize=bufferSize;
544 pCDH->fStatusMiniEventID|=0x1<<(AliHLTOUT::kCDHStatusFlagsOffset + AliHLTOUT::kCDHFlagsHLTPayload);
547 iResult=(int)bufferSize;
556 int AliHLTOUTComponent::WriteDigitArray(int hltddl, const AliHLTUInt8_t* pBuffer, unsigned int bufferSize)
558 // wite a buffer to the associated digit array
560 assert(hltddl<fNofDDLs);
561 if (hltddl>=fNofDDLs) return -ERANGE;
563 if (!fppDigitArrays) {
564 fppDigitArrays=new TArrayC*[fNofDDLs];
565 if (fppDigitArrays) {
566 for (int i=0; i<fNofDDLs; i++) {
567 fppDigitArrays[i]=new TArrayC(0);
571 if (fppDigitArrays && fppDigitArrays[hltddl]) {
572 fppDigitArrays[hltddl]->Set(bufferSize, reinterpret_cast<const Char_t*>(pBuffer));
579 int AliHLTOUTComponent::WriteDigits(int /*eventNo*/, AliRunLoader* /*runLoader*/)
581 // fill tree with digit arrays and write to file
582 // all links must be written, even in round robin mode, where all links but one
583 // do not contain any data blocks.
584 // This is a limitation of storing the links in a tree
587 fpDigitFile=new TFile(fDigitFileName, "RECREATE");
589 if (fpDigitFile && !fpDigitFile->IsZombie()) {
591 fpDigitTree=new TTree("rawhltout","HLTOUT raw data");
592 if (fpDigitTree && fppDigitArrays) {
593 for (int i=0; i<fNofDDLs; i++) {
594 const char* branchName=AliDAQ::DdlFileName("HLT", i);
595 if (fppDigitArrays[i]) fpDigitTree->Branch(branchName, "TArrayC", &fppDigitArrays[i], 32000/*just as the default*/, 0);
601 int res=fpDigitTree->Fill();
602 HLTDebug("writing digit tree: %d", res);
604 res=fpDigitTree->Write("",TObject::kOverwrite);
605 HLTDebug("writing digit tree: %d", res);
609 fpDigitTree->Write("",TObject::kOverwrite);
611 if (fppDigitArrays) for (int i=0; i<fNofDDLs; i++) {
612 if (fppDigitArrays[i]) fppDigitArrays[i]->Set(0);
616 const char* errorMsg="";
617 if (GetEventCount()==5) {
618 errorMsg=" (suppressing further error messages)";
620 if (GetEventCount()<5) {
621 HLTError("can not open HLT digit file %s%s", fDigitFileName.Data(), errorMsg);
628 int AliHLTOUTComponent::WriteRawFile(int eventNo, AliRunLoader* /*runLoader*/, int hltddl, const AliHLTUInt8_t* pBuffer, unsigned int bufferSize)
630 // write buffer to raw file in the current directory
631 // creates the event raw directories in the current directory
633 const char* fileName=AliDAQ::DdlFileName("HLT", hltddl);
634 assert(fileName!=NULL);
636 filePath.Form("raw%d/", eventNo);
637 if (gSystem->AccessPathName(filePath)!=0) {
638 // note: AccessPathName return 0 if the path is existing
639 TString command="mkdir "; command+=filePath;
640 gSystem->Exec(command);
644 ios::openmode filemode=(ios::openmode)0;
645 ofstream rawfile(filePath.Data(), filemode);
646 if (rawfile.good()) {
647 if (pBuffer && bufferSize>0) {
648 rawfile.write(reinterpret_cast<const char*>(pBuffer), bufferSize);
650 HLTWarning("writing zero length raw data file %s");
652 HLTDebug("wrote %d byte(s) to file %s", bufferSize, filePath.Data());
654 HLTError("can not open file %s for writing", filePath.Data());
662 void AliHLTOUTComponent::SetGlobalOption(unsigned int options)
664 // set the global options
668 void AliHLTOUTComponent::ClearGlobalOption(unsigned int options)
670 // reset the global options
674 bool AliHLTOUTComponent::TestGlobalOption(unsigned int option)
677 return (fgOptions&option)!=0;