cluster finder
[u/mrichter/AliRoot.git] / PMD / AliPMDClusterFinder.cxx
CommitLineData
01709453 1//-----------------------------------------------------//
2// //
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) //
7// //
8//-----------------------------------------------------//
9
10#include <Riostream.h>
11#include <TBRIK.h>
12#include <TNode.h>
13#include <TTree.h>
14#include <TGeometry.h>
15#include <TObjArray.h>
16#include <TClonesArray.h>
17#include <TFile.h>
18#include <TNtuple.h>
19#include <TParticle.h>
20
21#include "AliRun.h"
22#include "AliPMD.h"
23#include "AliDetector.h"
24#include "AliRunLoader.h"
25#include "AliLoader.h"
26#include "AliHeader.h"
27
28#include "AliPMDdigit.h"
29#include "AliPMDClusterFinder.h"
30#include "AliPMDClustering.h"
31#include "AliPMDContainer.h"
32#include "AliPMDcluster.h"
33#include "AliPMDrecpoint.h"
34
35
36ClassImp(AliPMDClusterFinder)
37//
38// Constructor
39//
40AliPMDClusterFinder::AliPMDClusterFinder()
41{
42 if (!fRecpoints) fRecpoints = new TClonesArray("AliPMDrecpoint", 1000);
43 fNpoint = 0;
44
45 for (Int_t i = 0; i < fTotSM; i++)
46 {
47 for (Int_t j = 0; j < fNCell; j++)
48 {
49 for (Int_t k = 0; k < fNCell; k++)
50 {
51 fCPV[i][j][k] = 0.;
52 fPMD[i][j][k] = 0.;
53 }
54 }
55 }
56
57}
58AliPMDClusterFinder::~AliPMDClusterFinder()
59{
60 delete fRecpoints;
61}
62//
63// Member functions
64//
65void AliPMDClusterFinder::OpengAliceFile(Char_t *file, Option_t *option)
66{
67
68 fRunLoader = AliRunLoader::Open(file,AliConfig::fgkDefaultEventFolderName,
69 "UPDATE");
70
71 if (!fRunLoader)
72 {
73 Error("Open","Can not open session for file %s.",file);
74 }
75
76 fRunLoader->LoadgAlice();
77 fRunLoader->LoadHeader();
78 fRunLoader->LoadKinematics();
79
80 gAlice = fRunLoader->GetAliRun();
81
82 if (gAlice)
83 {
84 printf("<AliPMDdigitizer::Open> ");
85 printf("AliRun object found on file.\n");
86 }
87 else
88 {
89 printf("<AliPMDdigitizer::Open> ");
90 printf("Could not find AliRun object.\n");
91 }
92 PMD = (AliPMD*)gAlice->GetDetector("PMD");
93 pmdloader = fRunLoader->GetLoader("PMDLoader");
94 if (pmdloader == 0x0)
95 {
96 cerr<<"OpengAlice : Can not find PMD or PMDLoader\n";
97 }
98
99 const char *cDR = strstr(option,"DR");
100
101 if (cDR)
102 {
103 pmdloader->LoadDigits("READ");
104 pmdloader->LoadRecPoints("recreate");
105 }
106}
107
108void AliPMDClusterFinder::Digits2RecPoints(Int_t ievt)
109{
110 Int_t det,smn;
111 Int_t cellno;
112 Int_t xpos,ypos;
113 Float_t adc;
114 Int_t isup, ix, iy;
115 Int_t idet;
116 Double_t d[72][72];
117 Float_t clusdata[7];
118 Int_t fMessage = 1;
119
120 fRunLoader->GetEvent(ievt);
121 //cout << " ***** Beginning::Digits2RecPoints *****" << endl;
122 treeD = pmdloader->TreeD();
123 if (treeD == 0x0)
124 {
125 cout << " Can not get TreeD" << endl;
126 }
127 AliPMDdigit *pmddigit;
128 TBranch *branch = treeD->GetBranch("PMDDigit");
129 branch->SetAddress(&fDigits);
130
131 ResetRecpoint();
132 treeR = pmdloader->TreeR();
133 if (treeR == 0x0)
134 {
135 pmdloader->MakeTree("R");
136 treeR = pmdloader->TreeR();
137 }
138
139 Int_t bufsize = 16000;
140 treeR->Branch("PMDRecpoint", &fRecpoints, bufsize);
141
142 Int_t nmodules = (Int_t) treeD->GetEntries();
143
144 for (Int_t imodule = 0; imodule < nmodules; imodule++)
145 {
146 treeD->GetEntry(imodule);
147 Int_t nentries = fDigits->GetLast();
148 for (Int_t ient = 0; ient < nentries+1; ient++)
149 {
150 pmddigit = (AliPMDdigit*)fDigits->UncheckedAt(ient);
151
152 det = pmddigit->GetDetector();
153 smn = pmddigit->GetSMNumber();
154 cellno = pmddigit->GetCellNumber();
155 adc = pmddigit->GetADC();
156
157 ypos = cellno/fNCell;
158 xpos = cellno - ypos*fNCell;
159
160 if (det == 0)
161 {
162 fPMD[smn][xpos][ypos] = adc;
163 }
164 else if(det == 1)
165 {
166 fCPV[smn][xpos][ypos] = adc;
167 }
168 }
169 } // modules
170
171 //
172 // Clustering started
173 //
174
175 // AliPMDContainer *pmdcont = new AliPMDContainer;
176
177 TObjArray *pmdcont = new TObjArray();
178
179 AliPMDcluster *pmdcl = new AliPMDcluster;
180
181 AliPMDClustering *pmdclust = new AliPMDClustering();
182 pmdclust->SetMessage(fMessage);
183
184 for (idet = 0; idet < 2; idet++)
185 {
186 for (isup = 0; isup < fTotSM; isup++)
187 {
188 for (ix = 0; ix < fNCell; ix++)
189 {
190 for (iy = 0; iy < fNCell; iy++)
191 {
192 if (idet == 0)
193 {
194 d[ix][iy] = (Double_t) fPMD[isup][ix][iy];
195 }
196 else if (idet == 1)
197 {
198 d[ix][iy] = (Double_t) fCPV[isup][ix][iy];
199 }
200
201 }
202 }
203 pmdclust->DoClust(idet,isup,d,pmdcont);
204
205 Int_t nentries = pmdcont->GetEntries();
206 cout << " nentries = " << nentries << endl;
207 for (Int_t ient = 0; ient < nentries; ient++)
208 {
209 clusdata[0] = (Float_t) idet;
210 clusdata[1] = (Float_t) isup;
211
212 pmdcl = (AliPMDcluster*)pmdcont->UncheckedAt(ient);
213
214 clusdata[2] = pmdcl->GetClusX();
215 clusdata[3] = pmdcl->GetClusY();
216 clusdata[4] = pmdcl->GetClusADC();
217 clusdata[5] = pmdcl->GetClusCells();
218 clusdata[6] = pmdcl->GetClusRadius();
219
220 AddRecPoint(clusdata);
221 }
222 pmdcont->Clear();
223
224 treeR->Fill();
225 ResetRecpoint();
226 } // SuperModule
227 } // Detector
228
229 ResetCellADC();
230
231 pmdloader->WriteRecPoints("OVERWRITE");
232
233 // delete the pointers
234 delete pmdclust;
235 delete pmdcont;
236
237 // cout << " ***** End::Digits2RecPoints *****" << endl;
238}
239
240
241void AliPMDClusterFinder::AddRecPoint(Float_t *clusdata)
242{
243 TClonesArray &lrecpoints = *fRecpoints;
244 AliPMDrecpoint *newrecpoint;
245 newrecpoint = new AliPMDrecpoint(clusdata);
246 new(lrecpoints[fNpoint++]) AliPMDrecpoint(newrecpoint);
247 delete newrecpoint;
248}
249void AliPMDClusterFinder::ResetCellADC()
250{
251 for (Int_t i = 0; i < fTotSM; i++)
252 {
253 for (Int_t j = 0; j < fNCell; j++)
254 {
255 for (Int_t k = 0; k < fNCell; k++)
256 {
257 fCPV[i][j][k] = 0.;
258 fPMD[i][j][k] = 0.;
259 }
260 }
261 }
262}
263
264void AliPMDClusterFinder::ResetRecpoint()
265{
266 fNpoint = 0;
267 if (fRecpoints) fRecpoints->Clear();
268}
269void AliPMDClusterFinder::UnLoad(Option_t *option)
270{
271 const char *cR = strstr(option,"R");
272
273 fRunLoader->UnloadgAlice();
274 fRunLoader->UnloadHeader();
275 fRunLoader->UnloadKinematics();
276
277 if (cR)
278 {
279 pmdloader->UnloadDigits();
280 }
281}