1 |
mocchiut |
1.1 |
/** |
2 |
|
|
* \file TrkLevel2.cpp |
3 |
|
|
* \author Elena Vannuccini |
4 |
|
|
*/ |
5 |
|
|
#include <TrkLevel2.h> |
6 |
|
|
#include <iostream> |
7 |
|
|
using namespace std; |
8 |
|
|
//...................................... |
9 |
|
|
// F77 routines |
10 |
|
|
//...................................... |
11 |
|
|
extern "C" { |
12 |
|
|
void dotrack_(int*, double*, double*, double*, double*, int*); |
13 |
pam-fi |
1.2 |
void dotrack2_(int*, double*, double*, double*, double*,double*, double*, double*,int*); |
14 |
mocchiut |
1.1 |
int readb_(const char*); |
15 |
|
|
} |
16 |
|
|
//-------------------------------------- |
17 |
|
|
// |
18 |
|
|
// |
19 |
|
|
//-------------------------------------- |
20 |
|
|
TrkTrack::TrkTrack(){ |
21 |
|
|
image = 0; |
22 |
|
|
chi2 = 0; |
23 |
|
|
for(int it1=0;it1<5;it1++){ |
24 |
|
|
al[it1] = 0; |
25 |
|
|
for(int it2=0;it2<5;it2++) |
26 |
|
|
coval[it1][it2] = 0; |
27 |
|
|
}; |
28 |
|
|
for(int ip=0;ip<6;ip++){ |
29 |
|
|
xgood[ip] = 0; |
30 |
|
|
ygood[ip] = 0; |
31 |
|
|
xm[ip] = 0; |
32 |
|
|
ym[ip] = 0; |
33 |
|
|
zm[ip] = 0; |
34 |
|
|
resx[ip] = 0; |
35 |
|
|
resy[ip] = 0; |
36 |
|
|
xv[ip] = 0; |
37 |
|
|
yv[ip] = 0; |
38 |
|
|
zv[ip] = 0; |
39 |
|
|
axv[ip] = 0; |
40 |
|
|
ayv[ip] = 0; |
41 |
|
|
dedx_x[ip] = 0; |
42 |
|
|
dedx_y[ip] = 0; |
43 |
|
|
}; |
44 |
|
|
}; |
45 |
|
|
//-------------------------------------- |
46 |
|
|
// |
47 |
|
|
// |
48 |
|
|
//-------------------------------------- |
49 |
|
|
TrkTrack::TrkTrack(const TrkTrack& t){ |
50 |
|
|
image = t.image; |
51 |
|
|
chi2 = t.chi2; |
52 |
|
|
for(int it1=0;it1<5;it1++){ |
53 |
|
|
al[it1] = t.al[it1]; |
54 |
|
|
for(int it2=0;it2<5;it2++) |
55 |
|
|
coval[it1][it2] = t.coval[it1][it2]; |
56 |
|
|
}; |
57 |
|
|
for(int ip=0;ip<6;ip++){ |
58 |
|
|
xgood[ip] = t.xgood[ip]; |
59 |
|
|
ygood[ip] = t.ygood[ip]; |
60 |
|
|
xm[ip] = t.xm[ip]; |
61 |
|
|
ym[ip] = t.ym[ip]; |
62 |
|
|
zm[ip] = t.zm[ip]; |
63 |
|
|
resx[ip] = t.resx[ip]; |
64 |
|
|
resy[ip] = t.resy[ip]; |
65 |
|
|
xv[ip] = t.xv[ip]; |
66 |
|
|
yv[ip] = t.yv[ip]; |
67 |
|
|
zv[ip] = t.zv[ip]; |
68 |
|
|
axv[ip] = t.axv[ip]; |
69 |
|
|
ayv[ip] = t.ayv[ip]; |
70 |
|
|
dedx_x[ip] = t.dedx_x[ip]; |
71 |
|
|
dedx_y[ip] = t.dedx_y[ip]; |
72 |
|
|
}; |
73 |
|
|
}; |
74 |
|
|
//-------------------------------------- |
75 |
|
|
// |
76 |
|
|
// |
77 |
|
|
//-------------------------------------- |
78 |
|
|
/** |
79 |
|
|
* Evaluates the trajectory in the apparatus associated to the track. |
80 |
|
|
* It integrates the equations of motion in the magnetic field. The magnetic field should be previously loaded ( by calling TrkLevel2::LoadField() ), otherwise an error message is returned. |
81 |
|
|
* @param t pointer to an object of the class Trajectory, |
82 |
|
|
* which z coordinates should be previously initialized by calling the proper constructor ( Trajectory::Trajectory(int n, float* zin) ). |
83 |
|
|
* @return error flag. |
84 |
|
|
*/ |
85 |
|
|
int TrkTrack::DoTrack(Trajectory* t){ |
86 |
|
|
|
87 |
|
|
double *dxout = new double[t->npoint]; |
88 |
|
|
double *dyout = new double[t->npoint]; |
89 |
|
|
double *dzin = new double[t->npoint]; |
90 |
|
|
double dal[5]; |
91 |
|
|
|
92 |
|
|
int ifail = 0; |
93 |
|
|
|
94 |
|
|
for (int i=0; i<5; i++) dal[i] = (double)al[i]; |
95 |
|
|
for (int i=0; i<t->npoint; i++) dzin[i] = (double)t->z[i]; |
96 |
|
|
|
97 |
|
|
dotrack_(&(t->npoint),dzin,dxout,dyout,dal,&ifail); |
98 |
|
|
|
99 |
|
|
for (int i=0; i<t->npoint; i++){ |
100 |
|
|
t->x[i] = (float)*dxout++; |
101 |
|
|
t->y[i] = (float)*dyout++; |
102 |
|
|
} |
103 |
|
|
|
104 |
|
|
// delete [] dxout; |
105 |
|
|
// delete [] dyout; |
106 |
|
|
// delete [] dzin; |
107 |
|
|
|
108 |
|
|
return ifail; |
109 |
|
|
}; |
110 |
|
|
//-------------------------------------- |
111 |
|
|
// |
112 |
|
|
// |
113 |
|
|
//-------------------------------------- |
114 |
pam-fi |
1.2 |
/** |
115 |
|
|
* Evaluates the trajectory in the apparatus associated to the track. |
116 |
|
|
* It integrates the equations of motion in the magnetic field. The magnetic field should be previously loaded ( by calling TrkLevel2::LoadField() ), otherwise an error message is returned. |
117 |
|
|
* @param t pointer to an object of the class Trajectory, |
118 |
|
|
* which z coordinates should be previously initialized by calling the proper constructor ( Trajectory::Trajectory(int n, float* zin) ). |
119 |
|
|
* @return error flag. |
120 |
|
|
*/ |
121 |
|
|
int TrkTrack::DoTrack2(Trajectory* t){ |
122 |
|
|
|
123 |
|
|
double *dxout = new double[t->npoint]; |
124 |
|
|
double *dyout = new double[t->npoint]; |
125 |
|
|
double *dthxout = new double[t->npoint]; |
126 |
|
|
double *dthyout = new double[t->npoint]; |
127 |
|
|
double *dtlout = new double[t->npoint]; |
128 |
|
|
double *dzin = new double[t->npoint]; |
129 |
|
|
double dal[5]; |
130 |
|
|
|
131 |
|
|
int ifail = 0; |
132 |
|
|
|
133 |
|
|
for (int i=0; i<5; i++) dal[i] = (double)al[i]; |
134 |
|
|
for (int i=0; i<t->npoint; i++) dzin[i] = (double)t->z[i]; |
135 |
|
|
|
136 |
|
|
dotrack2_(&(t->npoint),dzin,dxout,dyout,dthxout,dthyout,dtlout,dal,&ifail); |
137 |
|
|
|
138 |
|
|
for (int i=0; i<t->npoint; i++){ |
139 |
|
|
t->x[i] = (float)*dxout++; |
140 |
|
|
t->y[i] = (float)*dyout++; |
141 |
|
|
t->thx[i] = (float)*dthxout++; |
142 |
|
|
t->thy[i] = (float)*dthyout++; |
143 |
|
|
t->tl[i] = (float)*dtlout++; |
144 |
|
|
} |
145 |
|
|
|
146 |
|
|
// delete [] dxout; |
147 |
|
|
// delete [] dyout; |
148 |
|
|
// delete [] dzin; |
149 |
|
|
|
150 |
|
|
return ifail; |
151 |
|
|
}; |
152 |
|
|
//-------------------------------------- |
153 |
|
|
// |
154 |
|
|
// |
155 |
|
|
//-------------------------------------- |
156 |
mocchiut |
1.1 |
//float TrkTrack::BdL(){ |
157 |
|
|
//}; |
158 |
|
|
//-------------------------------------- |
159 |
|
|
// |
160 |
|
|
// |
161 |
|
|
//-------------------------------------- |
162 |
|
|
Float_t TrkTrack::GetRigidity(){ |
163 |
|
|
Float_t rig=0; |
164 |
|
|
if(chi2>0)rig=1./al[4]; |
165 |
|
|
if(rig<0) rig=-rig; |
166 |
|
|
return rig; |
167 |
|
|
}; |
168 |
|
|
// |
169 |
|
|
Float_t TrkTrack::GetDeflection(){ |
170 |
|
|
Float_t def=0; |
171 |
|
|
if(chi2>0)def=al[4]; |
172 |
|
|
return def; |
173 |
|
|
}; |
174 |
|
|
// |
175 |
|
|
Float_t TrkTrack::GetDEDX(){ |
176 |
|
|
Float_t dedx=0; |
177 |
|
|
for(Int_t i=0; i<6; i++)dedx+=dedx_x[i]*xgood[i]+dedx_y[i]*ygood[i]; |
178 |
|
|
dedx = dedx/(this->GetNX()+this->GetNY()); |
179 |
|
|
return dedx; |
180 |
|
|
}; |
181 |
|
|
|
182 |
|
|
//-------------------------------------- |
183 |
|
|
// |
184 |
|
|
// |
185 |
|
|
//-------------------------------------- |
186 |
|
|
void TrkTrack::Dump(){ |
187 |
|
|
cout << endl << "========== Track " ; |
188 |
|
|
cout << endl << "al : "; for(int i=0; i<5; i++)cout << al[i] << " "; |
189 |
|
|
cout << endl << "chi^2 : "<< chi2; |
190 |
|
|
cout << endl << "xgood : "; for(int i=0; i<6; i++)cout << xgood[i] ; |
191 |
|
|
cout << endl << "ygood : "; for(int i=0; i<6; i++)cout << ygood[i] ; |
192 |
|
|
cout << endl << "xm : "; for(int i=0; i<6; i++)cout << xm[i] << " "; |
193 |
|
|
cout << endl << "ym : "; for(int i=0; i<6; i++)cout << ym[i] << " "; |
194 |
|
|
cout << endl << "zm : "; for(int i=0; i<6; i++)cout << zm[i] << " "; |
195 |
|
|
cout << endl << "dedx_x : "; for(int i=0; i<6; i++)cout << dedx_x[i] << " "; |
196 |
|
|
cout << endl << "dedx_y : "; for(int i=0; i<6; i++)cout << dedx_y[i] << " "; |
197 |
|
|
} |
198 |
|
|
//-------------------------------------- |
199 |
|
|
// |
200 |
|
|
// |
201 |
|
|
//-------------------------------------- |
202 |
|
|
TrkSinglet::TrkSinglet(){ |
203 |
|
|
plane = 0; |
204 |
|
|
coord[0] = 0; |
205 |
|
|
coord[1] = 0; |
206 |
|
|
sgnl = 0; |
207 |
|
|
}; |
208 |
|
|
//-------------------------------------- |
209 |
|
|
// |
210 |
|
|
// |
211 |
|
|
//-------------------------------------- |
212 |
|
|
TrkSinglet::TrkSinglet(const TrkSinglet& s){ |
213 |
|
|
plane = s.plane; |
214 |
|
|
coord[0] = s.coord[0]; |
215 |
|
|
coord[1] = s.coord[1]; |
216 |
|
|
sgnl = s.sgnl; |
217 |
|
|
}; |
218 |
|
|
//-------------------------------------- |
219 |
|
|
// |
220 |
|
|
// |
221 |
|
|
//-------------------------------------- |
222 |
|
|
void TrkSinglet::Dump(){ |
223 |
|
|
int i=0; |
224 |
|
|
cout << endl << "========== Singlet " ; |
225 |
|
|
cout << endl << "plane : " << plane; |
226 |
|
|
cout << endl << "coord[2] : "; while( i<2 && cout << coord[i] << " ") i++; |
227 |
|
|
cout << endl << "sgnl : " << sgnl; |
228 |
|
|
} |
229 |
|
|
//-------------------------------------- |
230 |
|
|
// |
231 |
|
|
// |
232 |
|
|
//-------------------------------------- |
233 |
|
|
TrkLevel2::TrkLevel2(){ |
234 |
|
|
good2 = -1; |
235 |
|
|
for(Int_t i=0; i<12 ; i++){ |
236 |
|
|
crc[i] = -1; |
237 |
|
|
}; |
238 |
|
|
Track = new TClonesArray("TrkTrack"); |
239 |
|
|
SingletX = new TClonesArray("TrkSinglet"); |
240 |
|
|
SingletY = new TClonesArray("TrkSinglet"); |
241 |
|
|
} |
242 |
|
|
//-------------------------------------- |
243 |
|
|
// |
244 |
|
|
// |
245 |
|
|
//-------------------------------------- |
246 |
|
|
void TrkLevel2::Dump(){ |
247 |
|
|
TClonesArray &t = *Track; |
248 |
|
|
TClonesArray &sx = *SingletX; |
249 |
|
|
TClonesArray &sy = *SingletY; |
250 |
|
|
|
251 |
|
|
cout << endl << endl << "=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-"; |
252 |
|
|
cout << endl << "good2 : " << good2; |
253 |
|
|
cout << endl << "crc : "; for(int i=0; i<12; i++) cout << crc[i]; |
254 |
|
|
cout << endl << "ntrk() : " << this->ntrk() ; |
255 |
|
|
cout << endl << "nclsx() : " << this->nclsx(); |
256 |
|
|
cout << endl << "nclsy() : " << this->nclsy(); |
257 |
|
|
for(int i=0; i<this->ntrk(); i++) ((TrkTrack *)t[i])->Dump(); |
258 |
|
|
for(int i=0; i<this->nclsx(); i++) ((TrkSinglet *)sx[i])->Dump(); |
259 |
|
|
for(int i=0; i<this->nclsy(); i++) ((TrkSinglet *)sy[i])->Dump(); |
260 |
|
|
} |
261 |
|
|
//-------------------------------------- |
262 |
|
|
// |
263 |
|
|
// |
264 |
|
|
//-------------------------------------- |
265 |
|
|
/** |
266 |
|
|
* Fills a TrkLevel2 object with values from a struct cTrkLevel2 (to get data from F77 common). |
267 |
|
|
*/ |
268 |
|
|
void TrkLevel2::FillCommonVar(cTrkLevel2 *l2){ |
269 |
|
|
// temporary objects: |
270 |
|
|
TrkSinglet* t_singlet = new TrkSinglet(); |
271 |
|
|
TrkTrack* t_track = new TrkTrack(); |
272 |
|
|
// general variables |
273 |
|
|
good2 = l2->good2; |
274 |
|
|
for(Int_t i=0; i<12 ; i++){ |
275 |
|
|
crc[i] = l2->crc[i]; |
276 |
|
|
}; |
277 |
|
|
// *** TRACKS *** |
278 |
|
|
TClonesArray &t = *Track; |
279 |
|
|
for(int i=0; i<l2->ntrk; i++){ |
280 |
|
|
t_track->image = l2->image[i]-1; |
281 |
|
|
t_track->chi2 = l2->chi2_nt[i]; |
282 |
|
|
for(int it1=0;it1<5;it1++){ |
283 |
|
|
t_track->al[it1] = l2->al_nt[i][it1]; |
284 |
|
|
for(int it2=0;it2<5;it2++) |
285 |
|
|
t_track->coval[it1][it2] = l2->coval[i][it2][it1]; |
286 |
|
|
}; |
287 |
|
|
for(int ip=0;ip<6;ip++){ |
288 |
|
|
t_track->xgood[ip] = l2->xgood_nt[i][ip]; |
289 |
|
|
t_track->ygood[ip] = l2->ygood_nt[i][ip]; |
290 |
|
|
t_track->xm[ip] = l2->xm_nt[i][ip]; |
291 |
|
|
t_track->ym[ip] = l2->ym_nt[i][ip]; |
292 |
|
|
t_track->zm[ip] = l2->zm_nt[i][ip]; |
293 |
|
|
t_track->resx[ip] = l2->resx_nt[i][ip]; |
294 |
|
|
t_track->resy[ip] = l2->resy_nt[i][ip]; |
295 |
|
|
t_track->xv[ip] = l2->xv_nt[i][ip]; |
296 |
|
|
t_track->yv[ip] = l2->yv_nt[i][ip]; |
297 |
|
|
t_track->zv[ip] = l2->zv_nt[i][ip]; |
298 |
|
|
t_track->axv[ip] = l2->axv_nt[i][ip]; |
299 |
|
|
t_track->ayv[ip] = l2->ayv_nt[i][ip]; |
300 |
|
|
t_track->dedx_x[ip] = l2->dedx_x[i][ip]; |
301 |
|
|
t_track->dedx_y[ip] = l2->dedx_y[i][ip]; |
302 |
|
|
}; |
303 |
|
|
new(t[i]) TrkTrack(*t_track); |
304 |
|
|
t_track->Clear(); |
305 |
|
|
}; |
306 |
|
|
// *** SINGLETS *** |
307 |
|
|
TClonesArray &sx = *SingletX; |
308 |
|
|
for(int i=0; i<l2->nclsx; i++){ |
309 |
|
|
t_singlet->plane = l2->planex[i]; |
310 |
|
|
t_singlet->coord[0] = l2->xs[i][0]; |
311 |
|
|
t_singlet->coord[1] = l2->xs[i][1]; |
312 |
|
|
t_singlet->sgnl = l2->signlxs[i]; |
313 |
|
|
new(sx[i]) TrkSinglet(*t_singlet); |
314 |
|
|
t_singlet->Clear(); |
315 |
|
|
} |
316 |
|
|
TClonesArray &sy = *SingletY; |
317 |
|
|
for(int i=0; i<l2->nclsy; i++){ |
318 |
|
|
t_singlet->plane = l2->planey[i]; |
319 |
|
|
t_singlet->coord[0] = l2->ys[i][0]; |
320 |
|
|
t_singlet->coord[1] = l2->ys[i][1]; |
321 |
|
|
t_singlet->sgnl = l2->signlys[i]; |
322 |
|
|
new(sy[i]) TrkSinglet(*t_singlet); |
323 |
|
|
t_singlet->Clear(); |
324 |
|
|
}; |
325 |
|
|
} |
326 |
|
|
//-------------------------------------- |
327 |
|
|
// |
328 |
|
|
// |
329 |
|
|
//-------------------------------------- |
330 |
|
|
void TrkLevel2::Clear(){ |
331 |
|
|
good2 = -1; |
332 |
|
|
for(Int_t i=0; i<12 ; i++){ |
333 |
|
|
crc[i] = -1; |
334 |
|
|
}; |
335 |
|
|
Track->RemoveAll(); |
336 |
|
|
SingletX->RemoveAll(); |
337 |
|
|
SingletY->RemoveAll(); |
338 |
|
|
} |
339 |
|
|
//-------------------------------------- |
340 |
|
|
// |
341 |
|
|
// |
342 |
|
|
//-------------------------------------- |
343 |
|
|
/** |
344 |
|
|
* Sort physical tracks and stores them in a TObjectArray, ordering by increasing chi**2 value (in case of track image, it selects the one with lower chi**2). The total number of physical tracks is given by GetNTracks() and the it-th physical track can be retrieved by means of the method GetTrack(int it). |
345 |
|
|
* This method is overridden by PamLevel2::GetTracks(), where calorimeter and TOF information is used. |
346 |
|
|
*/ |
347 |
|
|
TClonesArray *TrkLevel2::GetTracks(){ |
348 |
|
|
|
349 |
|
|
TClonesArray *sorted = new TClonesArray("TrkTrack"); |
350 |
|
|
TClonesArray &t = *Track; |
351 |
|
|
TClonesArray &ts = *sorted; |
352 |
|
|
int N=this->ntrk(); |
353 |
|
|
vector<int> m(N); for(int i=0; i<N; i++)m[i]=1; |
354 |
|
|
|
355 |
|
|
int indo=0; |
356 |
|
|
int indi=0; |
357 |
|
|
while(N != 0){ |
358 |
|
|
float chi2ref=1000000; |
359 |
|
|
for(int i=0; i<this->ntrk(); i++){ |
360 |
|
|
if(((TrkTrack *)t[i])->chi2 < chi2ref && m[i]==1){ |
361 |
|
|
chi2ref = ((TrkTrack *)t[i])->chi2; |
362 |
|
|
indi = i; |
363 |
|
|
} |
364 |
|
|
} |
365 |
|
|
if( ((TrkTrack *)t[indi])->image != -1 ){ |
366 |
|
|
m[((TrkTrack *)t[indi])->image] = 0; |
367 |
|
|
N--; |
368 |
|
|
} |
369 |
|
|
new(ts[indo]) TrkTrack(*(TrkTrack*)t[indi]); |
370 |
|
|
m[indi] = 0; |
371 |
|
|
N--; |
372 |
|
|
indo++; |
373 |
|
|
} |
374 |
|
|
return sorted; |
375 |
|
|
} |
376 |
|
|
//-------------------------------------- |
377 |
|
|
// |
378 |
|
|
// |
379 |
|
|
//-------------------------------------- |
380 |
|
|
/** |
381 |
|
|
* Retrieves the is-th stored track. |
382 |
|
|
* @param it Track number, ranging from 0 to ntrk(). |
383 |
|
|
* Fitted tracks ( images included ) are stored in a TObjectArray ( TrkLevel2::Track ) in the same order they are returned by the F77 fitting routine. |
384 |
|
|
*/ |
385 |
|
|
TrkTrack *TrkLevel2::GetStoredTrack(int is){ |
386 |
|
|
|
387 |
|
|
if(is >= this->ntrk()){ |
388 |
|
|
cout << "** TrkLevel2 ** Track "<< is << "doen not exits! " << endl; |
389 |
|
|
cout << " Stored tracks ntrk() = "<< this->ntrk() << endl; |
390 |
|
|
return 0; |
391 |
|
|
} |
392 |
|
|
TClonesArray &t = *(Track); |
393 |
|
|
TrkTrack *track = (TrkTrack*)t[is]; |
394 |
|
|
return track; |
395 |
|
|
} |
396 |
|
|
//-------------------------------------- |
397 |
|
|
// |
398 |
|
|
// |
399 |
|
|
//-------------------------------------- |
400 |
|
|
/** |
401 |
|
|
* Retrieves the it-th "physical" track, sorted by the method GetNTracks(). |
402 |
|
|
* @param it Track number, ranging from 0 to GetNTracks(). |
403 |
|
|
*/ |
404 |
|
|
TrkTrack *TrkLevel2::GetTrack(int it){ |
405 |
|
|
|
406 |
|
|
if(it >= this->GetNTracks()){ |
407 |
|
|
cout << "** TrkLevel2 ** Track "<< it << "does not exits! " << endl; |
408 |
|
|
cout << " Physical tracks GetNTracks() = "<< this->ntrk() << endl; |
409 |
|
|
return 0; |
410 |
|
|
} |
411 |
|
|
TrkTrack *track = (TrkTrack*)(*(this->GetTracks()))[it]; |
412 |
|
|
return track; |
413 |
|
|
} |
414 |
|
|
//-------------------------------------- |
415 |
|
|
// |
416 |
|
|
// |
417 |
|
|
//-------------------------------------- |
418 |
|
|
/** |
419 |
|
|
* Retrieves (if present) the image of the it-th "physical" track, sorted by the method GetNTracks(). |
420 |
|
|
* @param it Track number, ranging from 0 to GetNTracks(). |
421 |
|
|
*/ |
422 |
|
|
TrkTrack *TrkLevel2::GetTrackImage(int it){ |
423 |
|
|
|
424 |
|
|
if(it >= this->GetNTracks()){ |
425 |
|
|
cout << "** TrkLevel2 ** Track "<< it << "does not exits! " << endl; |
426 |
|
|
cout << " Physical tracks GetNTracks() = "<< this->ntrk() << endl; |
427 |
|
|
return 0; |
428 |
|
|
} |
429 |
|
|
TrkTrack *track = (TrkTrack*)(*(this->GetTracks()))[it]; |
430 |
|
|
if(!track->HasImage()){ |
431 |
|
|
cout << "** TrkLevel2 ** Track "<< it << "does not have image! " << endl; |
432 |
|
|
return 0; |
433 |
|
|
} |
434 |
|
|
TrkTrack *image = (TrkTrack*)(*Track)[track->image]; |
435 |
|
|
return image; |
436 |
|
|
|
437 |
|
|
} |
438 |
|
|
//-------------------------------------- |
439 |
|
|
// |
440 |
|
|
// |
441 |
|
|
//-------------------------------------- |
442 |
|
|
/** |
443 |
|
|
* Loads the magnetic field. |
444 |
|
|
* @param s Path of the magnetic-field files. |
445 |
|
|
*/ |
446 |
|
|
void TrkLevel2::LoadField(TString s){ |
447 |
|
|
readb_(s.Data()); |
448 |
|
|
}; |
449 |
|
|
//-------------------------------------- |
450 |
|
|
// |
451 |
|
|
// |
452 |
|
|
//-------------------------------------- |
453 |
|
|
/** |
454 |
pam-fi |
1.2 |
* Trajectory default constructor. |
455 |
|
|
* (By default is created with z-coordinates inside the tracking volume) |
456 |
|
|
*/ |
457 |
|
|
Trajectory::Trajectory(){ |
458 |
|
|
npoint = 10; |
459 |
|
|
x = new float[npoint]; |
460 |
|
|
y = new float[npoint]; |
461 |
|
|
z = new float[npoint]; |
462 |
|
|
thx = new float[npoint]; |
463 |
|
|
thy = new float[npoint]; |
464 |
|
|
tl = new float[npoint]; |
465 |
|
|
float dz = ((ZTRKUP)-(ZTRKDW))/(npoint-1); |
466 |
|
|
for(int i=0; i<npoint; i++){ |
467 |
|
|
x[i] = 0; |
468 |
|
|
y[i] = 0; |
469 |
|
|
z[i] = (ZTRKUP) - i*dz; |
470 |
|
|
thx[i] = 0; |
471 |
|
|
thy[i] = 0; |
472 |
|
|
tl[i] = 0; |
473 |
|
|
} |
474 |
|
|
} |
475 |
|
|
//-------------------------------------- |
476 |
|
|
// |
477 |
|
|
// |
478 |
|
|
//-------------------------------------- |
479 |
|
|
/** |
480 |
mocchiut |
1.1 |
* Trajectory constructor. |
481 |
pam-fi |
1.2 |
* (By default is created with z-coordinates inside the tracking volume) |
482 |
mocchiut |
1.1 |
* \param n Number of points |
483 |
|
|
*/ |
484 |
|
|
Trajectory::Trajectory(int n){ |
485 |
pam-fi |
1.2 |
if(n<=0){ |
486 |
|
|
cout << "NB! Trajectory must have at least 1 point >>> created with 10 points" << endl; |
487 |
|
|
n=10; |
488 |
|
|
} |
489 |
mocchiut |
1.1 |
npoint = n; |
490 |
|
|
x = new float[npoint]; |
491 |
|
|
y = new float[npoint]; |
492 |
|
|
z = new float[npoint]; |
493 |
pam-fi |
1.2 |
thx = new float[npoint]; |
494 |
|
|
thy = new float[npoint]; |
495 |
|
|
tl = new float[npoint]; |
496 |
|
|
float dz = ((ZTRKUP)-(ZTRKDW))/(npoint-1); |
497 |
mocchiut |
1.1 |
for(int i=0; i<npoint; i++){ |
498 |
pam-fi |
1.2 |
x[i] = 0; |
499 |
mocchiut |
1.1 |
y[i] = 0; |
500 |
pam-fi |
1.2 |
z[i] = (ZTRKUP) - i*dz; |
501 |
|
|
thx[i] = 0; |
502 |
|
|
thy[i] = 0; |
503 |
|
|
tl[i] = 0; |
504 |
mocchiut |
1.1 |
} |
505 |
|
|
} |
506 |
|
|
//-------------------------------------- |
507 |
|
|
// |
508 |
|
|
// |
509 |
|
|
//-------------------------------------- |
510 |
|
|
/** |
511 |
|
|
* Trajectory constructor. |
512 |
|
|
* \param n Number of points |
513 |
|
|
* \param pz Pointer to float array, defining z coordinates |
514 |
|
|
*/ |
515 |
|
|
Trajectory::Trajectory(int n, float* zin){ |
516 |
pam-fi |
1.2 |
npoint = 10; |
517 |
|
|
if(n>0)npoint = n; |
518 |
mocchiut |
1.1 |
x = new float[npoint]; |
519 |
|
|
y = new float[npoint]; |
520 |
|
|
z = new float[npoint]; |
521 |
pam-fi |
1.2 |
thx = new float[npoint]; |
522 |
|
|
thy = new float[npoint]; |
523 |
|
|
tl = new float[npoint]; |
524 |
|
|
int i=0; |
525 |
|
|
do{ |
526 |
mocchiut |
1.1 |
x[i] = 0; |
527 |
|
|
y[i] = 0; |
528 |
|
|
z[i] = zin[i]; |
529 |
pam-fi |
1.2 |
thx[i] = 0; |
530 |
|
|
thy[i] = 0; |
531 |
|
|
tl[i] = 0; |
532 |
|
|
i++; |
533 |
|
|
}while(zin[i-1] > zin[i] && i < npoint); |
534 |
|
|
npoint=i; |
535 |
|
|
if(npoint != n)cout << "NB! Trajectory created with "<<npoint<<" points"<<endl; |
536 |
mocchiut |
1.1 |
} |
537 |
|
|
//-------------------------------------- |
538 |
|
|
// |
539 |
|
|
// |
540 |
|
|
//-------------------------------------- |
541 |
|
|
/** |
542 |
|
|
* Dump the trajectory coordinates. |
543 |
|
|
*/ |
544 |
|
|
void Trajectory::Dump(){ |
545 |
|
|
cout <<endl<< "Trajectory ========== "<<endl; |
546 |
|
|
for (int i=0; i<npoint; i++){ |
547 |
pam-fi |
1.2 |
cout << i <<" >> " << x[i] <<" "<< y[i] <<" "<< z[i] ; |
548 |
|
|
cout <<" -- " << thx[i] <<" "<< thy[i] ; |
549 |
|
|
cout <<" -- " << tl[i] << endl; |
550 |
mocchiut |
1.1 |
}; |
551 |
|
|
} |
552 |
pam-fi |
1.2 |
//-------------------------------------- |
553 |
|
|
// |
554 |
|
|
// |
555 |
|
|
//-------------------------------------- |
556 |
|
|
/** |
557 |
|
|
* Get trajectory length between two points |
558 |
|
|
* @param ifirst first point (default 0) |
559 |
|
|
* @param ilast last point (default npoint) |
560 |
|
|
*/ |
561 |
|
|
float Trajectory::GetLength(int ifirst, int ilast){ |
562 |
|
|
if( ifirst<0 ) ifirst = 0; |
563 |
|
|
if( ilast>=npoint) ilast = npoint-1; |
564 |
|
|
float l=0; |
565 |
|
|
for(int i=ifirst;i<=ilast;i++){ |
566 |
|
|
l=l+tl[i]; |
567 |
|
|
}; |
568 |
|
|
if(z[ilast] > ZINI)l=l-tl[ilast]; |
569 |
|
|
if(z[ifirst] < ZINI) l=l-tl[ifirst]; |
570 |
|
|
|
571 |
|
|
return l; |
572 |
mocchiut |
1.1 |
|
573 |
pam-fi |
1.2 |
} |
574 |
mocchiut |
1.1 |
ClassImp(TrkLevel2); |
575 |
|
|
ClassImp(TrkSinglet); |
576 |
|
|
ClassImp(TrkTrack); |
577 |
|
|
ClassImp(Trajectory); |