]> git.uio.no Git - u/mrichter/AliRoot.git/blame_incremental - RAW/AliMDC.cxx
SDD QA updated in order to deal with acquisition through HLT (M. Siciliano)
[u/mrichter/AliRoot.git] / RAW / AliMDC.cxx
... / ...
CommitLineData
1// @(#)alimdc:$Name: $:$Id$
2// Author: Fons Rademakers 26/11/99
3// Updated: Dario Favretto 15/04/2003
4
5/**************************************************************************
6 * Copyright(c) 1998-2003, ALICE Experiment at CERN, All rights reserved. *
7 * *
8 * Author: The ALICE Off-line Project. *
9 * Contributors are mentioned in the code where appropriate. *
10 * *
11 * Permission to use, copy, modify and distribute this software and its *
12 * documentation strictly for non-commercial purposes is hereby granted *
13 * without fee, provided that the above copyright notice appears in all *
14 * copies and that both the copyright notice and this permission notice *
15 * appear in the supporting documentation. The authors make no claims *
16 * about the suitability of this software for any purpose. It is *
17 * provided "as is" without express or implied warranty. *
18 **************************************************************************/
19
20/* $Id$ */
21
22//////////////////////////////////////////////////////////////////////////
23// //
24// AliMDC //
25// //
26// Set of classes defining the ALICE RAW event format. The AliRawEvent //
27// class defines a RAW event. It consists of an AliEventHeader object //
28// an AliEquipmentHeader object, an AliRawData object and an array of //
29// sub-events, themselves also being AliRawEvents. The number of //
30// sub-events depends on the number of DATE LDC's. //
31// The AliRawEvent objects are written to a ROOT file using different //
32// technologies, i.e. to local disk via AliRawDB or via rfiod using //
33// AliRawRFIODB or via rootd using AliRawRootdDB or to CASTOR via //
34// rootd using AliRawCastorDB (and for performance testing there is //
35// also AliRawNullDB). //
36// The AliStats class provides statics information that is added as //
37// a single keyed object to each raw file. //
38// The AliTagDB provides an interface to a TAG database. //
39// The AliMDC class is usid by the "alimdc" stand-alone program //
40// that reads data directly from DATE. //
41// //
42//////////////////////////////////////////////////////////////////////////
43
44#include <sys/types.h>
45#include <sys/stat.h>
46
47#include <errno.h>
48
49#include <TSystem.h>
50#include <TROOT.h>
51#include <TStopwatch.h>
52#include <TPluginManager.h>
53
54#include <sys/uio.h>
55#ifdef USE_EB
56#include "libDateEb.h"
57#endif
58
59#include "AliMDC.h"
60
61#include <AliLog.h>
62#include <AliESDEvent.h>
63
64#include "AliRawEvent.h"
65#include "AliRawEventHeaderBase.h"
66#include "AliRawEquipment.h"
67#include "AliRawEquipmentHeader.h"
68#include "AliRawData.h"
69#include "AliStats.h"
70#include "AliRawDB.h"
71#include "AliRawRFIODB.h"
72#include "AliRawCastorDB.h"
73#include "AliRawRootdDB.h"
74#include "AliRawNullDB.h"
75#include "AliTagDB.h"
76#include "AliRawEventTag.h"
77#include "AliFilter.h"
78
79
80
81ClassImp(AliMDC)
82
83
84// Filter names
85const char* const AliMDC::fgkFilterName[kNFilters] = {"AliHoughFilter"};
86
87//______________________________________________________________________________
88AliMDC::AliMDC(Int_t compress, Bool_t deleteFiles, EFilterMode filterMode,
89 Double_t maxSizeTagDB, const char* fileNameTagDB,
90 const char *guidFileFolder,
91 Int_t basketsize) :
92 fEvent(new AliRawEvent),
93 fESD(NULL),
94 fStats(NULL),
95 fRawDB(NULL),
96 fTagDB(NULL),
97 fEventTag(new AliRawEventTag),
98 fCompress(compress),
99 fBasketSize(basketsize),
100 fDeleteFiles(deleteFiles),
101 fFilterMode(filterMode),
102 fFilters(),
103 fStop(kFALSE),
104 fIsTagDBCreated(kFALSE),
105 fMaxSizeTagDB(maxSizeTagDB),
106 fFileNameTagDB(fileNameTagDB),
107 fGuidFileFolder(guidFileFolder)
108{
109 // Create MDC processor object.
110 // compress is the file compression mode.
111 // If deleteFiles is kTRUE the raw data files will be deleted after they
112 // were closed.
113 // If the filterMode if kFilterOff no filter algorithm will be, if it is
114 // kFilterTransparent the algorthims will be run but no events will be
115 // rejected, if it is kFilterOn the filters will be run and the event will
116 // be rejected if all filters return kFALSE.
117 // If maxSizeTagDB is greater than 0 it determines the maximal size of the
118 // tag DB and then fileNameTagDB is the directory name for the tag DB.
119 // Otherwise fileNameTagDB is the file name of the tag DB. If it is NULL
120 // no tag DB will be created.
121 // Optional 'guidFileFolder' specifies the folder in which *.guid files
122 // will be stored. In case this option is not given, the *.guid files
123 // will be written to the same folder as the raw-data files.
124
125 // Set the maximum tree size to 19GB
126 // in order to allow big raw data files
127 TTree::SetMaxTreeSize(20000000000LL);
128
129 // This line is needed in case of a stand-alone application w/o
130 // $ROOTSYS/etc/system.rootrc file
131 gROOT->GetPluginManager()->AddHandler("TVirtualStreamerInfo",
132 "*",
133 "TStreamerInfo",
134 "RIO",
135 "TStreamerInfo()");
136
137 if (fFilterMode != kFilterOff) {
138 fESD = new AliESDEvent();
139 }
140
141 // Create the guid files folder if it does not exist
142 if (!fGuidFileFolder.IsNull()) {
143 gSystem->ResetErrno();
144 gSystem->MakeDirectory(fGuidFileFolder.Data());
145 if (gSystem->GetErrno() && gSystem->GetErrno() != EEXIST) {
146 SysError("AliMDC", "mkdir %s", fGuidFileFolder.Data());
147 }
148 }
149
150 // install SIGUSR1 handler to allow clean interrupts
151 gSystem->AddSignalHandler(new AliMDCInterruptHandler(this));
152
153 // create the high level filters
154 if (fFilterMode != kFilterOff) {
155 for (Int_t iFilter = 0; iFilter < kNFilters; iFilter++) {
156 TClass* filterClass = gROOT->GetClass(fgkFilterName[iFilter]);
157 if (!filterClass) {
158 Warning("AliMDC", "no filter class %s found", fgkFilterName[iFilter]);
159 continue;
160 }
161 AliFilter* filter = (AliFilter*) filterClass->New();
162 if (!filter) {
163 Warning("AliMDC", "creation of filter %s failed", fgkFilterName[iFilter]);
164 continue;
165 }
166 fFilters.Add(filter);
167 }
168 }
169}
170
171//______________________________________________________________________________
172AliMDC::~AliMDC()
173{
174// destructor
175
176 fFilters.Delete();
177 if(fTagDB) delete fTagDB;
178 delete fRawDB;
179 delete fStats;
180 delete fESD;
181 delete fEvent;
182 delete fEventTag;
183}
184
185//______________________________________________________________________________
186Int_t AliMDC::Open(EWriteMode mode, const char* fileName,
187 Double_t maxFileSize,
188 const char* fs1, const char* fs2)
189{
190// open a new raw DB file
191
192 if (mode == kRFIO)
193 fRawDB = new AliRawRFIODB(fEvent, fESD, fCompress, fileName, fBasketSize);
194 else if (mode == kROOTD)
195 fRawDB = new AliRawRootdDB(fEvent, fESD, fCompress, fileName, fBasketSize);
196 else if (mode == kCASTOR)
197 fRawDB = new AliRawCastorDB(fEvent, fESD, fCompress, fileName, fBasketSize);
198 else if (mode == kDEVNULL)
199 fRawDB = new AliRawNullDB(fEvent, fESD, fCompress, fileName, fBasketSize);
200 else
201 fRawDB = new AliRawDB(fEvent, fESD, fCompress, fileName, fBasketSize);
202 fRawDB->SetDeleteFiles(fDeleteFiles);
203
204 if (fRawDB->IsZombie()) {
205 delete fRawDB;
206 fRawDB = NULL;
207 return -1;
208 }
209
210 if (fileName == NULL) {
211 fRawDB->SetMaxSize(maxFileSize);
212 fRawDB->SetFS(fs1, fs2);
213 if (!fRawDB->Create()) {
214 delete fRawDB;
215 fRawDB = NULL;
216 return -1;
217 }
218 }
219
220 if (!fRawDB->WriteGuidFile(fGuidFileFolder)) {
221 delete fRawDB;
222 fRawDB = NULL;
223 return -2;
224 }
225
226 Info("Open", "Filling raw DB %s\n", fRawDB->GetDBName());
227
228 // Create AliStats object
229 fStats = new AliStats(fRawDB->GetDBName(), fCompress,
230 fFilterMode != kFilterOff);
231 return 0;
232}
233
234//______________________________________________________________________________
235Int_t AliMDC::ProcessEvent(void* event, Bool_t isIovecArray)
236{
237// Convert the DATE event to an AliRawEvent object and store it in the raw DB,
238// optionally also run the filter.
239// event is either a pointer to the streamlined event
240// or, if isIovecArray is kTRUE, a pointer to an array of iovecs with one
241// iovec per subevent (used by the event builder).
242// The return value is the number of written bytes or an error code
243 const Long64_t kFileSizeErrorLevel = 19000000000LL;
244
245 Long64_t currentFileSize = GetTotalSize();
246 AliDebug(1,Form("current file size is %lld bytes",currentFileSize));
247 if(currentFileSize > kFileSizeErrorLevel) {
248 Error("ProcessEvent", "file size (%lu) exceeds the limit "
249 , currentFileSize);
250 return kErrFileSize;
251 }
252
253 Int_t status;
254 char* data = (char*) event;
255 if (isIovecArray) data = (char*) ((iovec*) event)[0].iov_base;
256
257 // Shortcut for easy header access
258 AliRawEventHeaderBase *header = fEvent->GetHeader(data);
259
260 // Read event header
261 if ((status = header->ReadHeader(data)) != (Int_t)header->GetHeadSize()) {
262 return kErrHeader;
263 }
264
265 if (AliDebugLevel() > 2) ToAliDebug(3, header->Dump(););
266
267 // Check event type and skip "Start of Run", "End of Run",
268 // "Start of Run Files" and "End of Run Files"
269 Int_t size = header->GetEventSize() - header->GetHeadSize();
270
271 AliDebug(1, Form("Processing %s (%d bytes)", header->GetTypeName(), size));
272
273 // Amount of data left to read for this event
274 Int_t toRead = size;
275
276 // StartOfRun, EndOfRun etc. events have no payload
277 // Nevertheless, store the event headers in the tree
278 if (toRead > 0) {
279
280 // If there is less data for this event than the next sub-event
281 // header, something is wrong. Skip to next event...
282 if (toRead < (Int_t)header->GetHeadSize()) {
283 Error("ProcessEvent", "header size (%d) exceeds number of bytes "
284 "to read (%d)", header->GetHeadSize(), toRead);
285 if (AliDebugLevel() > 0) ToAliDebug(1, header->Dump(););
286 return kErrHeaderSize;
287 }
288
289 // Loop over all sub-events... (LDCs)
290 Int_t nsub = 1;
291 while (toRead > 0) {
292 if (isIovecArray) data = (char*) ((iovec*) event)[nsub].iov_base;
293
294 AliDebug(1, Form("reading LDC %d", nsub));
295
296 AliRawEvent *subEvent = fEvent->NextSubEvent();
297
298 // Read sub-event header
299 AliRawEventHeaderBase *subHeader = subEvent->GetHeader(data);
300 if ((status = subHeader->ReadHeader(data)) != (Int_t)subHeader->GetHeadSize()) {
301 return kErrSubHeader;
302 }
303
304 if (AliDebugLevel() > 2) ToAliDebug(3, subHeader->Dump(););
305
306 toRead -= subHeader->GetHeadSize();
307
308 Int_t rawSize = subHeader->GetEventSize() - subHeader->GetHeadSize();
309
310 // Make sure raw data less than left over bytes for current event
311 if (rawSize > toRead) {
312 Warning("ProcessEvent", "raw data size (%d) exceeds number of "
313 "bytes to read (%d)\n", rawSize, toRead);
314 if (AliDebugLevel() > 0) ToAliDebug(1, subHeader->Dump(););
315 return kErrDataSize;
316 }
317
318 // Read Equipment Headers (in case of physics or calibration event)
319 if (header->Get("Type") == AliRawEventHeaderBase::kPhysicsEvent ||
320 header->Get("Type") == AliRawEventHeaderBase::kCalibrationEvent ||
321 header->Get("Type") == AliRawEventHeaderBase::kSystemSoftwareTriggerEvent ||
322 header->Get("Type") == AliRawEventHeaderBase::kDetectorSoftwareTriggerEvent ||
323 header->Get("Type") == AliRawEventHeaderBase::kStartOfData ||
324 header->Get("Type") == AliRawEventHeaderBase::kEndOfData) {
325 while (rawSize > 0) {
326 AliRawEquipment &equipment = *subEvent->NextEquipment();
327 AliRawEquipmentHeader &equipmentHeader =
328 *equipment.GetEquipmentHeader();
329 Int_t equipHeaderSize = equipmentHeader.HeaderSize();
330 if ((status = ReadEquipmentHeader(equipmentHeader, header->DataIsSwapped(),
331 data)) != equipHeaderSize) {
332 return kErrEquipmentHeader;
333 }
334
335 if (AliDebugLevel() > 2) ToAliDebug(3, equipmentHeader.Dump(););
336
337 toRead -= equipHeaderSize;
338 rawSize -= equipHeaderSize;
339
340 // Read equipment raw data
341 AliRawData &subRaw = *equipment.GetRawData();
342
343 Int_t eqSize = equipmentHeader.GetEquipmentSize() - equipHeaderSize;
344 if ((status = ReadRawData(subRaw, eqSize, data)) != eqSize) {
345 return kErrEquipment;
346 }
347 toRead -= eqSize;
348 rawSize -= eqSize;
349
350 }
351
352 } else { // Read only raw data but no equipment header
353 if (rawSize) {
354 AliRawEquipment &equipment = *subEvent->NextEquipment();
355 AliRawData &subRaw = *equipment.GetRawData();
356 if ((status = ReadRawData(subRaw, rawSize, data)) != rawSize) {
357 return kErrEquipment;
358 }
359 toRead -= rawSize;
360 }
361 }
362
363 nsub++;
364 }
365 }
366
367 // High Level Event Filter
368 if (fFilterMode != kFilterOff) {
369 if (header->Get("Type") == AliRawEventHeaderBase::kPhysicsEvent ||
370 header->Get("Type") == AliRawEventHeaderBase::kCalibrationEvent ||
371 header->Get("Type") == AliRawEventHeaderBase::kSystemSoftwareTriggerEvent ||
372 header->Get("Type") == AliRawEventHeaderBase::kDetectorSoftwareTriggerEvent ||
373 header->Get("Type") == AliRawEventHeaderBase::kStartOfData ||
374 header->Get("Type") == AliRawEventHeaderBase::kEndOfData) {
375 Bool_t result = kFALSE;
376 for (Int_t iFilter = 0; iFilter < fFilters.GetEntriesFast(); iFilter++) {
377 AliFilter* filter = (AliFilter*) fFilters[iFilter];
378 if (!filter) continue;
379 if (filter->Filter(fEvent, fESD)) result = kTRUE;
380 }
381 if ((fFilterMode == kFilterOn) && !result) return kFilterReject;
382 }
383 }
384
385 // Set stat info for first event of this file
386 if (fRawDB->GetEvents() == 0)
387 fStats->SetFirstId(header->Get("RunNb"), header->GetP("Id")[0]);
388
389 // Store raw event in tree
390 Int_t nBytes = fRawDB->Fill();
391
392 // Fill the event tag object
393 fEventTag->SetHeader(fEvent->GetHeader());
394 fEventTag->SetGUID(fRawDB->GetDB()->GetUUID().AsString());
395 fEventTag->SetEventNumber(fRawDB->GetEvents()-1);
396
397 // Create Tag DB here only after the raw data header
398 // version was already identified
399 if (!fIsTagDBCreated) {
400 if (!fFileNameTagDB.IsNull()) {
401 if (fMaxSizeTagDB > 0) {
402 fTagDB = new AliTagDB(fEventTag, NULL);
403 fTagDB->SetMaxSize(fMaxSizeTagDB);
404 fTagDB->SetFS(fFileNameTagDB.Data());
405 if (!fTagDB->Create()) return kErrTagFile;
406 } else {
407 fTagDB = new AliTagDB(fEventTag, fFileNameTagDB.Data());
408 if (fTagDB->IsZombie()) return kErrTagFile;
409 }
410 }
411 fIsTagDBCreated = kTRUE;
412 }
413
414 // Store event tag in tree
415 if (fTagDB) fTagDB->Fill();
416
417 // Make top event object ready for next event data
418 fEvent->Reset();
419 // Clean up HLT ESD for the next event
420 if (fESD) fESD->Reset();
421
422 if(nBytes >= 0)
423 return nBytes;
424 else
425 return kErrWriting;
426}
427
428//______________________________________________________________________________
429Long64_t AliMDC::GetTotalSize()
430{
431// return the total current raw DB file size
432
433 if (!fRawDB) return -1;
434
435 return fRawDB->GetTotalSize();
436}
437
438//______________________________________________________________________________
439Long64_t AliMDC::Close()
440{
441// close the current raw DB file
442
443 if (!fRawDB) return -1;
444
445 fRawDB->WriteStats(fStats);
446 Long64_t filesize = fRawDB->Close();
447 delete fRawDB;
448 fRawDB = NULL;
449 delete fStats;
450 fStats = NULL;
451 return filesize;
452}
453
454//______________________________________________________________________________
455Long64_t AliMDC::AutoSave()
456{
457 // Auto-save the raw-data
458 // and esd (if any) trees
459
460 if (!fRawDB) return -1;
461
462 return fRawDB->AutoSave();
463}
464
465//______________________________________________________________________________
466Int_t AliMDC::Run(const char* inputFile, Bool_t loop,
467 EWriteMode mode, Double_t maxFileSize,
468 const char* fs1, const char* fs2)
469{
470 // Run the MDC processor. Read from the input stream and only return
471 // when the input gave and EOF or a fatal error occured. On success 0
472 // is returned, 1 in case of a fatality.
473 // inputFile is the name of the DATE input file; if NULL the input will
474 // be taken from the event builder.
475 // If loop is set the same input file will be reused in an infinite loop.
476 // mode specifies the type of the raw DB.
477 // maxFileSize is the maximal size of the raw DB.
478 // fs1 and fs2 are the file system locations of the raw DB.
479
480 Info("Run", "input = %s, rawdb size = %f, filter = %s, "
481 "looping = %s, compression = %d, delete files = %s",
482 inputFile ? inputFile : "event builder", maxFileSize,
483 fFilterMode == kFilterOff ? "off" :
484 (fFilterMode == kFilterOn ? "on" : "transparent"),
485 loop ? "yes" : "no", fCompress, fDeleteFiles ? "yes" : "no");
486
487 // Open the input file
488 Int_t fd = -1;
489 if (inputFile) {
490 if ((fd = open(inputFile, O_RDONLY)) == -1) {
491 Error("Run", "cannot open input file %s", inputFile);
492 return 1;
493 }
494 }
495
496 // Used for statistics
497 TStopwatch timer;
498 timer.Start();
499 Double_t told = 0, tnew = 0;
500 Float_t chunkSize = maxFileSize/100, nextChunk = chunkSize;
501
502 // Create new raw DB.
503 if (fRawDB) Close();
504
505 if (Open(mode,NULL,maxFileSize,fs1,fs2) < 0) return 1;
506
507 // Process input stream
508#ifdef USE_EB
509 Int_t eorFlag = 0;
510#endif
511 char* event = NULL;
512 UInt_t eventSize = 0;
513 Int_t numEvents = 0;
514
515 AliRawEventHeaderBase header;
516
517 while (kTRUE) {
518
519 // If we were in looping mode stop directly after a SIGUSR1 signal
520 if (fStop) {
521 Info("Run", "Stopping loop, processed %d events", numEvents);
522 break;
523 }
524
525 if (!inputFile) { // get data from event builder
526#ifdef USE_EB
527 if ((eorFlag = ebEor())) break;
528 if ((event = (char*)ebGetNextEvent()) == (char*)-1) {
529 Error("Run", "error getting next event (%s)", ebGetLastError());
530 break;
531 }
532 if (event == 0) {
533 // no event, sleep for 1 second and try again
534 gSystem->Sleep(1000);
535 continue;
536 }
537#else
538 Error("Run", "AliMDC was compiled without event builder support");
539 delete fRawDB;
540 fRawDB = NULL;
541 delete fStats;
542 fStats = NULL;
543 return 1;
544#endif
545
546 } else { // get data from a file
547 {
548 Int_t nrecv;
549 if ((nrecv = Read(fd, header.HeaderBaseBegin(), header.HeaderBaseSize())) !=
550 header.HeaderBaseSize()) {
551 if (nrecv == 0) { // eof
552 if (loop) {
553 ::lseek(fd, 0, SEEK_SET);
554 continue;
555 } else {
556 break;
557 }
558 } else {
559 Error("Run", "error reading base header");
560 Close();
561 delete[] event;
562 return 1;
563 }
564 }
565 }
566 char *data = (char *)header.HeaderBaseBegin();
567 AliRawEventHeaderBase *hdr = AliRawEventHeaderBase::Create(data);
568 Int_t nrecv;
569 if ((nrecv = Read(fd, hdr->HeaderBegin(), hdr->HeaderSize())) !=
570 hdr->HeaderSize()) {
571 if (nrecv == 0) { // eof
572 if (loop) {
573 ::lseek(fd, 0, SEEK_SET);
574 delete hdr;
575 continue;
576 } else {
577 delete hdr;
578 break;
579 }
580 } else {
581 Error("Run", "error reading header");
582 Close();
583 delete[] event;
584 delete hdr;
585 return 1;
586 }
587 }
588 if (eventSize < hdr->GetEventSize()) {
589 delete[] event;
590 eventSize = 2 * hdr->GetEventSize();
591 event = new char[eventSize];
592 }
593 memcpy(event, header.HeaderBaseBegin(), header.HeaderBaseSize());
594 memcpy(event+hdr->HeaderBaseSize(), hdr->HeaderBegin(), hdr->HeaderSize());
595 if (hdr->GetExtendedDataSize() != 0)
596 memcpy(event+hdr->HeaderBaseSize()+hdr->HeaderSize(),
597 hdr->GetExtendedData(), hdr->GetExtendedDataSize());
598 Int_t size = hdr->GetEventSize() - hdr->GetHeadSize();
599 if (Read(fd, event + hdr->GetHeadSize(), size) != size) {
600 Error("Run", "error reading data");
601 Close();
602 delete[] event;
603 delete hdr;
604 return 1;
605 }
606 delete hdr;
607 }
608
609 Int_t result = ProcessEvent(event, !inputFile);
610 if(result < -1)
611 Error("Run", "error writing data. Error code: %d",result);
612
613 if (result >= 0) {
614 numEvents++;
615 if (!(numEvents%10))
616 printf("Processed event %d (%d)\n", numEvents, fRawDB->GetEvents());
617 }
618
619 if (result > 0) {
620 // Filling time statistics
621 if (fRawDB->GetBytesWritten() > nextChunk) {
622 tnew = timer.RealTime();
623 fStats->Fill(tnew-told);
624 told = tnew;
625 timer.Continue();
626 nextChunk += chunkSize;
627 }
628
629 // Check size of raw db. If bigger than maxFileSize, close file
630 // and continue with new file.
631 if (fRawDB->GetBytesWritten() > maxFileSize) {
632
633 printf("Written raw DB at a rate of %.1f MB/s\n",
634 fRawDB->GetBytesWritten() / timer.RealTime() / 1000000.);
635
636 // Write stats object to raw db, run db, MySQL and AliEn
637 fRawDB->WriteStats(fStats);
638 delete fStats;
639 fStats = NULL;
640
641 if (!fRawDB->NextFile()) {
642 Error("Run", "error opening next raw data file");
643 Close();
644 if (inputFile) delete[] event;
645 return 1;
646 }
647
648 printf("Filling raw DB %s\n", fRawDB->GetDBName());
649 fStats = new AliStats(fRawDB->GetDBName(), fCompress,
650 fFilterMode != kFilterOff);
651
652 timer.Start();
653 told = 0, tnew = 0;
654 nextChunk = chunkSize;
655 }
656
657 // Check size of tag db
658 if (fTagDB && fTagDB->FileFull()) {
659 if (!fTagDB->NextFile()) {
660 delete fTagDB;
661 fTagDB = 0;
662 } else {
663 printf("Filling tag DB %s\n", fTagDB->GetDBName());
664 }
665 }
666 }
667
668 // Make top event object ready for next event data
669 //printf("Event %d has %d sub-events\n", numEvents, fEvent->GetNSubEvents());
670 // fEvent->Reset();
671 // Clean up HLT ESD for the next event
672 // if (fESD) fESD->Reset();
673
674 if (!inputFile) {
675#ifdef USE_EB
676 if (!ebReleaseEvent((iovec*)event)) {
677 Error("Run", "problem releasing event (%s)", ebGetLastError());
678 break;
679 }
680#endif
681 }
682 }
683
684 printf("Written raw DB at a rate of %.1f MB/s\n",
685 fRawDB->GetBytesWritten() / timer.RealTime() / 1000000.);
686
687 // Write stats to raw db and run db and delete stats object
688 Close();
689
690 if (!inputFile) {
691#ifdef USE_EB
692 // Print eor flag
693 if (eorFlag) {
694 Info("Run", "event builder reported end of run (%d)", eorFlag);
695 }
696#endif
697 } else {
698 // Close input source
699 close(fd);
700 delete [] event;
701 }
702
703 return 0;
704}
705
706//______________________________________________________________________________
707Int_t AliMDC::Read(Int_t fd, void *buffer, Int_t length)
708{
709 // Read exactly length bytes into buffer. Returns number of bytes
710 // received, returns -1 in case of error and 0 for EOF.
711
712 errno = 0;
713
714 if (fd < 0) return -1;
715
716 Int_t n, nrecv = 0;
717 char *buf = (char *)buffer;
718
719 for (n = 0; n < length; n += nrecv) {
720 if ((nrecv = read(fd, buf+n, length-n)) <= 0) {
721 if (nrecv == 0)
722 break; // EOF
723 if (errno != EINTR)
724 SysError("Read", "read");
725 return -1;
726 }
727 }
728 return n;
729}
730
731//______________________________________________________________________________
732Int_t AliMDC::ReadEquipmentHeader(AliRawEquipmentHeader &header,
733 Bool_t isSwapped, char*& data)
734{
735 // Read equipment header info from DATE data stream. Returns bytes read
736 // (i.e. AliRawEquipmentHeader::HeaderSize()), -1 in case of error and
737 // 0 for EOF. If isSwapped is kTRUE the event data is byte swapped
738 // and we will swap the header to host format.
739
740 memcpy(header.HeaderBegin(), data, header.HeaderSize());
741 data += header.HeaderSize();
742
743 // Swap equipment header data if needed
744 if (isSwapped)
745 header.Swap();
746
747 if (header.GetEquipmentSize() < (UInt_t)header.HeaderSize()) {
748 Error("ReadEquipmentHeader", "invalid equipment header size");
749 // try recovery... how?
750 return -1;
751 }
752
753 return header.HeaderSize();
754}
755
756//______________________________________________________________________________
757Int_t AliMDC::ReadRawData(AliRawData &raw, Int_t size, char*& data)
758{
759 // Read raw data from DATE data stream. Returns bytes read (i.e.
760 // size), -1 in case of error and 0 for EOF.
761
762 raw.SetBuffer(data, size);
763 data += size;
764
765 return size;
766}
767
768//______________________________________________________________________________
769void AliMDC::Stop()
770{
771 // Stop the event loop
772
773 fStop = kTRUE;
774 if (fRawDB) fRawDB->Stop();
775}
776
777