cb876016f076a0aa7ac8ec1c40024e444125b80e
[u/mrichter/AliRoot.git] / ITS / AliITSMultReconstructor.cxx
1 /**************************************************************************
2  * Copyright(c) 2007-2009, 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 /* $Id$ */
17
18 //____________________________________________________________________
19 // 
20 // AliITSMultReconstructor - find clusters in the pixels (theta and
21 // phi) and tracklets.
22 // 
23 // These can be used to extract charged particles multiplcicity from the ITS.
24 //
25 // A tracklet consist of two ITS clusters, one in the first pixel
26 // layer and one in the second. The clusters are associates if the 
27 // differencies in Phi (azimuth) and Zeta (longitudinal) are inside 
28 // a fiducial volume. In case of multiple candidates it is selected the
29 // candidate with minimum distance in Phi. 
30 // The parameter AssociationChoice allows to control if two clusters
31 // in layer 2 can be associated to the same cluster in layer 1 or not.
32 // (TRUE means double associations exluded; default = TRUE)
33 //
34 // Two methods return the number of traklets and the number of clusters 
35 // in the first SPD layer (GetNTracklets GetNSingleClusters)
36 //
37 // -----------------------------------------------------------------
38 // 
39 // NOTE: The cuts on phi and zeta depend on the interacting system (p-p  
40 //  or Pb-Pb). Please, check the file AliITSMultReconstructor.h and be 
41 //  sure that SetPhiWindow and SetZetaWindow are defined accordingly.
42 // 
43 //  Author :  Tiziano Virgili 
44 //  
45 //  Recent updates (D. Elia, INFN Bari):
46 //     - multiple association forbidden (fOnlyOneTrackletPerC2 = kTRUE) 
47 //     - phi definition changed to ALICE convention (0,2*TMath::pi()) 
48 //     - cluster coordinates taken with GetGlobalXYZ()
49 //     - fGeometry removed
50 //     - number of fired chips on the two layers
51 //
52 //____________________________________________________________________
53
54 #include <TClonesArray.h>
55 #include <TH1F.h>
56 #include <TH2F.h>
57 #include <TTree.h>
58
59 #include "AliITSMultReconstructor.h"
60 #include "AliITSsegmentationSPD.h"
61 #include "AliITSRecPoint.h"
62 #include "AliITSgeom.h"
63 #include "AliLog.h"
64
65 //____________________________________________________________________
66 ClassImp(AliITSMultReconstructor)
67
68
69 //____________________________________________________________________
70 AliITSMultReconstructor::AliITSMultReconstructor():
71 fClustersLay1(0),
72 fClustersLay2(0),
73 fTracklets(0),
74 fSClusters(0),
75 fAssociationFlag(0),
76 fNClustersLay1(0),
77 fNClustersLay2(0),
78 fNTracklets(0),
79 fNSingleCluster(0),
80 fPhiWindow(0),
81 fZetaWindow(0),
82 fOnlyOneTrackletPerC2(0),
83 fHistOn(0),
84 fhClustersDPhiAcc(0),
85 fhClustersDThetaAcc(0),
86 fhClustersDZetaAcc(0),
87 fhClustersDPhiAll(0),
88 fhClustersDThetaAll(0),
89 fhClustersDZetaAll(0),
90 fhDPhiVsDThetaAll(0),
91 fhDPhiVsDThetaAcc(0),
92 fhDPhiVsDZetaAll(0),
93 fhDPhiVsDZetaAcc(0),
94 fhetaTracklets(0),
95 fhphiTracklets(0),
96 fhetaClustersLay1(0),
97 fhphiClustersLay1(0){
98
99   fNFiredChips[0] = 0;
100   fNFiredChips[1] = 0;
101
102   // Method to reconstruct the charged particles multiplicity with the 
103   // SPD (tracklets).
104
105
106   SetHistOn();
107   SetPhiWindow();
108   SetZetaWindow();
109   SetOnlyOneTrackletPerC2();
110
111   fClustersLay1       = new Float_t*[300000];
112   fClustersLay2       = new Float_t*[300000];
113   fTracklets          = new Float_t*[300000];
114   fSClusters          = new Float_t*[300000];
115   fAssociationFlag    = new Bool_t[300000];
116
117   for(Int_t i=0; i<300000; i++) {
118     fClustersLay1[i]       = new Float_t[6];
119     fClustersLay2[i]       = new Float_t[6];
120     fTracklets[i]          = new Float_t[5];
121     fSClusters[i]           = new Float_t[2];
122     fAssociationFlag[i]    = kFALSE;
123   }
124
125   // definition of histograms
126   fhClustersDPhiAcc   = new TH1F("dphiacc",  "dphi",  100,0.,0.1);
127   fhClustersDPhiAcc->SetDirectory(0);
128   fhClustersDThetaAcc = new TH1F("dthetaacc","dtheta",100,-0.1,0.1);
129   fhClustersDThetaAcc->SetDirectory(0);
130   fhClustersDZetaAcc = new TH1F("dzetaacc","dzeta",100,-1.,1.);
131   fhClustersDZetaAcc->SetDirectory(0);
132
133   fhDPhiVsDZetaAcc = new TH2F("dphiVsDzetaacc","",100,-1.,1.,100,0.,0.1);
134   fhDPhiVsDZetaAcc->SetDirectory(0);
135   fhDPhiVsDThetaAcc = new TH2F("dphiVsDthetaAcc","",100,-0.1,0.1,100,0.,0.1);
136   fhDPhiVsDThetaAcc->SetDirectory(0);
137
138   fhClustersDPhiAll   = new TH1F("dphiall",  "dphi",  100,0.0,0.5);
139   fhClustersDPhiAll->SetDirectory(0);
140   fhClustersDThetaAll = new TH1F("dthetaall","dtheta",100,-0.5,0.5);
141   fhClustersDThetaAll->SetDirectory(0);
142   fhClustersDZetaAll = new TH1F("dzetaall","dzeta",100,-5.,5.);
143   fhClustersDZetaAll->SetDirectory(0);
144
145   fhDPhiVsDZetaAll = new TH2F("dphiVsDzetaall","",100,-5.,5.,100,0.,0.5);
146   fhDPhiVsDZetaAll->SetDirectory(0);
147   fhDPhiVsDThetaAll = new TH2F("dphiVsDthetaAll","",100,-0.5,0.5,100,0.,0.5);
148   fhDPhiVsDThetaAll->SetDirectory(0);
149
150   fhetaTracklets  = new TH1F("etaTracklets",  "eta",  100,-2.,2.);
151   fhetaTracklets->SetDirectory(0);
152   fhphiTracklets  = new TH1F("phiTracklets",  "phi",  100, 0., 2*TMath::Pi());
153   fhphiTracklets->SetDirectory(0);
154   fhetaClustersLay1  = new TH1F("etaClustersLay1",  "etaCl1",  100,-2.,2.);
155   fhetaClustersLay1->SetDirectory(0);
156   fhphiClustersLay1  = new TH1F("phiClustersLay1", "phiCl1", 100, 0., 2*TMath::Pi());
157   fhphiClustersLay1->SetDirectory(0);
158 }
159
160 //______________________________________________________________________
161 AliITSMultReconstructor::AliITSMultReconstructor(const AliITSMultReconstructor &mr) : TObject(mr),
162 fClustersLay1(mr.fClustersLay1),
163 fClustersLay2(mr.fClustersLay2),
164 fTracklets(mr.fTracklets),
165 fSClusters(mr.fSClusters),
166 fAssociationFlag(mr.fAssociationFlag),
167 fNClustersLay1(mr.fNClustersLay1),
168 fNClustersLay2(mr.fNClustersLay2),
169 fNTracklets(mr.fNTracklets),
170 fNSingleCluster(mr.fNSingleCluster),
171 fPhiWindow(mr.fPhiWindow),
172 fZetaWindow(mr.fZetaWindow),
173 fOnlyOneTrackletPerC2(mr.fOnlyOneTrackletPerC2),
174 fHistOn(mr.fHistOn),
175 fhClustersDPhiAcc(mr.fhClustersDPhiAcc),
176 fhClustersDThetaAcc(mr.fhClustersDThetaAcc),
177 fhClustersDZetaAcc(mr.fhClustersDZetaAcc),
178 fhClustersDPhiAll(mr.fhClustersDPhiAll),
179 fhClustersDThetaAll(mr.fhClustersDThetaAll),
180 fhClustersDZetaAll(mr.fhClustersDZetaAll),
181 fhDPhiVsDThetaAll(mr.fhDPhiVsDThetaAll),
182 fhDPhiVsDThetaAcc(mr.fhDPhiVsDThetaAcc),
183 fhDPhiVsDZetaAll(mr.fhDPhiVsDZetaAll),
184 fhDPhiVsDZetaAcc(mr.fhDPhiVsDZetaAcc),
185 fhetaTracklets(mr.fhetaTracklets),
186 fhphiTracklets(mr.fhphiTracklets),
187 fhetaClustersLay1(mr.fhetaClustersLay1),
188 fhphiClustersLay1(mr.fhphiClustersLay1) {
189   // Copy constructor
190
191 }
192
193 //______________________________________________________________________
194 AliITSMultReconstructor& AliITSMultReconstructor::operator=(const AliITSMultReconstructor& mr){
195   // Assignment operator
196   this->~AliITSMultReconstructor();
197   new(this) AliITSMultReconstructor(mr);
198   return *this;
199 }
200
201 //______________________________________________________________________
202 AliITSMultReconstructor::~AliITSMultReconstructor(){
203   // Destructor
204
205   // delete histograms
206   delete fhClustersDPhiAcc;
207   delete fhClustersDThetaAcc;
208   delete fhClustersDZetaAcc;
209   delete fhClustersDPhiAll;
210   delete fhClustersDThetaAll;
211   delete fhClustersDZetaAll;
212   delete fhDPhiVsDThetaAll;
213   delete fhDPhiVsDThetaAcc;
214   delete fhDPhiVsDZetaAll;
215   delete fhDPhiVsDZetaAcc;
216   delete fhetaTracklets;
217   delete fhphiTracklets;
218   delete fhetaClustersLay1;
219   delete fhphiClustersLay1;
220
221   // delete arrays
222   for(Int_t i=0; i<300000; i++) {
223     delete [] fClustersLay1[i];
224     delete [] fClustersLay2[i];
225     delete [] fTracklets[i];
226     delete [] fSClusters[i];
227   }
228   delete [] fClustersLay1;
229   delete [] fClustersLay2;
230   delete [] fTracklets;
231   delete [] fSClusters;
232
233   delete [] fAssociationFlag;
234 }
235
236 //____________________________________________________________________
237 void
238 AliITSMultReconstructor::Reconstruct(TTree* clusterTree, Float_t* vtx, Float_t* /* vtxRes*/) {
239   //
240   // - calls LoadClusterArray that finds the position of the clusters
241   //   (in global coord) 
242   // - convert the cluster coordinates to theta, phi (seen from the
243   //   interaction vertex). The third coordinate is used for ....
244   // - makes an array of tracklets 
245   //   
246   // After this method has been called, the clusters of the two layers
247   // and the tracklets can be retrieved by calling the Get'er methods.
248
249   // reset counters
250   fNClustersLay1 = 0;
251   fNClustersLay2 = 0;
252   fNTracklets = 0; 
253   fNSingleCluster = 0; 
254   // loading the clusters 
255   LoadClusterArrays(clusterTree);
256
257   // find the tracklets
258   AliDebug(1,"Looking for tracklets... ");  
259
260   //###########################################################
261   // Loop on layer 1 : finding theta, phi and z 
262   for (Int_t iC1=0; iC1<fNClustersLay1; iC1++) {    
263     Float_t x = fClustersLay1[iC1][0] - vtx[0];
264     Float_t y = fClustersLay1[iC1][1] - vtx[1];
265     Float_t z = fClustersLay1[iC1][2] - vtx[2];
266
267     Float_t r    = TMath::Sqrt(TMath::Power(x,2) +
268                                TMath::Power(y,2) +
269                                TMath::Power(z,2));
270     
271     fClustersLay1[iC1][0] = TMath::ACos(z/r);                   // Store Theta
272     fClustersLay1[iC1][1] = TMath::Pi() + TMath::ATan2(-y,-x);  // Store Phi
273     fClustersLay1[iC1][2] = z/r;                                // Store scaled z
274     if (fHistOn) {
275       Float_t eta=fClustersLay1[iC1][0];
276       eta= TMath::Tan(eta/2.);
277       eta=-TMath::Log(eta);
278       fhetaClustersLay1->Fill(eta);    
279       fhphiClustersLay1->Fill(fClustersLay1[iC1][1]);
280     }      
281   }
282   
283   // Loop on layer 2 : finding theta, phi and r   
284   for (Int_t iC2=0; iC2<fNClustersLay2; iC2++) {    
285     Float_t x = fClustersLay2[iC2][0] - vtx[0];
286     Float_t y = fClustersLay2[iC2][1] - vtx[1];
287     Float_t z = fClustersLay2[iC2][2] - vtx[2];
288    
289     Float_t r    = TMath::Sqrt(TMath::Power(x,2) +
290                                TMath::Power(y,2) +
291                                TMath::Power(z,2));
292     
293     fClustersLay2[iC2][0] = TMath::ACos(z/r);                   // Store Theta
294     fClustersLay2[iC2][1] = TMath::Pi() + TMath::ATan2(-y,-x);  // Store Phi
295     fClustersLay2[iC2][2] = z;                                  // Store z
296
297  // this only needs to be initialized for the fNClustersLay2 first associations
298     fAssociationFlag[iC2] = kFALSE;
299   }  
300   
301   //###########################################################
302   // Loop on layer 1 
303   for (Int_t iC1=0; iC1<fNClustersLay1; iC1++) {    
304
305     // reset of variables for multiple candidates
306     Int_t  iC2WithBestDist = 0;     // reset 
307     Float_t distmin        = 100.;  // just to put a huge number! 
308     Float_t dPhimin        = 0.;  // Used for histograms only! 
309     Float_t dThetamin      = 0.;  // Used for histograms only! 
310     Float_t dZetamin       = 0.;  // Used for histograms only! 
311     
312     // Loop on layer 2 
313     for (Int_t iC2=0; iC2<fNClustersLay2; iC2++) {      
314       
315       // The following excludes double associations
316       if (!fAssociationFlag[iC2]) {
317         
318         // find the difference in angles
319         Float_t dTheta = fClustersLay2[iC2][0] - fClustersLay1[iC1][0];
320         Float_t dPhi   = TMath::Abs(fClustersLay2[iC2][1] - fClustersLay1[iC1][1]);
321         // take into account boundary condition
322         if (dPhi>TMath::Pi()) dPhi=2.*TMath::Pi()-dPhi; 
323
324         // find the difference in z (between linear projection from layer 1
325         // and the actual point: Dzeta= z1/r1*r2 -z2)   
326         Float_t r2   = fClustersLay2[iC2][2]/TMath::Cos(fClustersLay2[iC2][0]);
327         Float_t dZeta  = fClustersLay1[iC1][2]*r2 - fClustersLay2[iC2][2];
328
329         if (fHistOn) {
330           fhClustersDPhiAll->Fill(dPhi);    
331           fhClustersDThetaAll->Fill(dTheta);    
332           fhClustersDZetaAll->Fill(dZeta);    
333           fhDPhiVsDThetaAll->Fill(dTheta, dPhi);
334           fhDPhiVsDZetaAll->Fill(dZeta, dPhi);
335         }
336         // make "elliptical" cut in Phi and Zeta! 
337         Float_t d = TMath::Sqrt(TMath::Power(dPhi/fPhiWindow,2) + TMath::Power(dZeta/fZetaWindow,2));
338
339         if (d>1) continue;      
340         
341         //look for the minimum distance: the minimum is in iC2WithBestDist
342         if (TMath::Sqrt(dZeta*dZeta+(r2*dPhi*r2*dPhi)) < distmin ) {
343           distmin=TMath::Sqrt(dZeta*dZeta + (r2*dPhi*r2*dPhi));
344           dPhimin = dPhi;
345           dThetamin = dTheta;
346           dZetamin = dZeta; 
347           iC2WithBestDist = iC2;
348         }
349       } 
350     } // end of loop over clusters in layer 2 
351     
352     if (distmin<100) { // This means that a cluster in layer 2 was found that mathes with iC1
353
354       if (fHistOn) {
355         fhClustersDPhiAcc->Fill(dPhimin);
356         fhClustersDThetaAcc->Fill(dThetamin);    
357         fhClustersDZetaAcc->Fill(dZetamin);    
358         fhDPhiVsDThetaAcc->Fill(dThetamin, dPhimin);
359         fhDPhiVsDZetaAcc->Fill(dZetamin, dPhimin);
360       }
361       
362       if (fOnlyOneTrackletPerC2) fAssociationFlag[iC2WithBestDist] = kTRUE; // flag the association
363       
364       // store the tracklet
365       
366       // use the theta from the clusters in the first layer
367       fTracklets[fNTracklets][0] = fClustersLay1[iC1][0];
368       // use the phi from the clusters in the first layer
369       fTracklets[fNTracklets][1] = fClustersLay1[iC1][1];
370       // store the difference between phi1 and phi2
371       fTracklets[fNTracklets][2] = fClustersLay1[iC1][1] - fClustersLay2[iC2WithBestDist][1];
372
373       // define dphi in the range [0,pi] with proper sign (track charge correlated)
374       if (fTracklets[fNTracklets][2] > TMath::Pi()) 
375           fTracklets[fNTracklets][2] = fTracklets[fNTracklets][2]-2.*TMath::Pi();
376       if (fTracklets[fNTracklets][2] < -TMath::Pi()) 
377           fTracklets[fNTracklets][2] = fTracklets[fNTracklets][2]+2.*TMath::Pi();
378
379       // find label
380       // if equal label in both clusters found this label is assigned
381       // if no equal label can be found the first labels of the L1 AND L2 cluster are assigned
382       Int_t label1 = 0;
383       Int_t label2 = 0;
384       while (label2 < 3)
385       {
386         if ((Int_t) fClustersLay1[iC1][3+label1] != -2 && (Int_t) fClustersLay1[iC1][3+label1] == (Int_t) fClustersLay2[iC2WithBestDist][3+label2])
387           break;
388         label1++;
389         if (label1 == 3)
390         {
391           label1 = 0;
392           label2++;
393         }
394       }
395
396       if (label2 < 3)
397       {
398         AliDebug(AliLog::kDebug, Form("Found label %d == %d for tracklet candidate %d\n", (Int_t) fClustersLay1[iC1][3+label1], (Int_t) fClustersLay2[iC2WithBestDist][3+label2], fNTracklets));
399         fTracklets[fNTracklets][3] = fClustersLay1[iC1][3+label1];
400         fTracklets[fNTracklets][4] = fClustersLay2[iC2WithBestDist][3+label2];
401       }
402       else
403       {
404         AliDebug(AliLog::kDebug, Form("Did not find label %d %d %d %d %d %d for tracklet candidate %d\n", (Int_t) fClustersLay1[iC1][3], (Int_t) fClustersLay1[iC1][4], (Int_t) fClustersLay1[iC1][5], (Int_t) fClustersLay2[iC2WithBestDist][3], (Int_t) fClustersLay2[iC2WithBestDist][4], (Int_t) fClustersLay2[iC2WithBestDist][5], fNTracklets));
405         fTracklets[fNTracklets][3] = fClustersLay1[iC1][3];
406         fTracklets[fNTracklets][4] = fClustersLay2[iC2WithBestDist][3];
407       }
408
409       if (fHistOn) {
410         Float_t eta=fTracklets[fNTracklets][0];
411         eta= TMath::Tan(eta/2.);
412         eta=-TMath::Log(eta);
413         fhetaTracklets->Fill(eta);    
414         fhphiTracklets->Fill(fTracklets[fNTracklets][1]);    
415       }
416       
417       AliDebug(1,Form(" Adding tracklet candidate %d ", fNTracklets));
418       AliDebug(1,Form(" Cl. %d of Layer 1 and %d of Layer 2", iC1, 
419                       iC2WithBestDist));
420       fNTracklets++;
421     }
422
423     // Delete the following else if you do not want to save Clusters! 
424
425     else { // This means that the cluster has not been associated
426
427       // store the cluster
428        
429       fSClusters[fNSingleCluster][0] = fClustersLay1[iC1][0];
430       fSClusters[fNSingleCluster][1] = fClustersLay1[iC1][1];
431       AliDebug(1,Form(" Adding a single cluster %d (cluster %d  of layer 1)",
432                       fNSingleCluster, iC1));
433       fNSingleCluster++;
434     }
435
436   } // end of loop over clusters in layer 1
437   
438   AliDebug(1,Form("%d tracklets found", fNTracklets));
439 }
440
441 //____________________________________________________________________
442 void
443 AliITSMultReconstructor::LoadClusterArrays(TTree* itsClusterTree) {
444   // This method
445   // - gets the clusters from the cluster tree 
446   // - convert them into global coordinates 
447   // - store them in the internal arrays
448   // - count the number of cluster-fired chips
449   
450   AliDebug(1,"Loading clusters and cluster-fired chips ...");
451   
452   fNClustersLay1 = 0;
453   fNClustersLay2 = 0;
454   fNFiredChips[0] = 0;
455   fNFiredChips[1] = 0;
456   
457   AliITSsegmentationSPD *seg = new AliITSsegmentationSPD();
458
459   TClonesArray* itsClusters = new TClonesArray("AliITSRecPoint");
460   TBranch* itsClusterBranch=itsClusterTree->GetBranch("ITSRecPoints");
461
462   itsClusterBranch->SetAddress(&itsClusters);
463
464   Int_t nItsSubs = (Int_t)itsClusterTree->GetEntries();  
465   Float_t cluGlo[3]={0.,0.,0.};
466  
467   // loop over the its subdetectors
468   for (Int_t iIts=0; iIts < nItsSubs; iIts++) {
469     
470     if (!itsClusterTree->GetEvent(iIts)) 
471       continue;
472     
473     Int_t nClusters = itsClusters->GetEntriesFast();
474
475     // number of clusters in each chip of the current module
476     Int_t nClustersInChip[5] = {0,0,0,0,0};
477     Int_t layer = 0;
478     
479     // loop over clusters
480     while(nClusters--) {
481       AliITSRecPoint* cluster = (AliITSRecPoint*)itsClusters->UncheckedAt(nClusters);
482       
483       layer = cluster->GetLayer();
484       if (layer>1) continue;            
485       
486       cluster->GetGlobalXYZ(cluGlo);
487       Float_t x = cluGlo[0];
488       Float_t y = cluGlo[1];
489       Float_t z = cluGlo[2];      
490
491       // find the chip for the current cluster
492       Float_t locz = cluster->GetDetLocalZ();
493       Int_t iChip = seg->GetChipFromLocal(0,locz);
494       nClustersInChip[iChip]++; 
495       
496       if (layer==0) {
497         fClustersLay1[fNClustersLay1][0] = x;
498         fClustersLay1[fNClustersLay1][1] = y;
499         fClustersLay1[fNClustersLay1][2] = z;
500         for (Int_t i=0; i<3; i++)
501                 fClustersLay1[fNClustersLay1][3+i] = cluster->GetLabel(i);
502         fNClustersLay1++;
503       }
504       if (layer==1) {
505         fClustersLay2[fNClustersLay2][0] = x;
506         fClustersLay2[fNClustersLay2][1] = y;
507         fClustersLay2[fNClustersLay2][2] = z;
508         for (Int_t i=0; i<3; i++)
509                 fClustersLay2[fNClustersLay2][3+i] = cluster->GetLabel(i);
510         fNClustersLay2++;
511       }
512       
513     }// end of cluster loop
514
515     // get number of fired chips in the current module
516     if(layer<2)
517     for(Int_t ifChip=0; ifChip<5; ifChip++) {
518       if(nClustersInChip[ifChip] >= 1)  fNFiredChips[layer]++;
519     }
520
521   } // end of its "subdetector" loop  
522
523   if (itsClusters) {
524     itsClusters->Delete();
525     delete itsClusters;
526     delete seg;
527     itsClusters = 0;
528   }
529   AliDebug(1,Form("(clusters in layer 1 : %d,  layer 2: %d)",fNClustersLay1,fNClustersLay2));
530   AliDebug(1,Form("(cluster-fired chips in layer 1 : %d,  layer 2: %d)",fNFiredChips[0],fNFiredChips[1]));
531 }
532 //____________________________________________________________________
533 void
534 AliITSMultReconstructor::LoadClusterFiredChips(TTree* itsClusterTree) {
535   // This method
536   // - gets the clusters from the cluster tree 
537   // - counts the number of (cluster)fired chips
538   
539   AliDebug(1,"Loading cluster-fired chips ...");
540   
541   fNFiredChips[0] = 0;
542   fNFiredChips[1] = 0;
543   
544   AliITSsegmentationSPD *seg = new AliITSsegmentationSPD();
545
546   TClonesArray* itsClusters = new TClonesArray("AliITSRecPoint");
547   TBranch* itsClusterBranch=itsClusterTree->GetBranch("ITSRecPoints");
548
549   itsClusterBranch->SetAddress(&itsClusters);
550
551   Int_t nItsSubs = (Int_t)itsClusterTree->GetEntries();  
552  
553   // loop over the its subdetectors
554   for (Int_t iIts=0; iIts < nItsSubs; iIts++) {
555     
556     if (!itsClusterTree->GetEvent(iIts)) 
557       continue;
558     
559     Int_t nClusters = itsClusters->GetEntriesFast();
560
561     // number of clusters in each chip of the current module
562     Int_t nClustersInChip[5] = {0,0,0,0,0};
563     Int_t layer = 0;
564     
565     // loop over clusters
566     while(nClusters--) {
567       AliITSRecPoint* cluster = (AliITSRecPoint*)itsClusters->UncheckedAt(nClusters);
568       
569       layer = cluster->GetLayer();
570       if (layer>1) continue;            
571
572       // find the chip for the current cluster
573       Float_t locz = cluster->GetDetLocalZ();
574       Int_t iChip = seg->GetChipFromLocal(0,locz);
575       nClustersInChip[iChip]++; 
576       
577     }// end of cluster loop
578
579     // get number of fired chips in the current module
580     if(layer<2)
581     for(Int_t ifChip=0; ifChip<5; ifChip++) {
582       if(nClustersInChip[ifChip] >= 1)  fNFiredChips[layer]++;
583     }
584
585   } // end of its "subdetector" loop  
586   
587   if (itsClusters) {
588     itsClusters->Delete();
589     delete itsClusters;
590     delete seg;
591     itsClusters = 0;
592   }
593   AliDebug(1,Form("(cluster-fired chips in layer 1 : %d,  layer 2: %d)",fNFiredChips[0],fNFiredChips[1]));
594 }
595 //____________________________________________________________________
596 void
597 AliITSMultReconstructor::SaveHists() {
598   // This method save the histograms on the output file
599   // (only if fHistOn is TRUE). 
600   
601   if (!fHistOn)
602     return;
603
604   fhClustersDPhiAll->Write();
605   fhClustersDThetaAll->Write();
606   fhClustersDZetaAll->Write();
607   fhDPhiVsDThetaAll->Write();
608   fhDPhiVsDZetaAll->Write();
609
610   fhClustersDPhiAcc->Write();
611   fhClustersDThetaAcc->Write();
612   fhClustersDZetaAcc->Write();
613   fhDPhiVsDThetaAcc->Write();
614   fhDPhiVsDZetaAcc->Write();
615
616   fhetaTracklets->Write();
617   fhphiTracklets->Write();
618   fhetaClustersLay1->Write();
619   fhphiClustersLay1->Write();
620 }