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 //-----------------------------------------------------//
18 // Date : August 05 2003 //
19 // This reads the file PMD.digits.root(TreeD), //
20 // calls the Clustering algorithm and stores the //
21 // clustering output in PMD.RecPoints.root(TreeR) //
23 //-----------------------------------------------------//
25 #include <Riostream.h>
27 #include <TObjArray.h>
28 #include <TClonesArray.h>
30 #include "AliRunLoader.h"
31 #include "AliLoader.h"
32 #include "AliRawReader.h"
34 #include "AliPMDdigit.h"
35 #include "AliPMDClusterFinder.h"
36 #include "AliPMDClustering.h"
37 #include "AliPMDcluster.h"
38 #include "AliPMDrecpoint1.h"
39 #include "AliPMDRawStream.h"
41 ClassImp(AliPMDClusterFinder)
43 AliPMDClusterFinder::AliPMDClusterFinder(AliRunLoader* runLoader):
44 fRunLoader(runLoader),
45 fPMDLoader(runLoader->GetLoader("PMDLoader")),
48 fDigits(new TClonesArray("AliPMDdigit", 1000)),
49 fRecpoints(new TClonesArray("AliPMDrecpoint1", 1000)),
58 // ------------------------------------------------------------------------- //
59 AliPMDClusterFinder::~AliPMDClusterFinder()
75 // ------------------------------------------------------------------------- //
77 void AliPMDClusterFinder::Digits2RecPoints(Int_t ievt)
79 // Converts digits to recpoints after running clustering
80 // algorithm on CPV plane and PREshower plane
82 Int_t det = 0,smn = 0;
89 TObjArray *pmdcont = new TObjArray();
90 AliPMDClustering *pmdclust = new AliPMDClustering();
91 pmdclust->SetDebug(fDebug);
92 pmdclust->SetEdepCut(fEcut);
94 fRunLoader->GetEvent(ievt);
95 //cout << " ***** Beginning::Digits2RecPoints *****" << endl;
97 fTreeD = fPMDLoader->TreeD();
100 cout << " Can not get TreeD" << endl;
102 AliPMDdigit *pmddigit;
103 TBranch *branch = fTreeD->GetBranch("PMDDigit");
104 branch->SetAddress(&fDigits);
108 fTreeR = fPMDLoader->TreeR();
111 fPMDLoader->MakeTree("R");
112 fTreeR = fPMDLoader->TreeR();
115 Int_t bufsize = 16000;
116 fTreeR->Branch("PMDRecpoint", &fRecpoints, bufsize);
118 Int_t nmodules = (Int_t) fTreeD->GetEntries();
120 for (Int_t imodule = 0; imodule < nmodules; imodule++)
123 fTreeD->GetEntry(imodule);
124 Int_t nentries = fDigits->GetLast();
125 for (Int_t ient = 0; ient < nentries+1; ient++)
127 pmddigit = (AliPMDdigit*)fDigits->UncheckedAt(ient);
129 det = pmddigit->GetDetector();
130 smn = pmddigit->GetSMNumber();
131 xpos = pmddigit->GetRow();
132 ypos = pmddigit->GetColumn();
133 adc = pmddigit->GetADC();
134 //Int_t trno = pmddigit->GetTrackNumber();
135 fCellADC[xpos][ypos] = (Double_t) adc;
140 pmdclust->DoClust(idet,ismn,fCellADC,pmdcont);
142 Int_t nentries1 = pmdcont->GetEntries();
143 // cout << " nentries1 = " << nentries1 << endl;
144 for (Int_t ient1 = 0; ient1 < nentries1; ient1++)
146 AliPMDcluster *pmdcl = (AliPMDcluster*)pmdcont->UncheckedAt(ient1);
147 idet = pmdcl->GetDetector();
148 ismn = pmdcl->GetSMN();
149 clusdata[0] = pmdcl->GetClusX();
150 clusdata[1] = pmdcl->GetClusY();
151 clusdata[2] = pmdcl->GetClusADC();
152 clusdata[3] = pmdcl->GetClusCells();
153 clusdata[4] = pmdcl->GetClusRadius();
155 AddRecPoint(idet,ismn,clusdata);
165 fPMDLoader = fRunLoader->GetLoader("PMDLoader");
166 fPMDLoader->WriteRecPoints("OVERWRITE");
168 // delete the pointers
172 // cout << " ***** End::Digits2RecPoints *****" << endl;
174 // ------------------------------------------------------------------------- //
176 void AliPMDClusterFinder::Digits2RecPoints(Int_t ievt, AliRawReader *rawReader)
178 // Converts RAW data to recpoints after running clustering
179 // algorithm on CPV and PREshower plane
184 TObjArray *pmdcont = new TObjArray();
185 AliPMDClustering *pmdclust = new AliPMDClustering();
186 pmdclust->SetDebug(fDebug);
187 pmdclust->SetEdepCut(fEcut);
189 fRunLoader->GetEvent(ievt);
193 fTreeR = fPMDLoader->TreeR();
196 fPMDLoader->MakeTree("R");
197 fTreeR = fPMDLoader->TreeR();
199 Int_t bufsize = 16000;
200 fTreeR->Branch("PMDRecpoint", &fRecpoints, bufsize);
202 const Int_t kDDL = 6;
203 const Int_t kRow = 48;
204 const Int_t kCol = 96;
209 for (Int_t indexDDL = 0; indexDDL < kDDL; indexDDL++)
215 else if (indexDDL >= 4)
220 precpvADC = new int **[iSMN];
221 for (Int_t i=0; i<iSMN; i++) precpvADC[i] = new int *[kRow];
222 for (Int_t i=0; i<iSMN;i++)
224 for (Int_t j=0; j<kRow; j++) precpvADC[i][j] = new int [kCol];
226 for (Int_t i = 0; i < iSMN; i++)
228 for (Int_t j = 0; j < kRow; j++)
230 for (Int_t k = 0; k < kCol; k++)
232 precpvADC[i][j][k] = 0;
238 AliPMDRawStream pmdinput(rawReader);
239 rawReader->Select(12, indexDDL, indexDDL);
240 while(pmdinput.Next())
242 Int_t det = pmdinput.GetDetector();
243 Int_t smn = pmdinput.GetSMN();
244 //Int_t mcm = pmdinput.GetMCM();
245 //Int_t chno = pmdinput.GetChannel();
246 Int_t row = pmdinput.GetRow();
247 Int_t col = pmdinput.GetColumn();
248 Int_t sig = pmdinput.GetSignal();
255 printf(" *** DDL %d and Detector NUMBER %d NOT MATCHING *** ",
257 indexsmn = smn - indexDDL * 6;
259 else if (indexDDL == 4)
262 printf(" *** DDL %d and Detector NUMBER %d NOT MATCHING *** ",
268 else if (smn >= 12 && smn < 18)
273 else if (indexDDL == 5)
276 printf(" *** DDL %d and Detector NUMBER %d NOT MATCHING *** ",
278 if (smn >= 6 && smn < 12)
282 else if (smn >= 18 && smn < 24)
287 precpvADC[indexsmn][row][col] = sig;
291 for (Int_t indexsmn = 0; indexsmn < iSMN; indexsmn++)
294 for (Int_t irow = 0; irow < kRow; irow++)
296 for (Int_t icol = 0; icol < kCol; icol++)
298 fCellADC[irow][icol] =
299 (Double_t) precpvADC[indexsmn][irow][icol];
304 ismn = indexsmn + indexDDL * 6;
307 else if (indexDDL == 4)
313 else if (indexsmn >= 6 && indexsmn < 12)
319 else if (indexDDL == 5)
325 else if (indexsmn >= 6 && indexsmn < 12)
327 ismn = indexsmn + 12;
333 pmdclust->DoClust(idet,ismn,fCellADC,pmdcont);
334 Int_t nentries1 = pmdcont->GetEntries();
335 // cout << " nentries1 = " << nentries1 << endl;
336 for (Int_t ient1 = 0; ient1 < nentries1; ient1++)
338 AliPMDcluster *pmdcl =
339 (AliPMDcluster*)pmdcont->UncheckedAt(ient1);
340 idet = pmdcl->GetDetector();
341 ismn = pmdcl->GetSMN();
342 clusdata[0] = pmdcl->GetClusX();
343 clusdata[1] = pmdcl->GetClusY();
344 clusdata[2] = pmdcl->GetClusADC();
345 clusdata[3] = pmdcl->GetClusCells();
346 clusdata[4] = pmdcl->GetClusRadius();
348 AddRecPoint(idet,ismn,clusdata);
358 for (Int_t i=0; i<iSMN; i++)
360 for (Int_t j=0; j<kRow; j++) delete [] precpvADC[i][j];
362 for (Int_t i=0; i<iSMN; i++) delete [] precpvADC[i];
368 fPMDLoader = fRunLoader->GetLoader("PMDLoader");
369 fPMDLoader->WriteRecPoints("OVERWRITE");
371 // delete the pointers
375 // cout << " ***** End::Digits2RecPoints :: Raw *****" << endl;
377 // ------------------------------------------------------------------------- //
378 void AliPMDClusterFinder::SetCellEdepCut(Float_t ecut)
382 // ------------------------------------------------------------------------- //
383 void AliPMDClusterFinder::SetDebug(Int_t idebug)
387 // ------------------------------------------------------------------------- //
388 void AliPMDClusterFinder::AddRecPoint(Int_t idet,Int_t ismn,Float_t *clusdata)
390 // Add Reconstructed points
392 TClonesArray &lrecpoints = *fRecpoints;
393 AliPMDrecpoint1 *newrecpoint;
394 newrecpoint = new AliPMDrecpoint1(idet, ismn, clusdata);
395 new(lrecpoints[fNpoint++]) AliPMDrecpoint1(newrecpoint);
398 // ------------------------------------------------------------------------- //
399 void AliPMDClusterFinder::ResetCellADC()
401 // Reset the individual cell ADC value to zero
403 for(Int_t irow = 0; irow < fgkRow; irow++)
405 for(Int_t icol = 0; icol < fgkCol; icol++)
407 fCellADC[irow][icol] = 0.;
411 // ------------------------------------------------------------------------- //
413 void AliPMDClusterFinder::ResetRecpoint()
415 // Clear the list of reconstructed points
417 if (fRecpoints) fRecpoints->Clear();
419 // ------------------------------------------------------------------------- //
420 void AliPMDClusterFinder::Load()
422 // Load all the *.root files
424 fPMDLoader->LoadDigits("READ");
425 fPMDLoader->LoadRecPoints("recreate");
427 // ------------------------------------------------------------------------- //
428 void AliPMDClusterFinder::LoadClusters()
430 // Load all the *.root files
432 fPMDLoader->LoadRecPoints("recreate");
434 // ------------------------------------------------------------------------- //
435 void AliPMDClusterFinder::UnLoad()
437 // Unload all the *.root files
439 fPMDLoader->UnloadDigits();
440 fPMDLoader->UnloadRecPoints();
442 // ------------------------------------------------------------------------- //
443 void AliPMDClusterFinder::UnLoadClusters()
445 // Unload all the *.root files
447 fPMDLoader->UnloadRecPoints();
449 // ------------------------------------------------------------------------- //