12 |
extern "C" { |
extern "C" { |
13 |
void dotrack_(int*, double*, double*, double*, double*, int*); |
void dotrack_(int*, double*, double*, double*, double*, int*); |
14 |
void dotrack2_(int*, double*, double*, double*, double*,double*, double*, double*,int*); |
void dotrack2_(int*, double*, double*, double*, double*,double*, double*, double*,int*); |
|
// int readb_(const char*); |
|
|
int readb_(); |
|
15 |
void mini2_(int*,int*,int*); |
void mini2_(int*,int*,int*); |
16 |
void guess_(); |
void guess_(); |
17 |
|
void gufld_(float*, float*); |
18 |
|
float risxeta2_(float *); |
19 |
|
float risxeta3_(float *); |
20 |
|
float risxeta4_(float *); |
21 |
|
float risyeta2_(float *); |
22 |
} |
} |
23 |
|
|
24 |
//-------------------------------------- |
//-------------------------------------- |
25 |
// |
// |
26 |
// |
// |
27 |
//-------------------------------------- |
//-------------------------------------- |
28 |
TrkTrack::TrkTrack(){ |
TrkTrack::TrkTrack(){ |
29 |
|
// cout << "TrkTrack::TrkTrack()" << endl; |
30 |
seqno = -1; |
seqno = -1; |
31 |
image = -1; |
image = -1; |
32 |
chi2 = 0; |
chi2 = 0; |
36 |
for(int it2=0;it2<5;it2++)coval[it1][it2] = 0; |
for(int it2=0;it2<5;it2++)coval[it1][it2] = 0; |
37 |
}; |
}; |
38 |
for(int ip=0;ip<6;ip++){ |
for(int ip=0;ip<6;ip++){ |
39 |
xgood[ip] = 0; |
xgood[ip] = 0; |
40 |
ygood[ip] = 0; |
ygood[ip] = 0; |
41 |
xm[ip] = 0; |
xm[ip] = 0; |
42 |
ym[ip] = 0; |
ym[ip] = 0; |
43 |
zm[ip] = 0; |
zm[ip] = 0; |
44 |
resx[ip] = 0; |
resx[ip] = 0; |
45 |
resy[ip] = 0; |
resy[ip] = 0; |
46 |
xv[ip] = 0; |
tailx[ip] = 0; |
47 |
yv[ip] = 0; |
taily[ip] = 0; |
48 |
zv[ip] = 0; |
xv[ip] = 0; |
49 |
axv[ip] = 0; |
yv[ip] = 0; |
50 |
ayv[ip] = 0; |
zv[ip] = 0; |
51 |
dedx_x[ip] = 0; |
axv[ip] = 0; |
52 |
dedx_y[ip] = 0; |
ayv[ip] = 0; |
53 |
// clx[ip] = 0; |
dedx_x[ip] = 0; |
54 |
// cly[ip] = 0; |
dedx_y[ip] = 0; |
55 |
}; |
}; |
56 |
clx = new TRefArray(6,0); |
// clx = 0; |
57 |
cly = new TRefArray(6,0); |
// cly = 0; |
58 |
|
// clx = new TRefArray(6,0); //forse causa memory leak??? |
59 |
|
// cly = new TRefArray(6,0); //forse causa memory leak??? |
60 |
|
// clx = TRefArray(6,0); |
61 |
|
// cly = TRefArray(6,0); |
62 |
|
|
63 |
|
TrkParams::SetTrackingMode(); |
64 |
|
TrkParams::SetPrecisionFactor(); |
65 |
|
TrkParams::SetStepMin(); |
66 |
|
TrkParams::SetPFA(); |
67 |
|
|
68 |
}; |
}; |
69 |
//-------------------------------------- |
//-------------------------------------- |
70 |
// |
// |
80 |
for(int it2=0;it2<5;it2++)coval[it1][it2] = t.coval[it1][it2]; |
for(int it2=0;it2<5;it2++)coval[it1][it2] = t.coval[it1][it2]; |
81 |
}; |
}; |
82 |
for(int ip=0;ip<6;ip++){ |
for(int ip=0;ip<6;ip++){ |
83 |
xgood[ip] = t.xgood[ip]; |
xgood[ip] = t.xgood[ip]; |
84 |
ygood[ip] = t.ygood[ip]; |
ygood[ip] = t.ygood[ip]; |
85 |
xm[ip] = t.xm[ip]; |
xm[ip] = t.xm[ip]; |
86 |
ym[ip] = t.ym[ip]; |
ym[ip] = t.ym[ip]; |
87 |
zm[ip] = t.zm[ip]; |
zm[ip] = t.zm[ip]; |
88 |
resx[ip] = t.resx[ip]; |
resx[ip] = t.resx[ip]; |
89 |
resy[ip] = t.resy[ip]; |
resy[ip] = t.resy[ip]; |
90 |
xv[ip] = t.xv[ip]; |
tailx[ip] = t.tailx[ip]; |
91 |
yv[ip] = t.yv[ip]; |
taily[ip] = t.taily[ip]; |
92 |
zv[ip] = t.zv[ip]; |
xv[ip] = t.xv[ip]; |
93 |
axv[ip] = t.axv[ip]; |
yv[ip] = t.yv[ip]; |
94 |
ayv[ip] = t.ayv[ip]; |
zv[ip] = t.zv[ip]; |
95 |
dedx_x[ip] = t.dedx_x[ip]; |
axv[ip] = t.axv[ip]; |
96 |
dedx_y[ip] = t.dedx_y[ip]; |
ayv[ip] = t.ayv[ip]; |
97 |
// clx[ip] = 0;//<<<<pointer |
dedx_x[ip] = t.dedx_x[ip]; |
98 |
// cly[ip] = 0;//<<<<pointer |
dedx_y[ip] = t.dedx_y[ip]; |
99 |
}; |
}; |
100 |
clx = new TRefArray(*(t.clx)); |
// clx = 0; |
101 |
cly = new TRefArray(*(t.cly)); |
// cly = 0; |
102 |
|
// if(t.clx)clx = new TRefArray(*(t.clx)); |
103 |
|
// if(t.cly)cly = new TRefArray(*(t.cly)); |
104 |
|
// clx = TRefArray(t.clx); |
105 |
|
// cly = TRefArray(t.cly); |
106 |
|
|
107 |
|
TrkParams::SetTrackingMode(); |
108 |
|
TrkParams::SetPrecisionFactor(); |
109 |
|
TrkParams::SetStepMin(); |
110 |
|
TrkParams::SetPFA(); |
111 |
|
|
112 |
|
}; |
113 |
|
//-------------------------------------- |
114 |
|
// |
115 |
|
// |
116 |
|
//-------------------------------------- |
117 |
|
void TrkTrack::Copy(TrkTrack& t){ |
118 |
|
|
119 |
|
t.seqno = seqno; |
120 |
|
t.image = image; |
121 |
|
t.chi2 = chi2; |
122 |
|
t.nstep = nstep; |
123 |
|
for(int it1=0;it1<5;it1++){ |
124 |
|
t.al[it1] = al[it1]; |
125 |
|
for(int it2=0;it2<5;it2++)t.coval[it1][it2] = coval[it1][it2]; |
126 |
|
}; |
127 |
|
for(int ip=0;ip<6;ip++){ |
128 |
|
t.xgood[ip] = xgood[ip]; |
129 |
|
t.ygood[ip] = ygood[ip]; |
130 |
|
t.xm[ip] = xm[ip]; |
131 |
|
t.ym[ip] = ym[ip]; |
132 |
|
t.zm[ip] = zm[ip]; |
133 |
|
t.resx[ip] = resx[ip]; |
134 |
|
t.resy[ip] = resy[ip]; |
135 |
|
t.tailx[ip] = tailx[ip]; |
136 |
|
t.taily[ip] = taily[ip]; |
137 |
|
t.xv[ip] = xv[ip]; |
138 |
|
t.yv[ip] = yv[ip]; |
139 |
|
t.zv[ip] = zv[ip]; |
140 |
|
t.axv[ip] = axv[ip]; |
141 |
|
t.ayv[ip] = ayv[ip]; |
142 |
|
t.dedx_x[ip] = dedx_x[ip]; |
143 |
|
t.dedx_y[ip] = dedx_y[ip]; |
144 |
|
|
145 |
|
}; |
146 |
|
|
147 |
|
// t.clx = TRefArray(clx); |
148 |
|
// t.cly = TRefArray(cly); |
149 |
|
|
150 |
}; |
}; |
151 |
//-------------------------------------- |
//-------------------------------------- |
152 |
// |
// |
171 |
for (int i=0; i<5; i++) dal[i] = (double)al[i]; |
for (int i=0; i<5; i++) dal[i] = (double)al[i]; |
172 |
for (int i=0; i<t->npoint; i++) dzin[i] = (double)t->z[i]; |
for (int i=0; i<t->npoint; i++) dzin[i] = (double)t->z[i]; |
173 |
|
|
174 |
|
TrkParams::Load(1); |
175 |
|
if( !TrkParams::IsLoaded(1) ){ |
176 |
|
cout << "int TrkTrack::DoTrack(Trajectory* t) --- ERROR --- m.field not loaded"<<endl; |
177 |
|
return 0; |
178 |
|
} |
179 |
dotrack_(&(t->npoint),dzin,dxout,dyout,dal,&ifail); |
dotrack_(&(t->npoint),dzin,dxout,dyout,dal,&ifail); |
180 |
|
|
181 |
for (int i=0; i<t->npoint; i++){ |
for (int i=0; i<t->npoint; i++){ |
215 |
for (int i=0; i<5; i++) dal[i] = (double)al[i]; |
for (int i=0; i<5; i++) dal[i] = (double)al[i]; |
216 |
for (int i=0; i<t->npoint; i++) dzin[i] = (double)t->z[i]; |
for (int i=0; i<t->npoint; i++) dzin[i] = (double)t->z[i]; |
217 |
|
|
218 |
|
TrkParams::Load(1); |
219 |
|
if( !TrkParams::IsLoaded(1) ){ |
220 |
|
cout << "int TrkTrack::DoTrack2(Trajectory* t) --- ERROR --- m.field not loaded"<<endl; |
221 |
|
return 0; |
222 |
|
} |
223 |
dotrack2_(&(t->npoint),dzin,dxout,dyout,dthxout,dthyout,dtlout,dal,&ifail); |
dotrack2_(&(t->npoint),dzin,dxout,dyout,dthxout,dthyout,dtlout,dal,&ifail); |
224 |
|
|
225 |
for (int i=0; i<t->npoint; i++){ |
for (int i=0; i<t->npoint; i++){ |
259 |
return def; |
return def; |
260 |
}; |
}; |
261 |
// |
// |
262 |
Float_t TrkTrack::GetDEDX(){ |
/** |
263 |
Float_t dedx=0; |
* Method to retrieve the dE/dx measured on a tracker view. |
264 |
for(Int_t i=0; i<6; i++)dedx+=dedx_x[i]*xgood[i]+dedx_y[i]*ygood[i]; |
* @param ip plane (0-5) |
265 |
dedx = dedx/(this->GetNX()+this->GetNY()); |
* @param iv view (0=x 1=y) |
266 |
return dedx; |
*/ |
267 |
|
Float_t TrkTrack::GetDEDX(int ip, int iv){ |
268 |
|
if(iv==0 && ip>=0 && ip<6)return fabs(dedx_x[ip]); |
269 |
|
else if(iv==1 && ip>=0 && ip<6)return fabs(dedx_y[ip]); |
270 |
|
else { |
271 |
|
cout << "TrkTrack::GetDEDX(int ip, int iv) -- wrong input parameters "<<ip<<iv<<endl; |
272 |
|
return 0.; |
273 |
|
} |
274 |
|
} |
275 |
|
/** |
276 |
|
* Method to evaluate the dE/dx measured on a tracker plane. |
277 |
|
* The two measurements on x- and y-view are averaged. |
278 |
|
* @param ip plane (0-5) |
279 |
|
*/ |
280 |
|
Float_t TrkTrack::GetDEDX(int ip){ |
281 |
|
if( (Int_t)XGood(ip)+(Int_t)YGood(ip) == 0 ) return 0; |
282 |
|
return (GetDEDX(ip,0)+GetDEDX(ip,1))/((Int_t)XGood(ip)+(Int_t)YGood(ip)); |
283 |
}; |
}; |
284 |
|
|
285 |
|
/** |
286 |
|
* Method to evaluate the dE/dx averaged over all planes. |
287 |
|
*/ |
288 |
|
Float_t TrkTrack::GetDEDX(){ |
289 |
|
Float_t dedx=0; |
290 |
|
for(Int_t ip=0; ip<6; ip++)dedx+=GetDEDX(ip,0)*XGood(ip)+GetDEDX(ip,1)*YGood(ip); |
291 |
|
dedx = dedx/(GetNX()+GetNY()); |
292 |
|
return dedx; |
293 |
|
}; |
294 |
|
/** |
295 |
|
* Returns 1 if the cluster on a tracker view includes bad strips. |
296 |
|
* @param ip plane (0-5) |
297 |
|
* @param iv view (0=x 1=y) |
298 |
|
*/ |
299 |
|
Bool_t TrkTrack::IsBad(int ip,int iv){ |
300 |
|
if(iv==0 && ip>=0 && ip<6)return (xgood[ip]<0) ; |
301 |
|
else if(iv==1 && ip>=0 && ip<6)return (ygood[ip]<0) ; |
302 |
|
else { |
303 |
|
cout << "TrkTrack::IsBad(int ip, int iv) -- wrong input parameters "<<ip<<iv<<endl; |
304 |
|
return 0.; |
305 |
|
} |
306 |
|
}; |
307 |
|
/** |
308 |
|
* Returns 1 if the signal on a tracker view is saturated. |
309 |
|
* @param ip plane (0-5) |
310 |
|
* @param iv view (0=x 1=y) |
311 |
|
*/ |
312 |
|
Bool_t TrkTrack::IsSaturated(int ip,int iv){ |
313 |
|
if(iv==0 && ip>=0 && ip<6)return (dedx_x[ip]<0) ; |
314 |
|
else if(iv==1 && ip>=0 && ip<6)return (dedx_y[ip]<0) ; |
315 |
|
else { |
316 |
|
cout << "TrkTrack::IsSaturated(int ip, int iv) -- wrong input parameters "<<ip<<iv<<endl; |
317 |
|
return 0.; |
318 |
|
} |
319 |
|
}; |
320 |
|
/** |
321 |
|
* Returns 1 if either the x or the y signal on a tracker plane is saturated. |
322 |
|
* @param ip plane (0-5) |
323 |
|
*/ |
324 |
|
Bool_t TrkTrack::IsSaturated(int ip){ |
325 |
|
return (IsSaturated(ip,0)||IsSaturated(ip,1)); |
326 |
|
}; |
327 |
|
/** |
328 |
|
* Returns 1 if there is at least a saturated signal along the track. |
329 |
|
*/ |
330 |
|
Bool_t TrkTrack::IsSaturated(){ |
331 |
|
for(int ip=0; ip<6; ip++)for(int iv=0; iv<2; iv++)if(IsSaturated(ip,iv))return true; |
332 |
|
return false; |
333 |
|
} |
334 |
|
/** |
335 |
|
* Returns the track "lever-arm" on the x view, defined as the distance (in planes) between |
336 |
|
* the upper and lower x measurements (the maximum value of lever-arm is 6). |
337 |
|
*/ |
338 |
|
Int_t TrkTrack::GetLeverArmX(){ |
339 |
|
int first_plane = -1; |
340 |
|
int last_plane = -1; |
341 |
|
for(Int_t ip=0; ip<6; ip++){ |
342 |
|
if( XGood(ip) && first_plane == -1 )first_plane = ip; |
343 |
|
if( XGood(ip) && first_plane != -1 )last_plane = ip; |
344 |
|
} |
345 |
|
if( first_plane == -1 || last_plane == -1){ |
346 |
|
cout<< "Int_t TrkTrack::GetLeverArmX() -- XGood(ip) always false ??? "<<endl; |
347 |
|
return 0; |
348 |
|
} |
349 |
|
return (last_plane-first_plane+1); |
350 |
|
} |
351 |
|
/** |
352 |
|
* Returns the track "lever-arm" on the y view, defined as the distance (in planes) between |
353 |
|
* the upper and lower y measurements (the maximum value of lever-arm is 6). |
354 |
|
*/ |
355 |
|
Int_t TrkTrack::GetLeverArmY(){ |
356 |
|
int first_plane = -1; |
357 |
|
int last_plane = -1; |
358 |
|
for(Int_t ip=0; ip<6; ip++){ |
359 |
|
if( YGood(ip) && first_plane == -1 )first_plane = ip; |
360 |
|
if( YGood(ip) && first_plane != -1 )last_plane = ip; |
361 |
|
} |
362 |
|
if( first_plane == -1 || last_plane == -1){ |
363 |
|
cout<< "Int_t TrkTrack::GetLeverArmY() -- YGood(ip) always false ??? "<<endl; |
364 |
|
return 0; |
365 |
|
} |
366 |
|
return (last_plane-first_plane+1); |
367 |
|
} |
368 |
//-------------------------------------- |
//-------------------------------------- |
369 |
// |
// |
370 |
// |
// |
376 |
cout << endl << "al : "; for(int i=0; i<5; i++)cout << al[i] << " "; |
cout << endl << "al : "; for(int i=0; i<5; i++)cout << al[i] << " "; |
377 |
cout << endl << "chi^2 : "<< chi2; |
cout << endl << "chi^2 : "<< chi2; |
378 |
cout << endl << "n.step : "<< nstep; |
cout << endl << "n.step : "<< nstep; |
379 |
cout << endl << "xgood : "; for(int i=0; i<6; i++)cout << xgood[i] ; |
cout << endl << "xgood : "; for(int i=0; i<6; i++)cout << XGood(i) ; |
380 |
cout << endl << "ygood : "; for(int i=0; i<6; i++)cout << ygood[i] ; |
cout << endl << "ygood : "; for(int i=0; i<6; i++)cout << YGood(i) ; |
381 |
cout << endl << "xm : "; for(int i=0; i<6; i++)cout << xm[i] << " "; |
cout << endl << "xm : "; for(int i=0; i<6; i++)cout << xm[i] << " "; |
382 |
cout << endl << "ym : "; for(int i=0; i<6; i++)cout << ym[i] << " "; |
cout << endl << "ym : "; for(int i=0; i<6; i++)cout << ym[i] << " "; |
383 |
cout << endl << "zm : "; for(int i=0; i<6; i++)cout << zm[i] << " "; |
cout << endl << "zm : "; for(int i=0; i<6; i++)cout << zm[i] << " "; |
386 |
cout << endl << "zv : "; for(int i=0; i<6; i++)cout << zv[i] << " "; |
cout << endl << "zv : "; for(int i=0; i<6; i++)cout << zv[i] << " "; |
387 |
cout << endl << "resx : "; for(int i=0; i<6; i++)cout << resx[i] << " "; |
cout << endl << "resx : "; for(int i=0; i<6; i++)cout << resx[i] << " "; |
388 |
cout << endl << "resy : "; for(int i=0; i<6; i++)cout << resy[i] << " "; |
cout << endl << "resy : "; for(int i=0; i<6; i++)cout << resy[i] << " "; |
389 |
|
cout << endl << "tailx : "; for(int i=0; i<6; i++)cout << tailx[i] << " "; |
390 |
|
cout << endl << "taily : "; for(int i=0; i<6; i++)cout << taily[i] << " "; |
391 |
cout << endl << "coval : "; for(int i=0; i<5; i++)cout << coval[0][i]<<" "; |
cout << endl << "coval : "; for(int i=0; i<5; i++)cout << coval[0][i]<<" "; |
392 |
cout << endl << " "; for(int i=0; i<5; i++)cout << coval[1][i]<<" "; |
cout << endl << " "; for(int i=0; i<5; i++)cout << coval[1][i]<<" "; |
393 |
cout << endl << " "; for(int i=0; i<5; i++)cout << coval[2][i]<<" "; |
cout << endl << " "; for(int i=0; i<5; i++)cout << coval[2][i]<<" "; |
413 |
for(int i=0; i<6; i++) resy[i]=*ry++; |
for(int i=0; i<6; i++) resy[i]=*ry++; |
414 |
} |
} |
415 |
/** |
/** |
416 |
|
* Set the TrkTrack tails position resolution |
417 |
|
*/ |
418 |
|
void TrkTrack::SetTail(double *tx, double *ty, double factor){ |
419 |
|
for(int i=0; i<6; i++) tailx[i]=factor*(*tx++); |
420 |
|
for(int i=0; i<6; i++) taily[i]=factor*(*ty++); |
421 |
|
} |
422 |
|
/** |
423 |
|
* Set the TrkTrack Student parameter (resx,resy,tailx,taily) |
424 |
|
* from previous gausian fit |
425 |
|
*@param flag =0 standard, =1 with noise correction |
426 |
|
*/ |
427 |
|
void TrkTrack::SetStudentParam(int flag){ |
428 |
|
float sx[11]={0.000128242, |
429 |
|
0.000136942, |
430 |
|
0.000162718, |
431 |
|
0.000202644, |
432 |
|
0.00025597, |
433 |
|
0.000317456, |
434 |
|
0.000349048, |
435 |
|
0.000384638, |
436 |
|
0.000457295, |
437 |
|
0.000512319, |
438 |
|
0.000538573}; |
439 |
|
float tx[11]={1.79402, |
440 |
|
2.04876, |
441 |
|
2.88376, |
442 |
|
3.3, |
443 |
|
3.14084, |
444 |
|
4.07686, |
445 |
|
4.44736, |
446 |
|
3.5179, |
447 |
|
3.38697, |
448 |
|
3.45739, |
449 |
|
3.18627}; |
450 |
|
float sy[11]={0.000483075, |
451 |
|
0.000466925, |
452 |
|
0.000431658, |
453 |
|
0.000428317, |
454 |
|
0.000433854, |
455 |
|
0.000444044, |
456 |
|
0.000482098, |
457 |
|
0.000537579, |
458 |
|
0.000636279, |
459 |
|
0.000741998, |
460 |
|
0.000864261}; |
461 |
|
float ty[11]={0.997032, |
462 |
|
1.11147, |
463 |
|
1.18526, |
464 |
|
1.61404, |
465 |
|
2.21908, |
466 |
|
3.08959, |
467 |
|
4.48833, |
468 |
|
4.42687, |
469 |
|
4.65253, |
470 |
|
4.52043, |
471 |
|
4.29926}; |
472 |
|
int index; |
473 |
|
float fact; |
474 |
|
for(int i=0; i<6; i++) { |
475 |
|
index = int((fabs(axv[i])+1.)/2.); |
476 |
|
if(index>10) index=10; |
477 |
|
tailx[i]=tx[index]; |
478 |
|
if(flag==1) { |
479 |
|
if(fabs(axv[i])<=10.) fact = resx[i]/risxeta2_(&(axv[i])); |
480 |
|
if(fabs(axv[i])>10.&&fabs(axv[i])<=15.) fact = resx[i]/risxeta3_(&(axv[i])); |
481 |
|
if(fabs(axv[i])>15.) fact = resx[i]/risxeta4_(&(axv[i])); |
482 |
|
} else fact = 1.; |
483 |
|
resx[i] = sx[index]*fact; |
484 |
|
} |
485 |
|
for(int i=0; i<6; i++) { |
486 |
|
index = int((fabs(ayv[i])+1.)/2.); |
487 |
|
if(index>10) index=10; |
488 |
|
taily[i]=ty[index]; |
489 |
|
if(flag==1) fact = resy[i]/risyeta2_(&(ayv[i])); |
490 |
|
else fact = 1.; |
491 |
|
resy[i] = sy[index]*fact; |
492 |
|
} |
493 |
|
} |
494 |
|
/** |
495 |
* Set the TrkTrack good measurement |
* Set the TrkTrack good measurement |
496 |
*/ |
*/ |
497 |
void TrkTrack::SetGood(int *xg, int *yg){ |
void TrkTrack::SetGood(int *xg, int *yg){ |
498 |
|
|
499 |
for(int i=0; i<6; i++) xgood[i]=*xg++; |
for(int i=0; i<6; i++) xgood[i]=*xg++; |
500 |
for(int i=0; i<6; i++) ygood[i]=*yg++; |
for(int i=0; i<6; i++) ygood[i]=*yg++; |
501 |
} |
} |
505 |
*/ |
*/ |
506 |
void TrkTrack::LoadField(TString path){ |
void TrkTrack::LoadField(TString path){ |
507 |
|
|
508 |
strcpy(path_.path,path.Data()); |
// strcpy(path_.path,path.Data()); |
509 |
path_.pathlen = path.Length(); |
// path_.pathlen = path.Length(); |
510 |
path_.error = 0; |
// path_.error = 0; |
511 |
readb_(); |
// readb_(); |
512 |
|
|
513 |
|
TrkParams::SetTrackingMode(); |
514 |
|
TrkParams::SetPrecisionFactor(); |
515 |
|
TrkParams::SetStepMin(); |
516 |
|
|
517 |
|
TrkParams::Set(path,1); |
518 |
|
TrkParams::Load(1); |
519 |
|
|
520 |
}; |
}; |
521 |
|
|
522 |
|
|
523 |
/** |
/** |
524 |
* Method to fill minimization-routine common |
* Method to fill minimization-routine common |
525 |
*/ |
*/ |
527 |
|
|
528 |
for(int i=0; i<6; i++){ |
for(int i=0; i<6; i++){ |
529 |
|
|
530 |
track.xgood[i]=xgood[i]; |
// cout << i<<" - "<<xgood[i]<<" "<<XGood(i)<<endl; |
531 |
track.ygood[i]=ygood[i]; |
// cout << i<<" - "<<ygood[i]<<" "<<YGood(i)<<endl; |
532 |
|
track.xgood[i]=XGood(i); |
533 |
|
track.ygood[i]=YGood(i); |
534 |
|
|
535 |
track.xm[i]=xm[i]; |
track.xm[i]=xm[i]; |
536 |
track.ym[i]=ym[i]; |
track.ym[i]=ym[i]; |
543 |
track.xm_b[i]=xm[i]; |
track.xm_b[i]=xm[i]; |
544 |
track.ym_a[i]=ym[i]; |
track.ym_a[i]=ym[i]; |
545 |
track.ym_b[i]=ym[i]; |
track.ym_b[i]=ym[i]; |
546 |
if( xgood[i] && !ygood[i] ){ |
if( XGood(i) && !YGood(i) ){ |
547 |
track.ym_a[i] = track.ym_a[i]+segment; |
track.ym_a[i] = track.ym_a[i]+segment; |
548 |
track.ym_b[i] = track.ym_b[i]-segment; |
track.ym_b[i] = track.ym_b[i]-segment; |
549 |
}else if( !xgood[i] && ygood[i]){ |
}else if( !XGood(i) && YGood(i)){ |
550 |
track.xm_a[i] = track.xm_a[i]+segment; |
track.xm_a[i] = track.xm_a[i]+segment; |
551 |
track.xm_b[i] = track.xm_b[i]-segment; |
track.xm_b[i] = track.xm_b[i]-segment; |
552 |
} |
} |
554 |
|
|
555 |
track.resx[i]=resx[i]; |
track.resx[i]=resx[i]; |
556 |
track.resy[i]=resy[i]; |
track.resy[i]=resy[i]; |
557 |
|
track.tailx[i]=tailx[i]; |
558 |
|
track.taily[i]=taily[i]; |
559 |
} |
} |
560 |
|
|
561 |
for(int i=0; i<5; i++) track.al[i]=al[i]; |
for(int i=0; i<5; i++) track.al[i]=al[i]; |
584 |
axv[i] = track->axv[i]; |
axv[i] = track->axv[i]; |
585 |
ayv[i] = track->ayv[i]; |
ayv[i] = track->ayv[i]; |
586 |
} |
} |
587 |
|
|
588 |
|
} |
589 |
|
/** |
590 |
|
* \brief Method to re-evaluate coordinates of clusters associated with a track. |
591 |
|
* |
592 |
|
* The method can be applied only after recovering level1 information |
593 |
|
* (either by reprocessing single events from level0 or from |
594 |
|
* the TrkLevel1 branch, if present); it calls F77 subroutines that |
595 |
|
* read the level1 common and fill the minimization-routine common. |
596 |
|
* Some clusters can be excluded or added by means of the methods: |
597 |
|
* |
598 |
|
* TrkTrack::ResetXGood(int ip) |
599 |
|
* TrkTrack::ResetYGood(int ip) |
600 |
|
* TrkTrack::SetXGood(int ip, int cid, int is) |
601 |
|
* TrkTrack::SetYGood(int ip, int cid, int is) |
602 |
|
* |
603 |
|
* NB! The method TrkTrack::SetGood(int *xg, int *yg) set the plane-mask (0-1) |
604 |
|
* for the minimization-routine common. It deletes the cluster information |
605 |
|
* (at least for the moment...) thus cannot be applied before |
606 |
|
* TrkTrack::EvaluateClusterPositions(). |
607 |
|
* |
608 |
|
* Different p.f.a. can be applied by calling (once) the method: |
609 |
|
* |
610 |
|
* TrkParams::SetPFA(0); //Set ETA p.f.a. |
611 |
|
* |
612 |
|
* @see TrkParams::SetPFA(int) |
613 |
|
*/ |
614 |
|
void TrkTrack::EvaluateClusterPositions(){ |
615 |
|
|
616 |
|
// cout << "void TrkTrack::GetClusterPositions() "<<endl; |
617 |
|
|
618 |
|
TrkParams::Load( ); |
619 |
|
if( !TrkParams::IsLoaded() )return; |
620 |
|
|
621 |
|
for(int ip=0; ip<6; ip++){ |
622 |
|
// cout << ip<<" ** "<<xm[ip]<<" / "<<ym[ip]<<endl;; |
623 |
|
int icx = GetClusterX_ID(ip)+1; |
624 |
|
int icy = GetClusterY_ID(ip)+1; |
625 |
|
int sensor = GetSensor(ip)+1;//<< convenzione "Paolo" |
626 |
|
if(ip==5 && sensor!=0)sensor=3-sensor;//<< convenzione "Elena" |
627 |
|
int ladder = GetLadder(ip)+1; |
628 |
|
float ax = axv[ip]; |
629 |
|
float ay = ayv[ip]; |
630 |
|
float v[3]; |
631 |
|
v[0]=xv[ip]; |
632 |
|
v[1]=yv[ip]; |
633 |
|
v[2]=zv[ip]; |
634 |
|
float bfx = 10*TrkParams::GetBX(v);//Tesla |
635 |
|
float bfy = 10*TrkParams::GetBY(v);//Tesla |
636 |
|
int ipp=ip+1; |
637 |
|
xyzpam_(&ipp,&icx,&icy,&ladder,&sensor,&ax,&ay,&bfx,&bfy); |
638 |
|
} |
639 |
} |
} |
640 |
/** |
/** |
641 |
* Tracking method. It calls F77 mini routine. |
* \brief Tracking method. It calls F77 mini routine. |
642 |
|
* |
643 |
|
* @param pfixed Particle momentum. If pfixed=0 the momentum |
644 |
|
* is left as a free parameter, otherwise it is fixed to the input value. |
645 |
|
* @param fail Output flag (!=0 if the fit failed). |
646 |
|
* @param iprint Flag to set debug mode ( 0 = no output; 1 = verbose; 2 = debug). |
647 |
|
* @param froml1 Flag to re-evaluate positions (see TrkTrack::GetClusterPositions()). |
648 |
|
* |
649 |
|
* The option to re-evaluate positions can be used only after recovering |
650 |
|
* level1 information, eg. by reprocessing the single event. |
651 |
|
* |
652 |
|
* Example: |
653 |
|
* |
654 |
|
* if( !event->GetTrkLevel0() )return false; |
655 |
|
* event->GetTrkLevel0()->ProcessEvent(); // re-processing level0->level1 |
656 |
|
* int fail=0; |
657 |
|
* event->GetTrkLevel2()->GetTrack(0)->Fit(0.,fail,0,1); |
658 |
|
* |
659 |
|
* @see EvaluateClusterPositions() |
660 |
|
* |
661 |
|
* The fitting procedure can be varied by changing the tracking mode, |
662 |
|
* the fit-precision factor and the minimum number of step. |
663 |
|
* @see SetTrackingMode(int) |
664 |
|
* @see SetPrecisionFactor(double) |
665 |
|
* @see SetStepMin(int) |
666 |
*/ |
*/ |
667 |
void TrkTrack::Fit(double pfixed, int& fail, int iprint){ |
void TrkTrack::Fit(double pfixed, int& fail, int iprint, int froml1){ |
668 |
|
|
669 |
float al_ini[] = {0.,0.,0.,0.,0.}; |
float al_ini[] = {0.,0.,0.,0.,0.}; |
670 |
|
|
671 |
|
TrkParams::Load( ); |
672 |
|
if( !TrkParams::IsLoaded() )return; |
673 |
|
|
674 |
extern cMini2track track_; |
extern cMini2track track_; |
675 |
fail = 0; |
fail = 0; |
676 |
FillMiniStruct(track_); |
FillMiniStruct(track_); |
677 |
|
|
678 |
|
if(froml1!=0)EvaluateClusterPositions(); |
679 |
|
|
680 |
// if fit variables have been reset, evaluate the initial guess |
// if fit variables have been reset, evaluate the initial guess |
681 |
if(al[0]==-9999.&&al[1]==-9999.&&al[2]==-9999.&&al[3]==-9999.&&al[4]==-9999.)guess_(); |
if(al[0]==-9999.&&al[1]==-9999.&&al[2]==-9999.&&al[3]==-9999.&&al[4]==-9999.)guess_(); |
682 |
|
|
695 |
|
|
696 |
// ------------------------------------------ |
// ------------------------------------------ |
697 |
// call mini routine |
// call mini routine |
698 |
|
// TrkParams::Load(1); |
699 |
|
// if( !TrkParams::IsLoaded(1) ){ |
700 |
|
// cout << "void TrkTrack::Fit(double pfixed, int& fail, int iprint) --- ERROR --- m.field not loaded"<<endl; |
701 |
|
// return; |
702 |
|
// } |
703 |
int istep=0; |
int istep=0; |
704 |
int ifail=0; |
int ifail=0; |
705 |
mini2_(&istep,&ifail, &iprint); |
mini2_(&istep,&ifail, &iprint); |
710 |
// ------------------------------------------ |
// ------------------------------------------ |
711 |
|
|
712 |
SetFromMiniStruct(&track_); |
SetFromMiniStruct(&track_); |
|
// cout << endl << "eta ===> " << track_.al[4] << endl; |
|
713 |
|
|
714 |
// for(int i=0; i<5; i++) al[i]=track_.al[i]; |
if(fail){ |
715 |
// chi2=track_.chi2; |
if(iprint)cout << " >>>> fit failed "<<endl; |
716 |
// nstep=track_.nstep; |
for(int i=0; i<5; i++) al[i]=al_ini[i]; |
717 |
// for(int i=0; i<6; i++) xv[i]=track_.xv[i]; |
} |
|
// for(int i=0; i<6; i++) yv[i]=track_.yv[i]; |
|
|
// for(int i=0; i<6; i++) zv[i]=track_.zv[i]; |
|
|
// for(int i=0; i<6; i++) axv[i]=track_.axv[i]; |
|
|
// for(int i=0; i<6; i++) ayv[i]=track_.ayv[i]; |
|
|
// for(int i=0; i<5; i++) { |
|
|
// for(int j=0; j<5; j++) coval[i][j]=track_.cov[i][j]; |
|
|
// } |
|
|
|
|
|
// if(fail && iprint){ |
|
|
// cout << " >>>> fit failed >>>> drawing initial par"<<endl; |
|
|
// for(int i=0; i<5; i++) al[i]=al_ini[i]; |
|
|
// } |
|
718 |
|
|
719 |
}; |
}; |
720 |
/* |
/** |
721 |
* Reset the fit parameters |
* Reset the fit parameters |
722 |
*/ |
*/ |
723 |
void TrkTrack::FitReset(){ |
void TrkTrack::FitReset(){ |
724 |
for(int i=0; i<5; i++) al[i]=-9999.; |
for(int i=0; i<5; i++) al[i]=-9999.; |
725 |
chi2=0.; |
chi2=0.; |
726 |
nstep=0; |
nstep=0; |
727 |
for(int i=0; i<6; i++) xv[i]=0.; |
// for(int i=0; i<6; i++) xv[i]=0.; |
728 |
for(int i=0; i<6; i++) yv[i]=0.; |
// for(int i=0; i<6; i++) yv[i]=0.; |
729 |
for(int i=0; i<6; i++) zv[i]=0.; |
// for(int i=0; i<6; i++) zv[i]=0.; |
730 |
for(int i=0; i<6; i++) axv[i]=0.; |
// for(int i=0; i<6; i++) axv[i]=0.; |
731 |
for(int i=0; i<6; i++) ayv[i]=0.; |
// for(int i=0; i<6; i++) ayv[i]=0.; |
732 |
for(int i=0; i<5; i++) { |
for(int i=0; i<5; i++) { |
733 |
for(int j=0; j<5; j++) coval[i][j]=0.; |
for(int j=0; j<5; j++) coval[i][j]=0.; |
734 |
} |
} |
735 |
} |
} |
736 |
|
/** |
737 |
|
* Set the tracking mode |
738 |
|
*/ |
739 |
|
void TrkTrack::SetTrackingMode(int trackmode){ |
740 |
|
extern cMini2track track_; |
741 |
|
track_.trackmode = trackmode; |
742 |
|
} |
743 |
|
/** |
744 |
|
* Set the factor scale for tracking precision |
745 |
|
*/ |
746 |
|
void TrkTrack::SetPrecisionFactor(double fact){ |
747 |
|
extern cMini2track track_; |
748 |
|
track_.fact = fact; |
749 |
|
} |
750 |
|
/** |
751 |
|
* Set the minimum number of steps for tracking precision |
752 |
|
*/ |
753 |
|
void TrkTrack::SetStepMin(int istepmin){ |
754 |
|
extern cMini2track track_; |
755 |
|
track_.istepmin = istepmin; |
756 |
|
} |
757 |
|
|
758 |
|
|
759 |
|
/** |
760 |
|
* Method to retrieve ID (0,1,...) of x-cluster (if any) associated to this track. |
761 |
|
* If no cluster is associated, ID=-1. |
762 |
|
* @param ip Tracker plane (0-5) |
763 |
|
*/ |
764 |
|
Int_t TrkTrack::GetClusterX_ID(int ip){ |
765 |
|
return ((Int_t)fabs(xgood[ip]))%10000000-1; |
766 |
|
}; |
767 |
|
/** |
768 |
|
* Method to retrieve ID (0-xxx) of y-cluster (if any) associated to this track. |
769 |
|
* If no cluster is associated, ID=-1. |
770 |
|
* @param ip Tracker plane (0-5) |
771 |
|
*/ |
772 |
|
Int_t TrkTrack::GetClusterY_ID(int ip){ |
773 |
|
return ((Int_t)fabs(ygood[ip]))%10000000-1; |
774 |
|
}; |
775 |
|
/** |
776 |
|
* Method to retrieve the ladder (0-4, increasing x) traversed by the track on this plane. |
777 |
|
* If no ladder is traversed (dead area) the metod retuns -1. |
778 |
|
* @param ip Tracker plane (0-5) |
779 |
|
*/ |
780 |
|
Int_t TrkTrack::GetLadder(int ip){ |
781 |
|
if(XGood(ip))return (Int_t)fabs(xgood[ip]/100000000)-1; |
782 |
|
if(YGood(ip))return (Int_t)fabs(ygood[ip]/100000000)-1; |
783 |
|
return -1; |
784 |
|
}; |
785 |
|
/** |
786 |
|
* Method to retrieve the sensor (0-1, increasing y) traversed by the track on this plane. |
787 |
|
* If no sensor is traversed (dead area) the metod retuns -1. |
788 |
|
* @param ip Tracker plane (0-5) |
789 |
|
*/ |
790 |
|
Int_t TrkTrack::GetSensor(int ip){ |
791 |
|
if(XGood(ip))return (Int_t)((Int_t)fabs(xgood[ip]/10000000)%10)-1; |
792 |
|
if(YGood(ip))return (Int_t)((Int_t)fabs(ygood[ip]/10000000)%10)-1; |
793 |
|
return -1; |
794 |
|
}; |
795 |
|
|
796 |
|
/** |
797 |
|
* \brief Method to include a x-cluster to the track. |
798 |
|
* @param ip Tracker plane (0-5) |
799 |
|
* @param clid Cluster ID (0,1,...) |
800 |
|
* @param is Sensor (0-1, increasing y) |
801 |
|
* @see Fit(double pfixed, int& fail, int iprint, int froml1) |
802 |
|
*/ |
803 |
|
void TrkTrack::SetXGood(int ip, int clid, int is){ |
804 |
|
int il=0; //ladder (temporary) |
805 |
|
bool bad=false; //ladder (temporary) |
806 |
|
xgood[ip]=il*100000000+is*10000000+clid; |
807 |
|
if(bad)xgood[ip]=-xgood[ip]; |
808 |
|
}; |
809 |
|
/** |
810 |
|
* \brief Method to include a y-cluster to the track. |
811 |
|
* @param ip Tracker plane (0-5) |
812 |
|
* @param clid Cluster ID (0,1,...) |
813 |
|
* @param is Sensor (0-1) |
814 |
|
* @see Fit(double pfixed, int& fail, int iprint, int froml1) |
815 |
|
*/ |
816 |
|
void TrkTrack::SetYGood(int ip, int clid, int is){ |
817 |
|
int il=0; //ladder (temporary) |
818 |
|
bool bad=false; //ladder (temporary) |
819 |
|
ygood[ip]=il*100000000+is*10000000+clid; |
820 |
|
if(bad)ygood[ip]=-ygood[ip]; |
821 |
|
}; |
822 |
|
|
823 |
//-------------------------------------- |
//-------------------------------------- |
824 |
// |
// |
825 |
// |
// |
826 |
//-------------------------------------- |
//-------------------------------------- |
827 |
void TrkTrack::Clear(){ |
void TrkTrack::Clear(){ |
828 |
seqno = -1; |
// cout << "TrkTrack::Clear()"<<endl; |
829 |
image = -1; |
seqno = -1; |
830 |
chi2 = 0; |
image = -1; |
831 |
nstep = 0; |
chi2 = 0; |
832 |
for(int it1=0;it1<5;it1++){ |
nstep = 0; |
833 |
al[it1] = 0; |
for(int it1=0;it1<5;it1++){ |
834 |
for(int it2=0;it2<5;it2++)coval[it1][it2] = 0; |
al[it1] = 0; |
835 |
}; |
for(int it2=0;it2<5;it2++)coval[it1][it2] = 0; |
836 |
for(int ip=0;ip<6;ip++){ |
}; |
837 |
xgood[ip] = 0; |
for(int ip=0;ip<6;ip++){ |
838 |
ygood[ip] = 0; |
xgood[ip] = 0; |
839 |
xm[ip] = 0; |
ygood[ip] = 0; |
840 |
ym[ip] = 0; |
xm[ip] = 0; |
841 |
zm[ip] = 0; |
ym[ip] = 0; |
842 |
resx[ip] = 0; |
zm[ip] = 0; |
843 |
resy[ip] = 0; |
resx[ip] = 0; |
844 |
xv[ip] = 0; |
resy[ip] = 0; |
845 |
yv[ip] = 0; |
tailx[ip] = 0; |
846 |
zv[ip] = 0; |
taily[ip] = 0; |
847 |
axv[ip] = 0; |
xv[ip] = 0; |
848 |
ayv[ip] = 0; |
yv[ip] = 0; |
849 |
dedx_x[ip] = 0; |
zv[ip] = 0; |
850 |
dedx_y[ip] = 0; |
axv[ip] = 0; |
851 |
// clx[ip] = 0; |
ayv[ip] = 0; |
852 |
// cly[ip] = 0; |
dedx_x[ip] = 0; |
853 |
}; |
dedx_y[ip] = 0; |
854 |
clx->Clear(); |
|
855 |
cly->Clear(); |
}; |
856 |
|
// if(clx)clx->Clear(); |
857 |
|
// if(cly)cly->Clear(); |
858 |
|
// clx.Clear(); |
859 |
|
// cly.Clear(); |
860 |
}; |
}; |
861 |
//-------------------------------------- |
//-------------------------------------- |
862 |
// |
// |
863 |
// |
// |
864 |
//-------------------------------------- |
//-------------------------------------- |
865 |
void TrkTrack::Delete(){ |
void TrkTrack::Delete(){ |
866 |
Clear(); |
// cout << "TrkTrack::Delete()"<<endl; |
867 |
clx->Delete(); |
Clear(); |
868 |
cly->Delete(); |
// if(clx)delete clx; |
869 |
|
// if(cly)delete cly; |
870 |
}; |
}; |
871 |
//-------------------------------------- |
//-------------------------------------- |
872 |
// |
// |
873 |
// |
// |
874 |
//-------------------------------------- |
//-------------------------------------- |
878 |
// |
// |
879 |
//-------------------------------------- |
//-------------------------------------- |
880 |
TrkSinglet::TrkSinglet(){ |
TrkSinglet::TrkSinglet(){ |
881 |
|
// cout << "TrkSinglet::TrkSinglet() " << GetUniqueID()<<endl; |
882 |
plane = 0; |
plane = 0; |
883 |
coord[0] = 0; |
coord[0] = 0; |
884 |
coord[1] = 0; |
coord[1] = 0; |
885 |
sgnl = 0; |
sgnl = 0; |
886 |
cls = 0; |
// cls = 0; |
887 |
}; |
}; |
888 |
//-------------------------------------- |
//-------------------------------------- |
889 |
// |
// |
890 |
// |
// |
891 |
//-------------------------------------- |
//-------------------------------------- |
892 |
TrkSinglet::TrkSinglet(const TrkSinglet& s){ |
TrkSinglet::TrkSinglet(const TrkSinglet& s){ |
893 |
|
// cout << "TrkSinglet::TrkSinglet(const TrkSinglet& s) " << GetUniqueID()<<endl; |
894 |
plane = s.plane; |
plane = s.plane; |
895 |
coord[0] = s.coord[0]; |
coord[0] = s.coord[0]; |
896 |
coord[1] = s.coord[1]; |
coord[1] = s.coord[1]; |
897 |
sgnl = s.sgnl; |
sgnl = s.sgnl; |
898 |
// cls = 0;//<<<<pointer |
// cls = 0;//<<<<pointer |
899 |
cls = TRef(s.cls); |
// cls = TRef(s.cls); |
900 |
}; |
}; |
901 |
//-------------------------------------- |
//-------------------------------------- |
902 |
// |
// |
913 |
// |
// |
914 |
// |
// |
915 |
//-------------------------------------- |
//-------------------------------------- |
916 |
|
void TrkSinglet::Clear(){ |
917 |
|
// cout << "TrkSinglet::Clear() " << GetUniqueID()<<endl; |
918 |
|
// cls=0; |
919 |
|
plane=-1; |
920 |
|
coord[0]=-999; |
921 |
|
coord[1]=-999; |
922 |
|
sgnl=0; |
923 |
|
|
924 |
|
} |
925 |
|
//-------------------------------------- |
926 |
|
// |
927 |
|
// |
928 |
|
//-------------------------------------- |
929 |
TrkLevel2::TrkLevel2(){ |
TrkLevel2::TrkLevel2(){ |
930 |
// good2 = -1; |
// cout <<"TrkLevel2::TrkLevel2()"<<endl; |
931 |
for(Int_t i=0; i<12 ; i++){ |
for(Int_t i=0; i<12 ; i++){ |
932 |
// crc[i] = -1; |
good[i] = -1; |
933 |
good[i] = -1; |
VKmask[i] = 0; |
934 |
}; |
VKflag[i] = 0; |
935 |
Track = new TClonesArray("TrkTrack"); |
}; |
936 |
SingletX = new TClonesArray("TrkSinglet"); |
Track = 0; |
937 |
SingletY = new TClonesArray("TrkSinglet"); |
SingletX = 0; |
938 |
|
SingletY = 0; |
939 |
|
|
940 |
// PhysicalTrack = new TClonesArray("TrkTrack"); |
} |
941 |
//sostituire con TRefArray... appena ho capito come si usa |
//-------------------------------------- |
942 |
|
// |
943 |
|
// |
944 |
|
//-------------------------------------- |
945 |
|
void TrkLevel2::Set(){ |
946 |
|
if(!Track)Track = new TClonesArray("TrkTrack"); |
947 |
|
if(!SingletX)SingletX = new TClonesArray("TrkSinglet"); |
948 |
|
if(!SingletY)SingletY = new TClonesArray("TrkSinglet"); |
949 |
} |
} |
950 |
//-------------------------------------- |
//-------------------------------------- |
951 |
// |
// |
953 |
//-------------------------------------- |
//-------------------------------------- |
954 |
void TrkLevel2::Dump(){ |
void TrkLevel2::Dump(){ |
955 |
|
|
|
TClonesArray &t = *Track; |
|
|
TClonesArray &sx = *SingletX; |
|
|
TClonesArray &sy = *SingletY; |
|
956 |
// |
// |
957 |
cout << endl << endl << "=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-"; |
cout << endl << endl << "=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-"; |
958 |
cout << endl << "good : "; for(int i=0; i<12; i++) cout << good[i]<<" "; |
cout << endl << "good : "; for(int i=0; i<12; i++) cout << good[i]<<" "; |
959 |
cout << endl << "ntrk() : " << this->ntrk() ; |
cout << endl << "ntrk() : " << this->ntrk() ; |
960 |
cout << endl << "nclsx() : " << this->nclsx(); |
cout << endl << "nclsx() : " << this->nclsx(); |
961 |
cout << endl << "nclsy() : " << this->nclsy(); |
cout << endl << "nclsy() : " << this->nclsy(); |
962 |
for(int i=0; i<this->ntrk(); i++) ((TrkTrack *)t[i])->Dump(); |
if(Track){ |
963 |
for(int i=0; i<this->nclsx(); i++) ((TrkSinglet *)sx[i])->Dump(); |
TClonesArray &t = *Track; |
964 |
for(int i=0; i<this->nclsy(); i++) ((TrkSinglet *)sy[i])->Dump(); |
for(int i=0; i<ntrk(); i++) ((TrkTrack *)t[i])->Dump(); |
965 |
|
} |
966 |
|
if(SingletX){ |
967 |
|
TClonesArray &sx = *SingletX; |
968 |
|
for(int i=0; i<nclsx(); i++) ((TrkSinglet *)sx[i])->Dump(); |
969 |
|
} |
970 |
|
if(SingletY){ |
971 |
|
TClonesArray &sy = *SingletY; |
972 |
|
for(int i=0; i<nclsy(); i++) ((TrkSinglet *)sy[i])->Dump(); |
973 |
|
} |
974 |
} |
} |
975 |
//-------------------------------------- |
//-------------------------------------- |
976 |
// |
// |
977 |
// |
// |
978 |
//-------------------------------------- |
//-------------------------------------- |
979 |
/** |
/** |
980 |
* Fills a TrkLevel2 object with values from a struct cTrkLevel2 (to get data from F77 common). |
* The method returns false if the viking-chip was masked |
981 |
|
* either apriori ,on the basis of the mask read from the DB, |
982 |
|
* or run-by-run, on the basis of the calibration parameters) |
983 |
|
* @param iv Tracker view (0-11) |
984 |
|
* @param ivk Viking-chip number (0-23) |
985 |
|
*/ |
986 |
|
Bool_t TrkLevel2::GetVKMask(int iv, int ivk){ |
987 |
|
Int_t whichbit = (Int_t)pow(2,ivk); |
988 |
|
return (whichbit&VKmask[iv])!=0; |
989 |
|
} |
990 |
|
/** |
991 |
|
* The method returns false if the viking-chip was masked |
992 |
|
* for this event due to common-noise computation failure. |
993 |
|
* @param iv Tracker view (0-11) |
994 |
|
* @param ivk Viking-chip number (0-23) |
995 |
|
*/ |
996 |
|
Bool_t TrkLevel2::GetVKFlag(int iv, int ivk){ |
997 |
|
Int_t whichbit = (Int_t)pow(2,ivk); |
998 |
|
return (whichbit&VKflag[iv])!=0; |
999 |
|
} |
1000 |
|
/** |
1001 |
|
* The method returns true if the viking-chip was masked, either |
1002 |
|
* forced (see TrkLevel2::GetVKMask(int,int)) or |
1003 |
|
* for this event only (TrkLevel2::GetVKFlag(int,int)). |
1004 |
|
* @param iv Tracker view (0-11) |
1005 |
|
* @param ivk Viking-chip number (0-23) |
1006 |
|
*/ |
1007 |
|
Bool_t TrkLevel2::IsMaskedVK(int iv, int ivk){ |
1008 |
|
return !(GetVKMask(iv,ivk)&&GetVKFlag(iv,ivk) ); |
1009 |
|
}; |
1010 |
|
|
1011 |
|
//-------------------------------------- |
1012 |
|
// |
1013 |
|
// |
1014 |
|
//-------------------------------------- |
1015 |
|
/** |
1016 |
|
* Fills a TrkLevel2 object with values from a struct cTrkLevel2 (to get data from F77 common). |
1017 |
|
* Ref to Level1 data (clusters) is also set. If l1==NULL no references are set. |
1018 |
|
* (NB It make sense to set references only if events are stored in a tree that contains also the Level1 branch) |
1019 |
*/ |
*/ |
1020 |
void TrkLevel2::SetFromLevel2Struct(cTrkLevel2 *l2){ |
void TrkLevel2::SetFromLevel2Struct(cTrkLevel2 *l2, TrkLevel1 *l1){ |
1021 |
|
|
1022 |
|
// cout << "void TrkLevel2::SetFromLevel2Struct(cTrkLevel2 *l2, TrkLevel1 *l1)"<<endl; |
1023 |
|
Clear(); |
1024 |
|
|
1025 |
// temporary objects: |
// temporary objects: |
1026 |
TrkSinglet* t_singlet = new TrkSinglet(); |
TrkSinglet* t_singlet = new TrkSinglet(); |
1027 |
TrkTrack* t_track = new TrkTrack(); |
TrkTrack* t_track = new TrkTrack(); |
1028 |
|
|
1029 |
// **** general variables **** |
// ----------------- |
1030 |
// good2 = l2->good2; |
// general variables |
1031 |
|
// ----------------- |
1032 |
for(Int_t i=0; i<12 ; i++){ |
for(Int_t i=0; i<12 ; i++){ |
1033 |
// crc[i] = l2->crc[i]; |
good[i] = l2->good[i]; |
1034 |
good[i] = l2->good[i]; |
VKmask[i]=0; |
1035 |
|
VKflag[i]=0; |
1036 |
|
for(Int_t ii=0; ii<24 ; ii++){ |
1037 |
|
Int_t setbit = (Int_t)pow(2,ii); |
1038 |
|
if( l2->vkflag[ii][i]!=-1 )VKmask[i]=VKmask[i]|setbit; |
1039 |
|
if( l2->vkflag[ii][i]!=0 )VKflag[i]=VKflag[i]|setbit; |
1040 |
}; |
}; |
1041 |
// *** TRACKS *** |
}; |
1042 |
|
// -------------- |
1043 |
|
// *** TRACKS *** |
1044 |
|
// -------------- |
1045 |
|
if(!Track) Track = new TClonesArray("TrkTrack"); |
1046 |
TClonesArray &t = *Track; |
TClonesArray &t = *Track; |
1047 |
|
|
1048 |
for(int i=0; i<l2->ntrk; i++){ |
for(int i=0; i<l2->ntrk; i++){ |
1049 |
t_track->seqno = i;// NBNBNBNB deve sempre essere = i |
t_track->seqno = i;// NBNBNBNB deve sempre essere = i |
1050 |
t_track->image = l2->image[i]-1; |
t_track->image = l2->image[i]-1; |
1051 |
// cout << "track "<<i<<t_track->seqno << t_track->image<<endl; |
t_track->chi2 = l2->chi2_nt[i]; |
1052 |
t_track->chi2 = l2->chi2_nt[i]; |
t_track->nstep = l2->nstep_nt[i]; |
1053 |
t_track->nstep = l2->nstep_nt[i]; |
for(int it1=0;it1<5;it1++){ |
1054 |
for(int it1=0;it1<5;it1++){ |
t_track->al[it1] = l2->al_nt[i][it1]; |
1055 |
t_track->al[it1] = l2->al_nt[i][it1]; |
for(int it2=0;it2<5;it2++) |
1056 |
for(int it2=0;it2<5;it2++) |
t_track->coval[it1][it2] = l2->coval[i][it2][it1]; |
1057 |
t_track->coval[it1][it2] = l2->coval[i][it2][it1]; |
}; |
1058 |
}; |
for(int ip=0;ip<6;ip++){ |
1059 |
for(int ip=0;ip<6;ip++){ |
// --------------------------------- |
1060 |
t_track->xgood[ip] = l2->xgood_nt[i][ip]; |
// new implementation of xgood/ygood |
1061 |
t_track->ygood[ip] = l2->ygood_nt[i][ip]; |
// --------------------------------- |
1062 |
t_track->xm[ip] = l2->xm_nt[i][ip]; |
t_track->xgood[ip] = l2->cltrx[i][ip]; //cluster ID |
1063 |
t_track->ym[ip] = l2->ym_nt[i][ip]; |
t_track->ygood[ip] = l2->cltry[i][ip]; //cluster ID |
1064 |
t_track->zm[ip] = l2->zm_nt[i][ip]; |
t_track->xgood[ip] += 10000000*l2->ls[i][ip]; // ladder+sensor |
1065 |
t_track->resx[ip] = l2->resx_nt[i][ip]; |
t_track->ygood[ip] += 10000000*l2->ls[i][ip]; // ladder+sensor |
1066 |
t_track->resy[ip] = l2->resy_nt[i][ip]; |
if(l2->xbad[i][ip]>0)t_track->xgood[ip]=-t_track->xgood[ip]; |
1067 |
t_track->xv[ip] = l2->xv_nt[i][ip]; |
if(l2->ybad[i][ip]>0)t_track->ygood[ip]=-t_track->ygood[ip]; |
1068 |
t_track->yv[ip] = l2->yv_nt[i][ip]; |
// if(l2->xbad[i][ip]>0 || l2->ybad[i][ip]>0){ |
1069 |
t_track->zv[ip] = l2->zv_nt[i][ip]; |
// if(l2->dedx_x[i][ip]<0 || l2->dedx_y[i][ip]<0){ |
1070 |
t_track->axv[ip] = l2->axv_nt[i][ip]; |
// cout << ip << " - "<< l2->cltrx[i][ip] << " "<<l2->cltry[i][ip]<<" "<<l2->ls[i][ip]<<endl; |
1071 |
t_track->ayv[ip] = l2->ayv_nt[i][ip]; |
// cout << ip << " - "<<t_track->xgood[ip]<<" "<<t_track->ygood[ip]<<endl; |
1072 |
t_track->dedx_x[ip] = l2->dedx_x[i][ip]; |
// cout << ip << " - "<<t_track->GetClusterX_ID(ip)<<" "<<t_track->GetClusterY_ID(ip)<<" "<<t_track->GetLadder(ip)<<" "<<t_track->GetSensor(ip)<<endl; |
1073 |
t_track->dedx_y[ip] = l2->dedx_y[i][ip]; |
// cout << ip << " - "<<t_track->BadClusterX(ip)<<" "<<t_track->BadClusterY(ip)<<endl; |
1074 |
// t_track->clx[ip] = 0; |
// cout << ip << " - "<<t_track->SaturatedClusterX(ip)<<" "<<t_track->SaturatedClusterY(ip)<<endl; |
1075 |
// t_track->cly[ip] = 0; |
// } |
1076 |
}; |
t_track->xm[ip] = l2->xm_nt[i][ip]; |
1077 |
new(t[i]) TrkTrack(*t_track); |
t_track->ym[ip] = l2->ym_nt[i][ip]; |
1078 |
t_track->Clear(); |
t_track->zm[ip] = l2->zm_nt[i][ip]; |
1079 |
|
t_track->resx[ip] = l2->resx_nt[i][ip]; |
1080 |
|
t_track->resy[ip] = l2->resy_nt[i][ip]; |
1081 |
|
t_track->tailx[ip] = l2->tailx[i][ip]; |
1082 |
|
t_track->taily[ip] = l2->taily[i][ip]; |
1083 |
|
t_track->xv[ip] = l2->xv_nt[i][ip]; |
1084 |
|
t_track->yv[ip] = l2->yv_nt[i][ip]; |
1085 |
|
t_track->zv[ip] = l2->zv_nt[i][ip]; |
1086 |
|
t_track->axv[ip] = l2->axv_nt[i][ip]; |
1087 |
|
t_track->ayv[ip] = l2->ayv_nt[i][ip]; |
1088 |
|
t_track->dedx_x[ip] = l2->dedx_x[i][ip]; |
1089 |
|
t_track->dedx_y[ip] = l2->dedx_y[i][ip]; |
1090 |
|
//----------------------------------------------------- |
1091 |
|
//----------------------------------------------------- |
1092 |
|
//----------------------------------------------------- |
1093 |
|
//----------------------------------------------------- |
1094 |
|
}; |
1095 |
|
// if(t_track->IsSaturated())t_track->Dump(); |
1096 |
|
new(t[i]) TrkTrack(*t_track); |
1097 |
|
t_track->Clear(); |
1098 |
}; |
}; |
1099 |
|
|
1100 |
|
// ---------------- |
1101 |
// *** SINGLETS *** |
// *** SINGLETS *** |
1102 |
|
// ---------------- |
1103 |
|
if(!SingletX)SingletX = new TClonesArray("TrkSinglet"); |
1104 |
TClonesArray &sx = *SingletX; |
TClonesArray &sx = *SingletX; |
1105 |
for(int i=0; i<l2->nclsx; i++){ |
for(int i=0; i<l2->nclsx; i++){ |
1106 |
t_singlet->plane = l2->planex[i]; |
t_singlet->plane = l2->planex[i]; |
1107 |
t_singlet->coord[0] = l2->xs[i][0]; |
t_singlet->coord[0] = l2->xs[i][0]; |
1108 |
t_singlet->coord[1] = l2->xs[i][1]; |
t_singlet->coord[1] = l2->xs[i][1]; |
1109 |
t_singlet->sgnl = l2->signlxs[i]; |
t_singlet->sgnl = l2->signlxs[i]; |
1110 |
// t_singlet->cls = 0; |
//----------------------------------------------------- |
1111 |
new(sx[i]) TrkSinglet(*t_singlet); |
// if(l1) t_singlet->cls = l1->GetCluster(l2->clsx[i]-1); |
1112 |
t_singlet->Clear(); |
//----------------------------------------------------- |
1113 |
|
new(sx[i]) TrkSinglet(*t_singlet); |
1114 |
|
t_singlet->Clear(); |
1115 |
} |
} |
1116 |
|
if(!SingletY)SingletY = new TClonesArray("TrkSinglet"); |
1117 |
TClonesArray &sy = *SingletY; |
TClonesArray &sy = *SingletY; |
1118 |
for(int i=0; i<l2->nclsy; i++){ |
for(int i=0; i<l2->nclsy; i++){ |
1119 |
t_singlet->plane = l2->planey[i]; |
t_singlet->plane = l2->planey[i]; |
1120 |
t_singlet->coord[0] = l2->ys[i][0]; |
t_singlet->coord[0] = l2->ys[i][0]; |
1121 |
t_singlet->coord[1] = l2->ys[i][1]; |
t_singlet->coord[1] = l2->ys[i][1]; |
1122 |
t_singlet->sgnl = l2->signlys[i]; |
t_singlet->sgnl = l2->signlys[i]; |
1123 |
// t_singlet->cls = 0; |
//----------------------------------------------------- |
1124 |
new(sy[i]) TrkSinglet(*t_singlet); |
// if(l1) t_singlet->cls = l1->GetCluster(l2->clsy[i]-1); |
1125 |
t_singlet->Clear(); |
//----------------------------------------------------- |
1126 |
}; |
new(sy[i]) TrkSinglet(*t_singlet); |
1127 |
|
t_singlet->Clear(); |
1128 |
delete t_track; |
}; |
|
delete t_singlet; |
|
|
} |
|
|
//-------------------------------------- |
|
|
// |
|
|
// |
|
|
//-------------------------------------- |
|
|
/** |
|
|
* Fills a TrkLevel2 object with values from a struct cTrkLevel2 (to get data from F77 common). |
|
|
* Ref to Level1 data (clusters) is also set. |
|
|
*/ |
|
|
void TrkLevel2::SetFromLevel2Struct(cTrkLevel2 *l2, TrkLevel1 *l1){ |
|
|
|
|
|
// temporary objects: |
|
|
TrkSinglet* t_singlet = new TrkSinglet(); |
|
|
TrkTrack* t_track = new TrkTrack(); |
|
|
// general variables |
|
|
// good2 = l2->good2; |
|
|
for(Int_t i=0; i<12 ; i++){ |
|
|
// crc[i] = l2->crc[i]; |
|
|
good[i] = l2->good[i]; |
|
|
}; |
|
|
// *** TRACKS *** |
|
|
TClonesArray &t = *Track; |
|
|
for(int i=0; i<l2->ntrk; i++){ |
|
|
t_track->seqno = i;// NBNBNBNB deve sempre essere = i |
|
|
t_track->image = l2->image[i]-1; |
|
|
// cout << "track "<<i<<t_track->seqno << t_track->image<<endl; |
|
|
t_track->chi2 = l2->chi2_nt[i]; |
|
|
t_track->nstep = l2->nstep_nt[i]; |
|
|
for(int it1=0;it1<5;it1++){ |
|
|
t_track->al[it1] = l2->al_nt[i][it1]; |
|
|
for(int it2=0;it2<5;it2++) |
|
|
t_track->coval[it1][it2] = l2->coval[i][it2][it1]; |
|
|
}; |
|
|
for(int ip=0;ip<6;ip++){ |
|
|
t_track->xgood[ip] = l2->xgood_nt[i][ip]; |
|
|
t_track->ygood[ip] = l2->ygood_nt[i][ip]; |
|
|
t_track->xm[ip] = l2->xm_nt[i][ip]; |
|
|
t_track->ym[ip] = l2->ym_nt[i][ip]; |
|
|
t_track->zm[ip] = l2->zm_nt[i][ip]; |
|
|
t_track->resx[ip] = l2->resx_nt[i][ip]; |
|
|
t_track->resy[ip] = l2->resy_nt[i][ip]; |
|
|
t_track->xv[ip] = l2->xv_nt[i][ip]; |
|
|
t_track->yv[ip] = l2->yv_nt[i][ip]; |
|
|
t_track->zv[ip] = l2->zv_nt[i][ip]; |
|
|
t_track->axv[ip] = l2->axv_nt[i][ip]; |
|
|
t_track->ayv[ip] = l2->ayv_nt[i][ip]; |
|
|
t_track->dedx_x[ip] = l2->dedx_x[i][ip]; |
|
|
t_track->dedx_y[ip] = l2->dedx_y[i][ip]; |
|
|
// cout << "traccia "<<i<<" -- "<< ip << " "<< l2->cltrx[i][ip] <<" "<< l2->cltry[i][ip] <<" "<< t_track->xgood[ip] << t_track->ygood[ip]<<endl; |
|
|
//----------------------------------------------------- |
|
|
// t_track->clx[ip] = l1->GetCluster(l2->cltrx[i][ip]-1); |
|
|
// t_track->cly[ip] = l1->GetCluster(l2->cltry[i][ip]-1); |
|
|
if(t_track->xgood[ip])t_track->clx->AddAt(l1->GetCluster(l2->cltrx[i][ip]-1),ip); |
|
|
if(t_track->ygood[ip])t_track->cly->AddAt(l1->GetCluster(l2->cltry[i][ip]-1),ip); |
|
|
// if(t_track->ygood[ip])cout<<" i "<<i<<" ip "<<ip<<" l2->cltry[i][ip] "<<l2->cltry[i][ip]<< " l1->GetCluster(l2->cltry[i][ip]-1) "<<l1->GetCluster(l2->cltry[i][ip]-1)<<endl; |
|
|
// if(t_track->xgood[ip])cout<<" i "<<i<<" ip "<<ip<<" l2->cltrx[i][ip] "<<l2->cltrx[i][ip]<< " l1->GetCluster(l2->cltrx[i][ip]-1) "<<l1->GetCluster(l2->cltrx[i][ip]-1)<<endl; |
|
|
//----------------------------------------------------- |
|
|
}; |
|
|
new(t[i]) TrkTrack(*t_track); |
|
|
t_track->Clear(); |
|
|
}; |
|
|
// *** SINGLETS *** |
|
|
TClonesArray &sx = *SingletX; |
|
|
for(int i=0; i<l2->nclsx; i++){ |
|
|
t_singlet->plane = l2->planex[i]; |
|
|
t_singlet->coord[0] = l2->xs[i][0]; |
|
|
t_singlet->coord[1] = l2->xs[i][1]; |
|
|
t_singlet->sgnl = l2->signlxs[i]; |
|
|
//----------------------------------------------------- |
|
|
// cout << "singolo x "<<i<<" -- "<< l2->clsx[i] <<endl; |
|
|
t_singlet->cls = l1->GetCluster(l2->clsx[i]-1); |
|
|
// cout<<" i "<<i<<" l2->clsx[i] "<<l2->clsx[i]<< " l1->GetCluster(l2->clsx[i]-1) "<<l1->GetCluster(l2->clsx[i]-1)<<endl; |
|
|
//----------------------------------------------------- |
|
|
new(sx[i]) TrkSinglet(*t_singlet); |
|
|
t_singlet->Clear(); |
|
|
} |
|
|
TClonesArray &sy = *SingletY; |
|
|
for(int i=0; i<l2->nclsy; i++){ |
|
|
t_singlet->plane = l2->planey[i]; |
|
|
t_singlet->coord[0] = l2->ys[i][0]; |
|
|
t_singlet->coord[1] = l2->ys[i][1]; |
|
|
t_singlet->sgnl = l2->signlys[i]; |
|
|
//----------------------------------------------------- |
|
|
// cout << "singolo y "<<i<<" -- "<< l2->clsy[i] <<endl; |
|
|
t_singlet->cls = l1->GetCluster(l2->clsy[i]-1); |
|
|
// cout<<" i "<<i<<" l2->clsy[i] "<<l2->clsy[i]<< " l1->GetCluster(l2->clsy[i]-1) "<<l1->GetCluster(l2->clsy[i]-1)<<endl; |
|
|
//----------------------------------------------------- |
|
|
new(sy[i]) TrkSinglet(*t_singlet); |
|
|
t_singlet->Clear(); |
|
|
}; |
|
1129 |
|
|
1130 |
delete t_track; |
delete t_track; |
1131 |
delete t_singlet; |
delete t_singlet; |
1132 |
} |
} |
1133 |
/** |
/** |
1134 |
* Fills a struct cTrkLevel2 with values from a TrkLevel2 object (to put data into a F77 common). |
* Fills a struct cTrkLevel2 with values from a TrkLevel2 object (to put data into a F77 common). |
1144 |
}; |
}; |
1145 |
// *** TRACKS *** |
// *** TRACKS *** |
1146 |
|
|
1147 |
l2->ntrk = Track->GetEntries(); |
if(Track){ |
1148 |
for(Int_t i=0;i<l2->ntrk;i++){ |
l2->ntrk = Track->GetEntries(); |
1149 |
l2->image[i] = 1 + ((TrkTrack *)Track->At(i))->image; |
for(Int_t i=0;i<l2->ntrk;i++){ |
1150 |
l2->chi2_nt[i] = ((TrkTrack *)Track->At(i))->chi2; |
l2->image[i] = 1 + ((TrkTrack *)Track->At(i))->image; |
1151 |
l2->nstep_nt[i] = ((TrkTrack *)Track->At(i))->nstep; |
l2->chi2_nt[i] = ((TrkTrack *)Track->At(i))->chi2; |
1152 |
for(int it1=0;it1<5;it1++){ |
l2->nstep_nt[i] = ((TrkTrack *)Track->At(i))->nstep; |
1153 |
l2->al_nt[i][it1] = ((TrkTrack *)Track->At(i))->al[it1]; |
for(int it1=0;it1<5;it1++){ |
1154 |
for(int it2=0;it2<5;it2++) |
l2->al_nt[i][it1] = ((TrkTrack *)Track->At(i))->al[it1]; |
1155 |
l2->coval[i][it2][it1] = ((TrkTrack *)Track->At(i))->coval[it1][it2]; |
for(int it2=0;it2<5;it2++) |
1156 |
}; |
l2->coval[i][it2][it1] = ((TrkTrack *)Track->At(i))->coval[it1][it2]; |
1157 |
for(int ip=0;ip<6;ip++){ |
}; |
1158 |
l2->xgood_nt[i][ip] = ((TrkTrack *)Track->At(i))->xgood[ip]; |
for(int ip=0;ip<6;ip++){ |
1159 |
l2->ygood_nt[i][ip] = ((TrkTrack *)Track->At(i))->ygood[ip]; |
l2->xgood_nt[i][ip] = ((TrkTrack *)Track->At(i))->XGood(ip); |
1160 |
l2->xm_nt[i][ip] = ((TrkTrack *)Track->At(i))->xm[ip]; |
l2->ygood_nt[i][ip] = ((TrkTrack *)Track->At(i))->YGood(ip); |
1161 |
l2->ym_nt[i][ip] = ((TrkTrack *)Track->At(i))->ym[ip]; |
l2->xm_nt[i][ip] = ((TrkTrack *)Track->At(i))->xm[ip]; |
1162 |
l2->zm_nt[i][ip] = ((TrkTrack *)Track->At(i))->zm[ip]; |
l2->ym_nt[i][ip] = ((TrkTrack *)Track->At(i))->ym[ip]; |
1163 |
l2->resx_nt[i][ip] = ((TrkTrack *)Track->At(i))->resx[ip]; |
l2->zm_nt[i][ip] = ((TrkTrack *)Track->At(i))->zm[ip]; |
1164 |
l2->resy_nt[i][ip] = ((TrkTrack *)Track->At(i))->resy[ip]; |
l2->resx_nt[i][ip] = ((TrkTrack *)Track->At(i))->resx[ip]; |
1165 |
l2->xv_nt[i][ip] = ((TrkTrack *)Track->At(i))->xv[ip]; |
l2->resy_nt[i][ip] = ((TrkTrack *)Track->At(i))->resy[ip]; |
1166 |
l2->yv_nt[i][ip] = ((TrkTrack *)Track->At(i))->yv[ip]; |
l2->tailx[i][ip] = ((TrkTrack *)Track->At(i))->tailx[ip]; |
1167 |
l2->zv_nt[i][ip] = ((TrkTrack *)Track->At(i))->zv[ip]; |
l2->taily[i][ip] = ((TrkTrack *)Track->At(i))->taily[ip]; |
1168 |
l2->axv_nt[i][ip] = ((TrkTrack *)Track->At(i))->axv[ip]; |
l2->xv_nt[i][ip] = ((TrkTrack *)Track->At(i))->xv[ip]; |
1169 |
l2->ayv_nt[i][ip] = ((TrkTrack *)Track->At(i))->ayv[ip]; |
l2->yv_nt[i][ip] = ((TrkTrack *)Track->At(i))->yv[ip]; |
1170 |
l2->dedx_x[i][ip] = ((TrkTrack *)Track->At(i))->dedx_x[ip]; |
l2->zv_nt[i][ip] = ((TrkTrack *)Track->At(i))->zv[ip]; |
1171 |
l2->dedx_y[i][ip] = ((TrkTrack *)Track->At(i))->dedx_y[ip]; |
l2->axv_nt[i][ip] = ((TrkTrack *)Track->At(i))->axv[ip]; |
1172 |
}; |
l2->ayv_nt[i][ip] = ((TrkTrack *)Track->At(i))->ayv[ip]; |
1173 |
|
l2->dedx_x[i][ip] = ((TrkTrack *)Track->At(i))->dedx_x[ip]; |
1174 |
|
l2->dedx_y[i][ip] = ((TrkTrack *)Track->At(i))->dedx_y[ip]; |
1175 |
|
}; |
1176 |
|
} |
1177 |
} |
} |
|
|
|
1178 |
// *** SINGLETS *** |
// *** SINGLETS *** |
1179 |
l2->nclsx = SingletX->GetEntries(); |
if(SingletX){ |
1180 |
for(Int_t i=0;i<l2->nclsx;i++){ |
l2->nclsx = SingletX->GetEntries(); |
1181 |
l2->planex[i] = ((TrkSinglet *)SingletX->At(i))->plane; |
for(Int_t i=0;i<l2->nclsx;i++){ |
1182 |
l2->xs[i][0] = ((TrkSinglet *)SingletX->At(i))->coord[0]; |
l2->planex[i] = ((TrkSinglet *)SingletX->At(i))->plane; |
1183 |
l2->xs[i][1] = ((TrkSinglet *)SingletX->At(i))->coord[1]; |
l2->xs[i][0] = ((TrkSinglet *)SingletX->At(i))->coord[0]; |
1184 |
l2->signlxs[i] = ((TrkSinglet *)SingletX->At(i))->sgnl; |
l2->xs[i][1] = ((TrkSinglet *)SingletX->At(i))->coord[1]; |
1185 |
} |
l2->signlxs[i] = ((TrkSinglet *)SingletX->At(i))->sgnl; |
1186 |
l2->nclsy = SingletY->GetEntries(); |
} |
1187 |
for(Int_t i=0;i<l2->nclsy;i++){ |
} |
1188 |
l2->planey[i] = ((TrkSinglet *)SingletY->At(i))->plane; |
|
1189 |
l2->ys[i][0] = ((TrkSinglet *)SingletY->At(i))->coord[0]; |
if(SingletY){ |
1190 |
l2->ys[i][1] = ((TrkSinglet *)SingletY->At(i))->coord[1]; |
l2->nclsy = SingletY->GetEntries(); |
1191 |
l2->signlys[i] = ((TrkSinglet *)SingletY->At(i))->sgnl; |
for(Int_t i=0;i<l2->nclsy;i++){ |
1192 |
|
l2->planey[i] = ((TrkSinglet *)SingletY->At(i))->plane; |
1193 |
|
l2->ys[i][0] = ((TrkSinglet *)SingletY->At(i))->coord[0]; |
1194 |
|
l2->ys[i][1] = ((TrkSinglet *)SingletY->At(i))->coord[1]; |
1195 |
|
l2->signlys[i] = ((TrkSinglet *)SingletY->At(i))->sgnl; |
1196 |
|
} |
1197 |
} |
} |
1198 |
} |
} |
1199 |
//-------------------------------------- |
//-------------------------------------- |
1201 |
// |
// |
1202 |
//-------------------------------------- |
//-------------------------------------- |
1203 |
void TrkLevel2::Clear(){ |
void TrkLevel2::Clear(){ |
|
// good2 = -1; |
|
1204 |
for(Int_t i=0; i<12 ; i++){ |
for(Int_t i=0; i<12 ; i++){ |
1205 |
// crc[i] = -1; |
good[i] = -1; |
1206 |
good[i] = -1; |
VKflag[i] = 0; |
1207 |
}; |
VKmask[i] = 0; |
1208 |
/* Track->RemoveAll(); |
}; |
1209 |
SingletX->RemoveAll(); |
// if(Track)Track->Clear("C"); |
1210 |
SingletY->RemoveAll();*/ |
// if(SingletX)SingletX->Clear("C"); |
1211 |
// modify to avoid memory leakage |
// if(SingletY)SingletY->Clear("C"); |
1212 |
Track->Clear(); |
if(Track)Track->Delete(); |
1213 |
SingletX->Clear(); |
if(SingletX)SingletX->Delete(); |
1214 |
SingletY->Clear(); |
if(SingletY)SingletY->Delete(); |
1215 |
} |
} |
1216 |
//-------------------------------------- |
// //-------------------------------------- |
1217 |
// |
// // |
1218 |
// |
// // |
1219 |
//-------------------------------------- |
// //-------------------------------------- |
1220 |
void TrkLevel2::Delete(){ |
void TrkLevel2::Delete(){ |
1221 |
|
|
1222 |
Clear(); |
// cout << "void TrkLevel2::Delete()"<<endl; |
1223 |
Track->Delete(); |
Clear(); |
1224 |
SingletX->Delete(); |
if(Track)delete Track; |
1225 |
SingletY->Delete(); |
if(SingletX)delete SingletX; |
1226 |
|
if(SingletY)delete SingletY; |
1227 |
|
|
1228 |
} |
} |
1229 |
//-------------------------------------- |
//-------------------------------------- |
1230 |
// |
// |
1236 |
*/ |
*/ |
1237 |
TRefArray *TrkLevel2::GetTracks_NFitSorted(){ |
TRefArray *TrkLevel2::GetTracks_NFitSorted(){ |
1238 |
|
|
1239 |
TRefArray *sorted = new TRefArray(); |
if(!Track)return 0; |
1240 |
|
|
1241 |
|
TRefArray *sorted = new TRefArray(); |
1242 |
|
|
1243 |
TClonesArray &t = *Track; |
TClonesArray &t = *Track; |
1244 |
// TClonesArray &ts = *PhysicalTrack; |
// TClonesArray &ts = *PhysicalTrack; |
1245 |
int N = ntrk(); |
int N = ntrk(); |
1246 |
vector<int> m(N); for(int i=0; i<N; i++)m[i]=1; |
vector<int> m(N); for(int i=0; i<N; i++)m[i]=1; |
1247 |
// int m[50]; for(int i=0; i<N; i++)m[i]=1; |
// int m[50]; for(int i=0; i<N; i++)m[i]=1; |
1248 |
|
|
1249 |
int indo=0; |
int indo=0; |
1250 |
int indi=0; |
int indi=0; |
1251 |
while(N != 0){ |
while(N > 0){ |
1252 |
int nfit =0; |
// while(N != 0){ |
1253 |
float chi2ref = numeric_limits<float>::max(); |
int nfit =0; |
1254 |
|
float chi2ref = numeric_limits<float>::max(); |
1255 |
|
|
1256 |
// first loop to search maximum num. of fit points |
// first loop to search maximum num. of fit points |
1257 |
for(int i=0; i < ntrk(); i++){ |
for(int i=0; i < ntrk(); i++){ |
1258 |
if( ((TrkTrack *)t[i])->GetNtot() >= nfit && m[i]==1){ |
if( ((TrkTrack *)t[i])->GetNtot() >= nfit && m[i]==1){ |
1259 |
nfit = ((TrkTrack *)t[i])->GetNtot(); |
nfit = ((TrkTrack *)t[i])->GetNtot(); |
1260 |
} |
} |
1261 |
} |
} |
1262 |
//second loop to search minimum chi2 among selected |
//second loop to search minimum chi2 among selected |
1263 |
for(int i=0; i<this->ntrk(); i++){ |
for(int i=0; i<ntrk(); i++){ |
1264 |
Float_t chi2 = ((TrkTrack *)t[i])->chi2; |
Float_t chi2 = ((TrkTrack *)t[i])->chi2; |
1265 |
if(chi2 < 0) chi2 = chi2*1000; |
if(chi2 < 0) chi2 = -chi2*1000; |
1266 |
if( chi2 < chi2ref |
if( chi2 < chi2ref |
1267 |
&& ((TrkTrack *)t[i])->GetNtot() == nfit |
&& ((TrkTrack *)t[i])->GetNtot() == nfit |
1268 |
&& m[i]==1){ |
&& m[i]==1){ |
1269 |
chi2ref = ((TrkTrack *)t[i])->chi2; |
chi2ref = ((TrkTrack *)t[i])->chi2; |
1270 |
indi = i; |
indi = i; |
1271 |
}; |
}; |
1272 |
}; |
}; |
1273 |
if( ((TrkTrack *)t[indi])->HasImage() ){ |
if( ((TrkTrack *)t[indi])->HasImage() ){ |
1274 |
m[((TrkTrack *)t[indi])->image] = 0; |
m[((TrkTrack *)t[indi])->image] = 0; |
1275 |
N--; |
N--; |
1276 |
|
|
1277 |
// cout << "i** "<< ((TrkTrack *)t[indi])->image << " " << nfiti <<" "<<chi2i<<endl; |
// cout << "i** "<< ((TrkTrack *)t[indi])->image << " " << nfiti <<" "<<chi2i<<endl; |
1278 |
}; |
}; |
1279 |
sorted->Add( (TrkTrack*)t[indi] ); |
sorted->Add( (TrkTrack*)t[indi] ); |
1280 |
|
|
1281 |
m[indi] = 0; |
m[indi] = 0; |
1282 |
// cout << "SORTED "<< indo << " "<< indi << " "<< N << endl; |
// cout << "SORTED "<< indo << " "<< indi << " "<< N << " "<<((TrkTrack *)t[indi])->image<<" "<<chi2ref<<endl; |
1283 |
N--; |
N--; |
1284 |
indo++; |
indo++; |
1285 |
} |
} |
1286 |
m.clear(); |
m.clear(); |
1287 |
// cout << "GetTracks_NFitSorted(it): Done"<< endl; |
// cout << "GetTracks_NFitSorted(it): Done"<< endl; |
1288 |
|
|
1289 |
return sorted; |
return sorted; |
1290 |
// return PhysicalTrack; |
// return PhysicalTrack; |
1291 |
} |
} |
1292 |
//-------------------------------------- |
//-------------------------------------- |
1305 |
cout << " Stored tracks ntrk() = "<< this->ntrk() << endl; |
cout << " Stored tracks ntrk() = "<< this->ntrk() << endl; |
1306 |
return 0; |
return 0; |
1307 |
} |
} |
1308 |
|
if(!Track){ |
1309 |
|
cout << "TrkTrack *TrkLevel2::GetStoredTrack(int is) >> (TClonesArray*) Track ==0 "<<endl; |
1310 |
|
}; |
1311 |
TClonesArray &t = *(Track); |
TClonesArray &t = *(Track); |
1312 |
TrkTrack *track = (TrkTrack*)t[is]; |
TrkTrack *track = (TrkTrack*)t[is]; |
1313 |
return track; |
return track; |
1327 |
cout << " Stored x-singlets nclsx() = "<< this->nclsx() << endl; |
cout << " Stored x-singlets nclsx() = "<< this->nclsx() << endl; |
1328 |
return 0; |
return 0; |
1329 |
} |
} |
1330 |
|
if(!SingletX)return 0; |
1331 |
TClonesArray &t = *(SingletX); |
TClonesArray &t = *(SingletX); |
1332 |
TrkSinglet *singlet = (TrkSinglet*)t[is]; |
TrkSinglet *singlet = (TrkSinglet*)t[is]; |
1333 |
return singlet; |
return singlet; |
1347 |
cout << " Stored y-singlets nclsy() = "<< this->nclsx() << endl; |
cout << " Stored y-singlets nclsy() = "<< this->nclsx() << endl; |
1348 |
return 0; |
return 0; |
1349 |
} |
} |
1350 |
|
if(!SingletY)return 0; |
1351 |
TClonesArray &t = *(SingletY); |
TClonesArray &t = *(SingletY); |
1352 |
TrkSinglet *singlet = (TrkSinglet*)t[is]; |
TrkSinglet *singlet = (TrkSinglet*)t[is]; |
1353 |
return singlet; |
return singlet; |
1370 |
} |
} |
1371 |
|
|
1372 |
TRefArray *sorted = GetTracks(); //TEMPORANEO |
TRefArray *sorted = GetTracks(); //TEMPORANEO |
1373 |
|
if(!sorted)return 0; |
1374 |
TrkTrack *track = (TrkTrack*)sorted->At(it); |
TrkTrack *track = (TrkTrack*)sorted->At(it); |
1375 |
sorted->Delete(); |
sorted->Clear(); |
1376 |
|
delete sorted; |
1377 |
return track; |
return track; |
1378 |
} |
} |
1379 |
/** |
/** |
1382 |
Int_t TrkLevel2::GetNTracks(){ |
Int_t TrkLevel2::GetNTracks(){ |
1383 |
|
|
1384 |
Float_t ntot=0; |
Float_t ntot=0; |
1385 |
|
if(!Track)return 0; |
1386 |
TClonesArray &t = *Track; |
TClonesArray &t = *Track; |
1387 |
for(int i=0; i<ntrk(); i++) { |
for(int i=0; i<ntrk(); i++) { |
1388 |
if( ((TrkTrack *)t[i])->GetImageSeqNo() == -1 ) ntot+=1.; |
if( ((TrkTrack *)t[i])->GetImageSeqNo() == -1 ) ntot+=1.; |
1401 |
*/ |
*/ |
1402 |
TrkTrack *TrkLevel2::GetTrackImage(int it){ |
TrkTrack *TrkLevel2::GetTrackImage(int it){ |
1403 |
|
|
1404 |
if(it >= this->GetNTracks()){ |
if(it >= this->GetNTracks()){ |
1405 |
cout << "** TrkLevel2 ** Track "<< it << "does not exits! " << endl; |
cout << "** TrkLevel2 ** Track "<< it << "does not exits! " << endl; |
1406 |
cout << " Physical tracks GetNTracks() = "<< this->ntrk() << endl; |
cout << " Physical tracks GetNTracks() = "<< this->ntrk() << endl; |
1407 |
return 0; |
return 0; |
1408 |
} |
} |
1409 |
|
|
1410 |
TRefArray* sorted = GetTracks(); //TEMPORANEO |
TRefArray* sorted = GetTracks(); //TEMPORANEO |
1411 |
TrkTrack *track = (TrkTrack*)sorted->At(it); |
if(!sorted)return 0; |
1412 |
|
TrkTrack *track = (TrkTrack*)sorted->At(it); |
1413 |
|
|
1414 |
if(!track->HasImage()){ |
if(!track->HasImage()){ |
1415 |
cout << "** TrkLevel2 ** Track "<< it << "does not have image! " << endl; |
cout << "** TrkLevel2 ** Track "<< it << "does not have image! " << endl; |
1416 |
return 0; |
return 0; |
1417 |
} |
} |
1418 |
TrkTrack *image = (TrkTrack*)(*Track)[track->image]; |
if(!Track)return 0; |
1419 |
|
TrkTrack *image = (TrkTrack*)(*Track)[track->image]; |
1420 |
|
|
1421 |
sorted->Delete(); |
sorted->Delete(); |
1422 |
|
delete sorted; |
1423 |
return image; |
|
1424 |
|
return image; |
1425 |
|
|
1426 |
} |
} |
1427 |
//-------------------------------------- |
//-------------------------------------- |
1434 |
*/ |
*/ |
1435 |
void TrkLevel2::LoadField(TString path){ |
void TrkLevel2::LoadField(TString path){ |
1436 |
// |
// |
1437 |
strcpy(path_.path,path.Data()); |
// strcpy(path_.path,path.Data()); |
1438 |
path_.pathlen = path.Length(); |
// path_.pathlen = path.Length(); |
1439 |
path_.error = 0; |
// path_.error = 0; |
1440 |
readb_(); |
// readb_(); |
1441 |
|
|
1442 |
|
TrkParams::SetTrackingMode(); |
1443 |
|
TrkParams::SetPrecisionFactor(); |
1444 |
|
TrkParams::SetStepMin(); |
1445 |
|
|
1446 |
|
TrkParams::Set(path,1); |
1447 |
|
TrkParams::Load(1); |
1448 |
|
|
1449 |
// |
// |
1450 |
}; |
}; |
1451 |
|
// /** |
1452 |
|
// * Get BY (kGauss) |
1453 |
|
// * @param v (x,y,z) coordinates in cm |
1454 |
|
// */ |
1455 |
|
// float TrkLevel2::GetBX(float* v){ |
1456 |
|
// float b[3]; |
1457 |
|
// gufld_(v,b); |
1458 |
|
// return b[0]/10.; |
1459 |
|
// } |
1460 |
|
// /** |
1461 |
|
// * Get BY (kGauss) |
1462 |
|
// * @param v (x,y,z) coordinates in cm |
1463 |
|
// */ |
1464 |
|
// float TrkLevel2::GetBY(float* v){ |
1465 |
|
// float b[3]; |
1466 |
|
// gufld_(v,b); |
1467 |
|
// return b[1]/10.; |
1468 |
|
// } |
1469 |
|
// /** |
1470 |
|
// * Get BY (kGauss) |
1471 |
|
// * @param v (x,y,z) coordinates in cm |
1472 |
|
// */ |
1473 |
|
// float TrkLevel2::GetBZ(float* v){ |
1474 |
|
// float b[3]; |
1475 |
|
// gufld_(v,b); |
1476 |
|
// return b[2]/10.; |
1477 |
|
// } |
1478 |
//-------------------------------------- |
//-------------------------------------- |
1479 |
// |
// |
1480 |
// |
// |
1571 |
tl = new float[npoint]; |
tl = new float[npoint]; |
1572 |
int i=0; |
int i=0; |
1573 |
do{ |
do{ |
1574 |
x[i] = 0; |
x[i] = 0; |
1575 |
y[i] = 0; |
y[i] = 0; |
1576 |
z[i] = zin[i]; |
z[i] = zin[i]; |
1577 |
thx[i] = 0; |
thx[i] = 0; |
1578 |
thy[i] = 0; |
thy[i] = 0; |
1579 |
tl[i] = 0; |
tl[i] = 0; |
1580 |
i++; |
i++; |
1581 |
}while(zin[i-1] > zin[i] && i < npoint); |
}while(zin[i-1] > zin[i] && i < npoint); |
1582 |
npoint=i; |
npoint=i; |
1583 |
if(npoint != n)cout << "NB! Trajectory created with "<<npoint<<" points"<<endl; |
if(npoint != n)cout << "NB! Trajectory created with "<<npoint<<" points"<<endl; |
1584 |
} |
} |
1585 |
|
void Trajectory::Delete(){ |
1586 |
|
|
1587 |
|
if(x) delete [] x; |
1588 |
|
if(y) delete [] y; |
1589 |
|
if(z) delete [] z; |
1590 |
|
if(thx) delete [] thx; |
1591 |
|
if(thy) delete [] thy; |
1592 |
|
if(tl) delete [] tl; |
1593 |
|
|
1594 |
|
} |
1595 |
//-------------------------------------- |
//-------------------------------------- |
1596 |
// |
// |
1597 |
// |
// |
1652 |
for (int i=0; i<5; i++) dal[i] = (double)al[i]; |
for (int i=0; i<5; i++) dal[i] = (double)al[i]; |
1653 |
for (int i=0; i<npoint; i++) dzin[i] = (double)z[i]; |
for (int i=0; i<npoint; i++) dzin[i] = (double)z[i]; |
1654 |
|
|
1655 |
|
TrkParams::Load(1); |
1656 |
|
if( !TrkParams::IsLoaded(1) ){ |
1657 |
|
cout << "int Trajectory::DoTrack2(float* al) --- ERROR --- m.field not loaded"<<endl; |
1658 |
|
return 0; |
1659 |
|
} |
1660 |
dotrack2_(&(npoint),dzin,dxout,dyout,dthxout,dthyout,dtlout,dal,&ifail); |
dotrack2_(&(npoint),dzin,dxout,dyout,dthxout,dthyout,dtlout,dal,&ifail); |
1661 |
|
|
1662 |
for (int i=0; i<npoint; i++){ |
for (int i=0; i<npoint; i++){ |