Implement fixed due to Boris Batyunya. Including a fix for a bug with the
[u/mrichter/AliRoot.git] / ITS / AliITSClusterFinderSPD.cxx
CommitLineData
b0f5e3fc 1/**************************************************************************
2 * Copyright(c) 1998-1999, ALICE Experiment at CERN, All rights reserved. *
3 * *
4 * Author: The ALICE Off-line Project. *
5 * Contributors are mentioned in the code where appropriate. *
6 * *
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 **************************************************************************/
15
16
17#include "AliITSClusterFinderSPD.h"
e8189707 18#include "AliITSMapA1.h"
19#include "AliITS.h"
b0f5e3fc 20#include "AliRun.h"
21
22
23
24ClassImp(AliITSClusterFinderSPD)
25
26//----------------------------------------------------------
27AliITSClusterFinderSPD::AliITSClusterFinderSPD
28(AliITSsegmentation *seg, TClonesArray *digits, TClonesArray *recp)
29{
30 // constructor
31 fSegmentation=seg;
32 fDigits=digits;
33 fClusters=recp;
34 fNclusters= fClusters->GetEntriesFast();
35 SetDx();
36 SetDz();
37 SetMap();
38 SetNCells();
39}
40
41//_____________________________________________________________________________
42AliITSClusterFinderSPD::AliITSClusterFinderSPD()
43{
44 // constructor
45 fSegmentation=0;
46 fDigits=0;
47 fClusters=0;
48 fNclusters=0;
e8189707 49 fMap=0;
b0f5e3fc 50 SetDx();
51 SetDz();
b0f5e3fc 52 SetNCells();
53
54}
55
e8189707 56//_____________________________________________________________________________
57AliITSClusterFinderSPD::~AliITSClusterFinderSPD()
58{
59 // destructor
60 if (fMap) delete fMap;
61
62
63}
b0f5e3fc 64//__________________________________________________________________________
65AliITSClusterFinderSPD::AliITSClusterFinderSPD(const AliITSClusterFinderSPD &source){
66 // Copy Constructor
67 if(&source == this) return;
68 this->fClusters = source.fClusters ;
69 this->fNclusters = source.fNclusters ;
70 this->fMap = source.fMap ;
71 this->fDz = source.fDz ;
72 this->fDx = source.fDx ;
73 this->fMinNCells = source.fMinNCells ;
74 return;
75}
76
77//_________________________________________________________________________
78AliITSClusterFinderSPD&
79 AliITSClusterFinderSPD::operator=(const AliITSClusterFinderSPD &source) {
80 // Assignment operator
81 if(&source == this) return *this;
82 this->fClusters = source.fClusters ;
83 this->fNclusters = source.fNclusters ;
84 this->fMap = source.fMap ;
85 this->fDz = source.fDz ;
86 this->fDx = source.fDx ;
87 this->fMinNCells = source.fMinNCells ;
88 return *this;
89}
90
91//_____________________________________________________________________________
92void AliITSClusterFinderSPD::SetMap()
93{
94 // set map
a3e16987 95
b0f5e3fc 96 if(!fMap) fMap=new AliITSMapA1(fSegmentation,fDigits);
a3e16987 97
b0f5e3fc 98}
99
100//_____________________________________________________________________________
101
102void AliITSClusterFinderSPD::Find1DClusters()
103{
104 // Find one dimensional clusters, i.e.
105 // in r*phi(x) direction for each colunm in z direction
106
107 AliITS *iTS=(AliITS*)gAlice->GetModule("ITS");
108
109 // retrieve the parameters
110 Int_t fNofPixels = fSegmentation->Npz();
111 Int_t fMaxNofSamples = fSegmentation->Npx();
112
113 // read in digits -> do not apply threshold
114 // signal in fired pixels is always 1
a3e16987 115
b0f5e3fc 116 fMap->FillMap();
117
118 Int_t nofFoundClusters = 0;
119
120 Int_t k,it,m;
121 for(k=0;k<fNofPixels;k++) {
122
123 Int_t mmax = 10; // a size of the window for the cluster finding
124
125 for(it=0;it<fMaxNofSamples;it++) {
126
127 Int_t lclx = 0;
128 Int_t xstart = 0;
129 Int_t xstop = 0;
130 Int_t id = 0;
131 Int_t ilcl =0;
132
133 for(m=0;m<mmax;m++) { // find the cluster inside the window
134 id = it+m;
135 if(id >= fMaxNofSamples) break; // ! no possible for the fadc
136
e8189707 137 if(fMap->TestHit(k,id) == kUnused) { // start of the cluster
b0f5e3fc 138 lclx += 1;
139 if(lclx == 1) xstart = id;
140
141 }
142
e8189707 143 if(lclx > 0 && fMap->TestHit(k,id) == kEmpty) {
b0f5e3fc 144 // end of cluster if a gap exists
145 xstop = id-1;
146 ilcl = 1;
147 break;
148 }
149
150 } // end of m-loop
151
152 if(lclx == 0 && ilcl == 0) it = id; // no cluster in the window,
153 // continue the "it" loop
154
155 if(id >= fMaxNofSamples && lclx == 0) break; // the x row finished
156
157 if(id < fMaxNofSamples && ilcl == 0 && lclx > 0) {
e8189707 158 // cluster end is outside of the window,
b0f5e3fc 159 mmax += 5; // increase mmax and repeat the cluster
e8189707 160 // finding
b0f5e3fc 161 it -= 1;
162 }
163
164 if(id >= fMaxNofSamples && lclx > 0) { // the x row finished but
e8189707 165 xstop = fMaxNofSamples - 1; // the end cluster exists
b0f5e3fc 166 ilcl = 1;
167 }
168
169 // --- Calculate z and x coordinates for one dimensional clusters
170
171 if(ilcl == 1) { // new cluster exists
172 it = id;
173 mmax = 10;
174 nofFoundClusters++;
175 Float_t clusterCharge = 0.;
b0f5e3fc 176 Float_t zpitch = fSegmentation->Dpz(k+1);
177 Float_t clusterZ, dummyX;
178 Int_t dummy=0;
e8189707 179 fSegmentation->GetPadCxz(dummy,k,dummyX,clusterZ);
b0f5e3fc 180 Float_t zstart = clusterZ - 0.5*zpitch;
181 Float_t zstop = clusterZ + 0.5*zpitch;
182 Float_t clusterX = 0.;
183 Int_t xstartfull = xstart;
184 Int_t xstopfull = xstop;
185 Int_t clusterSizeX = lclx;
186 Int_t clusterSizeZ = 1;
187
188 Int_t its;
189 for(its=xstart; its<=xstop; its++) {
190 Int_t firedpixel=0;
e8189707 191 if (fMap->GetHitIndex(k,its)>=0) firedpixel=1;
b0f5e3fc 192 clusterCharge += firedpixel;
193 clusterX +=its + 0.5;
194 }
195 Float_t fRphiPitch = fSegmentation->Dpx(dummy);
196 clusterX /= (clusterSizeX/fRphiPitch); // center of gravity for x
197
198
199 //printf("ClusterZ ClusterX %f %f \n",clusterZ, clusterX);
200
201 // Int_t nclusters = fClusters->GetEntriesFast();
202 // cout << nclusters << " clusters" << endl;
203 // cout<< "Create point"<<endl;
204
205 // Write the points (coordinates and some cluster information) to the
206 // AliITSRawClusterSPD object
207
208 AliITSRawClusterSPD *clust = new AliITSRawClusterSPD(clusterZ,clusterX,clusterCharge,clusterSizeZ,clusterSizeX,xstart,xstop,xstartfull,xstopfull,zstart,zstop,k);
209 // fClusters->Add(point);
210 iTS->AddCluster(0,clust);
211 // cout << "Cluster at Ladder: " << fLadder << ", Detector: " <<fDetector<<endl;
212
213 // cout<<" end of cluster finding for Z pixel "<<endl;
214
215 } // new cluster (ilcl=1)
216 } // X direction loop (it)
217 } // Z direction loop (k)
e8189707 218
219 //fMap->ClearMap();
b0f5e3fc 220 return;
221
222}
223
224//_____________________________________________________________________________
225void AliITSClusterFinderSPD::GroupClusters()
226{
227 // Find two dimensional clusters, i.e. group one dimensional clusters
228 // into two dimensional ones (go both in x and z directions).
229
230 // get number of clusters for this module
231 Int_t nofClusters = fClusters->GetEntriesFast();
232 nofClusters -= fNclusters;
b0f5e3fc 233
234 AliITSRawClusterSPD *clusterI;
235 AliITSRawClusterSPD *clusterJ;
236
e8189707 237 Int_t *label=new Int_t[nofClusters];
b0f5e3fc 238 Int_t i,j;
239 for(i=0; i<nofClusters; i++) label[i] = 0;
240 for(i=0; i<nofClusters; i++) {
241 if(label[i] != 0) continue;
242 for(j=i+1; j<nofClusters; j++) {
243 if(label[j] != 0) continue;
244 clusterI = (AliITSRawClusterSPD*) fClusters->At(i);
245 clusterJ = (AliITSRawClusterSPD*) fClusters->At(j);
246 Bool_t pair = clusterI->Brother(clusterJ,fDz,fDx);
247 if(pair) {
248
249 // if((clusterI->XStop() == clusterJ->XStart()-1)||(clusterI->XStart()==clusterJ->XStop()+1)) cout<<"!! Diagonal cluster"<<endl;
250 /*
251 cout << "clusters " << i << "," << j << " before grouping" << endl;
e8189707 252 clusterI->PrintInfo();
253 clusterJ->PrintInfo();
b0f5e3fc 254 */
255 clusterI->Add(clusterJ);
256 // cout << "remove cluster " << j << endl;
257 label[j] = 1;
258 fClusters->RemoveAt(j);
259 /*
260 cout << "cluster " << i << " after grouping" << endl;
e8189707 261 clusterI->PrintInfo();
b0f5e3fc 262 */
263 } // pair
264 } // J clusters
265 label[i] = 1;
266 } // I clusters
267 fClusters->Compress();
a3e16987 268 // Int_t totalNofClusters = fClusters->GetEntriesFast();
269 // cout << " Nomber of clusters at the group end ="<< totalNofClusters<<endl;
b0f5e3fc 270
dccfe216 271 delete [] label;
e8189707 272
b0f5e3fc 273 return;
274
275
276}
277//_____________________________________________________________________________
278
a3e16987 279void AliITSClusterFinderSPD::TracksInCluster()
280{
281
282 // Find tracks creating one cluster
283
284 // get number of clusters for this module
285 Int_t nofClusters = fClusters->GetEntriesFast();
286 nofClusters -= fNclusters;
287
288 Int_t i, ix, iz, jx, jz, xstart, xstop, zstart, zstop, nclx, nclz;
289 // Int_t signal, track0, track1, track2;
290 Int_t trmax = 100;
291 Int_t cltracks[trmax], itr, tracki, ii, is, js, ie, ntr, tr0, tr1, tr2;
292
293 for(i=0; i<nofClusters; i++) {
294 ii = 0;
295 memset(cltracks,-1,sizeof(int)*trmax);
296 tr0=tr1=tr2=-1;
297
298 AliITSRawClusterSPD *clusterI = (AliITSRawClusterSPD*) fClusters->At(i);
299
300 nclx = clusterI->NclX();
301 nclz = clusterI->NclZ();
302 xstart = clusterI->XStartf();
303 xstop = clusterI->XStopf();
304 zstart = clusterI->Zend()-nclz+1;
305 zstop = clusterI->Zend();
306
307 Int_t ind;
308
309 for(iz=0; iz<nclz; iz++) {
310 jz = zstart + iz;
311 for(ix=0; ix<nclx; ix++) {
312 jx = xstart + ix;
313 ind = fMap->GetHitIndex(jz,jx);
314 if(ind == 0 && iz >= 0 && ix > 0) {
315 continue;
316 }
317 if(ind == 0 && iz > 0 && ix >= 0) {
318 continue;
319 }
320 if(ind == 0 && iz == 0 && ix == 0 && i > 0) {
321 continue;
322 }
323
324 AliITSdigitSPD *dig = (AliITSdigitSPD*)fMap->GetHit(jz,jx);
325 /*
326 signal=dig->fSignal;
327 track0=dig->fTracks[0];
328 track1=dig->fTracks[1];
329 track2=dig->fTracks[2];
330 */
331 for(itr=0; itr<3; itr++) {
332 tracki = dig->fTracks[itr];
333 if(tracki >= 0) {
334 ii += 1;
335 cltracks[ii-1] = tracki;
336 }
337 }
338 } // ix pixel
339 } // iz pixel
340
341 for(is=0; is<trmax; is++) {
342 if(cltracks[is]<0) continue;
343 for(js=is+1; js<trmax; js++) {
344 if(cltracks[js]<0) continue;
345 if(cltracks[js]==cltracks[is]) cltracks[js]=-5;
346 }
347 }
348
349 ntr = 0;
350 for(ie=0; ie<trmax; ie++) {
351 if(cltracks[ie] >= 0) {
352 ntr=ntr+1;
353 if(ntr==1) tr0=cltracks[ie];
354 if(ntr==2) tr1=cltracks[ie];
355 if(ntr==3) tr2=cltracks[ie];
356 }
357 }
358 // if delta ray only
359 if(ntr == 0) ntr = 1;
360
361 clusterI->SetNTracks(ntr);
362 clusterI->SetTracks(tr0,tr1,tr2);
363
364 } // I cluster
365
366}
367//_____________________________________________________________________________
368
b0f5e3fc 369void AliITSClusterFinderSPD::GetRecPoints()
370{
371 // get rec points
372 AliITS *iTS=(AliITS*)gAlice->GetModule("ITS");
373
374 // get number of clusters for this module
375 Int_t nofClusters = fClusters->GetEntriesFast();
376 nofClusters -= fNclusters;
b0f5e3fc 377 const Float_t kconv = 1.0e-4;
378 const Float_t kRMSx = 12.0*kconv; // microns -> cm ITS TDR Table 1.3
379 const Float_t kRMSz = 70.0*kconv; // microns -> cm ITS TDR Table 1.3
380
e8189707 381 Float_t spdLength = fSegmentation->Dz();
382 Float_t spdWidth = fSegmentation->Dx();
383
a3e16987 384 Int_t i;
385 Int_t track0, track1, track2;
386
b0f5e3fc 387 for(i=0; i<nofClusters; i++) {
a3e16987 388
b0f5e3fc 389 AliITSRawClusterSPD *clusterI = (AliITSRawClusterSPD*) fClusters->At(i);
a3e16987 390 clusterI->GetTracks(track0, track1, track2);
b0f5e3fc 391 AliITSRecPoint rnew;
a3e16987 392
e8189707 393 rnew.SetX((clusterI->X() - spdWidth/2)*kconv);
394 rnew.SetZ((clusterI->Z() - spdLength/2)*kconv);
b0f5e3fc 395 rnew.SetQ(1.);
396 rnew.SetdEdX(0.);
397 rnew.SetSigmaX2(kRMSx*kRMSx);
398 rnew.SetSigmaZ2(kRMSz*kRMSz);
a3e16987 399 rnew.fTracks[0]=track0;
400 rnew.fTracks[1]=track1;
401 rnew.fTracks[2]=track2;
b0f5e3fc 402 iTS->AddRecPoint(rnew);
403 } // I clusters
e8189707 404
405 fMap->ClearMap();
b0f5e3fc 406
407}
408//_____________________________________________________________________________
409
410void AliITSClusterFinderSPD::FindRawClusters()
411{
412 // find raw clusters
413 Find1DClusters();
414 GroupClusters();
a3e16987 415 TracksInCluster();
b0f5e3fc 416 GetRecPoints();
e8189707 417
b0f5e3fc 418}
419