/[PAMELA software]/DarthVader/TrackerLevel2/src/TrkLevel2.cpp
ViewVC logotype

Annotation of /DarthVader/TrackerLevel2/src/TrkLevel2.cpp

Parent Directory Parent Directory | Revision Log Revision Log


Revision 1.16 - (hide annotations) (download)
Wed Nov 8 16:42:28 2006 UTC (18 years, 1 month ago) by pam-fi
Branch: MAIN
Changes since 1.15: +21 -16 lines
fixed bug in readB

1 mocchiut 1.1 /**
2     * \file TrkLevel2.cpp
3     * \author Elena Vannuccini
4     */
5     #include <TrkLevel2.h>
6     #include <iostream>
7 pam-fi 1.14 #include <math.h>
8 mocchiut 1.1 using namespace std;
9     //......................................
10     // F77 routines
11     //......................................
12     extern "C" {
13     void dotrack_(int*, double*, double*, double*, double*, int*);
14 pam-fi 1.2 void dotrack2_(int*, double*, double*, double*, double*,double*, double*, double*,int*);
15 pam-fi 1.16 // int readb_(const char*);
16     int readb_();
17 pam-fi 1.13 void mini2_(int*,int*,int*);
18 pam-fi 1.15 void guess_();
19 mocchiut 1.1 }
20     //--------------------------------------
21     //
22     //
23     //--------------------------------------
24     TrkTrack::TrkTrack(){
25 pam-fi 1.3 seqno = -1;
26     image = -1;
27 mocchiut 1.1 chi2 = 0;
28 pam-fi 1.14 nstep = 0;
29     for(int it1=0;it1<5;it1++){
30     al[it1] = 0;
31     for(int it2=0;it2<5;it2++)coval[it1][it2] = 0;
32 mocchiut 1.1 };
33     for(int ip=0;ip<6;ip++){
34 pam-fi 1.9 xgood[ip] = 0;
35     ygood[ip] = 0;
36     xm[ip] = 0;
37     ym[ip] = 0;
38     zm[ip] = 0;
39     resx[ip] = 0;
40     resy[ip] = 0;
41     xv[ip] = 0;
42     yv[ip] = 0;
43     zv[ip] = 0;
44     axv[ip] = 0;
45     ayv[ip] = 0;
46     dedx_x[ip] = 0;
47     dedx_y[ip] = 0;
48     // clx[ip] = 0;
49     // cly[ip] = 0;
50     };
51     clx = new TRefArray(6,0);
52     cly = new TRefArray(6,0);
53 mocchiut 1.1 };
54     //--------------------------------------
55     //
56     //
57     //--------------------------------------
58     TrkTrack::TrkTrack(const TrkTrack& t){
59 pam-fi 1.3 seqno = t.seqno;
60 mocchiut 1.1 image = t.image;
61     chi2 = t.chi2;
62 pam-fi 1.14 nstep = t.nstep;
63     for(int it1=0;it1<5;it1++){
64     al[it1] = t.al[it1];
65     for(int it2=0;it2<5;it2++)coval[it1][it2] = t.coval[it1][it2];
66 mocchiut 1.1 };
67     for(int ip=0;ip<6;ip++){
68 pam-fi 1.9 xgood[ip] = t.xgood[ip];
69     ygood[ip] = t.ygood[ip];
70     xm[ip] = t.xm[ip];
71     ym[ip] = t.ym[ip];
72     zm[ip] = t.zm[ip];
73     resx[ip] = t.resx[ip];
74     resy[ip] = t.resy[ip];
75     xv[ip] = t.xv[ip];
76     yv[ip] = t.yv[ip];
77     zv[ip] = t.zv[ip];
78     axv[ip] = t.axv[ip];
79     ayv[ip] = t.ayv[ip];
80     dedx_x[ip] = t.dedx_x[ip];
81     dedx_y[ip] = t.dedx_y[ip];
82     // clx[ip] = 0;//<<<<pointer
83     // cly[ip] = 0;//<<<<pointer
84     };
85     clx = new TRefArray(*(t.clx));
86     cly = new TRefArray(*(t.cly));
87    
88 mocchiut 1.1 };
89     //--------------------------------------
90     //
91     //
92     //--------------------------------------
93     /**
94     * Evaluates the trajectory in the apparatus associated to the track.
95     * 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.
96     * @param t pointer to an object of the class Trajectory,
97     * which z coordinates should be previously initialized by calling the proper constructor ( Trajectory::Trajectory(int n, float* zin) ).
98     * @return error flag.
99     */
100     int TrkTrack::DoTrack(Trajectory* t){
101    
102     double *dxout = new double[t->npoint];
103     double *dyout = new double[t->npoint];
104     double *dzin = new double[t->npoint];
105     double dal[5];
106    
107     int ifail = 0;
108    
109     for (int i=0; i<5; i++) dal[i] = (double)al[i];
110     for (int i=0; i<t->npoint; i++) dzin[i] = (double)t->z[i];
111    
112     dotrack_(&(t->npoint),dzin,dxout,dyout,dal,&ifail);
113    
114     for (int i=0; i<t->npoint; i++){
115     t->x[i] = (float)*dxout++;
116     t->y[i] = (float)*dyout++;
117     }
118    
119     // delete [] dxout;
120     // delete [] dyout;
121     // delete [] dzin;
122    
123     return ifail;
124     };
125     //--------------------------------------
126     //
127     //
128     //--------------------------------------
129 pam-fi 1.2 /**
130     * Evaluates the trajectory in the apparatus associated to the track.
131     * 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.
132     * @param t pointer to an object of the class Trajectory,
133     * which z coordinates should be previously initialized by calling the proper constructor ( Trajectory::Trajectory(int n, float* zin) ).
134     * @return error flag.
135     */
136     int TrkTrack::DoTrack2(Trajectory* t){
137    
138     double *dxout = new double[t->npoint];
139     double *dyout = new double[t->npoint];
140     double *dthxout = new double[t->npoint];
141     double *dthyout = new double[t->npoint];
142     double *dtlout = new double[t->npoint];
143     double *dzin = new double[t->npoint];
144     double dal[5];
145    
146     int ifail = 0;
147    
148     for (int i=0; i<5; i++) dal[i] = (double)al[i];
149     for (int i=0; i<t->npoint; i++) dzin[i] = (double)t->z[i];
150    
151     dotrack2_(&(t->npoint),dzin,dxout,dyout,dthxout,dthyout,dtlout,dal,&ifail);
152    
153     for (int i=0; i<t->npoint; i++){
154     t->x[i] = (float)*dxout++;
155     t->y[i] = (float)*dyout++;
156     t->thx[i] = (float)*dthxout++;
157     t->thy[i] = (float)*dthyout++;
158     t->tl[i] = (float)*dtlout++;
159     }
160    
161     // delete [] dxout;
162     // delete [] dyout;
163     // delete [] dzin;
164    
165     return ifail;
166     };
167     //--------------------------------------
168     //
169     //
170     //--------------------------------------
171 mocchiut 1.1 //float TrkTrack::BdL(){
172     //};
173     //--------------------------------------
174     //
175     //
176     //--------------------------------------
177     Float_t TrkTrack::GetRigidity(){
178     Float_t rig=0;
179     if(chi2>0)rig=1./al[4];
180     if(rig<0) rig=-rig;
181     return rig;
182     };
183     //
184     Float_t TrkTrack::GetDeflection(){
185     Float_t def=0;
186     if(chi2>0)def=al[4];
187     return def;
188     };
189     //
190     Float_t TrkTrack::GetDEDX(){
191     Float_t dedx=0;
192     for(Int_t i=0; i<6; i++)dedx+=dedx_x[i]*xgood[i]+dedx_y[i]*ygood[i];
193     dedx = dedx/(this->GetNX()+this->GetNY());
194     return dedx;
195     };
196    
197     //--------------------------------------
198     //
199     //
200     //--------------------------------------
201     void TrkTrack::Dump(){
202     cout << endl << "========== Track " ;
203 pam-fi 1.13 cout << endl << "seq. n. : "<< seqno;
204     cout << endl << "image n. : "<< image;
205     cout << endl << "al : "; for(int i=0; i<5; i++)cout << al[i] << " ";
206 mocchiut 1.1 cout << endl << "chi^2 : "<< chi2;
207 pam-fi 1.13 cout << endl << "n.step : "<< nstep;
208     cout << endl << "xgood : "; for(int i=0; i<6; i++)cout << xgood[i] ;
209 mocchiut 1.1 cout << endl << "ygood : "; for(int i=0; i<6; i++)cout << ygood[i] ;
210     cout << endl << "xm : "; for(int i=0; i<6; i++)cout << xm[i] << " ";
211     cout << endl << "ym : "; for(int i=0; i<6; i++)cout << ym[i] << " ";
212     cout << endl << "zm : "; for(int i=0; i<6; i++)cout << zm[i] << " ";
213 pam-fi 1.13 cout << endl << "xv : "; for(int i=0; i<6; i++)cout << xv[i] << " ";
214     cout << endl << "yv : "; for(int i=0; i<6; i++)cout << yv[i] << " ";
215     cout << endl << "zv : "; for(int i=0; i<6; i++)cout << zv[i] << " ";
216     cout << endl << "resx : "; for(int i=0; i<6; i++)cout << resx[i] << " ";
217     cout << endl << "resy : "; for(int i=0; i<6; i++)cout << resy[i] << " ";
218     cout << endl << "coval : "; for(int i=0; i<5; i++)cout << coval[0][i]<<" ";
219     cout << endl << " "; for(int i=0; i<5; i++)cout << coval[1][i]<<" ";
220     cout << endl << " "; for(int i=0; i<5; i++)cout << coval[2][i]<<" ";
221     cout << endl << " "; for(int i=0; i<5; i++)cout << coval[3][i]<<" ";
222     cout << endl << " "; for(int i=0; i<5; i++)cout << coval[4][i]<<" ";
223 mocchiut 1.1 cout << endl << "dedx_x : "; for(int i=0; i<6; i++)cout << dedx_x[i] << " ";
224     cout << endl << "dedx_y : "; for(int i=0; i<6; i++)cout << dedx_y[i] << " ";
225 pam-fi 1.14 cout << endl;
226 mocchiut 1.1 }
227 pam-fi 1.13 /**
228     * Set the TrkTrack position measurements
229     */
230     void TrkTrack::SetMeasure(double *xmeas, double *ymeas, double *zmeas){
231     for(int i=0; i<6; i++) xm[i]=*xmeas++;
232     for(int i=0; i<6; i++) ym[i]=*ymeas++;
233     for(int i=0; i<6; i++) zm[i]=*zmeas++;
234     }
235     /**
236     * Set the TrkTrack position resolution
237     */
238     void TrkTrack::SetResolution(double *rx, double *ry){
239     for(int i=0; i<6; i++) resx[i]=*rx++;
240     for(int i=0; i<6; i++) resy[i]=*ry++;
241     }
242     /**
243     * Set the TrkTrack good measurement
244     */
245     void TrkTrack::SetGood(int *xg, int *yg){
246     for(int i=0; i<6; i++) xgood[i]=*xg++;
247     for(int i=0; i<6; i++) ygood[i]=*yg++;
248     }
249    
250     /**
251     * Load the magnetic field
252     */
253 pam-fi 1.16 void TrkTrack::LoadField(TString path){
254    
255     strcpy(path_.path,path.Data());
256     path_.pathlen = path.Length();
257     path_.error = 0;
258     readb_();
259    
260 pam-fi 1.13 };
261     /**
262     * Tracking method. It calls F77 mini routine.
263     */
264     void TrkTrack::Fit(double pfixed, int& fail, int iprint){
265 pam-fi 1.15
266     float al_ini[] = {0.,0.,0.,0.,0.};
267    
268 pam-fi 1.13 extern cMini2track track_;
269     fail = 0;
270 pam-fi 1.15
271 pam-fi 1.13 for(int i=0; i<6; i++) track_.xm[i]=xm[i];
272     for(int i=0; i<6; i++) track_.ym[i]=ym[i];
273     for(int i=0; i<6; i++) track_.zm[i]=zm[i];
274     for(int i=0; i<6; i++) track_.resx[i]=resx[i];
275     for(int i=0; i<6; i++) track_.resy[i]=resy[i];
276     for(int i=0; i<6; i++) track_.xgood[i]=xgood[i];
277     for(int i=0; i<6; i++) track_.ygood[i]=ygood[i];
278 pam-fi 1.15
279 pam-fi 1.13 // initial guess of "al" with linear fit
280 pam-fi 1.15 // if(al[0]==-9999.&&al[1]==-9999.&&al[2]==-9999.&&al[3]==-9999.&&al[4]==-9999.){
281     // cout << "initial guess "<<endl;
282     // double szz=0., szx=0., szy=0., ssx=0., ssy=0., sz=0., s1=0.;
283     // double det, ax, ay, bx, by;
284     // for(int i=0; i<NPLANE; i++) {
285     // szz=szz+zm[i]*zm[i];
286     // szx=szx+zm[i]*xm[i];
287     // szy=szy+zm[i]*ym[i];
288     // ssx=ssx+xm[i];
289     // ssy=ssy+ym[i];
290     // sz=sz+zm[i];
291     // s1=s1+1.;
292     // }
293     // det=szz*s1-sz*sz;
294     // ax=(szx*s1-sz*ssx)/det;
295     // bx=(szz*ssx-szx*sz)/det;
296     // ay=(szy*s1-sz*ssy)/det;
297     // by=(szz*ssy-szy*sz)/det;
298     // al[0]=ax*23.5+bx; // ZINI = 23.5 !!! it should be the same parameter in all codes
299     // al[1]=ay*23.5+by; // "
300     // al[2]= sqrt(pow(ax,2)+pow(ay,2))/ sqrt(pow(ax,2)+pow(ay,2)+1.);
301     // al[3]=0.;
302     // if( (ax!=0.)||(ay!=0.) ) {
303     // al[3]= asin(ay/ sqrt(pow(ax,2)+pow(ay,2)));
304     // if(ax<0.) al[3]=acos(-1.)-al[3];
305     // }
306     // al[4]=0.;
307     // }
308 pam-fi 1.13 // end guess
309 pam-fi 1.16
310 pam-fi 1.15 for(int i=0; i<5; i++) track_.al[i]=al[i];
311 pam-fi 1.16 track_.zini = 23.5;
312     // ZINI = 23.5 !!! it should be the same parameter in all codes
313 pam-fi 1.15
314     if(al[0]==-9999.&&al[1]==-9999.&&al[2]==-9999.&&al[3]==-9999.&&al[4]==-9999.)guess_();
315    
316 pam-fi 1.16 // --------------------- free momentum
317 pam-fi 1.14 if(pfixed==0.) {
318 pam-fi 1.15 // al[4]=0.; // free momentum
319 pam-fi 1.14 track_.pfixed=0.; // "
320     }
321 pam-fi 1.16 // --------------------- fixed momentum
322 pam-fi 1.13 if(pfixed!=0.) {
323     al[4]=1./pfixed; // to fix the momentum
324     track_.pfixed=pfixed; // "
325     }
326 pam-fi 1.15
327     // store temporarily the initial guess
328     for(int i=0; i<5; i++) al_ini[i]=track_.al[i];
329    
330    
331 pam-fi 1.13 int istep=0;
332     int ifail=0;
333     mini2_(&istep,&ifail, &iprint);
334     if(ifail!=0) {
335     if(iprint==1)cout << "ERROR: ifail= " << ifail << endl;
336     fail = 1;
337     // return;
338     }
339 pam-fi 1.15
340 pam-fi 1.13
341     // cout << endl << "eta ===> " << track_.al[4] << endl;
342    
343     for(int i=0; i<5; i++) al[i]=track_.al[i];
344     chi2=track_.chi2;
345     nstep=track_.nstep;
346     for(int i=0; i<6; i++) xv[i]=track_.xv[i];
347     for(int i=0; i<6; i++) yv[i]=track_.yv[i];
348     for(int i=0; i<6; i++) zv[i]=track_.zv[i];
349     for(int i=0; i<6; i++) axv[i]=track_.axv[i];
350     for(int i=0; i<6; i++) ayv[i]=track_.ayv[i];
351     for(int i=0; i<5; i++) {
352     for(int j=0; j<5; j++) coval[i][j]=track_.cov[i][j];
353     }
354 pam-fi 1.15
355     if(fail){
356     cout << " >>>> fit failed >>>> drawing initial par"<<endl;
357     for(int i=0; i<5; i++) al[i]=al_ini[i];
358     }
359    
360 pam-fi 1.13 };
361     /*
362 pam-fi 1.15 * Reset the fit parameters
363 pam-fi 1.13 */
364     void TrkTrack::FitReset(){
365     for(int i=0; i<5; i++) al[i]=-9999.;
366     chi2=0.;
367     nstep=0;
368     for(int i=0; i<6; i++) xv[i]=0.;
369     for(int i=0; i<6; i++) yv[i]=0.;
370     for(int i=0; i<6; i++) zv[i]=0.;
371     for(int i=0; i<6; i++) axv[i]=0.;
372     for(int i=0; i<6; i++) ayv[i]=0.;
373     for(int i=0; i<5; i++) {
374     for(int j=0; j<5; j++) coval[i][j]=0.;
375     }
376     }
377    
378 mocchiut 1.1 //--------------------------------------
379     //
380     //
381     //--------------------------------------
382 pam-fi 1.10 void TrkTrack::Clear(){
383     seqno = -1;
384     image = -1;
385     chi2 = 0;
386     nstep = 0;
387     for(int it1=0;it1<5;it1++){
388     al[it1] = 0;
389     for(int it2=0;it2<5;it2++)coval[it1][it2] = 0;
390     };
391     for(int ip=0;ip<6;ip++){
392     xgood[ip] = 0;
393     ygood[ip] = 0;
394     xm[ip] = 0;
395     ym[ip] = 0;
396     zm[ip] = 0;
397     resx[ip] = 0;
398     resy[ip] = 0;
399     xv[ip] = 0;
400     yv[ip] = 0;
401     zv[ip] = 0;
402     axv[ip] = 0;
403     ayv[ip] = 0;
404     dedx_x[ip] = 0;
405     dedx_y[ip] = 0;
406     // clx[ip] = 0;
407     // cly[ip] = 0;
408     };
409     clx->Clear();
410     cly->Clear();
411     };
412     //--------------------------------------
413     //
414     //
415     //--------------------------------------
416 pam-fi 1.11 void TrkTrack::Delete(){
417     Clear();
418     clx->Delete();
419     cly->Delete();
420     };
421     //--------------------------------------
422     //
423     //
424     //--------------------------------------
425 pam-fi 1.10
426     //--------------------------------------
427     //
428     //
429     //--------------------------------------
430 mocchiut 1.1 TrkSinglet::TrkSinglet(){
431     plane = 0;
432     coord[0] = 0;
433     coord[1] = 0;
434     sgnl = 0;
435 pam-fi 1.9 cls = 0;
436 mocchiut 1.1 };
437     //--------------------------------------
438     //
439     //
440     //--------------------------------------
441     TrkSinglet::TrkSinglet(const TrkSinglet& s){
442     plane = s.plane;
443     coord[0] = s.coord[0];
444     coord[1] = s.coord[1];
445     sgnl = s.sgnl;
446 pam-fi 1.9 // cls = 0;//<<<<pointer
447     cls = TRef(s.cls);
448 mocchiut 1.1 };
449     //--------------------------------------
450     //
451     //
452     //--------------------------------------
453     void TrkSinglet::Dump(){
454     int i=0;
455     cout << endl << "========== Singlet " ;
456     cout << endl << "plane : " << plane;
457     cout << endl << "coord[2] : "; while( i<2 && cout << coord[i] << " ") i++;
458     cout << endl << "sgnl : " << sgnl;
459     }
460     //--------------------------------------
461     //
462     //
463     //--------------------------------------
464     TrkLevel2::TrkLevel2(){
465 pam-fi 1.10 // good2 = -1;
466 mocchiut 1.1 for(Int_t i=0; i<12 ; i++){
467 pam-fi 1.10 // crc[i] = -1;
468     good[i] = -1;
469     };
470 mocchiut 1.1 Track = new TClonesArray("TrkTrack");
471     SingletX = new TClonesArray("TrkSinglet");
472     SingletY = new TClonesArray("TrkSinglet");
473 pam-fi 1.6
474     // PhysicalTrack = new TClonesArray("TrkTrack");
475     //sostituire con TRefArray... appena ho capito come si usa
476 mocchiut 1.1 }
477     //--------------------------------------
478     //
479     //
480     //--------------------------------------
481     void TrkLevel2::Dump(){
482 pam-fi 1.10
483 mocchiut 1.1 TClonesArray &t = *Track;
484     TClonesArray &sx = *SingletX;
485     TClonesArray &sy = *SingletY;
486 pam-fi 1.10 //
487 mocchiut 1.1 cout << endl << endl << "=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-=-";
488 pam-fi 1.10 cout << endl << "good : "; for(int i=0; i<12; i++) cout << good[i]<<" ";
489 mocchiut 1.1 cout << endl << "ntrk() : " << this->ntrk() ;
490     cout << endl << "nclsx() : " << this->nclsx();
491     cout << endl << "nclsy() : " << this->nclsy();
492     for(int i=0; i<this->ntrk(); i++) ((TrkTrack *)t[i])->Dump();
493     for(int i=0; i<this->nclsx(); i++) ((TrkSinglet *)sx[i])->Dump();
494     for(int i=0; i<this->nclsy(); i++) ((TrkSinglet *)sy[i])->Dump();
495     }
496     //--------------------------------------
497     //
498     //
499     //--------------------------------------
500     /**
501     * Fills a TrkLevel2 object with values from a struct cTrkLevel2 (to get data from F77 common).
502     */
503 pam-fi 1.7 void TrkLevel2::SetFromLevel2Struct(cTrkLevel2 *l2){
504 pam-fi 1.10
505     // temporary objects:
506 mocchiut 1.1 TrkSinglet* t_singlet = new TrkSinglet();
507     TrkTrack* t_track = new TrkTrack();
508 pam-fi 1.10
509     // **** general variables ****
510     // good2 = l2->good2;
511 mocchiut 1.1 for(Int_t i=0; i<12 ; i++){
512 pam-fi 1.10 // crc[i] = l2->crc[i];
513     good[i] = l2->good[i];
514     };
515     // *** TRACKS ***
516 mocchiut 1.1 TClonesArray &t = *Track;
517     for(int i=0; i<l2->ntrk; i++){
518 pam-fi 1.9 t_track->seqno = i;// NBNBNBNB deve sempre essere = i
519     t_track->image = l2->image[i]-1;
520     // cout << "track "<<i<<t_track->seqno << t_track->image<<endl;
521     t_track->chi2 = l2->chi2_nt[i];
522 pam-fi 1.10 t_track->nstep = l2->nstep_nt[i];
523 pam-fi 1.9 for(int it1=0;it1<5;it1++){
524     t_track->al[it1] = l2->al_nt[i][it1];
525     for(int it2=0;it2<5;it2++)
526     t_track->coval[it1][it2] = l2->coval[i][it2][it1];
527     };
528     for(int ip=0;ip<6;ip++){
529     t_track->xgood[ip] = l2->xgood_nt[i][ip];
530     t_track->ygood[ip] = l2->ygood_nt[i][ip];
531     t_track->xm[ip] = l2->xm_nt[i][ip];
532     t_track->ym[ip] = l2->ym_nt[i][ip];
533     t_track->zm[ip] = l2->zm_nt[i][ip];
534     t_track->resx[ip] = l2->resx_nt[i][ip];
535     t_track->resy[ip] = l2->resy_nt[i][ip];
536     t_track->xv[ip] = l2->xv_nt[i][ip];
537     t_track->yv[ip] = l2->yv_nt[i][ip];
538     t_track->zv[ip] = l2->zv_nt[i][ip];
539     t_track->axv[ip] = l2->axv_nt[i][ip];
540     t_track->ayv[ip] = l2->ayv_nt[i][ip];
541     t_track->dedx_x[ip] = l2->dedx_x[i][ip];
542     t_track->dedx_y[ip] = l2->dedx_y[i][ip];
543     // t_track->clx[ip] = 0;
544     // t_track->cly[ip] = 0;
545     };
546     new(t[i]) TrkTrack(*t_track);
547     t_track->Clear();
548 mocchiut 1.1 };
549     // *** SINGLETS ***
550     TClonesArray &sx = *SingletX;
551     for(int i=0; i<l2->nclsx; i++){
552 pam-fi 1.9 t_singlet->plane = l2->planex[i];
553     t_singlet->coord[0] = l2->xs[i][0];
554     t_singlet->coord[1] = l2->xs[i][1];
555     t_singlet->sgnl = l2->signlxs[i];
556     // t_singlet->cls = 0;
557     new(sx[i]) TrkSinglet(*t_singlet);
558     t_singlet->Clear();
559 mocchiut 1.1 }
560     TClonesArray &sy = *SingletY;
561     for(int i=0; i<l2->nclsy; i++){
562 pam-fi 1.9 t_singlet->plane = l2->planey[i];
563     t_singlet->coord[0] = l2->ys[i][0];
564     t_singlet->coord[1] = l2->ys[i][1];
565     t_singlet->sgnl = l2->signlys[i];
566     // t_singlet->cls = 0;
567     new(sy[i]) TrkSinglet(*t_singlet);
568     t_singlet->Clear();
569     };
570    
571     delete t_track;
572     delete t_singlet;
573     }
574     //--------------------------------------
575     //
576     //
577     //--------------------------------------
578     /**
579     * Fills a TrkLevel2 object with values from a struct cTrkLevel2 (to get data from F77 common).
580     * Ref to Level1 data (clusters) is also set.
581     */
582     void TrkLevel2::SetFromLevel2Struct(cTrkLevel2 *l2, TrkLevel1 *l1){
583    
584     // temporary objects:
585     TrkSinglet* t_singlet = new TrkSinglet();
586     TrkTrack* t_track = new TrkTrack();
587     // general variables
588 pam-fi 1.10 // good2 = l2->good2;
589 pam-fi 1.9 for(Int_t i=0; i<12 ; i++){
590 pam-fi 1.10 // crc[i] = l2->crc[i];
591     good[i] = l2->good[i];
592 pam-fi 1.9 };
593     // *** TRACKS ***
594     TClonesArray &t = *Track;
595     for(int i=0; i<l2->ntrk; i++){
596     t_track->seqno = i;// NBNBNBNB deve sempre essere = i
597     t_track->image = l2->image[i]-1;
598     // cout << "track "<<i<<t_track->seqno << t_track->image<<endl;
599     t_track->chi2 = l2->chi2_nt[i];
600 pam-fi 1.10 t_track->nstep = l2->nstep_nt[i];
601 pam-fi 1.9 for(int it1=0;it1<5;it1++){
602     t_track->al[it1] = l2->al_nt[i][it1];
603     for(int it2=0;it2<5;it2++)
604     t_track->coval[it1][it2] = l2->coval[i][it2][it1];
605     };
606     for(int ip=0;ip<6;ip++){
607     t_track->xgood[ip] = l2->xgood_nt[i][ip];
608     t_track->ygood[ip] = l2->ygood_nt[i][ip];
609     t_track->xm[ip] = l2->xm_nt[i][ip];
610     t_track->ym[ip] = l2->ym_nt[i][ip];
611     t_track->zm[ip] = l2->zm_nt[i][ip];
612     t_track->resx[ip] = l2->resx_nt[i][ip];
613     t_track->resy[ip] = l2->resy_nt[i][ip];
614     t_track->xv[ip] = l2->xv_nt[i][ip];
615     t_track->yv[ip] = l2->yv_nt[i][ip];
616     t_track->zv[ip] = l2->zv_nt[i][ip];
617     t_track->axv[ip] = l2->axv_nt[i][ip];
618     t_track->ayv[ip] = l2->ayv_nt[i][ip];
619     t_track->dedx_x[ip] = l2->dedx_x[i][ip];
620     t_track->dedx_y[ip] = l2->dedx_y[i][ip];
621     // cout << "traccia "<<i<<" -- "<< ip << " "<< l2->cltrx[i][ip] <<" "<< l2->cltry[i][ip] <<" "<< t_track->xgood[ip] << t_track->ygood[ip]<<endl;
622     //-----------------------------------------------------
623     // t_track->clx[ip] = l1->GetCluster(l2->cltrx[i][ip]-1);
624     // t_track->cly[ip] = l1->GetCluster(l2->cltry[i][ip]-1);
625     if(t_track->xgood[ip])t_track->clx->AddAt(l1->GetCluster(l2->cltrx[i][ip]-1),ip);
626     if(t_track->ygood[ip])t_track->cly->AddAt(l1->GetCluster(l2->cltry[i][ip]-1),ip);
627 pam-fi 1.10 // 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;
628     // 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;
629 pam-fi 1.9 //-----------------------------------------------------
630     };
631     new(t[i]) TrkTrack(*t_track);
632     t_track->Clear();
633     };
634     // *** SINGLETS ***
635     TClonesArray &sx = *SingletX;
636     for(int i=0; i<l2->nclsx; i++){
637     t_singlet->plane = l2->planex[i];
638     t_singlet->coord[0] = l2->xs[i][0];
639     t_singlet->coord[1] = l2->xs[i][1];
640     t_singlet->sgnl = l2->signlxs[i];
641     //-----------------------------------------------------
642     // cout << "singolo x "<<i<<" -- "<< l2->clsx[i] <<endl;
643     t_singlet->cls = l1->GetCluster(l2->clsx[i]-1);
644 pam-fi 1.10 // cout<<" i "<<i<<" l2->clsx[i] "<<l2->clsx[i]<< " l1->GetCluster(l2->clsx[i]-1) "<<l1->GetCluster(l2->clsx[i]-1)<<endl;
645 pam-fi 1.9 //-----------------------------------------------------
646     new(sx[i]) TrkSinglet(*t_singlet);
647     t_singlet->Clear();
648     }
649     TClonesArray &sy = *SingletY;
650     for(int i=0; i<l2->nclsy; i++){
651     t_singlet->plane = l2->planey[i];
652     t_singlet->coord[0] = l2->ys[i][0];
653     t_singlet->coord[1] = l2->ys[i][1];
654     t_singlet->sgnl = l2->signlys[i];
655     //-----------------------------------------------------
656     // cout << "singolo y "<<i<<" -- "<< l2->clsy[i] <<endl;
657     t_singlet->cls = l1->GetCluster(l2->clsy[i]-1);
658 pam-fi 1.10 // cout<<" i "<<i<<" l2->clsy[i] "<<l2->clsy[i]<< " l1->GetCluster(l2->clsy[i]-1) "<<l1->GetCluster(l2->clsy[i]-1)<<endl;
659 pam-fi 1.9 //-----------------------------------------------------
660     new(sy[i]) TrkSinglet(*t_singlet);
661     t_singlet->Clear();
662 mocchiut 1.1 };
663 pam-fi 1.5
664     delete t_track;
665     delete t_singlet;
666 mocchiut 1.1 }
667 pam-fi 1.7 /**
668     * Fills a struct cTrkLevel2 with values from a TrkLevel2 object (to put data into a F77 common).
669     */
670    
671     void TrkLevel2::GetLevel2Struct(cTrkLevel2 *l2) const {
672    
673     // general variables
674 pam-fi 1.10 // l2->good2 = good2 ;
675 pam-fi 1.7 for(Int_t i=0; i<12 ; i++){
676 pam-fi 1.10 // l2->crc[i] = crc[i];
677     l2->good[i] = good[i];
678 pam-fi 1.7 };
679     // *** TRACKS ***
680    
681     l2->ntrk = Track->GetEntries();
682     for(Int_t i=0;i<l2->ntrk;i++){
683     l2->image[i] = 1 + ((TrkTrack *)Track->At(i))->image;
684     l2->chi2_nt[i] = ((TrkTrack *)Track->At(i))->chi2;
685 pam-fi 1.10 l2->nstep_nt[i] = ((TrkTrack *)Track->At(i))->nstep;
686     for(int it1=0;it1<5;it1++){
687 pam-fi 1.7 l2->al_nt[i][it1] = ((TrkTrack *)Track->At(i))->al[it1];
688     for(int it2=0;it2<5;it2++)
689     l2->coval[i][it2][it1] = ((TrkTrack *)Track->At(i))->coval[it1][it2];
690     };
691     for(int ip=0;ip<6;ip++){
692     l2->xgood_nt[i][ip] = ((TrkTrack *)Track->At(i))->xgood[ip];
693     l2->ygood_nt[i][ip] = ((TrkTrack *)Track->At(i))->ygood[ip];
694     l2->xm_nt[i][ip] = ((TrkTrack *)Track->At(i))->xm[ip];
695     l2->ym_nt[i][ip] = ((TrkTrack *)Track->At(i))->ym[ip];
696     l2->zm_nt[i][ip] = ((TrkTrack *)Track->At(i))->zm[ip];
697     l2->resx_nt[i][ip] = ((TrkTrack *)Track->At(i))->resx[ip];
698     l2->resy_nt[i][ip] = ((TrkTrack *)Track->At(i))->resy[ip];
699     l2->xv_nt[i][ip] = ((TrkTrack *)Track->At(i))->xv[ip];
700     l2->yv_nt[i][ip] = ((TrkTrack *)Track->At(i))->yv[ip];
701     l2->zv_nt[i][ip] = ((TrkTrack *)Track->At(i))->zv[ip];
702     l2->axv_nt[i][ip] = ((TrkTrack *)Track->At(i))->axv[ip];
703     l2->ayv_nt[i][ip] = ((TrkTrack *)Track->At(i))->ayv[ip];
704     l2->dedx_x[i][ip] = ((TrkTrack *)Track->At(i))->dedx_x[ip];
705     l2->dedx_y[i][ip] = ((TrkTrack *)Track->At(i))->dedx_y[ip];
706     };
707     }
708    
709     // *** SINGLETS ***
710     l2->nclsx = SingletX->GetEntries();
711     for(Int_t i=0;i<l2->nclsx;i++){
712     l2->planex[i] = ((TrkSinglet *)SingletX->At(i))->plane;
713     l2->xs[i][0] = ((TrkSinglet *)SingletX->At(i))->coord[0];
714     l2->xs[i][1] = ((TrkSinglet *)SingletX->At(i))->coord[1];
715     l2->signlxs[i] = ((TrkSinglet *)SingletX->At(i))->sgnl;
716     }
717     l2->nclsy = SingletY->GetEntries();
718     for(Int_t i=0;i<l2->nclsy;i++){
719     l2->planey[i] = ((TrkSinglet *)SingletY->At(i))->plane;
720     l2->ys[i][0] = ((TrkSinglet *)SingletY->At(i))->coord[0];
721     l2->ys[i][1] = ((TrkSinglet *)SingletY->At(i))->coord[1];
722     l2->signlys[i] = ((TrkSinglet *)SingletY->At(i))->sgnl;
723     }
724     }
725 mocchiut 1.1 //--------------------------------------
726     //
727     //
728     //--------------------------------------
729     void TrkLevel2::Clear(){
730 pam-fi 1.10 // good2 = -1;
731 mocchiut 1.1 for(Int_t i=0; i<12 ; i++){
732 pam-fi 1.10 // crc[i] = -1;
733     good[i] = -1;
734     };
735 pam-fi 1.5 /* Track->RemoveAll();
736 mocchiut 1.1 SingletX->RemoveAll();
737 pam-fi 1.5 SingletY->RemoveAll();*/
738     // modify to avoid memory leakage
739     Track->Clear();
740     SingletX->Clear();
741     SingletY->Clear();
742 mocchiut 1.1 }
743     //--------------------------------------
744     //
745     //
746     //--------------------------------------
747 pam-fi 1.11 void TrkLevel2::Delete(){
748    
749     Clear();
750     Track->Delete();
751     SingletX->Delete();
752     SingletY->Delete();
753     }
754     //--------------------------------------
755     //
756     //
757     //--------------------------------------
758 mocchiut 1.1 /**
759     * Sort physical tracks and stores them in a TObjectArray, ordering by increasing chi**2 value (in case of track image, it selects the one with lower chi**2). The total number of physical tracks is given by GetNTracks() and the it-th physical track can be retrieved by means of the method GetTrack(int it).
760     * This method is overridden by PamLevel2::GetTracks(), where calorimeter and TOF information is used.
761     */
762 pam-fi 1.8 TRefArray *TrkLevel2::GetTracks_NFitSorted(){
763 pam-fi 1.3
764 pam-fi 1.8 TRefArray *sorted = new TRefArray();
765    
766     TClonesArray &t = *Track;
767     // TClonesArray &ts = *PhysicalTrack;
768     int N = ntrk();
769     vector<int> m(N); for(int i=0; i<N; i++)m[i]=1;
770     // int m[50]; for(int i=0; i<N; i++)m[i]=1;
771    
772     int indo=0;
773     int indi=0;
774     while(N != 0){
775     int nfit =0;
776     float chi2ref = numeric_limits<float>::max();
777 pam-fi 1.9
778 pam-fi 1.8 // first loop to search maximum num. of fit points
779     for(int i=0; i < ntrk(); i++){
780     if( ((TrkTrack *)t[i])->GetNtot() >= nfit && m[i]==1){
781     nfit = ((TrkTrack *)t[i])->GetNtot();
782     }
783     }
784     //second loop to search minimum chi2 among selected
785     for(int i=0; i<this->ntrk(); i++){
786 pam-fi 1.9 Float_t chi2 = ((TrkTrack *)t[i])->chi2;
787     if(chi2 < 0) chi2 = chi2*1000;
788     if( chi2 < chi2ref
789     && ((TrkTrack *)t[i])->GetNtot() == nfit
790     && m[i]==1){
791 pam-fi 1.8 chi2ref = ((TrkTrack *)t[i])->chi2;
792     indi = i;
793 pam-fi 1.9 };
794     };
795 pam-fi 1.8 if( ((TrkTrack *)t[indi])->HasImage() ){
796     m[((TrkTrack *)t[indi])->image] = 0;
797     N--;
798    
799     // cout << "i** "<< ((TrkTrack *)t[indi])->image << " " << nfiti <<" "<<chi2i<<endl;
800 pam-fi 1.9 };
801 pam-fi 1.8 sorted->Add( (TrkTrack*)t[indi] );
802 pam-fi 1.3
803 pam-fi 1.8 m[indi] = 0;
804     // cout << "SORTED "<< indo << " "<< indi << " "<< N << endl;
805     N--;
806     indo++;
807 pam-fi 1.3 }
808 pam-fi 1.8 m.clear();
809     // cout << "GetTracks_NFitSorted(it): Done"<< endl;
810    
811     return sorted;
812 pam-fi 1.6 // return PhysicalTrack;
813 pam-fi 1.3 }
814 mocchiut 1.1 //--------------------------------------
815     //
816     //
817     //--------------------------------------
818     /**
819     * Retrieves the is-th stored track.
820     * @param it Track number, ranging from 0 to ntrk().
821     * Fitted tracks ( images included ) are stored in a TObjectArray ( TrkLevel2::Track ) in the same order they are returned by the F77 fitting routine.
822     */
823     TrkTrack *TrkLevel2::GetStoredTrack(int is){
824    
825     if(is >= this->ntrk()){
826     cout << "** TrkLevel2 ** Track "<< is << "doen not exits! " << endl;
827     cout << " Stored tracks ntrk() = "<< this->ntrk() << endl;
828     return 0;
829     }
830     TClonesArray &t = *(Track);
831     TrkTrack *track = (TrkTrack*)t[is];
832     return track;
833     }
834     //--------------------------------------
835     //
836     //
837     //--------------------------------------
838     /**
839 pam-fi 1.6 * Retrieves the is-th stored X singlet.
840     * @param it Singlet number, ranging from 0 to nclsx().
841     */
842     TrkSinglet *TrkLevel2::GetSingletX(int is){
843    
844     if(is >= this->nclsx()){
845     cout << "** TrkLevel2 ** Singlet "<< is << "doen not exits! " << endl;
846     cout << " Stored x-singlets nclsx() = "<< this->nclsx() << endl;
847     return 0;
848     }
849     TClonesArray &t = *(SingletX);
850     TrkSinglet *singlet = (TrkSinglet*)t[is];
851     return singlet;
852     }
853     //--------------------------------------
854     //
855     //
856     //--------------------------------------
857     /**
858     * Retrieves the is-th stored Y singlet.
859     * @param it Singlet number, ranging from 0 to nclsx().
860     */
861     TrkSinglet *TrkLevel2::GetSingletY(int is){
862    
863     if(is >= this->nclsy()){
864     cout << "** TrkLevel2 ** Singlet "<< is << "doen not exits! " << endl;
865     cout << " Stored y-singlets nclsy() = "<< this->nclsx() << endl;
866     return 0;
867     }
868     TClonesArray &t = *(SingletY);
869     TrkSinglet *singlet = (TrkSinglet*)t[is];
870     return singlet;
871     }
872     //--------------------------------------
873     //
874     //
875     //--------------------------------------
876     /**
877 mocchiut 1.1 * Retrieves the it-th "physical" track, sorted by the method GetNTracks().
878     * @param it Track number, ranging from 0 to GetNTracks().
879     */
880 pam-fi 1.10
881 pam-fi 1.8 TrkTrack *TrkLevel2::GetTrack(int it){
882    
883     if(it >= this->GetNTracks()){
884     cout << "** TrkLevel2 ** Track "<< it << "does not exits! " << endl;
885     cout << " Physical tracks GetNTracks() = "<< this->ntrk() << endl;
886     return 0;
887     }
888    
889     TRefArray *sorted = GetTracks(); //TEMPORANEO
890     TrkTrack *track = (TrkTrack*)sorted->At(it);
891     sorted->Delete();
892     return track;
893 mocchiut 1.1 }
894 pam-fi 1.6 /**
895     * Give the number of "physical" tracks, sorted by the method GetTracks().
896     */
897 pam-fi 1.5 Int_t TrkLevel2::GetNTracks(){
898 pam-fi 1.8
899     Float_t ntot=0;
900     TClonesArray &t = *Track;
901 mocchiut 1.12 for(int i=0; i<ntrk(); i++) {
902 pam-fi 1.8 if( ((TrkTrack *)t[i])->GetImageSeqNo() == -1 ) ntot+=1.;
903     else ntot+=0.5;
904     }
905     return (Int_t)ntot;
906    
907 pam-fi 1.5 };
908 mocchiut 1.1 //--------------------------------------
909     //
910     //
911     //--------------------------------------
912     /**
913     * Retrieves (if present) the image of the it-th "physical" track, sorted by the method GetNTracks().
914     * @param it Track number, ranging from 0 to GetNTracks().
915     */
916 pam-fi 1.8 TrkTrack *TrkLevel2::GetTrackImage(int it){
917    
918     if(it >= this->GetNTracks()){
919     cout << "** TrkLevel2 ** Track "<< it << "does not exits! " << endl;
920     cout << " Physical tracks GetNTracks() = "<< this->ntrk() << endl;
921     return 0;
922     }
923    
924     TRefArray* sorted = GetTracks(); //TEMPORANEO
925     TrkTrack *track = (TrkTrack*)sorted->At(it);
926    
927     if(!track->HasImage()){
928     cout << "** TrkLevel2 ** Track "<< it << "does not have image! " << endl;
929     return 0;
930     }
931     TrkTrack *image = (TrkTrack*)(*Track)[track->image];
932    
933     sorted->Delete();
934    
935     return image;
936    
937 mocchiut 1.1 }
938     //--------------------------------------
939     //
940     //
941     //--------------------------------------
942     /**
943     * Loads the magnetic field.
944     * @param s Path of the magnetic-field files.
945     */
946 pam-fi 1.16 void TrkLevel2::LoadField(TString path){
947     //
948     strcpy(path_.path,path.Data());
949     path_.pathlen = path.Length();
950     path_.error = 0;
951     readb_();
952     //
953 mocchiut 1.1 };
954     //--------------------------------------
955     //
956     //
957     //--------------------------------------
958     /**
959 pam-fi 1.6 * Get tracker-plane (mechanical) z-coordinate
960     * @param plane_id plane index (1=TOP,2,3,4,5,6=BOTTOM)
961     */
962     Float_t TrkLevel2::GetZTrk(Int_t plane_id){
963     switch(plane_id){
964     case 1: return ZTRK1;
965     case 2: return ZTRK2;
966     case 3: return ZTRK3;
967     case 4: return ZTRK4;
968     case 5: return ZTRK5;
969     case 6: return ZTRK6;
970     default: return 0.;
971     };
972     };
973     //--------------------------------------
974     //
975     //
976     //--------------------------------------
977     /**
978 pam-fi 1.2 * Trajectory default constructor.
979     * (By default is created with z-coordinates inside the tracking volume)
980     */
981     Trajectory::Trajectory(){
982     npoint = 10;
983     x = new float[npoint];
984     y = new float[npoint];
985     z = new float[npoint];
986     thx = new float[npoint];
987     thy = new float[npoint];
988     tl = new float[npoint];
989 pam-fi 1.6 float dz = ((ZTRK1)-(ZTRK6))/(npoint-1);
990 pam-fi 1.2 for(int i=0; i<npoint; i++){
991     x[i] = 0;
992     y[i] = 0;
993 pam-fi 1.6 z[i] = (ZTRK1) - i*dz;
994 pam-fi 1.2 thx[i] = 0;
995     thy[i] = 0;
996     tl[i] = 0;
997     }
998     }
999     //--------------------------------------
1000     //
1001     //
1002     //--------------------------------------
1003     /**
1004 mocchiut 1.1 * Trajectory constructor.
1005 pam-fi 1.2 * (By default is created with z-coordinates inside the tracking volume)
1006 mocchiut 1.1 * \param n Number of points
1007     */
1008     Trajectory::Trajectory(int n){
1009 pam-fi 1.2 if(n<=0){
1010     cout << "NB! Trajectory must have at least 1 point >>> created with 10 points" << endl;
1011     n=10;
1012     }
1013 mocchiut 1.1 npoint = n;
1014     x = new float[npoint];
1015     y = new float[npoint];
1016     z = new float[npoint];
1017 pam-fi 1.2 thx = new float[npoint];
1018     thy = new float[npoint];
1019     tl = new float[npoint];
1020 pam-fi 1.6 float dz = ((ZTRK1)-(ZTRK6))/(npoint-1);
1021 mocchiut 1.1 for(int i=0; i<npoint; i++){
1022 pam-fi 1.2 x[i] = 0;
1023 mocchiut 1.1 y[i] = 0;
1024 pam-fi 1.6 z[i] = (ZTRK1) - i*dz;
1025 pam-fi 1.2 thx[i] = 0;
1026     thy[i] = 0;
1027     tl[i] = 0;
1028 mocchiut 1.1 }
1029     }
1030     //--------------------------------------
1031     //
1032     //
1033     //--------------------------------------
1034     /**
1035     * Trajectory constructor.
1036     * \param n Number of points
1037     * \param pz Pointer to float array, defining z coordinates
1038     */
1039     Trajectory::Trajectory(int n, float* zin){
1040 pam-fi 1.2 npoint = 10;
1041     if(n>0)npoint = n;
1042 mocchiut 1.1 x = new float[npoint];
1043     y = new float[npoint];
1044     z = new float[npoint];
1045 pam-fi 1.2 thx = new float[npoint];
1046     thy = new float[npoint];
1047     tl = new float[npoint];
1048     int i=0;
1049     do{
1050 pam-fi 1.9 x[i] = 0;
1051     y[i] = 0;
1052     z[i] = zin[i];
1053     thx[i] = 0;
1054     thy[i] = 0;
1055     tl[i] = 0;
1056     i++;
1057 pam-fi 1.2 }while(zin[i-1] > zin[i] && i < npoint);
1058     npoint=i;
1059     if(npoint != n)cout << "NB! Trajectory created with "<<npoint<<" points"<<endl;
1060 mocchiut 1.1 }
1061     //--------------------------------------
1062     //
1063     //
1064     //--------------------------------------
1065     /**
1066     * Dump the trajectory coordinates.
1067     */
1068     void Trajectory::Dump(){
1069     cout <<endl<< "Trajectory ========== "<<endl;
1070     for (int i=0; i<npoint; i++){
1071 pam-fi 1.2 cout << i <<" >> " << x[i] <<" "<< y[i] <<" "<< z[i] ;
1072     cout <<" -- " << thx[i] <<" "<< thy[i] ;
1073     cout <<" -- " << tl[i] << endl;
1074 mocchiut 1.1 };
1075     }
1076 pam-fi 1.2 //--------------------------------------
1077     //
1078     //
1079     //--------------------------------------
1080     /**
1081     * Get trajectory length between two points
1082     * @param ifirst first point (default 0)
1083     * @param ilast last point (default npoint)
1084     */
1085     float Trajectory::GetLength(int ifirst, int ilast){
1086     if( ifirst<0 ) ifirst = 0;
1087     if( ilast>=npoint) ilast = npoint-1;
1088     float l=0;
1089     for(int i=ifirst;i<=ilast;i++){
1090     l=l+tl[i];
1091     };
1092     if(z[ilast] > ZINI)l=l-tl[ilast];
1093     if(z[ifirst] < ZINI) l=l-tl[ifirst];
1094    
1095     return l;
1096 mocchiut 1.1
1097 pam-fi 1.2 }
1098 pam-fi 1.6
1099 pam-fi 1.10
1100 mocchiut 1.1 ClassImp(TrkLevel2);
1101     ClassImp(TrkSinglet);
1102     ClassImp(TrkTrack);
1103     ClassImp(Trajectory);

  ViewVC Help
Powered by ViewVC 1.1.23