Script supersceeded by AliForwarddNdetaTask.C and
[u/mrichter/AliRoot.git] / TPC / AliTPCpolyTrack.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 // this class is used for a simple
17 // track fitting using the polynom
18 // bending in both directions are allowed
19 //-----------------------------------------------------------------
20
21
22 #include "AliTPCpolyTrack.h"
23 #include "TMath.h"
24
25 ClassImp(AliTPCpolyTrack)
26
27
28 AliTPCpolyTrack::AliTPCpolyTrack()
29   :TObject(),
30    fA(0), 
31    fB(0), 
32    fC(0), 
33    fD(0), 
34    fE(0), 
35    fF(0), 
36    fMaxX(0), 
37    fMinX(0), 
38    fSumW(0),    
39    fSumX(0),    
40    fSumX2(0),   
41    fSumX3(0),    
42    fSumX4(0),   
43    fSumY(0),    
44    fSumYX(0),   
45    fSumYX2(0),  
46    fSumZ(0),     
47    fSumZX(0),   
48    fSumZX2(0), 
49    fNPoints(0)
50 {
51   Reset();
52 }
53
54 void   AliTPCpolyTrack::Reset()
55 {
56   //
57   // reset track
58   fSumX = fSumX2= fSumX3=fSumX4 = fSumY=fSumYX=fSumYX2=fSumZ=fSumZX=fSumZX2=fSumW =0;
59   fNPoints = 0;
60 }
61
62 void AliTPCpolyTrack::AddPoint(Double_t x, Double_t y, Double_t z,Double_t sy, Double_t sz)
63 {
64   //
65   //
66   if (fNPoints==0){
67     fMaxX = x;
68     fMinX = x;
69   }else{
70     if (x>fMaxX) fMaxX=x;
71     if (x<fMinX) fMinX=x;
72   }
73
74   Double_t x2 = x*x; 
75   Double_t w = 2./(sy+sz);
76   fSumW += w;
77   //
78   fSumX       += x*w;
79   fSumX2      += x2*w;
80   fSumX3      += x2*x*w;
81   fSumX4      += x2*x2*w;
82   //
83   fSumY       +=y*w;
84   fSumYX      +=y*x*w;
85   fSumYX2     +=y*x2*w;
86   //
87   fSumZ       +=z*w;
88   fSumZX      +=z*x*w;
89   fSumZX2     +=z*x2*w;
90   //
91   fX[fNPoints] = x;
92   fY[fNPoints] = y;
93   fZ[fNPoints] = z;  
94   fSY[fNPoints] = sy;
95   fSZ[fNPoints] = sz;  
96
97   fNPoints++;
98   
99 }
100
101 void  AliTPCpolyTrack::UpdateParameters()
102 {
103   //
104   //
105   //Update fit parameters
106   if (fNPoints>4){
107     Fit2(fA,fB,fC);
108     //    Fit2(fSumZ,fSumZX,fSumZX2,fSumX,fSumX2,fSumX3,fSumX4,fNPoints,fD,fE,fF);
109     Fit1(fD,fE,fF);
110   }
111   else
112     {
113       Fit1(fA,fB,fC);
114       Fit1(fD,fE,fF);
115     }
116 }
117
118 void  AliTPCpolyTrack::UpdateParameters(Int_t ny, Int_t nz)
119 {
120   //
121   //
122   //Update fit parameters
123   if (ny>1){
124     Fit2(fA,fB,fC);
125   }
126   else{ 
127     Fit1(fA,fB,fC);
128   }
129   if (nz>1){
130     Fit2(fD,fE,fF);
131   }
132   else{
133     Fit1(fD,fE,fF);
134   }
135 }
136
137
138
139 void  AliTPCpolyTrack::Fit2(
140             Double_t &a, Double_t &b, Double_t &c)
141 {
142   //fit of second order
143   Double_t det = 
144     fSumW* (fSumX2*fSumX4-fSumX3*fSumX3) -
145     fSumX*      (fSumX*fSumX4-fSumX3*fSumX2)+
146     fSumX2*     (fSumX*fSumX3-fSumX2*fSumX2);
147     
148   if (TMath::Abs(det)> 0.000000000000001) {    
149     a = 
150       (fSumY * (fSumX2*fSumX4-fSumX3*fSumX3)-
151        fSumX *(fSumYX*fSumX4-fSumYX2*fSumX3)+
152        fSumX2*(fSumYX*fSumX3-fSumYX2*fSumX2))/det; 
153     b=
154       (fSumW*(fSumYX*fSumX4-fSumX3*fSumYX2)-
155       fSumY*(fSumX*fSumX4-fSumX3*fSumX2)+
156       fSumX2*(fSumX*fSumYX2-fSumYX*fSumX2))/det;
157     c=
158       (fSumW*(fSumX2*fSumYX2-fSumYX*fSumX3)-
159        fSumX*(fSumX*fSumYX2-fSumYX*fSumX2)+
160        fSumY*(fSumX*fSumX3-fSumX2*fSumX2))/det;  
161   }
162 }
163
164 void  AliTPCpolyTrack::Fit1( Double_t &a, Double_t &b, Double_t &c)
165 {
166   //
167   //
168   //
169   Double_t det = fSumW*fSumX2-fSumX*fSumX;
170   if (TMath::Abs(det)> 0.000000000000001) { 
171     b = (fSumW*fSumYX-fSumX*fSumY)/det;
172     a = (fSumX2*fSumY-fSumX*fSumYX)/det;
173     c = 0;
174   }else{
175     a =fSumYX/fSumX;
176     b =0;
177     c =0;
178   }
179
180 }
181
182 void AliTPCpolyTrack::Refit(AliTPCpolyTrack &track, Double_t deltay, Double_t deltaz)
183 {
184   //
185   // refit with cut on distortion
186   //
187   track.Reset();
188   //first refit to temporary
189   AliTPCpolyTrack track0;
190   track0.Reset();
191   for (Int_t i=0;i<fNPoints;i++){
192     Double_t y,z;
193     GetFitPoint(fX[i],y,z);
194     if ( (TMath::Abs(y-fY[i])<deltay)&&(TMath::Abs(z-fZ[i])<deltaz)){
195       track0.AddPoint(fX[i],y,z);
196     }
197   }
198   if (track0.GetN()>2) 
199     track0.UpdateParameters();
200   else 
201     return;
202   //
203   for (Int_t i=0;i<fNPoints;i++){
204     Double_t y,z;
205     track0.GetFitPoint(fX[i],y,z);
206     if ( (TMath::Abs(y-fY[i])<deltay)&&(TMath::Abs(z-fZ[i])<deltaz)){
207       track.AddPoint(fX[i],y,z);
208     }
209   }
210   if (track.GetN()>2) 
211     track.UpdateParameters();
212
213 }
214
215 void AliTPCpolyTrack::Refit(AliTPCpolyTrack &track, Double_t deltay, Double_t deltaz, Int_t nfirst, Int_t ny, Int_t nz)
216 {
217   //
218   // refit with cut on distortion
219   //
220   track.Reset();
221   //first refit to temporary
222   AliTPCpolyTrack track0;
223   track0.Reset();
224   for (Int_t i=0;i<fNPoints;i++){
225     Double_t y,z;
226     GetFitPoint(fX[i],y,z);
227     if ( (TMath::Abs(y-fY[i])<deltay)&&(TMath::Abs(z-fZ[i])<deltaz)){
228       track0.AddPoint(fX[i],y,z);
229     }    
230   }
231   if (track0.GetN()>2){ 
232     if (track0.GetN()>nfirst)
233       track0.UpdateParameters(ny,nz);
234     else 
235       track0.UpdateParameters(1,1);
236   }
237   else 
238     return;
239   //
240   for (Int_t i=0;i<fNPoints;i++){
241     Double_t y,z;
242     track0.GetFitPoint(fX[i],y,z);
243     if ( (TMath::Abs(y-fY[i])<deltay)&&(TMath::Abs(z-fZ[i])<deltaz)){
244       track.AddPoint(fX[i],y,z);
245     }
246   }
247   if (track.GetN()>2) 
248     track.UpdateParameters(ny,nz);
249
250 }
251