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