Line data Source code
1 : /**************************************************************************
2 : * Copyright(c) 1998-2003, 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 : // Implementation of the ITS clusterer V2 class //
17 : // //
18 : // Origin: Iouri Belikov, CERN, Jouri.Belikov@cern.ch //
19 : // Unfolding switch from AliITSRecoParam: D. Elia, INFN Bari //
20 : // //
21 : ////////////////////////////////////////////////////////////////////////////
22 :
23 :
24 : #include <TGeoGlobalMagField.h>
25 : #include "AliITSCalibrationSPD.h"
26 : #include "AliITSClusterFinderV2SPD.h"
27 : #include "AliITSRecPoint.h"
28 : #include "AliITSgeomTGeo.h"
29 : #include "AliITSDetTypeRec.h"
30 : #include "AliITSReconstructor.h"
31 : #include "AliRawReader.h"
32 : #include "AliITSRawStreamSPD.h"
33 : #include <TClonesArray.h>
34 : #include "AliITSdigitSPD.h"
35 : #include "AliITSFOSignalsSPD.h"
36 : #include "AliITSRecPointContainer.h"
37 : #include "AliMagF.h"
38 : #include "AliITSsegmentationSPD.h"
39 :
40 118 : ClassImp(AliITSClusterFinderV2SPD)
41 :
42 : //__________________________________________________________________________
43 2 : AliITSClusterFinderV2SPD::AliITSClusterFinderV2SPD(AliITSDetTypeRec* dettyp):AliITSClusterFinder(dettyp),
44 4 : fLastSPD1(AliITSgeomTGeo::GetModuleIndex(2,1,1)-1),
45 2 : fNySPD(256),
46 2 : fNzSPD(160),
47 2 : fYpitchSPD(0.0050),
48 2 : fZ1pitchSPD(0.0425),
49 2 : fZ2pitchSPD(0.0625),
50 2 : fHwSPD(0.64),
51 12 : fHlSPD(3.48){
52 :
53 : //Default constructor
54 :
55 2 : fYSPD[0]=0.5*fYpitchSPD;
56 1024 : for (Int_t m=1; m<fNySPD; m++) fYSPD[m]=fYSPD[m-1]+fYpitchSPD;
57 2 : fZSPD[0]=fZ1pitchSPD;
58 640 : for (Int_t m=1; m<fNzSPD; m++) {
59 318 : Double_t dz=fZ1pitchSPD;
60 954 : if (m==31 || m==32 || m==63 || m==64 || m==95 || m==96 ||
61 652 : m==127 || m==128) dz=fZ2pitchSPD;
62 318 : fZSPD[m]=fZSPD[m-1]+dz;
63 : }
64 644 : for (Int_t m=0; m<fNzSPD; m++) {
65 320 : Double_t dz=0.5*fZ1pitchSPD;
66 960 : if (m==31 || m==32 || m==63 || m==64 || m==95 || m==96 ||
67 656 : m==127 || m==128) dz=0.5*fZ2pitchSPD;
68 320 : fZSPD[m]-=dz;
69 : }
70 :
71 4 : }
72 : //__________________________________________________________________________
73 : void AliITSClusterFinderV2SPD::FindRawClusters(Int_t mod){
74 : //Find clusters V2
75 1920 : SetModule(mod);
76 960 : FindClustersSPD(fDigits);
77 :
78 960 : }
79 : //__________________________________________________________________________
80 : void AliITSClusterFinderV2SPD::RawdataToClusters(AliRawReader* rawReader){
81 : //------------------------------------------------------------
82 : // This function creates ITS clusters from raw data
83 : //------------------------------------------------------------
84 8 : rawReader->Reset();
85 4 : fNClusters = 0; //RS
86 4 : AliITSRawStreamSPD inputSPD(rawReader);
87 4 : FindClustersSPD(&inputSPD);
88 :
89 4 : }
90 : //__________________________________________________________________________
91 : Int_t AliITSClusterFinderV2SPD::ClustersSPD(AliBin* bins, TClonesArray* digits,TClonesArray* clusters,Int_t maxBins,Int_t nzbins,Int_t iModule,Bool_t rawdata){
92 :
93 : //Cluster finder for SPD (from digits and from rawdata)
94 : const Double_t kmictocm = 1.0e-4; // convert microns to cm.
95 : const Double_t defaultField = 5.0; // default Bz value at which Tan(theta_Lorentz) is given in RecoParam
96 :
97 : static AliITSRecoParam *repa = NULL;
98 1030 : if(!repa){
99 2 : repa = (AliITSRecoParam*) AliITSReconstructor::GetRecoParam();
100 2 : if(!repa){
101 0 : repa = AliITSRecoParam::GetHighFluxParam();
102 0 : AliWarning("Using default AliITSRecoParam class");
103 0 : }
104 : }
105 1030 : const TGeoHMatrix *mT2L=AliITSgeomTGeo::GetTracking2LocalMatrix(iModule);
106 :
107 : // Lorentz angle correction
108 : Double_t tanLorentzAngle=0;
109 1030 : AliITSsegmentationSPD *seg = (AliITSsegmentationSPD*)(GetDetTypeRec()->GetSegmentationModel(0));
110 1030 : Double_t thick = 0.5*kmictocm*seg->Dy(); // Half Thickness in cm
111 1030 : if(repa->GetCorrectLorentzAngleSPD()) { // only if CorrectLorentzAngleSPD required
112 : // here retrieve the value of the field
113 3090 : AliMagF* field = dynamic_cast<AliMagF*>(TGeoGlobalMagField::Instance()->GetField());
114 1030 : if (field == 0) {
115 0 : AliError("Cannot get magnetic field from TGeoGlobalMagField");
116 0 : }
117 : else {
118 1030 : Float_t magField = field->SolenoidField();
119 1030 : tanLorentzAngle=repa->GetTanLorentzAngleHolesSPD() * magField / defaultField ;
120 : }
121 1030 : }
122 : //
123 :
124 1030 : if (repa->GetSPDRemoveNoisyFlag()) {
125 : // Loop on noisy pixels and reset them
126 : AliITSCalibrationSPD *cal =
127 1030 : (AliITSCalibrationSPD*) fDetTypeRec->GetCalibrationModel(iModule);
128 2060 : for(Int_t ipix = 0; ipix<cal->GetNrBad(); ipix++){
129 0 : Int_t row, col;
130 0 : cal->GetBadPixel(ipix,row,col);
131 : //PH printf(" module %d row %d col %d \n",iModule,row,col);
132 0 : Int_t index = (row+1) * nzbins + (col+1);
133 :
134 0 : bins[index].SetQ(0);
135 0 : bins[index].SetMask(0xFFFFFFFE);
136 0 : }
137 1030 : }
138 :
139 1030 : if (repa->GetSPDRemoveDeadFlag()) {
140 : // Loop on dead pixels and reset them
141 : AliITSCalibrationSPD *cal =
142 1030 : (AliITSCalibrationSPD*) fDetTypeRec->GetSPDDeadModel(iModule);
143 1131 : if (cal->IsBad()) return 0; // if all ladder is dead, return to save time
144 849570 : for(Int_t ipix = 0; ipix<cal->GetNrBad(); ipix++){
145 423856 : Int_t row, col;
146 423856 : cal->GetBadPixel(ipix,row,col);
147 423856 : Int_t index = (row+1) * nzbins + (col+1);
148 423856 : bins[index].SetQ(0);
149 423856 : bins[index].SetMask(0xFFFFFFFE);
150 423856 : }
151 929 : AliITSCalibrationSPD *calSparse = (AliITSCalibrationSPD*) fDetTypeRec->GetSPDSparseDeadModel(iModule);
152 1858 : for(Int_t ipix = 0; ipix<calSparse->GetNrBad(); ipix++){
153 0 : Int_t row, col;
154 0 : calSparse->GetBadPixel(ipix,row,col);
155 0 : Int_t index = (row+1) * nzbins + (col+1);
156 0 : bins[index].SetQ(0);
157 0 : bins[index].SetMask(0xFFFFFFFE);
158 0 : }
159 :
160 929 : }
161 :
162 : Int_t nclu=0;
163 77658826 : for(Int_t iBin =0; iBin < maxBins;iBin++){
164 38828484 : if(bins[iBin].IsUsed()) continue;
165 146 : Int_t nBins = 0;
166 146 : Int_t idxBins[200];
167 146 : FindCluster(iBin,nzbins,bins,nBins,idxBins);
168 146 : if (nBins == 200){
169 0 : Error("ClustersSPD","SPD Too big cluster !\n");
170 0 : continue;
171 : }
172 146 : Int_t milab[10];
173 3212 : for(Int_t ilab=0;ilab<10;ilab++){
174 1460 : milab[ilab]=-2;
175 : }
176 146 : if(rawdata){
177 73 : milab[3]=fNdet[iModule];
178 73 : }
179 : Int_t ymin,ymax,zmin,zmax;
180 146 : if(rawdata){
181 73 : ymin = (idxBins[0] / nzbins) - 1;
182 : ymax = ymin;
183 73 : zmin = (idxBins[0] % nzbins) - 1;
184 : zmax = zmin;
185 73 : }
186 : else{
187 73 : AliITSdigitSPD* dig = (AliITSdigitSPD*)digits->UncheckedAt(idxBins[0]);
188 73 : ymin=dig->GetCoord2();
189 : ymax=ymin;
190 73 : zmin=dig->GetCoord1();
191 : zmax=zmin;
192 : }
193 : //PH if(iModule == 24 || iModule == 25) printf("\n");
194 840 : for (Int_t idx = 0; idx < nBins; idx++) {
195 : Int_t iy;
196 : Int_t iz;
197 274 : if(rawdata){
198 137 : iy = (idxBins[idx] / nzbins) - 1;
199 137 : iz = (idxBins[idx] % nzbins) - 1;
200 137 : }
201 : else{
202 137 : AliITSdigitSPD* dig = (AliITSdigitSPD*)digits->UncheckedAt(idxBins[idx]);
203 137 : iy = dig->GetCoord2();
204 137 : iz = dig->GetCoord1();
205 : //if(iModule == 24 || iModule == 25) printf(" || iy %d iz %d in Module %d \n",iy,iz,iModule);
206 : }
207 274 : if (ymin > iy) ymin = iy;
208 382 : if (ymax < iy) ymax = iy;
209 274 : if (zmin > iz) zmin = iz;
210 288 : if (zmax < iz) zmax = iz;
211 :
212 : }
213 146 : if(!rawdata){
214 420 : for(Int_t l=0;l<nBins;l++){
215 137 : AliITSdigitSPD* dig = (AliITSdigitSPD*)digits->UncheckedAt(idxBins[l]);
216 3014 : for(Int_t dlab=0;dlab<10;dlab++){
217 1370 : Int_t digitlab = (dig->GetTracks())[dlab];
218 2551 : if(digitlab<0) continue;
219 189 : AddLabel(milab,digitlab);
220 189 : }
221 137 : if (milab[9]>0) CheckLabels2(milab);
222 : }
223 73 : CheckLabels2(milab);
224 73 : }
225 :
226 : Int_t idy =0; //max 2 clusters
227 146 : if((iModule <= fLastSPD1) &&idy<3) idy=3;
228 146 : if((iModule > fLastSPD1) &&idy<4) idy=4;
229 : Int_t idz=3;
230 :
231 : // Switch the unfolding OFF/ON
232 146 : if(!repa->GetUseUnfoldingInClusterFinderSPD()) {
233 146 : idy=ymax-ymin+1;
234 146 : idz=zmax-zmin+1;
235 146 : }
236 : //
237 584 : for(Int_t iiz=zmin; iiz<=zmax;iiz+=idz){
238 584 : for(Int_t iiy=ymin;iiy<=ymax;iiy+=idy){
239 :
240 : Int_t ndigits=0;
241 : Float_t y=0.,z=0.,q=0.;
242 840 : for(Int_t idx=0;idx<nBins;idx++){
243 : Int_t iy;
244 : Int_t iz;
245 274 : if(rawdata){
246 137 : iy = (idxBins[idx] / nzbins)-1;
247 137 : iz = (idxBins[idx] % nzbins)-1;
248 137 : }
249 : else{
250 137 : AliITSdigitSPD* dig = (AliITSdigitSPD*)digits->UncheckedAt(idxBins[idx]);
251 137 : iy = dig->GetCoord2();
252 137 : iz = dig->GetCoord1();
253 : }
254 548 : if(zmax-zmin>=idz || ymax-ymin>=idy){
255 0 : if(TMath::Abs(iy-iiy)>0.75*idy) continue;
256 0 : if(TMath::Abs(iz-iiz)>0.75*idz) continue;
257 : }
258 274 : ndigits++;
259 : Float_t qBin=0.;
260 274 : if(rawdata) {
261 137 : qBin = bins[idxBins[idx]].GetQ();
262 137 : if (fRawID2ClusID) { // RS: Register cluster id in raw words list
263 0 : int rwid = bins[idxBins[idx]].GetRawID();
264 0 : if (fRawID2ClusID->GetSize()<=rwid) fRawID2ClusID->Set( (rwid+10)<<1 );
265 0 : (*fRawID2ClusID)[rwid] = milab[0] = fNClusters+1; // RS: store clID+1 as a reference to the cluster
266 0 : }
267 : }
268 274 : if(!rawdata){
269 137 : AliITSdigitSPD* dig = (AliITSdigitSPD*)digits->UncheckedAt(idxBins[idx]);
270 137 : qBin = (Float_t)dig->GetSignal();
271 137 : }
272 274 : y+= qBin * fYSPD[iy];
273 274 : z+= qBin * fZSPD[iz];
274 274 : q+= qBin;
275 274 : }// for idx
276 146 : if(ndigits==0) continue;
277 :
278 146 : y /= q;
279 146 : z /= q;
280 146 : y -= fHwSPD;
281 146 : z -= fHlSPD;
282 :
283 : // Lorentz drift effect in local y
284 146 : y -= tanLorentzAngle*thick;
285 : //
286 :
287 146 : Float_t hit[6]; //y,z,sigma(y)^2, sigma(z)^2, charge
288 : {
289 146 : Double_t loc[3]={y,0.,z},trk[3]={0.,0.,0.};
290 146 : mT2L->MasterToLocal(loc,trk);
291 146 : hit[0]=trk[1];
292 146 : hit[1]=trk[2];
293 146 : }
294 146 : hit[2] = fYpitchSPD*fYpitchSPD/12.;
295 146 : hit[3] = fZ1pitchSPD*fZ1pitchSPD/12.;
296 146 : hit[4] = 1.;
297 146 : hit[5] = 0.;
298 :
299 219 : if(!rawdata) milab[3]=fNdet[iModule];
300 146 : Int_t info[3] = {ymax-ymin+1,zmax-zmin+1,fNlayer[iModule]};
301 146 : if(!rawdata){
302 73 : AliITSRecPoint cl(milab,hit,info);
303 73 : cl.SetType(nBins);
304 73 : fDetTypeRec->AddRecPoint(cl);
305 73 : }
306 : else{
307 73 : Int_t label[4]={milab[0],milab[1],milab[2],milab[3]};
308 73 : AliITSRecPoint cl(label, hit,info);
309 73 : cl.SetType(nBins);
310 219 : new (clusters->AddrAt(nclu))
311 73 : AliITSRecPoint(cl);
312 73 : }
313 146 : nclu++;
314 146 : fNClusters++;
315 146 : }// for iiy
316 : }// for iiz
317 292 : }//end for iBin
318 : return nclu;
319 :
320 1030 : }
321 : //__________________________________________________________________________
322 : void AliITSClusterFinderV2SPD::FindClustersSPD(AliITSRawStreamSPD* input)
323 : {
324 : //------------------------------------------------------------
325 : // SPD cluster finder for raw data (this method is called once per event)
326 : // Now also fills fast-or fired map
327 : //------------------------------------------------------------
328 :
329 8 : AliITSRecPointContainer* rpc = AliITSRecPointContainer::Instance();
330 : Int_t nClustersSPD = 0;
331 4 : Int_t kNzBins = fNzSPD + 2;
332 4 : Int_t kNyBins = fNySPD + 2;
333 4 : Int_t kMaxBin = kNzBins * kNyBins;
334 334376 : AliBin *binsSPD = new AliBin[kMaxBin];
335 334376 : AliBin *binsSPDInit = new AliBin[kMaxBin];
336 : AliBin* bins = NULL;
337 :
338 : // read raw data input stream
339 : int countRW = 0; //RS
340 4 : if (fRawID2ClusID) fRawID2ClusID->Reset(); //RS if array was provided, we shall store the rawID -> ClusterID
341 : //
342 : while (kTRUE) {
343 152 : Bool_t next = input->Next();
344 300 : if (!next || input->IsNewModule()) {
345 74 : Int_t iModule = input->GetPrevModuleID();
346 :
347 : // when all data from a module was read, search for clusters
348 74 : if (bins) {
349 70 : TClonesArray* clusters = rpc->UncheckedGetClusters(iModule);
350 70 : Int_t nClusters = ClustersSPD(bins,0,clusters,kMaxBin,kNzBins,iModule,kTRUE);
351 70 : nClustersSPD += nClusters;
352 : bins = NULL;
353 70 : }
354 :
355 78 : if (!next) break;
356 : bins = binsSPD;
357 70 : memcpy(binsSPD,binsSPDInit,sizeof(AliBin)*kMaxBin);
358 70 : }
359 :
360 148 : if (next && bins) {
361 : // fill the current digit into the bins array
362 148 : Int_t index = (input->GetCoord2()+1) * kNzBins + (input->GetCoord1()+1);
363 148 : bins[index].SetIndex(index);
364 148 : bins[index].SetQ(1);
365 148 : bins[index].SetMask(1);
366 148 : bins[index].SetRawID(countRW); //RS register raw id
367 148 : }
368 148 : countRW++; //RS
369 148 : }
370 :
371 8 : delete [] binsSPDInit;
372 8 : delete [] binsSPD;
373 :
374 12 : AliDebug(1,Form("found clusters in ITS SPD: %d", nClustersSPD));
375 :
376 :
377 : // Fill the FastOr fired map
378 168 : for (UInt_t eq=0; eq<20; eq++) {
379 1120 : for (UInt_t hs=0; hs<6; hs++) {
380 10560 : for (UInt_t chip=0; chip<10; chip++) {
381 4800 : if (input->GetFastOrSignal(eq,hs,chip)) {
382 74 : fDetTypeRec->SetFastOrFiredMapOnline(eq,hs,chip);
383 74 : }
384 : }
385 : }
386 : }
387 :
388 4 : }
389 : //__________________________________________________________________________
390 : void AliITSClusterFinderV2SPD::FindClustersSPD(TClonesArray *digits) {
391 : //------------------------------------------------------------
392 : // SPD cluster finder for digits (this method is called for each module)
393 : // Now also fills the fast-or fired map
394 : //------------------------------------------------------------
395 :
396 :
397 1920 : Int_t kNzBins = fNzSPD + 2;
398 960 : const Int_t kMAXBIN=kNzBins*(fNySPD+2);
399 :
400 960 : Int_t ndigits=digits->GetEntriesFast();
401 80250240 : AliBin *bins=new AliBin[kMAXBIN];
402 :
403 : Int_t idig;
404 : AliITSdigitSPD *digit=0;
405 2216 : for (idig=0; idig<ndigits; idig++) {
406 148 : digit=(AliITSdigitSPD*)digits->UncheckedAt(idig);
407 148 : Int_t i=digit->GetCoord2()+1; //y
408 148 : Int_t j=digit->GetCoord1()+1;
409 148 : Int_t index=i*kNzBins+j;
410 :
411 148 : bins[index].SetIndex(idig);
412 148 : bins[index].SetQ(1);
413 148 : bins[index].SetMask(1);
414 : }
415 :
416 :
417 960 : Int_t nClustersSPD = ClustersSPD(bins,digits,0,kMAXBIN,kNzBins,fModule,kFALSE);
418 1920 : delete [] bins;
419 :
420 2880 : AliDebug(1,Form("found clusters in ITS SPD: %d", nClustersSPD));
421 :
422 : // Fill the FastOr fired map
423 960 : AliITSFOSignalsSPD* foSignals = fDetTypeRec->GetFOSignals();
424 960 : if (foSignals) {
425 960 : UInt_t eq = AliITSRawStreamSPD::GetOnlineEqIdFromOffline(fModule);
426 960 : UInt_t hs = AliITSRawStreamSPD::GetOnlineHSFromOffline(fModule);
427 11520 : for(UInt_t ch=0; ch<5; ch++) {
428 4800 : UInt_t chip = AliITSRawStreamSPD::GetOnlineChipFromOffline(fModule,ch*32);
429 4800 : if (foSignals->GetSignal(eq,hs,chip)) {
430 74 : fDetTypeRec->SetFastOrFiredMapOnline(eq,hs,chip);
431 74 : }
432 : }
433 960 : }
434 :
435 960 : }
|