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