1 /**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
7 * Permission to use, copy, modify and distribute this software and its *
8 * documentation strictly for non-commercial purposes is hereby granted *
9 * without fee, provided that the above copyright notice appears in all *
10 * copies and that both the copyright notice and this permission notice *
11 * appear in the supporting documentation. The authors make no claims *
12 * about the suitability of this software for any purpose. It is *
13 * provided "as is" without express or implied warranty. *
14 **************************************************************************/
16 ////////////////////////////////////////////////////////////////////////////
18 // Decoding data from the TRD raw stream //
19 // and translation into ADC values, on-line tracklets and tracks //
21 // Author: J. Klein (jochen.klein@cern.ch) //
23 ////////////////////////////////////////////////////////////////////////////
28 #include "TClonesArray.h"
32 #include "AliRawReader.h"
33 #include "AliTRDdigitsManager.h"
34 #include "AliTRDdigitsParam.h"
35 #include "AliTRDcalibDB.h"
36 #include "AliTRDmcmSim.h"
37 #include "AliTRDtrapConfig.h"
38 #include "AliTRDarrayADC.h"
39 #include "AliTRDarrayDictionary.h"
40 #include "AliTRDSignalIndex.h"
41 #include "AliTRDtrackletWord.h"
42 #include "AliTRDtrackletMCM.h"
43 #include "AliESDTrdTrack.h"
44 #include "AliTreeLoader.h"
45 #include "AliLoader.h"
47 #include "AliTRDrawStream.h"
50 #include "AliRunLoader.h"
52 ClassImp(AliTRDrawStream)
54 // some static information
55 Int_t AliTRDrawStream::fgMcmOrder[] = {12, 13, 14, 15,
59 Int_t AliTRDrawStream::fgRobOrder [] = {0, 1, 2, 3};
60 const Int_t AliTRDrawStream::fgkNlinks = 12;
61 const Int_t AliTRDrawStream::fgkNstacks = 5;
62 const Int_t AliTRDrawStream::fgkNsectors = 18;
63 const Int_t AliTRDrawStream::fgkNtriggers = 12;
64 const UInt_t AliTRDrawStream::fgkDataEndmarker = 0x00000000;
65 const UInt_t AliTRDrawStream::fgkTrackletEndmarker = 0x10001000;
66 const UInt_t AliTRDrawStream::fgkStackEndmarker[] = { 0xe0d01000, 0xe0d10000 };
68 const char* AliTRDrawStream::fgkErrorMessages[] = {
70 "Link monitor active",
71 "Event counter mismatch",
72 "not a TRD equipment (1024-1041)",
73 "Invalid Stack header",
74 "Invalid detector number",
75 "No digits could be retrieved from the digitsmanager",
77 "HC check bits wrong",
78 "Unexpected position in readout stream",
79 "Invalid testpattern mode",
80 "Testpattern mismatch",
81 "Number of timebins changed",
82 "ADC mask inconsistent",
83 "ADC check bits invalid",
85 "Missing expected ADC channels",
86 "Missing MCM headers",
91 Int_t AliTRDrawStream::fgErrorDebugLevel[] = {
114 AliTRDrawStream::ErrorBehav_t AliTRDrawStream::fgErrorBehav[] = {
115 AliTRDrawStream::kTolerate,
116 AliTRDrawStream::kDiscardHC,
117 AliTRDrawStream::kTolerate,
118 AliTRDrawStream::kAbort,
119 AliTRDrawStream::kAbort,
120 AliTRDrawStream::kAbort,
121 AliTRDrawStream::kAbort,
122 AliTRDrawStream::kDiscardHC,
123 AliTRDrawStream::kDiscardHC,
124 AliTRDrawStream::kTolerate,
125 AliTRDrawStream::kTolerate,
126 AliTRDrawStream::kTolerate,
127 AliTRDrawStream::kTolerate,
128 AliTRDrawStream::kTolerate,
129 AliTRDrawStream::kTolerate,
130 AliTRDrawStream::kTolerate,
131 AliTRDrawStream::kTolerate,
132 AliTRDrawStream::kTolerate,
133 AliTRDrawStream::kTolerate,
134 AliTRDrawStream::kTolerate
137 AliTRDrawStream::AliTRDrawStream(AliRawReader *rawReader) :
138 fStoreError(&AliTRDrawStream::ForgetError),
139 fRawReader(rawReader),
156 fCurrSmHeaderSize(0),
157 fCurrSmHeaderVersion(0),
158 fCurrTrailerReadout(0),
159 fCurrTrgHeaderAvail(0),
160 fCurrTrgHeaderReadout(0),
161 fCurrTrkHeaderAvail(0),
162 fCurrStackEndmarkerAvail(0),
164 fCurrTriggerEnable(0),
165 fCurrTriggerFired(0),
167 fCurrTrackletEnable(0),
169 fCurrTrkHeaderIndexWord(0x0),
170 fCurrTrkHeaderSize(0x0),
172 fCurrTrgHeaderIndexWord(0x0),
173 fCurrTrgHeaderSize(0x0),
175 fCurrStackIndexWord(0x0),
176 fCurrStackHeaderSize(0x0),
177 fCurrStackHeaderVersion(0x0),
179 fCurrCleanCheckout(0x0),
183 fCurrLinkMonitorFlags(0x0),
184 fCurrLinkDataTypeFlags(0x0),
185 fCurrLinkDebugFlags(0x0),
186 fCurrMatchFlagsSRAM(0),
187 fCurrMatchFlagsPostBP(0),
188 fCurrChecksumStack(),
213 // default constructor
215 fCurrTrkHeaderIndexWord = new UInt_t[fgkNstacks];
216 fCurrTrkHeaderSize = new UInt_t[fgkNstacks];
217 fCurrTrkFlags = new ULong64_t[fgkNsectors*fgkNstacks];
218 fCurrTrgHeaderIndexWord = new UInt_t[fgkNtriggers];
219 fCurrTrgHeaderSize = new UInt_t[fgkNtriggers];
220 fCurrTrgFlags = new UInt_t[fgkNsectors];
221 fCurrStackIndexWord = new UInt_t[fgkNstacks];
222 fCurrStackHeaderSize = new UInt_t[fgkNstacks];
223 fCurrStackHeaderVersion = new UInt_t[fgkNstacks];
224 fCurrLinkMask = new UInt_t[fgkNstacks];
225 fCurrCleanCheckout = new UInt_t[fgkNstacks];
226 fCurrBoardId = new UInt_t[fgkNstacks];
227 fCurrHwRevTMU = new UInt_t[fgkNstacks];
228 fCurrLinkMonitorFlags = new UInt_t[fgkNstacks * fgkNlinks];
229 fCurrLinkDataTypeFlags = new UInt_t[fgkNstacks * fgkNlinks];
230 fCurrLinkDebugFlags = new UInt_t[fgkNstacks * fgkNlinks];
231 for (Int_t iSector = 0; iSector < fgkNsectors; iSector++)
232 fCurrTrgFlags[iSector] = 0;
233 for (Int_t i = 0; i < 100; i++)
236 // preparing TClonesArray
237 fTrackletArray = new TClonesArray("AliTRDtrackletWord", 256);
239 // setting up the error tree
240 fErrors = new TTree("errorStats", "Error statistics");
241 fErrors->SetDirectory(0x0);
242 fErrors->Branch("error", &fLastError);
243 fErrors->SetCircular(1000);
244 for (Int_t i = 0; i < 100; i++) {
250 AliTRDrawStream::~AliTRDrawStream()
256 delete [] fCurrTrkHeaderIndexWord;
257 delete [] fCurrTrkHeaderSize;
258 delete [] fCurrTrkFlags;
259 delete [] fCurrTrgHeaderIndexWord;
260 delete [] fCurrTrgHeaderSize;
261 delete [] fCurrTrgFlags;
262 delete [] fCurrStackIndexWord;
263 delete [] fCurrStackHeaderSize;
264 delete [] fCurrStackHeaderVersion;
265 delete [] fCurrLinkMask;
266 delete [] fCurrCleanCheckout;
267 delete [] fCurrBoardId;
268 delete [] fCurrHwRevTMU;
269 delete [] fCurrLinkMonitorFlags;
270 delete [] fCurrLinkDataTypeFlags;
271 delete [] fCurrLinkDebugFlags;
274 Bool_t AliTRDrawStream::ReadEvent(TTree *trackletTree)
276 // read the current event from the raw reader and fill it to the digits manager
279 AliError("No raw reader available");
284 ConnectTracklets(trackletTree);
289 // loop over all DDLs
290 // data starts with GTU payload, i.e. SM index word
291 UChar_t *buffer = 0x0;
293 while (fRawReader->ReadNextData(buffer)) {
295 fCurrEquipmentId = fRawReader->GetEquipmentId();
296 AliDebug(2, Form("equipment: %i", fCurrEquipmentId));
298 if (fCurrEquipmentId < kDDLOffset || fCurrEquipmentId > kDDLMax) {
299 EquipmentError(kNonTrdEq, "Skipping");
304 new ((*fMarkers)[fMarkers->GetEntriesFast()])
305 AliTRDrawStreamError(-kSecactive, fCurrEquipmentId - kDDLOffset);
307 ReadGTUHeaders((UInt_t*) buffer);
309 if (fCurrTrailerReadout)
312 // loop over all active links
313 AliDebug(2, Form("Stack mask 0x%02x", fCurrStackMask));
314 for (Int_t iStack = 0; iStack < fgkNstacks; iStack++) {
316 if ((fCurrStackMask & (1 << fCurrSlot)) == 0)
319 AliDebug(2, Form("Stack %i, Link mask: 0x%02x", fCurrSlot, fCurrLinkMask[fCurrSlot]));
320 for (Int_t iLink = 0; iLink < fgkNlinks; iLink++) {
322 fCurrHC = (fCurrEquipmentId - kDDLOffset) * fgkNstacks * fgkNlinks +
323 fCurrSlot * fgkNlinks + iLink;
324 if ((fCurrLinkMask[fCurrSlot] & (1 << fCurrLink)) == 0)
330 // check for link monitor error flag
331 if (fCurrLinkMonitorFlags[fCurrSlot*fgkNlinks + fCurrLink] != 0) {
332 LinkError(kLinkMonitor);
333 if (fgErrorBehav[kLinkMonitor] == kTolerate)
334 size += ReadLinkData();
337 // read the data from one HC
338 size += ReadLinkData();
340 // read all data endmarkers
341 size += SeekNextLink();
344 // continue with next stack
353 Bool_t AliTRDrawStream::NextDDL()
355 // continue reading with the next equipment
360 fCurrEquipmentId = 0;
364 UChar_t *buffer = 0x0;
366 while (fRawReader->ReadNextData(buffer)) {
368 fCurrEquipmentId = fRawReader->GetEquipmentId();
369 AliDebug(2, Form("equipment: %i", fCurrEquipmentId));
371 if (fCurrEquipmentId < kDDLOffset || fCurrEquipmentId > kDDLMax) {
372 EquipmentError(kNonTrdEq, "Skipping");
377 new ((*fMarkers)[fMarkers->GetEntriesFast()])
378 AliTRDrawStreamError(-kSecactive, fCurrEquipmentId - kDDLOffset);
380 ReadGTUHeaders((UInt_t*) buffer);
382 if (fCurrTrailerReadout)
392 Int_t AliTRDrawStream::NextChamber(AliTRDdigitsManager *digMgr)
394 // read the data for the next chamber
395 // in case you only want to read the data of a single chamber
396 // to read all data ReadEvent(...) is recommended
398 fDigitsManager = digMgr;
403 // tracklet output preparation
404 TTree *trklTree = 0x0;
405 AliRunLoader *rl = AliRunLoader::Instance();
406 AliLoader* trdLoader = rl ? rl->GetLoader("TRDLoader") : NULL;
407 AliDataLoader *trklLoader = trdLoader ? trdLoader->GetDataLoader("tracklets") : NULL;
409 AliTreeLoader *trklTreeLoader = (AliTreeLoader*) trklLoader->GetBaseLoader("tracklets-raw");
411 trklTree = trklTreeLoader->Tree();
413 trklTree = trklLoader->Tree();
416 if (fTrackletTree != trklTree)
417 ConnectTracklets(trklTree);
420 AliError("No raw reader available");
424 while (fCurrSlot < 0 || fCurrSlot >= fgkNstacks) {
429 while ((fCurrSlot < fgkNstacks) &&
430 (((fCurrStackMask & (1 << fCurrSlot)) == 0) ||
431 ((fCurrLinkMask[fCurrSlot] & (1 << fCurrLink))) == 0)) {
433 if (fCurrLink >= fgkNlinks) {
440 AliDebug(2, Form("Stack %i, Link %i, mask: 0x%02x", fCurrSlot, fCurrLink, fCurrLinkMask[fCurrSlot]));
441 fCurrHC = (fCurrEquipmentId - kDDLOffset) * fgkNlinks * fgkNstacks +
442 fCurrSlot * fgkNlinks + fCurrLink;
444 if (fCurrLinkMonitorFlags[fCurrSlot*fgkNlinks + fCurrLink] != 0) {
445 LinkError(kLinkMonitor);
446 if (fgErrorBehav[kLinkMonitor] == kTolerate)
450 // read the data from one HC
453 // read all data endmarkers
456 if (fCurrLink % 2 == 0) {
457 // if we just read the A-side HC then also check the B-side
460 if (fCurrLinkMask[fCurrSlot] & (1 << fCurrLink)) {
461 if (fCurrLinkMonitorFlags[fCurrSlot*fgkNlinks + fCurrLink] != 0) {
462 LinkError(kLinkMonitor);
463 if (fgErrorBehav[kLinkMonitor] == kTolerate)
474 if ((fCurrStackMask & (1 << fCurrSlot)) == 0) {
480 if (fCurrLink >= fgkNlinks) {
486 } while ((fCurrSlot < fgkNstacks) &&
487 (((fCurrStackMask & (1 << fCurrSlot)) == 0) ||
488 ((fCurrLinkMask[fCurrSlot] & (1 << fCurrLink))) == 0));
490 // return chamber information from HC if it is valid
491 // otherwise return information from link position
492 if (fCurrSm < 0 || fCurrSm >= fgkNsectors || fCurrStack < 0 || fCurrStack >= fgkNstacks || fCurrLayer < 0 || fCurrLayer >= fgkNlinks/2)
493 return ((fCurrEquipmentId-kDDLOffset) + fCurrSlot * fgkNlinks/2 + fCurrLink/2);
495 return (fCurrSm * fgkNstacks*fgkNlinks/2 + fCurrStack * fgkNlinks/2 + fCurrLayer);
499 Int_t AliTRDrawStream::ReadGTUHeaders(UInt_t *buffer)
501 // check the data source and read the headers
503 if (fCurrEquipmentId >= kDDLOffset && fCurrEquipmentId <= kDDLMax) {
506 // setting the pointer to data and current reading position
507 fPayloadCurr = fPayloadStart = buffer;
508 fPayloadSize = fRawReader->GetDataSize() / sizeof(UInt_t);
509 fStats.fStatsSector[fCurrEquipmentId - kDDLOffset].fBytes = fRawReader->GetDataSize();
510 AliDebug(2, Form("Read buffer of size: %i", fRawReader->GetDataSize()));
512 AliDebug(1, DumpRaw("raw data", fPayloadCurr, TMath::Min(fPayloadSize, 1000)));
515 if (ReadSmHeader() < 0) {
516 AliError(Form("Reading SM header failed, skipping this DDL %i", fCurrEquipmentId));
520 // read tracking headers (if available)
521 if (fCurrTrkHeaderAvail) {
522 for (Int_t iStack = 0; iStack < fgkNstacks; iStack++) {
523 if ((fCurrStackMask & (1 << iStack)) != 0)
524 ReadTrackingHeader(iStack);
528 // read trigger header(s) (if available)
529 if (fCurrTrgHeaderAvail)
530 ReadTriggerHeaders();
533 for (Int_t iStack = 0; iStack < fgkNstacks; iStack++) {
534 if ((fCurrStackMask & (1 << iStack)) != 0)
535 ReadStackHeader(iStack);
544 Int_t AliTRDrawStream::ReadSmHeader()
546 // read the SMU index header at the current reading position
547 // and store the information in the corresponding variables
549 if (fPayloadCurr - fPayloadStart >= fPayloadSize - 1) {
550 EquipmentError(kUnknown, "SM Header incomplete");
554 fCurrTrgFlags[fCurrEquipmentId-kDDLOffset] = 0;
556 fCurrSmHeaderSize = ((*fPayloadCurr) >> 16) & 0xffff;
557 fCurrSmHeaderVersion = ((*fPayloadCurr) >> 12) & 0xf;
558 fCurrTrackEnable = ((*fPayloadCurr) >> 6) & 0x1;
559 fCurrTrackletEnable = ((*fPayloadCurr) >> 5) & 0x1;
560 fCurrStackMask = ((*fPayloadCurr) ) & 0x1f;
561 fCurrHwRev = (fPayloadCurr[1] >> 12) & 0xffff;
562 fCurrStackEndmarkerAvail = 0;
564 switch (fCurrSmHeaderVersion) {
566 fCurrTrailerReadout = 0;
567 fCurrTrgHeaderAvail = 0;
569 fCurrTrkHeaderAvail = 0;
575 fCurrTrailerReadout = ((*fPayloadCurr) >> 10) & 0x1;
576 fCurrTrgHeaderAvail = 1;
577 fCurrTrgHeaderReadout = ((*fPayloadCurr) >> 9) & 0x1;
578 fCurrEvType = ((*fPayloadCurr) >> 7) & 0x3;
579 fCurrTrkHeaderAvail = fCurrTrackEnable;
580 fCurrTriggerEnable = (fPayloadCurr[2] >> 8) & 0xfff;
581 fCurrTriggerFired = (fPayloadCurr[2] >> 20) & 0xfff;
582 fCurrTrgFlags[fCurrEquipmentId-kDDLOffset] = fCurrTriggerFired;
586 fCurrTrailerReadout = ((*fPayloadCurr) >> 10) & 0x1;
587 fCurrTrgHeaderAvail = 1;
588 fCurrTrgHeaderReadout = ((*fPayloadCurr) >> 9) & 0x1;
589 fCurrEvType = ((*fPayloadCurr) >> 7) & 0x3;
590 fCurrTrkHeaderAvail = fCurrTrackEnable;
591 fCurrTriggerEnable = (fPayloadCurr[2] >> 8) & 0xfff;
592 fCurrTriggerFired = (fPayloadCurr[2] >> 20) & 0xfff;
593 fCurrTrgFlags[fCurrEquipmentId-kDDLOffset] = fCurrTriggerFired;
594 fCurrStackEndmarkerAvail = 1;
598 AliError(Form("unknown SM header version: 0x%x", fCurrSmHeaderVersion));
601 AliDebug(5, Form("SM header: size: %i, version: %i, track enable: %i, tracklet enable: %i, stack mask: %2x, trailer: %i, trgheader: %i, trkheader: %i",
603 fCurrSmHeaderVersion,
609 fCurrTrkHeaderAvail ));
611 // jump to the first word after the SM header
612 fPayloadCurr += fCurrSmHeaderSize + 1;
614 return fCurrSmHeaderSize + 1;
617 Int_t AliTRDrawStream::DecodeGTUtracks()
619 // decode GTU track words
620 // this depends on the hardware revision of the SMU
622 Int_t sector = fCurrEquipmentId-kDDLOffset;
624 if ((sector < 0) || (sector > 17)) {
625 AliError(Form("Invalid sector %i for GTU tracks", sector));
629 AliDebug(1, DumpRaw(Form("GTU tracks in sector %2i (hw rev %i)", sector, fCurrHwRev),
630 fPayloadCurr + 4, 10, 0xffe0ffff));
632 fCurrTrgFlags[sector] = 0;
634 if (fCurrHwRev < 1772) {
635 UInt_t fastWord; // fast trigger word
636 ULong64_t trackWord = 0; // extended track word
639 for (UInt_t iWord = 4; iWord < fCurrSmHeaderSize; iWord++) {
640 if (fPayloadCurr[iWord] == 0x10000000) { // stack boundary marker
646 ((fPayloadCurr[iWord] & 0xfffff0f0) == 0x13370000)) {
647 fastWord = fPayloadCurr[iWord];
648 fCurrTrgFlags[sector] |= 1 << (stack+11); // assume tracking done if fast word sent
649 AliDebug(1, Form("stack %i: fast trigger word: 0x%08x", stack, fastWord));
652 else if ((idx & 0x1) == 0x1) {
653 trackWord |= ((ULong64_t) fPayloadCurr[iWord]) << 32;
654 AliDebug(1,Form("track debug word: 0x%016llx", trackWord));
656 AliESDTrdTrack *trk = new ((*fTracks)[fTracks->GetEntriesFast()])
659 trk->SetSector(sector);
660 trk->SetStack((trackWord >> 60) & 0x7);
664 trk->SetLayerMask((trackWord >> 16) & 0x3f);
665 trk->SetTrackletIndex((trackWord >> 22) & 0x3f, 0);
666 trk->SetTrackletIndex((trackWord >> 28) & 0x3f, 1);
667 trk->SetTrackletIndex((trackWord >> 34) & 0x3f, 2);
668 trk->SetTrackletIndex((trackWord >> 40) & 0x3f, 3);
669 trk->SetTrackletIndex((trackWord >> 46) & 0x3f, 4);
670 trk->SetTrackletIndex((trackWord >> 52) & 0x3f, 5);
676 Float_t pt = (((Int_t) (trackWord & 0xffff) ^ 0x8000) - 0x8000)/128.;
677 if (TMath::Abs(pt) > 0.1) {
678 trk->SetA((Int_t) (0.15*51625./100./pt / 160e-4 * 2));
683 trackWord = fPayloadCurr[iWord];
689 else if (fCurrHwRev < 1804) {
690 UInt_t fastWord; // fast trigger word
691 ULong64_t trackWord = 0; // extended track word
694 for (UInt_t iWord = 4; iWord < fCurrSmHeaderSize; iWord++) {
695 if (fPayloadCurr[iWord] == 0xffe0ffff) { // stack boundary marker
701 ((fPayloadCurr[iWord] & 0xfffff0f0) == 0x13370000)) {
702 fastWord = fPayloadCurr[iWord];
703 fCurrTrgFlags[sector] |= 1 << (stack+11); // assume tracking done if fast word sent
704 AliDebug(1, Form("stack %i: fast trigger word: 0x%08x", stack, fastWord));
707 else if ((idx & 0x1) == 0x1) {
708 trackWord |= ((ULong64_t) fPayloadCurr[iWord]) << 32;
709 AliDebug(1, Form("track debug word: 0x%016llx", trackWord));
711 AliESDTrdTrack *trk = new ((*fTracks)[fTracks->GetEntriesFast()])
714 trk->SetSector(fCurrEquipmentId-kDDLOffset);
715 trk->SetStack((trackWord >> 60) & 0x7);
719 trk->SetLayerMask((trackWord >> 16) & 0x3f);
720 trk->SetTrackletIndex((trackWord >> 22) & 0x3f, 0);
721 trk->SetTrackletIndex((trackWord >> 28) & 0x3f, 1);
722 trk->SetTrackletIndex((trackWord >> 34) & 0x3f, 2);
723 trk->SetTrackletIndex((trackWord >> 40) & 0x3f, 3);
724 trk->SetTrackletIndex((trackWord >> 46) & 0x3f, 4);
725 trk->SetTrackletIndex((trackWord >> 52) & 0x3f, 5);
731 Float_t pt = (((Int_t) (trackWord & 0xffff) ^ 0x8000) - 0x8000)/128.;
732 if (TMath::Abs(pt) > 0.1) {
733 trk->SetA((Int_t) (-0.15*51625./100./pt / 160e-4 * 2));
738 trackWord = fPayloadCurr[iWord];
744 else if (fCurrHwRev < 1819) {
745 UInt_t fastWord; // fast trigger word
746 ULong64_t trackWord = 0; // extended track word
749 for (UInt_t iWord = 4; iWord < fCurrSmHeaderSize; iWord++) {
750 if (fPayloadCurr[iWord] == 0xffe0ffff) { // stack boundary marker
756 ((fPayloadCurr[iWord] & 0xfffff0f0) == 0x13370000)) {
757 fastWord = fPayloadCurr[iWord];
758 if (fastWord & (1 << 13))
759 fCurrTrgFlags[sector] |= 1 << (stack+11);
760 AliDebug(1, Form("stack %i: fast trigger word: 0x%08x", stack, fastWord));
763 else if ((idx & 0x1) == 0x1) {
764 trackWord |= ((ULong64_t) fPayloadCurr[iWord]) << 32;
765 AliDebug(1, Form("track debug word: 0x%016llx", trackWord));
768 AliESDTrdTrack *trk = new ((*fTracks)[fTracks->GetEntriesFast()])
771 trk->SetSector(fCurrEquipmentId-kDDLOffset);
772 trk->SetStack((trackWord >> 60) & 0x7);
775 // trk->SetPt(((trackWord & 0xffff) ^ 0x8000) - 0x8000);
777 trk->SetLayerMask((trackWord >> 16) & 0x3f);
778 trk->SetTrackletIndex((trackWord >> 22) & 0x3f, 0);
779 trk->SetTrackletIndex((trackWord >> 28) & 0x3f, 1);
780 trk->SetTrackletIndex((trackWord >> 34) & 0x3f, 2);
781 trk->SetTrackletIndex((trackWord >> 40) & 0x3f, 3);
782 trk->SetTrackletIndex((trackWord >> 46) & 0x3f, 4);
783 trk->SetTrackletIndex((trackWord >> 52) & 0x3f, 5);
789 Float_t pt = (((Int_t) (trackWord & 0xffff) ^ 0x8000) - 0x8000)/128.;
790 if (TMath::Abs(pt) > 0.1) {
791 trk->SetA((Int_t) (0.15*51625./100./trk->Pt() / 160e-4 * 2));
796 trackWord = fPayloadCurr[iWord];
802 else if (fCurrHwRev < 1860) {
803 UInt_t fastWord; // fast trigger word
804 ULong64_t trackWord = 0; // extended track word
807 Bool_t upperWord = kFALSE;
809 for (UInt_t iWord = 4; iWord < fCurrSmHeaderSize; iWord++) {
810 if (fPayloadCurr[iWord] == 0xffe0ffff) { // stack boundary marker
816 // assemble the 32-bit words out of 16-bit blocks
818 word |= (fPayloadCurr[iWord] & 0xffff0000);
822 // lower word is read first
823 word = (fPayloadCurr[iWord] & 0xffff0000) >> 16;
828 if ((word & 0xffff0008) == 0x13370008) {
830 AliDebug(1, Form("stack %i: fast track word: 0x%08x", stack, fastWord));
831 if (fastWord & (1 << 13))
832 fCurrTrgFlags[sector] |= 1 << (stack+11);
835 else if ((idx & 0x1) == 0x1) {
836 trackWord |= ((ULong64_t) word) << 32;
837 AliDebug(1, Form("track debug word: 0x%016llx", trackWord));
839 AliESDTrdTrack *trk = new ((*fTracks)[fTracks->GetEntriesFast()])
842 trk->SetSector(fCurrEquipmentId-kDDLOffset);
843 trk->SetStack((trackWord >> 60) & 0x7);
847 trk->SetLayerMask((trackWord >> 16) & 0x3f);
848 trk->SetTrackletIndex((trackWord >> 22) & 0x3f, 0);
849 trk->SetTrackletIndex((trackWord >> 28) & 0x3f, 1);
850 trk->SetTrackletIndex((trackWord >> 34) & 0x3f, 2);
851 trk->SetTrackletIndex((trackWord >> 40) & 0x3f, 3);
852 trk->SetTrackletIndex((trackWord >> 46) & 0x3f, 4);
853 trk->SetTrackletIndex((trackWord >> 52) & 0x3f, 5);
859 Float_t pt = (((Int_t) (trackWord & 0xffff) ^ 0x8000) - 0x8000)/128.;
860 if (TMath::Abs(pt) > 0.1) {
861 trk->SetA((Int_t) (0.15*51625./100./pt / 160e-4 * 2));
874 ULong64_t trackWord = 0; // this is the debug word
877 Bool_t upperWord = kFALSE;
879 for (UInt_t iWord = 4; iWord < fCurrSmHeaderSize; iWord++) {
880 if (fPayloadCurr[iWord] == 0xffe0ffff) {
886 // assemble the 32-bit words out of 16-bit blocks
888 word |= (fPayloadCurr[iWord] & 0xffff0000);
892 // lower word is read first
893 word = (fPayloadCurr[iWord] & 0xffff0000) >> 16;
898 if ((word & 0xffff0008) == 0x13370008) {
899 AliDebug(1, Form("stack %i: fast track word: 0x%08x", stack, word));
902 else if ((word & 0xffff0010) == 0x13370010) {
903 AliDebug(1, Form("stack %i: tracking done word: 0x%08x", stack, word));
904 fCurrTrgFlags[sector] |= 1 << (stack+11);
907 else if ((idx & 0x1) == 0x1) {
908 trackWord |= ((ULong64_t) word) << 32;
909 AliDebug(1, Form("track debug word: 0x%16llx", trackWord));
911 AliESDTrdTrack *trk = new ((*fTracks)[fTracks->GetEntriesFast()])
913 trk->SetSector(fCurrEquipmentId-kDDLOffset);
914 trk->SetStack((trackWord >> 60) & 0x7);
918 trk->SetLayerMask((trackWord >> 16) & 0x3f);
919 trk->SetTrackletIndex((trackWord >> 22) & 0x3f, 0);
920 trk->SetTrackletIndex((trackWord >> 28) & 0x3f, 1);
921 trk->SetTrackletIndex((trackWord >> 34) & 0x3f, 2);
922 trk->SetTrackletIndex((trackWord >> 40) & 0x3f, 3);
923 trk->SetTrackletIndex((trackWord >> 46) & 0x3f, 4);
924 trk->SetTrackletIndex((trackWord >> 52) & 0x3f, 5);
930 Float_t pt = (((Int_t) (trackWord & 0xffff) ^ 0x8000) - 0x8000)/128.;
931 if (TMath::Abs(pt) > 0.1) {
932 trk->SetA(-(Int_t) (0.15*51625./100./pt / 160e-4 * 2));
946 Int_t AliTRDrawStream::ReadTrackingHeader(Int_t stack)
948 // read the tracking information and store it for the given stack
952 fCurrTrkHeaderIndexWord[stack] = *fPayloadCurr;
953 fCurrTrkHeaderSize[stack] = ((*fPayloadCurr) >> 16) & 0x3ff;
955 AliDebug(1, Form("tracking header index word: 0x%08x, size: %i (hw rev: %i)",
956 fCurrTrkHeaderIndexWord[stack], fCurrTrkHeaderSize[stack], fCurrHwRev));
957 Int_t trackingTime = *fPayloadCurr & 0x3ff;
959 fCurrTrgFlags[fCurrEquipmentId-kDDLOffset] |= ((fCurrTrkHeaderIndexWord[stack] >> 10) & 0x1) << (22 + stack);
963 ULong64_t trackWord = 0;
965 Int_t trackIndex = fTracks ? fTracks->GetEntriesFast() : -1;
967 for (UInt_t iWord = 0; iWord < fCurrTrkHeaderSize[stack]; iWord++) {
970 // first part of 64-bit word
971 trackWord = fPayloadCurr[iWord];
974 trackWord |= ((ULong64_t) fPayloadCurr[iWord]) << 32;
976 if (trackWord & (1ul << 63)) {
977 if ((trackWord & (0x3ful << 56)) != 0) {
979 AliDebug(2, Form("track word: 0x%016llx", trackWord));
982 AliESDTrdTrack *trk = new ((*fTracks)[fTracks->GetEntriesFast()])
985 trk->SetSector(fCurrEquipmentId-kDDLOffset);
986 trk->SetLayerMask((trackWord >> 56) & 0x3f);
987 trk->SetA( (((trackWord >> 38) & 0x3ffff) ^ 0x20000) - 0x20000);
988 trk->SetB( (((trackWord >> 20) & 0x3ffff) ^ 0x20000) - 0x20000);
989 trk->SetC( (((trackWord >> 8) & 0xffff) ^ 0x8000) - 0x8000);
990 trk->SetPID((trackWord >> 0) & 0xff);
991 trk->SetStack(stack);
994 // now compare the track word with the one generated from the ESD information
995 if (trackWord != trk->GetTrackWord(0)) {
996 AliError(Form("track word 0x%016llx does not match the read one 0x%016llx",
997 trk->GetTrackWord(0), trackWord));
1002 // done marker (so far only used to set trigger flag)
1003 fCurrTrgFlags[fCurrEquipmentId-kDDLOffset] |= 1 << (27 + stack);
1004 fCurrTrkFlags[(fCurrEquipmentId-kDDLOffset)*fgkNstacks + stack] = trackWord;
1006 AliDebug(2, Form("tracking done marker: 0x%016llx, trigger flags: 0x%08x",
1007 trackWord, fCurrTrgFlags[fCurrEquipmentId-kDDLOffset]));
1008 AliDebug(2, Form("seg / stack / first / last / done / index : %i %i %lli %lli %lli %i",
1009 fCurrEquipmentId - kDDLOffset, stack,
1010 (trackWord >> 20) & 0x3ff,
1011 (trackWord >> 10) & 0x3ff,
1012 (trackWord >> 0) & 0x3ff,
1017 // extended track word
1018 AliDebug(2, Form("extended track word: 0x%016llx", trackWord));
1021 AliESDTrdTrack *trk = (AliESDTrdTrack*) (*fTracks)[trackIndex];
1023 trk->SetFlags((trackWord >> 52) & 0x7ff);
1024 trk->SetReserved((trackWord >> 49) & 0x7);
1025 trk->SetY((trackWord >> 36) & 0x1fff);
1026 trk->SetTrackletIndex((trackWord >> 0) & 0x3f, 0);
1027 trk->SetTrackletIndex((trackWord >> 6) & 0x3f, 1);
1028 trk->SetTrackletIndex((trackWord >> 12) & 0x3f, 2);
1029 trk->SetTrackletIndex((trackWord >> 18) & 0x3f, 3);
1030 trk->SetTrackletIndex((trackWord >> 24) & 0x3f, 4);
1031 trk->SetTrackletIndex((trackWord >> 30) & 0x3f, 5);
1033 if (trackWord != trk->GetExtendedTrackWord(0)) {
1034 AliError(Form("extended track word 0x%016llx does not match the read one 0x%016llx",
1035 trk->GetExtendedTrackWord(0), trackWord));
1045 fPayloadCurr += fCurrTrkHeaderSize[stack];
1047 return fCurrTrkHeaderSize[stack];
1050 Int_t AliTRDrawStream::ReadTriggerHeaders()
1052 // read all trigger headers present
1054 AliDebug(1, Form("trigger mask: 0x%03x, fired: 0x%03x\n",
1055 fCurrTriggerEnable, fCurrTriggerFired));
1056 // loop over potential trigger blocks
1057 for (Int_t iTrigger = 0; iTrigger < fgkNtriggers; iTrigger++) {
1058 // check for trigger enable
1059 if (fCurrTriggerEnable & (1 << iTrigger)) {
1060 // check for readout mode and trigger fired
1061 if ((fCurrTrgHeaderReadout == 0) || (fCurrTriggerFired & (1 << iTrigger))) {
1063 AliDebug(1, Form("trigger index word %i: 0x%08x\n", iTrigger, *fPayloadCurr));
1064 fCurrTrgHeaderIndexWord[iTrigger] = *fPayloadCurr;
1065 fCurrTrgHeaderSize[iTrigger] = ((*fPayloadCurr) >> 16) & 0x3ff;
1066 if (iTrigger == 7) {
1067 // timeout trigger, use to extract tracking time
1068 fCurrTrgFlags[fCurrEquipmentId-kDDLOffset] |= (*fPayloadCurr & 0x3ff) << 12;
1073 fPayloadCurr += fCurrTrgHeaderSize[iTrigger];
1081 Int_t AliTRDrawStream::ReadStackHeader(Int_t stack)
1083 // read the stack header
1084 // and store the information in the corresponding variables
1086 fCurrStackIndexWord[stack] = *fPayloadCurr;
1087 fCurrStackHeaderSize[stack] = (((*fPayloadCurr) >> 16) & 0xffff) + 1;
1088 fCurrStackHeaderVersion[stack] = ((*fPayloadCurr) >> 12) & 0xf;
1089 fCurrLinkMask[stack] = (*fPayloadCurr) & 0xfff;
1091 // dumping stack header
1092 AliDebug(1, DumpRaw(Form("stack %i header", stack), fPayloadCurr, fCurrStackHeaderSize[stack]));
1094 if (fPayloadCurr - fPayloadStart >= fPayloadSize - (Int_t) fCurrStackHeaderSize[stack]) {
1095 EquipmentError(kStackHeaderInvalid, "Stack index header %i incomplete", stack);
1096 // dumping stack header
1097 AliError(DumpRaw(Form("stack %i header", stack), fPayloadCurr, fCurrStackHeaderSize[stack]));
1102 switch (fCurrStackHeaderVersion[stack]) {
1104 if (fCurrStackHeaderSize[stack] < 8) {
1105 LinkError(kStackHeaderInvalid, "Stack header smaller than expected!");
1109 fCurrCleanCheckout[stack] = fPayloadCurr[1] & 0x1;
1110 fCurrBoardId[stack] = (fPayloadCurr[1] >> 8) & 0xff;
1111 fCurrHwRevTMU[stack] = (fPayloadCurr[1] >> 16) & 0xffff;
1113 for (Int_t iLayer = 0; iLayer < 6; iLayer++) {
1115 fCurrLinkMonitorFlags [stack*fgkNlinks + iLayer*2] = fPayloadCurr[iLayer+2] & 0xf;
1116 fCurrLinkDataTypeFlags [stack*fgkNlinks + iLayer*2] = (fPayloadCurr[iLayer+2] >> 4) & 0x3;
1117 fCurrLinkDebugFlags [stack*fgkNlinks + iLayer*2] = (fPayloadCurr[iLayer+2] >> 12) & 0xf;
1119 fCurrLinkMonitorFlags [stack*fgkNlinks + iLayer*2 + 1] = (fPayloadCurr[iLayer+2] >> 16) & 0xf;
1120 fCurrLinkDataTypeFlags [stack*fgkNlinks + iLayer*2 + 1] = (fPayloadCurr[iLayer+2] >> 20) & 0x3;
1121 fCurrLinkDebugFlags [stack*fgkNlinks + iLayer*2 + 1] = (fPayloadCurr[iLayer+2] >> 28) & 0xf;
1126 LinkError(kStackHeaderInvalid, "Invalid Stack Header version %x", fCurrStackHeaderVersion[stack]);
1129 fPayloadCurr += fCurrStackHeaderSize[stack];
1131 return fCurrStackHeaderSize[stack];
1134 Int_t AliTRDrawStream::ReadGTUTrailer()
1136 // read the SM trailer containing CRCs from various stages
1138 UInt_t* trailer = fPayloadStart + fPayloadSize -1;
1140 // look for the trailer index word from the end
1141 for (Int_t iWord = 0; iWord < fPayloadSize; iWord++) {
1142 if ((fPayloadStart[fPayloadSize-1-iWord] & 0xffff) == 0x1f51) {
1143 trailer = fPayloadStart + fPayloadSize - 1 - iWord;
1148 if (((*trailer) & 0xffff) == 0x1f51) {
1149 UInt_t trailerIndexWord = (*trailer);
1150 Int_t trailerSize = (trailerIndexWord >> 16) & 0xffff;
1151 AliDebug(2, DumpRaw("GTU trailer", trailer, trailerSize+1));
1152 // parse the trailer
1153 if (trailerSize >= 4) {
1154 // match flags from GTU
1155 fCurrMatchFlagsSRAM = (trailer[1] >> 0) & 0x1f;
1156 fCurrMatchFlagsPostBP = (trailer[1] >> 5) & 0x1f;
1157 // individual checksums
1158 fCurrChecksumStack[0] = (trailer[1] >> 16) & 0xffff;
1159 fCurrChecksumStack[1] = (trailer[2] >> 0) & 0xffff;
1160 fCurrChecksumStack[2] = (trailer[2] >> 16) & 0xffff;
1161 fCurrChecksumStack[3] = (trailer[3] >> 0) & 0xffff;
1162 fCurrChecksumStack[4] = (trailer[3] >> 16) & 0xffff;
1163 fCurrChecksumSIU = trailer[4];
1165 if ((fCurrMatchFlagsSRAM & fCurrStackMask) != fCurrStackMask)
1166 EquipmentError(kCRCmismatch, "CRC mismatch SRAM: 0x%02x", fCurrMatchFlagsSRAM);
1167 if ((fCurrMatchFlagsPostBP & fCurrStackMask) != fCurrStackMask)
1168 EquipmentError(kCRCmismatch, "CRC mismatch BP: 0x%02x", fCurrMatchFlagsPostBP);
1172 LinkError(kUnknown, "Invalid GTU trailer");
1176 EquipmentError(kUnknown, "trailer index marker mismatch");
1181 Int_t AliTRDrawStream::ReadLinkData()
1183 // read the data in one link (one HC) until the data endmarker is reached
1184 // returns the number of words read!
1187 UInt_t* startPosLink = fPayloadCurr;
1189 AliDebug(1, DumpRaw(Form("link data from seg %2i slot %i link %2i", fCurrEquipmentId-kDDLOffset, fCurrSlot, fCurrLink),
1190 fPayloadCurr, TMath::Min((Int_t) (fPayloadSize - (fPayloadCurr-fPayloadStart)), 100), 0x00000000));
1193 new ((*fMarkers)[fMarkers->GetEntriesFast()])
1194 AliTRDrawStreamError(-kHCactive, fCurrEquipmentId-kDDLOffset, fCurrSlot, fCurrLink);
1196 if (fErrorFlags & kDiscardHC)
1199 if (fCurrTrackletEnable) {
1200 count += ReadTracklets();
1201 if (fErrorFlags & kDiscardHC)
1205 AliDebug(1, DumpRaw("HC header", fPayloadCurr, 4, 0x00000000));
1206 count += ReadHcHeader();
1207 if (fErrorFlags & kDiscardHC)
1210 Int_t det = fCurrSm * 30 + fCurrStack * 6 + fCurrLayer;
1212 if (det > -1 && det < 540) {
1214 // ----- check which kind of data -----
1215 if (fCurrMajor & 0x40) {
1216 if ((fCurrMajor & 0x7) == 0x7) {
1217 AliDebug(1, "This is a config event");
1218 UInt_t *startPos = fPayloadCurr;
1219 while (fPayloadCurr - fPayloadStart < fPayloadSize &&
1220 *fPayloadCurr != fgkDataEndmarker)
1222 count += fPayloadCurr - startPos;
1224 // feeding TRAP config
1225 AliTRDtrapConfig *trapcfg = AliTRDcalibDB::Instance()->GetTrapConfig();
1226 AliTRDmcmSim::ReadPackedConfig(trapcfg, fCurrHC, startPos, fPayloadCurr - startPos);
1229 Int_t tpmode = fCurrMajor & 0x7;
1230 AliDebug(1, Form("Checking testpattern (mode %i) data", tpmode));
1231 count += ReadTPData(tpmode);
1235 // reading real data
1236 if (fDigitsManager) {
1237 if ((fAdcArray = fDigitsManager->GetDigits(det))) {
1238 //fAdcArray->Expand();
1239 if (fAdcArray->GetNtime() != fCurrNtimebins)
1240 fAdcArray->Allocate(16, 144, fCurrNtimebins);
1243 LinkError(kNoDigits);
1246 if (!fDigitsParam) {
1247 fDigitsParam = fDigitsManager->GetDigitsParam();
1250 fDigitsParam->SetPretriggerPhase(det, fCurrPtrgPhase);
1251 fDigitsParam->SetNTimeBins(det, fCurrNtimebins);
1252 fDigitsParam->SetADCbaseline(det, 10);
1255 if (fDigitsManager->UsesDictionaries()) {
1256 fDigitsManager->GetDictionary(det, 0)->Reset();
1257 fDigitsManager->GetDictionary(det, 1)->Reset();
1258 fDigitsManager->GetDictionary(det, 2)->Reset();
1261 if ((fSignalIndex = fDigitsManager->GetIndexes(det))) {
1262 fSignalIndex->SetSM(fCurrSm);
1263 fSignalIndex->SetStack(fCurrStack);
1264 fSignalIndex->SetLayer(fCurrLayer);
1265 fSignalIndex->SetDetNumber(det);
1266 if (!fSignalIndex->IsAllocated())
1267 fSignalIndex->Allocate(16, 144, fCurrNtimebins);
1270 if (fCurrMajor & 0x20) {
1271 AliDebug(1, "This is a zs event");
1272 count += ReadZSData();
1275 AliDebug(1, "This is a nozs event");
1276 count += ReadNonZSData();
1280 // just read until data endmarkers
1281 while (fPayloadCurr - fPayloadStart < fPayloadSize &&
1282 *fPayloadCurr != fgkDataEndmarker)
1288 LinkError(kInvalidDetector, "%i", det);
1289 while (fPayloadCurr - fPayloadStart < fPayloadSize &&
1290 *fPayloadCurr != fgkDataEndmarker)
1294 if (fCurrSm > -1 && fCurrSm < 18) {
1295 fStats.fStatsSector[fCurrSm].fStatsHC[fCurrHC%60].fBytes += (fPayloadCurr - startPosLink) * sizeof(UInt_t);
1296 fStats.fStatsSector[fCurrSm].fStatsHC[fCurrHC%60].fBytesRead += count * sizeof(UInt_t);
1297 fStats.fStatsSector[fCurrSm].fBytesRead += count * sizeof(UInt_t);
1298 fStats.fBytesRead += count * sizeof(UInt_t);
1304 Int_t AliTRDrawStream::ReadTracklets()
1306 // read the tracklets from one HC
1308 fTrackletArray->Clear();
1310 UInt_t *start = fPayloadCurr;
1311 while (*(fPayloadCurr) != fgkTrackletEndmarker &&
1312 fPayloadCurr - fPayloadStart < fPayloadSize) {
1313 new ((*fTrackletArray)[fTrackletArray->GetEntriesFast()]) AliTRDtrackletWord(*(fPayloadCurr), fCurrHC);
1318 if (fTrackletArray->GetEntriesFast() > 0) {
1319 AliDebug(1, Form("Found %i tracklets in %i %i %i (ev. %i)", fTrackletArray->GetEntriesFast(),
1320 (fCurrEquipmentId-kDDLOffset), fCurrSlot, fCurrLink, fRawReader->GetEventIndex()));
1321 if (fCurrSm > -1 && fCurrSm < 18) {
1322 fStats.fStatsSector[fCurrSm].fStatsHC[fCurrHC%60].fNTracklets += fTrackletArray->GetEntriesFast();
1323 fStats.fStatsSector[fCurrSm].fNTracklets += fTrackletArray->GetEntriesFast();
1326 fTrackletTree->Fill();
1328 for (Int_t iTracklet = 0; iTracklet < fTrackletArray->GetEntriesFast(); iTracklet++) {
1329 new ((*fTracklets)[fTracklets->GetEntriesFast()]) AliTRDtrackletWord(*((AliTRDtrackletWord*)(*fTrackletArray)[iTracklet]));
1333 // loop over remaining tracklet endmarkers
1334 while ((*(fPayloadCurr) == fgkTrackletEndmarker &&
1335 fPayloadCurr - fPayloadStart < fPayloadSize))
1338 return fPayloadCurr - start;
1341 Int_t AliTRDrawStream::ReadHcHeader()
1343 // read and parse the HC header of one HC
1344 // and store the information in the corresponding variables
1346 AliDebug(1, Form("HC header: 0x%08x", *fPayloadCurr));
1347 UInt_t *start = fPayloadCurr;
1348 // check not to be at the data endmarker
1349 if (*fPayloadCurr == fgkDataEndmarker)
1352 fCurrSpecial = (*fPayloadCurr >> 31) & 0x1;
1353 fCurrMajor = (*fPayloadCurr >> 24) & 0x7f;
1354 fCurrMinor = (*fPayloadCurr >> 17) & 0x7f;
1355 fCurrAddHcWords = (*fPayloadCurr >> 14) & 0x7;
1356 fCurrSm = (*fPayloadCurr >> 9) & 0x1f;
1357 fCurrLayer = (*fPayloadCurr >> 6) & 0x7;
1358 fCurrStack = (*fPayloadCurr >> 3) & 0x7;
1359 fCurrSide = (*fPayloadCurr >> 2) & 0x1;
1360 fCurrCheck = (*fPayloadCurr) & 0x3;
1362 if ((fCurrSm != (((Int_t) fCurrEquipmentId) - kDDLOffset)) ||
1363 (fCurrStack != fCurrSlot) ||
1364 (fCurrLayer != fCurrLink / 2) ||
1365 (fCurrSide != fCurrLink % 2)) {
1366 LinkError(kHCmismatch,
1367 "HC: %i, %i, %i, %i\n 0x%08x 0x%08x 0x%08x 0x%08x",
1368 fCurrSm, fCurrStack, fCurrLayer, fCurrSide,
1369 fPayloadCurr[0], fPayloadCurr[1], fPayloadCurr[2], fPayloadCurr[3]);
1371 if (fCurrCheck != 0x1) {
1372 LinkError(kHCcheckFailed);
1375 if (fCurrAddHcWords > 0) {
1376 fCurrNtimebins = (fPayloadCurr[1] >> 26) & 0x3f;
1377 fCurrBC = (fPayloadCurr[1] >> 10) & 0xffff;
1378 fCurrPtrgCnt = (fPayloadCurr[1] >> 6) & 0xf;
1379 fCurrPtrgPhase = (fPayloadCurr[1] >> 2) & 0xf;
1382 fPayloadCurr += 1 + fCurrAddHcWords;
1384 return (fPayloadCurr - start);
1387 Int_t AliTRDrawStream::ReadTPData(Int_t mode)
1389 // testing of testpattern 1 to 3 (hardcoded), 0 missing
1390 // evcnt checking missing
1392 Int_t cpufromchannel[] = {0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 2, 2, 2, 2, 2, 3, 3, 3, 3, 3, 3};
1396 Int_t mcmcount = -1;
1397 Int_t wordcount = 0;
1398 Int_t channelcount = 0;
1400 UInt_t expadcval = 0;
1402 Int_t lastmcmpos = -1;
1403 Int_t lastrobpos = -1;
1405 UInt_t* start = fPayloadCurr;
1407 while (*(fPayloadCurr) != fgkDataEndmarker &&
1408 fPayloadCurr - fPayloadStart < fPayloadSize - 1) {
1410 // ----- Checking MCM Header -----
1411 AliDebug(2, DumpMcmHeader("MCM header: ", *fPayloadCurr));
1412 UInt_t *startPosMCM = fPayloadCurr;
1415 // ----- checking for proper readout order - ROB -----
1416 fCurrRobPos = ROB(*fPayloadCurr);
1417 if (GetROBReadoutPos(ROB(*fPayloadCurr) / 2) >= lastrobpos) {
1418 if (GetROBReadoutPos(ROB(*fPayloadCurr) / 2) > lastrobpos)
1420 lastrobpos = GetROBReadoutPos(ROB(*fPayloadCurr) / 2);
1423 ROBError(kPosUnexp, Form("#%i after #%i in readout order", GetROBReadoutPos(ROB(*fPayloadCurr) / 2), lastrobpos));
1426 // ----- checking for proper readout order - MCM -----
1427 fCurrMcmPos = MCM(*fPayloadCurr);
1428 if (GetMCMReadoutPos(MCM(*fPayloadCurr)) > lastmcmpos) {
1429 lastmcmpos = GetMCMReadoutPos(MCM(*fPayloadCurr));
1432 MCMError(kPosUnexp, Form("#%i after #%i in readout order", GetMCMReadoutPos(MCM(*fPayloadCurr)), lastmcmpos));
1435 if (EvNo(*fPayloadCurr) != evno) {
1437 evno = EvNo(*fPayloadCurr);
1440 MCMError(kEvCntMismatch, "%i <-> %i", evno, EvNo(*fPayloadCurr));
1446 evcnt = 0x3f & *fPayloadCurr >> 26;
1449 while (channelcount < 21) {
1451 if (cpu != cpufromchannel[channelcount]) {
1452 cpu = cpufromchannel[channelcount];
1453 expadcval = (1 << 9) | (fCurrRobPos << 6) | (fCurrMcmPos << 2) | cpu;
1457 while (count < 10) {
1458 if (*fPayloadCurr == fgkDataEndmarker) {
1459 MCMError(kMissTpData);
1460 return (fPayloadCurr - start);
1463 if (channelcount % 2 == 0)
1470 expword |= expadcval << 2;
1471 expadcval = ( (expadcval << 1) | ( ( (expadcval >> 9) ^ (expadcval >> 6) ) & 1) ) & 0x3FF;
1472 expword |= expadcval << 12;
1473 expadcval = ( (expadcval << 1) | ( ( (expadcval >> 9) ^ (expadcval >> 6) ) & 1) ) & 0x3FF;
1474 expword |= expadcval << 22;
1475 expadcval = ( (expadcval << 1) | ( ( (expadcval >> 9) ^ (expadcval >> 6) ) & 1) ) & 0x3FF;
1477 else if (mode == 2) {
1478 // ----- TP 2 ------
1479 expword = ((0x3f & evcnt) << 26) | ((fCurrSm + 1) << 21) | ((fCurrLayer + 1) << 18) |
1480 ((fCurrStack + 1) << 15) |
1481 (fCurrRobPos << 12) | (fCurrMcmPos << 8) | (cpu << 6) | (wordcount + 1);
1483 else if (mode == 3) {
1485 expword = ((0xfff & evcnt) << 20) | (fCurrSm << 15) | (fCurrLink/2 << 12) | (fCurrStack << 9) |
1486 (fCurrRobPos << 6) | (fCurrMcmPos << 2) | (cpu << 0);
1490 LinkError(kTPmodeInvalid, "Just reading");
1493 diff = *fPayloadCurr ^ expword;
1494 AliDebug(11, Form("Comparing ch %2i, word %2i (cpu %i): 0x%08x <-> 0x%08x",
1495 channelcount, wordcount, cpu, *fPayloadCurr, expword));
1498 MCMError(kTPmismatch,
1499 "Seen 0x%08x, expected 0x%08x, diff: 0x%08x, 0x%04x, 0x%02x - word %2i (cpu %i, ch %i)",
1500 *fPayloadCurr, expword, diff,
1501 0xffff & (diff | diff >> 16),
1502 0xff & (diff | diff >> 8 | diff >> 16 | diff >> 24),
1503 wordcount, cpu, channelcount);;
1508 if (*fPayloadCurr == fgkDataEndmarker)
1509 return (fPayloadCurr - start);
1513 // continue with next MCM
1515 if (IsDumping() && DumpingMCM(fCurrHC/2, fCurrRobPos, fCurrMcmPos)) {
1516 AliInfo(DumpRaw(Form("Event %i: Det %3i ROB %i MCM %2i", fRawReader->GetEventIndex(), fCurrHC/2, fCurrRobPos, fCurrMcmPos),
1517 startPosMCM, fPayloadCurr - startPosMCM));
1521 return fPayloadCurr - start;
1525 Int_t AliTRDrawStream::ReadZSData()
1527 // read the zs data from one link from the current reading position
1529 UInt_t *start = fPayloadCurr;
1532 Int_t mcmcountExp = fCurrStack == 2 ? 48 : 64;
1533 Int_t channelcount = 0;
1534 Int_t channelcountExp = 0;
1535 Int_t channelcountMax = 0;
1537 Int_t currentTimebin = 0;
1540 Int_t lastmcmpos = -1;
1541 Int_t lastrobpos = -1;
1543 if (fCurrNtimebins != fNtimebins) {
1545 LinkError(kNtimebinsChanged,
1546 "No. of timebins changed from %i to %i", fNtimebins, fCurrNtimebins);
1547 fNtimebins = fCurrNtimebins;
1550 timebins = fNtimebins;
1552 while (*(fPayloadCurr) != fgkDataEndmarker &&
1553 fPayloadCurr - fPayloadStart < fPayloadSize) {
1555 // ----- Checking MCM Header -----
1556 AliDebug(2, DumpMcmHeader("MCM header: ", *fPayloadCurr));
1557 UInt_t *startPosMCM = fPayloadCurr;
1559 // ----- checking for proper readout order - ROB -----
1560 fCurrRobPos = ROB(*fPayloadCurr);
1561 if (GetROBReadoutPos(ROB(*fPayloadCurr) / 2) >= lastrobpos) {
1562 if (GetROBReadoutPos(ROB(*fPayloadCurr) / 2) > lastrobpos)
1564 lastrobpos = GetROBReadoutPos(ROB(*fPayloadCurr) / 2);
1567 ROBError(kPosUnexp, Form("#%i after #%i and #%i in readout order",
1568 GetROBReadoutPos(ROB(*fPayloadCurr) / 2), lastrobpos, GetROBReadoutPos(fCurrRobPos)));
1571 // ----- checking for proper readout order - MCM -----
1572 fCurrMcmPos = MCM(*fPayloadCurr);
1573 if (GetMCMReadoutPos(MCM(*fPayloadCurr)) > lastmcmpos) {
1574 lastmcmpos = GetMCMReadoutPos(MCM(*fPayloadCurr));
1577 MCMError(kPosUnexp, Form("#%i after #%i and #%i in readout order",
1578 GetMCMReadoutPos(MCM(*fPayloadCurr)), lastmcmpos, GetMCMReadoutPos(fCurrMcmPos)));
1581 if (EvNo(*fPayloadCurr) != evno) {
1583 evno = EvNo(*fPayloadCurr);
1585 MCMError(kEvCntMismatch, "%i <-> %i", evno, EvNo(*fPayloadCurr));
1588 Int_t adccoloff = AdcColOffset(*fPayloadCurr);
1589 Int_t padcoloff = PadColOffset(*fPayloadCurr);
1590 Int_t row = Row(*fPayloadCurr);
1593 if ((row > 11) && (fCurrStack == 2)) {
1594 MCMError(kUnknown, "Data in padrow > 11 for stack 2");
1597 // ----- Reading ADC channels -----
1598 AliDebug(2, DumpAdcMask("ADC mask: ", *fPayloadCurr));
1600 // ----- analysing the ADC mask -----
1602 channelcountExp = GetNActiveChannelsFromMask(*fPayloadCurr);
1603 channelcountMax = GetNActiveChannels(*fPayloadCurr);
1604 Int_t channelmask = GetActiveChannels(*fPayloadCurr);
1605 Int_t channelno = -1;
1608 if (channelcountExp != channelcountMax) {
1609 if (channelcountExp > channelcountMax) {
1610 Int_t temp = channelcountExp;
1611 channelcountExp = channelcountMax;
1612 channelcountMax = temp;
1614 while (channelcountExp < channelcountMax && channelcountExp < 21 &&
1615 fPayloadCurr - fPayloadStart < fPayloadSize - 10 * channelcountExp - 1) {
1616 MCMError(kAdcMaskInconsistent,
1617 "Possible MCM-H: 0x%08x, possible ADC-mask: 0x%08x",
1618 *(fPayloadCurr + 10 * channelcountExp),
1619 *(fPayloadCurr + 10 * channelcountExp + 1) );
1620 if (!CouldBeMCMhdr( *(fPayloadCurr + 10 * channelcountExp)) && !CouldBeADCmask( *(fPayloadCurr + 10 * channelcountExp + 1)))
1626 MCMError(kAdcMaskInconsistent,
1627 "Inconsistency in no. of active channels: Counter: %i, Mask: %i, chosen: %i!",
1628 GetNActiveChannels(fPayloadCurr[-1]), GetNActiveChannelsFromMask(fPayloadCurr[-1]), channelcountExp);
1630 AliDebug(2, Form("expecting %i active channels, %i timebins", channelcountExp, fCurrNtimebins));
1632 // ----- reading marked ADC channels -----
1633 while (channelcount < channelcountExp && *(fPayloadCurr) != fgkDataEndmarker) {
1636 while (channelno < 20 && (channelmask & 1 << channelno) == 0)
1639 if (fCurrNtimebins > 30) {
1640 currentTimebin = ((*fPayloadCurr >> 2) & 0x3f);
1641 timebins = ((*fPayloadCurr >> 8) & 0xf) * 3;
1648 Int_t nADCwords = (timebins + 2) / 3;
1649 AliDebug(3, Form("Now reading %i words for channel %2i", nADCwords, channelno));
1650 Int_t adccol = adccoloff - channelno;
1651 Int_t padcol = padcoloff - channelno;
1652 // if (adccol < 3 || adccol > 165)
1653 // AliInfo(Form("writing channel %i of det %3i %i:%2i to adcrow/-col: %i/%i padcol: %i",
1654 // channelno, fCurrHC/2, fCurrRobPos, fCurrMcmPos, row, adccol, padcol));
1656 while ((adcwc < nADCwords) &&
1657 (*(fPayloadCurr) != fgkDataEndmarker) &&
1658 (fPayloadCurr - fPayloadStart < fPayloadSize)) {
1659 int check = 0x3 & *fPayloadCurr;
1660 if (channelno % 2 != 0) { // odd channel
1661 if (check != 0x2 && channelno < 21) {
1662 MCMError(kAdcCheckInvalid,
1663 "%i for %2i. ADC word in odd channel %i",
1664 check, adcwc+1, channelno);
1667 else { // even channel
1668 if (check != 0x3 && channelno < 21) {
1669 MCMError(kAdcCheckInvalid,
1670 "%i for %2i. ADC word in even channel %i",
1671 check, adcwc+1, channelno);
1675 // filling the actual timebin data
1676 int tb2 = 0x3ff & *fPayloadCurr >> 22;
1677 int tb1 = 0x3ff & *fPayloadCurr >> 12;
1678 int tb0 = 0x3ff & *fPayloadCurr >> 2;
1679 if (adcwc != 0 || fCurrNtimebins <= 30)
1680 fAdcArray->SetDataByAdcCol(row, adccol, currentTimebin++, tb0);
1683 fAdcArray->SetDataByAdcCol(row, adccol, currentTimebin++, tb1);
1684 fAdcArray->SetDataByAdcCol(row, adccol, currentTimebin++, tb2);
1690 if (adcwc != nADCwords)
1691 MCMError(kAdcDataAbort);
1694 if (padcol > 0 && padcol < 144) {
1695 fSignalIndex->AddIndexRC(row, padcol);
1701 if (fCurrSm > -1 && fCurrSm < 18) {
1702 fStats.fStatsSector[fCurrSm].fStatsHC[fCurrHC%60].fNChannels += channelcount;
1703 fStats.fStatsSector[fCurrSm].fNChannels += channelcount;
1705 if (channelcount != channelcountExp)
1706 MCMError(kAdcChannelsMiss);
1709 if (fCurrSm > -1 && fCurrSm < 18) {
1710 fStats.fStatsSector[fCurrSm].fStatsHC[fCurrHC%60].fNMCMs++;
1711 fStats.fStatsSector[fCurrSm].fNMCMs++;
1714 if (IsDumping() && DumpingMCM(fCurrHC/2, fCurrRobPos, fCurrMcmPos)) {
1715 AliInfo(DumpRaw(Form("Event %i: Det %3i ROB %i MCM %2i", fRawReader->GetEventIndex(), fCurrHC/2, fCurrRobPos, fCurrMcmPos),
1716 startPosMCM, fPayloadCurr - startPosMCM));
1719 // continue with next MCM
1722 // check for missing MCMs (if header suppression is inactive)
1723 if (((fCurrMajor & 0x1) == 0) && (mcmcount != mcmcountExp)) {
1724 LinkError(kMissMcmHeaders,
1725 "No. of MCM headers %i not as expected: %i",
1726 mcmcount, mcmcountExp);
1729 return (fPayloadCurr - start);
1732 Int_t AliTRDrawStream::ReadNonZSData()
1734 // read the non-zs data from one link from the current reading position
1736 UInt_t *start = fPayloadCurr;
1739 Int_t mcmcountExp = fCurrStack == 2 ? 48 : 64;
1740 Int_t channelcount = 0;
1741 Int_t channelcountExp = 0;
1743 Int_t currentTimebin = 0;
1746 Int_t lastmcmpos = -1;
1747 Int_t lastrobpos = -1;
1749 if (fCurrNtimebins != fNtimebins) {
1751 LinkError(kNtimebinsChanged,
1752 "No. of timebins changed from %i to %i", fNtimebins, fCurrNtimebins);
1753 fNtimebins = fCurrNtimebins;
1756 timebins = fNtimebins;
1758 while (*(fPayloadCurr) != fgkDataEndmarker &&
1759 fPayloadCurr - fPayloadStart < fPayloadSize - 2) {
1761 // ----- Checking MCM Header -----
1762 AliDebug(2, Form("MCM header: 0x%08x", *fPayloadCurr));
1764 // ----- checking for proper readout order - ROB -----
1765 fCurrRobPos = ROB(*fPayloadCurr);
1766 if (GetROBReadoutPos(ROB(*fPayloadCurr) / 2) >= lastrobpos) {
1767 if (GetROBReadoutPos(ROB(*fPayloadCurr) / 2) > lastrobpos)
1769 lastrobpos = GetROBReadoutPos(ROB(*fPayloadCurr) / 2);
1772 ROBError(kPosUnexp, Form("#%i after #%i in readout order", GetROBReadoutPos(ROB(*fPayloadCurr) / 2), lastrobpos));
1775 // ----- checking for proper readout order - MCM -----
1776 fCurrMcmPos = MCM(*fPayloadCurr);
1777 if (GetMCMReadoutPos(MCM(*fPayloadCurr)) > lastmcmpos) {
1778 lastmcmpos = GetMCMReadoutPos(MCM(*fPayloadCurr));
1781 MCMError(kPosUnexp, Form("#%i after #%i in readout order", GetMCMReadoutPos(MCM(*fPayloadCurr)), lastmcmpos));
1784 if (EvNo(*fPayloadCurr) != evno) {
1786 evno = EvNo(*fPayloadCurr);
1788 MCMError(kEvCntMismatch, "%i <-> %i", evno, EvNo(*fPayloadCurr));
1793 channelcountExp = 21;
1796 Int_t adccoloff = AdcColOffset(*fPayloadCurr);
1797 Int_t padcoloff = PadColOffset(*fPayloadCurr);
1798 Int_t row = Row(*fPayloadCurr);
1802 // ----- reading marked ADC channels -----
1803 while (channelcount < channelcountExp &&
1804 *(fPayloadCurr) != fgkDataEndmarker) {
1811 Int_t nADCwords = (timebins + 2) / 3;
1812 AliDebug(2, Form("Now looking %i words", nADCwords));
1813 Int_t adccol = adccoloff - channelno;
1814 Int_t padcol = padcoloff - channelno;
1815 while ((adcwc < nADCwords) &&
1816 (*(fPayloadCurr) != fgkDataEndmarker) &&
1817 (fPayloadCurr - fPayloadStart < fPayloadSize)) {
1818 int check = 0x3 & *fPayloadCurr;
1819 if (channelno % 2 != 0) { // odd channel
1820 if (check != 0x2 && channelno < 21) {
1821 MCMError(kAdcCheckInvalid,
1822 "%i for %2i. ADC word in odd channel %i",
1823 check, adcwc+1, channelno);
1826 else { // even channel
1827 if (check != 0x3 && channelno < 21) {
1828 MCMError(kAdcCheckInvalid,
1829 "%i for %2i. ADC word in even channel %i",
1830 check, adcwc+1, channelno);
1834 // filling the actual timebin data
1835 int tb2 = 0x3ff & *fPayloadCurr >> 22;
1836 int tb1 = 0x3ff & *fPayloadCurr >> 12;
1837 int tb0 = 0x3ff & *fPayloadCurr >> 2;
1838 if (adcwc != 0 || fCurrNtimebins <= 30)
1839 fAdcArray->SetDataByAdcCol(row, adccol, currentTimebin++, tb0);
1842 fAdcArray->SetDataByAdcCol(row, adccol, currentTimebin++, tb1);
1843 fAdcArray->SetDataByAdcCol(row, adccol, currentTimebin++, tb2);
1849 if (adcwc != nADCwords)
1850 MCMError(kAdcDataAbort);
1853 if (padcol > 0 && padcol < 144) {
1854 fSignalIndex->AddIndexRC(row, padcol);
1860 if (channelcount != channelcountExp)
1861 MCMError(kAdcChannelsMiss);
1863 // continue with next MCM
1866 // check for missing MCMs (if header suppression is inactive)
1867 if (mcmcount != mcmcountExp) {
1868 LinkError(kMissMcmHeaders,
1869 "%i not as expected: %i", mcmcount, mcmcountExp);
1872 return (fPayloadCurr - start);
1875 Int_t AliTRDrawStream::SeekNextStack()
1877 // proceed in raw data stream till the next stack
1879 if (!fCurrStackEndmarkerAvail)
1882 UInt_t *start = fPayloadCurr;
1884 // read until data endmarkers
1885 while ((fPayloadCurr - fPayloadStart < fPayloadSize-1) &&
1886 ((fPayloadCurr[0] != fgkStackEndmarker[0]) ||
1887 (fPayloadCurr[1] != fgkStackEndmarker[1])))
1890 if ((fPayloadCurr - start) != 0)
1891 StackError(kUnknown, "skipped %i words to reach stack endmarker", fPayloadCurr - start);
1893 AliDebug(2, Form("stack endmarker: 0x%08x 0x%08x", fPayloadCurr[0], fPayloadCurr[1]));
1899 return (fPayloadCurr-start);
1902 Int_t AliTRDrawStream::SeekNextLink()
1904 // proceed in raw data stream till the next link
1906 UInt_t *start = fPayloadCurr;
1908 // read until data endmarkers
1909 while (fPayloadCurr - fPayloadStart < fPayloadSize &&
1910 *fPayloadCurr != fgkDataEndmarker)
1913 // read all data endmarkers
1914 while (fPayloadCurr - fPayloadStart < fPayloadSize &&
1915 *fPayloadCurr == fgkDataEndmarker)
1918 return (fPayloadCurr - start);
1921 Bool_t AliTRDrawStream::ConnectTracklets(TTree *trklTree)
1923 // connect the tracklet tree used to store the tracklet output
1925 fTrackletTree = trklTree;
1929 if (!fTrackletTree->GetBranch("hc"))
1930 fTrackletTree->Branch("hc", &fCurrHC, "hc/I");
1932 fTrackletTree->SetBranchAddress("hc", &fCurrHC);
1934 if (!fTrackletTree->GetBranch("trkl"))
1935 fTrackletTree->Branch("trkl", &fTrackletArray);
1937 fTrackletTree->SetBranchAddress("trkl", &fTrackletArray);
1943 void AliTRDrawStream::EquipmentError(ErrorCode_t err, const char *const msg, ...)
1945 // register error according to error code on equipment level
1946 // and return the corresponding error message
1948 fLastError.fSector = fCurrEquipmentId - kDDLOffset;
1949 fLastError.fStack = -1;
1950 fLastError.fLink = -1;
1951 fLastError.fRob = -1;
1952 fLastError.fMcm = -1;
1953 fLastError.fError = err;
1954 (this->*fStoreError)();
1957 if (fgErrorDebugLevel[err] > 10)
1958 AliDebug(fgErrorDebugLevel[err],
1959 Form("Event %6i: Eq. %2d - %s : %s",
1960 fRawReader->GetEventIndex(), fCurrEquipmentId, fgkErrorMessages[err],
1961 (va_start(ap, msg), vsprintf(fErrorBuffer, msg, ap), va_end(ap), fErrorBuffer) ));
1963 AliError(Form("Event %6i: Eq. %2d - %s : %s",
1964 fRawReader->GetEventIndex(), fCurrEquipmentId, fgkErrorMessages[err],
1965 (va_start(ap, msg), vsprintf(fErrorBuffer, msg, ap), va_end(ap), fErrorBuffer) ));
1966 fErrorFlags |= fgErrorBehav[err];
1970 void AliTRDrawStream::StackError(ErrorCode_t err, const char *const msg, ...)
1972 // register error according to error code on stack level
1973 // and return the corresponding error message
1975 fLastError.fSector = fCurrEquipmentId - kDDLOffset;
1976 fLastError.fStack = fCurrSlot;
1977 fLastError.fLink = -1;
1978 fLastError.fRob = -1;
1979 fLastError.fMcm = -1;
1980 fLastError.fError = err;
1981 (this->*fStoreError)();
1984 if (fgErrorDebugLevel[err] > 0)
1985 AliDebug(fgErrorDebugLevel[err],
1986 Form("Event %6i: Eq. %2d S %i - %s : %s",
1987 fRawReader->GetEventIndex(), fCurrEquipmentId, fCurrSlot, fgkErrorMessages[err],
1988 (va_start(ap, msg), vsprintf(fErrorBuffer, msg, ap), va_end(ap), fErrorBuffer) ));
1990 AliError(Form("Event %6i: Eq. %2d S %i - %s : %s",
1991 fRawReader->GetEventIndex(), fCurrEquipmentId, fCurrSlot, fgkErrorMessages[err],
1992 (va_start(ap, msg), vsprintf(fErrorBuffer, msg, ap), va_end(ap), fErrorBuffer) ));
1993 fErrorFlags |= fgErrorBehav[err];
1997 void AliTRDrawStream::LinkError(ErrorCode_t err, const char *const msg, ...)
1999 // register error according to error code on link level
2000 // and return the corresponding error message
2002 fLastError.fSector = fCurrEquipmentId - kDDLOffset;
2003 fLastError.fStack = fCurrSlot;
2004 fLastError.fLink = fCurrLink;
2005 fLastError.fRob = -1;
2006 fLastError.fMcm = -1;
2007 fLastError.fError = err;
2008 (this->*fStoreError)();
2011 if (fgErrorDebugLevel[err] > 0)
2012 AliDebug(fgErrorDebugLevel[err],
2013 Form("Event %6i: Eq. %2d S %i l %2i - %s : %s",
2014 fRawReader->GetEventIndex(), fCurrEquipmentId, fCurrSlot, fCurrLink, fgkErrorMessages[err],
2015 (va_start(ap, msg), vsprintf(fErrorBuffer, msg, ap), va_end(ap), fErrorBuffer) ));
2017 AliError(Form("Event %6i: Eq. %2d S %i l %2i - %s : %s",
2018 fRawReader->GetEventIndex(), fCurrEquipmentId, fCurrSlot, fCurrLink, fgkErrorMessages[err],
2019 (va_start(ap, msg), vsprintf(fErrorBuffer, msg, ap), va_end(ap), fErrorBuffer) ));
2020 fErrorFlags |= fgErrorBehav[err];
2024 void AliTRDrawStream::ROBError(ErrorCode_t err, const char *const msg, ...)
2026 // register error according to error code on ROB level
2027 // and return the corresponding error message
2029 fLastError.fSector = fCurrEquipmentId - kDDLOffset;
2030 fLastError.fStack = fCurrSlot;
2031 fLastError.fLink = fCurrLink;
2032 fLastError.fRob = fCurrRobPos;
2033 fLastError.fMcm = -1;
2034 fLastError.fError = err;
2035 (this->*fStoreError)();
2038 if (fgErrorDebugLevel[err] > 0)
2039 AliDebug(fgErrorDebugLevel[err],
2040 Form("Event %6i: Eq. %2d S %i l %2i ROB %i - %s : %s",
2041 fRawReader->GetEventIndex(), fCurrEquipmentId, fCurrSlot, fCurrLink, fCurrRobPos, fgkErrorMessages[err],
2042 (va_start(ap, msg), vsprintf(fErrorBuffer, msg, ap), va_end(ap), fErrorBuffer) ));
2044 AliError(Form("Event %6i: Eq. %2d S %i l %2i ROB %i - %s : %s",
2045 fRawReader->GetEventIndex(), fCurrEquipmentId, fCurrSlot, fCurrLink, fCurrRobPos, fgkErrorMessages[err],
2046 (va_start(ap, msg), vsprintf(fErrorBuffer, msg, ap), va_end(ap), fErrorBuffer) ));
2047 fErrorFlags |= fgErrorBehav[err];
2051 void AliTRDrawStream::MCMError(ErrorCode_t err, const char *const msg, ...)
2053 // register error according to error code on MCM level
2054 // and return the corresponding error message
2056 fLastError.fSector = fCurrEquipmentId - kDDLOffset;
2057 fLastError.fStack = fCurrSlot;
2058 fLastError.fLink = fCurrLink;
2059 fLastError.fRob = fCurrRobPos;
2060 fLastError.fMcm = fCurrMcmPos;
2061 fLastError.fError = err;
2062 (this->*fStoreError)();
2065 if (fgErrorDebugLevel[err] > 0)
2066 AliDebug(fgErrorDebugLevel[err],
2067 Form("Event %6i: Eq. %2d S %i l %2i ROB %i MCM %2i - %s : %s",
2068 fRawReader->GetEventIndex(), fCurrEquipmentId, fCurrSlot, fCurrLink, fCurrRobPos, fCurrMcmPos, fgkErrorMessages[err],
2069 (va_start(ap, msg), vsprintf(fErrorBuffer, msg, ap), va_end(ap), fErrorBuffer) ));
2071 AliError(Form("Event %6i: Eq. %2d S %i l %2i ROB %i MCM %2i - %s : %s",
2072 fRawReader->GetEventIndex(), fCurrEquipmentId, fCurrSlot, fCurrLink, fCurrRobPos, fCurrMcmPos, fgkErrorMessages[err],
2073 (va_start(ap, msg), vsprintf(fErrorBuffer, msg, ap), va_end(ap), fErrorBuffer) ));
2074 fErrorFlags |= fgErrorBehav[err];
2077 const char* AliTRDrawStream::GetErrorMessage(ErrorCode_t errCode)
2079 // return the error message for the given error code
2081 if (errCode > 0 && errCode < kLastErrorCode)
2082 return fgkErrorMessages[errCode];
2087 void AliTRDrawStream::AliTRDrawStats::ClearStats()
2089 // clear statistics (includes clearing sector-wise statistics)
2092 for (Int_t iSector = 0; iSector < 18; iSector++) {
2093 fStatsSector[iSector].ClearStats();
2098 void AliTRDrawStream::AliTRDrawStats::AliTRDrawStatsSector::ClearStats()
2100 // clear statistics (includes clearing HC-wise statistics)
2108 for (Int_t iHC = 0; iHC < 60; iHC++) {
2109 fStatsHC[iHC].ClearStats();
2113 void AliTRDrawStream::AliTRDrawStats::AliTRDrawStatsSector::AliTRDrawStatsHC::ClearStats()
2124 void AliTRDrawStream::SetDumpMCM(Int_t det, Int_t rob, Int_t mcm, Bool_t dump)
2126 // mark MCM for dumping of raw data
2129 fDumpMCM[fNDumpMCMs++] = (det << 7) | (rob << 4) | mcm;
2133 for (iMCM = 0; iMCM < fNDumpMCMs; iMCM++) {
2134 if (fDumpMCM[iMCM] == ((det << 7) | (rob << 4) | mcm)) {
2139 for ( ; iMCM < fNDumpMCMs; iMCM++) {
2140 fDumpMCM[iMCM] = fDumpMCM[iMCM+1];
2145 Bool_t AliTRDrawStream::DumpingMCM(Int_t det, Int_t rob, Int_t mcm) const
2147 // check if MCM data should be dumped
2149 for (Int_t iMCM = 0; iMCM < fNDumpMCMs; iMCM++) {
2150 if (fDumpMCM[iMCM] == ((det << 7) | (rob << 4) | mcm)) {
2157 TString AliTRDrawStream::DumpRaw(TString title, const UInt_t *start, Int_t length, UInt_t endmarker)
2162 for (Int_t pos = 0; pos < length; pos += 4) {
2163 if ((start[pos+0] != endmarker) && pos+0 < length)
2164 if ((start[pos+1] != endmarker && pos+1 < length))
2165 if ((start[pos+2] != endmarker && pos+2 < length))
2166 if ((start[pos+3] != endmarker && pos+3 < length))
2167 title += Form(" 0x%08x 0x%08x 0x%08x 0x%08x\n",
2168 start[pos+0], start[pos+1], start[pos+2], start[pos+3]);
2170 title += Form(" 0x%08x 0x%08x 0x%08x 0x%08x\n",
2171 start[pos+0], start[pos+1], start[pos+2], start[pos+3]);
2175 title += Form(" 0x%08x 0x%08x 0x%08x\n",
2176 start[pos+0], start[pos+1], start[pos+2]);
2180 title += Form(" 0x%08x 0x%08x\n",
2181 start[pos+0], start[pos+1]);
2185 title += Form(" 0x%08x\n",
2193 TString AliTRDrawStream::DumpMcmHeader(TString title, UInt_t word)
2195 title += Form("0x%08x -> ROB: %i, MCM: %2i",
2196 word, ROB(word), MCM(word));
2200 TString AliTRDrawStream::DumpAdcMask(TString title, UInt_t word)
2202 title += Form("0x%08x -> #ch : %2i, 0x%06x (%2i ch)",
2203 word, GetNActiveChannels(word), GetActiveChannels(word), GetNActiveChannelsFromMask(word));
2207 AliTRDrawStream::AliTRDrawStreamError::AliTRDrawStreamError(Int_t error, Int_t sector, Int_t stack, Int_t link, Int_t rob, Int_t mcm) :
2219 void AliTRDrawStream::SortTracklets(TClonesArray *trklArray, TList &sortedTracklets, Int_t *indices)
2221 // sort tracklets for referencing from GTU tracks
2223 Int_t nTracklets = trklArray->GetEntriesFast();
2226 for (Int_t iTracklet = 0; iTracklet < nTracklets; iTracklet++) {
2227 AliTRDtrackletBase *trkl = (AliTRDtrackletBase*) ((*trklArray)[iTracklet]);
2228 Int_t hc = trkl->GetHCId();
2229 if ((hc < 0) || (hc >= 1080)) {
2230 AliErrorClass(Form("HC for tracklet: 0x%08x out of range: %i", trkl->GetTrackletWord(), trkl->GetHCId()));
2233 AliDebugClass(5, Form("hc: %4i : 0x%08x z: %2i", hc, trkl->GetTrackletWord(), trkl->GetZbin()));
2235 AliDebugClass(2, Form("set tracklet index for HC %i to %i", hc, iTracklet));
2236 indices[hc] = iTracklet + 1;
2241 for (Int_t iDet = 0; iDet < 540; iDet++) {
2242 Int_t trklIndexA = indices[2*iDet + 0] - 1;
2243 Int_t trklIndexB = indices[2*iDet + 1] - 1;
2244 Int_t trklIndex = sortedTracklets.GetEntries();
2245 AliTRDtrackletBase *trklA = trklIndexA > -1 ? (AliTRDtrackletBase*) ((*trklArray)[trklIndexA]) : 0x0;
2246 AliTRDtrackletBase *trklB = trklIndexB > -1 ? (AliTRDtrackletBase*) ((*trklArray)[trklIndexB]) : 0x0;
2247 AliTRDtrackletBase *trklNext = 0x0;
2248 while (trklA != 0x0 || trklB != 0x0) {
2249 AliDebugClass(5, Form("det %i - A: %i/%i -> %p, B: %i/%i -> %p",
2250 iDet, trklIndexA, nTracklets, trklA, trklIndexB, nTracklets, trklB));
2254 trklB = trklIndexB < nTracklets ? (AliTRDtrackletBase*) ((*trklArray)[trklIndexB]) : 0x0;
2255 if (trklB && trklB->GetHCId() != 2*iDet + 1)
2258 else if (trklB == 0x0) {
2261 trklA = trklIndexA < nTracklets ? (AliTRDtrackletBase*) ((*trklArray)[trklIndexA]) : 0x0;
2262 if (trklA && trklA->GetHCId() != 2*iDet)
2266 if ((trklA->GetZbin() < trklB->GetZbin()) ||
2267 ((trklA->GetZbin() == trklB->GetZbin()) && (trklA->GetYbin() < trklB->GetYbin()))) {
2270 trklA = trklIndexA < nTracklets ? (AliTRDtrackletBase*) ((*trklArray)[trklIndexA]) : 0x0;
2271 if (trklA && trklA->GetHCId() != 2*iDet)
2277 trklB = trklIndexB < nTracklets ? (AliTRDtrackletBase*) ((*trklArray)[trklIndexB]) : 0x0;
2278 if (trklB && trklB->GetHCId() != 2*iDet + 1)
2283 Int_t label = -2; // mark raw tracklets with label -2
2284 if (AliTRDtrackletMCM *trklMCM = dynamic_cast<AliTRDtrackletMCM*> (trklNext))
2285 label = trklMCM->GetLabel();
2286 AliESDTrdTracklet *esdTracklet = new AliESDTrdTracklet(trklNext->GetTrackletWord(), trklNext->GetHCId(), label);
2287 sortedTracklets.Add(esdTracklet);
2290 // updating tracklet indices as in output
2291 if (sortedTracklets.GetEntries() != trklIndex) {
2292 indices[2*iDet + 0] = indices[2*iDet + 1] = trklIndex;
2295 indices[2*iDet + 0] = indices[2*iDet + 1] = -1;