Adding AliReconstruction class (T.Kuhr)
[u/mrichter/AliRoot.git] / ITS / AliITS.cxx
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 /* $Id$ */
17
18
19 ///////////////////////////////////////////////////////////////////////////////
20 //
21 //      An overview of the basic philosophy of the ITS code development
22 // and analysis is show in the figure below.
23 //Begin_Html
24 /*
25 <img src="picts/ITS/ITS_Analysis_schema.gif">
26 </pre>
27 <br clear=left>
28 <font size=+2 color=red>
29 <p>Roberto Barbera is in charge of the ITS Offline code (1999).
30 <a href="mailto:roberto.barbera@ct.infn.it">Roberto Barbera</a>.
31 </font>
32 <pre>
33 */
34 //End_Html
35 //
36 //  AliITS. Inner Traking System base class.
37 //  This class contains the base procedures for the Inner Tracking System
38 //
39 //Begin_Html
40 /*
41 <img src="picts/ITS/AliITS_Class_Diagram.gif">
42 </pre>
43 <br clear=left>
44 <font size=+2 color=red>
45 <p>This show the class diagram of the different elements that are part of
46 the AliITS class.
47 </font>
48 <pre>
49 */
50 //End_Html
51 //
52 // Version: 0
53 // Written by Rene Brun, Federico Carminati, and Roberto Barbera
54 //
55 // Version: 1
56 // Modified and documented by Bjorn S. Nilsen
57 // July 11 1999
58 //
59 // Version: 2
60 // Modified and documented by A. Bologna
61 // October 18 1999
62 //
63 // AliITS is the general base class for the ITS. Also see AliDetector for
64 // futher information.
65 //
66 ///////////////////////////////////////////////////////////////////////////////
67
68 #include <Riostream.h>
69 #include <stdlib.h>
70
71 #include <TBranch.h>
72 #include <TClonesArray.h>
73 #include <TFile.h>
74 #include <TMath.h>
75 #include <TObjectTable.h>
76 #include <TROOT.h>
77 #include <TRandom.h>
78 #include <TString.h>
79 #include <TTree.h>
80 #include <TVector.h>
81 #include <TVirtualMC.h>
82
83 #include "AliConfig.h"
84 #include "AliHeader.h"
85 #include "AliITS.h"
86 #include "AliITSClusterFinderSDD.h"
87 #include "AliITSClusterFinderSPD.h"
88 #include "AliITSClusterFinderSSD.h"
89 #include "AliITSDetType.h"
90 #include "AliITSLoader.h"
91 #include "AliITSRawClusterSPD.h"
92 #include "AliITSRawClusterSDD.h"
93 #include "AliITSRawClusterSSD.h"
94 #include "AliITSRecPoint.h"
95 #include "AliITSdigit.h"
96 #include "AliITSgeom.h"
97 #include "AliITShit.h"
98 #include "AliITSmodule.h"
99 #include "AliITSpList.h"
100 #include "AliITSresponseSDD.h"
101 #include "AliITSresponseSPD.h"
102 #include "AliITSresponseSSD.h"
103 #include "AliITSsegmentationSDD.h"
104 #include "AliITSsegmentationSPD.h"
105 #include "AliITSsegmentationSSD.h"
106 #include "AliITSsimulationSDD.h"
107 #include "AliITSsimulationSPD.h"
108 #include "AliITSsimulationSSD.h"
109 #include "AliMC.h"
110 #include "AliITSDigitizer.h"
111 #include "AliITSclustererV2.h"
112 #include "AliITStrackerV2.h"
113 #include "AliITSpidESD.h"
114 #include "AliV0vertexer.h"
115 #include "AliCascadeVertexer.h"
116 #include "AliESD.h"
117
118 ClassImp(AliITS)
119
120 //______________________________________________________________________
121 AliITS::AliITS() : AliDetector() {
122     // Default initializer for ITS
123     //      The default constructor of the AliITS class. In addition to
124     // creating the AliITS class it zeros the variables fIshunt (a member
125     // of AliDetector class), fEuclidOut, and fIdN, and zeros the pointers
126     // fITSpoints, fIdSens, and fIdName. The AliDetector default constructor
127     // is also called.
128     // Inputs:
129     //      none.
130     // Outputs:
131     //      none.
132     // Return:
133     //      Blank ITS class.
134
135     fIshunt     = 0;   // not zeroed in AliDetector.
136
137     // AliITS variables.
138     fEuclidOut  = 0;
139     fITSgeom    = 0;
140     fITSmodules = 0;
141     fOpt        = "All";
142 //    SetDetectors(); // default to fOpt="All". This variable not written out.
143
144     fIdN        = 0;
145     fIdName     = 0;
146     fIdSens     = 0;
147
148     fNDetTypes  = kNTYPES;
149     fDetTypes   = 0;
150
151     fSDigits    = 0;
152     fNSDigits   = 0;
153
154     fNdtype     = 0;
155     fDtype      = 0;
156
157     fCtype      = 0;
158     fNctype     = 0;
159
160     fRecPoints  = 0;
161     fNRecPoints = 0;
162
163     SetMarkerColor(kRed);
164 }
165 //______________________________________________________________________
166 AliITS::AliITS(const char *name, const char *title):AliDetector(name,title){
167     //     The standard Constructor for the ITS class. In addition to 
168     // creating the AliITS class, it allocates memory for the TClonesArrays 
169     // fHits, fSDigits, fDigits, fITSpoints, and the TObjArray of fCtype 
170     // (clusters). It also zeros the variables
171     // fIshunt (a member of AliDetector class), fEuclidOut, and fIdN, and zeros
172     // the pointers fIdSens and fIdName. To help in displaying hits via the
173     // ROOT macro display.C AliITS also sets the marker color to red. The
174     // variables passes with this constructor, const char *name and *title,
175     // are used by the constructor of AliDetector class. See AliDetector
176     // class for a description of these parameters and its constructor
177     // functions.
178     // Inputs:
179     //      const char *name      Detector name. Should always be "ITS"
180     //      const char *title     Detector title.
181     // Outputs:
182     //      none.
183     // Return:
184     //      An ITS class.
185
186     fIshunt     = 0;  // not zeroed in AliDetector
187     fHits       = new TClonesArray("AliITShit", 1560);//not done in AliDetector
188     gAlice->GetMCApp()->AddHitList(fHits);  // Not done in AliDetector.
189
190     fEuclidOut  = 0;
191     fITSgeom    = 0;
192     fITSmodules = 0;
193     SetDetectors(); // default to fOpt="All". This variable not written out.
194
195     fIdN        = 0;
196     fIdName     = 0;
197     fIdSens     = 0;
198
199     fNDetTypes  = kNTYPES;
200     fDetTypes   = new TObjArray(fNDetTypes);
201
202     fSDigits    = new TClonesArray("AliITSpListItem",1000);
203     fNSDigits   = 0;
204
205     fNdtype     = new Int_t[fNDetTypes];
206     fDtype      = new TObjArray(fNDetTypes);
207
208     fCtype      = new TObjArray(fNDetTypes);
209     fNctype     = new Int_t[fNDetTypes];
210
211     fRecPoints  = new TClonesArray("AliITSRecPoint",1000);
212     fNRecPoints = 0;
213
214     Int_t i;
215     for(i=0;i<fNDetTypes;i++) {
216         fDetTypes->AddAt(new AliITSDetType(),i); 
217         fNdtype[i] = 0;
218         fNctype[i] = 0;
219     } // end for i
220
221     SetMarkerColor(kRed);
222
223 }
224 //______________________________________________________________________
225 AliITS::~AliITS(){
226     // Default destructor for ITS.
227     //     The default destructor of the AliITS class. In addition to deleting
228     // the AliITS class it deletes the memory pointed to by the fHits, fDigits,
229     // fSDigits, fCtype, fITSmodules, fITSgeom, fRecPoints, fIdSens, fIdName, 
230     // fITSpoints, fDetType and it's contents.
231     // Inputs:
232     //      none.
233     // Outputs:
234     //      none.
235     // Return:
236     //      none.
237
238     if (fHits) {
239       fHits->Delete();
240       delete fHits;
241       fHits=0;
242     }
243     if (fSDigits) {
244       fSDigits->Delete();
245       delete fSDigits;
246       fSDigits=0;
247     }
248     if (fDigits) {
249       fDigits->Delete();
250       delete fDigits;
251       fDigits=0;
252     }
253     if (fRecPoints) {
254       fRecPoints->Delete();
255       delete fRecPoints;
256       fRecPoints=0;
257     }
258     delete[] fIdName;  // Array of TStrings
259     delete[] fIdSens;
260     if(fITSmodules) {
261         this->ClearModules();
262         delete fITSmodules;
263         fITSmodules = 0;
264     }// end if fITSmodules!=0
265
266     if(fDtype) {
267         fDtype->Delete();
268         delete fDtype;
269     } // end if fDtype
270     delete [] fNdtype;
271     if (fCtype) {
272         fCtype->Delete();
273         delete fCtype;
274         fCtype = 0;
275     } // end if fCtype
276     delete [] fNctype;
277
278     if (fDetTypes) {
279         fDetTypes->Delete();
280         delete fDetTypes;
281         fDetTypes = 0;
282     } // end if fDetTypes
283
284
285     if (fITSgeom) delete fITSgeom;
286 }
287 //______________________________________________________________________
288 AliITS::AliITS(const AliITS &source) : AliDetector(source){
289     // Copy constructor. This is a function which is not allowed to be
290     // done to the ITS. It exits with an error.
291     // Inputs:
292     //      AliITS &source  An AliITS class.
293     // Outputs:
294     //      none.
295     // Return:
296     //      none.
297
298     if(this==&source) return;
299     Error("Copy constructor",
300           "You are not allowed to make a copy of the AliITS");
301     exit(1);
302 }
303 //______________________________________________________________________
304 AliITS& AliITS::operator=(AliITS &source){
305     // Assignment operator. This is a function which is not allowed to be
306     // done to the ITS. It exits with an error.
307     // Inputs:
308     //      AliITS &source  An AliITS class.
309     // Outputs:
310     //      none.
311     // Return:
312     //      none.
313
314     if(this==&source) return *this;
315     Error("operator=","You are not allowed to make a copy of the AliITS");
316     exit(1);
317     return *this; //fake return
318 }
319 //______________________________________________________________________
320 Int_t AliITS::DistancetoPrimitive(Int_t,Int_t) const{
321     // Distance from mouse to ITS on the screen. Dummy routine
322     //     A dummy routine used by the ROOT macro display.C to allow for the
323     // use of the mouse (pointing device) in the macro. In general this should
324     // never be called. If it is it returns the number 9999 for any value of
325     // x and y.
326     // Inputs:
327     //      Int_t     Dummy screen coordinate.
328     //      Int_t     Dummy screen coordinate.
329     // Outputs:
330     //      none.
331     // Return:
332     //      Int_t     Dummy = 9999 distance to ITS.
333
334     return 9999;
335 }
336 //______________________________________________________________________
337 void AliITS::Init(){
338     // Initializer ITS after it has been built
339     //     This routine initializes the AliITS class. It is intended to be
340     // called from the Init function in AliITSv?. Besides displaying a banner
341     // indicating that it has been called it initializes the array fIdSens
342     // and sets the default segmentation, response, digit and raw cluster
343     // classes therefore it should be called after a call to CreateGeometry.
344     // Inputs:
345     //      none.
346     // Outputs:
347     //      none.
348     // Return:
349     //      none.
350     Int_t i;
351
352     SetDefaults();
353     // Array of TStrings
354     for(i=0;i<fIdN;i++) fIdSens[i] = gMC->VolId(fIdName[i]);
355 }
356 //______________________________________________________________________
357 void AliITS::SetDefaults(){
358     // sets the default segmentation, response, digit and raw cluster classes.
359     // Inputs:
360     //      none.
361     // Outputs:
362     //      none.
363     // Return:
364     //      none.
365
366     if(fDebug) printf("%s: SetDefaults\n",ClassName());
367
368     AliITSDetType *iDetType;
369
370     //SPD
371     iDetType=DetType(0); 
372     if (!iDetType->GetSegmentationModel()) {
373         AliITSsegmentationSPD *seg0=new AliITSsegmentationSPD(fITSgeom);
374         SetSegmentationModel(0,seg0); 
375     } // end if
376     if (!iDetType->GetResponseModel()) {
377         SetResponseModel(0,new AliITSresponseSPD()); 
378     } // end if
379     // set digit and raw cluster classes to be used
380
381     const char *kData0=(iDetType->GetResponseModel())->DataType();
382     if (strstr(kData0,"real")) {
383         iDetType->ClassNames("AliITSdigit","AliITSRawClusterSPD");
384     } else iDetType->ClassNames("AliITSdigitSPD","AliITSRawClusterSPD");
385
386     // SDD
387     iDetType=DetType(1); 
388     if (!iDetType->GetResponseModel()) {
389         SetResponseModel(1,new AliITSresponseSDD("simulated")); 
390     } // end if
391     AliITSresponse *resp1=iDetType->GetResponseModel();
392     if (!iDetType->GetSegmentationModel()) {
393         AliITSsegmentationSDD *seg1=new AliITSsegmentationSDD(fITSgeom,resp1);
394         SetSegmentationModel(1,seg1); 
395     } // end if
396     const char *kData1=(iDetType->GetResponseModel())->DataType();
397     const char *kopt=iDetType->GetResponseModel()->ZeroSuppOption();
398     if((!strstr(kopt,"2D"))&&(!strstr(kopt,"1D")) || strstr(kData1,"real") ){
399         iDetType->ClassNames("AliITSdigit","AliITSRawClusterSDD");
400     } else iDetType->ClassNames("AliITSdigitSDD","AliITSRawClusterSDD");
401
402     // SSD
403     iDetType=DetType(2); 
404     if (!iDetType->GetSegmentationModel()) {
405         AliITSsegmentationSSD *seg2=new AliITSsegmentationSSD(fITSgeom);
406         SetSegmentationModel(2,seg2); 
407     } // end if
408     if (!iDetType->GetResponseModel()) {
409         SetResponseModel(2,new AliITSresponseSSD("simulated"));
410     } // end if
411     const char *kData2=(iDetType->GetResponseModel())->DataType();
412     if (strstr(kData2,"real")) {
413         iDetType->ClassNames("AliITSdigit","AliITSRawClusterSSD");
414     } else iDetType->ClassNames("AliITSdigitSSD","AliITSRawClusterSSD");
415
416     if (kNTYPES>3) {
417         Warning("SetDefaults",
418                 "Only the three basic detector types are initialized!");
419     }  // end if
420 }
421 //______________________________________________________________________
422 void AliITS::SetDefaultSimulation(){
423     // sets the default simulation.
424     // Inputs:
425     //      none.
426     // Outputs:
427     //      none.
428     // Return:
429     //      none.
430
431     AliITSDetType *iDetType;
432     AliITSsimulation *sim;
433     iDetType=DetType(0);
434     sim = iDetType->GetSimulationModel();
435     if (!sim) {
436         AliITSsegmentation *seg0=
437             (AliITSsegmentation*)iDetType->GetSegmentationModel();
438         AliITSresponse *res0 = (AliITSresponse*)iDetType->GetResponseModel();
439         AliITSsimulationSPD *sim0=new AliITSsimulationSPD(seg0,res0);
440         SetSimulationModel(0,sim0);
441     }else{ // simulation exists, make sure it is set up properly.
442         ((AliITSsimulationSPD*)sim)->Init(
443             (AliITSsegmentationSPD*) iDetType->GetSegmentationModel(),
444             (AliITSresponseSPD*) iDetType->GetResponseModel());
445 //        if(sim->GetResponseModel()==0) sim->SetResponseModel(
446 //            (AliITSresponse*)iDetType->GetResponseModel());
447 //        if(sim->GetSegmentationModel()==0) sim->SetSegmentationModel(
448 //            (AliITSsegmentation*)iDetType->GetSegmentationModel());
449     } // end if
450     iDetType=DetType(1);
451     sim = iDetType->GetSimulationModel();
452     if (!sim) {
453         AliITSsegmentation *seg1=
454             (AliITSsegmentation*)iDetType->GetSegmentationModel();
455         AliITSresponse *res1 = (AliITSresponse*)iDetType->GetResponseModel();
456         AliITSsimulationSDD *sim1=new AliITSsimulationSDD(seg1,res1);
457         SetSimulationModel(1,sim1);
458     }else{ // simulation exists, make sure it is set up properly.
459         ((AliITSsimulationSDD*)sim)->Init(
460             (AliITSsegmentationSDD*) iDetType->GetSegmentationModel(),
461             (AliITSresponseSDD*) iDetType->GetResponseModel());
462 //        if(sim->GetResponseModel()==0) sim->SetResponseModel(
463 //            (AliITSresponse*)iDetType->GetResponseModel());
464 //        if(sim->GetSegmentationModel()==0) sim->SetSegmentationModel(
465 //            (AliITSsegmentation*)iDetType->GetSegmentationModel());
466     } //end if
467     iDetType=DetType(2);
468     sim = iDetType->GetSimulationModel();
469     if (!sim) {
470         AliITSsegmentation *seg2=
471             (AliITSsegmentation*)iDetType->GetSegmentationModel();
472         AliITSresponse *res2 = (AliITSresponse*)iDetType->GetResponseModel();
473         AliITSsimulationSSD *sim2=new AliITSsimulationSSD(seg2,res2);
474         SetSimulationModel(2,sim2);
475     }else{ // simulation exists, make sure it is set up properly.
476         ((AliITSsimulationSSD*)sim)->Init(
477             (AliITSsegmentationSSD*) iDetType->GetSegmentationModel(),
478             (AliITSresponseSSD*) iDetType->GetResponseModel());
479 //        if(sim->GetResponseModel()==0) sim->SetResponseModel(
480 //            (AliITSresponse*)iDetType->GetResponseModel());
481 //        if(sim->GetSegmentationModel()==0) sim->SetSegmentationModel(
482 //            (AliITSsegmentation*)iDetType->GetSegmentationModel());
483     } // end if
484 }
485 //______________________________________________________________________
486 void AliITS::SetDefaultClusterFinders(){
487     // Sets the default cluster finders. Used in finding RecPoints.
488     // Inputs:
489     //      none.
490     // Outputs:
491     //      none.
492     // Return:
493     //      none.
494
495     MakeTreeC();
496     AliITSDetType *iDetType;
497
498     // SPD
499     iDetType=DetType(0);
500     if (!iDetType->GetReconstructionModel()) {
501         AliITSsegmentation *seg0 =
502             (AliITSsegmentation*)iDetType->GetSegmentationModel();
503         TClonesArray *dig0=DigitsAddress(0);
504         TClonesArray *recp0=ClustersAddress(0);
505         AliITSClusterFinderSPD *rec0 = new AliITSClusterFinderSPD(seg0,dig0,
506                                                                   recp0);
507         SetReconstructionModel(0,rec0);
508     } // end if
509
510     // SDD
511     iDetType=DetType(1);
512     if (!iDetType->GetReconstructionModel()) {
513         AliITSsegmentation *seg1 =
514             (AliITSsegmentation*)iDetType->GetSegmentationModel();
515         AliITSresponse *res1 = (AliITSresponse*)iDetType->GetResponseModel();
516         TClonesArray *dig1=DigitsAddress(1);
517         TClonesArray *recp1=ClustersAddress(1);
518         AliITSClusterFinderSDD *rec1 =
519             new AliITSClusterFinderSDD(seg1,res1,dig1,recp1);
520       SetReconstructionModel(1,rec1);
521     } // end if
522
523     // SSD
524     iDetType=DetType(2);
525     if (!iDetType->GetReconstructionModel()) {
526         AliITSsegmentation *seg2=
527             (AliITSsegmentation*)iDetType->GetSegmentationModel();
528         TClonesArray *dig2=DigitsAddress(2);
529         AliITSClusterFinderSSD *rec2= new AliITSClusterFinderSSD(seg2,dig2);
530         SetReconstructionModel(2,rec2);
531     } // end if
532 }
533 //______________________________________________________________________
534 void AliITS::MakeBranch(Option_t* option){
535     // Creates Tree branches for the ITS.
536     // Inputs:
537     //      Option_t *option    String of Tree types S,D, and/or R.
538     //      const char *file    String of the file name where these branches
539     //                          are to be stored. If blank then these branches
540     //                          are written to the same tree as the Hits were
541     //                          read from.
542     // Outputs:
543     //      none.
544     // Return:
545     //      none.
546     Bool_t cH = (strstr(option,"H")!=0);
547     Bool_t cS = (strstr(option,"S")!=0);
548     Bool_t cD = (strstr(option,"D")!=0);
549     Bool_t cR = (strstr(option,"R")!=0);
550     Bool_t cRF = (strstr(option,"RF")!=0);
551     
552     if(cRF)cR = kFALSE;
553     if(cH && (fHits == 0x0)) fHits  = new TClonesArray("AliITShit", 1560);
554     
555     AliDetector::MakeBranch(option);
556
557     if(cS) MakeBranchS(0);
558     if(cD) MakeBranchD(0);
559     if(cR) MakeBranchR(0);
560     if(cRF) MakeBranchRF(0);
561 }
562 //______________________________________________________________________
563 void AliITS::SetTreeAddress(){
564     // Set branch address for the Trees.
565     // Inputs:
566     //      none.
567     // Outputs:
568     //      none.
569     // Return:
570     //      none.
571     TTree *treeS = fLoader->TreeS();
572     TTree *treeD = fLoader->TreeD();
573     TTree *treeR = fLoader->TreeR();
574     if (fLoader->TreeH() && (fHits == 0x0)) fHits = new TClonesArray("AliITShit", 1560);
575       
576     AliDetector::SetTreeAddress();
577
578     SetTreeAddressS(treeS);
579     SetTreeAddressD(treeD);
580     SetTreeAddressR(treeR);
581 }
582 //______________________________________________________________________
583 AliITSDetType* AliITS::DetType(Int_t id){
584     // Return pointer to id detector type.
585     // Inputs:
586     //      Int_t id   detector id number.
587     // Outputs:
588     //      none.
589     // Return:
590     //      returned, a pointer to a AliITSDetType.
591
592     return ((AliITSDetType*) fDetTypes->At(id));
593 }
594 //______________________________________________________________________
595 void AliITS::SetResponseModel(Int_t id, AliITSresponse *response){
596     // Set the response model for the id detector type.
597     // Inputs:
598     //      Int_t id        detector id number.
599     //      AliITSresponse* a pointer containing an instance of AliITSresponse
600     //                      to be stored/owned b y AliITSDetType.
601     // Outputs:
602     //      none.
603     // Return:
604     //      none.
605
606     ((AliITSDetType*) fDetTypes->At(id))->ResponseModel(response);
607 }
608 //______________________________________________________________________
609 void AliITS::SetSegmentationModel(Int_t id, AliITSsegmentation *seg){
610     // Set the segmentation model for the id detector type.
611     // Inputs:
612     //      Int_t id            detector id number.
613     //      AliITSsegmentation* a pointer containing an instance of 
614     //                          AliITSsegmentation to be stored/owned b y 
615     //                          AliITSDetType.
616     // Outputs:
617     //      none.
618     // Return:
619     //      none.
620
621     ((AliITSDetType*) fDetTypes->At(id))->SegmentationModel(seg);
622 }
623 //______________________________________________________________________
624 void AliITS::SetSimulationModel(Int_t id, AliITSsimulation *sim){
625     // Set the simulation model for the id detector type.
626     // Inputs:
627     //      Int_t id        detector id number.
628     //      AliITSresponse* a pointer containing an instance of AliITSresponse
629     //                      to be stored/owned b y AliITSDetType.
630     // Outputs:
631     //      none.
632     // Return:
633     //      none.
634
635    ((AliITSDetType*) fDetTypes->At(id))->SimulationModel(sim);
636
637 }
638 //______________________________________________________________________
639 void AliITS::SetReconstructionModel(Int_t id, AliITSClusterFinder *reconst){
640     // Set the cluster finder model for the id detector type.
641     // Inputs:
642     //      Int_t id             detector id number.
643     //      AliITSClusterFinder* a pointer containing an instance of 
644     //                           AliITSClusterFinder to be stored/owned b y 
645     //                           AliITSDetType.
646     // Outputs:
647     //      none.
648     // Return:
649     //      none.
650
651     ((AliITSDetType*) fDetTypes->At(id))->ReconstructionModel(reconst);
652 }
653 //______________________________________________________________________
654 void AliITS::SetClasses(Int_t id, const char *digit, const char *cluster){
655     // Set the digit and cluster classes name to be used for the id detector
656     // type.
657     // Inputs:
658     //      Int_t id            detector id number.
659     //      const char *digit   Digit class name for detector id.
660     //      const char *cluster Cluster class name for detector id.
661     // Outputs:
662     //      none.
663     // Return:
664     //      none.
665
666     ((AliITSDetType*) fDetTypes->At(id))->ClassNames(digit,cluster);
667 }
668 //______________________________________________________________________
669 void AliITS::AddHit(Int_t track, Int_t *vol, Float_t *hits){
670     // Add an ITS hit
671     //     The function to add information to the AliITShit class. See the
672     // AliITShit class for a full description. This function allocates the
673     // necessary new space for the hit information and passes the variable
674     // track, and the pointers *vol and *hits to the AliITShit constructor
675     // function.
676     // Inputs:
677     //      Int_t   track   Track number which produced this hit.
678     //      Int_t   *vol    Array of Integer Hit information. See AliITShit.h
679     //      Float_t *hits   Array of Floating Hit information.  see AliITShit.h
680     // Outputs:
681     //      none.
682     // Return:
683     //      none.
684
685     TClonesArray &lhits = *fHits;
686     new(lhits[fNhits++]) AliITShit(fIshunt,track,vol,hits);
687 }
688 //______________________________________________________________________
689 void AliITS::InitModules(Int_t size,Int_t &nmodules){
690     // Initialize the modules array.
691     // Inputs:
692     //      Int_t size  Size of array of the number of modules to be
693     //                  created. If size <=0 then the number of modules
694     //                  is gotten from AliITSgeom class kept in fITSgeom.
695     // Outputs:
696     //      Int_t &nmodules The number of modules existing.
697     // Return:
698     //      none.
699
700     if(fITSmodules){ 
701         fITSmodules->Delete();
702         delete fITSmodules;
703     } // end fir fITSmoudles
704
705     Int_t nl,indexMAX,index;
706
707     if(size<=0){ // default to using data stored in AliITSgeom
708         if(fITSgeom==0) {
709             Error("InitModules","fITSgeom not defined");
710             return;
711         } // end if fITSgeom==0
712         nl = fITSgeom->GetNlayers();
713         indexMAX = fITSgeom->GetModuleIndex(nl,fITSgeom->GetNladders(nl),
714                                             fITSgeom->GetNdetectors(nl))+1;
715         nmodules = indexMAX;
716         fITSmodules = new TObjArray(indexMAX);
717         for(index=0;index<indexMAX;index++){
718             fITSmodules->AddAt( new AliITSmodule(index),index);
719         } // end for index
720     }else{
721         fITSmodules = new TObjArray(size);
722         for(index=0;index<size;index++) {
723             fITSmodules->AddAt( new AliITSmodule(index),index);
724         } // end for index
725
726         nmodules = size;
727     } // end i size<=0
728 }
729 //______________________________________________________________________
730 void AliITS::FillModules(Int_t evnt,Int_t bgrev,Int_t nmodules,
731                          Option_t *option,Text_t *filename){
732     // fill the modules with the sorted by module hits; add hits from
733     // background if option=Add.
734     // Inputs:
735     //      Int_t evnt       Event to be processed.
736     //      Int_t bgrev      Background Hit tree number.
737     //      Int_t nmodules   Not used.
738     //      Option_t *option String indicating if merging hits or not. To
739     //                       merge hits set equal to "Add". Otherwise no
740     //                       background hits are considered.
741     //      Test_t *filename File name containing the background hits..
742     // Outputs:
743     //      none.
744     // Return:
745     //      none.
746     static TTree *trH1;                 //Tree with background hits
747     static Bool_t first=kTRUE;
748     static TFile *file;
749     const char *addBgr = strstr(option,"Add");
750
751     evnt = nmodules; // Dummy use of variables to remove warnings
752     if (addBgr ) {
753         if(first) {
754             file=new TFile(filename);
755         } // end if first
756         first=kFALSE;
757         file->cd();
758         file->ls();
759         // Get Hits Tree header from file
760         if(trH1) delete trH1;
761         trH1=0;
762
763         char treeName[20];
764         sprintf(treeName,"TreeH%d",bgrev);
765         trH1 = (TTree*)gDirectory->Get(treeName);
766         if (!trH1) {
767             Error("FillModules","cannot find Hits Tree for event:%d",bgrev);
768         } // end if !trH1
769         // Set branch addresses
770     } // end if addBgr
771
772     FillModules(fLoader->TreeH(),0); // fill from this file's tree.
773     
774     if (addBgr ) {
775         FillModules(trH1,10000000); // Default mask 10M.
776         TTree *fAli=fLoader->GetRunLoader()->TreeK();
777         TFile *fileAli=0;
778         if (fAli) fileAli =fAli->GetCurrentFile();
779         fileAli->cd();
780     } // end if add
781 }
782 //______________________________________________________________________
783 void AliITS::FillModules(TTree *treeH, Int_t mask) {
784     // fill the modules with the sorted by module hits; 
785     // can be called many times to do a merging
786     // Inputs:
787     //      TTree *treeH  The tree containing the hits to be copied into
788     //                    the modules.
789     //      Int_t mask    The track number mask to indecate which file
790     //                    this hits came from.
791     // Outputs:
792     //      none.
793     // Return:
794     //      none.
795
796     if (treeH == 0x0)
797      {
798        Error("FillModules","Tree is NULL");
799      }
800     Int_t lay,lad,det,index;
801     AliITShit *itsHit=0;
802     AliITSmodule *mod=0;
803     char branchname[20];
804     sprintf(branchname,"%s",GetName());
805     TBranch *branch = treeH->GetBranch(branchname);
806     if (!branch) {
807         Error("FillModules","%s branch in TreeH not found",branchname);
808         return;
809     } // end if !branch
810     branch->SetAddress(&fHits);
811     Int_t nTracks =(Int_t) treeH->GetEntries();
812     Int_t iPrimTrack,h;
813     for(iPrimTrack=0; iPrimTrack<nTracks; iPrimTrack++){
814         ResetHits();
815         Int_t nBytes = treeH->GetEvent(iPrimTrack);
816         if (nBytes <= 0) continue;
817         Int_t nHits = fHits->GetEntriesFast();
818         for(h=0; h<nHits; h++){
819             itsHit = (AliITShit *)fHits->UncheckedAt(h);
820             itsHit->GetDetectorID(lay,lad,det);
821             if (fITSgeom) {
822                 index = fITSgeom->GetModuleIndex(lay,lad,det);
823             } else {
824                 index=det-1; // This should not be used.
825             } // end if [You must have fITSgeom for this to work!]
826             mod = GetModule(index);
827             itsHit->SetTrack(itsHit->GetTrack()+mask); // Set track mask.
828             mod->AddHit(itsHit,iPrimTrack,h);
829         } // end loop over hits 
830     } // end loop over tracks
831 }
832 //______________________________________________________________________
833 void AliITS::ClearModules(){
834     // Clear the modules TObjArray.
835     // Inputs:
836     //      none.
837     // Outputs:
838     //      none.
839
840     if(fITSmodules) fITSmodules->Delete();
841 }
842 //______________________________________________________________________
843 void AliITS::MakeBranchS(const char *fl){
844     // Creates Tree Branch for the ITS summable digits.
845     // Inputs:
846     //      cont char *fl  File name where SDigits branch is to be written
847     //                     to. If blank it write the SDigits to the same
848     //                     file in which the Hits were found.
849     // Outputs:
850     //      none.
851     // Return:
852     //      none.
853     Int_t buffersize = 4000;
854     char branchname[30];
855
856     // only one branch for SDigits.
857     sprintf(branchname,"%s",GetName());
858     
859
860     if(fLoader->TreeS()){
861         if (fSDigits == 0x0)  fSDigits  = new TClonesArray("AliITSpListItem",1000);
862         MakeBranchInTree(fLoader->TreeS(),branchname,&fSDigits,buffersize,fl);
863     } // end if
864 }
865 //______________________________________________________________________
866 void AliITS::SetTreeAddressS(TTree *treeS){
867     // Set branch address for the ITS summable digits Trees.
868     // Inputs:
869     //      TTree *treeS   Tree containing the SDigits.
870     // Outputs:
871     //      none.
872     // Return:
873     //      none.
874     char branchname[30];
875
876     if(!treeS) return;
877     if (fSDigits == 0x0)  fSDigits = new TClonesArray("AliITSpListItem",1000);
878     TBranch *branch;
879     sprintf(branchname,"%s",GetName());
880     branch = treeS->GetBranch(branchname);
881     if (branch) branch->SetAddress(&fSDigits);
882 }
883 //______________________________________________________________________
884 void AliITS::MakeBranchInTreeD(TTree *treeD,const char *file){
885     // Creates Tree branches for the ITS.
886     // Inputs:
887     //      TTree     *treeD Pointer to the Digits Tree.
888     //      cont char *file  File name where Digits branch is to be written
889     //                       to. If blank it write the SDigits to the same
890     //                       file in which the Hits were found.
891     // Outputs:
892     //      none.
893     // Return:
894     //      none.
895     Int_t buffersize = 4000;
896     char branchname[30];
897
898     sprintf(branchname,"%s",GetName());
899     // one branch for digits per type of detector
900     const char *det[3] = {"SPD","SDD","SSD"};
901     char digclass[40];
902     char clclass[40];
903     Int_t i;
904     for (i=0; i<kNTYPES ;i++) {
905         DetType(i)->GetClassNames(digclass,clclass);
906         // digits
907         if (fDtype == 0x0) fDtype = new TObjArray(fNDetTypes);
908         if(!(fDtype->At(i))) fDtype->AddAt(new TClonesArray(digclass,1000),i);
909         else ResetDigits(i);
910     } // end for i
911     for (i=0; i<kNTYPES ;i++) {
912         if (kNTYPES==3) sprintf(branchname,"%sDigits%s",GetName(),det[i]);
913         else  sprintf(branchname,"%sDigits%d",GetName(),i+1);      
914         if (fDtype && treeD) {
915             MakeBranchInTree(treeD, branchname, &((*fDtype)[i]),buffersize,file);
916         } // end if
917     } // end for i
918 }
919 //______________________________________________________________________
920 void AliITS::SetTreeAddressD(TTree *treeD){
921     // Set branch address for the Trees.
922     // Inputs:
923     //      TTree *treeD   Tree containing the Digits.
924     // Outputs:
925     //      none.
926     // Return:
927     //      none.
928     char branchname[30];
929     const char *det[3] = {"SPD","SDD","SSD"};
930     TBranch *branch;
931     char digclass[40];
932     char clclass[40];
933     Int_t i;
934
935     if(!treeD) return;
936     for (i=0; i<kNTYPES; i++) {
937         DetType(i)->GetClassNames(digclass,clclass);
938         // digits
939         if (fDtype == 0x0) fDtype = new TObjArray(fNDetTypes);
940         if(!(fDtype->At(i))) fDtype->AddAt(new TClonesArray(digclass,1000),i);
941         else ResetDigits(i);
942         if (kNTYPES==3) sprintf(branchname,"%sDigits%s",GetName(),det[i]);
943         else  sprintf(branchname,"%sDigits%d",GetName(),i+1);
944         if (fDtype) {
945             branch = treeD->GetBranch(branchname);
946             if (branch) branch->SetAddress(&((*fDtype)[i]));
947         } // end if fDtype
948     } // end for i
949 }
950 //______________________________________________________________________
951 void AliITS::Hits2SDigits(){
952     // Standard Hits to summable Digits function.
953     // Inputs:
954     //      none.
955     // Outputs:
956     //      none.
957
958 //    return; // Using Hits in place of the larger sDigits.
959     fLoader->LoadHits("read");
960     fLoader->LoadSDigits("recreate");
961     AliRunLoader* rl = fLoader->GetRunLoader(); 
962
963     for (Int_t iEvent = 0; iEvent < rl->GetNumberOfEvents(); iEvent++) {
964     // Do the Hits to Digits operation. Use Standard input values.
965     // Event number from file, no background hit merging , use size from
966     // AliITSgeom class, option="All", input from this file only.
967       rl->GetEvent(iEvent);
968       if (!fLoader->TreeS()) fLoader->MakeTree("S");
969       MakeBranch("S");
970       SetTreeAddress();
971       HitsToSDigits(iEvent,0,-1," ",fOpt," ");
972     }
973     
974     fLoader->UnloadHits();
975     fLoader->UnloadSDigits();
976 }
977 //______________________________________________________________________
978 void AliITS::Hits2PreDigits(){
979     // Standard Hits to summable Digits function.
980     // Inputs:
981     //      none.
982     // Outputs:
983     //      none.
984
985     AliHeader *header=fLoader->GetRunLoader()->GetHeader(); // Get event number from this file.
986     // Do the Hits to Digits operation. Use Standard input values.
987     // Event number from file, no background hit merging , use size from
988     // AliITSgeom class, option="All", input from this file only.
989     HitsToPreDigits(header->GetEvent(),0,-1," ",fOpt," ");
990 }
991 //______________________________________________________________________
992 AliDigitizer* AliITS::CreateDigitizer(AliRunDigitizer* manager) const
993 {
994   return new AliITSDigitizer(manager);
995 }
996 //______________________________________________________________________
997 void AliITS::SDigitsToDigits(Option_t *opt){
998     // Standard Summable digits to Digits function.
999     // Inputs:
1000     //      none.
1001     // Outputs:
1002     //      none.
1003     char name[20] = {0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0};
1004
1005     if(!GetITSgeom()) return; // need transformations to do digitization.
1006     AliITSgeom *geom = GetITSgeom();
1007
1008     const char *all = strstr(opt,"All");
1009     const char *det[3] = {strstr(opt,"SPD"),strstr(opt,"SDD"),
1010                           strstr(opt,"SSD")};
1011     if( !det[0] && !det[1] && !det[2] ) all = "All";
1012     else all = 0;
1013     static Bool_t setDef=kTRUE;
1014     if (setDef) SetDefaultSimulation();
1015     setDef=kFALSE;
1016
1017     AliITSsimulation *sim      = 0;
1018     AliITSDetType    *iDetType = 0;
1019     TTree            *trees    = fLoader->TreeS();
1020     if( !(trees && this->GetSDigits()) ){
1021         Error("SDigits2Digits","Error: No trees or SDigits. Returning.");
1022         return;
1023     } // end if
1024     sprintf( name, "%s", this->GetName() );
1025     TBranch *brchSDigits = trees->GetBranch( name );
1026     
1027     Int_t id,module;
1028     for(module=0;module<geom->GetIndexMax();module++){
1029         id       = geom->GetModuleType(module);
1030         if (!all && !det[id]) continue;
1031         iDetType = DetType(id);
1032         sim      = (AliITSsimulation*)iDetType->GetSimulationModel();
1033         if (!sim) {
1034             Error("SDigit2Digits",
1035                   "The simulation class was not instanciated!");
1036             exit(1);
1037         } // end if !sim
1038         sim->InitSimulationModule(module,gAlice->GetEvNumber());
1039 //
1040         // add summable digits to module
1041         this->GetSDigits()->Clear();
1042         brchSDigits->GetEvent(module);
1043         sim->AddSDigitsToModule(GetSDigits(),0);
1044 //
1045         // Digitise current module sum(SDigits)->Digits
1046         sim->FinishSDigitiseModule();
1047
1048         // fills all branches - wasted disk space
1049         fLoader->TreeD()->Fill();
1050         this->ResetDigits();
1051     } // end for module
1052
1053     fLoader->TreeD()->GetEntries();
1054
1055     fLoader->TreeD()->AutoSave();
1056     // reset tree
1057     fLoader->TreeD()->Reset();
1058     
1059 }
1060 //______________________________________________________________________
1061 void AliITS::Hits2Digits(){
1062     // Standard Hits to Digits function.
1063     // Inputs:
1064     //      none.
1065     // Outputs:
1066     //      none.
1067
1068     fLoader->LoadHits("read");
1069     fLoader->LoadDigits("recreate");
1070     AliRunLoader* rl = fLoader->GetRunLoader(); 
1071
1072     for (Int_t iEvent = 0; iEvent < rl->GetNumberOfEvents(); iEvent++) {
1073     // Do the Hits to Digits operation. Use Standard input values.
1074     // Event number from file, no background hit merging , use size from
1075     // AliITSgeom class, option="All", input from this file only.
1076       rl->GetEvent(iEvent);
1077       if (!fLoader->TreeD()) fLoader->MakeTree("D");
1078       MakeBranch("D");
1079       SetTreeAddress();   
1080       HitsToDigits(iEvent,0,-1," ",fOpt," ");
1081     }
1082
1083     fLoader->UnloadHits();
1084     fLoader->UnloadDigits();
1085 }
1086 //______________________________________________________________________
1087 void AliITS::HitsToSDigits(Int_t evNumber,Int_t bgrev,Int_t size,
1088                           Option_t *option, Option_t *opt,Text_t *filename){
1089     // keep galice.root for signal and name differently the file for 
1090     // background when add! otherwise the track info for signal will be lost !
1091     // the condition below will disappear when the geom class will be
1092     // initialized for all versions - for the moment it is only for v5 !
1093     // 7 is the SDD beam test version. Dummy routine. Hits are ITS's Summable
1094     // Digits.
1095     // Inputs:
1096     //      Int_t evnt       Event to be processed.
1097     //      Int_t bgrev      Background Hit tree number.
1098     //      Int_t nmodules   Not used.
1099     //      Option_t *option String indicating if merging hits or not. To
1100     //                       merge hits set equal to "Add". Otherwise no
1101     //                       background hits are considered.
1102     //      Test_t *filename File name containing the background hits..
1103     // Outputs:
1104     //      none.
1105     // Return:
1106     //      none.
1107 //    return; // using Hits instead of the larger sdigits.
1108
1109     HitsToPreDigits(evNumber,bgrev,size,option,opt,filename);
1110 }
1111 //______________________________________________________________________
1112 void AliITS::HitsToPreDigits(Int_t evNumber,Int_t bgrev,Int_t size,
1113                           Option_t *option, Option_t *opt,Text_t *filename){
1114     //   Keep galice.root for signal and name differently the file for 
1115     // background when add! otherwise the track info for signal will be lost !
1116     // the condition below will disappear when the geom class will be
1117     // initialized for all versions - for the moment it is only for v5 !
1118     // 7 is the SDD beam test version.
1119     // Inputs:
1120     //      Int_t evnt       Event to be processed.
1121     //      Int_t bgrev      Background Hit tree number.
1122     //      Int_t nmodules   Not used.
1123     //      Option_t *option String indicating if merging hits or not. To
1124     //                       merge hits set equal to "Add". Otherwise no
1125     //                       background hits are considered.
1126     //      Test_t *filename File name containing the background hits..
1127     // Outputs:
1128     //      none.
1129     // Return:
1130     //      none.
1131
1132     if(!GetITSgeom()) return; // need transformations to do digitization.
1133     AliITSgeom *geom = GetITSgeom();
1134
1135     const char *all = strstr(opt,"All");
1136     const char *det[3] = {strstr(opt,"SPD"),strstr(opt,"SDD"),
1137                           strstr(opt,"SSD")};
1138     static Bool_t setDef=kTRUE;
1139     if (setDef) SetDefaultSimulation();
1140     setDef=kFALSE;
1141
1142     Int_t nmodules;
1143     InitModules(size,nmodules);
1144     FillModules(evNumber,bgrev,nmodules,option,filename);
1145
1146     AliITSsimulation *sim      = 0;
1147     AliITSDetType    *iDetType = 0;
1148     AliITSmodule     *mod      = 0;
1149     Int_t id,module;
1150     for(module=0;module<geom->GetIndexMax();module++){
1151         id       = geom->GetModuleType(module);
1152         if (!all && !det[id]) continue;
1153         iDetType = DetType(id);
1154         sim      = (AliITSsimulation*)iDetType->GetSimulationModel();
1155         if (!sim) {
1156             Error("HitsToSDigits",
1157                   "The simulation class was not instanciated!");
1158             exit(1);
1159         } // end if !sim
1160         mod      = (AliITSmodule *)fITSmodules->At(module);
1161         sim->SDigitiseModule(mod,module,evNumber);
1162         // fills all branches - wasted disk space
1163         fLoader->TreeS()->Fill(); 
1164         ResetSDigits();
1165     } // end for module
1166
1167     ClearModules();
1168
1169     fLoader->TreeS()->GetEntries();
1170     fLoader->TreeS()->AutoSave();
1171     fLoader->WriteSDigits("OVERWRITE");
1172     // reset tree
1173     fLoader->TreeS()->Reset();
1174 }
1175 //______________________________________________________________________
1176 void AliITS::HitsToDigits(Int_t evNumber,Int_t bgrev,Int_t size,
1177                           Option_t *option, Option_t *opt,Text_t *filename){
1178     //   Keep galice.root for signal and name differently the file for 
1179     // background when add! otherwise the track info for signal will be lost !
1180     // the condition below will disappear when the geom class will be
1181     // initialized for all versions - for the moment it is only for v5 !
1182     // 7 is the SDD beam test version.
1183     // Inputs:
1184     //      Int_t evnt       Event to be processed.
1185     //      Int_t bgrev      Background Hit tree number.
1186     //      Int_t nmodules   Not used.
1187     //      Option_t *option String indicating if merging hits or not. To
1188     //                       merge hits set equal to "Add". Otherwise no
1189     //                       background hits are considered.
1190     //      Test_t *filename File name containing the background hits..
1191     // Outputs:
1192     //      none.
1193     // Return:
1194     //      none.
1195
1196     if(!GetITSgeom()) return; // need transformations to do digitization.
1197     AliITSgeom *geom = GetITSgeom();
1198
1199     const char *all = strstr(opt,"All");
1200     const char *det[3] = {strstr(opt,"SPD"),strstr(opt,"SDD"),
1201                           strstr(opt,"SSD")};
1202     static Bool_t setDef=kTRUE;
1203     if (setDef) SetDefaultSimulation();
1204     setDef=kFALSE;
1205
1206     Int_t nmodules;
1207     InitModules(size,nmodules);
1208     FillModules(evNumber,bgrev,nmodules,option,filename);
1209
1210     AliITSsimulation *sim      = 0;
1211     AliITSDetType    *iDetType = 0;
1212     AliITSmodule     *mod      = 0;
1213     Int_t id,module;
1214     for(module=0;module<geom->GetIndexMax();module++){
1215         id       = geom->GetModuleType(module);
1216         if (!all && !det[id]) continue;
1217         iDetType = DetType(id);
1218         sim      = (AliITSsimulation*)iDetType->GetSimulationModel();
1219         if (!sim) {
1220             Error("HitsToDigits",
1221                   "The simulation class was not instanciated!");
1222             exit(1);
1223         } // end if !sim
1224         mod      = (AliITSmodule *)fITSmodules->At(module);
1225         sim->DigitiseModule(mod,module,evNumber);
1226         // fills all branches - wasted disk space
1227         fLoader->TreeD()->Fill(); 
1228         ResetDigits();
1229     } // end for module
1230
1231     ClearModules();
1232
1233     fLoader->TreeD()->GetEntries();
1234     fLoader->TreeD()->AutoSave();
1235     // reset tree
1236     fLoader->TreeD()->Reset();
1237 }
1238 //______________________________________________________________________
1239 void AliITS::ResetSDigits(){
1240     // Reset the Summable Digits array.
1241     // Inputs:
1242     //      none.
1243     // Outputs:
1244     //      none.
1245
1246     if (fSDigits) fSDigits->Clear();
1247     fNSDigits = 0;
1248 }
1249 //______________________________________________________________________
1250 void AliITS::ResetDigits(){
1251     // Reset number of digits and the digits array for the ITS detector.
1252     // Inputs:
1253     //      none.
1254     // Outputs:
1255     //      none.
1256
1257     if (!fDtype) return;
1258
1259     Int_t i;
1260     for (i=0;i<kNTYPES;i++ ) {
1261         if (fDtype->At(i))    ((TClonesArray*)fDtype->At(i))->Clear();
1262         if (fNdtype)  fNdtype[i]=0;
1263     } // end for i
1264 }
1265 //______________________________________________________________________
1266 void AliITS::ResetDigits(Int_t i){
1267     // Reset number of digits and the digits array for this branch.
1268     // Inputs:
1269     //      none.
1270     // Outputs:
1271     //      none.
1272
1273     if (fDtype->At(i))    ((TClonesArray*)fDtype->At(i))->Clear();
1274     if (fNdtype)  fNdtype[i]=0;
1275 }
1276 //______________________________________________________________________
1277 void AliITS::AddSumDigit(AliITSpListItem &sdig){
1278     // Adds the a module full of summable digits to the summable digits tree.
1279     // Inputs:
1280     //      AliITSpListItem &sdig   SDigit to be added to SDigits tree.
1281     // Outputs:
1282     //      none.
1283     // Return:
1284     //      none.
1285
1286     TClonesArray &lsdig = *fSDigits;
1287     new(lsdig[fNSDigits++]) AliITSpListItem(sdig);
1288 }
1289 //______________________________________________________________________
1290 void AliITS::AddRealDigit(Int_t id, Int_t *digits){
1291     //   Add a real digit - as coming from data.
1292     // Inputs:
1293     //      Int_t id        Detector type number.
1294     //      Int_t *digits   Integer array containing the digits info. See 
1295     //                      AliITSdigit.h
1296     // Outputs:
1297     //      none.
1298     // Return:
1299     //      none.
1300
1301     TClonesArray &ldigits = *((TClonesArray*)fDtype->At(id));
1302     new(ldigits[fNdtype[id]++]) AliITSdigit(digits);
1303 }
1304 //______________________________________________________________________
1305 void AliITS::AddSimDigit(Int_t id, AliITSdigit *d){
1306     //    Add a simulated digit.
1307     // Inputs:
1308     //      Int_t id        Detector type number.
1309     //      AliITSdigit *d  Digit to be added to the Digits Tree. See 
1310     //                      AliITSdigit.h
1311     // Outputs:
1312     //      none.
1313     // Return:
1314     //      none.
1315
1316     TClonesArray &ldigits = *((TClonesArray*)fDtype->At(id));
1317
1318     switch(id){
1319     case 0:
1320         new(ldigits[fNdtype[id]++]) AliITSdigitSPD(*((AliITSdigitSPD*)d));
1321         break;
1322     case 1:
1323         new(ldigits[fNdtype[id]++]) AliITSdigitSDD(*((AliITSdigitSDD*)d));
1324         break;
1325     case 2:
1326         new(ldigits[fNdtype[id]++]) AliITSdigitSSD(*((AliITSdigitSSD*)d));
1327         break;
1328     } // end switch id
1329 }
1330 //______________________________________________________________________
1331 void AliITS::AddSimDigit(Int_t id,Float_t phys,Int_t *digits,Int_t *tracks,
1332                          Int_t *hits,Float_t *charges){
1333     //   Add a simulated digit to the list.
1334     // Inputs:
1335     //      Int_t id        Detector type number.
1336     //      Float_t phys    Physics indicator. See AliITSdigits.h
1337     //      Int_t *digits   Integer array containing the digits info. See 
1338     //                      AliITSdigit.h
1339     //      Int_t *tracks   Integer array [AliITSdigitS?D::GetNTracks()] 
1340     //                      containing the track numbers that contributed to
1341     //                      this digit.
1342     //      Int_t *hits     Integer array [AliITSdigitS?D::GetNTracks()]
1343     //                      containing the hit numbers, from AliITSmodule, that
1344     //                      contributed to this digit.
1345     //      Float_t *charge Floating point array of the signals contributed
1346     //                      to this digit by each track.
1347     // Outputs:
1348     //      none.
1349     // Return:
1350     //      none.
1351
1352     TClonesArray &ldigits = *((TClonesArray*)fDtype->At(id));
1353     switch(id){
1354     case 0:
1355         new(ldigits[fNdtype[id]++]) AliITSdigitSPD(digits,tracks,hits);
1356         break;
1357     case 1:
1358         new(ldigits[fNdtype[id]++]) AliITSdigitSDD(phys,digits,tracks,
1359                                                    hits,charges);
1360         break;
1361     case 2:
1362         new(ldigits[fNdtype[id]++]) AliITSdigitSSD(digits,tracks,hits);
1363         break;
1364     } // end switch id
1365 }
1366 //______________________________________________________________________
1367 void AliITS::MakeTreeC(Option_t *option){
1368   //   Create a separate tree to store the clusters.
1369   // Inputs:
1370   //      Option_t *option  string which must contain "C" otherwise
1371   //                        no Cluster Tree is created.
1372   // Outputs:
1373   //      none.
1374   // Return:
1375   //      none.
1376
1377   AliITSLoader *pITSLoader = (AliITSLoader*)fLoader;    
1378     
1379   if (pITSLoader == 0x0) {
1380     Error("MakeTreeC","fLoader == 0x0 option=%s",option);
1381     return;
1382   }
1383   if (pITSLoader->TreeC() == 0x0) pITSLoader->MakeTree("C");
1384   MakeBranchC();
1385 }
1386
1387 void AliITS::MakeBranchC()
1388 {
1389 //Makes barnches in treeC
1390   AliITSLoader *pITSLoader = (AliITSLoader*)fLoader;    
1391   if (pITSLoader == 0x0) 
1392    {
1393     Error("MakeTreeC","fLoader == 0x0");
1394     return;
1395    }
1396   TTree * lTC = pITSLoader->TreeC();
1397   if (lTC == 0x0)
1398    {
1399      Error("MakeTreeC","Can not get TreeC from Loader");
1400      return;
1401    }
1402
1403   Int_t buffersize = 4000;
1404   char branchname[30];
1405   const char *det[3] = {"SPD","SDD","SSD"};
1406   char digclass[40];
1407   char clclass[40];
1408
1409     // one branch for Clusters per type of detector
1410   Int_t i;   
1411   for (i=0; i<kNTYPES ;i++) 
1412     {
1413         AliITSDetType *iDetType=DetType(i); 
1414         iDetType->GetClassNames(digclass,clclass);
1415         // clusters
1416         if (fCtype == 0x0) fCtype  = new TObjArray(fNDetTypes);
1417         if(!ClustersAddress(i))
1418          {
1419           fCtype->AddAt(new TClonesArray(clclass,1000),i);
1420          }
1421         if (kNTYPES==3) sprintf(branchname,"%sClusters%s",GetName(),det[i]);
1422         else  sprintf(branchname,"%sClusters%d",GetName(),i+1);
1423         if (fCtype  && lTC) 
1424          {
1425            if (lTC->GetBranch(branchname))
1426             {
1427               Warning("MakeBranchC","Branch %s alread exists in TreeC",branchname);
1428             }
1429            else
1430             {
1431               Info("MakeBranchC","Creating branch %s in TreeC",branchname);
1432               lTC->Branch(branchname,&((*fCtype)[i]), buffersize);
1433             }
1434          } // end if fCtype && lTC
1435   } // end for i
1436 }
1437
1438 //______________________________________________________________________
1439 void AliITS::GetTreeC(Int_t event){
1440   //    Get the clusters tree for this event and set the branch address.
1441   // Inputs:
1442   //      Int_t event    Event number for the cluster tree.
1443   // Outputs:
1444   //      none.
1445   // Return:
1446   //      none.
1447   char branchname[30];
1448   const char *det[3] = {"SPD","SDD","SSD"};
1449
1450   AliITSLoader *pITSLoader = (AliITSLoader*)fLoader;
1451   TTree * lTC = pITSLoader->TreeC();
1452
1453   ResetClusters();
1454   if (lTC) {
1455     pITSLoader->CleanRawClusters();
1456   } // end if TreeC()
1457
1458
1459   TBranch *branch;
1460
1461   if (lTC) {
1462     Int_t i;
1463         char digclass[40];
1464         char clclass[40];
1465         for (i=0; i<kNTYPES; i++) {
1466       AliITSDetType *iDetType=DetType(i); 
1467       iDetType->GetClassNames(digclass,clclass);
1468       // clusters
1469       if (fCtype == 0x0) fCtype  = new TObjArray(fNDetTypes);
1470       if(!fCtype->At(i)) fCtype->AddAt(new TClonesArray(clclass,1000),i);
1471       if(kNTYPES==3) sprintf(branchname,"%sClusters%s",GetName(),det[i]);
1472       else  sprintf(branchname,"%sClusters%d",GetName(),i+1);
1473       if (fCtype) {
1474                 branch = lTC->GetBranch(branchname);
1475         if (branch) branch->SetAddress(&((*fCtype)[i]));
1476       } // end if fCtype
1477         } // end for i
1478   } else {
1479         Error("GetTreeC","cannot find Clusters Tree for event:%d",event);
1480   } // end if lTC
1481 }
1482 //______________________________________________________________________
1483 void AliITS::AddCluster(Int_t id, AliITSRawCluster *c){
1484     //   Add a cluster to the list.
1485     // Inputs:
1486     //      Int_t id             Detector type number.
1487     //      AliITSRawCluster *c  Cluster class to be added to the tree of
1488     //                           clusters.
1489     // Outputs:
1490     //      none.
1491     // Return:
1492     //      none.
1493
1494     TClonesArray &lc = *((TClonesArray*)fCtype->At(id));
1495
1496     switch(id){
1497     case 0:
1498         new(lc[fNctype[id]++]) AliITSRawClusterSPD(*((AliITSRawClusterSPD*)c));
1499         break;
1500     case 1:
1501         new(lc[fNctype[id]++]) AliITSRawClusterSDD(*((AliITSRawClusterSDD*)c));
1502         break;
1503     case 2:
1504         new(lc[fNctype[id]++]) AliITSRawClusterSSD(*((AliITSRawClusterSSD*)c));
1505         break;
1506     } // end switch id
1507 }
1508 //______________________________________________________________________
1509 void AliITS::ResetClusters(){
1510     // Reset number of clusters and the clusters array for ITS.
1511     // Inputs:
1512     //      none.
1513     // Outputs:
1514     //      none.
1515
1516     Int_t i;
1517     for (i=0;i<kNTYPES;i++ ) ResetClusters(i);
1518 }
1519 //______________________________________________________________________
1520 void AliITS::ResetClusters(Int_t i){
1521     //    Reset number of clusters and the clusters array for this branch.
1522     // Inputs:
1523     //      Int_t i        Detector type number.
1524     // Outputs:
1525     //      none.
1526     // Return:
1527     //      none.
1528
1529     if (fCtype->At(i))    ((TClonesArray*)fCtype->At(i))->Clear();
1530     if (fNctype)  fNctype[i]=0;
1531 }
1532 //______________________________________________________________________
1533 void AliITS::MakeBranchR(const char *file, Option_t *opt){
1534     // Creates Tree branches for the ITS Reconstructed points.
1535     // Inputs:
1536     //      cont char *file  File name where RecPoints branch is to be written
1537     //                       to. If blank it write the SDigits to the same
1538     //                       file in which the Hits were found.
1539     // Outputs:
1540     //      none.
1541     // Return:
1542     //      none.
1543     Int_t buffsz = 4000;
1544     char branchname[30];
1545
1546     // only one branch for rec points for all detector types
1547     Bool_t oFast= (strstr(opt,"Fast")!=0);
1548     if(oFast){
1549       sprintf(branchname,"%sRecPointsF",GetName());
1550     } else {
1551       sprintf(branchname,"%sRecPoints",GetName());
1552     }
1553     
1554     
1555     if(!fRecPoints)fRecPoints = new TClonesArray("AliITSRecPoint",1000);
1556     if (fLoader->TreeR()) {
1557         if (fRecPoints == 0x0) fRecPoints = new TClonesArray("AliITSRecPoint",1000);
1558         MakeBranchInTree(fLoader->TreeR(),branchname,&fRecPoints,buffsz,file);
1559     } // end if
1560 }
1561 //______________________________________________________________________
1562 void AliITS::SetTreeAddressR(TTree *treeR){
1563     // Set branch address for the Reconstructed points Trees.
1564     // Inputs:
1565     //      TTree *treeR   Tree containing the RecPoints.
1566     // Outputs:
1567     //      none.
1568     // Return:
1569     //      none.
1570     char branchname[30];
1571
1572     if(!treeR) return;
1573     if (fRecPoints == 0x0) fRecPoints = new TClonesArray("AliITSRecPoint",1000);
1574     TBranch *branch;
1575     sprintf(branchname,"%sRecPoints",GetName());
1576     branch = treeR->GetBranch(branchname);
1577     if (branch) {
1578       branch->SetAddress(&fRecPoints);
1579     }
1580     else {
1581       sprintf(branchname,"%sRecPointsF",GetName());
1582       branch = treeR->GetBranch(branchname);
1583       if (branch) {
1584         branch->SetAddress(&fRecPoints);
1585       }
1586     }
1587 }
1588 //______________________________________________________________________
1589 void AliITS::AddRecPoint(const AliITSRecPoint &r){
1590     // Add a reconstructed space point to the list
1591     // Inputs:
1592     //      const AliITSRecPoint &r RecPoint class to be added to the tree
1593     //                              of reconstructed points TreeR.
1594     // Outputs:
1595     //      none.
1596     // Return:
1597     //      none.
1598
1599     TClonesArray &lrecp = *fRecPoints;
1600     new(lrecp[fNRecPoints++]) AliITSRecPoint(r);
1601 }
1602 //______________________________________________________________________
1603 void AliITS::HitsToFastRecPoints(Int_t evNumber,Int_t bgrev,Int_t size,
1604                                   Option_t *opt0,Option_t *opt1,Text_t *flnm){
1605     // keep galice.root for signal and name differently the file for 
1606     // background when add! otherwise the track info for signal will be lost !
1607     // the condition below will disappear when the geom class will be
1608     // initialized for all versions - for the moment it is only for v5 !
1609     // Inputs:
1610     //      Int_t evnt       Event to be processed.
1611     //      Int_t bgrev      Background Hit tree number.
1612     //      Int_t size       Size used by InitModules. See InitModules.
1613     //      Option_t *opt0   Option passed to FillModules. See FillModules.
1614     //      Option_t *opt1   String indicating if merging hits or not. To
1615     //                       merge hits set equal to "Add". Otherwise no
1616     //                       background hits are considered.
1617     //      Test_t *flnm     File name containing the background hits..
1618     // Outputs:
1619     //      none.
1620     // Return:
1621     //      none.
1622
1623     if(!GetITSgeom()) return;
1624     AliITSLoader *pITSloader = (AliITSLoader*)fLoader;
1625     AliITSgeom *geom = GetITSgeom();
1626
1627     const char *all = strstr(opt1,"All");
1628     const char *det[3] ={strstr(opt1,"SPD"),strstr(opt1,"SDD"),
1629                          strstr(opt1,"SSD")};
1630     Int_t nmodules;
1631     InitModules(size,nmodules);
1632     FillModules(evNumber,bgrev,nmodules,opt0,flnm);
1633
1634     AliITSsimulation *sim      = 0;
1635     AliITSDetType    *iDetType = 0;
1636     AliITSmodule     *mod      = 0;
1637     Int_t id,module;
1638
1639     //m.b. : this change is nothing but a nice way to make sure
1640     //the CPU goes up !
1641     
1642     cout<<"HitsToFastRecPoints: N mod = "<<geom->GetIndexMax()<<endl;
1643     
1644     for(module=0;module<geom->GetIndexMax();module++)
1645      {
1646         id       = geom->GetModuleType(module);
1647         if (!all && !det[id]) continue;
1648         iDetType = DetType(id);
1649         sim      = (AliITSsimulation*)iDetType->GetSimulationModel();
1650         if (!sim) 
1651          {
1652            Error("HitsToFastPoints","The simulation class was not instanciated!");
1653            exit(1);
1654          } // end if !sim
1655         mod      = (AliITSmodule *)fITSmodules->At(module);
1656         sim->CreateFastRecPoints(mod,module,gRandom);
1657         cout<<module<<"\r";fflush(0);
1658         //gAlice->TreeR()->Fill();
1659         TTree *lTR = pITSloader->TreeR();
1660         TBranch *br=lTR->GetBranch("ITSRecPointsF");
1661         br->Fill();
1662         ResetRecPoints();
1663     } // end for module
1664
1665     ClearModules();
1666     
1667     fLoader->WriteRecPoints("OVERWRITE");
1668 }
1669 //______________________________________________________________________
1670 void AliITS::Digits2Reco(){
1671     // Find clusters and reconstruct space points.
1672     // Inputs:
1673     //      none.
1674     // Outputs:
1675     //      none.
1676
1677     AliHeader *header=fLoader->GetRunLoader()->GetHeader();
1678     // to Digits to RecPoints for event in file, all digits in file, and
1679     // all ITS detectors.
1680     DigitsToRecPoints(header->GetEvent(),0,fOpt);
1681 }
1682 //______________________________________________________________________
1683 void AliITS::DigitsToRecPoints(Int_t evNumber,Int_t lastentry,Option_t *opt){
1684   // cluster finding and reconstruction of space points
1685   // the condition below will disappear when the geom class will be
1686   // initialized for all versions - for the moment it is only for v5 !
1687   // 7 is the SDD beam test version
1688   // Inputs:
1689   //      Int_t evNumber   Event number to be processed.
1690   //      Int_t lastentry  Offset for module when not all of the modules
1691   //                       are processed.
1692   //      Option_t *opt    String indicating which ITS sub-detectors should
1693   //                       be processed. If ="All" then all of the ITS
1694   //                       sub detectors are processed.
1695   // Outputs:
1696   //      none.
1697   // Return:
1698   //      none.
1699
1700   if(!GetITSgeom()) return;
1701   AliITSgeom *geom = GetITSgeom();
1702     
1703   const char *all = strstr(opt,"All");
1704   const char *det[3] = {strstr(opt,"SPD"),strstr(opt,"SDD"),
1705                         strstr(opt,"SSD")};
1706   static Bool_t setRec=kTRUE;
1707   if (setRec) SetDefaultClusterFinders();
1708   setRec=kFALSE;
1709
1710   AliITSLoader *pITSloader = (AliITSLoader*)fLoader;
1711   TTree *treeC=pITSloader->TreeC();
1712   AliITSClusterFinder *rec     = 0;
1713   AliITSDetType      *iDetType = 0;
1714   Int_t id,module,first=0;
1715   for(module=0;module<geom->GetIndexMax();module++){
1716       id       = geom->GetModuleType(module);
1717       if (!all && !det[id]) continue;
1718       if(det[id]) first = geom->GetStartDet(id);
1719       iDetType = DetType(id);
1720       rec = (AliITSClusterFinder*)iDetType->GetReconstructionModel();
1721       TClonesArray *itsDigits  = this->DigitsAddress(id);
1722       if (!rec) {
1723           Error("DigitsToRecPoints",
1724                 "The reconstruction class was not instanciated! event=%d",
1725                 evNumber);
1726           exit(1);
1727       } // end if !rec
1728       this->ResetDigits();
1729       TTree *lTD = pITSloader->TreeD();
1730       if (all) {
1731           lTD->GetEvent(lastentry+module);
1732       }else {
1733           lTD->GetEvent(lastentry+(module-first));
1734       }
1735       Int_t ndigits = itsDigits->GetEntriesFast();
1736       if (ndigits) rec->FindRawClusters(module);
1737       pITSloader->TreeR()->Fill(); 
1738       ResetRecPoints();
1739       treeC->Fill();
1740       ResetClusters();
1741   } // end for module
1742
1743
1744   pITSloader->WriteRecPoints("OVERWRITE");
1745   pITSloader->WriteRawClusters("OVERWRITE");
1746 }
1747 //______________________________________________________________________
1748 void AliITS::ResetRecPoints(){
1749     // Reset number of rec points and the rec points array.
1750     // Inputs:
1751     //      none.
1752     // Outputs:
1753     //      none.
1754
1755     if (fRecPoints) fRecPoints->Clear();
1756     fNRecPoints = 0;
1757 }
1758 //______________________________________________________________________
1759 AliLoader* AliITS::MakeLoader(const char* topfoldername)
1760
1761   //builds ITSgetter (AliLoader type)
1762   //if detector wants to use castomized getter, it must overload this method
1763
1764   Info("MakeLoader","Creating AliITSLoader. Top folder is %s.",topfoldername);
1765   fLoader = new AliITSLoader(GetName(),topfoldername);
1766   return fLoader;
1767 }
1768
1769
1770 //_____________________________________________________________________________
1771 void AliITS::Reconstruct() const
1772 {
1773 // reconstruct clusters
1774
1775   AliLoader* loader = GetLoader();
1776   loader->LoadRecPoints("recreate");
1777   loader->LoadDigits("read");
1778
1779   AliITSclustererV2 clusterer(GetITSgeom());
1780   AliRunLoader* runLoader = loader->GetRunLoader();
1781   Int_t nEvents = runLoader->GetNumberOfEvents();
1782
1783   for (Int_t iEvent = 0; iEvent < nEvents; iEvent++) {
1784     runLoader->GetEvent(iEvent);
1785
1786     TTree* treeClusters = loader->TreeR();
1787     if (!treeClusters) {
1788       loader->MakeTree("R");
1789       treeClusters = loader->TreeR();
1790     }
1791     TTree* treeDigits = loader->TreeD();
1792     if (!treeDigits) {
1793       Error("Reconstruct", "Can't get digits tree !");
1794       return;
1795     }
1796
1797     clusterer.Digits2Clusters(treeDigits, treeClusters);
1798          
1799     loader->WriteRecPoints("OVERWRITE");
1800   }
1801
1802   loader->UnloadRecPoints();
1803   loader->UnloadDigits();
1804 }
1805
1806 //_____________________________________________________________________________
1807 AliTracker* AliITS::CreateTracker() const
1808 {
1809 // create an ITS tracker
1810
1811   return new AliITStrackerV2(GetITSgeom());
1812 }
1813
1814 //_____________________________________________________________________________
1815 void AliITS::FillESD(AliESD* esd) const
1816 {
1817 // make PID, find V0s and cascades
1818
1819   Double_t parITS[] = {34., 0.15, 10.};
1820   AliITSpidESD itsPID(parITS);
1821   itsPID.MakePID(esd);
1822
1823   // V0 finding
1824   Double_t cuts[]={33,  // max. allowed chi2
1825                    0.16,// min. allowed negative daughter's impact parameter 
1826                    0.05,// min. allowed positive daughter's impact parameter 
1827                    0.080,// max. allowed DCA between the daughter tracks
1828                    0.998,// max. allowed cosine of V0's pointing angle
1829                    0.9,  // min. radius of the fiducial volume
1830                    2.9   // max. radius of the fiducial volume
1831   };
1832   AliV0vertexer vtxer(cuts);
1833   Double_t vtx[3], cvtx[6];
1834   esd->GetVertex(vtx, cvtx);
1835   vtxer.SetVertex(vtx);
1836   vtxer.Tracks2V0vertices(esd);
1837
1838   // cascade finding
1839   Double_t cts[]={33.,    // max. allowed chi2
1840                   0.05,   // min. allowed V0 impact parameter 
1841                   0.008,  // window around the Lambda mass 
1842                   0.035,  // min. allowed bachelor's impact parameter 
1843                   0.10,   // max. allowed DCA between a V0 and a track
1844                   0.9985, //max. allowed cosine of the cascade pointing angle
1845                   0.9,    // min. radius of the fiducial volume
1846                   2.9     // max. radius of the fiducial volume
1847   };
1848   AliCascadeVertexer cvtxer=AliCascadeVertexer(cts);
1849   cvtxer.V0sTracks2CascadeVertices(esd);
1850 }
1851