1 /**************************************************************************
2 * Copyright(c) 1998-2003, 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 **************************************************************************/
18 Revision 1.12 2006/08/22 13:29:42 arcelli
19 removal of effective c++ warnings (C.Zampolli)
21 Revision 1.11 2006/08/10 14:46:54 decaro
22 TOF raw data format: updated version
24 Revision 1.10.1 2006/06/28 A.De Caro
25 Update TOF raw data format
26 according to the final version
27 (see ALICE internal note in preparation
28 'ALICE TOF raw data format')
30 Revision 0.02 2005/7/25 A.De Caro
31 Update number of bits allocated for time-of-flight
32 and 'charge' measurements
34 Revision 0.01 2004/6/11 A.De Caro, S.B.Sellitto, R.Silvestri
35 First implementation: global methods RawDataTOF
39 ////////////////////////////////////////////////////////////////////
41 // This class contains the methods to create the Raw Data files //
42 // for the TOF detector starting from the Digits. //
43 // In this implementation, we defined the structure //
44 // of the ALICE-TOF raw data (according to the //
45 // ALICE technical note, in preparation) //
46 // starting from the TOF digit format. //
48 ////////////////////////////////////////////////////////////////////
50 #include "Riostream.h"
53 #include "TClonesArray.h"
57 #include "AliBitPacking.h"
60 //#include "AliRawDataHeader.h"
61 #include "AliRawDataHeaderSim.h"
63 #include "AliTOFDDLRawData.h"
64 #include "AliTOFDigitMap.h"
65 #include "AliTOFdigit.h"
66 #include "AliTOFGeometry.h"
67 #include "AliTOFRawStream.h"
69 extern TRandom *gRandom;
71 ClassImp(AliTOFDDLRawData)
73 //---------------------------------------------------------------------------
74 AliTOFDDLRawData::AliTOFDDLRawData():
77 fPackedAcquisition(kTRUE),
79 fTOFdigitMap(new AliTOFDigitMap()),
81 fTOFrawStream(new AliTOFRawStream())
86 //----------------------------------------------------------------------------
87 AliTOFDDLRawData::AliTOFDDLRawData(AliTOFGeometry *tofGeom):
90 fPackedAcquisition(kTRUE),
91 fTOFgeometry(tofGeom),
92 fTOFdigitMap(new AliTOFDigitMap()),
94 fTOFrawStream(new AliTOFRawStream())
100 //----------------------------------------------------------------------------
101 AliTOFDDLRawData::AliTOFDDLRawData(const AliTOFDDLRawData &source) :
105 fPackedAcquisition(kTRUE),
107 fTOFdigitMap(new AliTOFDigitMap()),
109 fTOFrawStream(new AliTOFRawStream())
112 this->fIndex=source.fIndex;
113 this->fVerbose=source.fVerbose;
114 this->fPackedAcquisition=source.fPackedAcquisition;
115 this->fTOFgeometry=source.fTOFgeometry;
116 this->fTOFdigitMap=source.fTOFdigitMap;
117 this->fTOFdigitArray=source.fTOFdigitArray;
118 this->fTOFrawStream=source.fTOFrawStream;
122 //----------------------------------------------------------------------------
123 AliTOFDDLRawData& AliTOFDDLRawData::operator=(const AliTOFDDLRawData &source) {
125 this->fIndex=source.fIndex;
126 this->fVerbose=source.fVerbose;
127 this->fPackedAcquisition=source.fPackedAcquisition;
128 this->fTOFgeometry=source.fTOFgeometry;
129 this->fTOFdigitMap=source.fTOFdigitMap;
130 this->fTOFdigitArray=source.fTOFdigitArray;
131 this->fTOFrawStream=source.fTOFrawStream;
135 //----------------------------------------------------------------------------
136 Int_t AliTOFDDLRawData::RawDataTOF(TBranch* branch)
139 // This method creates the Raw data files for TOF detector
142 const Int_t kSize = 5000; //max number of digits per DDL file times 2
148 fTOFdigitArray = * (TClonesArray**) branch->GetAddress();
151 ofstream outfile; // logical name of the output file
153 //AliRawDataHeader header;
154 AliRawDataHeaderSim header;
156 UInt_t sizeRawData = 0;
168 UInt_t nWordsPerTRM = 0;
170 //loop over TOF DDL files
171 for (nDDL=0; nDDL<AliDAQ::NumberOfDdls("TOF"); nDDL++) {
173 strcpy(fileName,AliDAQ::DdlFileName("TOF",nDDL)); //The name of the output file
175 outfile.open(fileName,ios::binary);
177 outfile.open(fileName);
180 iDDL = fTOFrawStream->GetDDLnumberPerSector(nDDL);
182 // write Dummy DATA HEADER
183 UInt_t dataHeaderPosition = outfile.tellp();
184 outfile.write((char*)(&header),sizeof(header));
186 // DRM section: trailer
191 buf[fIndex] = MakeFiller();
196 // loop on TRM number per DRM
197 for (nTRM=AliTOFGeometry::NTRM(); nTRM>=3; nTRM--) {
201 // the slot number 3 of the even DRM (i.e. right) doesn't contain TDC digit data
202 if (iDDL%2==0 && nTRM==3) continue;
204 // TRM global trailer
207 // loop on TRM chain number per TRM
208 for (iChain=AliTOFGeometry::NChain()-1; iChain>=0; iChain--) {
211 MakeTRMchainTrailer(iChain, buf);
215 MakeTDCdigits(nDDL, nTRM, iChain, buf, nWordsPerTRM);
218 MakeTRMchainHeader(nTRM, iChain, buf);
221 } // end loop on iChain
224 MakeTRMheader(nTRM, buf);
226 // TRM filler in case where TRM data number is odd
227 if (nWordsPerTRM%2!=0) MakeTRMfiller(buf, nWordsPerTRM);
229 } // end loop on nTRM
232 MakeDRMheader(nDDL, buf);
234 ReverseArray(buf, fIndex+1);
236 outfile.write((char *)buf,((fIndex+1)*sizeof(UInt_t)));
238 for (jj=0; jj<(fIndex+1); jj++) buf[jj]=0;
241 //Write REAL DATA HEADER
242 UInt_t currentFilePosition = outfile.tellp();
243 sizeRawData = currentFilePosition - dataHeaderPosition - sizeof(header);
244 header.fSize = currentFilePosition - dataHeaderPosition;
245 header.SetAttribute(0); // valid data
246 outfile.seekp(dataHeaderPosition);
247 outfile.write((char*)(&header),sizeof(header));
248 outfile.seekp(currentFilePosition);
252 } //end loop on DDL file number
258 //----------------------------------------------------------------------------
259 void AliTOFDDLRawData::GetDigits()
262 // Fill the TOF volumes' map with the TOF digit indices
265 Int_t vol[5] = {-1,-1,-1,-1,-1};
268 Int_t ndigits = fTOFdigitArray->GetEntries();
272 // loop on TOF digits
273 for (digit=0; digit<ndigits; digit++) {
274 digs = (AliTOFdigit*)fTOFdigitArray->UncheckedAt(digit);
276 vol[0] = digs->GetSector(); // Sector Number (0-17)
277 vol[1] = digs->GetPlate(); // Plate Number (0-4)
278 vol[2] = digs->GetStrip(); // Strip Number (0-14/18)
279 vol[3] = digs->GetPadx(); // Pad Number in x direction (0-47)
280 vol[4] = digs->GetPadz(); // Pad Number in z direction (0-1)
282 fTOFdigitMap->AddDigit(vol, digit);
284 } // close loop on digit del TOF
288 //----------------------------------------------------------------------------
289 void AliTOFDDLRawData::MakeDRMheader(Int_t nDDL, UInt_t *buf)
295 Int_t iDDL = fTOFrawStream->GetDDLnumberPerSector(nDDL);
297 Int_t iSector = fTOFrawStream->GetSectorNumber(nDDL);
304 word = 1; // 0001 -> DRM data are coming from the VME slot number 1
305 AliBitPacking::PackWord(word,baseWord, 0, 3);
306 word = 0; // event CRC
307 AliBitPacking::PackWord(word,baseWord, 4,19);
308 word = 0; // reserved for future use
309 AliBitPacking::PackWord(word,baseWord,20,27);
310 word = 4; // 0100 -> DRM header ID
311 AliBitPacking::PackWord(word,baseWord,28,31);
313 buf[fIndex]=baseWord;
315 // DRM status header 3
317 word = 1; // 0001 -> DRM data are coming from the VME slot number 1
318 AliBitPacking::PackWord(word,baseWord, 0, 3);
319 word = 0; // TTC event counter
320 AliBitPacking::PackWord(word,baseWord, 4,27);
321 word = 4; // 0100 -> DRM header ID
322 AliBitPacking::PackWord(word,baseWord,28,31);
324 buf[fIndex]=baseWord;
326 // DRM status header 2
328 word = 1; // 0001 -> DRM data are coming from the VME slot number 1
329 AliBitPacking::PackWord(word,baseWord, 0, 3);
332 word = 2047; // enable ID: [00000000000;11111111111] for odd
333 // (i.e. right) crates
334 AliBitPacking::PackWord(word,baseWord, 4,14);
336 word = 2045; // enable ID: [00000000000;11111111101] for even
337 // (i.e. left) crates
338 AliBitPacking::PackWord(word,baseWord, 4,14);
342 AliBitPacking::PackWord(word,baseWord,15,15);
343 word = 0; // fault ID
344 AliBitPacking::PackWord(word,baseWord,16,27);
345 word = 4; // 0100 -> DRM header ID
346 AliBitPacking::PackWord(word,baseWord,28,31);
348 buf[fIndex]=baseWord;
350 // DRM status header 1
352 word = 1; // 0001 -> DRM data are coming from the VME slot number 1
353 AliBitPacking::PackWord(word,baseWord, 0, 3);
356 word = 2047; // slot ID: [00000000000;11111111111] for odd
357 // (i.e. right) crates
358 AliBitPacking::PackWord(word,baseWord, 4,14);
360 word = 2045; // slot ID: [00000000000;11111111101] for even
361 // (i.e. left) crates
362 AliBitPacking::PackWord(word,baseWord, 4,14);
365 word = 1; // LHC clock status: 1/0
366 AliBitPacking::PackWord(word,baseWord,15,15);
367 word = 0; // reserved for future use
368 AliBitPacking::PackWord(word,baseWord,16,27);
369 word = 4; // 0100 -> DRM header ID
370 AliBitPacking::PackWord(word,baseWord,28,31);
372 buf[fIndex]=baseWord;
376 word = 1; // 0001 -> DRM data are coming from the VME slot number 1
377 AliBitPacking::PackWord(word,baseWord, 0, 3);
378 word = fIndex+1 + 1; // event words
379 AliBitPacking::PackWord(word,baseWord, 4,20);
380 word = iDDL; // crate ID [0;3]
381 AliBitPacking::PackWord(word,baseWord,21,22);
382 word = iSector; // sector ID [0;17]
383 AliBitPacking::PackWord(word,baseWord,23,27);
384 word = 4; // 0100 -> DRM header ID
385 AliBitPacking::PackWord(word,baseWord,28,31);
387 buf[fIndex]=baseWord;
391 //----------------------------------------------------------------------------
392 void AliTOFDDLRawData::MakeDRMtrailer(UInt_t *buf)
395 // DRM global trailer
402 word = 1; // 0001 -> DRM data are coming from the VME slot number 1
403 AliBitPacking::PackWord(word,baseWord, 0, 3);
404 word = 0; // local event counter
405 AliBitPacking::PackWord(word,baseWord, 4,15);
406 word = 0; // reserved for future use
407 AliBitPacking::PackWord(word,baseWord,16,27);
408 word = 5; // 0101 -> DRM trailer ID
409 AliBitPacking::PackWord(word,baseWord,28,31);
411 buf[fIndex]=baseWord;
415 //----------------------------------------------------------------------------
416 void AliTOFDDLRawData::MakeLTMheader(UInt_t *buf)
426 word = 2; // 0010 -> LTM data are coming from the VME slot number 2
427 AliBitPacking::PackWord(word,baseWord, 0, 3);
428 word = 35; // event words
429 AliBitPacking::PackWord(word,baseWord, 4,16);
430 word = 0; // crc error
431 AliBitPacking::PackWord(word,baseWord,17,17);
433 AliBitPacking::PackWord(word,baseWord,18,23);
435 AliBitPacking::PackWord(word,baseWord,24,27);
436 word = 4; // 0100 -> LTM header ID
437 AliBitPacking::PackWord(word,baseWord,28,31);
439 buf[fIndex]=baseWord;
443 //----------------------------------------------------------------------------
444 void AliTOFDDLRawData::MakeLTMdata(UInt_t *buf)
454 word = 100; // Local temperature in LTM5 -> 4 X 25 degree (environment temperature)
455 AliBitPacking::PackWord(word,baseWord, 0, 9);
456 word = 100; // Local temperature in LTM6 -> 4 X 25 degree (environment temperature)
457 AliBitPacking::PackWord(word,baseWord,10,19);
458 word = 100; // Local temperature in LTM7 -> 4 X 25 degree (environment temperature)
459 AliBitPacking::PackWord(word,baseWord,20,29);
461 AliBitPacking::PackWord(word,baseWord,30,31);
464 buf[fIndex]=baseWord;
466 // Local temperature in LTM2, LMT3, LTM4 -> 4 X 25 degree (environment temperature)
468 buf[fIndex]=baseWord;
470 // Local temperature in FEAC7, FEAC8, LTM1 -> 4 X 25 degree (environment temperature)
472 buf[fIndex]=baseWord;
474 // Local temperature in FEAC4, FEAC5, FEAC6 -> 4 X 25 degree (environment temperature)
476 buf[fIndex]=baseWord;
478 // Local temperature in FEAC1, FEAC2, FEAC3 -> 4 X 25 degree (environment temperature)
480 buf[fIndex]=baseWord;
483 word = 0; // GND-FEAC15 -> Voltage drop between GND and FEAC15
484 AliBitPacking::PackWord(word,baseWord, 0, 9);
485 word = 0; // VTH16 -> Thereshould voltage for FEAC16
486 AliBitPacking::PackWord(word,baseWord,10,19);
487 word = 0; // GND-FEAC16 -> Voltage drop between GND and FEAC16
488 AliBitPacking::PackWord(word,baseWord,20,29);
490 AliBitPacking::PackWord(word,baseWord,30,31);
493 buf[fIndex]=baseWord;
495 // VTH14 -> Thereshould voltage for FEAC14
496 // GND-FEAC14 -> Voltage drop between GND and FEAC14
497 // VTH15 -> Thereshould voltage for FEAC15
499 buf[fIndex]=baseWord;
501 // GND-FEAC12 -> Voltage drop between GND and FEAC12
502 // VTH13 -> Thereshould voltage for FEAC13
503 // GND-FEAC13 -> Voltage drop between GND and FEAC13
505 buf[fIndex]=baseWord;
507 // VTH11 -> Thereshould voltage for FEAC11
508 // GND-FEAC11 -> Voltage drop between GND and FEAC11
509 // VTH12 -> Thereshould voltage for FEAC12
511 buf[fIndex]=baseWord;
513 // GND-FEAC9 -> Voltage drop between GND and FEAC9
514 // VTH10 -> Thereshould voltage for FEAC10
515 // GND-FEAC10 -> Voltage drop between GND and FEAC10
517 buf[fIndex]=baseWord;
519 // VTH8 -> Thereshould voltage for FEAC8
520 // GND-FEAC8 -> Voltage drop between GND and FEAC8
521 // VTH9 -> Thereshould voltage for FEAC9
523 buf[fIndex]=baseWord;
525 // GND-FEAC6 -> Voltage drop between GND and FEAC6
526 // VTH7 -> Thereshould voltage for FEAC7
527 // GND-FEAC7 -> Voltage drop between GND and FEAC7
529 buf[fIndex]=baseWord;
531 // VTH5 -> Thereshould voltage for FEAC5
532 // GND-FEAC5 -> Voltage drop between GND and FEAC5
533 // VTH6 -> Thereshould voltage for FEAC6
535 buf[fIndex]=baseWord;
537 // GND-FEAC3 -> Voltage drop between GND and FEAC3
538 // VTH4 -> Thereshould voltage for FEAC4
539 // GND-FEAC4 -> Voltage drop between GND and FEAC4
541 buf[fIndex]=baseWord;
543 // VTH2 -> Thereshould voltage for FEAC2
544 // GND-FEAC2 -> Voltage drop between GND and FEAC2
545 // VTH3 -> Thereshould voltage for FEAC3
547 buf[fIndex]=baseWord;
549 // LV16 -> Low Voltage measured by FEAC16
550 // GND-FEAC1 -> Voltage drop between GND and FEAC1
551 // VTH1 -> Thereshould voltage for FEAC1
553 buf[fIndex]=baseWord;
555 // Low Voltage measured by FEAC13, FEAC14, FEAC15
557 buf[fIndex]=baseWord;
559 // Low Voltage measured by FEAC10, FEAC11, FEAC12
561 buf[fIndex]=baseWord;
563 // Low Voltage measured by FEAC7, FEAC8, FEAC9
565 buf[fIndex]=baseWord;
567 // Low Voltage measured by FEAC4, FEAC5, FEAC6
569 buf[fIndex]=baseWord;
571 // Low Voltage measured by FEAC1, FEAC2, FEAC3
573 buf[fIndex]=baseWord;
577 word = 0; // PDL45 -> Delay Line setting for PDL45
578 AliBitPacking::PackWord(word,baseWord, 0, 7);
579 word = 0; // PDL46 -> Delay Line setting for PDL46
580 AliBitPacking::PackWord(word,baseWord, 8,15);
581 word = 0; // PDL47 -> Delay Line setting for PDL47
582 AliBitPacking::PackWord(word,baseWord,16,23);
583 word = 0; // PDL48 -> Delay Line setting for PDL48
584 AliBitPacking::PackWord(word,baseWord,24,31);
586 buf[fIndex]=baseWord;
588 // Delay Line setting for PDL41, PDL42, PDL43, PDL44
590 buf[fIndex]=baseWord;
592 // Delay Line setting for PDL37, PDL38, PDL39, PDL40
594 buf[fIndex]=baseWord;
596 // Delay Line setting for PDL33, PDL34, PDL35, PDL36
598 buf[fIndex]=baseWord;
600 // Delay Line setting for PDL29, PDL30, PDL31, PDL32
602 buf[fIndex]=baseWord;
604 // Delay Line setting for PDL25, PDL26, PDL27, PDL28
606 buf[fIndex]=baseWord;
608 // Delay Line setting for PDL21, PDL22, PDL23, PDL24
610 buf[fIndex]=baseWord;
612 // Delay Line setting for PDL17, PDL18, PDL19, PDL20
614 buf[fIndex]=baseWord;
616 // Delay Line setting for PDL13, PDL14, PDL15, PDL16
618 buf[fIndex]=baseWord;
620 // Delay Line setting for PDL9, PDL10, PDL11, PDL12
622 buf[fIndex]=baseWord;
624 // Delay Line setting for PDL5, PDL6, PDL7, PDL8
626 buf[fIndex]=baseWord;
628 // Delay Line setting for PDL1, PDL2, PDL3, PDL4
630 buf[fIndex]=baseWord;
634 //----------------------------------------------------------------------------
635 void AliTOFDDLRawData::MakeLTMtrailer(UInt_t *buf)
645 word = 2; // 0010 -> LTM data are coming from the VME slot number 2
646 AliBitPacking::PackWord(word,baseWord, 0, 3);
647 word = 0; // event crc
648 AliBitPacking::PackWord(word,baseWord, 4,15);
649 word = 0; // event number
650 AliBitPacking::PackWord(word,baseWord,16,27);
651 word = 5; // 0101 -> LTM trailer ID
652 AliBitPacking::PackWord(word,baseWord,28,31);
654 buf[fIndex]=baseWord;
658 //----------------------------------------------------------------------------
659 void AliTOFDDLRawData::MakeTRMheader(Int_t nTRM, UInt_t *buf)
662 // TRM header for the TRM number nTRM [ 3;12]
665 if (nTRM<3 || nTRM>12) {
666 AliWarning(Form(" TRM number is out of the right range [3;12] (nTRM = %3i",nTRM));
674 word = nTRM; // TRM data coming from the VME slot number nTRM
675 AliBitPacking::PackWord(word,baseWord, 0, 3);
676 word = 0; // event words
677 AliBitPacking::PackWord(word,baseWord, 4,16);
679 if (fPackedAcquisition)
680 word = 0; // ACQuisition mode: [0;3] see document
682 word = 3; // ACQuisition mode: [0;3] see document
683 AliBitPacking::PackWord(word,baseWord,17,18);
684 word = 0; // description of a SEU inside LUT tables for INL compensation;
685 // the data are unaffected
686 AliBitPacking::PackWord(word,baseWord,19,19);
687 word = 0; // Must Be Zero (MBZ)
688 AliBitPacking::PackWord(word,baseWord,20,27);
689 word = 4; // 0100 -> TRM header ID
690 AliBitPacking::PackWord(word,baseWord,28,31);
692 buf[fIndex]=baseWord;
696 //----------------------------------------------------------------------------
697 void AliTOFDDLRawData::MakeTRMtrailer(UInt_t *buf)
707 word = 15; // 1111 -> TRM trailer ID 1
708 AliBitPacking::PackWord(word,baseWord, 0, 3);
709 word = 0; // event CRC
710 AliBitPacking::PackWord(word,baseWord, 4,15);
711 word = 0; // local event counter == DRM local event counter
712 AliBitPacking::PackWord(word,baseWord,16,27);
713 word = 5; // 0101 -> TRM trailer ID 2
714 AliBitPacking::PackWord(word,baseWord,28,31);
716 buf[fIndex]=baseWord;
720 //----------------------------------------------------------------------------
721 void AliTOFDDLRawData::MakeTRMchainHeader(Int_t nTRM, Int_t iChain,
731 if (nTRM<3 || nTRM>12) {
732 AliWarning(Form(" TRM number is out of the right range [3;12] (nTRM = %3i", nTRM));
736 if (iChain<0 || iChain>1) {
737 AliWarning(Form(" Chain number is out of the right range [0;1] (iChain = %3i", iChain));
742 word = nTRM; // TRM data coming from the VME slot ID nTRM
743 AliBitPacking::PackWord(word,baseWord, 0, 3);
744 word = 0; // bunch ID
745 AliBitPacking::PackWord(word,baseWord, 4,15);
746 word = 100; // PB24 temperature -> 4 X 25 degree (environment temperature)
747 AliBitPacking::PackWord(word,baseWord,16,23);
748 word = (Int_t)(5 * gRandom->Rndm()); // PB24 ID [0;4]
749 AliBitPacking::PackWord(word,baseWord,24,26);
751 AliBitPacking::PackWord(word,baseWord,27,27);
754 word = 0; // 0000 -> TRM chain 0 ID
757 word = 2; // 0010 -> TRM chain 1 ID
760 AliBitPacking::PackWord(word,baseWord,28,31);
762 buf[fIndex]=baseWord;
766 //----------------------------------------------------------------------------
767 void AliTOFDDLRawData::MakeTRMfiller(UInt_t *buf, UInt_t nWordsPerTRM)
776 for (jj=fIndex; jj>fIndex-(Int_t)nWordsPerTRM; jj--) {
780 buf[fIndex-nWordsPerTRM] = MakeFiller();
784 //----------------------------------------------------------------------------
785 UInt_t AliTOFDDLRawData::MakeFiller()
788 // Filler word definition: to make even the number of words per TRM/LTM
795 word = 0; // 0000 -> filler ID 1
796 AliBitPacking::PackWord(word,baseWord, 0, 3);
798 AliBitPacking::PackWord(word,baseWord, 4,27);
799 word = 7; // 0111 -> filler ID 2
800 AliBitPacking::PackWord(word,baseWord, 28,31);
806 //----------------------------------------------------------------------------
807 void AliTOFDDLRawData::MakeTRMchainTrailer(Int_t iChain, UInt_t *buf)
813 if (iChain<0 || iChain>1) {
814 AliWarning(Form(" Chain number is out of the right range [0;1] (iChain = %3i", iChain));
823 AliBitPacking::PackWord(word,baseWord, 0, 3);
825 AliBitPacking::PackWord(word,baseWord, 4,15);
826 word = 0; // event counter
827 AliBitPacking::PackWord(word,baseWord,16,27);
830 word = 1; // 0001 -> TRM chain 0 trailer ID
833 word = 3; // 0011 -> TRM chain 1 trailer ID
836 AliBitPacking::PackWord(word,baseWord,28,31);
838 buf[fIndex]=baseWord;
842 //----------------------------------------------------------------------------
843 void AliTOFDDLRawData::MakeTDCdigits(Int_t nDDL, Int_t nTRM, Int_t iChain,
844 UInt_t *buf, UInt_t &nWordsPerTRM)
850 const Double_t kOneMoreFilledCell = 1./(fTOFgeometry->NPadXSector()*fTOFgeometry->NSectors());
851 Double_t percentFilledCells = Double_t(fTOFdigitMap->GetFilledCellNumber())/(fTOFgeometry->NPadXSector()*fTOFgeometry->NSectors());
853 if (nDDL<0 || nDDL>71) {
854 AliWarning(Form(" DDL number is out of the right range [0;71] (nDDL = %3i", nDDL));
858 if (nTRM<3 || nTRM>12) {
859 AliWarning(Form(" TRM number is out of the right range [3;12] (nTRM = %3i", nTRM));
863 if (iChain<0 || iChain>1) {
864 AliWarning(Form(" Chain number is out of the right range [0;1] (iChain = %3i", iChain));
869 UInt_t localBuffer[1000];
870 Int_t localIndex = -1;
872 Int_t iDDL = nDDL%AliTOFGeometry::NDDL();
874 Int_t volume[5] = {-1, -1, -1, -1, -1};
875 Int_t indexDigit[3] = {-1, -1, -1};
877 Int_t totCharge = -1;
878 Int_t timeOfFlight = -1;
880 Int_t trailingSpurious = -1;
881 Int_t leadingSpurious = -1;
894 if (fVerbose==2) ftxt.open("TOFdigits.txt",ios::app);
896 for (jj=0; jj<5; jj++) volume[jj] = -1;
898 // loop on TDC number
899 for (nTDC=AliTOFGeometry::NTdc()-1; nTDC>=0; nTDC--) {
901 // the DRM odd (i.e. left) slot number 3 doesn't contain TDC digit data
902 // for TDC numbers 3-14
903 if (iDDL%2==1 && nTRM==3 && (Int_t)(nTDC/3.)!=0) continue;
905 // loop on TDC channel number
906 for (iCH=AliTOFGeometry::NCh()-1; iCH>=0; iCH--) {
908 fTOFrawStream->EquipmentId2VolumeId(nDDL, nTRM, iChain, nTDC, iCH, volume);
910 if (volume[0]==-1 || volume[1]==-1 || volume[2]==-1 ||
911 volume[3]==-1 || volume[4]==-1) continue;
913 for (jj=0; jj<3; jj++) indexDigit[jj] = -1;
915 fTOFdigitMap->GetDigitIndex(volume, indexDigit);
917 if (indexDigit[0]<0) {
919 trailingSpurious = Int_t(8192*gRandom->Rndm()) + Int_t(Int_t(256*gRandom->Rndm())*AliTOFGeometry::ToTBinWidth()/AliTOFGeometry::TdcBinWidth());
920 leadingSpurious = Int_t(8192*gRandom->Rndm());
922 if ( ( fPackedAcquisition && percentFilledCells<0.12 && gRandom->Rndm()<(0.12-percentFilledCells)) ||
923 (!fPackedAcquisition && percentFilledCells<0.24 && gRandom->Rndm()<(0.24-percentFilledCells)) ) {
925 percentFilledCells+=kOneMoreFilledCell;
930 word = trailingSpurious; // trailing edge measurement
934 word = leadingSpurious; // leading edge measurement
939 if (nDDL<10) ftxt << " " << nDDL;
940 else ftxt << " " << nDDL;
941 if (nTRM<10) ftxt << " " << nTRM;
942 else ftxt << " " << nTRM;
943 ftxt << " " << iChain;
944 if (nTDC<10) ftxt << " " << nTDC;
945 else ftxt << " " << nTDC;
947 if (volume[0]<10) ftxt << " -> " << volume[0];
948 else ftxt << " -> " << volume[0];
949 ftxt << " " << volume[1];
950 if (volume[2]<10) ftxt << " " << volume[2];
951 else ftxt << " " << volume[2];
952 ftxt << " " << volume[4];
953 if (volume[3]<10) ftxt << " " << volume[3];
954 else ftxt << " " << volume[3];
956 if (word<10) ftxt << " " << word;// << endl;
957 else if (word>=10 && word<100) ftxt << " " << word;// << endl;
958 else if (word>=100 && word<1000) ftxt << " " << word;// << endl;
959 else ftxt << " " << word;// << endl;
960 ftxt << " " << dummyPS << endl;
963 AliBitPacking::PackWord(word,baseWord, 0,20);
964 word = iCH; // TDC channel ID [0;7]
965 AliBitPacking::PackWord(word,baseWord,21,23);
966 word = nTDC; // TDC ID [0;14]
967 AliBitPacking::PackWord(word,baseWord,24,27);
968 word = 0; // error flag
969 AliBitPacking::PackWord(word,baseWord,28,28);
970 word = dummyPS; // Packing Status [0;3]
971 AliBitPacking::PackWord(word,baseWord,29,30);
972 word = 1; // TRM TDC digit ID
973 AliBitPacking::PackWord(word,baseWord,31,31);
976 localBuffer[localIndex]=baseWord;
977 psArray[localIndex]=dummyPS;
982 } // if (fPackedAcquisition && percentFilledCells<0.12 && gRandom->Rndm()<(0.12-percentFilledCells)) or ...
984 } // if (indexDigit[0]<0)
986 for (jj=0; jj<3;jj++) {
988 if (indexDigit[jj]<0) continue;
989 digs = (AliTOFdigit*)fTOFdigitArray->UncheckedAt(indexDigit[jj]);
991 if (digs->GetSector()!=volume[0] ||
992 digs->GetPlate() !=volume[1] ||
993 digs->GetStrip() !=volume[2] ||
994 digs->GetPadx() !=volume[3] ||
995 digs->GetPadz() !=volume[4]) AliWarning(" --- ERROR --- ");
997 timeOfFlight = (Int_t)(digs->GetTdc())%8192;
998 totCharge = (Int_t)(digs->GetToT());//digs->GetAdc();
1000 if (totCharge<0) totCharge = TMath::Abs(totCharge);
1001 if (totCharge>=256) totCharge = 255;
1003 if (fPackedAcquisition) {
1006 if (nDDL<10) ftxt << " " << nDDL;
1007 else ftxt << " " << nDDL;
1008 if (nTRM<10) ftxt << " " << nTRM;
1009 else ftxt << " " << nTRM;
1010 ftxt << " " << iChain;
1011 if (nTDC<10) ftxt << " " << nTDC;
1012 else ftxt << " " << nTDC;
1014 if (volume[0]<10) ftxt << " -> " << volume[0];
1015 else ftxt << " -> " << volume[0];
1016 ftxt << " " << volume[1];
1017 if (volume[2]<10) ftxt << " " << volume[2];
1018 else ftxt << " " << volume[2];
1019 ftxt << " " << volume[4];
1020 if (volume[3]<10) ftxt << " " << volume[3];
1021 else ftxt << " " << volume[3];
1022 if (totCharge<10) ftxt << " " << totCharge;
1023 else if (totCharge>=10 && totCharge<100) ftxt << " " << totCharge;
1024 else ftxt << " " << totCharge;
1025 if (timeOfFlight<10) ftxt << " " << timeOfFlight << endl;
1026 else if (timeOfFlight>=10 && timeOfFlight<100) ftxt << " " << timeOfFlight << endl;
1027 else if (timeOfFlight>=100 && timeOfFlight<1000) ftxt << " " << timeOfFlight << endl;
1028 else ftxt << " " << timeOfFlight << endl;
1031 word = timeOfFlight; // time-of-fligth measurement
1032 AliBitPacking::PackWord(word,baseWord, 0,12);
1034 word = totCharge; // time-over-threshould measurement
1035 AliBitPacking::PackWord(word,baseWord,13,20);
1037 word = iCH; // TDC channel ID [0;7]
1038 AliBitPacking::PackWord(word,baseWord,21,23);
1039 word = nTDC; // TDC ID [0;14]
1040 AliBitPacking::PackWord(word,baseWord,24,27);
1041 word = 0; // error flag
1042 AliBitPacking::PackWord(word,baseWord,28,28);
1043 word = 0; // Packing Status [0;3]
1044 AliBitPacking::PackWord(word,baseWord,29,30);
1045 word = 1; // TRM TDC digit ID
1046 AliBitPacking::PackWord(word,baseWord,31,31);
1049 localBuffer[localIndex]=baseWord;
1054 if (percentFilledCells<0.12 && gRandom->Rndm()<(0.12-percentFilledCells)) {
1056 percentFilledCells+=kOneMoreFilledCell;
1058 trailingSpurious = Int_t(8192*gRandom->Rndm()) + Int_t(Int_t(256*gRandom->Rndm())*AliTOFGeometry::ToTBinWidth()/AliTOFGeometry::TdcBinWidth());
1059 leadingSpurious = Int_t(8192*gRandom->Rndm());
1064 word = trailingSpurious; // trailing edge measurement
1068 word = leadingSpurious; // leading edge measurement
1073 if (nDDL<10) ftxt << " " << nDDL;
1074 else ftxt << " " << nDDL;
1075 if (nTRM<10) ftxt << " " << nTRM;
1076 else ftxt << " " << nTRM;
1077 ftxt << " " << iChain;
1078 if (nTDC<10) ftxt << " " << nTDC;
1079 else ftxt << " " << nTDC;
1081 if (volume[0]<10) ftxt << " -> " << volume[0];
1082 else ftxt << " -> " << volume[0];
1083 ftxt << " " << volume[1];
1084 if (volume[2]<10) ftxt << " " << volume[2];
1085 else ftxt << " " << volume[2];
1086 ftxt << " " << volume[4];
1087 if (volume[3]<10) ftxt << " " << volume[3];
1088 else ftxt << " " << volume[3];
1090 if (word<10) ftxt << " " << word;
1091 else if (word>=10 && word<100) ftxt << " " << word;
1092 else if (word>=100 && word<1000) ftxt << " " << word;
1093 else ftxt << " " << word;
1094 ftxt << " " << dummyPS << endl;
1097 AliBitPacking::PackWord(word,baseWord, 0,20);
1098 word = iCH; // TDC channel ID [0;7]
1099 AliBitPacking::PackWord(word,baseWord,21,23);
1100 word = nTDC; // TDC ID [0;14]
1101 AliBitPacking::PackWord(word,baseWord,24,27);
1102 word = 0; // error flag
1103 AliBitPacking::PackWord(word,baseWord,28,28);
1104 word = dummyPS; // Packing Status [0;3]
1105 AliBitPacking::PackWord(word,baseWord,29,30);
1106 word = 1; // TRM TDC digit ID
1107 AliBitPacking::PackWord(word,baseWord,31,31);
1110 localBuffer[localIndex]=baseWord;
1111 psArray[localIndex]=dummyPS;
1116 } // if (percentFilledCells<0.12 && gRandom->Rndm()<(0.12-percentFilledCells))
1119 } // if (fPackedAcquisition)
1120 else { // if (!fPackedAcquisition)
1122 if (percentFilledCells<0.24 && gRandom->Rndm()<(0.24-percentFilledCells) && HeadOrTail()) {
1124 percentFilledCells+=kOneMoreFilledCell;
1126 trailingSpurious = Int_t(8192*gRandom->Rndm()) + Int_t(Int_t(256*gRandom->Rndm())*AliTOFGeometry::ToTBinWidth()/AliTOFGeometry::TdcBinWidth());
1127 word = trailingSpurious;
1131 if (nDDL<10) ftxt << " " << nDDL;
1132 else ftxt << " " << nDDL;
1133 if (nTRM<10) ftxt << " " << nTRM;
1134 else ftxt << " " << nTRM;
1135 ftxt << " " << iChain;
1136 if (nTDC<10) ftxt << " " << nTDC;
1137 else ftxt << " " << nTDC;
1139 if (volume[0]<10) ftxt << " -> " << volume[0];
1140 else ftxt << " -> " << volume[0];
1141 ftxt << " " << volume[1];
1142 if (volume[2]<10) ftxt << " " << volume[2];
1143 else ftxt << " " << volume[2];
1144 ftxt << " " << volume[4];
1145 if (volume[3]<10) ftxt << " " << volume[3];
1146 else ftxt << " " << volume[3];
1148 if (word<10) ftxt << " " << word;
1149 else if (word>=10 && word<100) ftxt << " " << word;
1150 else if (word>=100 && word<1000) ftxt << " " << word;
1151 else ftxt << " " << word;
1152 ftxt << " " << dummyPS << endl;
1155 AliBitPacking::PackWord(word,baseWord, 0,20);
1156 word = iCH; // TDC channel ID [0;7]
1157 AliBitPacking::PackWord(word,baseWord,21,23);
1158 word = nTDC; // TDC ID [0;14]
1159 AliBitPacking::PackWord(word,baseWord,24,27);
1160 word = 0; // error flag
1161 AliBitPacking::PackWord(word,baseWord,28,28);
1162 word = dummyPS; // Packing Status [0;3]
1163 AliBitPacking::PackWord(word,baseWord,29,30);
1164 word = 1; // TRM TDC digit ID
1165 AliBitPacking::PackWord(word,baseWord,31,31);
1168 localBuffer[localIndex]=baseWord;
1169 psArray[localIndex]=dummyPS;
1174 } // if (percentFilledCells<0.24 && gRandom->Rndm()<(0.24-percentFilledCells))
1177 word = timeOfFlight + Int_t(totCharge*AliTOFGeometry::ToTBinWidth()/AliTOFGeometry::TdcBinWidth()); // trailing edge measurement
1180 if (nDDL<10) ftxt << " " << nDDL;
1181 else ftxt << " " << nDDL;
1182 if (nTRM<10) ftxt << " " << nTRM;
1183 else ftxt << " " << nTRM;
1184 ftxt << " " << iChain;
1185 if (nTDC<10) ftxt << " " << nTDC;
1186 else ftxt << " " << nTDC;
1188 if (volume[0]<10) ftxt << " -> " << volume[0];
1189 else ftxt << " -> " << volume[0];
1190 ftxt << " " << volume[1];
1191 if (volume[2]<10) ftxt << " " << volume[2];
1192 else ftxt << " " << volume[2];
1193 ftxt << " " << volume[4];
1194 if (volume[3]<10) ftxt << " " << volume[3];
1195 else ftxt << " " << volume[3];
1197 if (word<10) ftxt << " " << word;
1198 else if (word>=10 && word<100) ftxt << " " << word;
1199 else if (word>=100 && word<1000) ftxt << " " << word;
1200 else ftxt << " " << word;
1201 ftxt << " " << 2 << endl;
1204 AliBitPacking::PackWord(word,baseWord, 0,20);
1206 word = iCH; // TDC channel ID [0;7]
1207 AliBitPacking::PackWord(word,baseWord,21,23);
1208 word = nTDC; // TDC ID [0;14]
1209 AliBitPacking::PackWord(word,baseWord,24,27);
1210 word = 0; // error flag
1211 AliBitPacking::PackWord(word,baseWord,28,28);
1212 word = 2; // Packing Status [0;3]
1213 AliBitPacking::PackWord(word,baseWord,29,30);
1214 word = 1; // TRM TDC digit ID
1215 AliBitPacking::PackWord(word,baseWord,31,31);
1218 localBuffer[localIndex]=baseWord;
1219 psArray[localIndex]=2;
1226 word = timeOfFlight; // leading edge measurement
1229 if (nDDL<10) ftxt << " " << nDDL;
1230 else ftxt << " " << nDDL;
1231 if (nTRM<10) ftxt << " " << nTRM;
1232 else ftxt << " " << nTRM;
1233 ftxt << " " << iChain;
1234 if (nTDC<10) ftxt << " " << nTDC;
1235 else ftxt << " " << nTDC;
1237 if (volume[0]<10) ftxt << " -> " << volume[0];
1238 else ftxt << " -> " << volume[0];
1239 ftxt << " " << volume[1];
1240 if (volume[2]<10) ftxt << " " << volume[2];
1241 else ftxt << " " << volume[2];
1242 ftxt << " " << volume[4];
1243 if (volume[3]<10) ftxt << " " << volume[3];
1244 else ftxt << " " << volume[3];
1246 if (word<10) ftxt << " " << word;
1247 else if (word>=10 && word<100) ftxt << " " << word;
1248 else ftxt << " " << word;
1249 ftxt << " " << 1 << endl;
1252 AliBitPacking::PackWord(word,baseWord, 0,20);
1254 word = iCH; // TDC channel ID [0;7]
1255 AliBitPacking::PackWord(word,baseWord,21,23);
1256 word = nTDC; // TDC ID [0;14]
1257 AliBitPacking::PackWord(word,baseWord,24,27);
1258 word = 0; // error flag
1259 AliBitPacking::PackWord(word,baseWord,28,28);
1260 word = 1; // Packing Status [0;3]
1261 AliBitPacking::PackWord(word,baseWord,29,30);
1262 word = 1; // TRM TDC digit ID
1263 AliBitPacking::PackWord(word,baseWord,31,31);
1266 localBuffer[localIndex]=baseWord;
1267 psArray[localIndex]=1;
1273 if (percentFilledCells<0.24 && gRandom->Rndm()<(0.24-percentFilledCells) && !HeadOrTail()) {
1275 percentFilledCells+=kOneMoreFilledCell;
1277 leadingSpurious = Int_t(8192*gRandom->Rndm());
1278 word = leadingSpurious;
1282 if (nDDL<10) ftxt << " " << nDDL;
1283 else ftxt << " " << nDDL;
1284 if (nTRM<10) ftxt << " " << nTRM;
1285 else ftxt << " " << nTRM;
1286 ftxt << " " << iChain;
1287 if (nTDC<10) ftxt << " " << nTDC;
1288 else ftxt << " " << nTDC;
1290 if (volume[0]<10) ftxt << " -> " << volume[0];
1291 else ftxt << " -> " << volume[0];
1292 ftxt << " " << volume[1];
1293 if (volume[2]<10) ftxt << " " << volume[2];
1294 else ftxt << " " << volume[2];
1295 ftxt << " " << volume[4];
1296 if (volume[3]<10) ftxt << " " << volume[3];
1297 else ftxt << " " << volume[3];
1299 if (word<10) ftxt << " " << word;// << endl;
1300 else if (word>=10 && word<100) ftxt << " " << word;// << endl;
1301 else if (word>=100 && word<1000) ftxt << " " << word;// << endl;
1302 else ftxt << " " << word;// << endl;
1303 ftxt << " " << dummyPS << endl;
1306 AliBitPacking::PackWord(word,baseWord, 0,20);
1307 word = iCH; // TDC channel ID [0;7]
1308 AliBitPacking::PackWord(word,baseWord,21,23);
1309 word = nTDC; // TDC ID [0;14]
1310 AliBitPacking::PackWord(word,baseWord,24,27);
1311 word = 0; // error flag
1312 AliBitPacking::PackWord(word,baseWord,28,28);
1313 word = dummyPS; // Packing Status [0;3]
1314 AliBitPacking::PackWord(word,baseWord,29,30);
1315 word = 1; // TRM TDC digit ID
1316 AliBitPacking::PackWord(word,baseWord,31,31);
1319 localBuffer[localIndex]=baseWord;
1320 psArray[localIndex]=dummyPS;
1325 } // if (percentFilledCells<0.24 && gRandom->Rndm()<(0.24-percentFilledCells))
1328 } // if (!fPackedAcquisition)
1330 } //end loop on digits in the same volume
1332 } // end loop on TDC channel number
1334 } // end loop on TDC number
1336 if (fPackedAcquisition) {
1338 for (Int_t jj=0; jj<=localIndex; jj++) {
1340 buf[fIndex] = localBuffer[jj];
1346 for (Int_t jj=0; jj<=localIndex; jj++) {
1347 if (psArray[jj]==2) {
1349 buf[fIndex] = localBuffer[jj];
1352 for (Int_t jj=0; jj<=localIndex; jj++) {
1353 if (psArray[jj]==1) {
1355 buf[fIndex] = localBuffer[jj];
1361 if (fVerbose==2) ftxt.close();
1365 //----------------------------------------------------------------------------
1366 void AliTOFDDLRawData::ReverseArray(UInt_t a[], Int_t n) const
1369 // Reverses the n elements of array a
1374 for (ii=0; ii<n/2; ii++) {
1384 //----------------------------------------------------------------------------
1385 Bool_t AliTOFDDLRawData::HeadOrTail() const
1388 // Returns the result of a 'pitch and toss'
1391 Double_t dummy = gRandom->Rndm();
1393 if (dummy<0.5) return kFALSE;