Line data Source code
1 : #include "TMath.h"
2 : #include <TGeoGlobalMagField.h>
3 :
4 : #include "AliMFTCATrack.h"
5 : #include "AliMFTCACell.h"
6 : #include "AliMFTTrackExtrap.h"
7 :
8 : #include "AliLog.h"
9 :
10 12 : ClassImp(AliMFTCATrack)
11 :
12 :
13 :
14 : //___________________________________________________________________________
15 :
16 : AliMFTCATrack::AliMFTCATrack() :
17 0 : TObject(),
18 0 : fGID(-1),
19 0 : fNcells(0),
20 0 : fStartLayer(-1),
21 0 : fCellGIDarray(),
22 0 : fMCflag(0),
23 0 : fVertX(0.),
24 0 : fVertY(0.),
25 0 : fVertZ(0.),
26 0 : fTheta(0.),
27 0 : fPhi(0.),
28 0 : fChiSqX(0.),
29 0 : fChiSqY(0.),
30 0 : fMCindex(-1),
31 0 : fChargeSign(0.),
32 0 : fPt(0.)
33 0 : {
34 :
35 0 : fCells = new TClonesArray("AliMFTCACell", fNDetMax);
36 :
37 0 : }
38 :
39 : //___________________________________________________________________________
40 :
41 : AliMFTCATrack::AliMFTCATrack(const AliMFTCATrack &track) :
42 0 : TObject(track),
43 0 : fGID(track.fGID),
44 0 : fNcells(track.fNcells),
45 0 : fStartLayer(track.fStartLayer),
46 0 : fMCflag(track.fMCflag),
47 0 : fVertX(track.fVertX),
48 0 : fVertY(track.fVertY),
49 0 : fVertZ(track.fVertZ),
50 0 : fTheta(track.fTheta),
51 0 : fPhi(track.fPhi),
52 0 : fChiSqX(track.fChiSqX),
53 0 : fChiSqY(track.fChiSqY),
54 0 : fMCindex(track.fMCindex),
55 0 : fChargeSign(track.fChargeSign),
56 0 : fPt(track.fPt)
57 0 : {
58 :
59 : // copy constructor
60 :
61 0 : fCells = new TClonesArray("AliMFTCACell", fNDetMax);
62 :
63 0 : for (Int_t icell = 0; icell < track.fNcells; icell++)
64 0 : fCellGIDarray[icell] = track.fCellGIDarray[icell];
65 :
66 0 : }
67 :
68 : //___________________________________________________________________________
69 :
70 : AliMFTCATrack& AliMFTCATrack::operator=(const AliMFTCATrack& track)
71 : {
72 :
73 : // assignment operator
74 :
75 : // check assignement to self
76 0 : if (this == &track) return *this;
77 :
78 0 : TObject::operator=(track);
79 :
80 0 : fGID = track.fGID;
81 0 : fNcells = track.fNcells;
82 0 : fStartLayer = track.fStartLayer;
83 0 : fCells = new TClonesArray("AliMFTCACell", fNDetMax);
84 0 : for (Int_t icell = 0; icell < track.fNcells; icell++)
85 0 : fCellGIDarray[icell] = track.fCellGIDarray[icell];
86 0 : fMCflag = track.fMCflag;
87 0 : fVertX = track.fVertX;
88 0 : fVertY = track.fVertY;
89 0 : fVertZ = track.fVertZ;
90 0 : fTheta = track.fTheta;
91 0 : fPhi = track.fPhi;
92 0 : fChiSqX = track.fChiSqX;
93 0 : fChiSqY = track.fChiSqY;
94 0 : fMCindex = track.fMCindex;
95 0 : fChargeSign = track.fChargeSign;
96 0 : fPt = track.fPt;
97 :
98 0 : }
99 :
100 : //___________________________________________________________________________
101 :
102 : void AliMFTCATrack::Clear(Option_t *) {
103 :
104 0 : }
105 : //___________________________________________________________________________
106 :
107 : /// Estimate the charge sign
108 : void AliMFTCATrack::EvalSignedPt(){
109 : const Int_t nMaxh = 100;
110 0 : Double_t xTr[nMaxh], yTr[nMaxh], zTr[nMaxh];
111 0 : Double_t r[nMaxh], u[nMaxh] , v[nMaxh];
112 :
113 : AliMFTCACell *celltr;
114 : Int_t nptr = 0;
115 :
116 0 : for (Int_t iCell = 0; iCell < GetNcells(); iCell++) {
117 0 : celltr = (AliMFTCACell*)fCells->At(iCell);
118 : // extract hit x,y,z
119 0 : if (nptr == 0) {
120 0 : xTr[nptr] = celltr->GetHit2()[0];
121 0 : yTr[nptr] = celltr->GetHit2()[1];
122 0 : zTr[nptr] = celltr->GetHit2()[2];
123 0 : nptr++;
124 0 : xTr[nptr] = celltr->GetHit1()[0];
125 0 : yTr[nptr] = celltr->GetHit1()[1];
126 0 : zTr[nptr] = celltr->GetHit1()[2];
127 0 : nptr++;
128 0 : } else {
129 0 : xTr[nptr] = celltr->GetHit1()[0];
130 0 : yTr[nptr] = celltr->GetHit1()[1];
131 0 : zTr[nptr] = celltr->GetHit1()[2];
132 0 : nptr++;
133 : }
134 : } // END : cells loop
135 0 : for (Int_t i = 0; i< nptr; i++) {
136 0 : AliInfo(Form("x,y,z = %f %f %f",xTr[i],yTr[i],zTr[i]));
137 : }
138 : // Double_t distL = 0., q2=0.;
139 : // Double_t sagitta = AliMFTTrackExtrap::Sagitta(nptr, xTr, yTr, distL,q2);
140 : // AliInfo(Form("sagitta = %f",sagitta));;
141 : // AliInfo(Form("distL = %f",distL));;
142 : //
143 : // AliInfo(Form("q2 = %f => pt = %f ",q2, 0.01/q2/2.*0.3*b[2]/10.));
144 : //
145 : // //fPt = 0.3*TMath::Abs(b[2]/10.)*distL*distL/8./sagitta;
146 : // fPt = 0.3*TMath::Abs(b[2]/10.)*AliMFTTrackExtrap::CircleRegression(nptr, xTr, yTr)*TMath::Sign(1.,sagitta);
147 : // //fPt =0.01/q2/2.*0.3*b[2]/10.;
148 :
149 :
150 :
151 :
152 :
153 0 : for (Int_t iptr = 0; iptr < nptr; iptr++) {
154 :
155 0 : r[iptr] = TMath::Sqrt(yTr[iptr]*yTr[iptr]+xTr[iptr]*xTr[iptr]);
156 0 : u[iptr] = xTr[iptr]/r[iptr]/r[iptr];
157 0 : v[iptr] = yTr[iptr]/r[iptr]/r[iptr];
158 :
159 :
160 0 : AliInfo(Form("u,v,r = %f %f %f ",u[iptr],v[iptr],r[iptr]));
161 :
162 : }
163 0 : Double_t fdump, slopeX_Z, slopeY_Z;
164 0 : AliMFTTrackExtrap::LinearRegression(nptr, xTr,zTr, fdump, slopeX_Z);
165 0 : AliMFTTrackExtrap::LinearRegression(nptr, yTr,zTr, fdump, slopeY_Z);
166 :
167 0 : Double_t phi = TMath::ATan2(slopeY_Z,slopeX_Z);
168 0 : Double_t xS,x0;
169 0 : Double_t chi2 = AliMFTTrackExtrap::LinearRegression(nptr, u,v, x0, xS);
170 0 : AliInfo(Form("chi2 = %f x0, xS %f %f",chi2,x0,xS));;
171 :
172 0 : if (TMath::Abs(x0)<1.e-10) {
173 0 : fPt = 1.e6;
174 0 : return;
175 : }
176 :
177 0 : Double_t rx = -xS/(2.*x0);
178 0 : Double_t ry = 1./(2.*x0);
179 0 : Double_t rr = TMath::Sqrt(rx*rx+ry*ry);
180 0 : AliInfo(Form("Rx,Ry = %f %f",rx,ry));;
181 0 : AliInfo(Form("R = %f",rr));;
182 0 : AliInfo(Form("Phi = %f",phi));;
183 :
184 0 : Double_t zmean = 0.5 * (AliMFTConstants::DefaultPlaneZ(0) + AliMFTConstants::DefaultPlaneZ(9));
185 0 : const Double_t x[3] = {0.,0.,zmean};
186 0 : Double_t b[3] = {0.,0.,0.};
187 0 : TGeoGlobalMagField::Instance()->Field(x,b);
188 0 : AliInfo(Form("Field = %e %e %e",b[0],b[1],b[2]));
189 :
190 0 : fPt = 0.3*TMath::Abs(b[2]/10.)*rr/100.*TMath::Sign(1.,b[2]*rx*(-phi));
191 :
192 :
193 0 : AliInfo(Form("pt = %f",fPt));;
194 :
195 :
196 0 : }
197 :
198 : //___________________________________________________________________________
199 :
200 : void AliMFTCATrack::AddCell(AliMFTCACell *cell) {
201 :
202 0 : if (fNcells >= fNDetMax) { printf("Max number of cells in this track!\n"); return; }
203 0 : new ((*fCells)[fNcells]) AliMFTCACell(*cell);
204 0 : fCellGIDarray[fNcells] = cell->GetGID();
205 0 : fNcells++;
206 :
207 0 : }
208 :
209 : //___________________________________________________________________________
210 :
211 : Double_t AliMFTCATrack::AddCellToChiSq(AliMFTCACell *cell) {
212 :
213 : const Int_t nMaxh = 100;
214 0 : Double_t xTr[nMaxh], yTr[nMaxh], zTr[nMaxh];
215 : Double_t a, ae, b, be, x0, xS, y0, yS, chisqx, chisqy;
216 0 : Double_t xTrErrDet = 0.0025/TMath::Sqrt(12.);
217 0 : Double_t yTrErrDet = 0.0025/TMath::Sqrt(12.);
218 : Double_t xTrErrMS = 0.00055; // estimated at p = 5.5 GeV/c
219 : Double_t yTrErrMS = 0.00055; // estimated at p = 5.5 GeV/c
220 0 : Double_t xTrErr[nMaxh], yTrErr[nMaxh];
221 0 : for (Int_t i = 0; i < nMaxh; i++) {
222 0 : xTrErr[i] = TMath::Sqrt(xTrErrDet*xTrErrDet+xTrErrMS*xTrErrMS);
223 0 : yTrErr[i] = TMath::Sqrt(yTrErrDet*yTrErrDet+yTrErrMS*yTrErrMS);
224 : }
225 : Int_t cellGID, ndof, nptr = 0;
226 :
227 : AliMFTCACell *celltr;
228 :
229 0 : for (Int_t iCell = 0; iCell < GetNcells(); iCell++) {
230 0 : celltr = (AliMFTCACell*)fCells->At(iCell);
231 : // extract hit x,y,z
232 0 : if (nptr == 0) {
233 0 : xTr[nptr] = celltr->GetHit2()[0];
234 0 : yTr[nptr] = celltr->GetHit2()[1];
235 0 : zTr[nptr] = celltr->GetHit2()[2];
236 0 : nptr++;
237 0 : xTr[nptr] = celltr->GetHit1()[0];
238 0 : yTr[nptr] = celltr->GetHit1()[1];
239 0 : zTr[nptr] = celltr->GetHit1()[2];
240 0 : nptr++;
241 0 : } else {
242 0 : xTr[nptr] = celltr->GetHit1()[0];
243 0 : yTr[nptr] = celltr->GetHit1()[1];
244 0 : zTr[nptr] = celltr->GetHit1()[2];
245 0 : nptr++;
246 : }
247 : } // END : cells loop
248 : // the new cell
249 0 : xTr[nptr] = cell->GetHit1()[0];
250 0 : yTr[nptr] = cell->GetHit1()[1];
251 0 : zTr[nptr] = cell->GetHit1()[2];
252 0 : nptr++;
253 : // linear regression
254 : // if (LinFit(nptr,zTr,xTr,xTrErr,a,ae,b,be)) {
255 : // x0 = b; xS = a;
256 : // if (LinFit(nptr,zTr,yTr,yTrErr,a,ae,b,be)) {
257 : // y0 = b; yS = a;
258 : // chisqx = 0.;
259 : // chisqy = 0.;
260 : // for (Int_t iptr = 0; iptr < nptr; iptr++) {
261 : // //printf("%d %f %f %f \n",iptr,xTr[iptr],yTr[iptr],zTr[iptr]);
262 : // chisqx += (xTr[iptr]-(xS*zTr[iptr]+x0))*(xTr[iptr]-(xS*zTr[iptr]+x0))/(xTrErr[iptr]*xTrErr[iptr]);
263 : // chisqy += (yTr[iptr]-(yS*zTr[iptr]+y0))*(yTr[iptr]-(yS*zTr[iptr]+y0))/(yTrErr[iptr]*yTrErr[iptr]);
264 : // }
265 : // }
266 : // }
267 0 : ndof = 2*nptr-4;
268 :
269 0 : return (chisqx+chisqy)/(Double_t)ndof;
270 :
271 0 : }
|
