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>
22 #include "AliRawDataHeader.h"
23 #include "AliBitPacking.h"
24 #include "AliPMDdigit.h"
25 #include "AliPMDRawStream.h"
26 #include "AliPMDDDLRawData.h"
29 ClassImp(AliPMDDDLRawData)
31 AliPMDDDLRawData::AliPMDDDLRawData():
32 fDigits(new TClonesArray("AliPMDdigit", 1000))
34 // Default Constructor
38 //____________________________________________________________________________
40 AliPMDDDLRawData::~AliPMDDDLRawData()
47 //____________________________________________________________________________
48 void AliPMDDDLRawData::WritePMDRawData(TTree *treeD)
50 // write digits into raw data format
54 TBranch *branch = treeD->GetBranch("PMDDigit");
56 branch->SetAddress(&fDigits);
58 // Int_t nmodules = (Int_t) treeD->GetEntries();
59 // cout << " nmodules = " << nmodules << endl;
61 const Int_t kSize = 4608;
63 Int_t modulePerDDL = 0;
66 AliRawDataHeader header;
67 UInt_t sizeRawData = 0;
75 for(Int_t iddl = 0; iddl < kDDL; iddl++)
77 sprintf(filename,"PMD_%d.ddl",AliPMDRawStream::kDDLOffset+iddl);
79 outfile.open(filename,ios::binary);
81 outfile.open(filename);
101 // Write the Dummy Data Header into the file
102 Int_t bHPosition = outfile.tellp();
103 outfile.write((char*)(&header),sizeof(header));
106 for(Int_t ium = 0; ium < modulePerDDL; ium++)
109 if (iddl == 4 && ium == 6) mmodule = 36;
110 if (iddl == 5 && ium == 6) mmodule = 42;
112 for (Int_t i = 0; i < kSize; i++)
116 // Extract energy deposition per cell and pack it
117 // in a 32-bit word and returns all the total words
118 // per one unit-module
120 GetUMDigitsData(treeD, mmodule, ium, iddl, totword, buffer);
122 outfile.write((char*)buffer,totword*sizeof(UInt_t));
128 // Write real data header
129 // take the pointer to the beginning of the data header
130 // write the total number of words per ddl and bring the
131 // pointer to the current file position and close it
132 UInt_t cFPosition = outfile.tellp();
133 sizeRawData = cFPosition - bHPosition - sizeof(header);
134 header.fSize = cFPosition - bHPosition;
135 header.SetAttribute(0); // valid data
136 outfile.seekp(bHPosition);
137 outfile.write((char*)(&header),sizeof(header));
138 outfile.seekp(cFPosition);
145 //____________________________________________________________________________
146 void AliPMDDDLRawData::GetUMDigitsData(TTree *treeD, Int_t imodule, Int_t ium,
147 Int_t ddlno, Int_t & totword,
159 Int_t det, smn, irow, icol;
162 treeD->GetEntry(imodule);
163 Int_t nentries = fDigits->GetLast();
164 totword = nentries+1;
166 for (Int_t ient = 0; ient < totword; ient++)
168 fPMDdigit = (AliPMDdigit*)fDigits->UncheckedAt(ient);
170 det = fPMDdigit->GetDetector();
171 smn = fPMDdigit->GetSMNumber();
172 irow = fPMDdigit->GetRow();
173 icol = fPMDdigit->GetColumn();
174 adc = (UInt_t) fPMDdigit->GetADC();
176 // cout << " imodule = " << imodule << " smn = " << smn << endl;
183 else if( smn >= 12 && smn < 24)
191 irownew = 95 - irownew1;
194 else if(smn >= 6 && smn < 12)
197 icolnew = 47 - icolnew1;
199 else if(smn >= 12 && smn < 18)
201 irownew = 47 - irownew1;
204 else if(smn >= 18 && smn < 24)
207 icolnew = 95 - icolnew1;
212 GetMCMCh(ddlno, ium, irownew, icolnew, mcmno, chno);
215 AliBitPacking::PackWord(adc,baseword,0,11);
216 AliBitPacking::PackWord(chno,baseword,12,17);
217 AliBitPacking::PackWord(mcmno,baseword,18,28);
218 AliBitPacking::PackWord(0,baseword,29,30);
219 AliBitPacking::PackWord(1,baseword,31,31);
221 buffer[ient] = baseword;
227 //____________________________________________________________________________
228 void AliPMDDDLRawData::GetMCMCh(Int_t ddlno, Int_t um,
229 Int_t row, Int_t col,
230 UInt_t &mcmno, UInt_t &chno)
232 // This part will be modified once the final track layout on the PCB is
233 // designed. This will only change the coordinate of the individual cell
235 static const UInt_t kCh[16][4] = { {12, 13, 18, 19},
252 if (ddlno == 0 || ddlno == 1)
254 // PRE plane, SU Mod = 1, 2
255 Int_t irownew = row%16;
256 Int_t icolnew = col%4;
257 Int_t irowdiv = row/16;
258 Int_t icoldiv = col/4;
260 mcmno = 72*um + 12*irowdiv + icoldiv;
261 chno = kCh[irownew][icolnew];
263 else if (ddlno == 2 || ddlno == 3)
265 // PRE plane, SU Mod = 3, 4
266 Int_t irownew = row%16;
267 Int_t icolnew = col%4;
268 Int_t irowdiv = row/16;
269 Int_t icoldiv = col/4;
271 mcmno = 72*um + 24*irowdiv + icoldiv;
272 chno = kCh[irownew][icolnew];
274 else if (ddlno == 4 || ddlno == 5)
276 // CPV plane, SU Mod = 1, 3 : ddl = 4
277 // CPV plane, SU Mod = 2, 4 : ddl = 5
279 Int_t irownew = row%16;
280 Int_t icolnew = col%4;
281 Int_t irowdiv = row/16;
282 Int_t icoldiv = col/4;
286 mcmno = 72*um + 12*irowdiv + icoldiv;
290 mcmno = 72*um + 24*irowdiv + icoldiv;
292 chno = kCh[irownew][icolnew];
296 //____________________________________________________________________________