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 **************************************************************************/
18 #include <TObjArray.h>
24 #include <TDirectory.h>
26 #include "TTreeStream.h"
31 #include "AliRawReader.h"
32 #include "AliPMDRawStream.h"
33 #include "AliPMDddldata.h"
34 #include "AliBitPacking.h"
36 #include "AliPMDCalibPedestal.h"
39 ClassImp(AliPMDCalibPedestal)
42 AliPMDCalibPedestal::AliPMDCalibPedestal() :
48 // default constructor
51 for (int i = 0; i < kDet; i++)
53 for (int j = 0; j < kMaxSMN; j++)
55 for (int k = 0; k < kMaxRow; k++)
57 for (int l = 0; l < kMaxCol; l++)
59 fPedVal[i][j][k][l] = 0.;
60 fPedValSq[i][j][k][l] = 0.;
61 fPedCount[i][j][k][l] = 0.;
67 for (int i = 0; i < 6; i++)
69 for (int j = 0; j < 51; j++)
71 for (int k = 0; k < 25; k++)
73 for (int l = 0; l < 64; l++)
75 fPedChain[i][j][k][l] = -1;
83 //_____________________________________________________________________
84 AliPMDCalibPedestal::AliPMDCalibPedestal(const AliPMDCalibPedestal &ped) :
86 fRunNumber(ped.fRunNumber),
87 fEventNumber(ped.fEventNumber)
92 for (int i = 0; i < kDet; i++)
94 for (int j = 0; j < kMaxSMN; j++)
96 for (int k = 0; k < kMaxRow; k++)
98 for (int l = 0; l < kMaxCol; l++)
100 fPedVal[i][j][k][l] = ped.fPedVal[i][j][k][l];
101 fPedValSq[i][j][k][l] = ped.fPedValSq[i][j][k][l];
102 fPedCount[i][j][k][l] = ped.fPedCount[i][j][k][l];
108 for (int i = 0; i < 6; i++)
110 for (int j = 0; j < 51; j++)
112 for (int k = 0; k < 25; k++)
114 for (int l = 0; l < 64; l++)
116 fPedChain[i][j][k][l] = -1;
123 //_____________________________________________________________________
124 AliPMDCalibPedestal& AliPMDCalibPedestal::operator = (const AliPMDCalibPedestal &source)
127 // assignment operator
129 if (&source == this) return *this;
130 new (this) AliPMDCalibPedestal(source);
134 //_____________________________________________________________________
135 AliPMDCalibPedestal::~AliPMDCalibPedestal()
141 //_____________________________________________________________________
142 Bool_t AliPMDCalibPedestal::ProcessEvent(AliRawReader *rawReader, TObjArray *pmdddlcont)
145 // Event processing loop - AliRawReader
148 const Int_t kDDL = AliDAQ::NumberOfDdls("PMD");
151 UInt_t pbus, mcm, chno;
153 fRunNumber = rawReader->GetRunNumber();
155 AliPMDRawStream rawStream(rawReader);
160 Int_t numberofDDLs = 0;
162 while ((iddl = rawStream.DdlData(pmdddlcont)) >=0) {
164 Int_t ientries = pmdddlcont->GetEntries();
165 //printf("iddl = %d ientries = %d\n",iddl, ientries);
166 for (Int_t ient = 0; ient < ientries; ient++)
168 AliPMDddldata *pmdddl = (AliPMDddldata*)pmdddlcont->UncheckedAt(ient);
170 Int_t det = pmdddl->GetDetector();
171 Int_t smn = pmdddl->GetSMN();
172 Int_t row = pmdddl->GetRow();
173 Int_t col = pmdddl->GetColumn();
174 Float_t sig = (Float_t) pmdddl->GetSignal();
176 pbus = (UInt_t) pmdddl->GetPatchBusId();
177 mcm = (UInt_t) pmdddl->GetMCM();
178 chno = (UInt_t) pmdddl->GetChannel();
181 AliBitPacking::PackWord(det,detsmnrowcol,0,7);
182 AliBitPacking::PackWord(smn,detsmnrowcol,8,15);
183 AliBitPacking::PackWord(row,detsmnrowcol,16,23);
184 AliBitPacking::PackWord(col,detsmnrowcol,24,31);
186 if (fPedChain[iddl][pbus][mcm][chno] == -1)
187 fPedChain[iddl][pbus][mcm][chno] = (Int_t)detsmnrowcol;
190 fPedVal[det][smn][row][col] += sig;
191 fPedValSq[det][smn][row][col] += sig*sig;
192 fPedCount[det][smn][row][col]++;
194 pmdddlcont->Delete();
197 if (numberofDDLs < kDDL)
201 //_____________________________________________________________________
203 void AliPMDCalibPedestal::Analyse(TTree *pedtree)
206 // Calculate pedestal Mean and RMS
210 Int_t det, sm, row, col;
211 Int_t idet, ism, irow, icol;
213 Float_t meansq, diff;
215 FILE *fpw0 = fopen("pedestal2304.ped","w");
216 FILE *fpw1 = fopen("pedestal2305.ped","w");
217 FILE *fpw2 = fopen("pedestal2306.ped","w");
218 FILE *fpw3 = fopen("pedestal2307.ped","w");
219 FILE *fpw4 = fopen("pedestal2308.ped","w");
220 FILE *fpw5 = fopen("pedestal2309.ped","w");
222 fprintf(fpw0,"//=============================================\n");
223 fprintf(fpw0,"// Pedestal file Calculated by Online DA\n");
224 fprintf(fpw0,"//=============================================\n");
225 fprintf(fpw0,"// RUN :%d\n",fRunNumber);
226 fprintf(fpw0,"// Statistics :%d\n",fEventNumber);
227 fprintf(fpw0,"//---------------------------------------------\n");
228 fprintf(fpw0,"//format:CHAIN_NO FEE_ID CHANNEL MEAN SIGMA\n");
229 fprintf(fpw0,"//---------------------------------------------\n");
231 fprintf(fpw1,"//=============================================\n");
232 fprintf(fpw1,"// Pedestal file Calculated by Online DA\n");
233 fprintf(fpw1,"//=============================================\n");
234 fprintf(fpw1,"// RUN :%d\n",fRunNumber);
235 fprintf(fpw1,"// Statistics :%d\n",fEventNumber);
236 fprintf(fpw1,"//---------------------------------------------\n");
237 fprintf(fpw1,"//format:CHAIN_NO FEE_ID CHANNEL MEAN SIGMA\n");
239 fprintf(fpw2,"//=============================================\n");
240 fprintf(fpw2,"// Pedestal file Calculated by Online DA\n");
241 fprintf(fpw2,"//=============================================\n");
242 fprintf(fpw2,"// RUN :%d\n",fRunNumber);
243 fprintf(fpw2,"// Statistics :%d\n",fEventNumber);
244 fprintf(fpw2,"//---------------------------------------------\n");
245 fprintf(fpw2,"//format:CHAIN_NO FEE_ID CHANNEL MEAN SIGMA\n");
246 fprintf(fpw2,"//---------------------------------------------\n");
248 fprintf(fpw3,"//=============================================\n");
249 fprintf(fpw3,"// Pedestal file Calculated by Online DA\n");
250 fprintf(fpw3,"//=============================================\n");
251 fprintf(fpw3,"// RUN :%d\n",fRunNumber);
252 fprintf(fpw3,"// Statistics :%d\n",fEventNumber);
253 fprintf(fpw3,"//---------------------------------------------\n");
254 fprintf(fpw3,"//format:CHAIN_NO FEE_ID CHANNEL MEAN SIGMA\n");
255 fprintf(fpw3,"//---------------------------------------------\n");
257 fprintf(fpw4,"//=============================================\n");
258 fprintf(fpw4,"// Pedestal file Calculated by Online DA\n");
259 fprintf(fpw4,"//=============================================\n");
260 fprintf(fpw4,"// RUN :%d\n",fRunNumber);
261 fprintf(fpw4,"// Statistics :%d\n",fEventNumber);
262 fprintf(fpw4,"//---------------------------------------------\n");
263 fprintf(fpw4,"//format:CHAIN_NO FEE_ID CHANNEL MEAN SIGMA\n");
264 fprintf(fpw4,"//---------------------------------------------\n");
266 fprintf(fpw5,"//=============================================\n");
267 fprintf(fpw5,"// Pedestal file Calculated by Online DA\n");
268 fprintf(fpw5,"//=============================================\n");
269 fprintf(fpw5,"// RUN :%d\n",fRunNumber);
270 fprintf(fpw5,"// Statistics :%d\n",fEventNumber);
271 fprintf(fpw5,"//---------------------------------------------\n");
272 fprintf(fpw5,"//format:CHAIN_NO FEE_ID CHANNEL MEAN SIGMA\n");
273 fprintf(fpw5,"//---------------------------------------------\n");
276 for(Int_t iddl = 0; iddl < 6; iddl++)
278 for(Int_t ibus = 1; ibus < 51; ibus++)
280 for(Int_t imcm = 1; imcm < 25; imcm++)
282 for(Int_t ich = 0; ich < 64; ich++)
285 if (fPedChain[iddl][ibus][imcm][ich] != -1)
287 detsmnrowcol = (UInt_t)fPedChain[iddl][ibus][imcm][ich];
289 idet = detsmnrowcol & 0x00FF;
290 ism = (detsmnrowcol >> 8) & 0x00FF;
291 irow = (detsmnrowcol >> 16) & 0x00FF;
292 icol = (detsmnrowcol >> 24) & 0x00FF;
296 if (fPedCount[idet][ism][irow][icol] > 0)
298 mean = fPedVal[idet][ism][irow][icol]/fPedCount[idet][ism][irow][icol];
300 meansq = fPedValSq[idet][ism][irow][icol]/fPedCount[idet][ism][irow][icol];
302 diff = meansq - mean*mean;
315 fprintf(fpw0,"%d %d %d %f %f\n",
316 ibus, imcm, ich, mean, rms);
320 fprintf(fpw1,"%d %d %d %f %f\n",
321 ibus, imcm, ich, mean, rms);
325 fprintf(fpw2,"%d %d %d %f %f\n",
326 ibus, imcm, ich, mean, rms);
330 fprintf(fpw3,"%d %d %d %f %f\n",
331 ibus, imcm, ich, mean, rms);
335 fprintf(fpw4,"%d %d %d %f %f\n",
336 ibus, imcm, ich, mean, rms);
340 fprintf(fpw5,"%d %d %d %f %f\n",
341 ibus, imcm, ich, mean, rms);
358 pedtree->Branch("det",&det,"det/I");
359 pedtree->Branch("sm",&sm,"sm/I");
360 pedtree->Branch("row",&row,"row/I");
361 pedtree->Branch("col",&col,"col/I");
362 pedtree->Branch("mean",&mean,"mean/F");
363 pedtree->Branch("rms",&rms,"rms/F");
365 for (idet = 0; idet < kDet; idet++)
367 for (ism = 0; ism < kMaxSMN; ism++)
369 for (irow = 0; irow < kMaxRow; irow++)
371 for (icol = 0; icol < kMaxCol; icol++)
380 if (fPedCount[idet][ism][irow][icol] > 0)
382 mean = fPedVal[idet][ism][irow][icol]/fPedCount[idet][ism][irow][icol];
384 meansq = fPedValSq[idet][ism][irow][icol]/fPedCount[idet][ism][irow][icol];
386 diff = meansq - mean*mean;
403 // -------------------------------------------------------------------