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 #include <Riostream.h>
17 #include <TClonesArray.h>
25 #include "AliRawDataHeaderSim.h"
26 #include "AliBitPacking.h"
27 #include "AliPMDdigit.h"
28 #include "AliPMDBlockHeader.h"
29 #include "AliPMDDspHeader.h"
30 #include "AliPMDPatchBusHeader.h"
31 #include "AliPMDRawStream.h"
32 #include "AliPMDddlinfoData.h"
33 #include "AliPMDMappingData.h"
34 #include "AliPMDDDLRawData.h"
36 #include "AliFstream.h"
37 #include "AliCDBManager.h"
38 #include "AliCDBStorage.h"
39 #include "AliCDBEntry.h"
41 ClassImp(AliPMDDDLRawData)
43 AliPMDDDLRawData::AliPMDDDLRawData():
44 fDdlinfo(GetDdlinfoData()),
45 fMapData(GetMappingData()),
46 fDigits(new TClonesArray("AliPMDdigit", 1000))
48 // Default Constructor
52 //____________________________________________________________________________
53 AliPMDDDLRawData::AliPMDDDLRawData(const AliPMDDDLRawData& ddlraw):
55 fDdlinfo(ddlraw.fDdlinfo),
56 fMapData(ddlraw.fMapData),
57 fDigits(ddlraw.fDigits)
61 //____________________________________________________________________________
62 AliPMDDDLRawData & AliPMDDDLRawData::operator=(const AliPMDDDLRawData& ddlraw)
67 fDdlinfo = ddlraw.fDdlinfo;
68 fMapData = ddlraw.fMapData;
69 fDigits = ddlraw.fDigits;
73 //____________________________________________________________________________
75 AliPMDDDLRawData::~AliPMDDDLRawData()
82 //____________________________________________________________________________
83 void AliPMDDDLRawData::WritePMDRawData(TTree *treeD)
85 // write digits into raw data format
89 TBranch *branch = treeD->GetBranch("PMDDigit");
92 AliError("PMD Digit branch not found");
95 branch->SetAddress(&fDigits);
97 Int_t nmodules = (Int_t) treeD->GetEntries();
98 AliDebug(1,Form("Number of modules inside treeD = %d",nmodules));
100 const Int_t kDDL = AliDAQ::NumberOfDdls("PMD");
103 AliRawDataHeaderSim header;
104 UInt_t sizeRawData = 0;
106 const Int_t kbusSize = 51;
107 const Int_t kSize = 1536;
108 UInt_t buffer[kSize];
110 UInt_t busPatch[kbusSize][1536];
112 Int_t contentsBus[kbusSize];
116 Int_t modulePerDDL = 0;
119 Int_t modulenoddl[12] = {-1,-1,-1,-1,-1,-1,-1,-1,-1,-1,-1,-1};
121 for(Int_t iddl = 0; iddl < kDDL; iddl++)
124 modulePerDDL = fDdlinfo->GetNoOfModulePerDdl(iddl);
125 if (modulePerDDL == 0) continue;
126 for (Int_t im = 0; im < 12; im++)
128 modulenoddl[im] = fDdlinfo->GetModulesPerDdl(iddl,im);;
131 strncpy(filename,AliDAQ::DdlFileName("PMD",iddl),79);
133 outfile = new AliFstream(filename);
135 // Write the Dummy Data Header into the file
136 Int_t bHPosition = outfile->Tellp();
137 outfile->WriteBuffer((char*)(&header),sizeof(header));
139 for (Int_t ibus = 0; ibus < kbusSize; ibus++)
141 contentsBus[ibus] = 0;
142 for (Int_t ich = 0; ich < kSize; ich++)
144 busPatch[ibus][ich] = 0;
148 for(Int_t ium = 0; ium < 12; ium++)
150 // Extract energy deposition per cell and pack it
151 // in a 32-bit word and returns all the total words
152 // per one unit-module
154 mmodule = modulenoddl[ium];
155 if(mmodule == -1) continue;
156 GetUMDigitsData(treeD, mmodule, iddl, contentsBus, busPatch);
162 for (Int_t i = 0; i < 10; i++)
166 for (Int_t ibus=0; ibus < 5; ibus++)
169 if (contentsBus[ij] > 0)
171 dsp[i] += contentsBus[ij];
175 // Add the patch Bus header to the DSP contents
176 dsp[i] += 4*dspBus[i];
179 Int_t dspBlockARDL = 0;
180 Int_t dspBlockBRDL = 0;
184 for (Int_t i = 0; i < 5; i++)
187 Int_t iodd = 2*i + 1;
190 dspBlockARDL += dsp[ieven];
191 remainder = dsp[ieven]%2;
199 dspBlockBRDL += dsp[iodd];
200 remainder = dsp[iodd]%2;
211 // Start writing the DDL file
213 AliPMDBlockHeader blHeader;
214 AliPMDDspHeader dspHeader;
215 AliPMDPatchBusHeader pbusHeader;
217 const Int_t kblHLen = blHeader.GetHeaderLength();
218 const Int_t kdspHLen = dspHeader.GetHeaderLength();
219 const Int_t kpbusHLen = pbusHeader.GetHeaderLength();
222 UInt_t dspBlockHeaderWord[8];
223 UInt_t dspHeaderWord[10];
224 UInt_t patchBusHeaderWord[4];
226 UInt_t ddlEndWord[2] = {0xDEADFACE, 0xDEADFACE};
228 for (Int_t iblock = 0; iblock < 2; iblock++)
232 for (Int_t i=0; i<kblHLen; i++)
234 dspBlockHeaderWord[i] = 0;
238 dspBlockHeaderWord[1] = (UInt_t) (dspBlockARDL + kblHLen);
239 dspBlockHeaderWord[2] = (UInt_t) dspBlockARDL;
241 else if (iblock == 1)
243 dspBlockHeaderWord[1] = (UInt_t) (dspBlockBRDL + kblHLen);
244 dspBlockHeaderWord[2] = (UInt_t) dspBlockBRDL;
247 outfile->WriteBuffer((char*)dspBlockHeaderWord,kblHLen*sizeof(UInt_t));
257 else if (iblock == 1)
266 for (Int_t idsp = 0; idsp < 5; idsp++)
273 dspRDL = (UInt_t) dsp[dspno];
275 else if (iblock == 1)
278 dspRDL = (UInt_t) dsp[dspno];
281 for (Int_t i=0; i<kdspHLen; i++)
283 dspHeaderWord[i] = 0;
285 remainder = dspRDL%2;
286 if (remainder == 1) dspRDL++;
288 dspHeaderWord[1] = dspRDL + kdspHLen;
289 dspHeaderWord[2] = dspRDL;
290 dspHeaderWord[3] = dspno;
291 if (remainder == 1) dspHeaderWord[8] = 1; // setting the padding word
294 outfile->WriteBuffer((char*)dspHeaderWord,kdspHLen*sizeof(UInt_t));
296 for (Int_t ibus = 0; ibus < 5; ibus++)
300 Int_t busno = iskip[idsp] + ibus + 1;
301 Int_t patchbusRDL = contentsBus[busno];
305 patchBusHeaderWord[0] = 0;
306 patchBusHeaderWord[1] = (UInt_t) (patchbusRDL + kpbusHLen);
307 patchBusHeaderWord[2] = (UInt_t) patchbusRDL;
308 patchBusHeaderWord[3] = (UInt_t) busno;
310 else if (patchbusRDL == 0)
312 patchBusHeaderWord[0] = 0;
313 patchBusHeaderWord[1] = (UInt_t) kpbusHLen;
314 patchBusHeaderWord[2] = (UInt_t) 0;
315 patchBusHeaderWord[3] = (UInt_t) busno;
319 outfile->WriteBuffer((char*)patchBusHeaderWord,4*sizeof(UInt_t));
321 for (Int_t iword = 0; iword < patchbusRDL; iword++)
323 buffer[iword] = busPatch[busno][iword];
326 outfile->WriteBuffer((char*)buffer,patchbusRDL*sizeof(UInt_t));
328 } // End of patch bus loop
331 // Adding a padding word if the total words odd
334 UInt_t paddingWord = dspHeader.GetDefaultPaddingWord();
335 outfile->WriteBuffer((char*)(&paddingWord),sizeof(UInt_t));
340 // Write two extra word at the end of each DDL file
341 outfile->WriteBuffer((char*)ddlEndWord,2*sizeof(UInt_t));
343 // Write real data header
344 // take the pointer to the beginning of the data header
345 // write the total number of words per ddl and bring the
346 // pointer to the current file position and close it
347 UInt_t cFPosition = outfile->Tellp();
348 sizeRawData = cFPosition - bHPosition - sizeof(header);
350 header.fSize = cFPosition - bHPosition;
351 header.SetAttribute(0); // valid data
352 outfile->Seekp(bHPosition);
353 outfile->WriteBuffer((char*)(&header),sizeof(header));
354 outfile->Seekp(cFPosition);
360 //____________________________________________________________________________
362 void AliPMDDDLRawData::GetUMDigitsData(TTree *treeD, Int_t imodule,
363 Int_t ddlno, Int_t *contentsBus,
364 UInt_t busPatch[][1536])
366 // Retrieves digits data UnitModule by UnitModule
368 const Int_t kMaxBus = 51;
371 UInt_t mcmno = 0, chno = 0;
373 Int_t det = 0, smn = 0, irow = 0, icol = 0;
377 Int_t patchBusNo[kMaxBus], mcmperBus[kMaxBus];
378 Int_t startRowBus[kMaxBus], startColBus[kMaxBus];
379 Int_t endRowBus[kMaxBus], endColBus[kMaxBus];
381 Int_t beginPatchBus = -1;
382 Int_t endPatchBus = -1;
383 for(Int_t i = 0; i < kMaxBus; i++)
393 // Fetch the DDL mapping info from the mapping database
395 DdlMapping(ddlno, imodule, beginPatchBus, endPatchBus,
396 patchBusNo, mcmperBus, startRowBus, endRowBus,
397 startColBus, endColBus);
399 // Read if some chains are off from the ddlinfo database
401 Int_t srowoff1[2][24], erowoff1[2][24];
402 Int_t scoloff1[2][24], ecoloff1[2][24];
403 Int_t srowoff2[2][24], erowoff2[2][24];
404 Int_t scoloff2[2][24], ecoloff2[2][24];
406 for (Int_t idet = 0; idet < 2; idet++)
408 for (Int_t im = 0; im < 24; im++)
410 srowoff1[idet][im] = fDdlinfo->GetStartRowA(idet,im);
411 erowoff1[idet][im] = fDdlinfo->GetEndRowA(idet,im);
412 scoloff1[idet][im] = fDdlinfo->GetStartColA(idet,im);
413 ecoloff1[idet][im] = fDdlinfo->GetEndColA(idet,im);
414 srowoff2[idet][im] = fDdlinfo->GetStartRowB(idet,im);
415 erowoff2[idet][im] = fDdlinfo->GetEndRowB(idet,im);
416 scoloff2[idet][im] = fDdlinfo->GetStartColB(idet,im);
417 ecoloff2[idet][im] = fDdlinfo->GetEndColB(idet,im);
421 treeD->GetEntry(imodule);
422 Int_t nentries = fDigits->GetLast();
423 Int_t totword = nentries+1;
425 AliPMDdigit *pmddigit = 0x0;
427 for (Int_t ient = 0; ient < totword; ient++)
429 pmddigit = (AliPMDdigit*)fDigits->UncheckedAt(ient);
431 det = pmddigit->GetDetector();
432 smn = pmddigit->GetSMNumber();
433 irow = pmddigit->GetRow();
434 icol = pmddigit->GetColumn();
435 Float_t aadc = pmddigit->GetADC();
436 if (aadc < 0.) aadc = 0.;
439 TransformS2H(smn,irow,icol);
441 // remove the non-existence channels
443 //printf("%d %d %d %d\n",det,smn,irow,icol);
444 //printf("--- %d %d %d %d\n",srowoff[det][smn],erowoff[det][smn],
445 // scoloff[det][smn],ecoloff[det][smn]);
447 if (irow >= srowoff1[det][smn] && irow <= erowoff1[det][smn])
449 if (icol >= scoloff1[det][smn] && icol <= ecoloff1[det][smn])
454 if (irow >= srowoff2[det][smn] && irow <= erowoff2[det][smn])
456 if (icol >= scoloff2[det][smn] && icol <= ecoloff2[det][smn])
463 GetMCMCh(imodule, irow, icol, beginPatchBus, endPatchBus,
464 mcmperBus, startRowBus, startColBus,
465 endRowBus, endColBus, busno, mcmno, chno);
468 AliBitPacking::PackWord(adc,baseword,0,11);
469 AliBitPacking::PackWord(chno,baseword,12,17);
470 AliBitPacking::PackWord(mcmno,baseword,18,28);
471 AliBitPacking::PackWord(0,baseword,29,30);
472 parity = ComputeParity(baseword); // generate the parity bit
473 AliBitPacking::PackWord(parity,baseword,31,31);
475 Int_t jj = contentsBus[busno];
476 busPatch[busno][jj] = baseword;
478 contentsBus[busno]++;
483 //____________________________________________________________________________
484 void AliPMDDDLRawData::TransformS2H(Int_t smn, Int_t &irow, Int_t &icol)
486 // Does the Software to Hardware coordinate transformation
492 // First in digits we have all dimension 48x96
493 // Transform into the realistic one, i.e, For SM 0&1 96(row)x48(col)
494 // and for SM 2&3 48(row)x96(col)
501 else if( smn >= 12 && smn < 24)
513 //____________________________________________________________________________
515 void AliPMDDDLRawData::GetMCMCh(Int_t imodule, Int_t row, Int_t col,
516 Int_t beginPatchBus, Int_t endPatchBus,
518 Int_t *startRowBus, Int_t *startColBus,
519 Int_t *endRowBus, Int_t *endColBus,
520 Int_t & busno, UInt_t &mcmno, UInt_t &chno)
522 // This converts row col to hardware channel number
523 // This is the final version of mapping supplied by Mriganka
527 static const UInt_t kChDdl01[16][4] = { {6, 4, 5, 7},
545 static const UInt_t kChDdl23[16][4] = { {57, 59, 58, 56},
563 static const UInt_t kChDdl41[16][4] = { {56, 58, 59, 57},
581 static const UInt_t kChDdl42[16][4] = { {7, 5, 4, 6},
599 static const UInt_t kChDdl51[16][4] = { {7, 5, 4, 6},
618 static const UInt_t kChDdl52[16][4] = { {56, 58, 59, 57},
636 for (Int_t i = 0; i < 16; i++)
638 for (Int_t j = 0; j < 4; j++)
641 if(imodule < 6) iCh[i][j] = kChDdl01[i][j];
642 if(imodule >= 6 && imodule <= 11) iCh[i][j] = kChDdl01[i][j];
643 if(imodule >= 12 && imodule <= 17) iCh[i][j] = kChDdl23[i][j];
644 if(imodule >= 18 && imodule <= 23) iCh[i][j] = kChDdl23[i][j];
646 if(imodule >= 24 && imodule <= 29) iCh[i][j] = kChDdl41[i][j];
647 if(imodule >= 42 && imodule <= 47) iCh[i][j] = kChDdl42[i][j];
648 if(imodule >= 36 && imodule <= 41) iCh[i][j] = kChDdl51[i][j];
649 if(imodule >= 30 && imodule <= 35) iCh[i][j] = kChDdl52[i][j];
655 Int_t irownew = row%16;
656 Int_t icolnew = col%4;
658 chno = iCh[irownew][icolnew];
661 for (Int_t ibus = beginPatchBus; ibus <= endPatchBus; ibus++)
663 Int_t srow = startRowBus[ibus];
664 Int_t erow = endRowBus[ibus];
665 Int_t scol = startColBus[ibus];
666 Int_t ecol = endColBus[ibus];
667 Int_t tmcm = mcmperBus[ibus];
669 if ((row >= srow && row <= erow) && (col >= scol && col <= ecol))
673 // Find out the MCM Number
676 if (imodule < 6) mcmno = (col-scol)/4 + 1;
677 if (imodule >= 6 && imodule < 12) mcmno = (col-scol)/4 + 1;
679 if (imodule >= 12 && imodule < 18)
681 icolnew = (col - scol)/4;
682 mcmno = tmcm - icolnew;
684 if (imodule >= 18 && imodule < 24)
686 icolnew = (col - scol)/4;
687 mcmno = tmcm - icolnew;
691 if (imodule >= 24 && imodule < 30)
695 Int_t rowdiff = endRowBus[ibus] - startRowBus[ibus];
698 Int_t midrow = srow + 16;
699 if(row >= srow && row < midrow)
701 mcmno = 12 + (col-scol)/4 + 1;
703 else if(row >= midrow && row <= erow)
706 mcmno = (col-scol)/4 + 1;
709 else if (rowdiff < 16)
711 mcmno = (col-scol)/4 + 1;
715 if (imodule >= 42 && imodule < 48)
717 Int_t rowdiff = endRowBus[ibus] - startRowBus[ibus];
720 Int_t midrow = srow + 16;
721 if (row >= midrow && row <= erow)
723 mcmno = 12 + (ecol -col)/4 + 1;
725 else if (row >= srow && row < midrow)
727 mcmno = (ecol - col)/4 + 1;
730 else if (rowdiff < 16)
732 mcmno = (ecol - col)/4 + 1;
737 if (imodule >= 30 && imodule < 36)
739 // CPV plane, SU Mod = 1, 2 : ddl = 5
742 Int_t rowdiff = endRowBus[ibus] - startRowBus[ibus];
745 Int_t midrow = srow + 16;
746 if(row >= srow && row < midrow)
748 mcmno = 12 + (col-scol)/4 + 1;
750 else if(row >= midrow && row <= erow)
752 mcmno = (col-scol)/4 + 1;
755 else if(rowdiff < 16)
757 mcmno = (col-scol)/4 + 1;
761 if (imodule >= 36 && imodule < 42)
763 Int_t rowdiff = endRowBus[ibus] - startRowBus[ibus];
766 Int_t midrow = srow + 16;
767 if (row >= midrow && row <= erow)
769 mcmno = 12 + (ecol - col)/4 + 1;
771 else if (row >= srow && row < midrow)
773 mcmno = (ecol - col)/4 + 1;
776 else if (rowdiff < 16)
778 mcmno = (ecol - col)/4 + 1;
786 //____________________________________________________________________________
788 UInt_t AliPMDDDLRawData::ComputeParity(UInt_t baseword)
790 // Generate the parity bit
793 for(Int_t j=0; j<29; j++)
795 if (baseword & 0x01 ) count++;
798 UInt_t parity = count%2;
801 //____________________________________________________________________________
802 void AliPMDDDLRawData::DdlMapping(Int_t iddl, Int_t imodule,
803 Int_t &beginPatchBus, Int_t &endPatchBus,
804 Int_t patchBusNo[], Int_t mcmperBus[],
805 Int_t startRowBus[], Int_t endRowBus[],
806 Int_t startColBus[], Int_t endColBus[])
808 // DDL Mapping fetching from mapping database
810 beginPatchBus = fMapData->GetBeginPatchBus(iddl,imodule);
811 endPatchBus = fMapData->GetEndPatchBus(iddl,imodule);
813 for(Int_t ibus = beginPatchBus; ibus < endPatchBus+1; ibus++)
815 patchBusNo[ibus] = ibus;
816 mcmperBus[ibus] = fMapData->GetMcmperBus(iddl,ibus);
817 startRowBus[ibus] = fMapData->GetStartRowBus(iddl,ibus);
818 startColBus[ibus] = fMapData->GetStartColBus(iddl,ibus);
819 endRowBus[ibus] = fMapData->GetEndRowBus(iddl,ibus);
820 endColBus[ibus] = fMapData->GetEndColBus(iddl,ibus);
824 //____________________________________________________________________________
826 AliPMDddlinfoData* AliPMDDDLRawData::GetDdlinfoData() const
828 // The run number will be centralized in AliCDBManager,
829 // you don't need to set it here!
830 AliCDBEntry *entry = AliCDBManager::Instance()->Get("PMD/Calib/Ddlinfo");
832 if(!entry) AliFatal("ddlinfo object retrieval failed!");
834 AliPMDddlinfoData *ddlinfo = 0;
835 if (entry) ddlinfo = (AliPMDddlinfoData*) entry->GetObject();
837 if (!ddlinfo) AliFatal("No ddl info data from database !");
841 //____________________________________________________________________________
842 AliPMDMappingData* AliPMDDDLRawData::GetMappingData() const
844 // The run number will be centralized in AliCDBManager,
845 // you don't need to set it here!
846 AliCDBEntry *entry = AliCDBManager::Instance()->Get("PMD/Calib/Mapping");
848 if(!entry) AliFatal("Mapping object retrieval failed!");
850 AliPMDMappingData *mapda = 0;
851 if (entry) mapda = (AliPMDMappingData*) entry->GetObject();
853 if (!mapda) AliFatal("No mapping data from database !");
857 //____________________________________________________________________________