Fix to avoid divide by zero problem in MedianHitG and MedianHitL for track
[u/mrichter/AliRoot.git] / ITS / AliITSmodule.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 /*
17 $Log$
18 Revision 1.10  2002/03/15 17:21:54  nilsen
19 Removed zero-ing of fModules variable in constructors.
20
21 Revision 1.9  2000/10/04 19:46:39  barbera
22 Corrected by F. Carminati for v3.04
23
24 Revision 1.8  2000/10/02 16:32:57  barbera
25 Forward declarations added and formatting
26
27 Revision 1.3.4.8  2000/10/02 15:55:26  barbera
28 Forward declarations added and formatting
29
30 Revision 1.7  2000/09/22 12:36:38  nilsen
31 Minor changes to improve compilation and create less NOISE.
32
33 Revision 1.6  2000/07/10 16:07:18  fca
34 Release version of ITS code
35
36 Revision 1.3.4.2  2000/03/02 21:42:29  nilsen 
37 Linked AliDetector::fDigit to AliITSmodule::fDigitsM and AliITS::fITSRecPoints
38 to AliITSmodule::fRecPointsM. Renamed AliITSmodule::fPointsM to fRecPointsM.
39 Removed the deletion of fDigitsM from the distructor since it is only a copy
40 of what is in AliDetector. Fixed a bug in the functions LineSegmentL and
41 LineSegmentG. Added two new versions of LineSegmentL and LineSegmentG to
42 additionaly return track number from the hit. Removed FastPoint function,
43 haven't found anywere it was used, also it had very many problems, should
44 just call FastPointSPD... .
45
46 Revision 1.3.4.1  2000/01/12 19:03:32  nilsen
47 This is the version of the files after the merging done in December 1999.
48 See the ReadMe110100.txt file for details
49
50 Revision 1.3  1999/10/04 15:20:12  fca
51 Correct syntax accepted by g++ but not standard for static members, remove minor warnings
52
53 Revision 1.2  1999/09/29 09:24:20  fca
54 Introduction of the Copyright and cvs Log
55
56 */
57
58 #include <TArrayI.h>
59
60 #include <stdlib.h>
61
62 #include "AliRun.h"
63 #include "AliITS.h"
64 #include "AliITShit.h"
65 #include "AliITSmodule.h"
66 #include "AliITSgeom.h"
67
68 ClassImp(AliITSmodule)
69
70 //_______________________________________________________________________
71 //
72 // Impementation of class AliITSmodule
73 //
74 // created by: A. Bouchm, W. Peryt, S. Radomski, P. Skowronski 
75 //             R. Barbers, B. Batyunia, B. S. Nilsen
76 // ver 1.0     CERN 16.09.1999
77 //_______________________________________________________________________
78 //________________________________________________________________________
79 // 
80 // Constructors and deconstructor
81 //________________________________________________________________________
82 //
83 AliITSmodule::AliITSmodule() {
84     // constructor
85
86     fHitsM       = 0;
87     fTrackIndex  = 0;
88     fHitIndex    = 0;
89     fITS         = 0;
90
91 }
92 //_________________________________________________________________________
93 AliITSmodule::AliITSmodule(Int_t index) {
94   // constructor
95
96     fIndex      = index;
97     fHitsM      = new TObjArray();
98     fTrackIndex = new TArrayI(16);
99     fHitIndex   = new TArrayI(16);
100     fITS        = (AliITS*)(gAlice->GetDetector("ITS"));
101 }
102 //__________________________________________________________________________
103 AliITSmodule::~AliITSmodule() {
104     // The destructor for AliITSmodule. Before destoring AliITSmodule
105     // we must first destroy all of it's members.
106
107     if(fHitsM){
108         for(Int_t i=0;i<fHitsM->GetEntriesFast();i++) 
109             delete ((AliITShit *)(fHitsM->At(i)));
110         // must delete each object in the TObjArray.
111         delete fHitsM;
112     } // end if
113     delete fTrackIndex;
114     delete fHitIndex;
115     fITS = 0; // We don't delete this pointer since it is just a copy.
116 }
117 //____________________________________________________________________________
118 AliITSmodule::AliITSmodule(const AliITSmodule &source){
119 ////////////////////////////////////////////////////////////////////////
120 //     Copy Constructor 
121 ////////////////////////////////////////////////////////////////////////
122   printf("AliITSmodule error: AliITSmodule class has not to be copied! Exit.\n");
123   exit(1);
124 }
125
126 //_____________________________________________________________________________
127 AliITSmodule& AliITSmodule::operator=(const AliITSmodule &source){
128 ////////////////////////////////////////////////////////////////////////
129 //    Assignment operator 
130 ////////////////////////////////////////////////////////////////////////
131   printf("AliITSmodule error: AliITSmodule class has not to be copied! Exit.\n");
132   exit(1);
133   return *this; // fake return neded on Sun
134
135
136 //_________________________________________________________________________
137 // 
138 // Hits management
139 //__________________________________________________________________________
140 Int_t AliITSmodule::AddHit(AliITShit* hit,Int_t t,Int_t h) {
141 // Hits management
142
143   //printf("AddHit: beginning hit %p t h %d %d\n",hit,t,h);
144     fHitsM->AddLast(new AliITShit(*hit));
145     Int_t fNhitsM = fHitsM->GetEntriesFast();
146     if(fNhitsM-1>=fTrackIndex->GetSize()){ // need to expand the TArrayI
147         TArrayI *p = new TArrayI(fNhitsM+64);
148         for(Int_t i=0;i<fTrackIndex->GetSize();i++) 
149             (*p)[i] = fTrackIndex->At(i);
150         delete fTrackIndex;
151         fTrackIndex = p;
152     } // end if
153     if(fNhitsM-1>=fHitIndex->GetSize()){ // need to expand the TArrayI
154         TArrayI *p = new TArrayI(fNhitsM+64);
155         for(Int_t i=0;i<fHitIndex->GetSize();i++) 
156             (*p)[i] = fHitIndex->At(i);
157         delete fHitIndex;
158         fHitIndex = p;
159     } // end if
160     (*fTrackIndex)[fNhitsM-1] = t;
161     (*fHitIndex)[fNhitsM-1]   = h;
162     return fNhitsM;
163 }
164
165 //___________________________________________________________________________
166 void AliITSmodule::MedianHitG(Int_t index,
167                               Float_t hitx1,Float_t hity1,Float_t hitz1,
168                               Float_t hitx2,Float_t hity2,Float_t hitz2,
169                               Float_t &xMg, Float_t &yMg, Float_t &zMg){
170   // median hit
171    AliITSgeom *gm = fITS->GetITSgeom();
172    Float_t x1l,y1l,z1l;
173    Float_t x2l,y2l,z2l;
174    Float_t xMl,yMl=0,zMl;
175    Float_t l[3], g[3];
176
177    g[0] = hitx1;
178    g[1] = hity1;
179    g[2] = hitz1;
180    gm->GtoL(index,g,l);
181    x1l = l[0];
182    y1l = l[1];
183    z1l = l[2];
184
185    g[0] = hitx2;
186    g[1] = hity2;
187    g[2] = hitz2;
188    gm->GtoL(index,g,l);
189    x2l = l[0];
190    y2l = l[1];
191    z2l = l[2];
192
193    // Modified by N.Carrer. In very rare occasions the track may be just
194    // tangent to the module. Therefore the entrance and exit points have the
195    // same y.
196    if( (y2l-y1l) != 0.0 ) {
197      xMl = (-y1l / (y2l-y1l))*(x2l-x1l) + x1l;
198      zMl = (-y1l / (y2l-y1l))*(z2l-z1l) + z1l;
199    } else {
200      xMl = 0.5*(x1l+x2l);
201      zMl = 0.5*(z1l+z2l);
202    }
203
204    l[0] = xMl;
205    l[1] = yMl;
206    l[2] = zMl;
207    gm->LtoG(index,l,g);
208    xMg = g[0];
209    yMg = g[1];
210    zMg = g[2];
211 }
212 //___________________________________________________________________________
213 Double_t AliITSmodule::PathLength(Int_t index,AliITShit *itsHit1,
214                                   AliITShit *itsHit2){
215   // path lenght
216    Float_t  x1g,y1g,z1g;   
217    Float_t  x2g,y2g,z2g;
218    Double_t s;
219
220    itsHit1->GetPositionG(x1g,y1g,z1g);
221    itsHit2->GetPositionG(x2g,y2g,z2g);
222
223    s = TMath::Sqrt( ((Double_t)(x2g-x1g)*(Double_t)(x2g-x1g)) +
224                     ((Double_t)(y2g-y1g)*(Double_t)(y2g-y1g)) +
225                     ((Double_t)(z2g-z1g)*(Double_t)(z2g-z1g))  );
226    return s;
227 }
228 //___________________________________________________________________________
229 void AliITSmodule::PathLength(Int_t index,
230                               Float_t x,Float_t y,Float_t z,
231                               Int_t status,Int_t &nseg,
232                               Float_t &x1,Float_t &y1,Float_t &z1,
233                               Float_t &dx1,Float_t &dy1,Float_t &dz1,
234                               Int_t &flag){
235   // path length
236     static Float_t x0,y0,z0;
237
238     if (status == 66){ // entering
239         x0 = x;
240         y0 = y;
241         z0 = z;
242         nseg = 0;
243         flag = 1;
244     }else{
245         x1 = x0;
246         y1 = y0;
247         z1 = z0;
248         dx1 = x-x1;
249         dy1 = y-y1;
250         dz1 = z-z1;
251         nseg++;
252         if (status == 68) flag = 0; //exiting
253         else flag = 2; //inside
254         x0 = x;
255         y0 = y;
256         z0 = z;
257     } // end if
258 }
259 //___________________________________________________________________________
260 Bool_t AliITSmodule::LineSegmentL(Int_t hitindex,Double_t &a,Double_t &b,
261                                   Double_t &c,Double_t &d,
262                                   Double_t &e,Double_t &f,Double_t &de){
263   // line segment
264     static Int_t hitindex0;
265     AliITShit *h0,*h1;
266
267     if(hitindex>= fHitsM->GetEntriesFast()) return kFALSE;
268
269     h1 = (AliITShit *) (fHitsM->At(hitindex));
270     if(h1->StatusEntering()){ // if track entering volume, get index for next
271                               // step
272         hitindex0 = hitindex;
273         return kFALSE;
274     } // end if StatusEntering()
275     // else stepping
276     h0 = (AliITShit *) (fHitsM->At(hitindex0));
277     de = h1->GetIonization();
278     h0->GetPositionL(a,c,e);
279     h1->GetPositionL(b,d,f);
280     b = b - a;
281     d = d - c;
282     f = f - e;
283     hitindex0 = hitindex;
284     return kTRUE;
285 }
286 //___________________________________________________________________________
287 Bool_t AliITSmodule::LineSegmentG(Int_t hitindex,Double_t &a,Double_t &b,
288                                   Double_t &c,Double_t &d,
289                                   Double_t &e,Double_t &f,Double_t &de){
290   // line segment
291     static Int_t hitindex0;
292     AliITShit *h0,*h1;
293
294     if(hitindex>= fHitsM->GetEntriesFast()) return kFALSE;
295
296     h1 = (AliITShit *) (fHitsM->At(hitindex));
297     if(h1->StatusEntering()){ // if track entering volume, get index for next
298                               // step
299         hitindex0 = hitindex;
300         return kFALSE;
301     } // end if StatusEntering()
302     // else stepping
303     h0 = (AliITShit *) (fHitsM->At(hitindex0));
304     de = h1->GetIonization();
305     h0->GetPositionG(a,c,e);
306     h1->GetPositionG(b,d,f);
307     b = b - a;
308     d = d - c;
309     f = f - e;
310     hitindex0 = hitindex;
311     return kTRUE;
312 }
313 //___________________________________________________________________________
314 Bool_t AliITSmodule::LineSegmentL(Int_t hitindex,Double_t &a,Double_t &b,
315                                   Double_t &c,Double_t &d,
316                                   Double_t &e,Double_t &f,
317                                   Double_t &de,Int_t &track){
318   // line segmente
319     static Int_t hitindex0;
320     AliITShit *h0,*h1;
321
322     if(hitindex>= fHitsM->GetEntriesFast()) return kFALSE;
323
324     h1 = (AliITShit *) (fHitsM->At(hitindex));
325     if(h1->StatusEntering()){ // if track entering volume, get index for next
326                               // step
327         hitindex0 = hitindex;
328         track = h1->GetTrack();
329         return kFALSE;
330     } // end if StatusEntering()
331     // else stepping
332     h0 = (AliITShit *) (fHitsM->At(hitindex0));
333     de = h1->GetIonization();
334     h0->GetPositionL(a,c,e);
335     h1->GetPositionL(b,d,f);
336     b = b - a;
337     d = d - c;
338     f = f - e;
339     hitindex0 = hitindex;
340     track = h1->GetTrack();
341     return kTRUE;
342 }
343 //___________________________________________________________________________
344 Bool_t AliITSmodule::LineSegmentG(Int_t hitindex,Double_t &a,Double_t &b,
345                                   Double_t &c,Double_t &d,
346                                   Double_t &e,Double_t &f,
347                                   Double_t &de,Int_t &track){
348   // line segment
349     static Int_t hitindex0;
350     AliITShit *h0,*h1;
351
352     if(hitindex>= fHitsM->GetEntriesFast()) return kFALSE;
353
354     h1 = (AliITShit *) (fHitsM->At(hitindex));
355     if(h1->StatusEntering()){ // if track entering volume, get index for next
356                               // step
357         hitindex0 = hitindex;
358         track = h1->GetTrack();
359         return kFALSE;
360     } // end if StatusEntering()
361     // else stepping
362     h0 = (AliITShit *) (fHitsM->At(hitindex0));
363     de = h1->GetIonization();
364     h0->GetPositionG(a,c,e);
365     h1->GetPositionG(b,d,f);
366     b = b - a;
367     d = d - c;
368     f = f - e;
369     hitindex0 = hitindex;
370     track = h1->GetTrack();
371     return kTRUE;
372 }
373 //___________________________________________________________________________
374 void AliITSmodule::MedianHitL(Int_t index, 
375                              AliITShit *itsHit1, 
376                              AliITShit *itsHit2, 
377                              Float_t &xMl, Float_t &yMl, Float_t &zMl){
378   // median hit
379    Float_t x1l,y1l,z1l;
380    Float_t x2l,y2l,z2l;
381
382    itsHit1->GetPositionL(x1l,y1l,z1l);
383    itsHit2->GetPositionL(x2l,y2l,z2l);
384
385    yMl = 0.0;
386    // Modified by N.Carrer. In very rare occasions the track may be just
387    // tangent to the module. Therefore the entrance and exit points have the
388    // same y.
389    if( (y2l-y1l) != 0.0 ) {
390      xMl = (-y1l / (y2l-y1l))*(x2l-x1l) + x1l;
391      zMl = (-y1l / (y2l-y1l))*(z2l-z1l) + z1l;       
392    } else {
393      xMl = 0.5*(x1l+x2l);
394      zMl = 0.5*(z1l+z2l);
395    }
396 }
397 //___________________________________________________________________________
398 void AliITSmodule::MedianHit(Int_t index,
399                              Float_t xg,Float_t yg,Float_t zg,
400                              Int_t status,
401                              Float_t &xMg,Float_t &yMg,Float_t &zMg,
402                              Int_t &flag){
403   // median hit
404    static Float_t x1,y1,z1;
405
406    if (status == 66){ // entering
407        x1 = xg;
408        y1 = yg;
409        z1 = zg;
410        flag = 1;
411    } // end if
412    if (status == 68){ // exiting
413        MedianHitG(index,x1,y1,z1,xg,yg,zg,xMg,yMg,zMg);
414        flag = 0;
415    } // end if
416    if ((status != 66) && (status != 68)) flag = 1;
417 }
418 //___________________________________________________________________________
419 void AliITSmodule::GetID(Int_t &lay,Int_t &lad,Int_t &det){
420   // get ID
421         fITS->GetITSgeom()->GetModuleId(fIndex,lay,lad,det);
422         return ;
423 }
424