1 //-----------------------------------------------------//
3 // Date : August 05 2003 //
4 // This reads the file PMD.digits.root(TreeD), //
5 // calls the Clustering algorithm and stores the //
6 // clustering output in PMD.RecPoints.root(TreeR) //
8 //-----------------------------------------------------//
10 #include <Riostream.h>
14 #include <TGeometry.h>
15 #include <TObjArray.h>
16 #include <TClonesArray.h>
19 #include <TParticle.h>
23 #include "AliDetector.h"
24 #include "AliRunLoader.h"
25 #include "AliLoader.h"
26 #include "AliHeader.h"
28 #include "AliPMDdigit.h"
29 #include "AliPMDClusterFinder.h"
30 #include "AliPMDClustering.h"
31 #include "AliPMDcluster.h"
32 #include "AliPMDrecpoint1.h"
35 ClassImp(AliPMDClusterFinder)
39 AliPMDClusterFinder::AliPMDClusterFinder()
41 if (!fRecpoints) fRecpoints = new TClonesArray("AliPMDrecpoint", 1000);
44 for (Int_t i = 0; i < fTotSM; i++)
46 for (Int_t j = 0; j < fNCell; j++)
48 for (Int_t k = 0; k < fNCell; k++)
57 AliPMDClusterFinder::~AliPMDClusterFinder()
64 void AliPMDClusterFinder::OpengAliceFile(Char_t *file, Option_t *option)
67 fRunLoader = AliRunLoader::Open(file,AliConfig::fgkDefaultEventFolderName,
72 Error("Open","Can not open session for file %s.",file);
75 fRunLoader->LoadgAlice();
76 fRunLoader->LoadHeader();
77 fRunLoader->LoadKinematics();
79 gAlice = fRunLoader->GetAliRun();
83 printf("<AliPMDdigitizer::Open> ");
84 printf("AliRun object found on file.\n");
88 printf("<AliPMDdigitizer::Open> ");
89 printf("Could not find AliRun object.\n");
91 PMD = (AliPMD*)gAlice->GetDetector("PMD");
92 pmdloader = fRunLoader->GetLoader("PMDLoader");
95 cerr<<"OpengAlice : Can not find PMD or PMDLoader\n";
98 const char *cDR = strstr(option,"DR");
102 pmdloader->LoadDigits("READ");
103 pmdloader->LoadRecPoints("recreate");
107 void AliPMDClusterFinder::Digits2RecPoints(Int_t ievt)
119 fRunLoader->GetEvent(ievt);
120 //cout << " ***** Beginning::Digits2RecPoints *****" << endl;
121 treeD = pmdloader->TreeD();
124 cout << " Can not get TreeD" << endl;
126 AliPMDdigit *pmddigit;
127 TBranch *branch = treeD->GetBranch("PMDDigit");
128 branch->SetAddress(&fDigits);
131 treeR = pmdloader->TreeR();
134 pmdloader->MakeTree("R");
135 treeR = pmdloader->TreeR();
138 Int_t bufsize = 16000;
139 treeR->Branch("PMDRecpoint", &fRecpoints, bufsize);
141 Int_t nmodules = (Int_t) treeD->GetEntries();
143 for (Int_t imodule = 0; imodule < nmodules; imodule++)
145 treeD->GetEntry(imodule);
146 Int_t nentries = fDigits->GetLast();
147 for (Int_t ient = 0; ient < nentries+1; ient++)
149 pmddigit = (AliPMDdigit*)fDigits->UncheckedAt(ient);
151 det = pmddigit->GetDetector();
152 smn = pmddigit->GetSMNumber();
153 cellno = pmddigit->GetCellNumber();
154 adc = pmddigit->GetADC();
156 ypos = cellno/fNCell;
157 xpos = cellno - ypos*fNCell;
161 fPMD[smn][xpos][ypos] = adc;
165 fCPV[smn][xpos][ypos] = adc;
171 // Clustering started
174 TObjArray *pmdcont = new TObjArray();
176 AliPMDcluster *pmdcl = new AliPMDcluster;
178 AliPMDClustering *pmdclust = new AliPMDClustering();
179 pmdclust->SetMessage(fMessage);
181 for (idet = 0; idet < 2; idet++)
183 for (isup = 0; isup < fTotSM; isup++)
185 for (ix = 0; ix < fNCell; ix++)
187 for (iy = 0; iy < fNCell; iy++)
191 d[ix][iy] = (Double_t) fPMD[isup][ix][iy];
195 d[ix][iy] = (Double_t) fCPV[isup][ix][iy];
200 pmdclust->DoClust(idet,isup,d,pmdcont);
202 Int_t nentries = pmdcont->GetEntries();
203 cout << " nentries = " << nentries << endl;
204 for (Int_t ient = 0; ient < nentries; ient++)
206 clusdata[0] = (Float_t) idet;
207 clusdata[1] = (Float_t) isup;
209 pmdcl = (AliPMDcluster*)pmdcont->UncheckedAt(ient);
211 clusdata[2] = pmdcl->GetClusX();
212 clusdata[3] = pmdcl->GetClusY();
213 clusdata[4] = pmdcl->GetClusADC();
214 clusdata[5] = pmdcl->GetClusCells();
215 clusdata[6] = pmdcl->GetClusRadius();
217 AddRecPoint(clusdata);
228 pmdloader->WriteRecPoints("OVERWRITE");
230 // delete the pointers
234 // cout << " ***** End::Digits2RecPoints *****" << endl;
238 void AliPMDClusterFinder::AddRecPoint(Float_t *clusdata)
240 TClonesArray &lrecpoints = *fRecpoints;
241 AliPMDrecpoint *newrecpoint;
242 newrecpoint = new AliPMDrecpoint(clusdata);
243 new(lrecpoints[fNpoint++]) AliPMDrecpoint(newrecpoint);
246 void AliPMDClusterFinder::ResetCellADC()
248 for (Int_t i = 0; i < fTotSM; i++)
250 for (Int_t j = 0; j < fNCell; j++)
252 for (Int_t k = 0; k < fNCell; k++)
261 void AliPMDClusterFinder::ResetRecpoint()
264 if (fRecpoints) fRecpoints->Clear();
266 void AliPMDClusterFinder::UnLoad(Option_t *option)
268 const char *cR = strstr(option,"R");
270 fRunLoader->UnloadgAlice();
271 fRunLoader->UnloadHeader();
272 fRunLoader->UnloadKinematics();
276 pmdloader->UnloadDigits();