21 |
// |
// |
22 |
// |
// |
23 |
//-------------------------------------- |
//-------------------------------------- |
|
/** |
|
|
* Evaluates the trajectory in the apparatus associated to the track. |
|
|
* 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. |
|
|
* @param t pointer to an object of the class Trajectory, |
|
|
* which z coordinates should be previously initialized by calling the proper constructor ( Trajectory::Trajectory(int n, float* zin) ). |
|
|
* @return error flag. |
|
|
*/ |
|
|
int Trajectory::DoTrack2(float* al){ |
|
|
|
|
|
double *dxout = new double[t->npoint]; |
|
|
double *dyout = new double[t->npoint]; |
|
|
double *dthxout = new double[t->npoint]; |
|
|
double *dthyout = new double[t->npoint]; |
|
|
double *dtlout = new double[t->npoint]; |
|
|
double *dzin = new double[t->npoint]; |
|
|
double dal[5]; |
|
|
|
|
|
int ifail = 0; |
|
|
|
|
|
for (int i=0; i<5; i++) dal[i] = (double)al[i]; |
|
|
for (int i=0; i<t->npoint; i++) dzin[i] = (double)z[i]; |
|
|
|
|
|
dotrack2_(&(t->npoint),dzin,dxout,dyout,dthxout,dthyout,dtlout,dal,&ifail); |
|
|
|
|
|
for (int i=0; i<t->npoint; i++){ |
|
|
x[i] = (float)*dxout++; |
|
|
y[i] = (float)*dyout++; |
|
|
thx[i] = (float)*dthxout++; |
|
|
thy[i] = (float)*dthyout++; |
|
|
tl[i] = (float)*dtlout++; |
|
|
} |
|
|
|
|
|
return ifail; |
|
|
}; |
|
|
//--------------------------------------------- |
|
|
//--------------------------------------------- |
|
24 |
TrkTrack::TrkTrack(){ |
TrkTrack::TrkTrack(){ |
25 |
|
// cout << "TrkTrack::TrkTrack()" << endl; |
26 |
seqno = -1; |
seqno = -1; |
27 |
image = -1; |
image = -1; |
28 |
chi2 = 0; |
chi2 = 0; |
32 |
for(int it2=0;it2<5;it2++)coval[it1][it2] = 0; |
for(int it2=0;it2<5;it2++)coval[it1][it2] = 0; |
33 |
}; |
}; |
34 |
for(int ip=0;ip<6;ip++){ |
for(int ip=0;ip<6;ip++){ |
35 |
xgood[ip] = 0; |
xgood[ip] = 0; |
36 |
ygood[ip] = 0; |
ygood[ip] = 0; |
37 |
xm[ip] = 0; |
xm[ip] = 0; |
38 |
ym[ip] = 0; |
ym[ip] = 0; |
39 |
zm[ip] = 0; |
zm[ip] = 0; |
40 |
resx[ip] = 0; |
resx[ip] = 0; |
41 |
resy[ip] = 0; |
resy[ip] = 0; |
42 |
xv[ip] = 0; |
xv[ip] = 0; |
43 |
yv[ip] = 0; |
yv[ip] = 0; |
44 |
zv[ip] = 0; |
zv[ip] = 0; |
45 |
axv[ip] = 0; |
axv[ip] = 0; |
46 |
ayv[ip] = 0; |
ayv[ip] = 0; |
47 |
dedx_x[ip] = 0; |
dedx_x[ip] = 0; |
48 |
dedx_y[ip] = 0; |
dedx_y[ip] = 0; |
49 |
// clx[ip] = 0; |
}; |
50 |
// cly[ip] = 0; |
clx = 0; |
51 |
}; |
cly = 0; |
|
clx = new TRefArray(6,0); |
|
|
cly = new TRefArray(6,0); |
|
52 |
}; |
}; |
53 |
//-------------------------------------- |
//-------------------------------------- |
54 |
// |
// |
64 |
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]; |
65 |
}; |
}; |
66 |
for(int ip=0;ip<6;ip++){ |
for(int ip=0;ip<6;ip++){ |
67 |
xgood[ip] = t.xgood[ip]; |
xgood[ip] = t.xgood[ip]; |
68 |
ygood[ip] = t.ygood[ip]; |
ygood[ip] = t.ygood[ip]; |
69 |
xm[ip] = t.xm[ip]; |
xm[ip] = t.xm[ip]; |
70 |
ym[ip] = t.ym[ip]; |
ym[ip] = t.ym[ip]; |
71 |
zm[ip] = t.zm[ip]; |
zm[ip] = t.zm[ip]; |
72 |
resx[ip] = t.resx[ip]; |
resx[ip] = t.resx[ip]; |
73 |
resy[ip] = t.resy[ip]; |
resy[ip] = t.resy[ip]; |
74 |
xv[ip] = t.xv[ip]; |
xv[ip] = t.xv[ip]; |
75 |
yv[ip] = t.yv[ip]; |
yv[ip] = t.yv[ip]; |
76 |
zv[ip] = t.zv[ip]; |
zv[ip] = t.zv[ip]; |
77 |
axv[ip] = t.axv[ip]; |
axv[ip] = t.axv[ip]; |
78 |
ayv[ip] = t.ayv[ip]; |
ayv[ip] = t.ayv[ip]; |
79 |
dedx_x[ip] = t.dedx_x[ip]; |
dedx_x[ip] = t.dedx_x[ip]; |
80 |
dedx_y[ip] = t.dedx_y[ip]; |
dedx_y[ip] = t.dedx_y[ip]; |
81 |
// clx[ip] = 0;//<<<<pointer |
}; |
82 |
// cly[ip] = 0;//<<<<pointer |
clx = new TRefArray(*(t.clx)); |
83 |
}; |
cly = new TRefArray(*(t.cly)); |
|
clx = new TRefArray(*(t.clx)); |
|
|
cly = new TRefArray(*(t.cly)); |
|
84 |
|
|
85 |
}; |
}; |
86 |
//-------------------------------------- |
//-------------------------------------- |
87 |
// |
// |
88 |
// |
// |
89 |
//-------------------------------------- |
//-------------------------------------- |
90 |
|
void TrkTrack::Copy(TrkTrack& t){ |
91 |
|
|
92 |
|
t.seqno = seqno; |
93 |
|
t.image = image; |
94 |
|
t.chi2 = chi2; |
95 |
|
t.nstep = nstep; |
96 |
|
for(int it1=0;it1<5;it1++){ |
97 |
|
t.al[it1] = al[it1]; |
98 |
|
for(int it2=0;it2<5;it2++)t.coval[it1][it2] = coval[it1][it2]; |
99 |
|
}; |
100 |
|
for(int ip=0;ip<6;ip++){ |
101 |
|
t.xgood[ip] = xgood[ip]; |
102 |
|
t.ygood[ip] = ygood[ip]; |
103 |
|
t.xm[ip] = xm[ip]; |
104 |
|
t.ym[ip] = ym[ip]; |
105 |
|
t.zm[ip] = zm[ip]; |
106 |
|
t.resx[ip] = resx[ip]; |
107 |
|
t.resy[ip] = resy[ip]; |
108 |
|
t.xv[ip] = xv[ip]; |
109 |
|
t.yv[ip] = yv[ip]; |
110 |
|
t.zv[ip] = zv[ip]; |
111 |
|
t.axv[ip] = axv[ip]; |
112 |
|
t.ayv[ip] = ayv[ip]; |
113 |
|
t.dedx_x[ip] = dedx_x[ip]; |
114 |
|
t.dedx_y[ip] = dedx_y[ip]; |
115 |
|
|
116 |
|
}; |
117 |
|
|
118 |
|
}; |
119 |
|
//-------------------------------------- |
120 |
|
// |
121 |
|
// |
122 |
|
//-------------------------------------- |
123 |
/** |
/** |
124 |
* Evaluates the trajectory in the apparatus associated to the track. |
* Evaluates the trajectory in the apparatus associated to the track. |
125 |
* 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. |
* 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. |
288 |
readb_(); |
readb_(); |
289 |
|
|
290 |
}; |
}; |
291 |
|
|
292 |
|
/** |
293 |
|
* Method to fill minimization-routine common |
294 |
|
*/ |
295 |
|
void TrkTrack::FillMiniStruct(cMini2track& track){ |
296 |
|
|
297 |
|
for(int i=0; i<6; i++){ |
298 |
|
|
299 |
|
track.xgood[i]=xgood[i]; |
300 |
|
track.ygood[i]=ygood[i]; |
301 |
|
|
302 |
|
track.xm[i]=xm[i]; |
303 |
|
track.ym[i]=ym[i]; |
304 |
|
track.zm[i]=zm[i]; |
305 |
|
|
306 |
|
// --- temporaneo ---------------------------- |
307 |
|
// andrebbe inserita la dimensione del sensore |
308 |
|
float segment = 100.; |
309 |
|
track.xm_a[i]=xm[i]; |
310 |
|
track.xm_b[i]=xm[i]; |
311 |
|
track.ym_a[i]=ym[i]; |
312 |
|
track.ym_b[i]=ym[i]; |
313 |
|
if( xgood[i] && !ygood[i] ){ |
314 |
|
track.ym_a[i] = track.ym_a[i]+segment; |
315 |
|
track.ym_b[i] = track.ym_b[i]-segment; |
316 |
|
}else if( !xgood[i] && ygood[i]){ |
317 |
|
track.xm_a[i] = track.xm_a[i]+segment; |
318 |
|
track.xm_b[i] = track.xm_b[i]-segment; |
319 |
|
} |
320 |
|
// --- temporaneo ---------------------------- |
321 |
|
|
322 |
|
track.resx[i]=resx[i]; |
323 |
|
track.resy[i]=resy[i]; |
324 |
|
} |
325 |
|
|
326 |
|
for(int i=0; i<5; i++) track.al[i]=al[i]; |
327 |
|
track.zini = 23.5; |
328 |
|
// ZINI = 23.5 !!! it should be the same parameter in all codes |
329 |
|
|
330 |
|
} |
331 |
|
/** |
332 |
|
* Method to set values from minimization-routine common |
333 |
|
*/ |
334 |
|
void TrkTrack::SetFromMiniStruct(cMini2track *track){ |
335 |
|
|
336 |
|
for(int i=0; i<5; i++) { |
337 |
|
al[i]=track->al[i]; |
338 |
|
for(int j=0; j<5; j++) coval[i][j]=track->cov[i][j]; |
339 |
|
} |
340 |
|
chi2 = track->chi2; |
341 |
|
nstep = track->nstep; |
342 |
|
for(int i=0; i<6; i++){ |
343 |
|
xv[i] = track->xv[i]; |
344 |
|
yv[i] = track->yv[i]; |
345 |
|
zv[i] = track->zv[i]; |
346 |
|
xm[i] = track->xm[i]; |
347 |
|
ym[i] = track->ym[i]; |
348 |
|
zm[i] = track->zm[i]; |
349 |
|
axv[i] = track->axv[i]; |
350 |
|
ayv[i] = track->ayv[i]; |
351 |
|
} |
352 |
|
|
353 |
|
|
354 |
|
} |
355 |
/** |
/** |
356 |
* Tracking method. It calls F77 mini routine. |
* Tracking method. It calls F77 mini routine. |
357 |
*/ |
*/ |
361 |
|
|
362 |
extern cMini2track track_; |
extern cMini2track track_; |
363 |
fail = 0; |
fail = 0; |
364 |
|
FillMiniStruct(track_); |
365 |
|
|
366 |
for(int i=0; i<6; i++) track_.xm[i]=xm[i]; |
// if fit variables have been reset, evaluate the initial guess |
|
for(int i=0; i<6; i++) track_.ym[i]=ym[i]; |
|
|
for(int i=0; i<6; i++) track_.zm[i]=zm[i]; |
|
|
for(int i=0; i<6; i++) track_.resx[i]=resx[i]; |
|
|
for(int i=0; i<6; i++) track_.resy[i]=resy[i]; |
|
|
for(int i=0; i<6; i++) track_.xgood[i]=xgood[i]; |
|
|
for(int i=0; i<6; i++) track_.ygood[i]=ygood[i]; |
|
|
|
|
|
// initial guess of "al" with linear fit |
|
|
// if(al[0]==-9999.&&al[1]==-9999.&&al[2]==-9999.&&al[3]==-9999.&&al[4]==-9999.){ |
|
|
// cout << "initial guess "<<endl; |
|
|
// double szz=0., szx=0., szy=0., ssx=0., ssy=0., sz=0., s1=0.; |
|
|
// double det, ax, ay, bx, by; |
|
|
// for(int i=0; i<NPLANE; i++) { |
|
|
// szz=szz+zm[i]*zm[i]; |
|
|
// szx=szx+zm[i]*xm[i]; |
|
|
// szy=szy+zm[i]*ym[i]; |
|
|
// ssx=ssx+xm[i]; |
|
|
// ssy=ssy+ym[i]; |
|
|
// sz=sz+zm[i]; |
|
|
// s1=s1+1.; |
|
|
// } |
|
|
// det=szz*s1-sz*sz; |
|
|
// ax=(szx*s1-sz*ssx)/det; |
|
|
// bx=(szz*ssx-szx*sz)/det; |
|
|
// ay=(szy*s1-sz*ssy)/det; |
|
|
// by=(szz*ssy-szy*sz)/det; |
|
|
// al[0]=ax*23.5+bx; // ZINI = 23.5 !!! it should be the same parameter in all codes |
|
|
// al[1]=ay*23.5+by; // " |
|
|
// al[2]= sqrt(pow(ax,2)+pow(ay,2))/ sqrt(pow(ax,2)+pow(ay,2)+1.); |
|
|
// al[3]=0.; |
|
|
// if( (ax!=0.)||(ay!=0.) ) { |
|
|
// al[3]= asin(ay/ sqrt(pow(ax,2)+pow(ay,2))); |
|
|
// if(ax<0.) al[3]=acos(-1.)-al[3]; |
|
|
// } |
|
|
// al[4]=0.; |
|
|
// } |
|
|
// end guess |
|
|
|
|
|
for(int i=0; i<5; i++) track_.al[i]=al[i]; |
|
|
track_.zini = 23.5; |
|
|
// ZINI = 23.5 !!! it should be the same parameter in all codes |
|
|
|
|
367 |
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_(); |
368 |
|
|
369 |
// --------------------- free momentum |
// --------------------- free momentum |
370 |
if(pfixed==0.) { |
if(pfixed==0.) { |
371 |
// al[4]=0.; // free momentum |
track_.pfixed=0.; |
|
track_.pfixed=0.; // " |
|
372 |
} |
} |
373 |
// --------------------- fixed momentum |
// --------------------- fixed momentum |
374 |
if(pfixed!=0.) { |
if(pfixed!=0.) { |
375 |
al[4]=1./pfixed; // to fix the momentum |
al[4]=1./pfixed; |
376 |
track_.pfixed=pfixed; // " |
track_.pfixed=pfixed; |
377 |
} |
} |
378 |
|
|
379 |
// store temporarily the initial guess |
// store temporarily the initial guess |
380 |
for(int i=0; i<5; i++) al_ini[i]=track_.al[i]; |
for(int i=0; i<5; i++) al_ini[i]=track_.al[i]; |
381 |
|
|
382 |
|
// ------------------------------------------ |
383 |
|
// call mini routine |
384 |
int istep=0; |
int istep=0; |
385 |
int ifail=0; |
int ifail=0; |
386 |
mini2_(&istep,&ifail, &iprint); |
mini2_(&istep,&ifail, &iprint); |
387 |
if(ifail!=0) { |
if(ifail!=0) { |
388 |
if(iprint==1)cout << "ERROR: ifail= " << ifail << endl; |
if(iprint)cout << "ERROR: ifail= " << ifail << endl; |
389 |
fail = 1; |
fail = 1; |
|
// return; |
|
390 |
} |
} |
391 |
|
// ------------------------------------------ |
392 |
|
|
393 |
|
SetFromMiniStruct(&track_); |
394 |
// cout << endl << "eta ===> " << track_.al[4] << endl; |
// cout << endl << "eta ===> " << track_.al[4] << endl; |
395 |
|
|
396 |
for(int i=0; i<5; i++) al[i]=track_.al[i]; |
// for(int i=0; i<5; i++) al[i]=track_.al[i]; |
397 |
chi2=track_.chi2; |
// chi2=track_.chi2; |
398 |
nstep=track_.nstep; |
// nstep=track_.nstep; |
399 |
for(int i=0; i<6; i++) xv[i]=track_.xv[i]; |
// for(int i=0; i<6; i++) xv[i]=track_.xv[i]; |
400 |
for(int i=0; i<6; i++) yv[i]=track_.yv[i]; |
// for(int i=0; i<6; i++) yv[i]=track_.yv[i]; |
401 |
for(int i=0; i<6; i++) zv[i]=track_.zv[i]; |
// for(int i=0; i<6; i++) zv[i]=track_.zv[i]; |
402 |
for(int i=0; i<6; i++) axv[i]=track_.axv[i]; |
// for(int i=0; i<6; i++) axv[i]=track_.axv[i]; |
403 |
for(int i=0; i<6; i++) ayv[i]=track_.ayv[i]; |
// for(int i=0; i<6; i++) ayv[i]=track_.ayv[i]; |
404 |
for(int i=0; i<5; i++) { |
// for(int i=0; i<5; i++) { |
405 |
for(int j=0; j<5; j++) coval[i][j]=track_.cov[i][j]; |
// for(int j=0; j<5; j++) coval[i][j]=track_.cov[i][j]; |
406 |
} |
// } |
407 |
|
|
408 |
if(fail){ |
if(fail){ |
409 |
cout << " >>>> fit failed >>>> drawing initial par"<<endl; |
if(iprint)cout << " >>>> fit failed >>>> drawing initial par"<<endl; |
410 |
for(int i=0; i<5; i++) al[i]=al_ini[i]; |
for(int i=0; i<5; i++) al[i]=al_ini[i]; |
411 |
} |
} |
412 |
|
|
413 |
}; |
}; |
433 |
// |
// |
434 |
//-------------------------------------- |
//-------------------------------------- |
435 |
void TrkTrack::Clear(){ |
void TrkTrack::Clear(){ |
436 |
seqno = -1; |
// cout << "TrkTrack::Clear()"<<endl; |
437 |
image = -1; |
seqno = -1; |
438 |
chi2 = 0; |
image = -1; |
439 |
nstep = 0; |
chi2 = 0; |
440 |
for(int it1=0;it1<5;it1++){ |
nstep = 0; |
441 |
al[it1] = 0; |
for(int it1=0;it1<5;it1++){ |
442 |
for(int it2=0;it2<5;it2++)coval[it1][it2] = 0; |
al[it1] = 0; |
443 |
}; |
for(int it2=0;it2<5;it2++)coval[it1][it2] = 0; |
444 |
for(int ip=0;ip<6;ip++){ |
}; |
445 |
xgood[ip] = 0; |
for(int ip=0;ip<6;ip++){ |
446 |
ygood[ip] = 0; |
xgood[ip] = 0; |
447 |
xm[ip] = 0; |
ygood[ip] = 0; |
448 |
ym[ip] = 0; |
xm[ip] = 0; |
449 |
zm[ip] = 0; |
ym[ip] = 0; |
450 |
resx[ip] = 0; |
zm[ip] = 0; |
451 |
resy[ip] = 0; |
resx[ip] = 0; |
452 |
xv[ip] = 0; |
resy[ip] = 0; |
453 |
yv[ip] = 0; |
xv[ip] = 0; |
454 |
zv[ip] = 0; |
yv[ip] = 0; |
455 |
axv[ip] = 0; |
zv[ip] = 0; |
456 |
ayv[ip] = 0; |
axv[ip] = 0; |
457 |
dedx_x[ip] = 0; |
ayv[ip] = 0; |
458 |
dedx_y[ip] = 0; |
dedx_x[ip] = 0; |
459 |
// clx[ip] = 0; |
dedx_y[ip] = 0; |
460 |
// cly[ip] = 0; |
|
461 |
}; |
}; |
462 |
clx->Clear(); |
if(clx)clx->Clear(); |
463 |
cly->Clear(); |
if(cly)cly->Clear(); |
464 |
}; |
}; |
465 |
//-------------------------------------- |
//-------------------------------------- |
466 |
// |
// |
467 |
// |
// |
468 |
//-------------------------------------- |
//-------------------------------------- |
469 |
void TrkTrack::Delete(){ |
void TrkTrack::Delete(){ |
470 |
Clear(); |
// cout << "TrkTrack::Delete()"<<endl; |
471 |
clx->Delete(); |
// Clear(); |
472 |
cly->Delete(); |
if(clx)delete clx; |
473 |
|
if(cly)delete cly; |
474 |
}; |
}; |
475 |
//-------------------------------------- |
//-------------------------------------- |
476 |
// |
// |
477 |
// |
// |
478 |
//-------------------------------------- |
//-------------------------------------- |
482 |
// |
// |
483 |
//-------------------------------------- |
//-------------------------------------- |
484 |
TrkSinglet::TrkSinglet(){ |
TrkSinglet::TrkSinglet(){ |
485 |
|
// cout << "TrkSinglet::TrkSinglet() " << GetUniqueID()<<endl; |
486 |
plane = 0; |
plane = 0; |
487 |
coord[0] = 0; |
coord[0] = 0; |
488 |
coord[1] = 0; |
coord[1] = 0; |
489 |
sgnl = 0; |
sgnl = 0; |
490 |
cls = 0; |
// cls = 0; |
491 |
}; |
}; |
492 |
//-------------------------------------- |
//-------------------------------------- |
493 |
// |
// |
494 |
// |
// |
495 |
//-------------------------------------- |
//-------------------------------------- |
496 |
TrkSinglet::TrkSinglet(const TrkSinglet& s){ |
TrkSinglet::TrkSinglet(const TrkSinglet& s){ |
497 |
|
// cout << "TrkSinglet::TrkSinglet(const TrkSinglet& s) " << GetUniqueID()<<endl; |
498 |
plane = s.plane; |
plane = s.plane; |
499 |
coord[0] = s.coord[0]; |
coord[0] = s.coord[0]; |
500 |
coord[1] = s.coord[1]; |
coord[1] = s.coord[1]; |
501 |
sgnl = s.sgnl; |
sgnl = s.sgnl; |
502 |
// cls = 0;//<<<<pointer |
// cls = 0;//<<<<pointer |
503 |
cls = TRef(s.cls); |
cls = TRef(s.cls); |
504 |
}; |
}; |
505 |
//-------------------------------------- |
//-------------------------------------- |
506 |
// |
// |
517 |
// |
// |
518 |
// |
// |
519 |
//-------------------------------------- |
//-------------------------------------- |
520 |
|
void TrkSinglet::Clear(){ |
521 |
|
// cout << "TrkSinglet::Clear() " << GetUniqueID()<<endl; |
522 |
|
// cls=0; |
523 |
|
plane=-1; |
524 |
|
coord[0]=-999; |
525 |
|
coord[1]=-999; |
526 |
|
sgnl=0; |
527 |
|
|
528 |
|
} |
529 |
|
//-------------------------------------- |
530 |
|
// |
531 |
|
// |
532 |
|
//-------------------------------------- |
533 |
TrkLevel2::TrkLevel2(){ |
TrkLevel2::TrkLevel2(){ |
534 |
// good2 = -1; |
// cout <<"TrkLevel2::TrkLevel2()"<<endl; |
535 |
for(Int_t i=0; i<12 ; i++){ |
for(Int_t i=0; i<12 ; i++){ |
|
// crc[i] = -1; |
|
536 |
good[i] = -1; |
good[i] = -1; |
537 |
}; |
}; |
538 |
Track = new TClonesArray("TrkTrack"); |
// Track = new TClonesArray("TrkTrack"); |
539 |
SingletX = new TClonesArray("TrkSinglet"); |
// SingletX = new TClonesArray("TrkSinglet"); |
540 |
SingletY = new TClonesArray("TrkSinglet"); |
// SingletY = new TClonesArray("TrkSinglet"); |
541 |
|
Track = 0; |
542 |
|
SingletX = 0; |
543 |
|
SingletY = 0; |
544 |
|
|
|
// PhysicalTrack = new TClonesArray("TrkTrack"); |
|
|
//sostituire con TRefArray... appena ho capito come si usa |
|
545 |
} |
} |
546 |
//-------------------------------------- |
//-------------------------------------- |
547 |
// |
// |
549 |
//-------------------------------------- |
//-------------------------------------- |
550 |
void TrkLevel2::Dump(){ |
void TrkLevel2::Dump(){ |
551 |
|
|
|
TClonesArray &t = *Track; |
|
|
TClonesArray &sx = *SingletX; |
|
|
TClonesArray &sy = *SingletY; |
|
552 |
// |
// |
553 |
cout << endl << endl << "=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-"; |
cout << endl << endl << "=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-"; |
554 |
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]<<" "; |
555 |
cout << endl << "ntrk() : " << this->ntrk() ; |
cout << endl << "ntrk() : " << this->ntrk() ; |
556 |
cout << endl << "nclsx() : " << this->nclsx(); |
cout << endl << "nclsx() : " << this->nclsx(); |
557 |
cout << endl << "nclsy() : " << this->nclsy(); |
cout << endl << "nclsy() : " << this->nclsy(); |
558 |
for(int i=0; i<this->ntrk(); i++) ((TrkTrack *)t[i])->Dump(); |
// TClonesArray &t = *Track; |
559 |
for(int i=0; i<this->nclsx(); i++) ((TrkSinglet *)sx[i])->Dump(); |
// TClonesArray &sx = *SingletX; |
560 |
for(int i=0; i<this->nclsy(); i++) ((TrkSinglet *)sy[i])->Dump(); |
// TClonesArray &sy = *SingletY; |
561 |
|
// for(int i=0; i<ntrk(); i++) ((TrkTrack *)t[i])->Dump(); |
562 |
|
// for(int i=0; i<nclsx(); i++) ((TrkSinglet *)sx[i])->Dump(); |
563 |
|
// for(int i=0; i<nclsy(); i++) ((TrkSinglet *)sy[i])->Dump(); |
564 |
|
if(Track){ |
565 |
|
TClonesArray &t = *Track; |
566 |
|
for(int i=0; i<ntrk(); i++) ((TrkTrack *)t[i])->Dump(); |
567 |
|
} |
568 |
|
if(SingletX){ |
569 |
|
TClonesArray &sx = *SingletX; |
570 |
|
for(int i=0; i<nclsx(); i++) ((TrkSinglet *)sx[i])->Dump(); |
571 |
|
} |
572 |
|
if(SingletY){ |
573 |
|
TClonesArray &sy = *SingletY; |
574 |
|
for(int i=0; i<nclsy(); i++) ((TrkSinglet *)sy[i])->Dump(); |
575 |
|
} |
576 |
} |
} |
577 |
//-------------------------------------- |
//-------------------------------------- |
578 |
// |
// |
594 |
good[i] = l2->good[i]; |
good[i] = l2->good[i]; |
595 |
}; |
}; |
596 |
// *** TRACKS *** |
// *** TRACKS *** |
597 |
|
if(!Track) Track = new TClonesArray("TrkTrack"); |
598 |
TClonesArray &t = *Track; |
TClonesArray &t = *Track; |
599 |
for(int i=0; i<l2->ntrk; i++){ |
for(int i=0; i<l2->ntrk; i++){ |
600 |
t_track->seqno = i;// NBNBNBNB deve sempre essere = i |
t_track->seqno = i;// NBNBNBNB deve sempre essere = i |
601 |
t_track->image = l2->image[i]-1; |
t_track->image = l2->image[i]-1; |
602 |
// cout << "track "<<i<<t_track->seqno << t_track->image<<endl; |
// cout << "track "<<i<<t_track->seqno << t_track->image<<endl; |
603 |
t_track->chi2 = l2->chi2_nt[i]; |
t_track->chi2 = l2->chi2_nt[i]; |
604 |
t_track->nstep = l2->nstep_nt[i]; |
t_track->nstep = l2->nstep_nt[i]; |
605 |
for(int it1=0;it1<5;it1++){ |
for(int it1=0;it1<5;it1++){ |
606 |
t_track->al[it1] = l2->al_nt[i][it1]; |
t_track->al[it1] = l2->al_nt[i][it1]; |
607 |
for(int it2=0;it2<5;it2++) |
for(int it2=0;it2<5;it2++) |
608 |
t_track->coval[it1][it2] = l2->coval[i][it2][it1]; |
t_track->coval[it1][it2] = l2->coval[i][it2][it1]; |
609 |
}; |
}; |
610 |
for(int ip=0;ip<6;ip++){ |
for(int ip=0;ip<6;ip++){ |
611 |
t_track->xgood[ip] = l2->xgood_nt[i][ip]; |
t_track->xgood[ip] = l2->xgood_nt[i][ip]; |
612 |
t_track->ygood[ip] = l2->ygood_nt[i][ip]; |
t_track->ygood[ip] = l2->ygood_nt[i][ip]; |
613 |
t_track->xm[ip] = l2->xm_nt[i][ip]; |
t_track->xm[ip] = l2->xm_nt[i][ip]; |
614 |
t_track->ym[ip] = l2->ym_nt[i][ip]; |
t_track->ym[ip] = l2->ym_nt[i][ip]; |
615 |
t_track->zm[ip] = l2->zm_nt[i][ip]; |
t_track->zm[ip] = l2->zm_nt[i][ip]; |
616 |
t_track->resx[ip] = l2->resx_nt[i][ip]; |
t_track->resx[ip] = l2->resx_nt[i][ip]; |
617 |
t_track->resy[ip] = l2->resy_nt[i][ip]; |
t_track->resy[ip] = l2->resy_nt[i][ip]; |
618 |
t_track->xv[ip] = l2->xv_nt[i][ip]; |
t_track->xv[ip] = l2->xv_nt[i][ip]; |
619 |
t_track->yv[ip] = l2->yv_nt[i][ip]; |
t_track->yv[ip] = l2->yv_nt[i][ip]; |
620 |
t_track->zv[ip] = l2->zv_nt[i][ip]; |
t_track->zv[ip] = l2->zv_nt[i][ip]; |
621 |
t_track->axv[ip] = l2->axv_nt[i][ip]; |
t_track->axv[ip] = l2->axv_nt[i][ip]; |
622 |
t_track->ayv[ip] = l2->ayv_nt[i][ip]; |
t_track->ayv[ip] = l2->ayv_nt[i][ip]; |
623 |
t_track->dedx_x[ip] = l2->dedx_x[i][ip]; |
t_track->dedx_x[ip] = l2->dedx_x[i][ip]; |
624 |
t_track->dedx_y[ip] = l2->dedx_y[i][ip]; |
t_track->dedx_y[ip] = l2->dedx_y[i][ip]; |
625 |
// t_track->clx[ip] = 0; |
// t_track->clx[ip] = 0; |
626 |
// t_track->cly[ip] = 0; |
// t_track->cly[ip] = 0; |
627 |
}; |
}; |
628 |
new(t[i]) TrkTrack(*t_track); |
new(t[i]) TrkTrack(*t_track); |
629 |
t_track->Clear(); |
t_track->Clear(); |
630 |
}; |
}; |
631 |
// *** SINGLETS *** |
// *** SINGLETS *** |
632 |
|
if(!SingletX)SingletX = new TClonesArray("TrkSinglet"); |
633 |
TClonesArray &sx = *SingletX; |
TClonesArray &sx = *SingletX; |
634 |
for(int i=0; i<l2->nclsx; i++){ |
for(int i=0; i<l2->nclsx; i++){ |
635 |
t_singlet->plane = l2->planex[i]; |
t_singlet->plane = l2->planex[i]; |
636 |
t_singlet->coord[0] = l2->xs[i][0]; |
t_singlet->coord[0] = l2->xs[i][0]; |
637 |
t_singlet->coord[1] = l2->xs[i][1]; |
t_singlet->coord[1] = l2->xs[i][1]; |
638 |
t_singlet->sgnl = l2->signlxs[i]; |
t_singlet->sgnl = l2->signlxs[i]; |
639 |
// t_singlet->cls = 0; |
// t_singlet->cls = 0; |
640 |
new(sx[i]) TrkSinglet(*t_singlet); |
new(sx[i]) TrkSinglet(*t_singlet); |
641 |
t_singlet->Clear(); |
t_singlet->Clear(); |
642 |
} |
} |
643 |
|
if(!SingletY)SingletY = new TClonesArray("TrkSinglet"); |
644 |
TClonesArray &sy = *SingletY; |
TClonesArray &sy = *SingletY; |
645 |
for(int i=0; i<l2->nclsy; i++){ |
for(int i=0; i<l2->nclsy; i++){ |
646 |
t_singlet->plane = l2->planey[i]; |
t_singlet->plane = l2->planey[i]; |
647 |
t_singlet->coord[0] = l2->ys[i][0]; |
t_singlet->coord[0] = l2->ys[i][0]; |
648 |
t_singlet->coord[1] = l2->ys[i][1]; |
t_singlet->coord[1] = l2->ys[i][1]; |
649 |
t_singlet->sgnl = l2->signlys[i]; |
t_singlet->sgnl = l2->signlys[i]; |
650 |
// t_singlet->cls = 0; |
// t_singlet->cls = 0; |
651 |
new(sy[i]) TrkSinglet(*t_singlet); |
new(sy[i]) TrkSinglet(*t_singlet); |
652 |
t_singlet->Clear(); |
t_singlet->Clear(); |
653 |
}; |
}; |
654 |
|
|
655 |
delete t_track; |
delete t_track; |
656 |
delete t_singlet; |
delete t_singlet; |
657 |
} |
} |
658 |
//-------------------------------------- |
//-------------------------------------- |
659 |
// |
// |
675 |
good[i] = l2->good[i]; |
good[i] = l2->good[i]; |
676 |
}; |
}; |
677 |
// *** TRACKS *** |
// *** TRACKS *** |
678 |
|
if(!Track) Track = new TClonesArray("TrkTrack"); |
679 |
TClonesArray &t = *Track; |
TClonesArray &t = *Track; |
680 |
|
if(!t_track->clx)t_track->clx = new TRefArray(6,0); |
681 |
|
if(!t_track->cly)t_track->cly = new TRefArray(6,0); |
682 |
for(int i=0; i<l2->ntrk; i++){ |
for(int i=0; i<l2->ntrk; i++){ |
683 |
t_track->seqno = i;// NBNBNBNB deve sempre essere = i |
t_track->seqno = i;// NBNBNBNB deve sempre essere = i |
684 |
t_track->image = l2->image[i]-1; |
t_track->image = l2->image[i]-1; |
707 |
t_track->dedx_y[ip] = l2->dedx_y[i][ip]; |
t_track->dedx_y[ip] = l2->dedx_y[i][ip]; |
708 |
// cout << "traccia "<<i<<" -- "<< ip << " "<< l2->cltrx[i][ip] <<" "<< l2->cltry[i][ip] <<" "<< t_track->xgood[ip] << t_track->ygood[ip]<<endl; |
// cout << "traccia "<<i<<" -- "<< ip << " "<< l2->cltrx[i][ip] <<" "<< l2->cltry[i][ip] <<" "<< t_track->xgood[ip] << t_track->ygood[ip]<<endl; |
709 |
//----------------------------------------------------- |
//----------------------------------------------------- |
710 |
// t_track->clx[ip] = l1->GetCluster(l2->cltrx[i][ip]-1); |
//----------------------------------------------------- |
|
// t_track->cly[ip] = l1->GetCluster(l2->cltry[i][ip]-1); |
|
711 |
if(t_track->xgood[ip])t_track->clx->AddAt(l1->GetCluster(l2->cltrx[i][ip]-1),ip); |
if(t_track->xgood[ip])t_track->clx->AddAt(l1->GetCluster(l2->cltrx[i][ip]-1),ip); |
712 |
if(t_track->ygood[ip])t_track->cly->AddAt(l1->GetCluster(l2->cltry[i][ip]-1),ip); |
if(t_track->ygood[ip])t_track->cly->AddAt(l1->GetCluster(l2->cltry[i][ip]-1),ip); |
713 |
// 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->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; |
714 |
// 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; |
// 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; |
715 |
//----------------------------------------------------- |
//----------------------------------------------------- |
716 |
|
//----------------------------------------------------- |
717 |
}; |
}; |
718 |
new(t[i]) TrkTrack(*t_track); |
new(t[i]) TrkTrack(*t_track); |
719 |
t_track->Clear(); |
t_track->Clear(); |
720 |
}; |
}; |
721 |
// *** SINGLETS *** |
// *** SINGLETS *** |
722 |
|
if(!SingletX)SingletX = new TClonesArray("TrkSinglet"); |
723 |
TClonesArray &sx = *SingletX; |
TClonesArray &sx = *SingletX; |
724 |
for(int i=0; i<l2->nclsx; i++){ |
for(int i=0; i<l2->nclsx; i++){ |
725 |
t_singlet->plane = l2->planex[i]; |
t_singlet->plane = l2->planex[i]; |
734 |
new(sx[i]) TrkSinglet(*t_singlet); |
new(sx[i]) TrkSinglet(*t_singlet); |
735 |
t_singlet->Clear(); |
t_singlet->Clear(); |
736 |
} |
} |
737 |
|
if(!SingletY)SingletY = new TClonesArray("TrkSinglet"); |
738 |
TClonesArray &sy = *SingletY; |
TClonesArray &sy = *SingletY; |
739 |
for(int i=0; i<l2->nclsy; i++){ |
for(int i=0; i<l2->nclsy; i++){ |
740 |
t_singlet->plane = l2->planey[i]; |
t_singlet->plane = l2->planey[i]; |
767 |
}; |
}; |
768 |
// *** TRACKS *** |
// *** TRACKS *** |
769 |
|
|
770 |
l2->ntrk = Track->GetEntries(); |
if(Track){ |
771 |
for(Int_t i=0;i<l2->ntrk;i++){ |
l2->ntrk = Track->GetEntries(); |
772 |
l2->image[i] = 1 + ((TrkTrack *)Track->At(i))->image; |
for(Int_t i=0;i<l2->ntrk;i++){ |
773 |
l2->chi2_nt[i] = ((TrkTrack *)Track->At(i))->chi2; |
l2->image[i] = 1 + ((TrkTrack *)Track->At(i))->image; |
774 |
l2->nstep_nt[i] = ((TrkTrack *)Track->At(i))->nstep; |
l2->chi2_nt[i] = ((TrkTrack *)Track->At(i))->chi2; |
775 |
for(int it1=0;it1<5;it1++){ |
l2->nstep_nt[i] = ((TrkTrack *)Track->At(i))->nstep; |
776 |
l2->al_nt[i][it1] = ((TrkTrack *)Track->At(i))->al[it1]; |
for(int it1=0;it1<5;it1++){ |
777 |
for(int it2=0;it2<5;it2++) |
l2->al_nt[i][it1] = ((TrkTrack *)Track->At(i))->al[it1]; |
778 |
l2->coval[i][it2][it1] = ((TrkTrack *)Track->At(i))->coval[it1][it2]; |
for(int it2=0;it2<5;it2++) |
779 |
}; |
l2->coval[i][it2][it1] = ((TrkTrack *)Track->At(i))->coval[it1][it2]; |
780 |
for(int ip=0;ip<6;ip++){ |
}; |
781 |
l2->xgood_nt[i][ip] = ((TrkTrack *)Track->At(i))->xgood[ip]; |
for(int ip=0;ip<6;ip++){ |
782 |
l2->ygood_nt[i][ip] = ((TrkTrack *)Track->At(i))->ygood[ip]; |
l2->xgood_nt[i][ip] = ((TrkTrack *)Track->At(i))->xgood[ip]; |
783 |
l2->xm_nt[i][ip] = ((TrkTrack *)Track->At(i))->xm[ip]; |
l2->ygood_nt[i][ip] = ((TrkTrack *)Track->At(i))->ygood[ip]; |
784 |
l2->ym_nt[i][ip] = ((TrkTrack *)Track->At(i))->ym[ip]; |
l2->xm_nt[i][ip] = ((TrkTrack *)Track->At(i))->xm[ip]; |
785 |
l2->zm_nt[i][ip] = ((TrkTrack *)Track->At(i))->zm[ip]; |
l2->ym_nt[i][ip] = ((TrkTrack *)Track->At(i))->ym[ip]; |
786 |
l2->resx_nt[i][ip] = ((TrkTrack *)Track->At(i))->resx[ip]; |
l2->zm_nt[i][ip] = ((TrkTrack *)Track->At(i))->zm[ip]; |
787 |
l2->resy_nt[i][ip] = ((TrkTrack *)Track->At(i))->resy[ip]; |
l2->resx_nt[i][ip] = ((TrkTrack *)Track->At(i))->resx[ip]; |
788 |
l2->xv_nt[i][ip] = ((TrkTrack *)Track->At(i))->xv[ip]; |
l2->resy_nt[i][ip] = ((TrkTrack *)Track->At(i))->resy[ip]; |
789 |
l2->yv_nt[i][ip] = ((TrkTrack *)Track->At(i))->yv[ip]; |
l2->xv_nt[i][ip] = ((TrkTrack *)Track->At(i))->xv[ip]; |
790 |
l2->zv_nt[i][ip] = ((TrkTrack *)Track->At(i))->zv[ip]; |
l2->yv_nt[i][ip] = ((TrkTrack *)Track->At(i))->yv[ip]; |
791 |
l2->axv_nt[i][ip] = ((TrkTrack *)Track->At(i))->axv[ip]; |
l2->zv_nt[i][ip] = ((TrkTrack *)Track->At(i))->zv[ip]; |
792 |
l2->ayv_nt[i][ip] = ((TrkTrack *)Track->At(i))->ayv[ip]; |
l2->axv_nt[i][ip] = ((TrkTrack *)Track->At(i))->axv[ip]; |
793 |
l2->dedx_x[i][ip] = ((TrkTrack *)Track->At(i))->dedx_x[ip]; |
l2->ayv_nt[i][ip] = ((TrkTrack *)Track->At(i))->ayv[ip]; |
794 |
l2->dedx_y[i][ip] = ((TrkTrack *)Track->At(i))->dedx_y[ip]; |
l2->dedx_x[i][ip] = ((TrkTrack *)Track->At(i))->dedx_x[ip]; |
795 |
}; |
l2->dedx_y[i][ip] = ((TrkTrack *)Track->At(i))->dedx_y[ip]; |
796 |
|
}; |
797 |
|
} |
798 |
} |
} |
|
|
|
799 |
// *** SINGLETS *** |
// *** SINGLETS *** |
800 |
l2->nclsx = SingletX->GetEntries(); |
if(SingletX){ |
801 |
for(Int_t i=0;i<l2->nclsx;i++){ |
l2->nclsx = SingletX->GetEntries(); |
802 |
l2->planex[i] = ((TrkSinglet *)SingletX->At(i))->plane; |
for(Int_t i=0;i<l2->nclsx;i++){ |
803 |
l2->xs[i][0] = ((TrkSinglet *)SingletX->At(i))->coord[0]; |
l2->planex[i] = ((TrkSinglet *)SingletX->At(i))->plane; |
804 |
l2->xs[i][1] = ((TrkSinglet *)SingletX->At(i))->coord[1]; |
l2->xs[i][0] = ((TrkSinglet *)SingletX->At(i))->coord[0]; |
805 |
l2->signlxs[i] = ((TrkSinglet *)SingletX->At(i))->sgnl; |
l2->xs[i][1] = ((TrkSinglet *)SingletX->At(i))->coord[1]; |
806 |
} |
l2->signlxs[i] = ((TrkSinglet *)SingletX->At(i))->sgnl; |
807 |
l2->nclsy = SingletY->GetEntries(); |
} |
808 |
for(Int_t i=0;i<l2->nclsy;i++){ |
} |
809 |
l2->planey[i] = ((TrkSinglet *)SingletY->At(i))->plane; |
|
810 |
l2->ys[i][0] = ((TrkSinglet *)SingletY->At(i))->coord[0]; |
if(SingletY){ |
811 |
l2->ys[i][1] = ((TrkSinglet *)SingletY->At(i))->coord[1]; |
l2->nclsy = SingletY->GetEntries(); |
812 |
l2->signlys[i] = ((TrkSinglet *)SingletY->At(i))->sgnl; |
for(Int_t i=0;i<l2->nclsy;i++){ |
813 |
|
l2->planey[i] = ((TrkSinglet *)SingletY->At(i))->plane; |
814 |
|
l2->ys[i][0] = ((TrkSinglet *)SingletY->At(i))->coord[0]; |
815 |
|
l2->ys[i][1] = ((TrkSinglet *)SingletY->At(i))->coord[1]; |
816 |
|
l2->signlys[i] = ((TrkSinglet *)SingletY->At(i))->sgnl; |
817 |
|
} |
818 |
} |
} |
819 |
} |
} |
820 |
//-------------------------------------- |
//-------------------------------------- |
822 |
// |
// |
823 |
//-------------------------------------- |
//-------------------------------------- |
824 |
void TrkLevel2::Clear(){ |
void TrkLevel2::Clear(){ |
|
// good2 = -1; |
|
825 |
for(Int_t i=0; i<12 ; i++){ |
for(Int_t i=0; i<12 ; i++){ |
826 |
// crc[i] = -1; |
good[i] = -1; |
827 |
good[i] = -1; |
}; |
828 |
}; |
// if(Track)Track->Clear("C"); |
829 |
/* Track->RemoveAll(); |
// if(SingletX)SingletX->Clear("C"); |
830 |
SingletX->RemoveAll(); |
// if(SingletY)SingletY->Clear("C"); |
831 |
SingletY->RemoveAll();*/ |
if(Track)Track->Delete(); |
832 |
// modify to avoid memory leakage |
if(SingletX)SingletX->Delete(); |
833 |
Track->Clear(); |
if(SingletY)SingletY->Delete(); |
834 |
SingletX->Clear(); |
} |
835 |
SingletY->Clear(); |
// //-------------------------------------- |
836 |
} |
// // |
837 |
//-------------------------------------- |
// // |
838 |
// |
// //-------------------------------------- |
|
// |
|
|
//-------------------------------------- |
|
839 |
void TrkLevel2::Delete(){ |
void TrkLevel2::Delete(){ |
840 |
|
|
841 |
Clear(); |
// cout << "void TrkLevel2::Delete()"<<endl; |
842 |
Track->Delete(); |
Clear(); |
843 |
SingletX->Delete(); |
if(Track)delete Track; |
844 |
SingletY->Delete(); |
if(SingletX)delete SingletX; |
845 |
|
if(SingletY)delete SingletY; |
846 |
|
|
847 |
} |
} |
848 |
//-------------------------------------- |
//-------------------------------------- |
849 |
// |
// |
855 |
*/ |
*/ |
856 |
TRefArray *TrkLevel2::GetTracks_NFitSorted(){ |
TRefArray *TrkLevel2::GetTracks_NFitSorted(){ |
857 |
|
|
858 |
TRefArray *sorted = new TRefArray(); |
if(!Track)return 0; |
859 |
|
|
860 |
|
TRefArray *sorted = new TRefArray(); |
861 |
|
|
862 |
TClonesArray &t = *Track; |
TClonesArray &t = *Track; |
863 |
// TClonesArray &ts = *PhysicalTrack; |
// TClonesArray &ts = *PhysicalTrack; |
864 |
int N = ntrk(); |
int N = ntrk(); |
865 |
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; |
866 |
// 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; |
867 |
|
|
868 |
int indo=0; |
int indo=0; |
869 |
int indi=0; |
int indi=0; |
870 |
while(N != 0){ |
while(N != 0){ |
871 |
int nfit =0; |
int nfit =0; |
872 |
float chi2ref = numeric_limits<float>::max(); |
float chi2ref = numeric_limits<float>::max(); |
873 |
|
|
874 |
// first loop to search maximum num. of fit points |
// first loop to search maximum num. of fit points |
875 |
for(int i=0; i < ntrk(); i++){ |
for(int i=0; i < ntrk(); i++){ |
876 |
if( ((TrkTrack *)t[i])->GetNtot() >= nfit && m[i]==1){ |
if( ((TrkTrack *)t[i])->GetNtot() >= nfit && m[i]==1){ |
877 |
nfit = ((TrkTrack *)t[i])->GetNtot(); |
nfit = ((TrkTrack *)t[i])->GetNtot(); |
878 |
} |
} |
879 |
} |
} |
880 |
//second loop to search minimum chi2 among selected |
//second loop to search minimum chi2 among selected |
881 |
for(int i=0; i<this->ntrk(); i++){ |
for(int i=0; i<this->ntrk(); i++){ |
882 |
Float_t chi2 = ((TrkTrack *)t[i])->chi2; |
Float_t chi2 = ((TrkTrack *)t[i])->chi2; |
883 |
if(chi2 < 0) chi2 = chi2*1000; |
if(chi2 < 0) chi2 = chi2*1000; |
884 |
if( chi2 < chi2ref |
if( chi2 < chi2ref |
885 |
&& ((TrkTrack *)t[i])->GetNtot() == nfit |
&& ((TrkTrack *)t[i])->GetNtot() == nfit |
886 |
&& m[i]==1){ |
&& m[i]==1){ |
887 |
chi2ref = ((TrkTrack *)t[i])->chi2; |
chi2ref = ((TrkTrack *)t[i])->chi2; |
888 |
indi = i; |
indi = i; |
889 |
}; |
}; |
890 |
}; |
}; |
891 |
if( ((TrkTrack *)t[indi])->HasImage() ){ |
if( ((TrkTrack *)t[indi])->HasImage() ){ |
892 |
m[((TrkTrack *)t[indi])->image] = 0; |
m[((TrkTrack *)t[indi])->image] = 0; |
893 |
N--; |
N--; |
894 |
|
|
895 |
// cout << "i** "<< ((TrkTrack *)t[indi])->image << " " << nfiti <<" "<<chi2i<<endl; |
// cout << "i** "<< ((TrkTrack *)t[indi])->image << " " << nfiti <<" "<<chi2i<<endl; |
896 |
}; |
}; |
897 |
sorted->Add( (TrkTrack*)t[indi] ); |
sorted->Add( (TrkTrack*)t[indi] ); |
898 |
|
|
899 |
m[indi] = 0; |
m[indi] = 0; |
900 |
// cout << "SORTED "<< indo << " "<< indi << " "<< N << endl; |
// cout << "SORTED "<< indo << " "<< indi << " "<< N << endl; |
901 |
N--; |
N--; |
902 |
indo++; |
indo++; |
903 |
} |
} |
904 |
m.clear(); |
m.clear(); |
905 |
// cout << "GetTracks_NFitSorted(it): Done"<< endl; |
// cout << "GetTracks_NFitSorted(it): Done"<< endl; |
906 |
|
|
907 |
return sorted; |
return sorted; |
908 |
// return PhysicalTrack; |
// return PhysicalTrack; |
909 |
} |
} |
910 |
//-------------------------------------- |
//-------------------------------------- |
923 |
cout << " Stored tracks ntrk() = "<< this->ntrk() << endl; |
cout << " Stored tracks ntrk() = "<< this->ntrk() << endl; |
924 |
return 0; |
return 0; |
925 |
} |
} |
926 |
|
if(!Track){ |
927 |
|
cout << "TrkTrack *TrkLevel2::GetStoredTrack(int is) >> (TClonesArray*) Track ==0 "<<endl; |
928 |
|
}; |
929 |
TClonesArray &t = *(Track); |
TClonesArray &t = *(Track); |
930 |
TrkTrack *track = (TrkTrack*)t[is]; |
TrkTrack *track = (TrkTrack*)t[is]; |
931 |
return track; |
return track; |
945 |
cout << " Stored x-singlets nclsx() = "<< this->nclsx() << endl; |
cout << " Stored x-singlets nclsx() = "<< this->nclsx() << endl; |
946 |
return 0; |
return 0; |
947 |
} |
} |
948 |
|
if(!SingletX)return 0; |
949 |
TClonesArray &t = *(SingletX); |
TClonesArray &t = *(SingletX); |
950 |
TrkSinglet *singlet = (TrkSinglet*)t[is]; |
TrkSinglet *singlet = (TrkSinglet*)t[is]; |
951 |
return singlet; |
return singlet; |
965 |
cout << " Stored y-singlets nclsy() = "<< this->nclsx() << endl; |
cout << " Stored y-singlets nclsy() = "<< this->nclsx() << endl; |
966 |
return 0; |
return 0; |
967 |
} |
} |
968 |
|
if(!SingletY)return 0; |
969 |
TClonesArray &t = *(SingletY); |
TClonesArray &t = *(SingletY); |
970 |
TrkSinglet *singlet = (TrkSinglet*)t[is]; |
TrkSinglet *singlet = (TrkSinglet*)t[is]; |
971 |
return singlet; |
return singlet; |
988 |
} |
} |
989 |
|
|
990 |
TRefArray *sorted = GetTracks(); //TEMPORANEO |
TRefArray *sorted = GetTracks(); //TEMPORANEO |
991 |
|
if(!sorted)return 0; |
992 |
TrkTrack *track = (TrkTrack*)sorted->At(it); |
TrkTrack *track = (TrkTrack*)sorted->At(it); |
993 |
sorted->Delete(); |
sorted->Clear(); |
994 |
|
delete sorted; |
995 |
return track; |
return track; |
996 |
} |
} |
997 |
/** |
/** |
1000 |
Int_t TrkLevel2::GetNTracks(){ |
Int_t TrkLevel2::GetNTracks(){ |
1001 |
|
|
1002 |
Float_t ntot=0; |
Float_t ntot=0; |
1003 |
|
if(!Track)return 0; |
1004 |
TClonesArray &t = *Track; |
TClonesArray &t = *Track; |
1005 |
for(int i=0; i<ntrk(); i++) { |
for(int i=0; i<ntrk(); i++) { |
1006 |
if( ((TrkTrack *)t[i])->GetImageSeqNo() == -1 ) ntot+=1.; |
if( ((TrkTrack *)t[i])->GetImageSeqNo() == -1 ) ntot+=1.; |
1019 |
*/ |
*/ |
1020 |
TrkTrack *TrkLevel2::GetTrackImage(int it){ |
TrkTrack *TrkLevel2::GetTrackImage(int it){ |
1021 |
|
|
1022 |
if(it >= this->GetNTracks()){ |
if(it >= this->GetNTracks()){ |
1023 |
cout << "** TrkLevel2 ** Track "<< it << "does not exits! " << endl; |
cout << "** TrkLevel2 ** Track "<< it << "does not exits! " << endl; |
1024 |
cout << " Physical tracks GetNTracks() = "<< this->ntrk() << endl; |
cout << " Physical tracks GetNTracks() = "<< this->ntrk() << endl; |
1025 |
return 0; |
return 0; |
1026 |
} |
} |
1027 |
|
|
1028 |
TRefArray* sorted = GetTracks(); //TEMPORANEO |
TRefArray* sorted = GetTracks(); //TEMPORANEO |
1029 |
TrkTrack *track = (TrkTrack*)sorted->At(it); |
if(!sorted)return 0; |
1030 |
|
TrkTrack *track = (TrkTrack*)sorted->At(it); |
1031 |
|
|
1032 |
if(!track->HasImage()){ |
if(!track->HasImage()){ |
1033 |
cout << "** TrkLevel2 ** Track "<< it << "does not have image! " << endl; |
cout << "** TrkLevel2 ** Track "<< it << "does not have image! " << endl; |
1034 |
return 0; |
return 0; |
1035 |
} |
} |
1036 |
TrkTrack *image = (TrkTrack*)(*Track)[track->image]; |
if(!Track)return 0; |
1037 |
|
TrkTrack *image = (TrkTrack*)(*Track)[track->image]; |
1038 |
|
|
1039 |
sorted->Delete(); |
sorted->Delete(); |
1040 |
|
delete sorted; |
1041 |
return image; |
|
1042 |
|
return image; |
1043 |
|
|
1044 |
} |
} |
1045 |
//-------------------------------------- |
//-------------------------------------- |
1154 |
tl = new float[npoint]; |
tl = new float[npoint]; |
1155 |
int i=0; |
int i=0; |
1156 |
do{ |
do{ |
1157 |
x[i] = 0; |
x[i] = 0; |
1158 |
y[i] = 0; |
y[i] = 0; |
1159 |
z[i] = zin[i]; |
z[i] = zin[i]; |
1160 |
thx[i] = 0; |
thx[i] = 0; |
1161 |
thy[i] = 0; |
thy[i] = 0; |
1162 |
tl[i] = 0; |
tl[i] = 0; |
1163 |
i++; |
i++; |
1164 |
}while(zin[i-1] > zin[i] && i < npoint); |
}while(zin[i-1] > zin[i] && i < npoint); |
1165 |
npoint=i; |
npoint=i; |
1166 |
if(npoint != n)cout << "NB! Trajectory created with "<<npoint<<" points"<<endl; |
if(npoint != n)cout << "NB! Trajectory created with "<<npoint<<" points"<<endl; |
1167 |
} |
} |
1168 |
|
void Trajectory::Delete(){ |
1169 |
|
|
1170 |
|
if(x) delete [] x; |
1171 |
|
if(y) delete [] y; |
1172 |
|
if(z) delete [] z; |
1173 |
|
if(thx) delete [] thx; |
1174 |
|
if(thy) delete [] thy; |
1175 |
|
if(tl) delete [] tl; |
1176 |
|
|
1177 |
|
} |
1178 |
//-------------------------------------- |
//-------------------------------------- |
1179 |
// |
// |
1180 |
// |
// |
1213 |
|
|
1214 |
} |
} |
1215 |
|
|
1216 |
|
/** |
1217 |
|
* Evaluates the trajectory in the apparatus associated to the track. |
1218 |
|
* 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. |
1219 |
|
* @param t pointer to an object of the class Trajectory, |
1220 |
|
* which z coordinates should be previously initialized by calling the proper constructor ( Trajectory::Trajectory(int n, float* zin) ). |
1221 |
|
* @return error flag. |
1222 |
|
*/ |
1223 |
|
int Trajectory::DoTrack2(float* al){ |
1224 |
|
|
1225 |
|
double *dxout = new double[npoint]; |
1226 |
|
double *dyout = new double[npoint]; |
1227 |
|
double *dthxout = new double[npoint]; |
1228 |
|
double *dthyout = new double[npoint]; |
1229 |
|
double *dtlout = new double[npoint]; |
1230 |
|
double *dzin = new double[npoint]; |
1231 |
|
double dal[5]; |
1232 |
|
|
1233 |
|
int ifail = 0; |
1234 |
|
|
1235 |
|
for (int i=0; i<5; i++) dal[i] = (double)al[i]; |
1236 |
|
for (int i=0; i<npoint; i++) dzin[i] = (double)z[i]; |
1237 |
|
|
1238 |
|
dotrack2_(&(npoint),dzin,dxout,dyout,dthxout,dthyout,dtlout,dal,&ifail); |
1239 |
|
|
1240 |
|
for (int i=0; i<npoint; i++){ |
1241 |
|
x[i] = (float)*dxout++; |
1242 |
|
y[i] = (float)*dyout++; |
1243 |
|
thx[i] = (float)*dthxout++; |
1244 |
|
thy[i] = (float)*dthyout++; |
1245 |
|
tl[i] = (float)*dtlout++; |
1246 |
|
} |
1247 |
|
|
1248 |
|
return ifail; |
1249 |
|
}; |
1250 |
|
|
1251 |
ClassImp(TrkLevel2); |
ClassImp(TrkLevel2); |
1252 |
ClassImp(TrkSinglet); |
ClassImp(TrkSinglet); |