Add online raw reader by Jochen as default version
[u/mrichter/AliRoot.git] / TRD / AliTRDrawStream.cxx
1 #include "TClonesArray.h"
2 #include "TTree.h"
3
4 #include "AliLog.h"
5 #include "AliRawReader.h"
6 #include "AliTRDdigitsManager.h"
7 #include "AliTRDdigitsParam.h"
8 #include "AliTRDtrapConfig.h"
9 #include "AliTRDarrayADC.h"
10 #include "AliTRDarrayDictionary.h"
11 #include "AliTRDSignalIndex.h"
12 #include "AliTRDtrackletWord.h"
13
14 #include "AliTRDrawStream.h"
15
16 // temporary
17 #include "AliRunLoader.h"
18
19 ClassImp(AliTRDrawStream)
20
21 // some static information 
22 const Int_t AliTRDrawStream::fgkMcmOrder[] = {12, 13, 14, 15, 
23                                                     8, 9, 10, 11, 
24                                                     4, 5, 6, 7, 
25                                                     0, 1, 2, 3};
26 const Int_t  AliTRDrawStream::fgkRobOrder [] = {0, 1, 2, 3};
27 const Int_t  AliTRDrawStream::fgkNlinks = 12;
28 const Int_t  AliTRDrawStream::fgkNstacks = 5;
29 const UInt_t AliTRDrawStream::fgkDataEndmarker     = 0x00000000;
30 const UInt_t AliTRDrawStream::fgkTrackletEndmarker = 0x10001000;
31
32 char* AliTRDrawStream::fgErrorMessages[] = {
33   "Unknown error",
34   "Link monitor active",
35   "Pretrigger counter mismatch",
36   "not a TRD equipment (1024-1041)",
37   "Invalid Stack header",
38   "Invalid detector number",
39   "No digits could be retrieved from the digitsmanager",
40   "HC header mismatch", 
41   "HC check bits wrong",
42   "Unexpected position in readout stream",
43   "Invalid testpattern mode",
44   "Testpattern mismatch",
45   "Number of timebins changed",
46   "ADC mask inconsistent", 
47   "ADC check bits invalid", 
48   "Missing ADC data",
49   "Missing expected ADC channels",
50   "Missing MCM headers"
51 };
52
53 AliTRDrawStream::AliTRDrawStream(AliRawReader *rawReader) :
54   fRawReader(rawReader),
55   fDigitsManager(0x0),
56   fDigitsParam(0x0),
57   fErrors(0x0),
58   fLastError(),
59   fPayloadStart(0x0),
60   fPayloadCurr(0x0),
61   fPayloadSize(0),
62   fNtimebins(-1),
63   fLastEvId(-1),
64   fCurrSlot(-1),
65   fCurrLink(-1),
66   fCurrRobPos(-1),
67   fCurrMcmPos(-1),
68   fCurrEquipmentId(0),
69   fCurrSmuIndexHeaderSize(0),
70   fCurrSmuIndexHeaderVersion(0),
71   fCurrTrackEnable(0),
72   fCurrTrackletEnable(0),
73   fCurrStackMask(0),
74   fCurrStackIndexWord(0x0),
75   fCurrStackHeaderSize(0x0),
76   fCurrStackHeaderVersion(0x0),
77   fCurrLinkMask(0x0),
78   fCurrCleanCheckout(0x0),
79   fCurrBoardId(0x0),
80   fCurrHwRev(0x0),
81   fCurrLinkMonitorFlags(0x0),
82   fCurrLinkDataTypeFlags(0x0),
83   fCurrLinkDebugFlags(0x0),
84   fCurrSpecial(-1),
85   fCurrMajor(-1),
86   fCurrMinor(-1),
87   fCurrAddHcWords(-1),
88   fCurrSm(-1),
89   fCurrStack(-1),
90   fCurrLayer(-1),
91   fCurrSide(-1),
92   fCurrHC(-1),
93   fCurrCheck(-1),
94   fCurrNtimebins(-1),
95   fCurrBC(-1),
96   fCurrPtrgCnt(-1),
97   fCurrPtrgPhase(-1),
98   fTrackletArray(0x0),
99   fAdcArray(0x0),
100   fSignalIndex(0x0),
101   fTrackletTree(0x0)
102 {
103   // default constructor
104
105   fCurrStackIndexWord     = new UInt_t[fgkNstacks];      
106   fCurrStackHeaderSize    = new UInt_t[fgkNstacks];      
107   fCurrStackHeaderVersion = new UInt_t[fgkNstacks];
108   fCurrLinkMask           = new UInt_t[fgkNstacks];              
109   fCurrCleanCheckout      = new UInt_t[fgkNstacks];      
110   fCurrBoardId            = new UInt_t[fgkNstacks];              
111   fCurrHwRev              = new UInt_t[fgkNstacks];             
112   fCurrLinkMonitorFlags   = new UInt_t[fgkNstacks * fgkNlinks];
113   fCurrLinkDataTypeFlags  = new UInt_t[fgkNstacks * fgkNlinks];
114   fCurrLinkDebugFlags     = new UInt_t[fgkNstacks * fgkNlinks];
115
116   // preparing TClonesArray
117   fTrackletArray = new TClonesArray("AliTRDtrackletWord", 256);
118
119   // setting up the error tree
120   fErrors = new TTree("errorStats", "Error statistics");
121   fErrors->SetDirectory(0x0);
122   fErrors->Branch("error", &fLastError, "sector/I:stack:link:error:rob:mcm");
123   fErrors->SetCircular(1000);
124 }
125
126 AliTRDrawStream::~AliTRDrawStream()
127 {
128   // destructor
129
130   delete fErrors;
131
132   delete [] fCurrStackIndexWord;
133   delete [] fCurrStackHeaderSize;
134   delete [] fCurrStackHeaderVersion;
135   delete [] fCurrLinkMask;
136   delete [] fCurrCleanCheckout;
137   delete [] fCurrBoardId;
138   delete [] fCurrHwRev;
139   delete [] fCurrLinkMonitorFlags;
140   delete [] fCurrLinkDataTypeFlags;
141   delete [] fCurrLinkDebugFlags;
142 }
143
144 Bool_t AliTRDrawStream::ReadEvent(TTree *trackletTree)
145 {
146   // read the current event from the raw reader and fill it to the digits manager
147
148   if (!fRawReader) {
149     AliError("No raw reader available");
150     return kFALSE;
151   }
152
153   // tracklet output
154   fTrackletTree = trackletTree;
155   if (fTrackletTree) {
156     if (!fTrackletTree->GetBranch("trkl")) {
157       fTrackletTree->Branch("hc", &fCurrHC, "hc/I");
158       fTrackletTree->Branch("trkl", &fTrackletArray);
159     } 
160     else {
161       fTrackletTree->SetBranchAddress("hc", &fCurrHC);
162       fTrackletTree->SetBranchAddress("trkl", &fTrackletArray);
163     }
164   }
165
166   // some preparations
167   fDigitsParam = 0x0;
168
169   // loop over all DDLs
170   // data starts with GTU payload, i.e. SMU index word
171   UChar_t *buffer = 0x0;
172
173   while (fRawReader->ReadNextData(buffer)) {
174
175     fCurrEquipmentId = fRawReader->GetEquipmentId();
176     AliDebug(2, Form("equipment: %i", fCurrEquipmentId));
177
178     if (fCurrEquipmentId < 1024 || fCurrEquipmentId > 1041) {
179       AliError(EquipmentError(kNonTrdEq, "Skipping"));
180       continue;
181     }
182
183     // setting the pointer to data and current reading position
184     fPayloadCurr = fPayloadStart = (UInt_t*) (buffer);
185     fPayloadSize = fRawReader->GetDataSize() / sizeof(UInt_t);
186     AliDebug(2, Form("Read buffer of size: %i", fRawReader->GetDataSize()));
187
188     // read SMU index header
189     if (ReadSmHeader() < 0) {
190       AliError(Form("Reading SMU header failed, skipping this DDL %i", fCurrEquipmentId));
191       continue;
192     }
193
194     // read stack index header
195     for (Int_t iStack = 0; iStack < 5; iStack++) {
196       if (fCurrStackMask & (1 << iStack) != 0) 
197         ReadStackIndexHeader(iStack);
198     }
199
200     for (Int_t iStack = 0; iStack < 5; iStack++) {
201       fCurrSlot = iStack;
202       if (fCurrStackMask & (1 << fCurrSlot) == 0) 
203         continue;
204
205       AliDebug(2, Form("Stack %i, Link mask: 0x%02x", fCurrSlot, fCurrLinkMask[fCurrSlot]));
206       for (Int_t iLink = 0; iLink < 12; iLink++) {
207         fCurrLink = iLink;
208         fCurrHC   = fCurrSm * 60 + fCurrSlot * 12 + iLink;
209         if ((fCurrLinkMask[fCurrSlot] & (1 << fCurrLink)) == 0)
210           continue;
211         
212         if (fCurrLinkMonitorFlags[fCurrSlot*fgkNlinks + fCurrLink] != 0)
213           LinkError(kLinkMonitor);
214         AliDebug(2, Form("Payload for S%il%i", fCurrSlot, fCurrLink));
215         for (Int_t i = 0; i < 10; i++) //???
216           AliDebug(5, Form("%3i: 0x%08x 0x%08x 0x%08x 0x%08x", i*4,
217                  fPayloadCurr[4*i], fPayloadCurr[4*i+1], 
218                  fPayloadCurr[4*i+2], fPayloadCurr[4*i+3]));
219
220         // read the data from one HC
221         ReadLinkData();
222
223         // read all data endmarkers
224         while (fPayloadCurr - fPayloadStart < fPayloadSize &&
225                *fPayloadCurr == fgkDataEndmarker)
226           fPayloadCurr++;
227
228       }
229     }
230   }
231   return kTRUE;
232 }
233
234
235 Bool_t AliTRDrawStream::NextDDL()
236 {
237   if (!fRawReader)
238     return kFALSE;
239
240   fCurrEquipmentId = 0;
241   fCurrSlot = 0;
242   fCurrLink = 0;
243
244   UChar_t *buffer = 0x0;
245
246   while (fRawReader->ReadNextData(buffer)) {
247
248     fCurrEquipmentId = fRawReader->GetEquipmentId();
249     AliDebug(2, Form("equipment: %i", fCurrEquipmentId));
250     
251     if (fCurrEquipmentId < 1024 || fCurrEquipmentId > 1041) {
252       AliError(EquipmentError(kNonTrdEq, "Skipping"));
253       continue;
254     }
255
256     // setting the pointer to data and current reading position
257     fPayloadCurr = fPayloadStart = (UInt_t*) (buffer);
258     fPayloadSize = fRawReader->GetDataSize() / sizeof(UInt_t);
259     AliDebug(2, Form("Read buffer of size: %i", fRawReader->GetDataSize()));
260
261     // read SMU index header
262     if (ReadSmHeader() < 0) {
263       AliError(Form("Reading SMU header failed, skipping this DDL %i", fCurrEquipmentId));
264       continue;
265     }
266     
267     // read stack index header
268     for (Int_t iStack = 0; iStack < 5; iStack++) {
269       if (fCurrStackMask & (1 << iStack) != 0) {
270         ReadStackIndexHeader(iStack);
271       }
272     }
273     return kTRUE;
274   }
275
276   return kFALSE;
277 }
278
279
280 Int_t AliTRDrawStream::NextChamber(AliTRDdigitsManager *digMgr, UInt_t ** /* trackletContainer */, UShort_t ** /* errorContainer */)
281 {
282   fDigitsManager = digMgr; 
283   fDigitsParam   = 0x0;
284
285   // tracklet output preparation
286   fTrackletTree = 0x0;
287   //AliDataLoader *trklLoader = AliRunLoader::Instance()->GetLoader("TRDLoader")->GetDataLoader("tracklets");
288   //if (trklLoader) 
289   //  fTrackletTree = trklLoader->Tree();
290
291   if (!fRawReader) {
292     AliError("No raw reader available");
293     return -1;
294   }
295
296   if (fCurrSlot < 0 || fCurrSlot >= 5) {
297     if (!NextDDL()) {
298       fCurrSlot = -1;
299       return -1;
300     }
301   }
302
303   AliDebug(2, Form("Stack %i, Link %i, mask: 0x%02x", fCurrSlot, fCurrLink, fCurrLinkMask[fCurrSlot]));
304   fCurrHC   = (fCurrEquipmentId - 1024) * 60 + fCurrSlot * 12 + fCurrLink;
305
306   if (fCurrLinkMonitorFlags[fCurrSlot*fgkNlinks + fCurrLink] != 0)
307     LinkError(kLinkMonitor);
308
309   // read the data from one HC
310   ReadLinkData();
311   
312   // read all data endmarkers
313   while (fPayloadCurr - fPayloadStart < fPayloadSize && 
314          *fPayloadCurr == fgkDataEndmarker) 
315     fPayloadCurr++;
316
317   if (fCurrLink % 2 == 0) {
318     // if we just read the A-side HC then also check the B-side
319     fCurrLink++;
320     if (fCurrLinkMask[fCurrSlot] & (1 << fCurrLink)) {
321       ReadLinkData();
322       while (fPayloadCurr - fPayloadStart < fPayloadSize && 
323              *fPayloadCurr == fgkDataEndmarker) 
324         fPayloadCurr++;
325     }
326   }
327
328   //??? to check 
329   do {
330     fCurrLink++; 
331     if (fCurrLink > 11) {
332       fCurrLink = 0;
333       fCurrSlot++;
334     }
335   } while ((fCurrSlot < 5) && 
336            ((fCurrStackMask & (1 << fCurrSlot) == 0) || 
337             (fCurrLinkMask[fCurrSlot] & (1 << fCurrLink)) == 0));
338
339   return fCurrHC / 2;
340 }
341
342
343 Int_t AliTRDrawStream::ReadSmHeader()
344 {
345   // read the SMU index header at the current reading position 
346   // and store the information in the corresponding variables
347
348   if (fPayloadCurr - fPayloadStart >= fPayloadSize - 1) {
349     AliError(EquipmentError(kUnknown, "SM Header aborted"));
350     return -1;
351   }
352
353   fCurrSmuIndexHeaderSize     = ((*fPayloadCurr) >> 16) & 0xffff;
354   fCurrSmuIndexHeaderVersion  = ((*fPayloadCurr) >> 12) &    0xf;
355   fCurrTrackEnable            = ((*fPayloadCurr) >>  6) &    0x1;
356   fCurrTrackletEnable         = ((*fPayloadCurr) >>  5) &    0x1;
357   fCurrStackMask              = ((*fPayloadCurr)      ) &   0x1f;
358
359   AliDebug(5, Form("SMU header: size: %i, version: %i, track enable: %i, tracklet enable: %i, stack mask: %2x",
360                    fCurrSmuIndexHeaderSize, 
361                    fCurrSmuIndexHeaderVersion, 
362                    fCurrTrackEnable, 
363                    fCurrTrackletEnable,
364                    fCurrStackMask));
365   
366   fPayloadCurr += fCurrSmuIndexHeaderSize + 1;
367
368   return fCurrSmuIndexHeaderSize + 1;
369 }
370
371 Int_t AliTRDrawStream::ReadStackIndexHeader(Int_t stack)
372 {
373   // read the stack index header 
374   // and store the information in the corresponding variables
375
376   fCurrStackIndexWord[stack]     = *fPayloadCurr;
377   fCurrStackHeaderSize[stack]    = (((*fPayloadCurr) >> 16) & 0xffff) + 1;
378   fCurrStackHeaderVersion[stack] = ((*fPayloadCurr) >> 12) & 0xf;
379   fCurrLinkMask[stack]           = (*fPayloadCurr) & 0xfff;
380
381   if (fPayloadCurr - fPayloadStart >= fPayloadSize - fCurrStackHeaderSize[stack]) {
382     AliError(StackError(kStackHeaderInvalid, "Stack index header aborted"));
383     return -1;
384   }
385
386   switch (fCurrStackHeaderVersion[stack]) {
387   case 0xa: 
388     if (fCurrStackHeaderSize[stack] < 8) {
389       AliError(StackError(kStackHeaderInvalid, "Stack header smaller than expected!"));
390       return -1;
391     }
392     
393     fCurrCleanCheckout[stack] = fPayloadCurr[1] & 0x1;
394     fCurrBoardId[stack]       = (fPayloadCurr[1] >> 8) & 0xff;
395     fCurrHwRev[stack]         = (fPayloadCurr[1] >> 16) & 0xffff;
396     
397     for (Int_t iLayer = 0; iLayer < 6; iLayer++) {
398       // A side
399       fCurrLinkMonitorFlags  [stack*fgkNlinks + iLayer*2]      = fPayloadCurr[iLayer+2] & 0xf;
400       fCurrLinkDataTypeFlags [stack*fgkNlinks + iLayer*2]      = (fPayloadCurr[iLayer+2] >> 4) & 0x3;
401       fCurrLinkDebugFlags    [stack*fgkNlinks + iLayer*2]      = (fPayloadCurr[iLayer+2] >> 12) & 0xf;
402       // B side
403       fCurrLinkMonitorFlags  [stack*fgkNlinks + iLayer*2 + 1]  = (fPayloadCurr[iLayer+2] >> 16) & 0xf;
404       fCurrLinkDataTypeFlags [stack*fgkNlinks + iLayer*2 + 1]  = (fPayloadCurr[iLayer+2] >> 20) & 0x3;
405       fCurrLinkDebugFlags    [stack*fgkNlinks + iLayer*2 + 1]  = (fPayloadCurr[iLayer+2] >> 28) & 0xf;
406     }
407     break;
408     
409   default:
410     AliError(StackError(kStackHeaderInvalid, Form("Invalid Stack Index Header version %x", 
411                                                   fCurrStackHeaderVersion[stack])));
412   }
413   
414   fPayloadCurr += fCurrStackHeaderSize[stack];
415
416   return fCurrStackHeaderSize[stack];
417 }
418
419 Int_t AliTRDrawStream::ReadLinkData()
420 {
421   // read the data in one link (one HC) until the data endmarker is reached
422
423   Int_t count = 0;
424
425   count += ReadTracklets();
426
427   count += ReadHcHeader();
428
429   Int_t det = fCurrSm * 30 + fCurrStack * 6 + fCurrLayer;
430
431   if (det > -1 && det < 540) {
432     
433     if ((fAdcArray = fDigitsManager->GetDigits(det))) {
434       //fAdcArray->Expand();
435       if (fAdcArray->GetNtime() != fCurrNtimebins)
436         fAdcArray->Allocate(16, 144, fCurrNtimebins);
437     }
438     else {
439       AliError(LinkError(kNoDigits));
440     }
441     
442     if (!fDigitsParam) {
443       fDigitsParam = fDigitsManager->GetDigitsParam();
444     }
445     if (fDigitsParam) {
446       fDigitsParam->SetPretriggerPhase(det, fCurrPtrgPhase);
447       fDigitsParam->SetNTimeBins(det, fCurrNtimebins);
448       fDigitsParam->SetADCbaseline(det, 10);
449     }
450     
451     if (fDigitsManager->UsesDictionaries()) {
452       fDigitsManager->GetDictionary(det, 0)->Reset();
453       fDigitsManager->GetDictionary(det, 1)->Reset();
454       fDigitsManager->GetDictionary(det, 2)->Reset();
455     }
456
457     if ((fSignalIndex = fDigitsManager->GetIndexes(det))) {
458       fSignalIndex->SetSM(fCurrSm);
459       fSignalIndex->SetStack(fCurrStack);
460       fSignalIndex->SetLayer(fCurrLayer);
461       fSignalIndex->SetDetNumber(det);
462       if (!fSignalIndex->IsAllocated())
463         fSignalIndex->Allocate(16, 144, fCurrNtimebins);
464     }
465     
466     // ----- check which kind of data -----
467     if (fCurrMajor & 0x40) {
468       if ((fCurrMajor & 0x7) == 0x7) {
469         AliDebug(1, "This is a config event");
470         UInt_t *startPos = fPayloadCurr;
471         while (fPayloadCurr - fPayloadStart < fPayloadSize &&
472                *fPayloadCurr != fgkDataEndmarker)
473           fPayloadCurr++;
474         count += fPayloadCurr - startPos;
475         
476         // feeding TRAP config
477         AliTRDtrapConfig *trapcfg = AliTRDtrapConfig::Instance();
478         //      trapcfg->ReadPackedConfig(fCurrHC, startPos, fPayloadCurr - startPos);
479       }
480       else {
481         Int_t tpmode = fCurrMajor & 0x7;
482         AliDebug(1, Form("Checking testpattern (mode %i) data", tpmode));
483         ReadTPData(tpmode);
484       }
485     }
486     else if (fCurrMajor & 0x20) {
487       AliDebug(1, "This is a zs event");
488       count += ReadZSData();
489     }
490     else {
491       AliDebug(1, "This is a nozs event");
492       count += ReadNonZSData();
493     }
494   }
495   else {
496     AliError(LinkError(kInvalidDetector, Form("%i", det)));
497     UInt_t *startPos = fPayloadCurr;
498     while (fPayloadCurr - fPayloadStart < fPayloadSize &&
499            *fPayloadCurr != fgkDataEndmarker)
500       fPayloadCurr++;
501     count += fPayloadCurr - startPos;
502   }
503   
504   return count;
505 }
506
507 Int_t AliTRDrawStream::ReadTracklets()
508 {
509   // read the tracklets from one HC
510
511   fTrackletArray->Clear();
512
513   UInt_t *start = fPayloadCurr;
514   while (*(fPayloadCurr) != fgkTrackletEndmarker && 
515          fPayloadCurr - fPayloadStart < fPayloadSize) {
516
517     new ((*fTrackletArray)[fTrackletArray->GetEntriesFast()]) AliTRDtrackletWord(*(fPayloadCurr));
518
519     fPayloadCurr++;
520   }
521
522   if (fTrackletArray->GetEntriesFast() > 0) {
523     AliDebug(1, Form("Found %i tracklets in %i %i %i (ev. %i)", fTrackletArray->GetEntriesFast(), 
524                      fCurrSm, fCurrSlot, fCurrLink, fRawReader->GetEventIndex()));
525     if (fTrackletTree)
526       fTrackletTree->Fill();
527   }
528
529   // loop over remaining tracklet endmarkers
530   while ((*(fPayloadCurr) == fgkTrackletEndmarker && 
531           fPayloadCurr - fPayloadStart < fPayloadSize)) 
532     fPayloadCurr++;
533   
534   return (fPayloadCurr - start) / sizeof(UInt_t);
535 }
536
537 Int_t AliTRDrawStream::ReadHcHeader()
538 {
539   // read and parse the HC header of one HC
540   // and store the information in the corresponding variables
541
542   UInt_t *start = fPayloadCurr;
543   // check not to be at the data endmarker
544   if (*fPayloadCurr == fgkDataEndmarker)
545     return 0;
546
547   fCurrSpecial    = (*fPayloadCurr >> 31) & 0x1;
548   fCurrMajor      = (*fPayloadCurr >> 24) & 0x7f;
549   fCurrMinor      = (*fPayloadCurr >> 17) & 0x7f;
550   fCurrAddHcWords = (*fPayloadCurr >> 14) & 0x7;
551   fCurrSm         = (*fPayloadCurr >> 9) & 0x1f;
552   fCurrLayer      = (*fPayloadCurr >> 6) & 0x7;
553   fCurrStack      = (*fPayloadCurr >> 3) & 0x7;
554   fCurrSide       = (*fPayloadCurr >> 2) & 0x1;
555   fCurrCheck      = (*fPayloadCurr) & 0x3;
556
557   if (fCurrSm != (fCurrEquipmentId - 1024) || 
558       fCurrStack != fCurrSlot || 
559       fCurrLayer != fCurrLink / 2 || 
560       fCurrSide != fCurrLink % 2) {
561     AliError(LinkError(kHCmismatch,
562                        Form("HC: %i, %i, %i, %i\n 0x%08x 0x%08x 0x%08x 0x%08x", 
563                             fCurrSm, fCurrStack, fCurrLayer, fCurrSide,
564                             fPayloadCurr[0], fPayloadCurr[1], fPayloadCurr[2], fPayloadCurr[3]) ));
565   }
566   if (fCurrCheck != 0x1) {
567     AliError(LinkError(kHCcheckFailed));
568   }
569   
570   if (fCurrAddHcWords > 0) {
571     fCurrNtimebins = (fPayloadCurr[1] >> 26) & 0x3f;
572     fCurrBC = (fPayloadCurr[1] >> 10) & 0xffff;
573     fCurrPtrgCnt = (fPayloadCurr[1] >> 6) & 0xf;
574     fCurrPtrgPhase = (fPayloadCurr[1] >> 2) & 0xf;
575   }
576   
577   fPayloadCurr += 1 + fCurrAddHcWords;
578   
579   return (fPayloadCurr - start) / sizeof(UInt_t);
580 }
581
582 Int_t AliTRDrawStream::ReadTPData(Int_t mode)
583 {
584   // testing of testpattern 1 to 3 (hardcoded), 0 missing
585   // evcnt checking missing
586   Int_t cpu = 0;
587   Int_t cpufromchannel[] = {0, 0, 0, 0, 0,  1, 1, 1, 1, 1,  2, 2, 2, 2, 2,  3, 3, 3, 3, 3, 3};
588   Int_t evcnt = 0;
589   Int_t count = 0;
590   Int_t mcmcount = -1;
591   Int_t wordcount = 0;
592   Int_t channelcount = 0;
593   UInt_t expword = 0;
594   UInt_t expadcval = 0;
595   UInt_t diff = 0;
596   Int_t lastmcmpos = -1;
597   Int_t lastrobpos = -1;
598
599   UInt_t* start = fPayloadCurr;
600
601   while (*(fPayloadCurr) != fgkDataEndmarker && 
602          fPayloadCurr - fPayloadStart < fPayloadSize - 1) {
603
604     // ----- Checking MCM Header -----
605     AliDebug(2, Form("MCM header: 0x%08x", *fPayloadCurr));
606     mcmcount++;
607     
608     // ----- checking for proper readout order - ROB -----
609     if (GetROBReadoutPos(ROB(*fPayloadCurr) / 2) >= lastrobpos) {
610       lastrobpos = GetROBReadoutPos(ROB(*fPayloadCurr) / 2);
611     }
612     else {
613       AliError(ROBError(kPosUnexp));
614     }
615     fCurrRobPos = ROB(*fPayloadCurr);
616     
617     // ----- checking for proper readout order - MCM -----
618     if (GetMCMReadoutPos(MCM(*fPayloadCurr)) >= (lastmcmpos + 1) % 16) {
619       lastmcmpos = GetMCMReadoutPos(MCM(*fPayloadCurr));
620     }
621     else {
622       AliError(MCMError(kPosUnexp));
623     }
624     fCurrMcmPos = MCM(*fPayloadCurr);
625     
626
627     fPayloadCurr++;
628     
629     evcnt = 0x3f & *fPayloadCurr >> 26;
630     cpu = -1;
631     channelcount = 0;
632     while (channelcount < 21) {
633       count = 0;
634       if (cpu != cpufromchannel[channelcount]) {
635         cpu = cpufromchannel[channelcount];
636         expadcval = (1 << 9) | (fCurrRobPos << 6) | (fCurrMcmPos << 2) | cpu;
637         wordcount = 0;
638       }
639       
640       while (count < 10) {
641         if (channelcount % 2 == 0)
642           expword = 0x3;
643         else 
644           expword = 0x2;
645         
646         if (mode == 1) {
647           // ----- TP 1 -----
648           expword |= expadcval << 2;
649           expadcval = ( (expadcval << 1) | ( ( (expadcval >> 9) ^ (expadcval >> 6) ) & 1) ) & 0x3FF;
650           expword |= expadcval << 12;
651           expadcval = ( (expadcval << 1) | ( ( (expadcval >> 9) ^ (expadcval >> 6) ) & 1) ) & 0x3FF;
652           expword |= expadcval << 22;
653           expadcval = ( (expadcval << 1) | ( ( (expadcval >> 9) ^ (expadcval >> 6) ) & 1) ) & 0x3FF;
654         }
655         else if (mode == 2) {
656           // ----- TP 2 ------
657           expword = ((0x3f & evcnt) << 26) | ((fCurrSm + 1) << 21) | ((fCurrLayer + 1) << 18) | 
658             ((fCurrStack + 1) << 15) | 
659             (fCurrRobPos << 12) | (fCurrMcmPos << 8) | (cpu << 6) | (wordcount + 1); 
660         }
661         else if (mode == 3) {
662           // ----- TP 3 -----
663           expword = ((0xfff & evcnt) << 20) | (fCurrSm << 15) | (fCurrLink/2 << 12) | (fCurrStack << 9) | 
664             (fCurrRobPos << 6) | (fCurrMcmPos << 2) | (cpu << 0); 
665         }
666         else {
667           expword = 0;
668           AliError(LinkError(kTPmodeInvalid, "Just reading"));
669         }
670
671         diff = *fPayloadCurr ^ expword;
672         if (diff != 0) {
673           AliError(MCMError(kTPmismatch,
674                             Form("Seen 0x%08x, expected 0x%08x, diff: 0x%08x (0x%02x)", 
675                                  *fPayloadCurr, expword, diff, 0xff & (diff | diff >> 8 | diff >> 16 | diff >> 24)) ));
676         }
677         fPayloadCurr++;
678         count++;
679         wordcount++;
680       }
681       channelcount++;
682     }
683     // continue with next MCM
684   }
685   return fPayloadCurr - start; 
686 }
687
688
689 Int_t AliTRDrawStream::ReadZSData()
690 {
691   // read the zs data from one link from the current reading position
692   
693   UInt_t *start = fPayloadCurr;
694   
695   Int_t mcmcount = 0;
696   Int_t mcmcount_exp = fCurrStack == 2 ? 48 : 64;
697   Int_t channelcount = 0;
698   Int_t channelcount_exp = 0;
699   Int_t channelcount_max = 0;
700   Int_t timebins;
701   Int_t currentTimebin = 0;
702   Int_t adcwc = 0;
703   Int_t evno = -1;
704   Int_t lastmcmpos = -1;
705   Int_t lastrobpos = -1;
706
707   if (fCurrNtimebins != fNtimebins) {
708     if (fNtimebins > 0) 
709       AliError(LinkError(kNtimebinsChanged,
710                          Form("No. of timebins changed from %i to %i", fNtimebins, fCurrNtimebins) ));
711     fNtimebins = fCurrNtimebins;
712   }
713   
714   timebins = fNtimebins;
715   
716   while (*(fPayloadCurr) != fgkDataEndmarker && 
717          fPayloadCurr - fPayloadStart < fPayloadSize) {
718     
719     // ----- Checking MCM Header -----
720     AliDebug(2, Form("MCM header: 0x%08x", *fPayloadCurr));
721     
722     // ----- checking for proper readout order - ROB -----
723     if (GetROBReadoutPos(ROB(*fPayloadCurr) / 2) >= lastrobpos) {
724       lastrobpos = GetROBReadoutPos(ROB(*fPayloadCurr) / 2);
725     }
726     else {
727       AliError(ROBError(kPosUnexp));
728     }
729     fCurrRobPos = ROB(*fPayloadCurr);
730     
731     // ----- checking for proper readout order - MCM -----
732     if (GetMCMReadoutPos(MCM(*fPayloadCurr)) > lastmcmpos) {
733       lastmcmpos = GetMCMReadoutPos(lastmcmpos);
734     }
735     else {
736       AliError(MCMError(kPosUnexp));
737     }
738     fCurrMcmPos = MCM(*fPayloadCurr);
739     
740     if (EvNo(*fPayloadCurr) != evno) {
741       if (evno == -1)
742         evno = EvNo(*fPayloadCurr);
743       else {
744         AliError(MCMError(kPtrgCntMismatch,
745                           Form("%i <-> %i", evno, EvNo(*fPayloadCurr)) ));
746       }
747     }
748     Int_t adccoloff = AdcColOffset(*fPayloadCurr);
749     Int_t padcoloff = PadColOffset(*fPayloadCurr);
750     Int_t row = Row(*fPayloadCurr);
751     fPayloadCurr++;
752     
753     // ----- Reading ADC channels -----
754     AliDebug(2, Form("ADC mask: 0x%08x", *fPayloadCurr));
755     
756     // ----- analysing the ADC mask -----
757     channelcount = 0;
758     channelcount_exp = GetNActiveChannelsFromMask(*fPayloadCurr);
759     channelcount_max = GetNActiveChannels(*fPayloadCurr);
760     Int_t channelmask = GetActiveChannels(*fPayloadCurr);
761     Int_t channelno = -1;
762     fPayloadCurr++; 
763
764     if (channelcount_exp != channelcount_max) {
765       if (channelcount_exp > channelcount_max) {
766         Int_t temp = channelcount_exp;
767         channelcount_exp = channelcount_max;
768         channelcount_max = temp;
769       }
770       while (channelcount_exp < channelcount_max && channelcount_exp < 21 && 
771              fPayloadCurr - fPayloadStart < fPayloadSize - 10 * channelcount_exp - 1) {
772         AliError(MCMError(kAdcMaskInconsistent,
773                           Form("Possible MCM-H: 0x%08x, possible ADC-mask: 0x%08x", 
774                                *(fPayloadCurr + 10 * channelcount_exp), 
775                                *(fPayloadCurr + 10 * channelcount_exp + 1) ) ));
776         if (!CouldBeMCMhdr( *(fPayloadCurr + 10 * channelcount_exp)) && !CouldBeADCmask( *(fPayloadCurr + 10 * channelcount_exp + 1))) 
777           channelcount_exp++;
778         else {
779           break;
780         }
781       }
782       AliError(MCMError(kAdcMaskInconsistent,
783                         Form("Inconsistency in no. of active channels: Counter: %i, Mask: %i, chosen: %i!", 
784                              GetNActiveChannels(fPayloadCurr[-1]), GetNActiveChannelsFromMask(fPayloadCurr[-1]), channelcount_exp) ));
785     }
786     AliDebug(2, Form("expecting %i active channels, timebins: %i", channelcount_exp, fCurrNtimebins));
787     
788     // ----- reading marked ADC channels -----
789     while (channelcount < channelcount_exp && *(fPayloadCurr) != fgkDataEndmarker) {
790       if (channelno < 21)
791         channelno++;
792       while (channelno < 21 && (channelmask & 1 << channelno) == 0)
793         channelno++;
794       
795       if (fCurrNtimebins > 30) {
796         currentTimebin = ((*fPayloadCurr >> 2) & 0x3f);
797         timebins = ((*fPayloadCurr >> 8) & 0xf) * 3;
798       } 
799       else {
800         currentTimebin = 0;
801       }
802       
803       adcwc = 0;
804       AliDebug(2, Form("Now looking %i words", timebins / 3));
805       Int_t adccol = adccoloff - channelno;
806       Int_t padcol = padcoloff - channelno;
807       while (adcwc < timebins / 3 && 
808              *(fPayloadCurr) != fgkDataEndmarker && 
809              fPayloadCurr - fPayloadStart < fPayloadSize) {
810         int check = 0x3 & *fPayloadCurr;
811         if (channelno % 2 != 0) { // odd channel
812           if (check != 0x2 && channelno < 21) {
813             AliError(MCMError(kAdcCheckInvalid,
814                               Form("Invalid check bits: %i for %2i. ADC word in odd channel: %i", 
815                                    check, adcwc+1, channelno) ));
816           }
817         }
818         else {                  // even channel
819           if (check != 0x3 && channelno < 21) {
820             AliError(MCMError(kAdcCheckInvalid,
821                               Form("Invalid check bits: %i for %2i. ADC word in even channel: %i", 
822                                    check, adcwc+1, channelno) ));
823           }
824         }
825         
826         // filling the actual timebin data
827         int tb2 = 0x3ff & *fPayloadCurr >> 22;
828         int tb1 = 0x3ff & *fPayloadCurr >> 12;
829         int tb0 = 0x3ff & *fPayloadCurr >> 2;
830         if (adcwc != 0 || fCurrNtimebins <= 30) 
831           fAdcArray->SetDataByAdcCol(row, adccol, currentTimebin++, tb0);
832         else
833           tb0 = -1;
834         fAdcArray->SetDataByAdcCol(row, adccol, currentTimebin++, tb1);
835         fAdcArray->SetDataByAdcCol(row, adccol, currentTimebin++, tb2);
836         
837         adcwc++;
838         fPayloadCurr++;
839       }
840       
841       if (adcwc != timebins / 3) 
842         AliError(MCMError(kAdcDataAbort));
843       
844       // adding index 
845       if (padcol > 0 && padcol < 144) {
846         fSignalIndex->AddIndexRC(row, padcol);
847       }
848       
849       channelcount++;
850     }
851     
852     if (channelcount != channelcount_exp)
853       AliError(MCMError(kAdcChannelsMiss));
854     
855     mcmcount++;
856     // continue with next MCM
857   }
858
859   // check for missing MCMs (if header suppression is inactive)
860   if (fCurrMajor & 0x1 == 0 && mcmcount != mcmcount_exp) {
861     AliError(LinkError(kMissMcmHeaders,
862                        Form("No. of MCM headers %i not as expected: %i", 
863                             mcmcount, mcmcount_exp) ));
864   }
865
866   return (fPayloadCurr - start);
867 }
868
869 Int_t AliTRDrawStream::ReadNonZSData()
870 {
871   // read the non-zs data from one link from the current reading position
872   
873   UInt_t *start = fPayloadCurr;
874   
875   Int_t mcmcount = 0;
876   Int_t mcmcount_exp = fCurrStack == 2 ? 48 : 64;
877   Int_t channelcount = 0;
878   Int_t channelcount_exp = 0;
879   Int_t timebins;
880   Int_t currentTimebin = 0;
881   Int_t adcwc = 0;
882   Int_t evno = -1;
883   Int_t lastmcmpos = -1;
884   Int_t lastrobpos = -1;
885
886   if (fCurrNtimebins != fNtimebins) {
887     if (fNtimebins > 0) 
888       AliError(LinkError(kNtimebinsChanged,
889                          Form("No. of timebins changed from %i to %i", fNtimebins, fCurrNtimebins) ));
890     fNtimebins = fCurrNtimebins;
891   }
892   
893   timebins = fNtimebins;
894   
895   while (*(fPayloadCurr) != fgkDataEndmarker && 
896          fPayloadCurr - fPayloadStart < fPayloadSize - 2) {
897     
898     // ----- Checking MCM Header -----
899     AliDebug(2, Form("MCM header: 0x%08x", *fPayloadCurr));
900     
901     // ----- checking for proper readout order - ROB -----
902     if (GetROBReadoutPos(ROB(*fPayloadCurr) / 2) >= lastrobpos) {
903       lastrobpos = GetROBReadoutPos(ROB(*fPayloadCurr) / 2);
904     }
905     else {
906       AliError(ROBError(kPosUnexp));
907     }
908     fCurrRobPos = ROB(*fPayloadCurr);
909     
910     // ----- checking for proper readout order - MCM -----
911     if (GetMCMReadoutPos(MCM(*fPayloadCurr)) >= (lastmcmpos + 1) % 16) {
912       lastmcmpos = GetMCMReadoutPos(*fPayloadCurr);
913     }
914     else {
915       AliError(MCMError(kPosUnexp));
916     }
917     fCurrMcmPos = MCM(*fPayloadCurr);
918     
919     if (EvNo(*fPayloadCurr) != evno) {
920       if (evno == -1)
921         evno = EvNo(*fPayloadCurr);
922       else {
923         AliError(MCMError(kPtrgCntMismatch,
924                           Form("%i <-> %i", evno, EvNo(*fPayloadCurr)) ));
925       }
926     }
927     
928     channelcount = 0;
929     channelcount_exp = 21;
930     int channelno = -1;
931
932     Int_t adccoloff = AdcColOffset(*fPayloadCurr);
933     Int_t padcoloff = PadColOffset(*fPayloadCurr);
934     Int_t row = Row(*fPayloadCurr);
935
936     fPayloadCurr++;
937
938     // ----- reading marked ADC channels -----
939     while (channelcount < channelcount_exp && 
940            *(fPayloadCurr) != fgkDataEndmarker) {
941       if (channelno < 21)
942         channelno++;
943       
944       currentTimebin = 0;
945       
946       adcwc = 0;
947       AliDebug(2, Form("Now looking %i words", timebins / 3));
948       Int_t adccol = adccoloff - channelno;
949       Int_t padcol = padcoloff - channelno;
950       while (adcwc < timebins / 3 && 
951              *(fPayloadCurr) != fgkDataEndmarker && 
952              fPayloadCurr - fPayloadStart < fPayloadSize) {
953         int check = 0x3 & *fPayloadCurr;
954         if (channelno % 2 != 0) { // odd channel
955           if (check != 0x2 && channelno < 21) {
956             AliError(MCMError(kAdcCheckInvalid,
957                               Form("Invalid check bits: %i for %2i. ADC word in odd channel: %i", 
958                                    check, adcwc+1, channelno) ));
959           }
960         }
961         else {                  // even channel
962           if (check != 0x3 && channelno < 21) {
963             AliError(MCMError(kAdcCheckInvalid,
964                               Form("Invalid check bits: %i for %2i. ADC word in even channel: %i", 
965                                    check, adcwc+1, channelno) ));
966           }
967         }
968         
969         // filling the actual timebin data
970         int tb2 = 0x3ff & *fPayloadCurr >> 22;
971         int tb1 = 0x3ff & *fPayloadCurr >> 12;
972         int tb0 = 0x3ff & *fPayloadCurr >> 2;
973         if (adcwc != 0 || fCurrNtimebins <= 30) 
974           fAdcArray->SetDataByAdcCol(row, adccol, currentTimebin++, tb0);
975         else
976           tb0 = -1;
977         fAdcArray->SetDataByAdcCol(row, adccol, currentTimebin++, tb1);
978         fAdcArray->SetDataByAdcCol(row, adccol, currentTimebin++, tb2);
979
980         adcwc++;
981         fPayloadCurr++;
982       }
983
984       if (adcwc != timebins / 3) 
985         AliError(MCMError(kAdcDataAbort));
986       
987       // adding index 
988       if (padcol > 0 && padcol < 144) {
989         fSignalIndex->AddIndexRC(row, padcol);
990       }
991
992       channelcount++;
993     }
994
995     if (channelcount != channelcount_exp)
996       AliError(MCMError(kAdcChannelsMiss));
997     mcmcount++;
998     // continue with next MCM
999   }
1000
1001   // check for missing MCMs (if header suppression is inactive)
1002   if (mcmcount != mcmcount_exp) {
1003     AliError(LinkError(kMissMcmHeaders,
1004                        Form("%i not as expected: %i", mcmcount, mcmcount_exp) ));
1005   }
1006
1007   return (fPayloadCurr - start);
1008 }
1009
1010
1011 Int_t AliTRDrawStream::GetNActiveChannelsFromMask(UInt_t adcmask) 
1012 {
1013   adcmask = GetActiveChannels(adcmask);
1014   adcmask = adcmask - ((adcmask >> 1) & 0x55555555);
1015   adcmask = (adcmask & 0x33333333) + ((adcmask >> 2) & 0x33333333);
1016   return (((adcmask + (adcmask >> 4)) & 0xF0F0F0F) * 0x1010101) >> 24;
1017 }
1018
1019
1020 TString AliTRDrawStream::EquipmentError(ErrorCode_t err, TString msg)
1021
1022   fLastError.fSector = fCurrEquipmentId - 1024;
1023   fLastError.fStack  = -1;
1024   fLastError.fLink   = -1;
1025   fLastError.fRob    = -1;
1026   fLastError.fMcm    = -1;
1027   fLastError.fError  = err;
1028   fErrors->Fill();
1029
1030   return (Form("Event %6i: Eq. %2d - %s : ", 
1031                fRawReader->GetEventIndex(), fCurrEquipmentId, fgErrorMessages[err]) + msg); 
1032 }                                                                               
1033
1034
1035 TString AliTRDrawStream::StackError(ErrorCode_t err, TString msg)
1036
1037   fLastError.fSector = fCurrEquipmentId - 1024;
1038   fLastError.fStack  = fCurrSlot;
1039   fLastError.fLink   = -1;
1040   fLastError.fRob    = -1;
1041   fLastError.fMcm    = -1;
1042   fLastError.fError  = err;
1043   fErrors->Fill();
1044
1045   return (Form("Event %6i: Eq. %2d S %i - %s : ", 
1046                fRawReader->GetEventIndex(), fCurrEquipmentId, fCurrSlot, fgErrorMessages[err]) + msg); 
1047
1048
1049
1050 TString AliTRDrawStream::LinkError(ErrorCode_t err, TString msg)
1051
1052   fLastError.fSector = fCurrEquipmentId - 1024;
1053   fLastError.fStack  = fCurrSlot;
1054   fLastError.fLink   = fCurrLink;
1055   fLastError.fRob    = -1;
1056   fLastError.fMcm    = -1;
1057   fLastError.fError  = err;
1058   fErrors->Fill();
1059
1060   return (Form("Event %6i: Eq. %2d S %i l %2i - %s : ", 
1061                fRawReader->GetEventIndex(), fCurrEquipmentId, fCurrSlot, fCurrLink, fgErrorMessages[err]) + msg); 
1062
1063
1064
1065 TString AliTRDrawStream::ROBError(ErrorCode_t err, TString msg)
1066
1067   fLastError.fSector = fCurrEquipmentId - 1024;
1068   fLastError.fStack  = fCurrSlot;
1069   fLastError.fLink   = fCurrLink;
1070   fLastError.fRob    = fCurrRobPos;
1071   fLastError.fMcm    = -1;
1072   fLastError.fError  = err;
1073   fErrors->Fill();
1074
1075   return (Form("Event %6i: Eq. %2d S %i l %2i ROB %i - %s : ", 
1076                fRawReader->GetEventIndex(), fCurrEquipmentId, fCurrSlot, fCurrLink, fCurrRobPos, fgErrorMessages[err]) + msg); 
1077
1078
1079
1080 TString AliTRDrawStream::MCMError(ErrorCode_t err, TString msg)
1081
1082   fLastError.fSector = fCurrEquipmentId - 1024;
1083   fLastError.fStack  = fCurrSlot;
1084   fLastError.fLink   = fCurrLink;
1085   fLastError.fRob    = fCurrRobPos;
1086   fLastError.fMcm    = fCurrMcmPos;
1087   fLastError.fError  = err;
1088   fErrors->Fill();
1089
1090   return (Form("Event %6i: Eq. %2d S %i l %2i ROB %i MCM %2i - %s : ", 
1091                fRawReader->GetEventIndex(), fCurrEquipmentId, fCurrSlot, fCurrLink, fCurrRobPos, fCurrMcmPos, fgErrorMessages[err]) + msg); 
1092 }
1093
1094 const char* AliTRDrawStream::GetErrorMessage(ErrorCode_t errCode)
1095
1096   if (errCode > 0 && errCode < kLastErrorCode) 
1097     return fgErrorMessages[errCode];
1098   else 
1099     return ""; 
1100